Chapter 3

Constitutive Models $–$ Relations between Stress and Strain

3.5 Hyperelasticity $–$ time independent behavior of rubbers and foams subjected to large strains

Hyperelastic constitutive laws are used to model materials that respond elastically when subjected to very large strains. They account both for nonlinear material behavior and large shape changes.  The main applications of the theory are (i) to model the rubbery behavior of a polymeric material, and (ii) to model polymeric foams that can be subjected to large reversible shape changes (e.g. a sponge).

In general, the response of a typical polymer is strongly dependent on temperature, strain history and loading rate.  The behavior will be described in more detail in the next section, where we present the theory of viscoelasticity.  For now, we note that polymers have various regimes of mechanical behavior, referred to as glassy,’ viscoelastic’ and rubbery.’   The various regimes can be identified for a particular polymer by applying a sinusoidal variation of shear stress to the solid and measuring the resulting shear strain amplitude.  A typical result is illustrated in the figure, which shows the apparent shear modulus (ratio of stress amplitude to strain amplitude) as a function of temperature.

At a critical temperature known as the glass transition temperature, a polymeric material undergoes a dramatic change in mechanical response.  Below this temperature, it behaves like a glass, with a stiff response. Near the glass transition temperature, the stress depends strongly on the strain rate.  At the glass transition temperature, there is a dramatic drop in modulus.  Above this temperature, there is a regime where the polymer shows rubbery’ behavior $–$ the response is elastic; the stress does not depend strongly on strain rate or strain history, and the modulus increases with temperature.  All polymers show these general trends, but the extent of each regime, and the detailed behavior within each regime, depend on the solid’s molecular structure.  Heavily cross-linked polymers (elastomers) are the most likely to show ideal rubbery behavior.   Hyperelastic constitutive laws are intended to approximate this rubbery’ behavior.

Features of the behavior of a solid rubber:

1.      The material is close to ideally elastic. i.e. (i) when deformed at constant temperature or adiabatically, stress is a function only of current strain and independent of history or rate of loading, (ii) the behavior is reversible: no net work is done on the solid when subjected to a closed cycle of strain under adiabatic or isothermal conditions.

2.      The material strongly resists volume changes.  The bulk modulus (the ratio of volume change to hydrostatic component of stress) is comparable to that of metals or covalently bonded solids;

3.      The material is very compliant in shear $–$ shear modulus is of the order of ${10}^{-5}$ times that of most metals;

4.      The material is isotropic $–$ its stress-strain response is independent of material orientation.

5.      The shear modulus is temperature dependent: the material becomes stiffer as it is heated, in sharp contrast to metals;

6.      When stretched, the material gives off heat.

Polymeric foams (e.g. a sponge) share some of these properties:

1.      They are close to reversible, and show little rate or history dependence.

2.      In contrast to rubbers, most foams are highly compressible $–$ bulk and shear moduli are comparable.

3.      Foams have a complicated true stress-true strain response, generally resembling the figure to the right.  The finite strain response of the foam in compression is quite different to that in tension, because of buckling in the cell walls.

4.      Foams can be anisotropic, depending on their cell structure.   Foams with a random cell structure are isotropic.

The literature on stress-strain relations for finite elasticity can be hard to follow, partly because nearly every paper uses a different notation, and partly because there are many different ways to write down the same stress-strain law.   You should find that most of the published literature is consistent with the framework given below $–$ but it may take some work to show the equivalence.

All hyperelastic models are constructed as follows:

1.      Define the stress-strain relation for the solid by specifying its strain energy density W  as a function of deformation gradient tensor: W=W(F).    This ensures that the material is perfectly elastic, and also means that we only need to work with a scalar function.  The general form of the strain energy density is guided by experiment; and the formula for strain energy density always contains material properties that can be adjusted to describe a particular material.

2.      The undeformed material is usually assumed to be isotropic $–$ i.e the behavior of the material is independent of the initial orientation of the material with respect to the loading.  If the strain energy density is a function of the Left Cauchy-Green deformation tensor $B=F\cdot {F}^{T}$ the constitutive equation is automatically isotropic.   If B is used as the deformation measure, then the strain energy must be a function of the invariants of B to ensure that the constitutive equation is objective (recall that the invariants of a tensor remain constant under a change of basis)

3.      Formulas for stress in terms of strain are calculated by differentiating the strain energy density as outlined below.

3.5.1 Deformation Measures used in finite elasticity

Suppose that a solid is subjected to a displacement field ${u}_{i}\left({x}_{k}\right)$. Define

The deformation gradient and its Jacobian

${F}_{ij}={\delta }_{ij}+\frac{\partial {u}_{i}}{\partial {x}_{j}}$       $J=\mathrm{det}\left(F\right)$

