Prediction Accuracy of Hyperelastic Material Models for Rubber Bumper under Compressive Load

Different hyperelastic material models (Mooney-Rivlin, Yeoh, Gent, Arruda-Boyce and Ogden) are able to estimate Treloar’s test data series containing uniaxial and biaxial tension and pure shear stress-strain characteristics of rubber. If the rubber behaviour is only determined for the specific load of the product, which, in the case of rubber bumpers, is the compression, the time needed for the laboratory test can be significantly decreased. The stress-strain characteristics of the uniaxial compression test of rubber samples were used to fit hyperelastic material models. Laboratory and numerical tests of a rubber bumper with a given compound and complex geometry were used to determine the accuracy of the material models. Designing rubber products requires special consideration of the numerical discretization process due to the nonlinear behaviours (material nonlinearity, large deformation, connections, etc.). Modelling considerations were presented for the finite element analysis of the rubber bumper. The results showed that if only uniaxial compression test data are available for the curve fitting of the material model, the Yeoh model performs the best in predicting the rubber product material response under compressive load and complex strain state.


Introduction
Rubber-based machine elements, like vibration dampers, bumpers, and seals, are becoming increasingly common due to synthetic rubber compounds' development and shock-absorbing behaviour.In the automotive sector, they perform a crucial function in stabilizing the engine, transmission or cabin, as a chassis component, or as a jounce bumper that limits the last stage of movement.Even with the advancements in electromobility, many of these components are still required in large quantities during the manufacture of automobiles.However, because of the modified structural designs, products must be redesigned considering the loads resulting from the altered total vehicle weight.The only way to satisfy the resulting market demands in time involves innovative design processes that use optimization and the finite element method to determine the rubber product's design [1][2][3][4][5][6][7][8][9][10].Since a large portion of the literature [11][12][13][14][15] works on the numerical analysis of rubber bumpers, the research aimed to determine which considerations have to be taken to accurately model the behaviour of these rubber products under compressive load.
Rubber belongs to the group of elastomers, i.e., a polymer with elastic and viscous properties.When modelling rubber products, it is necessary to deal with nonlinearities due to material, geometric and contact effects.Numerical solvers for nonlinear equation systems can be found in commercial finite element software packages.Selecting a calibrated material model that describes the stress-strain relation is necessary for a structural analysis.The material response can be derived from either the phenomenological modelling of measurable macroscopic quantities like the relationship between stress and strain or the micromechanical modelling of local material behaviours and interactions.A hyperelastic constitutive model that describes the relationship between stress and strain using relevant continuum mechanics theories and the strain energy density function can predict the nonlinearly elastic and isotropic behaviour of rubber materials [16][17][18].In phenomenological hyperelastic material models, the material response is approximated by mathematical functions.Mooney-Rivlin [19,20], Ogden [21], Yeoh [22], Arruda-Boyce [23] and Gent [24] are a few examples of these model formulas; the latter two were derived using micromechanical model approximations.The load mode, the raw material, or its availability in the finite element software are the main factors that influence the choice of material model [25].Since the rubber compound used in the bumper is a trade secret, raw material measurements for various load cases are always required.Using these data, a process called curve fitting can be used to determine the values of the parameters included in the constitutive equations [26,27].
The combination of the uniaxial, biaxial or simple shear test will enhance the characterization accuracy of the hyperelastic behaviour of the rubber.Since only the uniaxial test is standardized, most material testing laboratories only have the equipment needed to perform this test.Treloar's natural rubber measurement data set [28] is used extensively in the literature [29][30][31][32][33] to evaluate the accuracy of material models.Examining twenty distinct material models was performed in one of the most often cited works.While the Mooney-Rivlin model performed well for a restricted strain of 250%, the Ogden model was among the best at describing the behaviour of the material over the whole deformation range [31].Forty-four different hyperelastic material models were fitted using a multi-objective optimization process based on a genetic algorithm.The number of material parameters was considered when selecting the best-performing model, in addition to the goodness of fit [32].Based on an evaluation system that considers the standard deviation of the approximation and the number of material parameters, a ranking of fifteen hyperelastic models was presented.The results show that the Ogden and Yeoh model is one of the best-performing models [33].In [34], 85 different isotropic strain energy function-based material models were fitted for Treolar's unfilled and Yeoh's filled rubber data set.The fitting results showed that the two data sets' best-performing material models are different.Numerical stress solution method was used to calibrate hyper-visco elastic solids subjected to various loading modes in [35][36][37], and prediction capability was tested by simulation and measurement of the radial force in an aged seal [38].Few studies have tested the accuracy of the fitted material model using rubber product finite element simulation [39][40][41].After a thorough analysis of 25 hyperelastic constitutive models, Steinmann and Hossain [29,30] concluded that when just one type of test data is available, the prediction of the natural rubber behaviour under the remaining two load modes deviates irregularly using the different constitutive models.In the lack of material test results, a satisfactory alternative could be the inverse material calibration.Nevertheless, it is a difficult and computationally demanding task, particularly when dealing with numerous material parameters [42][43][44].
Based on the reviewed literature, there is enough knowledge available about the behaviour and relative accuracy of hyperelastic material models fitted to the Treloar measurement data set.These showed that the Mooney-Rivlin, Ogden, Yeoh, Arruda-Boyce, and Gent were the most accurate phenomenological models.However, one of the engineering applicability problems of these results is that, in the case of rubber products under compression, the strain range is smaller than the tested one, and the Treloar data set does not contain uniaxial compression measurement.The rubber bumper material is a blend of natural rubber (NR) and styrene-butadiene rubber (SBR), which differ from the natural rubber used in the Treolar dataset and may have an impact on the prediction accuracy of the different hyperelastic material models.To numerically calculate the operating characteristics of a rubber product under compressive load, it is vital to select an appropriate hyperelastic material model and to be aware of its accuracy.It is often possible to identify the primary load and its extent based on the design requirements.In addition, a measure for the element strain state called biaxiality indication will be used to determine the load mode required for the curve fitting of the hyperelastic model.In the case of the investigated rubber bumper, this is the uniaxial compressive load up to a 30% drop in height.None of the reviewed studies examined the accuracy of predicting the behaviour of a rubber product under compressive load if only the uniaxial compression data set is used for the material model calibration.Consequently, this research aim was to fit various hyperelastic material models (Mooney-Rivlin, Ogden, Yeoh, Arruda-Boyce, Gent) for the stress-strain characteristics derived from the uniaxial compression test.Furthermore, the study sought to evaluate the accuracy of the fitted material model using the rubber bumper finite element simulation for compressive load.

