Comptes Rendus Mécanique

. This paper concerns numerical simulations of time-dependent problems in computational solid mechanics. A perturbation method, with the time as perturbation parameter, is proposed to solve two classical problems: an elastic bar excited by an end force and the dynamic buckling of a cylindrical panel. Specific quadratic recast of the equations is proposed to solve the nonlinear problems. Numerical results show that asymptotic time expansions is robust, e ffi cient and gives more accurate solutions than the ones obtained with classical time-integration schemes (implicit or explicit), even when the considered meshes are coarse.

This study is focused on the use of the Asymptotic Numerical Method for time integration of equations arising in structural dynamics.
Over the past decades, several time integration schemes were proposed to solve transient dynamics problems.Newmark methods were established to be applied directly for second order differential equations [1].Thus, the Hilber-Hugues-Taylor (HHT) method was established to enhance numerical properties, such as maintaining second order accuracy and controlling numerical damping [2].A two parameters scheme was proposed to generalize time solver for second order systems of structural mechanics, the so-called generalized-α method [3].The readers can refer to the monograph [4] for more details.
In non-linear structural mechanics, Energy-Momentum methods were designed to ensure unconditional stability [5].Main advantages of this class of algorithms are the well-known numerical damping features, allowing for larger time steps, and the use of adaptive time stepping schemes.
More recently, a serie of studies [6][7][8][9] were focused on implicit time integration for wave propagations in structural dynamics.Let cite the proposal of the ρ ∞ -Bathe method.This latter is based on a composite strategy for time integration containing as special cases the Bathe and Newmark methods.
This topic also emerges in non-smooth dynamics featured by a lack of regularity, such as contact or multibody problems.For this kind of problems, the previous time-stepping schemes are used [10].Let note that alternative approaches based on event-driven schemes may be implemented [11,12].
Previous approaches produce discrete time solutions.An alternative improvment is to compute a continuous solution in time.Thus, perturbation methods can be used and time variable becomes the perturbation parameter.
In this context, the Asymptotic Numerical Method is a semi-analytical method where the unknwons are expanded in power series with respect to the time variable.By grouping terms of the same power of the time variable, a set of linear systems are obtained with same matrices.
This strategy was employed to solve ordinary differential equations encountered in various fields of physics.Various problems have been solved, such as neutron kinetics equations [13,14], trajectory of a conservative non-linear pendulum [15] or the famous Duffing equation [16].
The Asymptotic Numerical Method was also employed to solve linear and non-linear structural problems in different ways.In [17], the ANM was used as an implicit time solver to study dynamic buckling of a cylindrical shell.Discrete solution was computed with a generalized αmethod as time solver.An homotopy technique was also used to keep valid for several time steps the Jacobian matrix in order to save computational time.The buckling instant was perfectly detected.Unfortunatly, this strategy failed to compute the entire response once the buckling appeared, the Jacobian matrix becoming singular.
Probably first attempts to use ANM as an explicit time solver [18,19], it was proposed to use this semi-analytical method to compute time-dependent displacement solution.Velocity and acceleration were easily computed by differentiation.Nevertheless, these studies were restricted to linear elastic structure excited by an harmonic external sollication.
Thus, it is proposed to apply the methodology presented in [16,19] to compute transient responses for structural dynamics.Problems of the elastic wave propagation in an elastic bar and the dynamic buckling of a cylindrical shell are adressed combining an ANM based time solver with finite element spatial discretizations.
Solutions of transient dynamic problems are not analytical in times, then it may seem surprising to use time asymptotic expansions to solve this kind of problems.Nevertheless, as underlying by Deeb et al. in [20], there are two possibilities when using time perturbation method for solving non linear transient problems.The first one is that the series is convergent and then a continuation technique based on the evaluation of the validity range of the asymptotic expansions is sufficient to get the whole non linear solution.Padé approximants [21], for example, can then be used to increase the validity range and consequently decrease the step numbers of the continuation technique.The second possibility is that this validity range becomes very small or worstly null.In this case, a resummation procedure must to be applied to the series expansion and then lead to a time-analytical representation of the solution.Several resummation methods exist in the literature : for example Borel-Padé-Laplace (BPL) [22], Inverse factorial series [23] or generalized hypergeometric Meijer G-function [24].Let us note that the BPL scheme seems to be today the only numerical method that can be used to solve transient linear or non linear problems.Aim of this study is the computation of semi-analytical continuous time solutions even for non linear problem without using any resummation technique.
Equations of motion are derived and the ANM algorithm is recalled in Section 1. Advantages and drawbacks of the method are exposed and influences of numerical parameters are evaluated on numerical precision and computation efficiency in Section 2. Finally, some recommandations are proposed in Section 3 to enhance the methodology.

