Search Continuum Mechanics Website

Keep in mind that auto-summation of duplicated subscripts will not apply to the equations on this page. Hopefully, the reasons for this will be obvious!

\[ \boldsymbol{\sigma}^\text{PK2} = \rho_o \, {\partial \Psi \over \partial {\bf E}^\text{el} } \]

The Helmholtz free energy contains thermal energy and mechanical strain energy. But in most every discussion of Mooney-Rivlin coefficients, the thermal part is neglected, leaving only the mechanical part, \(W\). (Actually, \(W\) is declared to represent \(\rho_o \Psi\), not just \(\Psi\)). Second, since all of the deformation of a hyperelastic material is elastic by definition, it is sufficient to write \({\bf E}^\text{el}\) simply as \({\bf E}\). This gives

\[ \boldsymbol{\sigma}^\text{PK2} = {\partial W \over \partial {\bf E} } \]

But there is a challenge with this general approach. It is the determination of off-diagonal (shear) terms. As with the shear terms in Hooke's Law, they are not independent of the normal terms, but must be consistent with coordinate transformations that transform normal components into shears and vice-versa. And as with Hooke's Law, the resolution is to define the material behavior for the principal values and rely on coordinate transformations to give the appropriate corresponding behavior of the shear terms.

\[ \sigma_i^\text{PK2} = {\partial W \over \partial E_i } \]

But alas, even this is not quite what is done in Mooney-Rivlin models. Instead, derivatives are taken with respect to

\[ \sigma_i^\text{PK2} = {\partial W \over \partial E_i } = {\partial W \over \partial \lambda_i } {\partial \lambda_i \over \partial E_i} \]

Recall that a principal Green strain equals

\[ E_i = {\Delta L_i \over L_o} + {1 \over 2} \left( {\Delta L_i \over L_o} \right)^2 \]

And insert the relationship \(\lambda_i - 1 = \Delta L_i / L_o\) to get

\[ E_i = (\lambda_i - 1) + {1 \over 2} (\lambda_i - 1)^2 \]

So \(\partial E_i / \partial \lambda_i = \lambda_i\) and therefore

\[ {\partial \lambda_i \over \partial E_i } = {1 \over \lambda_i} \]

Substituting this into the equation containing \(\sigma^\text{PK2}\) and \(W\) gives

\[ \lambda_i \sigma_i^\text{PK2} = {\partial W \over \partial \lambda_i } \]

This equation by itself is not helpful. But recall from the energetic conjugates page that \(\sigma_i = (1 + \epsilon_i)^2 \sigma_i^\text{PK2}\) and therefore \(\sigma_i = \lambda_i^2 \sigma_i^\text{PK2}\). So if we simply multiply both sides of the above equation by another \(\lambda_i\), we will have something much more useful.

\[ \sigma_i = \lambda_i^2 \sigma_i^\text{PK2} = \lambda_i {\partial W \over \partial \lambda_i } \]

So each principal Cauchy stress can be related directly to the partial derivative of strain energy through

\[ \sigma_i = \lambda_i {\partial W \over \partial \lambda_i } \]

In this case, the stretch ratio is in a principal strain direction, and will give a principal Cauchy stress in that same principal direction.

\[ W = {1 \over 2} E (\lambda - 1)^2 \]

Then the Cauchy stress is

\[ \begin{eqnarray} \sigma & = & \lambda {\partial \left[ {1 \over 2} E (\lambda - 1)^2 \right] \over \partial \lambda } \\ \\ \\ & = & \lambda E (\lambda - 1) \end{eqnarray} \]

Recall that \(\sigma^\text{Eng} = \sigma / \lambda\), so dividing through by \(\lambda\) gives

\[ \sigma^\text{Eng} = E (\lambda - 1) \]

and \(\lambda = 1 + \epsilon^\text{Eng}\), so \(\lambda -1 = \epsilon^\text{Eng}\). And this leads to

\[ \sigma^\text{Eng} = E \, \epsilon^\text{Eng} \]

