Large Strain Viscoelastic Material Model for Deformation, Stress and Strain Analysis of O-rings

Primary aim of this paper is to investigate the mechanical behavior of an O-ring seal made of NBR rubber material using finite element method. Further aim is to identify the O-ring ’ s critical locations where it can fail under normal operating conditions. For these a large-strain viscoelastic material model was set up. The recommended strain range where the material model should be sufficient is 0-100 % for normal conditions, but when the effect of overload or the temperature is the question, the strain range should be 0-150 %. From designer point of view, a large-strain viscoelastic material model with a 5 param eter Mooney-Rivlin term is appropriate. To make sure the mate rial model ′ s accurate response, uniaxial tension, pure shear and equibiaxial measurements are needed, completed with dynamic measurements with a DMTA machine. The results show two crit ical locations where the failure of the O-ring is expected. One is on the edge of the O-ring near the gap, when friction is present, the other one is under the edge near the gap.


Introduction
Because of their versatility and simple manufacturability, O-rings are the most widely used type of seals in the industry recently. They are used in both static and dynamic applications and play an important role in the automotive-, shipping-, agricultural machinery-, building machinery-, energy-, and healthcare industry as well as in household appliances. Static seals are used between two stationary surfaces. Dynamic seals are installed between surfaces, which are either rotating or reciprocating relative to each other. The failure of these seals can lead to malfunction of the equipment and cause serious consequences. Fig. 1 (a) shows the O-ring before installation. O-rings fit in the groove and the housing with overlap, which provides initial preload after assembly, Fig. 1 (b). This preload and the operating pressure act together creating the sealing effect. The O-ring is pressed to the side of the groove and as it deforms, the pressure is transferred to the surrounding surfaces Fig. 1 (c). During the assembly, the primary load of the O-ring is compression. After assembly as the operating pressure is applied, the O-ring deforms into the gap between the shaft and the housing, experiencing shear load and gets into a complex stress state.
To characterize the rubber′s response in a complex stress state the rubber had to be tested in different ways. For the numerical simulation uniaxial tension, pure shear, and equibiaxial tension tests were carried out, followed by the measurement of the dynamic properties using Dynamic Mechanical Thermal Analysis. The aim of the present study is, to develop FE models to evaluate the deformation, strain and stress conditions during operation of O-rings.

Short literature review of failure of O-rings
Many researchers have investigated the behavior and failure of O-rings. Karaszkiewicz studied and calculated analytically the deformed geometry, the arising contact pressure and the contact force of O-rings [2]. Zhang and Zhang carried out finite element calculations to study the sealing performance of both static and dynamic O-rings [3]. One typical failure of O-rings is caused by the too rough surfaces of the contacting elements, which leads to excessive wear as shown in Fig. 2 (a). Gawliński investigated the influence of friction and wear of static and dynamic elastomer seals [4]. Overpressure or too large gap cause extrusion of the material. Unequivocal sign of this failure mode is the sharp edges appearing on the O-ring as shown in Fig. 2 (b). Eshel studied the extrusion of O-rings and proposed a theoretical model to estimate the failure in this mode [5]. In Fig. 2 (c) circumferential cracks can be seen on the O-ring seal caused by excessive initial compression. There are several other factors that can lead to the deterioration of the O-ring [1,6,7]. Excessive remaining compression set, which is probably the most common failure type of O-rings, is caused by too high temperatures, low heat resistance of the material, poor compression-set properties of the elastomer, swelling of O-ring material, excessive initial compression or operating pressure. Chemical failure caused by the incompatibility of the contacting materials which causes change of the rubber′s physical properties. Damage during installation. Excessive temperature causes cracks in radial direction on the O-ring.

Material models
For the FE simulations a large-strain viscoelastic material model were used. The model is based on the extensions of the linear viscoelastic so called generalized Standard-Solid model [8]. The extensions required a 5 parameter Mooney-Rivlin hyperelastic term [8] for the material model.

