Modelling Strategy and Parametric Study of Metal Gaskets for Automotive
Applications

This paper is focused on finite element simulation of cylinder head gaskets. Finite element codes support several methodologies, each of which has its own strengths and weaknesses. One of the key points lies in the influence of the gasket geometry on its final behaviour. Such a contribution can come from the detailed modelling of the gasket or by defining a global non-linear behaviour in which material and geometry non-linearities are summarised. Two approaches were used to simulate the gasket behaviour. The first one consists in using a 2D approach, which allows to model through-thickness non-linear behaviour of gasket. The second one consists in using conventional 3D finite element modelling. The numerical methods have been discussed and compared in relation to the accordance with experimental data, amount of information supplied and computational time required. Finally, a parametric study shows how some geometric parameters influence the compressive load and the elastic recovery of a single-layer steel gasket.


Introduction
Automotive cylinder head gaskets must prevent leakage from the cylinder bore and the coolant and oil passages; in particular its main function is to concentrate the load coming from the cylinder head bolts on small areas around these openings in order to guarantee the sealing effect.
The most widespread gaskets in automotive applications are made of a number of metallic layers with different geometries and functions. In general, the sealing effect is obtained by the elastic action of beads, i.e., semicircular embossments, which are deformed by the flanges during clamping and act like springs against the mating flanges.
In a common design procedure, the overall clamping force is not adjustable: it comes from the layout of the joint, i.e., the position and dimensions of the bolts. The geometry of the gasket is then a key parameter to properly distribute the available load and to achieve the best sealing performance. Therefore the geometry of the bead must be accurately defined in relation to the clamping force. Moreover an optimum design of the gasket can minimise the contact pressure fluctuation and the amplitude of motion between cylinder head and engine block, due to the combustion pressure peak [3]. A variable cross section geometry of the gasket can achieve this objective, taking advantage of the recovery capability of the bead. Fatigue cracking can occur if the beads are required to compensate large displacements; in this case a stopper layer can act as a compression limiter [4].
Gaskets are not only present inside the engine, but in many applications, but the behaviour is more or less the same as in many other mechanical components gaskets do the same job, preventing leakages and contamination. The way to study them is almost the same, because the gaskets are not the components under investigation, but usually they are the interposed element to consider for the pressure distribution. In [5,6] the gaskets are modelled as a component in order to estimate their behaviour, by using finite element analysis (FEA). The goal was principally to evaluate pressure distribution.
The analysis of the force-displacement curve is a simple way to relate the compressive behaviour of the gasket to its geometry; such a curve can be obtained both by experimental tests and by FEA. Finite element codes supply two main methods for simulating the compressive behaviour of a gasket. When a forcedisplacement curve is available, the first one consists in simulating of the overall gasket behaviour by using 3D gasket elements [7][8][9][10]. In this case the geometric and material non-linearities are summarised in the so-called 'gasket behaviour' where a force-displacement relation describes the gasket behaviour. A simplified geometry consisting of a uniform thickness layer with an equivalent behaviour can be used [11]. Such approach is not suitable for a parametric study because the contribution of the geometry to the overall behaviour is not clearly evaluable. The second method consists in the detailed modelling of the gasket in which each layer is fully represented and discretised by 2D or 3D finite elements [12,13]. The elastic-plastic material behaviour does not depend on the geometry so a parametric analysis can be carried out.
Another important aspect to consider is the forming process of beads. In fact, to obtain a correct estimation of the clamping force, the contribution of the residual stresses to the overall behaviour of the gasket is not negligible and must be taken into account [14,15]. This solution is very attractive, especially during the design phase, because it can lead to a reduction in the development time and costs.
The present study shows the comparison of two modelling strategies and how some geometric parameters influence the compressive load and the elastic recovery of a single-layer steel gasket. Both approaches were used to simulate the gasket behaviour. In particular, the paper focuses on the forcedisplacement characteristic curve of a single-layer gasket for which compression tests provided the experimental force-displacement curve that is needed for validating the numerical models. The main simulation methods are discussed, and particular attention is given to the choice of the appropriate type of finite elements to model the effects of residual stresses. Difficulties arising from the two approaches are discussed, both considering the objective to achieve and the required computational time. A parametric investigation is also presented to show the influence of three geometric parameters upon the compressive load and the elastic recovery of the investigated gasket.

