Vibration analysis of variable thickness plates and shells by the Generalized Differential Quadrature method

Accepted Manuscript Vibration analysis of variable thickness plates and shells by the generalized differential quadrature method Michele Bacciocchi, Moshe Eisenberger, Nicholas Fantuzzi, Francesco Tornabene, Erasmo Viola PII: DOI: Reference:

S0263-8223(15)01081-8 http://dx.doi.org/10.1016/j.compstruct.2015.12.004 COST 7033

To appear in:

Composite Structures

Please cite this article as: Bacciocchi, M., Eisenberger, M., Fantuzzi, N., Tornabene, F., Viola, E., Vibration analysis of variable thickness plates and shells by the generalized differential quadrature method, Composite Structures (2015), doi: http://dx.doi.org/10.1016/j.compstruct.2015.12.004

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

VIBRATION ANALYSIS OF VARIABLE THICKNESS PLATES AND SHELLS BY THE GENERALIZED DIFFERENTIAL QUADRATURE METHOD

Michele Bacciocchi1, Moshe Eisenberger2, Nicholas Fantuzzi1, Francesco Tornabene*1, Erasmo Viola1

ABSTRACT. The main purpose of this work is to perform the free vibration analysis of several laminated composite doubly-curved shells, singly-curved shells and plates, characterized by a continuous thickness variation. Variable thickness could affect the design of shell structures since it allows to tailor the stiffness features in the most stressed areas within the domain, keeping the weight constant. As a consequence, an improved dynamic behavior may

be exhibited. The governing equations are solved numerically by the

Generalized Differential Quadrature (GDQ) method, which has proven to be an accurate, stable and reliable numerical tool. Its accuracy is tested by means of several comparisons with analytical and semi-analytical results available in the literature, and with the solutions obtained by a three-dimensional Finite Element (FE) model. The theoretical approach considered in the current paper is general and allows consideration of many higher-order structural theories in a unified manner, in which the order of the kinematic expansion can be chosen arbitrarily.

KEYWORDS: Variable thickness, Doubly-curved shells, Laminated composite structures, Natural frequencies, Generalized Differential Quadrature method.

1

DICAM - Department, School of Engineering and Architecture, University of Bologna, Italy. Department of Structural Engineering and Construction Management, Faculty of Civil and Environmental Engineering, Technion, Israel. * Corresponding author: [email protected]; web page: http://software.dicam.unibo.it/diqumaspab-project 2

1

1.

INTRODUCTION

The use of shell structures has grown hugely during the past five decades. Nowadays, in fact, it is easy to see examples of shell applications in many fields, such as aerospace, civil, mechanical and naval engineering. Thus, the need of more efficient structural elements has led to industrial advancement, and academic research towards great improvements of the performance of shells, due to their extraordinary features [1-10]. This aspect raised the interest of many researchers who aimed to understand in a clearer way the mechanical behavior of such structures. As a consequence, several structural models and theories have been derived. Since the development of thin shell theory by Love in 1944 [1] and enhanced subsequently through the introduction of differential geometry by means of Kraus’s work [2], many scientists aimed to improve these engineering theories in order to capture some effects that the original simpler models cannot predict. This gave rise to the development of several higher-order theories, as extensively stated in the books by Tornabene et al. [3, 4]. One of the most fascinating and investigated issue is the vibration analysis of shell elements since their dynamic characteristics affect the design such structures, as proven by the huge published literature on this topic, like the review type publications of Leissa [5], Liew et al. [6] and Qatu [7], with the most relevant publications since the early sixties are listed. The development of composite materials, such as fiber reinforced and laminated composites, and their consequent employment for shell manufacturing, allowed to highlight more and more the great features of shell structures, increasing for example the strength without a growth in terms of weight. Therefore, the advantages of composite materials coupled with the structural efficiency of shell elements gave start to an extensive use of such structures. The books by Reddy [8], Leissa and Qatu [9] and Tornabene and Fantuzzi [10] represent some of the fundamental contributions given to this topic. In the literature, it is possible to find lots of works that deal with several kinds of structural problems related to composite plates and shells. They aim to improve the mechanical behavior of these structures by means of a proper design of the lamination scheme or through the insertion of new classes of composites [11-40]. In addition, it is more and more common to read interesting papers about laminated shell structures containing functionally graded layers [41-57]. For the sake of completeness, it should be reminded that a functionally graded material (FGM) belongs to the class of granular composites and it is characterized by a continuous gradual variation of the mechanical properties of the material within the reference domain. The reason of their constant development is connected with the reduction of thermal stresses, residual stresses and stress concentrations which can affect a conventional laminate. Similar issues suggested the growth of nano-composites, such as carbon nanotube reinforced shells [58-63] and smart composites [64, 65].

2

The main purpose of the present paper is to investigate the free vibrations of laminated composite doubly-curved shells, singly-curved shells and plates, characterized by a continuous variation of the thickness profile. Structural elements with variable thickness, in fact, could exhibit better performance in bending, buckling and vibrations, when compared to structures with constant thickness. In addition, the use of variable thickness could lead to a shape optimization process which can allow to redistribute and increase the stiffness in the most stressed zones of the structure, while keeping the weight constant. According to the authors’ knowledge, several researchers dealt with this topic during the last years, focusing their works mainly on isotropic plates and singly-curved shells. In the following, a brief historical review is presented. An extensive analysis for vibration of isotropic rectangular plates with linearly variable thickness was carried out by Mizusawa [66], for several boundary conditions and geometric configurations. This work was expanded by Shufrin and Eisenberger [67], who presented a semi-analytical solution using the multi-term exatended Kantorovich method (MTEKM) for rectangular plates with variable thickness. The variation was assumed in polynomial form and, consequently, parabolic profiles are considered too. Dozio and Carrera [68] implemented the variable thickness in their numerical technique for vibration study of quadrilateral plates. Their solution was obtained combining Carrera’s unified formulation with a new version of the well-known Ritz method. In addition many analyses were performed on circular and annular plates. Here, the works of Eisenberger and Jabareen [69] and Efraim and Eisenberger [70] represent important benchmarks for this topic with the exact solutions for axisymmetric vibrations of circular and annular plates, made fromboth isotropic and functionally graded materials. Similarly, the natural frequency analysis of circular plates were the main topic of many other works [71-76], in which different numerical techniques were employed to solve the problem in hand. A three-dimensional method of analysis employing Ritz Method, was presented by Kang and Leissa for determining the free vibration frequencies and the corresponding mode shapes of spherical segments [77] and paraboloidal shells of revolution [78] with variable thickness. Another example is related to the free vibration analysis of a paraboloidal shell with variable thickness is the work of Leissa and Kang [79]. The same structure was considered also by Efraim and Eisenberger [80] and by Efraim [81], who employed the dynamic stiffness method for the exact frequency analysis. The presented high-precision numerical results allowed even the correction of some mistakes contained in the first paper by Kang and Leissa [77], who published as a consequence a corrigendum to their previous work about spherical shells [82]. The evaluation of natural frequencies of cylindrical shells with variable thickness was also performed by Duan and Koh [83], El-Kaabazi and Kennedy [84] and Ghannad et al. [85]. A deep

3

doubly-curved shell element was developed for free vibration analysis of general shells of revolution by Naghsh et al. [86], who presented many results for various structures with variable thickness. The free vibration problem of circular and annular membranes characterized by a variable thickness profile was solved exactly deriving the dynamic stiffness matrix by Jabareen and Eisenberger [87]. Lastly, the work of Jiang and Redekop [88] on the free vibration characteristics of linear elastic orthotropic toroidal shells with variable thickness is noted. According to the authors’ knowledge, this paper represents one of the few examples available in literature in which a sinusoidal thickness variation is applied to a composite structure made of orthotropic material. Having in mind all the works cited up to this point, it is the authors’ intention to highlight the main novelty introduced by the current paper. In particular, the evaluation of natural frequencies for several completely doubly-curved shells is performed taking into account different kinds of thickness variations that can be applied along both the two principal directions of the domain. In addition, the shell reference surface is described accurately by the application of the differential geometry. The solutions are obtained numerically by using the well-known Generalized Differential Quadrature (GDQ) method. This technique represents a valid tool for solving differential equations which leads to accurate results employing a small number of degrees of freedom. It should be noticed that the governing equations of motion are solved in their strong form, contrary to what usually finite element (FE) methods do by considering effectively the weak formulation of the same equations. A complete survey which is able to clarify all these aspects related to the numerical methods is represented in the work by Tornabene et al. [89]. An extensive account of the accuracy, stability and reliability characteristics of the GDQ method can be found in [90-97].

2.

SHELL THEORETICAL MODEL

A shell structure is a three-dimensional solid which can be studied by means of the classic theory of 3D elasticity [1-10]. Nevertheless, the calculation based on this approach turns out to be certainly burdensome, especially if the governing equations have to be solved numerically. In order to avoid this issue, the initial three-dimensional problem is converted into a two-dimensional problem by introducing a set of appropriate hypothesis [11, 12, 43]. According to the so-called Equivalent Single Layer (ESL) approach, the mechanical behavior of a shell structure is analyzed considering its middle surface as reference domain in which the fundamental equations, as well as each geometric and mechanical parameter, are defined. The mathematical description of the reference surface can be done conveniently by using the differential geometry [2, 3, 10]. With reference to Fig. 1, in which a

4

generic doubly-curved shell element with variable thickness is depicted, an arbitrary point P that lies inside the structure can be identified by the following position vector R (α 1 , α 2 , ζ ) = r (α 1 , α 2 ) +

h (α 1 , α 2 )

2

z n (α 1 , α 2 )

(1)

where r (α1 ,α 2 ) represents the position vector of the corresponding point P′ on the middle surface, n (α1 , α 2 ) is the outward unit normal vector and z = 2ζ h (α1 ,α 2 ) ∈ [ −1,1] . For the sake of completeness, a local reference system O ′α1α 2ζ has to be defined on the reference surface. In particular, α1 , α2 identify the coordinates upon the middle surface, whereas ζ is taken along the shell thickness. It should be underlined that in the present work the parametric lines of curvature α1 , α2 are orthogonal and principal, too. In addition, the coordinates α1 ,α 2 have to follow some restraints, since the shell three-dimensional domain is limited. Thus, the relation

α i ∈ α i0 , α i1  , for i = 1, 2 , is taken into account. Similarly, the parameter ζ cannot exceed the limits imposed by the shell thickness h (α1 , α 2 ) , which in general can vary along α1 , α2 , directions. It is important to emphasize that the thickness variation involves each layer if a laminated composite shell is considered (Fig. 2). In general, an arbitrary thickness profile is described through the following expression h (α1 , α 2 ) = h0 (1 + f (α1 , α 2 ) )

(2)

where f (α1 , α 2 ) can be any function of the dimensionless coordinates α1 , α 2

α1 =

α1 − α10 , α11 − α10

α2 =

α 2 − α 20 α 21 − α 20

(3)

In any case, the total thickness h (α1 , α 2 ) of a laminated composite doubly-curved shell made of l layers is given by l

h (α1 ,α 2 ) = ∑ hk (α1 , α 2 )

(4)

k =1

It should be underlined that the coordinates α1 , α2 are generic and assume a different meaning depending on the kind of structural element which is considered. For example, for a shell of revolution, they become α1 = ϕ ,α 2 = ϑ , respectively. Once the position vector r (α1 , α 2 ) is defined, by means of the differential geometry it is possible to evaluate the Lamé parameters A1 , A2 and the principal radii of curvature R1 , R 2 as

5

A1 (α1 , α 2 ) = r,1 ⋅ r,1 , R 1 (α 1 , α 2 ) = −

r,1 ⋅ r,1 r,11 ⋅ n

A2 (α1 , α 2 ) = r, 2 ⋅ r,2 R 2 (α 1 , α 2 ) = −

,

(5)

r,2 ⋅ r,2 r,22 ⋅ n

where r, i = ∂r ∂α i , for i = 1, 2 . Furthermore, the outward unit normal vector n (α1 , α 2 ) can be expressed as n=

r,1 × r,2

(6)

A1 A 2

For the sake of clarity, it should be specified that the symbols “ ⋅ ” and “ × ” employed in (5) and (6) denote the scalar product and the vector product, respectively. Considering the same unified formulation employed by the authors in their previous works [12, 23, 24, 28, 43], the three-dimensional displacement field can be written as U1   Fτ   N +1  U =  2  ∑ 0 U 3  τ = 0  0  

0

Fτ 0

(τ ) 0  u1   τ  0  u2( )    Fτ  u (τ )   3 

N +1

⇔

U = ∑ Fτ u (

τ)

in which U = U (α1 , α 2 , ζ , t ) the 3D displacement component vector, Fτ = Fτ (ζ matrix and u( ) = u( τ

τ)

(α1 , α 2 ,t )

(7)

τ =0

)

is the thickness function

is τ -th order generalized displacement component vector of points that belong

