- Email: [email protected]
Nonlinear coupled zigzag theory for buckling of hybrid piezoelectric plates
Composite Structures 74 (2006) 253–264 www.elsevier.com/locate/compstruct
Nonlinear coupled zigzag theory for buckling of hybrid piezoelectric plates S. Kapuria *, G.G.S. Achary Applied Mechanics Department, Indian Institute of Technology Delhi, Hauz Khas, New Delhi 110016, India Available online 6 June 2005
Abstract A coupled zigzag theory for hybrid piezoelectric plates, recently developed by the first author, is extended to include geometric nonlinearity in the Von Karman sense. In this theory, the potential field is approximated as piecewise linear across sublayers. The deflection approximation accounts for the transverse normal strain due to an electric field. The inplane displacements are assumed to follow a global third order variation across the thickness with a layerwise linear variation. The shear continuity conditions at the layer interfaces and the shear traction-free conditions at the top and bottom are enforced to formulate the theory in terms of only five primary displacement variables, independent of the number of layers. The coupled nonlinear equations of equilibrium and the boundary conditions are derived from a variational principle. The nonlinear theory is used to obtain the initial buckling response of symmetrically laminated hybrid plates under inplane electromechanical loading. Analytical solutions for buckling of simply supported plates under uniaxial and bi-axial inplane strains and electric potential are obtained for comparing the results with the available exact three-dimensional piezoelasticity solution. The comparison establishes that the present theory is very accurate for buckling response of hybrid plates even for highly inhomogeneous lay-ups. 2005 Elsevier Ltd. All rights reserved. Keywords: Buckling; Zigzag theory; Electromechanical coupling; Geometric nonlinearity; Hybrid plate
1. Introduction Laminated composite and sandwich structures offer superior replacement of their metallic counterparts in aerospace, aeronautical, automotive and other applications due to their superior strength to weight and stiffness to weight ratios, high fatigue life as well as high temperature resistance. Distributed piezoelectric sensors and actuators are surface-mounted or embedded in these laminates to introduce adaptive capabilities for active shape control, vibration suppression, acoustic control and buckling control, among others. A considerable amount of research has been dedicated to the study of buckling and post-buckling behaviour of composite *
Corresponding author. Tel.: +91 112 659 1218; fax: +91 112 658 1119. E-mail address: [email protected] (S. Kapuria). 0263-8223/$ - see front matter 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruct.2005.04.010
and sandwich plates [1,2]. Some of the recent works include the static and dynamic buckling analyses of composite plates using first order shear deformation theory (FSDT) [3–5] and refined third order theory (TOT) [5– 7], and of sandwich plates using layerwise models [8,9]. Three-dimensional (3D) elasticity solutions [10,11] of buckling of elastic laminated plates have been presented for simply supported boundary conditions. These 3D solutions serve as useful benchmarks for evaluation of various 2D plate theories. Study of buckling of hybrid laminated structures with surface-bonded or embedded piezoelectric layers and enhancement of their buckling capacity by proper actuation began in the last decade [12]. Batra and Geng [13] presented a numerical analysis for dynamic buckling of a plate with surface bonded piezocermic elements, based on 3D elasticity without considering two-way electromechanical coupling (direct piezoelectric effect). An exact 3D coupled piezoelasticity
254
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
solution for buckling of simply supported symmetrically laminated hybrid plates was presented by Kapuria and Achary [14], wherein the effect of piezoelectric coupling was demonstrated. 2D coupled classical plate theory [15] and first order theory [16] have been applied to study the buckling response of piezoelectric plates. Chandrashekhra and Bhatia [17] presented a finite element model based on uncoupled FSDT for active buckling control of hybrid composite plates. Varelis and Saravanos [18,19] developed an eight noded coupled nonlinear laminated piezoelectric plate element based on FSDT approximations for the displacement field and layerwise approximation for the electric potential and studied initial buckling and post-buckling behaviour of hybrid plates. FSDT requires arbitrary shear correction factors. Correia et al. [20] employed coupled TOT to study optimum design of hybrid composite plates for buckling. From the evaluation of 2D theories for hybrid plates in comparison with 3D exact solutions [21–23], it is now well understood that for accurate global and local lamina level response of these structures (1) the kinematic field approximations must account for the layerwise nature of distributions of the inplane displacements across the thickness, (2) the electric field must allow for the nonlinear distribution of electric potential across the piezoelectric layers and (3) the two-way electromechanical coupling should be considered. FSDT and TOT are equivalent single layer theories, which use the same global variations for displacements across the entire laminate thickness and cannot account for the zigzag nature of variation of the inplane displacements. The first author [23] has presented a coupled efficient layerwise (zigzag) theory for linear static analysis of hybrid plates, which considers layerwise variations for the displacements, but the number of primary displacement variables is reduced to only five (as in FSDT and TOT), by enforcing the shear traction-free conditions at the top and the bottom surfaces and the transverse shear stress continuity conditions at the layer interfaces. Comparison of results of this theory with the exact 3D solutions for simply supported hybrid plates of highly inhomogeneous lay-ups established the high accuracy of this theory for static electromechanical response. A nonlinear coupled 1D zigzag theory for buckling of hybrid beams has also been presented and tested for accuracy [24]. This work presents an efficient coupled geometrically nonlinear zigzag theory for hybrid plates by extending the linear zigzag theory of Kapuria [23]. The geometric nonlinearity is included due to deflection only in the sense of Von Karman. The potential field is approximated as piecewise linear across a number of subdivisions in the piezoelectric layers. Both the inplane and the transverse electric fields are considered. The deflection field explicitly accounts for the transverse normal strain induced by the electric field. The inplane displace-
ments are approximated as a combination of a global third order variation across the thickness and a layerwise linear variation. The displacement field is expressed in terms of five primary displacement variables and a set of electric potential variables by enforcing exactly the conditions of zero transverse shear stresses at the top and bottom and their continuity at the layer interfaces. The nonlinear coupled equilibrium equations and boundary conditions are derived using a variational principle. This nonlinear theory is used to obtain the initial buckling response of symmetrically laminated hybrid plates under inplane electromechanical loading. The accuracy of the new theory for buckling is tested by direct comparison with the available 3D exact solution [14] for simply supported hybrid plates of highly inhomogeneous, composite and sandwich laminates.
2. Zigzag theory approximations Consider a hybrid plate (Fig. 1) made of L perfectly bonded orthotropic laminas of total thickness h with the midplane chosen as the xy-plane z = 0. Some of the layers are of orthorhombic piezoelectric materials of class mm2 symmetry, with poling along z-direction. The plate is subjected to inplane and transverse mechanical loads and actuation potentials at some piezoelectric surfaces. The reference plane z = 0 either passes through or is the bottom surface of the k0th layer. The kth ply from the bottom has bottom surface at z = zk1. Considering the usual assumption of 2D plate theories, rz ’ 0 [23], the 3D constitutive equations of piezoelasticity relating stresses r, s and electric displacements Dx, Dy, Dz with strains e, c and electric field components Ex, Ey, Ez reduce to eT Ez ; r ¼ Qe 3
^ ^eE; s ¼ Qc
D ¼ ^eT c þ ^gE;
Dz ¼ e3 e þ g33 Ez ;
Fig. 1. Geometry of hybrid plate.
ð1Þ
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
where 2
3
rx 6 7 r ¼ 4 r y 5; " c¼
sxy czx cyz
" s¼
# E¼
# ;
syz "
;
szx
Ex Ey
" D¼
Dx
2
#
Dy
;
contribution of rx, ry via PoissonÕs effect compared to that due to the electric field:, i.e., ez = w,z = S13r1 + S23r2 + d33Ez ’ d33/,z )
3
ex 6e 7 e ¼ 4 y 5;
j
cxy
#
ð2Þ
;
and for general angle-ply laminates 2 12 Q 16 3 " Q11 Q 55 Q 6 7 ^¼ ¼ 4Q 12 Q 22 Q16 5; Q Q 45 Q 16 Q26 Q 66 Q " # " # e15 e14 g12 g11 ^e ¼ ; ^ g¼ ; e3 e25 e24 g12 g22
# 45 Q ; 44 Q ¼ ½ e31
þ z2 nðx; yÞ þ z3 gðx; yÞ;
e32
e36 .
ij , eij , Q gij are the reduced elastic stiffnesses, piezoelectric stress constants and electric permittivities. Partially geometrically nonlinear strain–displacement relations in the spirit of Von Karman are employed. The Lagrange strain–displacement relations including geometric nonlinearity due to deflection w0(x, y) = w(x, y, 0) of the midplane, and the electric field-potential relations are
ez ¼ w;z ; Ex ¼ /;x ;
cyz ¼ uy;z þ w;y ; Ey ¼ /;y ;
where " u¼
ux
#
"
w0;x
#
ð7Þ "
uk x
#
; w0d ¼ ; uk ¼ ; uk y uy w0;y " # " # " # wkx gx nx wk ¼ ; n¼ ; g¼ ; wky gy ny
ð8Þ
uk is the translation and wk is related to the shear rotation of the kth layer. Substituting ux, uy, w from Eqs. (7) and (6), and / from Eq. (5) into Eq. (4) and using Eq. (1) yield the transverse shear stresses s as: ^ kW j ðzÞ/j ; ^ k ½wk þ 2zn þ 3z2 g þ ½^ek Wj ðzÞ Q s¼Q d / /
ð9Þ
T
ð4Þ
czx ¼ ux;z þ w;x ; Ez ¼ /;z ;
where ux, uy; w denote the inplane and transverse displacements and / denotes the electric potential. The subscript comma denotes differentiation. The inplane electric field components Ex, Ey are not considered as zero, since these may be induced by the piezoelectric coupling or may be applied through actuation of segmented piezoelectric actuators. The potential / is approximated as piecewise linear between n/ points zj/ ; j ¼ 1; 2; . . . ; n/ , across the thickness: /ðx; y; zÞ ¼ Wj/ ðzÞ/j ðx; yÞ;
ðzÞ/j ðx; yÞ; wðx; y; zÞ ¼ w0 ðx; yÞ W ð6Þ / R j ðzÞ ¼ z d 33 Wj ðzÞ dz is a piecewise linear funcwhere W /;z / 0 tion. The inplane displacements ux, uy for the kth layer are assumed to follow a global third order variation across the thickness with layerwise linear variation: uðx; y; zÞ ¼ uk ðx; yÞ zw0d ðx; yÞ þ zwk ðx; yÞ
ð3Þ
1 1 ex ¼ ux;x þ w20;x ; ey ¼ uy;y þ w20;y ; 2 2 cxy ¼ ux;y þ uy;x þ w0;x w0;y ;
255
ð5Þ
where /j ðx; yÞ ¼ /ðx; y; zj/ Þ. Wj/ ðzÞ are linear interpolation functions and summation convention is used for indices j, j 0 (used later). This description, which was first introduced by Kapuria [25] and later adopted by Wang [26], allows the piezoelectric layers to be divided into a number of sublayers and a series of elastic layers to be combined into one, for effective discretisation of / across the thickness. The variation of deflection w is obtained by integrating the constitutive equation for ez by neglecting the
where /jd ¼ ½/j;x /j;y . For the mid-plane which passes through the k0th layer, denote u0 ðx; yÞ ¼ uk0 ðx; yÞ ¼ uðx; y; 0Þ, w0 ðx; yÞ ¼ wk0 ðx; yÞ. Thus u0 and w0 are the midplane inplane displacements and shear rotations. Using the 2(L 1) conditions each for the continuity of s and u at the layer interfaces and the four shear traction-free conditions s = 0 at z = z0, zL, the (4L + 4) variables uk, wk, n, g are expressed in terms of only four variables u0 and w0 as in [23]. This yields uðx; y; zÞ ¼ u0 ðx; yÞ zw0d ðx; yÞ þ Rk ðzÞw0 ðx; yÞ þ Rkj ðzÞ/jd ðx; yÞ;
ð10Þ
where Rk(z), Rkj(z) are 2 · 2 matrices of layerwise cubic functions of z whose coefficients are dependent on the material properties and the lay-ups. Thus, even though u, w have layerwise distributions, they are expressed in terms of only five displacement variables u0x , u0y , w0, w0x , w0y by Eqs. (6) and (10).
3. Coupled nonlinear equations of equilibrium The variational principle for the piezoelectric medium [27] R zk be expressed, using the notation h i ¼ PL can k¼1 zþ ð Þ dz for integration across the thickness, as k1
256
Z
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
½hrx dex þ ry dey þ sxy dcxy þ syz dcyz þ szx dczx þ Dx d/;x
A
þ Dy d/;y þ Dz d/;z i p1z dwðx; y; z0 Þ p2z dwðx; y; zL Þ þ Dz ðx; y; z0 Þd/1 Dz ðx; y; zL Þd/n/ qji d/ji dA Z hrn dun þ sns dus þ snz dw þ Dn d/i ds ¼ 0; CL
8 du0 ; dw0 ; dw0 ; d/j ;
ð11Þ
where A denotes the mid-plane surface area of the plate and CL is the boundary curve of the midplane of the plate with normal n and tangent s. p1z ; p2z are the normal forces per unit area on the bottom and top surfaces of the plate in direction z. qji is the the extraneous surface j charge density at the interface z ¼ z/i , where /ji is prescribed. The total number of such prescribed potentials is n/ . The above variational equation is expressed in terms of du0 ; dw0 ; dw0 ; d/j and stress and electric displacement resultants to yield the nonlinear equilibrium equations and the boundary conditions. The stress resulT T tants N ¼ ½N x N y N xy , M ¼ ½M x M y M xy , P ¼ ½P x P yx j ¼ ½Q j P xy P y T , S j ¼ ½S j S j S j S j T , Q ¼ ½Qx Qy T , Q x T
yx
xy
y
x
j T , V ¼ ½V x V y , V j ¼ ½V j V j T and the electric disQ / /x /y y j
½H jx
T H jy ; Gj
placement resultants H ¼ are defined by h iT F 1 ¼ N T M T P T S j T ¼ ½hf3T ri; h iT j Q j F 2 ¼ Qx Qy Q ¼ ½hf4T si x y ð12Þ j si; H j ¼ hWj ðzÞDi; V ¼ hsi; V j ¼ hW j
G ¼
/ hWj/;z ðzÞDz i;
where f3 ¼ I 3
zI 3
/
/
k
U
kj
U
; f4 ¼
h
Rk;z
Rkj ;z
j ðzÞI 2 W /
N xy;x þ N y;y ¼ 0;
M x;xx þ 2M xy;xy þ M y;yy þ ðN x w0;x þ N xy w0;y Þ;x þ ðN xy w0;x þ N y w0;y Þ;x þ F 3 ¼ 0; P x;x þ P yx;y Qx ¼ 0;
þ H jy;y Gj þ F j6 ¼ 0;
j ¼ 1; 2; . . . ; n/ ;
where the mechanical load F 3 ¼ p1z þ p2z and the electrical loads F j6 ¼ Dz ðx; y; zL Þdjn/ Dz ðx; y; z0 Þdj1 þ qji djji . The boundary conditions on CL are the prescribed values of one of the factors of each of the following products: u0n N n ; u0s N ns ; w0 ðV n þ M ns;s þ N n w0;n þ N ns w0;s Þ; w0;n M n ; w0n P n ; w0s P ns ; /j;n S jn ; /j ½H jn V j/n S jns;s and at corners si:
w0 ðsi ÞDM ns ðsi Þ;
i ;
Rkj 22
0
/j ðsi ÞDS jns ðsi Þ
with V n ¼ ðM x;x þ M xy;y Þnx þ ðM xy;x þ M y;y Þny ; j nx Q j ny . V j ¼ ðS j þ S j Þnx þ ðS j þ S j Þny Q x;x
yx;y
y;y
xy;x
bj ¼ hf3T ðzÞeT3 Wj/;z ðzÞi;
ð13Þ
Rkj 12
0
j ¼ hf T ðzÞ^eWj0 ðzÞi; b / 4
x
y
The resultants F1, F2, Hj, Gj defined in Eq. (12) are expressed in terms of the displacement and potential variables by substituting the expressions of r, s, D, Dz from Eq. (1) into Eq. (12) and using Eqs. (4)–(6) and (10): 0 0 0 1 j / j0 ; e2 þ b F 1 ¼ Ae1 þ bj /j þ A Uw w0d ; F 2 ¼ A d 2 T 0 T 0 1 0 j e2 E jj /jd0 ; Gj ¼ bj e1 Ejj /j þ bjw Uw w0d ; H j ¼ b 2 ð16Þ
h iT e1 ¼ u0x ;x u0y ;y u0x ;y þ u0y ;x w0;xx w0;yy 2w0;xy w0x ;x w0x ;y w0y ;x w0y ;y /j;xx /j;xy /j;yx /j;yy ; 2 3 0 w0;x h iT j j e2 ¼ w0x w0y /;x /;y ; Uw ¼ 4 0 w0;y 5; w0;y w0;x T ^ f 4 ðzÞi; bjw ¼ he3 Wj i; ½A; A ¼ hf3 ðzÞQ½f3 ðzÞ; I 3 i; A ¼ hf4T ðzÞQ /;z 0
ð14Þ
P xy;x þ P y;y Qy ¼ 0;
j þ Q j Sj Sj Sj Sj þ H j Q x;x x;xx xy;xy yx;xy y;yy x;x y;y
/n
Rk21 Rk11 Rk22 Rk12 2 kj 3 0 R11 0 Rkj 12 6 7 kj kj 7 Ukj ¼ 6 4 0 R21 0 R22 5. Rkj 11
N x;x þ N xy;y ¼ 0;
ð15Þ
In being a n · n identity matrix and 2 k 3 R11 0 Rk12 0 6 7 k k 7 Uk ¼ 6 4 0 R21 0 R22 5;
Rkj 21
It can be shown that elements of Rk, Rkj, N, M, P, Sj transform as second order tensors and elements of j ; H j transform as vectors for the coplanar V ; V j/ ; Q; Q axes x, y and n, s. Using expressions of /, w, u from Eqs. (5), (6), (10) and using definitions in Eq. (12), the area integral in Eq. (11) is expressed in terms of du0x , du0y , dw0, dw0x , dw0y , d/j, by using GreenÕs theorem if required, and the terms involving du0x ; du0y ; dw0x ; dw0y ; dw0;x ; dw0;y ; d/j;x ; d/j;y in the integrand of CL are expressed in terms of components n, s. Thus Eq. (11) yields coupled nonlinear field equations consisting of five equations of force equilibrium and n/ equations for charge equilibrium:
0
0
Ejj ¼ hg33 Wj/;z ðzÞWj/;z ðzÞi;
0
jj ¼ h^gWj ðzÞWj0 ðzÞi; E / /
ð17Þ
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
where 2
A11 A12 6 6A 6 21 A22 6 6 . .. 6 .. . 6 6 6 A10;1 A10;2 A¼6 6 j 6 A11;1 Aj11;2 6 6 j 6 A12;1 Aj12;2 6 6 j 6 A13;1 Aj13;2 4 Aj14;1 Aj14;2 2 j0 3 b1 6 j0 7 6 b2 7 6 7 6 . 7 6 . 7 6 . 7 6 j0 7 6b 7 0 j 10 7 b ¼6 6 jj0 7; 6 b11 7 6 7 6 jj0 7 6 b12 7 6 07 6 jj 7 4 b13 5
. . . A1;10
0 Aj1;11
0 Aj1;12
0 Aj1;13
0
0
0
0 Aj1;14
3
7 7 7 7 7 7 7 7 j0 A10;14 7 7 ¼ AT ; 7 jj0 A11;14 7 7 0 7 Ajj12;14 7 7 0 7 Ajj13;14 7 5 0 Ajj14;14 0
. . . A2;10 Aj2;11 Aj2;12 Aj2;13 Aj2;14 .. .. .. .. .. .. . . . . . . 0
0
0
. . . A10;10 Aj10;11 Aj10;12 Aj10;13 0
0
0
0
0
0
0
0
0
0
0
0
. . . Aj11;10 Ajj11;11 Ajj11;12 Ajj11;13 . . . Aj12;10 Ajj12;11 Ajj12;12 Ajj12;13 . . . Aj13;10 Ajj13;11 Ajj13;12 Ajj13;13 . . . Aj14;10 Ajj14;11 Ajj14;12 Ajj14;13
Taking into account the corresponding zero elements of the matrices defined in Eq. (18), the nonlinear terms are obtained as Ln U Þ ¼ 1 ½A11 w2 þ A12 w2 þ ½A33 w0;x w0;y ; ðLn U 0;x 0;y ;x 1 ;y 2 Þ ¼ 1 ½A21 w2 þ A22 w2 þ ½A33 w0;x w0;y ; ðLn U 0;x 0;y ;y 2 ;x 2 n ðL U Þ3 ¼ N x w0;xx þ N y w0;yy þ 2N xy w0;xy þ w0;x ðN x;x þ N xy;y Þ 1 þ w0;y ðN xy;x þ N y;y Þ þ ½A41 w20;x þ A42 w20;y ;xx 2 1 þ ½A51 w20;x þ A52 w20;y ;yy þ ½2A63 w0;x w0; y;xy ; 2 1 Þ ¼ ½A71 w2 þ A72 w2 þ ½A83 w0;x w0;y ; ðLn U 0;x 0;y ;x 4 ;y 2 1 Þ ¼ ½A93 w0;x w0;y þ ½A10;1 w2 þ A10;2 w2 ; ðLn U 0;x 0;y ;y 5 ;x 2 1 1 1 Þ ¼ bjw ðLn U w2 þ bjw w2 þ ½Aj w2 þ Aj11;2 w20;y ;xx 5þj 2 1 0;x 2 2 0;y 2 11;1 0;x þ ½ðAj12;3 þ Aj13;3 Þw0;x w0;y ;xy 1 þ ½Aj14;1 w20;x þ Aj14;2 w20;y ;yy ; j ¼ 1; . . . ; n/ . 2
0
bjj14 2 A11 6 A21 6 6 . 6 . 6 . 6 6 A10;1 6 A ¼6 j 6 A11;1 6 6 j 6 A12;1 6 6 Aj 4 13;1
ð21Þ
3
A12 A13 A22 A23 7 7 2 .. 7 .. 7 11 A . 7 . 7 6 6 A10;2 A10;3 7 21 7 6A 7; A ¼ 6 j j j 6A A11;2 A11;3 7 7 4 31 7 Aj12;2 Aj12;3 7 j A 7 41 j j 7 A13;2 A13;3 5
Aj14;1 Aj14;2 Aj14;3 0 0 3 j j b b 11 12 6 j0 " jj0 0 7 6b b j 7 0 E j0 jj 6 7 21 22 ¼6 0 ¼ 110 b ; E 7 0 jj 6b jj 7 jj E 4 31 b 21 32 5 0 0 jj b jj b
0
0
12 A j A j A 13 14
The nonzero elements of load vector P are P 3 ¼ F 3 ; P j6 ¼ F j6 . The elements of L are listed in Ref. [23].
3
0 0 7 22 A j A j 7 A 23 24 7 T; 7¼A j jj0 jj0 7 A32 A33 A34 5 0 0 j A jj A jj A
42
43
44
2
41
257
0
jj E 12 0
jj E 22
#
; bjw ¼ bjw . bjw bjw 1 2 3
42
ð18Þ
Substituting Eq. (16) for the resultants into Eq. (14) yields the nonlinear coupled equations of equilibrium in terms of the primary displacement and potential vari: ables, U ¼ P ; þ Ln U ð19Þ LU where h iT ¼ u0x u0y w0 w0 w0 /1 /2 . . . /n/ ; U x y n T P ¼ P 1 P 2 P 3 P 4 P 5 P 16 P 26 . . . P 6/ .
4. Buckling under uniform electromechanical load 4.1. Pre-buckling solution Consider a symmetrically laminated plate which is in equilibrium under applied uniform inplane normal strains e0x ; e0y , zero shear strain c0xy ¼ 0 and actuation potentials independent of x, y coordinates. Let this pre-buckling equilibrium state be denoted by superscript ( )0. For the symmetrically laminated plate under symmetrical loading about xy-plane,
0
0
0
0
N 0x ¼ A11 e0x þ A12 e0y þ bj1 /0j ; N 0y ¼ A12 e0x þ A22 e0y þ bj2 /0j ; 0j0
¼
Defining U ¼ /1
/2
...
/n/
h bj ¼ b1j
b2j
...
bj /
þ
bj2 e0y
F 0j 6 .
bj1 e0x
E /
L is a symmetric matrix of linear differential opera are the nonlinear terms due to geotors in x and y. Ln U 45 ¼ 0; metric nonlinearity. For cross-ply plates, Q Q16 ¼ Q26 ¼ 0; e14 ¼ e25 ¼ 0, b6 ¼ 0; g12 ¼ 0; e36 ¼ 0.
ð22Þ
Considering this, the plate constitutive Eq. (16) and equilibrium Eq. (14) yield
jj0
ð20Þ
w0x ¼ w0y ¼ 0.
w00 ¼ 0;
n
T iT
;
ð23Þ
258
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
and
F 6 ¼ F 16
F 26
...
n F 6/
T
buckling is obtained for simply supported rectangular plates of sides a, b along the axes x, y for the boundary conditions
:
Eq. (23)3 can be written in matrix form as EU0 ¼ b1 e0x þ b2 e0y F 06 .
ð24Þ
In order to consider open as well as closed circuit electric boundary conditions, U0 is partitioned into a set of unknown output voltages U0s at zj/ Õs where / is not prescribed and a set of known input actuation voltages U0a at the actuated surfaces. Accordingly, Eq. (24) can be partitioned and arranged as " # " 0 # b2s 0 b1s 0 Ess Esa F 06s Us e þ e ¼ . Eas Eaa U0a b1a x b2a y F 06a ð25Þ Eq. (25) is solved for U0s and the same is substituted into Eq. (23) to obtain the pre-buckling forces N 0x ; N 0y in terms of the known loading parameters T 1 0 0 N 0x ¼ ðA11 þ bT1s E1 ss b1s Þex þ ðA12 þ b1s Ess b2s Þey T 1 0 0 þ ðbT1a bT1s E1 ss Esa ÞUa b1s E ss F 6s ; T 1 0 0 N 0y ¼ ðA12 þ bT2s E1 ss b1s Þex þ ðA22 þ b2s Ess b2s Þey
ð26Þ
T 1 0 0 þ ðbT2a bT2s E1 ss Esa ÞUa b2s E ss F 6s .
4.2. Stability equations Let the solution for just after buckling be denoted by ð^Þ on the entities. The governing equations for buckling are derived using the approach of Noor and Burton [28] by introducing a small parameter for the buckling mode. is described by an The size of the buckling mode U ^ ¼U 0 þ U with U given arbitrary small parameter : U by Eq. (20)1. Substituting this solution into Eqs. (19), using Eq. (22) and considering upto the first order terms : in yield the following stability equations for U h iT þ 0 0 ðN 0x w0;xx þ N 0y w0;yy Þ 0 0 0 LU T ¼ P ¼ 0 0 0 0 0 F j6 ; ð27Þ ^ Þ ; ðLn U ^ Þ ; ðLn U ^ Þ ; ðLn U ^ Þ ; ðLn U ^ Þ are of secsince ðLn U 1 2 4 5 5þj 2 n^ 0 ond order and ðL U Þ3 ¼ ðN x w0;xx þ N 0y w0;yy Þ is of first j order. The load F 6 is the contribution of the incremental only. For a set of zero incremental potential solution U at the actuator locations and zero incremental electric displacement at the unknown potential locations, the j load F^ 6 is zero for index j corresponding to such surfaces. Note that whereas the equilibrium equations (14) are nonlinear, the stability equations (27) are linear in terms of pre-buckling loads N 0x ; N 0y . To assess the accuracy of the theory developed herein, by comparison with the exact 3D piezoelasticity solution [14], analytical Navier solution of Eq. (27) for
at x ¼ 0; a:
N x ; u0y ; w0 ; w0y ; M x ; P x ; /j ; S jx ¼ 0;
at y ¼ 0; b:
N y ; u0x ; w0 ; w0x ; M y ; P y ; /j ; S jy ¼ 0.
ð28Þ
The solution for the (m, n)th spatial mode of buckling is taken as: 2 3 2 3 sinðnyÞ ½ w0 /j mn sinðmxÞ w0 /j 6u 7 6 7 sinðnyÞ 5 4 0x w0x 5 ¼ 4 ½ u0x w0x mn cosðmxÞ u0y w0y cosðnyÞ ½ u0y w0y mn sinðmxÞ ð29Þ ¼ mp=a; n ¼ np=b. Substituting these into Eq. with m (27) yields ðK K G ÞU
mn
mn
¼ P ;
ð30Þ
KG is the geometric stiffness matrix with only one non 2 N 0x n2 N 0y . K is symmetric zero element K G ð3; 3Þ ¼ m stiffness matrix whose elements are not listed here for brevity. Partitioning electric potentials U into the unknown part Us and known zero incremental actuation potentials, Ua as before, Eq. (30) can be written for the T ~ mn ¼ u0x mn u0y mn w0mn w0x mn w0y nm Umn unknowns U s as: ~ K ~ G ÞU ~ mn ¼ ðK ~ kK ÞU ~ mn ¼ 0; ðK G
ð31Þ
~G ¼ where K Eq. (31) is a generalised linear eigenvalue problem and the eigenvalue k is the buckling load factor. The eigenvalues and eigenvectors are obtained by QR algorithm after reducing to Heissenberg form using NAG subroutines. kK G .
5. Numerical evaluation of the theory The accuracy of the present theory for buckling response is established by direct comparison with the exact 3D piezoelasticity [14] solution for simply supported cross-ply symmetrically laminated hybrid plates. Since TOT has the same global third order variation for the displacement field without the layerwise terms and has the same number of displacement variables as the present theory (ZIGT), present results are also compared with the coupled TOT [20,29] extended for the buckling case. This comparison will establish the effect of the layerwise terms for the displacements incorporated in the ZIGT. Results are presented for hybrid plates of three different laminate configurations (a), (b) and (c). The stacking order is mentioned from the bottom. The elastic substrate of plate (a) has five plys of equal thickness 0.16h of materials 2/3/1/3/2 with orientations of the principal material axis 1 as [90/0/0/0/90]. It is a good test case for assessing a 2D theory since the plies have highly
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
inhomogeneous stiffness in tension and shear. The substrate of plate (b) is a graphite–epoxy composite laminate of material 4 with four layers of equal thickness 0.2h with lay-up [0/90/90/0]. The substrate of plate (c) is a five-layer sandwich having graphite–epoxy faces [0/90] and a soft core with thicknesses .04h/.04h/.64h/ .04h/.04h. All the plates have two PZT-5A layers, each of thickness 0.1h, bonded to their elastic substrate on its top and bottom surfaces. The PZT-5A layers have poling in +z-direction. The top and the bottom of the substrate are grounded. The material properties are selected as [14]: [(Y1, Y2, Y3, G12, G23, G31), m12, m13, m23] = Material 1: [(6.9, 6.9, 6.9, 1.38, 1.38, 1.38) GPa, 0.25, 0.25, 0.25] Material 2: [(224.25, 6.9, 6.9, 56.58, 1.38, 56.58) GPa, 0.25, 0.25, 0.25] Material 3: [(172.5, 6.9, 6.9, 3.45, 1.38, 3.45) GPa, 0.25, 0.25, 0.25] Material 4: [(181, 10.3, 10.3, 7.17, 2.87, 7.17) GPa, 0.28, 0.28, 0.33] Face: [(131.1, 6.9, 6.9, 3.588, 2.3322, 3.588) GPa, 0.32, 0.32, 0.49] Core: [(0.2208, 0.2001, 2760, 16.56, 455.4, 545.1) MPa, 0.99, 3 · 105, 3 · 105] PZT-5A: [(61.0, 61.0, 53.2, 22.6, 21.1, 21.1) GPa, 0.35, 0.38, 0.38], and [(d31, d32, d33, d15, d24), (g11, g22, g33)] = [(171, 171, 374, 584, 584) · 1012 m/V, (1.53, 1.53, 1.5) · 108 F/m] where Yi, Gij, mij, dij and gij denote YoungÕs moduli, shear moduli, PoissonÕs ratios, piezoelectric strain constants and electric permittivities, respectively. Four pre-buckling load cases are considered: 1. Uniaxial inplane strain e0x ¼ e0 with e0y ¼ 0. The top and bottom surfaces are at closed circuit condition with zero potential ½/1 ¼ /n/ ¼ 0. 2. Uniaxial inplane strain e0y ¼ e0 with e0x ¼ 0 and the top and the bottom surfaces at zero potential.
259
3. Biaxial inplane strain e0x ¼ e0y ¼ e0 with top and bottom surfaces at zero potential. 4. Applied actuation potentials at top and bottom surfaces /n/ ¼ /1 ¼ /0 with the immovable edges (e0x ¼ e0y ¼ 0). For load cases 1 to 3, the lowest value of the strain e0 for buckling is denoted as ecr. For load case 4, the critical value of potential /0 for buckling is denoted as /cr. The results for these cases are nondimensionalised with S = a/h, d0 = 374 · 1012 CN1, Y0 = 6.9 GPa for laminates (a) and (c), and 10.3 GPa for laminate (b): ¼ / d 0 S 3 =a; ecr ¼ ecr S 2 ; / cr cr Þ ¼ ðSu; Sv; wÞ= maxðwÞ; ðu; v; w y Þ ¼ ðrx ; ry ÞS 2 h=Y 0 maxðwÞ; ð rx ; r ðszx ; syz Þ ¼ ðszx ; syz ÞS 3 h=Y 0 maxðwÞ; where max(w) denotes the largest value of w through the thickness. The critical buckling strains ecr obtained from the exact 3D piezoelasticity solution [14] and percentage (%) errors of the present theory and coupled TOT for square hybrid test plate (a) are presented in Table 1 for load cases 1, 2 and 3. Four values of the span-to-thickness ratio S = 5, 7.5, 10, 20 are considered for comparison. The values of (m, n) corresponding to the critical buckling mode are also given in the table. Similar results for composite plate (b) and sandwich plate (c) are presented in Tables 2 and 3, respectively. It is observed that the present theory yields highly accurate results for ecr for all laminate configurations and load cases with a maximum error of 1.5% even for thick plates with S = 5. In contrast, TOT results have errors upto 133.9%, 11.0% and 55.9% for plates (a), (b) and (c), respectively. Even for thin plates with S = 20, the error in TOT is as high as 15.4% for test plate (a) and 12.1% for sandwich plate (c). Even the critical mode (m, n) is wrongly estimated by TOT in some cases. The 3D exact results and % errors of (load ZIGT and TOT for the critical electric potential / cr case 4) are presented in Table 4 for the three square
Table 1 Exact 3D results and % errors of ZIGT and TOT for buckling inplane strain for plate (a) S
Entity
Load case 1 Exact
5 7.5 10 20
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
0.65980 (2, 1) 1.19904 (2, 1) 1.73289 (2, 1) 2.88725 (1, 1)
Load case 2 % Error
Exact
ZIGT
TOT
0.3 (2, 1) 1.0 (2, 1) 1.5 (2, 1) 0.0 (1, 1)
59.1 (3, 1) 56.2 (2, 1) 43.3 (2, 1) 15.2 (1, 1)
0.61949 (1, 2) 1.09362 (1, 2) 1.54017 (1, 1) 2.38148 (1, 1)
Load case 3 % Error
Exact
ZIGT
TOT
0.5 (1, 2) 0.2 (1, 2) 0.3 (1, 1) 0.2 (1, 1)
133.9 (1, 1) 81.4 (1, 1) 49.4 (1, 1) 15.4 (1, 1)
0.3744 (1, 1) 0.6280 (1, 1) 0.8449 (1, 1) 1.3048 (1, 1)
% Error ZIGT
TOT
0.6 (1, 1) 0.0 (1, 1) 0.1 (1, 1) 0.1 (1, 1)
111.9 (1, 1) 72.9 (1, 1) 49.1 (1, 1) 15.3 (1, 1)
260
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
Table 2 Exact 3D results and % errors of ZIGT and TOT for buckling inplane strain for plate (b) S
Entity
Load case 1
Load case 2
Exact
5 7.5 10 20
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
% Error
0.92953 (2, 1) 1.35619 (1, 1) 1.59769 (1, 1) 1.95146 (1, 1)
Exact
ZIGT
TOT
0.8 (2, 1) 0.1 (1, 1) 0.1 (1, 1) 0.1 (1, 1)
11.0 (2, 1) 6.3 (1, 1) 4.5 (1, 1) 1.4 (1, 1)
0.97300 (1, 1) 1.36354 (1, 1) 1.60444 (1, 1) 1.95454 (1, 1)
Load case 3 % Error
Exact
ZIGT
TOT
1.3 (1, 1) 0.6 (1, 1) 0.4 (1, 1) 0.1 (1, 1)
7.8 (1, 1) 5.7 (1, 1) 4.0 (1, 1) 1.3 (1, 1)
0.48543 (1, 1) 0.68011 (1, 1) 0.80061 (1, 1) 0.97651 (1, 1)
% Error ZIGT
TOT
1.0 (1, 1) 0.4 (1, 1) 0.2 (1, 1) 0.0 (1, 1)
8.0 (1, 1) 5.9 (1, 1) 4.2 (1, 1) 1.4 (1, 1)
Table 3 Exact 3D results and % errors of ZIGT and TOT for buckling inplane strain for plate (c) S
5 7.5 10 20
Entity
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
Load case 1
Load case 2
Exact
% Error ZIGT
TOT
0.71404 (2, 1) 1.39331 (2, 1) 2.08738 (1, 1) 3.59027 (1, 1)
0.1 (2, 1) 0.4 (2, 1) 0.2 (1, 1) 0.1 (1, 1)
50.9 (2, 1) 48.9 (1, 1) 31.8 (1, 1) 12.1 (1, 1)
Load case 3
Exact
% Error
Exact
ZIGT
TOT
0.67831 (1, 2) 1.28178 (1, 2) 1.98108 (1, 2) 3.59088 (1, 1)
0.2 (1, 2) 0.5 (1, 2) 0.4 (1, 2) 0.1 (1, 1)
53.3 (1, 2) 55.9 (1, 2) 38.9 (1, 1) 12.1 (1, 1)
0.40547 (1, 1) 0.73426 (1, 1) 1.04388 (1, 1) 1.79529 (1, 1)
% Error ZIGT
TOT
0.7 (1, 1) 0.5 (1, 1) 0.2 (1, 1) 0.1 (1, 1)
50.9 (1, 1) 41.3 (1, 1) 31.8 (1, 1) 12.1 (1, 1)
Table 4 Exact 3D results and % errors of ZIGT and TOT for critical electric potential for buckling S
Entity
Plate (a) Exact
5 7.5 10 20
/ cr (m, n) / cr (m, n) / cr (m, n) / cr (m, n)
Plate (b) % Error
0.36515 (1, 1) 0.62092 (1, 1) 0.84116 (1, 1) 1.31539 (1, 1)
Exact
ZIGT
TOT
3.5 (1, 1) 2.7 (1, 1) 2.1 (1, 1) 0.9 (1, 1)
116.2 (2, 1) 77.7 (1, 1) 52.2 (1, 1) 16.2 (1, 1)
0.52176 (1, 1) 0.74734 (1, 1) 0.89057 (1, 1) 1.10546 (1, 1)
Plate (c) % Error
Exact
ZIGT
TOT
5.1 (1, 1) 3.5 (1, 1) 2.4 (1, 1) 0.8 (1, 1)
14.7 (1, 1) 10.0 (1, 1) 6.9 (1, 1) 2.2 (1, 1)
0.14028 (1, 1) 0.25532 (1, 1) 0.36407 (1, 1) 0.62974 (1, 1)
% Error ZIGT
TOT
1.1 (1, 1) 0.8 (1, 1) 0.8 (1, 1) 0.5 (1, 1)
53.6 (1, 1) 43.1 (1, 1) 33.1 (1, 1) 12.6 (1, 1)
Table 5 Effect of b/a on the errors of ZIGT and TOT results for buckling parameters for plate (a) b/a
Entity
Load case 1 Exact
Entity % Error ZIGT
1 2 4 10
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
1.73289 (2, 1) 0.79996 (1, 1) 0.56101 (1, 1) 0.49803 (1, 1)
1.5 (2, 1) 1.2 (1, 1) 1.8 (1, 1) 2.1 (1, 1)
Load case 4 Exact
TOT 43.3 (2, 1) 17.3 (1, 1) 9.3 (1, 1) 7.1 (1, 1)
/ cr (m, n) / cr (m, n) / cr (m, n) / cr (m, n)
0.84116 (1, 1) 0.54542 (1, 1) 0.44376 (1, 1) 0.41278 (1, 1)
% Error ZIGT
TOT
2.1 (1, 1) 2.7 (1, 1) 3.1 (1, 1) 3.3 (1, 1)
52.2 (1, 1) 19.0 (1, 1) 10.6 (1, 1) 8.3 (1, 1)
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
261
y ; syz for critical buckling mode of square hybrid plate (a) under load case 1. Fig. 2. Distributions of v; r
x ; szx for critical buckling mode of square hybrid plate (b) under load case 1. Fig. 3. Distributions of u; r
plates. It is revealed that the maximum error in ZIGT for thick plates with S = 7.5 is 3.5% whereas the error in TOT is 77.7%, 10.0%, 43.1% for plates (a), (b), (c),
respectively. Even for thin plates with S = 20, the error is as high as 16.2% and 12.6% for plates in TOT for / cr (a) and (c). The errors in the 2D theories are dependent
262
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
x ; szx for critical buckling mode of square hybrid plate (c) under load case 1. Fig. 4. Distributions of u; r
x ; szx for critical buckling mode of square hybrid plate (a) under potential load (case 4). Fig. 5. Distributions of u; r
on the uniaxial or bi-axial nature of inplane loads as well as the laminate configurations. The extent of error in
TOT results for plate (a) confirms that this plate indeed offers a good test case for assessing the accuracy of 2D
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
263
y ; syz for critical buckling mode of plate (a) under load case 1. Fig. 6. Effect of b/a on distributions of v; r
theories for buckling. The present theory successfully passes this test. As is well known, equivalent single layer theories like TOT yield highly inaccurate results for sandwich plates, whereas ZIGT yields excellent results for these plates too. for modThe errors in ZIGT and TOT for ecr and / cr erately thick plate (a) with S = 10 are presented in Table 5 for four values of b/a ratio ranging from 1 to 10. It is observed that whereas the error in ZIGT generally shows a marginal increase, the error in TOT reduces with the increase in b/a. The through-the-thickness distributions of predominant modal inplane displacements u=v, normal stresses x = r ry and transverse shear stresses szx =syz for the critical buckling mode at (x, y) locations where they are maximum are compared in Figs. 2–4 for square hybrid plates (a), (b) and (c) for load case 1 for S = 5 and 10. The location parameters a; b are defined as a ¼ a=2m; b ¼ b=2n where (m, n) correspond to the critical buckling mode. It is observed that the present distributions are in excellent agreement with the exact modal distributions for all plates even for the thick ones with S = 5. In contrast, the distributions obtained from TOT have large error even for S = 10. The present theory is able to accurately capture the layerwise zigzag distribution of the modal inplane displacements obtained from the 3D solution, but the global distribution of the inplane displacements in TOT is unable to capture it. This explains why ZIGT yields accurate results and TOT does not.
The presented results cover a wide variety of possible cases as is evident from the qualitative distinct nature of the through-the-thickness distributions of the modal entities for the three plates. The distributions of the modal displacements and stresses across the thickness for square plate (a) under electric potential load of case 4 are compared in Fig. 5. For this case too, the present distributions are accurate. The effect of b/a ratio on the accuracy of the through-the-thickness distributions is illustrated for plate (a) in Fig. 6 for S = 10. It is seen that whereas ZIGT yields accurate distributions for all b/a ratios, the error in the TOT reduces with the increase in b/a.
6. Conclusions An efficient coupled zigzag theory incorporating geometric nonlinearity has been developed for hybrid piezoelectric plates. The theory considers layerwise distribution of the inplane displacements across the thickness and accounts for the transverse normal strain caused by the potential field, but the displacement field is expressed in terms of only five displacement variables as in TOT. The nonlinear theory is used to obtain the buckling response of symmetrically laminated plates due to inplane electromechanical loading. The accuracy of the theory is established directly by comparison with the exact 3D piezoelasticity solution for a simply
264
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
supported highly inhomogeneous hybrid test plate, a composite plate and a sandwich plate under uniaxial and bi-axial inplane strains and electric potential load, covering a wide range of possible cases. The comparisons establish that the present theory can accurately predict the zigzag distributions of inplane displacements as obtained from 3D exact solutions and hence yields accurate results for the buckling loads and modal distributions of the response entities. In contrast, TOT results may be highly inaccurate even for thin plates with S = 20. The present theory is also applicable for hybrid sandwich plates, for which TOT yields poor results. Acknowledgements The first author is thankful to Department of Science and Technology, Government of India, for providing financial assistance for this work under SERC scheme. References [1] Leissa AW. A review of laminated composite plate buckling. Appl Mech Rev 1987;40:575–91. [2] Chia CY. Geometrically nonlinear behavior of composite plates: A review. Appl Mech Rev 1988;41:439–51. [3] Kam TY, Chang PR. Buckling of shear deformable laminated composite plates. Compos Struct 1992;22(4):223–34. [4] Balamurugan M, Ganapathi M, Varadan TK. Nonlinear dynamic instability of laminated composite plates using finite element method. Compos Struct 1996;60(1):125–30. [5] Fares ME, Zankour AM. Buckling and free vibration of nonhomogeneous composite cross-ply laminated plates with various plate theories. Compos Struct 1999;44(4):279–87. [6] Chattopadhyay A, Radu AG. Dynamic instability of composite laminates using a higher order theory. Compos Struct 2000;77: 453–60. [7] Chakrabarti A, Sheikh AH. Buckling of laminated composite plates by a new element based on higher order shear deformation theory. Adv Mater Struct 2003;10(4):303–17. [8] Librescu L, Hause T, Jhonson T. Buckling and nonlinear response sandwich curved panels to combined mechanical loads-implications of face-sheet elastic tailoring. J Sandwich Struct Mater 2000;2:246–69. [9] Yuan WX, Dawe DJ. Free vibration and stability analysis of stiffened sandwich plates. Compos Struct 2004;63(1):123–37. [10] Srinivas S, Rao AK. Bending, vibration and buckling of simply supported thick orthotropic rectangular plates and laminates. Int J Solids Struct 1970;11:1463–81.
[11] Gu H, Chattopadhyay A. Three-dimensional elasticity solution for buckling of composite laminates. Compos Struct 1996;34: 292–9. [12] Thomson SP, Loughan J. The control of the post-buckling response in thin composite plates using smart technology. ThinWalled Struct 2000;36:231–63. [13] Batra RC, Geng TS. Enhancement of dynamic buckling load for a plate by using piezoceramic actuators. Smart Mater Struct 2001; 10:925–33. [14] Kapuria S, Achary GGS. Exact 3D piezoelasticity solution for buckling of hybrid cross-ply plates using transfer matrices. Acta Mech 2004;170(1–2):25–46. [15] Yang JS. Buckling of a piezoelectric plate. Int J Appl Electromagnetics Mech 1998;9:399–408. [16] Hu YT, Yang JS, Jiang Q. A model for electroelastic plates under biasing fields with applications in buckling analysis. Int J Solids Struct 2002;39:2629–42. [17] Chandrashekhara K, Bhatia K. Active buckling control of smart composite plates-finite-element analysis. Smart Mater Struct 1993; 2:31–9. [18] Varelis D, Saravanos DA. Nonlinear coupled mechanics and initial buckling of composite plates with piezoelectric actuators and sensors. Smart Mater Struct 2002;11:330–6. [19] Varelis D, Saravanos DA. Coupled buckling and postbuckling analysis of active laminated piezoelectric composite plates. Int J Solids Struct 2004;41(5–6):1519–38. [20] Correia VMF, Soares CMM, Soares CAM. Buckling optimization of composite laminated adaptive structures. Compos Struct 2003;62(3–4):315–21. [21] Yang JS. Equations for thick elastic plates with partially electroded piezoelectric actuators and higher order electric field. Smart Mater Struct 1999;8:73–82. [22] Bisegna P, Caruso G. Evaluation of higher-order theories of piezoelectric plates in bending and in stretching. Int J Solids Struct 2001;38:8805–30. [23] Kapuria S. A coupled zigzag third order theory for hybrid crossply plates. ASME J Appl Mech 2004;71(5):604–14. [24] Kapuria S, Alam N. Zigzag theory for buckling of hybrid piezoelectric beams under mechanical loads. Int J Mech Sci 2004; 46:1–25. [25] Kapuria S. An efficient coupled theory for multi-layered beams with embedded piezoelectric sensory and active layers. Int J Solids Struct 2001;38:9179–99. [26] Wang SY. A finite element model for the static and dynamic analysis of a piezoelectric bimorph. Int J Solids Struct 2004;41: 4075–96. [27] Tiersten HF. Linear piezoelectric plate vibrations. New York: Plenum Press; 1969. 25–37, 54, 55. [28] Noor AK, Burton WS. Three-dimensional solutions for thermal buckling of multilayered anisotropic plates. J Eng Mech 1992; 118(4):683–701. [29] Kapuria S, Achary GGS. A coupled consistent third order theory for hybrid piezoelectric plates. Compos Struct 2005;70(1):120–33.
Nonlinear coupled zigzag theory for buckling of hybrid piezoelectric plates S. Kapuria *, G.G.S. Achary Applied Mechanics Department, Indian Institute of Technology Delhi, Hauz Khas, New Delhi 110016, India Available online 6 June 2005
Abstract A coupled zigzag theory for hybrid piezoelectric plates, recently developed by the first author, is extended to include geometric nonlinearity in the Von Karman sense. In this theory, the potential field is approximated as piecewise linear across sublayers. The deflection approximation accounts for the transverse normal strain due to an electric field. The inplane displacements are assumed to follow a global third order variation across the thickness with a layerwise linear variation. The shear continuity conditions at the layer interfaces and the shear traction-free conditions at the top and bottom are enforced to formulate the theory in terms of only five primary displacement variables, independent of the number of layers. The coupled nonlinear equations of equilibrium and the boundary conditions are derived from a variational principle. The nonlinear theory is used to obtain the initial buckling response of symmetrically laminated hybrid plates under inplane electromechanical loading. Analytical solutions for buckling of simply supported plates under uniaxial and bi-axial inplane strains and electric potential are obtained for comparing the results with the available exact three-dimensional piezoelasticity solution. The comparison establishes that the present theory is very accurate for buckling response of hybrid plates even for highly inhomogeneous lay-ups. 2005 Elsevier Ltd. All rights reserved. Keywords: Buckling; Zigzag theory; Electromechanical coupling; Geometric nonlinearity; Hybrid plate
1. Introduction Laminated composite and sandwich structures offer superior replacement of their metallic counterparts in aerospace, aeronautical, automotive and other applications due to their superior strength to weight and stiffness to weight ratios, high fatigue life as well as high temperature resistance. Distributed piezoelectric sensors and actuators are surface-mounted or embedded in these laminates to introduce adaptive capabilities for active shape control, vibration suppression, acoustic control and buckling control, among others. A considerable amount of research has been dedicated to the study of buckling and post-buckling behaviour of composite *
Corresponding author. Tel.: +91 112 659 1218; fax: +91 112 658 1119. E-mail address: [email protected] (S. Kapuria). 0263-8223/$ - see front matter 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruct.2005.04.010
and sandwich plates [1,2]. Some of the recent works include the static and dynamic buckling analyses of composite plates using first order shear deformation theory (FSDT) [3–5] and refined third order theory (TOT) [5– 7], and of sandwich plates using layerwise models [8,9]. Three-dimensional (3D) elasticity solutions [10,11] of buckling of elastic laminated plates have been presented for simply supported boundary conditions. These 3D solutions serve as useful benchmarks for evaluation of various 2D plate theories. Study of buckling of hybrid laminated structures with surface-bonded or embedded piezoelectric layers and enhancement of their buckling capacity by proper actuation began in the last decade [12]. Batra and Geng [13] presented a numerical analysis for dynamic buckling of a plate with surface bonded piezocermic elements, based on 3D elasticity without considering two-way electromechanical coupling (direct piezoelectric effect). An exact 3D coupled piezoelasticity
254
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
solution for buckling of simply supported symmetrically laminated hybrid plates was presented by Kapuria and Achary [14], wherein the effect of piezoelectric coupling was demonstrated. 2D coupled classical plate theory [15] and first order theory [16] have been applied to study the buckling response of piezoelectric plates. Chandrashekhra and Bhatia [17] presented a finite element model based on uncoupled FSDT for active buckling control of hybrid composite plates. Varelis and Saravanos [18,19] developed an eight noded coupled nonlinear laminated piezoelectric plate element based on FSDT approximations for the displacement field and layerwise approximation for the electric potential and studied initial buckling and post-buckling behaviour of hybrid plates. FSDT requires arbitrary shear correction factors. Correia et al. [20] employed coupled TOT to study optimum design of hybrid composite plates for buckling. From the evaluation of 2D theories for hybrid plates in comparison with 3D exact solutions [21–23], it is now well understood that for accurate global and local lamina level response of these structures (1) the kinematic field approximations must account for the layerwise nature of distributions of the inplane displacements across the thickness, (2) the electric field must allow for the nonlinear distribution of electric potential across the piezoelectric layers and (3) the two-way electromechanical coupling should be considered. FSDT and TOT are equivalent single layer theories, which use the same global variations for displacements across the entire laminate thickness and cannot account for the zigzag nature of variation of the inplane displacements. The first author [23] has presented a coupled efficient layerwise (zigzag) theory for linear static analysis of hybrid plates, which considers layerwise variations for the displacements, but the number of primary displacement variables is reduced to only five (as in FSDT and TOT), by enforcing the shear traction-free conditions at the top and the bottom surfaces and the transverse shear stress continuity conditions at the layer interfaces. Comparison of results of this theory with the exact 3D solutions for simply supported hybrid plates of highly inhomogeneous lay-ups established the high accuracy of this theory for static electromechanical response. A nonlinear coupled 1D zigzag theory for buckling of hybrid beams has also been presented and tested for accuracy [24]. This work presents an efficient coupled geometrically nonlinear zigzag theory for hybrid plates by extending the linear zigzag theory of Kapuria [23]. The geometric nonlinearity is included due to deflection only in the sense of Von Karman. The potential field is approximated as piecewise linear across a number of subdivisions in the piezoelectric layers. Both the inplane and the transverse electric fields are considered. The deflection field explicitly accounts for the transverse normal strain induced by the electric field. The inplane displace-
ments are approximated as a combination of a global third order variation across the thickness and a layerwise linear variation. The displacement field is expressed in terms of five primary displacement variables and a set of electric potential variables by enforcing exactly the conditions of zero transverse shear stresses at the top and bottom and their continuity at the layer interfaces. The nonlinear coupled equilibrium equations and boundary conditions are derived using a variational principle. This nonlinear theory is used to obtain the initial buckling response of symmetrically laminated hybrid plates under inplane electromechanical loading. The accuracy of the new theory for buckling is tested by direct comparison with the available 3D exact solution [14] for simply supported hybrid plates of highly inhomogeneous, composite and sandwich laminates.
2. Zigzag theory approximations Consider a hybrid plate (Fig. 1) made of L perfectly bonded orthotropic laminas of total thickness h with the midplane chosen as the xy-plane z = 0. Some of the layers are of orthorhombic piezoelectric materials of class mm2 symmetry, with poling along z-direction. The plate is subjected to inplane and transverse mechanical loads and actuation potentials at some piezoelectric surfaces. The reference plane z = 0 either passes through or is the bottom surface of the k0th layer. The kth ply from the bottom has bottom surface at z = zk1. Considering the usual assumption of 2D plate theories, rz ’ 0 [23], the 3D constitutive equations of piezoelasticity relating stresses r, s and electric displacements Dx, Dy, Dz with strains e, c and electric field components Ex, Ey, Ez reduce to eT Ez ; r ¼ Qe 3
^ ^eE; s ¼ Qc
D ¼ ^eT c þ ^gE;
Dz ¼ e3 e þ g33 Ez ;
Fig. 1. Geometry of hybrid plate.
ð1Þ
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
where 2
3
rx 6 7 r ¼ 4 r y 5; " c¼
sxy czx cyz
" s¼
# E¼
# ;
syz "
;
szx
Ex Ey
" D¼
Dx
2
#
Dy
;
contribution of rx, ry via PoissonÕs effect compared to that due to the electric field:, i.e., ez = w,z = S13r1 + S23r2 + d33Ez ’ d33/,z )
3
ex 6e 7 e ¼ 4 y 5;
j
cxy
#
ð2Þ
;
and for general angle-ply laminates 2 12 Q 16 3 " Q11 Q 55 Q 6 7 ^¼ ¼ 4Q 12 Q 22 Q16 5; Q Q 45 Q 16 Q26 Q 66 Q " # " # e15 e14 g12 g11 ^e ¼ ; ^ g¼ ; e3 e25 e24 g12 g22
# 45 Q ; 44 Q ¼ ½ e31
þ z2 nðx; yÞ þ z3 gðx; yÞ;
e32
e36 .
ij , eij , Q gij are the reduced elastic stiffnesses, piezoelectric stress constants and electric permittivities. Partially geometrically nonlinear strain–displacement relations in the spirit of Von Karman are employed. The Lagrange strain–displacement relations including geometric nonlinearity due to deflection w0(x, y) = w(x, y, 0) of the midplane, and the electric field-potential relations are
ez ¼ w;z ; Ex ¼ /;x ;
cyz ¼ uy;z þ w;y ; Ey ¼ /;y ;
where " u¼
ux
#
"
w0;x
#
ð7Þ "
uk x
#
; w0d ¼ ; uk ¼ ; uk y uy w0;y " # " # " # wkx gx nx wk ¼ ; n¼ ; g¼ ; wky gy ny
ð8Þ
uk is the translation and wk is related to the shear rotation of the kth layer. Substituting ux, uy, w from Eqs. (7) and (6), and / from Eq. (5) into Eq. (4) and using Eq. (1) yield the transverse shear stresses s as: ^ kW j ðzÞ/j ; ^ k ½wk þ 2zn þ 3z2 g þ ½^ek Wj ðzÞ Q s¼Q d / /
ð9Þ
T
ð4Þ
czx ¼ ux;z þ w;x ; Ez ¼ /;z ;
where ux, uy; w denote the inplane and transverse displacements and / denotes the electric potential. The subscript comma denotes differentiation. The inplane electric field components Ex, Ey are not considered as zero, since these may be induced by the piezoelectric coupling or may be applied through actuation of segmented piezoelectric actuators. The potential / is approximated as piecewise linear between n/ points zj/ ; j ¼ 1; 2; . . . ; n/ , across the thickness: /ðx; y; zÞ ¼ Wj/ ðzÞ/j ðx; yÞ;
ðzÞ/j ðx; yÞ; wðx; y; zÞ ¼ w0 ðx; yÞ W ð6Þ / R j ðzÞ ¼ z d 33 Wj ðzÞ dz is a piecewise linear funcwhere W /;z / 0 tion. The inplane displacements ux, uy for the kth layer are assumed to follow a global third order variation across the thickness with layerwise linear variation: uðx; y; zÞ ¼ uk ðx; yÞ zw0d ðx; yÞ þ zwk ðx; yÞ
ð3Þ
1 1 ex ¼ ux;x þ w20;x ; ey ¼ uy;y þ w20;y ; 2 2 cxy ¼ ux;y þ uy;x þ w0;x w0;y ;
255
ð5Þ
where /j ðx; yÞ ¼ /ðx; y; zj/ Þ. Wj/ ðzÞ are linear interpolation functions and summation convention is used for indices j, j 0 (used later). This description, which was first introduced by Kapuria [25] and later adopted by Wang [26], allows the piezoelectric layers to be divided into a number of sublayers and a series of elastic layers to be combined into one, for effective discretisation of / across the thickness. The variation of deflection w is obtained by integrating the constitutive equation for ez by neglecting the
where /jd ¼ ½/j;x /j;y . For the mid-plane which passes through the k0th layer, denote u0 ðx; yÞ ¼ uk0 ðx; yÞ ¼ uðx; y; 0Þ, w0 ðx; yÞ ¼ wk0 ðx; yÞ. Thus u0 and w0 are the midplane inplane displacements and shear rotations. Using the 2(L 1) conditions each for the continuity of s and u at the layer interfaces and the four shear traction-free conditions s = 0 at z = z0, zL, the (4L + 4) variables uk, wk, n, g are expressed in terms of only four variables u0 and w0 as in [23]. This yields uðx; y; zÞ ¼ u0 ðx; yÞ zw0d ðx; yÞ þ Rk ðzÞw0 ðx; yÞ þ Rkj ðzÞ/jd ðx; yÞ;
ð10Þ
where Rk(z), Rkj(z) are 2 · 2 matrices of layerwise cubic functions of z whose coefficients are dependent on the material properties and the lay-ups. Thus, even though u, w have layerwise distributions, they are expressed in terms of only five displacement variables u0x , u0y , w0, w0x , w0y by Eqs. (6) and (10).
3. Coupled nonlinear equations of equilibrium The variational principle for the piezoelectric medium [27] R zk be expressed, using the notation h i ¼ PL can k¼1 zþ ð Þ dz for integration across the thickness, as k1
256
Z
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
½hrx dex þ ry dey þ sxy dcxy þ syz dcyz þ szx dczx þ Dx d/;x
A
þ Dy d/;y þ Dz d/;z i p1z dwðx; y; z0 Þ p2z dwðx; y; zL Þ þ Dz ðx; y; z0 Þd/1 Dz ðx; y; zL Þd/n/ qji d/ji dA Z hrn dun þ sns dus þ snz dw þ Dn d/i ds ¼ 0; CL
8 du0 ; dw0 ; dw0 ; d/j ;
ð11Þ
where A denotes the mid-plane surface area of the plate and CL is the boundary curve of the midplane of the plate with normal n and tangent s. p1z ; p2z are the normal forces per unit area on the bottom and top surfaces of the plate in direction z. qji is the the extraneous surface j charge density at the interface z ¼ z/i , where /ji is prescribed. The total number of such prescribed potentials is n/ . The above variational equation is expressed in terms of du0 ; dw0 ; dw0 ; d/j and stress and electric displacement resultants to yield the nonlinear equilibrium equations and the boundary conditions. The stress resulT T tants N ¼ ½N x N y N xy , M ¼ ½M x M y M xy , P ¼ ½P x P yx j ¼ ½Q j P xy P y T , S j ¼ ½S j S j S j S j T , Q ¼ ½Qx Qy T , Q x T
yx
xy
y
x
j T , V ¼ ½V x V y , V j ¼ ½V j V j T and the electric disQ / /x /y y j
½H jx
T H jy ; Gj
placement resultants H ¼ are defined by h iT F 1 ¼ N T M T P T S j T ¼ ½hf3T ri; h iT j Q j F 2 ¼ Qx Qy Q ¼ ½hf4T si x y ð12Þ j si; H j ¼ hWj ðzÞDi; V ¼ hsi; V j ¼ hW j
G ¼
/ hWj/;z ðzÞDz i;
where f3 ¼ I 3
zI 3
/
/
k
U
kj
U
; f4 ¼
h
Rk;z
Rkj ;z
j ðzÞI 2 W /
N xy;x þ N y;y ¼ 0;
M x;xx þ 2M xy;xy þ M y;yy þ ðN x w0;x þ N xy w0;y Þ;x þ ðN xy w0;x þ N y w0;y Þ;x þ F 3 ¼ 0; P x;x þ P yx;y Qx ¼ 0;
þ H jy;y Gj þ F j6 ¼ 0;
j ¼ 1; 2; . . . ; n/ ;
where the mechanical load F 3 ¼ p1z þ p2z and the electrical loads F j6 ¼ Dz ðx; y; zL Þdjn/ Dz ðx; y; z0 Þdj1 þ qji djji . The boundary conditions on CL are the prescribed values of one of the factors of each of the following products: u0n N n ; u0s N ns ; w0 ðV n þ M ns;s þ N n w0;n þ N ns w0;s Þ; w0;n M n ; w0n P n ; w0s P ns ; /j;n S jn ; /j ½H jn V j/n S jns;s and at corners si:
w0 ðsi ÞDM ns ðsi Þ;
i ;
Rkj 22
0
/j ðsi ÞDS jns ðsi Þ
with V n ¼ ðM x;x þ M xy;y Þnx þ ðM xy;x þ M y;y Þny ; j nx Q j ny . V j ¼ ðS j þ S j Þnx þ ðS j þ S j Þny Q x;x
yx;y
y;y
xy;x
bj ¼ hf3T ðzÞeT3 Wj/;z ðzÞi;
ð13Þ
Rkj 12
0
j ¼ hf T ðzÞ^eWj0 ðzÞi; b / 4
x
y
The resultants F1, F2, Hj, Gj defined in Eq. (12) are expressed in terms of the displacement and potential variables by substituting the expressions of r, s, D, Dz from Eq. (1) into Eq. (12) and using Eqs. (4)–(6) and (10): 0 0 0 1 j / j0 ; e2 þ b F 1 ¼ Ae1 þ bj /j þ A Uw w0d ; F 2 ¼ A d 2 T 0 T 0 1 0 j e2 E jj /jd0 ; Gj ¼ bj e1 Ejj /j þ bjw Uw w0d ; H j ¼ b 2 ð16Þ
h iT e1 ¼ u0x ;x u0y ;y u0x ;y þ u0y ;x w0;xx w0;yy 2w0;xy w0x ;x w0x ;y w0y ;x w0y ;y /j;xx /j;xy /j;yx /j;yy ; 2 3 0 w0;x h iT j j e2 ¼ w0x w0y /;x /;y ; Uw ¼ 4 0 w0;y 5; w0;y w0;x T ^ f 4 ðzÞi; bjw ¼ he3 Wj i; ½A; A ¼ hf3 ðzÞQ½f3 ðzÞ; I 3 i; A ¼ hf4T ðzÞQ /;z 0
ð14Þ
P xy;x þ P y;y Qy ¼ 0;
j þ Q j Sj Sj Sj Sj þ H j Q x;x x;xx xy;xy yx;xy y;yy x;x y;y
/n
Rk21 Rk11 Rk22 Rk12 2 kj 3 0 R11 0 Rkj 12 6 7 kj kj 7 Ukj ¼ 6 4 0 R21 0 R22 5. Rkj 11
N x;x þ N xy;y ¼ 0;
ð15Þ
In being a n · n identity matrix and 2 k 3 R11 0 Rk12 0 6 7 k k 7 Uk ¼ 6 4 0 R21 0 R22 5;
Rkj 21
It can be shown that elements of Rk, Rkj, N, M, P, Sj transform as second order tensors and elements of j ; H j transform as vectors for the coplanar V ; V j/ ; Q; Q axes x, y and n, s. Using expressions of /, w, u from Eqs. (5), (6), (10) and using definitions in Eq. (12), the area integral in Eq. (11) is expressed in terms of du0x , du0y , dw0, dw0x , dw0y , d/j, by using GreenÕs theorem if required, and the terms involving du0x ; du0y ; dw0x ; dw0y ; dw0;x ; dw0;y ; d/j;x ; d/j;y in the integrand of CL are expressed in terms of components n, s. Thus Eq. (11) yields coupled nonlinear field equations consisting of five equations of force equilibrium and n/ equations for charge equilibrium:
0
0
Ejj ¼ hg33 Wj/;z ðzÞWj/;z ðzÞi;
0
jj ¼ h^gWj ðzÞWj0 ðzÞi; E / /
ð17Þ
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
where 2
A11 A12 6 6A 6 21 A22 6 6 . .. 6 .. . 6 6 6 A10;1 A10;2 A¼6 6 j 6 A11;1 Aj11;2 6 6 j 6 A12;1 Aj12;2 6 6 j 6 A13;1 Aj13;2 4 Aj14;1 Aj14;2 2 j0 3 b1 6 j0 7 6 b2 7 6 7 6 . 7 6 . 7 6 . 7 6 j0 7 6b 7 0 j 10 7 b ¼6 6 jj0 7; 6 b11 7 6 7 6 jj0 7 6 b12 7 6 07 6 jj 7 4 b13 5
. . . A1;10
0 Aj1;11
0 Aj1;12
0 Aj1;13
0
0
0
0 Aj1;14
3
7 7 7 7 7 7 7 7 j0 A10;14 7 7 ¼ AT ; 7 jj0 A11;14 7 7 0 7 Ajj12;14 7 7 0 7 Ajj13;14 7 5 0 Ajj14;14 0
. . . A2;10 Aj2;11 Aj2;12 Aj2;13 Aj2;14 .. .. .. .. .. .. . . . . . . 0
0
0
. . . A10;10 Aj10;11 Aj10;12 Aj10;13 0
0
0
0
0
0
0
0
0
0
0
0
. . . Aj11;10 Ajj11;11 Ajj11;12 Ajj11;13 . . . Aj12;10 Ajj12;11 Ajj12;12 Ajj12;13 . . . Aj13;10 Ajj13;11 Ajj13;12 Ajj13;13 . . . Aj14;10 Ajj14;11 Ajj14;12 Ajj14;13
Taking into account the corresponding zero elements of the matrices defined in Eq. (18), the nonlinear terms are obtained as Ln U Þ ¼ 1 ½A11 w2 þ A12 w2 þ ½A33 w0;x w0;y ; ðLn U 0;x 0;y ;x 1 ;y 2 Þ ¼ 1 ½A21 w2 þ A22 w2 þ ½A33 w0;x w0;y ; ðLn U 0;x 0;y ;y 2 ;x 2 n ðL U Þ3 ¼ N x w0;xx þ N y w0;yy þ 2N xy w0;xy þ w0;x ðN x;x þ N xy;y Þ 1 þ w0;y ðN xy;x þ N y;y Þ þ ½A41 w20;x þ A42 w20;y ;xx 2 1 þ ½A51 w20;x þ A52 w20;y ;yy þ ½2A63 w0;x w0; y;xy ; 2 1 Þ ¼ ½A71 w2 þ A72 w2 þ ½A83 w0;x w0;y ; ðLn U 0;x 0;y ;x 4 ;y 2 1 Þ ¼ ½A93 w0;x w0;y þ ½A10;1 w2 þ A10;2 w2 ; ðLn U 0;x 0;y ;y 5 ;x 2 1 1 1 Þ ¼ bjw ðLn U w2 þ bjw w2 þ ½Aj w2 þ Aj11;2 w20;y ;xx 5þj 2 1 0;x 2 2 0;y 2 11;1 0;x þ ½ðAj12;3 þ Aj13;3 Þw0;x w0;y ;xy 1 þ ½Aj14;1 w20;x þ Aj14;2 w20;y ;yy ; j ¼ 1; . . . ; n/ . 2
0
bjj14 2 A11 6 A21 6 6 . 6 . 6 . 6 6 A10;1 6 A ¼6 j 6 A11;1 6 6 j 6 A12;1 6 6 Aj 4 13;1
ð21Þ
3
A12 A13 A22 A23 7 7 2 .. 7 .. 7 11 A . 7 . 7 6 6 A10;2 A10;3 7 21 7 6A 7; A ¼ 6 j j j 6A A11;2 A11;3 7 7 4 31 7 Aj12;2 Aj12;3 7 j A 7 41 j j 7 A13;2 A13;3 5
Aj14;1 Aj14;2 Aj14;3 0 0 3 j j b b 11 12 6 j0 " jj0 0 7 6b b j 7 0 E j0 jj 6 7 21 22 ¼6 0 ¼ 110 b ; E 7 0 jj 6b jj 7 jj E 4 31 b 21 32 5 0 0 jj b jj b
0
0
12 A j A j A 13 14
The nonzero elements of load vector P are P 3 ¼ F 3 ; P j6 ¼ F j6 . The elements of L are listed in Ref. [23].
3
0 0 7 22 A j A j 7 A 23 24 7 T; 7¼A j jj0 jj0 7 A32 A33 A34 5 0 0 j A jj A jj A
42
43
44
2
41
257
0
jj E 12 0
jj E 22
#
; bjw ¼ bjw . bjw bjw 1 2 3
42
ð18Þ
Substituting Eq. (16) for the resultants into Eq. (14) yields the nonlinear coupled equations of equilibrium in terms of the primary displacement and potential vari: ables, U ¼ P ; þ Ln U ð19Þ LU where h iT ¼ u0x u0y w0 w0 w0 /1 /2 . . . /n/ ; U x y n T P ¼ P 1 P 2 P 3 P 4 P 5 P 16 P 26 . . . P 6/ .
4. Buckling under uniform electromechanical load 4.1. Pre-buckling solution Consider a symmetrically laminated plate which is in equilibrium under applied uniform inplane normal strains e0x ; e0y , zero shear strain c0xy ¼ 0 and actuation potentials independent of x, y coordinates. Let this pre-buckling equilibrium state be denoted by superscript ( )0. For the symmetrically laminated plate under symmetrical loading about xy-plane,
0
0
0
0
N 0x ¼ A11 e0x þ A12 e0y þ bj1 /0j ; N 0y ¼ A12 e0x þ A22 e0y þ bj2 /0j ; 0j0
¼
Defining U ¼ /1
/2
...
/n/
h bj ¼ b1j
b2j
...
bj /
þ
bj2 e0y
F 0j 6 .
bj1 e0x
E /
L is a symmetric matrix of linear differential opera are the nonlinear terms due to geotors in x and y. Ln U 45 ¼ 0; metric nonlinearity. For cross-ply plates, Q Q16 ¼ Q26 ¼ 0; e14 ¼ e25 ¼ 0, b6 ¼ 0; g12 ¼ 0; e36 ¼ 0.
ð22Þ
Considering this, the plate constitutive Eq. (16) and equilibrium Eq. (14) yield
jj0
ð20Þ
w0x ¼ w0y ¼ 0.
w00 ¼ 0;
n
T iT
;
ð23Þ
258
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
and
F 6 ¼ F 16
F 26
...
n F 6/
T
buckling is obtained for simply supported rectangular plates of sides a, b along the axes x, y for the boundary conditions
:
Eq. (23)3 can be written in matrix form as EU0 ¼ b1 e0x þ b2 e0y F 06 .
ð24Þ
In order to consider open as well as closed circuit electric boundary conditions, U0 is partitioned into a set of unknown output voltages U0s at zj/ Õs where / is not prescribed and a set of known input actuation voltages U0a at the actuated surfaces. Accordingly, Eq. (24) can be partitioned and arranged as " # " 0 # b2s 0 b1s 0 Ess Esa F 06s Us e þ e ¼ . Eas Eaa U0a b1a x b2a y F 06a ð25Þ Eq. (25) is solved for U0s and the same is substituted into Eq. (23) to obtain the pre-buckling forces N 0x ; N 0y in terms of the known loading parameters T 1 0 0 N 0x ¼ ðA11 þ bT1s E1 ss b1s Þex þ ðA12 þ b1s Ess b2s Þey T 1 0 0 þ ðbT1a bT1s E1 ss Esa ÞUa b1s E ss F 6s ; T 1 0 0 N 0y ¼ ðA12 þ bT2s E1 ss b1s Þex þ ðA22 þ b2s Ess b2s Þey
ð26Þ
T 1 0 0 þ ðbT2a bT2s E1 ss Esa ÞUa b2s E ss F 6s .
4.2. Stability equations Let the solution for just after buckling be denoted by ð^Þ on the entities. The governing equations for buckling are derived using the approach of Noor and Burton [28] by introducing a small parameter for the buckling mode. is described by an The size of the buckling mode U ^ ¼U 0 þ U with U given arbitrary small parameter : U by Eq. (20)1. Substituting this solution into Eqs. (19), using Eq. (22) and considering upto the first order terms : in yield the following stability equations for U h iT þ 0 0 ðN 0x w0;xx þ N 0y w0;yy Þ 0 0 0 LU T ¼ P ¼ 0 0 0 0 0 F j6 ; ð27Þ ^ Þ ; ðLn U ^ Þ ; ðLn U ^ Þ ; ðLn U ^ Þ ; ðLn U ^ Þ are of secsince ðLn U 1 2 4 5 5þj 2 n^ 0 ond order and ðL U Þ3 ¼ ðN x w0;xx þ N 0y w0;yy Þ is of first j order. The load F 6 is the contribution of the incremental only. For a set of zero incremental potential solution U at the actuator locations and zero incremental electric displacement at the unknown potential locations, the j load F^ 6 is zero for index j corresponding to such surfaces. Note that whereas the equilibrium equations (14) are nonlinear, the stability equations (27) are linear in terms of pre-buckling loads N 0x ; N 0y . To assess the accuracy of the theory developed herein, by comparison with the exact 3D piezoelasticity solution [14], analytical Navier solution of Eq. (27) for
at x ¼ 0; a:
N x ; u0y ; w0 ; w0y ; M x ; P x ; /j ; S jx ¼ 0;
at y ¼ 0; b:
N y ; u0x ; w0 ; w0x ; M y ; P y ; /j ; S jy ¼ 0.
ð28Þ
The solution for the (m, n)th spatial mode of buckling is taken as: 2 3 2 3 sinðnyÞ ½ w0 /j mn sinðmxÞ w0 /j 6u 7 6 7 sinðnyÞ 5 4 0x w0x 5 ¼ 4 ½ u0x w0x mn cosðmxÞ u0y w0y cosðnyÞ ½ u0y w0y mn sinðmxÞ ð29Þ ¼ mp=a; n ¼ np=b. Substituting these into Eq. with m (27) yields ðK K G ÞU
mn
mn
¼ P ;
ð30Þ
KG is the geometric stiffness matrix with only one non 2 N 0x n2 N 0y . K is symmetric zero element K G ð3; 3Þ ¼ m stiffness matrix whose elements are not listed here for brevity. Partitioning electric potentials U into the unknown part Us and known zero incremental actuation potentials, Ua as before, Eq. (30) can be written for the T ~ mn ¼ u0x mn u0y mn w0mn w0x mn w0y nm Umn unknowns U s as: ~ K ~ G ÞU ~ mn ¼ ðK ~ kK ÞU ~ mn ¼ 0; ðK G
ð31Þ
~G ¼ where K Eq. (31) is a generalised linear eigenvalue problem and the eigenvalue k is the buckling load factor. The eigenvalues and eigenvectors are obtained by QR algorithm after reducing to Heissenberg form using NAG subroutines. kK G .
5. Numerical evaluation of the theory The accuracy of the present theory for buckling response is established by direct comparison with the exact 3D piezoelasticity [14] solution for simply supported cross-ply symmetrically laminated hybrid plates. Since TOT has the same global third order variation for the displacement field without the layerwise terms and has the same number of displacement variables as the present theory (ZIGT), present results are also compared with the coupled TOT [20,29] extended for the buckling case. This comparison will establish the effect of the layerwise terms for the displacements incorporated in the ZIGT. Results are presented for hybrid plates of three different laminate configurations (a), (b) and (c). The stacking order is mentioned from the bottom. The elastic substrate of plate (a) has five plys of equal thickness 0.16h of materials 2/3/1/3/2 with orientations of the principal material axis 1 as [90/0/0/0/90]. It is a good test case for assessing a 2D theory since the plies have highly
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
inhomogeneous stiffness in tension and shear. The substrate of plate (b) is a graphite–epoxy composite laminate of material 4 with four layers of equal thickness 0.2h with lay-up [0/90/90/0]. The substrate of plate (c) is a five-layer sandwich having graphite–epoxy faces [0/90] and a soft core with thicknesses .04h/.04h/.64h/ .04h/.04h. All the plates have two PZT-5A layers, each of thickness 0.1h, bonded to their elastic substrate on its top and bottom surfaces. The PZT-5A layers have poling in +z-direction. The top and the bottom of the substrate are grounded. The material properties are selected as [14]: [(Y1, Y2, Y3, G12, G23, G31), m12, m13, m23] = Material 1: [(6.9, 6.9, 6.9, 1.38, 1.38, 1.38) GPa, 0.25, 0.25, 0.25] Material 2: [(224.25, 6.9, 6.9, 56.58, 1.38, 56.58) GPa, 0.25, 0.25, 0.25] Material 3: [(172.5, 6.9, 6.9, 3.45, 1.38, 3.45) GPa, 0.25, 0.25, 0.25] Material 4: [(181, 10.3, 10.3, 7.17, 2.87, 7.17) GPa, 0.28, 0.28, 0.33] Face: [(131.1, 6.9, 6.9, 3.588, 2.3322, 3.588) GPa, 0.32, 0.32, 0.49] Core: [(0.2208, 0.2001, 2760, 16.56, 455.4, 545.1) MPa, 0.99, 3 · 105, 3 · 105] PZT-5A: [(61.0, 61.0, 53.2, 22.6, 21.1, 21.1) GPa, 0.35, 0.38, 0.38], and [(d31, d32, d33, d15, d24), (g11, g22, g33)] = [(171, 171, 374, 584, 584) · 1012 m/V, (1.53, 1.53, 1.5) · 108 F/m] where Yi, Gij, mij, dij and gij denote YoungÕs moduli, shear moduli, PoissonÕs ratios, piezoelectric strain constants and electric permittivities, respectively. Four pre-buckling load cases are considered: 1. Uniaxial inplane strain e0x ¼ e0 with e0y ¼ 0. The top and bottom surfaces are at closed circuit condition with zero potential ½/1 ¼ /n/ ¼ 0. 2. Uniaxial inplane strain e0y ¼ e0 with e0x ¼ 0 and the top and the bottom surfaces at zero potential.
259
3. Biaxial inplane strain e0x ¼ e0y ¼ e0 with top and bottom surfaces at zero potential. 4. Applied actuation potentials at top and bottom surfaces /n/ ¼ /1 ¼ /0 with the immovable edges (e0x ¼ e0y ¼ 0). For load cases 1 to 3, the lowest value of the strain e0 for buckling is denoted as ecr. For load case 4, the critical value of potential /0 for buckling is denoted as /cr. The results for these cases are nondimensionalised with S = a/h, d0 = 374 · 1012 CN1, Y0 = 6.9 GPa for laminates (a) and (c), and 10.3 GPa for laminate (b): ¼ / d 0 S 3 =a; ecr ¼ ecr S 2 ; / cr cr Þ ¼ ðSu; Sv; wÞ= maxðwÞ; ðu; v; w y Þ ¼ ðrx ; ry ÞS 2 h=Y 0 maxðwÞ; ð rx ; r ðszx ; syz Þ ¼ ðszx ; syz ÞS 3 h=Y 0 maxðwÞ; where max(w) denotes the largest value of w through the thickness. The critical buckling strains ecr obtained from the exact 3D piezoelasticity solution [14] and percentage (%) errors of the present theory and coupled TOT for square hybrid test plate (a) are presented in Table 1 for load cases 1, 2 and 3. Four values of the span-to-thickness ratio S = 5, 7.5, 10, 20 are considered for comparison. The values of (m, n) corresponding to the critical buckling mode are also given in the table. Similar results for composite plate (b) and sandwich plate (c) are presented in Tables 2 and 3, respectively. It is observed that the present theory yields highly accurate results for ecr for all laminate configurations and load cases with a maximum error of 1.5% even for thick plates with S = 5. In contrast, TOT results have errors upto 133.9%, 11.0% and 55.9% for plates (a), (b) and (c), respectively. Even for thin plates with S = 20, the error in TOT is as high as 15.4% for test plate (a) and 12.1% for sandwich plate (c). Even the critical mode (m, n) is wrongly estimated by TOT in some cases. The 3D exact results and % errors of (load ZIGT and TOT for the critical electric potential / cr case 4) are presented in Table 4 for the three square
Table 1 Exact 3D results and % errors of ZIGT and TOT for buckling inplane strain for plate (a) S
Entity
Load case 1 Exact
5 7.5 10 20
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
0.65980 (2, 1) 1.19904 (2, 1) 1.73289 (2, 1) 2.88725 (1, 1)
Load case 2 % Error
Exact
ZIGT
TOT
0.3 (2, 1) 1.0 (2, 1) 1.5 (2, 1) 0.0 (1, 1)
59.1 (3, 1) 56.2 (2, 1) 43.3 (2, 1) 15.2 (1, 1)
0.61949 (1, 2) 1.09362 (1, 2) 1.54017 (1, 1) 2.38148 (1, 1)
Load case 3 % Error
Exact
ZIGT
TOT
0.5 (1, 2) 0.2 (1, 2) 0.3 (1, 1) 0.2 (1, 1)
133.9 (1, 1) 81.4 (1, 1) 49.4 (1, 1) 15.4 (1, 1)
0.3744 (1, 1) 0.6280 (1, 1) 0.8449 (1, 1) 1.3048 (1, 1)
% Error ZIGT
TOT
0.6 (1, 1) 0.0 (1, 1) 0.1 (1, 1) 0.1 (1, 1)
111.9 (1, 1) 72.9 (1, 1) 49.1 (1, 1) 15.3 (1, 1)
260
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
Table 2 Exact 3D results and % errors of ZIGT and TOT for buckling inplane strain for plate (b) S
Entity
Load case 1
Load case 2
Exact
5 7.5 10 20
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
% Error
0.92953 (2, 1) 1.35619 (1, 1) 1.59769 (1, 1) 1.95146 (1, 1)
Exact
ZIGT
TOT
0.8 (2, 1) 0.1 (1, 1) 0.1 (1, 1) 0.1 (1, 1)
11.0 (2, 1) 6.3 (1, 1) 4.5 (1, 1) 1.4 (1, 1)
0.97300 (1, 1) 1.36354 (1, 1) 1.60444 (1, 1) 1.95454 (1, 1)
Load case 3 % Error
Exact
ZIGT
TOT
1.3 (1, 1) 0.6 (1, 1) 0.4 (1, 1) 0.1 (1, 1)
7.8 (1, 1) 5.7 (1, 1) 4.0 (1, 1) 1.3 (1, 1)
0.48543 (1, 1) 0.68011 (1, 1) 0.80061 (1, 1) 0.97651 (1, 1)
% Error ZIGT
TOT
1.0 (1, 1) 0.4 (1, 1) 0.2 (1, 1) 0.0 (1, 1)
8.0 (1, 1) 5.9 (1, 1) 4.2 (1, 1) 1.4 (1, 1)
Table 3 Exact 3D results and % errors of ZIGT and TOT for buckling inplane strain for plate (c) S
5 7.5 10 20
Entity
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
Load case 1
Load case 2
Exact
% Error ZIGT
TOT
0.71404 (2, 1) 1.39331 (2, 1) 2.08738 (1, 1) 3.59027 (1, 1)
0.1 (2, 1) 0.4 (2, 1) 0.2 (1, 1) 0.1 (1, 1)
50.9 (2, 1) 48.9 (1, 1) 31.8 (1, 1) 12.1 (1, 1)
Load case 3
Exact
% Error
Exact
ZIGT
TOT
0.67831 (1, 2) 1.28178 (1, 2) 1.98108 (1, 2) 3.59088 (1, 1)
0.2 (1, 2) 0.5 (1, 2) 0.4 (1, 2) 0.1 (1, 1)
53.3 (1, 2) 55.9 (1, 2) 38.9 (1, 1) 12.1 (1, 1)
0.40547 (1, 1) 0.73426 (1, 1) 1.04388 (1, 1) 1.79529 (1, 1)
% Error ZIGT
TOT
0.7 (1, 1) 0.5 (1, 1) 0.2 (1, 1) 0.1 (1, 1)
50.9 (1, 1) 41.3 (1, 1) 31.8 (1, 1) 12.1 (1, 1)
Table 4 Exact 3D results and % errors of ZIGT and TOT for critical electric potential for buckling S
Entity
Plate (a) Exact
5 7.5 10 20
/ cr (m, n) / cr (m, n) / cr (m, n) / cr (m, n)
Plate (b) % Error
0.36515 (1, 1) 0.62092 (1, 1) 0.84116 (1, 1) 1.31539 (1, 1)
Exact
ZIGT
TOT
3.5 (1, 1) 2.7 (1, 1) 2.1 (1, 1) 0.9 (1, 1)
116.2 (2, 1) 77.7 (1, 1) 52.2 (1, 1) 16.2 (1, 1)
0.52176 (1, 1) 0.74734 (1, 1) 0.89057 (1, 1) 1.10546 (1, 1)
Plate (c) % Error
Exact
ZIGT
TOT
5.1 (1, 1) 3.5 (1, 1) 2.4 (1, 1) 0.8 (1, 1)
14.7 (1, 1) 10.0 (1, 1) 6.9 (1, 1) 2.2 (1, 1)
0.14028 (1, 1) 0.25532 (1, 1) 0.36407 (1, 1) 0.62974 (1, 1)
% Error ZIGT
TOT
1.1 (1, 1) 0.8 (1, 1) 0.8 (1, 1) 0.5 (1, 1)
53.6 (1, 1) 43.1 (1, 1) 33.1 (1, 1) 12.6 (1, 1)
Table 5 Effect of b/a on the errors of ZIGT and TOT results for buckling parameters for plate (a) b/a
Entity
Load case 1 Exact
Entity % Error ZIGT
1 2 4 10
ecr (m, n) ecr (m, n) ecr (m, n) ecr (m, n)
1.73289 (2, 1) 0.79996 (1, 1) 0.56101 (1, 1) 0.49803 (1, 1)
1.5 (2, 1) 1.2 (1, 1) 1.8 (1, 1) 2.1 (1, 1)
Load case 4 Exact
TOT 43.3 (2, 1) 17.3 (1, 1) 9.3 (1, 1) 7.1 (1, 1)
/ cr (m, n) / cr (m, n) / cr (m, n) / cr (m, n)
0.84116 (1, 1) 0.54542 (1, 1) 0.44376 (1, 1) 0.41278 (1, 1)
% Error ZIGT
TOT
2.1 (1, 1) 2.7 (1, 1) 3.1 (1, 1) 3.3 (1, 1)
52.2 (1, 1) 19.0 (1, 1) 10.6 (1, 1) 8.3 (1, 1)
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
261
y ; syz for critical buckling mode of square hybrid plate (a) under load case 1. Fig. 2. Distributions of v; r
x ; szx for critical buckling mode of square hybrid plate (b) under load case 1. Fig. 3. Distributions of u; r
plates. It is revealed that the maximum error in ZIGT for thick plates with S = 7.5 is 3.5% whereas the error in TOT is 77.7%, 10.0%, 43.1% for plates (a), (b), (c),
respectively. Even for thin plates with S = 20, the error is as high as 16.2% and 12.6% for plates in TOT for / cr (a) and (c). The errors in the 2D theories are dependent
262
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
x ; szx for critical buckling mode of square hybrid plate (c) under load case 1. Fig. 4. Distributions of u; r
x ; szx for critical buckling mode of square hybrid plate (a) under potential load (case 4). Fig. 5. Distributions of u; r
on the uniaxial or bi-axial nature of inplane loads as well as the laminate configurations. The extent of error in
TOT results for plate (a) confirms that this plate indeed offers a good test case for assessing the accuracy of 2D
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
263
y ; syz for critical buckling mode of plate (a) under load case 1. Fig. 6. Effect of b/a on distributions of v; r
theories for buckling. The present theory successfully passes this test. As is well known, equivalent single layer theories like TOT yield highly inaccurate results for sandwich plates, whereas ZIGT yields excellent results for these plates too. for modThe errors in ZIGT and TOT for ecr and / cr erately thick plate (a) with S = 10 are presented in Table 5 for four values of b/a ratio ranging from 1 to 10. It is observed that whereas the error in ZIGT generally shows a marginal increase, the error in TOT reduces with the increase in b/a. The through-the-thickness distributions of predominant modal inplane displacements u=v, normal stresses x = r ry and transverse shear stresses szx =syz for the critical buckling mode at (x, y) locations where they are maximum are compared in Figs. 2–4 for square hybrid plates (a), (b) and (c) for load case 1 for S = 5 and 10. The location parameters a; b are defined as a ¼ a=2m; b ¼ b=2n where (m, n) correspond to the critical buckling mode. It is observed that the present distributions are in excellent agreement with the exact modal distributions for all plates even for the thick ones with S = 5. In contrast, the distributions obtained from TOT have large error even for S = 10. The present theory is able to accurately capture the layerwise zigzag distribution of the modal inplane displacements obtained from the 3D solution, but the global distribution of the inplane displacements in TOT is unable to capture it. This explains why ZIGT yields accurate results and TOT does not.
The presented results cover a wide variety of possible cases as is evident from the qualitative distinct nature of the through-the-thickness distributions of the modal entities for the three plates. The distributions of the modal displacements and stresses across the thickness for square plate (a) under electric potential load of case 4 are compared in Fig. 5. For this case too, the present distributions are accurate. The effect of b/a ratio on the accuracy of the through-the-thickness distributions is illustrated for plate (a) in Fig. 6 for S = 10. It is seen that whereas ZIGT yields accurate distributions for all b/a ratios, the error in the TOT reduces with the increase in b/a.
6. Conclusions An efficient coupled zigzag theory incorporating geometric nonlinearity has been developed for hybrid piezoelectric plates. The theory considers layerwise distribution of the inplane displacements across the thickness and accounts for the transverse normal strain caused by the potential field, but the displacement field is expressed in terms of only five displacement variables as in TOT. The nonlinear theory is used to obtain the buckling response of symmetrically laminated plates due to inplane electromechanical loading. The accuracy of the theory is established directly by comparison with the exact 3D piezoelasticity solution for a simply
264
S. Kapuria, G.G.S. Achary / Composite Structures 74 (2006) 253–264
supported highly inhomogeneous hybrid test plate, a composite plate and a sandwich plate under uniaxial and bi-axial inplane strains and electric potential load, covering a wide range of possible cases. The comparisons establish that the present theory can accurately predict the zigzag distributions of inplane displacements as obtained from 3D exact solutions and hence yields accurate results for the buckling loads and modal distributions of the response entities. In contrast, TOT results may be highly inaccurate even for thin plates with S = 20. The present theory is also applicable for hybrid sandwich plates, for which TOT yields poor results. Acknowledgements The first author is thankful to Department of Science and Technology, Government of India, for providing financial assistance for this work under SERC scheme. References [1] Leissa AW. A review of laminated composite plate buckling. Appl Mech Rev 1987;40:575–91. [2] Chia CY. Geometrically nonlinear behavior of composite plates: A review. Appl Mech Rev 1988;41:439–51. [3] Kam TY, Chang PR. Buckling of shear deformable laminated composite plates. Compos Struct 1992;22(4):223–34. [4] Balamurugan M, Ganapathi M, Varadan TK. Nonlinear dynamic instability of laminated composite plates using finite element method. Compos Struct 1996;60(1):125–30. [5] Fares ME, Zankour AM. Buckling and free vibration of nonhomogeneous composite cross-ply laminated plates with various plate theories. Compos Struct 1999;44(4):279–87. [6] Chattopadhyay A, Radu AG. Dynamic instability of composite laminates using a higher order theory. Compos Struct 2000;77: 453–60. [7] Chakrabarti A, Sheikh AH. Buckling of laminated composite plates by a new element based on higher order shear deformation theory. Adv Mater Struct 2003;10(4):303–17. [8] Librescu L, Hause T, Jhonson T. Buckling and nonlinear response sandwich curved panels to combined mechanical loads-implications of face-sheet elastic tailoring. J Sandwich Struct Mater 2000;2:246–69. [9] Yuan WX, Dawe DJ. Free vibration and stability analysis of stiffened sandwich plates. Compos Struct 2004;63(1):123–37. [10] Srinivas S, Rao AK. Bending, vibration and buckling of simply supported thick orthotropic rectangular plates and laminates. Int J Solids Struct 1970;11:1463–81.
[11] Gu H, Chattopadhyay A. Three-dimensional elasticity solution for buckling of composite laminates. Compos Struct 1996;34: 292–9. [12] Thomson SP, Loughan J. The control of the post-buckling response in thin composite plates using smart technology. ThinWalled Struct 2000;36:231–63. [13] Batra RC, Geng TS. Enhancement of dynamic buckling load for a plate by using piezoceramic actuators. Smart Mater Struct 2001; 10:925–33. [14] Kapuria S, Achary GGS. Exact 3D piezoelasticity solution for buckling of hybrid cross-ply plates using transfer matrices. Acta Mech 2004;170(1–2):25–46. [15] Yang JS. Buckling of a piezoelectric plate. Int J Appl Electromagnetics Mech 1998;9:399–408. [16] Hu YT, Yang JS, Jiang Q. A model for electroelastic plates under biasing fields with applications in buckling analysis. Int J Solids Struct 2002;39:2629–42. [17] Chandrashekhara K, Bhatia K. Active buckling control of smart composite plates-finite-element analysis. Smart Mater Struct 1993; 2:31–9. [18] Varelis D, Saravanos DA. Nonlinear coupled mechanics and initial buckling of composite plates with piezoelectric actuators and sensors. Smart Mater Struct 2002;11:330–6. [19] Varelis D, Saravanos DA. Coupled buckling and postbuckling analysis of active laminated piezoelectric composite plates. Int J Solids Struct 2004;41(5–6):1519–38. [20] Correia VMF, Soares CMM, Soares CAM. Buckling optimization of composite laminated adaptive structures. Compos Struct 2003;62(3–4):315–21. [21] Yang JS. Equations for thick elastic plates with partially electroded piezoelectric actuators and higher order electric field. Smart Mater Struct 1999;8:73–82. [22] Bisegna P, Caruso G. Evaluation of higher-order theories of piezoelectric plates in bending and in stretching. Int J Solids Struct 2001;38:8805–30. [23] Kapuria S. A coupled zigzag third order theory for hybrid crossply plates. ASME J Appl Mech 2004;71(5):604–14. [24] Kapuria S, Alam N. Zigzag theory for buckling of hybrid piezoelectric beams under mechanical loads. Int J Mech Sci 2004; 46:1–25. [25] Kapuria S. An efficient coupled theory for multi-layered beams with embedded piezoelectric sensory and active layers. Int J Solids Struct 2001;38:9179–99. [26] Wang SY. A finite element model for the static and dynamic analysis of a piezoelectric bimorph. Int J Solids Struct 2004;41: 4075–96. [27] Tiersten HF. Linear piezoelectric plate vibrations. New York: Plenum Press; 1969. 25–37, 54, 55. [28] Noor AK, Burton WS. Three-dimensional solutions for thermal buckling of multilayered anisotropic plates. J Eng Mech 1992; 118(4):683–701. [29] Kapuria S, Achary GGS. A coupled consistent third order theory for hybrid piezoelectric plates. Compos Struct 2005;70(1):120–33.