Free-Form Laminated Doubly-Curved Shells and Panels of Revolution Resting on Winkler-Pasternak Elastic Foundations : A 2-D GDQ Solution for Static and Free Vibration Analysis

This work presents the static and dynamic analyses of laminated doubly-curved shells and panels of revolution resting on Winkler-Pasternak elastic foundations using the Generalized Differential Quadrature (GDQ) method. The analyses are worked out considering the First-order Shear Deformation Theory (FSDT) for the above mentioned moderately thick structural elements. The effect of the shell curvatures is included from the beginning of the theory formulation in the kinematic model. The solutions are given in terms of generalized displacement components of points lying on the middle surface of the shell. Simple Rational Bézier curves are used to define the meridian curve of the revolution structures. The discretization of the system by means of the GDQ technique leads to a standard linear problem for the static analysis and to a standard linear eigenvalue problem for the dynamic analysis. Comparisons between the present formulation and the Reissner-Mindlin theory are presented. Furthermore, GDQ results are compared with those obtained by using commercial programs. Very good agreement is observed. Finally, new results are presented in order to investigate the effects of the Winkler modulus, the Pasternak modulus and the inertia of the elastic foundation on the behavior of laminated shells of revolution.


Introduction
During the last sixty years, two-dimensional linear theories of thin shells have been developed including important contributions by Timoshenko and Woinowsky-Krieger [1], Flügge [2], Gol'denveizer [3], Novozhilov [4], Vlasov [5], Ambartusumyan [6], Kraus [7], Leissa [8,9], Markuš [10], Ventsel and Krauthammer [11] and Soedel [12].All these contributions are based on the Kirchhoff-Love assumptions.The transverse shear deformation has been incorporated into shell theories by following the theory of Reissner-Mindlin [13], also named First-order Shear Deformation Theory (FSDT).Abandoning the assumption related to the preservation of the normals to the shell middle surface after the deformation, a comprehensive analysis for elastic isotropic shells was made by Kraus [7], Gould [14,15] and Qatu [16,17].The present work is just based on the FSDT.In order to include the effect of the initial curvature in the evaluation of the stress resultants a generalization of the Reissner-Mindlin (RM) theory has been proposed in literature by Kraus [7], Qatu [16,17] and Toorani and Lakis [18,19].There are three different ways to evaluate the engineering elastic constants in the study of curved shells.The first is the Reissner-Mindlin approach [7] that consists in neglecting the effect of curvatures.Using this approach the engineering elastic stiffnesses are constant and do not depend on curvatures.The second one, proposed by Kraus [7] and Toorani and Lakis [18], is based on the Taylor expansion, while the third one proposed by Qatu [16] consists in the exact integration of the elastic constants.As a consequence of the use of these considerations, the stress resultants directly depend on the geometry of the structure in terms of the curvature coefficients.In this latter case, the hypothesis of the symmetry of the in-plane shearing force resultants and the torsional couples de-clines.A further improvement of the previous theories of shells has been proposed by Toorani and Lakis [19].In the present work their kinematic model is used in order to include the effect of the curvature from the beginning of the shell formulation.In this way, the strain relationships have to change and, as a consequence, the equilibrium equations in terms of displacements have to be modified.In the present paper, the proposed shell theory, named General Shell Theory (GST), is considered and compared with the Reissner-Mindlin (RM) theory.Comparisons between these two different formulations are presented in this paper.Several studies dealing with the shells theory have been presented years before.The most popular numerical tool used to perform the static and dynamic analyses is currently the finite element method [14,15,20].The generalized collocation method based on the ring element method has also been applied.In this method, each static and kinematic variable is transformed into a theoretically infinite Fourier series of harmonic components, with respect to the circumferential co-ordinate [21,22].In other words, when dealing with a completely closed shell, the 2D problem can be reduced using standard Fourier decomposition.For a panel, however, it is not possible to perform such a reduction operation, and the two dimensional field must be directly dealt, as it will just be done in the present work.Furthermore, the system of second-order linear partial differential equations is solved, without resorting to the one-dimensional formulation of the equilibrium of the shell.Complete revolution shells are obtained as special cases of shell panels by satisfying the kinematical and physical compatibility at the common meridian with 0, 2π   . The excellent mathematical and computational algorithmic properties, combined with successful industrial applications, have contributed to the enormous popularity of the Rational Bézier and Non-Uniform Rational B-Splines (NURBS) curves [23][24][25].These curves allow to generalize the shape of the shell meridian and can be used for the optimization of the structure itself.By introducing the Differential Quadrature rule [26] and the simple mathematical formulation of the Rational Bézier and NURBS curves [23][24][25], it is possible to numerically evaluate the geometric parameters of a free-form shell of revolution.For the sake of simplicity and without loosing generality, only Rational Bézier curves are used in this study.Due to the increasing importance of the interaction of shells with the elastic medium, the Winkler-Pasternak foundation is introduced.Differently from papers presented in literature [27][28][29][30][31], all the effects of the foundation, except the damping, are separately considered.New results are presented in order to investigate the effects of the Winkler modulus, the Pasternak modulus and the inertia of the elastic foundation on the behavior of laminated shells of revolution.The mathematical fundamen-tals and recent developments of the GDQ method as well as its major applications in engineering are discussed in detail in the book by Shu [26].The interest of researches in this procedure is increasing due to its great simplicity and versatility.As shown in the literature [32], GDQ technique is a global method which can obtain very accurate numerical results by using a considerably small number of grid points.Therefore, this simple direct procedure has been applied in a large number of cases  to circumvent the difficulties of programming complex algorithms for the computer, as well as to reduce the computational time.In conclusion, the aim of the present paper is to demonstrate an efficient and accurate application of the Differential Quadrature approach, by solving the equations governing the static and the free vibration of laminated composite doubly-curved moderately thick shells and panels of revolution.Summarizing, this research deals with four aspects.The first is the improvement of the Reissner-Mindlin Theory using a different kinematical model.In this way the effect of the curvature of the shell structure is considered from the beginning of the theory derivation.The second is the generalization of the shape of the shell meridian.The Differential Quadrature rule is used to evaluate the geometric parameters needed to describe the geometry of the structure when a Rational Bézier meridian curve is assumed.The third is the investigation of the effects of Winkler-Pasternak foundations on the behavior of the shell structures in static and dynamic analyses.All the effects of the foundation are separately considered.The fourth is the use of the Generalized Differential Quadrature method to solve the governing shell equations.

