Isotropic Elastic Material

Two kinds of isotropic elastic materials are considered:

Linear elastic materials with Hooke’s law,
Nonlinear elastic materials with Ogden, Mooney-Rivlin and Arruda-Boyce laws.

Linear Elastic Material (LAW1)

This material law is used to model purely elastic materials, or materials that remain in the elastic range. The Hooke's law requires only two values to be defined; the Young's or elastic modulus $E$ $E$ , and Poisson's ratio, $υ$ $υ$ . The law represents a linear relation between stress and strain.

Ogden Materials (LAW42, LAW69 and LAW82)

Ogden's law is applied to slightly compressible materials as rubber or elastomer foams undergoing large deformation with an elastic behavior. The detailed theory for Odgen material models can be found in ¹. The strain energy

W

$W$ is expressed in a general form as a function of

W (λ_{1}, λ_{2}, λ_{3})

$W (λ_{1}, λ_{2}, λ_{3})$ :(1)

W (λ_{1}, λ_{2}, λ_{3}) = 5 \sum p = 1 \frac{μ_{p}}{α_{p}} ({¯ λ}_{1}_{p}^{α_{p}} + {¯ λ}_{2}_{p}^{α_{p}} + {¯ λ}_{3}_{p}^{α_{p}} - 3) + \frac{K}{2} {(J - 1)}^{2}

$W (λ_{1}, λ_{2}, λ_{3}) = \sum_{p = 1}^{5} \frac{μ_{p}}{α_{p}} ({\bar{λ}}_{1}^{α_{p}} + {\bar{λ}}_{2}^{α_{p}} + {\bar{λ}}_{3}^{α_{p}} - 3) + \frac{K}{2} {(J - 1)}^{2}$

Where, $λ_{1}$ $λ_{1}$ , i_th principal stretch

$λ_{i} = 1 + ε_{i}$ $λ_{i} = 1 + ε_{i}$ , with $ε_{i}$ $ε_{i}$ being the i_th principal engineering strain

J

$J$ is relative volume with:(2)

J = λ_{1} \cdot λ_{2} \cdot λ_{3} = \frac{ρ_{0}}{ρ}

$J = λ_{1} \cdot λ_{2} \cdot λ_{3} = \frac{ρ_{0}}{ρ}$

{¯ λ}_{i}

${\bar{λ}}_{i}$ is the deviatoric stretch(3)

{¯ λ}_{i} = J^{- \frac{1}{3}} λ_{i}

${\bar{λ}}_{i} = J^{- \frac{1}{3}} λ_{i}$

$α_{p}$ $α_{p}$ and $μ_{p}$ $μ_{p}$ are the material constants.

$p$ $p$ is order of Ogden model and defines the number of coefficients pairs $(α_{p}, μ_{p})$ $(α_{p}, μ_{p})$ .

This law is very general due to the choice of coefficient pair

(α_{p}, μ_{p})

$(α_{p}, μ_{p})$ .

If $p$ $p$ =1, then one pair $(α_{1}, μ_{1})$ $(α_{1}, μ_{1})$ of material constants is needed andin this case if $α_{1} = 2$ $α_{1} = 2$ then it becomes a Neo-hookean material model.
If $p$ $p$ =2 then two pairs $(α_{1}, μ_{1}), (α_{2}, μ_{2})$ $(α_{1}, μ_{1}), (α_{2}, μ_{2})$ of material constants are needed and in this case if $α_{1} = 2$ $α_{1} = 2$ and $α_{1} = -2$ $α_{1} = -2$ then it becomes a Mooney-Rivlin material model

For uniform dilitation:(4)

λ_{1} = λ_{2} = λ_{3} = λ

$λ_{1} = λ_{2} = λ_{3} = λ$

The strain energy function can be decomposed into deviatoric part

¯ ¯¯ ¯ W ({¯ λ}_{1}, {¯ λ}_{2}, {¯ λ}_{3})