Formulation of the problem
Let t be the time variable and u the displacement unknown.Dot symbol stands for time derivative d u d t = u.For an elastic structure with structural damping and submitted to some external sollicitation, equations of motion R(u, u, ü) = 0 read: where L and D are linear operators standing for elastic and damping properties and M is a constant mass operator.The non-linear operator Q stands for the geometrically non-linear effects.Operator F contains the forcing of the structure.Vectors u, u and ü stand, respectively, for displacement, velocity and acceleration fields of the structure.
For spatial discretization, a finite element method is applied to equation (1).Finite dimensional operators M, L, N and F are the discretized counterpart of previous operators in (1) and vector q is the vector of nodal unknowns associated with the variational formulation.Thus, the following set of ordinary differential equations R(q, q, q) = 0 is obtained: To solve this problem, it is proposed to apply the ANM on the variable t = t 0 + t , by considering t as the perturbation parameter and t 0 stands for the initial time.

A perturbation-continuation based solver
The use of the ANM to solve non-linear quasi-static problem in structural mechanics is well established.A high-order perturbation method is associated with a continuation technique in order to determine a semi-analytical solution [25,26].This procedure was applied in fluid mechanics to determine bifurcation points [27,28] and to determine global and local elastic instabilities [29][30][31].It is worth noticing that solution obtained with the ANM is independent from the discretization technique.It was shown that meshless method [32,33], finite difference method [34] and finite element method [29,35] can be used.Perturbation in time is introduced as t = t 0 + t .Vectors q and F are sought as power series of time t of order N and it is supposed that q 0 and F 0 at instant t 0 are known: Expressions (3) are introduced in equation (2).After identifying same power of t , a set of linear systems are obtained.A key point in applying ANM is the quadratic recast of eq. ( 2), specially the operator Q. Full details can be found in [36,37].Expansion of external force F requires a specific treatment with introduction of a regularization parameter (see Section 1.3).This procedure was defined in previous works for calculation of stress in plastic buckling [38,39].
Terms q 0 and q 1 supposed to be knowns, identification of equal power of t allows to define the 2nd order term: 2Mq 2 + L(q 0 ) + D(q 1 ) + Q(q 0 , q 0 ) − F 0 = 0 (4) Defining the operator L q 0 (q k ) = L(q k )+Q(q 0 , q k )+Q(q k , q 0 ), terms of higher order (k = 1, . . ., N −2) can be computed: Solving problem (2) consist in solving eq. ( 4)- (5).To progress in time, a continuation technique is applied.The validity of one step is evaluated according to the following criteria adapted from works on buckling analysis with the ANM [26] where δ is a user-defined numerical parameter: A new continuation step is initialized from solution obtained at the actual step according to the following procedure: Because solution q is expanded into power series, computation of velocity by differentiation is done very easily.Finally, let note that this process is initialized, for the very first step, by knowledge of the initial conditions q ini and qini : q 0 = q ini q 1 = qini (8)