Shell Fundamental Equations
The basic configuration of the problem herein considered is a laminated composite doubly-curved shell [83] as shown in Figure 1.The co-ordinates along the meridian and circumferential directions of the reference surface are  and s , respectively.The distance of each point from the shell mid-surface along the normal is  .It is considered a laminated composite shell made of laminae or plies, where the total thickness of the shell is defined as: in which is the thickness of the k-th lamina or ply.In this work, doubly-curved shells of revolution are considered.For this type of structures the analytical expressions of the meridian curve are reported in the work by Tornabene [78], so that no further considerations will be introduced.The angle formed by the extended normal n to the reference surface and the axis of rotation x , or the geometric axis 3 x of the meridian curve, is defined as the meridian angle  ; the angle between the radius of the parallel circle x axis is designated as the circumferential angle  , as shown in Figure 2.
A simple way to define a general meridian curve is to use the well-known Rational Bézier representation of a plane curve [24,25,80].In particular, it is possible to describe a Rational Bézier curve in the following manner:   s For these structures the parametric co-ordinates ,  define, respectively, the meridian curves and the parallel circles upon the middle surface of the shell.The curvilinear abscissa   s  of a generic parallel is related to the circumferential angle  by the relation

 
In this way, only the co-ordinates of the curve , are known in the co-ordinate reference system 1 3 .In order to solve the shell problem, it is important to express the horizontal radius   of a generic parallel and the radii of curvature  in the meridian and circumferential di- rections as functions of  .Based on the differential geometry [7,12,57,80,83], the radius of curvature of the meridian curve can be described as a function of 3 x using the following expression: It is worth noting that the derivatives of the meridian curve are not known a priori, so that a numeric method to evaluate the first and second derivatives of the meridian curve is required.The differential quadrature rule allows to approximate these derivatives using the following definition [26]: where   n ij  are the weighting coefficients of the n-th order derivative.By discretizing the domain using the Chebyshev-Gauss-Lobatto (C-G-L) grid distribution: and interpolating the 1 x  co-ordinates of the curve points derived by the Equations (2) using the previous calculated points (6), the general curve can be represented by the new co-ordinates points   ˆ, i i x x , for .Applying the differential quadrature definition (5), the expression (4) assumes the following discrete aspect: where As a results of the differential geometry [7,12,57,80,83], it is possible to introduce the following relation: By using the differential quadrature definition (5), the relation ( 7) can be expressed in the discrete form: By discretizing the domain and interpolating the 1 x and 3 x co-ordinates of the curve points using the calculated points (9), the general curve can be represented by the following new co-ordinate points   x x   , for .Thus, all the discrete points of the curve are determined in terms of the co-ordinates x x   and the angle i  .In the Figure 3, a Rational Bézier curve, its control points and the curve co-ordinates   x x   , evaluated as above exposed, are represented [80].The vectors of the control points and the weights used in Figure 3 are the following: 0.2 0.7 1.2 1.4 1.4 1.2 0 0.2 0.6 1 1.5 2 (10) Based on the previous considerations, the horizontal radius   0 R  of a shell of revolution assumes the following discrete form: , for 1, 2, , For doubly-curved revolution shells the Gaus da s-Cozzi relation can be expressed as follows: By using the differential quadrature definition (5), it is po  ssible to determine the radius of curvature   R   in meridian direction and its first and second der in discrete form: ,for 1, 2, , cos , its control points   and curve evaluated discrete points   . Figure 3.A Rational Bézier curve Finally, as a results of the differential geometry [7,12,57,80,83], the radius of curvature   s R  in circumferent a discrete form as follows: ial direction for a shell of revolutio be expressed n can in   0 , for , sin It is worth noting that, follo 1, 2, wing the previous considerations, all the useful geometric parameters the surface of revolution under consideration in discrete form (11)- (16).As shown, the differential qu sumed; 6) the rotary inertia and the initial curvature are also taken into account.Consistent with th assumptions f a mod above, the displacement field considered in this study follows the describing are known adrature rule (5) has been used to approximate the derivatives needed for the definition of the geometry of a shell of revolution.As concerns the shell theory, the present work is based on the following assumptions: 1) the transverse normal is inextensible so that the normal strain is equal to zero: 2) the transverse shear deformation is considered to influence the governing equations so that normal lines to the reference surface of the shell before deformation remain straight, but not necessarily norma relaxed Kirchhoff-Love hypothesis); 3) the shell deflections are small and the strains are infinitesimal; 4) the shell is moderately thick, therefore it is possible to assume that the thickness direction normal stress is negligible so that the in-plane assumption can be invoked:


he in-plane displacements U  and s U vary linearly through the thickness W remains independent of , while  .It should be also remarked that, differently from the previous works by Tornaben 65,77,7 , the displacement field has been improved taking into e [ 8] account the real geometry of the shell and in particular the curvature effect has been directly introduced into the kinematical model as proposed by Toorani and Lakis [19].
Due to the change of the kinematical model the relationships between strains and displacements along the shell reference surface   that are different from those presented in previous papers [77,78].In the above Equations ( 19), the first four strains     are the in-plane meridian, circumferential and shearing components, and 0 0 0 0 , , ,     are the analogous curvature changes.The last two components 0 0 , n sn    are the transverse shearing stra ition assumed in the following i composite linear elastic material.Accordingly, the following constitutive equations relate internal stress resul-ins.The shell s a laminated compos tants and internal couples with generalized strain components (18) on the middle surface: where the elastic engineering stiffnesses   q ijm A which depend on curvatures are defined as follows (see Appendix for more details): Several approaches can be found in literature to evaluate the engineering elastic constants   q ijm A [7,[16][17][18].It is worth noting that due to the elastic engineering stiffnesses fact that the   q ijm A rivatives respec depend on curvatures, the corresponding de t to the co-ordinates along the meridian  and circumferential s directions of the reference surface have to be evaluated.In order to erform this operation, the Differentia uadrature rule p l Q [26] is used.Thus, the derivatives of the elastic engineering stiffnesses   q ijm A are numerically eval d.
uate  is e shear corre th ction factor, which is usually taken equal to 5 6   , such as of in the present l the determination shear correction factor or com laminated structures is still an unresolved issue, because these factors depend on various parameters [18].In Equa-tions (19) where is the orientation angle of the principal material co-ordinate system Ô s    inate co-ordi of the k-th orthotropic ply with respect to the lam nate system Furthermore, the elastic constants in the terial co-ordinate system are e pressed as x follows: where 1 2 13 23   (23) where: mass inertias and is the mass density of the material per unit of volume of the k-th ply, while , ell can [12,79,83] the , circumferential s and normal  directions, while the last two are rotational equi- where , , , , , he eff [12,79,83].Furthermore, as additional hypothesis, t ect of the damping of the foundation is neglected.The three basic set of equations, namely the kinematic (18), constitutive (19) and motion (23) equations may be combined to give the fundamental system of equations, also known as the governing system of equations.By replacing the kinematic equations (18) into the constitutive equations (19) and the result of this substitution into the motion equations ( 23), the complete equations of motion in terms of displacements can be written as: where are the equilibrium operators and the tias are defined as follows: , , It is worth noting that, differently from previous works by Tornabene [77,78], the mass matrix and the equilibrium operators , introduced in Equation ( 27), have changed due to hoice of using the new kinematical model (18).Fu ore, the second derivative of the principal radius Soft simply supported edge boundary conditions (S) Free edge boundary conditions (F) 0 0 1 0 sin 0 at or , 0 where Physical compatibility conditions along the closing where 0 1 ˆ,0, ,0, In an analogous way, in order to consider a toroidal shell of revolution, it is necessary to implement the kinematic and physical compatibility conditions between the two computational parallels with 0 0   and with 1 2π   : Kinematic compatibility conditions along the closing parallel   0, 2π Ph mpatib itions along the closing paral

Num cal Implementation
The Differential Quadrature method used tize the derivatives in the governing equations s of displacements as well as boundary and co ence s will be mpatibility conditions.Since a review of the GDQ Method is presented in Tornabene [65,83], the same approach is used in the present work about the GDQ technique.The Chebyshev-Gauss-Lobatto (C-G-L) grid distribution is assumed.Since the co-ordinates of the grid points of the refer urface in the  direction are introduced in Equation ( 9), then the co-ordinates of the grid points in the s di ction are the following: re for 1, 2, , , for 0, with where M is the total number of sampling discretize the domain in points used to s direction of the doublycurv It has been proven that interpo polynomials, the Chebyshe samp guarantees convergence and effi-GDQ technique [57][58][59]6 nertias are set to zero, the GDQ pr m (26) in into a ed sum of node values of independent variables.Each approximated equation is valid sampling point.Thus, the whole system of differential equations has been discretized and the global assembling leads to the following set of linear alg ed shell.
In the above mentioned matrices and vectors, the partitioning is set forth by subscripts and ref to the system degrees of freedom and standin for boundary and domain, respectively.In this sense, -equations represent the discrete bound y cond which are valid only for the points lying on const of the shell; while -equations are the equ equations, assigned to the interior nodes.In order to e comp ore e The deflection of the structures considered can be determined by solving the linear algebraic problem (40).In particular, the solution procedure by means of the GDQ technique has been implemented in a MATLAB code.Otherwise, when the external forces , , , , , are set to zero, the free vibration of laminated composite doubly-curved shells and panels of revolution can be studied.Using the d of variabl on, i metho e bl s separati t is possie to seek solutions that are harmonic in time and whose frequency i 2π f   .The displace written as follows: where the vibration spatial amplitude values , , ,  the wh fulfill the fundamental differential system ole system of differential equations cretized and the global assembling leads to the following linear eigenvalue problem: . Thus, has been dis- In order to make the computation more efficient, kinematic condensation of non-domain degrees of freedom is performed: The natural frequencies of the structure considered can be determined by solving the standard eigenvalue problem (43).In particular, the solution procedure by means GDQ been implemented in a MA-TL ned using t e eigs function of MATLAB software.More details regarding the way to obtained the Equations ( 40) and ( 43) can be found in the previous works [57][58][59]63,[65][66][67][68]83].It is worth noting tha the present approach, differing from the finite element method, no integration occurs prior to the global assembl putational cost saving in favor of the Differential Quadrature technique.

