Left Ventricle Geometry and Shape Morphing
A 3D left ventricular myocardium at end-diastole (ED) was reconstructed from the cardiac magnetic resonance (CMR) imaging of a healthy adult human volunteer. The protocol was approved by the Surrey Research Ethics Committee (protocol 10/H0701/112), and informed consent were obtained from all participants. Segmentation and reconstruction were performed with VMTK (www.vmtk.org), while smoothing was performed with Geomagic Studio (Geomagic Inc., Morrisville, NC, USA). The healthy LV model was later morphed to one with concentric hypertrophy and one with dilated eccentric hypertrophy, to investigate biomechanics during such LV geometry alternations. Concentric hypertrophy was modelled as a 100% increase in wall thickness on the load-free geometry, via an offset of the epicardium outwards, while dilation was modelled as an 80% increase in end-diastolic volume (EDV) without a change in wall thickness on the load-free geometry, via an offset of both the epicardium and endocardium outwards (Fig. 2). These were in accordance with clinical measurements of wall thicknesses and EDV for HCM and DCM patients (Nielles-Vallespin et al. 2017b).
Computational Finite Element Modeling Framework
A computational finite element method was employed from the previous studies coupled with a LV close-loop Windkessel circulatory (Ong et al. 2021; Shavik et al. 2018). Simulations were conducted with FeniCS (www.fenicsproject.org), by minimizing the Lagrangian function detailed by (Shavik et al. 2018). FEM utilized a transversely isotropic Fung-type strain energy function to back compute the loading-free geometry with assigned end-diastolic pressure and myocardial stiffness information. The passive myocardial constitutive model, Fung-type strain energy function WLV, was given by:
Where Q was calculated as:
C in the Eq. (1a) was the global myocardial stiffness, while \({B}_{FF}\), \({B}_{SS}\), \({B}_{NN}\), \({B}_{NS}\), and \({B}_{FS}\) and \({B}_{FN}\) were the passive material parameters in various directions, where F, S, and N denoting the myocardial fiber, sheetlet, and normal directions, respectively. As shown in Fig. 1, F is the orientation of the myocyte, S is perpendicular to F in the sheetlet plane, while N is normal to the sheetlet plane. In particular, \({B}_{NS}\) denoted the shear stiffness in the sheet and normal directions, and was the parameter that was reduced to model sheetlet sliding. E was the Green-Lagrange strain tensor. For simulations of healthy LV, C was assumed to be 100 Pa, which was consistent with past simulation work (Rumindo et al. 2020; Shavik et al. 2018), but for HCM and DCM diseased conditions, they were increased to 200 Pa or 300 Pa, as informed by findings that diastolic dysfunction increases myocardial stiffness by 2–3 times (Klotz et al. 2005; Wang et al. 2018; Westermann et al. 2008).
The Guccione model was employed to simulate the active myocardial contractile mechanical behavior (Shavik et al. 2018), which is modelled as the maximum tension (\({T}_{max}\)) multiplied by a calcium activation curve over time, details of which are given in the Appendix. For healthy LVs, \({T}_{max}\) was assumed to be 150 kPa, but for diseased HCM and DCM LVs, they were assumed to be 105 or 75 kPa, in accordance to previous simulation work (Shavik et al. 2021).
The boundary conditions for the FEM simulations were a constraint at the LV base on out-of-plane motion and a low-stiffness (60 Pa) spring constraint on the entire epicardial surface to emulate interactions with surrounding tissues, and more details were given in the Appendix.
Before simulating LV systolic behavior, we first estimated the load-free geometry. To do this, we assume that the end diastolic pressure (EDP) of all LVs to be 5mmHg (Westermann et al. 2008), and we estimated the unloading deformation from the end-diastole state to the load-free state as the inverse of the loading deformation for the same pressure difference. Once the load-free geometry is obtained, we morphed the healthy load-free LV to HCM and DCM LV as described above. With the load-free geometry as the new starting point, we performed FEM simulation of the entire cycle of the LV. Simulations were conducted for 10 cycles to allow the lumped parameter model to converge.
A simplified Windkessel lumped-parameter model was coupled to the LV FEM, as shown in the Fig. 2D. It consists of peripheral vascular and venous resistances (Rper and Rven) and aortic valve and mitral valve resistances (Rao and Rmv), and arterial and venous compliances (Cart, and Cven). Initial volumes of arterial and venous were tuned together with values of resistances and compliances to obtain the expected pressure-volume loop for the normal LV. Thereafter, the same lumped parameter model was used for other LV geometries and cardiac contractilities and passive stiffnesses. Details of the lumped parameter model parameters are given in the Appendix.
Estimation of Myocardial Normal Stiffness and Shear Stiffness
To obtain the specification of the passive stiffness model, we performed numerical modelling of simple shear mechanical testing of a cuboid piece of myocardium using our passive stiffness model, to match data obtained by mechanical testing experiments (Dokos et al. 2002; Sommer et al. 2015). The relationship between myocardial passive stress (\({\sigma }_{p}\)) and Green-Lagrange strain (E) was modelled via finite strain theory as:
Where \({F}_{matrix}\) is the deformation gradient tensor and J is the Jacobian of the deformation gradient tensor \({F}_{matrix}\). According to Sommer et al., the myocardial stiffness in the fiber direction was about twice as stiff as in the cross-fiber direction from the biaxial extension testing (Sommer et al. 2015). Therefore, we specified \({B}_{FF}\) to be 29.8, twice that of \({B}_{SS}\) and \({B}_{NN}\), which were 14.9, in accordance to the previous FEM studies (Rumindo et al. 2020; Shavik et al. 2018).
To model the sheetlet sliding, we reduced \({B}_{SN}\) and \({B}_{NS}\) to be the same and reduced from other shear stiffness components \({B}_{FS}\), \({B}_{SF}\), \({B}_{NF}\) and \({B}_{FN}\). From Sommer et al (Sommer et al. 2015), the shear stress to achieve a shear deformation of 0.4 in the F-S or F-N direction was 4.6\(\pm\)1.0 kPa or 4.2\(\pm\)1.3 kPa. From our modelling of simple shear in the F-S or F-N direction, a \({B}_{FS}\) or \({B}_{FN}\) of 19.2 would match this behavior. From Sommer et al, the shear stress to achieve a shear deformation of 0.4 in the N-S or S-N direction was 2.2\(\pm\)0.8 kPa. From our modelling of simple shear in the N-S direction, a \({B}_{NS}\) of 15.4 would match this behavior. We performed FEM modelling with four values of \({B}_{NS}\), 9.3 12.1, 15.0 and 17.8, representing − 1.6, -1.0, -0.2, + 1.0 standard deviations from the mean N-S stress value from Sommer et al.’s data. These range of stiffness were within ranges investigated by previous FEM studies (Rumindo et al. 2020; Zhang et al. 2021).
We further tested reducing shear stiffness \({B}_{FS}\) and \({B}_{FN}\), on top of reducing \({B}_{SN}\) and \({B}_{NS}\), but found that \({B}_{FS}\) and \({B}_{FN}\) had very minimal influence on cardiac functions (Appendix table A.3), likely because cardiac deformations did not engage shear in these directions much. \({B}_{FS}\) and \({B}_{FN}\) were thus held constant value 19.2 in our simulations, and only reduced stiffness in \({B}_{SN}\) and \({B}_{NS}\) were used to model sheetlet sliding.
Assignments of Myocardial Orientations
Myocardial myocyte helix angle (HA) was defined as the angle between the projection of myocyte (F) direction onto the local longitudinal – circumferential plane and the circumferential axis (Fig. 3). Helix angle of the healthy geometry was set to vary linearly from + 60\(^\circ\) and − 60\(^\circ\) from the endocardium to the epicardium. We assumed that the healthy geometry could be transformed into diseased geometries via a homogeneous deformation, increased wall thickness for HCM and LV dilation for DCM. We further assumed that myocyte orientations would undergo the same transformation, to be realigned according to the deformation. Thus, a dilation in LV diameter would stretch the myocardium circumferentially and reduce the helix angle magnitudes via a cosine rule:
$$\begin{array}{c}\frac{\text{cos}\left({\theta }_{healthy}\right)}{\text{cos}\left({\theta }_{diseased}\right)}=\frac{{L}_{longi,healthy}/{L}_{circ, healthy}}{{L}_{longi, diseased}/{L}_{circ,diseased}}=\frac{{D}_{diseased}}{{D}_{healthy}}\#\left(3\right)\end{array}$$
Where\(\theta\) was the helix orientation, \(L\) was the projected length of the myofiber, projected to the longitudinal or circumferential direction (indicated by subscripts, longi. or circ.), and \(D\) was the diameter of the LV. Based on this, myocyte helix orientation of HCM was calculated to be + 60\(^\circ\) to -51.4\(^\circ\) from endocardium to epicardium, and for DCM was + 49.7\(^\circ\) and − 51.1\(^\circ\).
Myocardial sheetlet angle (SA) was the angle between the sheet (S) direction and its projection on the local longitudinal – circumferential plane (Fig. 3). We investigated three diastolic SA (DSA) in the healthy LV, 0\(^\circ\) and 18\(^\circ\) and 48\(^\circ\), to gauge the effects of the diastolic sheetlet angle on the LV functions. 18\(^\circ\) was the average diastolic sheetlet angle for healthy and DCM LVs while 48\(^\circ\) was that for HCM LV (Nielles-Vallespin et al. 2017b).
Myofiber transverse angle (TA) was the angle between the projection of myocyte (F) direction on the local radial – circumferential plane and the circumferential axis (Fig. 3). Vendelin et al. found that TA near to 10\(^\circ\) provided the best cardiac function efficiency (Vendelin et al. 2002). Here, we tested four transverse angles, 0\(^\circ\) and 10\(^\circ\) and 20\(^\circ\) and 30\(^\circ\). The spatial variation in transverse angles, \(\delta\), was modelled as follows, as proposed by Vendelin et al.,
$$\begin{array}{c}\delta =transverse angle*\omega \left(1-{\left(1-2\xi \right)}^{2}\right)\#\left(4\right)\end{array}$$
Where \(\xi\) was the normalized distance linearly ranging from − 1 in the endocardium surface to 1 in the epicardial surface, \(\omega\) was a linearly varied coefficient ranging from 0.5 at the base to -1 at the apex.