Issue 
Mechanics & Industry
Volume 25, 2024



Article Number  22  
Number of page(s)  20  
DOI  https://doi.org/10.1051/meca/2024018  
Published online  29 July 2024 
Original Article
Convergence analysis and mesh optimization of finite element analysis related to helical springs
ICA, Université de Toulouse, UPS, INSA, ISAESUPAERO, MINESALBI, CNRS, 3 rue Caroline Aigle, 31400 Toulouse, France
^{*} email: cadet@insatoulouse.fr
Received:
24
November
2023
Accepted:
9
June
2024
Helical springs are widely used in engineering applications. In order to reduce cost in “try and error” time consuming experimental campaigns, numerical simulations became an essential tool for engineers. Indeed, it saves considerable time in the ahead design phase of a project to ensure the feasibility of structures. However, these simulations run thanks to a lot of parameters, which all must be selected carefully to get access to reliable results. In this paper, ten main modeling parameters are presented. Thanks to a valuable literature statistical analysis, seven of them are settled. Three remain to be studied: the mesh density, the order of the elements and the integration method. Then, three convergence analyses are performed with ABAQUS about the circular geometry accuracy of the tessellated surface, the axial stiffness (and axial load) accuracy of the helical spring and the maximal Von Mises stress accuracy within the helical spring. The numerical campaign is led with 8 mesh densities along the circumference and 6 element types. After comparison, in order to get both fast and accurate results, a limited list of nearoptimal combination of density and element type are proposed. The users are free to use any of the presented solutions in function of the desired admissible accuracy of their model.The proposed meshing technique can be exploited for any helical structure with circular cross section, mainly loaded in torsion and shear, such as extension and compression springs.
Key words: Mesh sensitivity study / FEA uncertainty / mesh density / stress / stiffness / geometry
© G. Cadet and M. Paredes, Published by EDP Sciences 2024
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Helical springs are mechanical components, commonly exploited to store and restore mechanical energy. They can be made in various materials and sizes. Their external shape may be cylindrical, conical or more complex, and they can be used in compression, extension or torsion to solve problems in many applications [1]. The most frequently used one is certainly the cylindrical helical spring with constant pitch and circular wire, which is known for its linear forcelength relation [2,3]. Analytical approach for spring calculations has been already fully developed in the 1960's, mostly by A.M. Wahl [4]. However, further effort to bring about new calculation equations have been performed to increase the accuracy of these formulations, often base on FEM [1,5–15]. However, they are often limited due to the complex helical shape of the spring and the influence of the boundary conditions. Therefore, it seems that only numerical analysis is left to improve the accuracy of these equations. Spring calculations focus on maximal stress, axial load/length behavior, thermal response, natural frequency or mass/volume properties. Thus, FEA about springs are generally meant to calculate one or several of these five outputs. The FE model can either be built with 1Dbeam or with 3Dcells, both having their own strengths and weaknesses. The 1Dbeam model is very fast to compute due to the simplicity and the low density of the elements. It is very accurate to compute frequencies or stiffness [9,16–21]. However, it is complex to model contact or to visualize surface stress with this model. This is why 3Dcells models are more used in literature [22–71]. This type of model enables to obtain more reliable results with contacts between the wire and the support surfaces and between coils. Maximal stress value can also be easily computed with it.
The objective of this paper is to help designers and researchers model their springs. In order to achieve this, an important literature review is performed with a statistical analysis of a 50papers corpus studying FEA related to helical springs. Then, a 48models campaign is led with investigations about their stiffness, stress and geometry accuracy. Finally, an optimization process is realized and nearoptimal solutions are presented.
2 Literature review
Finite element analysis (FEA), also called the finite element method (FEM), is a method for numerical solution of field problems [72]. In mechanical engineering, we may seek the distribution of temperature, displacement, stress, pressure or load for example. Individual finite elements (FE) mesh the studied region in a finite number of small pieces, which are represented by a system of algebraic equations to be solved for unknowns at their nodes [72]. These sets of equations depend on the mesh geometry, size and resolution method. Thus, the structure field quantity is only approximated and accumulates the approximations, but can be optimized by selecting appropriate parameters. Unfortunately, this accuracy has often a cost: the central processing unit (CPU) time, also called solver time. In engineering and research fields, a compromise needs to be set between accuracy and calculation time. While it is possible to use FEA programs while having little knowledge of the analysis and procedure method used by the software, it is thus recommended to learn how elements work in order to perform reliable and fast calculus in specific applications [72]. The geometry, the material and the boundary conditions are often dependent in these applications. However, it is still possible to identify 10 major choices that can induce a significant impact on the result accuracy and calculation time:
The FEA Software (ANSYS/ABAQUS/...).
The dimension of the FE model space (2D/3D).
The procedure analysis (implicit or explicit).
The dimension of the elements (1D/2D/3D).
The shape of the elements (number of faces).
The quality of the elements (aspect ratio, skewness).
The density of the elements (size, uniformity and distribution).
The meshing technique (sweep/free/structured).
The order of the elements (linear/quadratic).
The integration method (hybrid/reduced/...).
In order to produce a reliable replicable scientific work, it is mandatory to specify each of these 10 selected options (or at least for some of them it must be visible on figures). The replicability of a scientific work is defined as the possibility to obtain similar results with a given accuracy by another research team, in a different laboratory, but with the same method [73–75]. Modern science relies on thirdparty verification of results. Any method must therefore be replicable. This allows future researchers to accurately examine the results and check the accuracy and the performance of the presented tool. In addition to these 10 items, an eleventh parameter can also be given if comparison with other works are performed in the paper: the number of CPU used for the calculation. An objective FEA comparison must be performed on similar computer power. In the worst cases, FEA are performed without mentioning or showing any of the 10 items above, or just giving the used FEA software [76–78]. It then becomes impossible to replicate these results.
Several FE commercial software are used in mechanical research and engineering with, among others, ABAQUS, ANSYS, NATRAN and HYPERMESH. Each of them has its own specifities, abilities and strengths. However, ANSYS and ABAQUS software are clearly ubiquitous in scientific literature. Both software seems to have similar performance capabilities [79]. M. Shimoseki et al. [9] states that, for linear analysis, commercial software are equivalent in quality, but a difference can be shown in their ability to treat nonlinear analysis. For the authors, ABAQUS is the best FEM software to deal with geometrical and boundary nonlinearities, compared with ANSYS.
Once the FEA software carefully selected, it is time to create the model. Because helical springs are helical structures, no symmetry simplifications are possible and the model must be threedimensional. The dimension of the solid, on the other hand, can either be made from 1Dbeams or 3Dsolid elements. The 1D element has the advantage to be solved really fast, but is not convenient to model complex contacts between coils and visualize localized stress along the wire surface. Indeed, G. Cadet et al. [5] highlighted that surface stress fields of helical springs admit significant gradient within the wire cross section due to curvature. 1D elements are not capable to predict such results. Threedimensional solids are thus more powerful, even if more time consuming. Solid (continuum) elements are the standard structural volume elements of ABAQUS. The 3DSolid elements in Abaqus can be used for complex linear and nonlinear analyses involving contact, plasticity, and large deformations. They enable to calculate stress and stiffness accurately [80].
Then, the procedure analysis of the software must be selected. The first step in solving a problem is to identify it [72]. Is the problem timeindependent or −dependent (i.e. is the problem static or dynamic)? Is nonlinearity involved so that iterative solution is necessary? What accuracy and what results are required? Answering to such questions influences how the problem is modeled and which procedure is adopted. Whereas implicit analysis (Standard in ABAQUS) must iterate to determine the solution to a nonlinear problem, explicit analysis determines the solution without iterating by explicitly advancing the kinematic state from the end of the previous increment [80]. Spring analysis can induce contacts between coils and “smooth” modification of geometry when loaded. Implicit analysis is dedicated to weak nonlinear problems, which makes it suitable for general spring FEA, even if the explicit procedure would be as efficient in some cases.
In general FEA commercial software, four solids are available for third dimensional elements: 4 faces tetrahedra (triangular based pyramids), 5 faces triangular prisms (wedges) and squarebased pyramids, and 6 faces hexahedra (bricks) [80]. Tetrahedral elements are used in many automatic meshing algorithms because they are geometrically versatile. It is very convenient to mesh a complex shape with tetrahedra. Indeed, they are almost nosensitive to the initial shape of the part. However, hexahedral elements usually provide more accurate results and have a better convergence rate than tetrahedra. In fact, tetrahedra are usually overly stiff [81] due to their triangular faces and extremely fine meshes are required to obtain accurate results. Wedges and squarebased pyramids also admit error due to the same problem. Only hexahedral elements are proof against this overstiffness thanks to their flexible rectangular faces. ABAQUS documentation [80] states that “these elements − tetrahedra and wedges − should not be used or avoided as much as possible, except as filler elements in noncritical areas [9,80]. If they are required, an extremely fine mesh may be needed to obtain results of sufficient accuracy. Therefore, try to use wellshaped elements −hexahedra in regions of interest”.
Then, we must choose between first and secondorder elements. In firstorder (linear) hexahedral solid elements, the strain operator provides constant volumetric strain throughout the element [80]. This constant strain prevents mesh “locking” in certain cases. Linear elements (C3D8) ensure faster results but can also induce major inaccuracy. Secondorder (quadratic) elements (C3D20) provide higher accuracy in ABAQUS/Standard than firstorder elements, in particular for problems that do not involve severe element distortions (like crimping, shock or explosion) [80]. They capture stress concentrations more efficiently and are better for modeling geometric curved features. Indeed, they can model a curved surface with fewer elements. Then, for helical springs, quadratic elements can be very effective. Thirdorder elements also exist in ANSYS, but are not available in ABAQUS, NASTRAN nor HYPERMESH. Because this paper aims to propose a general solution achievable on every simulation software, thirdorder elements are not studied.
In addition, reducedintegration, hybrid, improved surface stress and incompatible mode options are available in ABAQUS/Standard. Selecting the appropriate element type for a given simulation is critical for obtaining accurate results. For instance, using fullyintegrated hexahedral elements in bending results in shear locking (which artificially stiffens the model) [82]; similarly, using incompatible mode elements in compression can often lead to convergence difficulties. Let's detail each of these methods. Reduced integration uses a lowerorder integration to form the element stiffness. The mass matrix and distributed loading use full integration. Reduced integration reduces running time by removing integration points (27 to 8 for quadratic elements and 8 to 1 for linear elements) [80]. Therefore, reduced linear elements accuracy is largely dependent on the nature of the problem. The elements with reduced integration (C3D8R and C3D20R) are also referred to as uniform strain or centroid strain elements with hourglass control [83]. Hourglassing can be a problem with uniqueintegration point elements (firstorder, reducedintegration) in stress/displacement analyses [84]. It is possible for the element to distort in such a way that the strains calculated at the integration point are all zeros, which leads to uncontrolled distortion of the mesh. Hourglass control is made to avoid such problem. Fully integrated elements do not hourglass but may suffer from “locking” behavior: both shear and volumetric locking. Shear locking occurs in firstorder, fully integrated elements (C3D8) that are subjected to bending. The numerical formulation of the elements gives rise to parasitic shear strains. Therefore, these elements are too stiff in bending [82]. This parasitic shear stress can be avoided with the use of incompatible mode elements (C3D8I). These elements adds internally incompatible deformation modes to eliminate the parasitic shear stresses, which make them more expensive than the regular firstorder elements, but more economical than secondorder elements [80]. Volumetric locking occurs in fully integrated elements when the material behavior is nearly incompressible (ν ≈ 0.5). Additional pressure stress appears at the integration points, causing an element to behave too stiffly for deformations that should cause no volume changes. If nearincompressible behavior occurs, a very small change in displacement produces extremely large changes in pressure. Therefore, a purely displacementbased solution is too sensitive to be useful numerically. Hybrid formulation removes this singular behaviour from the system by treating the pressure stress independently. Hybrid elements (C3D8H and C3D20H) have more internal variables than their nonhybrid counterparts and are slightly more expensive [80]. Finally, the C3D8S elements have been developed to provide a superior stress visualization on the element surface by avoiding errors due to the extrapolation of stress components from the integration points to the nodes. These elements have the same degrees of freedom and use the same element linear interpolation as C3D8 but use a 27point integration scheme [80]. The presented configurations are summarized in Table 1.
Another mesh geometric parameter to be attentive to is the aspect ratio. Aspect ratio is a mesh quality metric [85,86] that measures the deviation of a mesh element. The aspect ratio of a hexahedral element is the ratio of its longest edge with its shortest one [86]. The ideal aspect ratio is 1 and it is its minimal value. A hexahedral element with aspect ratio of 1 is a perfect cube. Higher this ratio is, lower is the quality of the element. Large aspect ratios lead to poorly conditioned matrices, worsening the speed of convergence and accuracy of the linear solver [87]. The objective is thus to minimize the aspect ratio. One solution is to set for every edge identical mesh size (i.e. equidistributed mesh seeds). Another quality metric is the skewness of the element. The skewness is defined either as (1) the deviation between the optimal cell volume (of a cube for hexahedral elements) and the real cell volume, or (2) the deviation between the optimal faces angle (90° for hexahedral elements) and the real faces angle. In both cases, the skewness measures the distortion of the mesh element and how its faces angles are sharp. In general, sharpest angles deteriorate the accuracy of the results. Solid elements are more accurate if not distorted [80].
The next step after the selection of element and the determination of analyzing content is the discretization of domain into elements [9]. As previously mentioned, this is one of the most important step which have simultaneously a significant impact both in accuracy and solving time. Large elements (due to low density) will generally incorrectly predict stress and loads for a given displacement due to their constant strain and the tessellation. As previously mentioned, linear elements possess a constant strain operator, which provides quasinull stress gradient. If the stress gradient of the spring is sufficiently small within a broad mesh domain, a good approximation may be obtained by linear elements. However, when the stress gradient is relatively high compared to the mesh cell size, the supposed to be constant stress in the domain does not hold and the operator give poor results [9]. Tessellation is a mathematical concept consisting in representing a continuous complex surface with discrete simple shapes. Indeed, circles, as wires cross sections, are impossible to reproduce with perfect accuracy with linear elements. Thus, this tessellation degradation can induce multiple geometric errors impacting the stiffness but also contacts behaviour. Indeed, the gap between two coils is therefore increased if the wire is no longer a cylinder [88]. The compact length of the spring, for example, is impacted. In contrast, small elements (due to high density) will better approximate the geometric complex shapes. However, it is not advisable to divide the domain indiscriminately fine. Such dividing increases the number of elements, thus the number of unknowns to be solved, resulting in a lengthy and costly analysis. If the density is doubled, because the problem is 3dimensional and the solving time is mostly correlated to the number of unknowns, the solving time is approximately multiplied by 8 (2^{3}). No general standard exists regarding the ideal fineness of discretization. Convergence analysis or mesh optimization must be performed for each application [89]. For springs FEA, several convergence analyses were already performed but only for 1Dbeam helical elements [9,16–21]. Wire analyses are often realised for frequencies investigations because they are sufficiently accurate and very fast to compute. They all agree that a 50 nodes by coils density is optimal to obtain accurate results. For 3D elements, spring FEA convergences are not as much studied and only performed to solve particular engineering problems. Unfortunately, no general solution can be extracted from these analyses due to multiple reasoning errors [17,27,48,81]. For example, Kelly and Knight [17] compare only two mesh densities, which give a poorly admissible convergence proof. Mulla et al. [27] compare 6 models (2 element types × 3 mesh densities) with tetrahedral elements. Unfortunately, no convergence is proven and the authors take the model with the closest result too the theoretical value rather than the converged value. Mirzaeifar et al. [48] state they found an optimal 14 nodes in circumference mesh density with C3D20R elements but did not reveal how they obtain such results. Finally, Yazar [81] proves that hexahedral elements are faster and give similar results than tetrahedral elements and also that fine tetrahedra mesh gives divergent results due to overstiffening, but decides to chose anyway tetrahedra to mesh his spring.
An associated mesh quality metric [85,86] to mesh density is the smoothness, also known as the volume ratio. The smoothness is the ratio of the volumes of two adjacent cells. Any sudden jump in mesh size can result in an error. A good smoothness, with uniform seed distribution, is highly desirable in simulations, as it enables an accurate interpretation of the model.
The last parameter to select is the meshing technique, the structure of the elements. Since a spring is a helically swept circular (or elliptical) cross section, the best consistent meshing technique along the wire is by sweeping the cross section pattern. The cross section can be elliptical if the cross section is not perfectly orthogonal with its sweeping path. Therefore, we must choose carefully the pattern of the cross section, which will impact the model accuracy. As previously discussed, tetrahedral elements are nonoptimal to mesh solids in mechanical engineering. Thus tetrahedralmade cross section, free (Fig. 1a) or structured (Fig. 1c) will not be discussed here. For hexahedra meshing, three main options are available. The first one is to use free meshing, using no preestablished mesh pattern (Fig. 1b). This is convenient for complex geometric parts but can be not really convenient for swept profiles. Indeed, with free pattern, no wire section along the spring is meshed identically, which can be a problem when extracting the results. In addition, this technique does not ensure uniformity and absence of distorted elements within the part. Thus, structured patterns are more suitable for swept profiles. For circular cross sections, two structured patterns predominate. The first one is a structured “pie pattern”, also called “midpoint generated pattern” (Fig. 1d). This pattern partitions the section into n “pie slices", with n the number of nodes along the circumference of the cross section. This fast mesh generation has four major disadvantages: (1) the more a mesh is close to the centroid, the worst its aspect ratio and skewness are, (2) there is a significant disequilibrium in mesh size between external and internal mesh elements, (3) it forces the meshing algorithm to use triangular prism elements close to the centroid and (4) the central node is connected to a very important amount of elements (2 × n). Such a phenomenon is called a singularity. A singularity refers to an irregular node where more (or less) than the regular number of elements meet, 8 in the interior of a hexahedral grid [90–92]. Singularities are very critical in structured meshing because they induce a stiffening of the affected node and stress calculation will be nonconvergent. The more elements are connected, the sharper are the corners angles (the higher is the skewness), and the higher will be the calculated stress. Indeed, stresses in sharpest corners tend to infinity. Thus, it reveals that the best structured pattern to fill a circular shape with hexahedra is the “Ogrid” partitioning, also called “Butterfly grid” [90] (shown in Fig. 1e). The objective of this second structured mesh generation is to split the circular cross section into five subdomains: a squarelike shape in the middle of the section and four circular arcs around it. The main advantage of such submodelling technique is that each subdomain has a topological structure with rectangular features; hence, each block can be easily meshed with hexahedra [90]. To improve the grid quality, the central square can have slightly curved edges. Furthermore, this pattern offers a very low aspect ratio and avoids the presence of singularities mentioned previously [90–92]. At most, a node of this structure is linked to 8 element cells.
Following this important literature review, some parameters of the FEA of helical springs can already been selected. However, some of them remain unclear. Thus, the goal of this study is to propose a general optimal mesh for helical spring applications.
Comparison of the labels and configurations procedures of mesh geometric order and integration method for element in ABAQUS and ANSYS software packages.
Fig. 1 Examples of mesh patterns in circular cross section with 16 nodes along the circumference: free tetrahedra (a), free hexahedra (b), structured tetrahedra (c), structured “pie” hexahedra (d) and structured “curved Ogrid” hexahedra (e). 
3 Method
3.1 Corpus analysis
In order to have a global review of the scientific research already performed about FEA of helical springs, 50 papers from 2000 to 2023 were gathered into a corpus [22–71]. Statistics about this corpus can help understanding what is commonly used and if general solutions or errors can be extracted from it. The data from all the 50 papers of the corpus are available in Appendix A.
As mentioned previously, springs FEA are performed to compute more accurate results than standard formulation give access to. The corpus analysis shows that 84% of the papers perform FEA on springs to obtain the maximal stress value, 84% also to obtain the stiffness (or any loaddeflection relations), 16% to access the natural frequency of the spring, 18% to thermal responses and 10% on geometrical and mass properties. Thermal analysis is often conducted on similar specific applications, like shapememory alloys (SMA) springs, composite springs or fatigue creep. And frequency analyses are more accurately computed with faster 1Dbeam FEA. Therefore, thermal and frequency analyses are not further examined in this paper. However, stress, stiffness, and geometry are. To obtain these results, there are only two options to load the spring: by displacement or by exercising an axial force. Even if the second option is chosen 74% of the time (and displacement 26%), the closest boundary condition to an experimental setup is by displacement. Indeed, standards [2,3] advice to drive the spring by a certain deflection and then measure the associated force. Displacement seems more reliable, this is why it is selected in this paper, but further results in Section 4.4 will help settle this decision.
About the 10 proposed options, mandatory to present a replicable scientific work, only 6% of the papers mention or show all the required data [29,37,64]. Actually, among the corpus, 100% show the element shape and 98% state the used software, but only 48% declare the order of the element, 46% show the mesh structure in the cross section, 30% inform the used integration method (element type), 22% mention explicitly the mesh density in the manuscript and finally only 16% specify the selected procedure analysis. This is the proof that researchers often neglect the FEA options in their manuscript, but certainly also in the modeling process. This is why such work as performed here is significant for the scientific community.
Then, we can compare which software they use. 58% of researchers use ANSYS, 40% of them use ABAQUS and 2% use NASTRAN. Even if ANSYS is more used, this data alone is not necessarily relevant. Indeed, only ABAQUS works mention the selected procedure analysis and the element type. No found ANSYS works specify them. This is mainly due to the specificity which ABAQUS states the element type directly in the selected mesh name, where integration method in ANSYS only appears in options (Tab. 1). Also, 62% of ANSYS works use tetrahedral elements where only 25% of ABAQUS works do. Finally, 100% of very low mesh density n (number of nodes by circumference), with n < 8, are performed with ANSYS. Therefore, ABAQUS use seems lead to better habits, more relevant works and so more reliable results. Also, because our application relates to “smooth” nonlinear geometry, material and boundary conditions, the choice of the software package seems insignificant. However, we choose ABAQUS for its richer mesh element library. A comparison with ANSYS elements can be found in Table 1, in order to transfer these results to ANSYS software.
Because only 8 of the 50 papers specify the selected procedure analysis, the extracted ratio of Implicit/Explicit selected procedure is no more reliable for a statistical analysis. However, it is possible to compare the number of simulations of the corpus which use tetrahedral and hexahedral elements. Almost half of the works (44%) use tetrahedral elements everywhere in the helical geometry of the spring, often with a very low mesh density. These two poor meshing choices lead to not reliable results, as previously mentioned [9,80,81]. Thus, hexahedral elements are used 56% of the time, which proves that this is not common knowledge and a particular attention must be paid to element shape definition. But elements are not only defined by their shape, but also by their order and their integration method. As for procedure analysis, papers stating the used element type are not enough to conduct a reliable statistical analysis. However, the paper mentioning the elements order are. 48% of these works use quadratic element and 52% linear ones. No general solution can be thus extracted from this data. Both must be examined in the campaign.
Finally, among the corpus, the number of elements in the circumference n is shown in Figure 2. Because this mesh density has a significant impact on the accuracy of the model, this data should not be neglected when researchers write their manuscript. Quite low mesh densities are used, with a statistical median at n = 12 nodes. We can also remark that multipleof4 mesh densities are more frequently used due to their symmetry efficiency. Indeed, in order to mesh the wire section in a structured way, the most efficient technique is to perform orthogonal partitions, slicing the circular section into four equal disk quarters. But this requires a supplementary modeling step, which 45% of the papers in the corpus using hexahedral elements and showing the section did not invest (i.e. they use free pattern). 45% use the optimal “Ogrid” structured pattern and 10% use the less efficient “pie” pattern.
By comparing the common literature review about the judicious choices in FEA and the actual FEA works realised about helical springs, we reveal that multiple errors are regularly done by researchers. No significant convergence analysis is already published in the literature, which certainly leads to these several deviations in the performed works. Therefore, there is an important scientific need to do a proper mesh optimization for helical springs. An exhaustive numerical campaign must be conducted in this way.
Fig. 2 Recurrence of each mesh density among the corpus [22–71]. The sum of all the bars adds up to 54 because some papers use several mesh densities. 
3.2 Numerical campaign
The objective of this section is to establish the convergence analysis of 3D solid FE helical springs and optimize the mesh density and integration method in order to obtain the best compromise between fast and accurate results. Investigations about geometry, stress and stiffness are performed. Frequency analysis and thermal response are not studied here. However, because frequency depends mainly on mass (i.e. volume, geometry) and stiffness, if the present paper succeeds in finding parameters to obtain accurate results about both geometry and stiffness, it ensures satisfying results about natural frequency of the spring. In order to do that, a numerical campaign is led with the following fixed parameters:
The FEA Software is ABAQUS.
The FE model is threedimensional.
The procedure analysis is implicit (ABAQUS/Standard).
The elements are 3D solid continuum.
The mesh density is uniform all over the solid.
The elements are hexahedra (“bricks”) with low aspect ratio, skewness and smoothness.
The meshing technique is “sweep” along the wire and “Ogrid” structured in the cross section.
The material is purely elastic, isotropic, compressible and homogeneous.
All of these parameters are set thanks to the previously presented literature review. The remaining parameters are tested and compared during the campaign. The order of the elements can be linear (C3D8) or quadratic (C3D20). The integration method can be “fully integration”, “reduced integration” (R), “incompatible modes” (I) or “improved surface stress” (S). Therefore, the examined elements are: C3D8, C3D8R, C3D8I, C3D8S, C3D20 and C3D20R (with default options). Because the material is compressible, hybrid elements (H) offer insufficient advantage. Several models with hybrids elements were run. Results are identical than their nonhybrid counterparts and are slightly more time expensive. Hence, hybrid elements are not examined furthermore. The mesh density is also studied. The number of nodes n along the circumference of the wire can evolve with n ∈ [8, 12, 16, 24, 36, 50, 75, 100]. The associated cross sections are shown in Figure 3. The 8nodes meshed section is the worst density ABAQUS can generate for swept helical geometry. The only meshed section which is not 4symmetrical is the 12nodes section. ABAQUS enables to mesh the 12 nodes circumference 4symmetrically, even with judicious seeds. In order to check if this specificity could induce an orthotropic response, an analysis is performed with a straight beam, clamped at one extremity, and loaded in shear at the other one and meshed with C3D8 elements. Two models are made, one with the 6 mesh diameter horizontal and the other vertical. As a result, a divergence of almost 3% in the maximal stress value is observed and a significant difference in the stress distribution (Fig. 4). Even if this phenomenon is identified, this mesh density is kept in the campaign for three reasons. First, it is still a judicious median low meshing density between 8 and 16 nodes circumference, where the results gradients are significant. Then, this mesh density is the third more used in the corpus [22–71], thus it is important to also examine it. And finally, even subjected to this phenomenon, perhaps this mesh density can bring satisfactory results.
The reference model is made of one singular helical coil of mean diameter D = 40.0, wire diameter d = 10.0 mm and pitch m = 22.0 mm. Therefore, its index (D/d ratio) is the smallest one admissible by standards [2,3]. This choice is due too the need of significant stress gradient into the wire, available in springs with low index. Because the solution of this study needs to be general for all springs, it is important that the proposed mesh type, order and density are reliable for extreme cases. With this tested geometry, we guarantee that the results are applicable and pertinent for any geometry of springs. For large index springs, a degradation of the aspect ratio could be tolerated to obtain faster results without damaging the accuracy of the results. Also, in order to perform the convergence analysis, there is no need of model the entire spring. Only the geometry and the boundary conditions need to be similar. Simulating only one coil will also save time during the campaign. The material used has the following properties: Young's modulus E = 180.000MPa and Poisson's ratio v = 0.3. In order to simulate the deflection as in an entire spring, two reference points (RP) are created on the central coiling axis of the spring, at the same heights as the centroid of both cross section extremities. Then, these reference points are coupled with the respective cross sections. The lowest RP is fixed in translation and wire torsion momentum. The upper RP has similar boundary conditions but is driven in (compression) translation in the axial direction and does not have rotation clamping. The axial deflection is u = 4 mm. This model is shown in Figure 5. Because it is a unique coil and m − d > u, no contact is examined here. In experimental and numerical studies, the spring analysis is conducted between two operating lengths. These two lengths limit the range the spring is used for in the system it is integrated in. The lowest length (i.e. the highest deflection) is commonly set that no contact can occur between the coils within the operating range of the spring [2,3].
Fig. 3 Meshed sections of the realised campaign with 8nodes circumference (a), 12nodes (b), 16nodes (c), 24nodes (d), 36nodes (e), 50nodes (f), 75nodes (g) and 100nodes (h). 
Fig. 4 Von Mises stress distribution in the section of the 12nodes circumference sections, with the 6mesh diameter horizontal (a) and vertical (b). 
4 Results and discussion
Now, the CPU solving time, the surface area of the cross section, the axial force induced by the deflection (stiffness) and the maximal resulting Von Mises stress of the models are investigated. The results of all 48 simulations are available in Appendix B.
4.1 Geometry accuracy
One of the fastest way to compare linear and quadratic meshes is to compare their geometric accuracy. Indeed, the meshing procedure degrades the surface accuracy for curved elements and any circle becomes a polygon (a regular polygon if we consider a structured equidistributed meshing technique). The resulting polygon is inscribed in the original exact circle and the following parameters can be introduced: R the radius of the original circumcircle, n the number of sides (nodes) of the polygon, a the apothem (radius of the inscribed circle) and s the length of the sides. These parameters are shown in Figure 6. The area ${A}_{n}^{1}$ of a regular polygon with n sides inscribed in a circle of radius R can be written as:
$${A}_{n}^{1}=\frac{nsa}{2}=n{R}^{2}cos\frac{\pi}{n}sin\frac{\pi}{n}.$$(1)
A quadratic meshed circular cross section would have the initial linear npolygon surface area, in addition with the quadratic sectors from each side. These n sectors are the local integrations of the quadratic function y = c_{1}x^{2 }+ c_{2} connecting the 3 following nodes: N_{1} and N_{2} the two original nodes of the regular polygon and N_{3} the new median node, located on the associated perfect circle arc circumference. These principles are shown in Figure 6. Geometric resolution of the equations gives:
$${c}_{1}=\frac{4\left(aR\right)}{{s}^{2}}{c}_{2}=Ra.$$(2)
Thus, the new area ${A}_{n}^{2}$ of this quadratic “polygonal” shape becomes:
$${A}_{n}^{2}={A}_{n}^{1}+n{\displaystyle \underset{s/2}{\overset{s/2}{{\displaystyle \int}}}}y\left(x\right)dx={A}_{n}^{1}+\frac{2ns\left(Ra\right)}{3}={A}_{n}^{1}+\frac{4n{R}^{2}}{3}sin\frac{\pi}{n}\left(1cos\frac{\pi}{n}\right),$$(3)
We can now compute the error of both formulations compared to the perfect circle surface area A = πr^{2 }= 78.53982 mm^{2}, in function of the number of nodes (or sides) n. The results are shown in Table 2.
Unsurprisingly, a significant amount of nodes is required to get close to the exact value of the circle area when linear elements are used, whereas quadratic elements are unhesitatingly more efficient. We can note that 100 linear elements give an error with same magnitude that 8 quadratic elements. It is important to remind that this efficiency has a cost, with more than the triple of integration points. And because the number of integration points is correlated to the CPU time, no conclusion can be extracted yet from these data alone. We must also compare the results accuracy and the solver time for all mesh density and types.
Fig. 5 3D singular coil model of the campaign, with the two coupled reference points. 
Fig. 6 Geometrical principles of the tessellation of a circular cross section with n = 6. 
Comparison of meshed section surface area in function of the order of the elements and the number of nodes along the circumference, with associated error with geometric perfect circular cross section area with radius R = 5 mm: A = 78.53982 mm^{2}.
4.2 Stiffness accuracy
The objective of a spring is to store energy and restore it. Therefore, its main operating parameter is its stiffness, which is the ratio between the axial load and the associated deflection. The stiffness behaviour of a spring can be represented by a load/length or load/deflection curve. Because the deflection is fixed in the model, the missing unknown variable is the associated load. This load can be measured at the end of the simulation at the driven reference point. All element types convergence analyses about this axial load are visible in Figure 7. All curves converge towards an axial load F = 5321.15N. C3D20 and C3D20R curves are superimposed, as are C3D8, C3D8I and C3D8S curves. Quadratic elements admit a very small error of 0.3% in the worst case (the 8nodes C3D20R). Meanwhile, linear elements have more difficulty to well predict the axial load (and stiffness) of the coil. C3D8R curve admits the largest error, due to the previously mentioned reduced integration problem for linear elements. Also, the quadratic elements overestimate the axial load with coarse meshes while linear elements underestimate it.
Some additional analyses can be extracted with the help of the previous geometric investigation. For linear elements, the load error is directly proportional to the geometric error. For C3D8, C3D8I and C3D8S, the load error is twice (in average) the geometric error of the element. This phenomenon is consistent for every mesh density. For C3D8R, this ratio is raised to 3.5. For quadratic elements, the load error is not directly proportional to the geometric error. The 2ratio with C3D8, C3D8I and C3D8S elements will be explained in Section 4.4.
Also, this load error can be compared with the CPU time for each simulation. In Figure 8, each point is the result of a simulation. The color of the point corresponds to the value of n and the shape of the point corresponds to the element type. A simple observation is that the higher the mesh density is, the more accurate the result get. This accuracy is directly proportional to the CPU time. Each of the linear elements admits a constant slope m (in the loglog graph), which the accuracy can be reduced by for increasing the CPU time. This value is estimated to be around m = –2/3. Therefore, to divide the load error by 100 for linear elements, you must increase CPU time by 1000. For quadratic elements, this slope is degraded by the numerical resolution (i.e. the number of digits kept by the software to make calculation). However, if we consider that the initial slope (n < 36) is constant, the quadratic slope is approximately m = –2. Therefore, to divide the load error by 100, you must increase the CPU time by only 10. This is another example of the strength of quadratic elements compared to linear ones. For quadratic elements, with only n = 16, the error is bellow 0.01%, while less than 10 times slower than the fastest simulation. Then, even if C3D8, C3D8I and C3D8S give similar results for each mesh density, as already visible in Figure 6, thanks to this plot it is possible to remark that C3D8 are always faster, then comes the C3D8I, and finally the C3D8S. So for the load analysis, “incompatible modes” and “improved surface stress visualisation” are more expensive for identical results, then they should be avoided. However, they could still have a better impact on stress calculation.
Fig. 7 Convergence analysis on the axial load value (normalized) in function of the density and the element type. 
Fig. 8 The error on the axial load estimation in function of the CPU time, for all simulations of the campaign. 
4.3 Stress accuracy
Then, the spring design must ensure that the material stay in the elastic domain. In order to evaluate that numerically, we measure the equivalent Von Mises Stress (or Tresca stress). This third parameter will impact the safety factor, which is a significant system parameter. Indeed, springs are intentionally built much stronger than needed for usual range to allow for emergency situations as impacts absorption or any other unexpected loads. Thus, the maximal stress value within the spring must be known with very good accuracy. The security of the system and its users depend on it. All element types convergence curves about this maximal Von Mises stress value are visible in Figure 9. The C3D8R elements are undoubtedly the worst elements in this scenario, having more than 19.6% of error for n = 16 and do not even succeed to reach the converged value with n = 100, still admitting an error of 3.2%. All other curves converge towards a maximal Von Mises stress value σ = 1311.54 MPa. This convergence is reached with about n = 24, which can be remarkable for linear elements. Indeed, linear elements are more efficient with low n than the quadratic elements. This accuracy can be explained by the errors in the load and the geometry, which are both underestimated and almost compensate each other. This will be further discussed in the Section 4.4. For C3D8I elements with n = 12, an anomaly occurs, certainly due to the previously mentioned 12nodes shear phenomenon induced by the only 2symmetry planes in the section.
As in load investigation, stress error can be compared with the CPU time for each simulation. In Figure 10, each point is the result of a simulation. For each element type, it is possible to evaluate the slope m (in the loglog graph), as previously performed with load graph. In increasing order, C3D8's slope is approximately m = –0.29, C3D8R's m = –0.33, C3D8I's m = –0.4, C3D20R's m = –0.5, C3D20's m = –0.75 and C3D8S's m = –1. Higher the absolute value of the slope is, more interesting are these element types. Thus, we can see here the strength of “Improved surface stress visualisation” in stress evaluation, even if this element admits a limit in accuracy, where the last two mesh densities do not progress anymore.
Fig. 9 Convergence analysis on the Von Mises stress value (normalized) in function of the density and the element type. 
Fig. 10 The error on the Von Mises stress estimation in function of the CPU time, for all simulations of the campaign. 
4.4 Relative error propagation
In this section, are discussed the reasons why there is a proportional 2ratio between the load error and the geometric error for linear elements, and why this ratio is compensated for the stress error. Both phenomena are explainable with simple arithmetic. This section is only applicable to linear elements without reduced integration (C3D8, C3D8I and C3D8S).
Let the relative error ϵ_{X} of the measure X be the ratio between the uncertainty ΔX over its mean value X. When the result X involves the product of two observed uncertain quantities a and b, the relative error in the result ϵ_{X} is equal to the sum of the relative error in the observed quantities ϵ_{a} and ϵ_{b} (Eq. (4)). This approximation is possible because the product Δa × Δb is infinitesimal compared to Δa or Δb. If the uncertainties are negative, the associated relative error is impacted. By contrast, the division of a by b would result in ϵ_{X} ≈ ϵ_{a} − ϵ_{b}.
$$\begin{array}{l}X=a\times b\hfill \\ X\pm \mathrm{\Delta}X=\left(a\pm \mathrm{\Delta}a\right)\left(b\pm \mathrm{\Delta}b\right)\hfill \\ =ab\pm a\mathrm{\Delta}b\pm b\mathrm{\Delta}a\pm \mathrm{\Delta}b\mathrm{\Delta}a\hfill \\ \pm \frac{\mathrm{\Delta}X}{X}=\pm \frac{\mathrm{\Delta}a}{a}\pm \frac{\mathrm{\Delta}b}{b}\pm \frac{\mathrm{\Delta}a}{a}\frac{\mathrm{\Delta}b}{b}\hfill \\ {\u03f5}_{X}\approx {\u03f5}_{a}+{\u03f5}_{b}\hfill \end{array}$$(4)
For springs, standard formulation [2,3] states that the axial load F is directly correlated to the wire diameter d at the power of 4 and the displacement u (Eq. (5)). This power is equivalent to the surface area squared. Thus, the relative error of the load is impacted twice by the relative error of the surface area. Then, the stiffness K is directly proportional to the load F and the inverse of the displacement u (Eq. (6)). Finally, spring stress standard formulation is based on beam theory loaded in torsion, which is a good approximation for most spring geometries. In beam theory loaded in torsion, the stress σ is directly correlated to the load and to the inverse of the torsion constant J, which is itself correlated to the beam diameter to the power of 4 (Eq. (7)). However, unlike in the Equation (4), both values a and b are known underevaluated. Therefore, the error of the load (previously explained twice the error of the surface area) and the error of the torsion constant (also twice the error of the surface area) compensate each other in the division, resulting in error for C3D8 less than 1% for n = 8.
$$F=\frac{G{d}^{4}}{8{n}_{a}{D}^{3}}\left({L}_{0}L\right)\propto {A}^{2}u\to {\u03f5}_{F}\propto 2{\u03f5}_{A}+{\u03f5}_{u},$$(5)
$$K=\frac{F}{u}\propto F{u}^{1}\to {\u03f5}_{K}\propto {\u03f5}_{F}{\u03f5}_{u},$$(6)
$$\sigma =\frac{FD}{2J}\propto F{A}^{2}\to {\u03f5}_{\sigma}\propto {\u03f5}_{F}2{\u03f5}_{A}.$$(7)
Because the surface area relative error −ϵ_{A} is always present due to the meshing tessellation degradation, we can approximate all the relative error in function of ϵ_{A}, depending if the spring boundary condition (BC) is set by displacement or by force.
Therefore, the selected BC has a significant impact on the accuracy of the results. As shown in Table 3, the spring driven by the axial force cannot evaluate accurately any of the desired measures (geometry, stiffness and stress), while the spring driven by deflection can at least compute the stress correctly. Also, it is important to remind that driving the deflection of the spring by displacement is closer to an experimental setup and that is how the standard want spring manufacturers to test their products. Thus, displacement BC must be preferred rather than force BC.
Propagation of the surface area relative error ϵ_{A} due to tessellation on the aproximated relative error of the displacement u, the axial load F, the stiffness K and the stress σ, for linear elements without reduced integration, for springs driven by displacement and by axial force.
4.5 Global accuracy
Now all three investigations were conducted, it is important to know how to concatenate them or compare them. The proposed method is to define a fourth evaluation: the global admissible accuracy. The global admissible accuracy is calculated as the minimal value of all three accuracy for each model. In this way, we ensure a maximal admissible error of ϵ_{T} for any selected model. A Pareto front is now visible on Figure 11, representing this global accuracy in function of the CPU time factor. Because we have the objective to maximize the accuracy and minimize the CPU time, we are looking for solution on the top left corner of the plot. No fully optimal solution exists in this case, where a model is both the more accurate and the fastest one. However, some of the best nearoptimal solutions are in this region of the plot. We summarize them into Table 4.
The CPU time needed to reduce the relative error of the model increases exponentially. The proposed solution must be before this significant slope increasing. As previously mentioned, the orthotropic behaviour of the 12nodes meshing density could affect the shear behaviour in particular cases. That is why we suggest as the optimal solution the 16nodes in circumference with C3D20 elements (the yellow Xcross in Fig. 11). The global accuracy of this solution is 99.8%, with being less than 25 times slower than the 8nodes C3D8 solution. With this proposition, if the calculated stress is 1500 MPa, the accuracy of the 16nodes C3D20 solution gives an uncertainty of ±3 MPa, which is relatively acceptable. The 8nodes C3D20 would give an uncertainty for the same stress of ±39 MPa, which is not as much negligible. However, the reader is free to use any of the solution of the Table 4 or Appendix B, but with full knowledge of the degradation in the result or solving time that could occur.
Fig. 11 The global estimated accuracy in function of the CPU time factor, for all simulations of the campaign. 
Proposition of the best nearoptimal combinations in mesh density n along the circumference and integration method for a given admissible global accuracy. The time factor is calculated with the lowest quality proposed solution: 8nodes C3D8.
5 Conclusion
Helical springs are widely used in engineering applications. In order to reduce cost in ‘try and error’ time consuming experimental campaigns, numerical simulations became an essential tool for engineers. Indeed, it saves considerable time in the ahead design phase of a project to ensure the feasibility of structures. However, these simulations run thanks to a lot of parameters, which all must be selected carefully to get access to reliable results. Each error in the modeling phase can lead to inaccurate, long or incoherent results. Two major degrees of freedom in the modeling step are the mesh density and the element type. Mesh density impact the tessellation magnitude of the solid, which can induce an associated discretization error. This tessellation degradation can be compensated with quadratic elements, but they need more time to compute. The element type impact the mathematical resolution procedure of each element, which can induce model resolution error.
In order to get both fast and accurate results, three convergence analyses and a mesh optimization are conducted with ABAQUS in this paper. First, a valuable literature review is realised to set all fixed parameters. Then, a statistical analysis of a 50papers corpus using helical spring FEA is examined. This preliminary investigation shows that no such optimized mesh solution already exists and that scientific researchers often underestimate the impact of these choices on their results. Then, the three convergence analyses are performed about the circular geometry accuracy of tessellated surface, the axial stiffness (and axial load) accuracy of the helical spring and the maximal Von Mises stress accuracy within the helical spring. The numerical campaign is led with 8 mesh densities along the circumference and 6 element types. These 48 models are examined and compared with the most fine meshing technique realised: 100 nodes in the circumference and quadratic fully integrated hexahedral elements (C3D20). The main observations of the literature review and the numerical campaign about the optimized meshing technique applied to helical spring are summarized as follows:
Choose hexahedral elements, equidistributed along the circumference and the wire.
Select “Ogrid” structured meshing technique within the cross section and “sweep” along the wire, guaranteeing low aspect ratio, smoothness and skewness.
Avoid as much as possible linear elements with reduced integration.
“Incompatible modes” and “Improved surface stress visualisation” do not offer more accurate results than “Fully integrated” elements but are more timeexpensive.
No mesh density less than 16 nodes in circumference should be used in order to perform reliable simulations.
Prefer quadratic elements to model curved shapes like the spring's coiled wire.
Drive the spring by displacement and not by load.
Finally, the mesh optimization process considering all three investigations is performed. Eight nearoptimal meshing techniques (density and element type) are highlighted. Each solution is reliable in function of the desired global admissible accuracy of the model. Nevertheless, the authors advice the 16nodes in circumference density with quadratic fully integrated hexahedral elements (C3D20). This mesh density and element type combination offers the best compromise between high accuracy and low CPU time. The proposed meshing technique can be exploited for any helical structure with circular cross section, mainly loaded in torsion and shear, such as extension and compression springs. A future study can be led with similar goals under bending and radial loads for torsion springs.
Even if the discretization error and the resolution errors are considerably reduced, it must be reminded that FEA are still not errorproof. As a computer does arithmetic, it introduces numerical errors by using numbers of finite precision (number of saved digits). Also, because the actual physical problem can not be examined directly, analytical method is applied to a model problem. This model admits material hypothesis, simplifying assumptions on boundary conditions or geometric restrictions. Each of them induces therefore modelling errors. Thus, engineers and researchers must remain critical towards their FEA results.
Future researches must investigate the dynamic behavior of the spring, where as well as the meshing parameters, the boundary conditions (i.e. the seat conditions of the end coils and the load application hypotheses) impact significantly the results.
This optimized meshing technique allows researchers, engineers and designers to evaluate accurately stress, stiffness and geometry of helical springs, saving simultaneously calculation time. Beyond saving industrial resources, saving calculation time is incorporated in the framework of the climate change and the modern environmental stakes. The environmental impact of numerical technology does not end with the raw materials extraction and manufacturing processes. Computer cores, data banks and servers uses also contributes to the greenhouse effect, responsible for global warming and air pollution by requiring vast quantities of energy each year. Each second of calculation time has a carbon cost and it is our goal, researchers and engineers, to reduce this cost while increasing technology performances. Thus, this method foresee to reduce computational carbon footprint for many numerical analysis of helical springs.
Funding
This research received no external funding.
Conflicts of interest
The authors disclosed receipt of the following financial support for the research of this article: this work was supported by CGR International, a private company.
Data availability statement
This article has no associated data generated and/or analyzed.
Author contribution statement
Conceptualization, GC and MP; Methodology, GC and MP; Software, GC; Validation, GC and MP; Formal Analysis, GC and MP; Investigation, GC and MP; Data Curation, GC; Writing − Original Draft Preparation, GC; Writing − Review & Editing, GC and MP; Visualization, GC; Supervision, MP.
References
 M. Paredes, Enhanced formulae for determining both free length and rate of cylindrical compression springs, ASME. J. Mech. Des. 138, (2015) [Google Scholar]
 NormeNFEN139061, Ressorts hélicoïdaux cylindriques fabriqués à partir de fils ronds et de barres, calcul et conception, partie 1 : Ressorts de compression, AFNOR (2002) [Google Scholar]
 IST, Essential Spring Design Training Course (Institute of Spring Technology, Sheffield, United Kingdom, 1980−2005) [Google Scholar]
 A.M. Wahl, Mechanical Springs (McGrawHill Book Company, Second Edition, 1963) [Google Scholar]
 G. Cadet, M. Paredes, A new exhaustive semianalytical method to calculate stress distribution on the surface of a curved beam with circular cross section, with an application to helical compression springs, Eur. J. Mech. − A (2023) [Google Scholar]
 G. Cadet, M. Paredes, H. Orcière, Improved analytical model for cylindrical compression springs not ground considering end behavior of end coils, Mech. Ind. 22, (2021) [Google Scholar]
 R. Palaninathan, P.S. Chandrasekharan, Curved beam element stiffness matrix formulation, Comput. Struct. 21, 663–669 (1985) [CrossRef] [Google Scholar]
 C.L. Dym, Consistent derivations of spring rates for helical springs, ASME J. Mech. Des. 131, (2009) [Google Scholar]
 M. Shimoseki, T. Hamano, and T. Imaizumi, FEM for Springs (SpringerVerlag Berlin, 2003) [CrossRef] [Google Scholar]
 F. De Crescenzo and P. Salvini, Influence of coil contact on static behavior of helical compression springs, IOP Conf. Ser.: Mater. Sci. Eng. 1038 (2020) [Google Scholar]
 F. Dammak, M. Taktak, S. Abid, A. Dhieb, M. Haddar, Finite element method for the stress analysis of isotropic cylindrical helical spring, Eur. J. Mech. A 24, 1068–1078 (2005) [CrossRef] [Google Scholar]
 A.Y. Babenko, B. Soltannia, P.S. Mobarakeh, Solving geometrically nonlinear problem on deformation of a helical spring through variational methods, Int. J. Mech. Appl. 8, 21–24 (2018) [Google Scholar]
 Z. Gu, X. Hou, J. Ye, Design and analysis method of nonlinear helical springs using a combining technique: finite element analysis, constrained latin hypercube sampling and genetic programming, J. Mech. Eng. Sci. 235, 5917–5930 (2021) [CrossRef] [Google Scholar]
 M. Taktak, F. Dammak, S. Abid, M. Haddar, A mixedhybrid finite element for threedimensional isotropic helical beam analysis, Int. J. Mech. Sci. 47, 209–229 (2005) [CrossRef] [Google Scholar]
 M. Taktak, F. Dammak, S. Abid, M. Haddar, A finite element for dynamic analysis of a cylindrical isotropic helical spring, J. Mech. Mater. Struct. 3, 641–658 (2008) [CrossRef] [Google Scholar]
 R. Provasi, C. d.A. Martins, A threedimensional curved beam element for helical components modeling, J. Offshore Mech. Arctic Eng. 136 (2014) [CrossRef] [Google Scholar]
 A.D. Kelly, C.E. Knight, Helical coil suspension springs in finite element models of compressors, Int. Compress. Eng. Conf. 870 (1992) [Google Scholar]
 A.N. Chaudhury, D. Datta, Analysis of prismatic springs of noncircular coil shape and nonprismatic springs of circular coil shape by analytical and finite element methods, J. Comput. Des. Eng. 4, 178–191 (2017) [Google Scholar]
 Y. Zhuo, Z. Qi, J. Zhang, G. Wang, A geometrically nonlinear spring element for structural analysis of helical springs, Arch. Appl. Mech. 92, 1789–1821 (2022) [CrossRef] [Google Scholar]
 M. Ermis, M.H. Omurtag, Static and dynamic analysis of conical helices based on exact geometry via mixed fem, Int. J. Mech. Sci. 131, 296–304 (2017) [CrossRef] [Google Scholar]
 J. Lee, Free vibration analysis of cylindrical helical springs by the pseudospectral method, J. Sound Vibr. 302, 185–196 (2007) [CrossRef] [Google Scholar]
 A.R. Udhaya, B. Rajeswari, T. Mugilan, Static structural investigation of helical compression spring utilizing different materials for an automobile suspension system, Mater. Today: Proc. 80, 653–658 (2023) [CrossRef] [Google Scholar]
 H.B. Pawar, A.R. Patil, S.B. Zope, Design and analysis of a front suspension coil spring for three wheeler vehicle, Int. J. Innov. Eng. Res. Technol. 3 (2016) [Google Scholar]
 S. Kushwah, S. Parekh, M. Mangrola, Optimization of coil spring by finite element analysis method of automobile suspension system using different materials, Mater. Today: Proc. 42, 827–831 (2021) [CrossRef] [Google Scholar]
 R. Puff, R. Barbieri, Effect of nonmetallic inclusions on the fatigue strength of helical spring wire, Eng. Fail. Anal. 44, 441–454 (2014) [CrossRef] [Google Scholar]
 W.G. Jiang, J.L. Henshall, A novel finite element model for helical springs, Finite Elements i Anal. Des. 35, 363–377 (2000) [CrossRef] [Google Scholar]
 T.M. Mulla, S.J. Kadam, V.S. Kengar, Finite element analysis of helical coil compression spring for three wheeler automotive front suspension, Int. J. Mech. Ind. Eng. 2, 74–77 (2012) [Google Scholar]
 H.B. Pawar, D.D. Desale, Optimization of three wheeler front suspension coil spring, Proc. Manufactur. 20, 428–433 (2018) [CrossRef] [Google Scholar]
 H. Font, G. Cadet, M. Paredes, H. Orcière, Enhanced formulae for determining solid height of axially guided compression springs with closed and unground ends, Wire Forming Technol. 25, (2022) [Google Scholar]
 M. Bakhshesh, M. Bakhshesh, Optimization of steel helical spring by composite spring, Int. J. Multidiscipl. Sci. Eng. 3, 47–51 (2012) [Google Scholar]
 R.R.D.A. Andoko, Coil spring type analysis using the finite element method, IOP Conf. Ser.: Mater. Sci. Eng. (2021) [Google Scholar]
 T.A. Jadhav, M.P. Angaj, V.N. Kapatkar, Finite element analysis of helical coil spring with cushioning buffer, Int. J. Eng. Res. Technol. (2019) [Google Scholar]
 A. Banerjee, Design and analysis of helical spring profiles in an electric vehicle suspension system using finite element method, Int. J. Adv. Res. Ideas Innov. Technol. (2020) [Google Scholar]
 I. Pöllänen, H. Martikka, Optimal redesign of helical springs using fuzzy design and fem, Adv. Eng. Softw. 41, 410–414 (2010) [CrossRef] [Google Scholar]
 D. Čakmak, Z. Tomičević, H. Wolf, Ž. Božić, D. Semenski, I. Trapić, Vibration fatigue study of the helical spring in the baseexcited inerterbased isolation system, Eng. Fail. Anal. 103, 44–56 (2019) [CrossRef] [Google Scholar]
 A. Jain, S. Misra, A. Jindal, P. Lakhian, Structural analysis of compression helical spring used in suspension system, AIP Conf. Proc. (2017) [Google Scholar]
 A. Tiwari, K.K. Ray, B. Pyttel, Very high cycle fatigue behavior of helical compression springs: numerical and experimental analysis, Thesis (2012) [Google Scholar]
 Y. Wang, C. Soutis, M. Yar, X. Zhou, Modelling corrosion effect on stiffness of automotivesuspension springs, Mater. Des. Process. Commun. 1 (2018) [Google Scholar]
 Y. Wang, Q. Wang, Z. Su, Numerical studies on the stiffness of arc elliptical crosssection helical spring subjected to circumference force, Mechanika 27, 327−334 (2021) [CrossRef] [Google Scholar]
 J. Ke, Z.Y. Wu, Y.S. Liu, Z. Xiang, X.D. Hu, Design method, performance investigation and manufacturing process of composite helical springs: a review, Compos. Struct. 252 (2020) [Google Scholar]
 L. Wu, L. Chen, H. Fu, Q. Jiang, X. Wu, Y. Tang, Carbon fiber composite multistrand helical springs with adjustable spring constant: design and mechanism studies, J. Mater. Res. Technol. 9, 5067–5076 (2020) [CrossRef] [Google Scholar]
 Y. Zhang, C. Yu, D. Song, Y. Zhu, Q. Kan, G. Kang, Solidstate cooling with high elastocaloric strength and low driving force via niti shape memory alloy helical springs: experiment and theoretical model, Mech. Mater. 178 (2023) [Google Scholar]
 Q. Jiang, Y. Qiao, F. Zhao, Z. Pan, X. Wu, L. Wu, H. Fu, Composite helical spring with skincore structure:structural design and compression propertyevaluation, Soc. Plast. Eng. Polym. Compos. 42, 1292–1304 (2021) [Google Scholar]
 B. Pyttel, K.K. Ray, I. Brunner, A. Tiwari, S.A. Kaoua, Investigation of probable failure position in helical compression springs used in fuel injection system of diesel engines, IOSR J. Mech. Civil Eng. 2 (2012) [Google Scholar]
 M. Baghani, R. Naghdabadi, J. Arghavani, A semianalytical study on helical springs made of shape memory polymer, Smart Mater. Struct. 21 (2012) [Google Scholar]
 A.F. Saleeb, B. Dhakal, M.S. Hosseini, S.A. Padula II, Large scale simulation of niti helical spring actuators under repeated thermomechanical cycles, Smart Mater. Struct. 22 (2013) [Google Scholar]
 X. Nong, W. Feng, J. Gao, C. Shi, N. Zhao, Stress relaxation constitutive relations and finite element analysis of t9a helical compression spring, Mater. Trans. 62, 962–967 (2021) [CrossRef] [Google Scholar]
 R. Mirzaeifar, R. DesRoches, A. Yavari, A combined analytical, numerical, and experimental study of shapememoryalloy helical springs, Int. J. Solids Struct. 48, 611–624 (2011) [CrossRef] [Google Scholar]
 C. ElMtili, A. Khamlichi, L. Hessissen, H.M.W. Badar, Forcedisplacement relationships for niti alloy helical springs by using ansys: superelasticity and shape memory effect, Int. Rev. Appl. Sci. Eng. 13 (2022) [Google Scholar]
 M. Muralidharan, R. Aravinth, J. Gafferkhan, R. Gandhi, Comparative design and analysis of helical and wave spring, Int. J. Eng. Technol. 7, 353–356 (2018) [CrossRef] [Google Scholar]
 H.A. Rasol, M.R. Ismail, A.A. Najam, Study the possibility of using fiber and polymer composite materials in helical spring manufacturing, IOP Conf. Ser.: Mater. Sci. Eng. (2021) [Google Scholar]
 S.N. Khurd, P.P. Kulkarni, S.D. Katekar, A.M. Chavan, Analysis of two wheeler suspension spring by using fea for different materials, Int. Res. J. Eng. Technol. 3, 833–839 (2016) [Google Scholar]
 A.I. Razooqi, H.A. Ameen, K.K.M. Mashloosh, Compression and impact characterization of helical and slotted cylinder springs, Int. J. Eng. Technol. 3, 268–278 (2014) [CrossRef] [Google Scholar]
 H.B. Pawar, A.R. Patil, S.B. Zope, Analysis and optimization of a helical compression coil spring used for twv, Int. J. Adv. Res. Innov. Ideas Educ. (2016) [Google Scholar]
 L. DelLlanoVizcaya, C. RubioGonzález, G. Mesmacque, T. CervantesHernández, Multiaxial fatigue and failure analysis of helical compression springs, Eng. Fail. Anal. 13, 1303–1313 (2006) [CrossRef] [Google Scholar]
 C. Stephen, R. Selvam, S. Suranjan, A comparative study of steel and composite helical springs using finite element analysis, Adv. Sci. Eng. Technol. Int. Conf. (ASET) (2019) [Google Scholar]
 M.R. Khudhair, Failure analysis of compression helical spring used in the suspension system by FEA, Int. J. Mech. Product. Eng. Res. Dev. 9 (2019) [Google Scholar]
 R. Sreenivasulu, N.Y. Krishna, M. Sukumar, O.N.G. Basha, N. ArunKumar, K. Heamanth, M.V. Krishna, Modeling and analysis of helical springs using catiav5r19 and ansys 16.0, AKGEC Int. J. Technol. 11, 41–50 (2020) [Google Scholar]
 D. Čakmak, Ž. Božić, H. Wolf, N. Alujević, Simultaneous vibration and fatigue optimization of an inerterbased vibration isolation system, Engineering (2017) [Google Scholar]
 L. Hou, Y. Hu, Gursontvergraadneedleman modelbased damage analyses of stainless steel springs at high temperature, J. Phys.: Conf. Ser. (2023) [Google Scholar]
 P. Sedlák, M. Frost, A. Kruisová, K. Hiřmanová, L. Heller, P. Šittner, Simulations of mechanical response of superelastic niti helical spring and its relation to fatigue resistance, J. Mater. Eng. Perform. 23, 2591–2598 (2014) [CrossRef] [Google Scholar]
 Siddharth, D. Yadav, S. Lata, Design development and analysis of cylindrical spring with variable pitch for two wheelers, Mater. Today: Proc. 47, 3105–3111 (2021) [CrossRef] [Google Scholar]
 Y. Wang, C. Soutis, L. Gagliardi, A finite element and experimental analysis of durability tested springs, MATEC Web Conf. 165 (2018) [Google Scholar]
 R. Mehrabi, M.R.K. Ravari, Simulation of superelastic sma helical springs, Smart Struct. Syst. 16, 183–194 (2015) [Google Scholar]
 S. Lutz, Kennlinie und eigenfrequenzen von schraubendruckfedern, Dissertation TU Ilmenau (2000) [Google Scholar]
 U. Kletzin, H.J. Schorcht, K. Zimmermann, H.K.U. Liebers, FiniteElementebasiertes Entwurfssystem für Federn und Federanordnungen, Technische Universität Ilmenau, Institut für Maschinenelemente und Konstruktion (2000) [Google Scholar]
 A. Roychoudhury, A. Banerjee, S. Dutt, S. Sinha, Studying the effect of electroless nickel coating on helical compression springs by finite element analysis, Int. J. Eng. Res. Technol. 6 (2017) [Google Scholar]
 T.M. Mulla, Fatigue life estimation of helical coil compression spring used in front suspension of a three wheeler vehicle, Proc. Mod. Era Res. Mech. Eng. (2016) [Google Scholar]
 K. Sathishkumar, G. Dinesh, Design and material analysis of a suspension system in scooter by using finite element analysis method, Int. Res. J. Multidiscip. Technov. 1, 25–37 (2019) [Google Scholar]
 A.S. Karad, P.D. Sonawwanay, C.Y. Bachhav, Finite element analysis of coil spring by using carbon fibre, carbon steel and epoxy resin materials, Mater. Today: Proc. (2023) [Google Scholar]
 V. Fegade, U. Ragavendran, M. Ramachandran, Numerical investigation of hybrid helical spring for total deformation and von mises analysis, Int. J. Mech. Product. Eng. Res. Dev. 181–186 (2018) [Google Scholar]
 R.D. Cook, D.S. Malkus, M.E. Plesha, R.J. Witt, Concepts and applications of finite element analysis, John Wiley and Sons, INC. (1974) [Google Scholar]
 L.A. Barba, Terminologies for reproducible research, Comput. Sci. (2018) [Google Scholar]
 F.C.Y. Benureau, N.P. Rougier, Rerun, repeat, reproduce, reuse, replicate: transforming code into scientific contributions, Comput. Sci. (2018) [Google Scholar]
 V. Stodden, Reproducing statistical results, Annu. Rev. Stat. Appl. 2, 1–19 (2015) [CrossRef] [Google Scholar]
 S. Tejesh, T. Srinath, Design and analysis of helical compression spring, Int. J. Innov. Res. Adv. Stud. 9 (2022) [Google Scholar]
 K. Sataynarayana, T. Ugesh, B. Gowri, D. Sai Adhitya Ganesh, K. Bharath, Design and static analysis on a helical spring for two wheeler, J. Compos. Theory 13 (2020) [Google Scholar]
 P. Ravinder Reddy, V. Mukesh Reddy, Determination of buckling loads of wave spring using ansys, Int. J. Res. Eng. Sci. 3, 48–56 (2015) [Google Scholar]
 I.A. Magomedov, Z.S. Sebaeva, Comparative study of finite element analysis software packages, J. Phys.: Conf. Ser. 1515 (2020) [Google Scholar]
 Massachusetts Institute of Technology, Abaqus documentation 2017, https://abaqusdocs.mit.edu/2017/English/SIMACAEEXCRefMap/simaexccdocproc.htm, (2023) [Google Scholar]
 G. Yazar, Design and analysis of helical coil spring forms for independent suspensions of automobiles, PhD Thesis, Graduate School of Natural and Applied Sciences (2015) [Google Scholar]
 G. Prathap, The poor bending response of the fournode plane stress quadrilateral, Int. J. Numer. Methods Eng. 21, 825–835 (1985) [CrossRef] [Google Scholar]
 D.P. Flanagan, T. Belytschko, A uniform strain hexahedron and quadrilateral with orthogonal hourglass control, Int. J. Numer. Methods Eng. 17, 679–706 (1981) [CrossRef] [Google Scholar]
 T. Belytschko, J.S.J. Ong, W.K. Liu, J.M. Kennedy, Hourglass control in linear and nonlinear problems, Comput. Methods Appl. Mech. Eng. 43, 251–276 (1984) [CrossRef] [Google Scholar]
 W. Lowrie, V.S. Lukin, U. Shumlak, A priori mesh quality metric error analysis applied to a highorder finite element method, J. Comput. Phys. 230, 5564–5586 (2011) [CrossRef] [MathSciNet] [Google Scholar]
 P.M. Knupp, Remarks on mesh quality, American Institute of Aeronautics and Astronautics Paper − 45th Aerospace Sciences Meeting and Exhibit, Reno, NV (2007) [Google Scholar]
 J.R. Sack, J. Urrutia, Handbook of computational geometry (Elsevier Science B.V., NorthHolland, 2000) [Google Scholar]
 G. Cadet, M. Paredes, H. Orciere, Improved design of singlelayered wire strand for combined tensile and crimping application with meshing optimization, DYNA 98, 274–281 (2023) [CrossRef] [Google Scholar]
 O.C. Zienkiewicz, R.L. Taylor, The Finite Element Method, Volume 2: Solid Mechanics (Butterworth Heinemann, Oxford, 2000) [Google Scholar]
 W.H. Cai, J.M. Zhan, Y.Y. Luo, Userintervened structured meshing methods and applications for complex flow fields based on multiblock partitioning, J. Comput. Des. Eng. 8 (2020) [Google Scholar]
 C.G. Armstrong, H.J. Fogg, C.M. Tierney, T.T. Robinson, Common themes in multiblock structured quad/hex mesh generation, Proc. Eng. 124, 70–82 (2015) [CrossRef] [Google Scholar]
 H.J. Fogg, L. Sun, J.E. Makem, C.G. Armstrong, T.T. Robinson, Singularities in structured meshes and crossfields, Comput. Aided Des. 105, 11–25 (2018) [CrossRef] [MathSciNet] [Google Scholar]
