APM  Vol.8 No.7 , July 2018
Anisotropic Continuum Stored Energy Functional Solved by Lie Group and Differential Geometry
Author(s) Fuzhang Zhao
ABSTRACT
An anisotropic continuum stored energy (CSE), which is essentially composed of invariant component groups (ICGs), is postulated to be balanced with its stress work done, constructing a partial differential equation (PDE). The anisotropic CSE PDE is generally solved by the Lie group and the ICGs through curvatures of elasticity tensor are particularly grouped by differential geometry, representing three general deformations: preferred translational deformations, preferred rotational deformations, and preferred powers of ellipsoidal deformations. The anisotropic CSE constitutive models have been curve-fitted for uniaxial tension tests of rabbit abdominal skins and porcine liver tissues, and biaxial tension and triaxial shear tests of human ventricular myocardial tissues. With the newly defined second invariant component, the anisotropic CSE constitutive models capture the transverse effects in uniaxial tension deformations and the shear coupling effects in triaxial shear deformations.

1. Introduction

Soft biological tissues (SBTs) are recognized as anisotropic hyperelastic materials since they are naturally made of fiber reinforcements and a fluid-like matrix for supporting reversible finite deformations. Rubber-like materials are often treated as isotropic hyperelastic materials. The normalized nominal stress-stretch curves in uniaxial tension for SBTs in a fiber direction and rubber-like materials are shown in Figure 1. In uniaxial tension tests of SBTs, less stiffening at initial

Figure 1. Normalized stress-stretch curves of soft tissue and rubber in uniaxial tension tests.

stretches and more stiffening at final stretches have been observed, and the stress-stretch curve behaves more linearly at both initial and final stretches, making SBTs far more efficient as compared to rubber-like materials. Anisotropic constitutive modeling for finite deformations of SBTs requires experimental characterizations combined with theoretical predictions. As a common practice, theoretical models are fitted with experimental data tested in certain deformation modes and the fitted models predict deformations in untested modes. Thus, the theoretical development of constitutive models, the optimal design of experimental tests, and their numerical implementations into finite element methods are crucial in analyses and designs associated with SBTs.

Kinematics, conservation principles, and constitutive relations are the three pillars of continuum mechanics. Fundamental continuum theories for modeling SBTs as anisotropic hyperelastic materials were presented in the monographs by Fung (1993) [1] , Holzapfel (2000) [2] , and Cowin (2013) [3] . In constitutive relations, the second Piola-Kirchhoff stress tensor is generally related to the right Cauchy-Green tensor through a stored energy density functional or simply a stored energy functional, which needs to be rigorously determined. Thus, the determination of an anisotropic CSE functional is indispensable for constitutive modeling of anisotropic hyperelastic materials, including SBTs and fiber reinforced polymers.

A CSE functional is the key to both anisotropic and isotropic constitutive modelings for hyperelastic materials. For anisotropic constitutive modelings, strain invariant approaches dominate the established anisotropic models for SBTs since material properties must be invariant from all angles of observation. This objectivity ensures frame-indifference. Polynomial, power, exponential, and other functions of invariant components have been used to construct ICGs. Anisotropic constitutive models for SBTs have been reviewed through different perspectives by Humphrey (2003) [4] , Gasser, Ogden, and Holzapfel (2006) [5] , Holzapfel and Ogden (2010) [6] , Chagnon, Rebouah, and Favier (2015) [7] , and Lanir (2017) [8] .

Finite deformations can be measured through the right Cauchy-Green tensor [2]

C = F T F , (1)

where F is the deformation gradient tensor. The three invariants of right Cauchy-Green tensor, I 1 , I 2 , and I 3 , are defined by

I 1 = C : I = tr C , (2)

I 2 = tr ( adj C ) = 1 2 [ ( tr C ) 2 tr C 2 ] = I 3 tr C 1 , (3)

I 3 = det C = 1 3 [ tr C 3 3 2 tr C tr C 2 + 1 2 ( tr C ) 3 ] , (4)

where I is the second-order unit tensor and the matrix operators, “:”, “tr”, “adj”, and “det”, denote double contraction, trace, adjugate, and determinant operations, respectively.

The adjugate, cofactor, and determinant of C are related by the following equations

adj C = ( cof C ) T = det ( C ) C 1 , (5)

where the superscript T denotes the transpose operator of a matrix.

The first-order derivatives of invariants are given by

I 1 C = I , I 2 C = I 1 I C , I 3 C = I 3 C 1 . (6)

The second-order derivatives of invariants are given by

2 I 1 C 2 = O , (7)

2 I 2 C 2 = I I I I , (8)

2 I 3 C 2 = I 3 ( C 1 C 1 C 1 C 1 ) , (9)

where O is the fourth-order zero tensor. The derivatives of second-order symmetric tensors are defined as

C C = I I = ( I I ) i j k l = 1 2 ( δ i k δ j l + δ i l δ j k ) , (10)

in which δ i k , for example, is the Kronecker delta and

C 1 C = C 1 C 1 . (11)

Anisotropic finite deformations can be modeled through invariant components of the right Cauchy-Green tensor and structural tensors. According to Spencer (1971, 1984) [9] [10] , one invariant component in a preferred direction i, I 4, i , can be defined by

I 4 , i = C : A i = C : a 0 , i a 0 , i = tr ( C A i ) , ( i = 1 , 2 , , n ) (12)

where A i is a structural tensor, a 0 , i is a unit vector along a preferred direction i in reference configuration, and n is the total number of different preferred directions. The detailed definitions of structural tensors and the related theory of invariants have been reviewed by Zheng (1994) [11] .

Polyconvexity provides an alternative way to define invariant components. The polyconvex functions, F , adj F , and det F , describe deformations of distance, area, and volume elements since they map the corresponding elements from reference configuration to current configuration, respectively. Thus, the polyconvex functions play an important role in the definition of polyconvex invariant components. The invariant component I 4, i defined in (12) happens to be polyconvex. Following Schröder and Neff (2003) [12] , another polyconvex invariant component in a preferred direction i, I 5, i , can be defined as