The Left Cauchy-Green deformation tensor

$B=F\cdot {F}^{T}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{B}_{ij}={F}_{ik}{F}_{jk}$

Invariants of B (these are the conventional definitions)

$\begin{array}{l}{I}_{1}=\text{trace}\left(B\right)={B}_{kk}\\ {I}_{2}=\frac{1}{2}\left({I}_{1}^{2}-B\cdot \cdot B\right)=\frac{1}{2}\left({I}_{1}^{2}-{B}_{ik}{B}_{ki}\right)\\ {I}_{3}=\mathrm{det}B={J}^{2}\end{array}$

An alternative set of invariants of B (more convenient for models of nearly incompressible materials $–$ note that ${\overline{I}}_{1},{\overline{I}}_{2}$ remain constant under a pure volume change)

$\begin{array}{l}{\overline{I}}_{1}=\frac{{I}_{1}}{{J}^{2/3}}\text{=}\frac{{B}_{kk}}{{J}^{2/3}}\\ {\overline{I}}_{2}=\frac{{I}_{2}}{{J}^{4/3}}=\frac{1}{2}\left({\overline{I}}_{1}^{2}-\frac{B\cdot \cdot B}{{J}^{4/3}}\right)=\frac{1}{2}\left({\overline{I}}_{1}^{2}-\frac{{B}_{ik}{B}_{ki}}{{J}^{4/3}}\right)\\ J=\sqrt{\mathrm{det}B}\end{array}$

Principal stretches and principal stretch directions

1.      Let ${e}_{1},{e}_{2},{e}_{3}$ denote the three eigenvalues of B.  The principal stretches are

${\lambda }_{1}=\sqrt{{e}_{1}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\lambda }_{2}=\sqrt{{e}_{2}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\lambda }_{3}=\sqrt{{e}_{3}}$

2.      Let ${b}_{1},{b}_{2},{b}_{3}$ denote three, mutually perpendicular unit eigenvectors of B. These define the principal stretch directions.  (Note: since B is symmetric its eigenvectors are automatically mutually perpendicular as long as no two eigenvalues are the same.  If two, or all three eigenvalues are the same, the eignevectors are not uniquely defined $–$ in this case any convenient mutually perpendicular set of eigenvectors can be used).

3.      Recall that B can be expressed in terms of its eigenvectors and eigenvalues as $B={\lambda }_{1}^{2}{b}^{\left(1\right)}\otimes {b}^{\left(1\right)}+{\lambda }_{2}^{2}{b}^{\left(2\right)}\otimes {b}^{\left(2\right)}+{\lambda }_{3}^{2}{b}^{\left(3\right)}\otimes {b}^{\left(3\right)}$

3.5.2 Stress Measures used in finite elasticity

Usually stress-strain laws are given as equations relating Cauchy stress (true’ stress) ${\sigma }_{ij}$ to left Cauchy-Green deformation tensor.  For some computations it may be more convenient to use other stress measures.  They are defined below, for convenience.

The Cauchy (“true”) stress represents the force per unit deformed area in the solid and is defined by

${n}_{i}{\sigma }_{ij}=\underset{dA\to 0}{Lim}\frac{d{P}_{j}^{\left(n\right)}}{dA}$

Kirchhoff stress  $\tau =J\sigma \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\tau }_{ij}=J{\sigma }_{ij}$

Nominal (First Piola-Kirchhoff) stress  $S=J{F}^{-1}\cdot \sigma \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{S}_{ij}=J{F}_{ik}^{-1}{\sigma }_{kj}$

Material (Second Piola-Kirchhoff) stress  $\Sigma =J{F}^{-1}\cdot \sigma \cdot {F}^{-T}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\Sigma }_{ij}=J{F}_{ik}^{-1}{\sigma }_{kl}{F}_{jl}^{-1}$

3.5.3 Calculating stress-strain relations from the strain energy density

The constitutive law for an isotropic hyperelastic material is defined by an equation relating the strain energy density of the material to the deformation gradient, or, for an isotropic solid, to the three invariants of the strain tensor

$W\left(F\right)=U\left({I}_{1},{I}_{2},{I}_{3}\right)=\overline{U}\left({\overline{I}}_{1},{\overline{I}}_{2},J\right)=\stackrel{˜}{U}\left({\lambda }_{1},{\lambda }_{2},{\lambda }_{3}\right)$

The stress-strain law must then be deduced by differentiating the strain energy density.   This can involve some tedious algebra.  Formulas are listed below for the stress-strain relations for each choice of strain invariant.  The results are derived below

Strain energy density in terms of ${F}_{ij}$

${\sigma }_{ij}=\frac{1}{J}{F}_{ik}\frac{\partial W}{\partial {F}_{kj}}$