to the reference surface. The so-called thickness functions Fτ = Fτ (ζ ) in (7) can assume different meanings in order to define various higher-order kinematic models. The interested reader can find a list of the most common thickness functions (or shear functions) employed for structural purposes in the work [11, 12]. It should be noted that the kinematic expansion order, represented by the free parameter τ , is completely arbitrary and can be set properly to obtain the desired structural theory. For example, considering an expansion of third order ( N = 3) , the kinematic model (7) can be written in extended form as U1 = u1( ) + ζ u1( ) + ζ 2u1( ) + ζ 3u1( ) + ( −1) zk u1( 0

1

2

k

3

4)

U 2 = u2( ) + ζ u2( ) + ζ 2 u2( ) + ζ 3u2( ) + ( −1) zk u 2( U3 = u

0

1

( 0) 3

(1) 3

2

2 (2) 3

+ζ u +ζ u

3

3 ( 3) 3

+ζ u

k

4)

k

( 4) 3

+ ( −1) zk u

(8)

where the thickness functions Fτ = Fτ (ζ ) is chosen as ζ τ Fτ =  k (−1) zk

for τ = 0,1,..., N for τ = N + 1

6

(9)

It should be underlined that the thickness function related to the ( N + 1) -th order of expansion represents the Murakami’s function needed to consider the zig-zag effect along the thickness. The dimensionless parameter zk = zk ( ζ ) ∈ [ −1,1] is defined as zk =

2

ζ k +1 − ζ k

ζ−

ζ k +1 + ζ k ζ k +1 − ζ k

(10)

when a generic laminated composite doubly-curved shell is considered (Fig. 2). For conciseness purpose, the displacement field (7) can be indicated by using a compact notation when the thickness functions are chosen according to (9) ED −  1   ζ   ζ 2  ...  ζ N  ζ N +1 

(11)

where the letter “E” indicates that an ESL approach is taken into account, whereas “D” states that the governing equations are written in terms of generalized displacements. The most common higher-order theories are listed below ED −  1   ζ  = ED1 ED −  1   ζ   ζ 2  = ED2 ED −  1   ζ   ζ 2   ζ 3  = ED3

(12)

ED −  1   ζ   ζ 2   ζ 3   ζ 4  = ED4

If the Murakami’s function is introduced in the kinematic model, letter “Z” will be inserted in the compact notation. For example, a third-order theory able to capture the zig-zag effect will be indicated as EDZ3. It is important to specify that the present formulation is able to obtain the well-known Reissner-Mindlin theory. Once the kinematic model is defined, it is possible to get the relation between the displacement component and the three-dimensional strain components. If the parametric lines of curvature α1 , α 2 are orthogonal and principal, the kinematic equations take the following form ε = Dζ DΩ U

where ε (α1 , α 2 , ζ , t ) = [ε 1 ε 2

γ 12

γ 13 γ 23

T

ε3 ]

(13)

is the strain component vector. The first differential

operator Dζ contains the terms related to ζ and is defined as

7

 1 H  1   0    0 Dζ =    0    0    0 

0

0

0

0

0

0

0

1 H2

0

0

0

0

0

0

0

1 H1

1 H2

0

0

0

0

0

0

0

1 H1

0

∂ ∂ζ

0

0

0

0

0

1 H2

0

∂ ∂ζ

0

0

0

0

0

0

0

 0    0    0    0    0   ∂   ∂ζ 

(14)

in which H i = 1 + ζ Ri , for i = 1, 2 , whereas the second differential operator DΩ includes the terms linked to the curvature along α1 , α 2 coordinates and takes the aspect shown below  1 ∂   A1 ∂α1  1 ∂A 2   A1 A 2 ∂α1   − 1 ∂A1  A A ∂α 2  1 2  1 ∂ DΩ =  A 2 ∂α 2   1  − R1    0   1   0  0 

   1  R2   0     0   1 ∂   A1 ∂α1   1 ∂  A 2 ∂α 2   0  0   1 

1 ∂A1 A1 A 2 ∂α 2

1 R1

1 ∂ A 2 ∂α 2 1 ∂ A1 ∂α1 −

1 ∂A2 A1 A2 ∂α1 0 −

1 R2 0 1 0

(15)

The expression of the generalized strain defined on the shell reference surface can be deduced easily by inserting the general displacement field (7) in the relation (13) N +1

N +1

N +1

N +1

ε = ∑ Dζ DΩ Fτ u ( ) = ∑ Dζ Fτ DΩ u ( ) = ∑ Z ( ) DΩ u ( ) = ∑ Z ( ) ε ( τ =0

where ε(τ ) (α1 , α 2 , t ) = ε1(τ ) 

τ

τ

τ =0

τ

τ

τ =0

τ

τ)

(16)

τ =0

T

ε 2(τ ) γ 1(τ ) γ 2(τ ) γ 13(τ ) γ 23(τ ) ω13(τ ) ω23(τ ) ε 3(τ )  represents the τ -order generalized 

τ strain component vector, whereas Z( ) is the matrix which includes the terms related to the shell thickness

defined as

8

Z

(τ )

 Fτ   H1   0    0  =   0    0    0 

0

0

0

0

0

0

0

Fτ H2

0

0

0

0

0

0

0

Fτ H1

Fτ H2

0

0

0

0

0

∂Fτ ∂ζ

0

0

0

0

Fτ H1

0

0

0

0

Fτ H2

0

∂Fτ ∂ζ

0

0

0

0

0

0

0

 0    0    0     0    0   ∂Fτ   ∂ζ 

(17)

The relationships between generalized strains ε(τ ) and generalized displacements u(τ ) can be deduced easily from equation (16)

ε ( ) = DΩ u (

τ)

τ

for τ = 0,1, 2,… , N , N + 1

(18)

If a laminated composite shell is considered, each element of the strain component vector ε (α1 , α 2 , ζ ,t ) has to be referred to the k -th layer (or ply) and can be related to the corresponding stress through the constitutive equations, assuming that each layer is made of a linear elastic material

σ ( ) = C( ) ε ( k

where σ (

k)

k

k)

(α1 ,α 2 , ζ , t ) = σ1( k ) σ 2( k ) τ 12( k ) τ 13( k ) τ 23( k ) σ 3( k ) 

(19)

T

is the stress component vector. For example, if

the k -th layer is made of orthotropic material, the constitutive operator C(

k)

takes the following extended

representation

C(

k)

 E11( k )  k  E12( )  (k ) E =  16  0   0  E (k )  13

E12( )

E16( )

0

0

E22( )

E26( ) k

0

0

(k )

(k )

0

0

0

0

(k )

E44

E45( )

0

0

E45( )

E55( )

E23( )

E36( )

0

0

k

k

k

E26

E66

k

k

k

k k

E13( )   k E23( )   k E36( )   0   0  k E33( )  k

(20)

( ) It should be underlined that the elastic constants are denoted generically by the symbol Enm . This choice is k

justified by the fact that the elements of C( k ) can be set equal to the plane stress-reduced elastic coefficients ( ) ( ) Enm = Qnm k

k

9

(21)

or taken as the non-reduced ones ( ) ( ) Enm = Cnm k

k

(22)

Such distinction allows to consider those structural theories that neglect the stretching effect along the shell thickness and need the reduced elastic coefficients. This is the case of the Reissner-Mindlin theory, for example. Once the constitutive equations are established, it is possible to write the indefinite dynamic equilibrium equations for a doubly-curved laminated shell through the well-known Hamilton’s principle. In order to simplify the problem under consideration, it is advantageous to express the internal actions in terms of stress resultants for each order of kinematic expansion l

S( ) = ∑ τ

ζ k +1

k =1

where S (τ ) (α1 , α 2 , t ) =  N1(τ ) 

∫ (Z ζ

N 2(

τ)

(τ )

T

)

σ ( ) H1 H 2 d ζ k

(23)

for τ = 0,1, 2,..., N , N + 1

k

N12( ) τ

N 21( ) τ

T1(

τ)

T2(

τ)

P1(

τ)

P2(

τ)

τ S3( )  

T

is the τ -th order generalized

stress resultant vector. By substituting relations (19) and (16) into the previous definition, it is possible to obtain the expression of the stress resultants S (

τ)

N +1

as a function of the generalized displacements u (

τ)

S ( ) = ∑ A ( ) DΩ u( τ

τη

η)

for τ = 0,1, 2,..., N , N + 1

(24)

η =0

l

where A(τη ) = ∑ k =1

ζ k +1

∫ (Z ζ

(τ )

T

)

k η C( ) Z( ) H1 H 2 d ζ represents the stiffness matrix, for τ ,η = 0,1, 2,..., N , N + 1 . In the

k

recurring case of an orthotropic medium, the stiffness matrix (or constitutive operator) can be written in extended matrix form as

A(

τη )

 A(τη )  11( 20)  A(τη )  12(11)  A(τη )  16( 20)  (τη )  A16(11)  = 0   0   0   0   A(τη )  13(10)

A12( (11) ) τη

A16( ( 20) ) τη

A16( (11) ) τη

0

0

0

0

A22( (02) ) τη

( ) A26 (11) τη

( ) A26 ( 02 ) τη

0

0

0

0

τη A26( (11) )

τη A66( ( 20) )

(τη ) A66 (11)

0

0

0

0

A26( (02) ) τη

A66( (11) ) τη

A66( ( 02) ) τη

0

0

0

0

0

0

0

A44( ( 20) ) τη

A45( (11) ) τη

A44( (10) ) τη

A45( (10) )

0

(τη )

A45(11)

(τη )

A55( 02)

(τη )

A45( 01)

(τη )