That was fun because of its simplicity and the familiar result we obtained, but it is in fact not very useful. Why? Because we have over-simplified things so much that we have lost the ability to incorporate Poisson Effects. So our proposed function for the strain energy, \(W\), in terms of \(\lambda\) could not be used in, say, a general 3-D finite element program. We will see how to overcome these limitations in the following sections.

The invariants are of the product of the deformation gradient with its transpose, \({\bf F}\cdot {\bf F}^T\). This result is symmetric because the product of any matrix with its transpose is always symmetric. This makes it possible to calculate principal values, which will become important shortly. Recall from the polar decomposition page, that \({\bf F} \cdot {\bf F}^T = {\bf V} \cdot {\bf V}^T\) because

\[ {\bf F} \cdot {\bf F}^T \;\; = \;\; \left( {\bf V} \cdot {\bf R} \right) \cdot \left( {\bf V} \cdot {\bf R} \right)^T \;\; = \;\; {\bf V} \cdot {\bf R} \cdot {\bf R}^T \! \cdot {\bf V}^T \;\; = \;\; {\bf V} \cdot {\bf V}^T \]

And that the principal values of \({\bf V}\) are

\[ {\bf V}_\text{Pr} = \left[ \matrix{ (L_F / L_o)_1 & 0 & 0 \\ 0 & (L_F / L_o)_2 & 0 \\ 0 & 0 & (L_F / L_o)_3 } \right] \]

\[ {\bf V}_\text{Pr} = \left[ \matrix{ \lambda_1 & 0 & 0 \\ 0 & \lambda_2 & 0 \\ 0 & 0 & \lambda_3 } \right] \]

Therefore, the principal values of \({\bf V} \cdot {\bf V}^T\) are

\[ ( {\bf V} \cdot {\bf V}^T)_\text{Pr} = \left[ \matrix{ \lambda^2_1 & 0 & 0 \\ 0 & \lambda^2_2 & 0 \\ 0 & 0 & \lambda^2_3 } \right] \]

\({\bf V} \cdot {\bf V}^T\) is symmetric, so it does have three invariants. (Recall that invariants are independent of coordinate transforms.)

\[ \begin{eqnarray} I_1 & = & \lambda^2_1 + \lambda^2_2 + \lambda^2_3 \\ \\ \\ I_2 & = & \lambda^2_1 \, \lambda^2_2 + \lambda^2_2 \, \lambda^2_3 + \lambda^2_1 \, \lambda^2_3 \\ \\ I_3 & = & \, \text{det}(...\!) = \; \lambda^2_1 \, \lambda^2_2 \, \lambda^2_3 = \left( {V_F \over V_o} \right)^{\!2} = J^2 \end{eqnarray} \]

Keep in mind that these equations for the invariants are special cases because they are of a symmetric matrix

\(I_3\) merits a few additional comments. First, it is the square of the Jacobian, \(J\), because \(\lambda_1 \, \lambda_2 \, \lambda_3 = \left( {V_F \over V_o} \right)\), and therefore \(I_3\) is the square of the volume ratio. Also, for incompressible materials, \(I_3\) always equals 1, and therefore \(\lambda_1 \, \lambda_2 \, \lambda_3 = 1\).

The incompressibility relationship is often used to simplify the terms of the invariants, particularly \(I_2\). For example, note that the first term in \(I_2\), \(\lambda^2_1 \, \lambda^2_2\), is equal to \(1 / \lambda^2_3\). This can also be applied to the two remaining terms to give the following set of equations.

\[ \begin{eqnarray} I_1 & = & \lambda^2_1 + \lambda^2_2 + \lambda^2_3 \\ \\ \\ I_2 & = & {1 \over \lambda^2_1} + {1 \over \lambda^2_2} + {1 \over \lambda^2_3} \\ \\ \\ I_3 & = & ( \lambda_1 \, \lambda_2 \, \lambda_3 )^2 = 1 \end{eqnarray} \]

Finally, it is worth noting that the invariants do

\[ W = \sum_i \sum_j C_{ij} \left( I_1 - 3 \right)^i \left( I_2 - 3 \right)^j + D \, (J - 1)^2 \]