Strain energy density in terms of ${I}_{1},{I}_{2},{I}_{3}$

${\sigma }_{ij}=\frac{2}{\sqrt{{I}_{3}}}\left[\left(\frac{\partial U}{\partial {I}_{1}}+{I}_{1}\frac{\partial U}{\partial {I}_{2}}\right){B}_{ij}-\frac{\partial U}{\partial {I}_{2}}{B}_{ik}{B}_{kj}\right]+2\sqrt{{I}_{3}}\frac{\partial U}{\partial {I}_{3}}{\delta }_{ij}$

Strain energy density in terms of ${\overline{I}}_{1},{\overline{I}}_{2},J$

${\sigma }_{ij}=\frac{2}{J}\left[\frac{1}{{J}^{2/3}}\left(\frac{\partial \overline{U}}{\partial {\overline{I}}_{1}}+{\overline{I}}_{1}\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}\right){B}_{ij}-\left({\overline{I}}_{1}\frac{\partial \overline{U}}{\partial {\overline{I}}_{1}}+2{\overline{I}}_{2}\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}\right)\frac{{\delta }_{ij}}{3}-\frac{1}{{J}^{4/3}}\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}{B}_{ik}{B}_{kj}\right]+\frac{\partial \overline{U}}{\partial J}{\delta }_{ij}$

Strain energy density in terms of ${\lambda }_{1},{\lambda }_{2},{\lambda }_{3}$

${\sigma }_{ij}=\frac{{\lambda }_{1}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{1}}{b}_{i}^{\left(1\right)}{b}_{j}^{\left(1\right)}+\frac{{\lambda }_{2}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{2}}{b}_{i}^{\left(2\right)}{b}_{j}^{\left(2\right)}+\frac{{\lambda }_{3}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{3}}{b}_{i}^{\left(3\right)}{b}_{j}^{\left(3\right)}$

Derivations:   We start by deriving the general formula for stress in terms of $W\left(F\right)$:

1.      Note that, by definition, if the solid is subjected to some history of strain, the rate of change of the strain energy density W (F)  must equal the rate of mechanical work done on the material per unit reference volume.

2.      Recall that the rate of work done per unit undeformed volume by body forces and surface tractions is expressed in terms of the nominal stress ${S}_{ij}$ as ${S}_{ji}{\stackrel{˙}{F}}_{ij}$.

3.      Therefore, for any deformation gradient Fij,

$\frac{dW}{dt}=\frac{\partial W}{\partial {F}_{ij}}\frac{\partial {F}_{ij}}{\partial t}={S}_{ji}\frac{\partial {F}_{ij}}{\partial t}$

This must hold for all possible ${\stackrel{˙}{F}}_{ij}$,so that

$\frac{\partial W}{\partial {F}_{ij}}={S}_{ji}$

4.      Finally, the formula for Cauchy stress follows from the equation relating ${\sigma }_{ij}$ to ${S}_{ij}$

${\sigma }_{ij}=\frac{1}{J}{F}_{ik}{S}_{kj}=\frac{1}{J}{F}_{ik}\frac{\partial W}{\partial {F}_{kj}}$

For an isotropic material, it is necessary to find derivatives of the invariants with respect to the components of F in order to compute the stress-strain function for a given strain energy density.  It is straightforward, but somewhat tedious to show that:

$\frac{\partial {I}_{1}}{\partial {F}_{ij}}=2{F}_{ij},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\frac{\partial {I}_{2}}{\partial {F}_{ij}}=2\left({I}_{1}{F}_{ij}-{B}_{ik}{F}_{kj}\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\frac{\partial {I}_{3}}{\partial {F}_{ij}}=2{I}_{3}{F}_{ji}^{-1}$

Then,

$\frac{\partial W}{\partial {F}_{ij}}=\frac{\partial U}{\partial {I}_{1}}\frac{\partial {I}_{1}}{\partial {F}_{ij}}+\frac{\partial U}{\partial {\overline{I}}_{2}}\frac{\partial {I}_{2}}{\partial {F}_{ij}}+\frac{\partial U}{\partial {I}_{3}}\frac{\partial {I}_{3}}{\partial {F}_{ij}}=2\left(\frac{\partial U}{\partial {I}_{1}}+{I}_{1}\frac{\partial U}{\partial {I}_{2}}\right){F}_{ij}-2\frac{\partial U}{\partial {I}_{2}}{B}_{ik}{F}_{kj}2{I}_{3}\frac{\partial U}{\partial {I}_{3}}{F}_{ji}^{-1}$

and

${\sigma }_{ij}=\frac{1}{\sqrt{{I}_{3}}}{F}_{ik}\frac{\partial W}{\partial {F}_{kj}}=\frac{2}{\sqrt{{I}_{3}}}\left[\left(\frac{\partial U}{\partial {I}_{1}}+{I}_{1}\frac{\partial U}{\partial {I}_{2}}\right){B}_{ij}-\frac{\partial U}{\partial {I}_{2}}{B}_{ik}{B}_{kj}\right]+2\sqrt{{I}_{3}}\frac{\partial U}{\partial {I}_{3}}{\delta }_{ij}$

When using a strain energy density of the form $\overline{U}\left({\overline{I}}_{1},{\overline{I}}_{2},J\right)$,  we will have to compute the derivatives of the invariants  with respect to the components of F in order to find

$\frac{\partial W}{\partial {F}_{ij}}=\frac{\partial \overline{U}}{\partial {\overline{I}}_{1}}\frac{\partial {\overline{I}}_{1}}{\partial {F}_{ij}}+\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}\frac{\partial {\overline{I}}_{2}}{\partial {F}_{ij}}+\frac{\partial \overline{U}}{\partial J}\frac{\partial J}{\partial {F}_{ij}}$

