Solid Dynamic Models for Analysis of Stress and Strain in Human Hearts

Abstract

This paper proposes a solid model based on four-dimensional trivariate B-spline for strain and stress analysis of ventricular myocardium. With a series of processing steps in the four-dimensional medical images, the feature points of ventricular inner and outer wall are obtained. A B-spline surface is then used to build the dynamic deformation model of the myocardial walls. With such a surface model, a hexahedron control mesh can be constructed by sweeping the cloud data, and the ventricular solid model is built by fitting the trivariate B-spline parameters. Based on these models, a method of isogeometric analysis can be applied to calculate the stress and strain continuously distributed in the ventricle. The model is represented smoothly in the cylindrical coordinate system and is easy to measure myocardium dynamics for finding abnormal motion. Experiments are carried out for comparing the stress and strain distribution. It is found that the solid model can determine ventricular dynamics which can well reflect the deformation distribution in the heart and imply early clues of cardiac diseases.

1. Introduction

Cardiovascular diseases are currently the leading cause of death in the world, and the rate of death is increasing each year in many countries. Therefore, more and more physiologists and researchers make efforts to understand how the heart works and how to diagnosis and treat the heart diseases. Dynamics and kinetics of the left cardiac ventricle are the primary representation of the cardiac motion, as a series of systolic and diastolic motions of the left cardiac ventricle make the heart pump blood to circulate the whole body. Furthermore, the stress and strain express the characteristics of elasticity and motion of myocardial walls [1].

With assist of computer medical imaging technology, such as magnetic resonance imaging (MRI), computed tomography (CT), positron emission tomography, single photon emission computed tomography, ultrasound, and X-ray, these imaging techniques can give different representation of the heart. In this paper, the endocardium and epicardium shapes of a ventricle are constructed with a segmentation technique in MRI images. For constructing model, finite element (FE) models are the most commonly used in biomechanics. Especially, hexahedral and tetrahedral FE models are popular in the representation of cardiac model, for example, Figure 1 [2, 3]. Comparing with the FE model, left ventricular CAD (computer-aided design) model is more representative of the true ventricle in the shape, such as B-spline surface model. Its continuity and smoothness are better than the FE model. However, it is difficult to apply FE analysis directly to the CAD model [4].

Figure 1 

Finite element models of left ventricle in the literature [2, 3].

The appearance of isogeometric analysis makes it possible to carry out mechanical analysis directly by a CAD model. As a new computational technique, isogeometric analysis improves on and generalizes the standard FE method and has been proven to be a powerful method exceeding the FE method [5]. Inspired by isogeometric analysis, this paper constructs a solid model based on trivariate B-spline and calculates the strain and stress of a ventricle by this model. Finally, it shows the distribution of the strain and stress on the solid model.

2. Model Representation

2.1. B-Spline Curves and Surfaces

One can get a B-spline curve of order 𝑘. By [6]𝐶(𝑡)=^𝑛_𝑖=0𝑝_𝑖𝑁_𝑖,𝑘(𝑡),𝑡_𝑘≤𝑡<𝑡_𝑘+1, (1) where 𝑝₀,𝑝₁,…,𝑝_𝑛 are the control points of 𝐶(𝑡), the linear interpolation of control points is called control polygon. 𝑁_𝑖,𝑘(𝑢) is B-spline basis function of degree 𝑘, which is determined by a group of nondecreasing normalized sequence 𝑇: 𝑡₀≤𝑡₁≤⋯≤𝑡_𝑛+𝑘+1. The sequence can be determined by Riesenfeld method and Hartley-Judd method [7].