I 5 , i = cof C : A i = tr [ cof ( C ) A i ] = tr ( C 1 A i ) I 3 , ( i = 1 , 2 , , n ) . (13)

The polyconvex invariant components, I 4, i and I 5, i , along with the invariant I 3 , have been selected as the arguments of an anisotropic stored energy functional by Itskov, Ehret, and Mavrilas (2006) [13]

Ψ a = Ψ a ( I 4 , i , I 5 , i , I 3 ) , ( i = 1 , 2 , , n ) . (14)

Soft tissues can mainly be characterized as transversely isotropic and orthotropic hyperelastic materials. Transverse isotropy represents a material symmetry with respect to only one principal direction and the structural tensors can be expressed in the reference configuration as

A 1 = a 0 , 1 a 0 , 1 , A 2 = A 3 = 1 2 ( I A 1 ) , (15)

where the principal material direction is denoted by the subscript 1. Orthotropy is characterized by symmetry with respect to three mutually orthogonal planes and the structural tensors read

A i = a 0 , i a 0 , i , ( i = 1 , 2 , 3 ) . (16)

Many influential discoveries in constitutive modeling of SBTs, including the structural tensors by Spencer (1971) [9] , the exponential strain energy function by Fung (1973) [14] , the anisotropic constitutive model by Holzapfel, Gasser, and Ogden (2000) [15] , the extension of polyconvexity to invariant components by Schröder and Neff (2003) [12] , and the selection of invariant components by Itskov, Ehret, and Mavrilas (2006) [13] , have been achieved. Existing stored energy functionals for modeling SBTs are usually assumed to be divided into three parts: a purely volumetric part, a purely isochoric part due to isotropic deformations, and a purely isochoric part due to anisotropic deformations. Anisotropic constitutive models are essentially composed of ICGs. The ICGs of existing models, however, are phenomenologically determined and invariant components are incompletely defined and inconsistently selected. To overcome such shortcomings, a physically relevant and mathematically covariant anisotropic CSE functional is needed for modeling and predicting finite deformations of SBTs.

The CSE functional, originally developed for modeling isotropic polymers by Zhao (2016, 2017) [16] [17] , will be generalized to model anisotropic hyperelastic materials such as SBTs. The main objectives, therefore, are to rigorously formulate the anisotropic CSE functional with the ICGs generally grouped by the Lie group and particularly grouped by differential geometry, to newly define the second invariant component to capture anisotropic transverse effects and shear coupling effects, and to accurately fit anisotropic CSE models with experimental data for SBTs.

2. Anisotropic Constitutive Modeling

2.1. Anisotropic CSE Functional

An anisotropic CSE functional can be conceptually expressed as

Ψ a = Ψ a ( I 1 , i , I 2 , i , I 3 ) , ( i = 1 , 2 , , n ) , (17)

where the range for subscript i, ( i = 1 , 2 , , n ) , is omitted for the rest of the equations for simplification and the number of different preferred directions, n = 3, is used for both orthotropic and transversely isotropic materials. The first invariant component can be expressed as

I 1 , i = tr ( C A i ) , (18)

and the second invariant component is newly defined as

I 2 , i = 1 2 [ I 1 I 1 , i tr ( C 2 A i ) ] . (19)

For convenient subsequent derivations, one can additionally define

I 3 , i = I 3 . (20)

The first-order derivatives of invariant components with respect to the right Cauchy-Green tensor C are given below

I 1 , i C = A i , (21)

I 2 , i C = 1 2 ( I 1 A i + I 1 , i I C A i A i C ) , (22)

I 3 , i C = I 3 C 1 . (23)

The second-order derivatives of invariant components with respect to the right Cauchy-Green tensor C have been obtained as

2 I 1 , i C 2 = O , (24)

2 I 2 , i C 2 = 1 2 ( A i I + I A i A i I I A i ) , (25)

2 I 3, i C 2 = I 3 ( C 1 C 1 C 1 C 1 ) . (26)

For the anisotropic CSE functional (17), along with the invariant components (18) through (20), the second Piola-Kirchhoff stress tensor for anisotropic hyperelastic materials, S i , reads

S i = 2 Ψ a C = 2 j = 1 3 Ψ a I j , i I j , i C . (27)

Substituting the first order derivatives (21), (22), and (23) into (27)2 yields

S i = 2 Ψ a I 1 , i A i + Ψ a I 2 , i ( I 1 A i + I 1 , i I C A i A i C ) + 2 I 3 Ψ a I 3 C 1 , (28)

and the Kirchhoff stress tensor, τ i , can be converted from the second Piola-Kirchhoff stress tensor by the following push-forward operation

τ i = F S i F T . (29)

For isothermal processes, the general CSE functional for anisotropic hyperelastic materials at finite deformations, Ψ a , is postulated to be balanced with its stress work done as

Ψ a = S i : C 2 = τ i : I 2 , (30)

Substituting (28) into (30)1 or (29) into (30)2, simplifying, and rearranging produces the CSE PDE

Ψ a = I 1 , i Ψ a I 1 , i + 2 I 2 , i Ψ a I 2 , i + 3 I 3 Ψ a I 3 . (31)

With Lie group methods, the characteristic system for the PDE is

d I 1, i I 1, i = d I 2, i 2 I 2, i = d I 3 3 I 3 = d Ψ a Ψ a , (32)

and taking its three first-integrals, ψ 1, i = I 2, i / I 1, i 2 , ψ 2 , i = I 3 / I 1 , i 3 , and ψ 3 , i = Ψ a / I 1 , i , the general solution is obtained and written as

Ψ a = I 1 , i [ f 1 ( I 2 , i / I 1 , i 2 ) + f 2 ( I 3 / I 1 , i 3 ) ] , (33)

where f 1 and f 2 are two arbitrary functions and the general solution defines a group of anisotropic CSE functionals.

The general solution (33) has two arbitrary functions to be determined for practical applications. For translational and rotational deformations, the first arbitrary function, f 1 , is selected as

