Buckling analysis of laminated plate structures with elastic edges using a novel semi-analytical finite strip method

Composite Structures 152 (2016) 85–95

Contents lists available at ScienceDirect

Composite Structures journal homepage: www.elsevier.com/locate/compstruct

Buckling analysis of laminated plate structures with elastic edges using a novel semi-analytical finite strip method Qingyuan Chen a, Pizhong Qiao a,b,⇑ a State Key Laboratory of Ocean Engineering, Collaborative Innovation Center for Advanced Ship and Deep-Sea Exploration, School of Naval Architecture, Ocean and Civil Engineering, Shanghai Jiao Tong University, Shanghai 200240, PR China b Department of Civil and Environmental Engineering, Washington State University, WA 99164-2910, USA

a r t i c l e

i n f o

Article history: Received 10 February 2016 Revised 21 April 2016 Accepted 2 May 2016 Available online 3 May 2016 Keywords: Buckling Finite strip method End boundary conditions Semi-analytical solution

a b s t r a c t A novel semi-analytical finite strip method is presented for buckling analysis of composite plate structures with boundary edges elastically supported. A set of unique Fourier series functions is introduced to represent the longitudinal variation of deflection along a strip, and they are capable of handling elastic edges with translational and rotational spring supports. The proposed hybrid method overcomes limitation of classical finite strip method only capable of handling simple end boundary conditions of structures, and it avoids the ill-conditioning when a set of standard Fourier series functions is used for buckling analysis. Accuracy and validity of the proposed method are demonstrated by the convergence and comparative studies in comparison with the numerical finite element method. As an example, the present method is applied to buckling analysis of a composite Z-stiffened panel under pure shear, and its capability and efficiency of treating different edge conditions in the panel skin and stiffeners are illustrated. Ó 2016 Elsevier Ltd. All rights reserved.

1. Introduction The semi-analytical finite strip method has been proposed by Cheung [1,2] for several decades, and there have been a great amount of developments in the method itself and in its application [3–19]. However, this technique is fraught with limitation of only handling simple end boundary conditions of structures (i.e., both ends simply supported, both ends clamped, one end simply supported and the other end clamped, both ends free, or one end clamped and the other end free). The structures with the end either simply supported, clamped or free are extreme cases; while in reality, the end edges of structures are usually elastically supported or restrained by adjacent components. When it comes to these complicated and more realistic end boundary conditions, many scholars turn to the already developed spline finite strip method [2]. Many researchers have investigated the vibration of plates with elastically-restrained boundary condition edges using series solutions. For example, Li et al. [20] developed an analytical method for vibration analysis of rectangular isotropic plates with elastically-restrained edges, in which the displacement solution

⇑ Corresponding author at: Department of Engineering Mechanics, Shanghai Jiao Tong University, Shanghai 200240, PR China. E-mail addresses: [email protected], [email protected] (P. Qiao). http://dx.doi.org/10.1016/j.compstruct.2016.05.008 0263-8223/Ó 2016 Elsevier Ltd. All rights reserved.

is expressed as a combination of a standard Fourier cosine series and several auxiliary closed-form functions. Beslin and Nicolas [21] proposed a hierarchical set with trigonometric functions for flexural vibration of rectangular isotropic plates. Barrette [22] used the hierarchical trigonometric functions of Beslin and Nicolas [21] as local trial functions in prediction of stiffened plate vibration, and the local functions are defined on the plate domain present between consecutive stiffeners. Then, Dozio [23–25] used this set of trigonometric functions with the Ritz method for general vibration of rectangular isotropic plates and in-plane vibration analysis of isotropic and composite plates, and all the plates considered have the arbitrary elastic boundaries. However, the set of trigonometric functions proposed by Li et al. [20] and Beslin and Nicolas [21] cannot be used for buckling analysis of isotropic and composite rectangular plates with complex and restrained boundary conditions, because of ill-conditioning. While in the authors’ opinion, the incapability of considered trigonometric functions for buckling analysis is caused by many terms in the geometric matrices being zero. In this paper, a novel semi-analytical finite strip method is presented for buckling analysis of laminated composite structures under shear and compressive loading, in which the longitudinal series functions are replaced by the set of trigonometric functions as proposed by Beslin and Nicolas [21]. The proposed novel semianalytical finite strip method retains the advantages of the finite

86

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

strip method (FSM), and at the same time it overcomes the shortcoming of the existing FSM only capable of handling the classical and simple end boundary condition of structures; more important, the hybrid method allows the use of the set of trigonometric functions proposed by Beslin and Nicolas [21] without ill-conditioning. 2. Theoretical formulations Consider a typical finite strip element with the length L and the T

T

width Be (Fig. 1). ky0 and ky1 are the translational (vertical) spring stiffness coefficients along the ends of the finite strip element of R ky0

R ky1

y ¼ 0 and y ¼ L, respectively; and are the rotational spring stiffness coefficients along the ends of the finite strip element of y ¼ 0 and y ¼ L, respectively. Based on the Classical Laminated Plate Theory (CLPT), the displacements of the middle surface of the plate uðx; yÞ, v ðx; yÞ and wðx; yÞ are expressed by the interpolation polynomial function in x-direction and smooth series functions in y-direction:

u ¼

v

¼

w ¼

r X m¼1 r X m¼1 r X

e Ym u fC u gfdgm

e Ym v fC v gfdgm

ð1Þ

fC w g ¼ f 0 0 C 4

m¼1

v im

wim him ukm

v km

wkm hkm ujm

v jm

ð2Þ where him , hjm and hkm are the rotation parameters at the three longitudinal nodal lines, respectively, and they are defined as h ¼ @w=@x. fC u g, fC v g, fC w g are the transverse interpolation shape functions, and they are given by

C7

0 0 C8

;

C9 g

ð3Þ

