Hyperelastic Materials

Hyperelastic materials are the special class of materials which tends to respond elastically when they are subjected to very large strains. They shows both a nonlinear material behavior as well as large shape changes. They posses certain characteristics:

1. They can undergo large elastic deformations in order of around 100 to 700 % which is fully recoverable i.e. initial shape is recovered when load is removed.
2. They are also nearly incompressible. Which means that when they are subjected to huge loads, they change there shape but overall volume remains almost constant.
3. They shows a highly nonlinear stress-strain relation.
4. The material softens and then become stiffer again when applied to tension. Whereas, under compression they have a quite stiff response.

All the polymers tends to show the hyperelastic behaviour such as elastomers, rubbers, sponges and other soft flexible materials.
Hyperelastic materials are mostly used where high flexibility on a long run is required under large loads. The typical examples of there uses are as elastomeric pads in brigdes, rail pads, car door seal, car tyres etc.

Hyperelasticity

In finite element analysis, hyperelasticity theory is used to represent the non-linear response of hyperelastic materials at large strains. Hyperlasticity is popular due to it’s ease of use in finite element models. Normally stress-strain curve data from experiments is used to find the constants of theoretical models to fit the material response. There are some good available choices of hyperelasticity models which are available in SimScale platform including:

1. Neo Hookean (1^st oFrder reduced polynomial).
2. Mooney Rivlin (1^st order polynomial).
3. Signorini.
4. Yeoh (3^rd order reduced polynomial).
5. 2^nd order reduced polynomial.
6. 2^nd order polynomial.
7. Ogden up to order 3.
8. Arruda-Boyce.

Calculating stress-strain relation from strain energy density

The stress-strain relation for hyperelastic materials is normally calculated through strain energy density function. The formulation of which is shown below:
Suppose a solid is subjected to a displacement which can be represented by u_i(x_k). Than, the deformation gradient F and it’s Jacobian J (the total change in volume at point) can be represeted by:
$$\mathbf{F} = F_{ij} = \delta_{ij} + \frac {\partial u_i}{\partial x_j}\ , J = det(\mathbf{F})\quad$$
For simplicity,
$${\mathbf{\bar F}} = J^{\frac {-1}{3}} \mathbf{F}$$
From this, we have a left Cauchy-Green deformation tensor,
\(\bar B\) = \(\bar F\bar F\)^T B_ij = F_ikF_jk
Invariants of B are represented as:
$$\bar I_1 = trace \ \mathbf{\bar B} = \bar B_{kk} \\
\bar I_2 = \frac {1}{2} (\bar I_1^2 – trace \ (\mathbf{\bar B}.\mathbf{\bar B)}) = \frac {1}{2} (\bar I_1^2 – \bar B_{ik}\bar B_{ki})\\
\bar I_3 = det(\mathbf{B}) = J^2$$
The invariant \(\bar I_3\) is a volumetric invariant (zero for incompressible case).
Furthermore,
$$\delta {\bar I_1} = 2\mathbf{\bar B}\ :\ \delta \mathbf{\epsilon_v}\\
\delta {\bar I_2} = 2(\bar I_2\mathbf{\bar B} – \mathbf{\bar B}.\mathbf{\bar B})\ :\ \delta \mathbf{\epsilon_v}\\
J = J \delta \mathbf{\epsilon_v}\\ (where,\ \delta \mathbf{\epsilon_v}\ is\ a\ virtual\ deviatoric\ strain\ rate)$$
Variation of the strain energy density \(U = U(\bar I_1,\bar I_2,\bar J)\) for isotropic and compressible materials can be represented as,
$$\delta U = \frac {\partial U}{\partial \bar I_1} \delta \bar I_1 + \frac {\partial U}{\partial \bar I_2} \delta \bar I_2 + \frac {\partial U}{\partial \bar J} \delta \bar J$$
Using \(\bar I_1\), \(\bar I_2\) and J as defined above,
$$\delta U = 2 \left[\left(\frac {\partial U}{\partial \bar I_1} + \bar I_1 \frac {\partial U}{\partial \bar I_2}\right)\mathbf{\bar B}
– \frac {\partial U}{\partial \bar I_2} \mathbf{\bar B}.\mathbf{\bar B}\right]\ :\ \mathbf{\epsilon_v} + J\frac {\partial U}{\partial \bar J} \delta \epsilon^v \\
(where,\ \delta \epsilon^v\ is\ a\ virtual\ volumetric\ strain\ rate)$$
This variation of strain energy density is used to deduce the stress-strain law in case of hyperelastic materials, Therefore one can write:
$$\mathbf{\sigma} = 2\ DEV\ \left[\left(\frac {\partial U}{\partial \bar I_1} + \bar I_1 \frac {\partial U}{\partial \bar I_2}\right)\mathbf{\bar B}
– \frac {\partial U}{\partial \bar I_2} \mathbf{\bar B}.\mathbf{\bar B}\right] + J\frac {\partial U}{\partial \bar J} \mathbf{I}$$
For the incompressible cases (e.g. rubber), the last term containing J is neglected.