f 1 ( I 2 , i / I 1 , i 2 ) = c 1 , i + c 2 , i I 2 , i / I 1 , i 2 = c 1 , i + c 2 , i I 2 , i / I 1 , i . (34)

For different powers of ellipsoidal deformations, the inverse power function is required. Thus, the second arbitrary function, f 2 , is chosen as

f 2 ( I 3 / I 1 , i 3 ) = c 3 , i ( I 3 / I 1 , i 3 ) c 4 , i / 3 = c 3 , i I 1 , i c 4 , i I 3 c 4 , i / 3 . (35)

Substituting (34) and (35) into (33) gives the particular solution of the anisotropic CSE functional

Ψ a p = c 1 , i I 1 , i + c 2 , i I 2 , i + c 3 , i I 1 , i c 4 , i + 1 I 3 c 4 , i / 3 , (36)

where the four coefficients for a preferred fiber direction i, c 1, i , c 2, i , c 3, i , and c 4, i , are unknown constitutive constants to be determined by experimental tests of anisotropic hyperelastic materials.

2.2. Anisotropic Elasticity Tensors

Elasticity tensors are crucial in studying mathematical properties of the constitutive relations. The anisotropic fourth-order elasticity tensor in material description is defined as

i = 2 S i C . (37)

Substituting (27)2 into (37) and taking derivatives yields the fourth-order anisotropic elasticity tensor in material description

i = 4 j = 1 3 [ k = 1 3 2 Ψ a I j , i I k , i I j , i C I k , i C + Ψ a I j , i 2 I j , i C 2 ] . (38)

Substituting (21) through (26) into (38) yields

i = Δ 1 , i A i A i + Δ 2 , i ( A i I + I A i ) + Δ 3 , i ( A i B i + B i A i ) + Δ 4 , i ( A i C 1 + C 1 A i ) + Δ 5 , i I I + Δ 6 , i ( I B i + B i I ) + Δ 7 , i ( B i B i ) + Δ 8 , i ( I C 1 + C 1 I ) + Δ 9 , i ( C 1 B i + B i C 1 ) + Δ 10 , i C 1 C 1 + Δ 11 , i ( A i I + I A i ) + Δ 12 , i C 1 C 1 , (39)

where B i = 1 2 ( C A i + A i C ) and the twelve parameters are defined as

Δ 1 , i = 4 2 Ψ a I 1 , i 2 + 4 I 1 2 Ψ a I 1 , i I 2 , i + I 1 2 2 Ψ a I 2 , i 2 , (40)

Δ 2 , i = 2 I 1 , i 2 Ψ a I 1 , i I 2 , i + I 1 I 1 , i 2 Ψ a I 2 , i 2 + 2 Ψ a I 2 , i , (41)

Δ 3 , i = 4 2 Ψ a I 1 , i I 2 , i 2 I 1 2 Ψ a I 2 , i 2 , Δ 4 , i = 4 I 3 2 Ψ a I 3 I 1 , i + 2 I 3 I 1 2 Ψ a I 2 , i I 3 , (42)

Δ 5 , i = I 1 , i 2 2 Ψ a I 2 , i 2 , Δ 6 , i = 2 I 1 , i 2 Ψ a I 2 , i 2 , Δ 7 , i = 4 2 Ψ a I 2 , i 2 , (43)

Δ 8 , i = 2 I 3 I 1 , i 2 Ψ a I 2 , i I 3 , Δ 9 , i = 4 I 3 2 Ψ a I 2 , i I 3 , (44)

Δ 10 , i = 4 I 3 2 2 Ψ a I 3 2 + 4 I 3 Ψ a I 3 , Δ 11 , i = 2 Ψ a I 2 , i , Δ 12 , i = 4 I 3 Ψ a I 3 . (45)

The fourth-order anisotropic elasticity tensor in spatial description, similar to the isotropic elasticity tensor [18] , is defined as

i = F F T 2 : i : F T F 2 . (46)

With I 3 = 1 , the anisotropic elasticity tensors can be further simplified.

2.3. CSE Models for Different Deformation Modes

Commonly used experimental tests are uniaxial tension, biaxial tension, and triaxial shear tests. Nominal stress and stretch are usually recorded in experimental tests. With the CSE functional (36), nominal stress as a function of stretch can be established by

( P j k ) i = Ψ a I 1 , i I 1 , i λ j k + Ψ a I 2 , i I 2 , i λ j k + Ψ a I 3 I 3 λ j k , ( j , k = 1 , 2 , 3 ) . (47)

The three derivatives of the CSE functional (36) for incompressible SBTs are given below

Ψ a I 1 , i = c 1 , i + c 3 , i ( c 4 , i + 1 ) I 1 , i c 4 , i , Ψ a I 2 , i = c 2 , i 2 I 2 , i , Ψ a I 3 = 0. (48)

Substituting (48), I 1 , i = I 2 , i = I 3 = 1 , and derivatives of invariant components with respect to normal stretches into (47) yields the incompressible stress-free (ISF) condition in reference configuration

2 c 1 , i + 0.5 c 2 , i + 2 c 3 , i ( c 4 , i + 1 ) = 0. (49)

Thus, the CSE constitutive models in three deformation modes can be derived based on the equations (47) through (49).

2.3.1. Uniaxial Tension Mode

The deformation of uniaxial tension can be modeled as

x i = λ i X i , i = 1 , 2 , 3 , (50)

where X 1 , X 2 , X 3 and x 1 , x 2 , x 3 denote the Cartesian coordinates of a typical particle in reference and current configurations, respectively. With I 3 = 1 ,

the tensors for the uniaxial tension mode in (50) are F = diag [ λ 1 , λ 2 , ( λ 1 λ 2 ) 1 ] and C = diag [ λ 1 2 , λ 2 2 , ( λ 1 λ 2 ) 2 ] . The first principal invariant is obtained as I 1 = λ 1 2 + λ 2 2 + ( λ 1 λ 2 ) 2 . In the longitudinal uniaxial tension for orthotropic materials, the structural tensor is A 1 = diag [ 1 , 0 , 0 ] . The two invariant components are I 1,1 = λ 1 2 and I 2 , 1 = 1 2 ( λ 1 2 λ 2 2 + λ 2 2 ) . Substituting (48), invariant components, and their derivatives into (47) yields the longitudinal nominal stress, P u o 1 , in uniaxial tension tests for orthotropic materials