Hyperelasticity
If we ignore the viscoelastic behaviour of the rubber, then its stress state depends on the current state and any stress measure depends on the deformation gradient F. The strain energy density function is used by a hyperelastic constitutive model to describe the stressstrain relation [16,18].If the strain energy is only a function of the initial and instantaneous state, i.e., it does not depend on the deformation history, then the hyperelastic constitutive equation describing the nonlinear elastic and compressible properties of rubber [45,46] is where P is the I. Piola-Kirchhoff stress tensor, and W(F) is the strain energy density function.The latter can be decomposed into two terms [16,47] where J = detF, C is the right Couchy-Green deformation tensor and W V (J) is the volumet- ric, while W D (C) is the deviatoric (no change in the volume) strain energy density function.

Hyperelastic Material Models
Polynomial formulas for the deviatoric strain energy density function can be used to describe the phenomenological models using the I 1 , I 2 scalar invariants of the Cauchy-Green deformation tensor, as originally proposed by Rivlin [19].Without claiming to be exhaustive, the material models Mooney-Rivlin, Yeoh, Arruda-Boyce, and Gent are available within the polynomial formula that uses invariants.The polynomial model of the strain energy density function for incompressible material is [25] where c ij material constants need to be characterized.Another approach, proposed by Valanis [48] and Ogden [21], describes the deviatoric strain energy density function using λ 1 , λ 2 , λ 3 principal stretches.

Mooney-Rivlin Model
The Mooney-Rivlin formula has good curve fitting because of the various order models that are available when choosing how many terms to include in the polynomial formula.The values H = 2 and c 20 = c 02 = 0 could be inserted into Equation (3) to derive the three-term Mooney-Rivlin model [19,20] (4)

Yeoh Model
Numerous research papers [22,49,50] demonstrate that the deviatoric strain energy for most elastomers, including rubber, depends significantly less on the I 2 second scalar invariant of the Cauchy-Green deformation tensor than on the I 1 .Thus, the Yeoh model could be written using the third-order polynom [22] W D,Y3 I 1 = W D,Y3 (c 10 , c 20 , c 30 ) = c 10 (5)

Gent Model
Given that the deviatoric strain energy is a logarithmic function of I 1 , the Gent model is an extension of the well-known Neo-Hooke material model [24] where µ is the shear modulus and J m is the dimensionless parameter that controls the finite extensibility.This extension makes the Gent model more suitable for describing the responses of elastomer-like materials to large deformations.

Ogden Model
A more accurate approximation in large strain ranges is provided by Ogden [21] using the principal strains directly in the polynomial model where λ 1 , λ 2 , λ 3 are the principal stretches, µ k is the shear modulus and α k is the dimensionless parameter.During the investigations, the commonly used H = 3 formula was used.

Arruda-Boyce Model
Arruda and Boyce used a micromechanical model and statistical mechanics methods to predict the response of the elastomer [23] where µ is the shear modulus and λ L is the limit of stretch.

Modelling the Volumetric Strain Energy Density Function
The additions in the mixture cause rubber vulcanizates to behave as nearly incompressible material under hydrostatic pressure.The formula most commonly used to express the volumetric strain energy density function to model this behaviour is [51] where κ is the bulk modulus, which is a true material coefficient.It is recommended that during the finite element analysis of a rubber bumper subjected to compressive load, the value κ = 1000 [MPa] be used, based on my prior tests [52].

