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 VanOsta2024 is build following this standard. This model is constructed by the following code.

  1"""
  2VanOsta2024 model.
  3
  4This model is first used in Van Osta et al. under review
  5"""
  6
  7from circadapt.model import Model
  8from circadapt.adapt import ModelAdapt
  9from circadapt.plot import triseg2022
 10import matplotlib.pyplot as plt
 11import numpy as np
 12
 13# define colors for the plot function
 14colors = {
 15    'Lv': np.array([213/255, 96/255, 98/255]),
 16    'Sv': np.array([236/255, 195/255, 11/255]),
 17    'Rv': np.array([6/255, 123/255, 194/255]),
 18    'Ra': np.array([27/255, 231/255, 255/255]),
 19    'La': np.array([250/255, 175/255, 150/255]),
 20    'PuArt': np.array([5/255, 75/255, 125/255]),
 21    'SyArt': np.array([150/255, 25/255, 25/255]),
 22    'RaRv': np.array([0, 0, 0]),
 23    'RvPuArt': np.array([0, 0, 0]),
 24    'LaLv': np.array([0, 0, 0]),
 25    'LvSyArt': np.array([0, 0, 0]),
 26    }
 27colors['RaRv'] = colors['Ra'] #+ 0.5*colors['Rv']
 28colors['RvPuArt'] = colors['Rv'] #+ 0.5*colors['PuArt']
 29colors['LaLv'] = colors['La'] #+ 0.5*colors['Lv']
 30colors['LvSyArt'] = colors['Lv'] #+ 0.5*colors['SyArt']
 31colors['Sv'] = (0.5*colors['Lv'] + 0.5*colors['Rv'])**1.3
 32
 33
 34class VanOsta2024(Model, ModelAdapt):
 35    def __init__(self,
 36                 solver: str = None,
 37                 path_to_circadapt: str = None,
 38                 model_state: dict = None,
 39                 ):
 40        if solver is None:
 41            solver = 'adams_moulton'
 42
 43        self._local_save_reference = True
 44
 45        ModelAdapt.__init__(self)
 46        Model.__init__(self,
 47                       solver,
 48                       path_to_circadapt=path_to_circadapt,
 49                       model_state=model_state,
 50                       )
 51
 52    def _update_settings(self):
 53        # rename model components for easy input/output
 54        self._module_rename['PressureFlowControl'] = 'PFC'
 55
 56    def build(self):
 57        pass
 58        # Circulation
 59        self.add_smart_component('ArtVen', build='SystemicCirculation')
 60        self.add_smart_component('ArtVen', build='PulmonaryCirculation')
 61        self.add_smart_component('Heart', patch_type='Patch', valve_type='Valve')
 62        self.add_smart_component('Timings')
 63        self.add_smart_component('PressureFlowControl')
 64
 65        self.set_component('PFC.circulation_volume_object', 'SyArt')
 66        self.set_component('PFC.circulation_volume_object', 'SyVen')
 67        self.set_component('PFC.circulation_volume_object', 'PuArt')
 68        self.set_component('PFC.circulation_volume_object', 'PuVen')
 69        self.set_component('PFC.circulation_volume_object', 'Peri.La')
 70        self.set_component('PFC.circulation_volume_object', 'Peri.Ra')
 71        self.set_component('PFC.circulation_volume_object', 'Peri.TriSeg.cLv')
 72        self.set_component('PFC.circulation_volume_object', 'Peri.TriSeg.cRv')
 73
 74        # # manually set papillary muscles
 75        self.set_component("Peri.RaRv.wPapMus", "Peri.TriSeg.wRv")
 76        self.set_component("Peri.LaLv.wPapMus", "Peri.TriSeg.wLv")
 77
 78    def set_reference(self):
 79        self['Chamber']['buckling'] = True
 80        self['Valve']['soft_closure'] = True
 81        self['Valve']['papillary_muscles'] = False
 82
 83        self['Solver']['dt'] = 0.001
 84        self['Solver']['dt_export'] = 0.002
 85        self.set('Solver.order',  2)
 86
 87        self.set('Model.t_cycle', 0.85)
 88        self['ArtVen']['p0'] = np.array([6306.25832487, 1000.        ])
 89        self['ArtVen']['q0'] = np.array([4.5e-05, 4.5e-05])
 90        self['ArtVen']['k'] = np.array([1., 2.])
 91        self['Tube0D']['l'] = np.array([0.4, 0.4, 0.2, 0.2])
 92        self['Tube0D']['A_wall'] = np.array([1.12362733e-04, 6.57883944e-05, 9.45910889e-05, 8.22655361e-05])
 93        self['Tube0D']['k'] = np.array([1.66666667, 2.33333333, 1.66666667, 2.33333333])
 94        self['Tube0D']['p0'] = np.array([12162.50457811,   287.7083132 ,  2132.51755623,   830.54673184])
 95        self['Tube0D']['A0'] = np.array([0.0004983 , 0.00049909, 0.00047138, 0.00050803])
 96        self['Tube0D']['target_wall_stress'] = np.array([500000., 500000., 500000., 500000.])
 97        self['Tube0D']['target_mean_flow'] = np.array([0.17, 0.17, 0.17, 0.17])
 98        self['Bag']['k'] = np.array([10.])
 99        self['Bag']['p_ref'] = np.array([1000.])