A55(01)

 τη A45( ( 01) )

 ) (τη A44 ( 00 )

A45( ( 00) )

0

0

τη

0

0

0

 ) (τη A44 (10)

0

0

0

( ) A45 (10 )  τη

A55( ( 01) )  τη

A45( ( 00) )  τη

A55( ( 00) )

( ) A23 ( 01)

A36( (10) )

A36( ( 01) )

0

0

0

0

 τη

 τη

 τη

10

 τη  τη

A13( (10) )   (τη )  A23( 01)  τη A36( (10) )  τη  A36( ( 01) )   0   0   0   0    τη A33( ( 00) )   τη

(25)

The elements of the stiffness matrix (25) are defined through the following relations l

( ) Anm ( pq ) = ∑ τη

k =1

∫ ζ

l

ζ k +1

( ) Anm ( pq ) = ∑  τη

(τη ) Anm ( pq )

k =1 l

ζ k +1

=∑

( ) Bnm Fη Fτ k

k

∫ ζ

( ) Bnm Fη k

k

∫

(k )

Bnm

k =1 ζ k l

( ) Anm ( pq ) = ∑  τη

ζ k +1

k =1

ζ k +1

∫ ζ k

( ) Bnm k

H1 H 2 dζ H1p H 2q

∂Fτ H1 H 2 dζ ∂ζ H1p H 2q

for τ ,η = 0,1, 2,..., N , N + 1

∂Fη

H1 H 2 Fτ dζ ∂ζ H1p H 2q

for n, m = 1, 2, 3, 4, 5, 6 for p , q = 0,1, 2

(26)

∂Fη ∂Fτ H1 H 2 dζ ∂ζ ∂ζ H1p H 2q

in which the parameters τ ,η specify the thickness functions Fτ , Fη , whereas τ ,η stand for the corresponding derivatives ∂Fτ ∂ζ , ∂Fs ∂ζ , respectively. It should be underlined that a numerical technique has to be employed for the computation of such constants. In the following section, more details will be given in this regard. In addition, the integration limits in (26) vary in each point of the domain if a shell with variable (k ) thickness is considered. The symbols Bnm are introduced in (26) to include more easily in this treatise the shear

correction factor κ ( ζ ) , since they have the following meaning (k )  Enm (k ) Bnm = (k ) κ (ζ ) Enm

for n, m = 1, 2, 3, 6 for n, m = 4, 5

(27)

It should be remarked that the shear correction factor has to be used for those structural theories that do not consider a parabolic profile of shear stresses along the shell thickness. Classically, the shear correction factor

κ (ζ ) is taken equal to the constant value of 5 6 , even if many authors suggested using a different one, as it can be noticed in the work by Efraim and Eisenberger [80]. On the other hand, Tornabene and Reddy proposed employing a correction that is a quadratic function of ζ , instead of a constant value [41]. In any case, the correction factor κ (ζ ) will be specified whenever it is used. For example, the notation ED2κ =5 6 points out that a structural theory defined by a second order of kinematic expansion is considered and the constant value of 5 6 is employed for the shear correction factor. In the same way, it will be specified the introduction of the plane stress-reduced elastic coefficients. The abbreviation ED1RS means that the reduced coefficients are taken into account and no shear correction factor is considered. The combination of these notations will be used if both the shear correction factor and the plane stress-reduced elastic coefficients are employed, such as for the Reissnerκ =5 6 Mindlin theory ( FSDTRS ). Finally, the structural problem is concluded by the governing equations of motion,

11

and the corresponding boundary conditions, which can be obtained by Hamilton’s principle. The three motion equations for each order of kinematic expansion can be written in compact matrix form as N +1

τ τη (η )  D*Ω S ( ) = ∑ M ( ) u

for τ = 0,1, 2,..., N , N + 1

(28)

η =0

in which D*Ω represents the equilibrium operator  1   A1        1 D*Ω =  A2            

whereas M (

τη )

∂ ∂α 1 −

+

1 ∂A2 A1 A2 ∂α1

−

1 ∂A 2 A1 A2 ∂α1

1 ∂A1 A1 A2 ∂α 2

1 ∂ 1 ∂A1 + A2 ∂α 2 A1 A2 ∂α 2

1 ∂A1 A1 A2 ∂α 2

1 ∂ 1 ∂A2 + A1 ∂α1 A1 A2 ∂α1

∂ 1 ∂A1 + ∂α 2 A1 A2 ∂α 2

    1  −  R2   0    0   ∂ 1 ∂A2  +  ∂α1 A1 A2 ∂α1  ∂ 1 ∂A1   + ∂α 2 A1 A2 ∂α 2   0   0  −1  −

1 ∂A2 A1 A2 ∂α1

1 R1

0

1 A1

0

1 R2

1 A2

−1

0

0

−1

0

0

1 R1

T

(29)

is the inertia matrix

M

(τη )

 I (τη )  = 0   0

0   0  for τ ,η = 0,1, 2,..., N , N + 1  τη I( ) 

0 I

(τη )

0

(30)

The mass inertia terms included in the matrix (30) can be written as I(

τη )

l

=∑ k =1

where ρ (

k)

ζ k +1

∫ ζ

ρ ( k ) Fτ Fη H 1 H 2 d ζ

(31)

k

denotes the density of the material for the k -th layer. Finally, it is possible to describe the

mechanical behavior of a doubly-curved shell structure as a function of the generalized displacements u (τ ) only. Thus, inserting the kinematic equations (18) and the constitutive relations (24) in (28), the equations of motion for each order of kinematic expansion can be written in terms of the generalized displacements and are named fundamental equations. The set of fundamental equations takes the following aspect

12

N +1

L( ∑ η

τη )

N +1

η τη (η )  u( ) = ∑ M ( ) u

=0

in which L(

τη )

for τ = 0,1, 2,..., N , N + 1

(32)

η =0

represents the fundamental operator defined as

L(

τη )

 L(11τη )  τη τη = D*Ω A ( ) DΩ =  L(21 )  (τη )  L31

L(12 ) τη

L(22 ) τη

(τη )

L32

L(13 )   τη L(23 )   τη L(33 )   τη

(33)

τη ) The complete expression of each element L(nm , for n, m = 1, 2,3 and τ ,η = 0,1, 2,..., N , N + 1 , can be found in

[12]. It should be noted that in general the total number of equations which constitute the fundamental system is equal to 3 × ( N + 2 ) . Finally, proper boundary conditions allow to solve the differential system in hand. The different boundary conditions taken into account in the following are written below Clamped edge boundary conditions (C) u1( ) = u2( ) = u3( ) = 0

for τ = 0,1, 2,..., N , N + 1,

u1( ) = u2( ) = u3( ) = 0

for τ = 0,1, 2,..., N , N + 1,

τ

τ

τ

τ

τ

τ

at α1 = α10 or α1 = α11 , α 20 ≤ α 2 ≤ α 12 at α 2 = α 20 or α 2 = α 21 , α10 ≤ α1 ≤ α11

(34) (35)

Free edge boundary conditions (F) τ τ τ N1( ) = 0, N12( ) = 0, T1( ) = 0

for τ = 0,1, 2,..., N , N + 1,

at α1 = α10 or α1 = α11 , α 20 ≤ α 2 ≤ α 12

(36)

τ τ (τ ) = 0, N 2( ) = 0, T2( ) = 0 N 21

for τ = 0,1, 2,..., N , N + 1,

at α 2 = α 20 or α 2 = α 21 , α10 ≤ α1 ≤ α11

(37)

Simply-supported boundary conditions (S) τ τ τ N1( ) = 0, u2( ) = 0, u3( ) = 0

for τ = 0,1, 2,..., N , N + 1,

at α1 = α10 or α1 = α11 , α 20 ≤ α 2 ≤ α 21

(38)

τ τ τ u1( ) = 0, N 2( ) = 0, u3( ) = 0

for τ = 0,1, 2,..., N , N + 1,

at α 2 = α 20 or α 2 = α 12 , α10 ≤ α1 ≤ α11

(39)

In addition, if the shell structure under consideration is characterized by two coincident sides, such as a complete revolution shell or a toroid, another set of kinematic and physical compatibility conditions has to be enforced to satisfy the structural congruence, since these edges represent the same side from the physical point of view. In particular, if the common edges are identified by the coordinate α1 = α10 and α1 = α11 , the following conditions have to be added, for τ = 0,1, 2,..., N , N + 1 and α 20 ≤ α 2 ≤ α 12

13

(α N ( ) (α T ( ) (α N1(

τ)

τ

12

τ

1

0 1 0 1 0 1

) ) (α , α , t ) , , α , t ) = N ( ) (α , α , t ) , , α , t ) = T ( ) (α , α , t ) , ,α 2 , t = N1(

τ

τ

2

12

τ

2

1 1 1 1

1 1

1

(α u ( ) (α u ( ) (α u1(

2

2

τ)

0 1

τ

0 1

2

τ

2

3

0 1

) ) (α , α , t ) , α , t ) = u ( ) (α , α , t ) , α , t ) = u ( ) (α , α , t ) , α 2 , t = u1(

τ

1 1

2

τ

1 1

2

1 1

2

2

2

2

3

τ

(40)

On the other hand, the next conditions have to be written if the common sides are characterized by the coordinates α 2 = α 20 and α 2 = α 21 , for τ = 0,1, 2,..., N , N + 1 and α10 ≤ α1 ≤ α11

( ) ( ) ( ) ( ) N (α , α , t ) = N (α , α , t ) , T ( ) (α , α , t ) = T ( ) (α , α , t ) ,

(α , α , t ) = u ( ) ( α , α , t ) u ( ) (α , α , t ) = u ( ) ( α , α , t ) u ( ) (α , α , t ) = u ( ) ( α , α , t )

( ) N 21 α1 , α 20 , t = N 21( ) α1 ,α 12 , t , τ

τ

2

1

0 2

1

0 2

τ

2

u1(

τ)

τ

τ

2

1

τ

2

1

0 2

1

1

0 2

2

1

0 2

3

τ

1 2

2

τ

1 2

3

τ

1

1

1 2

1

1 2

1

1 2

τ

τ

(41)

With reference to Fig. 1, a convention should be introduced now to identify univocally the edge where the boundary conditions have to be enforced. In particular, when a shell structure has no common edges, that is the case of a doubly-curved panel, the following notation is used to specify the bounded side West edge (W) South edge (S) East edge (E) North edge (N)

α10 ≤ α1 ≤ α11 , α 2 = α 20

→ → → →

α1 = α11 , α 20 ≤ α 2 ≤ α 21 α10 ≤ α1 ≤ α11 , α 2 = α 21

(42)

α1 = α10 , α 20 ≤ α 2 ≤ α 21

In the numerical examples presented below, the symbolism WSEN is used to specify the boundary condition sequence for the panel at issue. For example, if a generic panel is clamped along the West (W), the South (S) and the East (E) edges and free along the North (N) side, the right sequence is CCCF. On the other hand, if a structure with two common edges is analyzed, only two boundary conditions are specified. For instance, if the common sides have coordinates α1 = α10 and α1 = α11 respectively, the acronym FC means that the Western edge is free and the Eastern side is clamped, whereas in the Southern and the Northern edges the compatibility conditions (40) have to be enforced. On the contrary, the notation FC specifies that the North edge is free and the Southern one is clamped, whereas the other two sides are involved by the compatibility conditions (41), if the common sides have coordinates α 2 = α 20 and α 2 = α 21 . Lastly, it should be remarked that the present theoretical model allows to analyze the mechanical behavior of thick and moderately thick shells, for which the following limits are valid

 h 1 h  1 ≤ max  , ≤ 100  Rmin Lmin  5

14

(43)

3.

NUMERICAL SOLUTION OF THE FREE VIBRATION PROBLEM

The GDQ method [89] is employed to obtain the numerical solution of the fundamental system of equations (32) . This numerical method, whose accuracy and excellent features have already been highlighted in many structural applications during the last years [11-28, 41-46, 89-97], approximates the n -th derivative at a generic point xi of a sufficiently smooth function f ( x ) as a weighted linear sum of the function values at some chosen grid points d n f ( x) dx n

IN

≅ ∑ ς ij( n ) f ( x j ) x = xi

i = 1, 2,..., I N

(44)

j =1

where I N is the total number of nodes (or grid points). From definition (44), it is clear that the core of the procedure is the computation of the weighting coefficients ς ij( n ) . It should be specified that expression (44) is related to a one-dimensional problem, but its extension to two-dimensional domains is rather immediate. In general, the n -th derivative with respect to x -coordinate, the m -th derivative with respect to y -coordinate and the mixed partial derivative of order ( n + m ) of a two-dimensional smooth function f ( x, y ) at a generic point

( xi , y j ) can be written respectively as ∂ n f ( x, y ) ∂x n

IN

≅ ∑ ς x( (ik) ) f ( xk , y j ) n

x = xi , y = y j

∂ m f ( x, y ) ∂y m

IM

≅ ∑ ς y( m( jl) ) f ( xi , yl ) x = xi , y = y j

∂ n + m f ( x, y ) n

∂x ∂y

i = 1, 2,..., I N , j = 1, 2,..., I M

k =1

n x = xi , y = y j

i = 1, 2,..., I N , j = 1, 2,..., I M

(45)

l =1

IN  IM  ≅ ∑ ς x( n(ik) )  ∑ ς y( m( jl) ) f ( xk , yl )  k =1  l =1 

i = 1, 2,..., I N , j = 1, 2,..., I M

where I N , I M are the effective number of grid points along x and y directions, whereas ς x( ( ik) ) and ς (y ( jl) ) n

m

represent the weighting coefficients for the derivatives. From expressions (45) we have that the derivatives with respect to x -coordinate contemplate every point of the domain along a line defined for a specific value of

y = y j . In addition,, the derivatives with respect to y -coordinate involve each node that lies along a line for x = xi fixed. Lastly, the mixed derivatives consider all the points of the two-dimensional domain since both indices k , l appear in the linear sum. In any case, it is clear that the derivative approximation can be done if a proper grid distribution is chosen. In the present work, the Chebyshev-Gauss-Lobatto grid distribution is used for all the next numerical applications due to its great efficiency and high reliability. Therefore, the points of the

15

reference surface are located according to the following expressions along the principal coordinates α1 , α 2 

 i −1

1 0   (α1 − α1 )

α1i = 1 − cos π    IN −1   

2

  j − 1   (α 2 − α 2 ) π   = 1 − cos + α 20 , I − 1 2  M    1

α2 j

+ α10 , i = 1, 2,..., I N ,

for α1 ∈ α10 , α11  (46)

0

j = 1, 2,..., I M ,

for α 2 ∈ α , α  0 2

1 2

The free vibrations of a moderately thick shell can be obtained from Equation (32). It is convenient at this point to write the system solution in the following form

u(

τ)

(α1 ,α 2 , t ) = U(τ ) (α1 , α2 ) eiωt

where U (τ ) = U 1(τ ) (α1 , α 2 ) U 2(τ ) (α1 , α 2 ) U 3(τ ) (α1 , α 2 )   

T

(47)

is the mode shape vector, whereas ω is the

corresponding circular frequency, from which it is possible to evaluate the natural frequency as f = ω 2π . By evaluating the time derivatives in Eqn. (47) and substituting the results in Eqn. (32), the fundamental system can be written as N +1

N +1

∑ L(τη ) U(η ) + ω η∑ M(τη ) U(η ) = 0 η 2

=0

for τ = 0,1, 2,..., N , N + 1

(48)

=0

Once the discrete form of Eqn. (48) is obtained by using the GDQ method, the governing equations take the following aspect

Kδ = ω 2 Mδ

(49)

in which K is the global stiffness matrix, M is the global mass matrix and δ is the mode shape vector. Finally, the numerical problem can be considerably simplified by separating the components linked to the boundary nodes (identified by the letter “ b ”) from those related to the internal points (specified by the letter “ d ”). Through the well-known kinematic condensation of non-domain degrees of freedom, the system (49) becomes

K bb K  db

K bd   δ b  0  δb  0 = ω2       K dd   δ d  0 M dd  δ d 

(50)

After performing some easy algebraic substitutions, it is possible to get the following result

(M (K −1 dd

dd

)

−1 − K db K bb K bd ) − ω 2 I δ d = 0

(51)

At this point, the numerical solution of the standard eigenvalue problem (51) allows to obtain the natural frequencies of the structure. Finally, it is worth spending some words about the numerical computation of the stiffness coefficients introduced

16

previously in Eqn. (26), since that issue has been only touched quickly. The exact solution for such integrals does not exist due to the curvature terms H1 , H 2 and the thickness functions that can be chosen arbitrarily. Thus, ( ) ( ) ( ) ( ) the elastic constants Anm ( pq ) , Anm ( pq ) , Anm ( pq ) , Anm ( pq ) have to be evaluated numerically. In the present paper, the τη

 τη

τη

 τη

Generalized Integral Quadrature (GIQ) method is employed for this purpose. This techniques is based on the fact that the integration is the inverse procedure of derivation. Thus, it can be obtained directly starting from the definition of derivative given by the GDQ method. Without addressing the complete demonstration, it is important to explain that the GIQ rule allows to compute the integral of a smooth function f ( x ) as a weighted linear sum of the function values in the whole domain by using the same concepts of GDQ. In general, the numerical integration of f ( x ) over a limited domain  xi , x j  can be written as xj

∫ xi

IT

f ( x ) dx = ∑ ( w jk − wik ) f ( xk )

(52)

k =1

in which w jk , wik represent the weighting coefficients for the integral that can be evaluated directly from the weighting coefficients of the first order derivative, whereas IT is the total number of grid points of the domain. In the present applications, the integration procedure takes place along the shell thickness, therefore a discretization along this direction has to be set. In particular, the Chebyshev-Gauss-Lobatto grid distribution illustrated in (46) is employed also in this circumstance, with IT = 51 for all the cases. The interested readers can find all the explanations they need about the GIQ rule in the recent review paper by Tornabene et al. [89].

4.

RESULTS AND DISCUSSION

The main aim of the present section is to show several numerical applications related to the free vibration analysis of plates and shells with variable thickness. The solutions are obtained by using a MATLAB code [98] that implements the GDQ technique to solve the strong form of the fundamental system (32). The accuracy of both the shell theoretical model and the numerical method presented above is validated through a wide comparison with the results that can be found in the literature and results from a 3D FE model.

5.1 Comparison study

The starting point of the current investigation is the comparison of the natural frequencies that can be evaluated for several isotropic plates. The reference solutions are shown in the work by Shufrin and Eisenberger [67] and

17

are based on the extended Kantorovich method. The first example is related to a SSSS isotropic square plate of sides Lx , Ly with a linear variation in x -direction (Fig 3a). Mathematically speaking, the thickness profile is given by the following relation  x  h ( x ) = h0  1 −  Lx  2 

(53)

in which h0 is the value of the thickness in the origin of the local coordinate reference system. The frequencies are expressed in terms of the dimensionless factor

λ=

2 f L2y

12 ρ (1 − ν 2 )

π

Eh02

(54)

where ρ is the material density, E is the Young’s modulus and ν is the Poisson’s ratio taken as 0.3. The mode shapes of vibration are defined by the parameters m and n , which represent the number of half-waves in the x and y directions, respectively. It should be noticed that the reference solutions are obtained by setting the value of the shear correction factor equal to κ = π 2 12 . Table 1 shows the natural frequency comparison for several values of the ratio h0 Ly . It is easy to see that the GDQ solutions, obtained for different structural theories indicated by the nomenclature specified above, are in good agreement with the reference ones. Since the Reissner-Mindlin theory gives the best results, it will be the only one considered in the following tests related to plates. This aspect was expected since the reference solutions are obtained using a first-order shear deformation theory. The second application deals with the comparison of the frequency factor λ for an isotropic CCCC square plate with linear thickness variation in both directions. The variable thickness is given by the following relation  x y xy  h ( x, y ) = h0  1 − β x − βy + βx β y   Lx Ly Lx Ly  

(55)

where β x , β y are the thickness increment factors along x and y directions, respectively. In this case, the adopted shear correction factor is equal to κ = 5 6 as in the reference solutions. Table 2 shows the excellent fit between the frequency parameters of the GDQ method and the ones presented by Shufrin and Eisenberger [67], for each value of the ratio h0 Ly and of the constants β x , β y . The frequency parameter values reported in Table 3 are referred to a parabolic variation of the plate thickness in x direction. In particular, three kinds of variation are considered and they are named respectively as “arched

18

form” (Fig. 3c), “concave form” (Fig. 3d) and “symmetric concave form” (Fig. 3e). Such profiles can be expresses by the following relations  1  x 2  h ( x ) = h0  1 −     2  Lx     2  1 x  x  h ( x ) = h0  1 +   −   2  Lx  Lx   

Arched form:

Concave form:

Symmetric concave form:

(56)

2   x  2x   h ( x ) = h0  1 + 2   −  Lx   Lx   

The comparison shown in Table 3 proves again the validity of the present approach, for every values of the ratio h0 Ly and for the three variations describe analytically in (56).

With the aim of demonstrating the effectiveness of the current technique also for shell structures, in the following the comparison is performed for two different spherical shells with linearly variable thickness. The reference solutions are the exact results from dynamic stiffness analysis taken from the paper by Efraim and Eisenberger [80]. The first application is related to the computation of the natural frequencies for an FF isotropic hemispherical annular shell of radius R shown in Fig. 4a. In the same picture is indicated also the reference surface position vector r (ϕ , ϑ ) and the geometric parameters needed for its representation, as well as the linear thickness variation along the shell meridian h (ϕ ) . The values of frequencies are given in Table 4 in a nondimensional form as Ω = 2π f R

ρ G

(57)

in which G is the shear modulus. The mode shapes of vibration are defined by the parameter n , which represent the number of half-waves along ϑ -direction. For n = 0 the axisymmetric modes (A) and the torsional modes (T) are obtained. It should be noted that in the spherical shell examples the shear correction factor is set equal to

κ = 5 ( 6 −ν ) , where ν = 0.3 is Poisson’s ratio. The comparison shown in Table 4 is performed with the results from FE analysis, with the results from a three-dimensional analysis by the Ritz method reported by Kang and Leissa in [82] and from the exact ones given by the dynamic stiffness method proposed by Efraim and Eisenberger [80]. The present solutions which are in good agreement with the references is obtained considering κ = 5 ( 6 −ν ) the FSDTRS as a structural theory. It should be noted that the number of grid point along ϑ -direction is

larger than the one along ϕ -direction. This choice is caused by the fact that the dimensions of the shell are quite

19

different along the two principal coordinates. The same considerations can be made for the FF isotropic barrel shell with linearly variable thickness shown in Fig. 4b, whose dimensionless frequencies are compared in Table 5. For the sake of completeness, Fig. 5 shows the first mode shape for each value of n for the FF barrel shell with variable thickness under consideration. The last comparison is based on the natural frequency computation for a FF isotropic paraboloidal shell with linearly variable thickness. Its geometry and thickness profile representations, as well as the reference surface position vector r (ϕ , ϑ ) , are shown in Fig. 4c. The reference solution is taken from Efraim’s thesis [81]. Even in this circumstance, with reference to Tab. 6, good agreement can be seen between the reference values and the GDQ solutions. It should be noted that the natural frequencies are presented in their dimensionless form, κ = 5 ( 6 −ν ) according to expression (57). As in the previous case, the adopted structural model is the FSDTRS , setting

the Poisson’s ratio equal to 0.3. Finally, Fig. 6 shows the first mode shape for each value of n for the FF paraboloidal shell with variable thickness under consideration. All things considered, both the GDQ method and the shell theoretical model presented above are validated in an excellent way through these comparisons.

5.2 Natural frequencies for singly-curved shells In this section , the free vibrations of two singly-curved panels are investigated. The first example is related to a FCCC conical panel with cubic thickness variation along y -coordinate. The structural element is shown in Fig. 7a, together with the position vector r (α1 , y ) for the description of the middle surface, the expression of the thickness variation h ( y ) , and the main geometric parameters. It is important to specify that the local coordinate reference system is always indicated in the representation of the geometry so that the identification of the direction affected by the thickness variation can be easily seen. For this purpose, it should be reminded that such reference system is right-handed in any case. Therefore, the other two directions are defined univocally following the principal coordinate lines. The structure is made of two isotropic layers with lamination scheme (Aluminum/Steel). The mechanical properties of the two materials are listed in Table 7. It should be noted that the thickness variation affects both layers, which are characterized by the same value of the initial value h1 = h2 = 0.05m , measured in the origin of the local coordinate reference system. The GDQ solutions are given

for seven structural theories, as it can be seen from Table 7, where the first ten natural frequencies are shown. In the same table, the results obtained by a 3D FEM model (61200 20-node brick elements of Abaqus) are reported

20

as reference. From the values shown in Table 7, it is possible to verify that each theory leads to similar results, which are in excellent agreement with the FE solution. The ED2κ =5 6 and the EDZ3 seem to be particularly accurate. In addition, it is important to notice that the structural theories with a kinematic expansion of second order should be employed with a shear correction factor. In fact, the ED2κ =5 6 gives better results than the corresponding ED2 , if compared with FE solutions. In general, the use of Murakami’s function for the zig-zag effect does not implicate significant improvements. Finally, the first nine mode shapes are depicted in Fig. 8. The second application deals with a FCFC cylinder of translation with elliptic profile made of two orthotropic layers of Graphite-Epoxy, whose mechanical properties are listed in Table 8. The lamination scheme is

( −45, 45) , as it can be seen from the same table. The geometry of the panel under consideration is depicted in Fig. 7b with the indication of the reference surface position vector r (ϕ , y ) and the main geometric features, as well as the variable thickness expression h (ϕ ) . Such variation involves ϕ -coordinate only and it is characterized by a sine wave shape. Table 8 shows the first ten natural frequencies for the translational panel obtained for several structural theories and for the 3D FEM model (100548 20-node brick elements of Abaqus). Good agreement can be observed if a comparison between the results for the two different numerical techniques is performed. In addition, it should be noted that even in this circumstance the structural theories with a first order kinematic expansion of should be used with the shear correction factor. In fact, the ED1κRS= 5 6 gives better results than the corresponding model without κ . Finally, Fig. 9 shows the first nine mode shapes for the current structure.

5.3 Natural frequencies for doubly-curved shells Here free vibrations analysis is performed for two laminated composite doubly-curved shells with variable thickness. Firstly, a FCFF panel of revolution with catenary shaped meridian made of three layers is considered. Its geometry is depicted in Fig. 7c, as well as the indication of the middle surface position vector r (ϕ , ϑ ) and the thickness variation expression h (ϕ , ϑ ) . It should be noticed that such variation affects both principal directions. In particular, a linear variation is considered along the ϕ -coordinate, whereas a sinusoidal variation is established along the ϑ -coordinate. The total thickness measured in the local reference system origin is equal to h0 = 0.045 m . This geometric configuration is analyzed for two different lamination schemes. Initially, three

21

isotropic layers are considered. In particular, the external plies are made of Aluminum, whereas the inner one is made of Zirconia. Their mechanical properties are listed in Table 9. In the upper part of the table, the first ten natural frequencies obtained for several structural models by using the GDQ method are compared to the 3D FEM solution (39420 20-node brick elements of Abaqus). An excellent agreement can be observed between the two different solutions, for each structural model. In particular, the EDZ3 theory leads to really close results. The same geometry and boundary conditions are employed for the next example. In this case, the three layers are made of orthotropic material (Graphite-Epoxy) and the lamination scheme is ( 30 / 65 / 45 ) . The first ten natural frequencies are shown in the lower part of Table 9, as well as the mechanical properties of the composite material. Even in this circumstance, the GDQ results are very similar to the solution given by the 3D FEM model (75336 20-node brick elements of Abaqus). The effect of the thickness variation on the mode shapes can be seen from the representations of Fig. 10. Finally, the doubly-curved translation panel made of two layers depicted in Fig. 7d is taken into account. Its shape is obtained by sliding a circumference over another circumference with the same geometric features. The reference surface position vector r (α1 , α 2 ) and the variable thickness expression h (α1 , α 2 ) are shown under the structure in Fig. 7d. The thickness profile is a sine wave variation applied in both the principal directions with the same characteristics. As in the previous case, the geometric configuration is kept for two different lamination schemes. The first analysis deals with the (Aluminum/Zirconia) stacking sequence. As usual, the mechanical properties of the constituents are listed in the upper part of Table 10, where the first ten natural frequencies for the completely clamped doubly-curved panel are reported too. The GDQ solutions are similar to the frequencies that can be achieved by a 3D FEM model (34496 20-node brick elements of Abaqus). The same considerations about the use of the shear correction factor illustrated before can be repeated also for this case, with connection to the structural theories defined by a first order kinematic expansion. Finally, the ( 0 / 90 ) lamination scheme is considered for the last analysis. Each ply is now made of Graphite-Epoxy, whose mechanical properties can be found in the lower part of Table 10. Even in this case, all the values obtained by the GDQ method for various structural theories are concordant to the 3D FEM ones (74880 20-node brick elements of Abaqus). The first nine mode shapes for a CCCC panel of translation with ( 0 / 90 ) as lamination scheme with the same sinusoidal thickness variation along the two principal directions are depicted in Fig. 11.

22

CONCLUSION All the numerical applications performed in the present work have proven that the GDQ method and different higher-order ESL structural theories give excellent results in terms of natural frequencies for plates, singlycurved shells and doubly-curved shells with variable thickness. The accuracy of the results has not been affected neither by any variable thickness profiles (such as linear, parabolic, cubic or sine wave) nor by the various lamination schemes. The outcome of the various analyses has been great for both isotropic and orthotropic layers. The validity of the present approach has been confirmed by means of the comparison with many results available in the literature and the solutions that were calculated by solving the dynamic problem through a threedimensional FE model. Thus, it is possible to conclude that the GDQ method is an efficient numerical tool which can deal easily with the free vibration analysis of laminated composite shells with variable thickness with a limited number of degrees of freedom if compared to 3D elasticity FEM.

ACKNOWLEDGMENTS The research topic is one of the subjects of the Centre of Study and Research for the Identification of Materials and Structures (CIMEST)-“M. Capurso” of the University of Bologna (Italy).

23

REFERENCES [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26]

Love A.E.H., A treatise on the mathematical theory of elasticity, Dover, 1944. Kraus H., Thin Elastic Shells, John Wiley & Sons, 1967. Tornabene F., Fantuzzi N., Bacciocchi M., Viola E., Strutture a Guscio in Materiale Composito. Geometria Differenziale. Teorie di Ordine Superiore, Esculapio, 2015. Tornabene F., Fantuzzi N., Bacciocchi M., Viola E., Strutture a Guscio in Materiale Composito. Quadtratura Differenziale e Integrale. Elementi Finiti in Forma Forte, Esculapio, 2015. Leissa A.W., Vibrations of shells, NASA SP-288, US Government Printing Office, 1973. Liew K.M., Wang C.M., Xiang Y., Kitipornchai S., Vibration of Mindlin Plates, Elsevier, 1998. Qatu M.S., Vibration of Laminated Shells and Plates, Elsevier, 2004. Reddy J.N., Mechanics of laminated composite plates and shells, CRC Press, 2004. Leissa A.W., Qatu M.S., Vibrations of Continuous Systems, McGraw-Hill, 2011. Tornabene F., Fantuzzi N., Mechanics of Laminated Composite Doubly-Curved Shell Structures. The Generalized Differential Quadrature Method and the Strong Formulation Finite Element Method, Esculapio, 2014. Viola E., Tornabene F., Fantuzzi N., General Higher-Order Shear Deformation Theories for the Free Vibration Analysis of Completely Doubly-Curved Laminated Shells and Panels, Compos. Struct., 2013, 95, 639-666. Tornabene F., Viola E., Fantuzzi N., General higher-order equivalent single layer theory for free vibrations of doubly-curved laminated composite shells and panels, Compos. Struct., 2013, 104, 94-117. Viola E., Tornabene F., Fantuzzi N., Generalized differential quadrature finite element method for cracked composite structures of arbitrary shape, Compos. Struct., 2013, 106, 815-834. Tornabene F., Fantuzzi N., Viola E., Ferreira A.J.M., Radial Basis Function method applied to doublycurved laminated composite shells and panels with a general higher-order equivalent single Layer theory, Compos. Part B-Eng., 2013, 55, 642-659. Fantuzzi N., Tornabene F., Strong formulation finite element method for arbitrarily shaped laminated plates - I. Theoretical analysis, Adv. Aircraft Space. Sci., 2014, 1, 124-142. Fantuzzi N., Tornabene F., Strong formulation finite element method for arbitrarily shaped laminated plates - II. Numerical analysis, Adv. Aircraft Space. Sci., 2014, 1, 143-173. Ferreira A.J.M., Carrera E., Cinefra M., Viola E., Tornabene F., Fantuzzi N., Zenkour A.M., Analysis of thick isotropic and cross-ply laminated plates by generalized differential quadrature method and a unified formulation, Compos. Part B-Eng., 2014, 58, 544-552. Tornabene F., Fantuzzi N., Bacciocchi M., The Local GDQ Method Applied to General Higher-Order Theories of Doubly-Curved Laminated Composite Shells and Panels: the Free Vibration Analysis, Compos. Struct., 2014, 116, 637-660. Fantuzzi N., Tornabene F., Viola E., Ferreira A.J.M., A Strong Formulation Finite Element Method (SFEM) based on RBF and GDQ techniques for the static and dynamic analyses of laminated plates of arbitrary shape, Meccanica, 2014, 49, 2503-2542. Fantuzzi N., Bacciocchi M., Tornabene F., Viola E., A.J.M. Ferreira, Radial Basis Functions Based on Differential Quadrature Method for the Free Vibration of Laminated Composite Arbitrary Shaped Plates, Compos. Part B-Eng., 2015, 78, 65-78. Tornabene F., Fantuzzi N., Bacciocchi M., Viola E., Accurate Inter-Laminar Recovery for Plates and Doubly-Curved Shells with Variable Radii of Curvature Using Layer-Wise Theories, Compos. Struct., 2015, 124, 368-393. Tornabene F., Fantuzzi N., Bacciocchi M., Viola E., A New Approach for Treating Concentrated Loads in Doubly-Curved Composite Deep Shells with Variable Radii of Curvature, Compos. Struct., 2015, 131, 433-452. Tornabene F., Fantuzzi N., Bacciocchi M., Dimitri R., Dynamic Analysis of Thick and Thin Elliptic Shell Structures Made of Laminated Composite Materials, Compos. Struct., 2015, 133, 278-299. Tornabene F., Fantuzzi N., Bacciocchi M., Dimitri R., Free Vibrations of Composite Oval and Elliptic Cylinders by the Generalized Differential Quadrature Method, Thin-Wall. Struct., 2015, 97, 114-129. Tornabene F., Brischetto S., Fantuzzi N., Viola E., Numerical and Exact Models for Free Vibration Analysis of Cylindrical and Spherical Shell Panels, Compos. Part B-Eng., 2015, 81, 231-250. Tornabene F., General Higher Order Layer-Wise Theory for Free Vibrations of Doubly-Curved Laminated Composite Shells and Panels, Mech. Adv. Mat. Struct., 2015. In Press.

24

[27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51]

Tornabene F., Fantuzzi N., Viola E., Inter-Laminar Stress Recovery Procedure for Doubly-Curved, Singly-Curved, Revolution Shells with Variable Radii of Curvature and Plates Using Generalized HigherOrder Theories and the Local GDQ Method, Mech. Adv. Mat. Struct., 2015. In Press. Tornabene F., Fantuzzi N., Bacciocchi M., Viola E., Higher-Order Theories for the Free Vibration of Doubly-Curved Laminated Panels with Curvilinear Reinforcing Fibers by Means of a Local Version of the GDQ Method, Compos. Part B-Eng., 2015, 81, 196-230. Eftekhari S.A., Jafari A.A., Accurate variational approach for free vibration of simply supported anisotropic rectangular plates, Arch. Appl. Mech., 2014, 84, 607-614. Ribeiro P., Akhavan H., Teter A., Warminski J., A review on the mechanical behaviour of curvilinear fibre composite laminated panels, J. Compos. Mater., 2014, 48, 2761-2777. Groh R.M.J., Weaver P.M., Buckling analysis of variable angle tow, variable thickness panels with transverse shear effects, Compos. Struct., 2014, 107, 482-493. Coburn B.H., Wu Z., Weaver P.M., Buckling analysis of stiffened variable angle tow panels, Compos. Struct., 2014, 111, 259-270. Hossein M.A., Dozio L., Prediction of natural frequencies of laminated curved panels using refined 2-D theories in the spectral collocation method, Curved Layer. Struct., 2015, 2, 1-14. Eisenträger J., Naumenko K., Altenbach H., Meenen J., A user-defined finite element for laminated glass panels and photovoltaic modules based on a layer-wise theory, Compos. Struct., 2015, 133, 265-277. Sarmila S., Free vibration behavior of laminated composite stiffened elliptic parabolic shell panel with cutout, Curved Layer. Struct., 2015, 2, 162-182. Yazdani S., Ribeiro P., A layerwise p-version finite element formulation for free vibration analysis of thick composite laminates with curvilinear fibres, Compos. Struct., 2015, 120, 531-542. Kurtaran H., Geometrically nonlinear transient analysis of moderately thick laminated composite shallow shells with generalized differential quadrature method, Compos. Struct., 2015, 125, 605-614. Raju G., Wu Z., Weaver P.M., Buckling and postbuckling of variable angle tow composite plates under in-plane shear loading, Int. J. Solids Struct., 2015, 58, 270-287. Sofiyev A.H., Buckling of heterogeneous orthotropic composite conical shells under external pressures within the shear deformation theory, Compos. Part B-Eng., 2016, 84, 175-187. Kulikov G.M., Mamontov A.A., Plotnikova S.V., Mamontov S.A., Exact geometry solid-shell element based on a sampling surface technique for 3D stress analysis of doubly-curved composite shells, Curved Layer. Struct., 2016, 3, 1-16. Tornabene F., Reddy J.N., FGM and laminated doubly-curved and degenerate shells resting on nonlinear elastic foundations: a GDQ solution for static analysis with a posteriori stress and strain recovery, J. Indian Inst. Sci., 2013, 93, 635-688. Tornabene F., Viola E., Static analysis of functionally graded doubly-curved shells and panels of revolution, Meccanica, 2013, 48, 901-930. Tornabene F., Fantuzzi N., Bacciocchi M., Free vibrations of free-form doubly-curved shells made of functionally graded materials using higher-order equivalent single layer theories, Compos. Part B Eng., 2014, 67, 490-509. Viola E., Rossetti L., Fantuzzi N., Tornabene F., Static analysis of functionally graded conical shells and panels using the generalized unconstrained third order theory coupled with the stress recovery, Compos. Struct., 2014, 112, 44-65. Tornabene F., Fantuzzi N., Viola E., Batra R.C., Stress and strain recovery for functionally graded freeform and doubly-curved sandwich shells using higher-order equivalent single layer theory, Compos. Struct., 2015, 119, 67-89. Fantuzzi N., Tornabene F., Viola E., Four-parameter functionally graded cracked plates of arbitrary shape: a GDQFEM solution for free vibrations, Mech. Adv. Mat. Struct., 2016, 23, 89-107. Quan T.Q., Tran P., Tuan N.D., Duc N.D., Nonlinear dynamic analysis and vibration of shear deformable eccentrically stiffened S-FGM cylindrical panels with metal-ceramic-metal layers resting on elastic foundations, Compos. Struct., 2015, 126, 16-33. Sofiyev A.H., Buckling analysis of freely-supported functionally graded truncated conical shells under external pressures, Compos. Struct., 2015, 132, 746-758. Sofiyev A.H., Influences of shear stresses on the dynamic instability of exponentially graded sandwich cylindrical shells, Compos. Part B-Eng., 2015, 77, 349-362. Sofiyev A.H., Kuruoglu N., Dynamic instability of three-layered cylindrical shells containing an FGM interlayer, Thin-Wall. Struct., 2015, 93, 10-21. Sofiyev A.H., Kuruoglu N., Parametric instability of shear deformable sandwich cylindrical shells containing an FGM core under static and time dependent periodic axial loads, Int. J. Mech. Sci., 2015,

25

[52] [53] [54] [55] [56]

[57] [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75] [76]

101-102, 114-123. Behravan Rad A., Shariyat M., Three-dimensional magneto-elastic analysis of asymmetric variable thickness porous FGM circular plates with non-uniform tractions and Kerr elastic foundations, Compos. Struct., 2015, 125, 558-574. Hosseini-Hashemi S., Abaei A.R., Ilkhani M.R., Free vibrations of functionally graded viscoelastic cylindrical panel under various boundary conditions, Compos. Struct., 2015, 126, 1-15. Alibeigloo A., Thermo elasticity solution of sandwich circular plate with functionally graded core using generalized differential quadrature method, Compos. Struct., 2016, 136, 229-240. Mantari J.L., Refined and generalized hybrid type quasi-3D shear deformation theory for the bending analysis of functionally graded shells, Compos. Part B-Eng., 2015, 83, 142-152. Sofiyev A.H., Hui D., Valiyev A.A., Kadioglu F., Turkaslan S., Yuan G.Q., Kalpakci V., Özdemir A., Effects of Shear Stresses and Rotary Inertia on the Stability and Vibration of Sandwich Cylindrical Shells with FGM Core Surrounded by Elastic Medium, Mech. Based Des. Struct., 2015, In Press. Doi: 10.1080/15397734.2015.1083870. Sofiyev A.H., Kuruoglu N., Domains of dynamic instability of FGM conical shells under time dependent periodic loads, Compos. Struct., 2016, 136, 139-148. He X.Q., Eisenberger M., Liew K.M., The effect of van der Waals interaction modeling on the vibration characteristics of multiwalled carbon nanotubes, J. Appl. Phys., 2006, 100, 124317. Cinefra M., Carrera E., Brischetto S., Refined shell models for the vibration analysis of multiwalled carbon nanotubes, Mech. Adv. Mater. Struc., 2011, 18, 476-483. Brischetto S., A continuum elastic three-dimensional model for natural frequencies of single-walled carbon nanotubes, Compos. Part B-Eng., 2014, 61, 222-228. Brischetto S., A continuum shell model including van der Waals interaction for free vibrations of doublewalled carbon nanotubes, CMES Comp. Model Eng., 2015, 104, 305-327. Liew K.M., Lei Z.X., Zhang L.W., Mechanical analysis of functionally graded carbon nanotube reinforced composites: a review, Compos. Struct., 2015, 120, 90-97. Zhang L.W., Liew K.M., Reddy J.N., Postbuckling of carbon nanotube reinforced functionally graded plateswith edges elastically restrained against translation and rotation underaxial compression, Comput. Method Appl. Mech. Engrg., 2016, 298, 1-28. Hadjiloizi D. A., Kalamkarov A.L., Metti Ch., Georgiades A. V., Analysis of Smart Piezo-MagnetoThermo-Elastic Composite and Reinforced Plates: Part I – Model Development, Curved Layer. Struct., 2014, 1, 11-31. Hadjiloizi D. A., Kalamkarov A.L., Metti Ch., Georgiades A. V., Analysis of Smart Piezo-MagnetoThermo-Elastic Composite and Reinforced Plates: Part II – Applications, Curved Layer. Struct., 2014, 1, 32-58. Mizusawa T., Vibration of rectangular Mindlin plates with tapered thickness by the spline strip method, Comput. Struct., 1993, 46, 451-463. Shufrin I., Eisenberger M., Vibration of shear deformable plates with variable thickness – first-order aand higher-order analyses, J. Sound Vib., 2006, 290, 465-489. Dozio L., Carrera E., A variable kinematic Ritz formulation for vibration study of quadrilateral plates with arbitrary thickness, J. Sound Vib., 2011, 330, 4611-4632. Eisenberger M., Jabareen M., Axisymmetric vibrations of circular and annular plates with variable thickness, Int. J. Struct. Stab. Dy., 2001, 1, 195-206. Efraim E., Eisenberger M., Exact vibration analysis of variable thickness thick annular isotropic and FGM plates, J. Sound Vib., 2007, 299, 720-738. Wu T.Y., Liu G.R., Free vibration analysis of circular plates with variable thickness by the generalized differential quadrature rule, Int. J. Solids Struct., 2001, 38, 7967-7980. Liang B., Zhang S.-F., Chen D.-Y., Natural frequencies of circular annular plates with variable thickness by a new method, Int. J. Pres. Ves. Pip., 2007, 84, 293-297. Xu Y., Zhou D., Three-dimensional elasticity solution of functionally graded rectangular plates with variable thickness, Compos. Struct., 2009, 91, 56-65. Hosseini-Hashemi S., Taher H.R.D., Akhavan H., Vibration analysis of radially FGM sectorial plates of variable thickness on elastic foundations, Compos. Struct, 2010, 1734-1743. Tajeddini V., Ohadi A., Sadighi M., Three-dimensional free vibration of variable thickness thick circular and annular isotropic and functionally graded plates on Pasternak foundation, Int. J. Mech. Sci., 2011, 53, 300-308. Lal R., Rani R., On radially symmetric vibrations of circular sandwich plates of non-uniform thickness, Int. J. Mech. Sci., 2015, 99, 29-39.

26

[77] [78] [79] [80] [81] [82] [83] [84] [85] [86] [87] [88] [89] [90] [91] [92] [93] [94] [95] [96] [97] [98]

Kang J.-H., Leissa A.W., Three-dimensional vibrations of thick spherical shell segments with variable thickness, Int. J. Solids Struct., 2000, 37, 4811-4823. Kang J.-H., Leissa A.W., Free vibration analysis of complete paraboloidal shells of revolution with variable thickness and solid paraboloids from a three-dimensional theory, Comput. Struct., 2005, 83, 2594-2608. Leissa A.W., Kang J.-H., Three-Dimensional Vibration Analysis of Paraboloidal Shells, JSME Int. J. C, 2002, 45, 2-7. Efraim E., Eisenberger M., Dynamic stiffness vibration analysis of thick spherical shell segments with variable thickness, J. Mech. Mater. Sruct., 2010, 5, 821-835. Efraim E., Free Vibrations of Thick Segmented Shells of Revolution, PhD Thesis, Techion - Israel Institute of Technology, Department of Civil and Environmental Engineering, 2006. Kang J.-H., Leissa A.W., Corrigendum to “Three-dimensional vibrations of thick spherical shell segments with variable thickness”, Int. J. Solids Struct., 2006, 43, 2848-2851. Duan W.H., Koh C.G., Axisymmetric transverse vibrations of circular cylindrical shells with variable thickness, J. Sound Vib., 2008, 317, 1035-1041. El-Kaabazi N., Kennedy D., Calculation of natural frequencies and vibration modes of variable thickness cylindrical shells using the Wittrick-Williams algorithm, Comput. Struct., 2012, 104-105, 4-12. Ghannad M., Rahimi G.H., Nejad M.Z., Elastic analysis of pressurized thick cylindrical shells with variable thickness made of functionally graded materials, Compos. Part B-Eng., 2013, 45, 388-396. Naghsh A., Saadatpour M.M., Azhari M., Free vibration analysis of stringer stiffened general shells of revolution using a meridional finite strip method, Thin-Wall. Struct., 2015, 94, 651-662. Jabareen M., Eisenberger M., Free vibrations of non-homogeneous circular and annular membranes, J. Sound Vib., 2001, 240, 409-429. Jiang W., Redekop D., Static and vibration analysis of orthotropic toroidal shells of variable thickness by differential quadrature, Thin-Wall. Struct.. 2003, 41, 461-478. Tornabene F., Fantuzzi N., Ubertini F., Viola E., Strong formulation finite element method based on differential quadrature: a survey, Appl. Mech. Rev., 2015, 67, 020801 (55 pages). Artioli E., Gould P., Viola E., A Differential Quadrature Method Solution for Shear-Deformable Shells of Revolution, Eng. Struct., 2005, 27, 1879-1892. Viola E., Tornabene F., Ferretti E., Fantuzzi N., Soft core plane state structures under static loads using GDQFEM and Cell method, CMES, 2013, 94, 301-329. Viola E., Tornabene F., Ferretti E., Fantuzzi N., GDQFEM numerical simulations of continuous media with cracks and discontinuities, CMES, 2013, 94, 331-368. Viola E., Tornabene F., Ferretti E., Fantuzzi N., On static analysis of plane state structures via GDQFEM and Cell method, CMES, 2013, 94, 421-458. Fantuzzi N., New insights into the strong formulation finite element method for solving elastostatic and elastodynamic problems, Curved Layer. Struct., 2014, 1, 93-126. Fantuzzi N., Tornabene F., Viola E., Generalized differential quadrature finite element method for vibration analysis of arbitrarily shaped membranes, Int. J. Mech. Sci., 2014, 79, 216-251. Tornabene F., Fantuzzi N., Bacciocchi M., The strong formulation finite element method: stability and accuracy, Fract. Struct. Integr., 2014, 29, 251-265. Viola E., Miniaci M., Fantuzzi N., Marzani A., Vibration Analysis of Multi-Stepped and Multi-Damaged Parabolic Arches Using GDQ, Curved Layer. Struct., 2015, 2, 28-49. Viola E., Tornabene F., Fantuzzi N., Bacciocchi M., DiQuMASPAB Software, DICAM Department, Alma Mater Studiorum - University of Bologna (http://software.dicam.unibo.it/diqumaspab-project).

27

Figure 1. Representation of a generic doubly-curved shell element with variable thickness.

Figure 2. Coordinate system and layer numbering and orientation employed for a laminated shell.

28

a)

b)

