Python controller scripts for Windkessel resistance tuning¶
Python control mechanism¶
The python controllers described in this document can be used to tune 3-element Windkessel model resistances dynamically during simulation runtime, thereby avoiding the tuning iterations typically required with manual parameter tuning. The controller adjusts the distal resistance of the Windkessel model to match the mean flow rate at the model outlet to a predefined mean target flow. A flow chart of the control mechanism is illustrated in the figure below
At the beginning of the simulation, the proximal and distal resistances are set to the initial values \(R_{p,i}^0\) and \(R_{d,i}^0\), where index \(i\) indicates the outlet. Initial resistance values can be estimated using the 1D theory. Since the controllers only allow for matching of the mean flow rate (cardiac-averaged), the controllers are activated after the simulation has passed a full cardiac cycle, \(T\) is the number of time steps per cycle. For each following simulation time step \(n > T\), a mean flow rate is calculated as \(\bar{Q}_i^n=∑_{k=n-T}^nQ_i^k /T\). The resistance is then changed proportionally to the difference between \(\bar{Q}_i^n\) and the mean target flow rate \(\bar{Q}_{target,i}\). Note that the controllers only adjust the distal resistance due to its high impact on the flow rate, whereas the proximal resistance is kept constant throughout the simulation. If desired, the ratio between the proximal and distal resistance can be adjusted after the controllers have converged. As long as the total resistance \(R_T=R_p+R_d\) remains constant, the mean flow rate is maintained.
Setting up python controllers¶
The following steps need to be repeated for each individual outlet that requires tuning
Create controller.py file¶
Make a copy of the [sample python controller] (https://www.dropbox.com/s/xvwta1ulkdw84kn/controller_sample.py?dl=0) file and make the following changes to it:
- Define the mean target flow rate
self.targetFlow
you want to match at the outlet - Define the initial estimates for the proximal and distal resistances
self.R_p
andself.R_d
- Define the number of timesteps per cardiac cycle
self.timeStepsCycle
- The values of distal resistance
R_d
and the error between the mean flow rate and the mean target flow ratehunger
are written out for each time step. Make sure to change theOUTLETNAME
in the filename
The controller script defines maximum limits for the hunger self.maximumAllowedNegativeHunger
and self.maximumAllowedPositiveHunger
to avoid extreme overshoot in the controller response in scenarios of large flow differences. These limits can be changed if the controller response is slow. However, this needs to be done carefully, since it may lead to an unstable controller behavior.
Create observer.py file¶
Make a copy of the [sample python observer] (https://www.dropbox.com/s/235h6occs8wimh4/observer_sample.py?dl=0) file and make the following changes to it:
- Define
OUTLETNAME
Set up Windkessel model in CRIMSON using netlist¶
In order to use the controller scripts, the Windkessel model needs to be defined in CRIMSON using the netlist module. Create a new BC and select Netlist. Then click on launch Editor
and create your Windkessel model as shown below
Under properties, define the initial Windkessel parameters and add the python controller scripts as following:
Observers.py: The observer is added under pressure P1
and “observes” the flow rate at the 3D model outlet. To add the script, change No Dynamic Control
to Custom Python Script
and specify the name of the observer.py name in the youPythonScriptNameHere.py
field.
Controller.py: The controller is added under the resistance R2
(distal resistance). Again, to add the script, change No Dynamic Control
to Custom Python Script
and specify the name of the controller.py name in the youPythonScriptNameHere.py
field.
Save the created Windkessel model in both .xml and .dat format. The .xml file can be opened in the editor to make changes later on and the .dat file will be used in the BC setup.
Next, the created files need to be loaded into the CRIMSON BC setup. Under circuit
load the Windkessel model in .dat format. Under Adjustment scripts
load both the observer.py and controller.py files
Running the python scripts¶
After writing out the simulation files, you will find a netlist_surfaces.dat
and a copy of the controller and observer scripts together with the rest of the simulation files. All these created files are necessary to run the simulation. After the simulation is started, the netlist_surfaces.dat
is converted to a netlist_surfaces.xml
. The .xml file can later be opened in any editor to make changes to the boundary conditions, avoiding the need to go back to the CRIMSON GUI.
The simulation should be run for at least 10 cardiac cycles to allow convergence of the controllers. In cases of large flow differences between the 3D CFD and the target values, more cycles might be necessary. To check if the controllers have converged, plot the R_d_history
. Additionally, the cardiac-averaged flow, obtained from FlowHist.dat
, and the mean target flow should be compared.