Control Design Aspects

[pattacini]

In this issue, I'd like to summarize ideas, approaches, results, and details that will come out from the control design and that might impact for example the FW development.

Current Model

Here's the latest model incorporating the sensing part, where HSC Pressure sensor and Zephyr airflow sensor from Honeywell have been employed. The model represents the use of a AP proportional inlet valve and a CFB On/Off outlet valve.

Current Controller

Below, you'll find a screenshot of the Volume/Pressure Controller:

Volume Control @ 1 kHz

Volume control, we know, is not critical:

Pressure Control @ 1 kHz

Pressure control, instead, requires feedback filtering and input shaping:

Gain Scheduling

However, it seems that the inlet valve with the whole circuit do respond more slowly while closing. In fact, we can observe nonnull steady-state error when we reach the final constant pressure setpoint (i.e. when input-shaping is done). Thus, scheduling the integral gain could improve performance:

As a result, the final pressure error gets reduced:

Pressure Control @ 100 Hz

We can lower the control rate and run at 100 Hz, as per this requirement:

Of course, performance gets degraded as for example we would need to curb filtering action (by increasing the cutoff frequency from 5 Hz to 10 Hz) in order not to lose too much in terms of control responsivity.

cc @maggia80 @simeonedussoni @andreamerello @mircodisalvo @RossiFederico @LucaBottazzi @marcoaccame @Nicogene @traversaro @fiorisi @RAR1989

[andreamerello]

I'm not sure about this, but..

.. By reading the sensors documentation I got some clues that make me suspecting that flow meters and pressure sensors perform free-running sampling, possibly around 1KHz, by themselves ( asynchronously wrt the I2C transactions, and thus wrt the control loop).

If case this turns out to be true, even if we could read sensors faster, I'm wondering whether the internal sensors sampling frequency and the 1KHz control loop frequency could badly "beat" one wrt the other, ending up in getting repeated samples or skipped samples in a random fashion, and whether this could adversely affect the control..

Maybe if the control loop runs at 100Hz this aspect would also become less relevant (?)

[pattacini]

@andreamerello what you report could be definitely happening. When I simulated the model at 1 kHz control rate, I assumed ideal conditions. As you said, at 100 Hz this behavior should be mitigated substantially.

[pattacini]

The model has been completed with the addition of Stage-4 flow-meter and pressure sensors plus the HEPA filters F2 and F3.

Shortly I'll be adding the quick-disconnections drops as well.

[pattacini]

Quick-disconnects drops added in https://github.com/icub-tech-iit/ventilator/commit/2243eaaa261026a9645fb0e64625b0761485e94f.

---

[pattacini]

As an example, https://github.com/icub-tech-iit/ventilator/commit/36a85202251dd2f3100e2dd9939be23a29d6e76c generates this code.

cc @marcoaccame

[marcoaccame]

ok @pattacini i will soon run it on the target to profile its exec duration.

[pattacini]

Great 👍 Remember to set the variable controller_U.enable to true.

[simeonedussoni]

@marcoaccame will u use the tool u showed us this morning?

[pattacini]

@marcoaccame following up this post, I've managed to get all the required files now. Here's the new archive.

[marcoaccame]

@marcoaccame will u use the tool u showed us this morning?

yes. i am going to release a new demo003 in branch code which have the capability to show a scope-like diagram of thread duration. but i will post the picture of the duration of the controller_Obj.step() as soon as i have it.

[pattacini]

@marcoaccame please deal with this controller.zip instead, which now makes use of single-precision floats.

Simulation runs just fine!

[marcoaccame]

result 1 on double precision:

the duration of the single step that i measured in two different ways is:

23 usec as printed by this code

duration = m0:u023

and 6 usec as displayed by the scope when i removed the trace print from code.

the time seems too small:

• the code inside the controller is actually executed because the debugger steps inside,
• the above two measure methods have proven to be pretty accurate in the past even if compared with a real scope.
• the differences on the two measures can be explained with the fact that writing to trace port slows things a bit.

conclusion

i want to double check this duration by feeding controller_Obj.step() with real input data from the @pattacini's simulation and comparing what the MPU produces with the associated output of simulation. i will confirm the duration only when the produced output will match the target one.

[marcoaccame]

result 2 on single precision:

the duration of the single step is:

10 usec (sometimes 9) as printed by this code

duration = m0:u010

yes @pattacini , it is faster.

i have also tested it by removing complier optimization w/ -O0 in the matlab generated cpp files: 14 usec

duration = m0:u014

nevertheless i want to do the same test as in my previous comment: use input data from simulation

[pattacini]

I've updated the model to incorporate the latest curves provided by Camozzi and system responses are very satisfactory with both inlet valve types.

AP Inlet Valve

CP Inlet Valve

👉 Importantly, it seems that we no longer need gain scheduling of the integral term!

