Week - 10 Hyperelastic Material Models
Objective:
To create a Hyperelastic material model with stress-strain curve data and extract Mooney Rivlin and Ogden constant parameters and Evaluation the material constant with actual material data.
Introduction
Hyper-elastic Material
Hyperelastic material models are regularly used to represent the large deformation behavior of materials with FEA. They are commonly used to model mechanical behaviors of unfilled/filled elastomers. In addition to elastomers, hyperelastic material models are also used to approximate the material behavior of biological tissues, polymeric foams, etc.
Linearly elastic materials are described through two material constants (like Young’s modulus and Poisson ratio). In contrast, hyperelastic materials are described through a strain-energy density function. The strain-energy density can be used to derive a nonlinear constitutive model (i.e., stresses as a function of large strain deformation measures like deformation) which are derived by researches and some acceptable models we are using popularly such as Neo-Hookean, Mooney-Rivlin and Signorini models, Ogden, Yeoh, Arruda-Boyce.
Mooney Rivlin
The strain energy density function for an incompressible Mooney–Rivlin material is defined as

Where `Cij`are empirically determined material constants, `I1`and `I2` are the first and the second invariant. By varying the Index values of `i` and `j` we can change the constitutive equations and use it appropriately for our component. Similarly, we will be having multiple coefficients C. This changing of the index will be basically changing the order of the polynomial of the energy density equation
We need to find out the coefficients of these material models for a different order of the polynomial and extract the best fit curve data from the D3HSp file and compare it with the material curve.
Procedure
Workflow:

- Hyper-elastic Material Card
As we are modeling a hyperelastic card, it is mandatory to use a material card that can access the necessary parameters pertaining to that of a hyperelastic material.
The following cards can be used to model a hyperelastic material,
- *MAT_007 *MAT_BLATZ-KO_RUBBER- This is Material Type 7. This one parameter material allows the modeling of nearly incompressible continuum rubber. The Poisson’s ratio is fixed to 0.463.
- *MAT_027 *MAT_MOONEY-RIVLIN_RUBBER-This is Material Type 27. A two-parametric material model for rubber can be defined.
- *MAT_077_h *MAT_HYPERELASTIC_RUBBER This is Material Type 77. This material model provides a general hyperelastic rubber model combined optionally with linear viscoelasticity
- *MAT_077_o *MAT_OGDEN_RUBBER This is also Material Type 77. This material model provides the Ogden rubber model combined optionally with linear viscoelasticity
Apart from these common hyperplastic materials for rubber there other material models too such as MAT_031, 038, 057, 063, 083, 092, 127, 129.
However, MAT_007 is a very simple rubber model with just one parameter to defined i.e. (Shear Modulus) we cant be used for our case and MAT_027 represents Mooney Rivlin constitutive model but this MAT_27 will only solve with First order of polynomial
So we will use MAT_077_h which represents Mooney Rivlin constitutive model and Different order of polynomial from 1 to 3 and find the Material constant for each order we can import a Load curve (experimental data value) and set SGL, SW, and ST parameters as 1 if the load curve is stress vs strain data.
- Defining The Experimental data in the Material model
The stress/strain for the material is now in the form of Engineering stress and Engineering strain which can be directly implemented on the MAT_77h material card without any conversion
These data point for the Stress/Strain curve are inserted *DEFINE_CURVE Card with LCID


Material Card
- MID Material identification. A unique number or label not exceeding 8 Characters must be specified
- RO Mass density= 0.00100
- PR Poisson’s ratio = 0.49100
- N Number of constants to solve for = 1-3 for each case
- LCSS Load curve ID defining Engineering stress versus Engineering strain LCID from *DEFINE_CURVE

- Tensile Test – Pre Processing
Procedure :
- Mesh model of Dogbone specimen is imported
- The material model we created for graphite iron casting is used here
- Section shell card is used with 16: Fully integrated shell element formulation and thickness as 1mm
Boundary Condition
Since we performing a Tensile Test where one end of the specimen is fixed and the displacement is added to the other end so we will use two types of boundary condition on the Specimen
- BOUNDARY_SPC_SET
- To fix one end of the specimen as test procedure we will constraint all the degree of freedom expect Y direction (because specimen will elongate in the X direction and simultaneously contract in Y) for edge node on one end of the specimen

- We will select the mid-node on edge of both ends of the specimen and constrain only the Y direction

- BOUNDARY_PRESCRIBED_MOTION_SET
- A set of a node is created on the edge specimen in which displacement is applied