The definition of 𝑁_𝑖,𝑘(𝑢) can be expressed as𝑁_𝑖,0(𝑡)=1,if𝑡_𝑖≤𝑡<𝑡_𝑖+1,𝑁0,else_𝑖,𝑘(𝑡)=𝑡−𝑡_𝑖𝑡_𝑖+𝑘−𝑡_𝑖𝑁_{𝑖,𝑘−1}(𝑡)+𝑡_𝑖+𝑘+1−𝑡𝑡_𝑖+𝑘+1−𝑡_𝑖+1𝑁_{𝑖+1,𝑘−1}(𝑡)0define0=0., (2) With a tensor product, given a (𝑛+1)×(𝑚+1) control mesh 𝑝_𝑖,𝑗 and knot vector 𝑇: 𝑡₀≤𝑡₁≤⋯≤𝑡_𝑛+𝑘+1, a B-spline surface of degree 𝑘₁ in the 𝑡 direction and degree 𝑘₂ in the 𝑤 direction is defined as𝑆(𝑡,𝑤)=^𝑛𝑚_𝑖=0𝑗=0𝑝_𝑖,𝑗𝑁_𝑖,𝑘1(𝑡)𝑁_𝑗,𝑘2(𝑤), (3) where 𝑆 is a point on the defined surface, and 𝑡 and 𝑤 are usually representing longitude and latitude; respectively, 𝑁_𝑖,𝑘1(𝑡) and 𝑁_𝑖,𝑘1(𝑤) are the basis functions in the 𝑡 and 𝑤 directions with a degree of 𝑘₁ and 𝑘₂.

2.2. Volumetric B-Splines

A volumetric B-spline can also be defined by𝑉(𝑡,𝑤,𝑢)=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘𝑝_{𝑖,𝑗,𝑙}𝑁_𝑖,𝑘1(𝑡)𝑁_𝑗,𝑘2(𝑤)𝑁_𝑘,𝑘3(𝑢). (4) Comparing with the model in Cartesian coordinates, a cylindrical coordinate model is better approximate to the shape of the heart [8]. Volumetric B-spline is the extension of B-spline surface, which can be given by the following two steps in cylindrical coordinates.

(1) Coordinate transform. A (n+1)×(m+1)×(l+1) control mesh 𝑝_{𝑖,𝑗,𝑙} with the points (𝑥_𝑖𝑗𝑘,𝑦_𝑖𝑗𝑘,𝑧_𝑖𝑗𝑘) in Cartesian coordinate is transformed into cylindrical coordinate by𝑟_𝑖𝑗𝑘=𝑥²_𝑖𝑗𝑘+𝑦²_𝑖𝑗𝑘,𝜃_𝑖𝑗𝑘=𝜋/2if𝑥_𝑖𝑗𝑘=0,𝑦_𝑖𝑗𝑘>0,3𝜋/2if𝑥_𝑖𝑗𝑘=0,𝑦_𝑖𝑗𝑘𝑦<0,𝑎tan^𝑖𝑗𝑘𝑥_𝑖𝑗𝑘if𝑥_𝑖𝑗𝑘>0,𝑦_𝑖𝑗𝑘𝑦≥0,𝜋+𝑎tan^𝑖𝑗𝑘𝑥_𝑖𝑗𝑘if𝑥_𝑖𝑗𝑘𝑦<0,2𝜋+𝑎tan^𝑖𝑗𝑘𝑥_𝑖𝑗𝑘if𝑥_𝑖𝑗𝑘>0,𝑦_𝑖𝑗𝑘𝑧<0,_𝑖𝑗𝑘=𝑧_𝑖𝑗𝑘. (5) (2) The representation of volumetric B-spline in cylindrical coordinate is𝑟(𝑡,𝑤,𝑢)=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘𝑟_{𝑖,𝑗,𝑘}𝑁_𝑖,𝑘1(𝑡)𝑁_𝑗,𝑘2(𝑤)𝑁_𝑘,𝑘3(𝑢),𝜃(𝑡,𝑤,𝑢)=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘𝜃_{𝑖,𝑗,𝑘}𝑁_𝑖,𝑘1(𝑡)𝑁_𝑗,𝑘2(𝑤)𝑁_𝑘,𝑘3(𝑢),𝑧(𝑡,𝑤,𝑢)=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘𝑧_{𝑖,𝑗,𝑘}𝑁_𝑖,𝑘1(𝑡)𝑁_𝑗,𝑘2(𝑤)𝑁_𝑘,𝑘3(𝑢). (6)

3. Determination of Stress and Strain