P u o 1 = 2 c 1 , 1 λ 1 + c 2 , 1 λ 1 λ 2 2 2 ( λ 1 2 λ 2 2 + λ 2 2 ) + 2 c 3 , 1 ( c 4 , 1 + 1 ) λ 1 2 c 4 , 1 + 1 , (51)

and by the same token for A 2 = diag [ 0 , 1 , 0 ] , the transverse nominal stress is written as

P u o 2 = 2 c 1 , 2 λ 2 + c 2 , 2 λ 2 λ 1 2 2 ( λ 1 2 λ 2 2 + λ 1 2 ) + 2 c 3 , 2 ( c 4 , 2 + 1 ) λ 2 2 c 4 , 2 + 1 . (52)

For transversely isotropic materials, the structural tensors are defined in (15). In longitudinal tension tests with cylindrical specimens, A 1 = diag [ 1 , 0 , 0 ] , the tensile stretch is λ 1 with isotropic contraction of λ 2 = λ 3 . With the equation (47), the longitudinal nominal stress is given by

P u t 1 = 2 c 1 , 1 λ 1 + c 2 , 1 2 λ 1 + 2 c 3 , 1 ( c 4 , 1 + 1 ) λ 1 2 c 4 , 1 + 1 , (53)

and in transverse tension tests, A 2 = diag [ 0.0 , 0.5 , 0.5 ] , the tensile stretch is λ 2 with anisotropic contraction and the shortening of fiber λ 1 is not necessarily identical to λ 3 . The transverse nominal stress is

P u t 2 = c 1 , 2 λ 2 + c 2 , 2 λ 2 ( λ 3 2 0.5 λ 2 4 ) 2 λ 2 2 λ 3 2 + λ 2 2 + λ 3 2 + c 3 , 2 ( c 4 , 2 + 1 ) λ 2 ( λ 2 2 + λ 3 2 2 ) c 4 , 2 . (54)

With the ISF condition (49), the CSE models (51), (52), (53), and (54) can be further simplified for subsequent curve fittings.

2.3.2. Biaxial Tension Mode

The deformation of biaxial tension can be generally modeled by

x 1 = λ 11 X 1 + λ 12 X 2 , x 2 = λ 21 X 1 + λ 22 X 2 , x 3 = λ 3 X 3 , (55)

where λ 12 or λ 21 is the amount of shear. The tensors F and C for cases in (55) are given by

F = ( λ 11 λ 12 0 λ 21 λ 22 0 0 0 λ 3 ) , (56)

C = ( λ 11 2 + λ 21 2 λ 11 λ 12 + λ 22 λ 21 0 λ 11 λ 12 + λ 22 λ 21 λ 22 2 + λ 12 2 0 0 0 λ 3 2 ) , (57)

where λ 3 can be determined by incompressible condition

λ 3 = ( λ 11 λ 22 λ 12 λ 21 ) 1 . (58)

For orthotropic materials with A 1 = diag [ 1 , 0 , 0 ] , I 1 , 1 = λ 11 2 + λ 21 2 , I 2 , 1 = 1 2 [ ( λ 11 2 + λ 21 2 ) λ 3 2 + λ 3 2 ] , and related derivatives, the longitudinal and transverse Cauchy normal stresses due to fiber 1 can be obtained as

( σ 11 ) 1 = 2 c 1 , 1 λ 11 2 + c 2 , 1 λ 11 2 λ 3 2 + λ 11 λ 22 λ 3 1 ( λ 11 2 + λ 21 2 ) λ 11 λ 22 λ 3 3 2 [ ( λ 11 2 + λ 21 2 ) λ 3 2 + ( λ 11 λ 22 λ 12 λ 21 ) 2 ] + 2 c 3 , 1 ( c 4 , 1 + 1 ) λ 11 2 ( λ 11 2 + λ 21 2 ) c 4 , 1 , (59)

( σ 22 ) 1 = c 2 , 1 λ 11 λ 22 λ 3 1 ( λ 11 2 + λ 21 2 ) λ 11 λ 22 λ 3 3 2 [ ( λ 11 2 + λ 21 2 ) λ 3 2 + ( λ 11 λ 22 λ 12 λ 21 ) 2 ] , (60)

and for orthotropic materials with A 2 = diag [ 0 , 1 , 0 ] , I 1 , 2 = λ 22 2 + λ 12 2 , I 2 , 2 = 1 2 [ ( λ 22 2 + λ 12 2 ) λ 3 2 + λ 3 2 ] , and related derivatives, similarly, the longitudinal and transverse Cauchy normal stresses due to fiber 2 turns out to be

( σ 11 ) 2 = c 2 , 2 λ 11 λ 22 λ 3 1 ( λ 22 2 + λ 12 2 ) λ 11 λ 22 λ 3 3 2 [ ( λ 22 2 + λ 12 2 ) λ 3 2 + ( λ 11 λ 22 λ 12 λ 21 ) 2 ] , (61)

( σ 22 ) 2 = 2 c 1,2 λ 22 2 + c 2,2 λ 22 2 λ 3 2 + λ 11 λ 22 λ 3 1 ( λ 22 2 + λ 12 2 ) λ 11 λ 22 λ 3 3 2 [ ( λ 22 2 + λ 12 2 ) λ 3 2 + ( λ 11 λ 22 λ 12 λ 21 ) 2 ] + 2 c 3,2 ( c 4,2 + 1 ) λ 22 2 ( λ 22 2 + λ 12 2 ) c 4,2 . (62)

Combining with longitudinal normal stress components ( σ 11 ) 1 and ( σ 11 ) 2 in (59) and (61) and summing up transverse normal stress components ( σ 22 ) 1 and ( σ 22 ) 2 in (60) and (62) yields the general CSE models for biaxial tension tests of orthotropic materials