Recast of the external forcing F
In transient dynamics, structures are submitted to external forcing F(t ).Models of forcing are generally combination of Dirac, Heaviside, geometric and hypergeometric functions.To set these kind of functions in a quadratic framework, specific recast is necessary.
In this study, a ramp function (see Figure 6b) is considered as forcing on cylindrical shell and can be defined by: where s = f c /t c , t c and f c stand, respectively, for the ramp slope, the duration and the maximum intensity of forcing.The symbol η stands for a regularization parameter.The quadratic recast is inspired from previous works on yield stress calculation for elastic/perfectly plastic material [38,39].After definition of auxiliary variable g (t ) = f (t ) 2 and parameters ξ = f 2 c , ζ = sξ and γ = ηξ, expression ( 9) is written: This quadratic recast allows to define power series expansion of f (t ) and g (t ) and to identify terms of identical order in t .Thus, the following sets of linear systems are obtained: Order 0: Order 1: These equations are invertible and require knowledge of f 0 , from which g 0 is easily deduced.
As for all power series expansion, a validity range must be defined according to the previous definition (6).Two validity durations t q max and t f k ,g k max can be defined from respective series and the lowest value must be considered.From a pratical point of view, the forcing signal has smoother temporal variations than the structural response.Thus, the continuation technique is performed on the lowest value which is, in general, t q max .

Numerical studies
In this section, two examples are considered.The first case-study is the well-known transient response of an elastic bar excited by an end force.Exact solution can be found in [4,40].The second case-study is related to transient dynamics of a cylindrical roof submitted to a point dead force.This problem can be considered as a reference one which was deeply examined [17,41].Numerical results are obtained with a spatial finite element discretization and the ANM based time integration procedure, previously introduced.

Elastic wave propagation in a rod
The model is considered as a rod of length L = 1 m which is clamped at its left hand-side on x = 0 m while its extremity on x = L m is submitted to an Heaviside force f (t ) = F 0 (for t ≥ 0).Force intensity and cross-section of the rod are set to F 0 = 1 N and A = 0.01 m 2 , respectively.A theoretical elastic material is defined according to Young modulus E = 100 Pa and mass density ρ = 1 kg.m −3 .Thus, longitudinal wave celerity c 0 = E /ρ = 10 m.s −1 is defined.Duration of simulation T max = 0.8 s corresponds to the time required for both incident and reflected waves to travel 4 round trips at speed c 0 .Four tracers are defined along the rod which abscissa are x 1 = 0.25 m, x 2 = 0.5 m, x 3 = 0.75 m and x 4 = 1.00 m (rod tip).The finite element model is composed of N F E = 20 elements which are based on linear Lagrange shape function.For this study, the lumped mass matrix is adopted.
For comparison purpose, following author's write convention [4], the ANM based solver is compared to the Central Difference scheme (CD scheme, γ = 1/2, β = 0), the Purely Explicit scheme (PE scheme, γ = β = 0) and the "Average constant acceleration" scheme from the Implicit Newmark schemes family (γ = 1/2, β = 1/4).Constant time-step ∆t is also adopted and calculated as a percentage of the critical time-step ∆t cr .This latter corresponds to the time required for a longitudinal wave to travel through an element of length L/N F E at the constant wave propagation speed c 0 .In this study, ∆t cr = 5 × 10 −3 s.
In Figure 1, transient responses at locations x 1 and x 4 and obtained with the AMN solver and the CD-scheme are compared with the exact solution.Phenomena of wave displacement, reflection and superposition are observed.Because the mass matrix is lumped, the CD scheme allows to recover the exact solution for ∆t = ∆t cr .In this figure, CD solution is computed with a ratio ∆t /∆t cr = 0.99 and exhibits the wave propagation.Nor attenuation neither periodicity error can be observed.Only sharp ripples of small amplitude can be noticed due to a time-step  different from the critical one.The ANM solution is obtained with a truncation order N = 10 and a continuation parameter δ = 10 −8 .Results show numerical solutions which preserve the wave propagation nature.But an attenuation of wave amplitudes, a periodicity error and oscillations are also noticed.
In order to better understand the ANM solver, comparisons are made with schemes relied on a Pure Explicit method and an implicit Newmark β-method.
In   In Figure 3, ANM solution at the location x 1 on the rod is compared to solutions obtained with the implicit Newmark β-method with parameter values set to β = 1/4 and γ = 1/2.Newmark solutions are computed for different time-step values.Because of the implicit nature of the solver, solutions are less influenced by the time-step value.In this case, time-step values correspond to 99%, 10%, 1%, 0.5% and 0.1% of the critical time-step.For values lesser than 1% of this later, no significant difference is visible.Only a very fine observation allows to distinguish numerical solutions.In this example, ANM solution is the one closest to the Newmark solution computed with the lowest time-step ∆t = 10 −3 × ∆t cr = 5 × 10 −6 s.Again, this value is of three orders magnitude lower than the value of < t max >.In previous simulations, the user-defined ANM parameter δ is set to 10 −8 .In Figure 4, transient responses measured at location x 1 for the values δ = 10 −10 , 10 −8 , 10 −6 , 10 −5 , 10 −4 and 10 −3 are shown.For value of δ ≥ 10 −4 , differences are noticeable.The wave propagation nature is still present but numerical oscillations are attenuated.For smaller values, no difference is observed.Thus, minimizing the δ value ensures a very high level of precision for the high-order prediction.The price to pay is decrease of the range of validity t max and increase of ANM continuation steps.In Table 1, influence of this parameter δ is highlighted.The lower the value of δ is, the larger the number of ANM steps and the mean value < t max > are.For this example, the calculation with δ = 10 −5 required 162 continuation steps with an average value < t max >= 4.9521 × 10 −3 instead of 579 ANM steps and a mean value < t max >= 1.3829 × 10 −3 s with δ = 10 −10 .This behaviour is in agreement with past observations for the non-linear buckling analysis with ANM [26,43,44].In the same way, influence of the truncation order N is studied.In Figure 5, transient responses measured at location x 1 are computed with the values N = 5, 10, 15, 20 and 30.In this figure, for convenience, only continuation points are plotted for the responses computed with N > 5.The wave propagation is well described and there is no observable difference between numerical predictions.Only very minor gaps can be detected: attenuation of the wave amplitude can be noticed with the smallest values of N .Data are summarized in Table 2.It is noticed that the higher truncation order is, the smaller number of continuation steps are and the greater mean value of < t max > is.
In the light of these previous analysis, it can be postulated that the explicit ANM solver allow to compute solution possessing the same properties of the Newmark solution.This later is unconditionally stable, without numerical damping and introducing periodicity error (as analysed in [4] for the Newmark integration methods).Some improvement can be reached with a pertinent choice of ANM parameter values.Increasing the truncation order and decreasing value of the δ parameter allow to decrease the number of continuation steps.Consequently, it is observed an increase of the mean value of < t max > which is greater than classical time-step value, outperforming the stability criterion associated with explicit time-integration schemes.