Experimental Tests
The selected gasket has an axisymmetric geometry consisting of a single layer with a single bead (Fig. 1). As the correct simulation of the forming process leads to a proper evaluation of the resultant shape of the bead, the dies and the reference gasket profiles were measured by a profilometer.
The gasket material is X5CrNi1810 steel (AISI 304) with yielding stress σ Y = 300 MPa, ultimate tensile stress UTS = 630 MPa, and elongation at break A = 50%. The geometry of the gasket is depicted in Fig. 2, where s is equal to 0.25 mm, w is 2.5 mm and the bead height h is 0.6 mm.
A Galdabini Sun 5 testing machine ( Fig. 3) (standard 5 kN load cell, accuracy of positioning repeatability 0.02 mm, internal stroke resolution 0.043 μm) was used for measuring the experimental curve force-displacement of the gasket; a 100 mm 2 punch compressed a part of the gasket according to the following procedure (graphical information in Fig. 3): Pre-load, a force of 73.5 N is kept constant to effectively establish contact (i.e., the zero displacement point); First compressive phase, mobile crosshead moves downward with a constant speed of 1 mm/min until a 0.6 mm displacement (corresponding to the nominal height of the bead) is reached; Second compressive phase, load increases with a constant speed of 50 N/s until a value of 5000 N (in this phase the displacement is very small because the gasket is almost fully compressed);   Holding phase, maximum load of 5000 N is hold for 40 s; Release phase, load decreases with a constant rate of 100 N/s.

Numerical Simulation
As previously stated, the actual non-linear behaviour of the gasket is due to two contributions: the material and the geometry. This investigation aims to outline the influence of some geometrical parameters on the overall compressive gasket behaviour, so the geometry must be defined in detail. Simulating methods are presented below dividing them into two groups: geometry-independent and geometry-dependent methods.
Geometry-independent methods do not directly depend on geometry, because its contribution to the overall gasket behaviour is not separated from the material contribution. The 'gasket' finite element (FE) of Abaqus code allows the definition of the through-thickness behaviour as a loading-unloading curve. The force-displacement relation can be obtained either from an experimental analysis or by a detailed simulation of the gasket cross section. This method has proved to be very strong with complex multilayer steel gaskets. In this case a single 3D gasket element layer can summarise the behaviour of several layers, each with different geometries and/or materials [11,18]. This allows great advantages when simulating large assemblies, in particular from a computational point of view. However the approach is not applicable during the first phase of design, when the objective consists in choosing a proper layout of the gasket, defining the number of layers and their geometry. A partial solution comes from the development of a simplified 1D model [16]. The non-linear behaviour of the cylinder head/gasket joint is simulated by springs. By changing few parameters the resultant load and deflection of each element changes, so an optimum gasket design can be quickly achieved. The approximation provides a near final design of the gasket that can be used for further 3D investigation.
Geometry-dependent methods directly depend on geometry: each gasket layer is fully represented and discretised by solid 2D or 3D elements. Material and geometry non-linearities have two distinct sources, so an optimisation process based on gasket layout variation is possible. The weakness of the method can be the computational time required. Several solid elements formulations were tested to evaluate their performances; material properties and boundary conditions were kept constant. Among these models, in the present paper, 2D and 3D modelling strategies are investigated.

Material
During compression, the material undergoes large plastic deformation. An isotropic elastic-plastic behaviour is considered, with yield surface defined by the law: (1) where σ eq is the equivalent stress, e pl 11 the plastic strain, and R the parameter that specifies the radius of the von Mises yield locus on the deviatoric plane.
The value of R is a function of the plastic strain e pl 11 and it can be obtained from the stress-strain curve of the material. As shown in Fig. 4, it results from the interpolation of two points: the yield and the ultimate material strengths.
According to the von Mises criterion, the equivalent stress σ eq is a function of the deviatoric stress tensor s:

Boundary Conditions
In the forming phase the gasket lies on a lower fixed die while a displacement condition is applied on the upper die; during the compression phase the gasket lies on a lower fixed flange while the punch moves and compresses it. The surface-to-surface formulation enforces contact between the master surfaces, i.e., the flanges, and the slave surface, i.e., the gasket.
A body weight force pushes the gasket against the flanges; the contact between them acts as a boundary condition for the gasket.
The tangential behaviour is governed by the basic Coulomb friction model while a hard pressureoverclosure relation is enforced as normal behaviour, i.e., the penetration of the slave surface into the master surface is minimised.
To facilitate convergence and control of rigid body motion before contact closure and friction are fully established an automatic stabilisation of the contact is applied, i.e., a small amount of viscous damping is activated at all slave nodes to prevent relative motion [17].

