Analysis of Sandwich Plates by Generalized Differential Quadrature Method

We combine a layer-wise formulation and a generalized differential quadrature technique for predicting the static deformations and free vibration behaviour of sandwich plates. Through numerical experiments, the capability and efficiency of this strong-form technique for static and vibration problems are demonstrated, and the numerical accuracy and convergence are thoughtfully examined.


Introduction
The two big groups of theories for the analysis of laminated plates are based on displacements only or on displacements and stresses at interfaces of the laminate.This paper deals with the first group, where again two approaches can be considered: the equivalent-single layer (ESL) theories (where all the layers are referred to the same degrees of freedom) or layer-wise (LW) theories (where specific degrees of freedom are linked to specific layers).
The classical laminated plate theory and the first-order shear deformation theory [1][2][3][4][5] are the simplest ESL theories, at the cost of introducing artificial shear correction factors.ESL higher-order theories usually consider higherorder expansion of displacements with respect to the thickness coordinate of the plate, with the further advantage of disregarding shear correction factors [6][7][8][9][10][11][12].Most of the previously mentioned ESL theories may ill-represent the transverse shear stresses for sandwich laminates, where the material properties of skins are typically much larger than those of the core.In such problems, LW theories that consider independent degrees of freedom for each layer should be considered instead.Several LW theories have been proposed in literature [13][14][15][16][17][18][19][20][21] usually considering translational degrees of freedom at the layer interfaces.Over the last decades, an automatic approach for generation of any kind of ESL or LW theories was developed by Carrera [18][19][20], which allows the user to freely develop new theories without much effort.Various examples of the use of Carrera's unified formulation with meshless methods by the first author can be found in [22][23][24][25][26][27].
Unlike LW theories that consider translations at each laminate interface, in this work we adopt a layer-wise theory that considers both translations and rotations, based on an expansion of Mindlin's first-order shear deformation theory in each layer.The displacement continuity at layer's interface is enforced.In each layer, the theory produces constant, very accurate transverse shear stress.
As is well known, the analysis of composite plates is typically performed by finite element method.On the contrary, in the present study, a layer-wise shear deformation theory is implemented within generalized differential quadrature technique.The latter meshfree method in combination with this layer-wise theory carries out an accurate calculation z  (3)   x  (2)   x  (1)   x z (3)   z (2)   z (1)   x Layers of the natural frequencies, displacements, and stresses of sandwich plates.This study appears for the first time and serves to fill the gap of knowledge in this research topic.