Non-linear elastic cylindrical shell submitted to a dead force
The second example concerns the Snap-through of a cylindrical shell whose geometry is given in Figure 6a.The material properties are E =2.1 10 11 N/m 2 , ν = 0.25 and ρ = 1000 kg/m 3 .This example is the same that in Refs.[41] and [17].
The time evolution of the load at shell center is shown in Figure 6b.Due to symmetries, only a quarter of the shell is discretized with 16 shell elements based on geometrically exact element [45].It is a quadrilateral eight nodes element with six degrees of freedom per node.Details relative to this finite element can be found in Ref. [45,46].This element is well-adapted to the perturbation method as shown in Ref. [39] and is the same element which was used in Ref. [17,41] for this numerical example.
Time evolution of the vertical displacement of point M (see Figure 6a) previously computed in Ref. [17,41] is plotted in Figure 7.In Ref. [41], Generalized Energy-Momentum Method was used to solve this non-linear dynamics shell problem with a time step equal to 1 ms.In [17], authors associated the Generalized Energy-Momentum Method to homotopy and perturbation methods with a time-step value set to 0.5 ms.The two computed time evolutions of the displacement were exactly the same as reported in Figure 7.   6a) according to solutions given in Baguet (2011) [17] and Kuhl (1999) [41].
Results obtained with the time perturbation method are now considered.For the following computations with ANM, the truncation order N and the prescribed tolerance δ (Eq.( 6)) are equal to 20 and 10 −6 respectively.These values have been chosen because on the one hand they lead to a good compromise between the step length (i.e.t max ) and the computation time of the right hand side and on the other hand they ensure a good accuracy of the solution (see Ref. [47]).The first difficulty with this approach is to correctly represent time evolution of the load which presents a non-smooth evolution at instant 0.2 s.Value of the regularization parameter η, introduced in Eq. ( 9), has to be well chosen to give accurate results.So, in Figure 8a, time evolution of the modified load (Eq.( 9)) are plotted for four values of η between 10 −8 and 10 −11 .According to these results, it seems that a value of this regularization parameter lower than 10 −10 should be chosen to correctly represent load evolution and its singularity at instant 0.2 s.   9)) on the applied load and on the response curve.The solution obtained with the perturbation method is carried out with the following parameters : N =20 and δ = 10 −6 .This is confirmed by results shown in Figure 8b where time evolutions of the vertical displacement at the centre of the shell are plotted for seven values of the regularization parameter η.In this figure, it is observed that values of parameter η greater than 10 −8 do not permit to compute a solution with the correct dynamic buckling of the shell.Additionally, it seems that only values lower than 10 −10 allow to compute the right instant for which buckling is triggered.
So in the sequel of the numerical study, value of this regularization parameter is set to η = 10 −15 .This extra-small value has a great influence on the validity range (i.e.t f k ,g k max ) of power series expansion of the load expression.Nevertheless, as explained above, step length of the continuation method is governed by two ranges of validity t max : one for the unknowns q, noted t q max , and an other one for the load, noted t f k ,g k max .As this latter, for all the numerical tests carried out in this study, is greater than the one for the unknowns (even in the case where η = 10 −15 ), this small value of the regularization parameter does not affect numerical simulations.
In Figure 8b, it is noticed that the time for which buckling appears obtained with the proposed method is lower than the ones obtained in Ref. [17,41].This difference is also illustrated in Figure 9, where solution obtained with Abaqus explicit software [48] is compared with solutions computed in Ref. [41] and by the proposed method.Computation with Abaqus is carried out by using 16 shell elements S4R (a 4 nodes doubly curved shell with reduced integration, hourglass control and finite membrane strains).In Figure 9, one can remark that the three numerical methods give three different instants for the panel buckling.6a) to solutions given in Ref. [41], with Abaqus Explicit software [48] and with the proposed method with η = 10 −15 , N=20 and δ = 10 −6 .For all the computations, 16 finite elements are considered.
To explain this difference, several computations have been performed with Abaqus for a different number of shell elements (varying from 25 to 506).Results are exposed in Figure 10.From this plot, it is observed that solutions computed with Abaqus converge towards the one obtained with the perturbation method when number of finite elements is greater than 300.So, the ANM based solution seems to be the correct one.Let recall that this solution is obtained with only 16 shell elements which are also used in Ref. [41].Consequently, the difference between solutions obtained with perturbation method and those of Ref. [41] is only due to the time integration method.One of advantages of a perturbation method is that the computed prediction is analytical, in time in this study.This property can explain the fact that, even with a small number of elements, the asymptotic expansions follow the correct solution emerging from the model.On the contrary, classical time integration schemes iteratively proceed with small timesteps introducing errors at each iteration.This can explain the difference on the buckling trigger instant between Abaqus predictions and ANM based ones.
Moreover, as shown in Figure 11, the solution obtained with the perturbation method is not sensitive to the number of finite elements considered.Indeed, the buckling trigger instant is always the same even if the number of elements is five times greater that the initial mesh.
It is known that classical time-integration methods also suffer from mesh dependency: CFL condition for explicit scheme (see [4,42]) and spectral conditioning number associated to linear systems after spatial and time discretization for implicit scheme see, for example, analysis [9] for hyperbolic problem and [49] for parabolic problems).
To illustrate this drawback, value of the Abaqus time increment obtained with two different meshes is added on Figure 12a.This time increment is approximately equal to 4.5.10 −5 s and 2.8.10 −5 s for 100 and 306 elements, respectively.On this figure, time evolution of the parameter t max (Eq.6), is also plotted for three distinct meshing, with an increasing number of element ranging from 16 to 100.Let recall this parameter t max evaluates the range of validity of asymptotic expansions.It is observed that t max value is approximately the same for the three meshing.It coincides with the time-step value 4.5.10 −5 s used in the Abaqus simulations with 100 elements.6a) versus time according to the solutions given in Ref. [41], with Abaqus Explicit software [48] with several number of elements and with the proposed method with η = 10 −15 , N=20 and δ = 10 −6 .Even if the value of parameter t max evolves a lot during computations, it seems do not depend on the mesh size as the time-step does with classical time integration schemes.
Finally, it is proposed to evaluate influence of the truncation order N of asymptotic expansions to get the whole solution up to time t = T max = 0.3 s.So, in Figure 12b, the number of continuation steps is plotted with respect to the parameter N .It is also reported the total number of time   6)) and number of ANM continuation steps with respect to the truncation order N .Continuation is performed until the end of simulation set to T max = 0.3 s.ANM parameter values are: η = 10 −15 and δ = 10 −6 .iterations required with Abaqus to perform the same simulation, with meshing involving 306 and 506 finite elements.From these results, it is noticed that an increase of the truncation order leads to a decrease of the number of continuation steps.This is mainly due to the fact the validity range of the time, parameter t max , depends on the truncation order.Greater is this latter, greater is the parameter t max .This has been highlighted in the past when using asymptotic expansions in computational solid mechanics [26] or fluid dynamics [27].
In this dynamic buckling analysis, the perturbation method, based on power series of time, is compared with classical time integration methods available in Abaqus software.Computations with the perturbation method are carried out with 16 shell finite elements.To get same coherent results with Abaqus, it is required more than 300 finite elements (see Figure 10).This implies very small time-step values and large number of iterations.In that case, truncation order N has to be greater than 10 to compute the solution with approximately the same number of continuation steps.As with a classical explicit time integration scheme, only a single matrix triangulation is needed for all the continuation step of the perturbation method.Nevertheless, each continuation step requires a number of backward/forward substitutions equals to N -1.The higher the truncation order is, the higher the computational time is.This can lead to a huge increase of the computational time especially for models involving a high number of degrees of freedom.One remedy, already used in computational mechanics, relies on the use of Padé approximants [21,50].In these past studies, a reduction of about approximately half the number of continuation step was observed.The Padé approximants have been tested in this study but results were not convincing.The reason seems mainly due to numerical instability issues.These instabilities appear in the computation of the Euclidian norm of vectors q i which dramatically increases with the truncation order N and leads to numerical difficulties.