Drucker's Stability
Only the nonlinear stress-strain characteristic satisfying the mathematical criteria is suitable for the material according to Drucker's stability [53,54].Examining Drucker's first stability criterion within the framework of this study is crucial.It stipulates that an increase in the external influence on the body (such as time, temperature, or load) must result in a positive change in the internal energy.Therefore, energy cannot be released, so the following criterion for the inner product of the tensors has to be fulfilled where T is the Cauchy stress tensor, while A t is the logarithmic strain tensor.Under tensile load, hyperelastic material models remain Drucker's stable even at large strains; however, in other load cases, they may become unstable even at small strains.Therefore, the stability of the material model response must be examined for the common load cases that are also used during material tests, to a pre-selected range of strain [25].To investigate the stability of the material models, the range of strain based on Table 1 was selected, considering the deformation range experienced during the rubber product operation.
Table 1.The selected range of strains for different modes of load to examine Drucker's stability.

Mode of Load Range of Strain
Uniaxial compression-tension

Compression Test of the Rubber Bumper
The compression test of the bumper was performed according to test conditions specified in Table 2 using an Instron 8801 (Norwood, MA, USA) uniaxial test machine.The measured compressive force-displacement characteristics in Figure 1 show the differences between the first and second cycles due to the viscous and Mullins effect, which no longer cause significant changes in the additional loading cycles.The rubber bumper working characteristic was determined using the data points of the fourth half cycle.

Compressive Stress-Strain Characteristics of the Rubber
The rubber bumper material is a blend of natural rubber (NR) and styrene-butadiene rubber (SBR); however, the blend ratio and the additional rubber compound elements are a trade secret.Therefore, rubber samples were machined out of the product to test the material behaviour according to ISO 23,529 using a milling machine.Figure 2 illustrates the location of sampling within the product and the prepared specimens, which have a diameter of 29 ± 0.5 (mm) and a height of 12 ± 0.5 (mm).As a reference to the rubber compound, the Shore hardness was measured according to the ISO 48-4 standard using a BAREISS Digi-Test II (Stouffville, ON, Canada) durometer, which is 77 Shore A for the test

Compressive Stress-Strain Characteristics of the Rubber
The rubber bumper material is a blend of natural rubber (NR) and styrene-butadiene rubber (SBR); however, the blend ratio and the additional rubber compound elements are a trade secret.Therefore, rubber samples were machined out of the product to test the material behaviour according to ISO 23529 using a milling machine.Figure 2 illustrates the location of sampling within the product and the prepared specimens, which have a diameter of 29 ± 0.5 (mm) and a height of 12 ± 0.5 (mm).As a reference to the rubber compound, the Shore hardness was measured according to the ISO 48-4 standard using a BAREISS Digi-Test II (Stouffville, ON, Canada) durometer, which is 77 Shore A for the test specimens and 78 Shore A for the product.

Compressive Stress-Strain Characteristics of the Rubber
The rubber bumper material is a blend of natural rubber (NR) and styrene-butadiene rubber (SBR); however, the blend ratio and the additional rubber compound elements are a trade secret.Therefore, rubber samples were machined out of the product to test the material behaviour according to ISO 23,529 using a milling machine.Figure 2 illustrates the location of sampling within the product and the prepared specimens, which have a diameter of 29 ± 0.5 (mm) and a height of 12 ± 0.5 (mm).As a reference to the rubber compound, the Shore hardness was measured according to the ISO 48-4 standard using a BAREISS Digi-Test II (Stouffville, ON, Canada) durometer, which is 77 Shore A for the test specimens and 78 Shore A for the product.The compression test of the rubber specimens was performed according to the ISO 7743 'A' method and the test conditions specified in Table 3.The tests were performed using an Instron 68TM-10 (Norwood, MA, USA) uniaxial test machine and the Instron AVE 2 non-contacting video extensometer.The operating rate of the traverse was calculated based on the strain rate used for the compression test of the rubber bumper.

