Introduction

FEA software development

Voxel based meshing

Voxel based 3D thinning

Material modeling of MC

Motion tracking

Stereo reconstruction

Current Project

The equivalent mechanical properties of the trabeculated myocardium (MC) are very difficult to be defined experimentally due to the prevailing boundary effects. In this work, we develop a numerical homogenization procedure for deriving the constitutive relations of the trabeculated myocardium on the basis of voxel modeling, Finite Element (FE) analysis, and nonlinear data regression. The new constitutive relations allow us to represent the effect of trabeculation in a continuous FE model, in which the trabeculated wall is replaced with homogeneous material. The result is essential for establishing an appropriate computational model of the entire trabeculated heart, and for future theoretical investigations.

Basic approach

1. Propose possible forms of constitutive relations defined in terms of strain energy density (SED) functions based on the preliminary studies of the general characteristics of the trabeculated myocardium.
2. Obtain three-dimensional stress-strain relations through FE simulations of a series of stretch experiments designed to generate tri-axial stress states on a representative volume element (RVE) of the trabeculated myocardium.
3. Use a modified multi-dimensional Levenberg-Marquardt procedure to compute the best-fit coefficients in the SED functions based on the obtained stress-strain relations.

back to top

Strain energy density functions
In our computational model, the material properties of the myocardial core in the trabeculae are defined with an exponential SED function

W = A [exp(B(I1-3)) - 1.0] + C Eff

The first term represents the passive constitutive relation of the incompressible homogeneous matrix. The second term is the additional muscle fiber contribution introduced during activation.

Preliminary FE testing indicates that the trabeculated myocardium behaves qualitatively like the homogeneous myocardium. Based on this finding, we introduce a new exponential SED function for the passive component of the trabeculated myocardium in the form

Wp = a1 [exp(a2 Exx*Exx + a3 Eyy*Eyy + a4 Ezz*Ezz + 2 a5 Exx*Eyy + 2 a6 Eyy*Ezz + 2 a7 Ezz*Exx) - 1.0]

, and a polynomial SED function for myofibrils

Wa = (I1 - 3)(a8 Exx*Exx + a9 Eyy*Eyy + a10 Ezz*Ezz + 2 a11 Exx*Eyy + 2 a12 Eyy*Ezz + 2 a13 Ezz*Exx)

back to top

Numerial simulations
We perform nonlinear FE analyses to simulate a series of stretch tests on the RVE. Boundary conditions are imposed such that fully tri-axial stress states are generated.

The loading protocol adopted for the numerical simulations: Exx = Eyy > 0, Ezz = 0 Eyy = Ezz > 0, Exx = 0 Ezz = Exx > 0, Eyy = 0
Figure 1. Boundary settings on a RVE. All six faces of the hexahedral RVE are attached to straight rigid plates. These plates are either fixed or move under a prescribed displacement.

Three analyses are performed for each numerical stretch simulation. Each analysis is completed in 10 load steps. Therefore, a total of 30 pairs of averaged stress-strain correspondences are obtained for each simualtion. Two sets of stress-strain data are obtained separately for passive and active myocardium.

Figure 2. Distribution of computed Cauchy stress at local myofibril direction on the deformed mesh (active analysis). The undeformed mesh is shown in outline.

back to top

Parameter estimation
We implement a Levenberg-Marquardt nonlinear regression procedure, with appropriate modifications for multi-dimensional regression, to compute the best-fit parameters in the proposed SED functions by minimizing the combined square error of three normal stress components at the same time. A two-step fitting procedure is used to compute all parameters in both passive and active SED functions. The seven parameters in Wp are determined first based on the passive stress-strain data. Then the six parameters in Wa are determined using active stress-strain data by fitting W = Wp+Wa, while maintaining the already computed seven parameters in Wp constant. We computed the parameter correlation matrix to check if a proposed SED function contains redundant parameters. We also verified the general applicability of the proposed SED functions.

back to top