Simulation Procedure
The simulation consists in two steps: the forming process and the compression. These two steps are necessary in order to understand the influence of the forming process on the mechanical behaviour of the gasket. So all models were run with both the steps or just the compression one. In all the cases the gasket is the only deformable component, the dies and the clamping flanges are considered rigid. The sequence of operations is exactly the same of the actual procedure. To define in detail the gasket geometry, 2D and 3D solid finite elements were used; in particular four element types were considered, for details on the elements, refer to [17]: 4-node linear plane stress quadrilateral, reduced integration (4pσ); 4-node linear plane strain quadrilateral, reduced integration (4pε); 4-node linear axisymmetric quadrilateral, reduced integration (4ax); 8-node linear brick, reduced integration (8l).
As depicted in Fig. 5, the 2D model presents just the mesh of the gasket cross section and the dies are modelled as rigid lines. The mesh quality is necessarily high, with a dense discretization. In Fig. 6, the 3D model is shown. It has a coarser mesh because of the very high number of nodes needed to describe the entire gasket, but at least 3 elements are used to discretize the thickness of the gasket.
Forming and compression simulations involve complex, changing contact conditions, so a fine mesh of linear reduced integration elements is needed. Only linear elements were taken into account because they provide better results when complex contact conditions are involved [17]. The parametric study takes into account three characteristic dimensions of the gasket: bead width w, bead height h, and metal layer thickness s (Fig. 2).

Results and Discussion
Model validation was performed comparing the computed load-displacement curves with the experimental one. Simulations without considering the gasket forming process leads to an underestimation of the force-displacement curve (Fig. 7), while those obtained taken into account the forming process (Fig. 8) show a good correlation with the experimental data, in particular for what regards the maximum compressive force and the elastic recovery phase. The region around 350 μm of deformation highlights the difference between the models and the experimental test. By taking into account the forming process, the shape of the curves of the numerical result approaches the experimental one more than the results without forming process. This region is characterised by a sudden stiffness  Tab. 1 reports the internal energies associated with deformation, i.e., the integrals of each curve. Regarding the simulations without accounting the gasket forming process, the discrepancy between numerical and experimental values is evident: the average internal energy of the three simulations spans between −47% and −35% of the experimental value. Also the stiffest finite element (FE), the axisymmetric one, provides a very rough approximation (−35%) of the experimental data. Better results are obtained when the forming process is taken into account. All the numerical curves shift upward (Fig. 8) and the deformation energies increase: the average internal energy of the three simulations (4pσ, 4ax, and 4pε) spans now between −20% and −7% of the experimental value; the gap reaches the minimum with the axisymmetric FE (−7%). The 3D FE (8l) overestimates the internal energy by +12% but outputs the best fitting curve (Fig. 8, black curve); this is due to the fact that the simulation accurately reproduces the experimental conditions.  The 2D elements allow a fine discretisation of the cross-section, i.e., a good evaluation of the stresses caused by the forming process, with a very low computational time. The 3D FE leads to a dramatic increase of the computational time, even if the number of elements in the cross section is much lower. Moreover, the improvement in the force-displacement curve does not match an equivalent improvement in the internal energy. Considering both the deformation energy and the maximum compressive force, the axisymmetric FE model represents the best compromise; it also allows a good computational efficiency.
Another aspect to evaluate is the shape of the gasket at the end of the forming process. Numerical and experimental profiles, for the axisymmetric FE model, are almost coincident, as shown in Fig. 10. This correlation means that the number of elements in the gasket width is sufficient to properly describe the formed shape. The axisymmetric FE model has then proven to be a good method to obtain reliable results because it properly describes both the forming and the clamping phases. In Fig. 9 the maximum in-plane principal and the Von Mises stresses during the forming and squeezing process are recorded. It can be noted that the residual stresses are important, around 250 MPa.
The low computational time of axisymmetric FE model makes it attractive also for a parametric study of the gasket. Starting from the reference gasket of Fig. 1, the inner and outer diameters were kept constant while the height h and the width w of the bead varied from the reference values. To obtain different beads, the dimensions of the dies were modified, as shown in Fig. 11. This map was repeated for three thickness values: 80%, 100% (reference value) and 120%.   For each thickness s, results are plotted as a percentage of the maximum value. Fig. 12 shows that the clamping force decreases as the bead width w increases and as the bead height h decreases. This trend is monotonic except for the 80% thick layer, where a local maximum is shown for w = 110% and h = 67%.  As the height h decreases, the influence of the width w becomes less relevant, i.e., the difference between the maximum and the lower values of w (at that height) decreases. The only exception occurs for s = 80%, where the influence of width w is maximum for intermediate values of height h. Intermediate values of w show the maximum influence of height h for thicknesses of 80% and 100%. Moreover this maximum is emphasised for a thickness of 100%. Regarding the 120% thickness gaskets, the maximum influence of height h is shown for the minimum width w. Finally, the difference between maximum and minimum value is almost the same for each thickness s (about 61%).
Recovery is defined as the increase in gasket thickness after the complete removal of the clamping force. In practical applications, however, it is relevant to understand the decrease rate of the clamping force, i.e., the displacement at which the force reaches a specific threshold. The latter is function of the characteristics of the joint, i.e., the geometry of the flanges, the pressure of the confined fluid, the operative temperature, and represents the minimum value that guarantees the functionality of the joint. In this paper that value was set to 2000 N; therefore the recovery is the upward displacement of the flanges that makes the force drop to 2000 N. For each thickness s, in Fig. 13 the results are normalised with respect to the maximum value. The greatest influence of width w is shown for s = 80%; for each value of height h the difference between maximum and minimum values is not lower than 23%. Moreover the surface monotonically decreases when the width w increases. For other thicknesses s the trend is not so clear and depends on the value of height h. Considering the influence of the height h, the monotonically decreasing trend of Fig. 13 gradually vanishes in Fig. 13 (middle) and Fig. 12 (bottom) and a maximum appears for h = 67%. The difference between maximum and minimum values decreases from 72% to 55% as the thickness s increases from 80% to 120%. Fig. 14 shows the influence of thickness s on clamping force and elastic recovery: values are normalised with respect to the absolute maximum. It is clear that the clamping force increases with the thickness; moreover this increase is almost uniform: the increase amounts to 80-100% and 100-120% (for specified w and h) and does not differ more than 5% of the maximum load. The only exception occurs for w = 100%, h = 34% where a small increase (4%) from 80% to 100% is followed by a great increase (17%) from 100% to 120%. As the thickness increases, the recovery increases. This increase, however, is not uniform: depending on the values h and w, it is greater between 80% and 100% or between 100% and 120%. For example, if w = 90%, h = 67% an increase of 36% between s = 80% and s = 100% is followed by an increase of 2% between s = 100% and s = 120%. Therefore the thickness increase does not cause any relevant improvements of the recovery. On the other hand, if w = 110%, h = 34% the total increase is evenly distributed between the two steps. Figure 14: Influence of thickness s on closing force (upper) and elastic recovery (lower)