where g ¼ y=L, and the admissible series function is proposed by Beslin and Nicolas [21] as:

/m ðgÞ ¼ sinðam g þ bm Þ sinðcm g þ dm Þ

ð4Þ

where the coefficients am , bm , cm and dm are listed in Table 1, and the first eight series functions /m ðgÞ defined by Eq. (4) are reported in Table 2. Note that the series functions when m > 4 have zero displacement and zero slope at the ends of the plate strip element; the free deflection and rotation at the end of y ¼ 0 are controlled by the first function /1 ðgÞ and the second function /2 ðgÞ; and the non-zero displacement and slope at the end of y ¼ L are dominated by the third function /3 ðgÞ and the fourth function /4 ðgÞ. The plate stiffness equations are expressed as:

fNg

"

 ¼

½A ½BT ½B

#

feg

 ð5Þ

fjg

½D

with

fNg ¼ fNx Ny Nxy gT ;

T

e0y c0xy g fMg ¼ fM x M y Mxy gT ; fjg ¼ fjx jy jxy gT 2

wjm hjm gT

0 0 C6

m m Ym u ¼ Y v ¼ Y w ¼ /m ðgÞ

fMg

where fdgem is a vector representing the mth term nodal displacement parameters at the nodal lines of the finite strip element. For the low order finite strip with three nodal lines (LO3, see Fig. 1), the following expression is held:

C5

where, C 1 ¼ 1  3x þ 2x2 , C 2 ¼ 4x  4 x2 , C 3 ¼ x þ 2x2 , C 4 ¼ 1 x3  68 x4 þ 24 x5 , C 5 ¼  x  6x2 þ 13x3  12x4 þ 4x5 , C 6 ¼ 23x2 þ 66 16 x2  32 x3 þ 16 x4 , C 7 ¼ 8x2 þ 32x3  40 x4 þ 16x5 , C 8 ¼ 7x2  3 4 5 2 3 4        34x þ 52x  24x , C 9 ¼ x þ 5x  8x þ 4x5 , x ¼ x=Be. m m The longitudinal series functions Y m u , Y v , Y w are defined as:



e Ym w fC w gfdgm

fdgem ¼ fuim

fC u g ¼ f C 1 0 0 0 C 2 0 0 0 C 3 0 0 0 g fC v g ¼ f 0 C 1 0 0 0 C 2 0 0 0 C 3 0 0 g

fe0 g ¼ fe0x

3 A11 A12 A16 6 7 ½A ¼ 4 A12 A22 A26 5; A16 A26 A66 2 3 D11 D12 D16 6 7 ½D ¼ 4 D12 D22 D26 5 D16 D26 D66

2

B11 6 ½B ¼ 4 B12 B16

B12

B16

3

B22

7 B26 5;

B26

B66

ð6Þ

where Nx, Ny and Nxy are, respectively, the membrane transverse, longitudinal and in-plane shear forces per unit length, and e0x , e0y and c0xy are the membrane strains; Mx, My and Mxy are, respectively, the transverse and longitudinal bending and twisting moments per unit length, and jx , jy and jxy are the flexural strains, known as the curvatures. Aij, Bij and Dij are, respectively, the extensional, bendingextension, and bending stiffness coefficients. The strain energy of the plate finite strip element can be expressed as

   ZZ   ½A ½B feg 1 dxdy fegT fjgT 2 fjg ½B ½D ZZ  1 1 ¼ fegT ½Afeg þ jT ½Bfeg þ fjgT ½Dfjg dxdy 2 2

Ue ¼

ð7Þ

Substituting the related variables in Eqs. (1)–(6) into Eq. (7), the strain energy equation is written as Table 1 Coefficients of the series function.

Fig. 1. LO3 plate finite strip element.

m

am

bm

cm

dm

1 2 3 4 >4

p/2 p/2 p/2 p/2

3p/2 p/2 0 0 (m4) p

p/2 p p/2 p p

3p/2 p 0 0

(m4) p

p

87

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

Performing the principle of minimum total potential energy, the governing equation of the plate strip element can be formulated as:

Table 2 The first eight series functions. m

Plot

m

1

Plot

2

@

e Y e

@fdg

¼

@ðU e þ U eS þ V e Þ ¼ ðK e  kP e Þfdge ¼ f0g @fdge

ð13Þ

In Eq. (13), K e and P e are the stiffness and geometric matrices, respectively, and they are defined as: 3

4

K e ¼ ½K ip emn þ ½K io emn þ ½K op emn þ ½K S emn e

e

ð14Þ

e

P e ¼ ½kg ip1 mn þ ½kg ip2 mn þ ½kg op mn T

5

6

7

8

U e ¼ U eip þ U eio þ U eop

where ½K s emn ¼ ½Rs emn þ ½Rs emn , and the other expression of the matrices are given in Appendix A. From the governing equation (Eq. (13)) of the plate strip element, the governing equation of the whole plate structure can be obtained by the assembly of corresponding strip matrices in the standard fashion. 3. Results and discussion

ð8Þ

where the expressions of U eip , U eio and U eop are given in Appendix A. The strain energy stored by the elastic edge restraints is given

1 2

U eS ¼

þ

Z

1 2

T

ky0 ðwjy¼0 Þ2 dx þ Z

R

ky0

1 2

Z

T

ky1 ðwjy¼L Þ2 dx

!2

!2 Z @w

1 @w

R k dx þ dx y1 @y y¼0 2 @y y¼L

r X r X

fdgem ½RS emn fdgen

ð10Þ

Z 





1 T 1 T T m

n

n

k fC w gT Y m ¼ w y¼0 Y w y¼0 fC w g þ ky0 fC w g Y w y¼L Y w y¼L fC w g 2 y0 2



1 R



Y nw;y fC w g þ ky0 fC w gT Y m w;y

y¼0 y¼0 2