Cite this article as: G. Cadet, M. Paredes, Convergence analysis and mesh optimization of finite element analysis related to helical springs, 25, 22 (2024)
All Tables
Comparison of the labels and configurations procedures of mesh geometric order and integration method for element in ABAQUS and ANSYS software packages.
Comparison of meshed section surface area in function of the order of the elements and the number of nodes along the circumference, with associated error with geometric perfect circular cross section area with radius R = 5 mm: A = 78.53982 mm^{2}.
Propagation of the surface area relative error ϵ_{A} due to tessellation on the aproximated relative error of the displacement u, the axial load F, the stiffness K and the stress σ, for linear elements without reduced integration, for springs driven by displacement and by axial force.
Proposition of the best nearoptimal combinations in mesh density n along the circumference and integration method for a given admissible global accuracy. The time factor is calculated with the lowest quality proposed solution: 8nodes C3D8.
All Figures
Fig. 1 Examples of mesh patterns in circular cross section with 16 nodes along the circumference: free tetrahedra (a), free hexahedra (b), structured tetrahedra (c), structured “pie” hexahedra (d) and structured “curved Ogrid” hexahedra (e). 

In the text 
Fig. 2 Recurrence of each mesh density among the corpus [22–71]. The sum of all the bars adds up to 54 because some papers use several mesh densities. 

In the text 
Fig. 3 Meshed sections of the realised campaign with 8nodes circumference (a), 12nodes (b), 16nodes (c), 24nodes (d), 36nodes (e), 50nodes (f), 75nodes (g) and 100nodes (h). 

In the text 
Fig. 4 Von Mises stress distribution in the section of the 12nodes circumference sections, with the 6mesh diameter horizontal (a) and vertical (b). 

In the text 
Fig. 5 3D singular coil model of the campaign, with the two coupled reference points. 

In the text 
Fig. 6 Geometrical principles of the tessellation of a circular cross section with n = 6. 

In the text 
Fig. 7 Convergence analysis on the axial load value (normalized) in function of the density and the element type. 

In the text 
Fig. 8 The error on the axial load estimation in function of the CPU time, for all simulations of the campaign. 

In the text 
Fig. 9 Convergence analysis on the Von Mises stress value (normalized) in function of the density and the element type. 

In the text 
Fig. 10 The error on the Von Mises stress estimation in function of the CPU time, for all simulations of the campaign. 

In the text 
Fig. 11 The global estimated accuracy in function of the CPU time factor, for all simulations of the campaign. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.