An improved C0 FE model for the analysis of laminated sandwich plate with soft core

Finite Elements in Analysis and Design 56 (2012) 20–31

Contents lists available at SciVerse ScienceDirect

Finite Elements in Analysis and Design journal homepage: www.elsevier.com/locate/finel

An improved C0 FE model for the analysis of laminated sandwich plate with soft core H.D. Chalak a,n, Anupam Chakrabarti a, Mohd. Ashraf Iqbal a, Abdul Hamid Sheikh b a b

Department of Civil Engineering, Indian Institute of Technology Roorkee, Roorkee 247667, India School of Civil, Environment and Mining Engineering, University of Adelaide, North Terrace, Adelaide, SA 5005, Australia

a r t i c l e i n f o

a b s t r a c t

Article history: Received 16 September 2011 Received in revised form 13 February 2012 Accepted 18 February 2012 Available online 22 March 2012

An improved C0 two dimensional finite element model based on higher order zigzag plate theory (HOZT) is developed and applied to the analysis of laminated composite and sandwich plates under different situations to study the performance of the model. In the proposed model, the in-plane displacements variation is considered to be cubic for both the face sheets and the core, while the transverse displacement is assumed to vary quadratically within the core and remains constant in the faces beyond the core. It satisfies the conditions of transverse shear stress continuity at the layer interfaces as well as satisfies the zero transverse shear stress condition at the top and bottom of the plate. The well-known problem of continuity requirement of the derivatives of transverse displacements is overcome by choosing the nodal field variables in an efficient manner. A nine-node C0 quadratic plate finite element is implemented to model the HOZT for the present analysis. Numerical examples covering different features of laminated composite and sandwich plates are presented to illustrate the accuracy of the present model. & 2012 Elsevier B.V. All rights reserved.

Keywords: Sandwich plate Soft core Zigzag theory Finite element Stress continuity

1. Introduction A laminated composite/sandwich structure has a layered construction, which consists of a number of lamina or ply of orthotropic materials stacked one over the other and bonded together to act as an integral structural element. The individual layers of the laminate may have different orientations which enables the structural designer to achieve required strength in the preferred direction. It also shows superior properties such as high strength/stiffness to weight ratio and greater resistance to environmental degradation compared to conventional metallic materials. Due to all these merits, the fiber reinforced laminated composite/sandwich is gaining wide acceptance in various structural applications. In order to fulfill the requirement of weight minimization in a more efficient manner, a sandwich construction having low strength core and high strength face sheets is used. The role of transverse shear deformation is very important in laminated composites, as the material is weak in shear due to its low shear modulus to extensional rigidity. Due to the large variation of material properties across the thickness, the behavior of laminated sandwich plates becomes very complex. A proper understating on the response of these structures under different

n

Corresponding author. Tel.: þ91 1332 285844; fax: þ91 1332 275568. E-mail address: [email protected] (H.D. Chalak).

0168-874X/$ - see front matter & 2012 Elsevier B.V. All rights reserved. doi:10.1016/j.finel.2012.02.005

loading conditions is extremely important for their safe design. In this context a number of plate theories have been developed where the major emphasis is to model the shear deformation in refined manner. These plate theories can be broadly divided into two categories based on their assumed displacement fields: (1) Single layer theory and (2) Layer-wise theory. In single layer theory, the deformation of plate is expressed in terms of unknown parameters of the reference plane, i.e. middle plane. In this theory the transverse shear strain is assumed to be uniform over the entire plate thickness and it is known as Reissner–Mindlin’s plate theory which is also known as the first order shear deformation theory (FSDT). Goyal and Kapania [1] developed a five node beam FE model based on FSDT. Moderately thick rectangular laminated composite plate was analyzed by Ferreira [2] using multiquadric radial basis function method (i.e. mesh free collocation method) based on FSDT. However, this theory (FSDT) requires a shear correction factor to compensate for the actual parabolic variation of the shear stress. The higher order shear deformation theories (HSDT) have been developed with the aim to avoid the use of shear correction factors by including the actual cross sectional warping and to get the realistic variation of the transverse shear strains and stresses throughout the plate thickness [3]. Kant [4] derived the complete set of equations for the analysis of thick elastic plates with the help of third order refined shear deformation theory (HSDT). The three-dimensional Hooke’s laws was also used for plate material

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

in the theory [4] which gave a more realistic quadratic variation of the transverse shearing strains and linear variation of the transverse normal strain through the plate thickness. Later, Kant and Swaminathan [5] have reported analytical solutions based on the higher order shear deformation theory (HSDT). Using Reddy’s displacement field [3] for third order shear deformation theory, a set of dynamic equations was derived for modeling of laminated structures by Aagaah et al. [6]. The theories proposed by Murthy et al. [7] and Subramanian [8] for the analysis of laminated beams, are also based on higher order shear deformation theory (HSDT). Pervez et al. [9] presented a two dimensional serendipity model based on HSDT for static analysis of laminated composite plate. Based on global local higher order shear deformation theory, Wu et al. [10] presented 4-node quadrilateral element and 3-node triangular element (assuming 13 field variables per node) satisfying weak continuity conditions for the static analysis of plates. Recently, Wu et al. [11] described the short coming of zigzag theories, i.e. requirement of C1 continuity condition and presented a six node triangular C0 FE by taking out the first derivatives of transverse displacements from the in-plane displacement fields for the static analysis of laminated sandwich plate. Kulkarni and Kapuria [12] proposed a new discrete Kirchhoff 4-node quadrilateral element (having 7 field variables per node) based on Reddy’s HSDT [3]. The C1 continuity requirement [12] is overcame by defining the derivatives of transverse displacement in terms of separate field variables. Aydogdu [13] presented HSDT for the static, vibration and buckling analysis of laminated composite plates where the shear deformation function was chosen according to 3-D analysis by using inverse method. A nine-node rectangular element with nine field variables at each node was developed by Tu et al. [14] based on HSDT for the bending and vibration analysis of laminated composite and sandwich plates. Ferreira et al. [15] presented the radial basis function collocation method for the static and vibration analysis of thick plates using FSDT and HSDT given by Kant [4]. For the analysis of thin and thick composite plates, Roque et al. [16] used higher order shear deformation theory. Due to different values of shear rigidity at the adjacent layers, HSDT shows discontinuity in the shear stress distribution at the layer interfaces with continuous variation of the transverse shear strain across the thickness. But the actual behavior of a composite laminate is opposite, i.e., the transverse shear stress must be continuous at the layer interface and the corresponding strain may be discontinuous [17]. In order to overcome the above disparity, the layer-wise theories are developed. The layer-wise theories may be further classified into discrete layer plate theory and refined plate theory. In discrete layer plate theory, unknown displacement components are taken at all the layer interfaces. Discrete layer theories proposed by Robbins and Reddy [18], Toledano and Murakami [19], Lu and Liu [20], Reddy [21] and many others assume unique displacement field in each layer and displacement continuity across the layers. Tahani [22] presented the analytical solution for laminated beams by using two theories based on layer-wise displacement fields. Ramesh et al. [23] presented a 45-node triangular element with 7 and 3 field variables at each node, based on the HSDT and layer-wise plate theory of Reddy respectively for the static analysis of the laminated composite plate. The performance of this plate theory is good but it required huge computational involvement as the number of unknowns increases directly with the increase in the number of layers. To solve the above problem, the unknowns at different interfaces are defined in terms of those at the reference plane in refined plate theories (also known as zigzag theories). In this theory, the in plane displacements have piecewise variation across the plate thickness and the number of unknowns are made independent of the number of layers by equating the transverse

21

