Two multilayered plate models with transverse shear warping functions issued from three dimensional elasticity equations
Abstract
A multilayered plate theory which uses transverse shear warping functions is presented. Two methods to obtain the transverse shear warping functions from threedimensional elasticity equations are proposed. The warping functions are issued from the variations of transverse shear stresses computed at specific points of a simply supported plate. The first method considers an exact 3D solution of the problem. The second method uses the solution provided by the model itself: the transverse shear stresses are computed integrating equilibrium equations. Hence, an iterative process is applied, the model is updated with the new warping functions, and so on. Once the sets of warping functions are obtained, the stiffness and mass matrices of the models are computed. These two models are compared to other models and to analytical solutions for the bending of simply supported plates. Four different laminates and a sandwich plate are considered. Their lengthtothickness ratios vary from 2 to 100. An additional analytical solution that simulates the behavior of laminates under the plane stress hypothesis – shared by all the considered models – is computed. Both presented models give results very close to this exact solution, for all laminates and all lengthtothickness ratios.
drive]A. Loredo drive]A. Castel
Two multilayered plate models with transverse shear warping functions issued from three dimensional elasticity equations
^{a}DRIVE, Université de Bourgogne, 49 rue Mlle Bourgeois, 58027 Nevers, France
Keywords: Plate theory, warping function, laminate, multilayered, composite, sandwich
1 Introduction
In many humanbuilt structures, plates and shells are present. These particular structures are distinguished from others because a dimension – the transverse dimension – is much smaller than the others. Hence, although it is always possible, their representation through a threedimensional domain is not the best way to study them. To understand and predict their mechanical behavior, plate models have been developed. These models permit to study plates and shells through a twodimensional domain while allowing at least membrane and bending deformations. What happens in the third direction is not ignored, it is precisely the purpose of the plate model to integrate the transverse behavior into its equations. The better this behavior is described, the more accurate the model will be. History of plate models probably begins with works of Kirchhoff, Love, and Rayleigh Kirchhoff (1850); Love (1888); Rayleigh (1944) leading to the LoveKirchhoff model in which no shear deformation is allowed. Because of the limitation of this model to thin plates, authors like Reissner, Hencky, Bolle, Uflyand, and MindlinReissner (1945); Hencky (1947); Bolle (1947); Uflyand (1948); Mindlin (1951) have proposed to integrate the shear phenomenon into their formulation. In particular, Reissner made the assumption that shear stresses have a parabolic distribution and consider the normal stress in his model which is derived from a complementary energy. On the other hand, Hencky and Mindlin considered a linear variation of the displacements and and no transverse strain. Mindlin and Reisnner theories are often associated but it is incorrect as demonstrated in reference Wang et al. (2001). The LoveKirchhoff and the Hencky–Mindlin models were first proposed for homogeneous plates but their pendent for laminated structures have been later presented and are called the Classical Laminated Plate Theory (CLPT) and the Firstorder Shear Deformation laminated plate Theory (FoSDT). The last one has been improved by the use of shear correction factors Whitney (1973); Noor and Scott Burton (1989); Pai (1995). The need to improve this accuracy has been motivated by the study of thick plates and laminated plates. In both cases, the early plate models fail to give accurate results. It is even worse if some layers have low mechanical properties compared to others, which is the case for sandwich structures or when viscoelastic layers are used to improve the damping. To overcome these difficulties, specific models have been proposed. However, a universal model that can manage plates of various lengthtothickness ratios, with any lamination scheme and various materials including functionally graded materials, with good accuracy for the static, dynamic, and damped dynamic behaviors is still a challenge.
Dealing with multilayered plates has given rise to another class of theories, called the LayerWise (LW) models, in which the number of unknowns depends on the number of layers, by opposition to the Equivalent Single Layer (ESL) family of models in which the number of unknowns is independent of the number of layers. Obviously, LW models are expected to be more accurate than ESL models but they are less easy to implement for complex structures and require more computational resources than ESL models. This has motivated researchers to propose ESL models in which the transverse shear is taken into account in a more precise way than for the first models. With the help of the classification given by Carrera Carrera (2002) among other review papers, we can distinguish different approaches that have given interesting ESL models:

Higherorder (and also nonpolynomial) shear deformation theories: the Vlasov–Levinson–Reddy Vlasov (1957); Levinson (1980); Reddy (1984) theory, also called Third order Shear Deformation Theory (ToSDT), among other third order theories, proposes a kinematic field with a thirdorder polynomial dependence on , motivated by the respect of the nullity of transverse shear at top and bottom faces of the plate. Higher order (order greater than 3) theories have been proposed. Nonpolynomial functions have also been used to take into account the shear phenomenon: for example, Touratier in reference Touratier (1991), uses a function, Soldatos Soldatos (1992) uses a function, and Thai & al. Thai et al. (2014) uses a function to integrate the transverse shear in the kinematic field. Note that, for the study of laminated plates, these models do not integrate information about the lamination sequence in their kinematic field, and do not verify the transverse stress continuity at interfaces.

Interlaminar continuous ESL models: some a priori LW models reduce to ESL models with the help of assumptions between the fields in each layer. ZigZag (ZZ) models enter in this category. Pioneering works of Lekhnitskii Lekhnitskii (1935) and Ambartsumyan Ambartsumyan (1958) have been classified as such by Carrera Carrera (2002) who also shows that other authors have integrated the multilayer structure in their model Whitney (1969); Sun and Whitney (1973); Cho and Parmerter (1993) in a very similar manner. The main idea of ZZ models is to let inplane displacements vary with according to the superposition of a zigzag law to a global law – cubic for example. With these models, shear stresses can satisfy both continuity at interfaces and null (or prescribed) values at top and bottom faces of the plate. An enhancement of the second model presented in reference Sun and Whitney (1973) for general lamination sequences is proposed in references Woodcock (2008); Loredo and Castel (2013). Models of references Loredo and Castel (2013); Cho and Parmerter (1993) are used for comparison in this study, they are presented at section 5.2. Other works can be cited. For example, Pai Pai (1995) presents a method to determine four warping functions of cubic order in each layer, and also a generic method to compute shear correction factors. In ref. Arya et al. (2002), the author proposes a warping function for beam problems which is a linear ZZ function superimposed to an overall sine function. In ref. Kim and Cho (2006), the authors develop an enhancement of the model given in ref. Cho and Parmerter (1993) for general lamination sequences. Note that references Loredo and Castel (2013); Cho and Parmerter (1993); Pai (1995); Kim and Cho (2006) propose models with four WF, which is necessary for the study of general lamination sequences.

Direct approaches and micropolar plate theories: based on the works of the Cosserat brothers Cosserat and Cosserat (1909) on generalized continuum mechanics and on deformable lines and surfaces. The micropolar theory of elasticity, which considers (microscopic) continuous rotations and couples in addition to the classical displacements and forces, has been applied to rods and shells by Ericksen Ericksen and Truesdell (1957). Later, it has been used in numerous works, including plate models, as it can be seen in the review paper Altenbach et al. (2010).