Characterizing the rubber behavior
In order to characterize the rubber material′s stress-strain response, tensile tests were carried out in uniaxial tension, pure shear and equibiaxial tension mode by a universal tensile tester. The tension test was carried out according to ISO 37 [1,9]. For the pure shear test a custom made specimen was used. The equibiaxial test was carried out by the Institute of Polymer Product Engineering at Johannes Kepler University in Linz, with a custom made biaxial device [10]. The biaxial test was used instead of the uniaxial compression test, because the results of the uniaxial compression tests greatly depend on the coefficient of friction between the rubber sample and the specimen grips [11]. The tensile speed in all three modes was 100 mm/s. Fig. 3 shows the dimensions and the shape of the specimens used for the tests. Fig. 4 shows the engineering stress-engineering strain curves obtained by the tensile tests.
To capture the NBR 70 rubber's viscoelastic response Dynamic Mechanical Thermal Analysis (DMTA) was carried out in tension mode on a GABO EPLEXOR 500N DMTA machine. During this test a specimen was stretched by a static load and then a sinusoidal dynamic load was superimposed on it. The static and dynamic load as an excitation is not in phase with the reaction force developed in the rubber [12]. This phase shift characterizes the viscoelastic response of the rubber. During this analysis the machine measures the storage modulus (E′), the loss modulus (E″), the complex modulus (E*) and the loss factor (tanδ) of the material. The DMTA test was carried out from -100 ºC to +80 ºC with 4 ºC temperature steps. At each temperature the specimen was excited with (b) (c) (a) 15 different frequencies from 0.2 Hz to 60 Hz. The rectangular cross section of the specimen was 5.96 mm by 1.79 mm. The distance between the specimen grips was 20 mm. Before each measurement the test specimen was preloaded with 2 N force, and the distance between the tensile specimen grips was recorded. These recorded data then served as the basis of the applied loads. The static load was 1% the dynamic load was 0.1% strain of the measured distances between the grips at each measurement points. The reason to use preload was to compensate the deformation of the specimen caused by the grips and to compensate the material's thermal expansion, thus to enhance the precision of the measurement. At each temperature the soak time was 600 s. Fig. 5 shows the E′ storage modulus and tanδ loss factor as a results of the measurement.   Fig. 6 shows the small-strain viscoelastic generalized Standard-Solid model, which is the extension of the Standard-Solid model used to describe stress-relaxation in weakly crosslinked elastomers [13]. It consists of a spring and series of a spring and a dashpot element connected in parallel with each other. The constitutive model is given by [14]:

Parameter identification
where σ is the Cauchy stress, e is the deviatoric strain, Δ is the volumetric strain, t is the past time, I is the identity tensor, G(t) and K(t) are the Prony series shear and bulk-relaxation moduli. The first part of (1) is the deviatoric part. It gives the stress arises from the shape change of the material. The second part is the volumetric part that gives the stress arises from volume change of the material. The rubber material was assumed to be incompressible therefore the volumetric response of the material was neglected. Because only dynamic frequency data is available for the rubber material the Prony series shear relaxation modulus can be determined from the following implicit equations [15,16]: where G′(ω) is the shear storage modulus, G″(ω) is the shear loss modulus, G* is the complex shear relaxation modulus, G 0 is the glassy shear modulus, g i are the relative shear moduli, τ i are the relaxation times for a given Prony series and N is the number of Prony terms. Assuming that the material is incompressible, the results of the DMTA measurement in tensile mode were converted into shear moduli: where G* is the complex shear relaxation modulus, E* is the measured complex tensile relaxation modulus and υ is the Poisson's ratio, which was assumed to be 0.5. To create the relaxation master curve the measured moduli had to be plotted in function of the frequency, and the measured values had to be (1) shifted at each temperature according to William-Landel-Ferry (WLF) theory. The shift factors can be calculated as: where C 1 and C 2 are the WLF constants, T ref is the reference temperature and T is the temperature which results is shifted. The reference temperature was set to 20 °C. Fig. 7 shows the measured storage modulus in function of the frequency. The constructed master curve can be seen in Fig. 8.    9 shows the measured and the calculated shift factors. The calculated factors show good correlation with the measured data. Moderate difference can be seen between the two curves below the material's glass transition temperature. Because rubber materials become rigid below the glass transition temperature (Tg) -their storage modulus can increase up to three orders of magnitude -they are not to be operated below Tg. In this case Tg≈ -60 °C. Thus the calculated shift functions can be accepted. The resulting difference could be reduced by using Arrhenius supplement of the WLF theory.
For the parameter fitting the shear moduli were calculated using (5) and the parameters were obtained by using the least squares method. Fig. 10 shows the shear storage and loss moduli and the loss factor obtained from the DMTA measurement, represented with continuous lines, and the calculated moduli and loss factor using (2) and (3) represented with dashed lines. All in all, a 55 term generalized Standard-Solid model was determined. It can be seen that the model shows a perfect match with the measured data in terms of the shear storage modulus. The model underestimates the loss modulus at low and at very high frequencies. Therefore, the ratio of the loss and storage moduli, the loss factor is slightly underestimated at low and at high frequencies, but follows the curve of the measured data properly. The small-strain viscoelastic model was then extended with a hyperelastic term. A 5 parameter Mooney-Rivlin term was used in order to best fit the results of the static tests. The resulting large-strain viscoelastic model is the modification of the model originally proposed by Simo [14,17].