shear stresses at the layer interfaces of the laminate. In some improved version of these theories, the condition of zero transverse shear stresses at the plate/beam top and bottom was also satisfied. The theories developed by Murakami [24], Di Sciuva [25], Lee et al. [26], Cho and Parmerter [27], Cho and Averill [28] and many other fall under this category. Carrera [29] presented a historical review on the zigzag theories used for the analysis of multilayered laminate plates and shells, in which the three basic theories have been discussed, namely: Lekhnitskii Multilayered Theory (LMT), Ambartsumian Multilayered Theory (AMT) and Reissner Multilayered Theory (RMT). A triangular element was presented by Chakrabarti and Sheikh [30] based on zigzag theory, which shows excellent performance though the element is does not satisfy the normal slope continuity requirement. Akhras and Li [31] developed a spline finite strip method based on higher order zigzag theory for the static analysis of the plate. Kapuria and Kulkarni [32] presented a four node quadrilateral element based on third order zigzag theory for the analysis of the laminates. The C1 continuity requirement is circumvented by using discrete Kirchoff constraint approach, where the derivaties of transverse displacement are replaced by rotational variables. Fares and Elmarghany [33] presented a first order zigzag theory of composite plates using Ressiner’s mixed variational formula. Recently, Ferreira et al. [34] and Rodrigues et al. [35] presented radial basis functions, finite differences collocation and unified formulation for the analysis of laminated plates based on zigzag theory. Zigzag models for laminated composite structures were developed by using trigonometric terms to represent the linear displacement field, transverse shear strains and stresses [36–38]. These theories (Zigzag) provide a very accurate approximation of the structural behavior even for lower span to thickness ratio. However, the zigzag theory has a problem in its finite element implementation as it requires C1 continuity of the transverse displacement at the nodes. To combine the benefits of the discrete layer-wise and higher order zigzag theories, Icardi [39,40], Yip and Averill [41] and many others developed theories which are known as sub-laminated models. Cho and Averill [42] presented an improved sublaminate model with first order zigzag approximation of displacement within each sub-laminate, which contains an eight node C0 FE having five displacement field variables at each node for each sub-laminate. Averill [43] developed a C0 finite element based on first order zigzag theory and overcome the C1 continuity requirement by incorporating the concepts of independent interpolations and penalty functions. Hermitian functions were used by Di Sciuva [44,45] to approximate the transverse displacement in his formulations. Carrera [46] used two different fields along the laminate thickness direction for transverse displacement and transverse shear stress respectively for formulation. Averill and Yip [47] developed a C0 finite element based on cubic zigzag theory, using interdependent interpolations for transverse displacement and rotations and penalty function concepts. Aitharaju and Averill [48] developed a new C0 FE based on a quadratic zigzag layer-wise theory. For eliminating shear locking phenomenon, the shear strain field is also made field consistent. The transverse normal stress was assumed to be constant through the thickness of the laminate. The new FE was applied to model the beam as combination of different sub-laminates. A C0 plate model based on enhanced first order theory (EFSDT) was presented by Kim and Cho [49,50], where it was shown that the displacements, in-plain strains and stresses can be approximated to those of the three dimensional theory or higher order theory, in the least square sense. Recently the authors [51], have developed a C0 model using the EFSDT based on mixed variational theorem, which also satisfy the lateral conditions at the top and bottom surfaces of the plate. The mixed FE approach was

22

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

presented by Rao et al. [52], Ramtekkar et al. [53,54], Bambole and Desai [55] and many others, where the stress components were assumed as field variables at interface nodes along with displacement field variables for the accurate analysis of the stresses. Laminated soft core sandwich structure is extensively used mainly in weight minimization applications. The transverse deformation is very significant in such sandwich structure as there is abrupt change in the values of transverse shear rigidity and thickness of face sheet and the core. As such more attention must be given for the accurate modeling of the variation of transverse displacement across the depth of a sandwich structure having soft core. Frostig [56] has presented the classical and the higher order computational models of unidirectional sandwich panels with incompressible and compressible cores to demonstrate the differences in overall response of the panels as well as in the vicinity of the localized loads and supports. Givil et al. [57] have presented the dynamic model based on higher order sandwich panel theory to study the behavior of soft core sandwich panel under dynamic loading. It is well known that the sandwich structures having low strength honeycomb or foam cores are subjected to complex stress distributions under different kinds of transverse loadings due to their typical construction and material characteristics. In addition to the effects of the transverse shear deformation, the bending response of the sandwich laminate is also greatly dependent on the transverse normal deformation of the flexible core. The situation becomes critical when the plate is thick and the difference in the properties between the core and the face sheets increases. So it is essential to introduce unknowns in the transverse displacement field across the depth in addition to that in the reference plane to represent the variation of transverse displacement in a soft core sandwich laminate. This can be done by using sub-laminate plate theories but the number of unknowns will increase with the increase in the number of sublaminates. On the other hand, introduction of additional unknowns in the transverse displacement field invites additional C1 continuity requirements in its finite element implementation using the zigzag theory as mentioned earlier. This requirement of C1 continuity in the HOZT imposes restriction in the finite element implementation. There are very few FE elements based on HOZT which can overcome the above C1 continuity requirement. Moreover, most of these elements are not attractive for practical application due to the presence of higher order derivative in the field variables and some other disadvantages. Pandit et al. [58–61] proposed a higher order zigzag theory for the analysis of sandwich plates with soft compressible core. To overcome the above problem of C1 continuity they have used separate shape functions to define the derivatives of transverse displacements and developed a C0 finite element model for the implementation of the theory. However, it has imposed some constrains, which are enforced variationally through penalty approach. The selection of suitable value for the penalty stiffness multiplier is quite arbitrary and is a well known problem in the finite element method. Recently, an improved 1-D C0 beam finite element model was also developed by the authors [62] based on higher order zigzag theory for the static analysis of laminated beams. In view of the above discussion, a new C0 two dimensional plate finite element model has been proposed in this paper for the analysis of laminated sandwich plates having soft compressible core using higher order zigzag theory (HOZT). The major contribution of the proposed FE model is threefold: (1) it overcomes the problem of the C1 continuity associated with the HOZT to implementing a C0 formulation, (2) it includes the effect of core compressibility in the formulation and (3) it eliminates the requirement of using penalty multiplier in the stiffness matrix

formation. In this model, the in-plane displacement fields are assumed as a combination of a linear zigzag function with different slopes at each layer and a cubically varying function over the entire thickness. The transverse displacement is considered to be quadratic within the core and constant in the face sheets. The model satisfies the transverse shear stress continuity conditions at the layer interfaces and the conditions of zero transverse shear stress at the top and bottom of the plate. The isoparametric quadratic plate element used has nine nodes with 11 field variables (i.e., in-plane displacements and transverse displacement at the reference midsurface, at the top and at the bottom of the plate along with rotational field variable only at the reference midsurface) at each node. In this formulation the displacement fields are chosen in such a manner that there is no need to define any new field variable [32] and no need to impose any penalty stiffness as observed in the works by Pandit et al. [58–61]. The element may also be matched quite conveniently with other C0 elements. The present FE model is primarily developed to address the problem associated with the modeling of laminated sandwich plates having soft core. However, the proposed model is also capable to predict relatively simpler behavior of ordinary laminated composite plates quite efficiently. Therefore, the present plate FE model is used to solve many problems on laminated composites and sandwich structures having soft compressible core for different loadings, geometry boundary conditions and others.

2. Mathematical formulations To ensure a piecewise parabolic variation of transverse shear strains across the thickness with discontinuity at the layer interface as expected in a layered plate [63], the in-plane displacement fields [27] (Fig. 1) are chosen as follows: U ¼ u0 þ zyx þ

nX u 1

ðzzui ÞHðzzui Þaixu þ

i¼1

n l 1 X

ðzzlj ÞHðz þ zlj Þajxl

j¼1

þ bx z2 þ Zx z3 V ¼ v0 þzyy þ

ð1Þ

nX u 1

n l 1 X

i¼1

j¼1

ðzzui ÞHðzzui Þaiyu þ

ðzzlj ÞHðz þzlj Þajyl

þ by z2 þ Zy z3

ð2Þ

where u0 and v0 denote the in-plane displacements of any point on the midsurface, yx and yy are the rotations of the normal to the middle plane about the y- and x-axis respectively, nu and nl are number of upper and lower layers respectively, bx, by, Zx and Zy

Fig. 1. General lamination scheme and displacement configuration.

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

Fig. 2. Variation of transverse displacement (w) through the thickness of laminated sandwich plate.

are the higher order unknown, aixu , aiyu , ajxl and ajyl are the slopes of ith/jth layer corresponding to upper and lower layers respectively and Hðzzui Þ and Hðz þ zlj Þ are the unit step functions. The transverse displacement is assumed to vary quadratically through the core thickness and constant over the face sheets (as shown in Fig. 2) and it may be expressed as W ¼ l1 wu þl2 w0 þ l3 wl for core ¼ wu for upper face layers ¼ wl for lower face layers

ð3Þ

where wu, w0 and wl are the values of the transverse displacement at the top layer, middle layer and bottom layer of the core, respectively, and l1, l2 and l3 are Lagrangian interpolation functions in the thickness co-ordinate as defined in Appendix A. The stress–strain relationship of an orthotropic layer/lamina (say kth layer) having any fiber orientation with respect to structural axes system (x–y–z) may be expressed as 9 2 3 8 8 exx > s 9 Q Q 12 Q 13 Q 14 0 0 > > > > > > 6 11 > > xx > 7 > > > > > 6Q > eyy > > > > s Q 22 Q 23 Q 24 0 0 7 > yy > 21 > > > > > > 6 > > 7 > > > > = < < szz = 6 Q 7 e zz 0 0 7 6 31 Q 32 Q 33 Q 34 ¼6 7 g t > > > > xy 6 7 xy > Q Q 42 Q 43 Q 44 0 0 > > > > 6 41 > > 7 > > > > > > > > > 6 7 g t > > > Q 55 Q 56 5 > 0 0 0 > > 4 0 > xz > > xz > > > > ; > : t ; : g yz 0 0 0 0 Q 65 Q 66 yz K