In the tensile test, the strain value as to be 100%. So we need to have the displacement same as the length of the specimen. The specimen length is 145 mm. So we need to have the displacement also as the 145mm.
- Displacement curve is created using Define_curve Card

- DOF in BOUNDARY_PRESCRIBED_MOTION_SET is set to be 1 which means displacement takes place in the X direction and VAD is set as 2 Displacement flag

Control card
1.*CONTROL_IMPLICIT_GENERAL
Activate implicit analysis and define associated control parameters. This keyword is required for all implicit analyses
- IMFLAG - Implicit/Explicit analysis type flag: 1- implicit analysis
- DT0 – Initial time step size for implicit analysis: 0.0050000
- IMFORM - Element formulation flag for seamless spring back: 2 retain original element formulation (default)
2.*CONTROL_IMPLICIT_SOLUTION
Use of these cards to specify whether a linear or nonlinear solution is desired
- NSOLVER – Solution method for implicit analysis: 12 For Non-linear problem
- ILIMIT-Iteration limit between automatic stiffness reformations: 11
- MAXREF- Stiffness reformation limit per time step: 15
- DCTOL- Displacement relative convergence tolerance: 001
- DNORM- Displacement norm for convergence test: 2 (Increment as a function of total displacement)
3.*CONTROL_IMPLICIT_SOLVER
The linear equation solver performs the CPU-intensive stiffness matrix inversion
- LSOLVR- Linear equation solver method: 4 (uses solve the simulation with single precision which consequently saves run time rather than going for a double-precision solver)
4.*CONTROL_IMPLICIT_AUTO card has been defined to adjusts the time step size
- IAUTOset to 1 which automatically adjusts the timestep size.
- ITEOPToptimizes the operation count step:20
- ITEWINis the iteration window: 5
- DTMIN Minimum allowable time step size: 1.0000e-05
- DTMAXMaximum allowable time step size:0.2500
Database Creation
- *BINARY_D3PLOT - It defines the frequency at which the animation file is to be created and is set to 5.000e-05 ms
- *DATABASE_HISTORY_NODE A nodes on Accel Block of head were defined to measure the acceleration and calculated HIC value
- *DATABASE_EXTEND_BINARY - this writes the strain tensor to d3plot when STRFLG is set to 1
ASCII_option keyword
The output request in ASCII format, The following keyword are Activated
- GLSTAT- Global Data
- MATSUM- Material Energies and
- ELOUT- Element Output Data
- Extracting Material constant and Fit Curve data for each order of polynomial from D3HSp file
The Mooney-Rivlin constants for each case are stored in the d3hsp file
1. Case 1: N = 1 - First order polynomial
- Mooney constant C01 = 0.177
- Mooney constant C10 = 0.147
Final Fit Data from d3hsp file are generally in Stretch vs True Stress but We have to Covert them into Engineering Strain vs Engineering Stress
Engineering Strain = Stretch -1
Engineering Stress = True Stress/Stretch

2. Case 2: N = 2 - Second order polynomial
Mooney constant C10 = 0.259
Mooney constant C01 = 0.4594E-01
Mooney constant C11 = -0.2092E-01
Mooney constant C20 = 0.2913E-02
Mooney constant C02 = 0.1016E-02

3. Case 3: N = 3 - Thrid order polynomial
Mooney constant C10 = 0.2611E+00
Mooney constant C01 = 0.4420E-01
Mooney constant C11 = -0.6076E-02
Mooney constant C20 = -0.2667E-02
Mooney constant C02 = -0.9966E-02
Mooney constant C30 = 0.1565E-03

- By Comparing, the 2nd & 3rd Order of polynomial Final Fit data matches the Experimental data accurately compared to the 1st Order of polynomial.
Now we will compare the Simulated Stress/Strain curve from the Tensile Test for Each order of the polynomial.
- Comparison of Simulation Stress/Strain with Material experimental Stress/Strain data
1. Case 1: N = 1 - First-order polynomial
Effective Stress
Max Stress: 1.657MPa

Mean Ipt X strain
Max Strain: 0.749

Stress and Strain curve is plotted for the element which is there in the mid-portion of the specimen

Since these are True stress Strain values we can to convert it to Engineering Stress and Strain so that we can compare with the Material experimental data

2. Case 2: N = 2 - Second-order polynomial
Effective Stress
Max Stress: 1.691MPa

Mean Ipt X strain
Max Strain: 0.7563

Stress and Strain curve is plotted for the element which is there in the mid-portion of the specimen

Since these are True stress Strain values we can to convert it to Engineering Stress and Strain so that we can compare with the Material experimental data