1 R

T n

Y fC w g dx þ ky0 fC w g Y m w;y

w;y

y¼L y¼L 2

The loss of the potential by the external load is:

where, V eip1 , V eip2 and V eop are also given in Appendix A.

K Ry0 ¼

ky0 L ; D22

T

K Ty1 ¼

ky1 L3 D22

ð15Þ

R

K Ry1 ¼

ky1 L D22

For the sake of convenience, the dimensionless parameter K Ty ¼ K Ty0 ¼ K Ty1 is introduced if the non-dimensional translational tional spring stiffness K Ty1 (at y = L); and similarly, the dimensionless parameter K Ry ¼ K Ry0 ¼ K Ry1 is introduced if the non-dimensional rotational spring constant K Ry0

(at y = 0) equals the non-

dimensional rotational spring constant K Ry1 (at y = L). Several elastic end boundaries are defined in Table 3. Note that the elastic edge (end) with the dimensionless translational (vertical) stiffness constant of 107 can be treated as the one which cannot move along the normal direction of the plate finite strip element; while the

Table 3 The definition of several elastic end boundaries.

ð11Þ

Substituting the related variables in Eq. (1) into Eq. (11), the loss of the potential is written as

V e ¼ V eip1  V eip2  V eop

ky0 L3 ; D22

spring stiffness K Ty0 (at y = 0) equals the non-dimensional transla-

where,

ZZ " " 2  2  2 # 1 @u @v @w Ve ¼  h r0xx þ þ 2 @x @x @x "  2  2 # 2 @u @v @w þ þ þr0yy @y @y @y   @u @u @ v @ v @w @w þ dxdy þ þ2s0xy @x @y @x @y @x @y

K Ty0 ¼

R

m¼1 n¼1

½RS emn

spring stiffness coefficients K Ty0 , K Ty1 , K Ry0 and K Ry1 are introduced in the present study for the buckling behavior of elasticallysupported (both the translational and rotational supports) composite plate structures, and they are given by: T

ð9Þ

and it is can be expressed in the matrix form

U eS ¼

In this section, buckling analysis of laminated composite plate structures under compression and shear is conducted using the proposed semi-analytical finite strip method, and their results are all compared with the numerical finite element method (FEM) to illustrate their accuracy. The finite element analysis as a comparison and validation tool to the proposed FSM is performed using the commercial software ABAQUS, and the shell element S4R is used. The non-dimensional translational (vertical) and rotational

ð12Þ

Symbol

Case

y¼0

y¼L

S

Simply supported

K Ty0 ¼ 1, K Ry0 ¼ 0

K Ty1 ¼ 1, K Ry1 ¼ 0

C R

Clamped

K Ty1 ¼ 1, K Ry1 ¼ 1

V

Vertical spring

E

Elastic support

F

Free

K Ty0 ¼ 1, K Ry0 ¼ 1 K Ty0 ¼ 1, K Ry0 –0 K Ty0 –0, K Ry0 ¼ 0 K Ty0 –0, K Ry0 –0 K Ty0 ¼ 0, K Ry0 ¼ 0

Rotational spring

End of the plate strip element

K Ty1 ¼ 1, K Ry1 –0 K Ty1 –0, K Ry1 ¼ 0 K Ty1 –0, K Ry1 –0 K Ty1 ¼ 0, K Ry1 ¼ 0

Note: the bar under S or C represents the edge at either y = 0 or y = L.

88

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

end with the dimensionless translational (vertical) stiffness constant of 107 combined with the dimensionless rotational stiffness constant of 104 can be treated as the clamped support. These extreme boundary conditions will be demonstrated in the following examples and discussions. 3.1. Convergence and comparison study for isotropic plates In order to illustrate feasibility of the present semi-analytical finite strip method, the convergence and comparison studies for the buckling behavior of isotropic plates are first conducted. The square aluminum alloy plate with all the four edges either simply supported or clamped are considered, and the material and geometric properties of aluminum alloy are: the modulus of elasticity E = 69 GPa, the Poisson’s ratio m = 0.33, and the width and thickness of the aluminum alloy plate are B = 0.2 m and t = 0.002 m, respectively. In the present semi-analytical finite strip analysis, the simply supported edges can be conceptually considered as the special case when the stiffness constants for the translational springs at the ends of plate become infinitely large (K Ty0 ¼ K Ty1 ¼ 1, and in this study, the finite values of vertical spring stiffness are assigned K Ty0 ¼ K Ty1 ¼ 107 to represent the case of simply supported edges), while the stiffness constants for the rotational springs at the ends of the finite strip element are zero (K Ry0 ¼ K Ry1 ¼ 0). The clamped edges can be considered as another special case in which both the stiffness constants for the translational and the rotationally springs are infinitely large (K Ty0 ¼ K Ty1 ¼ 1 and K Ry0 ¼ K Ry1 ¼ 1; and in this study, the finite values are assigned K Ty0 ¼ K Ty1 ¼ 107 and K Ry0 ¼ K Ry1 ¼ 104 to represent the case of clamped edges). The buckling loads of the square aluminum alloy plate calculated using different number of strips and different number of terms in the series expansion are shown in Tables 4 and 5, where the results of exact solution and numerical finite element method are also given. The exact solution of the isotropic plates is given in the form

Ncr ¼ k

p2 D b

ð16Þ

2