Verification of the material model
To verify the behavior of the material model, FE simulatios were performed in uniaxial tension, pure shear and equibiaxial tension modes. The modelled specimen's dimensions are in accordance with the ones seen in Fig. 3. For the finite element simulation 1/8-th of the test specimens were modeled to reduce computational time. The tensile speed in all three cases was 100 mm/s. In Fig. 11 one can see the modelled specimens.  Fig. 12 shows the response of the model compared to the measurements. The material model supplemented with the 5 parameter Mooney-Rivlin term was fitted in the region of 0 -150 % strain and gives a good correlation with the measured data. The model follows the inflexion of the measured curves though it slightly overrates the rubbers response at moderate strains.

Finite element model of the O-ring developed
For the finite element simulation of the O-ring an axisymmetric model was built. The model consists of three parts namely the O-ring, the housing and the shaft, which was modeled as a tube shaft. This pipe-like shaft geometry helped to reduce the number of elements used in the FE model and hence the simulation time. The nominal dimensions of the parts and their tolerances can be seen in Fig. 13 and Fig. 14. First, a FE simulation based on the nominal dimensions was carried out. Secondly, the tolerances of the dimensions were taken into account, the worst and the best case scenarios was modelled from the point of view of the arising sealing pressure. Table 1 summarizes the main dimensions, the corresponding tolerances and the extreme values of the FE model. The tolerance of the inside diameter of the O-ring (d1) was calculated according to ISO 3601-1 [18]. The dimensions marked with blue correspond to the maximal gap (minimal sealing pressure) and the ones marked with red to the minimal gap (maximal sealing pressure) case.
The O-ring was chosen to have 101.19 mm inside diameter and 3.53 mm cross sectional diameter, which is an O-ring with commonly used dimensions. The dimensions of the housing were chosen so that the static pre-stress of the O-ring is 15 % with nominal dimensions. This 15 % is the lower limit of the pre-stress recommended by the manufacturers for O-rings made of NBR material with 3.53 mm cross sectional diameter [1,7]. Considering the tolerances the minimal pre-stress is 10.2 % and the maximal is 16.9 %. The resulting gap between the housing and the shaft is 0.08 mm with nominal dimensions, which is a recommended value for this O-ring type [7]. With the tolerances taken into account the gap is 0.1425 mm and 0.098 mm for the worst and best case scenarios.      Fig. 16 shows the time curve of the FE simulation. First, the O-ring is installed into its groove. ANSYS automatically softens the normal contact stiffness between the O-ring and the contacting elements. As a result, the software eliminates the overlap and the O-ring deforms into the groove by the time of 1 s. The operating pressure is gradually applied on the lower semicircle of the O-ring by the end of 2 s, then the pressure acts for 1.5 hours, hence the end of the simulation is at t=5402 s in order to model the time-dependent behavior.

Finite element results of the O-ring in operation 6.1 Results for nominal case
The maximal displacement at t=1 s happens to be at the side of the O-ring where it touches the shaft and the housing, as seen in Fig. 17. This is because of the initial compression of the sealing. After applying the operating pressure, the maximal displacement occurs at the upper left and right corners of the O-ring. This is because the O-ring fills out the groove and touches the groove's upper side and its fillet. The displacement field is more or less symmetric, but the effect of the gap between the housing and the shaft distorts the symmetry. After 1.5 hours the total displacement of the O-ring is significantly greater close to the gap than in the upper right corner. The cause of this phenomenon is that the rubber material shows viscoelastic behavior and penetrates into the gap under the influence of constant pressure. Because of the great displacement and the presence of friction, high strain values are expected in this area, which can be an indicator where damage of the O-ring may occur. High displacement values can also be seen at the side of the O-ring where it is connecting to the shaft. This indicates that the contact pressure is higher in this region. Fig. 18 shows the equivalent true strain of the O-ring. It can be seen, that the maximal equivalent strain occurs to be on the edges of the seal near the gap. This is because the friction between the rubber O-ring and the housing or the shaft try to prevent the seal to penetrate into the gap. Because of the friction and the constant applied pressure, the material elongates along the contacting edges and high strain and stress occur. The maximum values is 0.488 mm/mm. The extension of the high strain area is local and caused by the friction acting between the O-ring and the housing or the shaft, and the gap between the housing and the shaft. Fig. 19 shows the maximum principal true strain of the O-ring. One can also see similar behavior compared to the equivalent true strain results. In fact, for maximal principal true strain, shear intensity, octahedral shear strain, and strain energy density the results showed similar distribution. The above mentioned factors differed only in values from each other. Thus, only the maximal principal true strain distribution is shown. Fig. 20 shows the detailed view of the O-ring in the vicinity of the gap for the results presented above. The effect of friction can clearly be seen in Fig. 20 (b) and (c). Fig. 21 shows the contact pressure distribution at t=1, t=2 and t=5402 s. Fig. 21 (a) illustrates, that because of the initial compression the O-ring functions as a seal before the operating pressure is applied. In Fig. 21 (b) one can see that the operating pressure is superimposed on the initial pressure, because the pressure propagates the material without weakening as it is nearly incompressible, although as the material shows Fig. 16 Time curve of the installation and the pressure build-up procedure viscoelastic behavior and relaxes the maximal contact pressure will be slightly less than the sum of the two pressures. Between t=2 s and t=5402 s as the material relaxes, the maximum values of the contact pressure are decreasing (Fig. 21 (c)). The overall contact area increases as the O-ring deforms. The maximum value of the contact pressure is at the side of the O-ring where it touches the shaft. Table 2 summarizes the results.