$\bar{W} ({\bar{λ}}_{1}, {\bar{λ}}_{2}, {\bar{λ}}_{3})$ and spherical part

U (J)

$U (J)$ :(5)

W = ¯ ¯¯ ¯ W ({¯ λ}_{1}, {¯ λ}_{2}, {¯ λ}_{3}) + U (J)

$W = \bar{W} ({\bar{λ}}_{1}, {\bar{λ}}_{2}, {\bar{λ}}_{3}) + U (J)$

With:

$\begin{matrix} \bar{W} ({\bar{λ}}_{1}, {\bar{λ}}_{2}, {\bar{λ}}_{3}) = \sum_{p} \frac{μ_{p}}{α_{p}} ({\bar{λ}}_{1}^{α_{p}} + {\bar{λ}}_{2}^{α_{p}} + {\bar{λ}}_{3}^{α_{p}} - 3) \\ U (J) = \frac{K}{2} {(J - 1)}^{2} \end{matrix}$

The stress

$σ_{i}$ corresponding to this strain energy is given by:(6)

$σ_{i} = \frac{λ_{i}}{J} \frac{\partial W}{\partial λ_{i}}$

which can be written as:(7)

$σ_{i} = \frac{λ_{i}}{J} \frac{\partial W}{\partial λ_{i}} = \frac{λ_{i}}{J} (\sum_{j = 1}^{3} \frac{\partial \bar{W}}{\partial {\bar{λ}}_{j}} \frac{\partial {\bar{λ}}_{j}}{\partial λ_{i}} + \frac{\partial U}{\partial J} \frac{\partial J}{\partial λ_{i}})$

Since $λ_{1} \frac{\partial J}{\partial λ_{i}} = J$ and $\frac{\partial {\bar{λ}}_{j}}{\partial λ_{i}} = \frac{2}{3} J^{- \frac{1}{3}}$ for i=j and $\frac{\partial {\bar{λ}}_{j}}{\partial λ_{i}} = \frac{1}{3} J^{- \frac{1}{3}} \frac{λ_{j}}{λ_{i}}$ for i≠j

Equation 7 is simplified to:(8)

$σ_{i} = \frac{1}{J} ({\bar{λ}}_{i} \frac{\partial \bar{W}}{\partial {\bar{λ}}_{i}} - (\frac{1}{3} \sum_{j = 1}^{3} {\bar{λ}}_{j} \frac{\partial \bar{W}}{\partial {\bar{λ}}_{j}} - J \frac{\partial U}{\partial J}))$

For which the deviator of the Cauchy stress tensor

$s_{i}$ , and the pressure

$P$ would be:(9)

$s_{i} = \frac{1}{J} ({\bar{λ}}_{i} \frac{\partial \bar{W}}{\partial {\bar{λ}}_{i}} - \frac{1}{3} \sum_{j = 1}^{3} {\bar{λ}}_{j} \frac{\partial \bar{W}}{\partial {\bar{λ}}_{j}})$

(10)

$p = - \frac{1}{3} \sum_{j = 1}^{3} σ_{j} = - \frac{\partial U}{\partial J}$

Only the deviatoric stress above is retained, and the pressure is computed independently:(11)

$P = K \cdot F s c a l e_{b l k} \cdot f_{b l k} (J) \cdot (J - 1)$

Where,

$f_{b l k} (J)$ a user-defined function related to the bulk modulus

$K$ in LAW42 and LAW69:(12)

$K = μ \cdot \frac{2 (1 + ν)}{3 (1 - 2 ν)}$

For an imcompressible material

$(ν \approx 0.5)$ ,

$J = 1$ and no pressure in material.(13)

$μ = \frac{\sum_{p} μ_{p} \cdot α_{p}}{2}$

With

$μ$ being the initial shear modulus, and

$υ$ the Poisson's ratio.

Note: For an incompressible material you have

$υ \approx 0.5$ . However,

$υ \approx 0.495$ is a good compromise to avoid too small time steps in explicit codes.