c)

d)

e) Figure 3. GDQ discretization and local co-ordinate system representation for several plate with different thickness profiles: a) linear variation in x-direction; b) linear variation in both directions; c) parabolic arched form; d) parabolic concave form; e) parabolic symmetric concave form.

29

a) Hemispherical annular shell

b) Spherical barrel shell

r (ϕ , ϑ ) = R sin ϕ cos ϑ e1 − R sin ϕ sin ϑ e 2 + R (1 − cos ϕ ) e3

r (ϕ , ϑ ) = R sin ϕ cos ϑ e1 − R sin ϕ sin ϑ e 2 + R (1 − cos ϕ ) e3

ϕ ∈ [ϕ0 , ϕ1 ] = [π 6, π 2 ] , θ ∈ [θ 0 , θ1 ] = [0, 2π ]

ϕ ∈ [ϕ0 , ϕ1 ] = [π 4 , 3π 4 ] , θ ∈ [θ0 , θ1 ] = [ 0, 2π ]

  ϕ − ϕ0   h (ϕ ) = h0  1 + 2      ϕ1 − ϕ0   

  ϕ − ϕ0   h (ϕ ) = h0  1 + 2      ϕ1 − ϕ 0    c) Paraboloidal shell

r (ϕ , ϑ ) = 2 R tan ϕ cosϑ e1 − 2 R tan ϕ sin ϑ e 2 + R tan ϕ e 3