or fsg ¼ ½Q K feg

ð4Þ

where fsg,feg and ½Q K  are the stress vector, the strain vector and the transformed rigidity matrix of kth lamina, respectively. Utilizing the conditions of zero transverse shear stress at the top and bottom surfaces of the plate and imposing the conditions of the transverse shear stress continuity at the interfaces between the layers along with the conditions, u ¼uu and v¼vu at the top and u ¼ul and v ¼vl at the bottom of the plate, bx, Zx, by, Zy, aixu , aixl ,aiyu , aiyl , (qwu/qx),(qwl/qx)(qwu/qy) and (qwl/qy) may be expressed in terms of the displacements u0, v0, yx,, yy, uu, ul, vu and vl as fBg ¼ ½Afag

ð5Þ

where

23

coordinates, unit step functions and material properties as defined in Appendix B. By imposing four additional conditions obtained by satisfying the in-plane displacements at the top and bottom of the plate, four first order derivatives terms of transverse displacements are replaced in terms of nodal field variables (C0) in Eq. (5). Thus the in-plane displacement fields expressed in Eqs. (6) and (7) do not contain any first order derivatives of the transverse displacement and therefore the requirement of Cl continuity of HOZT has been avoided very efficiently without defining new field variables [32] and without using any penalty method [58–61]. The generalized displacement vector {d} for the present plate model can now be written with the help of Eqs. (3), (6) and (7) as fdg ¼ fu0 v0 w0 yx yy uu vu wu ul vl wl g

T

Using linear strain–displacement relation and Eqs. (1)–(5), the strain field may be expressed in terms of unknowns (for the structural deformation) as   @U @V @W @U @V @U @W @V @W þ þ þ or feg ¼ ½Hfeg feg ¼ ð8Þ @x @y @z @x @y @z @x @z @x where feg ¼ ½u0 v0 w0 yx yy uu vu wu ul vl wl ð@u0 =@xÞ ð@u0 =@yÞ ð@v0 =@xÞ ð@v0 =@yÞ ð@w0 =@xÞ ð@w0 =@yÞ ð@yx =@xÞ ð@yx =@yÞ ð@yy =@xÞ ð@yy =@yÞ ð@uu =@xÞ ð@uu =@yÞ ð@vu =@xÞ ð@vu =@yÞ ð@wu =@xÞ ð@wu =@yÞ ð@ul =@xÞ ð@ul =@yÞ ð@vl =@xÞ ð@vl =@yÞ ð@wl =@xÞ ð@wl =@yÞ and the elements of [H] are functions of z and unit step functions, as given in Appendix C. With the quantities found in the above equations, the total potential energy of the system under the action of transverse load may be expressed as

Pe ¼ U s 2W ext

ð9Þ

where Us is the strain energy and Wext is the energy due to the external transverse static load. Using Eqs. (3) and (6), the strain energy (Us) is given by ZZ n ZZZ 1X 1 fegT ½Dfegdx dy ð10Þ fegT ½Q k fegdx dy dz ¼ Us ¼ 2k¼1 2 where ½D ¼

n Z X

½HT ½Q k ½H dz

ð11Þ

k¼1

and the energy due to externally applied distributed transverse static load of intensity q(x,y) can be calculated as ZZ ð12Þ W ext: ¼ wq dxdy

fBg ¼ fbx Zx by Zy a1xu a2xu . . .anu1 a1xl a2xl . . .anl1 a1yu a2yu . . .anu1 a1yl a2yl . . .anl1 xu yu xl yl ð@wu =@xÞð@wu =@yÞð@wl =@xÞð@wl =@yÞgT , T

fag ¼ fu0 v0 yx yy uu vu ul vl g

and the elements of [A] are dependent on material properties. It is to be noted that last four entries of the vector {B} helps to define the derivatives of transverse displacement at the top and bottom faces of the plate in terms of the displacements u0, v0, yx, yy, uu, vu, ul and vl to overcome the problem of C1 continuity as mentioned before. Using the above equations, the in-plane displacement fields as given in Eqs. (1) and (2) may be expressed as U ¼ b1 u0 þ b2 v0 þb3 yx þb4 yy þb5 uu þb6 vu þ b7 ul þ b8 vl

ð6Þ

V ¼ c1 u0 þc2 v0 þ c3 yx þ c4 yy þc5 uu þc6 vu þc7 ul þ c8 vl

ð7Þ

where the coefficients bi’s and ci’s are function of thickness

In the present problem, a nine-node quadratic element with 11 field variables (u0, v0, w0, yx, yy, uu, vu, wu, ul, vl and wl) per node is employed. Using finite element method the generalized displacement vector {d} at any point may be expressed as fdg ¼

n X

N i fdgi

ð13Þ

i¼1

where {d}¼{u0 v0 w0 yx yy uu vu wu ul vl wl}T as defined earlier, {d}i is the displacement vector corresponding to node i, Ni is the shape function associated with the node i and n is the number of nodes per element, which is nine in the present study. With the help of Eq. (13), the strain vector {e} that appeared in Eq. (8) may be expressed in terms of unknowns (for the structural deformation) as feg ¼ ½Bfdg

ð14Þ

24

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

where [B] is the strain–displacement matrix in the Cartesian coordinate system. The elemental potential energy as given in Eq. (9) may be rewritten with the help of Eqs. (10)–(14) as ZZ ZZ   1 1 T T fdg ½BT ½D½B d dxdy fdg ½BT ½Nw T q dx dy &e ¼ 2 2   1 1 T T ð15Þ ¼ fdg ½K e  d  fdg fP e g 2 2 where ZZ ½K e  ¼ ½BT ½D½Bdx dy ZZ fPe g ¼

boundary, the field variables w0, wu and wl are restrained while u0, yx, uu and ul are unrestrained. 2. Clamped boundary condition: All the nodal field variables at the boundary are fully restrained. 3. Free boundary condition: All the nodal field variables at the boundary are unrestrained. The following non-dimensional quantities are used to show different results in this paper: Non-dimensional in-plane stresses, 2 ðsx , sy , sxy Þ ¼ ðh =q0 a2 Þðsx , sx , sxy Þ Non-dimensional transverse shear stresses, ðtxz , tyz Þ ¼ ðh=q0 aÞðtxz , tyz Þ Non-dimensional transverse displacement, 3 w ¼ ð100ETf h w=a4 q0 Þ,

ð16Þ

½Nw T q dxdy

ð17Þ

where [Nw] is the shape function like matrix with non-zero terms associated only with the corresponding transverse nodal displacements. The equilibrium equation can be obtained by minimizing &e as given in Eq. (15) with respect to {d} as ½K e fdg ¼ fP e g

ð18Þ

where [Ke] is the element stiffness matrix and {Pe} is the nodal load vector. The global stiffness matrix and global load vector for the whole plate is then formed by taking the contribution of all the plate elements. Finally, the global linear simultaneous equations are formed and solved for the problem of the sandwich plate after incorporation of appropriate boundary conditions. The stresses are calculated with the help of constitutive relationship by using the condition of stress continuity as in Eq. (5). A numerical code is developed to implement the above mentioned operations involved in the proposed FE model for calculating deflections and stresses in composites and sandwich laminates. The skyline technique has been used to store the global stiffness matrix in a single array and Gaussian decomposition scheme is adopted for the solution.

3. Numerical results To check the accuracy and applicability of the present C0 2Dplate FE model, many problems of laminated composites/sandwich plates are solved under static loading. The results obtained are compared with different results reported in the literatures and 3D elasticity results [63]. The following different boundary conditions are used: 1. Simply supported boundary condition: (Parallel to x-axis) The field variables u0, v0, w0, yx, uu, wu, ul, wl are restrained while yy, vu and vl are unrestrained in one boundary. In the other boundary, the field variables w0, wu and wl are restrained while v0, yy, vu and vl are unrestrained. (Parallel to y-axis) The field variables u0, v0, w0, yy, vu, wu, vl, wl are restrained while yx, uu and ul are unrestrained in one boundary. In the other