We find that

$\begin{array}{c}\frac{\partial J}{\partial {F}_{ij}}=J{F}_{ji}^{-1}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\frac{\partial {\overline{I}}_{1}}{\partial {F}_{ij}}=\frac{1}{{J}^{2/3}}\frac{\partial {I}_{1}}{\partial {F}_{ij}}-\frac{2{I}_{1}}{3{J}^{5/3}}\frac{\partial J}{\partial {F}_{ij}}\text{\hspace{0.17em}}=\frac{2}{{J}^{2/3}}\left({F}_{ij}-\frac{{I}_{1}}{3}\text{\hspace{0.17em}}{F}_{ji}^{-1}\right)=\frac{2}{{J}^{2/3}}{F}_{ij}-\frac{2}{3}{\overline{I}}_{1}\text{\hspace{0.17em}}{F}_{ji}^{-1}\\ \frac{\partial {\overline{I}}_{2}}{\partial {F}_{ij}}=\frac{1}{{J}^{4/3}}\frac{\partial {I}_{2}}{\partial {F}_{ij}}-\frac{4{I}_{2}}{3{J}^{7/3}}\frac{\partial J}{\partial {F}_{ij}}=\frac{2}{{J}^{4/3}}\left({I}_{1}{F}_{ij}-{B}_{ik}{F}_{kj}-\frac{2{I}_{2}}{3}{F}_{ji}^{-1}\right)=\frac{2}{{J}^{2/3}}{\overline{I}}_{1}{F}_{ij}-\frac{2}{{J}^{4/3}}{B}_{ik}{F}_{kj}-\frac{4{\overline{I}}_{2}}{3}{F}_{ji}^{-1}\end{array}$

Thus,

${\sigma }_{ij}=\frac{1}{J}{F}_{ik}\frac{\partial W}{\partial {F}_{kj}}=\frac{2}{{J}^{5/3}}\left(\frac{\partial \overline{U}}{\partial {\overline{I}}_{1}}+{\overline{I}}_{1}\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}\right){B}_{ij}-\frac{2}{3J}\left({\overline{I}}_{1}\frac{\partial \overline{U}}{\partial {\overline{I}}_{1}}+2{\overline{I}}_{2}\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}\right){\delta }_{ij}-\frac{2}{{J}^{7/3}}\frac{\partial \overline{U}}{\partial {\overline{I}}_{2}}{B}_{ik}{B}_{kj}+\frac{\partial \overline{U}}{\partial J}{\delta }_{ij}$

Next, we derive the stress-strain relation in terms of a strain energy density $\stackrel{˜}{U}\left({\lambda }_{1},{\lambda }_{2},{\lambda }_{3}\right)$ that is expressed as a function of the principal  strains.  Note first that

$\stackrel{˜}{U}\left({\lambda }_{1},{\lambda }_{2},{\lambda }_{3}\right)=U\left({I}_{1},{I}_{2},{I}_{3}\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}{I}_{1}={\lambda }_{1}^{2}+{\lambda }_{2}^{2}+{\lambda }_{3}^{2},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{I}_{2}={\lambda }_{1}^{2}{\lambda }_{3}^{2}+{\lambda }_{2}^{2}{\lambda }_{3}^{2}+{\lambda }_{1}^{2}{\lambda }_{3}^{2},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{I}_{3}={\lambda }_{1}^{2}{\lambda }_{2}^{2}{\lambda }_{3}^{2}$

so that the chain rule gives