ϕ ∈ [ϕ0 , ϕ1 ] = [π 6, π 3] , θ ∈ [θ 0 , θ1 ] = [0, 2π ]  3  ϕ − ϕ0   h (ϕ ) = h0 1 +   2 ϕ − ϕ    1 0   Figure 4. GDQ discretization and local co-ordinate system representation for three different doubly-curved shells with linearly variable thickness: a) hemispherical annular shell; b) spherical barrel shell; c) paraboloidal shell. The corresponding variable sections are shown on the right of the shells. For each structure is shown the position vector r (ϕ , ϑ ) which describes the reference surface and the thickness variation analytical expression.

30

Ω = 0.3614

Ω = 0.9175

Ω = 1.3927

Ω = 1.5417

n=2

n=3

n =1

n = 0 (A)

Ω = 1.5617

Ω = 2.2744

Ω = 2.4712

n=4

n=5

n = 0 (T)

Figure 5. First mode shape for each value of n of a FF barrel shell with linearly variable thickness.

Ω = 0.0672

Ω = 0.1808

Ω = 0.3308

Ω = 0.4427

n=2

n=3

n=4

n = 0 (A)

Ω = 0.4756

Ω = 0.5107

Ω = 1.1314

n =1

n=5

n = 0 (T)

Figure 6. First mode shape for each value of n of a completely free paraboloidal shell with linearly variable thickness.