Conclusion
In this study, it was proposed a specific time-integration method in order to solve transient problems encountered in structural dynamics.
The new method relies on a high-order perturbation method in which the time variable is chosen as the perturbation parameter.Based on the use of a regularization parameter, it was shown how to recast equations in order to deal with arbitrary time variable external forcing.Two classical problems were addressed.The first problem was the wave propagation in an elastic rod submitted to a constant end force.The second problem was the dynamic buckling of a cylindrical panel.
For both problems, results obtained with the ANM were analysed in the light of available results.For the first problem, these later were obtained with classical time integration schemes (Pure explicit method, Newmark β-method and Central Difference method).For the second problem, comparisons were carried out with results available in the literature and computed with the commercial Abaqus software.
For both problems, it was shown that explicit ANM predictions were very coherent.For the first problem, the wave propagation nature of the solution is preserved.For same meshing, the ANM solution is close to the implicit Newmark based solution computed with the smallest time-step value.
For the second problem, influence of the regularization parameter in the external forcing was established and it was demonstrated that a very small value was necessary.ANM solution also described the non-linear dynamics of the panel.Trigger instant of the buckling and postbifurcated displacement were correctly obtained.Surprisingly, this study allowed for questioning some previous well-established results.More importantly, it was also demonstrated that the ANM solution is very similar to the Abaqus one which was based on a very refined meshing.
Finally, it was demonstrated that the ANM can be used as a pure explicit integration scheme without suffering from a stability criterion of CFL kind.Under the hypothesis of a correct choice of the continuation tolerance δ and regularization parameter η, it has seemed that the ANM allowed to intrinsically compute the right solution.To get solutions with the same level of precision, the alternative methods have been used with very fine meshing and very small timesteps.
However, in order the ANM based solver to become competitive, further developments are still necessary.As explained, the ANM scheme required high value of truncation order.But, the higher truncation order is, the higher computational time is, becoming harmful for models involving great number of unknowns.The use of a convergence accelerator could be a traditional solution.Unfortunately first attempt with Padé approximants failed.So, alternative solutions have to be sought.In order to enhance the range of validity of series and to decrease number of continuation steps, it could be envisaged use of the generalized factorial series or the Borel resummation technique [16,23,51].Alternatively, use of special functions, as generalized hypergeometric Meijer G-function [24], could also be tested.

