Research Article  Open Access
Fatima Boussaoui, Hassane Lahmam, Bouazza Braikat, "Numerical HighOrder Model for the Nonlinear Elastic Computation of Helical Structures", Modelling and Simulation in Engineering, vol. 2021, Article ID 6655909, 15 pages, 2021. https://doi.org/10.1155/2021/6655909
Numerical HighOrder Model for the Nonlinear Elastic Computation of Helical Structures
Abstract
In this work, we propose a highorder algorithm based on the asymptotic numerical method (ANM) for the nonlinear elastic computation of helical structures without neglecting any nonlinear term. The nonlinearity considered in the following study will be a geometric type, and the kinematics adopted in this numerical modeling takes into account the hypotheses of Timoshenko and de SaintVenant. The finite element used in the discretization of the middle line of this structure is curvilinear with twelve degrees of freedom. Using a simple example, we show the efficiency of the algorithm which was carried out in this context and which resides in the reduction of the number of inversions of the tangent matrix compared to the incremental iterative algorithm of NewtonRaphson.
1. Introduction
Helical structures, particularly springs, are considered to be the most used mechanical elements in several industries such as aeronautics, automotive, and civil engineering. Due to their reversible elastic deformability, springs can store or release energy when they are subjected to different loads such as normal force, torque, or both [1, 2].
Numerical modeling of the linear elastic behavior of helical springs has been studied by many searchers. Taktak et al. [3] proposed a mixed hybrid formulation, including shear effects, based on the use of a curvilinear finite element with two nodes and six degrees of freedom per node. Dammak et al. [4] have shown that it is possible to obtain the distribution of the various generalized stresses along the middle line of the helical spring with great precision by using a single curvilinear finite element with twelve degrees of freedom. Taktak et al. [5, 6] also used this same finite element for the calculation of the natural frequencies and the dynamic response of a simple or assembled helical spring. Other researchers [7–11] have used theoretical and numerical methods for the dynamic behavior study as well as the buckling of coil springs. However, most of these studies are limited to the linear elastic case.
The objective of this work is to propose a highorder algorithm based on the asymptotic numerical method (ANM) [12–17] to calculate the displacements and rotations at a point of the middle line of a helical structure. The originality of our numerical approach is to calculate these unknowns without neglecting any nonlinear term of the strain. The finite element used in our numerical approach is of curvilinear type, defined along the middle line of the structure considered [18, 19]. The kinematics adopted in our theoretical formulation takes into account the hypotheses of Timoshenko and de SaintVenant [3–6]. The efficiency of the proposed algorithm lies in the reduction of the number of inversions of the tangent matrix compared to the incremental iterative algorithm of NewtonRaphson [20, 21].
2. Problem Modeling
We consider a helical structure of a circular section having a middle radius , a constant pitch , and a length as shown in Figure 1, where and is an arbitrary value.
The position vector of each point of the middle line is written, in the Cartesian coordinate system base, as follows:where is the polar angle.
The direct orthonormal basis of the Frenet coordinate system , defined at point , is expressed in the Cartesian basis as follows:where .
The pass matrix from the Cartesian coordinate system to the Frenet coordinate system is given by
The derivatives of these vectors with respect to the curvilinear abscissa are given bywhere and designate, respectively, the radius of curvature and the radius of torsion, defined bywith .
Using the Timoshenko and de SaintVenant hypotheses [3–6], the components of the displacement vector at a point in the cross section are written in the Frenet coordinate system as follows:where and are the Cartesian coordinates in the cross section. , , , , , and are the displacements and rotations at a point in the middle line as shown in Figure 2.
The gradient of the displacements at each point of the cross section is given by
The GreenLagrange strain vector is given by the following formula:where , , and , with , , and designating the membrane strains; , , , and are the strains in the plane ; and , , , , and are, respectively, the strains in the planes and; and and are the shear strains (see Appendix A).
The stationarity of the potential energy of the studied structure results inwhere represents the elastic strain energy, is the virtual work of the external forces, and is the loading parameter.
The expression of the elastic strain energy is given bywhere is the second PiolaKirchhoff tensor.
The elementary volume at any point of the considered structure is written as
The variation of the elastic strain energy is given by the following expression:where and are, respectively, the vectors of the generalized stresses and strains given bywhere is the normal effort, and are the shearing forces, and are the bending moments along the transversal axes, are quantities equivalent to the moment, and are additional quantities (see Appendix B).
The gradient vector of generalized displacements is written in the following form:
Taking into account the expression (14), the vector of generalized strains can also be written in the following form:where is a constant matrix which depends on the radius of curvature and is a matrix which depends on the gradient of the unknown of the considered problem (see Appendix C).
The generalized behavior law in nonlinear elasticity is given bywhere is the elastic constant matrix which also depends on the radius of curvature of the studied structure (see Appendix D).
By virtue of the previous relations, problem (9) is written:
3. Resolution Strategy
The solution of the nonlinear problem (17) is obtained using a highorder algorithm based on the discretization of the study domain into finite elements, the development of the different unknowns in the Taylor series, and the continuation technique.
3.1. Finite Element Discretization
The discrete form of the problem (17) is obtained by applying the finite element method [18, 19]. The degrees of freedom of each curvilinear element as shown in Figure 3 are represented by the displacement vector as follows:
The interpolation of the displacement vector in each element is given bywhere the interpolation functions are given bywhere is the length of the finite element “” and is the number of nodes.
The elementary displacement vector and elementary displacement gradient are given bywhere and are, respectively, the matrix of interpolation functions and the matrix of their derivatives with respect to the curvilinear abscissa (see Appendix E).
The variations of elementary vectors are written:where is the variation of the nodal displacement vector.
Thus, the discretized form of the global problem (17) is written as follows:where and is a matrix which depends on the unknown vector .
3.2. Taylor Series Expansion
In this stage, we write the unknown vectors, as well as the loading parameter in the form of Taylor series developments truncated to the order , that is to say,where , , and are known solutions at the initial point and “” is the parameter of the asymptotic development defined by [15]
By injecting the developments (24) into (23) and (25) and by identifying the coefficients according to the different powers of the parameter “,” we obtain a succession of linear problems given as follows:
Problem at order :
Problem at order :where the tangent stiffness matrix of the discretized structure, evaluated at the initial point, is given by
is a vector which depends on solutions to the previous orders, given byand is the matrix of initial stresses (see Appendix F).
Thus, the computation of each asymptotic branch of the nonlinear problem (23) requires a single inversion of the tangent matrix. The validity range of the series (24) is defined by the following criterion [17]:where is the tolerance parameter.
3.3. Continuation Technique
For the computation of the entire solution of the nonlinear problem (23), we apply the continuation technique [17] which consists in updating the tangent matrix (28) at the new initial point defined by , , and .
4. Results and Discussions
In this part, we consider a helical structure whose mechanical and geometric characteristics are mentioned in Table 1. One of the ends of this structure is embedded, the other end is subjected to external actions represented by the torsor as shown in Figure 4.