σ b o 11 = ( σ 11 ) 1 + ( σ 11 ) 2 , (63)

and

σ b o 22 = ( σ 22 ) 1 + ( σ 22 ) 2 . (64)

Without shear deformations, λ 12 = λ 21 = 0 , in-plane normal stretches become in-plane principal stretches and normal stresses become principal stresses. Ignoring both shear deformations and shear coupling effects produces the simplified equations. For the structure tensor A 1 = diag [ 1 , 0 , 0 ] , the simplified longitudinal Cauchy principal stress for equibiaxial tension tests is

σ b o 1 = 2 c 1 , 1 λ 1 2 + c 2 , 1 λ 1 4 2 ( λ 1 4 + λ 1 2 ) + 2 c 3 , 1 ( c 4 , 1 + 1 ) λ 1 2 c 4 , 1 + 2 , (65)

and for the structure tensor A 2 = diag [ 0 , 1 , 0 ] , the simplified transverse Cauchy principal stress is

σ b o 2 = 2 c 1 , 2 λ 2 2 + c 2 , 2 λ 2 4 2 ( λ 2 4 + λ 2 2 ) + 2 c 3 , 2 ( c 4 , 2 + 1 ) λ 2 2 c 4 , 2 + 2 . (66)

With the ISF condition (49), the CSE models (63), (64), (65), and (66) can be further simplified in subsequent curve fitting processes.

2.3.3. Triaxial Shear Mode

The deformation of triaxial shear can be modeled as

x 1 = X 1 + λ 12 X 2 , x 2 = X 2 , x 3 = X 3 , (67)

where the angle of shear is t a n 1 λ 12 . The tensors F , C , and C 2 for cases in (67) are worked out as

F = ( 1 λ 12 0 0 1 0 0 0 1 ) , C = ( 1 λ 12 0 λ 12 λ 12 2 + 1 0 0 0 1 ) , (68)

C 2 = ( λ 12 2 + 1 λ 12 3 + 2 λ 12 0 λ 12 3 + 2 λ 12 λ 12 4 + 3 λ 12 2 + 1 0 0 0 1 ) , (69)

and with the right Cauchy-Green tensor (68)2 the three principal invariants are obtained as I 1 = I 2 = λ 12 2 + 3 and I 3 = 1 . For orthotropic materials, the structural tensor is defined in (16).

For the longitudinal shear, the structural tensor is A 1 = diag [ 1 , 0 , 0 ] . The two invariant components are I 1 , 1 = I 2 , 1 = 1 . With all three zero derivatives in (47), the longitudinal shear stress component ( P 12 ) 1 is

( P 12 ) 1 = 0. (70)

For the perpendicular shear, the structural tensor is A 2 = diag [ 0 , 1 , 0 ] . The two invariant components are I 1 , 2 = λ 12 2 + 1 and I 2 , 2 = λ 12 2 / 2 + 1 . The two derivatives are I 1 , 2 / λ 12 = 2 λ 12 and I 2 , 2 / λ 12 = λ 12 . Substituting the derivatives (48) for i = 2, the two invariant components I 1,2 and I 2,2 , and their derivatives into the stress Equation (47) produces the perpendicular shear stress component ( P 12 ) 2

( P 12 ) 2 = 2 c 1 , 2 λ 12 + c 2 , 2 λ 12 2 λ 12 2 + 4 + 2 c 3 , 2 ( c 4 , 2 + 1 ) λ 12 ( λ 12 2 + 1 ) c 4 , 2 . (71)

For the transverse shear, the structural tensor is A 3 = diag [ 0 , 0 , 1 ] . The two invariant components are I 1 , 3 = 1 , I 2 , 3 = λ 12 2 / 2 + 1 , and their derivatives are I 1 , 3 / λ 12 = 0 and I 2 , 3 / λ 12 = λ 12 . Substituting the derivatives (48) for i = 3, I 1,3 and I 2,3 , and their derivatives into the stress Equation (47) gives the transverse shear stress component ( P 12 ) 3

( P 12 ) 3 = c 2 , 3 λ 12 2 λ 12 2 + 4 . (72)

Combining shear stress components in longitudinal shear (70), perpendicular shear (71), and transverse shear (72) results in the CSE model for triaxial shear tests

P t o 12 = 2 c 1 , 2 λ 12 + ( c 2 , 2 + c 2 , 3 ) λ 12 2 λ 12 2 + 4 + 2 c 3 , 2 ( c 4 , 2 + 1 ) λ 12 ( λ 12 2 + 1 ) c 4 , 2 . (73)

3. Fitting CSE Models to Testing Data

A self-developed graphics digitizer with MATLAB has been used to read out pairs of stress values and corresponding stretches from stress-stretch curves cited in following subsections.

3.1. Modeling for Rabbit Abdominal Skins in Uniaxial Tension

Rabbit abdominal skins are usually treated as incompressible orthotropic hyperelastic materials and they have been tested in constrained uniaxial tension and biaxial tension. In the constrained uniaxial tension tests, the longitudinal tension is conducted under the fixed transverse stretch while the transverse tension is conducted under the fixed longitudinal stretch by Lanir and Fung (1974) [19] . The norminal stress-stretch experimental data in longitudinal and transverse directions has been calculated from the corresponding force-stretch curves and initial areas for sample 36 provided by Tong and Fung (1976) [20] . The uniaxial tension experimental data of rabbit abdominal skins in both longitudinal or 1 and transverse or 2 directions have been used to fit the uniaxial tension models (51) and (52) and the two sets of constitutive constants have been extracted by an iterative least square method. The comparison between the CSE model and the uniaxial tension test data of rabbit abdominal skins is shown in Figure 2.

Figure 2. Comparison between CSE Models and uniaxial tension tests for rabbit abdominal skins.

3.2. Modeling for Porcine Liver Tissues in Uniaxial Tension

