Hyperelastic materials are the special class of materials which respond elastically even when they are subjected to large deformations. They show both a nonlinear material behavior as well as large shape changes. They are characterized by:
Many polymers show hyperelastic behavior, such as elastomers, rubbers, and other similar soft flexible materials.
Hyperelastic materials are mostly used in applications where high flexibility, in the long run, is required, under the presence of high loads. Some typical examples of their use are as elastomeric pads in bridges, rail pads, car door seal, car tires, and fluid seals.
In finite element analysis, the hyperelasticity theory is used to represent the non-linear response of hyperelastic materials at large deformations. Hyperelasticity is popular due to its ease of use in finite element models. Usually, stress-strain curve data from experimental tests is used to fit the constants of theoretical models, thus approximating the material response.
The choices of hyperelasticity models which are available in SimScale platform are:
The stress-strain relation for hyperelastic materials is normally calculated with a strain energy density function. A brief theoretical description follows.
Consider a solid body subjected to a great deformation. A point inside the body with position \(X\) with respect to the original configuration is displaced to the position \(x\) in the final configuration through the displacement vector \(u\):
$$ x = X + u $$
The tensor of the gradient of the deformation \(F\) relating the initial configuration to the final configuration is given by:
$$ F = \frac{ \partial x }{ \partial X } = I + \nabla_{X} u $$
The local change of volume \(J\) is given by its determinant \(J = det(F)\). This tensor turns out to not be the best way to describe large deformations, so the right Cauchy-Green tensor is used instead:
$$ C = F^T F $$
This tensor is symmetrical, with its invariants given by:
$$ \mathrm{I}_C = tr(C) $$
$$ \mathrm{II}_C = \frac{1}{2} \Big( tr(C)^2 – tr(C^2) \Big) $$
$$ \mathrm{III}_C = det(C) $$
The third invariant \(\mathrm{III}_C\) relates to the change in volume, and can also be written as:
$$ \mathrm{III}_C = det(F)^2 = J^2 $$
In the case of an incompressible material, \(J = 1\).
Of high interest for our purpose is to express the right Cauchy-Green tensor and its invariants in terms of the principal stretches \( \lambda_i \):
$$ \lambda_i = 1 + \varepsilon_i $$
$$ C_{ii} = \lambda_i^2 $$
$$ C_{ij} = 0, i \neq j $$
$$ \mathrm{I}_C = \lambda_1^2 + \lambda_2^2 + \lambda_3^2 $$
$$ \mathrm{II}_C = \lambda_1^2 \lambda_2^2 + \lambda_2^2 \lambda_3^2 + \lambda_3^2 \lambda_1^2 $$
$$ \mathrm{III}_C = \lambda_1^2 \lambda_2^2 \lambda_3^2 $$
where \( \varepsilon_i \) are the nominal strains in the principal directions (measured in the original configuration). Notice that the stretches \( \lambda_i \) are the components of the gradient of deformation tensor \(F\).
It is assumed that a strain energy density function \(U\) exists such that the hyperelastic (second) Piola-Kirchhoff stress tensor \(S\) in the material can be related to the right Cauchy-Green deformation tensor:
$$ S = 2 \frac{ \partial U }{ \partial C } $$
The true stress tensor on the material is related to the second Piola-Kirchhoff stress tensor. It can be expressed (after some algebra, not shown here) in terms of the right Cauchy-Green tensor and the strain energy density function:
$$ \sigma_{ij} = -p\delta_{ij} + 2 \frac{ \partial U }{ \partial I_1 } C_{ij} – 2 \frac{ \partial U }{ \partial I_2 } C_{ij}^{-1} $$
Here, \(p\) is the hydrostatic external pressure, which causes pure volumetric change, and \(delta_{ij}\) is the kronecker delta. Notice that due to the incompressibility assumption, the \(J\) term is dropped.
For the material models available in SimScale, the strain energy density function is given by:
$$ U = C_{10}(I_1 – 3) + C_{01}(I_2 – 3) + C_{20}(I_1 -3)^2 + \frac{1}{2} K (J – 1)^2 $$
$$ I_1 = \mathrm{I}_c J^{-\frac{2}{3}} $$
$$ I_2 = \mathrm{II}_c J^{-\frac{4}{3}} $$
$$ J = \mathrm{III}_c^{\frac{1}{2}} $$
You will find that in the user interface, the compressibility is controlled by the \(D_1\) constant, given by:
$$ D_1 = \frac{2}{K} $$
Before giving the appropriate material parameters to define specific hyperelastic materials, one should know the strain energy density forms of the hyperelasticity models. Following are the strain energy density forms of all the available hyperelasticity models on the SimScale platform (as mentioned above).
With \( C_{01} = C_{20} = 0 \) in the above formulation, one obtains the most basic form known as neo-Hookean, given by:
$$ U = C_{10}(I_1 – 3) + \frac{1}{D_1} (J – 1)^2 $$
With \( C_{20} = 0 \) in the above formulation, one obtains the enhancement to neo-Hookean form known as Mooney-Rivlin, given by:
$$ U = C_{10}(I_1 – 3) + C_{01}(I_2 – 3) + \frac{1}{D_1} (J – 1)^2 $$
Strain energy density function for Signorini is represented as:
$$ U = C_{10}(I_1 – 3) + C_{01}(I_2 – 3) + C_{20}(I_1 -3)^2 +\frac{1}{D_1} (J – 1)^2 $$
The final piece of the puzzle to be able to fully characterize the hyperelastic material is the definition of the stretches \( \lambda_i \), which are determined by the loading state. The same material will have different behavior depending on how it is loaded, as shown in the following plot (taken from Axel Product, Inc.). In a typical test, the specimen is subjected to the desired load state, and the stress-strain relation is measured.
For the case of uniaxial tension loading of an incompressible material, the specimen is loaded along one axis. The stretches are given by:
$$ \lambda_1 = 1 + \varepsilon $$
$$ \lambda_2 = \lambda_3 = \lambda_1^{ -1/2 } $$
For the case of equibiaxial tension loading of an incompressible material, the specimen is loaded along two perpendicular axes with the same magnitude. The stretches are given by:
$$ \lambda_1 = \lambda_2 = 1 + \varepsilon $$
$$ \lambda_3 = \lambda_1^{ -2 } $$
For the case of pure shear loading of an incompressible material, the specimen is loaded in two perpendicular directions, one in tension and the other in compression. Plane strain condition is assumed, so there is no deformation in the unloaded direction. The stretches are given by:
$$ \lambda_1 = 1 + \varepsilon $$
$$ \lambda_2 = 1 $$
$$ \lambda_3 = \lambda_1^{ -1 } $$
In order to fit the material properties to the experimental stress-strain curve, one has to take into account all the relations presented above. Combining the equations, the stress is expressed in terms of the strain, and the constants of the model are determined through a least-squares fit of that function to the test data. 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 Signorini model:
$$ U = C_{10}(I_1 – 3) + C_{01}(I_2 – 3) + C_{20}(I_1 -3)^2 +\frac{1}{D_1} (J – 1)^2 $$
Using the given definition for uniaxial stretches, the principal true stress can be written as:
$$ \sigma_{11} = -p + 2 \frac{\partial U}{\partial I_1} \lambda^2 – 2 \frac{\partial U}{\partial I_2} \lambda^{-2} $$
And:
$$ \sigma_{22} = -p + 2 \frac{\partial U}{\partial I_1} \lambda^{-1} – 2 \frac{\partial U}{\partial I_2} \lambda = 0 $$
After subtracting to get rid of the pressure term, factoring out, and performing the derivatives we get:
$$ \sigma_{11} = 2 ( \lambda^2 – \lambda^{-1} ) \Big[ \frac{\partial U}{\partial I_1} + \lambda^{-1} \frac{\partial U}{\partial I_2} \Big] $$
$$ \frac{\partial U}{\partial I_1} = C_{10} + 2 C_{20} (I_1 – 3) $$
$$ \frac{\partial U}{\partial I_2} = C_{01} $$
This completes the stress-strain function that should be input to the least-square fit algorithm to determine the \( C_{10}, C_{01}, C_{20} \) constants.
Important
Notice that when using the equations as defined in this procedure for the relations among stretch, strain and stress, the input quantities for the parameter curve fitting should be nominal (engineering) strain and true stress.
Similar expressions can also be derived for the equibiaxial tension case:
$$ \sigma_{11} = 2 ( \lambda^2 – \lambda^{-4} ) \Big[ \frac{\partial U}{\partial I_1} + \lambda^{2} \frac{\partial U}{\partial I_2} \Big] $$
And for pure shear case:
$$ \sigma_{11} = 2 ( \lambda^2 – \lambda^{-2} ) \Big[ \frac{\partial U}{\partial I_1} + \frac{\partial U}{\partial I_2} \Big] $$
Notice that due to the incompressibility condition, the \(D_1\) parameter can not be estimated with this procedure. The recommendation is to give an approximate value using the elasticity relations:
$$ D_1 = \frac{ 1 – 2 \nu }{ C_{01} + C_{10} } $$
For an almost incompressible material, the value of \( \nu \) is close to 0.5.
Once the material constants of a specific model are obtained from a good fit to experimental data, one can input these constants to the specified material model in SimScale under material properties:
Last updated: March 4th, 2021
We appreciate and value your feedback.