Temperature and relative humidity of the test 25 (°C); 50%
The number, magnitude and speed of the preload cycles 3 cycles; 6.4 (mm); 6.76 (mm•min −1 ) The number, magnitude and speed of the measured cycle 4th cycle; 6.4 (mm); 6.76 (mm•min −1 ) Figure 2. The sampled specimens and their location within the rubber product.
The compression test of the rubber specimens was performed according to the ISO 7743 'A' method and the test conditions specified in Table 3.The tests were performed using an Instron 68TM-10 (Norwood, MA, USA) uniaxial test machine and the Instron AVE 2 non-contacting video extensometer.The operating rate of the traverse was calculated based on the strain rate used for the compression test of the rubber bumper.The measurement setup is shown in Figure 3 where to reduce friction, the compression platens were prepared by mechanical polishing and PTFE (polytetrafluoroethylene) oil lubricant.Figure 3 shows the deformation state of the specimen under maximum load, based on which it can be inferred that the barrelling effect as a result of friction during the compression test is minor.Test setup Polished and lubricated compression platens The measurement setup is shown in Figure 3 where to reduce friction, the compression platens were prepared by mechanical polishing and PTFE (polytetrafluoroethylene) oil lubricant.Figure 3 shows the deformation state of the specimen under maximum load, based on which it can be inferred that the barrelling effect as a result of friction during the compression test is minor.In this instance, there is an almost linear relationship between the displacement and force up to approximately 15% compressive deformation due to the specimen geometry and the test conditions can both be considered ideal.This behaviour was taken into account to correct data points that differ from linear caused by geometric errors.Assuming frictionless sliding between the rubber specimen and compression platens, the stress is In this instance, there is an almost linear relationship between the displacement and force up to approximately 15% compressive deformation due to the specimen geometry and the test conditions can both be considered ideal.This behaviour was taken into account to correct data points that differ from linear caused by geometric errors.Assuming frictionless sliding between the rubber specimen and compression platens, the stress is homogeneous, hence the engineering stress where F is the compressive force, A o is the original cross-sectional area.The engineering strain is expressed as the ratio of the change in dimension ∆L to the original dimension L o in the direction of the applied compressive load The true strain, which is the change in instantaneous dimension L in the direction of the applied compressive load over time The rubber behaves as a nearly incompressible material; thus, λ 1 λ 2 λ 3 ≈ 1.It allows the expression of the instantaneous cross-section A, and thus the true stress can be written Using Equations ( 11)-( 14), Figure 4 shows the calculated characteristics for specimen number 1, where an inflexion point can be found in the true stress-strain curve.Therefore, a higher-order polynomial form was selected for the hyperelastic material model wherever it was available.
Polymers 2024, 16, x FOR PEER REVIEW 8 of 19 a higher-order polynomial form was selected for the hyperelastic material model wherever it was available.The stresses were calculated for the strain  = {0: 0.01: 0.45} by linear interpolation between the measured stress values.This was undertaken to be able to calculate the average engineering stress-strain characteristic of the four rubber specimens and to reduce the points in the data set.The result is shown in Figure 5 and Table A1.The stresses were calculated for the strain ε = {0 : 0.01 : 0.45} by linear interpolation between the measured stress values.This was undertaken to be able to calculate the average engineering stress-strain characteristic of the four rubber specimens and to reduce the points in the data set.The result is shown in Figure 5 and Table A1.The stresses were calculated for the strain  = {0: 0.01: 0.45} by linear interpolation between the measured stress values.This was undertaken to be able to calculate the average engineering stress-strain characteristic of the four rubber specimens and to reduce the points in the data set.The result is shown in Figure 5 and Table A1.

Fitting the Material Parameters of the Phenomenological Hyperelastic Models Describing the Deviatoric Strain Energy Density Function
The selection of the material parameters of the hyperelastic model is an optimization task, where the objective function is given as the difference between the measured stressstrain characteristic and the one predicted by the material model.Let  be the set of the  vector of variable material coefficients and  , the measured stress for the i-th deformation state, while the () , the stress predicted by the hyperelastic model.Minimizing the objective function, which is calculated as a normalized mean absolute difference (NMAD), is one way to find the   , the vector of optimal material parameters

Fitting the Material Parameters of the Phenomenological Hyperelastic Models Describing the Deviatoric Strain Energy Density Function
The selection of the material parameters of the hyperelastic model is an optimization task, where the objective function is given as the difference between the measured stressstrain characteristic and the one predicted by the material model.Let C be the set of the c vector of variable material coefficients and σ i,e the measured stress for the i-th deformation state, while the σ(c) i,r the stress predicted by the hyperelastic model.Minimizing the objective function, which is calculated as a normalized mean absolute difference (NMAD), is one way to find the c opt , the vector of optimal material parameters where every pair of the stress-strain data set is taken into account with the same weight under the curve fitting process.The process shown in Figure 6 is called the curve-fitting method and is used to find hyperelastic material parameters.The average engineering stress-strain data set from the uniaxial compression test, shown in Figure 5 and Table A1, was used as input for the curve-fitting process.
where every pair of the stress-strain data set is taken into account with the same weight under the curve fitting process.The process shown in Figure 6 is called the curve-fitting method and is used to find hyperelastic material parameters.The average engineering stress-strain data set from the uniaxial compression test, shown in Figure 5 and Table A1, was used as input for the curve-fitting process.Table 4 lists the material parameters that were determined by fitting the phenomenological material models Mooney-Rivlin, Yeoh, Gent, Arruda-Boyce, and Ogden.Drucker's stability test of the material models was performed for the load modes given in Table 1; only the Ogden shows unstable behaviour when the biaxial compressive load mode exceeds the   = −0.42strain value.Table 4 lists the material parameters that were determined by fitting the phenomenological material models Mooney-Rivlin, Yeoh, Gent, Arruda-Boyce, and Ogden.Drucker's stability test of the material models was performed for the load modes given in Table 1; only the Ogden shows unstable behaviour when the biaxial compressive load mode exceeds the ε bia = −0.42strain value.