$\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{i}}=2{\lambda }_{i}\left(\frac{\partial U}{\partial {I}_{1}}+\left({I}_{1}-{\lambda }_{i}^{2}\right)\frac{\partial U}{\partial {I}_{2}}+\frac{{I}_{3}}{{\lambda }_{i}^{2}}\frac{\partial U}{\partial {I}_{3}}\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\left(i=1,2,3\right)$

Using this and the expression that relates the stress components to the derivatives of U,

${\sigma }_{ij}=\frac{1}{\sqrt{{I}_{3}}}{F}_{ik}\frac{\partial W}{\partial {F}_{kj}}=\frac{2}{\sqrt{{I}_{3}}}\left[\left(\frac{\partial U}{\partial {I}_{1}}+{I}_{1}\frac{\partial U}{\partial {I}_{2}}\right){B}_{ij}-\frac{\partial U}{\partial {I}_{2}}{B}_{ik}{B}_{kj}\right]+2\sqrt{{I}_{3}}\frac{\partial U}{\partial {I}_{3}}{\delta }_{ij}$

we find that the principal stresses ${\sigma }_{1},{\sigma }_{2},{\sigma }_{3}$ are related to the corresponding principal stretches ${\lambda }_{1},{\lambda }_{2},{\lambda }_{3}$ (square-roots of the eigenvalues of B) through

${\sigma }_{i}=\frac{{\lambda }_{i}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{i}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\left(i=1,2,3\right)$

The spectral decomposition for B in terms of its eigenvalues ${\lambda }_{1}^{2},{\lambda }_{2}^{2},{\lambda }_{3}^{2}$ and eigenvectors ${b}^{\left(1\right)},{b}^{\left(2\right)},{b}^{\left(3\right)}$:

${B}_{ij}={\lambda }_{1}^{2}{b}_{i}^{\left(1\right)}{b}_{j}^{\left(1\right)}+{\lambda }_{2}^{2}{b}_{i}^{\left(2\right)}{b}_{j}^{\left(2\right)}+{\lambda }_{3}^{2}{b}_{i}^{\left(3\right)}{b}_{j}^{\left(3\right)}$ now allows the stress tensor to be written as

${\sigma }_{ij}=\frac{{\lambda }_{1}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{1}}{b}_{i}^{\left(1\right)}{b}_{j}^{\left(1\right)}+\frac{{\lambda }_{2}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{2}}{b}_{i}^{\left(2\right)}{b}_{j}^{\left(2\right)}+\frac{{\lambda }_{3}}{{\lambda }_{1}{\lambda }_{2}{\lambda }_{3}}\frac{\partial \stackrel{˜}{U}}{\partial {\lambda }_{3}}{b}_{i}^{\left(3\right)}{b}_{j}^{\left(3\right)}$

3.5.4 A note on perfectly incompressible materials

The preceding formulas assume that the material has some (perhaps small) compressibility $–$ that is to say, if you load it with hydrostatic pressure, its volume will change by a measurable amount.   Most rubbers strongly resist volume changes, and in hand calculations it is sometimes convenient to approximate them as perfectly incompressible.   The material model for incompressible materials is specified as follows:

The deformation must satisfy J=1 to preserve volume.

The strain energy density is therefore only a function of two invariants; furthermore, both sets of invariants defined above are identical.  We can use a strain energy density of the form $U\left({I}_{1},{I}_{2}\right)$.

Because you can apply any pressure to an incompressible solid without changing its shape, the stress cannot be uniquely determined from the strains.   Consequently, the stress-strain law only specifies the deviatoric stress ${\overline{\sigma }}_{ij}={\sigma }_{ij}-{\sigma }_{kk}{\delta }_{ij}/3$.  In problems involving quasi-static loading, the hydrostatic stress $p={\sigma }_{kk}/3$ can usually be calculated, by solving the equilibrium equations (together with appropriate boundary conditions).   Incompressible materials should not be used in a dynamic analysis, because the speed of elastic pressure waves is infinite.

The formula for stress in terms of $U\left({I}_{1},{I}_{2}\right)$ has the form

${\sigma }_{ij}=2\left[\left(\frac{\partial U}{\partial {I}_{1}}+{I}_{1}\frac{\partial U}{\partial {I}_{2}}\right){B}_{ij}-\left({I}_{1}\frac{\partial U}{\partial {I}_{1}}+2{I}_{2}\frac{\partial U}{\partial {I}_{2}}\right)\frac{{\delta }_{ij}}{3}-\frac{\partial U}{\partial {I}_{2}}{B}_{ik}{B}_{kj}\right]+p{\delta }_{ij}$

The hydrostatic stress p is an unknown variable, which must be calculated by solving the boundary value problem.

3.5.5 Specific forms of the strain energy density

Generalized Neo-Hookean solid  (Adapted from Treloar, Proc Phys Soc 60 135-44 1948)

$\overline{U}=\frac{{\mu }_{1}}{2}\left({\overline{I}}_{1}-3\right)+\frac{{K}_{1}}{2}{\left(J-1\right)}^{2}$

where ${\mu }_{1}$ and ${K}_{1}$ are material properties (for small deformations, ${\mu }_{1}$ and ${K}_{1}$ are the shear modulus and bulk modulus of the solid, respectively). Elementary statistical mechanics treatments predict that ${\mu }_{1}=NkT$, where N is the number of polymer chains per unit volume, k is the Boltzmann constant, and T is temperature.  This is a rubber elasticity model, for rubbers with very limited compressibility, and should be used with ${K}_{1}>>{\mu }_{1}$.  The stress-strain relation follows as

${\sigma }_{ij}=\frac{{\mu }_{1}}{{J}^{5/3}}\left({B}_{ij}-\frac{1}{3}{B}_{kk}{\delta }_{ij}\right)+{K}_{1}\left(J-1\right){\delta }_{ij}$

The fully incompressible limit can be obtained by setting ${K}_{1}\left(J-1\right)=p/3$ in the stress-strain law.

Generalized Mooney-Rivlin solid (Adapted from Mooney, J Appl Phys 11 582 1940)

$\overline{U}=\frac{{\mu }_{1}}{2}\left({\overline{I}}_{1}-3\right)+\frac{{\mu }_{2}}{2}\left({\overline{I}}_{2}-3\right)+\frac{{K}_{1}}{2}{\left(J-1\right)}^{2}$

where ${\mu }_{1},{\mu }_{2}$ and ${K}_{1}$ are material properties.  For small deformations, the shear modulus and bulk modulus of the solid are $\mu ={\mu }_{1}+{\mu }_{2}$ and $K={K}_{1}$.  This is a rubber elasticity model, and should be used with ${K}_{1}>>{\mu }_{1}$. The stress-strain relation follows as

${\sigma }_{ij}=\frac{{\mu }_{1}}{{J}^{5/3}}\left({B}_{ij}-\frac{1}{3}{B}_{kk}{\delta }_{ij}\right)+\frac{{\mu }_{2}}{{J}^{7/3}}\left({B}_{kk}{B}_{ij}-\frac{1}{3}{\left[{B}_{kk}^{}\right]}^{2}{\delta }_{ij}-{B}_{ik}{B}_{kj}+\frac{1}{3}{B}_{kn}{B}_{nk}{\delta }_{ij}\right)+{K}_{1}\left(J-1\right){\delta }_{ij}$

Generalized polynomial rubber elasticity potential

$\overline{U}=\sum _{i+j=1}^{N}{C}_{ij}{\left({\overline{I}}_{1}-3\right)}^{i}{\left({\overline{I}}_{2}-3\right)}^{j}+\sum _{i=1}^{N}\frac{{K}_{i}}{2}{\left(J-1\right)}^{2i}$

where ${C}_{ij}$ and ${K}_{i}$ are material properties.  For small strains the shear modulus and bulk modulus follow as $\mu =2\left({C}_{01}+{C}_{10}\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}K=2{K}_{1}$. This model is implemented in many finite element codes.  Both the neo-Hookean solid and the Mooney-Rivlin solid are special cases of the law (with N=1 and appropriate choices of ${C}_{ij}$ ).  Values of $N>2$ are rarely used, because it is difficult to fit such a large number of material properties to experimental data.

Ogden model (Ogden, Proc R Soc Lond A326, 565-84 (1972), ibid A328 567-83 (1972))

$\stackrel{˜}{U}=\sum _{i=1}^{N}\frac{2{\mu }_{i}}{{\alpha }_{i}^{2}}\left({\overline{\lambda }}_{1}^{{\alpha }_{i}}+{\overline{\lambda }}_{2}^{{\alpha }_{i}}+{\overline{\lambda }}_{3}^{{\alpha }_{i}}-3\right)+\frac{{K}_{1}}{2}{\left(J-1\right)}^{2}$

where ${\overline{\lambda }}_{i}={\lambda }_{i}/{J}^{1/3}$ , and ${\mu }_{i},{\alpha }_{i},K$ are material properties.  For small strains the shear modulus and bulk modulus follow as $\mu ={\sum }_{i=1}^{N}{\mu }_{i},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}K={K}_{1}$. This is a rubber elasticity model, and is intended to be used with ${K}_{1}>>{\mu }_{i}$.  The stress can be computed using the formulas in 3.4.3, but are too lengthy to write out in full here.

Arruda-Boyce 8 chain model (J. Mech. Phys. Solids, 41, (2) 389-412, 1992)

$\overline{U}=\mu \left\{\frac{1}{2}\left({\overline{I}}_{1}^{}-3\right)+\frac{1}{20{\beta }^{2}}\left({\overline{I}}_{1}^{2}-9\right)+\frac{11}{1050{\beta }^{4}}\left({\overline{I}}_{1}^{3}-27\right)+...\right\}+\frac{K}{2}{\left(J-1\right)}^{2}$

where $\mu ,\beta ,K$ are material properties.  For small deformations $\mu ,K$ are the shear and bulk modulus, respectively. This is a rubber elasticity model, so $K>>\mu$.    The potential was derived by calculating the entropy of a simple network of long-chain molecules, and the series is the result of a Taylor expansion of an inverse Langevin function.  The reference provided lists more terms if you need them.  The stress-strain law is

${\sigma }_{ij}=\frac{\mu }{{J}^{5/3}}\left(1+\frac{{B}_{kk}}{5{J}^{2/3}{\beta }^{2}}+\frac{33{\left({B}_{kk}\right)}^{2}}{525{\beta }^{4}{J}^{4/3}}+...\right)\left({B}_{ij}-\frac{{B}_{kk}}{3}{\delta }_{ij}\right)+K\left(J-1\right){\delta }_{ij}$

Ogden-Storakers hyperelastic foam

$\stackrel{˜}{U}=\sum _{i=1}^{N}\frac{2{\mu }_{i}}{{\alpha }_{i}^{2}}\left({\lambda }_{1}^{{\alpha }_{i}}+{\lambda }_{2}^{{\alpha }_{i}}+{\lambda }_{3}^{{\alpha }_{i}}-3+\frac{1}{{\beta }_{i}}\left({J}^{-{\alpha }_{i}{\beta }_{i}}-1\right)\right)$

where ${\mu }_{i},{\alpha }_{i},{\beta }_{i}$ are material properties.   For small strains the shear modulus and bulk modulus follow as $\mu ={\sum }_{i=1}^{N}{\mu }_{i},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}K={\sum }_{i=1}^{N}2{\mu }_{i}\left({\beta }_{i}+1/3\right)$.   This is a foam model, and can model highly compressible materials.  The shear and compression responses are coupled.