A particular case of the Ogden material model is the Mooney-Rivlin material law which has two basic assumptions:

The rubber is incompressible and isotropic in unstrained state
The strain energy expression depends on the invariants of Cauchy tensor

The three invariants of the Cauchy-Green tensor are:(14)

$I_{1} = λ_{1}^{2} + λ_{2}^{2} + λ_{3}^{2}$

(15)

$I_{2} = λ_{1}^{2} λ_{2}^{2} + λ_{2}^{2} λ_{3}^{2} + λ_{3}^{2} λ_{1}^{2}$

For incompressible material:(16)

$I_{3} = λ_{1}^{2} λ_{2}^{2} λ_{3}^{2} = 1$

The Mooney-Rivlin law gives the closed expression of strain energy as:(17)

$W = C_{10} (I_{1} - 3) + C_{01} (I_{2} - 3)$

with:(18)

$μ_{1} = 2 \cdot C_{10}$

$μ_{2} = - 2 \cdot C_{01}$
$α_{1} = 2$
$α_{2} = - 2$

The model can be generalized for a compressible material.

Viscous Effects in LAW42

Viscous effects are modeled through the Maxwell model:

Figure 1. Maxwell Model

Where, the shear modulus of the hyper-elastic law

$μ$ is exactly the long-term shear modulus

$G_{\infty}$ .(19)

$μ = \frac{\sum_{p} μ_{p} \cdot α_{p}}{2} = G_{\infty}$

$τ_{i}$ are relaxation times: $τ_{i} = \frac{η_{i}}{G_{i}}$

Rate effects are modeled through visco-elasticity using convolution integral using Prony series. This corresponds to extension of small deformation theory to finite deformation.

This viscous stress is added to the elastic one.

The visco-Kirchoff stress is given by:(20)

$τ^{v} = {\sum_{i = 1}^{M} G_{i} \int_{0}^{t} e}^{- \frac{t - s}{τ_{i}}} \frac{d}{d s} [dev (\bar{F} {\bar{F}}^{T})] d s$

Where,

$M$: Order of the Maxwell model
$F$: Deformation gradient matrix
$\bar{F} = J^{- \frac{1}{3}} F$
$dev (\bar{F} {\bar{F}}^{T})$: Denotes the deviatoric part of tensor $\bar{F} {\bar{F}}^{T}$

The viscous-Cauchy stress is written as:(21)

$σ^{v} (t) = \frac{1}{J} τ^{v} (t)$

LAW69, Ogden Material Law (Using Test Data as Input)

This law, like /MAT/LAW42 (OGDEN) defines a hyperelastic and incompressible material specified using the Ogden or Mooney-Rivlin material models. Unlike LAW42 where the material parameters are input this law computes the material parameters from an input engineering stress-strain curve from a uniaxial tension and compression tests. This material can be used with shell and solid elements.

The strain energy density formulation used depends on the law_ID.

law_ID =1, Ogden law (Same as LAW42):: $W (λ_{1}, λ_{2}, λ_{3}) = \sum_{p = 1}^{5} \frac{μ_{p}}{α_{p}} ({\bar{λ}}_{1}^{α_{p}} + {\bar{λ}}_{2}^{α_{p}} + {\bar{λ}}_{3}^{α_{p}} - 3) + \frac{K}{2} {(J - 1)}^{2}$
law_ID =2, Mooney-Rivlin law: $W = C_{10} (I_{1} - 3) + C_{01} (I_{2} - 3)$

Curve Fitting

After reading the stress-strain curve (fct_ID₁), Radioss calculates the corresponding material parameter pairs using a nonlinear least-square fitting algorithm. For classic Ogden law, (law_ID =1), the calculated material parameter pairs are $μ_{p}$ and $α_{p}$ where the value of $p$ is defined via the N_pair input. The maximum value of N_pair = 5 with a default value of 2.