Verification of Material Models Using the Finite Element Model of the Compression Test
The finite element simulation of the uniaxial compression test with the material parameters in Table 4 was used to evaluate the accuracy of the fitted material models.Axisymmetric linear quadrilateral element (PLANE183) with the settings shown in Table 5 was chosen for the finite element discretization due to the isotropic material, the axisymmetric geometry and boundary conditions.It should be noted that during the compression tests, it is assumed that the stress distribution in the rubber specimen is homogeneous.Consequently, the stress response can be obtained precisely with just one element.4.,κ = 1000 (MPa) [52]   At nodes 1 and 2 on the upper edge of the test specimen, the prescribed displacement of 5.625 (mm) (45% compression) in direction −y is applied as a load; see UY in Figure 7. Furthermore, to model the frictionless connection between the specimen and the compression platen, the displacement is constrained in the y direction at nodes 3 and 4 on the lower edge; see TY in Figure 7.
Post-processing the biaxiality indication could be one way to interpret the strain state of a rubber product under compressive load.The ratio of the smaller to larger principal stress-ignoring the principal stress closest to zero-defines the biaxiality indication.Hence, the locations of the uniaxial stress state are reported by a value of zero, pure shear by minus one, and biaxial by one.According to Figure 7, the numerical test results of the compression test show nearly zero biaxiality indication to the whole specimen, which was expected under the uniaxial load case.
of a rubber product under compressive load.The ratio of the smaller to larger principal stress-ignoring the principal stress closest to zero-defines the biaxiality indication.Hence, the locations of the uniaxial stress state are reported by a value of zero, pure shear by minus one, and biaxial by one.According to Figure 7, the numerical test results of the compression test show nearly zero biaxiality indication to the whole specimen, which was expected under the uniaxial load case.The ANSYS Mechanical solver and Newton-Raphson numerical technique were used to solve the nonlinear problem.The number of sub-steps was selected in the analysis setup to match the number of data points on the average engineering σ-ε characteristic.The engineering stress-strain characteristics were computed numerically using the fitted hyperelastic models; see Figure 8 for the results.The ANSYS Mechanical solver and Newton-Raphson numerical technique were used to solve the nonlinear problem.The number of sub-steps was selected in the analysis setup to match the number of data points on the average engineering σ-ε characteristic.The engineering stress-strain characteristics were computed numerically using the fitted hyperelastic models; see Figure 8 for the results.Figure 9 shows the errors of the stresses predicted to the specific elongation values  = {0: 0.01: 0.45} of the fitted hyperelastic models that were determined relative to the stresses measured on the laboratory test.Based on these results, in the range of strains above 10%, the Yeoh and Mooney-Rivlin material models approximate the measured stress-strain characteristic of the rubber within the error of less than 5%.The prediction of the Mooney-Rivlin and Yeoh material models differ significantly from rubber's linear behaviour in the strain region below 10%.This difference can be attributed to the polynomial function approximation used to approximate the deformation energy density.The reason for this phenomenon can be caused by the use of the polynomial form to approximate the deviatoric strain energy density function.Figure 9 shows the errors of the stresses predicted to the specific elongation values ε = {0 : 0.01 : 0.45} of the fitted hyperelastic models that were determined relative to the stresses measured on the laboratory test.Based on these results, in the range of strains above 10%, the Yeoh and Mooney-Rivlin material models approximate the measured stress-strain characteristic of the rubber within the error of less than 5%.The prediction of the Mooney-Rivlin and Yeoh material models differ significantly from rubber's linear behaviour in the strain region below 10%.This difference can be attributed to the polynomial function approximation used to approximate the deformation energy density.The reason for this phenomenon can be caused by the use of the polynomial form to approximate the deviatoric strain energy density function.
stress-strain characteristic of the rubber within the error of less than 5%.The prediction of the Mooney-Rivlin and Yeoh material models differ significantly from rubber's linear behaviour in the strain region below 10%.This difference can be attributed to the polynomial function approximation used to approximate the deformation energy density.The reason for this phenomenon can be caused by the use of the polynomial form to approximate the deviatoric strain energy density function.Figure 10 shows the mean relative absolute error (MRAE) calculated using the values of the relative errors previously determined in Figure 9 can be used to determine generally the goodness of the material model's prediction capability.The Gent, Arruda-Boyce and Ogden material models have close to 10% MRAE, mostly because the predicted stress values in the upper strain range are more accurate.Thus, it is not recommended to use these models to predict the material stress-strain behaviour over the whole range of a given Figure 10 shows the mean relative absolute error (MRAE) calculated using the values of the relative errors previously determined in Figure 9 can be used to determine generally the goodness of the material model's prediction capability.The Gent, Arruda-Boyce and Ogden material models have close to 10% MRAE, mostly because the predicted stress values in the upper strain range are more accurate.Thus, it is not recommended to use these models to predict the material stress-strain behaviour over the whole range of a given strain.The results showed that the Mooney-Rivlin and Yeoh material models can predict the material response to the compressive load with less than 2.5% MRAE.strain.The results showed that the Mooney-Rivlin and Yeoh material models can predict the material response to the compressive load with less than 2.5% MRAE.

