Build your own model class
Todo
Explain how to create a custom model class
When the model is created, it must be parameterized properly. It is explained here how to create your own reference model state.
The model VanOsta2023 is build following this standard. This model is constructed by the following code.
1"""
2VanOsta2023 model.
3
4"""
5
6from circadapt.model import Model
7from circadapt.adapt import ModelAdapt
8from circadapt.plot import triseg2022
9import matplotlib.pyplot as plt
10import numpy as np
11
12
13class VanOsta2023(Model, ModelAdapt):
14 def __init__(self,
15 solver: str = None,
16 path_to_circadapt: str = None,
17 model_state: dict = None,
18 ):
19 if solver is None:
20 solver = 'backward_differential'
21
22 self._local_save_reference = True
23
24 ModelAdapt.__init__(self)
25 Model.__init__(self,
26 solver,
27 path_to_circadapt=path_to_circadapt,
28 model_state=model_state,
29 )
30
31 def _update_settings(self):
32 # rename model components for easy input/output
33 self._module_rename['Wall2022'] = 'Wall'
34 self._module_rename['Patch2022'] = 'Patch'
35 self._module_rename['Chamber2022'] = 'Chamber'
36 self._module_rename['Valve2022'] = 'Valve'
37 self._module_rename['TriSeg2022'] = 'TriSeg'
38 self._module_rename['PressureFlowControl'] = 'PFC'
39
40 def build(self):
41 # Circulation
42 self.add_smart_component('ArtVen', build='SystemicCirculation')
43 self.add_smart_component('ArtVen', build='PulmonaryCirculation')
44 self.add_smart_component('Heart', patch_type='Patch2022')
45 self.add_smart_component('Timings')
46 self.add_smart_component('PressureFlowControl')
47
48 # manually set papillary muscles
49 self.set_component("Peri.RaRv.wPapMus", "Peri.TriSeg.wRv")
50 self.set_component("Peri.LaLv.wPapMus", "Peri.TriSeg.wLv")
51
52 def set_reference(self):
53 self['Chamber']['buckling'] = False
54 self['Valve']['soft_closure'] = True
55 self['Valve']['papillary_muscles'] = True
56
57 self['Solver']['dt'] = 0.001
58 self['Solver']['dt_export'] = 0.002
59 self.set('Solver.order', 2)
60
61 self.set('Model.t_cycle', 0.85)
62
63 self.set('Model.PFC.fac', 0.5)
64 self.set('Model.PFC.epsilon', 0.1)
65 self.set('Model.PFC.fac_pfc', 1)
66 self.set('Model.PFC.stable_threshold', 0.001)
67
68 # Set parameters Circulation
69 self['Tube0D']['A0'] = [0.00049828, 0.00049899, 0.00046929, 0.00050574]
70 self['Tube0D']['A_wall'] = [1.13187792e-04, 3.72289316e-05, 7.27489286e-05, 4.31617949e-05]
71 self['Tube0D']['k'] = [1.66666667, 2.33333333, 1.66666667, 2.33333333]
72 self['Tube0D']['l'] = [0.4, 0.4, 0.19, 0.19]
73 self['Tube0D']['p0'] = [12155.25602777, 297.44066943, 1912.6684696 , 580.51425824]
74
75 self['ArtVen']['p0'] = [6299.86481149, 1000. ]
76 self['ArtVen']['q0'] = [4.5e-05, 4.5e-05]
77 self['ArtVen']['k'] = [1. , 1.88]
78
79 # Set volume state variables
80 # not needed to reset, only if initial model crashes
81 if True:
82 self.set('Model.SyArt.V', 2e-4)
83 self.set('Model.PuArt.V', 1e-4)
84 self.set('Model.SyVen.V', 4e-4)
85 self.set('Model.PuVen.V', 2e-4)
86
87 self.set('Model.Peri.TriSeg.V', 56e-6)
88 self.set('Model.Peri.TriSeg.Y', 36e-3)
89
90 self.set('Model.Peri.TriSeg.cLv.V', 150e-6)
91 self.set('Model.Peri.TriSeg.cRv.V', 110e-6)
92 self.set('Model.Peri.La.V', 25e-6)
93 self.set('Model.Peri.Ra.V', 150e-6)
94
95 # all wall parameters
96 self['Patch']['l_se0'] = 0.04
97 self['Patch']['l_s_ref'] = 2.0
98 self['Patch']['l_s0'] = [1.79, 1.79, 1.8 , 1.8 , 1.8 ]
99 self['Patch']['dl_s_pas'] = 0.6
100 self['Patch']['k1'] = 10
101 self['Patch']['dt'] = 0
102 self['Patch']['C_rest'] = 0
103 self['Patch']['l_si0'] = 1.51
104 self['Patch']['LDAD'] = [1.057, 1.057, 0.74 , 0.74 , 0.74 ]
105 self['Patch']['ADO'] = [0.081, 0.081, 0.75 , 0.75 , 0.75 ]
106 self['Patch']['LDCC'] = [4. , 4. , 3., 3., 3.]
107 self['Patch']['v_max'] = 7.
108 self['Patch']['v_max'][:2] = 14.
109
110 self['Patch']['tr'] = [0.4 , 0.4 , 0.21, 0.21, 0.21]
111 self['Patch']['td'] = [0.4 , 0.4 , 0.24, 0.24, 0.24]
112
113 # wall specific
114 self['Patch']['Sf_act'] = [100000., 100000., 120000., 120000., 120000.]
115 self['Patch']['Am_ref'] = [0.00400721, 0.00311263, 0.00918281, 0.00518683, 0.01321999]
116 self['Patch']['V_wall'] = [1.59935455e-05, 7.72140495e-06, 1.16441768e-04, 4.51241560e-05,
117 6.88196851e-05]
118 self['Patch']['Sf_pas'] = [ 9.51064014, 10.3647024 , 424.3537027 , 417.74300178,
119 439.21574206]
120
121 # Valves
122 self['Valve']['adaptation_A_open_fac'] = [1. , 1.12, 1. , 1. , 1.12, 1. ]
123 self['Valve']['A_open'] = [0.00050044, 0.00052598, 0.00046962, 0.00050745, 0.0005583 ,
124 0.00049849]
125 self['Valve']['A_leak'] = [2.64705882e-04, 2.64705882e-10, 2.64705882e-10, 2.64705882e-04,
126 2.64705882e-10, 2.64705882e-10]
127 self['Valve']['l'] = 0.01626978
128 self['Valve']['rho_b'] = 1050
129 self['Valve']['papillary_muscles'] = True
130 self['Valve']['papillary_muscles_slope'] = 100
131 self['Valve']['papillary_muscles_min'] = 0.1
132 self['Valve']['papillary_muscles_A_open_fac'] = 0.1
133 self['Valve']['soft_closure'] = True
134
135 self['Bag']['k'] = 10
136 self['Bag']['V_ref'] = [0.00057937]
137 self['Bag']['p_ref'] = 100
138
139 self.set('Model.Peri.TriSeg.Y', 0.035)
140
141 self['Timings']['law_tau_av'] = 2
142 self['Timings']['c_tau_av0'] = 0.189
143 self['Timings']['c_tau_av1'] = -0.029
144
145 # Adaptation
146 self['Patch']['adapt_gamma'] = 0.5
147 self.set('Solver.store_beats', 1)
148
149 self['PFC']['stable_threshold'] = 1e-3
150
151 self['Patch']['SfPasMaxT'] = [3600., 3600., 4200., 4200., 4200.]
152 self['Patch']['FacSfActT'] = [0.28, 0.28, 0.69, 0.69, 0.69]
153 self['Patch']['SfPasActT'] = [2800., 2800., 6600., 6600., 6600.]
154 self['Patch']['LsPasActT'] = [3. , 3. , 2.31, 2.31, 2.31]
155
156 self.run(stable=True)
157
158 return
159
160
161 def get_unittest_targets(self):
162 """Hardcoded results after initializing and running 1 beat."""
163 return {
164 'LVEDV': 130.9,
165 'LVESV': 58.5,
166 }
167
168 def get_unittest_results(self, model):
169 """Real-time results after initializing and running 1 beat."""
170 LVEDV = np.max(model['Cavity']['V'][:, 'cLv'])*1e6
171 LVESV = np.min(model['Cavity']['V'][:, 'cLv'])*1e6
172 return {
173 'LVEDV': LVEDV,
174 'LVESV': LVESV,
175 }
176
177 def plot(self, fig=None):
178 # TODO
179 self.plot_extended(fig)
180
181 def plot_extended(self, fig=None):
182 if fig is None:
183 fig = 1
184 if isinstance(fig, int):
185 fig = plt.figure(fig, clear=True, figsize=(12, 8))
186
187 # Settings
188 grid_size = [32, 32]
189
190 def get_lim(module, signal, locs=slice(None, None, None)):
191 signal = self[module][signal][:, locs]
192 lim = np.array([np.min(signal), np.max(signal)])
193 lim += np.array([-1, 1]) * 0.1*np.diff(lim)
194 return lim
195
196 lim_V = get_lim('Cavity', 'V', ['cRv', 'Ra', 'La', 'cLv']) * 1e6
197 lim_V[0] = 0
198 lim_p = get_lim('Cavity', 'p', ['cRv', 'Ra', 'La', 'cLv']) / 133
199 lim_p[0] = np.min([lim_p[0], 0])
200 lim_Ls = get_lim('Patch', 'l_s')
201 lim_Sf = get_lim('Patch', 'Sf') * 1e-3
202 lim_q = get_lim('Valve', 'q',
203 ['LaLv', 'RaRv', 'LvSyArt', 'RaPuArt']) * 1e6
204
205 all_lim = [lim_V, lim_p, lim_Ls, lim_Sf, lim_q]
206 if (np.any(np.isnan(all_lim)) or np.any(np.isinf(all_lim))):
207 lim_V = [0, 200]
208 lim_p = [0, 150]
209 lim_Ls = [1.5, 2.0]
210 lim_Sf = [0, 100]
211 lim_q = [-1e-3, 1e-3]
212
213 # Pressure Volume plot
214 axPV = plt.subplot2grid(grid_size, (0, 17), rowspan=15, colspan=15, fig=fig)
215 axPV.plot(self['Cavity']['V'][:, 'cLv']*1e6, self['Cavity']['p'][:, 'cLv']/133)
216 axPV.plot(self['Cavity']['V'][:, 'cRv']*1e6, self['Cavity']['p'][:, 'cRv']/133)
217 axPV.plot(self['Cavity']['V'][:, 'La']*1e6, self['Cavity']['p'][:, 'La']/133)
218 axPV.plot(self['Cavity']['V'][:, 'Ra']*1e6, self['Cavity']['p'][:, 'Ra']/133)
219 axPV.spines[['top', 'right']].set_visible(False)
220 axPV.set_title('Pressure-Volume loop', weight='bold')
221 axPV.set_xlabel('Volume [mL]')
222 axPV.set_ylabel('Pressure [mmHg]')
223 axPV.spines[['bottom', 'left']].set_position(('outward', 5))
224
225 ylabel_x_left = -0.25
226 ylabel_x_right = 1.25
227
228 # Volumes
229 t = self['Solver']['t']*1e3
230
231 axVRv = plt.subplot2grid(grid_size, (0, 0), rowspan=8, colspan=6, fig=fig)
232 axVRv.plot(t, self['Cavity']['V'][:, 'cRv']*1e6)
233 axVRv.plot(t, self['Cavity']['V'][:, 'Ra']*1e6)
234 axVRv.set_ylim(lim_V)
235 axVRv.set_ylabel('Volume\n[mL]')
236 axVRv.spines[['top', 'right']].set_visible(False)
237 axVRv.set_title('Right Heart')
238 # axVRv.set_xticks([])
239 axVRv.tick_params(axis='both', direction='in')
240 axVRv.yaxis.set_label_coords(ylabel_x_left, 0.5)
241
242
243 axVLv = plt.subplot2grid(grid_size, (0, 6), rowspan=8, colspan=6, fig=fig)
244 axVLv.plot(t, self['Cavity']['V'][:, 'cLv']*1e6)
245 axVLv.plot(t, self['Cavity']['V'][:, 'La']*1e6)
246 axVLv.set_ylabel('Volume\n[mL]')
247 axVLv.set_ylim(lim_V)
248 axVLv.yaxis.set_ticks_position('right')
249 axVLv.yaxis.set_label_position('right')
250 axVLv.spines['right'].set_position(('outward', 0))
251 axVLv.spines[['top', 'left']].set_visible(False)
252 axVLv.set_title('Left Heart')
253 # axVLv.set_xticks([])
254 axVLv.tick_params(axis='both', direction='in')
255 axVLv.yaxis.set_label_coords(ylabel_x_right, 0.5)
256
257 # Pressures
258 axpRv = plt.subplot2grid(grid_size, (8, 0), rowspan=8, colspan=6, fig=fig)
259 axpRv.plot(t, self['Cavity']['p'][:, 'cRv']/133)
260 axpRv.plot(t, self['Cavity']['p'][:, 'Ra']/133)
261 axpRv.plot(t, self['Cavity']['p'][:, 'PuArt']/133)
262 axpRv.spines[['top', 'right']].set_visible(False)
263 # axpRv.set_xticks([])
264 axpRv.tick_params(axis='both', direction='in')
265 axpRv.set_ylim(lim_p)
266 axpRv.set_ylabel('Pressure\n[mmHg]')
267 axpRv.yaxis.set_label_coords(ylabel_x_left, 0.5)
268
269 axpLv = plt.subplot2grid(grid_size, (8, 6), rowspan=8, colspan=6, fig=fig)
270 axpLv.plot(t, self['Cavity']['p'][:, 'cLv']/133)
271 axpLv.plot(t, self['Cavity']['p'][:, 'La']/133)
272 axpLv.plot(t, self['Cavity']['p'][:, 'SyArt']/133)
273 axpLv.yaxis.set_ticks_position('right')
274 axpLv.yaxis.set_label_position('right')
275 axpLv.spines['right'].set_position(('outward', 0))
276 axpLv.spines[['top', 'left']].set_visible(False)
277 # axpLv.set_xticks([])
278 axpLv.tick_params(axis='both', direction='in')
279 axpLv.set_ylim(lim_p)
280 axpLv.set_ylabel('Pressure\n[mmHg]')
281 axpLv.yaxis.set_label_coords(ylabel_x_right, 0.5)
282
283 # Valves
284 ax = plt.subplot2grid(grid_size, (16, 0), rowspan=6, colspan=6, fig=fig)
285 ax.plot(t, self['Valve']['q'][:, 'RaRv']*1e6)
286 ax.plot(t, self['Valve']['q'][:, 'RvPuArt']*1e6)
287 ax.spines[['top', 'right']].set_visible(False)
288 ax.set_ylim(lim_q)
289 # ax.set_xticks([])
290 ax.set_ylabel('Flow\n[mL/s]')
291 ax.yaxis.set_label_coords(ylabel_x_left, 0.5)
292
293 ax = plt.subplot2grid(grid_size, (16, 6), rowspan=6, colspan=6, fig=fig)
294 ax.plot(t, self['Valve']['q'][:, 'LaLv']*1e6)
295 ax.plot(t, self['Valve']['q'][:, 'LvSyArt']*1e6)
296 ax.spines[['top', 'left']].set_visible(False)
297 ax.set_ylim(lim_q)
298 ax.yaxis.set_ticks_position('right')
299 ax.yaxis.set_label_position('right')
300 # ax.set_xticks([])
301 ax.set_ylabel('Flow\n[mL/s]')
302 ax.yaxis.set_label_coords(ylabel_x_right, 0.5)
303
304 # Stress
305 ax = plt.subplot2grid(grid_size, (22, 0), rowspan=4, colspan=6, fig=fig)
306 ax.plot(t, self['Patch']['Sf'][:, 'pRv0']*1e-3)
307 ax.plot(t, self['Patch']['Sf'][:, 'pRa0']*1e-3)
308 ax.spines[['top', 'right']].set_visible(False)
309 # ax.set_xticks([])
310 ax.set_ylim(lim_Sf)
311 ax.set_ylabel('Total\nstress [kPa]')
312 ax.yaxis.set_label_coords(ylabel_x_left, 0.5)
313
314 ax = plt.subplot2grid(grid_size, (22, 6), rowspan=4, colspan=6, fig=fig)
315 ax.plot(t, self['Patch']['Sf'][:, 'pLv0']*1e-3)
316 ax.plot(t, self['Patch']['Sf'][:, 'pSv0']*1e-3)
317 ax.plot(t, self['Patch']['Sf'][:, 'pLa0']*1e-3)
318 ax.spines[['top', 'left']].set_visible(False)
319 ax.yaxis.set_ticks_position('right')
320 ax.yaxis.set_label_position('right')
321 # ax.set_xticks([])
322 ax.set_ylim(lim_Sf)
323 ax.set_ylabel('Total\nstress [kPa]')
324 ax.yaxis.set_label_coords(ylabel_x_right, 0.5)
325
326 # Sarcomere Length
327 ax = plt.subplot2grid(grid_size, (26, 0), rowspan=6, colspan=6, fig=fig)
328 ax.plot(t, self['Patch']['l_s'][:, 'pRv0'])
329 ax.plot(t, self['Patch']['l_s'][:, 'pRa0'])
330 ax.spines[['top', 'right']].set_visible(False)
331 ax.spines['bottom'].set_position(('outward', 5))
332 ax.set_ylim(lim_Ls)
333 ax.set_ylabel('Sarcomere\nlength [$\mu$m]')
334 ax.yaxis.set_label_coords(ylabel_x_left, 0.5)
335
336 ax = plt.subplot2grid(grid_size, (26, 6), rowspan=6, colspan=6, fig=fig)
337 ax.plot(t, self['Patch']['l_s'][:, 'pLv0'])
338 ax.plot(t, self['Patch']['l_s'][:, 'pSv0'])
339 ax.plot(t, self['Patch']['l_s'][:, 'pLa0'])
340 ax.spines[['top', 'left']].set_visible(False)
341 ax.spines['bottom'].set_position(('outward', 5))
342 ax.yaxis.set_ticks_position('right')
343 ax.yaxis.set_label_position('right')
344 ax.set_ylim(lim_Ls)
345 ax.set_ylabel('Sarcomere\nlength [$\mu$m]')
346 ax.yaxis.set_label_coords(ylabel_x_right, 0.5)
347
348 # ax.set_xlabel('Time [ms]')
349 # ax.xaxis.set_label_coords(0, -0.3)
350
351
352 # Plot TriSeg
353 titles = ['Pre-A', 'Onset QRS', 'Peak LV \n pressure', 'AV close']
354 idx = [0,
355 np.argmax(np.diff(self['Patch']['C'][:, 'pLv0'])>0),
356 np.argmax(self['Cavity']['p'][:, 'cLv']),
357 len(t) - 1 - np.argmax(
358 np.diff(self['Valve']['q'][:, 'LvSyArt'][::-1])>0)
359 ]
360
361 for i in range(4):
362 ax = plt.subplot2grid(grid_size, (26, 16+4*i),
363 rowspan=5, colspan=4, fig=fig)
364 triseg2022(self, ax, idx[i])
365 ax.spines[['top', 'right', 'bottom', 'left']].set_visible(False)
366 ax.set_xticks([])
367 ax.set_yticks([])
368 plt.xlabel(titles[i], fontsize=12)
369
370 # Plot settings
371 plt.subplots_adjust(
372 top=0.96,
373 bottom=0.05,
374 left=0.075,
375 right=0.98,
376 hspace=5,
377 wspace=0.5)
378 plt.draw()
379
380 # Stress strain
381 ax_right_stress_strain = plt.subplot2grid(grid_size, (18, 17),
382 rowspan=7, colspan=7, fig=fig)
383 ax_left_stress_strain = plt.subplot2grid(grid_size, (18, 24),
384 rowspan=7, colspan=7, fig=fig)
385
386 ax_right_stress_strain.plot(
387 self['Patch']['Ef'][:, [1, 4]],
388 self['Patch']['Sf'][:, [1, 4]]*1e-3,
389 )
390 ax_left_stress_strain.plot(
391 self['Patch']['Ef'][:, [0, 2, 3]],
392 self['Patch']['Sf'][:, [0, 2, 3]]*1e-3,
393 )
394
395 ylim = [np.min(self['Patch']['Sf'])*1e-3,
396 np.max(self['Patch']['Sf'])*1e-3]
397 ylim_ptp = np.ptp(ylim)*0.025
398 ylim[0] -= ylim_ptp
399 ylim[1] += ylim_ptp
400 ax_right_stress_strain.set_ylim(ylim)
401 ax_left_stress_strain.set_ylim(ylim)
402 xlim = [np.min(self['Patch']['Ef']),
403 np.max(self['Patch']['Ef'])]
404 xlim_ptp = np.ptp(xlim)*0.025
405 xlim[0] -= xlim_ptp
406 xlim[1] += xlim_ptp
407 ax_right_stress_strain.set_xlim(xlim)
408 ax_left_stress_strain.set_xlim(xlim)
409
410 ax_right_stress_strain.spines[['left', 'bottom']].set_position(('outward', 5))
411 ax_left_stress_strain.spines[['right', 'bottom']].set_position(('outward', 5))
412
413 ax_right_stress_strain.spines[['top', 'right']].set_visible(False)
414 ax_right_stress_strain.set_ylabel('[kPa]', rotation=0, va='bottom', ha='right')
415 ax_right_stress_strain.yaxis.set_label_coords(0., 1.05)
416 ax_right_stress_strain.set_xlabel('Strain [-]')
417 ax_left_stress_strain.spines[['top', 'left']].set_visible(False)
418 ax_left_stress_strain.yaxis.set_ticks_position('right')
419 ax_left_stress_strain.yaxis.set_label_position('right')
420 ax_left_stress_strain.set_ylabel('Stress\n[kPa]', rotation=0, va='bottom', ha='left')
421 ax_left_stress_strain.yaxis.set_label_coords(1., 1.05)
422 ax_left_stress_strain.set_xlabel('Strain [-]')
423
424
425class VanOsta2023_2310(VanOsta2023):
426 def set_reference(self):
427 self['Chamber']['buckling'] = False
428 self['Valve']['soft_closure'] = True
429 self['Valve']['papillary_muscles'] = True
430
431 self['Solver']['dt'] = 0.001
432 self['Solver']['dt_export'] = 0.002
433 self.set('Solver.order', 2)
434
435 self.set('Model.t_cycle', 0.85)
436
437 self.set('Model.PFC.fac', 0.5)
438 self.set('Model.PFC.epsilon', 0.1)
439 self.set('Model.PFC.fac_pfc', 1)
440 self.set('Model.PFC.stable_threshold', 0.001)
441
442 # Set parameters Circulation
443 self['Tube0D']['A0'] = [0.00049828, 0.00049899, 0.00046929, 0.00050574]
444 self['Tube0D']['A_wall'] = [1.13187792e-04, 3.72289316e-05, 7.27489286e-05, 4.31617949e-05]
445 self['Tube0D']['k'] = [1.66666667, 2.33333333, 1.66666667, 2.33333333]
446 self['Tube0D']['l'] = [0.4, 0.4, 0.19, 0.19]
447 self['Tube0D']['p0'] = [12155.25602777, 297.44066943, 1912.6684696 , 580.51425824]
448
449 self['ArtVen']['p0'] = [6299.86481149, 1000. ]
450 self['ArtVen']['q0'] = [4.5e-05, 4.5e-05]
451 self['ArtVen']['k'] = [1. , 1.88]
452
453 # Set volume state variables
454 # not needed to reset, only if initial model crashes
455 if True:
456 self.set('Model.SyArt.V', 2e-4)
457 self.set('Model.PuArt.V', 1e-4)
458 self.set('Model.SyVen.V', 4e-4)
459 self.set('Model.PuVen.V', 2e-4)
460
461 self.set('Model.Peri.TriSeg.V', 56e-6)
462 self.set('Model.Peri.TriSeg.Y', 36e-3)
463
464 self.set('Model.Peri.TriSeg.cLv.V', 150e-6)
465 self.set('Model.Peri.TriSeg.cRv.V', 110e-6)
466 self.set('Model.Peri.La.V', 25e-6)
467 self.set('Model.Peri.Ra.V', 150e-6)
468
469 # all wall parameters
470 self['Patch']['l_se0'] = 0.04
471 self['Patch']['l_s_ref'] = 2.0
472 self['Patch']['l_s0'] = [1.79, 1.79, 1.8 , 1.8 , 1.8 ]
473 self['Patch']['dl_s_pas'] = 0.6
474 self['Patch']['k1'] = 10
475 self['Patch']['dt'] = 0
476 self['Patch']['C_rest'] = 0
477 self['Patch']['l_si0'] = 1.51
478 self['Patch']['LDAD'] = [1.057, 1.057, 0.74 , 0.74 , 0.74 ]
479 self['Patch']['ADO'] = [0.081, 0.081, 0.75 , 0.75 , 0.75 ]
480 self['Patch']['LDCC'] = [4. , 4. , 3., 3., 3.]
481 self['Patch']['v_max'] = 7.
482 self['Patch']['v_max'][:2] = 14.
483
484 self['Patch']['tr'] = [0.4 , 0.4 , 0.21, 0.21, 0.21]
485 self['Patch']['td'] = [0.4 , 0.4 , 0.24, 0.24, 0.24]
486
487 # wall specific
488 self['Patch']['Sf_act'] = [100000., 100000., 120000., 120000., 120000.]
489 self['Patch']['Am_ref'] = [0.00400721, 0.00311263, 0.00918281, 0.00518683, 0.01321999]
490 self['Patch']['V_wall'] = [1.59935455e-05, 7.72140495e-06, 1.16441768e-04, 4.51241560e-05,
491 6.88196851e-05]
492 self['Patch']['Sf_pas'] = [ 9.51064014, 10.3647024 , 424.3537027 , 417.74300178,
493 439.21574206]
494
495 # Valves
496 self['Valve']['adaptation_A_open_fac'] = [1. , 1.12, 1. , 1. , 1.12, 1. ]
497 self['Valve']['A_open'] = [0.00050044, 0.00052598, 0.00046962, 0.00050745, 0.0005583 ,
498 0.00049849]
499 self['Valve']['A_leak'] = [2.64705882e-04, 2.64705882e-10, 2.64705882e-10, 2.64705882e-04,
500 2.64705882e-10, 2.64705882e-10]
501 self['Valve']['l'] = 0.01626978
502 self['Valve']['rho_b'] = 1050
503 self['Valve']['papillary_muscles'] = True
504 self['Valve']['papillary_muscles_slope'] = 100
505 self['Valve']['papillary_muscles_min'] = 0.1
506 self['Valve']['papillary_muscles_A_open_fac'] = 0.1
507 self['Valve']['soft_closure'] = True
508
509 self['Bag']['k'] = 10
510 self['Bag']['V_ref'] = [0.00057937]
511 self['Bag']['p_ref'] = 100
512
513 self.set('Model.Peri.TriSeg.Y', 0.035)
514
515 self['Timings']['law_tau_av'] = 2
516 self['Timings']['c_tau_av0'] = 0.22
517 self['Timings']['c_tau_av1'] = -0.03
518
519 # Adaptation
520 self['Patch']['adapt_gamma'] = 0.5
521 self.set('Solver.store_beats', 1)
522
523 self['PFC']['stable_threshold'] = 1e-3
524
525 self['Patch']['SfPasMaxT'] = [3600., 3600., 4200., 4200., 4200.]
526 self['Patch']['FacSfActT'] = [0.28, 0.28, 0.69, 0.69, 0.69]
527 self['Patch']['SfPasActT'] = [2800., 2800., 6600., 6600., 6600.]
528 self['Patch']['LsPasActT'] = [3. , 3. , 2.31, 2.31, 2.31]
529
530 self.run(stable=True)
531
532 return
533
534 def get_unittest_targets(self):
535 """Hardcoded results after initializing and running 1 beat."""
536 return {
537 'LVEDV': 121.75,
538 'LVESV': 48.0,
539 }
540
541
542class VanOsta2023_2306(VanOsta2023):
543 def set_reference(self):
544 self['Chamber']['buckling'] = False
545 self['Valve']['soft_closure'] = True
546
547 self['Solver']['dt'] = 0.0005
548 self['Solver']['dt_export'] = self['Solver']['dt']
549
550 self.set('Model.t_cycle', 0.85)
551
552 self.set('Model.PFC.fac', 0.5)
553 self.set('Model.PFC.epsilon', 0.1)
554 self.set('Model.PFC.fac_pfc', 1)
555 self.set('Model.PFC.stable_threshold', 0.001)
556
557 # Set parameters Circulation
558 self['Tube0D']['A0'] = [
559 0.00049833,
560 0.00050034,
561 0.0004684 ,
562 0.00051849,
563 ]
564 self['Tube0D']['A_wall'] = [
565 1.14408188e-04,
566 4.81082964e-05,
567 5.19472921e-05,
568 5.79933632e-05,
569 ]
570 self['Tube0D']['k'] = np.array([ 8., 10., 8., 10.]) / 3 - 1
571 self['Tube0D']['l'] = [0.4, 0.4, 0.2, 0.2]
572 self['Tube0D']['p0'] = [
573 12161.21763497,
574 112.09788908,
575 1899.95009955,
576 572.39790382,
577 ]
578
579 self['ArtVen']['p0'] = [6395.43973044, 500.]
580 self['ArtVen']['q0'] = [4.5e-05, 4.5e-05]
581 self['ArtVen']['k'] = [1, 2]
582
583 # Set volume state variables
584 # not needed to reset, only if initial model crashes
585 if True:
586 self.set('Model.SyArt.V', 2e-4)
587 self.set('Model.PuArt.V', 1e-4)
588 self.set('Model.SyVen.V', 4e-4)
589 self.set('Model.PuVen.V', 2e-4)
590
591 self.set('Model.Peri.TriSeg.V', 44e-6)
592 self.set('Model.Peri.TriSeg.Y', 34.6e-3)
593
594 self.set('Model.Peri.TriSeg.cLv.V', 150e-6)
595 self.set('Model.Peri.TriSeg.cRv.V', 100e-6)
596 self.set('Model.Peri.La.V', 100e-6)
597 self.set('Model.Peri.Ra.V', 50e-6)
598
599 # all wall parameters
600 self['Patch']['l_se'] = 0.04
601 self['Patch']['l_s_ref'] = 2.0
602 self['Patch']['l_s0'] = 1.8
603 self['Patch']['dl_s_pas'] = 0.6
604 self['Patch']['k1'] = 10
605 self['Patch']['dt'] = 0
606 self['Patch']['C_rest'] = 0
607 self['Patch']['l_si0'] = 1.51
608 self['Patch']['LDAD'] = 1.057
609 self['Patch']['ADO'] = 0.65
610 self['Patch']['LDCC'] = 4.
611 self['Patch']['v_max'] = 7.
612 self['Patch']['v_max'][:2] = 14.
613
614 self['Patch']['tr'] = 0.25
615 self['Patch']['td'] = 0.25
616 self['Patch']['tr'][:2] = 0.4
617 self['Patch']['td'][:2] = 0.4
618
619 # wall specific
620 self['Patch']['Am_ref'] = [
621 0.00706159,
622 0.00603949,
623 0.0093191 ,
624 0.00507053,
625 0.01313792,
626 ]
627 self['Patch']['V_wall'] = [
628 1.86518494e-05,
629 6.12511827e-06,
630 9.51121820e-05,
631 3.74835076e-05,
632 4.97478018e-05,
633 ]
634 self['Patch']['Sf_pas'] = np.array([
635 51591.64007141,
636 52277.48834349,
637 22370.46712926,
638 22051.80366169,
639 23278.22158555,
640 ]) * 0.0349
641 self['Patch']['Sf_act'] = [84e3, 84e3, 120e3, 120e3, 120e3]
642
643 # Valves
644 self['Valve']['adaptation_A_open_fac'] = [
645 1.,
646 1.5,
647 1.,
648 1.,
649 1.5,
650 1.]
651 self['Valve']['A_open'] = [
652 0.00050078,
653 0.0007021,
654 0.00046806,
655 0.00051776,
656 0.00074755,
657 0.00049837]
658 self['Valve']['A_leak'] = [
659 2.64705882e-04,
660 2.64705882e-10,
661 2.64705882e-10,
662 2.64705882e-04,
663 2.64705882e-10,
664 2.64705882e-10]
665 self['Valve']['l'] = 0.01626978
666 self['Valve']['rho_b'] = 1050
667 self['Valve']['papillary_muscles'] = False
668 self['Valve']['papillary_muscles_slope'] = 100
669 self['Valve']['papillary_muscles_min'] = 0.1
670 self['Valve']['papillary_muscles_A_open_fac'] = 0.1
671 self['Valve']['soft_closure'] = True
672
673 self['Bag']['k'] = 10
674 self['Bag']['V_ref'] = 0.0005
675 self['Bag']['p_ref'] = 100
676
677 self.set('Model.Peri.TriSeg.Y', 0.035)
678
679 self['Timings']['law_tau_av'] = 2
680 self['Timings']['c_tau_av0'] = 0.2
681 self['Timings']['c_tau_av1'] = -0.485*60e-3
682
683 # Adaptation
684 self['Patch']['adapt_gamma'] = 0.75
685
686 self.run(stable=True)
687
688 def get_unittest_targets(self):
689 """Hardcoded results after initializing and running 1 beat."""
690 return {
691 'LVEDV': 141.6,
692 'LVESV': 69.4,
693 }
694
695
696class VanOsta2023_2308(VanOsta2023):
697 def set_reference(self):
698 self['Chamber']['buckling'] = False
699 self['Valve']['soft_closure'] = True
700 self['Valve']['papillary_muscles'] = True
701
702 self['Solver']['dt'] = 0.001
703 self['Solver']['dt_export'] = 0.002
704 self.set('Solver.order', 2)
705
706 self.set('Model.t_cycle', 0.85)
707
708 self.set('Model.PFC.fac', 0.5)
709 self.set('Model.PFC.epsilon', 0.1)
710 self.set('Model.PFC.fac_pfc', 1)
711 self.set('Model.PFC.stable_threshold', 0.001)
712
713 # Set parameters Circulation
714 self['Tube0D']['A0'] = [0.0004982 , 0.00049959, 0.00045184, 0.00051768]
715 self['Tube0D']['A_wall'] = [1.13597208e-04, 3.79106863e-05, 8.91394588e-05, 4.26698697e-05]
716 self['Tube0D']['k'] = [1.66666667, 2.33333333, 1. , 2.33333333]
717 self['Tube0D']['l'] = [0.4, 0.4, 0.1, 0.1]
718 self['Tube0D']['p0'] = [12154.79898845, 213.26733196, 1913.05480512, 600.72294959]
719
720 self['ArtVen']['p0'] = [6345.26731406, 950. ]
721 self['ArtVen']['q0'] = [4.5e-05, 4.5e-05]
722 self['ArtVen']['k'] = [1. , 1.72]
723
724 # Set volume state variables
725 # not needed to reset, only if initial model crashes
726 if True:
727 self.set('Model.SyArt.V', 2e-4)
728 self.set('Model.PuArt.V', 1e-4)
729 self.set('Model.SyVen.V', 4e-4)
730 self.set('Model.PuVen.V', 2e-4)
731
732 self.set('Model.Peri.TriSeg.V', 56e-6)
733 self.set('Model.Peri.TriSeg.Y', 36e-3)
734
735 self.set('Model.Peri.TriSeg.cLv.V', 150e-6)
736 self.set('Model.Peri.TriSeg.cRv.V', 110e-6)
737 self.set('Model.Peri.La.V', 25e-6)
738 self.set('Model.Peri.Ra.V', 150e-6)
739
740 # all wall parameters
741 self['Patch']['l_se'] = 0.04
742 self['Patch']['l_s_ref'] = 2.0
743 self['Patch']['l_s0'] = 1.8
744 self['Patch']['dl_s_pas'] = 0.6
745 self['Patch']['k1'] = 10
746 self['Patch']['dt'] = 0
747 self['Patch']['C_rest'] = 0
748 self['Patch']['l_si0'] = 1.51
749 self['Patch']['LDAD'] = [1.057, 1.057, 0.64 , 0.64 , 0.64 ]
750 self['Patch']['ADO'] = [0.65, 0.65, 0.75, 0.75, 0.75]
751 self['Patch']['LDCC'] = [4. , 4. , 3.2, 3.2, 3.2]
752 self['Patch']['v_max'] = 7.
753 self['Patch']['v_max'][:2] = 14.
754
755 self['Patch']['tr'] = [0.4 , 0.4 , 0.24, 0.24, 0.24]
756 self['Patch']['td'] = [0.4 , 0.4 , 0.23, 0.23, 0.23]
757
758 # wall specific
759 self['Patch']['Sf_act'] = [80e3, 80e3, 120e3, 120e3, 120e3]
760 self['Patch']['Am_ref'] = [0.00458522, 0.00381765, 0.00806997, 0.00454307, 0.01176514]
761 self['Patch']['V_wall'] = [2.75050615e-05, 1.25959326e-05, 8.80910385e-05, 3.46281364e-05,
762 5.16086110e-05]
763 self['Patch']['Sf_pas'] = [ 16.51117993, 17.0639141 , 580.39235927, 569.67392472,
764 606.50514861]
765
766 # Valves
767 self['Valve']['adaptation_A_open_fac'] = [1. , 1.11, 1. , 1. , 1.11, 1. ]
768 self['Valve']['A_open'] = [0.00050042, 0.00050203, 0.00045228, 0.00051764, 0.00055321,
769 0.00049838]
770 self['Valve']['A_leak'] = [2.64705882e-04, 2.64705882e-10, 2.64705882e-10, 2.64705882e-04,
771 2.64705882e-10, 2.64705882e-10]
772 self['Valve']['l'] = 0.01626978
773 self['Valve']['rho_b'] = 1050
774 self['Valve']['papillary_muscles'] = True
775 self['Valve']['papillary_muscles_slope'] = 100
776 self['Valve']['papillary_muscles_min'] = 0.1
777 self['Valve']['papillary_muscles_A_open_fac'] = 0.1
778 self['Valve']['soft_closure'] = True
779
780 self['Bag']['k'] = 10
781 self['Bag']['V_ref'] = [0.00051968]
782 self['Bag']['p_ref'] = 100
783
784 self.set('Model.Peri.TriSeg.Y', 0.035)
785
786 self['Timings']['law_tau_av'] = 2
787 self['Timings']['c_tau_av0'] = 0.172
788 self['Timings']['c_tau_av1'] = -0.485*60e-3
789
790 # Adaptation
791 self['Patch']['adapt_gamma'] = 0.5
792 self.set('Solver.store_beats', 1)
793
794 self['PressureFlowControl']['stable_threshold'] = 1e-3
795
796 self['Patch']['SfPasMaxT'] = [6400., 6400., 6600., 6600., 6600.]
797 self['Patch']['FacSfActT'] = [0.44, 0.44, 0.61, 0.61, 0.61]
798 self['Patch']['SfPasActT'] = [4800., 4800., 6600., 6600., 6600.]
799 self['Patch']['LsPasActT'] = [3. , 3. , 2.23, 2.23, 2.23]
800
801 self.run(stable=True)
802
803 return
804
805 def get_unittest_targets(self):
806 """Hardcoded results after initializing and running 1 beat."""
807 return {
808 'LVEDV': 119.7,
809 'LVESV': 47.2,
810 }
811
812
813class Walmsley2015(VanOsta2023):
814 def set_reference(self):
815
816 # change model parameters here
817 self['Chamber']['buckling'] = False
818 self['Valve']['soft_closure'] = True
819 self['Valve']['papillary_muscles'] = False
820
821 self['Solver']['dt'] = 0.001
822 self['Solver']['dt_export'] = 0.002
823 self.set('Solver.order', 2)
824
825 self.set('Model.t_cycle', 0.85)
826
827 self.set('Model.PFC.fac', 0.5)
828 self.set('Model.PFC.epsilon', 0.1)
829 self.set('Model.PFC.fac_pfc', 1)
830 self.set('Model.PFC.stable_threshold', 0.001)
831
832 # Set parameters Circulation
833 self['Tube0D']['A0'] = [4.9720e-04, 4.9830e-04, 4.6850e-04, 5.0910e-04]
834 self['Tube0D']['A_wall'] = [1.1370e-04, 6.410e-05, 9.670e-05, 5.950e-05]
835 self['Tube0D']['k'] = [1.66666667, 2.33333333, 1.66666667, 2.33333333]
836 # self['Tube0D']['k'] = [8, 10, 8, 10]
837 self['Tube0D']['l'] = [0.4, 0.4, 0.2, 0.2]
838 self['Tube0D']['p0'] = [12164, 228, 2630 , 678]
839
840 self['ArtVen']['p0'] = [6396, 1500. ]
841 self['ArtVen']['q0'] = [4.5e-05, 4.5e-05]
842 self['ArtVen']['k'] = [1. , 2]
843
844 # Set volume state variables
845 # not needed to reset, only if initial model crashes
846 if True:
847 self.set('Model.SyArt.V', 2e-4)
848 self.set('Model.PuArt.V', 1e-4)
849 self.set('Model.SyVen.V', 4e-4)
850 self.set('Model.PuVen.V', 2e-4)
851
852 self.set('Model.Peri.TriSeg.V', 56e-6)
853 self.set('Model.Peri.TriSeg.Y', 36e-3)
854
855 self.set('Model.Peri.TriSeg.cLv.V', 150e-6)
856 self.set('Model.Peri.TriSeg.cRv.V', 110e-6)
857 self.set('Model.Peri.La.V', 25e-6)
858 self.set('Model.Peri.Ra.V', 150e-6)
859
860 # all wall parameters
861 self['Patch']['l_se0'] = 0.04
862 self['Patch']['l_s_ref'] = 2.0
863 self['Patch']['l_s0'] = [1.8, 1.8, 1.8 , 1.8 , 1.8 ]
864 self['Patch']['dl_s_pas'] = 0.6
865 self['Patch']['k1'] = 10
866 self['Patch']['dt'] = 0
867 self['Patch']['C_rest'] = 0
868 self['Patch']['l_si0'] = 1.51
869 self['Patch']['LDAD'] = [1.057, 1.057, 0.74 , 0.74 , 0.74 ]
870 self['Patch']['ADO'] = [0.081, 0.081, 0.75 , 0.75 , 0.75 ]
871 self['Patch']['LDCC'] = [4. , 4. , 3., 3., 3.]
872 self['Patch']['v_max'] = 7.
873 self['Patch']['v_max'][:2] = 14.
874
875 self['Patch']['tr'] = [0.4 , 0.4 , 0.25, 0.25, 0.25]
876 self['Patch']['td'] = [0.4 , 0.4 , 0.25, 0.25, 0.25]
877
878 # wall specific
879 self['Patch']['Sf_act'] = [60000., 60000., 100000., 100000., 100000.]
880 self['Patch']['Am_ref'] = [0.0071, 0.0058, 0.0099, 0.0052, 0.0135]
881 self['Patch']['V_wall'] = [1.3377e-05, 5.4588e-06, 1.1382e-04, 3.8320e-05,
882 8.0093e-05]
883 self['Patch']['Sf_pas'] = [ 2.8942e3, 3.6579e3 , 794.1495 , 790.6246,
884 801.1993]
885
886 # Valves
887 self['Valve']['adaptation_A_open_fac'] = [1. , 1.12, 1. , 1. , 1.12, 1. ]
888 self['Valve']['A_open'] = [0.00049811, 0.00070281, 0.00046854, 0.00050903, 0.00074578 ,
889 0.00049719]
890 self['Valve']['A_leak'] = [2.6471e-04, 2.6471e-10, 2.6471e-10, 2.6471e-04,
891 2.6471e-10, 2.6471e-10]
892 self['Valve']['l'] = 0.0163
893 self['Valve']['rho_b'] = 1050
894 self['Valve']['papillary_muscles'] = False
895 self['Valve']['papillary_muscles_slope'] = 100
896 self['Valve']['papillary_muscles_min'] = 0.1
897 self['Valve']['papillary_muscles_A_open_fac'] = 0.1
898 self['Valve']['soft_closure'] = True
899
900 self['Bag']['k'] = 10
901 self['Bag']['V_ref'] = [0.00072379]
902 self['Bag']['p_ref'] = 100
903
904 self.set('Model.Peri.TriSeg.Y', 0.035)
905
906 # self['Timings']['law_tauAv'] = 1
907 # self['Timings']['c_tauAv0'] = 0.22
908 # self['Timings']['c_tauAv1'] = -0.03
909
910 # Adaptation
911 # self['Patch2022']['adapt_gamma'] = 0.5
912 self.set('Solver.store_beats', 1)
913
914 self['PFC']['stable_threshold'] = 1e-3
915
916 # self['Patch']['SfPasMaxT'] = [3600., 3600., 4200., 4200., 4200.]
917 # self['Patch']['FacSfActT'] = [0.28, 0.28, 0.69, 0.69, 0.69]
918 # self['Patch']['SfPasActT'] = [2800., 2800., 6600., 6600., 6600.]
919 # self['Patch']['LsPasActT'] = [3. , 3. , 2.31, 2.31, 2.31]
920
921 self.run(stable=True)