Computer Programs in Biomedicine 7 (1977) 219-231

219

© Elsevier/North-Holland Biomedical Press

FINITE ELEMENT STRESS ANALYSIS OF THE HUMAN LEFT VENTRICLE WHOSE IRREGULAR SHAPE IS DEVELOPED FROM SINGLE PLANE CINEANGIOCARDIOGRAM Dhanjoo N. GHISTA Biomedical Research Division, NASA, Ames Research Center. Moffett Field, Calif. 94035, USA

and M. Sahul HAMID Applied Mathematics Department, Indian Institute of Technology, Hauz Khas, New DelhL India

The three-dimensional left ventricular chamber geometrical model is developed from single plane cineangiocardiogram. This left ventricular model is loaded by an internal pressure monitored by cardiac catheter~ation. The resulting stresses in the left ventricular model chamber's wall are determined by computerized finite element procedure. For the discretization of this left ventricular model structure, a 20-node, isoparametric finite element is employed. The analysis and formulation of the computerised procedure is presented in the paper, along with the detailed algorithms and computer programs, The procedure is applied to determine the stresses in a left ventricle at an instant, during systole. Next, a portion (represented by a finite element) of this left ventricular chamber is simulated as being infarcted by making its active-state modulus value equal to its passive-state value; the neighbouring elements are shown to relieve the 'infarcted' element of stress by themselves taking on more stress. Left Ventricle

Stress

Analysis

Finite Element

1. Introduction The human left ventricular chamber o f the heart caters to and responds to the requirement o f the systemic system. The determination of the in vivo stress distribution in the wall of the human left ventricle (LV) is o f interest for a variety o f reasons: to determine the peak wall stress, to employ the in vivo stress state as an index of compensation o f the left ventricle due to chronic hypertension and dilatation, to determine the stresses in the infarcted wall segments and predict the development o f anneurysms. Previous studies o f stress states o f the left ventricular myocardium have adopted idealized geometries of the LV as spheres and ellipsoids and these analyses have employed thin- and thick-walled theories [ 1 , 2 ] and an elasticity solution [3], with non-homogeneity and anisotropy [4]. Attempts are being made to increasingly incorporate the irregular shape o f the LV by means o f finite element analysis [5,6]. In this study, three-dimensional isoparametric finite

Cineangiocardiogram

elements are employed for discretization o f the left ventricular chamber so as to incorporate the irregularity o f its non-symmetric geometry. This paper presents analyses and results for the stress distributions in the left ventricular model wall (during systole) for: (i) An idealized ellipsoidal left ventricular geometrical model, in order to verify the accuracy o f the 3dimensional finite element model o f the LV (ii) The actual left ventricular chamber geometry (at an instant during systole) generated from cineangiocardiography in anterioposterior projection (iii) The above LV with an infarcted segment (to determine the stresses in the infarcted segments and in their neighbourhood) in order to study the mechanisms o f load sharing of an infarcted LV. The available medical data consists of instantaneous outlines of the LV, obtained from single plane cineangiocardiography, following injection of an opaque contrast medium into the ventricular cavity by transeptal and/or retrogate catheterization. Figure 1 shows the outline o f the LV o f a human subject, trace~ from cine at mid-

220

D.N. Ghista, M.S. HamM, Finite element stress analysis o~ the human left ventricle

Fig. 1. Cine of a human left ventricle at mid-ejection.

ejection; this outline will be utilised for stress analysis. An important aspect of this paper is to provide the analysis, the associated computer solution procedure and the programs, which can be employed (with small laboratory-type computers) for continuing study of the functional mechanisms of the normal pathological and compensated LV(s).

2. Stress analysis of the LV The following considerations are adopted in the stress analysis of the LV:

1. The left ventricular chamber is loaded internally by forces equal to the left ventricular blood pressure monitored by cardiac catheterization. . The deformations of the aortic and mitral rings (due to internal pressure, at the instant of cine) are neglected. 3. The instantaneous Young's modulus of the myocardium employed in the analysis is calculated based on the method of Ghista et al. [7] and equals 2850 g/cm 2 for the instant corresponding to the outline shown in fig. 1.

221