Figure 1 .
Figure 1.Computed displacement u(t ) with ANM solver.Comparisons with the CD solver based solution and theoretical one [4, 42].Study of the transient rod response.

Figure 2 ,
solutions measured in x 1 = 0.25 m and obtained with the ANM solver are compared with the ones computed with the PE solver.Five different values of the time-step ∆t varying from 20% to 0.1% of the critical time-step ∆t cr are used to compute PE solutions.The lower time-step is, the lesser oscillatory the solutions are.Interestingly, the ANM solution tends to the PE one obtained with the lowest value of ∆t .The ANM solution also presents the smallest amplitude of oscillations.In this simulation, values of the ANM parameters are identical: δ = 10 −8 and N = 10 leading to 348 continuation steps and a mean value < t max >= 2.3044 × 10 −3 s for the range of validity.This value can be compared to the value ∆t = 5.000 × 10 −6 s used for the finest PE simulation.The difference is of three orders magnitude.

( a )
Displacement at location x 1 = 0.25 m.(b) Enlarged vision on the displacement at location x 1 = 0.25 m.

Figure 2 .
Figure 2. Computed displacement u(t ) with ANM solver.Comparisons with the closedform solution [4, 42] and PE solver based solutions computed with 5 different values of time-step ∆t .Study of the transient rod response.