Verification of Material Models Using the Finite Element Model of the Rubber Bumper Compression Test
The finite element analysis of the compression test of the rubber bumper was used to evaluate the prediction capability of the hyperelastic models, fitted for the compression test, for different strain states.The inhomogeneous strain state within the rubber product under compressive load is due to the complex geometry, large deformation and connections.This investigation requires post-processing of the numerical analysis results for each of the material models to determine the predicted working characteristics of the product.These results must then be compared to the laboratory test results shown in Figure 1. Figure 11 shows the geometric dimensions and the meridian section of the compression platens and the product.

Verification of Material Models Using the Finite Element Model of the Rubber Bumper Compression Test
The finite element analysis of the compression test of the rubber bumper was used to evaluate the prediction capability of the hyperelastic models, fitted for the compression test, for different strain states.The inhomogeneous strain state within the rubber product under compressive load is due to the complex geometry, large deformation and connections.This investigation requires post-processing of the numerical analysis results for each of the material models to determine the predicted working characteristics of the product.These results must then be compared to the laboratory test results shown in Figure 1. Figure 11 shows the geometric dimensions and the meridian section of the compression platens and the product.
The finite element analysis of the compression test of the rubber bumper was used to evaluate the prediction capability of the hyperelastic models, fitted for the compression test, for different strain states.The inhomogeneous strain state within the rubber product under compressive load is due to the complex geometry, large deformation and connections.This investigation requires post-processing of the numerical analysis results for each of the material models to determine the predicted working characteristics of the product.These results must then be compared to the laboratory test results shown in Figure 1. Figure 11 shows the geometric dimensions and the meridian section of the compression platens and the product.Mesh locking under large deformations is a phenomenon that can be observed in meshes discretized with nearly incompressible material.In such a case, for the meshing, linear quadrilateral elements that are less prone to mesh locking must be used.To reduce the numerical error in the determination of the operating characteristic of the rubber bumper a force convergence test was carried out while changing the mesh density.Investigating Figure A1, by choosing an element size of 1 (mm), the numerical solution error is under 0.5%.Table 6 contains the settings for the finite element discretization, while the resulting finite element mesh for the bumper is seen in Figure 12.  4., κ = 1000 (MPa) [52]   The rubber bumper contacts the upper and lower compression platens as it operates while resulting in a 30% reduction in height.Hence, A and B frictional edge-to-edge connections were defined, where the contact elements were created on the rubber bumper, while the target elements were created on the compression platens according to Figure 12.The static friction coefficient was selected as s = 0.6 [47]; thus, the state of contact between the contacting elements is mostly sticking under the compression test.At the edge of the upper compression platen, the prescribed displacement of 28 (mm) in direction −y is applied as a load; see UY in Figure 12.Furthermore, the displacement is constrained in the x and y directions at the edge of the lower compression platen; see TX, TY.
while the target elements were created on the compression platens according to Figure 12.The static friction coefficient was selected as µ  = 0.6 [47]; thus, the state of contact between the contacting elements is mostly sticking under the compression test.At the edge of the upper compression platen, the prescribed displacement of 28 (mm) in direction -y is applied as a load; see UY in Figure 12.Furthermore, the displacement is constrained in the x and y directions at the edge of the lower compression platen; see TX, TY.After the converged nonlinear numerical analysis, Figure 12 shows the deformation state under maximum compressive displacement, allowing the appropriate operation of the discretized model and connections to be verified.Additionally, the stress biaxiality After the converged nonlinear numerical analysis, Figure 12 shows the deformation state under maximum compressive displacement, allowing the appropriate operation of the discretized model and connections to be verified.Additionally, the stress biaxiality indication in the contour plot shows, in line with previous intuition, that the stress state of the rubber bumper varies significantly by location.In a significant part of the inside of the rubber bumper, the state of stress is nearly uniaxial; however, close to the compression platen where the contact is sticking, it is biaxial, while shear also occurs near the diameter.In locations where the biaxiality indicator deviates significantly from zero, the prediction of the material models may differ significantly from the material real behaviour, as the fitted hyperelastic material models lack information for load cases other than the uniaxial compression.After running the finite element analysis with the different hyperelastic models, Figure 13 shows to what extent the predicted characteristics differ from the measured one.indication in the contour plot shows, in line with previous intuition, that the stress state of the rubber bumper varies significantly by location.In a significant part of the inside of the rubber bumper, the state of stress is nearly uniaxial; however, close to the compression platen where the contact is sticking, it is biaxial, while shear also occurs near the diameter.
In locations where the biaxiality indicator deviates significantly from zero, the prediction of the material models may differ significantly from the material real behaviour, as the fitted hyperelastic material models lack information for load cases other than the uniaxial compression.After running the finite element analysis with the different hyperelastic models, Figure 13 shows to what extent the predicted characteristics differ from the measured one.Figure 14 shows the error of the reaction forces predicted by the hyperelastic models for ten uniformly distributed compression levels relative to compressive force from the laboratory test.The Gent, Arruda-Boyce and Ogden material models have large deviations in the simulated compressive force under almost the whole compressive displace-   Figure 14 shows the error of the reaction forces predicted by the hyperelastic models for ten uniformly distributed compression levels relative to compressive force from the laboratory test.The Gent, Arruda-Boyce and Ogden material models have large deviations in the simulated compressive force under almost the whole compressive displacement.In comparison to the characteristics measured and simulated by the Yeoh material model, the numerical result analysed with the Mooney-Rivlin material model shows a stiffer behaviour in almost the whole compressive displacement.This is due to the strain energy based hyperelastic model fitted to the uniaxial compressive load case, which cannot accurately predict material responses different from the known strain state.A deviation of more than 10% can be observed when investigating the MRAE shown in Figure 15 for the Gent, Arruda-Boyce, and Ogden models.In comparison to the values observed during the finite element run of the compression test, the MRAE of the Yeoh material model increased by 3% and the Mooney-Rivlin increased by 7% based on simulation results shown in Figure 10.Nevertheless, the material response to an inhomogeneous strain state can be estimated by the Mooney-Rivlin and Yeoh material models with an error significantly below 10%, which is acceptable in engineering simulation problems.A deviation of more than 10% can be observed when investigating the MRAE shown in Figure 15 for the Gent, Arruda-Boyce, and Ogden models.In comparison to the values observed during the finite element run of the compression test, the MRAE of the Yeoh material model increased by 3% and the Mooney-Rivlin increased by 7% based on simulation results shown in Figure 10.Nevertheless, the material response to an inhomogeneous strain state can be estimated by the Mooney-Rivlin and Yeoh material models with an error significantly below 10%, which is acceptable in engineering simulation problems.