Porcine liver tissues, composed of liver lobules and connective tissues, are transversely isotropic with the principal axis along the direction of the lobule. Uniaxial tension and compression tests of porcine liver tissues, using cylindrical specimens, have been conducted by Chui et al. (2007) [21] . Uniaxial tension experimental data averaging over five porcine liver tissue specimens in both longitudinal and transverse directions have been used to fit the uniaxial tension models (53) and (54). Two sets of constitutive constants have been solved by an iterative least square method. The comparison between the anisotropic CSE model and the uniaxial tension test data of porcine liver tissues is shown in Figure 3.

3.3. Modeling for Human Myocardial Tissues in Biaxial Tension

Heart walls consist of three distinct layers: the endocardium, the myocardium, and the epicardium. The ventricular myocardium is the main functional tissue of the heart wall. Biaxial loading is closer to the actual loading condition of passive ventricular myocardial tissues than uniaxial loading. Biaxial tension tests of human passive ventricular myocardial tissues have been conducted by Sommer et al. (2015) [22] . The Cauchy stress-stretch equibiaxial tension experimental data in the mean-fiber direction and the cross-fiber direction has been averaged among 26 specimens, respectively. The averaged equibiaxial tension experimental data have been used to fit the biaxial tension models (65) and (66) and the two

Figure 3. Comparison between CSE models and uniaxial tension tests for porcine liver tissues.

sets of constitutive constants for both directions have been obtained by an iterative least square method. The comparison between the anisotropic CSE models and the averaged equibiaxial tension test data of human passive ventricular myocardial tissues is shown in Figure 4.

3.4. Modeling for Human Myocardial Tissues in Triaxial Shear

Triaxial shear tests for porcine myocardial tissues demonstrate the orthotropic behavior of myocardial tissues by Dokos et al. (2002) [23] . The passive myocardial tissue can therefore be modeled as an orthotropic material with three orthogonal directions: fiber, sheet, and normal. The fiber refers to muscle fiber, the sheet covers a network of collagen fibers, and the normal direction is perpendicular to both muscular and collagen fiber directions [24] . The triaxial shear tests for human myocardial tissues have also been conducted by Sommer et al. (2015) [22] . Human myocardial tissues have also been modeled as orthotropic materials with three orthogonal directions: fiber, sheet, and normal. The data for triaxial perpendicular shear (TP-Shear) to fiber direction are averaged between those of fiber-sheet mode and fiber-normal mode, the data for TP-Shear to sheet direction are averaged between those of sheet-fiber mode and sheet-normal mode, and the data for TP-Shear to normal direction are averaged between those of normal-fiber mode and normal-sheet mode. The triaxial shear experimental data of ventricular myocardial tissues in fiber, sheet, and normal directions have been used to fit the CSE model for triaxial shear (73) and three sets of constitutive constants have been resolved by an iterative least square

Figure 4. Comparison between CSE models and biaxial tension tests for human myocardium.

method. The comparison between the CSE model and the triaxial shear test data of human passive myocardial tissues is shown in Figure 5.

3.5. Constitutive Constants

In anisotropic CSE constitutive models, four constitutive constants are needed to describe anisotropic finite deformations in a preferred direction. Based on the ISF condition (49), one of the four constitutive constants can be treated as a dependent constant. The first constant, c 1, i , is treated as a dependent constant since the last three constants, c 2, i , c 3, i , and c 4, i , are contained in the anisotropic elasticity tensors (39) and (46).

The ISF condition may be applied during extraction of constitutive constants for curve fittings of uniaxial tension and biaxial tension tests. There are two methods for curve fitting experimental data with the ISF condition. In the first method, three unknowns, c 2, i , c 3, i , and c 4, i , are to be solved by the iterative least square method and c 1, i is eliminated from a model by the ISF condition (49). In the second method, four unknowns, c 1, i , c 2, i , c 3, i , and c 4, i , are to be solved by the iterative least square method and the ISF condition (49) is simply subtracted from a model.

The constitutive constants for modeling rabbit abdominal skins (RAS) and porcine liver tissues (PLT) in uniaxial tension, human myocardial tissues (HMT) in biaxial tension and triaxial shear are collected in Table 1. Extraction of constitutive constants for uniaxial tensions in longitudinal (L) and transverse (T) directions, biaxial tensions in the mean-fiber (MF) and the cross-fiber (CF)

Figure 5. Comparison between CSE model and triaxial shear tests for human myocardium.

Table 1. Constitutive constants for different SBTs.

directions was conducted by the first method. Extraction of constitutive constants for triaxial shear tests in fiber (F), sheet (S), and normal (N) directions was conducted directly from the model since the ISF condition is automatically satisfied.

In the triaxial shear model (73), the extracted second constitutive constants, however, are the sums of and. For the three orthogonal directions of fiber, sheet, and normal, they are, , and, respectively. The data for triaxial perpendicular shear of fiber are averaged between those of fiber-sheet mode and fiber-normal mode. As an example for the fiber direction, we have and. Similarly, the two other equations can be established for the sheet and normal directions. Combining equations together for three orthogonal directions yields

(74)

Solving linear Equation (74) simultaneously gives

(75)

Note that the values of, , and are used to plot Figure 5 while the values of for triaxial shear tests listed in Table 1 are the final solution of equation (75).

4. Discussion

4.1. Anisotropic CSE Functional

In the anisotropic CSE functional, and instead of and are named and used albeit since the summation of invariant components or results in the corresponding invariant or

(76)

Moreover, the use of defined in (19) rather than in (13) captures anisotropic transverse deformations and avoids the unnecessary calculation of inverse right Cauchy-Green tensor.

In the anisotropic CSE functional (36), the first term, , represents the work done of normal stress and translational deformation. The second term, , describes the work done of shear stress and rotational deformation. The third term, , captures the work done of stress for finite deformations with different anisotropic degrees of ellipsoidal deformations.

With the fourth-order anisotropic elasticity tensor (38), the most general fourth-order elasticity tensor for isotropic hyperelastic materials can be recovered by the mappings of, , , …, and

(77)

Substituting the first and second order derivatives of invariants (6), (7), (8), and (9) into (77), simplifying, and rearranging produces