Sandwich theories: sandwich panels have been early considered apart from laminated composite because of their particular structure. They are made of two face layers (themselves made of one or several layers) and a core layer. The thickness of the faces is small compared to the thickness of the core (typically 10 times less) and Young modulus of the faces is high compared to Young modulus of the core (typically more than 1000 times higher). For these reasons, specific models have been proposed, for example by Reissner Reissner (1947), followed by many others, as it can be seen in the review papers Noor et al. (1996); Kreja (2011).
Additional references can be found in recent reviews on the subject Wanji and Zhen (2008); Kreja (2011); Khandan et al. (2012).
When theories are issued from assumptions on displacement or stress fields, there is no guarantee for these theories to be consistent. This means that, with respect to the required order of the considered theory, some terms of the threedimensional strain energy may appear with an erroneous coefficient, or may even be omitted. The consistency of plate theories has been discussed by many authors, one can for example refer to recent works Schneider and Kienzler (2011); Schneider et al. (2014). Precise rules for the consistency are formulated for homogeneous materials, but it seems difficult to extrapolate them for inhomogeneous materials like laminates and sandwiches, especially when advanced kinematic and/or static assumptions are made. Although this subject is not treated in this study, it can be a significant aspect as it is mentioned in the conclusion.
In this work, the kinematic assumptions are similar of those taken in Ref. Kim and Cho (2006) but they differ because the choice of the functions describing the transverse behavior is not made a priori. These functions, called warping functions (WF) are the core of the model. The model has been entirely formulated in Loredo and Castel (2013) and it has been shown that, according to specific choices of the WF, it can also represent most of the classical models (CLPT, FoSDT, ToSDT) and others, as it will appear in the following. However, and it is precisely the subject of this article, it is possible to choose and adapt the WF in a completely free manner.
This article presents two different ways of obtaining new sets of WF issued from three dimensional elasticity laws. The first way consists in building the WF from threedimensional solutions. Threedimensional solutions for the bending of laminates have been first obtained by Pagano and by Srinivas & al. Pagano (1970); Srinivas et al. (1970) for crossply laminates, by Noor & Burton Noor and Burton (1990) for antisymmetric angleply laminates, and have been recently obtained for general lamination schemes Loredo (2014). These solutions are achieved for particular boundary conditions and load, which can be seen as the principal limitation of their use in the present model. The second way is to derive the WF from the equilibrium equations. This leads to an iterative model: starting with “classical” WF, for example Reddy’s formula, WF are issued from the equilibrium equations and then integrated to the model, and so on until no significant changes are detected.
2 Considered plate theory
In this section, we recall the main components of the theory presented in Ref. Loredo and Castel (2013). It is a plate theory based on the use of transverse shear WF. It allows the simulation of multilayer laminates made of orthotropic plies using different sets of transverse shear WF. The purpose of this paper is to propose enhanced WF issued from 3D elasticity equations (see section 3), but several plate models (among ESL and ZZ models) issued from the literature can also be formulated in terms of transverse shear WF, as it was shown in Ref. Loredo and Castel (2013). This is of practical interest when comparing results issued from different models, because these models can be implemented in a similar manner.
2.1 Laminate definition and index convention
The laminate, of height , is composed of layers. All the quantities will be related to those of the middle plane^{1}^{1}1The reference plane can be arbitrarily chosen in the laminate assuming that the corresponding transverse shear WF are adapted in consequence. which is placed at ; they are marked with the superscript . In the following, Greek subscripts take values or and Latin subscripts take values , or . The Einstein’s summation convention is used for subscripts only. The comma used as a subscript index means the partial derivative with respect to the following index(ices).
2.2 Displacement field
The kinematic assumptions of the present theory are:
(1) 
where , and , are respectively the inplane displacements, the deflection and the engineering transverse shear strains evaluated at the middle plane. The are the four WF. The associated strain field is derived from equation (1):
(2a)  
(2b)  
(2c) 
which, with the use of Hooke’s law, leads to the following stress field:
(3a)  
(3b)  
(3c) 
where are the generalized plane stress stiffnesses and are the components of Hooke’s tensor corresponding to the transverse shear stiffnesses.
2.3 Static laminate behavior
The model requires introduction of generalized forces:
(4a)  
(4b) 
They are then set, by type, into vectors:
(5) 
and the same is done for the corresponding generalized strains:
(6) 
Generalized forces are linked with the generalized strains by the and following stiffness matrices:
(7) 
with the following definitions:
(8a)  
(8b) 
2.4 Laminate equations of motion
Weighted integration of equilibrium equations leads to
(9a)  
(9b)  
(9c) 
where is the value of the transverse loading (no tangential forces are applied on the top and bottom of the plate), and:
(10) 
3 Warping functions issued from transverse shear stress analysis
We introduce two different ways to obtain WF from transverse shear stress analysis: a first set of WF is issued from an analytical solution and a second set is issued from an iterative process using the integration of equilibrium equations. First, we shall examine a way to link the WF to the transverse shear stresses.
3.1 From transverse shear stresses to WF
Considering equation (3b), we see that the are directly linked to the . Introducing the transverse shear stresses at into this equation leads to
(11) 
where are components of the compliance tensor.
This can be written
(12) 
where:
(13) 
The cannot be directly issued from equation (12) because there are four functions to determine from the variations of two transverse shear stresses, leading to infinitely many solutions. The main idea is to consider a simply supported plate subjected to a doublesine load, and to issue the four functions from the variations of the transverse shear stresses at two separate points of the plate. Since the deformation of the plate is of the form (19), the transverse shear stresses at the reference plane are of the form:
(14a)  
(14b) 
These shear stresses are evaluated at points and :