Particular strain energy density forms of the available SimScale hyperelasticity models

Before giving the appropriate material parameters to define a specific hyperelastic materials, one should know the strain energy density forms of the hyperelasticity models. In this section, we have mentioned the strain energy density forms of all the available hyperelasticity models on SimScale platform (already mentioned in Hyperelasticity).

Reduced polynomial

The basic strain energy density form of reduced polynomial is represented as:
$$U = \displaystyle\sum_{i=1}^{N} C_{i0}(\bar I_1 – 3)^i + \displaystyle\sum_{i=1}^{N} \frac {1}{D_i}(J – 1)^{2i}$$
where, C_i0(i = 1, 2, 3) are the constants that needs to be deduced in order to get the appropriate material behaviour. For the compressibility of a material, one have to give the value of D_i(i = 1, 2, 3). As a first approximation, one can use a basic formulation for elastic case:
$$D_i = \frac {3(1 – 2 \nu)}{G_o (1 + \nu)}\quad with,\ G_o = 2 C_{10}$$
If material is incompressible, than one can take the value of ν nearest to 5.0 i.e. ν = 0.499.

Neo Hookean (reduced polynomial of order 1)

For N = 1 in above formulation, one obtains the most basic form known as Neo Hookean, represented as:
$$U = C_{10}(\bar I_1 – 3) + \frac {1}{D_1}(J – 1)^2$$

Reduced polynomial of order 2

For N = 2 in above formulation, one obtains:
$$U = C_{10}(\bar I_1 – 3) + C_{20}(\bar I_1 – 3)^2 + \frac {1}{D_1}(J – 1)^2 + \frac {1}{D_2}(J – 1)^4$$

Yeoh (reduced polynomial of order 3)

For N = 3 in above formulation, one obtains Yeoh, represented as:
$$U = C_{10}(\bar I_1 – 3) + C_{20}(\bar I_1 – 3)^2 + C_{30}(\bar I_1 – 3)^3 + \frac {1}{D_1}(J – 1)^2 + \frac {1}{D_2}(J – 1)^4 + \frac {1}{D_3}(J – 1)^6$$

Polynomial

The basic strain energy density form of polynomial is represented as:
$$U = \displaystyle\sum_{i+j=1}^{N}C_{ij}(\bar I_1 – 3)^i(\bar I_2 – 3)^j + \displaystyle\sum_{i=1}^{N} \frac {1}{D_i}(J – 1)^{2i}$$
As a first approximation for D_i, same eq. eq2 can be used with G_o = 2(C₁₀ + C₀₁).

Mooney Rivlin (polynomial of order 1)

For N = 1 in above formulation, one obtains the enhancement to Neo Hookean form known as Mooney Rivlin, represented as:
$$U = C_{10}(\bar I_1 – 3) + C_{01}(\bar I_2 – 3) + \frac {1}{D_1}(J – 1)^2$$

Polynomial of order 2

For N = 2 in above formulation, one obtains:
$$U = C_{10}(\bar I_1 – 3) + C_{01}(\bar I_2 – 3) + C_{20}(\bar I_1 – 3)^2 + C_{02}(\bar I_2 – 3)^2 + C_{11}(\bar I_1 – 3)(\bar I_2 – 3)
+ \frac {1}{D_1}(J – 1)^2 + \frac {1}{D_2}(J – 1)^4$$

Signorini

Strain energy density form for Signorini is represented as:
$$U = C_{10}(\bar I_1 – 3) + C_{01}(\bar I_2 – 3) + C_{20}(\bar I_1 – 3)^2 + \frac {1}{D_1}(J – 1)^2$$

Ogden

Compared to the previous forms, Ogden is represented by the principal stretches \(\bar \lambda_i (i = 1,2,3)\) compared to the deviatoric strain invariants \(\bar I_i (i = 1,2,3)\). The basic strain energy density form of Ogden is represented as:
$$U = \displaystyle\sum_{i=1}^{N} \frac {2 \mu_i}{\alpha_{i}^{2}}(\bar \lambda_{1}^{\alpha_i} + \bar \lambda_{2}^{\alpha_i} + \bar \lambda_{3}^{\alpha_i} – 3)
+ \displaystyle\sum_{i=1}^{N} \frac {1}{D_i}(J – 1)^{2i}$$
where, μ_i and α_i for (i = 1, 2, 3) are the constants that needs to be deduced in order to get the appropriate material behaviour. For the compressibility of a material, eq. eq2 can be used as a first approximation for D_i with \(G_o = \frac {1}{2} \displaystyle\sum_{i=1}^{N} \mu_i \alpha_i\)

