Selection of a Hyper-Elastic Material Model – a Case Study for a Polyurethane Component

The aim of this study was to select an appropriate material model for a polyurethane rail pad. The considered material was assumed to be hyper-elastic, and model selection was limited only to the Mooney–Rivlin (M-R) material models available in the MSC.Marc/Mentat system. The selection of the most appropriate material model was based on the most accurate approximation of the experimental results, considering also CPU time, the stability of the material model, and convergence of the FE analyses. The current study was limited to static analysis without cyclical deformation modes. The conclusions reached and the observations from the conducted analyses should not be generalized to the whole group of hyper-elastic materials, as they relate only to the specific material under consideration. However, the results may provide a foundation for researchers involved in the analysis of such types of engineering materials.


INTRODUCTION
Finite element (FE) analysis is widely used to predict deformations and stress distributions, including for hyperelastic materials.There are many material models that describe these materials, for which Hooke's law does not apply.Moreover, proper analysis of rubber-like components requires nonlinear FE analysis tools, which differ greatly from those used for metallic components.The unique properties of rubber-like materials include the following (MSC 2010): • an ability to undergo large deformations under load, sustaining strains of up to 500% in engineering applications; • markedly non-linear load-extension behaviour; and • near incompressibility, with no appreciable change in volume in response to stress.
The successful modelling and design of rubber components depends on the selection of an appropriate strain energy function and accurate determination of the material constants in this function.However, the testing of hyperelastic materials for the purpose of defining material models is often very difficult, with experimental requirements that differ from those for most standardized test methods.Moreover, national and international standards organizations have yet to clearly define appropriate experiments (MSC 2010).
The aim of this study was to select an appropriate material model for a polyurethane rail pad.The main role of the elastic rail pads applied in fastening systems is to reduce the dynamic effect of a rolling railway vehicle on concrete sleepers.The pads also provide electrical insulation between the rail foot and the sleeper.The stiffness of a rail pad depends on the track class and the maximum track capacity (Szurgott et al. 2012).The rail pad considered in this study was designed using computer-aided engineering (CAE) systems as a result of the modernization of existing solutions and was optimized mainly in terms of ensuring the appropriate stiffness.The basic pads were made of polyurethane without any specified material parameters, and unspecified material was applied for the modernized pad as well.The conducted studies were limited to static analysis, and therefore cyclical deformation modes were omitted.
When selecting a constitutive model, the degree of complexity is an important factor, and compromise is often required.More complex models will better approximate the experimental data but may result in instability of the analysis and increased central processing unit (CPU) time.Among constitutive models, hyper-elastic Mooney-Rivlin (M-R) based material models have been applied in numerical simulations of various problems.
For example, Martins et al. (2006) considered seven material models, including some described in the current paper, i.e.Neo-Hookean, M-R and Yeoh.Material constants were obtained via the inverse method, and the fitting accuracy of the experimental stress-strain relation was evaluated.The correlation coefficient exceeded 0.995 for the M-R material model and was close to 1 for the Yeoh model.Uniaxial stretching and equibiaxial tests were proposed by Sasso et al. (2008) to fit the model and extract material parameters.Numerical simulations of the test procedures gave a useful comparison between the numerical and experimental data.Four material models from the M-R group of models were taken into consideration, and the second order M-R model was found to be the most accurate.Berselli et al. (2011) presented accurate, hyper-elastic constitutive models describing photopolymer behaviour.The authors fit the numerical data to the experimental data by applying the Neo-Hookean and 2-term M-R models.Uniaxial compression and tension tests were taken into account.The Neo-Hookean model could not adequately capture the material behaviour within the stretch range under consideration, whereas the M-R model provided acceptable fitting accuracy within 10% relative error.
With respect to rubbery materials, Hamzah and Razaq (2013) took into account the complex behaviour of these materials in theoretical and experimental works in order to attain an appropriate constitutive model of rubber components for a specific application.Uniaxial tensile tests were carried out to obtain stress-stretch curves.The 2term M-R model adequately described the hyper-elasticity of the material under consideration, and the numerical results exhibited good agreement with the experimental data.Jahani and Mahmoodzade (2014) identified dynamic material constants of the 2-term M-R constitutive model for elastomeric components.The material constants were obtained by curve fitting of a uniaxial stress-strain curve obtained from an experimental tension test.The selection of hyper-elastic material models for a natural rubber reinforced with carbon black was proposed by Shahzad et al. (2015).Four experimental tests were performed: uniaxial, biaxial, planar shear and volumetric.The results showed that the Yeoh model was the most suitable choice, taking into account all four tests.Kumar and Rao (2016) performed experimental curve fitting using four M-R material models with 2, 3, 5 or 9 C ij.The 5-term model yielded the most accurate approximation while maintaining the correct behaviour of the.The 9-parameter M-R model was used to describe high-temperature silicone rubber in (Lobdell and Croop, 2016).Stress-strain curves were obtained from the uniaxial tension/compression and planar tension tests, and a biaxial curve was selected on the basis of the uniaxial compression data.Abdelsalam et al. (2017) selected material models to appropriately describe styrene butadiene rubber based on a uniaxial tensile test only and various elastomeric material models available in MSC.Marc software.The relative error was calculated to determine the accuracy of data fitting.The third order deformation model proved to be the best because it generated the smallest relative percentage error in experimental data fitting.Finally, Łagan and Liber-Kneć (2017) applied an M-R material model for swine skin tissue modelling and obtained good matching to the experimental data from uniaxial tensile tests for the 2-parameter M-R and Yeoh models.
To obtain reliable results, the constitutive model must be selected carefully, as demonstrated by Tobajas et al. (2018).The authors considered six material models widely cited in the literature, including the Neo-Hookean, 2-term M-R and Yeoh models.To determine the quality of each of the models, the R 2 correlation coefficient between the experimental and calculated values was employed.The 2-parameter M-R model was most capable of reproducing the real behaviour of the material under consideration.Frequently, when obtaining the material constants for an M-R model, only a uniaxial data set is used.Keerthiwansa et al. (2018) described the risk of using such an approach for fitting the material model, considering the 2-and 3-parameter M-R models and the Yeoh model.The authors discussed the variations related to three basic deformation modes (uniaxial, equibiaxial, pure-shear) and stated that using only the uniaxial set of data may be insufficient.
Based on this analysis of the literature, it is clear that relatively simple M-R models-those with two or three Cij constants-dominate numerical analyses of hyper-elastic materials.Applying more complex models with a large number of Cij constants does not guarantee better correlation and higher accuracy.The selection of material constants is usually based on the correct approximation of experimental results, and other evaluation criteria are rarely considered.Moreover, uniaxial tests are often the only available experimental data.To address these shortcomings, the current paper presents a broader view of the selection process for hyper-elastic material constants.