cc @maggia80 @simeonedussoni @RossiFederico @LucaBottazzi

@marcoaccame next I'll be gathering the signals you'd need to feed the software.

[pattacini]

The model has been slightly refactored – here's the latest snapshot:

Further, @marcoaccame you'll find in the Easter egg the archive containing the latest code along with a testing dataset, as promised 🥚

[pattacini]

I've integrated the latest indications from MVM about the pressure input profile plus an initial estimate of the tidal volume as measured from the integration of the inlet and outlet flows.

This is the MVM pressure profile with a snap onset and a gentle increase toward the target pressure:

---

Here's how the Simscape model copies it out:

cc @maggia80

PS: en passant, what's the use of S4 as it'll be always measuring PEEP? Safety?

[pattacini]

With the new pressure input profiles, it seems that we don't have enough margin w/ AP as we're very close to the 50 L/min flow limit.

In fact, with a faster request, we get saturated:

---

By contrast, CP behaves just fine − as expected − since it doesn't suffer from the same flow limitation:

[maggia80]

PS: en passant, what's the use of S4 as it'll be always measuring PEEP? Safety?

yes, indeed we moved the S4 before the valve in the schema V9.

[pattacini]

Here's the analysis w/:

• CP inlet valve along with its full dynamics (time constant of 0.5/4 s)
• FFW+PID fine-tuned
• 2 sets of input shaping parameters

As it is shown, we can now relax specifications on feedback filtering.

The behavior doesn't seem dramatic, although the CP works wide open (= 100% DC) for a small amount of time at the onset in order to compensate for the delay in the flow delivery. Could this be somehow critical?

Clearly, reducing the CP full-dynamics time constant could improve performance.

cc @maggia80 @simeonedussoni @traversaro

[marcoaccame]

Hi @pattacini, i opened my easter egg.

it contained proof that the step() of the controller executes in 10usec in the evaluation board with compiler optimization and in 14usec w/out optimization. or maybe slightly less because the retrieval of the delta time affects the measure.

i have fed the controller on the evaluation board with the input data you gave me and the controller on the board computed the correct output of the simulation.

hooray!

here is the trace that the board has printed:

->
... omissis 
@ S0:m114:u016 tCTRL executes the control
duration of controller_Obj.step() = m0:u010
verification controller_Obj.step() @ step # 10 is OK: IV_P_fbk = 2211.840088 -> IV_dc = 39.747360 [%FS] IV_P_ref = 5.421200 [cmH2O]  P_hsc = 6.934096 [cmH2O]  P_hsc_filt = 6.754559 [cmH2O]
... omissis 
<-

->
... omissis 
@ S0:m124:u016 tCTRL executes the control
duration of controller_Obj.step() = m0:u010
verification controller_Obj.step() @ step # 11 is OK: IV_P_fbk = 2228.224121 -> IV_dc = 39.456169 [%FS] IV_P_ref = 5.539822 [cmH2O]  P_hsc = 7.138040 [cmH2O]  P_hsc_filt = 6.778560 [cmH2O]
... omissis 
<-

and here is the code which i ran (soon to be committed in a new test004):

const AllActuators& compute(const AllSensors &s)
{
    static uint32_t step = 0;
    vnt::tools::dataFmt d;
    fdr->get(d); // retrieve the IV_P_fbk from the stored simulation input in the vnt::tools::Feeder

    controllerModelClass::ExtU_controller_T inputs;
    inputs.enable = true;
    inputs.IV_Vdot_fbk = 0.F;
    inputs.IV_P_fbk = (float)d; // assign the IV_P_fbk as input of the controller 
    inputs.params[0] = 0.F;
    inputs.params[1] = 160.F;
    controller_Obj.setExternalInputs(&inputs);        

    volatile vnt::core::Time t0 = vnt::core::now();
    controller_Obj.step();
//  vnt::core::delay(10*vnt::core::time1microsec); // used to verify the time w/ a proven delay function
    volatile vnt::core::Time t1 = vnt::core::now();
    vnt::bsp::trace::puts("duration of controller_Obj.step() = " + vnt::core::TimeFormatter(t1-t0).to_string(vnt::core::TimeFormatter::Mode::milli));

    const controllerModelClass::ExtY_controller_T& outputs = controller_Obj.getExternalOutputs();
    vnt::tools::Result r {outputs.IV_dc, outputs.IV_P_ref, outputs.signals[1], outputs.signals[2]}; // P_hsc_filt(i)

    volatile bool res = fdr->check(r); // compare the result w/ the stored simulation result in the vnt::tools::Feeder 

    vnt::bsp::trace::puts("verification controller_Obj.step() @ step # " + std::to_string(step) + std::string(res? " is OK" : " is KO") + ": IV_P_fbk = " + std::to_string(d) + " ->" + 
        " IV_dc = " + std::to_string(r.v0) + " [%FS]" + 
        " IV_P_ref = " + std::to_string(r.v1) + " [cmH2O]" +
        "  P_hsc = " + std::to_string(r.v2) +  " [cmH2O]" +
        "  P_hsc_filt = " + std::to_string(r.v3) +  " [cmH2O]"
    );

    step++;

    status = s.s1.pressure + s.s2.pressure  + s.s3.pressure;       
    ac.v1 = status + 1;
    ac.v2 = status + 2;
    return ac;
}          