3. Case 3: N = 3 - Third-order polynomial
Effective Stress
Max Stress: 1.680MPa

Mean Ipt X strain
Max Strain: 0.7553

Stress and Strain curve is plotted for the element which is there in the mid-portion of the specimen

Since these are True stress Strain values we can to convert it to Engineering Stress and Strain so that we can compare with the Material experimental data

Overall Comparison

As we can see the First Order Polynomial N=1 is closer to the Experimental data So we can consider this as the best fit for our case and we can use the Mooney constant from First Order Polynomial N=1 in the Material Deck.
- *MAT_077_o *MAT_OGDEN_RUBBER

Ogden Hyper-elastic model is the expressed the strain energy function W in terms of principal stretches `λ1`, `λ2`, and `λ3` Where N is the order of the model, usually, N takes a number between 1–8. `µ_i` and `α_1` are material constants, where the unit of `µ_i` is pressure, `α_1` is a dimensionless quantity.
The same procedure is followed for building the Ogden Material Deck Since there 8 order of polynomial in this material model we will use 1,3,5 order of polynomial for this case
Density, Passion ratio, SGL, SW, and ST parameters are set with the same value we used for MAT_077_h
- Extracting Material constant and Fit Curve data for each order of polynomial from D3HSp file
The Ogden constants for each case are stored in the d3hsp file
1. Case 1: N = 1 - First order polynomial
`µ_(1)` = 0.910
`α_(1)` = 1.34

2. Case 2: N =3 – Third order polynomial
`µ_(1)`= -3.18E-02
`µ_(2)` = 2.61E-01
`µ_(3)` = 6.82E-01
`α_(1)` = 1.93E+00
`α_(2)` = 1.55E+00
`α_(3)` = 1.29E+00

3. Case 3: N =5 – Fifth order polynomial
`µ_(1)`= -3.18E-02
`µ_(2)` = 2.61E-01
`µ_(3)` = 6.82E-01
`µ_(4)`= -3.18E-02
`µ_(5)` = 2.61E-01
`α_(1)` = 1.93E+00
`α_(2)` = 1.55E+00
`α_(3)` = 1.29E+00
`α_(4)` = 1.93E+00
`α_(5)` = 1.55E+00

By Comparing, the 1st, 3rd, and 5th Order of polynomial Final Fit data they are very similar but they are very deviated from the Experimental data. so we will compare the Simulated Stress/Strain curve from the Tensile Test for Each order of the polynomial and check which the Order of polynomial is better.

- Comparison of Simulation Stress/Strain with Material experimental Stress/Strain data
1. Case 1: N = 1 - First-order polynomial
Effective Stress
Max Stress: 1.694MPa

Mean Ipt X strain
Max Strain: 0.7559

Stress and Strain curve is plotted for the element which is there in the mid-portion of the specimen

Since these are True stress Strain values we can to convert it to Engineering Stress and Strain so that we can compare with the Material experimental data

2. Case 2: N = 3 - Third-order polynomial
Effective Stress
Max Stress: 1.693MPa

Mean Ipt X strain
Max Strain: 0.7560

Stress and Strain curve is plotted for the element which is there in the mid-portion of the specimen

Since these are True stress Strain values we can to convert it to Engineering Stress and Strain so that we can compare with the Material experimental data

2. Case 3: N = 5 - Fifth-order polynomial
Effective Stress
Max Stress: 1.694MPa

Mean Ipt X strain
Max Strain: 0.7561

Stress and Strain curve is plotted for the element which is there in the mid-portion of the specimen

Since these are True stress Strain values we can to convert it to Engineering Stress and Strain so that we can compare with the Material experimental data

As we can see all the order of the polynomial looks same but if we look at the data points values more closer we can say N=3 is better than the other two order of the polynomial
Conclusion :
We created two different hyperelastic material models and understood how each constitutive model works We extracted the material constant for Mooney–Rivlin and Odgen material for each order of the polynomial.
The solver tries to perform curve fitting for each order of polynomial with the actual material data And we compared which order tries to fit the material curve closer and as conclusion 2nd and 3rd order of polynomial has closer value with material data and
We perform the simulation and compare the simulation stress and strain curve with the actual material curve but the Simulation stress and Strain curve data depends on the material properties like Poisson ratio and density and since We have approximately defined these like Poisson ratio and density value so we can’t expect them to match the material curve exactly.
After comparison, we can Say N=2 is better For Mooney–Rivlin Hyperelastic material
and N=3 is better for Odgen Hyperelastic Material Its because they fit the material original Stress and Strain curve closer than other order of polynomial