Anisotropic hyperelastic material behavior

Hyperelastic materials are described in terms of a “strain energy potential,” which defines the strain energy stored in the material per unit of reference volume (volume in the initial configuration) as a function of the deformation at that point in the material. Two distinct formulations, strain-based and invariant-based, are used for the representation of the strain energy potential of anisotropic hyperelastic materials.

The following topics are discussed:

Strain-based formulation
Invariant-based formulation
Particular forms of the strain energy potential

Strain-based formulation

In this case the strain energy function is expressed directly in terms of the components of a suitable strain tensor, such as the Green strain tensor (see Strain measures):

U = U (ε^{G}),

where $ε^{G} = \frac{1}{2} (C - I)$ is Green's strain, $C = F^{T} \cdot F$ is the right Cauchy-Green strain tensor, $F$ is the deformation gradient, and $I$ is the identity matrix. Without loss of generality, the strain energy function can be written in the form

\begin{matrix} (1) & U = U ({\bar{ε}}^{G}, J), \end{matrix}

where ${\bar{ε}}^{G} = \frac{1}{2} (\bar{C} - I)$ is the modified Green strain tensor, $\bar{C} = J^{- \frac{2}{3}} C$ is the distortional part of the right Cauchy-Green strain, and $J = det (F)$ is the volume change.

The underlying assumption in models based on the strain-based formulation is that the preferred materials directions are initially aligned with an orthogonal coordinate system in the reference (stress-free) configuration. These directions may become nonorthogonal only after deformation. Examples of this form of strain energy function include the generalized Fung-type form (see Generalized Fung form below).

From Equation 1 the variation of $U$ is given as

\begin{matrix} (2) & δ U = \frac{\partial U}{\partial {\bar{ε}}^{G}} : δ {\bar{ε}}^{G} + \frac{\partial U}{\partial J} δ J = \bar{F} \cdot \frac{\partial U}{\partial {\bar{ε}}^{G}} \cdot {\bar{F}}^{T} : δ e + \frac{\partial U}{\partial J} δ J . \end{matrix}

Using the principle of virtual work, the variation of the strain energy potential can be written as

\begin{matrix} (3) & δ W_{I} = \int_{V^{0}} J (S : δ e - p δ ε^{vol}) d V^{0} = \int_{V^{0}} δ U d V^{0} \end{matrix}

(see Equation 7).

For a compressible material the strain variations are arbitrary, so this equation defines the stress components for such a material as

S = \frac{1}{J} DEV [\bar{F} \cdot \frac{\partial U}{\partial {\bar{ε}}^{G}} \cdot {\bar{F}}^{T}]

and

p = - \frac{\partial U}{\partial J} .

When the material response is almost incompressible, the pure displacement formulation, in which the strain invariants are computed from the kinematic variables of the finite element model, can behave poorly. One difficulty is that from a numerical point of view the stiffness matrix is almost singular because the effective bulk modulus of the material is so large compared to its effective shear modulus, thus causing difficulties with the solution of the discretized equilibrium equations. Similarly, in Abaqus/Explicit the high bulk modulus increases the dilatational wave speed, thus reducing the stable time increment substantially. To avoid such problems, Abaqus/Standard offers a “mixed” formulation for such cases (refer to Hyperelastic material behavior).

Invariant-based formulation

Using the continuum theory of fiber-reinforced composites (Spencer, 1984), the strain energy function can be expressed directly in terms of the invariants of the deformation tensor and fiber directions. For example, consider a composite material that consists of an isotropic hyperelastic matrix reinforced with $N$ families of fibers. The directions of the fibers in the reference configuration are characterized by a set of unit vectors $A_{α}$ , ( $α = 1, \dots, N$ ). Assuming that the strain energy depends not only on deformation, but also on the fiber directions, the following form is postulated:

U = U (C, A_{α}); α = 1, \dots, N .

The strain energy of the material must remain unchanged if both matrix and fibers in the reference configuration undergo a rigid body rotation. Then, following Spencer (1984), the strain energy can be expressed as an isotropic function of an irreducible set of scalar invariants that form the integrity basis of the tensor $C$ and the vectors $A_{α}$ :

\begin{matrix} (4) & U = U ({\bar{I}}_{1}, {\bar{I}}_{2}, J, {\bar{I}}_{4 (α β)}, {\bar{I}}_{5 (α β)}, ζ_{α β}), α = 1, \dots, N; β = 1, \dots, α, \end{matrix}

where ${\bar{I}}_{1}$ and ${\bar{I}}_{2}$ are the first and second deviatoric strain invariants; $J$ is the volume ratio (or third strain invariant); and ${\bar{I}}_{4 (α β)}$ and ${\bar{I}}_{5 (α β)}$ are the pseudo-invariants of $\bar{C}$ , $A_{α}$ , and $A_{β}$ , defined as