Blatz-Ko foam rubber

$U\left({I}_{1},{I}_{2},{I}_{3}\right)=\frac{\mu }{2}\left(\frac{{I}_{1}}{{I}_{2}}+2\sqrt{{I}_{3}}\right)$

where $\mu$ is a material parameter corresponding to the shear modulus at infinitesimal strains. Poisson’s ratio for such a material is 0.25.

3.5.6 Calibrating nonlinear elasticity models

To use any of these constitutive relations, you will need to determine values for the material constants.  In some cases this is quite simple (the incompressible neo-Hookean material only has 1 constant!); for models like the generalized polynomial or Ogden’s it is considerably more involved.

Conceptually, however, the procedure is straightforward.  You can perform various types of test on a sample of the material, including simple tension, pure shear, equibiaxial tension, or volumetric compression. It is straightforward to calculate the predicted stress-strain behavior for the specimen for each constitutive law.  The parameters can then be chosen to give the best fit to experimental behavior.

Here are some guidelines on how best to do this:

1.      When modeling the behavior of rubber under ambient pressure, you can usually assume that the material is nearly incompressible, and don’t need to characterize response to volumetric compression in detail.  For the rubber elasticity models listed above, you can take ${K}_{1}\approx {10}^{5}$ MPa. To fit the remaining parameters, you can assume the material is perfectly incompressible.