Arruda-Boyce

The strain energy density of Arruda-Boyce has the following form:
$$U = \mu \displaystyle\sum_{i=1}^{5} \frac {C_i}{\lambda_{m}^{2i-2}} \left(\bar I_{1}^{i} – 3 \right) + \frac {1}{D} \left(\frac {J^2 – 1}{2} – lnJ \right)$$
where, \(C_1 = \frac {1}{2}\), \(C_2 = \frac {1}{20}\), \(C_3 = \frac {11}{1050}\), \(C_4 = \frac {19}{7000}\) and \(C_5 = \frac {519}{673750}\).
For the compressibility of a material, eq. eq2 can be used as a first approximation for D_i with \(G_o = \mu (1 + \frac {3}{5 \lambda_{m}^{2}} + \frac {99}{175 \lambda_{m}^{4}} + \frac {513}{875 \lambda_{m}^{6}} + \frac {42039}{67375 \lambda_{m}^{8}})\)

Getting material parameters from experimental data

In experiments, the response of a material for a specific load is normally deduced by the stress-strain relationship. As an example, a basic stress-strain curve of an elastomer for Uniaxial (Tension), Biaxial and Shear load (provided by Axel Product, Inc.) is shown below:

In order to fit the material properties to the experimental stress-strain curve, one has to take in to account eq. eq1 above. This equation shows the stress as a function of invariants containing principal stretches (1 + “strain values”) and strain energy density (of any of the above model). Using this equation, one can calculate the constants of the undertaken model through numerical scheme to get an accurate fit to the stress-strain material curve. In most of the cases, the accuracy of the fit increases with the increasing order of the model.
As an example, let us consider a stress-strain relation of an incompressible uniaxial case for a general polynomial model:
The principal stretches of uniaxial case are:
$$\lambda_1 = \lambda, \quad \lambda_2 = \lambda_3 = \lambda^\frac{-1}{2} \\
(where, \lambda = 1 + \epsilon)$$
and the deviatoric strain invariants are:
$$\bar I_1 = \lambda^2 + 2 \lambda^{-1}, \quad \bar I_2 = \lambda^{-2} + 2 \lambda$$
From principle of virtual work, it follows that:
δU = σδλ
Than stress can be represented as:
$$\sigma = \frac {\partial U}{\partial \lambda} = \frac {\partial U}{\partial \bar I_1} \frac {\partial \bar I_1}{\partial \lambda}
+ \frac {\partial U}{\partial \bar I_2} \frac {\partial \bar I_2}{\partial \lambda}\\
= 2(\lambda – \lambda^{-2}) \left(\frac {\partial U}{\partial \bar I_1} + \frac {\partial U}{\partial \bar I_2} \right)$$
Therefore, for the general polynomial case:
$$\sigma = 2(\lambda – \lambda^{-2}) \displaystyle\sum_{i+j=1}^{N}[iC_{ij}(\bar I_1 – 3)^i(\bar I_2 – 3)^j + jC_{ij}(\bar I_1 – 3)^i(\bar I_2 – 3)^j]$$
By using the experimental stress-strain data in the above equation, one can calculate the material constants for a specific order of polynomial through a numerical scheme e.g. Least square method.

Giving material parameters in SimScale

Once the material constants of a specific model are obtained for a good fit to experimental data, one can input these constants to the specified material model in SimScale under material properties. In case a polynomial of order 2 is used, the following can be interpreted:

• The provided constants in this case are C₁₀, C₀₁, C₂₀, C₁₁ and C₀₂ as obtained from fitting.
• It is observed that in some cases, if any of the obtained constant value is comparatively higher from the other constants, than the results obtained after the simulation are not accurate.
• For the incompressible case (as mentioned here), one can:
  • Use eq. eq2 for calculating values of D₁ and D₂ with G_o = 2(C₁₀ + C₀₁) and ν = 0.499 as a first approximation.
  • Get more accurate results by giving more lesser value to D₁ and D₂.
  • Avoid convergence issues by not giving much lesser value to D₁ and D₂, specially in case of Code_Aster.
  • Give value to D₁ and D₂ atmost up to 1e^− 10.

Same cases are applied to all other models.