31

a) Conical panel

b) Cylinder with elliptic profile

(singly-curved panel of revolution)

(singly-curved panel of translation)

(

)

r (α1 , y ) = a cos α1 + y sin α sin ( arctan ( − cot α1 ) ) e1 +

  1  e3 e1 − y e 2 + b  1 − 2 2   1 + k tan ϕ 1 + k tan ϕ   ϕ ∈ [ϕ0 , ϕ1 ] = [ −π 4 , π 4] , y ∈ [ y0 , y0 ] = 0, Ly  , k = a b

r (ϕ , y ) =

− y cos α e 2 +

(

)

+ a sin α1 − y sin α cos ( arctan ( − cot α1 )) e3

α1 ∈ α , α  = [0, 2π 3] , y ∈ [ y0 , y1 ] = 0, Ly  a = 0.5 m, Ly = 2 m, α = π 6 0 1

2

2

a = 1m, b = 0.5 m, Ly = 2 m

1 1

  ϕ − ϕ0   π   , h0 = 0.1m h (ϕ ) = h0  1 + sin   ϕ1 − ϕ 0   

  y − y 0 3  h ( y ) = h0  1 +    , h0 = 0.1m   y1 − y0    

c) Panel with catenary meridian (doubly-curved panel of revolution)

d) Doubly-curved translational panel (a circumference slides over another circumference)

r (α1 , α 2 ) = ( Rα0 1 (α1 ) − x3α2 (α 2 ) sin α1 ) e1 +

r (ϕ , ϑ ) = d arcsinh ( tan ϕ ) cos ϑ e1 +

− d arcsinh ( tan ϕ ) sin ϑ e 2 +

(

ak tan ϕ

− Rα0 2 (α 2 ) e 2 + ( x3α1 (α1 ) + x3α2 (α 2 ) cos α1 ) e 3

)

+ d cosh ( arcsinh ( tan ϕ )) − 1 e3

R 0 (α1 ) = a sin α1 ,

x3α1 (α1 ) = a (1 − cos α1 )

ϕ ∈ [ϕ0 , ϕ1 ] = [π 12 , π 4] , θ ∈ [θ 0 ,θ1 ] = [0, π 2] ,

Rα0 2 (α 2 ) = a sin α 2 ,

x3α 2 (α 2 ) = a (1 − cos α 2 ) , a = 1m

 3  ϕ − ϕ0  1  θ − θ0   h (ϕ , ϑ ) = h0  1 +  π    − sin  2 − 2 ϕ ϕ  1 0  θ1 − θ0    d = 1m, h0 = 0.045 m

 π π  π π α1 ∈ α10 , α11  = − ,  , α 2 ∈ α 20 , α 21  = − ,   4 4  4 4

α1

 1  α −α0 π h (α1 , α 2 ) = h0 1 − sin  2 11 10 π +  + 4 2  α1 − α1  1  α − α 20 π  − sin  2 21 π +   , h0 = 0.1m 4  α 2 − α 20 2 

Figure 7. Four panel structures with variable thickness representation and GDQ grid definition. For each structure is shown the position vector r (α1 , α 2 ) , which describes the reference surface, the thickness variation analytical expression and the local coordinates reference system.

32

1st mode

2nd mode

3rd mode

4th mode

5th mode

6th mode

7th mode

8th mode

9th mode