A Layer-Wise Theory
The layer-wise theory proposed in this paper is based on the assumption of a first-order shear deformation theory in each layer and the imposition of displacement continuity at layer's interfaces.Given that, a layer-wise formulation links the degrees of freedom at each layer.Nevertheless, the problem may become large if a large number of layers are considered, although the present formulation is adequate for a general laminate.In this paper, we concentrate on sandwich plates built from two external, stiff layers and one central, soft core.The three-layer laminate is shown schematically in Figure 1.
In order to simplify the imposition of displacement continuity at the core-skin interfaces, we first define the displacement field for the middle layer (the core in the sandwich plate) as  (2) (, , ) =  0 (, ) +  (2)  (2)   (, ) , where  and V are the in-plane displacements at any point (, , ),  0 and V 0 are the in-plane displacements of the points (, , 0) on the midplane,  is the deflection, and  (2)    and  (2)   are the rotations of the normals to the midplane about the and -axes, respectively.
The present layer-wise formulation neglects the transverse normal deformations, by assuming a constant transverse displacement in every layer.Neglecting  ()   for each orthotropic layer, the stress-strain relations in the fiber local coordinate system can be expressed as where subscripts 1 and 2 are, respectively, the principal middle surface directions, 3 is the direction normal to the middle plane of the plate, and the reduced stiffness components,  ()  , are given by in which  () 1 ,  () 2 , ] () 12 ,  () 12 ,  () 23 , and  () 31 are material properties of the lamina .Note that shear correction factors are not considered in the stress-strain laws, as they must be done in the first-order shear deformation theories [21].
The stress-strain relations in the global -- coordinate system can be obtained by coordinate transformation as By considering  as the angle between -axis and 1-axis, with 1-axis being the first principal material axis, connected usually with fiber direction, the components can be computed by coordinate transformation (see [13] for details).
The equations of motion of this layer-wise theory are derived from the dynamic version of the principle of virtual displacements.In the present work, only symmetric laminates are considered; therefore,  0 , V 0 , and related stress resultants can be discarded.The virtual strain energy (), the virtual kinetic energy (), and the virtual work done by applied forces (), assuming a three-layer laminate, are given by It should be noted that the integrals in the thickness direction are now in the layer  domain, from −ℎ  /2 to ℎ  /2.Consider where boldface quantities are vectors and (⋅) denotes time differentiation.Finally, the virtual work is given as where  quantities are virtual, Ω 0 denotes the midplane of the laminate,  is the external distributed load, and where ,  take the symbols , .Substituting , , and  into the virtual work statement, using the fundamental lemma of the calculus of variations, and integrating by parts to relieve from any derivatives of the generalized displacements (see [13] for details), the equations of motion with respect to 7 degrees of freedom ( 0 ,  (1)   ,  (1)   ,  (2)   ,  (2)   ,  (3)   , and  (3)  ) are obtained.Consider : +  (1)    =  (1)  0 ( (1)   : +  (1)    =  (1)  0 ( : : : where where   , the specific mass of the material of the  is layer and ℎ  is the thickness of the th layer.

The Generalized Differential Quadrature Method
The differential quadrature (DQ) method was presented by Bellman and his associates in the early 1970s [28].It is a numerical discretization technique for the approximation of derivatives.The DQ method was initiated from the idea of conventional integral quadrature.In fact, one problem which arises frequently in structural mechanics and in many other engineering problems is the evaluation of the integral ∫   ().If  is a function such that / = , then the value of the given integral is () − ().Unfortunately, in practical problems, it is extremely difficult, and most of the times impossible, to obtain an explicit expression for .The values of , perhaps, can be known at a discrete set of points, and, in this situation, a numerical approach is essential.Following the idea of integral quadrature, Bellman et al. [28] suggested that the first-order derivative of the function () with respect to  at a grid point   is approximated by a linear sum of all the functional values in the whole domain as where  , (1)    represent the weighting coefficients,  is the total number of grid points  1 ,  2 , . . .,   in the whole domain, and (  ) is the calculated value of () at the point  =   .However, the weighting coefficients which were used depend, in the Bellman approach, on the number of grid points and the type of discretization that was used.Some years later the DQ approach has been extended by Shu [29], who introduced the generalized differential quadrature (GDQ) method.In [29], the weighting coefficients calculation is performed through recursive formulae, and they are independent of the number of grid points  and the discretization type.The interested reader can find a brief review on GDQ applications in [30][31][32][33][34][35][36][37][38][39][40].
To compute the first-order weighting coefficients  ,(1)  , the following algebraic formulae are derived: where  (1) So any derivative of an order greater than the first can be written as for  = 1, 2, . . ., ,  = 2, 3, . . .,  − 1. ( The previous derivative discretization was presented for the one-dimension problems only.Nevertheless, the simple onedimensional case can be directly extended to the multidimensional one for any regular shape, such as a rectangular or circular.Given a two-dimensional physical system, in Cartesian coordinates, described by a function (, ), the problem values depend on the nodal coordinates  and .Since GDQ works with regular grid along the main axes,  and  will indicate the points along  and , respectively.Using the discretization of the derivatives in the one-dimensional case (23), a generic higher-order derivative along  and  for the two-dimensional case can be reported.Consider where  ,()  and

𝑦,(𝑚) 𝑗𝑙
are the weighting coefficients of order  and  along  and , respectively.Moreover, the mixed derivative can be written in the same way as (24).Consider where  ,()  and

Grid Distributions.
As it has been widely proven [34][35][36][37], the GDQ numerical accuracy strongly depends on the grid distribution choice for a certain numerical problem.On the one hand, the natural and simplest choice of uniform grids does not lead to stable and accurate results for any GDQ computation; on the other hand, the Chebyshev-Gauss-Lobatto (C-G-L) grid point distribution gives accurate results in most of cases.In the numerical results proposed in this paper a C-G-L grid distribution is considered for both directions  and .Consider Mathematical Problems in Engineering 7 where ,  are the total number of points in the given directions , , respectively.

The Static
Problem.The partial differential system of equations developed in Section 2 has been discretised and numerically solved using the GDQ method.Since GDQ allows turning any derivative in an algebraic expression, the fundamental system of equations is solved in its strong form directly, within its boundary conditions.A static boundary problem can be indicated as follows: where ℓ is a partial differential operator and contains all the spatial derivatives of the problem with respect to a Cartesian coordinate system x.The model degrees of freedom are (x), and the external forces acting on the domain Ω are indicated by (x).It is obvious that the differential problem cannot be solved without its boundary conditions, where ℓ  is the differential operator related to the boundary conditions of the given system of equations and g are the boundary loads.
Once the boundary-value problem ( 27) is discretised by GDQ the following set of algebraic equations is obtained: where U is a vector representing a set of unknown functional values at all the interior points, L is a matrix carried out from the GDQ method, and Q, G are the known vectors arising from the functions (x) and g.
In order to obtain the global stiffness matrix of the physical problem under study, the discretised boundary conditions must be embedded into the discrete fundamental system of equations.For the sake of conciseness, the global static system (28) is reported in matrix form as follows: where Q is the external boundary and domain forces and K is the global stiffness matrix which can be inverted leading to the problem solution.The interested reader can refer to the paper by Viola et al. [39] for further details about static problem resolutions of plates and shells.

The Eigenproblem.
Given a generic set of equations of motion, the linear eigenvalue problem can be achieved eliminating the external forces, as it is shown in the following: where the right-hand term represents the inertia term, given by the product of density  and the acceleration term ü (x, ).
It is noted that the symbol ̈stands for the second derivative respect to time of a given variable.In general, an eigenproblem seeks particular numbers, called eigenvalues (), and certain vectors, called eigenvectors (x).Mathematically speaking, the eigenvalue problem is obtained from (30) setting as follows: Equation (31) indicates that variable separation technique is employed.In fact, the problem degree of freedom  which is function of space and time is given by the product of two functions, which are functions of space and time separately.Hence, substituting relation ( 31) into (30), a generalized eigenvalue problem can be achieved.For the sake of completeness, it is presented with its boundary conditions as Applying GDQ rule, the eigenvalue problem ( 32) can be written in discretised form as follows: where M indicates the mass matrix of the structural system under study.Solving the discretised system (33), the eigenvalues  and eigenvectors U of the given system can be obtained.Finally, the reader can refer to [38] for solutions of eigenvalue problems regarding plates and shells.The plate is subjected to a sinusoidal vertical pressure of the form   =  sin(/) sin(/) with the origin of the coordinate system located at the lower left corner on the midplane and  the maximum load (at the plate center).
The in-plane displacements, the transverse displacements, the normal stresses, and the in-plane and transverse shear stresses are presented in normalized form as In Table 1 for the present layer-wise theory, using 7 × 7 up to 19 × 19 grid points, the results with higher-order solutions by Akhras et al. [41] and Reddy [11], FSDT solutions by Reddy and Chao [42], and an exact solution by Pagano [43] are compared.Furthermore, the numerical computations by the authors using radial basis functions (RBFs) collocation, by the authors with Reddy's theory [44], and a layer-wise theory [45] are compared too.The present layer-wise displacements are in excellent agreement for thinner or thicker plates.Highly accurate normal stresses and transverse shear stresses are Transverse displacement and stresses for a sandwich plate are indicated in Tables 2, 3, and 4 and compared with various sources, for three  ratios.This present approach results are in excellent agreement both in transverse displacement and transverse shear stresses.In particular for higher values of , the present formulation computes more accurate results than those of Pandya and kant [47].Some other sources such as laminated shell finite elements [48], multiquadrics [21], and third-order formulation presented by Ferreira et al. [44] are compared.In all cases, the lesson learned is that layer-wise approach should be used instead of ESL theories in the analysis of sandwich soft-core plates, where the skin properties are much higher than core properties.The present technique shows that the thickness-stretching approach is important for thicker plates.

Free Vibration Problems of Cross-Ply Laminated Plates.
In this example, all layers of the laminate are assumed to be of the same thickness and density and made of the same linearly elastic composite material.The following material parameters of a layer are used:  1 / 2 = 10, 20, 30, or 40;  12 =  13 = 0.6 2 ;  3 = 0.5 2 ; ] 12 = 0.25.The subscripts 1 and 2 denote the directions normal and transverse to the fiber direction in a lamina, which may be oriented at an angle to the plate axes.The ply angle of each layer is measured from the global -axis to the fiber direction.The example under consideration is a simply supported square plate with a cross-ply lamination scheme [0 ∘ /90 ∘ /90 ∘ /0 ∘ ].The thickness and length of the plate are denoted by ℎ and , respectively.The thickness-to-span ratio ℎ/ = 0.2 is employed in the computation.Table 5 lists the fundamental frequency of the simply supported laminate made of various modulus ratios of  1 / 2 .Figure 2 illustrates the modes of vibration for  1 / 2 = 40 with 15 × 15 grid points.It is found that the present meshless results are in very close agreement with the values of [49] and the meshfree results of Liew et al. [50] based on the FSDT.The small differences may be due to the consideration of the throughthe-thickness deformations in the present formulation.

Conclusions
The first-order and the third-order shear deformation theories are laminate-wise, with laminate degrees of freedom, where all layers have the same rotations.Layer-wise formulations can better accommodate the deformation behaviour of some laminates, in particular for sandwich laminates, where core and skin materials are so different.
In this paper the static and free vibration analysis of composite laminated plates by the use of generalized differential quadrature and using a layer-wise theory with independent rotations in each layer is performed.
The equations of motion were derived and interpolated; moreover, boundary conditions interpolation was schematically formulated.
Finally, composite laminated plate and sandwich plate were considered for testing the present methodology.The obtained results show an excellent accuracy for all cases.

Table 2 :
Square laminated plate under uniform load for  = 5.

Table 3 :
Square laminated plate under uniform load for  = 10.

Table 4 :
Square laminated plate under uniform load for  = 15.
[46]Three-Layer Square Sandwich Plate in Bending, under Uniform Load.A simply supported sandwich square plate, under a uniform transverse load, is considered.This is the classical sandwich example of Srinivas[46].The material