where a, h are the dimension of plate along x- direction and z-direction respectively and ETf is the transverse modulus of elasticity of face layer. 3.1. Simply supported laminated sandwich plate (0/90/0) under sinusoidal loading The problem of a three ply (0/90/0) square laminate subjected to sinusoidal loading in both the directions is considered for the analysis. The layers are of equal thickness. The material properties used in this problem and all subsequent problems of laminated composite plates are as shown in Table 1. Considering different thickness ratio (h/a) ranging from 0.01 to 0.25, the results for the non-dimensional displacements (transverse) and stresses (transverse shear and in-plane normal) are presented in Table 2 mainly to study the rate of convergence and validation of the present results. The full plate is analyzed with different mesh sizes. It may be observed in Table 2 that the displacements are converged at mesh size (6  6). However, more elements are required for the convergence of the stresses as expected in a displacement based formulation. As such a mesh size of (12  12) is taken for all subsequent analysis to get sufficiently accurate results for displacements as well as stresses. For the comparison of the present results a computer code is also developed to generate results based on 3-D elasticity solution [63]. The present FE results are also compared with the results obtained by Ramesh et al. [23] using 45-node triangular elements with three field variables and seven field variables at each node based on LT and TSDT respectively. The total number of unknowns per element is 99 in the present FE as compared to 315 (TSDT) and 135 (LT) in Ramesh et al. [23] case. The present FE model gives better results than the results reported by Ramesh et al. [23] with less number of unknowns per element. The present 2D results are also compared with the results given by Wu et al. [10] obtained from two finite elements (4-node plate and 3-node triangular with 13 field variables per node) developed by them satisfying weak-continuity conditions based on HSDT.

Table 1 Material properties for laminated plates. Example

Layer/sheet

Material properties E1 (psi)

E2 (psi)

E3 (psi)

G12 (psi)

G13 (psi)

G23 (psi)

u12 ¼u12 ¼ u23

Composite plates

All layers

25.0E06

1.0E06

1.0E06

0.5E06

0.5E06

0.2E06

0.25

Sandwich plates

Face Core

25.0E06 4.0E06

1.0E06 5.0E04

1.0E06 5.0E04

0.5E06 16.0E04

0.5E06 6.0E04

0.2E06 6.0E04

0.25 0.25

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

25

Table 2 Convergence and comparison of non-dimensional deflection and stresses of a simply supported square laminate (0/90/0) under sinusoidal load. h/a

Ref.

w

sx

txz

txy

0.01

Present (4  4)a Present (8  8) Present (12  12) Pagano [63] % errorb Ramesh et al.c [23] Ramesh et al.d [23] Aydogdu [13] Kulkarni and Kapuria [12] Wu et al. [10] Chakrabarti and Sheikh [30] Sheikh and Chakrabarti [17] Kant and Swaminathan [5]

0.4300 0.4298 0.4298 0.4347 (  1.12) 0.4349 0.4345 0.435 0.4349 0.4350 0.4358 0.4350 0.4343

0.5627 0.5439 0.5401 0.5392 (0.17) 0.5395 0.5394 0.5389 0.5403 0.5238 – 0.5496 0.5392

0.4192 0.4040 0.3874 0.3947 (  1.85) 0.3921 0.3952 0.3003 0.2592 0.3853 – 0.2401 –

0.0221 0.0214 0.0212 0.0213 ( 0.469) 0.0214 0.0214 0.0214 0.0214 0.0213 – 0.0215 0.0214

0.02

Present Pagano [63] % error Aydogdu [13] Wu et al. [10] Chakrabarti and Sheikh [30] Kant and Swaminathan [5]

0.4401 0.4451 (  1.12) 0.4448 0.4461 0.4462 0.4432

0.5418 0.5410 (0.15) 0.540 0.5267 – 0.5406

0.3864 0.3934 (  1.78) 0.299 0.3887 – –

0.0215 0.0216 ( 0.46) 0.0216 0.0213 – 0.0216

0.05

Present Pagano [63] % error Ramesh et al.c [23] Ramesh et al.d [23] Aydogdu [13] Wu et al. [10] Sheikh and Chakrabarti [30] Kant and Swaminathan [5] Liou and Sun [64]

0.5113 0.5164 (  0.99) 0.5166 0.5060 0.511 0.520 0.5066 0.5053 0.5170

0.5536 0.5525 (0.20) 0.5528 0.5509 0.548 0.5457 – 0.5504 0.553

0.3777 0.3846 (  1.79) 0.3866 0.3871 0.295 0.3826 – – 0.395

0.0232 0.0234 ( 0.86) 0.0231 0.0231 0.0232 0.0234 – 0.0231 0.0233

0.10

Present Pagano [63] % error Ramesh et al.c [23] Ramesh et al.d [23] Aydogdu [13] Kulkarni and Kapuria [22] Wu et al. [10] Chakrabarti and Sheikh [30] Sheikh and Chakrabarti [17] Kant and Swaminathan [5] Liou and Sun [64] Reddy [3]

0.7480 0.7530 (  0.66) 0.7535 0.7178 0.7336 0.7136 0.7637 0.7522 0.7140 0.7151 0.7546 0.713

0.5918 0.5908 (0.17) 0.5910 0.5850 0.578 0.5696 0.5864 – 0.5806 0.5836 0.590 0.5684

0.3513 0.3573 (  1.68) 0.3576 0.3671 0.282 0.2453 0.3561 – 0.2437 – 0.357 –

0.0287 0.0290 ( 1.03) 0.0289 0.0281 0.0284 0.0277 0.0290 – 0.0279 0.0279 0.0289 –

0.25

Present Pagano [63] % error Ramesh et al.c [23] Ramesh et al.d [23] Aydogdu [13] Kulkarni and Kapuria [12] Wu et al. [10] Chakrabarti and Sheikh [30] Sheikh and Chakrabarti [17] Kant and Swaminathan [5] Liou and Sun [64] Reddy [3]

2.0151 2.0059 (0.46) 1.9927 1.9136 1.9856 1.9248 2.0724 1.9502 1.9230 1.8948 2.020 1.922

0.7633 0.7548 (1.13) 0.8014 0.7672 0.755 0.7357 0.7669 – 0.7500 0.7648 0.717 0.7345

0.2519 0.2559 (  1.56) 0.2562 0.2809 0.226 0.2029 0.2576 – 0.2023 – 0.263 –

0.0490 0.0505 ( 2.97) 0.0511 0.0500 0.0524 0.0498 0.0523 – 0.0499 0.0487 0.0476 –

a

Mess size. Percentage error¼ [(present result  3D elasticity result)/3D elasticity result]  100. c Layer-wise theory. d TSDT. b

A 4-node quadratic plate element (7 field variable per node) is used by Kulkarni and Kapuria [12], where the problem of C1 continuity is overcome by defining the derivatives of transverse displacements in terms of separate field variables. The results reported by Sheikh and Chakrabarti [17] based on HSDT using 6-node non-conforming triangular element are also shown in Table 2 for comparison. In addition to this, the analytical results

reported by Reddy [3], Kant and Swaminathan [5] and Aydogdu [13] are also used to show the accuracy of the present model. From Table 2, it may be observed that the performance of the present FE model is quite good compared to the other models especially at lower thickness ratio (h/a) and the present model also uses less number of element field variables. The percentage error (with respect to the 3D results) in predicting

26

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

the non-dimensional transverse deformations is less than 1.2 while it lies around 1 for in-plane normal stresses. The maximum percentage error in calculating the transverse shear stress is 1.85 and for that of transverse in-plane shear stress is 2.97. These values should be considered quite accurate in predicting results of a complex 3D structure by using only a 2D model. The plate is also analyzed for different aspect ratios (b/a ¼3) by taking thickness ratio as 0.01, 0.02, 0.05, 0.10 and 0.25. The results obtained for the non-dimensional transverse displacements and stresses are presented in Table 3 along with many other published results based on different shear deformation theories [12,13,17,31,33]. The present 2D results are in good agreement with 3D results compared to the other results reported in the literature [12,13,17,31,33] as well as the results based on 3D hybrid stress finite element method [64]. The maximum percentage errors in predicting the in-plane stresses (sy ) and transverse shear stress (txz ) are 6.5 and 2.4 respectively and for the other non-dimensional results, it is around 1.0.

The present 2D results are compared with 3D results reported by Piskunov et al. [65] in Table 4. The present results based on two dimensional theory are reasonably close to the three dimensional analysis results reported by Piskunov et al. [65].

3.3. Symmetric laminated sandwich plate (f/c/f) In this example, a simply supported symmetric laminated square sandwich plate (0/C/0) is considered for the analysis. The Table 4 Deflection and stresses of a simply supported un-symmetric square angle ply laminate (  y/y/  y/y) under sinusoidal load. h/a

Angle

Ref.

0.10

15

Present Piskunov et al. [66] % errora Present Piskunov et al. [66] % error Present Piskunov et al. [66] % error