(78)

where the eight parameters are defined by

(79)

(80)

(81)

(82)

where the isotropic CSE functional is given by

(83)

The three constitutive constants, , and are generally used for modeling finite deformations of isotropic hyperelastic materials. The two constitutive constants, and, determine the isotropic elasticity tensor (78) with the eight parameters (79) through (82).

4.2. Comments on Different Experimental Tests

In uniaxial tension tests for SBTs, not only does the stress-stretch data in the tensile direction need to be captured but also at least one other orthogonal stretch must be measured for incompressible SBTs since needs to be evaluated for (19) in the CSE model, and transverse deformation effects on longitudinal deformations are indispensable in uniaxial tension tests. The transverse deformation effects in modeling SBTs has been elaborated on by Latorre, Romero, and Montáns (2016) [25] . With uniaxial tension measurements in the direction, for example, we also need to capture the corresponding transverse stretch or and use the incompressible condition to determine the other stretch or, respectively since and are generally not equivalent due to orthotropy. In constitutive modeling of the transverse uniaxial tension test of porcine liver tissues, the tensile stretch is and the shortenings of fibers and cross-fibers are and, respectively. Considering the anisotropic compression test results, the anisotropic contraction is simulated with and the constitutive constants have been obtained and listed in the PLT-T row in Table 1. If the isotropic contraction of were assumed, the constitutive constants would be and in which they are quite different from the PLT-T row of results in Table 1. Thus, it is essential to measure both the longitudinal and transverse stretches in uniaxial tension tests for incompressible SBTs.

In biaxial tension tests, it is very difficult to maintain both uniform force distribution and uniform normal deformations. In general, the stress-stretch curves in biaxial tension tests from cruciform specimens to square specimens are not accurate. For biaxial tension tests of square specimens, the stress as a function of stretch is generally overestimated. The overestimation, the correction factor, and the inverse finite element method regarding biaxial tension tests have been studied by Nolan and McGarry (2016) [26] . In the CSE models for biaxial tension tests, shear coupling effects do exist, making curve fittings harder. Without simplifications, at least four arguments, , , , and, have to be measured, causing material characterizations to be more complex.

In triaxial shear tests, according to different orthogonal fiber reinforcement orientations, shear deformations have been classified as longitudinal shear, perpendicular shear, and transverse shear by Destrade, Horgan, and Murphy (2015) [27] . The three orthogonal directions, fiber, sheet, and normal, for myocardial tissues generate six deformation modes: fiber-sheet mode, fibernormal mode, sheet-normal mode, sheet-fiber mode, normal-fiber mode, and normal-sheet mode. In the CSE model for triaxial shear tests, the longitudinal shear does not contribute any anisotropy, only the perpendicular shear elongates fiber reinforcements, and the perpendicular shear and transverse shear are coupled together. Without shear coupling effects, there would be three different experimental curves in the six deformation modes. The triaxial shear tests in six deformation modes for both porcine myocardial tissues by Dokos et al. (2002) [23] and human myocardial tissues by Sommer et al. (2015) [22] have produced more than three different curves, indicating shear coupling effects as well.

5. Conclusions

An anisotropic CSE functional, for isothermal processes, is postulated to be balanced with its stress work done, constructing a PDE. The anisotropic CSE PDE with respect to the three invariant components, , , and, is generally solved by Lie group methods. A three-term particular solution, which is essentially composed of ICGs, is particularly grouped by differential geometry to capture the three fundamental deformations. In a preferred direction i, the

term monitors translational deformations, the term captures

curvatures of rotational deformations, and the term describes curvatures of different powers of ellipsoidal deformations.

The anisotropic CSE constitutive models for uniaxial tension, biaxial tension, and triaxial shear tests have been derived. The constitutive constants have been solved by an iterative least square method for the uniaxial tension tests of rabbit abdominal skins and porcine liver tissues, the biaxial tension, and the triaxial shear tests of human ventricular myocardial tissues. With the newly defined second invariant component, the anisotropic CSE models capture the transverse effects in uniaxial tension deformations and the coupling effects between the perpendicular shear and transverse shear in triaxial shear deformations.

For anisotropic CSE models, the first constitutive constant, , can be treated as a dependent constant due to the ISF condition if necessary and the last three constitutive constants, , , and, are independent constants, defining the anisotropic elasticity tensor. For isotropic CSE models, the three constitutive constants, , , and, are needed for modeling finite deformations of isotropic hyperelastic materials. The two constitutive constants, and, are required for the isotropic elasticity tensor.

Acknowledgements

The author’s heartfelt thanks go to Jim S.J. Chen and Kenneth B. Margulies for their collaboration in earlier studies of rat hearts through ex-vivo tests of diastolic pressure as a function of balloon volume. He would also like to express his gratitude to Jianming and Jiesi Zhao for their discussion and assistance.

Cite this paper
Zhao, F. (2018) Anisotropic Continuum Stored Energy Functional Solved by Lie Group and Differential Geometry. Advances in Pure Mathematics, 8, 631-651. doi: 10.4236/apm.2018.87037.
References
[1]   Fung, Y.C. (1993) Biomechanics. Springer-Verlag, New York.
https://doi.org/10.1007/978-1-4757-2257-4

[2]   Holzapfel, G.A. (2000) Nonlinear Solid Mechanics: A Continuum Approach for Engineering. John Wiley & Sons, Chichester.

[3]   Cowin, S.C. (2013) Continuum Mechanics of Anisotropic Materials. Springer, New York.
https://doi.org/10.1007/978-1-4614-5025-2

[4]   Humphrey, J.D. (2003) Continuum Biomechanics of Soft Biological Tissues. Proceedings of the Royal Society of London A, 459, 3-46.
https://doi.org/10.1098/rspa.2002.1060

[5]   Gasser, T.C., Ogden, R.W. and Holzapfel, G.A. (2006) Hyperelastic Modelling of Arterial Layers with Distributed Collagen Fibre Orientations. Journal of the Royal Society Interface, 3, 15-35.
https://doi.org/10.1098/rsif.2005.0073