MODELLING OF HYPER-ELASTIC MATERIALS
Numerical simulations were carried out using MSC.Marc/Mentat computer code.Therefore, the terminology used here corresponds to the MSC software nomenclature.Selection of the material model was limited only to the M-R material models available in the mentioned software.Other models such as Ogden (principal stretch model) or Arruda-Boyce and Gent (micro-mechanical models) were not taken into account at this stage.
The Mooney-Rivlin constitutive models are based on the strain energy function W. This function, disregarding terms associated with volumetric deformation and temperature changes, can be written in the following polynomial form (MSC 2010): where I1 and I2 are deformation invariants and Cij are constants defined for a particular material.The deformation invariants are given as follows: where λ1, λ2, λ3 are the elongations of the considered element along the respective directions.M-R constitutive models differ from each other in the degree of complexity of the strain energy function resulting from the number of Cij constants.Older versions of the MSC.Marc/Mentat software gave the possibility of using one of seven different models due to the combination of five Cij constants.In the most recent versions, the family of polynomial strain energy functions can be generalized to a complete 5th order as a combination of constants Cij, where i = 1, 2, 3, 4, 5 and j = 1, 2, 3, 4, 5.However, the current paper will show that a higher-order fit does not necessarily mean a more correct model.
The Cij constants expressed in Eq. ( 1) can be determined on the basis of the stress-strain curves obtained from experimental tension/compression tests.In the current study, the Cij constants were determined using the Experimental Data Fitting algorithm included in the MSC.Mentat preprocessor, which is based on the method of least squares (MSC 2010).Since uniaxial tension/ compression test results were available, it was necessary to verify the model behaviour in other states, e.g.biaxial, planar shear etc. (Table 1).Moreover, the predicted response of the model beyond the range of strains experienced during the experimental test was checked.
Table 1 Basic deformation modes and relationships between respective elongations (MSC 2010).
Assuming incompressibility of the material: and particular elongations λ1, λ2, λ3 as a function of elongation λ = 1 + ε, the deformation invariants can be expressed as follows: • for the uniaxial deformation mode (5) • for the planar shear deformation mode Having defined the strain energy function W, one can determine the stress components along the direction of interest as a derivative σ(ε) = ∂W/∂ε, where σ is the engineering stress and ε is the engineering strain.For each of the M-R models, the strain energy function W was converted to a form including the Cij constants and elongation λ.The σ(ε) relation was not derived explicitly due to the complexity of the mathematical transformations.The σ(λ) relations were limited to the uniaxial deformation mode only (MSC 2010): 1. Neo-Hookean model 5. Second order invariant model 6. Third order deformation model 7. Yeoh model The number of components (terms) in the σ(λ) relation increases with the number of constants Cij.The polynomial form of the stress-strain relations (7)-( 13) would be even more complex after expansion of the (1 + ε) n expression.