where k is the weight coefficient and D = Et3/(1m2). A total of 1600 equal length shell elements of S4R are used in the finite element method (FEM)-based analysis. Based on the results shown in Tables 4 and 5, the conclusion can be drawn that nine and ten approximate terms of the series expansion in two LO3 strips for the whole plate should be used in the present FSM for buckling of SSSS and CCCC (where the under-bar represents the boundary condition at either y = 0 or y = L) isotropic plates, respectively, under uniaxial compression or shear. 3.2. Convergence and comparison study for laminates The symmetrically laminated anisotropic plates are then considered for the convergence and compassion studies, and all the laminates are constructed of the carbon/epoxy T800-3900-2. The material properties of the lamina are [26]: longitudinal modulus EL = 155.8 GPa, transverse modulus ET = 8.89 GPa, shear modulus GLT = 5.14 GPa, and Poisson’s ratio m = 0.3. The width of laminates is 0.5 m, and the ply thickness of the lamina is 0.125 mm. The shell element with the size of 0.0125 m  0.0125 m is used for the results obtained by employing the commercial FEM software ABAQUS. The convergence study of clamped laminates with different pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi material orthotropic parameters 4 D11 =D22 are investigated. When pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi the material orthotropic parameter 4 D11 =D22 is small (e.g., for the ½ð45 Þ3 =ð0 Þ2 S laminate), the square symmetrically laminated plate buckles in a symmetric mode under shear (see Fig. 2); when pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi the material orthotropic parameter 4 D11 =D22 is large enough (e.g., for the ½ð15 Þ3 =ð0 Þ2 S laminate), the square symmetrically laminated plate buckles in an antisymmetric mode under shear (see Fig. 3). The clamped edges are regarded as that both the stiffness constants for the translational spring and rotational springs on

Table 4 Convergence study of buckling loads for a square thin plate with four edges simply supported (SSSS). Terms

8 9 10 11 12 20

Compression (KN/m) Number of strips

Shear (KN/m) Number of strips

1

2

3

4

1

2

3

4

50.97 50.96 50.96 50.96 50.96 50.96

50.97 50.95 50.95 50.95 50.95 50.95

50.97 50.95 50.95 50.95 50.95 50.95

50.97 50.95 50.95 50.95 50.95 50.95

122.6 122.3 122.3 122.2 122.2 122.1

119.3 119.0 119.0 118.9 118.9 118.8

119.3 119.0 119.0 118.9 118.9 118.8

119.3 119.0 119.0 118.9 118.9 118.8

50.95 (4.00 p2D/b2) 50.79

Exact FEM

119.0 (9.34 p2D/b2) 11.88

Table 5 Convergence study of buckling loads for a square thin plate with four edges clamped (CCCC). Terms

8 9 10 11 12 20 Exact FEM

Compression (KN/m) Number of strips

Shear (KN/m) Number of strips

1

2

3

4

1

2

3

4

129.0 128.8 128.8 128.7 128.7 128.6

128.8 128.5 128.5 128.4 128.4 128.3

128.8 128.5 128.5 128.4 128.4 128.3

128.8 128.5 128.5 128.4 128.4 128.3

213.2 212.6 211.6 211.4 211.0 210.6

189.3 188.0 187.4 187.1 186.9 186.6

189.2 187.9 187.4 187.1 186.8 186.5

189.2 187.9 187.4 187.0 186.9 186.5

128.3 (10.07 p2D/b2) 128.6

187.4 (14.71 p2D/b2) 187.5

89

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

Fig. 2. The buckling modes of CCCC ½ð45 Þ3 =ð0 Þ2 S under (a) transverse compression; (b) longitudinal compression; (c) positive shear; (d) negative shear.

Fig. 3. The buckling mode of CCCC ½ð15 Þ3 =ð0 Þ2 S under (a) transverse compression; (b) longitudinal compression; (c) positive shear; (d) negative shear.

Table 6 Convergence study of buckling loads for the ½ð45 Þ3 =ð0 Þ2 S laminate with four edges clamped (CCCC) under compression. Terms

8 9 10 11 12 20 FEM

Longitudinal compression (KN/m) Number of strips

Transverse compression (KN/m) Number of strips

2

3

4

2

3

4

16.87 16.83 16.83 16.82 16.82 16.81

16.86 16.82 16.82 16.80 16.80 16.79

16.86 16.82 16.81 16.80 16.80 16.79

17.13 16.99 16.97 16.92 16.91 16.86

17.14 16.99 16.97 16.92 16.91 16.85

17.13 16.99 16.97 16.92 16.91 16.85

16.84

Table 7 Convergence study of buckling loads for the ½ð45 Þ3 =ð0 Þ2 S laminate with four edges clamped (CCCC) under shear. Terms

8 9 10 11 12 20 FEM

Positive shear (KN/m) Number of strips

Negative shear (KN/m) Number of strips

2

3

4

2

3

4

24.21 23.97 23.89 23.82 23.79 23.70

24.20 23.96 23.88 23.81 23.78 23.69

24.20 23.96 23.88 23.81 23.77 23.69

30.51 30.10 30.22 30.19 30.12 30.03

30.50 30.39 30.20 30.18 30.11 30.02

30.50 30.39 30.20 30.18 30.10 30.01

23.83

30.17

16.89

the ends of the strip are infinitely large (K Ty0 ¼ K Ty1 ¼ 1 and K Ry0 ¼ K Ry1 ¼ 1; and in this study, the finite values are assigned K Ty0 ¼ K Ty1 ¼ 107 and K Ry0 ¼ K Ry1 ¼ 104 ) in the present finite strip method to represent the clamped conditions. The buckling modes of the ½ð45 Þ3 =ð0 Þ2 S and ½ð15 Þ3 =ð0 Þ2 S laminates are shown in Figs. 2 and 3, respectively. The buckling loads of the ½ð45 Þ3 =ð0 Þ2 S laminate are for all the edges clamped under compression and shear, respectively. The buckling modes of the ½ð45 Þ3 =ð0 Þ2 S and ½ð15 Þ3 =ð0 Þ2 S laminates under transverse compression (compression along x-axis) are the same; while the buckling behavior of laminates with more than one half waves along y-axis under compression is shown for the ½ð15 Þ3 =ð0 Þ2 S laminate, due to its large material orthotropic parameter p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4 D11 =D22 . The convergence of the buckling loads for the ½ð15 Þ3 =ð0 Þ2 S laminate under shear are presented in Table 6. As shown in Table 6, it can be seen that 9 and 11 approximate terms of the series expansion with two LO3 strips should be used in the present FSM for the case of all the edges clamped under compression along y-axis and x-axis, respectively. As demonstrated in Tables 7, 8 and 11 approximate terms of the series expansion with