Note that the series is not a function of \(I_3\) because it is a constant value, 1. The constants, \(C_{ij}\) and \(D\), will be determined by curve-fitting measured stress-strain curves to the derivative of the equation. The number of terms in the expansion is determined by the application's accuracy requirements.

As an example, the first few terms of the series are

\[ W = C_{10} \left( I_1 - 3 \right) + C_{01} \left( I_2 - 3 \right) + C_{11} \left( I_1 - 3 \right) \left( I_2 - 3 \right) + C_{20} \left( I_1 - 3 \right)^2 + \, ... \, + D \, (J - 1)^2 \]

\[ \sigma_1 = \lambda_1 {\partial W \over \partial \lambda_1} = \lambda_1 \left( {\partial W \over \partial I_1} {\partial I_1 \over \partial \lambda_1} + {\partial W \over \partial I_2} {\partial I_2 \over \partial \lambda_1} + {\partial W \over \partial J} {\partial J \over \partial \lambda_1} \right) \]

The derivatives of the strain energy with respect to the invariants, and \(J\), are

\[ {\partial W \over \partial I_1 } = C_{10} + C_{11} \left( I_2 - 3 \right) + 2 C_{20} \left( I_1 - 3 \right) + \, ... \qquad \qquad {\partial W \over \partial I_2 } = C_{01} + C_{11} \left( I_1 - 3 \right) + \, ... \]

\[ {\partial W \over \partial J } = 2 \, D \, (J-1) \]

And the derivatives of the invariants, and \(J\), with respect to \(\lambda_1\) are

\[ {\partial I_1 \over \partial \lambda_1} = 2 \lambda_1 \qquad \qquad {\partial I_2 \over \partial \lambda_1} = - { 2 \over \lambda^3_1 } \qquad \qquad {\partial J \over \partial \lambda_1} = \lambda_2 \, \lambda_3 \]

All of these terms can be combined to give polynomials relating stretch ratios to principal stresses, with coefficients such as \(C_{10}\), \(C_{01}\), \(C_{11}\), and \(C_{20}\) that are determined from curve-fitting these equations to experimental data.

\[ W = C_{10} ( I_1 - 1 ) + D (J - 1)^2 \]

Taking derivatives with respect to \(\lambda_1\) to obtain \(\sigma_1\) gives

\[ \begin{eqnarray} \sigma_1 & = & \lambda_1 \left( {\partial W \over \partial I_1} {\partial I_1 \over \partial \lambda_1} + {\partial W \over \partial J} {\partial J \over \partial \lambda_1} \right) \\ \\ \\ \\ & = & \lambda_1 \big[ C_{10} ( 2 \lambda_1 ) + 2 \, D \, (J - 1) (\lambda_2 \lambda_3) \big] \\ \\ \\ \\ & = & 2 \, C_{10} \lambda_1^2 + 2 \, D \, (J - 1) (\lambda_1 \lambda_2 \lambda_3) \end{eqnarray} \]

Once again, the \(D\) term requires some explanation. First, since \(J = V_F / V_o\), then \((J-1) = \Delta V / V_o = \epsilon_\text{Vol}\). So we will replace \((J-1)\) with the volumetric strain. Second, the product \((\lambda_1 \lambda_2 \lambda_3)\) is also \(J\) or \(V_F / V_o\). We will assume that this value is so close to 1 for an incompressible material that it can be dropped.

An alert reader will notice an apparent contradiction between the dropping of \(J\) and the retaining of \((J-1)\) as \(\epsilon_\text{Vol}\). It is trully suspicious at first glance. The justification is that \((J-1)\) represents the difference of two relative large numbers to give a small result that is nevertheless critical to the build-up of hydrostatic stresses. This difference may be small compared to 1, but it is not small when compared to a reference volumetric strain of zero. In contrast, once the small value is obtained, the effect of multiplying it by \(J\) again is negligible.

So the equation can be written as

\[ \sigma_1 = 2 \, C_{10} \lambda_1^2 + 2 \, D \, \epsilon_\text{Vol} \]

Likewise, the equations for the other two principal Cauchy stresses are