EXPERIMENTAL TESTS
Experimental tests were carried out to determine the engineering stress-strain characteristic in the uniaxial tension/compression state for the polyurethane material under consideration.These tests included tensile tests of 5 dumbbell specimens and compression tests of 10 cylindrical specimens.Both tests were carried out using an INSTRON 8802 universal testing machine as described in detail previously (Szurgott et al. 2012), (MUT 2010).The tension tests were conducted according to the ASTM (2002) standard.The tension force and vertical displacement of the machine crosshead were registered during the test.Moreover, the longitudinal strain was measured using an extensometer attached to the specimen.Each test was conducted once the longitudinal strain reached 200%.The compression test was carried out according to the Polish standard (PCS 1954).The diameter-to-height ratio for the specimen corresponded to the values recommended in (ASTM 2001).
The results of the tension and compression tests were averaged and combined into the final curve depicted in Figure 1.This curve was used to determine the material constants in the constitutive models.The number of σengεeng characteristic points was limited to 53 with intervals of 0.05 for the strain.This set of points was intended for approximation using the EXPERIMENTAL DATA FIT module in the MSC.Mentat pre-processor.The set of 53 points used for data fitting is plotted.

Experimental data fitting
The Cij constants were determined using the EXPERIMENTAL DATA FIT algorithm included in the MSC.Mentat preprocessor, which is based on least squares fitting (MSC 2010).The least squares error is given by either where ΔR and ΔA are the relative and absolute least squares errors, respectively (the relative error is the default in MSC.Mentat), n is the total number of data points, σc is the calculated stress, and σm is the measured engineering stress (taken from the experimental data).Figures 2a to 8a (Jarzębski 2015) show the σengεeng curves obtained as a result of data fitting for each material model.The results relate to the uniaxial mode only since these data were available from the experimental test.The graphs are presented for a deliberately increased range of strain to permit an assessment of the behaviour of the model beyond the range of strain received during the experimental test.
The material Cij constants in the constitutive models under consideration determined using the EXPERIMENTAL DATA FIT algorithm are provided in Tables 2 and 3 for the relative and absolute error criteria, respectively.As the material Cij constants and, in turn, the complexity of the constitutive model increase, the relative and absolute errors decrease.The approximation of the experimental data is more accurate, but the behaviour of the model is not appropriate in some cases.Several material models are unstable since they do not meet Drucker's stability postulate (Drucker 1957): where dσ is the stress increment tensor associated with the strain increment tensor dε through the constitutive relation.Relation ( 16) states that the incremental internal energy of a material can only increase.Stability was assessed based on the waveforms from Figures 2a-8a.The approximations of the experimental values by the different material models according to the model stability criterion are summarized in Table 4.Only two material models-the Neo-Hookean and 2-term Mooney-are stable in the uniaxial tension/compression state regardless of the accepted error criterion.The next two models-the Yeoh and 3rd order deformation-are stable during compression and tension only when one error criterion is used.The remaining models are generally unstable in or beyond the range of strain received during the experimental test.

Approximation error analysis
Error analysis is performed for each model to more clearly evaluate the match between the calculated and experimental values, as depicted in Figures 2b-8b.The normalized error (%) is defined as follows: where σm(ε) and σc(ε) are the measured and calculated engineering stress, respectively.The normalized error for the strain ε = 0 is omitted in the graphs.Each graph includes curves for the relative and absolute error criteria.The values correspond to the engineering stress-engineering strain (σeng -εeng)curves from Figures 2a-8a.
The values of the normalized errors are averaged according to the following relation: for compression and tension independently (Figure 9a and 9b, respectively).A better approximation is reached for tension.The average of the normalized errors for the whole range of strain (Figure 9c) is proposed to qualitatively compare the approximations obtained using the different material models and assess the best compatibility between the experimental and calculated results.
Latin American Journal of Solids and Structures, 2019, 16(5), e191 9/16 Based on the comparison depicted in Figure 9, a better approximation is obtained using the default relative error criterion.The weighted average does not exceed 3.7% except for the Neo-Hookean model (9%).Once the absolute error criterion is applied, the average of the normalized errors increases and varies between -12% and 11.4%.Therefore, the 3-term, Signiorini, 2nd and 3rd order material models approximate the experimental values most accurately.However, some models have been eliminated due to instability.

