DOI 10.1007/s11071-011-0056-9 O R I G I N A L PA P E R
Dynamic buckling analysis of composite cylindrical shells
using a finite element based perturbation method
T. Rahman· E.L. Jansen · Z. Gürdal
Received: 19 September 2010 / Accepted: 13 April 2011 / Published online: 13 May 2011 © Springer Science+Business Media B.V. 2011
Abstract In this paper a finite element formulation of a reduction method for dynamic buckling analysis of imperfection-sensitive shell structures is presented. The reduction method makes use of a perturbation ap-proach, initially developed for static buckling and later extended to dynamic buckling analysis. The imple-mentation of a single-mode dynamic buckling analysis in a general purpose finite element code is described. The effectiveness of the approach is illustrated by ap-plication to the dynamic buckling of composite cylin-drical shells under axial and radial step loads. Results of the reduction method are compared with results available in the literature. The results are also com-pared with full model finite element explicit dynamic analysis, and a reasonable agreement is obtained.
T. Rahman
TNO DIANA BV, Delftechpark 19a, 2624 XJ Delft, The Netherlands
e-mail:T.Rahman@tnodiana.com E.L. Jansen (
)Institut für Statik und Dynamik, Leibniz Universität Hannover, Appelstrasse 9A, 30167 Hannover, Germany e-mail:E.Jansen@isd.uni-hannover.de
Z. Gürdal
Faculty of Aerospace Engineering, Delft University of Technology, Kluyverweg 1, 2629 HS Delft, The Netherlands
e-mail:Z.Gurdal@tudelft.nl
Keywords Dynamic buckling· Thin-walled structures· Finite elements · Perturbation method · Reduction method
1 Introduction
Shell structures are important structural components in various branches of engineering. These thin-walled structures are prone to buckling instabilities under static and dynamic compressive loading. Two main classes of dynamic stability can be distinguished [1,2], parametric excitation (vibration induced by pul-sating parametric loading), and dynamic buckling un-der step loading or impulsive loading. Cylindrical shells have attracted special attention because of their theoretical and practical significance. Investigations on parametric excitation of cylindrical shells in recent years can be found in e.g. [3–7], which also include literature overviews on this topic. The present work is concerned with the second class of dynamic-stability problems (dynamic buckling under step loading or im-pulsive loading), and more specifically, considers the case of dynamic buckling under step loading. Shell structures typically exhibit unstable post-buckling be-havior which makes them highly sensitive to geomet-ric or load imperfections. For imperfection-sensitive structures, under step loading buckling may occur at a load level that is lower than the corresponding static buckling load. Several criteria that can be used to esti-mate the critical load in the dynamic buckling case are
discussed in [2]. Following Budiansky and Roth [8], in the present work the critical load is defined as the step load level at which, when the amplitude of the step load is gradually increased in small increments, a distinct jump in the maximum response occurs.
The dynamic buckling of discrete systems under step loading was investigated by Kounadis [9]. Im-portant early studies on dynamic buckling of shells were done by Budiansky and Roth [8] and by Roth and Klosner [10]. Tamura and Babcock [11] included parametric excitation effects in their study of dynamic buckling of cylindrical shells under axial step loading. Finite element studies on the dynamic buckling anal-ysis of shells include the investigations by Saigal et al. [12], Ganapathi et al. [13], Yaffe and Abramovich [14], and Bisagni [15]. The latter study points out the important influence of the load duration on the dy-namic buckling load level.
The possibility to carry out a nonlinear transient analysis for shell structures is standardly available within finite element codes. However, finite element based transient analysis is computationally expensive and not suitable for the repeated runs necessary for a design and optimization process, and often it is dif-ficult to interpret the result and to judge its correct-ness. Reduction methods for general transient analysis of structures have therefore received considerable at-tention in recent years [16–19].
The reduction method proposed in the present pa-per is based on a finite element based pa-perturbation approach that has been successfully applied earlier by various researchers as the foundation of a reduc-tion method for static buckling analysis, e.g. [20–23]. In earlier studies using this perturbation approach for static buckling, often problems have been considered under the assumption of a linear pre-buckling state. For shell structures, in relevant cases the inclusion of pre-buckling nonlinearity is essential [24,25]. Koiter’s perturbation approach was presented using a func-tional notation by Budiansky [26] and by Byskov and Hutchinson [27].
Budiansky [26] also proposed an extension of the perturbation approach to dynamic buckling analysis. On the basis of a perturbation approach, Virgin and Chen [28] studied post-buckling dynamics of plates. Recently Teter [29] studied the dynamic buckling of long, prismatic columns consisting of rectangu-lar plates using the dynamic buckling approach in-troduced by Budiansky. The model was based on the
differential equations governing the dynamics of thin plates. Budiansky’s approach was applied by Schokker et al. [30] for the dynamic buckling analysis of com-posite cylindrical shells under hydrostatic pressure us-ing the p-version of the finite element method.
Objective of the present work is provide a reduction method for dynamic buckling analysis of general shell structures on the basis of the approach introduced by Budiansky. The development environment of the gen-eral purpose finite element code DIANA [31] is used as implementation platform. DIANA’s original imple-mentation of the perturbation approach for static post-buckling analysis and its recent extension to include the effect of pre-buckling nonlinearity by the present authors [25] is used as a basis, and will be further extended in the present work to include the effect of inertia. The effectiveness of the approach will be il-lustrated by application to the dynamic buckling of composite cylindrical shells under axial and radial step loads.
2 Perturbation method
In the present work the earlier implementation of a perturbation approach for static buckling analysis [25] will be extended to account for inertial forces in order to be able to analyse dynamic buckling problems. The extension is done along the lines proposed by Budi-ansky [26]. A short description of the perturbation ap-proach, both for static buckling and for dynamic buck-ling, will be given in the following sections. For clarity and conciseness, the approach will be presented for a linear pre-buckling state, and the modifications neces-sary in order to include a nonlinear pre-buckling state will be inserted in the final equations.
2.1 Static analysis
Introducing the generalized displacement u, strain , load f, and stress σ , it is assumed that the nonlinear strain-displacement relation and the linear elastic con-stitutive relation of a structure can be written as
= L1(u)+ 1
2L2(u) (1)
where L1 and H are linear functionals and L2 is a quadratic functional. The equilibrium equation in vari-ational form is written as
σ· δ − f · δu = 0 (3)
Here σ· δ and f · δu denote, respectively, the internal virtual work of the stress σ through the strain variation
δ, and the external virtual work of the load f through the displacement variation δu, both integrated over the entire structure. Further, if the bilinear functional L11 is defined such that
L2(u+ v) = L2(u)+ 2L11(u, v)+ L2(v) (4) then it follows from (1) that the first order strain vari-ation δ produced by δu can be written as
δ= L1(δu)+ L11(u, δu) (5)
It will also be assumed that the reciprocity relation
σi· j= σj· i (i, j= 1, 2 . . .) (6) holds. In this study proportional loading will be con-sidered, i.e. f= λf0, where λ is a normalized load parameter. Now the variables (u, , σ ) of the post-buckling equilibrium state can be expanded in the fol-lowing perturbation series:
u= λu0+ u1ξ+ u2ξ2+ u3ξ3+ · · ·
= λ0+ 1ξ+ 2ξ2+ 3ξ3+ · · · (7) σ= λσ0+ σ1ξ+ σ2ξ2+ σ3ξ3+ · · ·
where ξ is a normalized mode amplitude. The pertur-bation expansions in (7) are assumed to be asymptot-ically valid in the neighborhood of the pre-buckling equilibrium state corresponding to a bifurcation buck-ling load λ= λc. In the following it will be assumed that this buckling load corresponds to the lowest buck-ling load. Substituting (7) into (1), (2) and (3), tak-ing the limit ξ→ 0, and with some further manipula-tions one obtains the equamanipula-tions for the (lowest) buck-ling load λcand the corresponding buckling mode u1,
1= L1(u1) (8)
σ1= H (1) (9)
σ1· L11(δu)+ λcσ0· L11(u1, δu)= 0 (10)
Assuming that the pre-buckling solution is linear the following relation holds:
L11(u0, δu)= 0 (11) In addition, it will be assumed that (λ− λc)admits the asymptotic perturbation expansion
λ− λc= aλcξ+ bλcξ2+ · · · (12) If a plot of load parameter λ versus the mode ampli-tude ξ is made, then in view of (12) the a and b co-efficients indicate the slope and curvature of the post-buckling curve, respectively. In the present work we consider “symmetric” structures with a post-buckling slope a= 0, and typically a negative post-buckling curvature b < 0, indicating unstable post-buckling be-havior.
Inserting (12) together with (7) into (1), (2), and (3), and equating the coefficients of ξ2 with the as-sumption of a= 0 (“symmetric” structures) one ob-tains the equations for the determination of the second order mode u2,
2= L1(u2)+ 1 2L2(u1) (13) σ2= H (2) (14) σ2· L11(δu)+ λcσ0· L11(u2, δu) + σ1· L11(u1, δu)= 0 (15)
For the second order mode u2, the following orthogo-nality condition is imposed:
σ0· L11(u1,u2)= 0 (16)
In order to obtain the expression for the b coefficient we set δu= u1in (10) and (15) and make use of the reciprocity relation (6). This gives
b=2σ1· L11(u1,u2)+ σ2· L2(u1)
σ1· 1
(17)
The behavior of the imperfect structure will be related to the properties of the perfect structure. If the initial geometric imperfection is denoted by ˆu = ¯ξu1, where ¯ξ is the normalized imperfection amplitude and u1is a geometric imperfection pattern in the shape of the first buckling mode, then the strain-displacement equation, (1), can be modified to
= L1(u)+ 1
Further, the asymptotic expansion as defined by (12) is modified to
ξ(λ− λc)= aλcξ2+ bλcξ3− λ¯ξ + · · · (19) Equation (19) can be rearranged as
1− λ λc ξ+ aξ2+ bξ3= λ λc ¯ξ (20)
Equation (20) can be seen as the reduced-order model for the static analysis.
2.2 Dynamic analysis
In order to analyse dynamic buckling, in this section the perturbation method discussed in the previous sec-tion will be extended to account for inertial forces. In-troducing inertial forces, (3) takes the form
σ· δ − f · δu + M ∂2u ∂t2 · δu = 0 (21) where M(∂2u ∂t2), which is linear in ∂2u ∂t2, represents the inertial loading. The following reciprocity relation, cf. (6) for the static analysis, holds for the dynamic case:
M(u)· v = M(v) · u (22)
The dynamic loading is assumed to have the form f=
λF (t )f0, where the time variation F (t) is normalized such that its maximum value is unity. The typical case that will be considered in the present paper is the case of a step loading, with F (t) corresponding to the unit step function. The dynamic counterpart of (7) can now be written as
u= λF (t)u0+ u1ξ(t )+ u2ξ2(t )+ u3ξ3(t )+ · · · = λF (t)0+ 1ξ(t )+ 2ξ2(t )+ 3ξ3(t )+ · · · (23) σ= λF (t)σ0+ σ1ξ(t )+ σ2ξ2(t )+ σ3ξ3(t )+ · · · By repeating the same procedure as for the static case, and neglecting the inertial forces associated with the pre-buckling displacements, the dynamic counterpart of (20) is obtained, 1 ω12 ξ(t )+ 1−λF (t ) λc ξ(t )+ aξ2(t )+ bξ3(t ) = λF (t ) λc ¯ξ (24)
where a prime represent differentiation with respect to time and where ω21is defined by
ω21= σ1· 1 M(u1)· u1
(25)
If u1happens to be a natural vibration mode then ω1 is its natural radial frequency, otherwise ω21 can be interpretated as a Rayleigh quotient corresponding to the radial frequency squared, based on the buckling mode u1.
For clarity and conciseness the analysis so far has been outlined for a linear pre-buckling state. In order to account for pre-buckling nonlinearity, (17) should be modified to b= −2σ1· L11(u1,u2)+ σ2· L2(u1) λcΔˆ (26) where ˆ Δ= 2σ1· L11(˙uc,u1)+ ˙σc· L2(u1) (27) Here the subscript ( )cdenotes pre-buckling quantities evaluated at λ= λcand the dots in this case represent differentiation with respect to λ, ˙( )=∂λ∂ ( ). It is noted that as compared with (17), the quantity σ1· 1is re-placed by −λcΔˆ in (26) to account for pre-buckling nonlinearity. Similarly (25) is adjusted to
ω21= −λcΔˆ M(u1)· u1
(28)
Finally, in order to include the effect of pre-buckling nonlinearity, (24), the reduced-order model for the dy-namic state, is modified to
1 ω21 ξ(t )+ 1−λF (t ) λc ξ(t )+ aξ2(t )+ bξ3(t ) = α ¯ξ − β 1−λF (t ) λc ¯ξ (29)
where the coefficients α and β are known as the first and second imperfection form factors, respectively [25], and are given by
α= (1/λcΔ)ˆ
σ1· L11(ˆu, uc)+ σc· L11(ˆu, u1)
(30)
β= (1/ ˆΔ)σ1· L11(ˆu, ˙uc)+ ˙σc· L11(ˆu, u1) + HL11(˙uc,u1) · L11(ˆu, uc) − αλc σ1· L11(¨uc,u1)+ (1/2) ¨σc· L11(u1,u1)
+ HL11(˙uc,u1)
· L11(˙uc,u1)
(31) Equation (29) can be regarded as the reduced-order model for the dynamic analysis. In the case of a lin-ear pre-buckling state both α and β are unity and (29) becomes identical to (24). In the case of a nonlinear pre-buckling state the bifurcation buckling analysis is carried out at a fundamental state close to the bifurca-tion buckling load.
On the basis of the previous discussion, the
reduc-tion method for dynamic buckling analysis proposed
in the present work is described by the following step-by-step procedure:
– Establish the reduced-order model for static buck-ling analysis, (20):
1. Computation of the pre-buckling state u0. 2. Computation of the buckling load λcand the
cor-responding buckling mode u1.
3. Computation of the second order mode u2. 4. Computation of the b coefficient (a= 0 for the
“symmetric” structures considered in the present study).
– Establish the reduced-order model for dynamic buckling analysis, (24).
– Carry out the reduced analysis:
1. Computation of the normalized mode amplitude
ξ(t ) as function of time by solving (24) by nu-merical time integration for specified load levels, increasing the amplitude of the applied load with small increments.
2. Recovering the displacement, stress, and strain by substituting the given loads and calculated amplitudes in (23).
3. Identification of the dynamic buckling load level
λ= λd, the level at which the maximum of ξ(t) shows a sharp rise.
3 Finite element implementation
The finite element implementation of the perturbation approach discussed in the previous section has been done in the development environment of the general purpose finite element code DIANA using a layered, eight-node quadrilateral iso-parametric curved shell element, the CQ40L element. A description of this el-ement is available in the DIANA manual [31]. In or-der to use this element in the perturbation approach
described, a modification of the element formulation is not required. However, it is necessary to construct BN L, the nonlinear part of the strain-displacement ma-trix. In finite element notation the strain-displacement relation is given as
= BLq+1
2BN L(q)q (32)
where BLand BN L, defined at each integration point, correspond to the L1 and L2 functionals in (1), re-spectively, and q is the vector of nodal displacements at each element corresponding to the displacement field u from the functional notation. In a similar way
L11(u, v) is translated to finite element notation as BN L(q1)q2, where the nodal displacement vectors q1 and q2 correspond to the displacement fields u and v in the functional notation, respectively. In the follow-ing, we will discuss the finite element implementation briefly. Details of the implementation are available in earlier work of the present authors [25].
The static pre-bucking state q0is obtained from a linear static analysis,
K0q0= f0 (33)
where K0 is the global stiffness matrix at the unde-formed state of the structure.
The buckling problem, in functional notation given by (8), (9), and (10), corresponds to the following eigenvalue problem:
[K0+ λcKG]q1= 0 (34)
where KG is the geometric stiffness matrix. Solution of (34) gives the lowest bifurcation buckling load λc and the corresponding buckling mode q1. In the case of follower loads the loading direction can change during the deformation, and the follower load gives rise to an additional stiffness term, known as the load stiffness. For instance, for fluid pressure loading, the pressure remains normal to the shell surface. In the present study, the contribution of the load stiffness is accounted for as discussed in [32]. The eigenvalue problem in (34) is then modified to
K0+ λc(KG+ KL)
q1= 0 (35)
The translation to finite element notation of the sec-ond order state problem, in functional notation repre-sented by (13), (14), and (15) gives
[K0+ φλcKG]q2= g (36)
while the orthogonality constraint given by (16) trans-lates to
qT1KGq2= 0 (37)
In (36) the right-hand side force vector g is given by
g= −1 2 BTLHBN L(q1)q1+ 2BTN L(q1)HBLq1 (38)
and φ is a factor such that φ≈ 1, but φ < 1. Usually
φ= 0.99 is applied, if φ = 1 then (36) becomes sin-gular and cannot be solved. When pre-buckling non-linearity is considered, as discussed in [25], (36) and (37) include additional terms.
The b coefficient is obtained from (17) or, if pre-buckling nonlinearity is considered, from (26), by computing the post-buckling stresses and strains at each integration point, and summing them up over all the elements in the structure.
4 Results and discussion
In this section, results of the reduction method for the dynamic buckling analysis are presented for unstiff-ened and ring-stiffunstiff-ened composite cylindrical shells. Firstly, a study of the dynamic buckling of compos-ite cylindrical shell structures under external pressure, based on a static buckling mode obtained using a linear pre-buckling state, is presented. Additionally results of the dynamic buckling analysis of a composite cylindri-cal shell under axial loading, based on a static buckling mode obtained using a nonlinear pre-buckling state, are presented. In two additional cases a comparison is made with earlier work from other investigators on the dynamic buckling of an unstiffened and a ring-stiffened composite shell. The dynamic buckling loads in the reduced model are computed in the present study by numerically solving the reduced (24) by a standard Runge-Kutta scheme of the 4th order. In all the numer-ical examples the buckling mode is normalized such that the maximum out-of-plane displacement is equal to the shell wall thickness.
4.1 Anisotropic composite shell under external pressure step loading
The first numerical example for composite shells is the dynamic buckling of Booton’s anisotropic cylin-drical shell [35] under external pressure. This shell was also used earlier in static stability investigations [33]. The load has been applied as a pressure with a direction that remains radial throughout the defor-mation process (as opposed to fluid pressure load-ing, where the loading direction remains normal to the shell surface). Referring to the discussion earlier on the eigenvalue problems in (34) and (35), the load stiffness term is therefore not included in the buck-ling analysis. The edges have classical ‘simply sup-ported’ boundary conditions (SS-3) [33]. In order to apply SS-3 boundary conditions in the finite element model, both radial and circumferential displacements at one end of the cylinder are restrained and also one node is fixed in the axial direction to suppress the rigid body motion in that direction. At the other end the radial displacements are restrained, and also the circumferential displacements relative to the dis-placement of the end plane are restrained, by specify-ing kinematic constraints in the finite element model. The material and geometric properties of the shell are given in Tables1and2. The results of the perturbation method for the static buckling analysis obtained from the finite element approach are compared with results from a semi-analytical method based on Donnell-type governing differential equations, which accurately ac-counts for the effect of boundary conditions at the shell edges [33]. In Table3, the first (i.e. lowest) bifurca-tion buckling load and the corresponding b coefficient
Table 1 Material properties of Booton’s shell
Young’s modulus, E1 5.83× 106psi
Young’s modulus, E2 2.42× 106psi
Poisson’s ratio, ν12 0.36
Shear modulus, G12 6.68× 105psi
Density, ρ 2.6× 10−4lb.s2/ in4
Table 2 Geometric properties of Booton’s shell
Radius, R 2.67 in
Thickness, t 0.0267 in
Length, L 3.776 in
Table 3 Normalized lowest buckling load pcr/pcl and b
co-efficient under external dead load pressure, and lowest natu-ral frequency f1for Booton’s anisotropic shell, SS-3 boundary
conditions. Comparison between present finite element results and semi-analytical results. Number of circumferential waves between parentheses
Normalized buckling load b coefficient Natural frequency (Hz)
Semi-analytical [33] Present Semi-analytical [33] Present Semi-analytical [34] Present 5.6569× 10−1(8) 5.7× 10−1(8) −6.3289 × 10−2 −6.2664 × 10−2 1.1595× 103(6) 1.129× 103(6)
Fig. 1 Buckling mode and second order mode of Booton’s shell under static external dead load pressure
are compared between the present study and the semi-analytical tool ANILISA [33]. The dead load buckling pressure p (lb/ in2) and its critical value pcrare nor-malized with respect to a reference buckling pressure
pcl = (E1t2)/(cR2), where c=
3(1− ν122). A rea-sonable agreement is obtained, a difference of approx-imately one percent both for the buckling pressure and for the b-factor. As an additional verification case for the finite element model used, the first (i.e. lowest) natural frequency f1 obtained with the present finite element model is compared with the results of a semi-analytical method for vibrations of composite cylindri-cal shells [34]. In Fig.1the first bifurcation buckling mode and the corresponding second order mode are depicted by means of deformed mesh plots.
By establishing the characteristics of the static be-havior of the perfect structure, the parameters for the reduded order model for the dynamic buckling analy-sis, (24), have been determined. The dynamic buckling load under step loading for the imperfect structure has been evaluated using this equation. The imperfection used in the analysis has the shape of the first buck-ling mode and its maximum out-of-plane amplitude is taken as 10% of the shell thickness. With the reduced analysis the dynamic response of the structure in terms of the modal amplitude ξ(t) is computed for increas-ing load levels in order to identify the load level at
which the maximum of ξ(t) shows a sharp rise or be-comes unbounded. This load level is considered to be the dynamic buckling load λd. With the reduced anal-ysis the calculation time is very low as compared with a dynamic analysis of the full finite element model be-cause each response analysis in the reduced analysis corresponds to solving only one ordinary differential equation, (24). In Fig.2 the response of the structure is shown in terms of the out-of-plane (radial) deflec-tion of one of the nodes at the cylinder’s midlength, just below and just above the dynamic buckling load level. The dashed line shows the structural response as a function of time at load level λ/λc= 0.82, computed using the reduced analysis. At load level λ/λc= 0.83 the reduced analysis (dash-dotted line) yields an un-bounded response, and the dynamic buckling load lies between λ/λc= 0.82 and λ/λc= 0.83. A dynamic analysis under step loading was also carried out using a full model explicit dynamic analysis with ABAQUS [36], using a laminated shell element similar to the ele-ment used in DIANA. It can be observed that with full model analysis the maximum response takes a sharp rise at the same load level as in the reduced analy-sis. The reduced-order model, corresponding to a one-degree-of-freedom system with a cubic nonlinearity, is not able to represent stabilizing effects for larger amplitudes, and its response tends to infinity, whereas
Table 4 Normalized lowest buckling load N0cr/Ncland b
co-efficient under static axial loading, and lowest natural frequency
f1 for Booton’s anisotropic shell, SS-4 boundary conditions.
Comparison between present finite element results and semi-analytical results. Number of circumferential waves between parentheses
Normalized buckling load b coefficient Natural frequency (Hz)
Semi-analytical [33] Present Semi-analytical [33] Present Semi-analytical [34] Present 0.4001 (8) 0.3942 (8) −0.3548 −0.3658 1.4224× 103(7) 1.3936× 103(7)
Fig. 2 Comparison of modal response between full model ex-plicit dynamic analysis and reduced-order model analysis, for load levels just below (λ/λc= 0.82) and just above the dynamic
buckling load (λ/λc= 0.83). Booton’s shell under external dead
load pressure step loading. Imperfection amplitude ¯ξ= 0.1
the nonlinearities in the multi-degree-of-freedom full model will limit the deflections at finite amplitudes. The characteristics of the response near the dynamic buckling load, for times that are sufficiently small and for amplitudes that are sufficiently small are captured by the reduced-order model, and the prediction of the dynamic buckling load with the reduced analysis is in good agreement with the result of the full model analysis (in the present case the buckling load of the reduced-order model and the full model occur at the same load level).
4.2 Anisotropic composite shell under axial step loading
In this example the same cylindrical shell as in the previous example is considered (Booton’s anisotropic shell), but instead of external pressure, an axial load
is applied, evenly distributed along the circumference at the shell edge. Another type of simply supported boundary condition is used, namely SS4. The SS-4 condition is similar to SS-3, but in the SS-4 case also the axial nodal displacements, relative to the edge plane, are restrained at both shell edges by specifying kinematic constraints in the finite element model.
In Table4, as for the previous example, the results of the perturbation method for the static buckling anal-ysis (the first bifurcation buckling loads and the b co-efficients) obtained from the finite element approach are compared with results from the semi-analytical method from [33]. The compressive buckling stress re-sultant N0(lb/ in) and its critical value N0cr are
nor-malized with respect to the reference buckling stress resultant Ncl= (E1t2)/(cR), where c=
3(1− ν122). Also in this case a reasonable agreement is obtained, a difference of about 1.5 percent for the buckling load, and of about 3 percent for the b-factors. As an addi-tional verification case for the finite element model used, the first natural frequency f1obtained with the present finite element model is again compared with the results from [34]. In Fig. 3 the first bifurcation buckling mode and the corresponding second order mode are shown by means of deformed mesh plots.
The dynamic buckling load under step loading for the imperfect structure has been evaluated using the reduced-order model, (24). The imperfection used in the analysis has, as in the previous example, the shape of the first buckling mode and its maximum out-of-plane amplitude is taken as 10% of the shell thickness. With the reduced analysis, dynamic buckling was de-tected around the load level λ/λc= 0.768, while with full model explicit dynamic analysis dynamic buck-ling was found at about λ/λc= 0.718. While the load levels from the full model analysis and the reduced analysis in the case of external pressure showed a very good agreement, under axial compressive loading a clear difference in the load levels from the two anal-yses can be observed (a difference of about 7 percent).
Fig. 3 Pre-buckling shape, buckling mode, and second order mode of Booton’s shell under static axial loading
The buckling behavior of cylindrical shells under axial loading is considerably more complicated than the be-havior under external pressure, since the pre-buckling path is nonlinear and furthermore the buckling loads can be very closely spaced. In the present work, a single-mode study has been performed. A multi-mode reduced analysis could improve the reduced analysis results in the present example of an axially loaded shell.
4.3 Composite shell under external fluid pressure step loading
The present example provides a comparison with the work of Schokker et al. [30]. They investigated unstiff-ened and ring-stiffunstiff-ened composite cylindrical shells under external pressure step loading using also the approach followed in the present study (a reduction method based on Budiansky’s approach), but within the framework of the p-version of the finite element method. In this example an unstiffened composite cylindrical shell under external pressure with simply supported (SS-3) boundary conditions will be anal-ysed. The material and geometric properties of this cylinder are given in Tables5and6. Schokker et al. ap-plied the external pressure as a fluid pressure loading. Therefore in the present analysis fluid pressure loading was used (contrary to the earlier “dead load” pressure used in the first analysis example). Now in the bifur-cation buckling (eigenvalue) problem the contribution of the load stiffness matrix is accounted for by using (35).
In Table 7 for fluid pressure, the buckling load is compared with results obtained using the general purpose finite element code ABAQUS [36], and a
Table 5 Material properties of Schokker’s unstiffened compos-ite shell
Young’s modulus, E11 16.0× 106psi
Young’s modulus, E22 1.48× 106psi
Young’s modulus, E33 1.48× 106psi
Shear modulus, G12 0.76× 106psi
Shear modulus, G23 0.51× 106psi
Shear modulus, G13 0.76× 106psi
Poisson’s ratio, ν12 0.33
Poisson’s ratio, ν23 0.45
Poisson’s ratio, ν13 0.33
Density, ρ 1.47× 10−4lb.s2/ in4
Table 6 Geometric properties of Schokker’s unstiffened com-posite shell
Radius, R 50 in
Ratio of radius to thickness, R/t 50 Ratio of length to radius, L/R 4
Stacking sequence [90/0]S
good agreement is observed (the buckling loads ob-tained are virtually equal). The pressure q (lb/ in2) and its critical value qcrhave been normalized by us-ing (q/E11)× 105. As a further verification of the model used, in this table also the lowest natural radial frequencies are compared with the results from [30] and with results obtained using ABAQUS. The lowest natural radial frequency ω1 (rad/s) is normalized by using ω1L2(ρ/E2t2)1/2. Also for the natural frequen-cies a good agreement is obtained (the frequenfrequen-cies ob-tained are virtually equal).
In Table8the normalized static and dynamic buck-ling pressures under step loading are compared for
Table 7 Normalized lowest buckling loads under static exter-nal fluid pressure (qcr/E11)×105and normalized natural radial
frequencies ω1L2(ρ/E2t2)1/2for Schokker’s unstiffened
com-posite shell, SS-3 boundary conditions. Comparison between present results, ABAQUS, and [30]. Number of circumferential waves between parentheses
Normalized buckling load Normalized natural frequency
Present ABAQUS Ref. [30] Present ABAQUS
0.8504 (3) 0.8503 (3) 144.7 (3) 144.89 (3) 144.85 (3)
Table 8 Comparison of normalized buckling load static and dynamic buckling load (qcr/E11)× 105for
Schokker’s unstiffened composite shell under external fluid pressure for varying imperfection amplitudes
¯ξ Normalized buckling load Normalized buckling load
Static Dynamic
Ref. [30] Present Ref. [30] Present
0.05 0.774 0.8025 0.756 0.7906
0.25 0.720 0.7194 0.655 0.6899
0.5 0.678 0.6550 0.579 0.6144
varying imperfection amplitude. In the presence of im-perfections, the static bifurcation buckling load which occurs for a perfect shell is replaced by buckling at a limit point at load level λs, which is evaluated us-ing (7). The dynamic buckling load is computed by solving (23) by means of numerical integration in time for increasing loading amplitude. A reasonable agree-ment between the results of the two investigations can be observed (for an imperfection amplitude ¯ξ of 0.5, a difference of about 3.5 percent for the static case, and of about 6 percent for the dynamic case).
In Fig.4the modal amplitude ξ is plotted against time at a load level far below the dynamic buckling load level, at λ/λc= 0.356, and at load levels just be-low and just above the dynamic buckling load level (at
λ/λc= 0.722 and λ/λc= 0.723, respectively). The re-sults correspond well with a similar plot shown in the paper by Schokker et al. [30], where dynamic buckling was detected between λ/λc=0.713 and 0.714. 4.4 Ring-stiffened composite shell under external
fluid pressure step loading
The present example provides a further comparison with the work of Schokker et al. [30]. In the present example, a ring-stiffened composite shell is consid-ered made of the same material as in the previous example. The boundary conditions and loading con-ditions are also identical. Geometric properties of the ring-stiffened shell are given in Table9. In Fig.5the
Fig. 4 Modal response for three characteristic values of step load amplitude λ/λc: Far below (λ/λc = 0.356), just below
(λ/λc = 0.722) and just above the dynamic buckling load
(λ/λc = 0.723). Schokker’s unstiffened composite shell
un-der external fluid pressure step loading. Imperfection amplitude
¯ξ = 0.5
first buckling mode and the corresponding second or-der mode are shown. The first buckling mode cor-responds to a global pattern with three circumferen-tial waves (n= 3), involving buckling of the ring-stiffeners. In Table 10 the static and dynamic buck-ling fluid pressures for varying imperfection ampli-tudes are reported. The pressure q (lb/ in2) and its critical value qcrhave again been normalized by using
Fig. 5 Buckling mode and second order mode of Schokker’s ring-stiffened composite shell under static external fluid pressure
Table 9 Geometric properties of Schokker’s ring-stiffened composite shell
Number of stiffeners, Ns 10
Ratio of depth of stiffener to stiffener thickness, d/ts 2.711
Ratio of distance between stiffeners to cylinder 0.22 radius, Sp/R
Radius, R 50
Ratio of cylinder length to radius, L/R 2.62 Ratio of cylinder radius to wall thickness, R/t 100 Ratio of stiffener thickness to cylinder wall 2.0 thickness, ts/t
Table 10 Comparison of normalized static buckling load and dynamic buckling load under step loading (qcr/E11)× 105for
Schokker’s ring-stiffened composite shell under external fluid pressure
¯ξ Normalized buckling load Normalized buckling load
Static Dynamic
Ref. [30] Present Ref. [30] Present
0.1t 1.80 1.8337 1.72 1.8103
0.3t 1.71 1.7394 1.56 1.6946
0.5t 1.65 1.6700 1.45 1.6107
1.0t 1.54 1.5400 1.25 1.4569
(q/E11)× 105. Comparing the present results with the results from Schokker et al. [30], a reasonable agree-ment for the static buckling loads is observed (for an imperfection amplitude ¯ξ of 0.1, a difference of about 2 percent, and smaller differences for the other im-perfection amplitudes investigated), but the dynamic buckling loads obtained in the present study are clearly higher than the buckling loads from [30] (for an imfection amplitude ¯ξ of 0.1, a difference of about 5 per-cent, and for an imperfection amplitude ¯ξ of 0.1, a dif-ference of about 17 percent). In the current implemen-tation, the dynamic buckling loads are computed by numerically solving the reduced-order model, (24), by a standard Runge–Kutta scheme of the 4th order. How-ever, to shed some more light on this discrepancy, it is noted that in the case of step loading this equation can be solved analytically and Budiansky [26] presented for this case the following relation between the static and dynamic buckling load (for a= 0 and b < 0): 1− (λD/λS)(λS/λC) 1− (λS/λC) 3 2 =√2 λD λS (39)
where λC is the static bifurcation buckling load, λS is the static limit-point buckling load, and λD is the dynamic buckling load. For the present example λC= 1.9275, and the results obtained through solution of (39) correspond closely to the dynamic buckling loads
of the present study reported in Table10, which seems to indicate that the current results are reasonable.
5 Concluding remarks
A perturbation approach has been used as the basis of a reduction method for the finite element dynamic buckling analysis of general shell structures. The im-plementation has been done within a general pur-pose finite element code. Single-mode dynamic buck-ing analysis was carried out for unstiffened and rbuck-ing- ring-stiffened composite cylindrical shells.
The main advantage of the approach presented lies in the quick dynamic buckling load estimates that can be obtained. Instead of carrying out a transient analysis with a (large) finite element model, one solves a single ordinary differential equation through numerical time integration.
The approach gives reasonable estimates of the dy-namic buckling load under step loading. The approach is based on the assumption that the dynamic buckling mode is identical to the static buckling mode. The ap-proach can therefore be expected to give appropriate estimates for the dynamic buckling load in cases in which the dynamic buckling mode resembles the static buckling mode. An extension of the approach to a multi-mode perturbation analysis will improve the dy-namic buckling estimates in more complicated cases in which modal interaction in the structural response is essential.
Acknowledgements This research was done at the Aerospace Structures Group, Faculty of Aerospace Engineering, Delft Uni-versity of Technology. The research was financially supported by the Dutch Technology Foundation STW, applied science di-vision of NWO and the Technology Program of the Ministry of Economic Affairs (STW project 06613). The authors would like to acknowledge the help of TNO DIANA BV with the finite element implementation.
References
1. Lindberg, H.E., Florence, A.L.: Dynamic Pulse Buckling. Nijhoff, Dordrecht (1987)
2. Simitses, G.J.: Dynamic Stability of Suddenly Loaded Structures. Springer, Berlin (1990)
3. Pellicano, F., Amabili, M.: Dynamic instability and chaos of empty and fluid-filled circular cylindrical shells under periodic axial loads. J. Sound Vib. 293, 227–252 (2006) 4. Amabili, M.: Nonlinear Vibrations and Stability of Shells
and Plates. Cambridge University Press, Cambridge (2008)
5. Gonçalves, P.B., Del Prado, Z.J.G.N.: Low-dimensional Galerkin models for nonlinear vibration and instability analysis of cylindrical shells. Nonlinear Dyn. 41, 129–145 (2005)
6. Mallon, N.J., Fey, R.H.B., Nijmeijer, H.: Dynamic stabil-ity of a thin cylindrical shell with top mass subjected to harmonic base-acceleration. Int. J. Solids Struct. 45, 1587– 1613 (2008)
7. Jansen, E.L.: Dynamic stability problems of anisotropic cylindrical shells via a simplified analysis. Nonlinear Dyn. 39, 349–367 (2005)
8. Budiansky, B., Roth, R.S.: Axisymmetric dynamic buck-ling of clamped shallow spherical shells. Technical report TN D-1510, NASA, 1962
9. Kounadis, A.: Nonlinear dynamic buckling of discrete dissipative or nondissipative systems under step loading. AIAA J. 29, 280–289 (1991)
10. Roth, R., Klosner, J.: Nonlinear response of cylindrical shells subjected to dynamic axial loads. AIAA J. 2, 1788– 1794 (1964)
11. Tamura, Y.S., Babcock, C.D.: Dynamic stability of cylin-drical shells under step loading. J. Appl. Mech. 42, 190– 194 (1975)
12. Saigal, S., Yang, T.Y., Kapania, R.K.: Dynamic buckling of imperfection sensitive shell structures. J. Aircr. 24(10), 718–725 (1987)
13. Ganapathi, M., Gupta, S.S., Patel, B.P.: Nonlinear axisym-metric dynamic buckling of laminated angle-ply composite spherical caps. Compos. Struct. 59, 89–97 (2003) 14. Yaffe, R., Abramovich, H.: Dynamic buckling of
cylindri-cal stringer stiffened shells. Comput. Struct. 81, 1031–1039 (2003)
15. Bisagni, C.: Dynamic buckling of fiber composite shells un-der impulsive axial compression. Thin-Walled Struct. 43, 499–514 (2005)
16. Spottswood, S.M., Hollkamp, J.J., Eason, T.G.: On the use of reduced-order models for a shallow curved beam under combined loading. In: Proceedings of the 49th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dy-namics and Materials Conference, Schaumburg, IL (2008) 17. Gordon, R.W., Hollkamp, J.J.: Reduced-order
model-ing of the random response of curved beams us-ing implicit condensation. In: Proceedus-ings of the 47th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dy-namics and Materials Conference, Newport, RI (2006) 18. Przekop, A., Rizzi, S.A.: Nonlinear reduced order
ran-dom response analysis of structures with shallow curvature. AIAA J. 44(8), 1767–1778 (2006)
19. Tiso, P., Jansen, E.L., Abdalla, M.M.: A reduction method for finite element nonlinear dynamic analysis of shells. In: Proceedings of the 47th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics and Materials Conference, Newport, RI (2006)
20. Olesen, J.F., Byskov, E.: Accurate determination of asymp-totic postbuckling stresses by the finite element method. Comput. Struct. 15, 157–163 (1982)
21. Peek, R., Kheyrkhahan, M.: Postbuckling behavior and im-perfection sensitivity of elastic structures by the Lyapunov– Schmidt–Koiter approach. Comput. Methods Appl. Mech. Eng. 108, 261–279 (1993)
22. Garcea, G., Casciaro, R., Attanasio, G., Giordano, F.: Per-turbation approach to elastic post-buckling analysis. Com-put. Struct. 66, 585–595 (1998)
23. Menken, C.M., Schreppers, G.M.A., Groot, W.J., Petterson, R.: Analyzing buckling mode interactions in elastic struc-tures using an asymptotic approach; theory and experiment. Comput. Struct. 64, 473–480 (1997)
24. Kheyrkhahan, M., Peek, R.: Postbuckling analysis and im-perfection sensitivity of general shells by the finite element method. Int. J. Solids Struct. 36, 2641–2681 (1999) 25. Rahman, T., Jansen, E.L.: Finite element based coupled
mode initial post-buckling analysis of a composite cylin-drical shell. Thin-Walled Struct. 48, 25–32 (2010) 26. Budiansky, B.: Dynamic buckling of elastic structures:
Cri-teria and estimates. In: Dynamic Stability of Structures. Pergamon, Elmsford (1967)
27. Byskov, E., Hutchinson, J.W.: Mode interaction in axially stiffened cylindrical shells. AIAA J. 15, 941–948 (1977) 28. Chen, H., Virgin, L.N.: Finite element analysis of
post-buckling dynamics in plates—Part I: An asymptotic ap-proach. Int. J. Solids Struct. 43, 3983–4007 (2006)
29. Teter, A.: Dynamic, multimode buckling of thin-walled columns subjected to in-plane pulse loading. Int. J. Non-Linear Mech. 45, 207–218 (2010)
30. Schokker, A., Sridharan, S., Kasagi, A.: Dynamic buckling of composite shells. Comput. Struct. 59, 43–53 (1996) 31. DIANA User’s Manual—Release 9.2. Delft, The
Nether-lands (2007)
32. Schweizerhof, K., Ramm, E.: Displacement dependent pressure loads in nonlinear finite element analyses. Com-put. Struct. 18(6), 1099–1114 (1984)
33. Arbocz, J., Hol, J.M.A.M.: Koiter’s stability theory in a computer-aided engineering (CAE) environment. Int. J. Solids Struct. 26, 945–975 (1990)
34. Jansen, E.L.: A perturbation method for nonlinear vibra-tions of imperfect structures: Application to cylindrical shell vibrations. Int. J. Solids Struct. 45, 1124–1145 (2008) 35. Booton, M.: Buckling of imperfect anisotropic cylinders
under combined loading. Technical report 203, UTIAS (1976)