Using this example, we solve the nonlinear problem (23) by the proposed highorder algorithm and the NewtonRaphson algorithm. The calculated displacements are expressed in the Cartesian coordinate system and the rotations in the Frenet coordinate system. In this application, we first study the influence of the mesh and the truncation order of the Taylor series on the solution quality of the problem (23) expressed here by the decimal logarithm of the norm of the residual vector (see Figures 5 and 6). For that, we have chosen in the first case a truncation order and in the second case a mesh with 59 elements and this by using the same tolerance parameter .
As shown in Figure 5, the optimal mesh corresponds to which ensures a good quality of the solution of the considered problem, expressed here by . Figure 6 shows that if we use a mesh with 59 elements in the proposed algorithm, then the optimal truncation order corresponds well to the value which ensures the same quality of the solution represented in Figure 5.
In the following, we are interested in the comparison of computation times of CPU used by the highorder algorithm and the iterative incremental algorithm of NewtonRaphson and this for various meshes as shown in Table 2. In this numerical experiment, we have fixed the number of ANM steps, that is to say, .

As shown in this table, the use of the NewtonRaphson algorithm requires for each mesh a high number of inversions of the tangent matrix in comparison with the proposed algorithm.
In this application, we also studied the variation of the degrees of freedom , , , , , and of the loaded node according to the parameter of loading by using the following data: , , , and for the proposed algorithm and for the NewtonRaphson algorithm.
From Figures 7 and 8, we note the good concordance between the numerical solutions obtained by our highorder modeling and those calculated by the iterative NewtonRaphson algorithm.
In Figure 9, we represent the initial state and the middle line deformation of the considered structure, and in Figure 10, we give the evolution of the solution quality of the solved problem according to the loading parameter by using both algorithms.
From Figure 10, we note that our numerical modeling based on the ANM allows us to calculate the displacements and rotations of the middle line of the considered structure with a better quality compared to the results obtained by the NewtonRaphson algorithm.
5. Conclusion
In this work, we propose a highorder numerical modeling for the nonlinear elastic computation of a thin structure of helical type without neglecting any nonlinear term for the first time. The algorithm realized in this context is based on the asymptotic numerical method steps. The finite element used in the discretization of the middle line of the structure considered is of curvilinear type with twelve degrees of freedom. Using a simple example, we have shown the efficiency of our numerical approach in comparison with the NewtonRaphson incremental iterative algorithm. First, we discussed the optimum choice of the mesh and the truncation order of the series allowing having the solution of the chosen example with a better quality. In the second stage, we presented a numerical analysis showing that the proposed highorder algorithm significantly reduces the number of inversions of the tangent matrix, and thus the CPU computation time, to have the entire solution of the chosen problem and this in comparison with the reference algorithm mentioned previously.
Appendix
A. Generalized Deformations
B. Generalized Stresses
C. Matrices and
with
D. Matrix
wherewith
E. Matrices and
F. Matrix
with
Data Availability
No data were used to support this study.
Consent
No individual participants were included in the study so there are no subjects for which informed consent requirements arise.
Conflicts of Interest
The authors declare that they have no conflict of interest.
References
 A. M. Wahl, Mechanical Springs, McGrawHill, New York, 2nd edition, 1963.
 S. P. Timoshenko, Mechanics of Materials, Chapman & Hall, 1991.
 M. Taktak, F. Dammak, S. Abid, and M. Haddar, “A mixedhybrid finite element for threedimensional isotropic helical beam analysis,” International Journal of Mechanical Sciences, vol. 47, no. 2, pp. 209–229, 2005. View at: Publisher Site  Google Scholar
 D. Fakhreddine, T. Mohamed, A. Said, D. Abderrazek, and H. Mohamed, “Finite element method for the stress analysis of isotropic cylindrical helical spring,” European Journal of MechanicsA/Solids, vol. 24, no. 6, pp. 1068–1078, 2005. View at: Publisher Site  Google Scholar
 M. Taktak, F. Dammak, S. Abid, and M. Haddar, “A finite element for dynamic analysis of a cylindrical isotropic helical spring,” Journal of Mechanics of Materials and Structures, vol. 3, no. 4, pp. 641–658, 2008. View at: Publisher Site  Google Scholar
 M. Taktak, K. Omheni, A. Aloui, F. Dammak, and M. Haddar, “Dynamic optimization design of a cylindrical helical spring,” AppliedAcoustics, vol. 77, pp. 178–183, 2014. View at: Publisher Site  Google Scholar
 Y. Lin and A. P. Pisano, “General dynamic equations of helical springs with static solution and experimental verification,” Journalof Applied Mechanics, vol. 54, no. 4, pp. 910–917, 1987. View at: Publisher Site  Google Scholar
 Y. Xiong and B. Tabarrok, “A finite element model for the vibration of spatial rods under various applied loads,” International Journal of Mechanical Sciences, vol. 34, no. 1, pp. 41–51, 1992. View at: Publisher Site  Google Scholar
 V. Yıldırım, “Numerical buckling analysis of cylindrical helical coil springs in a dynamic manner,” International Journal of Engineering & Applied Sciences, vol. 1, no. 1, pp. 20–32, 2009. View at: Google Scholar
 V. Yıldırım, “On the linearized disturbance dynamic equations for buckling and free vibration of cylindrical helical coil springs under combined compression and torsion,” Meccanica, vol. 47, no. 4, pp. 1015–1033, 2012. View at: Publisher Site  Google Scholar
 I. Kacar and V. Yildirim, “Free vibration/buckling analyses of noncylindrical initially compressed helical composite springs,” Mechanics Based Design of Structures and Machines, vol. 44, no. 4, pp. 340–353, 2015. View at: Publisher Site  Google Scholar
 N. Damil and M. PotierFerry, “A new method to compute perturbed bifurcations: application to the buckling of imperfect elastic structures,” International Journal of Engineering Science, vol. 28, no. 9, pp. 943–957, 1990. View at: Publisher Site  Google Scholar
 M. PotierFerry, N. Damil, B. Braikat et al., “Treatment of strong nonlinearities by the asymptotic numerical method,” Comptes Rendus de l'Académie des Sciences  Series IIB  MechanicsPhysicsChemistryAstronomy, vol. 324, no. 3, pp. 171–177, 1997. View at: Publisher Site  Google Scholar
 B. Cochelin, N. Damil, and M. PotierFerry, Méthode asymptotique numérique, Hermés Lavoisier, 2007.
 H. Mottaqui, B. Braikat, and N. Damil, “Discussion about parameterization in the asymptotic numerical method: application to nonlinear elastic shells,” Computer Methods in Applied Mechanics and Engineering, vol. 199, no. 25–28, pp. 1701–1709, 2010. View at: Publisher Site  Google Scholar
 H. Mottaqui, B. Braikat, and N. Damil, “Local parameterization and the asymptotic numerical method,” Mathematical Modelling of Natural Phenomena, vol. 5, no. 7, pp. 16–22, 2010. View at: Publisher Site  Google Scholar
 B. Cochelin, “A pathfollowing technique via an asymptoticnumerical method,” Computers & Structures, vol. 53, no. 5, pp. 1181–1192, 1994. View at: Publisher Site  Google Scholar
 K. J. Bathe, Finite Element Procedure in Engineering Analysis, Pren ticeHall, Englewood Cliffs, 1982.
 O. C. Zienkiewicz and R. L. Taylor, The Finite Element Method, Berk Shire: Mc GrawHill, 1991.
 O. C. Zienkiewicz, “Incremental displacement in nonlinear analysis,” International Journal for Numerical Methods in Engineering, vol. 3, no. 4, pp. 587588, 1971. View at: Publisher Site  Google Scholar
 M. A. Crisfield, “A faster modified newtonraphson iteration,” Computer Methods in Applied Mechanics and Engineering, vol. 20, no. 3, pp. 267–278, 1979. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2021 Fatima Boussaoui et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.