D.N. Ghista, M.S. Hamid, Finite element stress analysis o f the human left ventricle 3A

Z,w

10

I

9

t13

,0% 8 q

--

S6

.

"k'-

20

~=-1

,~ 77~ LOCAL SYSTEM X Y Z GLOBAL SYSTEM ~r Y,v

Fig. 2.20-Node isoparametric element. 20

4. The material is taken to be isotropic, with Poisson's ratio of 0.45 (nearly incompressible).

y= ~ &yi= {Ni}r o, i)

5. The external-wall pressure and also the tribological pressures arising from movements of the LV chamber in the pericardial sac are neglected.

Z = ~ ]Via i = { N i } T ( z i } i=1

6. The cross-sections of the LV normal to anterioposterior project plane of the cine picture are taken to be ellipses; the associated ratio of the major axes to minor axes equals 0.9, based on the studies by Sandier et al. [8], of 1204 observations that A - P and lateral diameters of the left ventricular chambers are closely related. 2.1. Isoparametric e l e m e n t properties

The curved 20-node, 60 degrees of freedom, isoparametric curved element used in this analysis, for left ventricular discretization is shown in fig. 2. The relationship between the global coordinates (x, y, z) and the local curvilinear coordinates (~, r/, ~') is given by

i=l 20

(1)

where N i are the interpolation functions in local coordinates and x i, Yi, zi are the nodal coordinates. The interpolation functions for the 20 node elements are as follows: (i) For corner nodes (~i = +1, r/i = -+1, ~i = ±1) N i = 8!(1 + ~i)(1 + r/r/i)(1 + ~'i) "

"(~i + r/r/i + ~'~'i- 2)

(2)

,(ii) For typical mid-side node (~i = O, r/i = -+1, ~'i = +1) N i = 41(1 - ~2) (1 + r/r/i) (1 + ~'~'i) •

(3)

In the isoparametric case, the same shape functions are used to delineate the displacement field (u, o, w) in terms of the nodal displacements (ui, oi, wi) i.e.,

20

x = ~ Nixi = {Ni)r{xi) i=l

u = (Ni}T{ui}

O= {Ni}T(oi}

w = { N i } T ( w i } (4)

222

D.N. Ghista, M.S. Hamid, Finite element stress analysis o f the human left ventricle

ON/

2.2. Development of the elemental stiffness matrix

0x The strain displacement is written for the element as 6x ey Cz

ON; By

= [J]

Ou/Ox

aN,.

ov/Oy

ow/Oz

0z = [B1, B2 ..... B2o] {3 }e

Ou/Oy + Ov/Ox Ov/Oz + Ow/Oy Ow/Ox + Ou/az

3'xy 7yz

1%x

where from, the derivatives ONi/Ox, ON`./Oy and ONi/Oz are obtained as follows:

(5)

ION`. ~x

wherein

~N`. = [ j ] - I -ON`. ax

0

0

ON/

The stress-strain law for the element is written as (6)

0

0

0z

0N`.

aN`.

0y

0x

0

ON,. Oz

ON`. Oy

0

ONi Ox,

0x

Oy

0

Oz

{o) =

aN`. _ 0N`. ax + ON`.__0y + 0N`. az Ox Ot

a), at

0z 0t

"ryz Tgx

at

Oy at

Ozat

I Ox an

ay On

Oz On

~

ay

Oz

I

Ox

I

i

ONi

a~

aN`. I O~

ax

where D, the elasticity matrix is given by [DI -

E

(10a)

(1 + v ) ( l - 2 v )

-(1 - v)

v

u

0

0

0

(1 - v )

v

0

0

0

(1 - v)

o

o

o

½(1 2v)

0

0

½(1 - 2 v )

0

(7)

The differentials ON`./O~,ONi/O~ and ONi/O~ can hence be written in matrix form as

'ON~

(10)

= [D] (e) "l'xy

In the above equation, terms like ON`./Ox, ONi/Oy etc can be obtained in terms of the local derivatives of N i (ONi/O~ etc) by invoking rules of partial differentiation:

at

aN`.

~z

0y

aNi Oz

(9)

0

0N`. [Bi] :

at ONi

ay

ON,.

0

(8a)

×

Symmetric

°x I 0N~

½(1-2v)

(8)

The stiffness matrix of the element is derived from virtual work or energy principle as

ON,.

[K] e =f[B] T[D] [B] dx dy dz

0z

The global stiffness matrix is expressed in terms of the

(ll)

D.N. Ghista, M.S. Hamid, Finite element stress analysis of the human left ventricle

A

B 1

l

,,3

/

223

- ,-r

l

1 r

r

P4

Fig. 3. Distributed surface-loadpattern.

local coordinates (~, r/, ~') as 1 1 1

[g] e = f f f [B] T[D] [B] detlJI d~ dr/d~"

and (12)

-1-1-1

The evaluation of eq. (12)by closed form integration is difficult because of the complexity of the form of the integrand. Hence numerical integration is resorted to and the stiffness matrix is written in Gaussian quadrature scheme as n [g] e =~

~

n i n j n m ( [ B ] T[D ] [B] )~i, ~j,~rn

X det IJ[

(13)

wherein Hi, Hi and Hm are the weighted coefficients (Zienkiewicz, 1971 [9]).

2.3. Consistent element nodal loads Since the 'lumping' of loads at various nodes is no longer correct for elements of higher order, it is necessary to obtain the joint 'equivalent' forces from virtual work principle. Equating the virtual works done by the nodal loads Pi and by the surface loads of intensity b (fig. 3a), we put down

~lgiPi = f f i

(~u)T (/3) dy dz

(14)

Area

The load intensity/3 and the virtual displacement 8u can be expressed in terms of the nodal values as

/3 = (wi}r (ffi}

On employing pertinent coordinate transformation (dy dz = IJi dr7 d~') and the relation (15), the eq. (14) becomes

{~Ui }T {Pi } = {6ui }T f l f 1 {Ni ) (Ni }Z Pi IJi dr/d~" -1 1

(16)

where IJI is the determinant of the Jacobian matrix. From the above expression, the equivalent (consistent) loads can be written as

n n

i j m

~u = {Ni )T (~ui}

(15)

+1 +1

(Pi)=f f {Ni)(Ni)T(/3i)iJldndf

(17)

-1 -1

For a uniformly distributed load the consistent nodal loads are given in fig. 3b.

2.4. Assembly and solution With the establishment of element stiffness and joint loads, the next step is the assembly of the individual elements by establishing equilibrium of forces at the nodes. By substituting for the nodal forces in terms of the nodal displacements and surface loads/3, we obtain the global stiffness matrix equation. Gaussian elimination procedure is used for its solution, by forming the equations in increasing nodal order and eliminating the rows formed immediately. After the displacements have been determined, the elemental stresses are calculated.

D.N. Ghista,M.S. Hamid, Finite element stress analysis of the human left ventricle

224

3. Computer programs The computer program formulation, employed here, entails compilation of a series of common modules, of a master program with subroutines for element stiffness, formation of overall stiffness, equation solving routines and stress calculations. The subroutines are generalized in order to handle problems of varying sizes.

3.1. Solution of the global stiJjhess equation The assembled global stiffness equation is of the form [K] {q} = {F}

(18)

wherein [K] is the global stiffness (banded) matrix, {q}, is the nodal unknown and {F} is the force acting at the nodes. In the Gaussian elimination procedure, the first equation of eq. (18) is written as

K l l q l +K12qz+K13q3 +...+Klmqm = F 1

(19)

where 'm' is the semibandwidth which is the distance from the diagonal term and the last term in any row. Solving for ql in the above equation, and substituting into the second equation yields (K22 - K21K-l~ K l l)qz

+

(K23 K21KlllK13)q3 +"" -

are eliminated as soon as they are constituted. The equations for a particular node are constituted by incorporating only those coefficients that are associated with that node, starting with the first non-zero coefficient; zero coefficients are not processed in the elimination procedure (thereby minimising the number of operations) and only those elements, which contain the particular node contribute to the associated stiffness coefficients. Hence, it is sufficient to compute only the pertinent rows, in the element stiffness matrices of elements, containing this node. Thus the element stiffness matrix need not be fully computed for each element; ifp be the number of nodes in an element, the storage required is only 1/p-times the full element size. After elimination of all the equations and determination of qn, the next step is the back substitution to determine the remaining nodal unknown displacements. For these operations, we need to store the complete upper diagonal matrix, not written on auxiliary units. The unused core storage is utilized for temporarily storing these equations in a vector form and then written on disk or tape as a block when this vector is full. After solving for nodal displacements, the load vectors for these displacements are calculated from the original equations stored in tape or disc for equilibrium check. Then stresses at various points are calculated from the subroutine (STRESS).

+ (K2m-K21KlllKlm)qm = F 2 - K 2 1 K l l F 1 (20) 4. Flow charts or

K~2q2 + K~3q3 + ... + K~rnqm = F~

(21)

wherein the second unknown q2 is obtained as q2

= / / ' * - - I ~tT*

* --1

'x22 12 - K 2 2

+ K~rnqm )

*

(K23q3 +K~'4q4 +.-(22)

From this operation, one can note that at each line of elimination (in order to formulate the updated set of equations) one needs operations like K21K~d K12 [eq. (22)]. This way the forward elimination process can be continued until the last equation is formed to directly yield the value of qn. Thereafter a back substitution procedure will result in the solution of the vector (q }. The above-described procedure can be efficiently programmed so that the stiffness equations for a node

Main program In this program, the necessary input data is created for nodal coordinates, element characteristics, boundary conditions, and loads. This routine (presented in fig. 4) calls different subroutines for the analytical operations described earlier. The notations employed therein are given below: Variable

Definition

MPOIN NELEM NB NDF

Total number of nodal points Number of elements Number of boundary points Number degrees of freedom per node

225

D.N. Ghista, M.S. Hamid, Finite element stress analysis of the human left ventricle

NBC NCN NCOLN NOD XC U BVV, BNN, NSS

READ NUMBER OF PROBLEMS

I I LOOPON NUMBER OF PROBLEMS

I

I

READAND PRINT INPUT DATA

I

CALL INITIALIZATION

[ NBAND ST

[ LOOPON NUMBER POINTS(NN)

-•

LOOPON NUMBER OF ELEMENTS I

NO

ALT SMTA, SSI, SZTA, H LLREC

I

The various subroutines used in the master program are delineated herewith:

CALCULATEMIN. AND MAX. NODE NUMBERS CALL ELEMENT STIFFNESS

I

--~

Boundary node numbers Number of nodes per element Number of loading cases Element connection array Nodal point coordinate array Load vector Variables used in defining boundary cor~ditions Semi-bandwidth Rectangular matrix for equations of the order (NDF, NBAND) Nodal equations Integration abscissae and weight coefficients of the Gaussian quadrature formula Index in equation formation

4.1. Subroutine initialization (INIT)

CALL FORM EQUATIONS

END LOOPON ELEMENTS I IS PRE.ELIMINATION NECESSARY I ~ YES I CALL PRE-ELIMINATION

This subroutine sets various indices used in the subroutines BSUB, SOLVE, STORE, REABA, etc and the vector X for temporary storing of the modified equations used in back-substitution process. The Flow Chart for INIT is given in fig. 5. 4.2. Subroutine element stiffness (FEM)

I CALL SOLUTION

I END LOOPON POINTS

I

CALL BACK SUBSTITUTION

I

CALL RESIDUALS

I

CALL STRESS

I

I

END LOOPON PROBLEMS

I

[

STOP

I

In the element stiffness matrix, k(i,/) (i = 1, NDF and j = 1 ..... NCN) are the submatrices of order equal to the number of degrees of freedom (NDF) and 'r' is

I

FLOW CHART (MAIN PROGRAM)

I

START

I

SET CONTROL VARIABLES

J

I

I ''CUL T 'E O "OFV TO" I

I

I

I

SUBROUTINE INITIALIZATION (INIT)

Fig. 4. Flow chart (main program).

Fig. 5. Subroutine initialization (INIT).

226

D.N. Ghista, M.S. Hamid, Finite element stress analysis o f the human left ventricle

the number of nodes per element (NCN). Suppose that a particular node point (NN) involved in the formation process is the pth-node in the elemental nodal connection array, it is necessary to calculate only the rows k(p, m). The required indices for this purpose are set in the main program itself and stiffness routine is called. In the computation of the stiffness matrix, in eq. (11), the elements in the matrix (B) corresponding to the required row (p) are chosen and the coefficients are evaluated. The flow chart for a typical 3-dimensional, 20-node isoparametric element for the numerical evaluation of eq. (11) is given in fig. 6 wherein

START

]

I

I

~oo~o~o~

I

I I ~,~o,.~o.~oo,~,o~ I I ~oo ~,~ ~o~,~,~ I

I

~oo~oo~o~o~

I I

I

I SUBROUTINE FORM EQUATIONS (ADDT)

Variable

Definition

D AJAC B AEST

Elasticity matrix Jacobian matrix Strain-displacement matrix Rectangular element stiffness coefficients corresponding to

,p,

START

I

Fig. 7. Subroutine form equation (ADDT).

4. 3. Subroutine form equations (ADDT) ADDT performs the operation of addition of stiff. ness coefficients from the element stiffness matrices that contribute to the nodal point. This subroutine (fig. 7) and the element stiffness routine result in the formation of the global stiffness equations.

]

I I

CALCULATE MATRIX D

START

I

LOOPON INTEGRATING GAUSSIANPOINTS(IG)

I

I

I

.oo~o~..~,.~,o. I

I

GENERATESTRAIN DISPLACEMENT MATRIX B

READ FROMTAPE ELEMENTS FOR ELIMINATION

I

I

I

l

I

I

I

I

SHIFT POSITIONTO STARTING I POINT

CALCULATESTIFFNESSCOEFFICIENTS AEST END LOOPON IG

I

I

I I

I I ~o~oo~o~ ~,,~,~,,o~

I

I

RETURN ] SUBROUTINEELEMENTSTIFFNESS(FEM)

Fig. 6. Subroutine element stiffness (FEM).

SUBROUTINEPRE-ELIMINATION(PRELIM) Fig. 8. Subroutine pre-elimination (PRELIM).

227

D.N. Ghista, M.S. Hamid, Finite element stress analysis of the human left ventricle 4.4. Subroutine pre-elimination (PRELIM)

I START I

I

PRELIM performs the elimination of the columns upto the diagonal point of global stiffness matrix, using the modified equations (of type A.3). The subroutine is given in fig. 8.

[

4.5. Subroutine solution (SOL VE)

WRITESTOREDDATA INTO TAPE

i

I

I

I

SUBROUTINESTOREEQUATIONS(STORE) Fig. 10. Subroutine store equations (STORE).

I LOOP ON DEGREESNDF OF FREEDOM

I

STOREEQUATIONSAND BOUNDARY I VARIABLES IN 'X'

I START I

I

I

I

SOLVE performs the process of line by line elimination of the particular nodal global equations, as per eq. (20), while incorporating the necessary boundary conditions. The subroutine is presented in fig. 9. l

I

I NO

I

START [ NO

I I LOOP ON NUMBER OF EQUATIONS NPZ I

I INCORPORATE PRESCRIBEDVALUE

I

I

I MODIFY CURRENT EQUATION

I

I

STORE EQUATIONS

I

I

ELIMINATE AS INDICATED

SHIFT POSITIONSOF RESULTING ] MATRICES BACK TO STARTING POINT I

END LOOPON NDF

I

I

I I

I

I

YES CALCULATEREACTIONS

1 I

I

READEQUATIONSAND BOUNDARYVARIABLES

II

I

SUBSTITUTEIN COMPUTED DISPLACEMENTVECTOR

I

SHIFTTHE KNOWNDISPLACEMENT MATRIX ONEPOSITIONDOWN

I I RETURN I

SUBROUTINE SOLUTION (SOLVE) Fig. 9. Subroutine solution (SOLVE).

SUBROUTINE BACK-SUBSTITUTION (BSUB)

Fig. 11. Subroutine back-substitution (BSUB).

D.N. Ghista, M.S. HamM, Finite element stress analysis of the human left ventricle

228

4. 7. Subroutine back-substitution (BSUB)

START

I

BSUB calculates the displacement vector by backsubstitution procedure and also computes the reactions at the boundary points. The subroiatine is given in fig.

IS RF~ORD STILL IN X

]NO I

11.

BACK SPACE TAPE AND READ IN NEW BLOCK

I I

4.8. Subroutine read equations and boundary variables (REABA)

RESET INDICES

I I PICK OUT NECESSARY EQUATION AND BOUNDARY VARIABLES

I

I

I

REABA reads back the equations and boundary variables from the vector X for back-substitution process and sets the necessary indices. This subroutine is given in fig. 12. 4. 9. Subroutine residuals (RESID)

SUBROUTINE REABA

RES1D calculates the load vector from the computed displacements for equilibrium check. This subroutine is presented in fig. 13.

Fig. 12. Subroutine REABA.

4. 6. Subroutine store equations (STORE) STORE stores the equations and the boundary values in the buffer area (vector X) in sequential order and adjust the indices. This subroutine is given in fig. 10,

I

-t i

I

I --I LOOPON NOMBEROF PO,NTS,NN' J

I

I LOOPON ELEMENTS

I START

LOOPON POINTS(NN)

I

IS 'NN' A NODAL NUMBER IN THE ELEMENT I YES CALCULATE B

I EINDJol ASPEREGN, ACCUMULATESTRESS

]

END LOOPON ELEMENTS

]

FIND AVERAGESTRESS

I

I

READ FROM TAPE ORIGINAL STIFFNESS EQUATIONS

I

I

I

PERFORM [K]lql NN = IF}NN

I

I

I

CALCULATE PRINCIPALSTRESS I

PRINT LOAD

I

END LOOPON POINTS PRINT PRINCIPA'LSTRESS

I

I RETURNI SUBROUTINE RESIDUALS(RESID) Fig. 13. Subroutine residuals (RESID).

SUBROUTINESTRESS Fig. 14. Subroutine stress.

t ~ [

D.N. Ghista, M.S. Hamid, Finite element stress analysis of the human left ventricle 4.10. Subroutine stress The stresses (o) at different nodal points are calculated from the strains {o} and elasticity matrix [D] as follows: (o) = [O] {e} = [O] [B] {~i}e where (5 }e is the nodal computed displacement vector for an element. The subroutine is presented in fig. 14.

5. Applications of the analysis The analysis, outlined above, is now applied to determine the (instantaneous) stress distributions in various left ventricular geometry simulating models. The results and their implications are presented in this following section.

5.1. Problems, solutions, results Problem I

229

The instantaneous left ventricular chamber geometry is idealized as an ellipsoid of revolution, whose dimensions are obtained by having its volume match that of the left ventricular chamber [8,10]. This problem is solved to check the accuracy of the 3-dimensional element by comparing the results (nature and variation of stresses) with that of the elasticity model of Ghista et al. [3]. The results are" shown in fig. 15. Problem 2 Stress analysis o f the irregular (asymmetric) left ventricular chambers The chamber geometry (at an instant) is developed from the AP-projection cine film. The left ventricutar chamber is obtained by constructing ellipses, normal to the AP-project outline, whose major axes equal the lengths of the corresponding segments in the AP-plane and whose minor axes lengths are 0.9-times the major axes length, based on anatomical observations of the left ventricular contour (Sand-

Stress analysis o f the idealized left ventricular ellipsoidal model LEFT VENTRICULAR FINITE ELEMENT DISCRETIZATION MODEL

------

FINE

O RADIAL

ELEMENT

ANALYTICAL

0

HOOP

I I I I I

3.0--

2.4

[3 M E R I D I O N A L

D

2.0--

1.5

--

1.0

--

P1 ~ 0

--

:1 -1.0

-I

Fig. 15. Comparison of the wall stress distribution results in an ellipsoidal model of the LV of the finite element and elasticity analyses.

2S5

APEX ELEMENT

Fig. 16. Left ventricular finite element discretization.

230

D.N. Ghista, M.S. Hamid, Finite element stress analysis o f the human left ventricle

ler et al. [8] ). The left ventricular chamber's reconstruction from the A cine film (in fig. 1, taken at mid-ejection) and the discretization is demonstrated in fig. 16. The chamber is discretized into 35 elements, with 768 degrees of freedom, resulting in a semi-bandwidth of 123. The solution of the global stiffness matrix is obtained by the Gaussian elimination procedure; its formulation is given in Section 3. The results, presented in fig. 17, depict the distribution of maximum principal stresses at the nodal points. Problem 3 Stress analysis of above 'constructed'

DISTRIBUTION OF MAXIMUM PRINCIPAL STRESS IN A PARTIALLY INFARCTED LV

o., zl~ ,,,

i

=3 ~2

,.8,

0.9

0

left ventricular chamber, with infarcted wall portion designated by element number 16 (fig. 18) For the passive infarcted portion, the Young's modulus is determined to be 2850 g/cm 2 the effective modulus of the LV at diastole. Figure 18 depicts the distributions of maximum tensile principal stresses at the nodal joints.

DISTRIBUTION OF MAXIMUM PRINCIPAL STRESS IN NORMAL LV WALL

A

Fig. 18. Distribution of m a x i m u m principal stresses in the wall of a partially infarcted LV.

(0.6)

(0.5)

0.8

/

[1.2)

13 ~ Io81 ~.

I0.05~

'roT'/ ~

/i

(8,51m ~ - , ~

(0.3)

0.4 ~"

~i. (1.3,

0"~ll,/"

I

T 03

(0.40)L 050

I (0.55) 0.5 Fig. 17. Distribution of m a x i m u m principal stresses in the left ventricular wall.

I

I

+,1-2;, I

1.4

A

,'"

(o7~

/

,'~"

/

(0.4) 0.6

/

/ _ /

/

i

~__~

//

(0.2) 0.5

1.0

,,

.... ( ) = INFARCTED

Fig. 19. Comparison of the m a x i m u m shear stress distributions in normal and infarcted elements.

D.N. Ghista, M.S. Hamid, Finite element stress analysis o f the human left ventricle

The magnitude of the maximum shear stress at each node equals half the difference of the maximum and minimum principal stresses. Figure 19 shows the maximum shear stresses at the junction of normal and infarcted areas. 5.2. Discussion o f problems' solutions

It is seen, from fig. 15, that the finite element analysis of the ellipsoid chamber yields stress distribution (for radial, meridonial, and hoop stresses) that match those obtained analytically by Ghista et al. [3]. This gives us a measure of confidence into the capability of the 20-node, 3-dimensional isoparametric element (fig. 2) to represent the thick-walled chamber for evaluating its stress response. For the finite element analysis of the irregular shaped LV (Problem 2, mentioned in the previous section) the distribution of the maximum principal stresses (in the AP-projection plane) provides the peak-stress concentration levels. The level of peak stress as well as the overall stress distribution, when compared with the results of relatively idealised geometry analyses [ 1,3, 5] indicates that the incorporation of asymmetricity of the LV-chamber geometry does make a significant differ-

231

ence in aftbrding us a more realistic representation of the left ventricular stress state. It is significant to note that whereas the stress within the infarcted element (fig. 18) is less than obtained for the same element, the 'normal' ventricle (fig. 17), the stress around the element is higher. This suggests that the neighbouring myocardium assists the belaboring wall element by partaking of its load. References [1] I. Mirsky, Bull. Math. Biophys. 32 (1971) 197. [2] A.Y.K. Wong and P.M. Rautaharju, Am. Heart J. 75 (1968) 649. [3] D.N. Ghista and H. Sandier, J. Biomech. 2 (1969) 35. [4] I. Mirsky, Bull. Math. Biophys. 32 (1971) 197. [5] P. Gould, D.N. Ghista, L. Brombolich and I. Mirsky, J. Biomecti. 4 (1972) 521. [6] Y.C. Pao, E.L. Ritman and E.H. Wood, J. Biomech. 6 (1974) 469. [7] D.N. Ghista, H. Sandier and W.H. Vayo, Med. Biol. Eng. 13 (1975) 151. [8] H. Sandier, H.T. Dodge, Am. Heart J. 75 (1968) 649. [9] O.C. Zienkiewicz, The finite element method in engineering science (McGraw-Hill, 1971) p. 147. [10] J.E. Davila and M.E. Summereo, Am. J. Cardiol. 18 (1966) 31.