30

3.2. Angle-ply laminated plate under sinusoidal loading An example of a simply supported angle-ply laminated plate ( y/y/ y/y) subjected to sinusoidal loading is considered to study the performance of the present 2D FE model in solving the problems having different angle ply orientations. The material properties are same as in the previous example. The layers are of equal thickness.

45

a Percentage result]  100.

error¼ [(Present

0.5951 0.6150 (3.24) 0.5339 0.5619 (4.98) 0.5207 0.5430 (4.11)

result  3D

sx

sy

txy

0.4650 0.4453 (4.42) 0.3033 0.2833 (7.06) 0.1809 0.1749 (3.43)

0.6500 0.6153 (5.64) 0.1192 0.1115 (6.91) 0.1809 0.1749 (3.43)

0.6581 0.6739 (2.35) 1.2283 1.3310 (7.72) 1.5634 1.6420 (4.84)

elasticity

result)/3D

elasticity

Table 3 Non-dimensional deflection and stresses of a simply supported rectangular (b/a ¼3) cross ply laminate (0/90/0) under sinusoidal load. h/a

Ref.

w

sx

sy

txz

tyz

txy

0.01

Present Pagano [63] % errora Aydogdu [13] Fares and Elmarghany [33] Akhras and Li [31] Kulkarni and Kapuria [12] Sheikh and Chakrabarti [17]

0.5052 0.5030 (0.44) 0.508 0.508 0.508 0.5078 0.5097

0.6276 0.6240 (0.58) 0.624 0.624 0.632 0.6255 0.6457

0.0253 0.0253 (0.0) 0.0253 0.0253 0.0256 0.0253 0.0253

0.4337 0.4390 (  1.21) 0.335 – 0.449 0.2893 0.2847

0.0108 0.0108 (0.0) 0.0134 – 0.0106 0.0130 0.0129

0.0083 0.0083 (0.0) 0.0083 0.0083 0.0083 0.0083 0.0084

0.02

Present Pagano [63] % error Aydogdu [13] Fares and Elmarghany [33]

0.5181 0.5205 ( 0.46) 0.520 0.520

0.6310 0.6276 (0.54) 0.626 0.627

0.0258 0.0260 ( 0.77) 0.0258 0.0258

0.4332 0.4387 (  1.25) 0.335 –

0.0110 0.0110 (0.0) 0.0135 –

0.0085 0.0085 (0.0) 0.0084 0.0084

0.05

Present Pagano [63] % error Aydogdu [13] Fares and Elmarghany [33] Sheikh and Chakrabarti [17] Liou and Sun [64]

0.6074 0.6100 ( 0.43) 0.602 0.609 0.5965 0.611

0.6537 0.6500 (0.57) 0.644 0.644 0.6634 0.653

0.0291 0.0299 ( 2.68) 0.0293 0.0294 0.0274 0.0298

0.4291 0.4340 (  1.13) 0.334 – 0.2859 0.450

0.0118 0.0119 (  0.84) 0.0145 – 0.0135 0.0118

0.0093 0.0093 (0.0) 0.0092 0.0092 0.0092 0.0093

0.10

Present Pagano [63] % error Aydogdu [13] Fares and Elmarghany [33] Akhras and Li [31] Kulkarni and Kapuria [12] Sheikh and Chakrabarti [17] Liou and Sun [64]

0.9181 0.9190 ( 0.10) 0.891 0.918 0.920 0.8636 0.8649 0.921

0.7297 0.7250 (0.65) 0.706 0.703 0.740 0.6939 0.7164 0.709

0.0405 0.0435 ( 4.90) 0.0418 0.0417 0.0425 0.0399 0.0383 0.0429

0.4151 0.4200 (  1.17) 0.330 – 0.436 0.2866 0.2851 0.428

0.0146 0.0152 (  3.95) 0.0118 – 0.0142 0.0170 0.0106 0.0151

0.0122 0.0123 ( 0.81) 0.0179 0.0120 0.0122 0.0115 0.0117 0.0121

0.25

Present Pagano [63] % error Aydogdu [13] Fares and Elmarghany [33] Akhras and Li [31] Kulkarni and Kapuria [12] Sheikh and Chakrabarti [17] Liou and Sun [64]

2.8373 2.8200 (0.61) 2.756 2.843 2.757 2.6452 2.6437 2.828

1.1045 1.1000 (0.41) 1.112 0.984 1.220 1.0377 1.0650 1.077

0.1113 0.1190 (6.47) 0.118 0.1103 0.109 0.1030 0.1209 0.108

0.3777 0.3870 (  2.40) 0.308 – 0.377 0.2730 0.2723 0.360

0.0329 0.0334 (  1.50) 0.0373 – 0.0277 0.0349 0.0320 0.0326

0.0276 0.0281 ( 1.78) 0.0278 0.0266 0.109 0.0263 0.0264 0.0269

a

Percentage error¼ [(Present result  3D elasticity result)/3D elasticity result]  100.

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

material properties used in this problem are as shown in Table 1 for sandwich plates. The thickness distribution amongst the layers (f/c/f) is (0.10h/0.8h /0.1h), where h ( ¼1 in.) is the overall thickness of the plate. The plate is analyzed for sinusoidal load. The non-dimensional transverse displacement and non-dimensional stresses for different thickness ratio are reported in Table 5. The results reported by the Pandit et al. [58] and Singh et al. [66] based on zigzag theory incorporating penalty approach as well as the results reported by Ramtekkar et al. [54] using hybrid formulation based on layerwise theory are used to check the accuracy of present model. It may be observed in Table 6 that the present results are sufficiently close to the results based on 3-D elasticity solution [63] as compared to the other models [54,58,61]. This shows the better ability of present FE model in predicting the response of soft core sandwich structure. The difference of the present results compared to the 3D elasticity results [63] is below 1% for transverse displacements and in-plane stresses and it is 7% for transverse shear stresses.

3.4. Un-symmetric laminated sandwich plate with different boundary conditions In this example an un-symmetrical simply supported laminated sandwich plate (0/90/C/0/90) is considered for the analysis under the sinusoidal loading. The core has a thickness of 0.8h while the two laminated faces are of 0.1h each, where h ( ¼1 in.) is the overall thickness of the plate. Material properties are as used in the previous example of sandwich plate (Table 1). The values of non-dimensional transverse displacements and the stresses are presented in Table 6 for different thickness ratio (h/a) ranging from 0.01 to 0.50. The percentage increase in transverse

27

deformation for the higher thickness ratios is much higher than those of lower thickness ratios for both the boundary conditions. This is due to the effects of transverse flexibility of the core and shear deformation, which are more pronounced at higher thickness ratios. Also in this case the thickness of the core is more (0.8h) so the compressibility of the core plays a very important role and it is reflected in Table 6. The results are in good agreement with the similar results reported by Pandit et al. [58] and Singh et al. [66]. The variations of the in-plane normal stress and the in-plane shear stress across the depth are also plotted in Figs. 3 and 4 with those obtained by the 3-D elasticity solution [63] for simply supported boundary conditions. The thickness ratio (h/a) is considered as 4. The results of these two models are found to match quite well in both the cases.

4. Conclusions An improved C0 two dimensional plate finite element (FE) model has been developed in this paper for the static analysis of laminated sandwich plate having soft core. The higher order zigzag shear deformation theory (HOZT) is used to define the inplane displacement fields whereas a quadratic displacement field is used to define the transverse displacement in the proposed model. The continuity requirement of the derivates of transverse displacement is circumvented by effectively choosing the nodal field variables. There is no need to use any penalty multiplier in the formulation as used by many previous researchers. A nine node C0 quadratic plate finite element is successfully implemented in the present model. Many numerical examples are solved for different problems of laminated composite and sandwich plates.

Table 5 Non-dimensional deflection and stresses at important points of a simply supported square sandwich plate (0/C/0) under sinusoidal loading. h/a

Ref.

w

sx

sy

txz

tyz

txy

0.01

Present Pagano [63] % errora Pandit et al.[58] Tu et al.[14] Singh et al.[66]

0.8814 0.8917 (  1.16) 0.8917 0.8919 0.9017

1.0982 1.0980 (0.02) 1.1093 1.1069 1.1020

0.0592 0.0550 (7.64) 0.0547 0.0573 –

0.3426 0.3240 (5.74) 0.3412 0.3312 0.4079

0.0322 0.0297 (8.42) 0.0324 0.0337 –

0.0433 0.0433 (0.0) 0.0434 0.0432 0.0453

0.02

Present Pagano [63] % error Pandit et al.[58] Singh et al. [66]

0.9234 0.9348 (  1.23) 0.9341 0.9458

1.0997 1.0990 (0.06) 1.0948 1.1050

0.0611 0.0569 (7.38) 0.0566 –

0.3300 0.3230 (2.17) 0.3403 0.3617

0.0321 0.0306 (4.90) 0.0333 –

0.0443 0.0446 ( 0.67) 0.0445 0.0465

0.05

Present Pagano [63] % error Pandit et al. [58] Singh et al. [66] Ramtekkar et al. [54]

1.2121 1.2264 (  1.17) 1.2252 1.2424 –

1.1103 1.1161 (  0.52) 1.1055 1.1100 1.110

0.0742 0.0700 (6.0) 0.0694 – 0.070

0.3272 0.3174 (3.09) 0.3342 0.3429 0.317

0.0399 0.0361 (1.79) 0.0392 – 0.036

0.0508 0.0511 ( 0.58) 0.0.509 0.0536 0.051

0.10

Present Pagano [63] % error Pandit et al. [58] Tu et al. [14] Singh et al. [66] Ramtekkar et al. [54]

2.1775 2.2004 (  1.04) 2.0020 2.2027 2.2389 –

1.1528 1.1530 (  0.02) 1.1483 1.1466 1.1594 1.159

0.1143 0.1104 (3.53) 0.1086 0.1105 – 0.111

0.3058 0.3000 (1.93) 0.3158 0.3181 0.3287 0.303

0.0575 0.0557 (2.28) 0.0570 0.0532 – 0.055

0.0705 0.0707 ( 0.28) 0.0709 0.0715 0.0707 0.071

0.25

Present Pagano [63] % error Pandit et al. [58] Tu et al. [14] Singh et al. [66] Ramtekkar et al. [54]

7.5822 7.5962 (  0.18) 7.6552 7.5610 7.8556 –

1.5306 1.5560 (  0.02) 1.5218 1.5518 1.5488 1.570

0.2581 0.2595 (  0.54) 0.2506 0.2483 – 0.260

0.2436 0.2390 (1.93) 0.2520 0.2447 0.2611 0.240

0.1147 0.1072 (7.0) 0.1156 0.1184 – 0.108

0.1445 0.1437 ( 0.56) 0.1468 0.1459 0.1671 0.145

a

Percentage error ¼[(Present result  3D elasticity result)/3D elasticity result]  100.

28

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

Table 6 Non-dimensional deflection and stresses at important points anti-symmetric square sandwich plate (0/90/core/0/90) under sinusoidal loading for different boundary conditions. h/a

B.C.

Ref.

w

sx

sy

txz

tyz

txy

0.01

SCSC

Present Pandit et al. [58] Singh et al. [66] Present Pandit et al. [58] Singh et al. [66]

0.3430 0.3453 0.3920 0.2267 0.2286 0.2260

0.4250 0.4077 0.5986 0.4371 0.4270 0.4283

0.0366 0.0326 – 0.0259 0.0228 –

0.0775 0.0778 0.0944 0.2171 0.2189 0.2348

0.3311 0.3086 – 0.2327 0.2189 –

0.0028 0.0036 – 0.0007 0.0005 –

Present Pandit et al. [58] Singh et al. [66] Present Pandit et al. [58] Singh et al. [66]

0.6022 0.6052 0.6080 0.4283 0.4296 0.4442

0.6046 0.5850 0.6138 0.4372 0.4275 0.4293

0.0378 0.0334 – 0.0268 0.0236 –

0.1012 0.1061 0.1542 0.1568 0.1828 0.2004

0.3209 0.2527 – 0.2319 0.1828 –

0.0084 0.0107 – 0.0018 0.0018 –

Present Pandit et al. [58] Singh et al. [66] Present Pandit et al. [58] Singh et al. [66]

1.2994 1.3026 1.3092 1.0484 1.0489 1.0213

0.8566 0.8310 0.7392 0.4708 0.4597 0.4621

0.0426 0.0372 – 0.0315 0.0279 –

0.1366 0.1418 0.1523 0.1308 0.1587 0.1651

0.2696 0.1967 – 0.2175 0.1586 –

0.0161 0.0189 – 0.0033 0.0036 –

Present Pandit et al. [58] Singh et al. [66] Present Pandit et al. [58] Singh et al. [66]

5.6453 3.8087 3.8500 5.2305 3.4521 3.3421

1.3176 1.1415 1.0189 0.7432 0.6170 0.6022

0.0786 0.0561 – 0.0699 0.0470 –

0.1667 0.1683 0.1620 0.1157 0.1396 0.1422

0.1907 0.1539 – 0.1766 0.1394 –

0.0272 0.0292 – 0.0085 0.0074 –

Present Pandit et al. [58] Singh et al. [66] Present Pandit et al. [58] Singh et al. [66]

20.1918 19.5512 19.580 19.0444 18.3454 18.345

2.1923 2.3071 2.4028 1.6291 1.8156 1.8150

0.1995 0.2047 – 0.1823 0.1902 –

0.1659 0.1691 0.1721 0.1182 0.1227 0.1325

0.1512 0.1296 – 0.1427 0.1217 –

0.0482 0.0667 – 0.0167 0.0215 –

CCCC

0.05

SCSC

CCCC

0.01

SCSC

CCCC

4

SCSC

CCCC

2

SCSC

0.6

0.75

0.4

0.50 Normalized depth

Normalized depth

CCCC

0.2 Present 0.0

Pagano[63]

-0.2

Present Pagano[63]

0.00 -0.25 -0.50

-0.4 -0.6 -1.5 -1.2 -0.9 -0.6 -0.3 0.0 0.3 0.6 0.9 Non-dimensional in-plane normal stress

0.25

1.2

-0.75 -0.08

1.5

Fig. 3. Variation of in-plane normal stress across the depth of un-symmetric sandwich plate (0/90/C/0/90) (h/a ¼0.25).

The results obtained by using the present FE model are successfully compared with many published results. The numerical results show that the performance of the present finite element model is excellent in predicting the response of thin and thick laminated composites as well as soft core sandwich structures as the percentage error with respect to the 3D elasticity solution is quite low. The present FE model may, therefore, be recommended for the accurate analysis of laminated sandwich plates with/ without soft cores.

-0.06

-0.04 -0.02 0.00 0.02 0.04 0.06 Non-dimensional in-plane shear stress

0.08

Fig. 4. Variation of in-plane shear stress across the depth of un-symmetric sandwich plate (0/90/C/0/90) (h/a¼ 0.25).

Appendix A l1 ¼

zðz þhl Þ hu ðhu þhl Þ

l2 ¼

ðhl þ zÞðhu zÞ hu hl

l3 ¼

zðhu zÞ hl ðhu þ hl Þ

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

Appendix B þ 2

3

b1 ¼ 1 þ A11 z þ A31 z þ

nu1 X

þ

Aði þ 4Þ,1 ðzziu ÞHðzziu Þ c5 ¼ A25 z2 þA45 z3 þ

Aðnu þ 3 þ jÞ,1 ðzzjl ÞHðz þzjl Þ þ 3

b2 ¼ A12 z þA32 z þ

nu1 X

þ

c6 ¼ A26 z2 þA46 z3 þ

Aðnu þ 3 þ jÞ,2 ðzzjl ÞHðz þzjl Þ þ nu1 X

Aði þ 4Þ,3 ðzziu ÞHðzziu Þ

n l1 X

Aðnu þ nl þ 2 þ iÞ,6 ðzziu ÞHðzziu Þ

Að2nu þ nl þ 1 þ jÞ,6 ðzzjl ÞHðz þ zjl Þ

c7 ¼ A27 z2 þA47 z3 þ

nu1 X

Aðnu þ nl þ 2 þ iÞ,7 ðzziu ÞHðzziu Þ

i¼1

Aðnu þ 3 þ jÞ,3 ðzzjl ÞHðz þzjl Þ þ

j¼1

b4 ¼ A14 z2 þA34 z3 þ

nu1 X

n l1 X

Að2nu þ nl þ 1 þ jÞ,7 ðzzjl ÞHðz þ zjl Þ

j¼1

Aði þ 4Þ,4 ðzziu ÞHðzziu Þ

i¼1

þ

nu1 X

j¼1

i¼1

n l1 X

Að2nu þ nl þ 1 þ jÞ,5 ðzzjl ÞHðz þ zjl Þ

i¼1

b3 ¼ z þ A13 z2 þ A33 z3 þ

þ

nl1 X

Aði þ 4Þ,2 ðzziu ÞHðzziu Þ

j¼1

n l1 X

Aðnu þ nl þ 2 þ iÞ,5 ðzziu ÞHðzziu Þ

j¼1

i¼1 n l1 X

nu1 X i¼1

j¼1 2

Að2nu þ nl þ 1 þ jÞ,4 ðzzjl ÞHðz þ zjl Þ

j¼1

i¼1 n l1 X

n l1 X

29

c8 ¼ A28 z2 þA48 z3 þ

nu1 X

Aðnu þ nl þ 2 þ iÞ,8 ðzziu ÞHðzziu Þ

i¼1

Aðnu þ 3 þ jÞ,4 ðzzjl ÞHðz þzjl Þ

þ

j¼1

b5 ¼ A15 z2 þA35 z3 þ

nu1 X

n l1 X

Að2nu þ nl þ 1 þ jÞ,8 ðzzjl ÞHðz þ zjl Þ

j¼1

Aði þ 4Þ,5 ðzziu ÞHðzziu Þ

i¼1

þ

n l1 X

Aðnu þ 3 þ jÞ,5 ðzzjl ÞHðz þzjl Þ

Appendix C

j¼1

b6 ¼ A16 z2 þA36 z3 þ

nu1 X

Aði þ 4Þ,6 ðzziu ÞHðzziu Þ

i¼1

þ

nl1 X

Aðnu þ 3 þ jÞ,6 ðzzjl ÞHðz þzjl Þ

j¼1

b7 ¼ A17 z2 þA37 z3 þ

nu1 X

Aði þ 4Þ,7 ðzziu ÞHðzziu Þ

Elements of matrix H 2

0 6 60 6 60 6 ½H ¼ 6 60 6 6 a1 4 e1

0

0

0

0

0

0

0

0

0

0

b1

0

b2

0

0

0 0 0 a2

0 d2 0 0

0 0 0 a3

0 0 0 a4

0 0 0 a5

0 0 0 a6

0 d1 0 0

0 0 0 a7

0 0 0 a8

0 d3 0 0

0 0 c1 0

c1 0 b1 0

0 0 c2 0

c2 0 b2 0

0 0 0 l2

e2

0

e3

e4

e5

e6

0

e7

e8

0

0

0

0

0

0

0

i¼1

þ

n l1 X

Aðnu þ 3 þ jÞ,7 ðzzjl ÞHðz þzjl Þ

j¼1

b8 ¼ A18 z2 þA38 z3 þ

nu1 X

Aði þ 4Þ,8 ðzziu ÞHðzziu Þ

i¼1

þ

n l1 X

Aðnu þ 3 þ jÞ,8 ðzzjl ÞHðz þzjl Þ

j¼1

c1 ¼ A21 z2 þ A41 z3 þ

nu1 X

n l1 X

b3

0

b4

0

b5

0

b6

0

0

0

b7

0

b8

0

0

0 0

0 0

c3 0

0 0

c4 0

0 0

c5 0

0 0

c6 0

0 0

0 0

0 0

c7 0

0 0

c8 0

0 0

0

c3

b3

c4

b4

c5

b5

c6

b6

0

0

c7

b7

c8

b8

0

0

0

0

0

0

0

0

0

0

l1

0

0

0

0

0

l3

l2

0

0

0

0

0

0

0

0

0

l1

0

0

0

0

0

where a1 ¼

Aðnu þ nl þ 2 þ iÞ,1 ðzziu ÞHðzziu Þ

i¼1

þ

0

Að2nu þ nl þ 1 þ jÞ,1 ðzzjl ÞHðz þzjl Þ

a2 ¼ a3 ¼

j¼1

c2 ¼ 1 þ A22 z2 þ A42 z3 þ

nu1 X

a4 ¼ Aðnu þ nl þ 2 þ iÞ,2 ðzziu ÞHðzziu Þ a5 ¼

i¼1

þ

n l1 X

Að2nu þ nl þ 1 þ jÞ,2 ðzzjl ÞHðz þzjl Þ

j¼1

c3 ¼ A23 z2 þ A43 z3 þ

nu1 X

Aðnu þ nl þ 2 þ iÞ,3 ðzziu ÞHðzziu Þ

i¼1

þ

nl1 X

a6 ¼ a7 ¼ a8 ¼

@b1 @z @b2 @z @b3 @z @b4 @z @b5 @z @b6 @z @b7 @z @b8 @z

Að2nu þ nl þ 1 þ jÞ,3 ðzzjl ÞHðz þzjl Þ

j¼1

c4 ¼ z þ A24 z2 þA44 z3 þ

nu1 X i¼1

@l1 @z @l2 d2 ¼ @z

d1 ¼ Aðnu þ nl þ 2 þ iÞ,4 ðzziu ÞHðzziu Þ

3

7 07 7 07 7 7 07 7 07 5 l3

30

d3 ¼ e1 ¼ e2 ¼ e3 ¼ e4 ¼ e5 ¼ e6 ¼ e7 ¼ e8 ¼

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

@l3 @z @c1 @z @c2 @z @c3 @z @c4 @z @c5 @z @c6 @z @c7 @z @c8 @z

References [1] V.K. Goyal, R.K. Kapania, A shear deformeable beam element for analysis of laminated composites, Finite Elem. Anal. Des. 43 (2007) 463–477. [2] A.J.M. Ferreira, A formulation of the multiquadratic radial basis function method for the analysis of laminated composite plates, Compos. Struct. 59 (2003) 385–392. [3] J.N. Reddy, A simple higher order theory for laminated composite plates, J. Appl. Mech. 51 (1984) 745–752. [4] Tarun Kant, Numerical analysis of thick plates, Comput. Methods Appl. Mech. Eng. 44 (4) (1982) 1–18. [5] T. Kant, A. Swaminathan, Analytical solutions for the static analysis of laminated composite and sandwich plates based on refined higher order shear deformation theory, Compos. Struct. 56 (2002) 329–344. [6] M.R. Aagaah, M. Mahinfalah, G.N. Jazar, Linear static analysis and finite element modeling for laminated composite plates using third order shear deformation theory, Compos. Struct. 62 (2003) 27–39. [7] M.V.V.S. Murthy, D.R. Mahapatra, K. Badarinarayana, S. Gopalakrishnan, A refined higher order finite element for asymmetric composite beams, Compos. Struct. 67 (2005) 27–35. [8] P. Subramanian, Dynamic analysis of laminated composite beams using higher order theories and finite elements, Compos. Struct. 73 (2006) 342–353. [9] T. Pervez, A.C. Seibi, F.K.S. Al-Jahwari, Analysis of thick orthotropic laminated composite plates based on higher order shear deformation theory, Compos. Struct. 71 (2005) 414–422. [10] Z. Wu, R. Chen, W. Chen, Refined laminated plate element based on global local higher order shear deformation theory, Compos. Struct. 70 (2005) 135–152. [11] Z. Wu, S.H. Lo, K.Y. Sze, W. Chen, A higher order finite element including transverse normal strain including for linear elastic composite plates with general lamination configurations, Finite Elem. Anal. Des. 48 (2012) 1346–1357. [12] S.D. Kulkarni, S. Kapuria, A new discrete Kirchhoff quadrilateral element based on the third order theory for composite plates, Comput. Mech. 39 (2007) 237–246. [13] M. Aydogdu, A new shear deformation theory for laminated composite plates, Compos. Struct. 89 (2009) 94–101. [14] T.M. Tu, L.N. Thach, T.H. Quoc, Finite element modeling for bending and vibration analysis of laminated and sandwich composite plates based on higher order theory, Comput. Mater. Sci. 49 (2010) 390–394. [15] A.J.M. Ferreira, C.M.C. Roque, E. Carrera, M. Cinefra, Analysis of thick isotropic and cross ply laminated plates by radial basis function and a unified formulation, J. Sound Vib. 330 (2011) 771–787. [16] C.M.C. Roque, D. Cunha, C. Shu, A.J.M. Ferreira, A local radial basis functionsfinite differences techniques for the analysis of composite plates, Eng. Anal. Bound. Elem. 35 (2011) 363–374. [17] A.H. Sheikh, A. Chakrabarti, A new llate bending element based on higher order shear deformation theory for the analysis of composite plates, Finite Elem. Anal. Des. 39 (2003) 883–903. [18] D.H. Robbins Jr., J.N. Reddy, Theories and computational models for composite laminates, Allied Mech. Rev. 49 (1996) 155–199. [19] A. Toledano, H. Murakami, Composite plate theory for arbitrary laminate configuration, J. Appl. Mech. 54 (1987) 181–189. [20] X. Lu, D. Liu, An interlaminar shear stress continuity theory for both thin and thick laminates, ASME J. Appl. Mech. 59 (1992) 502–509. [21] J.N. Reddy, A generalization of two dimensional theories of laminated composite plates, Commun. Appl. Numer. Methods 3 (1987) 173–180.

[22] M. Tahani, Analysis of laminated composite beams using layer-wise displacement theories, Compos. Struct. 79 (2007) 535–547. [23] S.S. Ramesh, C.M. Wang, J.N. Reddy, K.K. Ang, A higher order plate element for the accurate prediction of interlaminar stresses in laminated composite plates, Compos. Struct. 91 (2009) 337–357. [24] H. Murakami, Laminated composite plate theory with improved in-plane responses, J. Appl. Mech. 53 (1986) 661–666. [25] M. Di Scuiva, An improved shear deformation theory for moderately thick multilayered anisotropic shells and plates, J. Appl. Mech. 54 (1987) 589–596. [26] S.P. Lee, N.R. Senthilnathan, S.P. Lim, S.T. Chow, An improved zigzag model for the bending of laminated composite plates, Compos. Struct. 15 (1990) 137–148. [27] M. Cho, R.R. Parmerter, Efficient higher order composite plate theory for general lamination configurations, AIAA J. 31 (7) (1993) 1299–1306. [28] Y.B. Cho, R.C. Averill, An improved theory and finite element model for laminated beams using first order zigzag sub-laminate approximations, Compos. Struct. 37 (1997) 281–298. [29] E. Carrera, Historical review of zig-zag theories for multilayered plates and shells, Appl. Mech. Rev. 56 (3) (2003) 287–308. [30] A. Chakrabarti, A.H. Sheikh, A new triangular element to model inter-laminar shear stress continuous plate theory, Int. J. Numer. Methods Eng. 60 (2004) 1237–1257. [31] G. Akhras, W. Li, Spline finite strip analysis of composite plates based on higher order zigzag composite plate theory, Compos. Struct. 78 (2007) 112–118. [32] S. Kapuria, S.D. Kulkarni, An improved discrete Kirchhoff element based on third order zigzag theory for static analysis of composite and sandwich plates, Int. J. Numer. Methods Eng. 69 (2007) 1948–1981. [33] M.E. Fares, M.Kh. Elmarghany, A refined zigzag non-linear first order shear deformation theory of composite laminated plates, Compos. Struct. 82 (2008) 71–83. [34] A.J.M. Ferreira, C.M.C. Roque, E. Carrera, M. Cinefra, O. Polit, Radial basis functions collocation and a unified formulation for bending, vibration and buckling analysis of laminated plates, according to variation of Murakami’s zigzag theory, Eur. J. Mech. 30 (4) (2011) 559–570. [35] J.D. Rodrigues, C.M.C. Roque, A.J.M. Ferreira, E. Carrera, M. Cinefra, Radial basis functions-finite differences collocation and a unified formulation for bending, vibration and buckling analysis of laminated plates, according to Murakami’s zigzag theory, Compos. Struct. 93 (7) (2011) 1613–1620. [36] R.P. Shimpi, Y.M. Ghugal, A new layerwise trigonometric shear deformation theory for two layered cross-ply beams, Compos. Sci. Technol. 61 (2001) 1271–1283. [37] R.P. Shimpi, A.V. Ainapure, A beam finite element based on layerwise trigonometric shear deformation theory, Compos. Struct. 53 (2001) 153–162. [38] H. Arya, R.P. Shimpi, N.K. Naik, A zigzag model for laminated composite beams, Compos. Struct. 56 (1) (2002) 21–24. [39] U. Icardi, A three dimensional zigzag theory for analysis of thick laminated beams, Compos. Struct. 53 (2001) 123–135. [40] U. Icardi, Applications of zigzag theories to sandwich beams, Mech. Adv. Mater. Struct. 10 (2003) 77–97. [41] Y.C. Yip, R.C. Averill, Thick beam theory and finite element model with zigzag sub-laminate approximations, J. AIAA 34 (1996) 1627–1632. [42] Y.B. Cho, R.C. Averill, First order zigzag sub-laminate plate theory and finite element model for laminated composite and sandwich panels, Compos. Struct. 50 (2000) 1–15. [43] R.C. Averill, Static and dynamic response of moderately thick laminated beams with damage, Compos. Eng. J. 4 (1994) 381–395. [44] M. Di Sciuva, Development of anisotropic multilayered shear deformable rectangular plate element, Comput. Struct. 21 (1985) 789–796. [45] M. Di Sciuva, A general quadrilateral multilayered plate element with continuous interlaminar stresses, Comput. Struct. 47 (1993) 91–105. [46] E. Carrera, C0 Reissner–Mindlin miltilayered plate element including zigzag and interlaminar stress continuity, Int. J. Numer. Methods Eng. 39 (1996) 1797–1820. [47] R.C. Averill, Y.C. Yip, Development of simple, robust finite elements based on refined theories for thick laminated beams, Compos. Struct. 59 (1996) 661–666. [48] V.R. Aitharaju, R.C. Averill, C0 Zigzag kinematic displacement models for the analysis of laminated composites, Mech. Compos. Mater. Struct. 6 (1999) 31–56. [49] J.S. Kim, M. Cho, Enhanced first-order shear deformation theory for laminated and sandwich plates, J. Appl. Mech. 72 (2005) 809–817. [50] J.S. Kim, M. Cho, Enhanced modeling of laminated and sandwich plates via strain energy transformation, Compos. Sci. Technol. 66 (2006) 1575–1587. [51] J.S. Kim, M. Cho, Enhanced first-order theory based on mixed formulation and transverse normal effect, Int. J. Solids Struct. 44 (3–4) (2007) 1256–1276. [52] M.K. Rao, Y.M. Desai, M.R. Chitnis, Free vibrations of laminated beams using mixed theory, Compos. Struct. 52 (2001) 149–160. [53] G.S. Ramtekkar, Y.M. Desai, Natural vibrations of laminated composite beams by using mixed finite element modeling, J. Sound Vib. 257 (4) (2002) 635–651. [54] G.S. Ramtekkar, Y.M. Desai, A.H. Shah, Application of a three dimensional mixed finite element model to the flexure of sandwich plate, Comput. Struct. 81 (2003) 2383–2398. [55] A.N. Bambole, Y.M. Desai, Hybrid-interface finite element for laminated composite and sandwich beams, Finite Elem. Anal. Des. 43 (2007) 1023–1036.

H.D. Chalak et al. / Finite Elements in Analysis and Design 56 (2012) 20–31

[56] Y. Frostig, Classical and high order computational models in the analysis of modern sandwich panels, Composites: Part B 34 (2003) 83–100. [57] H.S. Givil, O. Rabinovitch, Y. Frostig, High-order non-linear contact effects in the dynamic behavior of delaminated sandwich panel with a flexible core, Int. J. Solids Struct. 44 (2007) 77–99. [58] M.K. Pandit, A.H. Sheikh, B.N. Singh, An improved higher order zigzag theory for the static analysis of laminated sandwich plate with soft-core, Finite Elem. Anal. Des. 44 (2008) 602–610. [59] M.K. Pandit, A.H. Sheikh, B.N. Singh, Buckling of laminated sandwich plates with soft core based on an improved higher order zigzag theory, Thin-Walled Struct. 46 (2008) 1183–1191. [60] M.K. Pandit, A.H. Sheikh, B.N. Singh, Stochastic perturbation based finite element for deflection statistics of soft core sandwich plate with random material properties, Int. J. Mech. Sci. 51 (5) (2010) 14–23. [61] M.K. Pandit, A.H. Sheikh, B.N. Singh, Stochastic free vibration response of soft core sandwich plates using an improved higher order zigzag theory, J. Aerosp. Eng. 23 (2010) 14–23.

31

[62] A. Chakrabarti, H.D. Chalak, A. Iqbal, A.H. Sheikh, A new FE model based on higher order zigzag theory for the analysis of laminated sandwich beam with soft core, Compos. Struct. 93 (2011) 271–279. [63] N.J. Pagano, Exact solutions for rectangular bidirectional composites and sandwich plates, Compos. Mater. J. 4 (1970) 20–35. [64] W.J. Liou, C.T. Sun, A three dimensional hybrid stress isoparametric element for analysis of laminated composite plates, Comput. Struct. 25 (2) (1987) 241–244. [65] V.G. Piskunov, V.E. Verijenko, S. Adali, P.Y. Tabakov, V.K. Prisyazhnyouk, S.G. Buryhin, Rational transverse shear deformation higher order theory of anisotropic laminated plates and shells, Int. J. Solids Struct. 38 (2001) 6491–6523. [66] S.K. Singh, A. Chakrabarti, P. Bera, J.S.D. Sony, An efficient C0 FE model for the analysis of composites and sandwich laminates with general layup, Latin Am. J. Solids Struct. 8 (2011) 197–212.