{\bar{I}}_{4 (α β)} = A_{α} \cdot \bar{C} \cdot A_{β}, and {\bar{I}}_{5 (α β)} = A_{α} \cdot {\bar{C}}^{2} \cdot A_{β}, α = 1, \dots, N; β = 1, \dots, α .

The terms $ζ_{α β}$ are geometrical constants (independent of deformation) equal to the cosine of the angle between the directions of any two families of fibers in the reference configuration,

ζ_{α β} = A_{α} \cdot A_{β}, α = 1, \dots, N; β = 1, \dots, α .

Unlike in the case of the strain-based formulation, in the invariant-based formulation the fiber directions need not be orthogonal in the initial configuration. An example of the invariant-based energy function is the form proposed by Holzapfel, Gasser, and Ogden (2000) for arterial walls (see Holzapfel-Gasser-Ogden form below).

From Equation 4 the variation of $U$ is given as

\begin{matrix} (5) & δ U = \frac{\partial U}{\partial {\bar{I}}_{1}} δ {\bar{I}}_{1} + \frac{\partial U}{\partial {\bar{I}}_{2}} δ {\bar{I}}_{2} + \frac{\partial U}{\partial J} δ J + \sum_{α = 1}^{N} \sum_{β = 1}^{α} (\frac{\partial U}{\partial {\bar{I}}_{4 (α β)}} δ {\bar{I}}_{4 (α β)} + \frac{\partial U}{\partial {\bar{I}}_{5 (α β)}} δ {\bar{I}}_{5 (α β)}) . \end{matrix}

Using the principle of virtual work (Equation 3) and after some lengthy derivations, the stress components for a compressible material are found to be given as

\begin{matrix} S & ​ = \frac{2}{J} DEV [(\frac{\partial U}{\partial {\bar{I}}_{1}} + {\bar{I}}_{1} \frac{\partial U}{\partial {\bar{I}}_{2}}) \bar{B} - \frac{\partial U}{\partial {\bar{I}}_{2}} \bar{B} \cdot \bar{B}] \\ + \overset{N}{\sum_{α = 1}} \overset{α}{\sum_{β = 1}} \frac{\partial U}{\partial {\bar{I}}_{4 (α β)}} DEV ({\bar{a}}_{α} {\bar{a}}_{β} + {\bar{a}}_{β} {\bar{a}}_{α}) \\ + \overset{N}{\sum_{α = 1}} \overset{α}{\sum_{β = 1}} \frac{\partial U}{\partial {\bar{I}}_{5 (α β)}} DEV ({\bar{a}}_{α} {\bar{a}}_{β}^{​^{'}} + {\bar{a}}_{α}^{​^{'}} {\bar{a}}_{β} + {\bar{a}}_{β} {\bar{a}}_{α}^{​^{'}} + {\bar{a}}_{β}^{​^{'}} {\bar{a}}_{α}), \end{matrix}

and

p = - \frac{\partial U}{\partial J};