Conclusions
Designing rubber products requires special consideration of the numerical discretization process due to the nonlinear behaviours (material nonlinearity, large deformation, connections).Modelling considerations were presented for the finite element analysis of

Conclusions
Designing rubber products requires special consideration of the numerical discretization process due to the nonlinear behaviours (material nonlinearity, large deformation, connections).Modelling considerations were presented for the finite element analysis of the rubber bumper.If the rubber behaviour is only determined for the specific load of the product, which, in the case of rubber bumpers, is the compression, the time needed for the laboratory test can be significantly decreased.
Based on the foregoing, my research aimed at fitting and selecting a hyperelastic material model suitable for the finite element modelling of the working characteristic of a rubber product under compressive load.The Mullins effect as well as the viscoelastic and hysteretic material behaviours were not modelled because the time dependence was not taken into consideration during the product testing of the rubber bumpers.
Rubber samples for compression tests were machined out of the product because the content of the rubber compounds is a trade secret.The stress-strain characteristics of the uniaxial compression test of rubber samples were used to fit the hyperelastic material models (Mooney-Rivlin, Yeoh, Gent, Arruda-Boyce, and Ogden).The results showed that the Mooney-Rivlin and Yeoh material models fitted for the uniaxial compression data set can predict the material response up to 45% uniaxial compressive deformation with less than 2.5% MRAE if the stress state is homogeneous.Nevertheless, the Gent, Arruda-Boyce and Ogden material models are unable to properly predict material response even in a homogeneous stress state.The resulting 10% MRAE is mostly because the predicted stress values in the upper strain range are more accurate; therefore, the use of these material models in the investigated strain range is not recommended in the case of compressive load.The laboratory and numerical tests of the rubber bumper were used to determine the accuracy of the material models if the stress state is changing.The distribution of the stress biaxiality in the contour plot is a good measure to determine which load mode is required for the curve fitting of the hyperelastic model.The results showed that if only uniaxial compression test data are available for the curve fitting of the material model, the Yeoh model performs the best in predicting how the rubber product behaves under compressive load and complex strain state.
The research pointed out that the constitutive modelling error in complex strain state could be significantly below 10% when using only the uniaxial compression data set to calibrate the Yeoh or Mooney Rivlin material models.This is acceptable in nonlinear engineering simulation problems.Non-standard mechanical tests like simple shear and biaxial could therefore be eliminated, saving a significant amount of time and money, which would help increase market competitiveness.
The developed process opens up numerous new research possibilities.One area is the investigation of other load modes or different strain ranges, for which the product can be selected using the introduced biaxiality measure.Other hyperelastic materials like 3Dprintable thermoplastic polyurethane blends could be investigated to shorten the specimen preparation and the product manufacturing time.Using the introduced biaxiality measure, an extension of the research could be to develop a strain distribution-based material model calibration method.

Temperature and Relative Humidity of the Test 25 (
• C); 50% The number, magnitude and speed of the preload cycles 3 cycles; 35 mm; 50 (mm•min −1 ) The number, magnitude and speed of the measured cycle 4th cycle; 35 (mm); 50 (mm•min −1 ) Test setup Polished compression platen without lubrication

Figure 1 .
Figure 1.Uniaxial compression test setup and the measured characteristics of the rubber bumper.

Figure 1 .
Figure 1.Uniaxial compression test setup and the measured characteristics of the rubber bumper.

Figure 1 .
Figure 1.Uniaxial compression test setup and the measured characteristics of the rubber bumper.

Figure 2 .
Figure 2. The sampled specimens and their location within the rubber product.