3.1. Elasticity

Stress and strain actually reflect elasticity of ventricular myocardial walls [9]. Calculations of the stress and strain can be inspired from the idea of elasticity theory. In 3D space, set u, v, w as the displacements in x, y, z direction. Strain vector and stress vector are 𝜎=[𝜎_𝑥,𝜎_𝑦,𝜎_𝑧,𝜏_𝑥𝑦,𝜏_𝑦𝑧,𝜏_𝑥𝑧]^𝑇 and 𝜀=[𝜀_𝑥,𝜀_𝑦,𝜀_𝑧,𝛾_𝑥𝑦,𝛾_𝑦𝑧,𝛾_𝑥𝑧]^𝑇, respectively, as shown in Figure 2.

Figure 2 

Three-dimensional stress and strain components.

The relation between strain and displacement is𝜀𝜀=^𝑥=𝜕𝑢,𝜀𝜕𝑥_𝑦=𝜕𝑣,𝜀𝜕𝑦_𝑧=𝜕𝑤,𝛾𝜕𝑧_𝑥𝑦=𝜕𝑢+𝜕𝑦𝜕𝑣,𝛾𝜕𝑥_𝑦𝑧=𝜕𝑣+𝜕𝑧𝜕𝑤,𝛾𝜕𝑦_𝑥𝑧=𝜕𝑢+𝜕𝑧𝜕𝑤.𝜕𝑥 (7) Physical relationship or the relation between strain and stress is𝜎=𝐷𝜀, (8) where 𝐷 is the elastic matrix defined as𝐸𝐷=1−𝜇×1𝜇1+𝜇1−2𝜇𝜇1−𝜇𝜇1−𝜇0001𝜇1−𝜇𝜇1−𝜇000𝜇1−𝜇1000000000𝐻000𝐻000𝐻1−𝜇000, (9) where 𝐸 is elasticity modulus, 𝜇 is Poisson’s ratio, and 𝐻=(1−2𝜇)/2(1−𝜇). With (9), strain and stress can be determined.

3.2. Isogeometric Analysis

We regard B-spline basis functions as the displacement function and also the basis function of the ventricular shape, to calculate the stress and strain. Solid models at adjacent time points can be reconstructed by (4) with corresponding control points. The continuous displacements of the solid model correspond the displacements of control points, that is, Disp𝑡_𝑐,𝑡_𝑙,𝑡_𝑟=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘=0𝑑_{𝑖,𝑗,𝑙}𝑁_𝑖,𝑘1𝑡_𝑐𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑘,𝑘3𝑡_𝑟, (10) where 𝑑_{𝑖,𝑗,𝑙} is the displacements of control points of Disp(𝑡_𝑐,𝑡_𝑙,𝑡_𝑟). 𝑛,𝑚,𝑙 define the control net, 𝑡_𝑐, 𝑡_𝑙, and 𝑡_𝑟 are the knot vectors in the directions of circumference, long axis and radius, respectively. For X, Y, Z directions, the displacements DX, DY, DZ are defined as𝐷𝑋𝑡_𝑐,𝑡_𝑙,𝑡_𝑟=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘=0𝑑^𝑥_{𝑖,𝑗,𝑙}𝑁_𝑖,𝑘1𝑡_𝑐𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑘,𝑘3𝑡_𝑟,𝐷𝑌𝑡_𝑐,𝑡_𝑙,𝑡_𝑟=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘=0𝑑^𝑦_{𝑖,𝑗,𝑙}𝑁_𝑖,𝑘1𝑡_𝑐𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑘,𝑘3𝑡_𝑟,𝐷𝑍𝑡_𝑐,𝑡_𝑙,𝑡_𝑟=^𝑛𝑚_𝑖=0^𝑙_𝑗=0𝑘=0𝑑^𝑧_{𝑖,𝑗,𝑙}𝑁_𝑖,𝑘1𝑡_𝑐𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑘,𝑘3𝑡_𝑟, (11) where 𝑑_{𝑥𝑖,𝑗,𝑙}, 𝑑_{𝑦𝑖,𝑗,𝑙}, 𝑑_{𝑧𝑖,𝑗,𝑙} are the displacements of control points in X, Y, Z directions.