Table 8 Convergence study of buckling loads for the square laminate ½ð15 Þ3 =ð0 Þ2 S with four edges clamped (CCCC) under shear. Terms

8 9 10 11 12 20 FEM

Positive shear (KN/m) Number of strips

Negative shear (KN/m) Number of strips

2

3

4

2

3

4

18.41 18.37 18.26 18.25 18.22 18.17

18.40 18.36 18.26 18.24 18.21 18.17

18.40 18.36 18.26 18.24 18.21 18.17

21.78 21.69 21.61 21.58 21.55 21.51

21.77 21.68 21.60 21.57 21.55 21.50

21.77 21.68 21.60 21.57 21.54 21.50

18.32

21.66

two LO3 strips should be used in the present FSM for the case of all the edges clamped under shear. The next convergence study is about different aspect ratio laminates with all the edges clamped. Like the above analysis, the clamped edges are considered when both the stiffness constants for the translational and the rotational springs on the ends of plates T

T

R

R

are infinitely large (ky0 ¼ ky1 ¼ 1 and ky0 ¼ ky1 ¼ 1; in this study,

90

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

the finite values are assigned K Ty0 ¼ K Ty1 ¼ 107 and K Ry0 ¼ K Ry1 ¼ 104 for the case of clamped edges) in the present semi-analytical finite strip method to represent the clamped boundary condition. The buckling modes of different aspect ratios of the ½ð45 Þ3 =ð0 Þ2 S laminate with all the edges clamped under uniaxial compression along y-axis and shear are shown in Figs. 4 and 5, respectively.

Based on the above analysis, it is concluded that two LO3 strips is adequate for the buckling analysis of different buckling modes under compression or shear. In Tables 9 and 10, the buckling loads for different aspect ratios of the ½ð45 Þ3 =ð0 Þ2 S laminate under the longitudinal compression and shear are presented against the number of terms in the series expansion. It is shown that the

Fig. 4. The buckling mode of CCCC ½ð45 Þ3 =ð0 Þ2 S with various aspect ratios under longitudinal compression.

Fig. 5. The buckling modes of CCCC ½ð45 Þ3 =ð0 Þ2 S with various aspect ratios under pure shear.

91

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95 Table 9 Convergence study of buckling loads for different aspect ratios of the ½ð45 Þ3 =ð0 Þ2 S laminate with four edges clamped (CCCC) under longitudinal compression. Terms

Longitudinal compression (KN/m) Aspect ratio 2/3

1/2

1/3

2/7

1/4

1/5

2/11

1/6

8 9 10 11 12 13 14 15 16 18 20

14.61 14.60 14.59 14.59 14.58 14.58 14.58 14.58 14.58 14.58 14.58

13.98 13.97 13.96 13.96 13.95 13.95 13.95 13.95 13.95 13.95 13.95

13.55 13.48 13.40 13.40 13.38 13.38 13.38 13.38 13.38 13.38 13.38

– 13.35 13.28 13.25 13.25 13.25 13.25 13.24 1.3.24 13.24 13.24

– 13.54 13.23 13.19 13.18 13.17 13.17 13.17 13.17 13.17 13.17

– – – – – 13.09 13.09 13.08 13.08 13.08 13.08

– – – – 13.11 13.07 13.05 13.05 13.04 13.04 13.04

– – – – – 13.07 13.03 13.02 13.02 13.02 13.02

14.59

13.96

13.39

13.25

13.18

13.09

13.05

13.03

FEM

Table 10 Convergence study of buckling loads for different aspect ratios of the ½ð45 Þ3 =ð0 Þ2 S laminate with four edges clamped (CCCC) under shear. (a) Positive shear Terms

8 9 10 11 12 13 14 15 16 18 20 FEM

Positive shear (KN/m) Aspect ratio 2

3

4

5

6

7

16.62 16.55 16.46 16.43 16.40 16.39 16.38 16.37 16.36 16.35 16.35

– 15.18 15.15 15.10 15.09 15.07 15.05 15.05 15.05 15.04 15.04

– – 14.66 14.65 14.62 14.62 14.61 14.60 14.60 14.59 14.59

– – – – – 14.41 14.40 14.39 14.39 14.39 14.38

– – – – – – 14.29 14.29 14.28 14.27 14.27

– – – – – – – 14.22 14.22 14.21 14.21

16.42

15.11

14.66

14.45

14.34

14.28

(b) Negative shear Terms

8 9 10 11 12 13 14 15 16 18 20 FEM

Negative shear (KN/m) Aspect ratio 2

3

4

5

6

7

21.07 20.88 20.85 20.77 20.76 20.73 20.73 20.71 20.70 20.69 20.69

– 19.07 19.00 18.99 18.97 18.97 18.96 18.96 18.95 18.95 18.94

– – 18.41 18.37 18.36 18.35 18.35 18.34 18.34 18.34 18.33

– – – – – 18.08 18.07 18.07 18.07 18.06 18.06

– – – – – – 17.93 17.92 17.92 17.92 17.92

– – – – – – – 17.84 17.84 17.83 17.83

20.77

19.03

18.42

18.14

18.00

17.91