[pattacini]

That's a great news @marcoaccame 🚀

[simeonedussoni]

@pattacini in the plots reported above seems there is an effect of decreasing measured tidal volume which drops below zero, do you have any idea on its origin?

[pattacini]

@simeonedussoni don't get scared, the red curve just represents the very first attempt at gauging the tidal volume by integrating the measured flows. I'd need some time to look deeper at it.

[simeonedussoni]

ok @pattacini this was my guess thanks for confirmation

[pattacini]

@marcoaccame, while talking with @maggia80, we thought it might be worth plotting a graph containing both the expected and the actual outputs. This way, we could be 100% certain that we're not underestimating specific details.

Is this doable?

[marcoaccame]

yes it is. i assume that:

• expected outputs = the values you gave me in the excel file
• actual outputs = the values computed by C++ code given the input.

i can anticipate that 95% of them return true when compared wit operator ==, and the remaining are very similar.

i will soon publish the values / plots.

[marcoaccame]

hi, I have collected the data and produced the graphs. they are in the attached excel file.

results.xlsx

the output values from matlab are exactly equal to those computed by the c++ code run in the stm32f7 board. their difference is equal to 0.

i attach the plots anyway.

input: IV_P_fbk

output: IV_dc

output: IV_P_ref

output: signals[1]

output: signals[2]

[pattacini]

Thanks @marcoaccame super 👍

[pattacini]

Here are the latest simulation results w/ CP controlled by Series 130 evo:

Notably, the CP valve no longer needs to operate fully open at the onset.

cc @marcoreds81 @maggia80 @simeonedussoni

[simeonedussoni]

feeco thus we have margin.

[pattacini]

I've made some progress on the part of assisted ventilation.

Contrary to what was my original intuition, triggering upon S5 (outlet flow) is not doable as the flowmeter can sense the flow only in one direction and the respiratory request may occur when we're close to null flow. As a result, we're forced to proceed w/ S3 sensor (mask pressure).

As pressure sensors are particularly noisy, we need to hive off a dedicated filter from the pressure filtering bank working at 25 Hz in order to be able to tune up the cut off frequency much lower, i.e. at 2.5 Hz. This way, the resulting signal is smooth enough to be processed by the derivative filter running with a time constant of 0.1 s (required to be one decade lower than the sample rate).

Here's below you'll find the outcome where the autonomous respiratory requests of the patient are reported in purple:

Triggering is elicited at a threshold of -10 cmH20/s:

As it can be easily figured out, we don't have much margin to play with on the threshold to be safe against false positives. Consider that the triggering peaks are then cut out but they are basically only slightly more negative than -10 cmH2O/s.

[pattacini]

Interestingly, the inlet valve (IV) attempts to attain full closure as the pressure is constantly slightly higher than the reference. In this case, we've succeeded to keep the set-point, but once the IV is fully closed, our only chance to still modulate the pressure is operating the OV (in PWM as it's an On/Off valve), which is currently kept closed while delivering flow.

[pattacini]

Triggering upon the pressure measured by S1 (when it goes under a given threshold - here set to 3 cmH2O) instead of its derivative is definitely doable, as illustrated by the graph below:

Moreover, triggering appears to be robust also upon an asynchronous request occurring when we're far from the threshold.

[pattacini]

FMEA: PR1 breakage

PRVC

The system turns out to be robust upon breakage of PR1 regulator, although the controller's performances in subsequent events are quite poor:

To get around this, one possibility is to proceed with gain scheduling:

VCV

The situation in Volume Control mode is pretty much similar, barring that the controller's gains need to be scheduled at once:

[simeonedussoni]

@pattacini great! I also noticed that you reset the output volume calculation upon the beginning of a new respiratory cycle. this looks good but I would like

1. to measure this volume in the applied case to see if we can do better
2. to provide an alarm if the calculated volume falls too far from zero, in order to have a diagnostic. this can be made in the firmware outside of the control loop, right?

[pattacini]

@simeonedussoni I'm still not focusing on tidal volume measurements and what you see from time to time are mere random tryouts. That said, I'm afraid the problem is intrinsically difficult to tackle. Nonetheless, yes we can certainly disable resetting while hoping that in real conditions those measurements will make more sense to us.

Regarding alarms, I can confirm that they can be profitably and more conveniently handled outside the controller.

[pattacini]

Let's close this issue as it's very likely that we'll be adapting the present structures in order to let it fit the real settings. No time for further developments now.