[6]   Holzapfel, G.A. and Ogden, R.W. (2010) Constitutive Modelling of Arteries. Proceedings of the Royal Society A, 466, 1551-1597.
https://doi.org/10.1098/rspa.2010.0058

[7]   Chagnon, G., Rebouah, M. and Favier, D. (2015) Hyperelastic Energy Densities for Soft Biological Tissues: A Review. Journal of Elasticity, 120, 129-160.
https://doi.org/10.1007/s10659-014-9508-z

[8]   Lanir, Y. (2017) Multi-Scale Structural Modeling of Soft Tissues Mechanics and Mechanobiology. Journal of Elasticity, 129, 7-48.
https://doi.org/10.1007/s10659-016-9607-0

[9]   Spencer, A.J.M. (1971) Theory of Invariants. In: Eringen, A.C., Ed., Continuum Physics, Vol. 1, Academic Press, New York, 239-353.
https://doi.org/10.1016/B978-0-12-240801-4.50008-X

[10]   Spencer, A.J.M. (1984) Continuum Theory for Strongly Anisotropic Solids. In: Spencer, A.J.M., Ed., Continuum Theory of the Mechanics of Fibre-Reinforced Composites, CISM Courses and Lectures, No. 282, International Centre for Mechanical Sciences, Springer-Verlag, Wien, 1-32.

[11]   Zheng, Q.-S. (1994) Theory of Representations for Tensor Functions—A Unified Invariant Approach to Constitutive Equations. Applied Mechanics Reviews, 47, 545-587.
https://doi.org/10.1115/1.3111066

[12]   Schroder, J. and Neff, P. (2003) Invariant Formulation of Hyperelastic Transverse Isotropy Based on Polyconvex Free Energy Functions. International Journal of Solids and Structures, 40, 401-445.
https://doi.org/10.1016/S0020-7683(02)00458-4

[13]   Itskov, M., Ehret, A.E. and Mavrilas, D. (2006) A Polyconvex Anisotropic Strain-Energy Function for Soft Collagenous Tissues. Biomechanics and Modeling in Mechanobiology, 5, 17-26.
https://doi.org/10.1007/s10237-005-0006-x

[14]   Fung, Y.C. (1973) Biorheology of Soft Tissues. Biorheology, 10, 139-155.
https://doi.org/10.3233/BIR-1973-10208

[15]   Holzapfel, G.A., Gasser, T.C. and Ogden, R.W. (2000) A New Framework for Arterial Wall Mechanics and a Comparative Study of Material Models. Journal of Elasticity, 61, 1-48.
https://doi.org/10.1023/A:1010835316564

[16]   Zhao, F.Z. (2016) Continuum Constitutive Modeling for Isotropic Hyperelastic Materials. Advances in Pure Mathematics, 6, 571-582.
https://doi.org/10.4236/apm.2016.69046

[17]   Zhao, F.Z. (2017) The Continuum Stored Energy for Constitutive Modeling Finite Deformations of Polymeric Materials. Advances in Pure Mathematics, 7, 597-613.
https://doi.org/10.4236/apm.2017.710036

[18]   Itscov, M. (2000) On the Theory of Fourth-Order Tensors and Their Applications in Computational Mechanics. Computer Methods in Applied Mechanics and Engineering, 189, 419-438.
https://doi.org/10.1016/S0045-7825(99)00472-7

[19]   Lanir, Y. and Fung, Y.C. (1974) Two-Dimensional Mechanical Properties of Rabbit Skin. II. Experimental Results. Journal of Biomechanics, 7, 171-182.
https://doi.org/10.1016/0021-9290(74)90058-X

[20]   Tong, P. and Fung, Y.C. (1976) The Stress-Strain Relationship for the Skin. Journal of Biomechanics, 9, 649-657.
https://doi.org/10.1016/0021-9290(76)90107-X

[21]   Chui, C., Kobayashi, E., Chen, X., Hisada, T. and Sakuma, I. (2007) Transversely Isotropic Properties of Porcine Liver Tissue: Experiments and Constitutive Modelling. Medical & Biological Engineering & Computing, 45, 99-106.
https://doi.org/10.1007/s11517-006-0137-y

[22]   Sommer, G., Schriefl, A.J., Andra, M., Sacherer, M., Viertler, C., Wolinski, H. and Holzapfel, G.A. (2015) Biomechanical Properties and Microstructure of Human Ventricular Myocardium. Acta Biomaterialia, 24, 172-192.
https://doi.org/10.1016/j.actbio.2015.06.031

[23]   Dokos, S., Smaill, B.H., Young, A.A. and LeGrice, I.J. (2002) Shear Properties of Passive Ventricular Myocardium. American Journal of Physiology-Heart and Circulatory Physiology, 283, H2650-H2659.
https://doi.org/10.1152/ajpheart.00111.2002

[24]   Holzapfel, G.A. and Ogden, R.W. (2009) Constitutive Modelling of Passive Myocardium: A Structurally Based Framework for Material Characterization. Philosophical Transactions of the Royal Society A, 367, 3445-3475.
https://doi.org/10.1098/rsta.2009.0091

[25]   Latorre, M., Romero, X. and Montáns, F.J. (2016) The Relevance of Transverse Deformation Effects in Modeling Soft Biological Tissues. International Journal of Solids and Structures, 99, 57-70.
https://doi.org/10.1016/j.ijsolstr.2016.08.006

[26]   Nolan, D.R. and McGarry, J.P. (2016) On the Correct Interpretation of Measured Force and Calculation of Material Stress in Biaxial Tests. Journal of the Mechanical Behavior of Biomedical Materials, 53, 187-199.
https://doi.org/10.1016/j.jmbbm.2015.08.019

[27]   Destrade, M., Horgan, C.O. and Murphy, J.G. (2015) Dominant Negative Poynting Effect in Simple Shearing of Soft Tissues. Journal of Engineering Mathematics, 95, 87-98.
https://doi.org/10.1007/s10665-014-9706-5

 
 
Top