convergence values are obtained with 10 terms in the series expansion when the number of half waves is less than 5, and the convergence values can be obtained with 11, 12, 13, 14, and 14 terms when the number of half waves is, respectively, 5, 6, 7, 8, and 9 for the considered laminate under the longitudinal compression. While for the laminates under shear, the convergence values are obtained with 11 terms in the series expansion when the number of half waves is less than 5, and the convergence values can be obtained with 13, 14, and 15 terms when the number of half waves is 5, 6 and 7, respectively.

Finally, the effects of translational (vertical) and rotational spring (restraint) stiffness on the shear and compressive buckling behavior of the restrained square symmetrically-laminated anisotropic plates with the ends of plate restrained against the respective translation and rotation are presented in Figs. 6–11. The effect of translational (vertical) restraint stiffness (K Ty ) on the shear and compressive buckling behavior of laminate ½ð45 Þ3 =ð0 Þ2 S with two opposite simply supported edges in the y-axis and the other two opposite ends vertically elastically restrained (SVSV) under pure shear and uniaxial compression along x-axis are shown in Figs. 6 and 7, respectively, where the shear or compressive buckling loads are plotted against the base 10 logarithm of the dimensionless vertical restrained stiffness K Ty . As shown in Figs. 6 and 7, it can be seen that both the shear and transverse compression buckling loads are proportional to the dimensionless vertical restraint stiffness K Ty , and there are almost no changes about the shear and compressive loads when K Ty is in the range of 106.5–107 (i.e., it approaches to the case of the simply supported boundary condition). Thus, the ends with the dimensionless vertical restraint stiffness of 107 or higher can be treated as the simply supported ends. In addition, it can also be seen that the difference between the positive and the negative shear buckling loads becomes larger with increase of vertical restraint stiffness K Ty . The effect of vertical restraint stiffness (K Ty1 ) on the shear and compressive buckling behavior of ½ð45 Þ3 =ð0 Þ2 S laminate with two opposite simply supported edges in the y-axis and the end at y ¼ 0 simply supported and the other end at y ¼ L vertically elastically restrained (SSSV) under pure shear and uniaxial compression along x-axis is evaluated (Figs. 8 and 9). The same conclusions are reached as those from Figs. 6 and 7. In Chen and Qiao [27], the shear buckling of composite laminated plates with all four edges elasticallyrestrained against rotation was analyzed, and there are almost no changes about the shear buckling load when the rotational restraint stiffness K Ry is in the range of 103.5–104. Thus, the edges with the dimensionless rotational restraint stiffness of 104 or higher can be treated as the clamped support. In this study, the effect of the rotational restraint stiffness (K Ry1 ) on the shear and compressive buckling behavior of the ½ð45 Þ3 =ð0 Þ2 S laminate with two opposite clamped edges in the y-axis and the end at y ¼ 0 clamped and the other end rotationally elastically restrained (R at y = L) (i.e., CCCR) under pure shear and uniaxial compression along x-axis is investigated, and the results are shown in Figs. 10 and 11. It can be seen that the same conclusions are reached as

92

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

Fig. 6. Buckling load vs. the vertical restraint stiffness parameters of K Ty for symmetrical anisotropic SVSV laminates under shear.

Fig. 7. Buckling load vs. the vertical restraint stiffness parameters of K Ty for symmetrical anisotropic SVSV laminates under compression.

Fig. 8. Buckling load vs. the rotational restraint stiffness parameters of K Ty1 for symmetrical anisotropic SSSV laminates under shear.

Fig. 9. Buckling load vs. the rotational restraint stiffness parameters of K Ty1 for symmetrical anisotropic SSSV laminates under compression.

Fig. 10. Buckling load vs. the rotational restraint stiffness parameters of K Ry1 for symmetrical anisotropic CCCR laminates under shear.

Fig. 11. Buckling load vs. the rotational restraint stiffness parameters of K Ry1 for symmetrical anisotropic CCCR laminates under compression.

93

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

those given in [27], but the difference between the positive and negative shear buckling loads is not affected much with the increase of the rotational restraint stiffness K Ry1 (see Fig. 10). 3.3. Application to the stiffened laminated panels As an application, the present hybrid finite strip method is used to analyze buckling of a stiffened laminated panel with different end boundary conditions. The mixed boundary conditions at the ends of the composite stiffened panel are considered, and the boundary conditions at the ends of the plate structure can be controlled by the stiffness constants K Ty0 , K Ty1 , K Ry0 and K Ry1 of every plate strip element. The composite Z-stiffened panel considered have a square planform and six equally spaced stiffeners, and the crosssection and a repeating element of the composite Z-stiffened panel are shown in Fig. 12. All the plate elements are constructed of the

graphite-epoxy, and the material properties of the lamina are [28]: longitudinal modulus EL = 131 GPa, transverse modulus ET = 13.0 GPa, shear modulus GLT = 6.41 GPa, Poisson’s ratio m = 0.38. The plate elements with numbers of 1 and 2 are the ½45 t1 =  45 t1 = 45 t1 =45 t1 =0 t2 =90 t3 S laminate, where t 1 ¼ 0:1618 mm, t2 ¼ 0:6325 mm and t3 ¼ 1:0566 mm; the plate elements with numbers of 3 and 4 are the ½45 t3 =  45 t3 =  45 t3 =45 t3 =0 t4 S laminate, where t3 ¼ 0:209 mm and t4 ¼ 1:7145 mm. The longitudinal edges of the Z-stiffened panel are simply supported, and there are two scenarios for the boundary condition at its transverse edges: one is with the simply supported boundary condition to the stiffeners as well as to the skin (Case A); and the other is with the simply supported boundary condition only to the skin and the free boundary condition on the stiffeners (Case B). When the shear load is applied on the Z-stiffened panel, the shear load is carried only by the main (skin) panel for both of the cases. The buckling modes of both the cases are shown in Fig. 13, and the shear buckling loads of both the cases are illustrated in Table 11 along with the FEM results. Note that the shell element size of approximate 0.016 is used in the finite element analysis. As expected, the ends of the stiffeners are distorted in the buckling mode of the composite Z-stiffened panel in Case B; while the ends of the stiffeners maintain the status quo in Case A. As expected, the shear buckling load of the former (Case A) is greater than that of the latter (Case B).

