Jtam-A4.dvi JOURNAL OF THEORETICAL AND APPLIED MECHANICS 51, 4, pp. 813-826, Warsaw 2013 BUCKLING ANALYSIS OF THREE-LAYERED RECTANGULAR PLATE WITH PIEZOELECTRIC LAYERS Hessameddin Yaghoobi Young Researchers and Elite Club, Central Tehran Branch, Islamic Azad University, Tehran, Iran e-mail: yaghoobi.hessam@gmail.com Iraj Rajabi K.N. Toosi University of Technology, Faculty of Mechanical Engineering, Tehran, Iran This paper employs an analytical method to analyze the buckling of piezoelectric coupled plates with different boundary conditions on the basis of the first order shear deformation plate theory.The structure is composedofahost isotropicplateand twobondedpiezoelectric layers. Convergence study is performed in order to verify the numerical stability of the presented method. Also, the present analysis is validated by comparing the results with those in the literature, and then the critical buckling load of the piezoelectric coupled plates is presented in tabular and graphical forms for different plate aspect ratios, thickness of the piezoelectric, actuator voltage and boundary conditions. Key words: three-layered rectangular plate, piezoelectric, critical buckling load 1. Introduction Piezoelectric materials have strongly attracted the attention of many research groups due to their unique electromechanical coupling characteristics which producemechanical deformations under application of electrical loads and electrical fields under application of mechanical loads. The main advantages of these types of smart materials are high precision, low weight and high sensibility. Smart structures (e.g., piezoelectric coupled plates) may be used as sensors or/and actuators in various engineering applications including vibration control, acoustic noise suppression, active damping and so on. They are commonly used as an embedded layer on host structures. As a result, a thorough understanding of the interaction between the host structure and piezoelectric layer is helpful in order to effectively utilize this combination in different applications. Shape or vibration control of laminated plates with integrated piezoelectric sensors and ac- tuators has been identified as an important field of study in recent years. However, relatively fewworks have been done on the compressive and/or thermal buckling of plates containing pie- zoelectric layers. Oh et al. (2000) investigated postbuckling and vibration analysis considering large thermopiezoelastic deflections for fully symmetric and partially eccentric piezolaminated composite plates. The non-linear finite element equations based on the layerwise displacement theorywere formulated for piezolaminated plates subject to thermal andpiezoelectric loads.The Newton-Raphson iteration method was used to solve the non-linear equation. Shen (2001a,b) analyzed compressive and thermal postbuckling of shear deformable laminated plates with fully covered or embedded piezoelectric actuators subjected to combined mechanical, electrical and thermal loads. A higher order shear deformation plate theory was adopted and the initial geo- metric imperfection of plates was accounted. It was found that the control voltage had a small effect on the postbuckling load-deflection relationship of shear deformable piezolaminated pla- tes with immovable unloaded edges, and almost no effect on the postbuckling load-deflection relationships of the same plate with movable edges. 814 H. Yaghoobi, I. Rajabi Coupledmulti-field generalized non-linearmechanics together with an associated plate finite element for analyzing the buckling andpostbuckling responseof active and sensorypiezoelectric- composite laminatedplates includingnon-linear effects dueto large rotations andstress stiffening were presented byVarelis andSaravanos (2004). The discrete coupled equations ofmotion of the smart structure were finally linearized and solved using an incremental-iterative method based on the Newton-Raphson technique. Kapuria and Achary (2006, 2008) employed 3D elasticity and zigzag theory for linear compressive and thermal buckling of laminated plates containing piezoelectric layers. Akhras and Li (2008) extended the finite layer method to thermal buckling analysis of rectangular simply supported symmetrical cross-ply piezoelectric composite plates. Using thismethod, the three-dimensional analysiswas transformed intoone-dimensional analysis by virtue of the orthogonal properties of trigonometric interpolation functions. Buckling optimization of laminated plates with integrated piezoelectric actuators can be fo- und in Correia et al. (2003). Assessment of third order smeared and zigzag theories for buckling and vibration of symmetrically laminated hybrid angle-ply plates containing piezoelectric layers can be found inDumir et al. (2009). Shariyat (2009) studied dynamic buckling of laminated pla- tes with piezoelectric sensors and actuators under thermo-electro-mechanical loads using a finite element formulation based on ahigher-order shear deformation theory.Anine-node second order Lagrangian element, an efficient numerical algorithm for solving the resulted highly non-linear governing equations, and an instability criterion already proposed by the authorwere employed. Shen and Zhu (2011) investigated compressive postbuckling under thermal environments and thermal postbuckling due to a uniform temperature rise for a shear deformable laminated plate with piezoelectric fiber reinforced composite (PFRC) actuators based on a higher order shear deformation plate theory that includes thermo-piezoelectric effects. In the present research, buckling analysis of a three-layered rectangular plate with piezoelec- tric layers is investigated.Basedon thefirst order sheardeformationplate theory, the equilibrium and stability equations are obtained. Introducing a new analytical method, the coupled stability equations are converted into independent partial differential equations. It is assumed that the plate is simply supported on two opposite edges and has arbitrary boundary conditions along the other edges. By using the Levy solution, these equations are converted into two ordinary differential equations, one of which has variable coefficients. For solving the equations accura- tely, the power series method of Frobenius (see Wylie and Barrett, 1951) is used. To examine accuracy of the present formulation and procedure, several convergence and comparison studies are investigated. Also, the effects of some parameters of the plate and piezoelectric layers on the critical buckling load are studied. 2. Stability equations Consider a three-layered rectangular plate, made of an isotropic substrate of thickness h and piezoelectric films of thickness ha that are perfectly bonded on its top and bottom surfaces as actuators, as shown in Fig. 1. The length and the width of the plate are a and b, respectively. Rectangular Cartesian coordinates (x,y,z) are assumed for derivations in this study. The displacement components of the plate based on the first-order shear deformation plate theory are considered as (Reddy, 1984, 2004) u(x,y,z) = u0(x,y)+zψx(x,y) v(x,y,z) = v0(x,y)+zψy(x,y) w(x,y,z) = w0(x,y) (2.1) where u, v and w represent the displacement of the plate in the x, y and z directions, respecti- vely, u0 and v0 are the displacements of the mid-plane, w0 is the transverse displacement, and ψx and ψy show rotation terms about y and x axes, respectively. The parameters u0, v0, w0, Buckling analysis of three-layered rectangular plate ... 815 Fig. 1. Coordinate system and geometry of a rectangular plate integrated with piezoelectric layers ψx and ψy are all functions of x and y variables. In this theory, the transverse normal do not remain perpendicular to themid-surface after deformation. This amounts to including transver- se shear strains in the theory. The inextensibility of the transverse normal requires that w not be a function of the thickness coordinate z (Reddy, 2004). Using the non-linear form of strain-displacement relations, the following strain components are obtained (Reddy, 2004)              εxx εyy γyz γxz γxy              =              u,x+w 2 ,x/2 v,y +w 2 ,y/2 ψy +w,y ψx+w,x u,y+v,x+w,xw,y              =                  ε (0) xx ε (0) yy γ (0) yz γ (0) xz γ (0) xy                  +z                ε (1) xx ε (1) yy 0 0 γ (1) xy                (2.2) where                  ε (0) xx ε (0) yy γ (0) yz γ (0) xz γ (0) xy                  =              u0,x+w 2 0,x/2 v0,y+w 2 0,y/2 ψy +w0,y ψx+w0,x u0,y +v0,x+w0,xw0,y              and        ε (1) xx ε (1) yy γ (1) xy        =      ψx,x ψy,y ψx,y +ψy,x      (2.3) Subscript ( , ) denotes derivationwith respect to the coordinates. The constitutive law for hybrid rectangular plates, taking into account the piezoelectric effect, is given by Liew et al. (2003)              σxx σyy τyz τxz τxy              =        Q11 Q12 0 0 0 Q21 Q22 0 0 0 0 0 Q44 0 0 0 0 0 Q55 0 0 0 0 0 Q66                     εxx εyy γyz γxz γxy              −        0 0 e31 0 0 e32 0 0 0 0 e24 0 e15 0 0             Ex Ey Ez      (2.4) where Qij (i,j =1, . . . ,6) is the elastic stiffness of the layers given by Q11 = Q22 = E 1−υ2 Q12 = Q21 = υE 1−υ2 Q44 = Q55 = Q66 = E 2(1+υ) (2.5) The piezoelectric stiffness e31, e32, e15 and e24 can be expressed in terms of the dielectric constants d31, d32, d15 and d24. The elastic stiffness Q a ij (i,j = 1, . . . ,6) of the piezoelectric actuator layers as e31 =(d31Q a 11+d32Q a 21) e32 =(d31Q a 12+d32Q a 22) e24 = d24Q a 44 e15 = d15Q a 55 (2.6) 816 H. Yaghoobi, I. Rajabi As only the transverse electric field component Ez is dominant in the plate type piezoelectric material, it is assumed that      Ex Ey Ez      =      0 0 Va/ha      (2.7) where Va is the voltage applied to the actuators in the thickness direction. The equilibrium equations undermechanical load may be derived on the basis of the statio- nary potential energy. The equilibrium equations of a plate can be obtained using the principle of minimum total potential energy (Reddy, 1984) Nxx,x+Nxy,y =0 Nxy,x+Nyy,y =0 Mxx,x+Mxy,y−Qxz =0 Mxy,x+Myy,y −Qyz =0 Nxxw,xx+2Nxyw,xy+Nyyw,yy+Qxz,x+Qyz,y =0 (2.8) Equations (2.8) are non-linear equilibrium equations based on the first-order shear deformation plate theory. In Eqs. (2.8), the terms N, Q and M are stress andmoment resultants. Using the constitutive relations, the stress andmoment resultants are defined as Nxx = h 2 +ha ∫ − h 2 −ha σxx dz = A1ε (0) xx +A3ε (0) yy +B1ψx,x+B3ψy,y−2Va(d31Q a 11+d32Q a 21) Nyy = h 2 +ha ∫ − h 2 −ha σyy dz = A3ε (0) xx +A1ε (0) yy +B3ψx,x+B1ψy,y−2Va(d31Q a 12+d32Q a 22) Nxy = h 2 +ha ∫ − h 2 −ha σxy dz =A2γ (0) xy +B2γ (1) xy Qxz = k 2 h 2 +ha ∫ − h 2 −ha σxz dz = Cγ (0) xz Qyz = k 2 h 2 +ha ∫ − h 2 −ha σyz dz = Cγ (0) yz Mxx = h 2 +ha ∫ − h 2 −ha σxxz dz = B1ε (0) xx +B3ε (0) yy +D1ψx,x+D3ψy,y Myy = h 2 +ha ∫ − h 2 −ha σyyz dz = B3ε (0) xx +B1ε (0) yy +D3ψx,x+D1ψy,y Mxy = h 2 +ha ∫ − h 2 −ha σxyz dz = B2γ (0) xy +D2γ (1) xy (2.9) where the constants Ai, Bi, Di, Ti, C in Eqs. (2.9) are Buckling analysis of three-layered rectangular plate ... 817 (A1,B1,D1)= h 2 +ha ∫ − h 2 −ha Q11(1,z,z 2) dz (A2,B2,D2)= h 2 +ha ∫ − h 2 −ha Q44(1,z,z 2) dz (A3,B3,D3)= h 2 +ha ∫ − h 2 −ha Q12(1,z,z 2) dz C = k2 h 2 +ha ∫ − h 2 −ha Q44(1,z,z 2) dz (2.10) In Eq. (2.10)4, k 2 is the shear correction factor. Equations (2.8) are five coupled equilibrium equations, which are non-linear in terms of the displacement components. In order to obtain the stability equations, the adjacent equilibrium criterion is used (Brush and Almroth, 1975), and the stability equations are obtained as N1xx,x+N 1 xy,y =0 N 1 xy,x+N 1 yy,y =0 M1xx,x+M 1 xy,y−Q 1 xz =0 M 1 xy,x+M 1 yy,y −Q 1 yz =0 N0xxw 1 ,xx+2N 0 xyw 1 ,xy+N 0 yyw 1 ,yy+Q 1 xz,x+Q 1 yz,y =0 (2.11) where N0xx, N 0 xy and N 0 yy are the pre-buckling force resultants. 3. Decoupling the stability equations In order to obtain the governing equations, the equivalent form of Eqs. (2.9) in terms of neigh- boring state displacements is substituted into Eqs. (2.11). Therefore, the stability equations in terms of the displacement components are obtained as follows A2(u 1 0,yy +v 1 0,xy)+B2(ψ 1 x,yy +ψ 1 y,xy)+A3v 1 0,xy+A1u 1 0,xx+B3ψ 1 y,xy+B1ψ 1 x,xx =0 A2(u 1 0,xy +v 1 0,xx)+B2(ψ 1 x,xy+ψ 1 y,xx)+A3u 1 0,xy+A1v 1 0,yy +B3ψ 1 x,xy+B1ψ 1 y,yy =0 B1u 1 0,xx+B3v 1 0,xy+D1ψ 1 x,xx+D3ψ 1 y,xy+B2(u 1 0,yy +v 1 0,xy)+D2(ψ 1 x,yy +ψ 1 y,xy) −C(ψ1x+w 1 0,x)= 0 B2(u 1 0,xy +v 1 0,xx)+D2(ψ 1 x,xy+ψ 1 y,xx)+B3u 1 0,xy+B1v 1 0,yy +D3ψ 1 x,xy+D1ψ 1 y,yy −C(ψ1y +w 1 0,y)= 0 N0xxw 1 0,xx+2N 0 xyw 1 0,xy +N 0 yyw 1 0,yy +C(ψ 1 x,x+w 1 0,xx)+C(ψ 1 y,y+w 1 0,yy)= 0 (3.1) Based on Eqs. (2.10), the coefficients A3, B3 and D3 can be rewritten as (A3,B3,D3)= (A1,B1,D1)−2(A2,B2,D2) (3.2) Equations (3.1) are five coupled equations in terms of the neighboringdisplacement components. To decouple governing stability equations (3.1), four new functions are introduced as ϕ1 =u 1 0,x+v 1 0,y ϕ2 = u 1 0,y−v 1 0,x ϕ3 = ψ 1 x,x+ψ 1 y,y ϕ4 = ψ 1 x,y−ψ 1 y,x (3.3) Using the functions introduced in relations (3.3), the stability Eqs. (3.1) can be expressed as A1ϕ1,x+B1ϕ3,x+A2ϕ2,y+B2ϕ4,y =0 A1ϕ1,y+B1ϕ3,y−A2ϕ2,x−B2ϕ4,x =0 (3.4) and B1ϕ1,x+D1ϕ3,x+B2ϕ2,y +D2ϕ4,y −k 2A2(ψ 1 x+w 1 0,x)= 0 B1ϕ1,y +D1ϕ3,y −B2ϕ2,x−D2ϕ4,x−k 2A2(ψ 1 y +w 1 0,y)= 0 (3.5) 818 H. Yaghoobi, I. Rajabi and k2A2(ϕ3+∇ 2w10)+N 0 xxw 1 0,xx+2N 0 xyw 1 0,xy+N 0 yyw 1 0,yy =0 (3.6) where ∇2 is the two-dimensional Laplace operator. It should be pointed out that the coefficients B1 and B2 are exactly equal to zero when the laminated plate is symmetric. By differentiation of Eqs. (3.5) with respect to the variables x and y, respectively, and simplifying the resulting equations, the function ϕ3 is related to the transverse displacement w as follows D1∇ 2ϕ3−k 2A2(ϕ3+∇ 2w10)= 0 (3.7) From Eq. (3.6), ϕ3 can be obtained as ϕ3 =− 1 k2A2 (N0xxw 1 0,xx+2N 0 xyw 1 0,xy+N 0 yyw 1 0,yy)−∇ 2w10 (3.8) Substituting Eq. (3.8) into Eq. (3.7) yields − D1 k2A2 ∇ 2(N0xxw 1 0,xx+2N 0 xyw 1 0,xy +N 0 yyw 1 0,yy)−D1∇ 2 ∇ 2w10 +(N0xxw 1 0,xx+2N 0 xyw 1 0,xy +N 0 yyw 1 0,yy)= 0 (3.9) Following the same procedure as was done to formulate Eq. (3.9), the following equation can be obtained in terms of function ϕ4 as follows D2∇ 2ϕ4−k 2A2ϕ4 =0 (3.10) Equations (3.9) and (3.10) are two decoupled equations in terms of the transverse displace- ment w0 and function ϕ4, respectively. Using Eqs. (3.5), and (3.8), the rotation functions ψx and ψy can be expressed in terms of w0 and ϕ4 as ψ1x = D1 k2A2 [ − k2A2 D1 w10 − 1 k2A2 (N0xxw 1 0,xx+2N 0 xyw 1 0,xy+N 0 yyw 1 0,yy)−∇ 2w10 ] ,x + D2 k2A2 ϕ4,y ψ1y = D1 k2A2 [ − k2A2 D1 w10 − 1 k2A2 (N0xxw 1 0,xx+2N 0 xyw 1 0,xy +N 0 yyw 1 0,yy)−∇ 2w10 ] ,y − D2 k2A2 ϕ4,x (3.11) 4. Boundary conditions It is assumed that two opposite edges of the plate at x =0 and x = a are simply supported (S) and have arbitrary boundary conditions at the other two edges. The arbitrary boundary con- ditions along the other edges, y = 0 and y = b can be clamped-clamped (CC), free-free (FF), simply supported-simply supported (SS), free-clamped (FC), free-simply supported (FS) and clamped-simply supported (CS). Each boundarymay have the following conditions: — Simply Supported v0 = w0 = ψx =0 (4.1) —Clamped w0 = ψx = ψy =0 (4.2) Buckling analysis of three-layered rectangular plate ... 819 —Free Myy = Mxy = Qyz +N 0 xyw 1 0,x+N 0 yyw 1 0,y =0 (4.3) The stress andmoment resultants Qyz, Myy and Mxy can be defined as Myy = D1(ψx,x+ψy,y)−2D2ψx,x Mxy = D2(ψx,y +ψy,x) Qyz = k 2A2(ψy +w0,y) (4.4) which are functions of the rotation functions ψx and ψy and the transverse displacement w0. 5. Buckling analysis To find the critical buckling load, the pre-buckling forces should be found. Thus, using the same procedure developed byDuc andTung (2010), the pre-buckling force resultants are found to be N0xx =− Px b N0xy =0 N0yy =− A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] (5.1) where Px is the uniformly distributed load along the edges x =0,a. Substituting relations (5.1) into Eq. (3.9), yields −D1∇ 2 ∇ 2w10 − D1 k2A2 ∇ 2 [ − Px b w10,xx + ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] ) w10,yy ] + [ − Px b w10,xx+ ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] ) w10,yy ] =0 (5.2) To analyze the buckling behavior, decoupled stability equations (3.10) and (5.2) should be solved. Asmentioned before, the edges of the plate in the x direction are assumed to be simply supported. Using the series solutions in the x direction, the functions w10 and ϕ4 are expressed as w10 = ∞ ∑ m=1 f(y)sin mπx a ϕ4 = ∞ ∑ m=1 g(y)cos mπx a (5.3) where m is the number of half-waves in the x direction. Series solutions (5.3) satisfy the simply supported boundary conditions in the x direction. Substituting relation (5.3) into Eqs. (3.10) and (5.2), yields two ordinary differential equations in terms of the functions f(y) and g(y) as follows D1 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)]+k 2A2 )d4f(y) dy4 + { −2D1 (mπ a )2[ − 1 2 Px b +k2A2+ 1 2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] )] −k2A2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] )}d2f(y) dy2 + [ D1 (mπ a )4( − Px b +k2A2 ) − (mπ a )2 k2A2 Px b ] f(y)= 0 D2 [ −g(y) (mπ a )2 + d2g(y) dy2 ] −k2A2g(y) = 0 (5.4) 820 H. Yaghoobi, I. Rajabi Equation (5.4)1 is an ordinary differential equation with variable coefficients. In order to solve this equation, thepower series solutionmethodofFrobenius (Wylie andBarrett, 1951) isutilized. To this end, the function f(y) is written in the following form f(y)= ∞ ∑ j=0 Cjy j (5.5) where Cj are arbitrary constant coefficients. Substitutingproposed solution (5.5) intoEq. (5.4)1, and shifting the indices, yields D1 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)]+k 2A2 ) · ∞ ∑ j=0 [(j +4)(j +3)(j +2)(j +1)Cj+4y j]+ { −2D1 (mπ a )2[ − 1 2 Px b +k2A2 + 1 2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] )] −k2A2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] )} ∞ ∑ j=0 [(j +2)(j +1)Cj+2y j] + [ D1 (mπ a )4( − Px b +k2A2 ) − (mπ a )2 k2A2 Px b ] ∞ ∑ j=0 (Cjy j)=0 (5.6) Collecting the coefficients of similar powers of j in Eq. (5.6), from the coefficient of y0, it can be obtained C4 =− 1 24D1 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)]+k 2A2 ) · { { −2D1 (mπ a )2[ − 1 2 Px b + 1 2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] ) +k2A2 ] −k2A2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] )} C2 (5.7) + [ D1 (mπ a )4( − Px b +k2A2 ) − (mπ a )2 k2A2 Px b ] C0 } Also, the coefficient of yj gives Cj+4 =− 1 D1 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)]+k 2A2 ) (j+4)(j+3)(j+2)(j+1) · { { −2D1 (mπ a )2( − 1 2 Px b + 1 2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] ) +k2A2 ] −k2A2 ( − A3 A1 Px b + A1−A3 A1 [−2Va(d31Q a 12+d32Q a 22)] )} (j +2)(j +1)Cj+2 (5.8) + [ D1 (mπ a )4( − Px b +k2A2 ) − (mπ a )2 k2A2 Px b ] Cj } Equations (5.7) and (5.8) are recursion relationships, and relation (5.8) is valid for j ­ 0. It should be noted that the coefficients Ci (i = 0,1,2,3) are arbitrary coefficients on account of which the other coefficients Cj (j ­ 4) would be obtained just from recurrence formulas expressed in terms of them selves. Buckling analysis of three-layered rectangular plate ... 821 For solvingEq. (3.10), substitutingproposed series solution (5.3)2 intoEq. (3.10), and solving the resulted ordinary differential equation, yields g(y) = C −1 sinh(λy)+C−2cosh(λy) (5.9) where λ = √ (mπ a )2 + k2A2 D2 (5.10) Imposing the boundary conditions at the edges of the plate in the y direction, a system of six homogeneous algebraic equations is obtained. Setting the determinant of the coefficients equal to zero, the buckling load of the hybrid plate is determined. Needless to say that the lowest value among all of these Px for each m is known as the critical buckling load Pxcr. 6. Results and discussion Numerical results for buckling analysis of a three-layered rectangular plate with piezoelectric layers for different boundary conditions are computed. The material properties are shown in Table 1. This table shows the characteristics of PZT-5A as the piezoelectric and aluminum as the host plate. Moreover, for all numerical results which are reported here, the following values of variables are used unless otherwise indicated by tables or graphs a b =1 b t =100 ha =0.001m h =0.01m Va =500V k 2 = 5 6 Table 1.Material properties of the aluminum and PZT-5A layers Aluminum PZT-5A Elastic modulus E [GPa] 70 63 Poisson’s ratio ν 0.3 0.3 Piezoelectric constant d31 [10 −10m/V] – 2.54 Piezoelectric constant d32 [10 −10m/V] – 2.54 In this Section, firstly, the convergence rate of the power series is checked. Secondly, compa- rison with the previously published related article is employed in order to verify the accuracy of the proposed method. Finally, the critical buckling load of the piezoelectric coupled pla- tes are presented in tabular and graphical forms for different plate aspect ratios, thickness of piezoelectric, actuator voltage and boundary conditions. To guarantee the accuracy of the buckling load obtained by the procedure described above, it is necessary to conduct a convergence study to determine the number of terms required in the power series solution. Since in real calculations a series solution will have to be truncated somewhere according to a pre-determined error bound, an exact solution really implies that the results can be obtained to any desired degree of accuracy. Therefore, the series expansion, Eq. (5.5), will have to be truncated in numerical calculations. Accordingly, to calculate a sufficient number of terms (N), a special case for all kinds of boundary conditions was studied. Table 2 shows the convergence of Pxcr for six different boundary conditions. From this table, it is clearly visible that for the SSSS case, more than 18 terms are needed to obtain the value of Pxcr, accurately to six significant digits. Also, it is seen that if the SCSC is chosen as a boundarycondition, at least 24 terms are required to obtain an extremely accurate value of Pxcr. The bold numbers in the table are those beyond which the sixth digit does not change as N 822 H. Yaghoobi, I. Rajabi Table 2.Convergence test of the critical buckling load, Pxcr [KN], with different combinations of boundary conditions N Boundary conditions SCSC SSSC SSSS SFSC SSSF SFSF 10 423.4500 346.6509 251.0307 118.7528 108.3371 77.6299 12 551.1659 357.9866 253.8557 120.6816 106.5636 81.2094 14 524.5552 355.9793 253.5745 120.9698 106.2536 82.1316 16 530.1927 356.1081 253.5927 120.9983 106.2160 82.2540 18 529.9826 356.0962 253.5918 121.0004 106.2127 82.2646 20 530.0266 356.0965 253.5918 121.0005 106.2125 82.2653 22 530.0258 356.0965 253.5918 121.0005 106.2125 82.2653 24 530.0259 356.0965 253.5918 121.0005 106.2125 82.2653 26 530.0259 356.0965 253.5918 121.0005 106.2125 82.2653 28 530.0259 356.0965 253.5918 121.0005 106.2125 82.2653 30 530.0259 356.0965 253.5918 121.0005 106.2125 82.2653 increases. As more terms are taken, Pxcr converges to its exact value. Therefore, the numerical results from the power series approach which are presented in the calculations were obtained by taking sufficient terms N to converge to the number of digits shown in the tables. In order to verify the accuracy of the present formulations, the buckling load obtained from the present method is compared with those available in the literature. In Table 3, comparison of the non-dimensional critical buckling loads for isotropic plates is made between the results obtained by the present method and those reported by Hosseini-Hashemi et al. (2008). For all boundary conditions, good agreements can be observed, and it is concluded that our formulation is completely trustful. After verifying the merit and accuracy of the present accurate solution, the following new results for the three-layered rectangular plate with piezoelectric actuators can be used as the benchmark for future research studies. Table 3.Comparison of non-dimensional critical buckling loads (Pcr = Pxa 2/D1) for an isotro- pic rectangular plate with different boundary conditions for uniaxial compressive loading in the x direction a b h a Boundary conditions SCSC SSSC SSSS SFSC SSSF SFSF 0.5 0.1 [6] 18.055467 16.247894 14.915722 10.642566 10.408425 9.325573 Pr. 18.05546694 16.24789354 14.91572216 10.64256680 10.40842450 9.32557275 0.2 [6] 15.851026 14.573567 13.580179 9.769217 9.590606 8.632796 Pr. 15.85102641 14.57356705 13.58017889 9.76921715 9.59060641 8.63279575 1 0.1 [6] 63.404106 51.727083 37.447690 15.394207 13.257101 9.112050 Pr. 63.40410565 51.72708294 37.44768979 15.39420741 13.25710052 9.11205022 0.2 [6] 43.567623 41.394657 32.441432 13.627190 12.054586 8.434126 Pr. 43.56762327 41.39465731 32.44143157 13.62719035 12.05458573 8.43412622 2 0.1 [6] 168.416063 151.127293 129.765726 46.387846 24.457914 8.891002 Pr. 168.41606196 151.12729228 129.7657263 46.38784629 24.45791485 8.89100235 0.2 [6] 80.032333 78.916068 76.902078 36.692284 21.430636 8.248819 Pr. 80.03233149 78.91606516 76.90207522 36.69228453 21.43063569 8.24881884 [6] – Hosseini-Hashemi et al. (2008); Pr. – present The variation of the critical buckling load versus the plate aspect ratio for three various voltage actuators are shown in Table 4 and Fig. 2. The primary conclusion tabulated from Buckling analysis of three-layered rectangular plate ... 823 Table 4 is that the critical buckling load diminishes as the plate aspect ratio increases.Moreover, the percentage decrease is about 89% for the SFSFplate and about 15% for the SCSC one from a/b =0.5 to a/b =1.5 under the same actuator voltage Va =−500. It is worthmentioning that increasing the constraints on boundary conditions results in an increase in the critical buckling load, i.e. for a fixed value of variables, the SCSC and SFSF have the highest and lowest Pxcr, respectively. Figure 2 illustrates the effect of a/b for three various voltage actuators, i.e. 500, 0 and −500 on the Pxcr for the SFSC plate. It is apparent from this figure that the Pxcr can be increased by applying a negative voltage on the actuator layers, and the effect of Va becomes more significant at higher plate aspect ratios. Table 4. Effect of the plate aspect ratio on the critical buckling load for different boundary conditions Va a/b Boundary conditions SCSC SSSC SSSS SFSC SSSF SFSF −500 0.5 620.6445 556.4839 511.1295 379.6773 371.7482 339.1097 1.0 564.8202 389.0031 283.1340 135.4563 117.3006 82.9526 1.5 525.5415 408.0511 270.2732 101.6814 72.9200 36.4562 0 0.5 615.4263 551.5939 506.6638 377.8521 370.1353 338.4924 1.0 547.4352 372.5638 268.3629 128.3953 111.8539 82.6416 1.5 514.6062 397.7806 244.4790 86.9342 61.9961 36.2055 +500 0.5 610.2019 546.6993 502.1981 375.9378 368.4899 337.9531 1.0 530.0259 356.0965 253.5918 121.0005 106.2125 82.2653 1.5 503.6551 387.4945 218.6847 71.3912 50.6116 35.8746 Fig. 2. Effect of the plate aspect ratio on the critical buckling load for the SFSC boundary condition Figures 3a,b show the critical buckling load for hybrid laminated plates with different bo- undary conditions subjected to various actuator voltages. The results presented herein reveal that theminus actuator voltages increase the buckling load, whereas the plus actuator voltages decrease the buckling load at the same condition. Very high voltages will be able to influence the buckling response of the hybrid laminated plate. However, such high voltages cannot be applied, because they lead to breakdown in thematerial properties. It can also be seen fromFig. 3b that when the SFSF is chosen as a boundary condition, the effect of voltage actuator on the critical buckling load is very small. In Table 5, the effect of the ratio of the piezoelectric layer thickness to thickness of the host layer on the critical buckling load at different boundary conditions is tabulated. Also, in Fig. 4, this effect for the SSSS boundary condition is depicted. It is seen that with an increase in the piezo-to-host thickness ratio, the Pxcr increases. 824 H. Yaghoobi, I. Rajabi Fig. 3. Effect of the actuator voltage on the critical buckling load for the SCSC, SSSC and SSSS boundary conditions (a) and for the SFSC, SSSF and SFSF boundary conditions (b) Table 5.Effect of the ratio of piezoelectric layer thickness to thickness of the host layer on the critical buckling load for different boundary conditions ha/h Boundary conditions SCSC SSSC SSSS SFSC SSSF SFSF 0 382.3593 256.3732 182.2478 86.8944 76.3741 59.5945 0.1 530.0259 356.0965 253.5918 121.0005 106.2125 82.2653 0.2 708.1540 476.5520 339.8699 162.2539 142.2768 109.5356 0.3 916.0362 617.2578 440.7347 210.4871 184.4213 141.2979 0.4 1153.2794 777.9461 555.9935 265.6071 232.5654 177.4930 Fig. 4. Effect of the ratio of the piezoelectric layer thickness to thickness of the host layer on the critical buckling load for the SSSS boundary condition 7. Conclusion In this article, mechanical buckling analysis has been presented for a three-layered rectangular plate with piezoelectric actuators subjected to the combined action of mechanical and electric loads. The derivations were based on the first-order plate theory and by employing an analy- tical approach, the five coupled governing stability equations are converted into two decoupled partial differential equations. By using the Levy solution, these equations are converted into two independent ordinary differential equations, and the power series method of Frobenius is used for solving these equations accurately. Extensive parametric studies for this structureunder different sets of electric loading and boundary conditions have been carried out. The following conclusions, from the numerical computations were drawn. Buckling analysis of three-layered rectangular plate ... 825 • The buckling load is decreased by increasing the plate aspect ratio in both negative and positive actuator voltages and all boundary conditions. • The application of negative voltage on the actuator layers can improve the mechanical buckling strength, and the critical buckling load can be controlled by applying a suitable voltage on the actuator layers. • The critical buckling load increases with the increase of ratio of the piezoelectric layer thickness to the thickness of the host layer. • Increasing the constraints on boundary conditions results in an increase in the critical buckling load. References 1. Akhras G., Li W.C., 2008, Three-dimensional thermal buckling analysis of piezoelectric compo- site plates using the finite layer method, Smart Materials and Structures, 17, 1-8 2. BrushD.O.,AlmrothB.O., 1975,Buckling of Bars, Plates, and Shells,McGraw-Hill,NewYork 3. Correia V.M.F., Soares C.M.M., Soares C.A.M., 2003, Buckling optimization of composite laminated adaptive structures,Composite Structures, 62, 315-321 4. DucN.D.,TingH.V., 2010,Mechanical and thermal postbuckling of shear-deformable fgmplates with temperature-dependent properties,Mechanics of Composite Materials, 46, 461-476 5. Dumir P.C., Kumari P., Kapuria S., 2009, Assessment of third order smeared and zigzag the- ories for buckling and vibration of flat angle-ply hybrid piezoelectric panels,Composite Structures, 90, 346-362 6. Hosseini-Hashemi Sh., Khorshidi K., Amabili M., 2008, Exact solution for linear buckling of rectangularMindlin plates, Journal of Sound and Vibration, 315, 318-342 7. Kapuria S., Achary G.G.S., 2004, Exact 3-D piezoelasticity solution for buckling of hybrid cross-ply plates using transfer matrices,Acta Mechanica, 170, 25-45 8. KapuriaS.,AcharyG.G.S.,2006,Nonlinear zigzag theory for electrothermomechanicalbuckling of piezoelectric composite and sandwich plates,Acta Mechanica, 184, 61-76 9. Liew K.M., Yang J., Kitipornchai S., 2003, Postbuckling of piezoelectric FGM plates sub- ject to thermo-electro mechanical loading, International Journal of Solids and Structures, 40, 15, 3869-3892 10. Oh K., Han J.H., Lee I., 2000, Postbuckling and vibration characteristics of piezolaminated composite plate subject to thermopiezoelectric loads, Journal of Sound and Vibration, 233, 1, 19-40 11. Reddy J.N., 1984,Energy Principles and Variational Methods in AppliedMechanics, JohnWiley, NewYork 12. Reddy J.N., 2004, Mechanics of Laminated Composite Plates and Shels: Theory and Analysis, CRCPress LLC, USA 13. Shariyat M., 2009, Dynamic buckling of imperfect laminated plates with piezoelectric sen- sors and actuators subjected to thermo-electro-mechanical loadings, considering the temperature- dependency of the material properties,Composite Structures, 88, 228-239 14. ShenH.-S., 2001a,Postbuckling of shear deformable laminated plateswith piezoelectric actuators under complex loading conditions, International Journal of Solids and Structures, 38, 7703-7721 15. Shen H.-S., 2001b, Thermal postbuckling of shear deformable laminated plates with piezoelectric actuators,Composites Science and Technology, 61, 1931-1943 826 H. Yaghoobi, I. Rajabi 16. Shen H.-S., Zhu Z.H., 2011, Compressive and thermal postbuckling behaviors of laminated pla- tes with piezoelectric fiber reinforced composite actuators, Applied Mathematical Modelling, 35, 1829-1845 17. Varelis D., Saravanos D.A., 2004, Coupled buckling and postbuckling analysis of active lami- nated piezoelectric composite plates, International Journal of Solids and Structures,41, 1519-1538 18. Wylie C.R., Barrett L.C., 1951,Advanced EngineeringMathematics, NewYork,McGraw-Hill Manuscript received September 14, 2012; accepted for print December 7, 2012