\[ \sigma_2 = 2 \, C_{10} \lambda_2^2 + 2 \, D \, \epsilon_\text{Vol} \] \[ \sigma_3 = 2 \, C_{10} \lambda_3^2 + 2 \, D \, \epsilon_\text{Vol} \]

The constants are determined by curve-fitting the equations to laboratory test data. They can then be input into finite element programs to describe the material behavior.

Note that \(2D\) serves as the bulk modulus, \(K\), in the equations. Therefore, the proper value for \(D\) is \(K/2\).

\[ \begin{eqnarray} \sigma_1 & = & \lambda_1 \Big\{ \big[ C_{10} + 2 \, C_{20} (I_1 - 3) \big] \big[ 2 \lambda_1 \big] + \big[ C_{01} \big] \big[ {-2 / \lambda_1^3 } \big] \Big\} + 2 \, D \, \epsilon_\text{Vol} \\ \\ \\ & = & 2 \, C_{10} \lambda_1^2 - 2 \, C_{01} {1 \over \lambda_1^2 } + 4 \, C_{20} \lambda_1^2 (\lambda_1^2 + \lambda_2^2 + \lambda_3^2 - 3) + 2 \, D \, \epsilon_\text{Vol} \end{eqnarray} \]

Likewise, the equations for \(\sigma_2\) and \(\sigma_3\) are

\[ \begin{eqnarray} \sigma_2 & = & 2 \, C_{10} \lambda_2^2 - 2 \, C_{01} {1 \over \lambda_2^2 } + 4 \, C_{20} \lambda_2^2 (\lambda_1^2 + \lambda_2^2 + \lambda_3^2 - 3) + 2 \, D \, \epsilon_\text{Vol} \\ \\ \\ \sigma_3 & = & 2 \, C_{10} \lambda_3^2 - 2 \, C_{01} {1 \over \lambda_3^2 } + 4 \, C_{20} \lambda_3^2 (\lambda_1^2 + \lambda_2^2 + \lambda_3^2 - 3) + 2 \, D \, \epsilon_\text{Vol} \\ \end{eqnarray} \]

Don't forget that the stresses here are Cauchy stresses, not engineering stresses.

Keep in mind that there is no physical insight built into these equations! They are merely a collection of polynomials that happen to curve this way and that, in ways that are similar to the stress-strain curves one gets from actual laboratory tests on rubber test specimens. The equations also possess one additional key critical property. As will be seen shortly, they predict higher principal stresses for equibiaxial tension than plane strain tension, and then higher again against uniaxial tension.

Developing Mooney-Rivlin equations for these special cases appears to be a very popular thing to do. There are a couple reasons for this. First, the deformation modes are (relatively) easily accomplished in lab tests, at least the uniaxial tension case is. The second reason is that the equations can be manipulated to eliminate the bulk modulus, which makes life a little easier. This will be demonstrated next.

The key to eliminating the bulk modulus for these three special cases is to recognize that \(\sigma_3 = 0\) in each. Therefore, one can subtract two equations from each other to accomplish \(\sigma_1 - \sigma_3\) as follows.

\[ \sigma_1 - \sigma_3 = 2 \, C_{10} (\lambda_1^2 - \lambda_3^2) + 2 \, C_{01} \left( {1 \over \lambda_3^2 } - {1 \over \lambda_1^2} \right) + 4 \, C_{20} (\lambda_1^2 - \lambda_3^2) (\lambda_1^2 + \lambda_2^2 + \lambda_3^2 - 3) \]

But \(\sigma_3 = 0\), so

\[ \sigma_1 = 2 \, C_{10} (\lambda_1^2 - \lambda_3^2) + 2 \, C_{01} \left( {1 \over \lambda_3^2 } - {1 \over \lambda_1^2} \right) + 4 \, C_{20} (\lambda_1^2 - \lambda_3^2) (\lambda_1^2 + \lambda_2^2 + \lambda_3^2 - 3) \]

An alternative popular form is obtained by dividing through by \(\lambda_1\) to obtain engineering stress. It is this form that is often used to determine the Mooney-Rivlin coefficients by fitting the equation(s) to experimental data.

