Nonlinear Plate Theory for Postbuckling Behaviour of Thin-Walled Structures Under Static and Dynamic Load

1. Introduction

A thin plate or thin-walled constructions are used in the sports industry, automotive, aerospace and civil engineering. As an example of such structural elements snowboard, skis, poles may be mentioned, as well as all kinds of crane girders, structural components of automobiles (car body sheathing or all longitudinal members), aircraft fuselages and wings, supporting structures of the walls and roofs of large halls and warehouses. All the above structures, as well as many others which can be regarded as a thin, exhaust carrying capacity not by exceeding the allowable stresses but by the stability loss. Therefore, not only critical load but also the postbuckling behaviour of thin-walled structures subjected to static and dynamic load is essential knowledge for designers. The use of more accurate mathematical models allows to explore the phenomena occurring after the loss of stability and to describe more precisely their behaviour. Engineers and designers need guidelines to construct as well as quick and easy software to use for analyse the behaviour of thin-walled structures. Therefore, the author of this chapter decided to explore this issue, propose a mathematical model and the method of analysis of orthotropic thin-walled structures subjected to static and dynamic load.

1.1. Static buckling

The buckling and postbuckling of thin-walled structures subjected to static load have been investigated by many authors for more than one hundred years. To the group of precursors of the investigation on the stability of thin-walled structures problem should be included following scientists: Euler [1], Timoshenko [2] and Volmir [3].

This chapter considered the thin plate or thin-walled structures composed of flat plates. Such structures have a various of buckling modes which can differ from one another both in quantitative (e.g., by the number of half-waves) and in qualitative (e.g., by global and local buckling) way.

The stability loss or buckling is a system transition from one equilibrium to another (the bifurcation point), or jump from the stable to the unstable equilibrium path (the limit point). Load resulting in the loss of stability is called the critical load. The behaviour of the structure subjected to load higher than the critical one can be described by a stable (the grow of displacement is caused by increased load) or unstable (displacements grow with decreasing load) postbuckling equilibrium path.

The postbuckling behaviour of the structures depends on their type. For example, the cylindrical shells subjected to axial compression change their equilibrium stage (buckling) by unstable bifurcation point or limit point. Long rods or columns subjected to axial compression have usually a sudden global buckling (bifurcation point of passage to the unstable postbuckling equilibrium path). Thin plates supported on all edges lose their stability having the local buckling mode and the stable postbuckling equilibrium path. Mentioned above type of buckling and postbuckling behaviour for given thin-walled structures are the same for ideal structures as for structures with geometrical imperfection. Columns made ​​of thin prismatic plates could have the local buckling mode, global (flexural, torsional or distorsional) one or coupled. The structures after local buckling are able to sustain further load, because increasing the displacement is only possible by increasing the load value (stable postbuckling equilibrium path), further increasing the load leads to plasticity or reaching the new, this time unstable bifurcation point (global buckling). The dangerous form of stability loss is the interactive buckling (coupled buckling), which usually causes the structure transition to the unstable equilibrium path what leads to the destruction of the structure with load lower than the critical load corresponding to each mode separately. The interaction of different buckling modes occurs when the critical loads corresponding to the different buckling modes are close to each other.

A more comprehensive review of the literature concerning the interactive buckling analysis of an isotropic structure can be found for example in Ali and Sridharan [3], Benito and Sridharan [5], Byskov [6], Koiter and Pignataro [7], Kolakowski [8–10], Manevich [11], Moellmann and Goltermann [12], Pignataro et al. [13], Pignataro and Luongo [14, 15], Sridharan and Ali [16, 17]. The interactive buckling of orthotropic structures can be found for example in [18, 19].

1.2. Dynamic buckling

In literature a quantity of ‘‘pulse intensity’’ [20] or ‘‘pulse velocity’’ [21] is introduced. The analysis of dynamic stability of plates under in-plane pulse loading can be divided into three categories depending on pulse duration and magnitude of its amplitude. For pulses of high intensity the impact phenomenon is observed whereas for pulses of low intensity the problem becomes quasi-static. The phenomenon of dynamic stability and dynamic buckling are often confused with each other. In this chapter the dynamic buckling phenomenon is examined but the concept of dynamic stability is broader and applies also to the stability of motion, which for thin-walled structures can be found for example in [22, 23]. The dynamic buckling occurs when the loading process is of intermediate amplitude and the pulse duration is close to the period of fundamental natural flexural vibrations (in range of milliseconds). In such case the effects of dumping are neglected [24]. Damping neglecting is only possible for problem solved in elastic range [25].

It should be noted that dynamic stability loss may occur only for structures with initial geometric imperfections; therefore the dynamic bifurcation load does not exist. For the ideal structures (without geometrical imperfection) the critical buckling amplitude of pulse loading tends to infinity [26]. The dynamic buckling load should be defined on the basis of the assumed buckling criterion.

The precise mathematical criteria were formulated for structures having unstable postcritical equilibrium path or having limit point [26, 27]. But for the structures having stable postbuckling equilibrium path (thin plate, thin-walled beam-columns with minimal critical load corresponding to local buckling) the precise mathematical criterion have not been defined till now.

Therefore Simitses [27] suggested not to define the dynamic buckling for the structures with stable postbuckling behaviour, but rather it should be defined as a dynamic response to pulse loads.

It is a reason why in world literature a lot of criteria can be found. In the sixties of the twenty century Volmir [28] proposed a criterion for plates subjected to in-plane pulse loading. The Volmir criterion - considered the easiest to use - states that the dynamic critical load corresponds to the amplitude of pulse force (of constant duration) at which the maximum plate deflection is equal to some constant value k (k - one half or one plate thickness) [28].

In many publications the dynamic buckling load is determined on the basis of stability criterion of Budiansky and Hutchinson [26, 29, 30]. However, this criterion was formulated for shell structures but also it can be used for the plate structures [31-34]. Budiansky and Hutchinson noticed that in some range of the amplitude value, the deflection of structures grows more rapidly than in other. Budiansky and Hutchinson formulated the following criterion: Dynamic stability loss occurs when the maximum deflection grows rapidly with the small variation of the load amplitude [26].

In the end of 90’s Ari-Gur and Simonetta [20] analysed laminated plates behaviour under impulse loading and formulated four own criteria of dynamic buckling, two of them of collapse-type conditions. One of them states: Dynamic buckling occurs when a small increase in the pulse intensity causes a decrease in the peak lateral deflection [20].

The failure criterion was proposed by Petry and Fahlbush [34], who suggest that for structures with stable postbuckling equilibrium path the Budiansky-Hutchinson criterion is conservative because it does not take into account load carrying- capacity of the structure.

Based on examples [35] it was noticed that for the thin-walled structures subjected to pulse loading, which lose their stability according to Budiansky-Hutchinson criterion or Volmir criterion, the maximal radius r_max calculated from characteristic root χ=a+jb (where j=$$\sqrt{-1}$$) of Jacoby matrix is equal or greater than unity in complex plane.

Therefore the criterion for thin-walled structures proposed by author [35] can be formulated as follows: Thin-walled structures subjected to pulse loading of finite duration lose their stability even if one characteristic root χ=a+jb of Jacoby matrix find for every time moment from 0 to 1.5T_p lies in the complex plane outside the circle with radius equal to unity.

Teter [36, 37] in his works analysed the long columns with longitudinal stiffeners and basing on the phase portrait for dynamic response of these structures defined the following criterion: The dynamic buckling load for the tracing time of solutions has been defined as the minimum value of the pulse load such that phase portrait is an open curve.