where ${\bar{a}}_{α} \overset{def}{=} \bar{F} \cdot A_{α}$ and ${\bar{a}}_{α}^{^{'}} \overset{def}{=} \bar{B} \cdot a_{α}$ .

Particular forms of the strain energy potential

Several particular forms of the strain energy potential are available in Abaqus. The incompressible or almost incompressible models make up:

the generalized Fung form and
the Holzapfel-Gasser-Ogden form.

In addition, Abaqus provides a general capability to support user-defined forms of the strain energy potential via two sets of user subroutines: one for strain-based and one for invariant-based formulations.

Generalized Fung form

The generalized Fung strain energy potential in Abaqus is based on the two-dimensional exponential form proposed by Fung et al. (1979), suitably generalized to arbitrary three-dimensional states following Humphrey (1995). It has the following form:

U = \frac{c}{2} (e^{Q} - 1) + \frac{1}{D} (\frac{{(J^{e ℓ})}^{2} - 1}{2} - \ln J^{e ℓ}),

where U is the strain energy per unit of reference volume, c and D are temperature-dependent material parameters, $J^{e ℓ}$ is the elastic volume ratio, and $Q$ is defined as

Q = {\bar{ε}}^{G} : b : {\bar{ε}}^{G} = {\bar{ε}}_{i j}^{G} b_{i j k l} {\bar{ε}}_{k l}^{G},

where $b_{i j l k}$ is a dimensionless symmetric fourth-order tensor of anisotropic material constants that can be temperature dependent and ${\bar{ε}}_{i j}^{G}$ are the components of the modified Green strain tensor.

The elastic volume ratio, $J^{e ℓ}$ , relates the total volume ratio, J, and the thermal volume ratio, $J^{t h}$ :

J^{e ℓ} = \frac{J}{J^{t h}} .

$J^{t h}$ is given by

J^{t h} = (1 + ε_{1}^{t h}) (1 + ε_{2}^{t h}) (1 + ε_{3}^{t h}),

where $ε_{i}^{t h}$ are the principal thermal expansion strains that are obtained from the temperature and the thermal expansion coefficients.

The initial deviatoric elasticity tensor, ${\bar{D}}_{0}$ , and bulk modulus, $K_{0}$ , are given by

{\bar{D}}_{0} = c b, K_{0} = \frac{2}{D} .

Abaqus supports two forms of the generalized Fung model: anisotropic and orthotropic. The number of independent components $b_{i j l k}$ that must be specified depends on the level of anisotropy of the material: 21 for the fully anisotropic case and 9 for the orthotropic case.

Holzapfel-Gasser-Ogden form

The form of the strain energy potential is based on that proposed by Holzapfel, Gasser, and Ogden (2000) and Gasser, Ogden, and Holzapfel (2006) for modeling arterial layers with distributed collagen fiber orientations:

U = C_{10} ({\bar{I}}_{1} - 3) + \frac{1}{D} (\frac{{(J^{e ℓ})}^{2} - 1}{2} - \ln J^{e ℓ}) + \frac{k_{1}}{2 k_{2}} \overset{N}{\sum_{α = 1}} {\exp [k_{2} {⟨ {\bar{E}}_{α} ⟩}^{2}] - 1},

with

{\bar{E}}_{α} \overset{def}{=} κ ({\bar{I}}_{1} - 3) + (1 - 3 κ) ({\bar{I}}_{4 (α α)} - 1),

where U is the strain energy per unit of reference volume; $C_{10}$ , D, $k_{1}$ , $k_{2}$ , and $κ$ are temperature-dependent material parameters; $N$ is the number of families of fibers ( $N \leq 3$ ); ${\bar{I}}_{1}$ is the first invariant of $\bar{C}$ ; $J^{e ℓ}$ is the elastic volume ratio, as defined previously; and ${\bar{I}}_{4 (α α)}$ are pseudo-invariants of $\bar{C}$ and $A_{α}$ .

The model assumes that the directions of the collagen fibers within each family are dispersed (with rotational symmetry) about a mean preferred direction. The parameter $κ$ ( $0 \leq κ \leq \frac{1}{3}$ ) describes the level of dispersion in the fiber directions. If $ρ (Θ)$ is the orientation density function that characterizes the distribution (it represents the normalized number of fibers with orientations in the range $[Θ, Θ + d Θ]$ with respect to the mean direction), the parameter $κ$ is defined as

κ = \frac{1}{4} \int_{0}^{π} ρ (Θ) \sin^{3} Θ d Θ .

It is also assumed that all families of fibers have the same mechanical properties and the same dispersion. When $κ = 0,$ the fibers are perfectly aligned (no dispersion). When $κ = \frac{1}{3},$ the fibers are randomly distributed and the material becomes isotropic; this corresponds to a spherical orientation density function.

The strain-like quantity ${\bar{E}}_{α}$ characterizes the deformation of the family of fibers with mean direction $A_{α}$ . ${\bar{E}}_{α} = {\bar{I}}_{4 (α α)} - 1$ for perfectly aligned fibers ( $κ = 0$ ), and ${\bar{E}}_{α} = ({\bar{I}}_{1} - 3) / 3$ for randomly distributed fibers ( $κ = \frac{1}{3}$ ).

The first two terms in the expression of the strain energy function represent the distortional and volumetric contributions of the noncollagenous isotropic ground material; and the third term represents the contributions from the different families of collagen fibers, taking into account the effects of dispersion. A basic assumption of the model is that collagen fibers can support only tension, as they would buckle under compressive loading. Thus, the anisotropic contribution in the strain energy function appears only when the strain of the fibers is positive or, equivalently, when ${\bar{E}}_{α} > 0$ . This condition is enforced by the term $⟨ {\bar{E}}_{α} ⟩$ , where the operator $⟨ \cdot ⟩$ stands for the Macauley bracket and is defined as $⟨ x ⟩ = \frac{1}{2} (| x | + x)$ .

The initial deviatoric elasticity tensor, ${\bar{D}}_{0}$ , and bulk modulus, $K_{0}$ , are given by

{\bar{D}}_{0} = 4 C_{10} ℑ + 2 {(1 - 3 κ)}^{2} k_{1} \sum_{α = 1}^{N} H ({\bar{E}}_{α}) A_{α} A_{α} A_{α} A_{α}, K_{0} = \frac{2}{D},

where $ℑ$ is the fourth-order unit tensor and $H (x)$ is the Heaviside unit step function.

User-defined form: strain-based

Abaqus also allows other forms of strain-based energy potentials to be defined via user subroutines UANISOHYPER_STRAIN in Abaqus/Standard and VUANISOHYPER_STRAIN in Abaqus/Explicit by programming the first and second derivatives of the strain energy potential with respect to the components of the modified Green strain and the elastic volume ratio, $J^{e ℓ}$ .

User-defined form: invariant-based

Abaqus also allows other forms of invariant-based energy potentials to be defined via user subroutines UANISOHYPER_INV in Abaqus/Standard and VUANISOHYPER_INV in Abaqus/Explicit by programming the first and second derivatives of the strain energy potential with respect to each invariant.