By elasticity theory, the strain of one point in the model can be calculated𝜀_𝑥=𝜕𝑢=𝜕𝑥𝜕𝐷𝑋=𝜕𝑥𝜕𝑁𝑑𝜕𝑥^𝑥,𝜀_𝑦=𝜕𝑣=𝜕𝑦𝜕𝐷𝑌=𝜕𝑦𝜕𝑁𝑑𝜕𝑦^𝑦,𝜀_𝑧=𝜕𝑤=𝜕𝑧𝜕𝐷𝑍=𝜕𝑧𝜕𝑁𝑑𝜕𝑧^𝑧,𝛾_𝑥𝑦=𝜕𝑢+𝜕𝑦𝜕𝑣=𝜕𝑥𝜕𝐷𝑋+𝜕𝑦𝜕𝐷𝑌=𝜕𝑥𝜕𝑁𝑑𝜕𝑦^𝑥+𝜕𝑁𝑑𝜕𝑥^𝑦,𝛾_𝑦𝑧=𝜕𝑣+𝜕𝑧𝜕𝑤=𝜕𝑦𝜕𝐷𝑌+𝜕𝑧𝜕𝐷𝑍=𝜕𝑦𝜕𝑁𝑑𝜕𝑧^𝑦+𝜕𝑁𝑑𝜕𝑦^𝑧,𝛾_𝑥𝑧=𝜕𝑢+𝜕𝑧𝜕𝑤=𝜕𝑥𝜕𝐷𝑋+𝜕𝑧𝜕𝐷𝑍=𝜕𝑥𝜕𝑁𝑑𝜕𝑧^𝑥+𝜕𝑁𝑑𝜕𝑥^𝑧, (12) where 𝑁=𝑁_𝑖,𝑘1(𝑡_𝑐)𝑁_𝑗,𝑘2(𝑡_𝑙)𝑁_𝑘,𝑘3(𝑡_𝑟) is the B-spline basis function, and dx, dy, and dz are displacement values of the corresponding points.

Then the derivatives of B-spline basis functions in X, Y, Z directions can be derived in each parameter direction𝜕𝑁𝜕𝑡_𝑐=𝑘₁𝑁−1^{𝑖,𝑘1−1}𝑡_𝑐𝑡_{𝑐,𝑖+𝑘1−1}−𝑡_𝑐,𝑖−𝑁_{𝑖+1,𝑘1−1}𝑡_𝑐𝑡_{𝑐,𝑖+𝑘1}−𝑡_𝑐,𝑖+1×𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑘,𝑘3𝑡_𝑟,𝜕𝑁𝜕𝑡_𝑙=𝑘₂𝑁−1^{𝑗,𝑘2−1}𝑡_𝑙𝑡_{𝑙,𝑗+𝑘2−1}−𝑡_𝑙,𝑗−𝑁_{𝑗+1,𝑘2−1}𝑡_𝑙𝑡_{𝑙,𝑗+𝑘2}−𝑡_𝑙,𝑗+1×𝑁_𝑖,𝑘1𝑡_𝑐𝑁_𝑘,𝑘3𝑡_𝑟,𝜕𝑁𝜕𝑡_𝑟=𝑘₃𝑁−1^{𝑘,𝑘3−1}𝑡_𝑟𝑡_{𝑟,𝑘+𝑘3−1}−𝑡_𝑟,𝑖−𝑁_{𝑘+1,𝑘3−1}𝑡_𝑟𝑡_{𝑟,𝑘+𝑘3}−𝑡_𝑟,𝑘+1×𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑖,𝑘1𝑡_𝑟. (13) To get the derivatives of B-spline basis functions in X, Y, Z directions, the partial derivative transform of isoparametric principles [10] is used𝜕𝑁𝜕𝑥𝜕𝑁𝜕𝑦𝜕𝑁𝜕𝑧=𝐽⁻¹𝜕𝑁𝜕𝑡_𝑐𝜕𝑁𝜕𝑡_𝑙𝜕𝑁𝜕𝑡_𝑟, (14) where 𝐽 is Jacobian matrix as follows:𝐽=𝜕𝑥𝜕𝑡_𝑐𝜕𝑥𝜕𝑡_𝑙𝜕𝑥𝜕𝑡_𝑟𝜕𝑦𝜕𝑡_𝑐𝜕𝑦𝜕𝑡_𝑙𝜕𝑦𝜕𝑡_𝑟𝜕𝑧𝜕𝑡_𝑐𝜕𝑧𝜕𝑡_𝑙𝜕𝑧𝜕𝑡_𝑟. (15) Here, take 𝜕𝑥/𝜕𝑡_𝑐 as an example. 𝑝^𝑥_{𝑖,𝑗,𝑙} is the control point coordinate in the 𝑋 direction.𝜕𝑥𝜕𝑡_𝑐=^𝑚𝑙_𝑖=0𝑗=0𝑁_𝑗,𝑘2𝑡_𝑙𝑁_𝑘,𝑘3𝑡_𝑟^𝑛_𝑖=0𝑁_{𝑖,𝑘1−1}𝑡_𝑐𝑝^𝑥_{𝑖,𝑗,𝑙}−𝑝^𝑥_{𝑖−1,𝑗,𝑙}𝑡_{𝑐,𝑖+𝑘1−1}−𝑡_𝑐,𝑖. (16) Therefore, the strain can be calculated in this way. By (8), the stress condition can also be calculated.