Numerical Results
In the present paragraph, some results and considerations about the static analysis and the free vibration problem of laminated composite doubly-curved shells and panels revolution are presented.The analysis has been c out by means of numerical procedures illustrated above.One of the aims of this paper is to compare results obtained through the GDQ method with the ones obtained through finite element techniques.In order to verify the accuracy of the present method, some comparisons and tests have been performed.Extensive attempts to validate de for the isotropic in the Ph.D.   the numerical procedure have been ma and anisotropic cases and can be found sis by Tornabene [57] and in the book by Tornabene [83].
In this work, the static deflection and the frequency parameters evaluated by the present formulation are in good agreement with the results obtained with the finite ele-e ur g the elas rical and material properties are the same reported in int Table 1 2. First ten frequencies for a SSSS spherical panel resting on elastic foundation.
Foundation properties:   0 90 0 Table 3.First ten frequencies for a SSSS spherical panel resting on elastic foundation.
Foundation properties:  T ta FC laminated composi o stic foun dation acting on the bottom surface.Also for these cases, the effect f the initial curvature is considered by com ral Shell Theory (G T) and the ry.T us, as r arding the influence of the initial curvature, it is worth noting that the difference between the two theories considered is low structures are reported the first six mode shapes for some of the structures considered above.

Conclusion
The atic v an f l d -curved shells and panels of revolution resting on have been preusing e GDQ m thod.A dation, except the damping, are separately introduced.problem for the static analysis and a standard linear eigenvalue problem for the dynamic analysis.Numerical solutions have been compared with the ones obtained using commercial programs.The comparisons conducted with FEM codes confirm how the GDQ simple numerical method provides accurate and computationally low cost results for all the structures considered.Furthermore, discretizing and programming procedures are quite easy.The GDQ results show to be precise and reliable.The numerical tests dem confirm the favorable sion of the Generalized Differential Quadrature Method.
where r is the maximum order of the series expa nsions.By neglecting the series terms (45) with order higher than r and by introducing the relations (45) into the integrals (44), it is obtained: Thus, the engineering elastic stiffnesses A can be evaluated in the following manner: The same results can be exactly obtained using the exressions proposed by Qatu [16], by approximating the logarithmic function using a Taylor series expansion.Furthermore, the relations (47) can be defined in the form: p here: w Finally, the exact integration of the expressions (47) assume the following aspect: henever the reduced elastic constants w do not deend on p  co-ordinate.The engineering elastic stiffess efinitions proposed in the present work rep-resent a generalization of the ones proposed in [7,16-18, 77-79,85].n Copyright © 2013 SciRes.
the linear elastic behavior of anisotropic materials is as-First-order Shear Deformation Theory and it can be put in the following form:


are the mass density of the material per unit of volume a ckness of the elastic foundation at the top and the bottom surface of the shell, respectively.librium equations about the s and  directions, respectively.Furthermore, the generalized external forces , , , , to the external forces and the Winkler-Pasternak elastic foundation acting v on the top and f the sh be e aluated using the principle and can be written rence surface of doubly-curved shell as bottom surfaces o static equivalence on the refe follows: and shell, re the bottom surface of the , the top and the bottom surface of the shell.It is worth noting that the formulation of the Winkler-Pasternak foundation is based on the first-order approximated assumption that the foundation is a homogeneous material of uniform thickness the shear modulus of the Pa k elastic , evaluated, as it can be d from the explicit form of the equilibrium .Three typ s of boundary conditions are co ely the lly clamped edge boundary condition (C), soft simply supported edge boundary conditions (S) and the free edge boundary condition (F).The equations describing the boundary conditions can be written as follows:Clamped edge boundary conditions (C) ].For the staocedure enables to write the equations of equilibriu discrete form, transforming each space derivative condensation of non-domain degrees of freedom is performed: an excellent agreement for all the cases considered.GDQ results are compared with the FEM results obtained with Straus 7 commercial software using 8 node parabolic shell elements.In order to illustrate the effect of the foundation Fig the stress resultants for a at the top surface.

Figure 4
illustrates the stress resultants obtained without considering the foundation, whil Fig e5presents the same quantities obtained considerin tic foundation.The geomet- . Static deflection for a SSSS spherical panel at the po

Figure 4
Figure 4. tress result S ants for a   30 45 sph with y dist at the top surf e and wit found

4 ,
the results, in ten frequend by the GDQ Me Ge el pared with the FEM o d software.The details regarding etry e m a GST) prese above with an without ela .For the present GDQ results, the grid distributions (9) and (38) with 31 N M   have been considered.Tables 2-4 present the first ten frequencies for a SSSS spherical panel characterized by   0 90 ,   0 90 0 and   30 45 70 lamination scheme, respectively.As can be seen, the numerical results from the GDQ methodology are very close to those obtained by the commercial program and show an excellent agreement.Tables 5-7 pre-

Figure 6 .
Figure 6.Mode shapes for the CC shell of Table5.

Figu
Figu 7. Mode pes for t C shell able 7. a e et al. [ , respec ly.In b cases th isan between ults obtained is closely zero.
st and free ibration alyses o aminated oubly Winkler-Pasternak elastic foundations sented th e ll the effects of the foun ed e curvature effect.be -of re lution.S ple Rati al Bézie lamina schem with dif

Table 1 .
As

Table 6 . First ten frequencies for a  0 CCC panel resting on Winkler-Pasternak elastic foundation 30 C .
F

shell resting on Winkler-Pasternak elastic foundation. Table 7. First ten frequencies for a 30
*Control points and weights of the Bézier curve: F