100        self['Bag']['V_ref'] = np.array([0.00054267])
101        self['Patch']['Am_ref'] = np.array([0.00425687, 0.00401573, 0.00966859, 0.00289936, 0.01084227])
102        self['Patch']['V_wall'] = np.array([4.46069398e-06, 2.14548521e-06, 7.35720515e-05, 1.88904978e-05, 3.67720116e-05,])
103        self['Patch']['v_max'] = np.array([14., 14.,  7.,  7.,  7.])
104        self['Patch']['l_se0'] = np.array([0.04, 0.04, 0.04, 0.04, 0.04])
105        self['Patch']['l_s0'] = np.array([1.8, 1.8, 1.8, 1.8, 1.8])
106        self['Patch']['l_s_ref'] = np.array([2., 2., 2., 2., 2.])
107        self['Patch']['dl_s_pas'] = np.array([0.6, 0.6, 0.6, 0.6, 0.6])
108        self['Patch']['Sf_pas'] = np.array([2248.53598   , 2684.76100348,  731.24545453,  729.06063771, 749.47694522,])
109        self['Patch']['tr'] = np.array([0.4 , 0.4 , 0.25, 0.25, 0.25])
110        self['Patch']['td'] = np.array([0.4 , 0.4 , 0.25, 0.25, 0.25])
111        self['Patch']['time_act'] = np.array([0.15 , 0.15 , 0.425, 0.425, 0.425])
112        self['Patch']['Sf_act'] = np.array([ 84000.,  84000., 120000., 120000., 120000.])
113        self['Patch']['fac_Sf_tit'] = np.array([0.01, 0.01, 0.01, 0.01, 0.01])
114        self['Patch']['k1'] = np.array([10., 10., 10., 10., 10.])
115        self['Patch']['dt'] = np.array([0., 0., 0., 0., 0.])
116        self['Patch']['C_rest'] = np.array([0., 0., 0., 0., 0.])
117        self['Patch']['l_si0'] = np.array([1.51, 1.51, 1.51, 1.51, 1.51])
118        self['Patch']['LDAD'] = np.array([1.057, 1.057, 1.057, 1.057, 1.057])
119        self['Patch']['ADO'] = np.array([0.65, 0.65, 0.65, 0.65, 0.65])
120        self['Patch']['LDCC'] = np.array([4., 4., 4., 4., 4.])
121        self['Patch']['SfPasMaxT'] = np.array([320000., 320000.,  16000.,  16000.,  16000.])
122        self['Patch']['SfPasActT'] = np.array([40000., 40000., 20000., 20000., 20000.])
123        self['Patch']['FacSfActT'] = np.array([0.8, 0.8, 1. , 1. , 1. ])
124        self['Patch']['LsPasActT'] = np.array([2.42, 2.42, 2.42, 2.42, 2.42])
125        self['Patch']['adapt_gamma'] = np.array([0.5, 0.5, 0.5, 0.5, 0.5])
126        self['Patch']['transmat00'] = np.array([-0.5751, -0.5751, -0.5751, -0.5751, -0.5751])
127        self['Patch']['transmat01'] = np.array([-0.7851, -0.7851, -0.7851, -0.7851, -0.7851])
128        self['Patch']['transmat02'] = np.array([0.6063, 0.6063, 0.6063, 0.6063, 0.6063])
129        self['Patch']['transmat03'] = np.array([-0.5565, -0.5565, -0.5565, -0.5565, -0.5565])
130        self['Patch']['transmat10'] = np.array([-0.1279, -0.1279, -0.1279, -0.1279, -0.1279])
131        self['Patch']['transmat11'] = np.array([0.0999, 0.0999, 0.0999, 0.0999, 0.0999])
132        self['Patch']['transmat12'] = np.array([0.2066, 0.2066, 0.2066, 0.2066, 0.2066])
133        self['Patch']['transmat13'] = np.array([-1.8441, -1.8441, -1.8441, -1.8441, -1.8441])
134        self['Patch']['transmat20'] = np.array([-0.1865, -0.1865, -0.1865, -0.1865, -0.1865])
135        self['Patch']['transmat21'] = np.array([-0.201, -0.201, -0.201, -0.201, -0.201])
136        self['Patch']['transmat22'] = np.array([1.3195, 1.3195, 1.3195, 1.3195, 1.3195])
137        self['Patch']['transmat23'] = np.array([-11.8745, -11.8745, -11.8745, -11.8745, -11.8745])
138        self['Valve']['adaptation_A_open_fac'] = np.array([1., 1., 1., 1., 1., 1.])
139        self['Valve']['A_open'] = np.array([0.00049916, 0.00047155, 0.00047155, 0.00050805, 0.00049835,
140               0.00049835])
141        self['Valve']['A_leak'] = np.array([0.00049916, 1.e-09, 1.e-09, 0.00050805, 1.e-09, 1.e-09])
142        self['Valve']['l'] = np.array([0.01260512, 0.01225144, 0.01225144, 0.01271679, 0.01259479,
143               0.01259479])
144        self['Valve']['L_fac_prox'] = np.array([0.75, 0.75, 0.75, 0.75, 0.75, 0.75])
145        self['Valve']['L_fac_dist'] = np.array([0.75, 0.75, 0.75, 0.75, 0.75, 0.75])
146        self['Valve']['L_fac_valve'] = np.array([1.5, 1.5, 1.5, 1.5, 1.5, 1.5])
147        self['Valve']['rho_b'] = np.array([1050., 1050., 1050., 1050., 1050., 1050.])
148        self['Valve']['papillary_muscles'] = np.array([False, False, False, False, False, False])
149        self['Valve']['papillary_muscles_slope'] = np.array([100., 100., 100., 100., 100., 100.])
150        self['Valve']['papillary_muscles_min'] = np.array([0.1, 0.1, 0.1, 0.1, 0.1, 0.1])
151        self['Valve']['papillary_muscles_A_open_fac'] = np.array([0.1, 0.1, 0.1, 0.1, 0.1, 0.1])
152        self['Valve']['soft_closure'] = np.array([ True,  True,  True,  True,  True,  True])
153        self['Valve']['fraction_A_open_Aext'] = np.array([0.9, 0.9, 0.9, 0.9, 0.9, 0.9])
154        self['Timings']['time_fac'] = np.array([1.])
155        self['Timings']['tau_av'] = np.array([0.150025])
156        self['Timings']['dtau_av'] = np.array([0.])
157        self['Timings']['law_tau_av'] = np.array([1])
158        self['Timings']['law_Ra2La'] = np.array([1])
159        self['Timings']['law_ta'] = np.array([1])
160        self['Timings']['law_tv'] = np.array([1])
161        self['Timings']['c_tau_av0'] = np.array([0.])
162        self['Timings']['c_tau_av1'] = np.array([0.1765])
163        self['Timings']['c_ta_rest'] = np.array([0.])
164        self['Timings']['c_ta_tcycle'] = np.array([0.17647059])
165        self['Timings']['c_tv_rest'] = np.array([0.1])
166        self['Timings']['c_tv_tcycle'] = np.array([0.4])
167        self['PFC']['p0'] = np.array([12200.])
168        self['PFC']['q0'] = np.array([8.5e-05])
169        self['PFC']['stable_threshold'] = np.array([0.0001])
170        self['PFC']['is_active'] = np.array([ True])
171        self['PFC']['fac'] = np.array([1.])
172        self['PFC']['fac_pfc'] = np.array([1.])
173        self['PFC']['epsilon'] = np.array([0.4])
174        self['TriSeg']['tau'] = np.array([2.])
175        self['TriSeg']['ratio_septal_LV_Am'] = np.array([0.3])
176        self['TriSeg']['max_number_of_iterations'] = np.array([100.])
177        self['TriSeg']['thresh_F'] = np.array([0.001])
178        self['TriSeg']['thresh_dV'] = np.array([1.e-09])
179        self['TriSeg']['thresh_dY'] = np.array([1.e-06])
180
181        self.run(stable=True)
182        # self.adapt()
183
184        # set target volume
185        self['PFC']['target_volume'] = self['PFC']['circulation_volume']
186
187        return
188
189
190    def get_unittest_targets(self):
191        """Hardcoded results after initializing and running 1 beat."""
192        return {
193            'LVEDV': 120.5,
194            'LVESV':  48.3,
195            }
196
197    def get_unittest_results(self, model):
198        """Real-time results after initializing and running 1 beat."""
199        LVEDV = np.max(model['Cavity']['V'][:, 'cLv'])*1e6
200        LVESV = np.min(model['Cavity']['V'][:, 'cLv'])*1e6
201        return {
202            'LVEDV': LVEDV,
203            'LVESV': LVESV,
204            }
205
206
207    def plot(self, fig=None):
208        # TODO
209        self.plot_extended(fig)
210
211    def plot_extended(self, fig=None):
212        if len(self['Solver']['t'])==0:
213            raise RuntimeError('Can not plot the model, as no data is available. Did you simulate the model?')
214
215        if fig is None:
216            fig = 1
217        if isinstance(fig, int):
218            fig = plt.figure(fig, clear=True, figsize=(12, 8))
219
220        # Settings
221        grid_size = [32, 32]
222
223        def get_lim(module, signal, locs=slice(None, None, None)):
224            signal = self[module][signal][:, locs]
225            lim = np.array([np.min(signal), np.max(signal)])
226            lim += np.array([-1, 1]) * 0.1*np.diff(lim)
227            return lim
228
229        lim_V = get_lim('Cavity', 'V', ['cRv', 'Ra', 'La', 'cLv']) * 1e6
230        lim_V[0] = 0
231        lim_p = get_lim('Cavity', 'p', ['cRv', 'Ra', 'La', 'cLv']) / 133
232        lim_p[0] = np.min([lim_p[0], 0])
233        lim_Ls = get_lim('Patch', 'l_s')
234        lim_Sf = get_lim('Patch', 'Sf') * 1e-3
235        lim_q = get_lim('Valve', 'q',
236                        ['LaLv', 'RaRv', 'LvSyArt', 'RaPuArt']) * 1e6
237
238        all_lim = [lim_V, lim_p, lim_Ls, lim_Sf, lim_q]
239        if (np.any(np.isnan(all_lim)) or np.any(np.isinf(all_lim))):
240            lim_V = [0, 200]
241            lim_p = [0, 150]
242            lim_Ls = [1.5, 2.0]
243            lim_Sf = [0, 100]
244            lim_q = [-1e-3, 1e-3]
245
246        # Pressure Volume plot
247        axPV = plt.subplot2grid(grid_size, (0, 17), rowspan=15, colspan=15, fig=fig)
248        axPV.plot(self['Cavity']['V'][:, 'cLv']*1e6,
249                  self['Cavity']['p'][:, 'cLv']/133,
250                  c=colors['Lv'],
251                  )
252        axPV.plot(self['Cavity']['V'][:, 'cRv']*1e6,
253                  self['Cavity']['p'][:, 'cRv']/133,
254                  c=colors['Rv'],
255                  )
256        axPV.plot(self['Cavity']['V'][:, 'La']*1e6,
257                  self['Cavity']['p'][:, 'La']/133,
258                  c=colors['La'],
259                  )
260        axPV.plot(self['Cavity']['V'][:, 'Ra']*1e6,
261                  self['Cavity']['p'][:, 'Ra']/133,
262                  c=colors['Ra'],
263                  )
264        axPV.spines[['top', 'right']].set_visible(False)
265        axPV.set_title('Pressure-Volume loop', weight='bold')
266        axPV.set_xlabel('Volume [mL]')
267        axPV.set_ylabel('Pressure [mmHg]')
268        axPV.spines[['bottom', 'left']].set_position(('outward', 5))
269
270        ylabel_x_left = -0.25
271        ylabel_x_right = 1.25
272
273        # Volumes
274        t = self['Solver']['t']*1e3
275
276        axVRv = plt.subplot2grid(grid_size, (0, 0), rowspan=8, colspan=6, fig=fig)
277        axVRv.plot(t, self['Cavity']['V'][:, 'cRv']*1e6,
278                   c=colors['Rv'],
279                   )
280        axVRv.plot(t, self['Cavity']['V'][:, 'Ra']*1e6,
281                   c=colors['Ra'],
282                   )
283        axVRv.set_ylim(lim_V)
284        axVRv.set_ylabel('Volume\n[mL]')
285        axVRv.spines[['top', 'right']].set_visible(False)
286        axVRv.set_title('Right Heart', weight='bold')
287        # axVRv.set_xticks([])
288        axVRv.tick_params(axis='both', direction='in')
289        axVRv.yaxis.set_label_coords(ylabel_x_left, 0.5)
290
291
292        axVLv = plt.subplot2grid(grid_size, (0, 6), rowspan=8, colspan=6, fig=fig)
293        axVLv.plot(t, self['Cavity']['V'][:, 'cLv']*1e6,
294                   c=colors['Lv'],
295                   )
296        axVLv.plot(t, self['Cavity']['V'][:, 'La']*1e6,
297                   c=colors['La'],
298                   )
299        axVLv.set_ylabel('Volume\n[mL]')
300        axVLv.set_ylim(lim_V)
301        axVLv.yaxis.set_ticks_position('right')
302        axVLv.yaxis.set_label_position('right')
303        axVLv.spines['right'].set_position(('outward', 0))
304        axVLv.spines[['top', 'left']].set_visible(False)
305        axVLv.set_title('Left Heart', weight='bold')
306        # axVLv.set_xticks([])
307        axVLv.tick_params(axis='both', direction='in')
308        axVLv.yaxis.set_label_coords(ylabel_x_right, 0.5)
309
310        # Pressures
311        axpRv = plt.subplot2grid(
312            grid_size, (8, 0), rowspan=8, colspan=6, fig=fig, sharex=axVRv)
313        axpRv.plot(t, self['Cavity']['p'][:, 'cRv']/133,
314                   c=colors['Rv'],
315                   )
316        axpRv.plot(t, self['Cavity']['p'][:, 'Ra']/133,
317                   c=colors['Ra'],
318                   )
319        axpRv.plot(t, self['Cavity']['p'][:, 'PuArt']/133,
320                   c=colors['PuArt'],
321                   )
322        axpRv.spines[['top', 'right']].set_visible(False)
323        # axpRv.set_xticks([])
324        axpRv.tick_params(axis='both', direction='in')
325        axpRv.set_ylim(lim_p)
326        axpRv.set_ylabel('Pressure\n[mmHg]')
327        axpRv.yaxis.set_label_coords(ylabel_x_left, 0.5)
328
329        axpLv = plt.subplot2grid(grid_size, (8, 6), rowspan=8, colspan=6, fig=fig, sharex=axVLv)
330        axpLv.plot(t, self['Cavity']['p'][:, 'cLv']/133,
331                   c=colors['Lv'],
332                   )
333        axpLv.plot(t, self['Cavity']['p'][:, 'La']/133,
334                   c=colors['La'],
335                   )
336        axpLv.plot(t, self['Cavity']['p'][:, 'SyArt']/133,
337                   c=colors['SyArt'],
338                   )
339        axpLv.yaxis.set_ticks_position('right')
340        axpLv.yaxis.set_label_position('right')
341        axpLv.spines['right'].set_position(('outward', 0))
342        axpLv.spines[['top', 'left']].set_visible(False)
343        # axpLv.set_xticks([])
344        axpLv.tick_params(axis='both', direction='in')
345        axpLv.set_ylim(lim_p)
346        axpLv.set_ylabel('Pressure\n[mmHg]')
347        axpLv.yaxis.set_label_coords(ylabel_x_right, 0.5)
348
349        # Valves
350        ax = plt.subplot2grid(grid_size, (16, 0), rowspan=6, colspan=6, fig=fig, sharex=axVRv)
351        ax.plot(t, self['Valve']['q'][:, 'RaRv']*1e6,
352                c=colors['RaRv'],
353                )
354        ax.plot(t, self['Valve']['q'][:, 'RvPuArt']*1e6,
355                c=colors['RvPuArt'],
356                )
357        ax.spines[['top', 'right']].set_visible(False)
358        ax.set_ylim(lim_q)
359        # ax.set_xticks([])
360        ax.set_ylabel('Flow\n[mL/s]')
361        ax.yaxis.set_label_coords(ylabel_x_left, 0.5)
362
363        ax = plt.subplot2grid(grid_size, (16, 6), rowspan=6, colspan=6, fig=fig, sharex=axVLv)
364        ax.plot(t, self['Valve']['q'][:, 'LaLv']*1e6,
365                   c=colors['LaLv'],
366                   )
367        ax.plot(t, self['Valve']['q'][:, 'LvSyArt']*1e6,
368                   c=colors['LvSyArt'],
369                   )
370        ax.spines[['top', 'left']].set_visible(False)
371        ax.set_ylim(lim_q)
372        ax.yaxis.set_ticks_position('right')
373        ax.yaxis.set_label_position('right')
374        # ax.set_xticks([])
375        ax.set_ylabel('Flow\n[mL/s]')
376        ax.yaxis.set_label_coords(ylabel_x_right, 0.5)
377
378        # Stress
379        ax = plt.subplot2grid(grid_size, (22, 0), rowspan=4, colspan=6, fig=fig, sharex=axVRv)
380        ax.plot(t, self['Patch']['Sf'][:, 'pRv0']*1e-3,
381                   c=colors['Rv'],
382                   )
383        ax.plot(t, self['Patch']['Sf'][:, 'pRa0']*1e-3,
384                   c=colors['Ra'],
385                   )
386        ax.spines[['top', 'right']].set_visible(False)
387        # ax.set_xticks([])
388        ax.set_ylim(lim_Sf)
389        ax.set_ylabel('Total\nstress [kPa]')
390        ax.yaxis.set_label_coords(ylabel_x_left, 0.5)
391
392        ax = plt.subplot2grid(grid_size, (22, 6), rowspan=4, colspan=6, fig=fig, sharex=axVLv)
393        ax.plot(t, self['Patch']['Sf'][:, 'pLv0']*1e-3,
394                   c=colors['Lv'],
395                   )
396        ax.plot(t, self['Patch']['Sf'][:, 'pSv0']*1e-3,
397                   c=colors['Sv'],
398                   )
399        ax.plot(t, self['Patch']['Sf'][:, 'pLa0']*1e-3,
400                   c=colors['La'],
401                   )
402        ax.spines[['top', 'left']].set_visible(False)
403        ax.yaxis.set_ticks_position('right')
404        ax.yaxis.set_label_position('right')
405        # ax.set_xticks([])
406        ax.set_ylim(lim_Sf)
407        ax.set_ylabel('Total\nstress [kPa]')
408        ax.yaxis.set_label_coords(ylabel_x_right, 0.5)
409
410        # Sarcomere Length
411        ax = plt.subplot2grid(grid_size, (26, 0), rowspan=6, colspan=6, fig=fig, sharex=axVRv)
412        ax.plot(t, self['Patch']['l_s'][:, 'pRv0'],
413                   c=colors['Rv'],
414                   )
415        ax.plot(t, self['Patch']['l_s'][:, 'pRa0'],
416                   c=colors['Ra'],
417                   )
418        ax.spines[['top', 'right']].set_visible(False)
419        ax.spines['bottom'].set_position(('outward', 5))
420        ax.set_ylim(lim_Ls)
421        ax.set_ylabel('Sarcomere\nlength [$\\mu$m]')
422        ax.yaxis.set_label_coords(ylabel_x_left, 0.5)
423
424        ax = plt.subplot2grid(grid_size, (26, 6), rowspan=6, colspan=6, fig=fig, sharex=axVLv)
425        ax.plot(t, self['Patch']['l_s'][:, 'pLv0'],
426                   c=colors['Lv'],
427                   )
428        ax.plot(t, self['Patch']['l_s'][:, 'pSv0'],
429                   c=colors['Sv'],
430                   )
431        ax.plot(t, self['Patch']['l_s'][:, 'pLa0'],
432                   c=colors['La'],
433                   )
434        ax.spines[['top', 'left']].set_visible(False)
435        ax.spines['bottom'].set_position(('outward', 5))
436        ax.yaxis.set_ticks_position('right')
437        ax.yaxis.set_label_position('right')
438        ax.set_ylim(lim_Ls)
439        ax.set_ylabel('Sarcomere\nlength [$\\mu$m]')
440        ax.yaxis.set_label_coords(ylabel_x_right, 0.5)
441
442        # ax.set_xlabel('Time [ms]')
443        # ax.xaxis.set_label_coords(0, -0.3)
444
445
446        # Plot TriSeg
447        titles = ['Pre-A', 'Onset QRS', 'Peak LV \n pressure', 'AV close']
448        idx = [0,
449               np.argmax(np.diff(self['Patch']['C'][:, 'pLv0'])>0),
450               np.argmax(self['Cavity']['p'][:, 'cLv']),
451               len(t) - 1 - np.argmax(
452                   np.diff(self['Valve']['q'][:, 'LvSyArt'][::-1])>0)
453               ]
454
455        for i in range(4):
456            ax = plt.subplot2grid(grid_size, (26, 16+4*i),
457                                  rowspan=5, colspan=4, fig=fig)
458            triseg2022(self, ax, idx[i],
459                       colors=[colors['Lv'], colors['Sv'], colors['Rv']])
460            ax.spines[['top', 'right', 'bottom', 'left']].set_visible(False)
461            ax.set_xticks([])
462            ax.set_yticks([])
463            plt.xlabel(titles[i], fontsize=12)
464
465        # Plot settings
466        plt.subplots_adjust(
467            top=0.96,
468            bottom=0.05,
469            left=0.075,
470            right=0.98,
471            hspace=5,
472            wspace=0.5)
473        plt.draw()
474
475        # Stress strain
476        ax_right_stress_strain = plt.subplot2grid(grid_size, (18, 17),
477                              rowspan=7, colspan=7, fig=fig)
478        ax_left_stress_strain = plt.subplot2grid(grid_size, (18, 24),
479                              rowspan=7, colspan=7, fig=fig)
480
481
482        ax_right_stress_strain.plot(
483            self['Patch']['Ef'][:, 'pRa0'],
484            self['Patch']['Sf'][:, 'pRa0']*1e-3,
485            c=colors['Ra'],
486            )
487        ax_right_stress_strain.plot(
488            self['Patch']['Ef'][:, 'pRv0'],
489            self['Patch']['Sf'][:, 'pRv0']*1e-3,
490            c=colors['Rv'],
491            )
492        ax_left_stress_strain.plot(
493            self['Patch']['Ef'][:, 'pLa0'],
494            self['Patch']['Sf'][:, 'pLa0']*1e-3,
495            c=colors['La'],
496            )
497        ax_left_stress_strain.plot(
498            self['Patch']['Ef'][:, 'pLv0'],
499            self['Patch']['Sf'][:, 'pLv0']*1e-3,
500            c=colors['Lv'],
501            )
502        ax_left_stress_strain.plot(
503            self['Patch']['Ef'][:, 'pSv0'],
504            self['Patch']['Sf'][:, 'pSv0']*1e-3,
505            c=colors['Sv'],
506            )
507
508        ylim = [np.min(self['Patch']['Sf'])*1e-3,
509                np.max(self['Patch']['Sf'])*1e-3]
510        ylim_ptp = np.ptp(ylim)*0.025
511        ylim[0] -= ylim_ptp
512        ylim[1] += ylim_ptp
513        ax_right_stress_strain.set_ylim(ylim)
514        ax_left_stress_strain.set_ylim(ylim)
515        xlim = [np.min(self['Patch']['Ef']),
516                np.max(self['Patch']['Ef'])]
517        xlim_ptp = np.ptp(xlim)*0.025
518        xlim[0] -= xlim_ptp
519        xlim[1] += xlim_ptp
520        ax_right_stress_strain.set_xlim(xlim)
521        ax_left_stress_strain.set_xlim(xlim)
522
523        ax_right_stress_strain.spines[['left', 'bottom']].set_position(('outward', 5))
524        ax_left_stress_strain.spines[['right', 'bottom']].set_position(('outward', 5))
525
526        ax_right_stress_strain.spines[['top', 'right']].set_visible(False)
527        ax_right_stress_strain.set_ylabel('[kPa]', rotation=0, va='bottom', ha='right')
528        ax_right_stress_strain.yaxis.set_label_coords(0., 1.05)
529        ax_right_stress_strain.set_xlabel('Strain [-]')
530        ax_left_stress_strain.spines[['top', 'left']].set_visible(False)
531        ax_left_stress_strain.yaxis.set_ticks_position('right')
532        ax_left_stress_strain.yaxis.set_label_position('right')
533        ax_left_stress_strain.set_ylabel('Stress\n[kPa]', rotation=0, va='bottom', ha='left')
534        ax_left_stress_strain.yaxis.set_label_coords(1., 1.05)
535        ax_left_stress_strain.set_xlabel('Strain [-]')
536