Figure 8. First nine mode shapes of a FCCC conical panel (Aluminum/Steel) with cubic thickness variation along y-direction.

33

1st mode

2nd mode

3rd mode

4th mode

5th mode

6th mode

7th mode

8th mode

9th mode

Figure 9. First nine mode shapes of a FCFC cylinder of translation

( −45 / 45 )

defined by an elliptic profile with sinusoidal thickness

variation along ϕ -direction.

1st mode

2nd mode

3rd mode

4th mode

5th mode

6th mode

7th mode

8th mode

9th mode

Figure 10. First nine mode shapes of a FCFF panel of revolution ( 30 / 65 / 45 ) with catenary meridian with linear thickness variation along

ϕ -direction and sinusoidal thickness variation along ϑ -direction.

34

1st mode

2nd mode

3rd mode

4th mode

5th mode

6th mode

7th mode

8th mode

9th mode

Figure 11. First nine mode shapes of a CCCC panel of translation ( 0 / 90 ) obtained by sliding a circumference over another circumference

with the same sinusoidal thickness variation along both directions.

35

h0 Ly

0.1

Theory

1,2

2,1

2,2

1,3

3,1

3,2

2,3

Shufrin and Eisenberger [67]

1.4504

3.4743

3.5058

5.4838

6.5345

6.7038

8.5303

8.5921

Mizusawa [66]

1.4504

3.4743

3.5058

5.4840

6.5347

6.7039

8.5302

-

κ =5 6

1.4507

3.4758

3.5074

5.4876

6.5394

6.7093

8.5390

8.6009

κ =π 2 12

1.4504

3.4743

3.5058

5.4838

6.5345

6.7038

8.5303

8.5921

GDQ - ED2

1.4505

3.4744

3.5063

5.4846

6.5352

6.7054

8.5324

8.5937

GDQ - ED3

1.4519

3.4820

3.5135

5.5028

6.5596

6.7310

8.5741

8.6375

GDQ - FSDTRS

GDQ - FSDTRS

κ =π 2 12

0.2

Shufrin and Eisenberger [67]

1.3738

3.1096

3.1276

4.6613

5.4883

5.5657

6.8435

6.8725

Mizusawa [66]

1.3738

3.1096

3.1276

4.6613

5.4881

5.5656

6.8437

6.8726

κ =5 6 GDQ - FSDTRS

1.3747

3.1139

3.1320

4.6703

5.4996

5.5778

6.8607

6.8900

1.3738

3.1096

3.1276

4.6613

5.4883

5.5657

6.8435

6.8725

GDQ - ED2

1.3740

3.1100

3.1290

4.6631

5.4896

5.5688

6.8469

6.8751

GDQ - ED3

1.3787

3.1324

3.1500

4.6923

5.5509

5.6310

6.9400

6.9721

Shufrin and Eisenberger [67]

1.1664

2.3603

2.3637

3.2845

3.7942

3.8050

4.5043

4.5105

Mizusawa [66]

1.1665

2.3603

2.3637

3.2845

3.7942

3.8050

4.5043

4.5105

1.1687

2.3675

2.3710

3.2969

3.8093

3.8204

4.5241

4.5304

1.1664

2.3603

2.3637

3.2845

3.7942

3.8050

4.5043

4.5105

GDQ - ED2

1.1669

2.3607

2.3657

3.2861

3.7952

3.8077

4.5065

4.5120

GDQ - ED3

1.1788

2.4040

2.4058

3.3619

3.8488

3.8926

4.6229

4.6368

κ =π 2 12

GDQ - FSDTRS

κ =π 2 12

0.4

Mode ( m, n ) 1,1

κ =5 6 GDQ - FSDTRS κ =π GDQ - FSDTRS

2

κ =π 2 12

12

Table 1. Comparison of frequency factor λ = 2π −1 fLy 12 ρ (1 −ν 2 ) Eh02 for a SSSS square plate with linearly variable thickness in xdirection. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = I M = 25 for the GDQ solution.

36

Var.

h0 Ly

1,1

1,2

2,1

2,2

1,3

3,1

3,2

2,3

3,3

Shufrin and Eisenberger [67]

2.2208

4.3495

7.3887

4.3771

6.3217

9.3160

7.5108

9.2594

11.9744

κ =5 6 GDQ - FSDTRS

2.2189

4.3368

7.3818

4.3836

6.3135

9.3277

7.5108

9.2281

11.9580

Shufrin and Eisenberger [67]

1.9951

3.6870

5.9211

3.6943

5.1208

7.1690

5.9556

7.1420

8.9131

GDQ - FSDTRS

1.9929

3.6747

5.8991

3.6996

5.1127

7.1708

5.9687

7.1192

8.8992

Shufrin and Eisenberger [67]

1.7458

3.0609

4.7226

3.0587

4.1284

5.6032

4.7272

5.5906

6.8489

GDQ - FSDTRS

1.7443

3.0640

4.6906

3.0495

4.1217

5.6039

4.7523

5.5742

6.8402

0.1

Shufrin and Eisenberger [67]

1.8870

3.7253

6.3773

3.7253

5.4755

8.1202

6.3773

8.1202

10.5425

κ =5 6 GDQ - FSDTRS

1.8849

3.6823

6.3660

3.7605

5.4588

8.1949

6.3841

8.0143

10.5039

0.2

Shufrin and Eisenberger [67]

1.7410

3.2764

5.3441

3.2764

4.6231

6.5491

5.3441

6.5491

8.2249

GDQ - FSDTRS

1.7385

3.2446

5.3194

3.2983

4.6084

6.5872

5.3586

6.4783

8.1958

Shufrin and Eisenberger [67]

1.5654

2.8060

4.3987

2.8060

3.8389

5.2674

4.3987

5.2674

6.4883

1.5633

2.7835

4.3696

2.8193

3.8266

5.2891

4.4184

5.2192

6.4685

0.1 βx=0.25, βy=0.50

Mode ( m, n )

Theory

0.2

κ =5 6

0.3

βx =0.50, βy =0.50

κ =5 6

0.3

κ =5 6

κ =5 6

GDQ - FSDTRS

Table 2. Comparison of the frequency factor λ = 2π −1 fLy 12 ρ (1 −ν 2 ) Eh02 for a CCCC square plate with linear thickness variation in both directions. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = I M = 25 for the GDQ solution.

37

Var.

h0 Ly

1,1

1,2

2,1

2,2

1,3

3,1

3,2

2,3

3,3

Shufrin and Eisenberger [67]

2.7420

5.2909

8.9859

5.4058

7.6223

10.9760

9.0766

11.1258

14.0809

κ =5 6 GDQ - FSDTRS

2.7392

5.2879

8.9747

5.4007

7.6143

10.9679

9.0768

11.1120

14.0674

Shufrin and Eisenberger [67]

2.3434

4.1952

6.6405

4.2697

5.7661

7.9029

6.7123

7.9491

9.7694

GDQ - FSDTRS

2.3410

4.1934

6.6183

4.2640

5.7588

7.8964

6.7243

7.9360

9.7593

Shufrin and Eisenberger [67]

1.9623

3.3237

5.0699

3.3712

4.4605

5.9663

5.1128

5.9835

7.2572

GDQ - FSDTRS

1.9609

3.3228

5.0446

3.3665

4.4551

5.9617

5.1322

5.9756

7.2518

0.1

Shufrin and Eisenberger [67]

2.2835

4.5231

7.7170

4.3327

6.4168

9.3996

7.2995

9.3726

12.1072

κ =5 6 GDQ - FSDTRS

2.2818

4.5207

7.7157

4.3298

6.4107

9.3922

7.2945

9.3631

12.0957

0.2

Shufrin and Eisenberger [67]

2.0494

3.7933

6.0592

3.6893

5.1880

7.2173

5.8872

7.2188

8.9821

GDQ - FSDTRS

2.0476

3.7915

6.0624

3.6851

5.1813

7.2066

5.8766

7.2109

8.9709

Shufrin and Eisenberger [67]

1.7915

3.1256

4.7835

3.0742

4.1753

5.6344

4.7155

5.6393

6.8881

GDQ - FSDTRS

1.7902

3.1245

4.7925

3.0699

4.1693

5.6296

4.7006

5.6301

6.8812

0.1

Shufrin and Eisenberger [67]

2.4367

4.7909

7.9973

4.2182

6.4917

9.5571

6.9687

9.2326

12.1207

κ =5 6 GDQ - FSDTRS

2.4356

4.7894

7.9971

4.2163

6.4882

9.5497

6.9645

9.2266

12.1120

2.1852

3.9874

6.2143

3.6437

5.2693

7.3133

5.7265

7.2047

9.0165

0.2

Shufrin and Eisenberger [67] GDQ - FSDTRS

2.1842

3.9864

6.2154

3.6414

5.2665

7.3096

5.7211

7.1992

9.0183

Shufrin and Eisenberger [67]

1.9087

3.2658

4.8745

3.0746

4.2482

5.6962

4.6488

5.6635

6.9251

1.9081

3.2654

4.8762

3.0724

4.2461

5.6940

4.6645

5.6600

6.9213

0.1 Arched form

Mode ( m, n )

Theory

0.2

κ =5 6

0.3

Concave form

κ =5 6

Symmetric concave form

0.3

0.3

κ =5 6

κ =5 6

κ =5 6

κ =5 6

GDQ - FSDTRS

(

Table 3. Comparison of the frequency factor λ = 2π −1 fLy 12 ρ 1 −ν 2

)

Eh02 for a CCCC square plate with parabolic thickness variation in

x-directions. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = I M = 25 for the GDQ solution.

38

n

0 (A)

0 (T)

1

2

3

4

5

2D FE 1600 S4R shell elements 10080 DOF [80] 1.5357 2.1152 2.3992 4.3117 5.9909 3.7144 6.5151 9.5088 1.5666 2.4016 2.4824 4.4149 4.6268 0.3495 0.7135 1.8514 2.8071 3.6503 0.9382 1.6815 2.6617 3.3486 4.8367 1.6823 2.6672 3.6246 4.0553 5.8216 2.5319 3.6775 4.5940 4.9368 6.6094

3D FE 2095 C3D20R solid elements 43020 DOF [80] 1.5361 2.0862 2.4216 4.3143 5.8511 3.6775 6.3843 9.1976 1.5687 2.3730 2.4819 4.4072 4.5782 0.3495 0.7031 1.8563 2.8044 3.5714 0.9368 1.6615 2.6422 3.3397 4.7495 1.6757 2.6368 3.5577 4.0540 5.8013 2.5140 3.6304 4.4686 4.9436 6.6536

Kang and Leissa [82] 3D, Ritz

DSM with 2 exact elements 15 DOF [80]

Present solution GDQ κ =5 ( 6−ν ) FSDTRS

1.536 2.090 2.423 4.322 5.844 3.669 6.385 9.200 1.568 2.374 2.485 4.414 4.572 0.349 0.705 1.857 2.809 3.572 0.934 1.666 2.646 3.346 4.746 1.671 2.643 3.564 4.060 5.793 2.507 3.637 4.475 4.951 6.665

1.5347 2.1138 2.4057 4.3041 5.9780 3.7053 6.4586 9.3248 1.5645 2.3990 2.4845 4.4101 4.6131 0.3494 0.7136 1.8443 2.8093 3.6372 0.9364 1.6794 2.6436 3.3494 4.8121 1.6753 2.6577 3.5878 4.0547 5.8063 2.5144 3.6558 4.5297 4.9250 6.6166

1.5298 2.1012 2.4116 4.3249 5.8810 3.6740 6.3950 9.2147 1.5685 2.3837 2.4786 4.4135 4.5802 0.3491 0.7041 1.8622 2.8173 3.5730 0.9364 1.6649 2.6540 3.3556 4.7492 1.6758 2.6451 3.5755 4.0681 5.8018 2.5150 3.6439 4.4912 4.9541 6.6585

Table 4. Comparison of the non-dimensional natural frequencies Ω = 2π fR ρ G for the FF hemispherical annular shell with linearly

variable thickness shown in Fig. 4a. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = 31, I M = 51 for the GDQ solution.

39

n

0 (A)

0 (T)

1

2

3

4

5

3D FE 2095 C3D20R solid elements 43020 DOF [80]

Kang and Leissa [82] 3D, Ritz

DSM with 3 exact elements 20 DOF [80]

Present solution GDQ κ =5 ( 6 −ν ) FSDTRS

1.5345 1.9099 1.9632 2.4772 3.8292 2.4774 4.2877 6.1689 1.3915 1.7690 2.3793 2.5836 3.4497 0.3593 0.6065 1.6842 2.0238 2.8902 0.9085 1.4032 2.2845 2.6510 3.3737 1.5436 2.1947 2.9024 3.5713 4.0506 2.2465 2.9445 3.6991 4.4644 4.9457

1.535 1.911 1.963 2.487 3.844 2.470 4.290 6.172 8.093 1.390 1.771 2.379 2.594 3.445 0.362 0.604 1.690 2.027 2.901 0.918 1.400 2.298 2.651 3.385 1.561 2.196 2.911 3.579 4.059 2.271 2.954 3.707 4.477 4.952