The dynamic buckling problem has been well known in the literature for over 50 years and was the subject of numerous works [20, 24, 26-34]. The extensive list of work dealing with dynamic buckling can be found for example in the book edited by Kowal-Michalska [38] or written by Simitses [39] or Grybos [40]. It seems that the analysis of dynamic buckling of thin-walled structures, especially structures with flat walls is not sufficiently investigated. There is a lack of both single-and multimodal analysis of dynamic buckling of columns with complex cross-sections made ​​of thin flat walls. The author of this paper decided to fill this gap presenting a method for the analysis of the local (single mode) and interactive (coupled mode – local and global) buckling of thin-walled structures subjected to pulse loading.

It should be mentioned that the presented method can be used only if the structures are in the elastic range. The case of dynamic buckling in elasto-plastic range including the viscoplastic effect has been investigated by Mania and Kowal-Michalska [25, 41-43]. In world literature it is also possible to find the paper dealing with the dynamic buckling of thin-walled structures subjected to combined load [44]. Czechowski [45] modelled the girder subjected to twist and bending considering only one plate subjected to shear and compression. The general summary showing which parameters have an influence on dynamic buckling of plated structures can be found in [46, 47].

2. Thin orthotropic plate theory

The thin isotropic or orthotropic plates with constant or widthwise variable material properties are considered. The thin-walled beam-columns or girders composed of mentioned above plates are also analysed. In order to taken into account all buckling modes (global, local and their interaction) the plate model was adopted to the analysed structures.

2.1. Basic assumptions

The basic assumption for thin plate are given by Kirchhoff for linear and by von Kármán and Marquerre for nonlinear thin plate theory. They made their assumption for isotropic material; lots of authors extended these assumptions for orthotropic or even for orthotropic multilayer thin plate [18, 48]. The assumptions are as follows:

• the plate is homogeneous (for example orthotropic homogenisation is made for fibre composite – resin matrix and fibre-reinforcement)

• the plate is thin – other dimensions (length and width) are at least 10 times higher than plate thickness;

• the material of the plate subjects to Hooke’s law;

• the plane stress state is considered for the plate – stress acting in the plate plane dominates the plate behaviour, stress acting in normal to plate plane direction are assumed to be zero;

• all strains (normal and shear) in plate plane are small compared to unity and they are linear;

• the strains of the plate to its normal direction are neglected (thickness of the plate do not change after deformation) – this assumption are made according to the Kirchhoff-Love hypothesis;

• straight lines normal to the mid-surface of the plate remain straight and normal to the mid-surface after deformation

• there is no interaction in normal direction between layers parallel to middle surface;

• deflections of the plate can be considerable in terms of nonlinear geometrical relations;

Additionally, it is assumed that principal axes of orthotropy are parallel to the edges of analysed structures (plate, beam, column, beam-column or girder).

2.2. Geometrical equations for thin plate

A plate model has been assumed for a thin plates and thin-walled beam-columns or girders. For easier explanation the plate (Figure 1a) or each i-th strip (Figure 1b) of the plate (or wall of the girder) or each i-th wall of the girder (Figure 1c) are called plate.

To describe the middle surface strains for each plate the following strain tensor have been assumed:

$$\begin{array}{l}{\epsilon }_{ix}^m=u_{i,x}+\frac{1}{2}(w_{i,x}^2+u_{i,x}^2+v_{i,x}^2),\\ {\epsilon }_{iy}^m=v_{i,y}+\frac{1}{2}(w_{i,y}^2+u_{i,y}^2+v_{i,y}^2),\\ {\gamma }_{ixy}^m=u_{i,y}+v_{i,x}+w_{i,x}{w}_{i,y}+u_{i,x}{u}_{i,y}+v_{i,x}{v}_{i,y},\end{array}$$E1

where: u_i, v_i, w_i - displacements parallel to the respective axes x_i, y_i, z_i of the local Cartesian system of co-ordinates, whose plane x_iy_i coincides with the middle surface of the i-th plate before its buckling (Figure 1).

In the majority of publications devoted to structure stability, the terms$$(u_{i,x}^2+v_{i,x}^2)$$, $$(u_{i,y}^2+v_{i,y}^2)$$and $$(u_{i,x}{u}_{i,y}+v_{i,x}{v}_{i,y})$$ are in general neglected for $${\epsilon }_{ix}^m,{\epsilon }_{iy}^m,{\gamma }_{ixy}^m$$ correspondingly, in (1) in the strain tensor components.

The change of the bending and twisting curvatures of the middle surface are assumed according to [48, 49] as follows:

$$\begin{array}{l}{\kappa }_{ix}=-w_{i,xx}\text{,}\\ {\kappa }_{iy}=-w_{i,yy}\text{,}\\ {\kappa }_{ixy}=-w_{i,xy}.\end{array}$$E2

The geometrical relationship given by equations (1) and (2) allow to consider both out-of-plane and in-plane bending of the plate.

Figure 1.

Possible models: plates, strips or wall with assumed dimension, coordinate systems and direction of deflections

2.3. Constitutive equations for orthotropy

Let’s consider orthotropic plate with principal axes of ortothropy 1 and 2 parallel to plate edges (Figure 2).

As same as in previous paragraph let’s consider i-th plate or strip of structures under analysis. The stress – strain relationship for orthotropic plate can be written in following form:

$$\left\{\begin{array}{c}{\sigma }_{i1}\\ {\sigma }_{i2}\\ {\tau }_{i12}\end{array}\right\}=\left[\begin{array}{ccc}{Q}_{i11}& Q_{i12}& 0\\ Q_{i21}& Q_{i22}& 0\\ 0& 0& Q_{66}\end{array}\right]\left\{\begin{array}{c}{\epsilon }_{i1}\\ {\epsilon }_{i2}\\ {\gamma }_{i12}\end{array}\right\},$$E3

where:

$$\left\{\begin{array}{c}{\epsilon }_{i1}\\ {\epsilon }_{i2}\\ {\gamma }_{i12}\end{array}\right\}=\left\{\begin{array}{c}{\epsilon }_{i1}^m\\ {\epsilon }_{i2}^m\\ {\gamma }_{i12}^m\end{array}\right\}+z\left\{\begin{array}{c}{\kappa }_{i1}\\ {\kappa }_{i2}\\ 2{\kappa }_{i12}\end{array}\right\}$$E4

$$\begin{array}{l}{Q}_{i11}=\frac{E_{i1}}{1-{\nu }_{i12}{\nu }_{i21}},\\ Q_{i12}=Q_{21}={\nu }_{i21}\frac{E_{i1}}{1-{\nu }_{i12}{\nu }_{i21}}={\nu }_{i12}\frac{E_{i2}}{1-{\nu }_{i12}{\nu }_{i21}},\\ Q_{i22}=\frac{E_{i2}}{1-{\nu }_{i12}{\nu }_{i21}},\\ Q_{i66}=G_{12}.\end{array}$$E5

and E_i1, E_i2 are Young modulus in longitudinal 1 and transverse 2 direction respectively, ν_i12 is a Poisson ratio for which strains are in longitudinal direction 1 and stress in transverse direction 2, G_i12 is a shear modulus (Kirchhoff modulus) in 12 plane.

Figure 2.

Plates or walls with principal axes of orthotropy

Young modulus and Poisson ratio occurring in (5) according to Betty-Maxwell theorem or according to symmetry condition of stress tensor should fulfil following relation:

$$E_{i1}{\nu }_{i21}=E_{i2}{\nu }_{i12}$$E6

.

For isotropic plate (wall of beam-columns) the constitutive equations are as follows:

$$\left\{\begin{array}{c}{\sigma }_{ix}\\ {\sigma }_{iy}\\ {\tau }_{ixy}\end{array}\right\}=\frac{E_i}{1-{\nu }_i^2}\left[\begin{array}{ccc}1& {\nu }_i& 0\\ {\nu }_i& 1& 0\\ 0& 0& \frac{1-{\nu }_i}{2}\end{array}\right]\left\{\left\{\begin{array}{c}{\epsilon }_{ix}^m\\ {\epsilon }_{iy}^m\\ {\gamma }_{ixy}^m\end{array}\right\}+z\left\{\begin{array}{c}{\kappa }_{ix}\\ {\kappa }_{iy}\\ 2{\kappa }_{ixy}\end{array}\right\}\right\}.$$E7

2.4. Generalized sectional forces

Substituting stress-strain relation from previous subchapter, the sectional moments and forces:

• for i-th isotropic plate or wall of beam-column are expressed by:

$$\begin{array}{l}\left\{\begin{array}{c}{N}_{ix}\\ N_{iy}\\ N_{ixy}\end{array}\right\}=\frac{E_i{h}_i}{1-{\nu }_i^2}\left[\begin{array}{ccc}1& {\nu }_i& 0\\ {\nu }_i& 1& 0\\ 0& 0& \frac{1-{\nu }_i}{2}\end{array}\right]\left\{\begin{array}{c}{\epsilon }_{ix}^m\\ {\epsilon }_{iy}^m\\ {\gamma }_{ixy}^m\end{array}\right\},\\ \left\{\begin{array}{c}{M}_{ix}\\ M_{iy}\\ M_{ixy}\end{array}\right\}=D_i\left[\begin{array}{ccc}1& {\nu }_i& 0\\ {\nu }_i& 1& 0\\ 0& 0& 1-{\nu }_i\end{array}\right]\left\{\begin{array}{c}{\kappa }_{ix}\\ {\kappa }_{iy}\\ {\kappa }_{ixy}\end{array}\right\},\end{array}$$E8

where:

$$D_i=\frac{E_i{h}_i^3}{12\left(1-{\nu }_i^2\right)}$$E9

• for i-th orthotropic strip or wall are:

$$\begin{array}{l}\left\{\begin{array}{c}{N}_{ix}\\ N_{iy}\\ N_{ixy}\end{array}\right\}=\frac{h_i}{1-{\nu }_{ixy}{\nu }_{iyx}}\left[\begin{array}{ccc}{E}_{ix}& {\nu }_{iyx}{E}_{ix}& 0\\ {\nu }_{ixy}{E}_{iy}& E_{iy}& 0\\ 0& 0& \left(1-{\nu }_{ixy}{\nu }_{iyx}\right)G_{ixy}\end{array}\right]\left\{\begin{array}{c}{\epsilon }_{ix}^m\\ {\epsilon }_{iy}^m\\ {\gamma }_{ixy}^m\end{array}\right\},\\ \left\{\begin{array}{c}{M}_{ix}\\ M_{iy}\\ M_{ixy}\end{array}\right\}=\left[\begin{array}{ccc}{D}_{ix}& {\nu }_{iyx}{D}_{ix}& 0\\ {\nu }_{ixy}{D}_{iy}& D_{iy}& 0\\ 0& 0& D_{ixy}\end{array}\right]\left\{\begin{array}{c}{\kappa }_{ix}\\ {\kappa }_{iy}\\ {\kappa }_{ixy}\end{array}\right\},\end{array}$$E10

where:

$$D_{ix}=\frac{E_{ix}{h}_i^3}{12\left(1-{\nu }_{ixy}{\nu }_{iyx}\right)},D_{iy}=\frac{E_{iy}{h}_i^3}{12\left(1-{\nu }_{ixy}{\nu }_{iyx}\right)},D_{ixy}=\frac{G_{ixy}{h}_i^3}{6}.$$E11

2.5. Dynamic equations of stability for thin plate

Differential equations of motion of the plate were derived basing on Hamilton’s principle. It states that the dynamics of a physical system is determined by a variation problem for a functional based on a single function, the Lagrangian, which contains all physical information concerning the system and the forces acting on it. In dynamic buckling problem the motion should be understand as the time dependent deflection.

The Hamilton’s principles for conservative systems states that the true evolution (compatible with constrains) of the system between two specific states in specific time range (t₀, t₁) is a stationary point (a point where the variation is zero) of the action functional Ψ. Action functional Ψ for i-th plate is described by following equation:

$$\Psi ={\displaystyle \underset{t_0}{\overset{t_1}{\int }}\Lambda dt={\displaystyle \underset{t_0}{\overset{t_1}{\int }}(K}-\Pi )dt}$$E12

where Λ is the Lagrangian function for the system, K is a kinetic energy of the system and Π is a total potential energy of the system.

The subscript i denoting i-th plate or strip in all equations in this subchapter is omitted – all equations are presented for one plate, which could be i-th plate, wall or strip of considered plate, beam-columns or girder (Figure 1).

Taking the action functional Ψ in form (9) the Hamilton’s principle can be written as:

$$\delta \Psi =\delta {\displaystyle \underset{t_0}{\overset{t_1}{\int }}\Lambda }dt=\delta {\displaystyle \underset{t_0}{\overset{t_1}{\int }}(K}-\Pi )dt=0$$E13

.

The total potential energy variation δΠ for i-th thin plate (or strip) can be written in form:

$$\delta \Pi =\delta Q-\delta W$$E14

where δQ is a variation of internal elastic strain energy:

$$\delta Q={\displaystyle \underset{\Omega }{\int }{({\sigma }_x\delta {\epsilon }_x+{\sigma }_y\delta {\epsilon }_y+{\tau }_{xy}\delta {\gamma }_{xy})}_{}}d\Omega$$E15

,

and Ω is the volume of the plate and S is its area, the volume can be expressed as l b h or S h.

The variation of internal elastic strain energy for i-th plate or strip could be expressed by strain and sectional forces and moments in a following way:

$$\begin{array}{c}\delta Q=\delta Q^m+\delta Q^b=\\ ={\displaystyle \underset{S}{\int }(N_x\delta {\epsilon }_x^m+N_y\delta {\epsilon }_y^m+N_{xy}\delta {\gamma }_{xy}^m)}dS-{\displaystyle \underset{S}{\int }(M_x\delta w_{,xx}+M_y\delta w_{,yy}+2M_{xy}\delta w_{,xy})}dS.\end{array}$$E16

The work W of external forces (neglecting the out-of plane load) done on i-th plate can be written as follows:

$$W={\displaystyle \underset{0}{\overset{b}{\int }}h[p_{}^0(y)u+{\tau }_{xy}^0(y)}v]dy+{\displaystyle \underset{o}{\overset{\ell }{\int }}h[p_{}^0(x)v}+{\tau }_{xy}^0(x)u]dx,$$E17

where: p⁰(x), p⁰(y), τ⁰_xy(x), τ⁰_xy(y) are the prebuckling load applied to the middle surface of the considered plate (wall or strip)

For thin plates, it is assumed that the displacements u and v do not depend on rotation w_,x and w_,y and therefore do not depend on the coordinate z. This approach results in exclusion of rotational inertia [50] in the equation for kinetic energy, which for the i-th thin plate (strip) can be written as:

$$K=\frac{1}{2}\rho {\displaystyle \underset{\Omega }{\int }\left({(\stackrel{.}{u})}^2+{(\stackrel{.}{v})}^2+{(\stackrel{.}{w})}^2\right)d\Omega }$$E18

The Hamilton’s principle, it is the variation of the action functional δΨ (10) for i-th thin plate (strip or wall) which after taking into consideration equations from (11) to (15) can be written as:

$$\delta \Psi ={\displaystyle \underset{t_0}{\overset{t_1}{\int }}\left(\delta K-\delta Q^m-\delta Q^b+\delta W\right)dt}=0$$E19

The Lagrangian function for the whole system is equal to the sum of the Lagrangian functions of all n plates of which the system was composed. To determine the variation of action δΨ for i-th plate, the following identity:

$$X\delta Y=\delta (XY)-Y\delta X$$E20

was used.

In the obtained equation, terms with the same variations were grouped, and then each of the obtained groups of terms (due to the mutual independence of variations) were equated to zero, giving:

• equilibrium equations:

$$\begin{array}{l}{\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\{[N_{x,x}+N_{xy,y}+{(N_x{u}_{,x})}_{,x}+{(N_y{u}_{,y})}_{,y}+{(N_{xy}{u}_{,x})}_{,y}+{(N_{xy}{u}_{,y})}_{,x}]}}-h\rho \stackrel{\mathrm{..}}{u}\}\delta udSdt=0\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\{[N_{xy,x}+N_{y,y}+{(N_x{v}_{,x})}_{,x}+{(N_y{v}_{,y})}_{,y}+{(N_{xy}{v}_{,x})}_{,y}+{(N_{xy}{v}_{,y})}_{,x}]}}-h\rho \stackrel{\mathrm{..}}{v}\}\delta vdSdt=0\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\{[}{M}_{x,xx}+M_{y,yy}+2M_{xy,xy}+{(N_x{w}_{,x})}_{,x}+{(N_y{w}_{,y})}_{,y}+{(N_{xy}{w}_{,x})}_{,y}+{(N_{xy}{w}_{,y})}_{,x}]+}\\ -h\rho \stackrel{\mathrm{..}}{w}\}\delta wdSdt=0\end{array}$$E21

• boundary conditions for lateral edges of the plate (x = const):

$$\begin{array}{l}{\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}[N_x+N_x{u}_{,x}+N_{xy}{u}_{,y}-h{p}^0(y)]\delta udydt|{}_{x=const}=0}}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}[N_{xy}+N_x{v}_{,x}+N_{xy}{v}_{,y}-h{\tau }_{xy}^0(y)]\delta vdydt|{}_{x=const}=0}}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}{M}_x}\delta w_{,x}dydt|{}_{x=const}=0}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}(M_{x,x}+2M_{xy,y}+N_x{w}_{,x}+N_{xy}{w}_{,y})\delta wdydt|{}_{x=const}=0}}\end{array}$$E22

• boundary conditions for longitudinal edges of the plate (y = const):

$$\begin{array}{l}{\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{\ell }{\int }}[N_y+N_y{v}_{,y}+N_{xy}{v}_{,x}-h{p}^0(x)]\delta vdxdt|{}_{y=const}=0}}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{\ell }{\int }}[N_{xy}+N_y{u}_{,y}+N_{xy}{u}_{,x}-h{\tau }_{xy}^0(x)]\delta udxdt|{}_{y=const}=0}}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{\ell }{\int }}{M}_y}\delta w_{,y}dxdt|{}_{y=const}=0}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{0}{\overset{\ell }{\int }}(M_{y,y}+2M_{xy,x}+N_y{w}_{,y}+N_{xy}{w}_{,x})\delta wdxdt|{}_{y=const}=0}}\end{array}$$E23

• boundary condition for the plate corners (x = const and y = const):

$${\displaystyle \underset{t_0}{\overset{t_1}{\int }}2M_{xy}\delta wdt|{}_{x=const}|{}_{y=const}=0}$$E24

• initial conditions for t = const:

$$\begin{array}{l}{\displaystyle \underset{S}{\int }h\rho \stackrel{.}{u}\delta udS|{}_{t=const}=0}\\ {\displaystyle \underset{S}{\int }h\rho \stackrel{.}{v}\delta vdS|{}_{t=const}=0}\\ {\displaystyle \underset{S}{\int }h\rho \stackrel{.}{w}\delta wdS|{}_{t=const}=0}\end{array}$$E25

Above conditions are fulfilled for the entire structure, so if one apply the restrictions in moment of the initial t₀ and in moment of the final t₁ that the displacement variations are zero at all points of the structure. Then the system of equations (23) vanishes.

• already used the relationship between deformations and internal forces and moments (8) or (9):

$$\begin{array}{l}{\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\left(E_{x}h{\epsilon }_x-N_x+{\nu }_{xy}{N}_y\right)}\delta N_{x}dSdt=0}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\left(E_{y}h{\epsilon }_y+{\nu }_{yx}{N}_x-N_y\right)\delta N_{y}dSdt=0}}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\left(2Gh{\epsilon }_{xy}-N_{xy}\right)\delta N_{xy}dSdt=0}}\end{array}$$E26

$$\begin{array}{l}{\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\left(\frac{E_x{h}^3}{12}{\kappa }_x-M_x+{\nu }_{xy}{M}_y\right)\delta M_{x}dSdt=0}}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\left(\frac{E_y{h}^3}{12}{\kappa }_y+{\nu }_{yx}{M}_x-M_y\right)\delta M_y}dSdt=0}\\ {\displaystyle \underset{t_0}{\overset{t_1}{\int }}{\displaystyle \underset{S}{\int }\left(\frac{G{h}^3}{6}{\kappa }_{xy}-M_{xy}\right)\delta M_{xy}dSdt=0}}\end{array}$$E27

(25)

3. Solution method

To determine the critical loads, natural frequencies and the coefficients of the equation describing the postbuckling equilibrium path, the analytical-numerical method has been employed. The proposed method also allows analysing dynamic response of the structure subjected to pulse loading. Taking the time courses of deflections and applying the relevant dynamic buckling criteria it is possible to determine the dynamic critical load.

3.1. Equilibrium equations

The differential equations of equilibrium for orthotropic plate or strip directly from the equations (19) can be derived and have the form:

$$\begin{array}{l}{N}_{x,x}+N_{xy,y}+\{{(}_{N_x}+{(}_{N_y}+{(}_{N_{xy}}+{(}_{N_{xy}}\}-h\rho \stackrel{\mathrm{..}}{u}=0\\ N_{xy,x}+N_{y,y}+\{{(}_{N_x}+{(}_{N_y}+{(}_{N_{xy}}+{(}_{N_{xy}}\}-h\rho \stackrel{\mathrm{..}}{v}=0\\ M_{x,xx}+M_{y,yy}+2M_{xy,xy}+{(}_{N_x}+{(}_{N_y}+{(}_{N_{xy}}+{(}_{N_{xy}}-h\rho \stackrel{\mathrm{..}}{w}=0\end{array}$$E28

Above equilibrium equations after omitting the inertia forces$$h\rho \stackrel{\mathrm{..}}{u}$$, $$h\rho \stackrel{\mathrm{..}}{v}$$and $$h\rho \stackrel{\mathrm{..}}{w}$$becomes the equilibrium equations for thin plates allowing analysis of both local and global buckling mode.

3.2. Boundary and initial condition

As the wave propagation effects have been neglected, the boundary conditions referring to the simply supported columns at their both ends, i.e. x = 0 and x = l, according to (20), are assumed to be:

$$\begin{array}{l}\frac{1}{b_i}{\displaystyle \int N_{ix}\left(x_i=0,y_i,t\right)d{y}_i}=\frac{1}{b_i}{\displaystyle \int N_{ix}\left(x_i=l,y_i,t\right)}d{y}_i=N_{ix}^{(0)},\\ v_i\left(x_i=0,y_i,t\right)=v_i\left(x_i=l,y_i,t\right)=0,\\ w_i\left(x_i=0,y_i,t\right)=w_i\left(x_i=l,y_i,t\right)=0,\\ M_{ix}\left(x_i=0,y_i,t\right)=M_{iy}\left(x_i=l,y_i,t\right)=0,\end{array}$$E29

The condition written as a first of equations of (27) is satisfied for the prebuckling state and first-order approximation, the condition for deflection v (27) is satisfied for the first and second order of approximations, while the other two conditions are met for prebuckling state as well as for the first and second order of approximation. The condition of displacement in the y direction in the prebuckling state can be found for example in [51]. This approach allows to take into account the Poisson effect on the edges of the walls of the column. The boundary conditions described by equations (27) assume the lack of displacement possibility of points lying at the loaded edges in the transverse v and normal w directions to the surface in a wall or column. Furthermore, it is assumed that the moments M_ix (as a vector parallel to the edge of the plate or end edge of the column walls) are zero.

For structures with material properties varying widthwise the strip model was adopted what forces the boundary conditions modification in the second order approximations [51]. Modification consists of changing the first condition of (27) onto the following form:

$${\displaystyle \sum _{i=1}^J\frac{1}{b_i}{\displaystyle \underset{0}{\overset{b_i}{\int }}{N_{ix}^{(2)}d{y}_i|}_{x=0;l}}}=0$$E30

Summation is performed only for the J number of the strips, between which the angle φ_i,i+1 (Figure 3) is equal to zero.

To determine the boundary conditions on the longitudinal edges of plates or free edges of columns with open cross-sections the equations (20) were used. Whereas, directly from equations (23) result the following initial conditions:

$$\begin{array}{l}\stackrel{.}{u_i}(x_i,y_i,t=t_0)=\widetilde{u_i}(x_i,y_i)\text{ and }{u}_i(x_i,y_i,t=t_0)=\overline{u_i}(x_i,y_i),\\ \stackrel{.}{v_i}(x_i,y_i,t=t_0)=\widetilde{v_i}(x_i,y_i)\text{ and }{v}_i(x_i,y_i,t=t_0)=\overline{v_i}(x_i,y_i),\\ \stackrel{.}{w_i}(x_i,y_i,t=t_0)=\widetilde{w_i}(x_i,y_i)\text{ and }{w}_i(x_i,y_i,t=t_0)=\overline{w_i}(x_i,y_i),\end{array}$$E31

where the following functions $$\overline{u_i},\overline{v_i},\overline{w_i},\stackrel{\text{~}}{u_i}\text{ ,}\stackrel{\text{~}}{v_i}\text{, }\stackrel{\text{~}}{w_i}$$ are given for the initial moment t = t₀.

3.3. Interaction condition between adjacent plates

Static and kinematic junction conditions on the longitudinal edges of adjacent plates (Figure 3), according to (21), can be written as:

$$\begin{array}{l}{u_{i+1}|}^{-}={u_i|}^{+},\\ {w_{i+1}|}^{-}={w_i|}^{+}\cos (\varphi )-{v_i|}^{+}\sin (\varphi ),\\ {v_{i+1}|}^{-}={w_i|}^{+}\sin (\varphi )+{v_i|}^{+}\cos (\varphi ),\\ {w_{i+1,y}|}^{-}={w_{i,y}|}^{+},\\ {M_{(i+1)y}|}^{-}={M_{iy}|}^{+},\\ {N_{(i+1)y}^{*}|}^{-}-{N_{iy}^{*}|}^{+}\cos (\varphi )-{Q_{iy}^{*}|}^{+}\sin (\varphi )=0,\\ {Q_{(i+1)y}^{*}|}^{-}+{N_{iy}^{*}|}^{+}\sin (\varphi )-{Q_{iy}^{*}|}^{+}\cos (\varphi )=0,\\ {N_{(i+1)xy}^{*}|}^{-}={N_{ixy}^{*}|}^{+},\end{array}$$E32

where:

$$\begin{array}{l}{N}_{ixy}^{*}=N_{iy}+N_{iy}{v}_{i,y}+N_{ixy}{v}_{i,x}\\ N_{ixy}^{*}=N_{ixy}+N_{ixy}{u}_{i,x}+N_{iy}{u}_{i,y}\\ M_{iy}=-{\eta }_i{D}_i(w_{i,yy}+{\nu }_i{w}_{i,xx}),\\ Q_{iy}^{*}=-{\eta }_i{D}_i{w}_{i,yyy}-({\nu }_i{\eta }_i{D}_i+2D_{1i})w_{i,xxy}+N_{iy}{w}_{i,y}+N_{ixy}{w}_{i,x},\\ \varphi ={\varphi }_{i;i+1}.\end{array}$$E33

Figure 3.

The geometrical dimensions and local coordinate systems adjacent plates

3.4. Buckling and postbuckling equilibrium paths

A non-linear stability problem has been solved by means of the Koiter’s asymptotic theory. The displacement field$$\overline{U}$$, and sectional force field $$\overline{N}$$ have been expanded into the power series with respect to the parameter ξ, - the buckling linear eigenvector amplitude (normalised with the equality condition between the maximum deflection and the thickness of the first plate h₁).

$$\begin{array}{l}\overline{U}=\lambda {\overline{U}}_{}^{(0)}+{\xi }_i{\overline{U}}_{}^{(i)}+{\xi }_i{\xi }_j{\overline{U}}_{}^{(ij)}+\mathrm{...}\\ \overline{N}=\lambda {\overline{N}}_{}^{(0)}+{\xi }_i{\overline{N}}_{}^{(i)}+{\xi }_i{\xi }_j{\overline{N}}_{}^{(ij)}+\mathrm{...}\end{array}$$E34

It was assumed that the dimensionless amplitude of the initial deflections (imperfections) correspond to the considered buckling mode (for s-th buckling mode) is:

$$\overline{U}={\xi }_s^{*}{\overline{U}}_{}^{(i)}.$$E35

By substituting expansions (32) into equations of equilibrium (26) with neglected inertia terms (static buckling problem), junction conditions (30) and boundary conditions (27), the boundary problem of the zero (superscript ⁽⁰⁾ in Equations (32) and further), first (superscript ⁽ⁱ⁾) and second (superscript ^(ij)) order has been obtained [18, 50, 52, 53]. The zero approximation describes the prebuckling state, whereas the first order approximation allows for determination of critical loads and the buckling modes corresponding to them, taking into account minimisation with respect to the number of half-waves m in the lengthwise direction. The second order approximation is reduced to a linear system of differential heterogeneous equations, which right-hand sides depend on the force field and the first order displacements only.

The most important advantage of this method is that it enables us to describe a complete range of behaviour of thin-walled structures from all global (i.e. flexural, flexural–torsional, lateral, distortional buckling and their combinations) to the local dynamic stability. In the solution obtained, the shear lag phenomenon, the effect of cross-sectional distortions and also the interaction between all the walls of structures are included.

Having found the solutions to the first and second order of the boundary problem, the coefficients a_ijs, b_ijks have been determined [18, 50, 52, 53]:

$$\begin{array}{l}{a}_{ijs}=\frac{{\sigma }^{(i)}*L_{11}(U^{(j)},U^{(s)})+0.5{\sigma }^{(s)}*L_{11}(U^{(i)},U^{(j)})}{-{\lambda }_s{\sigma }^{(0)}*L_2(U^{(s)})},\\ \\ b_{ijks}=\frac{2{\sigma }^{(i)}*L_{11}(U^{(jk)},U^{(s)})+{\sigma }^{(ij)}*L_{11}(U^{(k)},U^{(s)})}{-{\lambda }_s{\sigma }^{(0)}*L_2(U^{(s)})},\end{array}$$E36

where: λ_s – is the critical load corresponding to the s-th mode, L₁₁ is the bilinear operator, L₂ is the quadratic operator and σ⁽ⁱ⁾, σ^(ij) are the stress field tensors in the first and second order.

The postbuckling static equilibrium paths for coupled buckling can be described by the equation:

$$\left(1-\frac{\lambda }{{\lambda }_s}\right){\xi }_s+a_{ijs}{\xi }_i{\xi }_j+b_{ijks}{\xi }_i{\xi }_j{\xi }_k=\frac{\lambda }{{\lambda }_s}{\xi }_s^{*};\left(\text{s }=\text{ 1},\dots ,N\right),$$E37

which for the uncoupled problem have the form:

$$\left(1-\frac{\lambda }{{\lambda }_{cr}}\right)\xi +a_{111}{\xi }^2+b_{1111}{\xi }^3={\xi }^{*}\frac{\lambda }{{\lambda }_{cr}}$$E38

where λ_cr is the critical load value.

In a special case, i.e. for the so-called ideal structure without initial imperfections (ξ*=0) and when the equilibrium path (a₁₁₁) is symmetrical, the postbuckling equilibrium path is defined by the equation:

$$\frac{\lambda }{{\lambda }_{cr}}=1+b_{1111}{\xi }^2$$E39

(37)

3.5. Natural frequencies

Determination of the natural frequencies is similar to the determination of critical buckling load and the natural frequencies are found by solving the eigenvalue problem.

Natural frequencies of thin-walled structures were determined by solving a dynamic problem, which uses the approach proposed by Koiter in his asymptotic stability theory of conservative systems in the first-order approximation [52].

To determine the natural frequencies [55] of the structure the adopted equilibrium equations (26) contain cross-sectional inertia forces acting in the direction normal to the middle surface of the plate (column wall) and in the middle plane of plate (i.e. $$h\rho \stackrel{\mathrm{..}}{u}$$≠ 0 and $$h\rho \stackrel{\mathrm{..}}{v}$$≠ 0).

3.6. Lagrange equations

In the dynamic analysis (while finding the frequency of natural vibrations [55]), the independent non-dimensional displacement ξ and the load factor λ become a function dependent on time, and dynamic terms were added to equations describing postbuckling equilibrium path. Neglecting the forces associated with the inertia terms of prebuckling state and the second-order approximations, and taking into account the orthogonality conditions for the displacement field in the first $${\overline{U}}_{}^{(i)}$$and second-order approximation$${\overline{U}}_{}^{(ij)}$$, the Lagrange equations can be written as [56]:

$$\frac{1}{{\omega }_s^2}{\stackrel{\mathrm{..}}{\xi }}_s+\left(1-\frac{\lambda }{{\lambda }_s}\right){\xi }_s+a_{ijs}{\xi }_i{\xi }_j+b_{ijks}{\xi }_i{\xi }_j{\xi }_k={\xi }_s^{\ast }\frac{\lambda }{{\lambda }_s};\left(s=\text{1},\text{2},\dots ,N\right)$$E40

where ω_s is a natural frequency with mode corresponding to buckling mode; a_ijs and b_ijks are the coefficients (34) describing the postbuckling behaviour of the structure (independent of time); however the parameters of load λ and the displacement ξ are the functions of time t.

For the uncoupled buckling, i.e. the single-mode buckling (where index s = N = 1), the equations of motion may be written in the form:

$$\frac{1}{{\omega }_1^2}\stackrel{\mathrm{..}}{{\xi }_1}+\left(1-\frac{\lambda }{{\lambda }_1}\right){\xi }_1+a_{111}{\xi }_1^2+b_{1111}{\xi }_1^3={\xi }_1^{*}\frac{\lambda }{{\lambda }_1};$$E41

(39)

It is assumed that in the initial moment of time t = 0 the non-dimensional displacement ξ, as well as the velocity of displacement are equal to zero, i.e.:

$$\xi \left(t=0\right)=0 \text{and} \stackrel{.}{\xi }(t=0)=0$$E42

.

The Runge-Kutta method [57] for solving the equation (39) requires the following substitutions:

$$\begin{array}{l}\stackrel{.}{\xi }=\Gamma (t),\\ \stackrel{.}{\Gamma }=-{\omega }_1^2\left(1-\frac{\lambda (t)}{{\lambda }_1}\right)\xi -{\omega }_1^2{b}_{111}{\xi }_1^2-{\omega }_1^2{b}_{1111}{\xi }_1^3+{\omega }_1^2\frac{\lambda (t)}{{\lambda }_1}{\xi }^{*},\end{array}$$E43

which lead to the system of two differential equations. ‘‘Complete’’ equations of motion (41) are solved with the numerical Runge–Kutta method of order 8 (5,3), thanks to Dormand and Price (with step-size control and density output).

4. Exemplary results of calculations

The exemplary results of numerical calculation are presented in this sub-chapter. All results are obtained using explained above proposed analytical-numerical method (ANM) based on the nonlinear orthotropic plate theory.

The material properties (E – Young modulus, ν – Poisson ratio, G=E/[2(1+ν)] – Kirchhoff modulus; ρ – density) for materials taken into account are presented in Table 1.

material type:	E; [GPa]	ν	ρ; [kg/m^3]
steel	200	0.3	7850
aluminium	70	0.33	2950
epoxy resin	3.5	0.33	1249
glass fibre	71	0.22	2450

Table 1.

Assumed material properties

The fibre composite material was modelled as orthotropic but for components (resin and fibre) the isotropic material properties (Table 1) was assumed. Necessary equations for material properties homogenization based on theory of mixture [57, 58] are as follows:

$$\begin{array}{l}{E}_x=E_m\left(1-f\right)+E_{f}f,\\ E_y=E_m\frac{E_m\left(1-\sqrt{f}\right)+E_f\sqrt{f}}{E_m{\left[1-\sqrt{f}\left(1-\sqrt{f}\right)\right]}^{}+E_f\sqrt{f}\left(1-\sqrt{f}\right)},\\ {\nu }_{yx}={\nu }_m\left(1-\sqrt{f}\right)+{\nu }_f\sqrt{f},\\ G=G_m\frac{G_m\sqrt{f}\left(1-\sqrt{f}\right)+G_f\left[1-\sqrt{f}\left(1-\sqrt{f}\right)\right]}{G_m\sqrt{f}+G_f\left(1-\sqrt{f}\right)}.\end{array}$$E44

where E_m and E_f are the Young’s modulus of elasticity for matrix and fibre, respectively, G_m and G_f are the shear modulus for matrix (subscript m) and fibre (subscript f), ν_m and ν_f are the Poisson’s ratios for matrix and fibre and f = V_f /(V_m + V_f) is the fibre volume fraction.

For static buckling the critical buckling load and corresponding modes are presented as well as the postbuckling equilibrium paths.

For dynamic buckling the proposed by Budiansky and Hutchinson parameter called Dynamic Load Factor DLF is introduced. The DLF is defined as a ratio of pulse loading amplitude to static buckling load. The results are presented of nondimensional deflection ξ versus DLF. The critical dynamic load factor DLF_cr corresponding to dynamic buckling has been estimated using different criteria – the obtained results were compared.

For the proposed method the validation of the results was made by comparison with the other Authors [34] calculations (Figure 4) or with the results obtained with FEM [38]. The results presented in Figure 4 were obtained for thin (ratio length to thickness equals 200) aluminium square plate simply supported at all edges and subjected to sinusoidal pulse load. The time of pulse duration was equal to the period of natural vibration of the plate. The considered plate has a geometrical imperfection corresponding to buckling mode with amplitude equal to 0.05 of the plate thickness.

Figure 4.

The results of different calculation comparison

4.1. Plates

The rectangular thin plates simply supported on loaded edges with different boundary conditions along the unloaded ones were considered (Figure 5). On the longitudinal edges five different boundary condition cases were taken into account. Following notations is used in Figure 5: s – simply supported edge, c – clamped edge, e – free edge.

Figure 5.

Analysed plates with different boundary conditions

material:	boundary condition	Pcr [kN]		ω [rad/s]
ANM	FEM	ANM	FEM
steel	ss	7.23	7.24	3016	3010
cc	15.6 (m = 1); 13.9 (m = 2)	15.7; 14.0	4423 (m = 1); 8363 (m = 2)	4423; 8344
se	2.53	2.54	1784	1784
ce	2.99	2.99	1935	1935
sc	10.38	10.41	3613	3607
composite; f = 0.5	ss	0.54	0.54	1703	1709
cc	0.93	0.95	2231	2237
se	0.34	0.35	1351	1351
ce	0.36	0.37	1389	1389
sc	0.69	0.70	1916	1923

Table 2.

Critical load P_cr and natural frequencies ω for analysed plates

Exemplary results were calculated for steel and epoxy glass composite (fibre volume factor f = 0.5) square plates subjected to rectangular compressive pulse loading. The buckling load for plate under analysis is presented in Table 2. The pulse duration T_p was equal to the period of natural vibration with mode corresponding to the buckling mode.

The dimensions of analysed plates were assumed as follows: the length (width) a= b= 100 mm and thickness h = 1 mm.

The geometrical imperfection was assumed in the shape corresponding to the buckling mode with amplitude ξ^* = 0.01, where ξ^* is an amplitude of deflection divided by the plate thickness. The Figure 6 presents postbuckling equilibrium paths for ideal flat composite plate (Figure 6a) and for plate with geometrical imperfection with amplitude ξ^* =0.01 (Figure 6b).

Figure 6.

Postbuckling equilibrium paths for square ideal plates (a) and plates with imperfection (b) with different boundary conditions on non-loaded edges

In the dynamic buckling case the results are shown as graphs presenting nondimensional deflection ξ or radius r calculated from real and imaginary part of maximal characteristic root of Jacoby matrix as a function of dynamic load factor DLF. The graphs mentioned above allow to find critical amplitude of pulse loading using the proposed criterion (PC) [35] and to compare the obtained results with Budiansky-Hutchinson (B-H) or Volmir (V) criteria. In brackets the notation used in Figures and Tables is given. The critical deflection according to Volmir criterion was assumed as ξ_cr= 1.

Figure 7.

Nondimensional deflection ξ (a) and maximum radius r_max (b) vs. DLF for square plates with different boundary conditions on non-loaded edges [35]

Basing on curves presented in Figure 7 the critical value of dynamic load factor can be found. The comparison of obtained critical DLF values using different criteria is presented in Table 3. All critical DLF values except the case denoted as se obtained from the proposed criterion (PC) are in the range obtained from Budiansky-Hutchinson criterion (B-H). For the case denoted as se the greatest differences between critical dynamic load factors from the proposed and Budiansky-Hutchinson criterion were obtained but these differences are less than 10%.

boundary conditions for unloaded edges	mode	B-H	V	PC
se	m=1	1.51.6	1.55	1.35
ce	m=1	1.31.4	1.58	1.35
ss	m=1	1.51.6	1.43	1.51
sc	m=1	1.41.8	1.46	1.58
cc	m=2	1.41.5	1.51	1.39

Table 3.

Comparison of DLF_cr obtained from different criteria for compressed plate

Figure 8.

Nondimensional deflection ξ vs. DLF for simply supported square plate

The comparison for postbuckling behaviour of the rectangular plate simply supported at all edges subjected to static and dynamic load are obtained using proposed analytical-numerical method and presented in Figure 8.

4.2. Segments of the girders

As a next example the static and dynamic buckling of composite (epoxy glass composite with different volume fibre fraction f) girders with open cross-section (Figure 9) is presented. The assumed boundary conditions on loaded edges correspond to simply support. The calculation was carried out for short segment of girder with length to web width ratio l/b₁ = 1 and for the following dimensions of the cross-section: b₁/h = 50, b₂/h = 25 and b₃/h = 12.5.

Figure 9.

Cross-sections of analysed segment of the girders

The geometrical imperfection was assumed in the shape corresponding to the buckling mode with amplitude ξ^* = 0.01. The static buckling load and fundamental flexural natural frequency obtained with analytical-numerical method for girders made of composite with different fibre fraction are presented in Tables 4. The static critical buckling loads are presented in Table 4 and postbuckling equilibrium paths are presented in Figure 10.

volume fibre fraction f: cross-section	critical load; Pcr [N]		natural frequency; ω [rad/s]
0.4	0.6	0.4	0.6
channel (Figure 9a)	1526	2281	1076	1076
channel with inner stiffeners (Figure 9c)	2821	4217	1308	1308
omega (Figure 9b)	2819	4214	1308	1308

Table 4.

Critical load and natural frequencies for analysed girder’s segments

Figure 10.

Postbuckling equilibrium paths for segment of girders

Buckling load and natural frequency for girder segment with omega and stiffened channel cross-section are similar – it is true only for local buckling case. For girder with channel cross-section the buckling was caused by flanges – this is a reason why for this cross-section the buckling load and natural frequency are smaller than for two others analysed cross-sections. Looking at obtained results (Table 4) it can be seen that increasing the volume fibre fraction f leads to an increasing the buckling loads as well as the natural frequencies. The postbuckling equilibrium paths for stiffened cross-section (channel with inner stiffeners and omega) overlap. The postbuckling path for channel cross-section lies below the equilibrium paths of girders with stiffened cross-section – it is obvious because the girders with stiffened cross-section have similar stiffness and the girder with channel cross-section is more flexible.

In dynamic buckling case the time of pulse duration T_p was assumed as the period of natural fundamental flexural vibration corresponding to the local buckling mode. Considered shapes of pulse loading are presented in Figure 11.

Figure 11.

Shapes of pulse loading: three triangular impulses T1 (a), T2 (b), T3 (c), rectangular R (d) and sinusoidal S (e)

Figure 12.

Dimensionless deflection vs. DLF for channel girder – results comparison obtained with analytical-numerical method ANM and finite element method FEM [56]

The results presented in Figure 12 were obtained with the proposed analytical-numerical method (ANM) and compared with FEM computations. They are similar for both assumed shapes of pulse loading (Figure 11): rectangular (R) and sinusoidal (S). Nevertheless, for rectangular shape pulse loading some small differences in deflection are visible for DLF greater than 2. It should be noted that obtained curves (Figure 12) from both methods allow to find the same critical dynamic load factor DLF_cr using Budiansky-Hutchinson or Volmir criterion – for rectangular pulse loading DLF_cr ≈1.4 (both criteria) and for sinusoidal pulse loading DLF_cr ≈2.1 (Budiansky-Hutchinson) or 2.3 (Volmir).

In Figure 13 the dynamic response comparison of girders made of composite material (f = 0.5) with different cross-section subjected to triangular T3 pulse was presented. The curves for omega cross-section and channel cross-section with inner stiffeners cover each other’s.

Figure 13.

Nondimensional deflection ξ as a DLF function for girder with different channel cross-section subjected to T3 pulse loading [56]

Dynamic responses for girder with channel cross-section with inner stiffeners for different pulse loading are presented in Figure 14. The curves denoted by S=R were obtained for sinusoidal pulse loading with the same area as rectangular pulse loading (for the same pulse duration the amplitude was higher for sinusoidal pulse). The highest deflection was obtained for pulse denoted by S=R because this pulse has the highest amplitude. The rest of compared pulses have the same amplitude and the same duration. Analysing the curves ξ(DLF) for rectangular, sinusoidal and three triangular impulses it can be seen that the highest increment of deflection for the smallest DLF takes place for rectangular pulse loading.

Figure 14.

Dimensionless deflection vs. dynamic load factor for different shapes of pulse loading - channel with inner stiffeners, composite material f = 0.7 [56]

The comparison of DLF_cr obtained using Budiansky-Hutchinson criterion for girders with different cross-section made of composite materials (f = 0.5) are presented in Table 5. Only the average values from the obtained critical ranges are presented. The dynamic load factors for different cross-sections are in the same relation as buckling loads (Table 4) – the same or similar DLF_cr for cross-sections with stiffeners (omega and channel with inner stiffeners).

4.3. Columns

Next, the exemplary results for the dynamic interactive buckling of channel cross-section the columns are presented. The columns subjected to rectangular compressive pulse loading were analysed. The calculation was carried out for columns with various length to web width ratio l/ b₁ = 4; 6 and 8 and for the following dimensions of the cross-section: width of the web to its thickness b₁/h = 100, width of the flange to its thickness b₂/h = 50.

analysed cross-section: type of pulse	channel	omega	channel with inner stiffeners
S; R; S=R; T1; T2; T3	2.0; 1.6; 1.4; 3.1; 2.3; 2.4	2.1; 1.6; 1.4; 3.1; 2.5; 2.5	2.1; 1.6; 1.4; 3.1; 2.5; 2.5

Table 5.

The DLF_cr value for different shape of applied pulses

Figure 15.

Nondimensional deflection ξ (a) and maximum radius r_max (b) as a function of DLF for channel beam-columns with different length ratio l/ b₁ and pulse duration T_p [35]

The interaction between global buckling mode m = 1 and the first local buckling mode m > 1 was considered. The geometrical imperfection were assumed in the shape corresponding to the buckling mode with amplitude ξ^* ≡ ξ₂^*equals 1/100 wall thickness for local mode and ξ₁^* equal to length to one thousand wall thickness (l/1000h) for global mode. Time of pulse duration T_p was assumed as T₁ equal to the period of natural fundamental flexural vibration or T_m (where m is a number of half waves of local buckling mode) equal to the period of natural vibration with mode corresponding to the local buckling mode.

Figure 15 presents dimensionless deflection ξ as a function of dynamic load factor and maximal radius r_max calculated for maximal characteristic root of Jacoby matrix as a function of DLF.

From curves presented in Figure 15a the critical value of dynamic load factor based on Volmir (V) or Budiansky-Hutchinson (B-H) criterion can be found. The curves presented in Figure 15b help to find critical DLF value based on proposed criterion (PC) [35]. The obtained critical DLF’s according to mentioned above criteria are presented in Table 6.

columns length ratio	pulse duration; Tp [ms]	buckling modes	B-H	V	PC
l/b1 = 4	T1= 1.6	m=1; 3 local	1.0$÷$1.15	1.07	1.16
l/ b1 = 6	T1= 1.9	m=1; 5 local	0.95$÷$1.0	1.01	1.08
l/ b1 = 8	T1= 2.7	m=1; 6 local	1.1$÷$1.15	1.09	1.09
l/ b1 = 4	T3= 0.8	m=1; 3 local	1.6$÷$1.75	1.40	1.58
l/ b1 = 6	T5= 0.7	m=1; 5 local	1.6$÷$1.75	1.39	1.59
l/ b1 = 8	T6= 0.8	m=1; 6 local	1.6$÷$2.05	1.43	1.51

Table 6.

Comparison of DLF_cr obtained from different criteria for interactive buckling

The comparison of the obtained results shows that they are in good agreement. In all cases with pulse duration equal to the period of natural vibration of the form corresponding to local buckling mode the results obtained from the proposed criterion (PC) are between the results obtained from Volmir (V) and Budiansky-Hutchinson (B-H) criteria. For loading with time of pulse duration equal to the period of natural fundamental vibration T₁ the critical DLF values obtained using the proposed criterion (PC) are equal or a bit greater (about 6%) than the critical dynamic load factors from Budiansky-Hutchinson (B-H) or Volmir (V) criteria.

Should be pointed out that in the dynamic buckling problem also for the short columns the multimodal buckling analysis should be carried out. It has been proven on exemplary channel columns with following dimensions: b₁/h = 100, b₂/h = 50, b₃/h = 25 and l/b₁ = 4.

The problem has been calculated with the analytical-numerical method and the finite element method [56].

The dimensionless deflection ξ as a function of dimensionless time (time divided by pulse duration) for channel column is presented in Figure 16. The characteristic points are located in the middle cross-section of the columns and in the middle of the web (point 1 – Figure 16) and on the edge between the web and the flange (point 2 – Figure 16).

Figure 16.

Deflection in time for channel columns subjected to rectangular pulse loading with the DLF = 1.6

Analysing the results presented in Figure 16, it can be said that the global buckling appears for channel cross-section columns. The column edge deflections are greater than deflections of the middle part of the web. The FEM results of calculations presented in Figure 16 have initiated the need for a multimodal buckling analysis also for short columns subjected to pulse loading.

Results for linear buckling and modal analyses obtained with proposed analytical-numerical method are presented in Table 7. As it will be presented below (see Figure 17) the four modes should take into consideration in ANM to obtain similar results to this obtained with FEM. The finite element method gives results (global mode) even in the case when only one buckling mode as the initial imperfection (for example, the local buckling mode m=3) has been taken into account [56].

Mode	σcr [MPa]	n [Hz]
local mode m = 3	53	614
primary local mode m = 1	123	312
secondary local mode m = 1’	972	880
global mode m = 1’’	5122	2001

Table 7.

Buckling stress and natural vibration for the channel column

A comparison between the results obtained with the analytical-numerical method and the finite element method on plots presenting a dimensionless deflection vs. a dynamic load factor for columns with channel cross-sections are shown in Figure 17. Some differences appear because the analytical-numerical model has only a few degrees of freedom in contrary to FE model, which has thousands DoF. However the curves presented in Figure 17 are different the critical DLF values estimated using the Budiasky-Hutchinson criterion is similar. From the ANM, the DLF_cr = 2.7 and from the FEM, the DLF_cr = 2.6.

Good agreement between the results obtained with ANM and FEM is possible because the interactive dynamic buckling problem has been solved in the analytical-numerical method. Four modes have been taken into account. The buckling stress and the natural frequency obtained with the analytical-numerical method for all the modes assumed in the multimodal analysis are listed in Table 7. The buckling modes taken into consideration in the multimodal analysis are presented in Figure 18, correspondingly.

Figure 17.

Dimensionless edge deflection ξ vs. the DLF for channel columns subjected to rectangular pulse loading T_p = T₃ = 1.6 ms [56]

Figure 18.

Buckling modes for the channel cross-section column

5. Conclusion

Taken into consideration the nonlinear thin plate theory for orthotropic material allows, as is presented in exemplary results of calculation, to analyse thin-walled structures composed of flat plates and subjected to static and dynamic load. The nonlinear orthotropic plate theory is the base for the proposed analytical-numerical method which allows to find buckling load with corresponding buckling mode, natural frequencies with corresponding modes and to analyse the postbuckling behaviour – drawing the postbuckling equilibrium paths for plate, segment of girders or columns made of isotropic, orthotropic or even composite materials. As it was shown not only static load can be considered but also dynamic load with intermediate velocity – the dynamic buckling can be analysed using assumed plate theory and the proposed method of solution.

The proposed analytical–numerical method gives almost the same results for eigenvalue problem (buckling loads, natural frequencies with corresponding modes) and similar results for dynamic buckling as the finite element method. However the dimensionless deflection versus dynamic load factor relation obtained with both (proposed and FEM) methods are not identical (especially for higher DLF value). These relations allow to find similar critical value of DLF taking into consideration one of the well-known criterion. The differences in the dimensionless deflection ξ appear because the numerical model in the FEM has more degrees of freedom than the model in the analytical–numerical method, but the results from the ANM are obtained in a significantly faster way than those from the finite element method.

Acknowledgement

This publication is a result of the research work carried out within the project subsidized over the years 2009-2012 from the state funds designated for scientific research (MNiSW - N N501 113636).