and implies and ,

and implies and .
Setting these local values into formula (12) leads to the following system:
(15) 
The are obtained from the resolution of this system; are then obtained using the reciprocal of equation (13):
(16) 
Then, integrating the so that gives the four WF .
3.2 WF issued from exact 3D solutions
Exact 3D solutions of simply supported plates subjected to a doublesine load are known for crossply and antisymmetric angleply laminates since the works of Pagano Pagano (1970) and Noor Noor and Burton (1990). Further works have shown that they can be obtained by several ways. Solution for general lamination have been proposed recently in Ref. Loredo (2014). The corresponding plate problem is solved with the appropriate method, and the transverse shear stresses are computed at points A and B. Then, the procedure described in the previous section is applied.
3.3 WF issued from iterative equilibrium equation integration
Since transverse shear stresses can be obtained from the equilibrium equations, it is also possible to get the warping functions following an iterative process^{2}^{2}2This iterative process, although based on a simpler formulation, has been proposed in the 1989 unpublished reference Loredo (1989). The process, described in the algorithm 1, starts with any known WF, says Reddy’s formula for example, and, at each iteration, the model is updated with WF issued from the transverse stresses of the previous iteration. The procedure is also based on the simply supported plate bending problem with a doublesine loading. Let us establish the needed formulas, starting from the equilibrium conditions within a solid, without body forces:
(17a)  
(17b) 
The transverse shear stresses, for the static case, are therefore computed using:
(18) 
Then, the spatial derivatives of the generalized displacements are eliminated accounting to the specific nature of the chosen functions in formulas (19).
3.4 Discussion
The two ways to obtain WF from transverse shear stresses described in the above sections are based on the study of a simply supported plate subjected to a doublesine load. As one shall see in the following, the two methods give very similar results. This proves that the model, which is strongly implicated in the iterative process, is able to fit the transverse shear stresses of the 3D solution with good agreement. Of course, there is no guarantee at this time that these WF will be the best candidates if another plate problem is studied, with different boundary conditions and/or different loading. It is precisely the reason why the iterative process is interesting because it might be adapted to a local strategy.
4 Solving method by a Navierlike procedure
A Navierlike procedure is implemented to solve both static and dynamic problems for a simply supported plate. For the static case, in order to respect the simply supported boundary condition for laminates which are not of crossply nor antisymmetrical angleply types, a specific loading is applied. This is done with the help of a Lagrange multiplier, as explained below. The dynamic study is restricted to the search of the natural frequencies.
The Fourier series is limited to one term, hence the generalized displacement field is
(19) 
with
where and are the length of the sides of the plate, and are wavenumbers, set to for static analysis or to arbitrary values for the dynamic study of the corresponding mode. Then for a given and , the motion equations of section 2.4 give a stiffness and a mass matrix, respectively and , related to the vector . The static case is treated by solving the linear system , where is a force vector containing for its third component (generally set to one). Solving the dynamic case consists in searching the generalized eigenvalues for matrices and . For crossply and antisymmetric angle ply, respects the simply supported conditions, i. e. . For the general laminates, the deflection under a doublesine loading gives a . Since we choose to keep simply supported boundary conditions, may be set to zero if a bicosine term is added to the loading. The amplitude of the bicosine term is obtained using a Lagrange multiplier. The stiffness matrix is then of size .
with being a vector with a one on it’s eighth component. For general laminates, this augmented system remains regular but no longer positive definite, hence the solving procedure must be chosen adequately. This is not detrimental because the system is of very small size. For the dynamic case, the matrix is augmented with a line and a column of zeros so it becomes a matrix. Detailed formulation of matrices and are given in Appendix A.1. Note that it is also possible to keep the loading of the form , then the simply supported boundary condition is no longer respected for general laminates.
5 Reference models
The model using WF issued from transverse shear stresses of analytical solutions, presented in section 3.2, will be denoted 3DWF. The model using WF obtained by the iterative process, presented in section 3.3, will be denoted IteWF. The results obtained with 3DWF and IteWF models are compared to those obtained with different models issued from the literature and with exact analytical solutions. These reference models are presented below.
5.1 Exact analytical solutions
5.1.1 The Exa solution
Each studied case is solved by a statespace algorithm described in reference Loredo (2014). The algorithm is a generalization for general lamination sequences of existing algorithms for crossply and antisymmetric angleply lamination sequences. In the reference Loredo (2014), the algorithm has been tested: exact solutions of previous works have been reproduced, and for general lamination sequences, the solution has been confronted with success to an accurate finite element computation. In this work, the algorithm is used to obtain deflections, stresses and natural frequencies taken as reference for comparisons, but it is also used to create a set of WF for the 3DWF model, as explained in section 3. The corresponding solution is denoted Exa in the following text and in tables. For all the studied configurations, the loading is divided into two equal parts which are applied to the top and bottom faces.
5.1.2 The Exa solution
All models involved in this study consider the generalized plane stress assumption which leads to the use of the reduced stiffnesses . Further, these models do not consider the normal strain in their formulation. It is well known that these assumptions are correct when the lengthtothickness ratio is high and when layers are made of similar materials. This hypothesis is no longer valid for low lengthtothickness ratio (say 2 and 4) and when materials with very different behavior are used, which is typical of sandwiches. Hence, it is not possible to draw a conclusion from a comparison of models with an analytical solution they cannot approximate. For this reason, another exact solution is computed for a virtual laminate which has modified stiffness values in order to, firstly, fit the generalized plane stress assumption and secondly, force . It is done by setting and for each material. This is equivalent to the replacement of the by the and forces to take very small values. This exact solution, designated by the Exa symbol, is computed with the same algorithm than the Exa solution. For this case too, the loading is divided into two equal parts which are applied to the top and bottom faces.
5.2 Other models
Some more or less classical models have been chosen for comparison matters. They offer the advantage to be easily simulated by the present model when appropriate sets of WF are selected:

Firstorder Shear Deformation Theory (denoted FoSDT in tables and figures): Often called Mindlin plate theory, it can be formulated setting in the generic model the following warping functions,
(20) where is Kronecker’s delta. Note that this theory is generally used with shear correction factors, often the factor which corresponds to an homogeneous plate. As there are several ways to compute these correction factors in the general case, we chose in this study not to use them. This is of course a serious penalty for this model.

Thirdorder Shear Deformation Theory (denoted ToSDT): often called Reddy’s third order theory, verifies that transverse shear stresses are null at the top and bottom faces of the plate. It is simulated using the following warping functions:
(21) 
Firstorder ZigZag model with 4 WF (denoted FoZZ4): This model verifies the continuity of transverse shear stresses at the layers’ interfaces. It was first presented by Sun & Whitney Sun and Whitney (1973) for crossply laminates and generalized for general lamination sequences by Woodcock Woodcock (2008). It is possible to formulate this model with the following WF as shown in Loredo and Castel (2013):
(22) 
Thirdorder ZigZag model with 4 WF: (denoted ToZZ4): This formulation, presented by Cho^{3}^{3}3In reference Carrera (2002), Carrera wrote that this model was a rediscovery of previous works, and called the model the Ambartsumyan–Whitney–Rath–Das theory. Cho and Parmerter (1993); Kim and Cho (2007) consists in superimposing a cubic displacement field, which permits the transverse shear stresses to be null at the top and bottom faces of the laminate, to a zigzag displacement field issued from the continuity of the transverse shear stresses at layer interfaces. The corresponding WF, which are polynomials of third order in , are not detailed. Indeed, as their computation involves the resolution of a system of equations, it is difficult to give here an explicit form. Note also that for coupled laminates, an extension of this model Kim and Cho (2006) has been proposed. This extension has not been implemented in this study, it could have given different results for the studied case (which is the only coupled laminate considered in this study).
As these models can be implemented by means of WF, the solving process described in section 4 can be applied to all models.
6 Numerical results
This section proposes the study of five configurations including four laminates and a sandwich plate. Only two materials are involved: an orthotropic composite material used in all laminates and an honeycombtype material used for the core of the sandwich panel. All the material properties are given in table 1.
Composite ply (p)  