Influence of the tolerances
The influence of the tolerances was also investigated for only the material model with the 5 parameter Mooney-Rivlin term. The resulting gap between the shaft and the housing is 0.1425 mm for the maximal and 0.098 mm for the minimal gap case instead of the 0.08 mm, which occurs with nominal dimensions. This means that because of the material's viscoelastic behavior and the larger gap, the O-ring is expected to penetrate more into the gap thus higher strain levels are expected near the gap on the edge of the O-ring. For maximal gap case, the maximal equivalent true strain occurred to be higher as expected, namely 0.532 mm/mm. For minimal gap case, the maximal equivalent strain value is 0.514 mm/mm compared to the value of 0.488 mm/mm for nominal case. The maximal contact pressure values increased around 2% for minimal gap case and decreased approximately 5% for maximal gap case. In Table 3, one can see the maximal values of the different parameters for nominal, minimal and maximal gap case.

Considering frictionless case
There are many situations where the friction between the O-rings and the contacting elements varies between 0 and its' maximum. These situations develop because of good lubrication, dynamic applications and vibrations from the surrounding machine elements. To determine the effect of friction the FE model with nominal dimensions and with the 5 parameter Mooney-Rivlin large-strain viscoelastic material model was investigated at 50 bars. To represent the frictionless case the applied coefficient of friction was µ=0.001. The reason for this was to enhance convergence during the FE simulation.
Compared to the model with µ=0.2 the lack of friction lets the O-ring penetrate more into the gap. Examining the strains on the other hand show that the O-Ring does not experiences high strain values on its' circumference, the maximum values there, are a lot lower than the case with friction ( Table 3). The area with the highest strain values is near the gap but not on the edge but under the edge of the O-ring. In Fig. 22 one can see that maximum value of equivalent true strain is 0.46 mm/mm instead of 0.488 mm/mm when friction is acting.
This shows the importance of lubrication during operation. Fig. 23 presents the influence of the gap for frictionless case. It is seen, that for maximal gap size the maximum value of the equivalent true strain increased about 4 % compared to the nominal case. The sealing pressure increased 3-5% compared to frictional case. Table 4 summarizes the occurring strain values in the different zones of the O-ring according to Fig. 23.  (Fig. 20 (a)). The reason of this behavior is the fact that the seal penetrates into the gap and fills it out. This is a place where high stress and strain values were found hence failure of the O-ring could be initiated from here ( Fig. 18 (c)). Another area of high strains was found on the edges of the O-ring near the gap (Fig. 20 (b)). These high strains caused by the friction between the O-ring and the contacting elements. Friction tries to prevent the O-ring to penetrate into the gap and stretches the material. This means that friction pays an important role in the arising strains and stresses. The failure of the rubber is expected near the gap on the edge of the O-ring, or if the friction is not significant, close to the gap, in the material. Summarizing the results it is clear that the cause of dynamic failure of O-rings is the presence of the gap between the shaft and the housing because it alters the strain distribution of the rubber. The influencing factors are the tolerances of the dimensions of shaft and the housing that are related to the gap size, the surface roughness of the housing and the shaft, the lubrication of the O-ring as they effect the friction conditions, the operating pressure, the temperature and the material of the O-ring. 5 Normally O-rings experience volumetric compression load as they are forced in the groove, where they are pushed to the side when pressure is applied. Because there is no room to deform the primary load is volumetric compression. On the other hand near the gap the rubber experiences high strains and shear stresses as it penetrates into the gap. 6 Considering the frictionless case, the O-ring can deeper penetrate into the gap and the strain is much lower in Zone A and B due to the missing local frictional force.