Correlation between measured and calculated values
Another determinant used in this matching process was the Pearson correlation coefficient between the measured and calculated values of engineering stress.It can be defined as follows (Cohen 1998): where σm, σc are the measured and calculated engineering stress, respectively, and m σ , c σ are the means of σm, σc.As mentioned before, the number of σ(ε) characteristic points was limited to 53.
It can be assumed that the correlation between the two considered variables is large if the correlation coefficient ranges between 0.5 and 1.0 (Cohen 1998).In the current study, the correlation is very large since the obtained values of the correlation coefficient (Table 5) are close to unity.However, the correlation for the Neo-Hookean model clearly stands out from the rest of the models as the lowest.The highest correlation is observed for the 3rd order deformation model when the absolute error criterion is taken into account.

Verification of material models in states other than uniaxial tension/compression
Only the uniaxial tension/compression test results were taken into consideration above, and the model behaviour in other states should be assessed.The MSC.Mentat pre-processor permits the verification of the material model behaviour under biaxial and planar shear (Figure 1) for the assumed Cij constants.The waveforms shown in Figures 10-16 (Jarzębski 2015) should be considered only in qualitative terms because neither biaxial nor planar shear tests were conducted.Nevertheless, the waveforms permit an assessment of the model behaviour in such states.The non-Latin American Journal of Solids and Structures, 2019, 16(5), e191 10/16 physical behaviour of the higher-order models with a larger number of Cij constants confirms the previous analyses.
Table 6 provides information on the computational stability of the material model-understood as a lack of non-physical behaviour-in all three deformation modes.
The results presented in Table 6 clearly show that the simplest models-the Neo-Hookean and the 2-term M-Rmost accurately describe the polyurethane under consideration.Furthermore, there was no non-physical behaviour of the 3rd order deformation model with Cij constants determined on the basis of the smallest absolute error criterion.The above-mentioned material models can thus be considered stable in all three deformation modes during both tension and compression.

RESULTS OF FINITE ELEMENT ANALYSES
The numerical analyses in the current study were limited only to the static compression test of the pad.An FE model of a quarter of the pad was developed due to the symmetry of the problem.The pad was modelled as a deformable body, whereas the top surface of the sleeper and the bottom face of the rail were treated as rigid bodies (Figure 17).The results of the maximum equivalent stress in the last three calculation steps are summarized in Figure 25.In general, the stress values increase in subsequent steps except for the Mooney(3) and 2nd order invariant material models with the relative error criterion.
The results obtained for the relative error criterion (Figure 25a) vary more greatly than those for the absolute error criterion (Figure 25b); the latter yields similar values regardless of the applied material model.Based on the results presented in Figure 25, the absolute error criterion seems to be a better criterion for selecting Cij constants in material models due to the full convergence in each of the seven cases and the similar values of equivalent stress obtained.
The CPU time was the last criterion for material model selection.All FE analyses were carried out on the same average-class computer.The shortest calculation time, 397 s, was obtained with the 2nd order invariant material model, whereas the longest analysis-excluding non-convergent cases-was obtained with the Neo-Hookean material model, with a time of 626 s.In both cases, the absolute error criterion for selecting Cij constants was assumed.

CONCLUSIONS
The selection of a material model on the basis of the most accurate approximation of the experimental results is an important and very difficult problem.Ultimately, more complex material models are not necessarily preferable.Increasing the number of criteria for model selection allows researchers to make the right choice.
The current study takes into account the following selection criteria: • the most accurate approximation of the experimental results in the uniaxial tension/compression state, including approximation error analysis and the correlation between the calculated and experimental values; • the complexity of the material model, including the number of Cij constants and the complexity of the stressstrain relation; • verification of the material models in states other than uniaxial tension/compression; • the results obtained from preliminary FE analyses of a simple component made of the hyper-elastic material, including the convergence of these analyses and CPU time.
Taking into account all of these factors allows appropriate selection of the material model.Based on the proposed comprehensive selection of a material model, the following conclusions can be drawn.Material models with a larger number of Cij constants approximate the stress-strain relation more accurately.However, they are less computationally stable and may lead to non-convergence in analysis.Application of the absolute error criterion during the determination of the Cij constants seems to be more reasonable because similar values are obtained regardless of the applied material model.A large number of Cij constants does not increase the CPU time, which remains at a similar level for most material models.
Observations during analyses should not be generalized to the whole group of hyper-elastic materials since they relate only to the specific material under consideration.However such observations may provide a foundation for researchers involved in the analysis of these types of engineering materials.In conclusion, the material model describing the polyurethane under consideration in the most appropriate manner-satisfying the above conditions and guaranteeing the stability of the model and convergence of the analysis-is the Mooney(2) model with only two constants, C10 and C01.