For the Mooney-Rivlin law (law_ID =2), the material parameter $C_{10}$ and $C_{01}$ are calculated remembering that $μ_{p}$ and $α_{p}$ for the LAW42 Ogden law can be calculated using this conversion.

$μ_{1} = 2 \cdot C_{10}$ , $μ_{2} = - 2 \cdot C_{01}$ , $α_{1} = 2$ and $α_{2} = - 2$ .

The minimum test data input should be a uniaxial tension engineering stress strain curve. If uniaxial compression data is available, the engineering strain should increate monotonically from a negative value in compression to a positive value in tension. In compression, the engineering strain should not be less than -1.0 since -100% strain is physically not possible.

This material law is stable when (with $p$ =1,…5) is satisfied for parameter pairs for all loading conditions. By default, Radioss tries to fit the curve by accounting for these conditions (I_check=2). If a proper fit cannot be found, then Radioss uses a weaker condition (I_check=1: ) which ensures that the initial shear hyperelastic modulus ( $μ$ ) is positive.

Once the material parameters are calculated by the Radioss Starter in LAW69, all the calculations done by LAW69 in the simulation are the same as LAW42.

LAW82

The Ogden model used in LAW82 is:(22)

$W = \sum_{i = 1}^{N} \frac{2 μ_{i}}{α_{i}^{2}} ({\bar{λ}}_{1}^{α_{i}} + {\bar{λ}}_{2}^{α_{i}} + {\bar{λ}}_{3}^{α_{i}} - 3) + \sum_{i = 1}^{N} \frac{1}{D_{i}} {(J - 1)}^{2 i}$

The Bulk Modulus is calculated as

$K = \frac{2}{D_{1}}$ based on these rules:

If $ν = 0$ , $D_{1}$ should be entered.
If $ν \neq 0$ , $D_{1}$ input is ignored and will be recalculated and output in the Starter output using:(23) $D_{1} = \frac{3 (1 - 2 v)}{μ (1 + v)}$
If $ν = 0$ and $D_{1}$ =0, a default value of $ν = 0.475$ is used and $D_{1}$ is calculated using Equation 23

LAW88, A simplified hyperelastic material with strain rate effects

This law utilizes tabulated uniaxial tension and compression engineering stress and strain test data at different strain rates to model incompressible materials. It is only compatible with solid elements. The material is based on Ogden’s strain energy density function but does not require curve fitting to extract material constants like most other hyperelastic material models. Strain rate effects can be modeled by including engineering stress strain test data at different strain rates. This can be easier than calculating viscous parameters for traditional hyperelastic material models. The following Ogden strain energy density function is used but instead of extracting material constants via curve fitting this law determines the Ogden function directly from the uniaxial engineering stress strain curve tabulated data. ⁵

(24)

$W = \underset{d e v i a t o r i c p a r t}{\underset{︸}{\sum_{i = 1}^{3} \sum_{j = 1}^{m} \frac{μ_{j}}{α_{j}} ({\bar{λ}}_{i}^{α_{j}} - 1)}} + \underset{s p h e r i c a l p a r t}{\underset{︸}{K (J - 1 - \ln J)}}$

Unloading can be represented using an unloading function, Fscale_unL, or by providing hysteresis, Hys and shape factor, Shape, inputs to a damage model based on energy.

When using the damage model, the loading curves are used for both loading and unloading and the unloading stress tensor is reduced by:(25)

$σ = (1 - D) σ$

(26)

$D = (1 - H y s) (1 - {(\frac{W_{c u r}}{W_{\max}})}^{S h a p e})$

If the unloading function, Fscale_unL, is entered, unloading is defined based on the unloading flag, Tension and the damage model is not used.

Arruda-Boyce Material (LAW92)

LAW92 describes the Arruda-Boyce material model, which can be used to model hyperelastic behavior. The Arruda-Boyce model is based on the statistical mechanics of a material with a cubic representative volume element containing eight chains along the diagonal directions. It assumes that the chain molecules are located on the average along the diagonals of the cubic in principal stretch space.