\[ \sigma_1^\text{Eng} = 2 \, C_{10} \left( \lambda_1 - {\lambda_3^2 \over \lambda_1} \right) + 2 \, C_{01} \left( {1 \over \lambda_1 \lambda_3^2 } - {1 \over \lambda_1^3} \right) + 4 \, C_{20} \left( \lambda_1 - {\lambda_3^2 \over \lambda_1} \right) (\lambda_1^2 + \lambda_2^2 + \lambda_3^2 - 3) \]

\[ \lambda_1 = \lambda \qquad \qquad \lambda_2 = \lambda_3 = {1 \over \sqrt{\lambda} } \]

Substituting these into the equation for engineering stress gives

\[ \sigma_\text{Uniaxial}^\text{Eng} = 2 C_{10} \left( \lambda - {1 \over \lambda^2 } \right) + 2 C_{01} \left( 1 - {1 \over \lambda^3} \right) + 4 \, C_{20} \left( \lambda - {1 \over \lambda^2 } \right) \left( \lambda^2 + {2 \over \lambda} - 3 \right) \]

\[ \lambda_1 = \lambda \qquad \qquad \lambda_2 = 1 \qquad \qquad \lambda_3 = {1 \over \lambda } \]

And the equation for stress is

\[ \sigma_{\begin{eqnarray}\text{Plane}\;\;\, \\ \text{Strain}\;\; \\ \text{Tension}\end{eqnarray}}^\text{Eng} = \; 2 \, C_{10} \left( \lambda - {1 \over \lambda^3 } \right) + 2 \, C_{01} \left( \lambda - {1 \over \lambda^3} \right) + 4 \, C_{20} \left( \lambda - {1 \over \lambda^3 } \right) \left( \lambda^2 + { 1 \over \lambda^2 } - 2 \right) \]

\[ \lambda_1 \; = \; \lambda_2 \; = \; \lambda \qquad \qquad \lambda_3 \; = \; {1 \over \lambda^2 } \]

And the stress is

\[ \sigma_\text{Equibiaxial}^\text{Eng} = 2 C_{10} \left( \lambda - {1 \over \lambda^5 } \right) + 2 C_{01} \left( \lambda^3 - {1 \over \lambda^3} \right) + 4 \, C_{20} \left( \lambda - {1 \over \lambda^5 } \right) \left( 2 \lambda^2 + {1 \over \lambda^4} - 3 \right) \]

For the case where \(C_{01} = 1\) and all other constants are zero...

For the case where \(C_{20} = 1\) and all other constants are zero...

\[ \sigma^\text{Eng} = 2 \, C_{10} \left( \lambda - {1 \over \lambda^2 } \right) + 2 \, C_{01} \left( 1 - {1 \over \lambda^3} \right) + 4 \, C_{20} \left( \lambda - {1 \over \lambda^2 } \right) \left( \lambda^2 + {2 \over \lambda} - 3 \right) \]

Insert \(\lambda = 1.10\) into the equation to get

\[ \sigma^\text{Eng}_\text{10%} = 0.547 \, C_{10} + 0.497 \, C_{01} + 0.031 C_{20} \]

And divide by 0.10 to get MA10.

\[ \text{MA10} \; = \; 5.47 \, C_{10} + 4.97 \, C_{01} + 0.31 \, C_{20} \]

And MA100 is

\[ \text{MA100} \; = \; 3.5 \, C_{10} + 1.75 \, C_{01} + 14 \, C_{20} \]

For reference, the slope at zero strain is the following

\[ E_{\epsilon=0} \; = \; \left. {d \sigma \over d \lambda}\right|_{\lambda=1} \; = \; 6 \, ( C_{10} + C_{01} ) \]

The \(C_{20}\) term has zero slope at \(\epsilon = 0\).

Hello,

I want to use Mooney-Rivlin parametrs from the paper publish in the internet - they received parameters :

C10: 2.7507e-2

C01: 1.3077e-1

C11: 6.8027e-2

D: 4.8694e-2

My problem is :

Authors of these paper didn't mention what units of stress and strain they used to obtained the C10, C01, C11, D parameters. Is the units of input data of stress-strain important?