Figure 1
Figure 1 Averaged engineering stress-engineering strain curve for the material under consideration.The set of 53 points used for data fitting is plotted.

Figure 2
Figure 2 Stress-strain curves (a) obtained as a result of data fitting and normalized error (b) between the experimental and calculated results; Neo-Hookean material model.

Figure 3 16 Figure 4
Figure 3 Stress-strain curves (a) obtained as a result of data fitting and normalized error (b) between the experimental and calculated results; Mooney(2) material model.

Figure 5
Figure 5 Stress-strain curves (a) obtained as a result of data fitting and normalized error (b)between the experimental and calculated results; Signiorini material model.

Figure 6
Figure 6 Stress-strain curves (a) obtained as a result of data fitting and normalized error (b)between the experimental and calculated results; 2nd order invariant material model.

Figure 7 16 Figure 8
Figure 7 Stress-strain curves (a) obtained as a result of data fitting and normalized error (b) between the experimental and calculated results; 3rd order deformation material model.

Figure 9
Figure 9 Averaged normalized error between the experimental and calculated results determined for compression (a) and tension (b) independently, and for tension/compression combined (c).

Figure 10
Figure 10 Engineering stress-engineering strain curves expected for the Neo-Hookean material model in biaxial (a) and planar shear (b) deformation modes.

Figure 11
Figure 11 Engineering stress-engineering strain curves expected for the Mooney(2) material model in biaxial (a) and planar shear (b) deformation modes.

Figure 12 16 Figure 13
Figure 12 Engineering stress-engineering strain curves expected for the Mooney(3) material model in biaxial (a) and planar shear (b) deformation modes.

Figure 14
Figure 14 Engineering stress-engineering strain curves expected for the 2nd order invariant material model in biaxial (a) and planar shear (b) deformation modes.

Figure 15
Figure 15 Engineering stress-engineering strain curves expected for the 3rd order deformation material model in biaxial (a) and planar shear (b) deformation modes.

Figure 16
Figure 16 Engineering stress-engineering strain curves expected for the Yeoh material model in biaxial (a) and planar shear (b) deformation modes.

16 Figure 18
Figure 18 Contours of equivalent stress (MPa) in the compressed rail pad.Neo-Hookean material model with relative (a) and absolute (b) error criteria.

Figure 19
Figure 19 Contours of equivalent stress (MPa) in the compressed rail pad.Mooney(2) material model with relative (a) and absolute (b) error criteria.

Figure 20
Figure 20 Contours of equivalent stress (MPa) in the compressed rail pad.Mooney(3) material model with relative (a) and absolute (b) error criteria.

Figure 21 16 Figure 22
Figure 21 Contours of equivalent stress (MPa) in the compressed rail pad.Signiorini material model with relative (a) and absolute (b) error criteria.

Figure 23
Figure 23 Contours of equivalent stress (MPa) in the compressed rail pad.3rd order deformation material model with relative (a) and absolute (b) error criteria.

Figure 24
Figure 24 Contours of equivalent stress (MPa) in the compressed rail pad.Yeoh material model with relative (a) and absolute (b) error criteria.

Figure 25
Figure 25 CPU time and maximum values of equivalent stress in the last three calculation steps registered for all material models with the relative (a) and absolute (b) error criteria.

Table 2
Cij constants (MPa) determined using the EXPERIMENTAL DATA FIT algorithm for different material models; relative error criterion according to Eq. (14).

Table 3
Cij constants (MPa) determined using the EXPERIMENTAL DATA FIT algorithm for different material models; absolute error criterion according to Eq. (15).

Table 4
Approximation of the experimental values by different material models according to the model stability criterion.

Table 5
Correlation coefficient between measured and calculated values of engineering stress.