1.5299 1.9268 1.9651 2.4846 3.8355 2.4938 4.3343 6.2411 8.1862 1.3898 1.7688 2.4037 2.5915 3.4897 0.3613 0.6076 1.6977 2.0281 2.9006 0.9169 1.3997 2.3155 2.6687 3.3824 1.5597 2.1917 2.9192 3.6300 4.0507 2.2696 2.9478 3.7014 4.5276 4.9700

1.5417 1.9086 1.9549 2.4947 3.8505 2.4712 4.2929 6.1767 8.0934 1.3927 1.7778 2.3769 2.6030 3.4454 0.3614 0.6034 1.6950 2.0378 2.9121 0.9175 1.4011 2.3062 2.6620 3.3957 1.5617 2.2014 2.9177 3.5908 4.0663 2.2744 2.9608 3.7095 4.4896 4.9563

Table 5. Comparison of the non-dimensional natural frequencies Ω = 2π fR ρ G for the FF barrel shell with linearly variable thickness

shown in Fig. 4b. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = 31, I M = 51 for the GDQ solution.

40

n

0 (A)

0 (T)

1

2

3

4

5

2D FE 1800 S8R shell elements 3348 DOF [81]

3D FE 1440 C3D20 solid elements 31590 DOF [81]

DSM with 3 exact elements 20 DOF [81]

Present solution GDQ κ =5 ( 6 −ν ) FSDTRS

0.4408 0.6214 0.8359 1.0170 1.1308 0.4744 0.6844 0.7547 1.0492 1.5066 0.0670 0.2149 0.6010 0.7170 1.1217 0.1800 0.4775 0.8272 0.9051 1.2833 0.3290 0.7178 1.0697 1.1942 1.4978 0.5073 0.9597 1.3633 1.4729 1.7816

0.4433 0.6185 0.8346 1.0152 1.1312 0.4753 0.6823 0.7537 1.0480 1.5041 0.0671 0.2156 0.6018 0.7165 1.1173 0.1807 0.4789 0.8299 0.9071 1.2875 0.3305 0.7207 1.0749 1.2015 1.5063 0.5100 0.9655 1.3720 1.4883 1.7950

0.4420 0.6193 0.8345 1.0149 1.1338 0.4746 0.6827 0.7547 1.0477 1.5087 0.0669 0.2154 0.6005 0.7183 1.1197 0.1805 0.4792 0.8296 0.9074 1.2866 0.3302 0.7215 1.0745 1.2008 1.5042 0.5095 0.9661 1.3717 1.4886 1.7914

0.4427 0.6196 0.8337 1.0153 1.1314 0.4756 0.6831 0.7535 1.0483 1.5047 0.0672 0.2153 0.6025 0.7177 1.1181 0.1808 0.4793 0.8304 0.9078 1.2876 0.3307 0.7220 1.0749 1.2020 1.5054 0.5106 0.9671 1.3722 1.4889 1.7931

Table 6. Comparison of the non-dimensional natural frequencies Ω = 2π fR ρ G for the FF paraboloidal shell with linearly variable

thickness shown in Fig. 4c. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = 31, I M = 51 for the GDQ solution.

41

Lamination scheme: (Aluminum/Steel), h1 = 0.05 m , h2 = 0.05 m . Aluminum: E = 70 GPa, ν = 0.3,

ρ = 2707 kg m 3 .

Steel: E = 210 GPa, ν = 0.3, ρ = 7800 kg m 3 . f [Hz]

κ =5 6 FSDTRS

FSDTZκRS=5 6

ED2κ =5 6

ED2

ED3

EDZ3

ED4

1 2 3 4 5 6 7 8 9 10

505.739 589.018 754.753 775.822 820.092 998.533 1012.434 1049.787 1076.337 1218.420

503.244 586.586 749.709 771.959 815.865 992.905 1003.778 1037.191 1067.710 1210.838

505.840 588.449 753.689 775.186 818.961 996.592 1009.339 1044.084 1072.800 1214.998

508.879 591.508 760.043 779.984 824.389 1003.676 1020.366 1060.288 1084.075 1224.910

507.985 590.666 758.166 778.681 822.870 1001.664 1017.057 1055.429 1080.726 1221.958

505.400 588.483 753.895 775.590 818.995 997.054 1009.557 1044.217 1073.153 1215.983

507.127 589.932 756.440 777.438 821.523 999.982 1014.263 1051.351 1077.891 1219.681

3D FEM 838359 DOF (C3D20 Abaqus) 506.010 588.860 753.951 775.643 819.484 997.245 1009.849 1044.877 1073.425 1215.772

Table 7. First ten natural frequencies for a FCCC conical panel with cubic thickness variation along y-direction, for several ESL theories. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = I M = 31 .

Lamination scheme: ( −45 / 45 ) , h1 = 0.05 m , h2 = 0.05 m . Graphite-Epoxy: E1 = 137.9 GPa , E2 = E3 = 8.96 GPa , G12 = G13 = 7.1 GPa , G23 = 6.21 GPa , ν 12 = ν 13 = 0.3 , ν 23 = 0.49 ,

ρ = 1450 kg m 3 . f [Hz]

κ =5 6 FSDTRS

FSDTZκRS=5 6

ED1κRS=5 6

ED1RS

ED3

EDZ3

ED4

1 2 3 4 5 6 7 8 9 10

345.199 361.898 372.516 489.917 530.994 686.366 736.640 805.493 837.194 979.873

343.414 358.891 368.220 484.690 524.192 677.076 723.464 796.377 824.092 957.569

342.280 356.426 371.811 482.462 521.177 679.189 718.688 803.927 830.684 964.076

344.486 358.839 378.306 485.993 527.672 687.061 729.298 819.391 846.933 982.569

343.773 359.390 369.613 485.520 526.052 682.516 731.204 803.370 831.471 973.639

340.042 355.627 367.607 481.341 521.989 678.579 725.433 798.852 826.157 962.845

343.219 358.417 367.678 483.631 523.587 679.179 726.477 799.709 826.785 965.476

3D FEM 1372977 DOF (C3D20 Abaqus) 342.036 357.363 366.794 479.987 523.961 675.414 724.020 795.709 824.051 958.537

Table 8. First ten natural frequencies for a FCFC cylinder of translation defined by an elliptic profile with sinusoidal thickness variation along ϕ -direction, for several ESL theories. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = I M = 31 .

42

Lamination scheme: ( Aluminum/Zirconia/Aluminum ) , h1 = 0.015 m , h2 = 0.015 m , h3 = 0.015 m . Aluminum: E = 70 GPa, ν = 0.3,

ρ = 2707 kg m 3 .

Zirconia: E = 168GPa, ν = 0.3, ρ = 5700 kg m3 . f [Hz]

κ =5 6 FSDTRS

FSDTZκRS=5 6

ED2κ =5 6

EDZ2κ =5 6

ED3

1 360.124 360.111 360.857 360.827 361.262 2 394.739 394.685 395.555 395.465 396.972 3 783.840 783.682 785.576 785.253 789.855 4 922.800 922.685 923.464 923.200 926.748 5 949.631 949.405 951.985 951.504 958.301 6 1491.039 1490.670 1491.103 1490.144 1501.818 7 1511.516 1511.084 1513.018 1512.076 1525.356 8 1714.465 1714.384 1716.696 1716.499 1719.222 9 1732.301 1731.759 1736.198 1734.782 1751.955 10 1917.820 1917.164 1922.526 1920.880 1942.105 Lamination scheme: ( 30 / 65 / 45 ) , h1 = 0.015 m , h2 = 0.015 m , h3 = 0.015 m .

EDZ3

ED4

361.100 396.228 787.130 925.182 954.870 1495.210 1518.142 1717.875 1742.946 1930.556

361.228 396.908 789.708 926.631 958.071 1501.653 1525.109 1719.121 1751.537 1941.657

3D FEM 524088 DOF (C3D20 Abaqus) 361.040 396.098 787.182 924.872 954.581 1495.359 1517.971 1717.930 1742.743 1930.544

Graphite-Epoxy: E1 = 137.9 GPa , E2 = E3 = 8.96 GPa , G12 = G13 = 7.1 GPa , G23 = 6.21 GPa , ν 12 = ν 13 = 0.3 , ν 23 = 0.49 ,

ρ = 1450 kg m 3 . f [Hz]

κ =5 6 FSDTRS

FSDTZκRS=5 6

ED2κ =5 6

EDZ2κ =5 6

ED3

EDZ3

ED4

1 2 3 4 5 6 7 8 9 10

324.775 431.345 742.367 783.448 1030.797 1213.818 1346.235 1497.123 1626.781 1708.800

324.061 429.149 737.569 777.585 1021.910 1206.911 1341.382 1483.561 1609.517 1690.239

322.505 431.044 740.806 782.225 1028.081 1209.866 1347.377 1494.923 1623.728 1705.827

321.669 428.832 736.324 776.250 1019.720 1203.585 1342.766 1482.304 1608.222 1688.969

323.233 432.425 744.046 784.678 1033.748 1215.289 1350.431 1502.478 1633.487 1714.925

321.822 428.995 742.634 777.454 1027.722 1214.116 1349.827 1496.585 1629.190 1705.260

322.743 431.953 743.116 783.903 1032.852 1213.724 1348.546 1501.049 1632.196 1713.430

3D FEM 978735 DOF (C3D20 Abaqus) 322.241 431.778 742.541 783.569 1031.940 1213.181 1348.310 1500.121 1630.542 1712.089

Table 9. First ten natural frequencies for a FCFF panel of revolution with catenary meridian with linear thickness variation along ϕ -direction and sinusoidal thickness variation along ϑ -direction, for several ESL theories. The Chebyshev-Gauss-Lobatto grid

distributions is employed with I N = I M = 31 .

43

Lamination scheme: ( Aluminum/Zirconia ) , h1 = 0.5 m , h2 = 0.05 m . Aluminum: E = 70 GPa, ν = 0.3,

ρ = 2707 kg m 3 .

Zirconia: E = 168GPa, ν = 0.3, ρ = 5700 kg m3 . f [Hz]

κ =5 6 FSDTRS

FSDTZκRS=5 6

1 948.621 946.455 2 954.979 952.647 3 1032.569 1030.623 4 1286.812 1281.922 5 1416.846 1407.872 6 1572.215 1564.999 7 1618.172 1608.297 8 1637.866 1628.892 9 1779.694 1768.247 10 1797.125 1788.569 Lamination scheme: ( 0 / 90 ) , h1 = 0.5 m ,

ED2κ =5 6

EDZ2κ =5 6

ED3

EDZ3

ED4

946.735 953.134 1029.975 1282.750 1412.008 1571.066 1611.792 1632.253 1773.668 1793.317 h2 = 0.05 m .

945.044 951.272 1025.983 1279.331 1407.080 1566.451 1607.365 1627.430 1766.289 1785.438

948.656 955.176 1031.831 1287.530 1420.985 1578.098 1621.998 1641.347 1785.838 1802.322

945.834 952.059 1026.476 1282.459 1411.613 1570.386 1613.451 1632.502 1773.141 1790.202

947.756 954.208 1031.232 1285.694 1417.515 1575.319 1618.416 1638.205 1781.731 1799.418

3D FEM 463185 DOF (C3D20 Abaqus) 946.793 953.244 1030.497 1283.536 1413.336 1572.083 1613.960 1634.313 1776.407 1795.716

Graphite-Epoxy: E1 = 137.9 GPa , E2 = E3 = 8.96 GPa , G12 = G13 = 7.1 GPa , G23 = 6.21 GPa , ν 12 = ν 13 = 0.3 , ν 23 = 0.49 ,

ρ = 1450 kg m 3 . f [Hz]

κ =5 6 FSDTRS

FSDTZκRS=5 6

ED1κRS=5 6

ED1RS

ED3

EDZ3

ED4

1 2 3 4 5 6 7 8 9 10

864.252 982.909 1013.509 1085.913 1357.628 1376.021 1420.052 1525.651 1608.613 1637.363

860.491 976.899 1008.057 1078.487 1347.831 1362.105 1403.402 1509.219 1580.294 1615.188

862.875 960.210 1002.574 1068.549 1340.688 1342.628 1393.042 1496.162 1605.174 1634.851

868.049 967.033 1009.881 1077.194 1356.002 1356.825 1411.160 1513.169 1636.286 1666.770

864.887 975.796 1012.130 1082.592 1355.038 1365.787 1416.267 1527.844 1607.351 1643.673

861.051 971.898 1009.911 1076.303 1348.348 1356.078 1413.344 1522.609 1599.294 1637.549

863.467 973.402 1009.869 1080.017 1351.158 1360.546 1409.986 1521.584 1597.148 1635.278

3D FEM 982911 DOF (C3D20 Abaqus) 862.743 972.360 1008.838 1078.818 1349.474 1358.461 1407.265 1518.824 1593.092 1631.854

Table 10. First ten natural frequencies for a CCCC panel of translation obtained by sliding a circumference over another circumference with the same sinusoidal thickness variation along both directions, for several ESL theories. The Chebyshev-Gauss-Lobatto grid distributions is employed with I N = I M = 31 .

44