The strain energy density function is:(27)

$W = \underset{W ({\bar{I}}_{1})}{\underset{︸}{μ \sum_{i = 1}^{5} \frac{c_{i}}{{(λ_{m})}^{2 i - 2}} ({\bar{I}}_{1}^{i} - 3^{i})}} + \underset{U (J)}{\underset{︸}{\frac{1}{D} (\frac{J^{2} - 1}{2} + \ln (J))}}$

Where, Material constant $c_{i}$ are:

$c_{1} = \frac{1}{2}, c_{2} = \frac{1}{20}, c_{3} = \frac{11}{1050}, c_{4} = \frac{19}{7000}, c_{5} = \frac{519}{673750}$

$μ$: Shear modulus
$μ_{0}$: Initial shear modulus

(28)

$μ_{0} = μ (1 + \frac{3}{5 λ_{m}^{2}} + \frac{99}{175 λ_{m}^{4}} + \frac{513}{875 λ_{m}^{6}} + \frac{42039}{67375 λ_{m}^{8}})$

$λ_{m}$ is the limit of stretch which describes the beginning of hardening phase in tension (locking strain in tension) and so it is also called the locking stretch.

Arruda-Boyce is always stable if positive values of the shear modulus, $μ$ , and the locking stretch, $λ_{m}$ are used.

${\bar{I}}_{1}$ is deviatoric part of first strain invarient

$I_{1}$

(29)

${\bar{I}}_{1} = {\bar{λ}}_{1}^{2} + {\bar{λ}}_{2}^{2} + {\bar{λ}}_{3}^{2} = J^{- 2 / 3} I_{1}$

with ${\bar{λ}}_{i} = J^{- \frac{1}{3}} λ_{i}$

$D$ is a material parameter for the bulk modulus computation given as:(30)

$D = \frac{2}{K}$

The Cauchy stress corresponding to above strain energy is:(31)

$σ_{i} = \frac{λ_{i}}{J} \frac{\partial W}{\partial λ_{i}}$

For incompressible materials, the Cauchy stress is then given by:

Uniaxial test(32) $σ = λ \frac{\partial W}{\partial λ} = 2 μ (λ^{2} - \frac{1}{λ}) \sum_{i = 1}^{5} \frac{i \cdot c_{i}}{λ_{m}^{2 i - 2}} {({\bar{λ}}_{1}^{2} + {\bar{λ}}_{2}^{2} + {\bar{λ}}_{3}^{2})}^{i - 1}$
with $λ_{1} = λ$ and $λ_{2} = λ_{3} = λ^{- \frac{1}{2}}$ , then ${\bar{I}}_{1} = λ^{2} + \frac{2}{λ}$

and nominal stress is:(33) $N^{t h} = \frac{\partial W}{\partial λ} = 2 μ (λ - λ^{- 2}) \sum_{i = 1}^{5} \frac{i \cdot c_{i}}{{(λ_{m})}^{2 i - 2}} {(λ^{2} + \frac{2}{λ})}^{i - 1}$
Equibiaxial test(34) $σ = λ \frac{\partial W}{\partial λ} = 2 μ (λ^{2} - \frac{1}{λ^{4}}) \sum_{i = 1}^{5} \frac{i \cdot c_{i}}{λ_{m}^{2 i - 2}} {({\bar{λ}}_{1}^{2} + {\bar{λ}}_{2}^{2} + {\bar{λ}}_{3}^{2})}^{i - 1}$
with $λ_{1} = λ_{2} = λ$ and $λ_{3} = λ^{- 2}$ , then ${\bar{I}}_{1} = 2 λ^{2} + \frac{1}{λ^{4}}$