Core material (c) 
Three nondimensionalized quantities are considered and compared to those issued from analytical solutions:

Deflections are nondimensionalized using equation
(23) 
First natural frequencies are nondimensionalized using equation
(24) 
Shear stresses are nondimensionalized using equation
(25)
where and are taken as values of the core material for the sandwich and as values of the composite ply for laminates.
6.1 Rectangular crossply composite plate
This plate is made of three composite plies of equal thickness, with a stacking sequence and . Results are presented in table 2. Note that Reddy’s model (ToSDT) – with relatively simple WF – gives quite good results for this laminate. Cho’s model (ToZZ4) – with more sophisticated WF – gives better values than the previous, except for the lengthtothickness ratio. Results also show that the two proposed models, 3DWF and IteWF, give a very satisfying accuracy for all lengthtothickness ratios. Compared to other models with reference to the exact solution, values for the deflection are among the best results, values for the transverse stresses at points A and B, and for the first natural frequency, are the best. However, we can note a weakness for the prediction of , which is probably due to a sandwichlike behavior of this laminate according to the direction. The FoZZ4 model, which is accurate with sandwiches, tends to confirm this point. Both proposed models, even though warping functions are generated with two very different methods, give almost identical results. Compared to the Exa plane stress exact solution, 3DWF and IteWF models give the best results. The model IteWF gives, for this problem, the same result as the exact solution. This point will be discussed later on section 6.6.
Figure 1 shows the corresponding WF for , for all plate models except the FoSDT one. Figure 2 presents the variations of the transverse shear stresses obtained by integration of the equilibrium equations for all models, compared to the exact solution, in the case.
a/h  Model  %  %  %  %  %  %  %  %  

2  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref.  
4  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref.  
10  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref.  
100  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref. 
6.2 Square sandwich plate
In order to study the behavior of a structure exhibiting a high variation of stiffness through the thickness, we propose to study a square sandwich plate with ply thicknesses and . The face sheets are made of one ply of unidirectional composite and the core is constituted of a honeycombtype material. Material properties are presented in table 1. Results presented in table 3 show that, this time, Reddy’s model (ToSDT) is not as accurate as in the previous case, although Cho’s model (ToZZ4) obtains better values. This is due to the particular nature of sandwich materials which gives typically zigzag variations for displacements through the thickness and then typically zigzag WF. Cho’s model (ToZZ4) is able to fit this kind of variation although Reddy’s (ToSDT) model is not. The Sun & Withney model (FoZZ4) has also been proved to be very efficient for sandwiches, which can be verified in this table. Note that the two proposed models, 3DWF and IteWF, globally give the best results. Comparison with the Exa model is discussed in section 6.6.
Figure 3 shows the corresponding WF for , for all plate models except the FoSDT one. Figure 4 presents the variations of the transverse shear stresses obtained by integration of the equilibrium equations for all models, compared to the exact solution, in the case.
a/h  Model  %  %  %  %  %  %  %  %  

2  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref.  
4  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref.  
10  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref.  
100  FoSDT  
ToSDT  
FoZZ4  
ToZZ4  
3DWF  
IteWF  
Exa  ref.  ref.  ref.  ref.  
Exa  ref.  ref.  ref.  ref. 
6.3 Square antisymmetric angleply composite plate
As the two previous cases enter in the crossply family, let us now consider an antisymmetric angleply square plate with two layers of equal thickness and a