( a )
Displacement at location x 1 = 0.25 m.(b) Enlarged vision on the displacement at location x 1 = 0.25 m.(c) Enlarged vision on the displacement at location x 1 = 0.25 m.

( a )
Displacement at location x 1 = 0.25 m.(b) Enlarged vision on the displacement at location x 1 = 0.25 m.

( a )
Displacement at location x 1 = 0.25 m.(b) Enlarged vision on the displacement at location x 1 = 0.25 m.

Figure 5 .
Figure 5. Computed displacement u(t ) with ANM solver.Parameter δ set to δ = 10 −5 .Influence of the truncation order N (for clarity purpose, for N > 5, only continuation points are plotted).Study of the transient rod response.

( a )
Geometry definition.(b) Time evolution of load.

Figure 6 .Figure 7 .
Figure 6.Geometry and load definitions of the cylindrical shell.Only 1/4 of the roof is discretized, with ρ=5m, L=2.5m, h=0.1m and θ = 30 • .The nodal force is applied at the center of the cylindrical shell.

Figure 8 .
Figure 8. Influence of the regularization parameter η (Eq.(9)) on the applied load and on the response curve.The solution obtained with the perturbation method is carried out with the following parameters : N =20 and δ = 10 −6 .

Figure 10 .
Figure10.Comparison of the evolution of the displacement (Point M, Figure6a) versus time according to the solutions given in Ref.[41], with Abaqus Explicit software[48] with several number of elements and with the proposed method with η = 10 −15 , N=20 and δ = 10 −6 .

Figure 11 .
Figure 11.Comparison of the evolution of the displacement (Point M, Figure 6a) versus time according to the solutions given with the Perturbation method with several number of elements with η = 10 −15 , N=20 and δ = 10 −6 .
Evolution of the parameter t max (Eq.(6)), N =20.Evolution of the number of steps.

Figure 12 .
Figure12.Evolutions of parameter t max (Eq.(6)) and number of ANM continuation steps with respect to the truncation order N .Continuation is performed until the end of simulation set to T max = 0.3 s.ANM parameter values are: η = 10 −15 and δ = 10 −6 .

Table 2 .
Influence of the truncation order N on the ANM precision with value δ = 10 −5 .Study of the transient rod response.