and the nominal stress is:(35) $N^{t h} = \frac{\partial W}{\partial λ} = 2 μ (λ - λ^{- 5}) \sum_{i = 1}^{5} \frac{i \cdot c_{i}}{{(λ_{m})}^{2 i - 2}} {(2 λ^{2} + \frac{1}{λ^{4}})}^{i - 1}$
Planar test(36) $σ = λ \frac{\partial W}{\partial λ} = 2 μ (λ^{2} - \frac{1}{λ^{2}}) \sum_{i = 1}^{5} \frac{i \cdot c_{i}}{λ_{m}^{2 i - 2}} {({\bar{λ}}_{1}^{2} + {\bar{λ}}_{2}^{2} + {\bar{λ}}_{3}^{2})}^{i - 1}$
with $λ_{1} = λ, λ_{3} = 1$ and $λ_{2} = λ^{- 1} $ , then ${\bar{I}}_{1} = λ^{2} + 1 + λ^{- 2}$

and nominal stress is:(37) $N^{t h} = \frac{\partial W}{\partial λ} = 2 μ (λ - λ^{- 3}) \sum_{i = 1}^{5} \frac{i c_{i}}{{(λ_{m})}^{2 i - 2}} {(λ^{2} + 1 + λ^{- 2})}^{i - 1}$

Additional information about Arruda-Boyce model. ² ³

Yeoh Material (LAW94)

The Yeoh model (LAW94) ⁴is a hyperelastic material model that can be used to describe incompressible materials. The strain energy density function of LAW94 only depends on the first strain invariant and is computed as:(38)

$W = \sum_{i = 1}^{3} [\underset{W ({\bar{I}}_{1})}{\underset{︸}{C_{i 0} {({\bar{I}}_{1} - 3)}^{i}}} + \underset{U (J)}{\underset{︸}{\frac{1}{D_{i}} {(J - 1)}^{2 i}}}]$

Where,

${\bar{I}}_{1} = {\bar{λ}}_{1}^{2} + {\bar{λ}}_{2}^{2} + {\bar{λ}}_{3}^{2}$: First strain invariant
${\bar{λ}}_{i} = J^{- \frac{1}{3}} λ_{i}$: Deviatoric stretch

The Cauchy stress is computed as:(39)

$σ_{i} = \frac{λ_{i}}{J} \frac{\partial W}{\partial λ_{i}}$

For incompressible materials with $i$ =1 only and $D_{1}$ are input and the Yeoh model is reduced to a Neo-Hookean model.

The material constant specify the deviatoric part (shape change) of the material and parameters $D_{1}$ , $D_{2}$ , $D_{3}$ specify the volumetric change of the material. These six material constants need to be calculated by curve fitting material test data. RD-E: 5600 Hyperelastic Material with Curve Input includes a Yeoh fitting Compose script for uniaxial test data. The Yeoh material model has been shown to model all deformation models even if the curve fit was obtained using only uniaxial test data.

The initial shear modulus and the bulk modulus are computed as:

$μ = 2 \cdot C_{10}$ and $K = \frac{2}{D_{1}}$

LAW94 is available only as an incompressible material model.

$D_{1}$ =0, an incompressible material is considered, where

$ν = 0.495$ and

$D_{1}$ is calculated as:(40)

$D_{1} = \frac{3 (1 - 2 v)}{μ (1 + v)}$

Ogden R.W., “Nonlinear Elastic Deformations”, Ellis Horwood, 1984.

Arruda, E.M. and Boyce, M.C., “A three-dimensional model for the large stretch behavior of rubber elastic materials”, J. Mech. Phys. Solids, 41(2), pp. 389–412, 1993.

Jörgen Bergström, “Mechanics of solid polymers: theory and computational modeling”, pp. 250-254, 2015.

Yeoh, O. H., “Some forms of the strain energy function for rubber”, Rubber Chemistry and Technology, Vol. 66, Issue 5, pp. 754-771, November 1993.

Kolling S., Du Bois P.A., Benson D.J., and Feng W.W., "A tabulated formulation of hyperelasticity with rate effects and damage." Computational Mechanics 40, no. 5 (2007).