Figure 3 .
Figure 3. Uniaxial compression test setup and the deformation state of the specimen under maximum load.

Figure 3 .
Figure 3. Uniaxial compression test setup and the deformation state of the specimen under maximum load.

Figure 4 .
Figure 4.The engineering and true σ-ε characteristics under uniaxial compression test for specimen number 1.

Figure 4 .
Figure 4.The engineering and true σ-ε characteristics under uniaxial compression test for specimen number 1.

Figure 5 .
Figure 5.The average engineering σ-ε characteristic under uniaxial compression test for the four rubber specimens.

Figure 5 .
Figure 5.The average engineering σ-ε characteristic under uniaxial compression test for the four rubber specimens.

Figure 6 .
Figure 6.The process of curve-fitting method to find the hyperelastic material parameters.

Figure 6 .
Figure 6.The process of curve-fitting method to find the hyperelastic material parameters.

Figure 7 .
Figure 7.The 2D axisymmetric discretized model of the compression specimen indicating the boundary conditions (prescribed displacement UY in direction -y, roller support TY) and the postprocessed deformation state under 45% prescribed compressive load with the distribution of the stress biaxiality in the contour plot.

Figure 7 .
Figure 7.The 2D axisymmetric discretized model of the compression specimen indicating the boundary conditions (prescribed displacement UY in direction −y, roller support TY) and the postprocessed deformation state under 45% prescribed compressive load with the distribution of the stress biaxiality in the contour plot.

19 Figure 8 .
Figure 8.Comparison of the different hyperelastic models' prediction with the laboratory test results of the average engineering σ-ε characteristic.

Figure 8 .
Figure 8.Comparison of the different hyperelastic models' prediction with the laboratory test results of the average engineering σ-ε characteristic.

Figure 9 .
Figure 9.The variation in the relative error of the different hyperelastic models compared to the average σ-ε characteristic measured on the specimens.

Figure 9 .
Figure 9.The variation in the relative error of the different hyperelastic models compared to the average σ-ε characteristic measured on the specimens.

Figure 10 .
Figure 10.The mean relative absolute error of the hyperelastic model prediction capability relative to the average σ-ε characteristic.

Figure 10 .
Figure 10.The mean relative absolute error of the hyperelastic model prediction capability relative to the average σ-ε characteristic.

Figure 11 .
Figure 11.Geometric dimensions and the meridian section of the compression platens and the rubber product.

Figure 11 .
Figure 11.Geometric dimensions and the meridian section of the compression platens and the rubber product.

Figure 12 .
Figure 12.The 2D axisymmetric discretized model of the rubber bumper indicating the boundary conditions and contact regions and the post-processed deformation state under 30% prescribed compressive load with the distribution of the stress biaxiality in the contour plot.

Figure 12 .
Figure 12.The 2D axisymmetric discretized model of the rubber bumper indicating the boundary conditions and contact regions and the post-processed deformation state under 30% prescribed compressive load with the distribution of the stress biaxiality in the contour plot.

Figure 13 .
Figure 13.The laboratory and the numerical test results of the rubber bumper working characteristic.

Figure 13 .
Figure 13.The laboratory and the numerical test results of the rubber bumper working characteristic.

Figure 14
Figure14shows the error of the reaction forces predicted by the hyperelastic models for ten uniformly distributed compression levels relative to compressive force from the laboratory test.The Gent, Arruda-Boyce and Ogden material models have large deviations in the simulated compressive force under almost the whole compressive displacement.In comparison to the characteristics measured and simulated by the Yeoh material model, the numerical result analysed with the Mooney-Rivlin material model shows a stiffer behaviour in almost the whole compressive displacement.This is due to the strain energy based hyperelastic model fitted to the uniaxial compressive load case, which cannot accurately predict material responses different from the known strain state.

Figure 13 .
Figure 13.The laboratory and the numerical test results of the rubber bumper working characteristic.

Figure 14 .
Figure 14.The variation in the relative error of the different hyperelastic models compared to the characteristic measured on the rubber bumper.

Figure 14 .
Figure 14.The variation in the relative error of the different hyperelastic models compared to the characteristic measured on the rubber bumper.

Figure 15 .
Figure 15.The mean relative absolute error of the hyperelastic models compared to the characteristic of the rubber product.

Figure 15 .
Figure 15.The mean relative absolute error of the hyperelastic models compared to the characteristic of the rubber product.

Funding:
This research was supported by the ÚNKP-23-4-II.New National Excellence Program of the Ministry for Culture and Innovation from the source of the National Research, Development and Innovation fund.Institutional Review Board Statement: Not applicable.

Table 2 .
Test conditions of the uniaxial compression of the rubber bumper.

Table 3 .
Test conditions of the uniaxial compression of the cylindrical rubber specimens.

Table 3 .
Test conditions of the uniaxial compression of the cylindrical rubber specimens.

Table 4 .
Hyperelastic material parameters found by the curve fitting method.

Table 4 .
Hyperelastic material parameters found by the curve fitting method.

Table 5 .
Settings for the finite element discretization.

Table 6 .
Settings for the finite element discretization.