4. Experiments and Results

4.1. Construction of Ventricular Models

Both inside and outside data points of ventricular myocardial walls can be obtained from 3D medical images [11, 12]. Figure 3(a) shows the points obtained by model-based segmentation [13]. Figure 3(b) is the corresponding B-spline surface model.

(a)

(b)

Figure 3 

(a) Ventricular point distribution model (b) B-spline surface model.

A sweeping method [14] can be used to obtain the control hexahedral, for example, the one in Figure 4, where Figure 4(a) is the quadrilateral mesh of inside wall, and Figure 4(b) is the hexahedral mesh. While control hexahedral mesh is obtained, the corresponding ventricular B-spline solid model can be reconstructed by (4) and (5).

(a)

(b)

Figure 4 

(a) Left ventricular inside wall quadrilateral mesh, (b) hexahedral control mesh.

4.2. Stress and Strain

Stress and strain of a ventricle are calculated based on the steps described in the above sections. Here, for calculation of the stress results, we set myocardial elastic modulus 11 Kpa and Poisson’s ratio of 0.49. Figure 5 shows the principal strain distribution in 𝑋 direction, and Figure 6 shows the principal stress distribution. In these figures, the left one shows the result by the finite element model, and the right shows that by the B-spline solid model. The color represents the change in stress or strain distribution, with the specific reference of the color bar. From the results, it can be seen that the stress and strain of the left ventricular model show overall nonuniform distribution [15, 16], which is consistent with the results by other researchers [17].

(a)

(b)

Figure 5 

Principal strain distribution in X direction: (a) the results by finite element model, (b) result by B-spline solid model.

(a)

(b)

Figure 6 

Principal stress distribution in X direction: (a) the results by finite element model, (b) result by B-spline solid model.

5. Conclusion

With the situation that traditional finite element methods are difficult for direct use in mechanics analysis which has continuous distribution in space, this paper proposes a solid B-spline model to construct a continuous ventricular mechanical model and applies isogeometric analysis. Stress and strain calculative formulas are derived. The proposed model features continuous, smooth, and inseparable. According to a set of ventricular hexahedral solid B-spline models sampled at different times in a cardiac cycle, strain and stress are determined for medical analysis.

Acknowledgments

Thanks are due to Z. Xu, T. T. Jiang, Y. Q. Xu, Y. H. Du, et al. for their helpful work and laboratory implementation in the project during their study in the university. This work was supported in part by the National Natural Science Foundation of China (NSFC nos. 61103140 and 61105073) and Zhejiang Provincial S&T Department (2010R10006 and 2010C33095).