Table 11 Buckling loads of composite Z-stiffened panels under pure shear. Case

Fig. 12. The cross-section and details of the repeating elements of Z-stiffened panel (dimensions in millimeters).

A B

Shear buckling loads (KN/m) Present approach

FEM

1387 1315

1423 1343

Fig. 13. Buckling mode shapes for the composite Z-stiffened panel under pure shear (a) case A: the simply supported boundary conditions are applied to the stiffeners as well as to the skin; (b) case B: the simply supported boundary conditions are only applied to the skin, and the free boundary conditions are with the stiffeners.

94

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95

4. Conclusions

V eip1

In this study, a novel semi-analytical finite strip method is presented, in which the longitudinal strip series functions is introduced by a set of admissible trigonometric functions to deal with the complex end boundary conditions of composite plate structures, and the ill-conditioning existed in the Ritz method is avoided. It can be seen that the bucking loads of composite plate structures with elastic edge supports obtained from the proposed hybrid method agree closely with those of exact and FEM solutions. Also, as demonstrated in the convergence study, only a limited number of finite strips and lower number of approximate terms in the expansion series are needed to provide accurate and efficient buckling analysis of isotropic and anisotropic laminated plates under compression or shear loading. As illustrated in the application example to the composite Z-stiffened panels with different edge boundary conditions in skin and stiffeners, the present hybrid finite strip method is capable of effectively and efficiently analyzing the buckling behavior of stiffened composite panels with different edge boundary conditions (including the rotational and translational restraints at edges) and subject to shear or compression.

ZZ "

1 ¼ h 2 ¼

0 xx

r

r X r X fdgem ½Rg ip1 mn fdgen

ZZ

r X r X 1 fegT ½Afegdxdy ¼ fdgem ½Rip emn fdgen 2 m¼1 n¼1

U eip ¼

ðA:1Þ

where,

½Rg ip1 emn ¼

ZZ

r0xx fC u gT;x Y mu Y nu fC u g;x þ r0yy fC u gT Y mu;y Y nu;y fC u g

V eip2

ZZ "

1 ¼ h 2 ¼

0 xx

r

#  2  2 @v @v 0 0 @v @v dxdy þ ryy þ 2sxy @x @y @x @y

r X r X fdgem ½Rg ip2 mn fdgen

ðA:5Þ

m¼1 n¼1

where,

½Rg ip2 emn ¼

1 h 2

ZZ

r0xx fC v gT;x Y mv Y nv fC v g;x þ r0yy fC v gT Y mv ;y Y nv ;y fC v g

n þ2r0xy fC v gT;x Y m v Y v ;y fC v g dxdy:

V eop ¼ ¼

1 h 2

ZZ "

r0xx

r X r X

#  2  2 @w @w @w @w dxdy þ r0yy þ 2s0xy @x @y @x @y

fdgem ½Rg op emn fdgen

ðA:6Þ

m¼1 n¼1

where,

½Rip emn

½Rg op emn ¼

ZZ  1 1 T m n n A11 fC u gT;x Y m ¼ u Y u fC u g;x þ A22 fC v g Y v ;y Y v ;y fC v g 2 2 T m n n þA12 fC u gT;x Y m u Y v ;y fC v g þ A66 fC u g Y u;y Y v fC v g;x

1 T m n n þ A66 fC u gT Y m u;y Y u;y fC u g þ fC v g;y Y v Y v fC v g;x 2 T m n n þA16 ðfC u gT;x Y m u Y u;y fC u g þ fC u g;x Y u Y v fC v g;x Þ

T m n n þA26 fC v gT Y m v ;y Y u;y fC u g þ fC v g Y v ;y Y v fC v g;x dxdy

1 h 2

ZZ

r0xx fC w gT;x Y mw Y nw fC w g;x þ r0yy fC w gT Y mw;y Y nw;y fC w g

n þ2r0xy fC w gT;x Y m w Y w;y fC w g dxdy:

2. The parts of stiffness and geometric matrices T

½K ip emn ¼ ½Rip emn þ ½Rip emn T

½K io emn ¼ ½Rio emn þ ½Rio emn T

ZZ

r X r X fjg ½Bfegdxdy ¼ fdgem ½Rio emn fdgen T

¼

½K op emn ¼ ½Rop emn þ ½Rop emn ðA:2Þ

m¼1 n¼1

where,

1 h 2

n þ2r0xy fC u gT;x Y m u Y u;y fC u g dxdy:

where,

U eio

ðA:4Þ

m¼1 n¼1

Appendix A 1. The expressions about the parts of the strain energy and the loss of the potential of the plate finite strip element

#  2  2 @u @u 0 0 @u @u dxdy þ ryy þ 2rxy @x @y @x @y

T

½Kg ip1 emn ¼ ½Rg ip1 emn þ ½Rg ip1 emn

ðA:7Þ

T

½Kg ip2 emn ¼ ½Rg ip2 emn þ ½Rg ip2 emn ZZ h T m n n B11 fC w gT;x Y m w Y u fC u g;x þ B22 fC w g Y w;y Y v ;y fC v g

T n n m þB12 fC w gT;x Y m w Y v ;y fC v g þ fC w g Y w;y Y u fC u g;x

T

½Kg op emn ¼ ½Rop emn þ ½Rop emn

½Rio emn ¼ 

References

T n m þB66 ð2fC w gT;x Y nw;y Y m u;y fC u g þ 2fC w g;x Y w;y Y v fC v g;x Þ T n m þB16 2fC w gT;x Y nw;y Y m u fC u g;x þ fC w g;x Y w Y v fC v g;x