2.      If rubber is subjected to large hydrostatic stress (>100 MPa) its volumetric and shear responses are strongly coupled. Compression increases the shear modulus, and high enough pressure can even induce a glass transition (see, e.g. D.L. Quested, K.D. Pae, J.L. Sheinbein and B.A. Newman, J. Appl. Phys, 52, (10) 5977 (1981)).  To account for this, you would have to use one of the foam models: in the rubber models the volumetric and shear responses are decoupled. You would also have to determine the material constants by testing the material under combined hydrostatic and shear loading.

3.      For the simpler material models, (e.g. the neo-Hookean solid, the Mooney-Rivlin material, or the Arruda-Boyce model, which contain only two material parameters in addition to the bulk modulus) you can estimate material parameters by fitting to the results of a uniaxial tension test.  There are various ways to actually do the fit $–$ you could match the small-strain shear modulus to experiment, and then select the remaining parameter to fit the stress-strain curve at a larger stretch.  Least-squared fits are also often used.  However, models calibrated in this way do not always predict material behavior under multiaxial loading accurately.

4.      A more accurate description of material response to multiaxial loading can be obtained by fitting the material parameters to multiaxial tests.  To help in this exercise, the nominal stress (i.e. force/unit undeformed area) $–$v- extension predicted by several constitutive laws are listed in the table below (assuming perfectly incompressible behavior, as suggested in 1.)

 Uniaxial Tension Biaxial Tension Pure Shear Invariants $\begin{array}{l}{l}_{1}/{L}_{1}=\lambda \\ {l}_{2}/{L}_{2}={l}_{3}/{L}_{3}={\lambda }^{-1/2}\\ {I}_{1}={\lambda }^{2}+2{\lambda }^{-1}\end{array}$ $\begin{array}{l}{l}_{1}/{L}_{1}={l}_{2}/{L}_{2}=\lambda \\ {l}_{3}/{L}_{3}={\lambda }^{-2}\\ {I}_{1}=2{\lambda }^{2}+{\lambda }^{-4}\end{array}$ $\begin{array}{l}{l}_{1}/{L}_{1}=\lambda \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{l}_{2}/{L}_{2}=1\\ {l}_{3}/{L}_{3}={\lambda }^{-1}\\ {I}_{1}=1+{\lambda }^{2}+{\lambda }^{-2}\end{array}$ Neo-Hookean $\begin{array}{l}{S}_{1}={\mu }_{1}\left(\lambda -{\lambda }^{-2}\right)\\ {S}_{2}={S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}={S}_{2}={\mu }_{1}\left(\lambda -{\lambda }^{-5}\right)\\ {S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}={\mu }_{1}\left(\lambda -{\lambda }^{-3}\right)\\ {S}_{2}={\mu }_{1}\left({\lambda }^{-1}-{\lambda }^{-3}\right)\\ {S}_{3}=0\end{array}$ Mooney-Rivlin $\begin{array}{l}{S}_{1}={\mu }_{1}\left(\lambda -{\lambda }^{-2}\right)\text{\hspace{0.17em}}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}+{\mu }_{2}\left(1-{\lambda }^{-3}\right)\\ {S}_{2}={S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}={S}_{2}={\mu }_{1}\left(\lambda -{\lambda }^{-5}\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}+{\mu }_{2}\left({\lambda }^{3}-{\lambda }^{-3}\right)\\ {S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}=\left({\mu }_{1}+{\mu }_{2}\right)\left(\lambda -{\lambda }^{-3}\right)\\ {S}_{2}={\mu }_{1}\left(1-{\lambda }^{-1}\right)+{\mu }_{2}\left({\lambda }^{2}-1\right)\\ {S}_{3}=0\end{array}$ Arruda-Boyce $\begin{array}{l}{S}_{1}=C\left(\lambda -{\lambda }^{-2}\right)\\ {S}_{2}={S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}={S}_{2}=C\left(\lambda -{\lambda }^{-5}\right)\\ {S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}=C\left(\lambda -{\lambda }^{-3}\right)\\ {S}_{2}=C\left(1-{\lambda }^{-2}\right)\end{array}$ $C=\mu \left(1+\frac{{I}_{1}}{5{\beta }^{2}}+\frac{33{I}_{1}^{2}}{525{\beta }^{4}}\right)$ Ogden $\begin{array}{l}{S}_{1}=\sum _{n}{\mu }_{n}\left({\lambda }^{{\alpha }_{n}}-{\lambda }^{-{\alpha }_{n}/2}\right)/\lambda \\ {S}_{2}={S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}={S}_{2}=\sum _{n}{\mu }_{n}\left({\lambda }^{{\alpha }_{n}}-{\lambda }^{-2{\alpha }_{n}}\right)/\lambda \\ {S}_{3}=0\end{array}$ $\begin{array}{l}{S}_{1}=\sum _{n}{\mu }_{n}\left({\lambda }^{{\alpha }_{n}}-{\lambda }^{-{\alpha }_{n}}\right)/\lambda \\ {S}_{2}=\sum _{n}{\mu }_{n}\left(1-{\lambda }^{-{\alpha }_{n}}\right)\end{array}$

3.5.7 Representative values of material properties for rubbers

The properties of rubber are strongly sensitive to its molecular structure, and for accurate predictions you will need to obtain experimental data for the particular material you plan to use.    As a rough guide, the experimental data of Treloar  (Trans. Faraday Soc. 40, 59.1944) for the behavior of vulcanized rubber under uniaxial tension, biaxial tension, and pure shear is shown in the picture.  The solid lines in the figure show the predictions of the Ogden model (which gives the best fit to the data).

Material parameters fit to this data for several constitutive laws are listed below.

 Neo-Hookean ${\mu }_{1}=0.4$ MPa Mooney-Rivlin ${\mu }_{1}=0.39$ MPa,  ${\mu }_{2}=0.015$ MPa Arruda-Boyce ${\mu }_{1}=0.4$ MPa, $\beta =10$ Ogden ${\mu }_{1}=0.62$ MPa, ${\alpha }_{1}=1.3$ ${\mu }_{2}=0.00118$ MPa, ${\alpha }_{2}=5$ ${\mu }_{3}=-0.00981$ MPa, ${\alpha }_{3}=-2$