Conclusions
The paper focuses on a single beaded metal gasket. Experimental tests were carried out on a reference gasket in order to obtain its compressive behaviour (compression test) and its profile. Since the forming process deeply influences the final behaviour of the gasket, information about the geometry of the dies was acquired by profilometer. The available numerical approaches are distinguished into two classes, according to their relation with the actual geometry of the gasket.
The geometry-dependent methods were investigated; in particular four element types were evaluated in relation to the compliance with experimental data and to the computational time.
The axisymmetric FE proved to be the best compromise between accuracy and computational time. Therefore this type of FE was chosen for the parametric study aimed to outline the influence of three geometric parameters (gasket height h, width w and thickness s) upon the clamping force and the recovery.

Several conclusions arise:
Stress field due to the forming process produces a non negligible increment in the compressive force; 2D elements are computationally efficient and allow a fine discretisation of the cross section; Axisymmetric elements represent a good compromise between computational time and compliance with experimental data; 3D elements accurately reproduce experimental conditions but they lead to a dramatic increase in the computational time; moreover, such an increase does not match a comparable improvement of accuracy; Clamping force increases if the height h of the bead increases and decreases if the width w increases; Clamping force increases almost evenly with thickness s; As thickness increases, the recovery is maximum for intermediate values of height h; Recovery increases unevenly with thickness s; the degree of recovery improvement due to a thickness increase strictly depends on the other two dimensions w and h.
In conclusion, in the selected gasket dimension range, a maximum clamping force does not match a maximum recovery. The design phase should then accurately evaluate which of these properties is essential to the functionality of the joint. Finally, a reliable method to obtain the force-displacement curve represents a solid starting point for the early definition of a complex multi-layer-metal gasket layout.
Funding Statement: The authors received no specific funding for this study.
Conflicts of Interest: The authors declare that they have no conflicts of interest to report regarding the present study.