T n m þfC w gT;x Y nw Y m u;y fC u g þ B26 2fC w g;x Y w;y Y v ;y fC v g

i T n m þfC w gT Y nw;y Y m v fC v g;x þ fC w g Y w;y Y u;y fC u g dxdy:

ZZ U eop ¼ where,

½Rop emn ¼

r X r X 1 fjgT ½Dfjgdxdy ¼ fdgem ½Rop mn fdgen 2 m¼1 n¼1

ðA:3Þ

ZZ  1 1 T m n n D11 fC w gT;x Y m w Y w fC w g;x þ D22 fC w g Y w;y Y w;y fC w g 2 2

T m n n þD12 fC w gT;x Y m w Y w;y fC w g þ 2D66 fC w g;x Y w;y Y w;y fC w g;x

T m n n þ2D16 fC w gT;x Y m w;y Y w fC w g;x þ 2D26 fC w g Y w;y Y w;y fC w g;x dxdy

[1] Cheung YK. Finite strip method in structural analysis. Publication of Pergamon Press Incorporated; 1976. [2] Cheung YK, Tham L. The finite strip method. CRC Press; 1997. [3] Dawe D, Lam S, Azizian Z. Non-linear finite strip analysis of rectangular laminates under end shortening, using classical plate theory. Int J Numer Methods Eng 1992;35:1087–110. [4] Lam S, Dawe D, Azizian Z. Non-linear analysis of rectangular laminates under end shortening, using shear deformation plate theory. Int J Numer Methods Eng 1993;36:1045–64. [5] Dawe DJ. Use of the finite strip method in predicting the behaviour of composite laminated structures. Compos Struct 2002;57:11–36. [6] Ovesy HR, GhannadPour SAM. Geometric non-linear analysis of imperfect composite laminated plates, under end shortening and pressure loading, using finite strip method. Compos Struct 2006;75:100–5. [7] Amoushahi H, Azhari M. Buckling of composite FRP structural plates using the complex finite strip method. Compos Struct 2009;90:92–9. [8] Milašinovic´ DD. Geometric non-linear analysis of thin plate structures using the harmonic coupled finite strip method. Thin-Walled Struct 2011;49:280–90. [9] Ovesy H, Assaee H, Hajikazemi M. Post-buckling of thick symmetric laminated plates under end-shortening and normal pressure using semi-energy finite strip method. Comput Struct 2011;89:724–32.

Q. Chen, P. Qiao / Composite Structures 152 (2016) 85–95 [10] Assaee H, Hajikazemi M, Ovesy H. The effect of anisotropy on post-buckling behavior of laminated plates using semi-energy finite strip method. Compos Struct 2012;94:1880–5. [11] Ghannadpour S, Ovesy H, Nassirnia M. Buckling analysis of functionally graded plates under thermal loadings using the finite strip method. Comput Struct 2012;108:93–9. [12] Ovesy H, Fazilati J. Buckling and free vibration finite strip analysis of composite plates with cutout based on two different modeling approaches. Compos Struct 2012;94:1250–8. [13] Bui HC. Semi-analytical finite strip method based on the shallow shell theory in buckling analysis of cold-formed sections. Thin-Walled Struct 2012;50:141–6. [14] Hancock GJ, Pham CH. Shear buckling of channel sections with simply supported ends using the semi-analytical finite strip method. Thin-Walled Struct 2013;71:72–80. [15] Chen QY, Qiao PZ. Post-buckling analysis of composite plates under combined compression and shear loading using finite strip method. Finite Elem Anal Des 2014;83:33–42. [16] Sarrami-Foroushani S, Azhari M. Nonlocal vibration and buckling analysis of single and multi-layered graphene sheets using finite strip method including van der Waals effects. Physica E 2014;57:83–95. [17] Ghannadpour SAM, Ovesy HR, Zia-Dehkordi E. Buckling and post-buckling behaviour of moderately thick plates using an exact finite strip. Comput Struct 2015;147:172–80. [18] Ovesy HR, Totounferoush A, Ghannadpour SAM. Dynamic buckling analysis of delaminated composite plates using semi-analytical finite strip method. J Sound Vib 2015;343:131–43.

95

[19] Hancock GJ, Cao HP. Buckling analysis of thin-walled sections under localised loading using the semi-analytical finite strip method. Thin-Walled Struct 2015;86:35–46. [20] Li WL, Zhang XF, Du JT, Liu ZG. An exact series solution for the transverse vibration of rectangular plates with general elastic boundary supports. J Sound Vib 2009;321:254–69. [21] Beslin O, Nicolas J. A hierarchical functions set for predicting very high order plate bending modes with any boundary conditions. J Sound Vib 1997;202:633–55. [22] Barrette M, Berry A, Beslin O. Vibration of stiffened plates using hierarchical trigonometric functions. J Sound Vib 2000;235:727–47. [23] Dozio L. Free in-plane vibration analysis of rectangular plates with arbitrary elastic boundaries. Mech Res Commun 2010;37:627–35. [24] Dozio L. In-plane free vibrations of single-layer and symmetrically laminated rectangular composite plates. Compos Struct 2011;93:1787–800. [25] Dozio L. On the use of the Trigonometric Ritz method for general vibration analysis of rectangular Kirchhoff plates. Thin-Walled Struct 2011;49:129–44. [26] Barbero EJ. Introduction to composite materials design. CRC Press; 2010. [27] Chen QY, Qiao PZ. Shear buckling of rotationally-restrained composite laminated plates. Thin-Walled Struct 2015;94:147–54. [28] Stroud WJ, Greene WH, Anderson MS. Buckling loads of stiffened panels subjected to combined longitudinal compression and shear: results obtained with PASCO, EAL and STAGS computer programs. TP 2215, NASA; 1984.