Vibration characteristics of rotating pretwisted composite tapered blade with graphene coating layers

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.1 (1-20)

Aerospace Science and Technology ••• (••••) ••••••

1

Contents lists available at ScienceDirect

1 2

2 3

Aerospace Science and Technology

4

3 4 5

5

6

6

www.elsevier.com/locate/aescte

7

7

8

8

9

9

10

10

11 12 13

Vibration characteristics of rotating pretwisted composite tapered blade with graphene coating layers

16 17 18 19 20

W. Zhang

a,∗

a

, Y. Niu , K. Behdinan

b

a

Beijing Key Laboratory of Nonlinear Vibrations and Strength of Mechanical Structures, College of Mechanical Engineering, Beijing University of Technology, Beijing 100124, PR China b Advanced Research Laboratory for Multifunctional Light Weight Structures, Department of Mechanical and Industrial Engineering, University of Toronto, 5 King’s College Road, Toronto, ON M5S 3G8, Canada

a r t i c l e

i n f o

a b s t r a c t

25 26 27

Article history: Received 27 October 2019 Received in revised form 2 December 2019 Accepted 11 December 2019 Available online xxxx

28 29 30 31 32

15 16 17 18 19 20 22 23

23 24

13

21

21 22

12 14

14 15

11

Keywords: Pretwisted tapered cylindrical panel Graphene coating layer Rotating speed Frequency veering

33 34 35 36 37 38 39

A new dynamic model of the rotating tapered cantilever cylindrical panel with the graphene coating layers is developed to investigate the vibration characteristics of the rotating pretwisted tapered blade. It is assumed that the graphene platelets (GPLs) are randomly oriented and uniformly dispersed in the top layer and the bottom layer of the rotating pretwisted composite tapered blade. The modified Halpin-Tsai model is used to estimate the effective Young’s modulus. The rule of the mixture is used to calculate the effective Poisson’s ratio and mass density. Based on the Green strain tensor, an accurate straindisplacement relationship is acquired. The effects of the centrifugal force and Coriolis force are considered in the formulation. The Chebyshev-Ritz method is utilized to obtain the natural frequencies and mode shapes of the rotating pretwisted composite tapered blade with the graphene coating layers. The accuracy of the proposed model is validated through several comparison studies with the results of the present literatures and ANSYS. The free vibration characteristics are analyzed by considering different material and geometry parameters of the rotating pretwisted composite tapered cantilever cylindrical panel with the graphene coating layers, such as the graphene platelet (GPL) geometry, GPL weight fraction, taper ratio, length-to-radius ratio, pretwist angle, presetting angle and rotating speed. The frequency veering and the mode shape shift phenomena are found in the rotating pretwisted tapered cantilever cylindrical panel with the graphene coating layers. © 2019 Elsevier Masson SAS. All rights reserved.

24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39

40

40

41

41

42

42

43 44

1. Introduction

47 48 49 50 51 52 53 54 55 56 57 58 59

Rotating blades are the critical components in many practical engineering, such as aero-engines, turbines and helicopters. The aeroengines rely on the numerous blades to compress the airflow, and then generate the strong thrust to achieve the efficient work. During the aforementioned work, the rotating blades will produce the strong vibration, which may lead to the structural failure. Thus, the vibration characteristics of the rotating blades should be fully investigated and accurately estimated. The beam models are usually used to simplify the rotating blades since the beam models are the simplest model. Song and Librescu [1] provided a refined theory to investigate the free vibration characteristics of the thin-walled composite beam. Considering the effects of the gyroscopic and centrifugal forces, Yang et al. [2,3] investigated the free vibrations of the rotating Bernoulli-Euler and Timoshenko beams. Wang and Zhang [4] analyzed the stability of the rotating blade for both linear and nonlinear dynamic models. Considering the harmonic gas pressure and pre-deformation, Zhang and Li [5] discussed the nonlinear internal resonances of the rotating pretwisted blade. Yao et al. [6,7] studied the nonlinear dynamics and the steady-state vibration responses of the pre-twisted, presetting, thin-walled rotating cantilever beam with varying rotating speed. Chen et al. [8] developed an analytical approach to analyze the free vibration characteristics and damping effect of the blisk with the damping hard coating on the blades. Kandil et al. [9,10] utilized the positive position feedback and the proportional-derivative controllers to suppress the vibrations of the rotating blade. Rafiee et al. [11] presented a review of the vibration and control for the rotating beams and blades.

64 65 66

47 48 49 50 51 52 53 54 55 56 57 58 59 61

61 63

46

60

60 62

44 45

45 46

43

*

Corresponding author. E-mail addresses: [email protected] (W. Zhang), [email protected] (Y. Niu), [email protected] (K. Behdinan).

https://doi.org/10.1016/j.ast.2019.105644 1270-9638/© 2019 Elsevier Masson SAS. All rights reserved.

62 63 64 65 66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.2 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

2

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 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43

Since the beam models are the one-dimensional model, they are unable to predict the chordwise bending vibration and the coupled spanwise bending-chordwise bending vibration [12], which are important in the rotating blades with small aspect ratio. Therefore, the plate models are developed to establish the analytical dynamic models of the rotating blades. Yoo et al. [13,14] analyzed the modal characteristics of the rotating composite cantilever plate and discussed the frequency loci veering and crossing. Cao et al. [15] developed a dynamic model of the sandwich plate with the thermal barrier coating layers to investigate the aero-engine turbine blade. Yao et al. [16] focused on the nonlinear dynamic behaviors of the rotating plate considering the warping effect of the cross-section. Considering the hollow circular cross-section, Dong et al. [17] used the assumed modes method and Lagrange’s equations to derive the dynamic equations of motion for the rotating functionally graded tapered cantilever beams and analyzed the vibration characteristics of the system. The aforementioned literatures are related to the plate models. Considering the pretwist angle of the rotating blades, the panel models are developed. Sun et al. [18] established a shell model to investigate the vibration for the rotating blade with the pretwist angle and stagger angle. Sinha and Zylka [19] studied the free vibrations of the rotating pretwisted airfoil by considering both the stress-softening and stress-stiffening effects. Liu et al. [20] utilized the Rayleigh-Ritz method to investigate the free vibration of the rotating functionally graded pretwisted sandwich plate. Gu et al. [21] discussed the free vibration of the rotating pretwisted panel with the initial geometric imperfection. Chen and Li [22] focused on dynamic model characteristics of the rotating composite laminated blade. The configuration of the rotating blades is complex and is usually cambered along the chordwise direction, pretwisted along the spanwise direction and tapered along the thickness direction. Hence, the cylindrical panel and conical panel models are more accurate to describe the dynamics of the rotating blades. Hu et al. [23,24] used the Rayleigh-Ritz method to investigate the free vibrations of the rotating pretwisted cylindrical panel and conical panel models. Yao et al. [25] analyzed the nonlinear dynamics of the rotating pretwisted cylindrical shell subjected to different loadings. In order to promote the property and efficiency of the rotating blades, advanced materials and novel structures should be used. The coating layers are widely used in the rotating blades to improve the operating temperature and corrosion resistance. Thermal barrier coating (TBC) materials [15,26] were the most common coating layers, which can resist the high temperature. As a new two dimensional carbon nanomaterial with the highest strength and good toughness, the graphene had been widely used to reinforce the matrix due to its excellent mechanical, thermal, electronic, physical and chemical properties [27,28]. In recent years, many researchers have found that a low content of graphene or its derivatives can promote the mechanical properties of the composites remarkably. Rafiee et al. [29] carried out the experimental study and reported that the Young’s modulus of the graphene reinforced composites is about 31% greater than the pure epoxy with the addition of 0.1% weight fraction (wt%) of the graphene platelets (GPLs). Das et al. [30] utilized the nano-indentation technique to demonstrate that the elastic modulus and hardness increase significantly when 0.6 wt% of the graphene is added. In addition to the experimental studies, the extensive theoretical studies have also illustrated that the graphene has a significant enhancement effect on the elastic modulus and hardness of the structures. Lin et al. [31] utilized the molecular dynamic simulation to perform the mechanical properties of nanocomposites reinforced with the graphene. Moradi-Dastjerdi and Behdinan [32] discuss the thermoelastic static and free vibrations of the graphene-reinforced nanocomposite cylinders. Many scholars have endeavored to investigate dynamic characteristics of the graphene reinforced nanocomposite structures, including bending [33], buckling [34] and postbuckling [35], nonlinear vibrations [36], transient analysis [37] and nonlinear internal resonances [38,39] of the beam [34], plate [40], trapezoidal plate [33], cylindrical panel [41] and cylindrical shell [42]. Due to the complexity of the configuration for the rotating blades, few studies have been focused on the rotating pretwisted tapered cylindrical panel. To the author’s knowledge, it is the first time to investigate the vibration characteristics of a rotating pretwisted tapered cylindrical panel with the graphene coating layers. In this paper, a new dynamic model of the rotating pretwisted tapered cylindrical panel with the graphene coating layers is developed to investigate the free vibrations of the rotating composite tapered blade. The ChebyshevRitz method is used to obtain the natural frequencies and mode shapes of the rotating composite tapered cylindrical panel. The accuracy of the proposed model is validated. Then, the effects of the GPL geometry, GPL weight fraction, taper ratio, length-to-radius ratio, pretwisted angle, presetting angle and rotating speed on the vibration characteristics of the rotating composite tapered cylindrical panel are discussed through the comprehensive parametric studies. The frequency veering and mode shape shift phenomena are found in the system.

48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 45

2. Theoretical formulation

46

46 47

2

44

44 45

1

Fig. 1 indicates a rotating blade of the aero-engine, which is simplified to a rotating pretwisted composite tapered cantilever cylindrical panel. The blade is mounted on a rigid rotating cylindrical hub with the radius R 0 , which rotates with a varying rotating speed (t ) = 0 + f cos 1 t. 0 and f cos 1 t represent the steady-state and the periodic perturbation rotating speeds, respectively. The presetting angle is φ and the pretwist angle at the tip end of the blade is K = kL, where k is the uniform pretwist ratio. The blade is assumed to be the length L, the radius R, the subtended angle  and the thickness h. Fig. 2 portrays the cross section of the rotating composite tapered blade along the blade length. The thickness at an arbitrary cross section is expressed by

 h(x) = h0 1 − ch

x L



47 48 49 50 51 52 53 54

,

(1)

where h0 is the thickness of the blade root and ch is the thickness taper ratio. In order to establish the governing equations of motion for the rotating pretwisted composite tapered cylindrical panel, the following coordinate systems are introduced. The inertial reference Cartesian coordinate system ( X 1 , Y 1 , Z 1 ) is located at the center of the rotating hub O 0 . The rotating coordinate system (x0 , y 0 , z0 ) is located at the blade root O with the unit vectors (i1 , i2 , i3 ), where x0 -axis is along the spanwise direction, y 0 -axis is along the chordwise direction and z0 -axis is along the radial direction, respectively. The curvilinear coordinate system (x, s, z) is defined on the mid-surface of the cylindrical panel with the center O 1 , where x-axis, s-axis and z-axis are along the spanwise, circumferential and radial directions, respectively. θ represents the angle measured from the z axis. e gives the distance between the center of the cylindrical panel O 1 and the blade root O . The rotating pretwisted composite tapered cantilever cylindrical panel is composed of the core, the top and bottom layers, which are the graphene platelet reinforced composite (GPLRC) layers. The thickness of the core, top layer and bottom layer at the blade root are hc0 ,

55 56 57 58 59 60 61 62 63 64 65 66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.3 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

3

1

1

2

2

3

3

4

4

5

5

6

6

7

7

8

8

9

9

10

10

11

11

12

12

13

13

14

14

15

15

16

16

17

17

18

18

19

19

20

20

21

21

22

22

23 24

Fig. 1. The configuration of a rotating pretwisted composite tapered cantilever cylindrical panel is portrayed: (a) the model of the rotating pretwisted composite tapered blade, (b) the coordinate systems of the rotating pretwisted composite tapered blade, (c) the cross section of the rotating pretwisted composite tapered blade.

23 24

25

25

26

26

27

27

28

28

29

29

30

30

31

31

Fig. 2. The cross sections of the rotating pretwisted tapered cylindrical panel with different taper ratios are given.

32

32 33 34 35 36 37

h G P L0 and h G P L0 , respectively. The GPLRCs are considered as isotropic. It is assumed that the GPLs are randomly oriented and uniformly dispersed in the top layer and the bottom layer of the rotating pretwisted composite tapered cylindrical panel. The GPL volume fraction in each layer for the rotating pretwisted composite tapered cylindrical panel are described by the following components. At the bottom layer, we have the GPL volume fraction

38

∗ V GPL = V GPL . (1 )

39 40 41

48

∗ V GPL =

51 52 53

56 57 58

61

64 65 66

(3)

,

EM +

5 1 + ξ W η W V GPL 8

1 − η W V GPL

50 51 52 53 54

EM,

(4)

55 56 57 58

ηL =

E GPL −1 EM E GPL + ξL EM



62 63

W GPL + (1 − W GPL )(ρGPL /ρ M )

8 1 − η L V GPL

46 48

where

59 60

42

49

W GPL

3 1 + ξ L η L V GPL

E eff =

55

39

47

where W GPL and ρGPL are the weight fraction and the mass density of the GPLs, respectively, ρ M is the mass density of the matrix. Based on the modified Halpin-Tsai model, the effective Young’s modulus of the GPLRC is calculated by

54

38

45

(2c)

∗ is the volume fraction of the GPLs where V GPL

50

37

44

(3 )

49

36

43

∗ V GPL = V GPL ,

46 47

(2b)

At the top layer, we give the GPL volume fraction

45

35

41

(2 )

V GPL = 0.

43

34

40

At the core, we obtain the GPL volume fraction

42 44

(2a)

33

ξL = 2

lGPL hGPL

,

ηW =



,

ξW

E GPL EM

−1

E GPL EM

+ ξW   w GPL , =2 hGPL

59

,

(5a)

60 61 62

(5b)

where E GPL , hGPL , lGPL and w GPL are the Young’s modulus, thickness, length and width of the GPLs, respectively, E M is the Young’s modulus of the matrix.

63 64 65 66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.4 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

4

1

The effective Poisson’s ratio

νeff and the effective mass density ρeff are expressed as follows

1 2

2 3 4 5 6 7 8 9

νeff = νGPL V GPL + νM V M ,

(6a)

ρeff = ρGPL V GPL + ρM V M ,

(6b)

where νGPL and ρGPL are the Poisson’s ratio and mass density of the GPLs, ν M and ρ M are the Poisson’s ratio and mass density of the matrix, respectively. The position vector of an arbitrary point on the mid-surface of the rotating pretwisted composite tapered cylindrical panel before the deformation is expressed as follows ( 0)

r0 = xi1 + R sin θ i2 + ( R cos θ − e )i3 .

11 12 14

di2

16 18

dx di3

19

dx

17

21 22

= −ki2 .

(8b)

27 28 29 30

( 0)

39

45

g2 =

47 48 50

53

62 63 64 65 66

19 20 21 22 24

27 29

z B



sin θ i2 +

R cos θ − e +

z B

 cos θ i3 .

1

1

B

B

33 34 35 36 37

(11)

38 39 40 41 42 43

(12a)

44 45 46

(12b)

47 48 49 50

sin θ i2 +

1 B

51

cos θ i3 ,

(12c)

52 53 54

a11 =

57

61

R sin θ +



where

55

60

B



ke sin θ i1 +

= − ke sin θ i1 +

52

59

z



∂ r(0) = g 21 i1 + g 22 i2 + g 23 i3 R ∂θ = za21 i1 + (cos θ + za22 )i2 − (sin θ + za23 )i3 ,

51

58

(10)

∂ r(0) = g 31 i1 + g 32 i2 + g 33 i3 g3 = ∂z

49

18

30

∂ r(0) g1 = = g 11 i1 + g 12 i2 + g 13 i3 ∂x = i1 − k( R cos θ − e + za11 )i2 + k( R sin θ + za12 )i3 ,

46

16

28

(9c)

The base vectors g1 , g2 and g3 are respectively denoted by the partial derivatives of r(0) with respect to the x, s and z

44

14

32

r(0) = r0 + za3 = x −

43

13

26

(9b)

1 + (ke sin θ)2 .



42

11

25

The position vector of an arbitrary point outside the mid-surface of the rotating pretwisted composite tapered cylindrical panel before the deformation is expressed by the following equation

38

56

(9a)



37

54

9

31

B=

33

41

8

23

where

32

40

7

17

∂r a1 = 0 = i1 + k(e − R cos θ)i2 + kR sin θ i3 , ∂x ∂ r(0) a2 = 0 = cos θ i2 − sin θ i3 , R ∂θ a1 × a2 1 a3 = = (−ke sin θ i1 + sin θ i2 + cos θ i3 ), |a1 × a2 | B

26

36

(8a)

( 0)

25

35

6

15

= ki3 ,

The corresponding base vectors a1 and a2 are denoted by the partial derivatives of r0 with respect to x and s, respectively. Thus, the unit vector a3 is perpendicular to a1 and a2 . The aforementioned base vectors a1 , a2 and a3 are derived by

24

34

5

12

(0)

23

31

(7)

Since the rotating pretwisted composite tapered cylindrical panel has a uniform pretwist rate k along the x axis, i2 and i3 are the rotating vectors around the x axis. Frenet-Serret formula is used to obtain the following expressions

15

20

4

10

10

13

3

a22 =

1

k cos θ,

B 1

BR

cos θ −

a12 =

1 B

k sin θ,

k2 e 2 sin2 θ cos θ B3 R

,

a21 = −

1

ke cos θ +

k3 e 3 sin2 θ cos θ

BR B3 R 2 2 1 k e sin θ cos2 θ a23 = . sin θ + BR B3 R

55

,

56 57

(13)

58 59 60

The displacement vector of the rotating pretwisted composite tapered cylindrical panel is described as follows

61

(14)

62

where the u, v and w denote the displacement components along the directions of a1 , a2 and a3 , respectively. The position vector of an arbitrary point outside the mid-surface of the rotating pretwisted composite tapered cylindrical panel after the deformation is given as follows

64

U = ua1 + va2 + wa3 ,

63 65 66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.5 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

r = R 0 i1 + r(0) + U

1



2

z

4 5 6 7 8

B

9 10

13 14 15 16 17 18

2 f i j = Gi G j − gi g j ,

21

24

27 28 29 30

2e i j =

34 35

39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61

3  3  k =1 l =1

36 38

5 6 7

(15)

∂ αk ∂ αl 2 f kl , ∂β i ∂β j

10 11

(16a) (16b)

13 14 15 17 18 19

(17)

20 21 22 23 24 25

(18a)

26

(18b)

28

(18c)

27 29 30 31 32 33

(i , j = 1, 2, 3),

(19)

34 35 36

where

37

k

3 

∂α = g kj (ji · g j ), ∂β i j =1 ⎤−1 ⎡ g1 · g1 g1 · g2 g1 · g3 g kj = ⎣ g1 · g2 g2 · g2 g2 · g3 ⎦ . g1 · g3 g2 · g3 g3 · g3

εi j (i , j = ξ, η, ζ ) are determined by ⎧ ⎫ ⎧ ⎫ εξ ξ ⎪ ⎪ e11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ εηη ⎪ ⎬ ⎪ ⎬ ⎨ e 22 ⎪ γξ η = 2e12 , ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 2e 13 ⎪ ⎪ ⎪ ⎩ γξ ζ ⎪ ⎭ ⎪ ⎭ ⎩

38

(20a)

41 42

(20b)

43 44 45

The strain components

γη ζ

39 40

46 47 48

(21)

49 50 51

2e 23

52

where (εξ ξ , εηη , γξ η ) and (γξ ζ , γηζ ) are the membrane strains and transverse strains, respectively, and the strain components are given in the Appendix B. Compared to the classical theory, the first-order shear deformation theory, which introduces the shear deformation, has been used to investigate the vibration of the structures. Moreover, many studies have been focused on the higher-order shear deformation theory, which considers the shear deformation and the changes in thickness direction of the plate and shell. However, the difference in the frequencies between the first-order shear deformation theory and higher-order shear deformation theory is insignificant when the length to thickness ratio up to 5 [40]. Therefore, in this paper, we use the first-order shear deformation theory to investigate the vibration of the rotating pretwisted composite tapered cylindrical panel. The displacement field of the rotating pretwisted composite tapered cylindrical panel is considered as follows

53 54 55 56 57 58 59 60 61 62

63

u = u 0 + z ϕx ,

(22a)

64

v = v 0 + zϕθ ,

(22b)

w = w0,

(22c)

66

12

16

(16c)

62

65

8 9

According to a transformation law [25], the strain tensors e i j in the new coordinate system (ξ, η, ζ ) are calculated as follows

33

37

(i , j = 1, 2, 3),

 1 j1 = i1 + ke sin2 θ i2 + ke sin θ cos θ i3 , B j2 = cos θ i2 − sin θ i3 , 1 j3 = (−ke sin θ i1 + sin θ i2 + cos θ i3 ). B

26

32

4

where the Green strain tensors are given in the Appendix A. A local orthogonal coordinate system (ξ, η, ζ ) issuing from the point of the rotating pretwisted composite tapered cylindrical panel before the deformation is established. The corresponding unit vector is calculated as follows

25

31

3

The Green strain tensors in the coordinate system (x, s, z) are determined as

20

23

B

∂r G1 = = g1 + G 11 i1 + G 12 i2 + G 13 i3 , ∂x ∂r G2 = = g2 + G 21 i1 + G 22 i2 + G 23 i3 , R ∂θ ∂r G3 = = g3 + G 31 i1 + G 32 i2 + G 33 i3 . ∂z

12

22

2

The corresponding base vectors G1 , G2 and G3 are derived by

11

19

w

1



= R 0 + x − ke sin θ + u − ke sin θ i1 B B   z w + R sin θ + sin θ + uk(e − R cos θ) + v R cos θ + sin θ i2 B B   z w + R cos θ − e + cos θ + ukR sin θ − v R sin θ + cos θ i3 .

3

5

63 64 65 66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.6 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

6

1 2 3 4

where u 0 , v 0 and w 0 are the displacements of any point on the mid-surface of the rotating pretwisted composite tapered cylindrical panel along the x-axis, s-axis and z-axis, respectively, ϕx and ϕθ denote the rotations of the transverse normal about the s-axis and x-axis, respectively. Substituting equation (22) into equation (21), the strain-displacement relationships are represented as follows

⎧ ( 0) εξ ξ ⎪ ⎪ ⎪ ⎪ ( 0) ⎪ ⎪ ⎪ ⎨ εηη

5

⎧ ⎫ εξ ξ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ εηη ⎪ ⎬

6 7 8

γξ η = γξ(η0) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ γ ( 0) ⎩ γξ ζ ⎪ ⎭ ⎪ ⎪ ⎪ ξζ γη ζ ⎪ ⎩ ( 0) γη ζ

9 10 11 12 13 14

18 19 20 21 22 23 24 25 26 27 28 29 30 31

35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55

58

65 66

7 8

(23)

,

9 10 11 12

15 16 17

(24a)

18 19

(24b)

20 21 22 23 24 25

(24c)

26 27

(24d) (24e)

28 29 30 31 32

∂ ϕx ke 2 ∂ ϕx 1 − + 2 k2 e 2 sin θ cos θ ϕθ , ∂x R ∂θ B R ke 2 ∂ ϕx 1 ∂ ϕθ (1 ) εηη = + , R ∂θ R ∂θ   ke 2 ∂ ϕx 1 ∂ ϕθ 1 1 B 2 ∂ ϕx ke 2 ∂ ϕθ ke sin θ − e 2k2 cos2 θ + 2k2 e cos θ − k2 R + γξ(η1) = + + ϕθ , − − B ∂x B ∂x B R R ∂θ B R ∂θ BR ke cos θ γξ(ζ1) = 2 ϕθ ,

εξ(1ξ) =

B R

γη(ζ1)

=

k

1

ϕx −

ϕθ .

B BR The kth (k = 1, 2, 3) layer constitutive equation of the rotating pretwisted composite tapered cylindrical panel is written by

⎧ ⎫ σξ ξ ⎪(k) ⎪ ⎪ ⎪ ⎪ ⎨ σηη ⎪ ⎬

σξ η ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ σξ ζ ⎪ ⎭ ση ζ

⎡ ⎢ ⎢ =⎢ ⎣

Q 11 Q 12 0 0 0

Q 12 Q 22 0 0 0

(k) (k) (k) where (σξ ξ , σηη , σξ η ) and

0 0 Q 66 0 0

33

(25a) (25b)

0 0 ⎥ ⎥ 0 ⎥ ⎦ 0 k0 Q 44

γξ η ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ γξ ζ ⎪ ⎭ γη ζ

(k)

(k)

Q 11 = Q 22 =

1−ν

2 eff

,

(k)

Q 12 =

νeff E eff

(25d)

1−ν

2 eff

,

(k)

(25e)

44 46 48

(26)

,

49 50 51 52

(k)

Q 44 = Q 55 = Q 66 =

E eff

2(1 + νeff )

z L  2  3 k+1 1

(k) 2 ρeff r˙ MdzRdθ dx,

0 −  k=1 zk 2

where the detailed expression and coefficients are given in the Appendix C.

53 54 55 56

.

(27)

57 58 59 60 61



2

43 45

The kinetic energy of the rotating pretwisted composite tapered cylindrical panel is calculated by

K=

41 42

(k)

(k)

39

47

(k)

E eff

37

40

(k) (k) (σξ ζ , σηζ ) are the membrane stresses and transverse shear stresses, respectively. k0 denotes the shear cor(k)

(k)

36 38

(25c)

⎫ ⎤(k) ⎧ εξ ξ ⎪(k) ⎪ ⎪ ⎪ ⎪ ⎨ εηη ⎪ ⎬

0 0 0 k0 Q 55 0

34 35

rection coefficient and is taken to be 5/6, Q i j (i , j = 1, 2, 4, 5, 6) are the stiffness coefficients

62 64

6

and

61 63

4

14

59 60

3

13

56 57

⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭

2

5

∂ u 0 ke 2 ∂ u 0 1 1 = − + 2 k2 e 2 sin θ cos θ v 0 + 3 k2 e 2 e cos θ w 0 ∂x R ∂θ B R B R 2  2  1 1 ∂ w0 ke 2 ∂ w 0 ∂ w 0 1 k2 e 22 ∂ w 0 + − , + 2 B2 ∂x 2 B2 R2 ∂θ B 2 R ∂ x ∂θ 2  ke 2 ∂ u 0 1 ∂ v0 1 1 1 ∂ w0 ( 0) εηη = , w0 + + + R ∂θ R ∂θ BR 2 R2 ∂θ   ke 2 ∂ u 0 1 ∂ v0 1 1 B 2 ∂ u0 − e 2k2 cos2 θ + 2k2 e cos θ − k2 R + γξ(η0) = + + B ∂x B ∂x B R R ∂θ 2  ke 2 ∂ v 0 ke sin θ 2ke cos θ 1 ∂ w0 ∂ w0 ke 2 ∂ w 0 − + , − − v0 − w 0 B R ∂θ BR B R ∂ x ∂θ ∂θ B2 R B R2 1 ∂ w0 ke 2 ∂ w 0 ke cos θ v 0 + B ϕx , γξ(ζ0) = − + B ∂x B R ∂θ B2 R 1 ∂ w0 k 1 γη(ζ0) = + u0 − v 0 + ke 2 ϕx + ϕθ , R ∂θ B BR

17

34

(1 ) + z γξ η ⎪ ⎪ ⎪ ⎪ (1 ) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ γξ ζ ⎪ ⎪ ⎭ ⎩ (1 ) γη ζ

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬

εξ(0ξ)

16

33

⎧ (1 ) εξ ξ ⎪ ⎪ ⎪ ⎪ (1 ) ⎪ ⎪ ⎪ ⎨ εηη

where

15

32

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬

1

62

(28)

63 64 65 66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.7 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

1

The terms of the mass inertia are obtained by

2 3

(I0, I1, I2) =

4 6

U s1 =

14 15

1 2

16

8

(30)

( Ai j , B i j , D i j ) =



23

+ σηη εηη + σξ(kη) γξ η

+ σξ(kζ ) γξ ζ

 + ση(kζ) γηζ MdzRdθ dx,

13

(31)

(k) 

Q ij



17 18 20

1, z, z2 dz,

(i , j = 1, 2, 4, 5, 6).

(32)

29

B

30 31

23

F 2 = F · j2 = ρeff 

32 33

F 3 = F · j3 =

34

z





1

ρeff 

B

B



z



B

z B

24



25

cos(φ − θ) − e cos φ

26 27 28

(33a)

 

 +

ke sin θ

2

F=

ρeff 

z

R0 + x −

B





40

+ ρeff 2

R+

42

R+

z B

(33b)





 cos(φ − θ) − e cos φ cos(φ − θ) ,

 ke sin θ



 i1 +

2

ρeff 

R+







z

z





37 38 39

(34)

51 

L N1 =

N2 =

F 1 dx, x

1 =



1 2



∂ w0 ∂x

2

 +

52

2

F 2 dθ, θ

∂ v0 ∂x

(36a)

2  ,

2 =

1 2



∂ w0 ∂θ

2

 +

∂ u0 ∂θ

65 66

54 56

.

(36b)

57 58 59

3. Solution procedure

60 61

61

64

53 55

2 

59

63

48 50

56

62

47 49

where

55

60

44 46

(35)

0 −  k=1 zk 2

52

58

43 45

49

57

41 42



48

54

34

40

cos(φ − θ) − e cos φ cos φ i3 .

B



cos(φ − θ) − e cos φ sin φ i2

B

z L  2  3 k+1 U s2 = ( N 1 1 + N 2 2 ) MdzRdθ dx,

47

53

32

35

The potential energy caused by the centrifugal force for the rotating pretwisted composite tapered cylindrical panel is obtained as follows

46

51

31 33

(33c)

45

50

29 30

cos(φ − θ) − e cos φ sin(φ − θ)],

−ke sin θ R 0 + x −





37

41

R+

ke sin θ



36

where

39

z

2

35

38

R0 + x −

 2

21 22

+ ke sin2 θ sin φ R + B B B B    1 z cos(φ − θ) − e cos φ , + ke sin θ cos θ cos φ R +

28

14

19

1

F 1 = F · j1 = ρeff 2

11

16

 

27

44

(k)

The centrifugal force components are calculated by

26

10

15

z 3 k+1 

25

43

σξ(kξ ) εξ ξ

k=1 zk

22

9

12

zk+1 

where the detailed expression is given in the Appendix D. The extensional stiffness A i j , the bending-extensional coupling stiffness B i j and the bending stiffness D i j are given by

21

36

7

0 −  k=1 zk 2

20

4 6

L  2  3

19

24

(29)

5



13

18

3

(k) ρeff dz, (k = 1, 2, 3).

The potential energy induced by the deformation for the rotating pretwisted composite tapered cylindrical panel is given as follows

12

17

1, z 1 , z

U s = U s1 + U s2 .

9 11

2

 2

The total potential energy of the rotating pretwisted composite tapered cylindrical panel is expressed by

8 10

z 3 k+1  

1

k=1 zk

5 7

7

In this section, the Ritz method is utilized to calculate the natural frequencies and the corresponding mode shapes of the rotating pretwisted composite tapered cantilever cylindrical panel. The displacements components are expanded to the product of the mode shape functions with the exponent functions

u 0 (x, θ, t ) = U (x, θ)ei ωt ,

62 63 64 65

(37a)

66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.8 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

8

1 2 3 4 5 6 7 8 9 10 11 12

1

Table 1 The boundary functions of different boundary conditions are listed.

F-F F-S S-F S-S F-C C-F S-C C-S C-C

2

F u (¯x)

F v (¯x)

F w (¯x)

F ϕx (¯x)

F ϕθ (¯x)

¯ F u (θ)

F v (θ¯ )

F w (θ¯ )

¯ F ϕx (θ)

¯ F ϕθ (θ)

1 1 1 1 1 − x¯ 1 + x¯ 1 − x¯ 1 + x¯ 1 − x¯ 2

1 1 − x¯ 1 + x¯ 1 − x¯ 2 1 − x¯ 1 + x¯ 1 − x¯ 2 1 − x¯ 2 1 − x¯ 2

1 1 − x¯ 1 + x¯ 1 − x¯ 2 1 − x¯ 1 + x¯ 1 − x¯ 2 1 − x¯ 2 1 − x¯ 2

1 1 1 1 1 − x¯ 1 + x¯ 1 − x¯ 1 + x¯ 1 − x¯ 2

1 1 − x¯ 1 + x¯ 1 − x¯ 2 1 − x¯ 1 + x¯ 1 − x¯ 2 1 − x¯ 2 1 − x¯ 2

1 1 − θ¯ 1 + θ¯ 1 − θ¯ 2 1 − θ¯ 1 + θ¯ 1 − θ¯ 2 1 − θ¯ 2 1 − θ¯ 2

1 1 1 1 1 − θ¯ 1 + θ¯ 1 − θ¯ 1 + θ¯ 1 − θ¯ 2

1 1 − θ¯ 1 + θ¯ 1 − θ¯ 2 1 − θ¯ 1 + θ¯ 1 − θ¯ 2 1 − θ¯ 2 1 − θ¯ 2

1 1 − θ¯ 1 + θ¯ 1 − θ¯ 2 1 − θ¯ 1 + θ¯ 1 − θ¯ 2 1 − θ¯ 2 1 − θ¯ 2

1 1 1 1 1 − θ¯ 1 + θ¯ 1 − θ¯ 1 + θ¯ 1 − θ¯ 2

15 16 17 18 19 20 21 22 23 24 25

32

35

i ωt

17

ϕθ (x, θ, t ) = θ (x, θ)eiωt ,

(37e)

where ω represents the natural frequency of the rotating pretwisted composite tapered cylindrical panel. The Chebyshev polynomials are selected to establish the mode shape functions. The mode shape functions are expressed as the double series of the Chebyshev polynomials multiplied by the corresponding boundary functions

U (¯x, θ¯ ) = F u (¯x, θ¯ )

N M  

U mn P m (¯x) P n (θ¯ ),

N M  

V (¯x, θ¯ ) = F v (¯x, θ¯ )

W (¯x, θ¯ ) = F (¯x, θ¯ )

44 45 46 47 48 49 50 51 52 53 54

¯ W mn P m (¯x) P n (θ),

(38c)

¯ xmn P m (¯x) P n (θ),

59 60 61 62 63 64 65 66

23 25

28

31 32 34

(38d)

35 36

θ mn P m (¯x) P n (θ¯ ),

37

(38e)

38 39

where U mn , V mn , W mn , xmn and θ mn are the coefficients to be determined, P m (¯x) and P n (θ¯ ) are the mth-order and nth-order Chebyshev polynomials of the first kind



 P m (¯x) = cos (m − 1) arccos(¯x) ,   P n (θ¯ ) = cos (n − 1) arccos(θ¯ ) .

(39a) (39b)

The Chebyshev polynomials are a set of orthogonal series in the interval [−1, 1]. Therefore, we introduce the dimensionless coordinates ¯ ( = u , v , w , ϕx , ϕθ ) are the boundary functions corresponding to the x¯ = 2x − 1 and θ¯ = 2θ . In addition, in equation (38), F  (¯x, θ) L boundary conditions

( = u , v , w , ϕx , ϕθ ),

(40)

where F  (¯x) and F  (θ¯ ) are the boundary functions along the x¯ and θ¯ directions, respectively, which are shown in Table 1 [43]. Substituting equation (37) into equations (28) and (30), and neglecting the nonlinear terms, the energy function  of the rotating pretwisted composite tapered cantilever cylindrical panel is derived by

40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55

 = U max − K max .

(41)

56 57

57 58

22

33

55 56

21

30

N M  

F  (¯x, θ¯ ) = F  (¯x) F  (θ¯ ),

20

29

m =1 n =1

39

43

(38b)

m =1 n =1

θ (¯x, θ¯ ) = F ϕθ (¯x, θ¯ )

19

27

V mn P m (¯x) P n (θ¯ ),

N M  

x (¯x, θ¯ ) = F ϕx (¯x, θ¯ )

18

26

m =1 n =1 w

15

24

(38a)

m =1 n =1

N M  

42

11

(37d)

37

41

10

,

ϕx (x, θ, t ) = x (x, θ)e

m =1 n =1

40

9

16

i ωt

36 38

8

(37c)

,

w 0 (x, θ, t ) = W (x, θ)ei ωt ,

33 34

7

14

30 31

6

(37b)

v 0 (x, θ, t ) = V (x, θ)e

27 29

5

13

26 28

4

12

Note: F, free edge; S, simply-supported edge; C, clamped edge.

13 14

3

Based on the Ritz method, we obtain

∂ = 0, ∂ U mn ∂ = 0, ∂ xmn

∂ = 0, ∂ V mn ∂ = 0, ∂ θ mn

58

∂ ∂ W mn

59

= 0,

(m = 1, ..., M ; n = 1, ..., N ).

60 61

(42)

63 64

The eigenvalue equation is rewritten in the following matrix form





K − ω 2 M d = 0,

62

65

(43)

66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.9 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

1 2 3 4

Source

2 3

Symbol

Value

Unit

Description

L R h0 h G P L0 R0 ch

0.1 0.5 0.002 0.0002 0.3 0.5 60 30 30

m m m m m

Length of the rotating blade Radius of the rotating blade Thickness of the rotating blade root Thickness of the graphene coating layer of the rotating blade root Radius of the rotating hub Taper ratio of the rotating blade Subtended angle of the rotating blade Presetting angle of the rotating blade Pretwist angle of the rotating blade

2.5 1.5 1.5 1 1.01 1060 0.186

μm μm nm % TPa kg/m3

3 × 109 1200 0.34

Pa kg/m3

4

Geometric parameter

6 7 8 9 10 11

 φ

12

K

◦ ◦ ◦

lGPL w GPL hGPL W GPL E GPL

15 16 17

ρGPL νGPL

18 19

Length of the graphene platelets Width of the graphene platelets Thickness of the graphene platelets Weight fraction of the graphene platelets Young’s modulus of the graphene platelets Density of the graphene platelets Poisson’s ratio of the graphene platelets

EM

22

ρM νM

27 28 29

b/h b/h = 10

b/h = 100

37

18 19 21 22 23

ν = 0.3 and

1.0

1.5

2.0

Shen and Wang [44] Kobayashi and Leissa [45] Chern and Chao [46] Present

1.31532 1.3360 1.31742 1.31489

0.55240 0.5563 0.55049 0.55197

0.40227 0.4044 0.39987 0.40181

0.34887 0.3505 0.34612 0.34842

Shen and Wang [44] Kobayashi and Leissa [45] Chern and Chao [46] Present

0.16149 0.16150 0.16066 0.16146

0.07429 0.07429 0.07368 0.07425

0.05053 0.05053 0.04912 0.05049

0.04039 0.04039 0.03925 0.04035

30

45

50 51 52 53

58 59 60 61 62 63 64 65 66

36 37

41

d = {U 11 , ..., U M N , V 11 , ..., V M N , W 11 , ..., W M N , x11 , ..., xM N , θ 11 , ..., θ M N }.

(44)

The natural frequencies and the corresponding mode shapes of the rotating composite tapered blade are obtained by solving equation (43).

42 43 44 45 46

4. Numerical results and discussions

47 48

The purpose of this study is to discuss the free vibration characteristics of the rotating pretwisted composite tapered cantilever cylindrical panel with the graphene coating layers. At the beginning, a series of comparison studies are performed. Then, the parameter studies are conducted to illustrate the effects of the involved parameters on the free vibration of the rotating pretwisted composite tapered cantilever cylindrical panel with the graphene coating layers. Unless otherwise stated, the geometry parameters and material parameters are given in Table 2 and the rotating speed is selected to be 3000 rpm.

49 50 51 52 53 54

4.1. Comparison studies

55 56

56 57

35

40

54 55

33

39

where K is the stiffness matrix and M is the mass matrix, respectively, d is the displacement coefficients vector

48 49

32

38

46 47

31

34

41

44

27 29

0.5

39

43

26 28

38

42

12

17

L /b

36

40

11

25

33 35

10

16

Young’s modulus of the matrix Density of the matrix Poisson’s ratio of the matrix

Table 3  ¯ = ωb (1 − ν 2 )ρ / E of the simply supported cylindrical panel is performed when R /b = 10, E = 1, The comparison of the dimensionless natural frequencies ω ρ = 1.

32 34

9

24

30 31

8

15

24 26

7

20

21

25

6

14

Epoxy matrix

23

5

13

Graphene platelets

14

20

1

Table 2 The geometric parameters and material parameters of the rotating pretwisted composite tapered blade are given.

5

13

9

¯ = The accuracy of our present model and the solution method is verified by comparing the dimensionless natural frequencies ω  ωb (1 − ν 2 )ρ / E with the results of Shen and Wang [44], Kobayashi and Leissa [45], Chern and Chao [46] for a simply supported cylindrical panel with R /b = 10, E = 1, ν = 0.3 and ρ = 1, as shown in Table 3. The accuracy of our present model is illustrated according to

57

the good agreement of our present results with the results in the references. √ ¯ = ω( R 2 /h0 ) ρ / E and the corresponding mode shapes The next study gives the comparison of the dimensionless natural frequencies ω for the first five modes of a pretwisted tapered cylindrical panel with the results of ANSYS. The geometric parameters and the material parameters are selected as L = 0.1 m, R = 0.05 m,  = 60◦ , φ = 30◦ , h0 = 0.002 m, ch = 0.5, E = 1.01 × 105 MPa, ρ = 4480 kg/m3 and ν = 0.32. Table 4 gives the comparison results of the dimensionless natural frequencies for the first five modes. Figs. 3 and 4 depict the corresponding mode shapes for the first five modes of a pretwisted tapered cantilever cylindrical panel. The good agreements of the first five dimensionless natural frequencies and mode shapes are found.

60

58 59 61 62 63 64 65 66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.10 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

10

1

Table 4 √ ¯ = ω( R 2 /h0 ) ρ / E for the first five modes of the The comparison of the dimensionless natural frequencies ω pretwisted tapered cylindrical panel is performed when L = 0.1 m, R = 0.05 m,  = 60◦ , φ = 30◦ , h0 = 0.002 m and ch = 0.5.

2 3 4

Mode

5 6 7

1 2 3 4 5

8 9 10

K =0

K = π /12

1 2 3 4 5

ANSYS

Present

Error

ANSYS

Present

Error

6

0.9937 1.0681 3.0070 4.0534 4.2642

1.0483 1.0698 3.0617 4.3392 4.4824

5.49% 0.16% 1.82% 7.05% 5.12%

0.7707 1.2984 2.8463 4.1735 4.3449

0.8253 1.3337 3.0271 4.3806 4.5067

7.08% 2.72% 6.35% 4.96% 3.72%

7 8 9 10

11

11

12

12

13

13

14

14

15

15

16

16

17

17

18

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29

30

30

31

31

32

32

33

33

34

34

35

35

36

36

37

37

38

38

39

39

40

40

41

41

42

42

43

43

44

44

45

45

46

46

47

47

48

48

49

49

50

50

51

51

52

52

53

53 54

54 55 56

Fig. 3. The mode shapes of the pretwisted tapered cylindrical panel with K = 0◦ are depicted.

Fig. 4. The mode shapes of the pretwisted tapered cylindrical panel with K = 15◦ are depicted.

4.2. Parameter studies

61 62 63 64 65 66

58 59

59 60

56 57

57 58

55

After conducting a series of comparison studies to verify the accuracy of our present model, the free vibration characteristics of the rotating pretwisted composite tapered cantilever cylindrical panel with the graphene coating layers are discussed in detail. It is worth √ ¯ = ω( R 2 /h0 ) ρ M / E M . noticing that the dimensionless natural frequencies are defined by ω We focus on the GPL parameters of the graphene coating layers at the beginning. Table 5 provides the dimensionless natural frequencies of the rotating pretwisted composite tapered blade with different GPL dimensions. The length of the graphene platelet is selected as 2.5 μm. With the increase of the lGPL /hGPL , the first five dimensionless natural frequencies increase. Besides, the increasing rate of the dimensionless natural frequencies becomes smaller when the lGPL /hGPL reaches to 1(×1000). Moreover, with the decrease of the lGPL / w GPL ,

60 61 62 63 64 65 66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.11 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

1 2 3 4

Table 5 √ ¯ = ω( R 2 /h0 ) ρ M / E M of the rotating pretwisted composite tapered blade are Dimensionless natural frequencies ω obtained with different GPL dimensions. lGPL / w GPL

lGPL /hGPL (×1000)

2

7 8 9 10 11 12 13 14

1

Mode

1 2 3 4

1

2

3

4

5

5

0.01 0.1 1.0 3.0 5.0

0.94914 1.02870 1.17041 1.20296 1.21073

1.76656 1.95626 2.29730 2.37559 2.39427

3.41952 3.85317 4.60231 4.77092 4.81102

4.63563 5.15358 6.04119 6.24371 6.29200

4.86295 5.65281 7.01808 7.32147 7.39340

6

0.01 0.1 1.0 3.0 5.0

0.95503 1.05581 1.18763 1.21037 1.21543

1.78051 2.02138 2.33872 2.39339 2.40556

3.45206 3.99860 4.69164 4.80912 4.83522

4.67627 5.32453 6.14836 6.28972 6.32119

4.92049 5.91990 7.17899 7.39001 7.43680

11

5 6

11

7 8 9 10 12 13 14

15

15

16

16

17

17

18

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29

30

30

31

31 32

32 33 34

Fig. 5. The effect of the GPL weight fraction on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade is given.

Fig. 6. The effect of the length-to-radius ratio on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade is given.

36

36 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

34 35

35 37

33

the w GPL increases, and then the dimensionless natural frequencies increase. This implies that the large contact area between the GPLs and the matrix enhances the stiffness of the structure more efficiently. Fig. 5 portrays the effect of the GPL weight fraction on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade. By dispersing a small amount of GPLs into the matrix, the stiffness of the structure is strengthened. Therefore, the dimensionless natural frequencies increase with the increase of the GPL weight fraction. In order to emphasize the enhancement effect of the graphene clearly, we define the relative frequency increase ω and give the following expression ω = (ωC − ω M )/ω M , where ω M and ωC are the natural frequencies of the rotating pretwisted tapered blade without and with the GPLs, respectively. When the GPL weight fraction W GPL is 1.2%, the relative frequency increases ω are 32.3%, 41.9%, 48.1%, 42.5% and 61.6% for the first five modes, respectively. We investigate the geometric parameters of the rotating pretwisted composite tapered blade. The effect of the length-to-radius ratio on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade is investigated, as shown in Fig. 6. It is observed from Fig. 6 that the length-to-radius ratio plays an important role to the dimensionless natural frequencies, especially for the higher modes. For the first five modes, as the length-to-radius ratio increases, the dimensionless natural frequencies decrease rapidly at first, and then decrease slowly for the rotating pretwisted composite tapered blade. Fig. 7 reveals the effect of the taper ratio on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade. With the increase of the taper ratio, the dimensionless natural frequencies of the first-order mode and the second-order mode vibrations increase, which means that the decrease of the blade thickness does not always lead to the decrease of the stiffness for the rotating pretwisted composite tapered blade. The dimensionless natural frequency of the third-order mode decreases with the increase of the taper ratio. The effect of the taper ratio is complex on the higher mode. It is noticed that the frequency veering phenomenon which was analyzed by Li and Zhang [47,48] for the rotating beam and plate is observed in Fig. 7. The frequency veering phenomenon associated with the fourth-order mode and the fifth-order mode occurs when the taper ratio is ch = 0.7. Fig. 8 demonstrates the corresponding fourth-order and fifth-order mode shapes of the rotating pretwisted composite tapered blade with the taper ratios ch = 0.6 and ch = 0.8, respectively. The shift phenomenon of the mode shape is observed. The fourth-order and fifth-order mode shapes exchange when the taper ratio changes from ch = 0.6 to ch = 0.8. Fig. 9 indicates the effects of the taper ratio and pretwist angle on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade when the taper ratio changes from ch = 0 to ch = 1, and the pretwist angle is selected as K = 0◦ , K = 15◦ and K = 30◦ , respectively. The frequency veering and mode shape shift phenomena are detected in the rotating pretwisted composite tapered blade. Three frequency veering regions depicted in Fig. 9(a) are associated with the first-order and second-order modes when ch = 0.7, with the fourth-order and fifth-order modes when ch = 0.6, and with the third-order and fourth-order modes when ch = 0.9. For the rotating pretwisted composite tapered blade with the pretwist angle K = 15◦ and K = 30◦ , it is observed from Figs. 9(b) and 9(c) that

37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

JID:AESCTE

12

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.12 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

1

1

2

2

3

3

4

4

5

5

6

6

7

7

8

8

9

9

10

10

11

11

12

12

13

13

14

14 15

15 16

Fig. 7. The effect of the taper ratio on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade is given.

16

17

17

18

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29

30

30

31

31

32

32

33

33

34

34

35

35 36

36 37

Fig. 8. The mode shapes of the rotating pretwisted composite tapered blade are depicted.

37

38

38

39

39

40

40

41

41

42

42

43

43

44

44

45

45

46

46

47

47

48

48

49

49

50

50

51

51

52

52

53

53

54

54

55

55

56

56

57

57

58

58

59

59

60

60

61

61

62

62

63

63

64

64 65

65 66

Fig. 9. The effects of the taper ratio and pretwist angle on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade are given.

66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.13 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

13

1

1

2

2

3

3

4

4

5

5

6

6

7

7

8

8

9

9

10

10

11

11

12

12

13

13

14

14

15

15

16

16 17

17 18

Fig. 10. The mode shapes of the rotating pretwisted composite tapered cantilever cylindrical panel are depicted with K = 0◦ .

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29

30

30

31

31

32

32

33

33

34

34

35

35

36

36

37

37

38

38

39

39

40

40

41

41 42

42 43

Fig. 11. The mode shapes of the rotating pretwisted composite tapered cantilever cylindrical panel are depicted with K = 0◦ .

45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

43 44

44

the frequency veering phenomena are found between the fourth-order and fifth-order modes under different taper ratios in the rotating pretwisted composite tapered blade. Figs. 10 and 11 portray the mode shape shift phenomena corresponding to Fig. 9(a) in the rotating pretwisted composite tapered blade. Figs. 10(a) and 10(b) give the first-order and second-order mode shapes versus ch = 0.6. Figs. 10(c) and 10(d) demonstrate the first-order and second-order mode shapes versus ch = 0.8. It is obviously known that when 0 < ch < 0.7, the first-order and second-order mode shapes exhibit the bending and torsional vibrations of the rotating pretwisted composite tapered blade, respectively. When 0.7 < ch < 1.0, the first-order and second-order mode shapes become the torsional and bending vibrations of the rotating pretwisted composite tapered blade, respectively. Figs. 11(a)-11(c), 11(d)-11(f) and 11(g)-11(i) depict the third-order, fourth-order and fifth-order mode shapes of the rotating pretwisted composite tapered blade when ch = 0.5, ch = 0.7 and ch = 1.0, respectively. It is observed that the fourth-order mode shape exchanges with the fifth-order mode shape at first, and then exchanges with the third-order mode shape. Fig. 12 demonstrates the first five dimensionless natural frequencies of the rotating pretwisted composite tapered blade under different pretwist angles. At the beginning, when the pretwist angle is K = 0◦ , the first two dimensionless natural frequencies are close. With the increase of the pretwist angles, the dimensionless natural frequencies of the first-order and second-order modes decrease and increase obviously, respectively. It is observed that the pretwist angle has a little effect on the dimensionless natural frequencies of the higher modes. The mode shapes of the rotating pretwisted composite tapered blade with different pretwist angles corresponding to Fig. 12 are given, as shown in Fig. 13. It is observed that the first-order mode shape of the rotating pretwisted composite tapered blade under different pretwist angles is the first-order bending mode. In other words, the increase of the pretwist angles can decrease the first-order bending stiffness. The second-order mode shape of the rotating pretwisted composite tapered blade under different pretwist angles is the first-order torsional mode. The increase of the pretwist angles leads to the increase of the first-order torsional stiffness. Fig. 14 illustrates the effect of the presetting angle on the first five dimensionless natural frequencies of the rotating pretwisted composite tapered blade. Compared to the pretwist angle, the presetting angle does not have much effect on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade.

45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

JID:AESCTE

14

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.14 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

1

1

2

2

3

3

4

4

5

5

6

6

7

7

8

8

9

9

10

10

11

11

12

12

13

13

14

14 15

15 16

Fig. 12. The effect of the pretwist angle on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade is obtained.

16

17

17

18

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29

30

30

31

31

32

32

33

33

34

34

35

35

36

36

37

37

38

38

39

39

40

40

41

41

42

42

43

43

44

44

45

45 46

46 47

Fig. 13. The mode shapes of the rotating pretwisted composite tapered blade with different pretwist angles are obtained.

47

48

48

49

49

50

50

51

51

52

52

53

53

54

54

55

55

56

56

57

57

58

58

59

59

60

60

61

61

62

62

63

63

64

64 65

65 66

Fig. 14. The effect of the presetting angle on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade is obtained.

66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.15 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

15

1

1

2

2

3

3

4

4

5

5

6

6

7

7

8

8

9

9

10

10

11

11 12

12

Fig. 15. The Campbell diagrams of the rotating pretwisted composite tapered blade are given.

13

13

14

14

15

15

16

16

17

17

18

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29 30

30

Fig. 16. The mode shapes of the rotating pretwisted composite tapered blade are depicted when ch = 0.3.

31 33 34 35 36 37 38 39 40 41 42 43 44 45 46

Fig. 15 depicts the Campbell diagrams of the rotating pretwisted composite tapered blade when the rotating speed changes from 0 to 20000 rpm, where Figs. 15(a) and 15(b) represent the Campbell diagrams when the taper ratio is selected as ch = 0.3 and ch = 0.5, respectively. With the increase of the rotating speed, the first five dimensionless natural frequencies increase for the rotating pretwisted composite tapered blade. This is due to the fact that the increase of the rotating speed leads to the increase of the stiffness. The Campbell diagram also exhibits the mode shape shift phenomenon. When the taper ratio is ch = 0.3, it is observed from Fig. 15(a) that the fourthorder mode shape exchanges with the fifth-order mode shape at 0 = 16000 rpm. In Fig. 15(b), when the taper ratio is ch = 0.5, the fourth-order mode shape firstly exchanges with the fifth-order mode shape at 0 = 10000 rpm, and then exchanges with the third-order mode shape at 0 = 14000 rpm. Fig. 16 depicts the mode shape shift phenomenon corresponding to Fig. 15(a). Figs. 16(a) and 16(b) are the fourth-order and fifthorder mode shapes when 0 = 14000 rpm. Figs. 16(c) and 16(d) demonstrate the fourth-order and fifth-order mode shapes when 0 = 18000 rpm. It is obviously known that the fourth-order mode shape exchanges with the fifth-order mode shape. 5. Conclusions

49 50 51 52 53 54 55 56 57 58 59 60 61 62

This paper develops a new dynamic model of the rotating pretwisted composite tapered cantilever cylindrical panel with graphene coating layers to investigate the free vibration characteristics. The Chebyshev-Ritz method is used to obtain the natural frequencies and mode shapes of the rotating pretwisted composite tapered blade. The accuracy of the proposed model is validated. The effects of the GPL geometry, GPL weight fraction, taper ratio, length-to-radius ratio, pretwist angle, presetting angle and rotating speed on the free vibration characteristics of the rotating pretwisted composite tapered blade are discussed in detail through the comprehensive parametric studies. The frequency veering and mode shape shift phenomena are found in the rotating pretwisted composite tapered blade. The numerical results demonstrate that by dispersing a small amount of GPLs into the matrix, the stiffness of the structure is strengthened. The large contact area between the GPLs and the matrix enhances the stiffness of the structure more efficiently. The taper ratio and pretwist angle both have a complex effect on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade. The presetting angle does not have much effect on the dimensionless natural frequencies of the rotating pretwisted composite tapered blade. The Campbell diagram illustrates that the high rotating speed results in the large stiffness of the rotating pretwisted composite tapered blade due to the effect of the centrifugal force. Based on the investigation of this paper about the free vibrations of the rotating pretwisted composite tapered cantilever cylindrical panel with the graphene coating layers, the further investigations on the nonlinear dynamics and internal resonances of the rotating pretwisted composite tapered cantilever cylindrical panel with the graphene coating layers will be conducted. Declaration of competing interest

35 36 37 38 39 40 41 42 43 44 45 46 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 64 65

65 66

34

63

63 64

33

47

47 48

31 32

32

The authors declare that there is no conflict of interest regarding the publication of this paper.

66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.16 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

16

1

Acknowledgements

1 2

2 3 4 5

The authors gratefully acknowledge the support of National Natural Science Foundation of China (NNSFC) through grant Nos. 11832002, 11290152 and 11427801, the Funding Project for Academic Human Resources Development in Institutions of Higher Learning under the Jurisdiction of Beijing Municipality (PHRIHLB). The authors also gratefully acknowledge the support of China Scholarship Council (CSC). Appendix A

10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

8

The Green strain tensors in equation (17) are represented by

51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

9



 ∂u 2 f 11 = 2 1 + R 2 k2 sin2 θ + e 21 k2 + z(−a11 e 1 k + a12 Rk sin θ) ∂x  ∂v 2 + 2 e 2k − z(a11 cos θ + a12 sin θ) + 2za12k eu ∂x  2  + 2 Rk sin θ cos θ + e 1k2 sin θ + z(−a11k sin θ + a12k cos θ) v    + 2 a12 Rk sin θ − a11 e 1 k + z a211 + a212 w + G 211 + G 212 + G 213 ,    ∂u 1 1 1 2 f 22 = 2 ke 2 + z a21 + a22 e 1 k − a23 k sin θ R R R ∂θ    1 1 1 ∂v +z + 2z(a22k sin θ − a23k cos θ)u a22 cos θ + a23 sin θ +2 R R R ∂θ      1 1 + 2z − a22 sin θ + a23 cos θ v + 2 a22 cos θ + a23 sin θ + z a221 + a222 + a223 w R

10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

R

+ + + ∂w 2 f 33 = 2 + G 231 + G 232 + G 233 , ∂z   ∂u  ∂v 2 f 12 = ke 2 + z(a21 + a22 e 1 k − a23 Rk sin θ) + 1 + z(a22 cos θ + a23 sin θ) ∂x ∂x    1 2 1 2 2 1 ∂u + B + e 2 k + z a12 k sin θ − a11 e 1 k R R R ∂θ    1 1 ∂v a11 cos θ + a12 sin θ + ke 2 − z R R ∂θ   2 + z −a11k sin θ + a12k cos θ − a22 Rk sin θ − a23 e 1k2 u    1 1 1 + − ek sin θ + z a22k sin θ − a23k cos θ + a11 sin θ − a12 cos θ v R R R   2k + − − 2z(a11a22 + a12a23 ) w + G 11 G 21 + G 12 G 22 + G 13 G 23 ,

25

 ∂u  ∂v ∂w  2 2 f 13 = + B + k2 e 22 + z(a12 Rk sin θ − a11 e 1k) + ke 2 − z(a11 cos θ + a12 sin θ) ∂x ∂z ∂z

42

G 221

G 222

G 223 ,

26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41

B

+

1

B ∂w

e 2 k2 u +

1 B

43 44 45

kv + G 11 G 31 + G 12 G 32 + G 13 G 33 ,

46

 ∂u  ∂v 2 f 23 = + 1 + z(a22 cos θ + a23 sin θ) + ke 2 + z(a21 + a22 e 1k − a23 Rk sin θ) R ∂θ ∂z ∂z

49 50

5 7

8 9

4 6

6 7

3

+

k B



u−

1

47 48 49

v + G 21 G 31 + G 22 G 32 + G 23 G 33 ,

BR ∂u ∂ w 1 ∂u ∂u ∂ w 1 ∂w 1 − G 21 = G 31 = − G 11 = − ke sin θ, ke sin θ + wa21 , ke sin θ, ∂x ∂x B R ∂θ ∂θ B R ∂z ∂z B ∂u ∂v ∂w 1 G 12 = k(e − R cos θ) + cos θ + sin θ − uk2 R sin θ + vk sin θ − wa11 , ∂x ∂x ∂x B ∂u ∂v ∂w 1 G 13 = kR sin θ − sin θ + cos θ + uk2 (e − R cos θ) + vk cos θ + wa12 , ∂x ∂x ∂x B ∂u ∂v ∂w 1 v G 22 = k(e − R cos θ) + uk sin θ + cos θ − sin θ + sin θ + wa22 , R ∂θ R ∂θ R ∂θ B R ∂u ∂v ∂w 1 v G 23 = k sin θ + uk cos θ − sin θ − cos θ + cos θ − wa23 , ∂θ R ∂θ R ∂θ B R ∂u ∂v ∂w 1 G 32 = k(e − R cos θ) + cos θ + sin θ, ∂z ∂z ∂z B ∂u ∂v ∂w 1 G 33 = kR sin θ − sin θ + cos θ, ∂z ∂z ∂z B e 2 = e cos θ − R . e 1 = e − R cos θ,

(A.1)

50 51 52 53 54 55 56 57 58 59 60 61 62 63

(A.2) (A.3)

64 65 66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.17 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

1

17

Appendix B

1 2

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

  ∂u 1 1 2 2 1 ∂u εξ ξ = B − Bke 2 + k e sin θ cos θ v + 2 k2 e 2 e cos θ w M ∂x R ∂θ BR B R 2    2  1 1 ∂w 1 ∂w ∂w 1 1 2 2 ∂w + 2 , − ke2 k e + 2 R ∂ x ∂θ 2 R2 ∂θ M 2 ∂x 2    1 1 ∂u B ∂v 1 1 1 ∂w εηη = , Bke 2 + + w + M R ∂θ R ∂θ R 2 M 2 R 2 ∂θ    1 ∂u ∂ v 1 B2 ∂u ke 2 ∂ v γξ η = ke 2 + + − e 2k2 cos2 θ + 2k2 e cos θ − k2 R + − M ∂x ∂x R R ∂θ R ∂θ 2     ke sin θ 2ke cos θ 1 B ∂w ∂w kBe 2 ∂ w , − − 2 v− w + 2 R BR R ∂ x ∂θ ∂θ M R   1 ∂w ke 2 ∂ w ke cos θ ∂u , γξ ζ = − + v + B2 M ∂x R ∂θ BR ∂z   1 B ∂w v ∂u ∂v , γη ζ = +B + ku − + Bke 2 M R ∂θ R ∂z ∂z   k2 ee 1 + 1 . M = B 1+z 3 1



B R

26 27

3

The strain components in equation (21) are represented by

30

33

K=

1 2

40 41

46 47 48 49 50

53

56 57 58

61

66

12 13 14 15

(B.3)

16 17 18

(B.4)

19 20 21

(B.5)

22 23 24

(B.6)

25

33 35



n3 I 2 x2

ϕ + 2I 1 u 0 ϕx +

I 0 u 20



1



2

+ n4 I 2 ϕθ + 2I 1 v 0 ϕθ + 2

I 0 v 20



1 2

+ +

41 43

1

45

2 1

2

2

46

44



2

1



I 2 x2

2

 ∂ ϕx ∂ u0 ∂ ϕx ∂ u0 ϕθ + I 1 ϕθ + I 1 n15 I 2 v0 + I0 v0 2 ∂t ∂t ∂t ∂t   1 ∂ ϕx ∂ u0 ∂ ϕx ∂ u0 + I1 + I1 n16 I 2 w0 + I0 w0 2 ∂t ∂t ∂t ∂t   1 ∂ ϕθ ∂ ϕθ ∂ v0 ∂ v0 ϕx + I 1 ϕx + I 0 n17 I 2 u0 + I 1 u0 2 ∂t ∂t ∂t ∂t   1 ∂ ϕθ ∂ v0 ∂ ϕθ ∂ v0 + I1 + I1 n18 I 2 w0 + I0 w0 2 ∂t ∂t ∂t ∂t     1 ∂ w0 ∂ w0 1 ∂ w0 ∂ w0 ϕx + I 0 ϕθ + I 0 n19 I 1 u 0 + n20 I 1 v0 2 ∂t ∂t 2 ∂t ∂t      1 ∂ ϕx ∂ u0 1 ∂ ϕθ ∂ v0 1 ∂ w0 + n22 I 1 + n23 I 0 + I0 + I0  n21 I 1 2 ∂t ∂t 2 ∂t ∂t 2 ∂t



+

40

2 1

 1  + n12 I 0 ( R 0 + x + u 0 ) + ϕ + 2I 1 ϕx ( R 0 + x + u 0 ) + n13 I 0 w 0 ( R 0 + x + u 0 ) 2 2   1 2 + I 2 ϕx + I 1 ( R 0 + x + u 0 ) + I 1 w 0 ϕx + n14 I 0 

+

38

2 1

+ n9 ( I 2 ϕθ + I 1 w 0 ϕθ + I 1 v 0 + I 0 v 0 w 0 ) + n10 ( I 1 ϕθ + I 0 v 0 ) + n11 ( I 0 w 0 + I 1 )

+

37

42

1

+ n7 ( I 2 ϕx + I 1 u 0 + I 1 w 0 ϕx + I 0 u 0 w 0 ) + n8 ( I 1 ϕx + I 0 u 0 )

+

36

39



2 1

64 65

11

34



1

62 63

(B.2)

32

ρe(kf )f r˙ 2 MdzRdθ dx

+ n5 I 0 w 20 + 2I 1 w 0 + I 2 + n6 ( I 2 ϕx ϕθ + I 1 u 0 ϕθ + I 1 v 0 ϕx + I 0 u 0 v 0 )

59 60

10

31

0 − 2

54 55

1 2

51 52

9

30

zk+1 

0 −  k=1 zk 2

=

44 45

8

29



42 43

2  3

L  2 

36

39



L

35

38

7

(B.1)

28

The detailed expression and coefficients in equation (28) are given as follows

34

37

6

27

Appendix C

31 32

5

26

28 29

4

47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

JID:AESCTE

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.18 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

18

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

 

 2  ∂ u 0 ∂ ϕx ∂ u0 + n1 I 2 + 2I 1 + I0 ∂t ∂t ∂t   2  2   2 ∂ ϕθ ∂ v0 ∂ w0 1 ∂ v 0 ∂ ϕθ 1 + B R I0 + I0 + 2I 1 B R I2 2 ∂t ∂t ∂t ∂t 2 ∂t   ∂ ϕx ∂ ϕθ ∂ u 0 ∂ ϕθ ∂ v 0 ∂ ϕx ∂ u0 ∂ v 0 + n2 I 2 + I1 + I1 + I0 dθ dx, ∂t ∂t ∂t ∂t ∂t ∂t ∂t ∂t   1 B R + k2 B R R 2 + e 2 − 2e R cos θ , n2 = kB R (e cos θ − R ), n1 = 2   2 n3 = k2 B R R cos( − θ) − e cos  , n4 = B R cos( − θ)2 ,   R R n5 = sin( − θ)2 + k2 e 2 sin2 θ, n6 = −2kB R R cos( − θ) − e cos  cos( − θ), B B   n7 = 2kR R cos( − θ) − e cos  sin( − θ),    n8 = 2kB R R cos( − θ) − e cos  R sin( − θ) − e sin  , n9 = −2R sin( − θ) cos( − θ),     n10 = −2B R cos( − θ) R sin( − θ) − e sin  , n11 = 2R sin( − θ) R sin( − θ) − e sin  ,  2 n12 = B R , n13 = −2ke R sin θ, n14 = B R R sin( − θ) − e sin  , n15 = −2B R cos( − θ),   2 2 2 n16 = 2R sin( − θ) − k e sin θ cos  + k e R sin θ cos( − θ) , n17 = 2B R cos( − θ),   2 n18 = −2ke R sin θ cos( − θ), n19 = −2k e R sin θ R cos( − θ) − e cos  − 2R sin( − θ), ∂ ϕx ∂t

2

n20 = 2ke R sin θ cos( − θ),

n22 = 2B R ( R 0 + x) cos( − θ),



28

 = φ + kx.

43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29

33



U s1 =

z L  2  3 k+1  1

2

σξ(kξ ) εξ ξ

(k)

+ σηη εηη + σξ(kη) γξ η

+ σξ(kζ ) γξ ζ

 + ση(kζ) γηζ MdzRdθ dx

0 −  k=1 zk 2 

L  2 

39

42

8

32

The detailed expression in equation (31) is given as follows

38

41

7

31

37

40

6

(C.1)

30

Appendix D

34 36

5

(C.2)

33 35

4



31 32

3





n23 = −2ke R sin θ R sin( − θ) − e sin  − 2R ( R 0 + x) sin( − θ),

30

2

n21 = B R 2R sin( − θ) − 2e sin  − 2kR ( R 0 + x) cos( − θ) + 2ke ( R 0 + x) cos  ,

27 29

1

=

1 2

0 − 2

34 35 36 37

 B R D 11 m11

∂ ϕx ∂ ϕx + m12 ϕθ + ∂θ ∂x

2

2  ∂ ϕθ 1 ∂ ϕx + B R D 22 −m11 + 2 ∂θ R ∂θ

     ∂ ϕx ∂ ϕθ ∂ ϕx ∂ ϕx 1 k ϕx ϕθ 2 + B Rk0 D 44 + B R D 12 −m11 m11 − + + m12 ϕθ + ∂θ R ∂θ ∂θ ∂x 2 B BR  2 1 1 ∂ ϕx ∂ ϕθ ∂ ϕx ∂ ϕθ + B Rk0 D 55m235 ϕθ2 + B R D 66 m31 + + m32 + m33 + m34 ϕθ 2 2 ∂x B∂ x ∂θ ∂θ 2 2   2   ∂ u0 ∂ w0 ∂ w0 ∂ u0 1 ∂ w0 ∂ w0 + B R B 11 + m11 + m14 + m12 v 0 + m13 w 0 + 2 + m15 ∂x ∂θ ∂x ∂ x ∂θ ∂θ 2B 2     ∂ ϕx ∂ v0 ∂ u0 w0 1 ∂ ϕx ∂ w0 + B R B 12 −m11 + m11 + + + m12 ϕθ + ∂θ R ∂θ BR ∂θ ∂θ ∂x 2R 2 2    2   ∂ u0 1 ∂ w0 ∂ w0 ∂ u0 ∂ w0 ∂ w0 + m11 + m14 + B R B 12 + m12 v 0 + m13 w 0 + 2 + m15 ∂x ∂θ ∂x ∂ x ∂θ ∂θ 2B 2     ∂ v0 ∂ ϕθ ∂ u0 w0 1 ∂ ϕx ∂ w0 −m11 + + B R B 22 −m11 + + + ∂θ R ∂θ BR ∂θ ∂θ R ∂θ 2R 2    ∂ w 0 ku 0 k ϕx ϕθ v0 + B Rk0 B 44 − − + ke 2 ϕx + ϕθ − B R B 55k0m35 ϕθ + B BR R ∂θ B BR    ∂ v0 ∂ w0 ∂ u0 ∂ u0 ∂ v0 ∂ w0 + m33 + + m32 − m35 v 0 + B ϕx + B R B 66 m31 + m33 + m34 v 0 B∂ x ∂θ ∂x B∂ x ∂θ ∂θ 2    1 ∂ w0 ∂ w0 ∂ ϕx ∂ ϕθ ∂ ϕx ∂ ϕθ ∂ w0 + 2m35 w 0 + m31 + + m32 + m36 + m33 + m34 ϕθ B R ∂ x ∂θ ∂θ ∂x B∂ x ∂θ ∂θ

38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

JID:AESCTE AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.19 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

 2 2  2  ∂ u0 1 ∂ w0 ∂ w0 ∂ u0 ∂ w0 ∂ w0 + B R A 11 + m11 + m14 + m12 v 0 + m13 w 0 + 2 + m15 2 ∂x ∂θ ∂x ∂ x ∂θ ∂θ 2B 2    ∂ w0 ∂ u0 ∂ v0 ∂ u0 w0 1 ∂ u0 ∂ w0 ∂ w0 + B R A 12 −m11 + + m11 + + + m12 v 0 + m14 2 ∂θ R ∂θ BR ∂θ ∂x ∂θ ∂ x ∂θ 2R 2  2 2  2    ∂ w0 ∂ w0 ∂ w0 ∂ v0 1 B R A 22 ∂ u0 w0 1 + m13 w 0 + 2 + −m11 + m15 + + + ∂x ∂θ 2 ∂θ R ∂θ BR ∂θ 2B 2R 2  2  2 1 k 1 1 ∂ w0 ∂ w0 ∂ w0 + B Rk0 A 44 + B Rk0 A 55 + m33 + u0 − − m35 v 0 + B ϕx v 0 + ke 2 ϕx + ϕθ 2 R ∂θ B BR 2 B∂ x ∂θ  2 2    ∂ v0 1 ∂ u0 ∂ u0 ∂ v0 1 ∂ w0 ∂ w0 ∂ w0 + B R A 66 m31 + + m32 dθ dx. + m33 + m34 v 0 + 2m35 w 0 + + m36 2 ∂x B∂ x ∂θ ∂θ B R ∂ x ∂θ ∂θ 1

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

19



1 2 3 4 5 6 7 8 9 10 11

(D.1)

References

18 19 20 21 22 23 24 25 26 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 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

13 14 15 16

16 17

12

[1] O. Song, L. Librescu, Structural modeling and free vibration analysis of rotating composite thin-walled beams, J. Am. Helicopter Soc. 42 (1997) 358–369. [2] X.D. Yang, S.W. Wang, W. Zhang, T.Z. Yang, C.W. Lim, Model formulation and modal analysis of a rotating elastic uniform Timoshenko beam with setting angle, Eur. J. Mech. A, Solids 72 (2018) 209–222. [3] X.D. Yang, Z. Li, W. Zhang, T.Z. Yang, C.W. Lim, On the gyroscopic and centrifugal effects in the free vibration of rotating beams, J. Vib. Control 25 (2019) 219–227. [4] F.X. Wang, W. Zhang, Stability analysis of a nonlinear rotating blade with torsional vibrations, J. Sound Vib. 331 (2012) 5755–5773. [5] B. Zhang, Y.M. Li, Nonlinear vibration of rotating pre-deformed blade with thermal gradient, Nonlinear Dyn. 86 (2016) 459–478. [6] M.H. Yao, Y.P. Chen, W. Zhang, Nonlinear vibrations of blade with varying rotating speed, Nonlinear Dyn. 68 (2012) 487–504. [7] M.H. Yao, W. Zhang, Y.P. Chen, Analysis on nonlinear oscillations and resonant responses of a compressor blade, Acta Mech. 225 (2014) 3483–3510. [8] Y.G. Chen, J.Y. Zhai, Q.K. Han, Vibration and damping analysis of the bladed disk with damping hard coating on blades, Aerosp. Sci. Technol. 58 (2016) 248–257. [9] A. Kandil, W.A. El-Ganaini, Investigation of the time delay effect on the control of rotating blade vibrations, Eur. J. Mech. A, Solids 72 (2018) 16–40. [10] A. Kandil, H. El-Gohary, Investigating the performance of a time delayed proportional-derivative controller for rotating blade vibrations, Nonlinear Dyn. 91 (2018) 2631–2649. [11] M. Rafiee, F. Nitzsche, M. Labrosse, Dynamics, vibration and control of rotating composite beams and blades: a critical review, Thin-Walled Struct. 119 (2017) 795–819. [12] S.K. Sinha, K.E. Turner, Natural frequencies of a pre-twisted blade in a centrifugal force field, J. Sound Vib. 330 (2011) 2655–2681. [13] H.H. Yoo, S.K. Kim, D.J. Inman, Modal analysis of rotating composite cantilever plates, J. Sound Vib. 258 (2002) 233–246. [14] H.H. Yoo, C. Pierre, Modal characteristic of a rotating rectangular cantilever plate, J. Sound Vib. 259 (2003) 81–96. [15] D.X. Cao, B.Y. Liu, M.H. Yao, W. Zhang, Free vibration analysis of a pre-twisted sandwich blade with thermal barrier coatings layers, Sci. China, Technol. Sci. 60 (2017) 1747–1761. [16] M.H. Yao, L. Ma, W. Zhang, Nonlinear dynamics of the high-speed rotating plate, Int. J. Aerosp. Eng. (2018) 5610915. [17] S.P. Dong, L. Li, D.G. Zhang, Vibration analysis of rotating functionally graded tapered beams with hollow circular cross-section, Aerosp. Sci. Technol. 95 (2019) 105476. [18] J. Sun, I.L. Arteaga, L. Kari, General shell model for a rotating pretwisted blade, J. Sound Vib. 332 (2013) 5804–5820. [19] S.K. Sinha, R.P. Zylka, Vibration analysis of composite airfoil blade using orthotropic thin shell bending theory, Int. J. Mech. Sci. 121 (2017) 90–105. [20] L.T. Liu, Y.X. Hao, W. Zhang, J. Chen, Free vibration analysis of rotating pretwisted functionally graded sandwich blades, Int. J. Aerosp. Eng. (2018) 2727452. [21] X.J. Gu, Y.X. Hao, W. Zhang, L.T. Liu, J. Chen, Free vibration of rotating cantilever pre-twisted panel with initial exponential function type geometric imperfection, Appl. Math. Model. 68 (2019) 327–352. [22] J. Chen, Q.S. Li, Vibration characteristics of a rotating pre-twisted composite laminated blade, Compos. Struct. 208 (2019) 78–90. [23] X.X. Hu, T. Tsuiji, Free vibration analysis of rotating twisted cylindrical thin panels, J. Sound Vib. 222 (1999) 209–224. [24] X.X. Hu, T. Sakiyama, H. Matsuda, C. Morita, Fundamental vibration of rotating cantilever blades with pre-twist, J. Sound Vib. 271 (2004) 47–66. [25] M.H. Yao, Y. Niu, Y.X. Hao, Nonlinear dynamic responses of rotating pretwisted cylindrical shells, Nonlinear Dyn. 95 (2019) 151–174. [26] L. Yang, Q.X. Liu, Y.C. Zhou, W.G. Mao, C. Lu, Finite element simulation on thermal fatigue of a turbine blade with thermal barrier coatings, J. Mater. Sci. Technol. 30 (2014) 371–380. [27] C. Lee, X. Wei, J.W. Kysar, J. Hone, Measurement of the elastic properties and intrinsic strength of monolayer graphene, Science 321 (2008) 385–388. [28] S. Ghosh, W. Bao, D.L. Nika, S. Subrina, E.P. Pokatilov, C.N. Lau, A.A. Balabdin, Dimensional crossover of thermal transport in few-layer graphene, Nat. Mater. 9 (2010) 555–558. [29] M.A. Rafiee, J. Rafiee, Z. Wang, H. Song, Z.Z. Yu, N. Koratkar, Enhanced mechanical properties of nanocomposites at low graphene content, ACS Nano 3 (2009) 3884–3890. [30] B. Das, K.E. Prasad, U. Ramamurty, C.N.R. Rao, Nano-indentation studies on polymer matrix composites reinforced by few-layer graphene, Nanotechnology 20 (2009) 125705. [31] F. Lin, Y. Xiang, H.S. Shen, Temperature dependent mechanical properties of graphene reinforced polymer nanocomposites-a molecular dynamics simulation, Composites, Part B, Eng. 111 (2017) 261–269. [32] R. Moradi-Dastjerdi, K. Behdinan, Thermoelastic static and vibrational behaviors of nanocomposite thick cylinders reinforced with graphene, Steel Compos. Struct. 31 (2019) 529–539. [33] Z. Zhao, C. Feng, Y. Wang, J. Yang, Bending and vibration analysis of functionally graded trapezoidal nanocomposite plates reinforced with graphene nanoplatelets (GPLs), Compos. Struct. 180 (2017) 799–808. [34] Y. Wang, C. Feng, C. Santiuste, Z. Zhao, J. Yang, Buckling and postbuckling of dielectric composite beam reinforced with Graphene Platelets (GPLs), Aerosp. Sci. Technol. 91 (2019) 208–218. [35] J.J. Mao, W. Zhang, Buckling and post-buckling analyses of functionally graded graphene reinforced piezoelectric plate subjected to electric potential and axial forces, Compos. Struct. 216 (2019) 392–405. [36] J.J. Mao, W. Zhang, Linear and nonlinear free and forced vibrations of graphene reinforced piezoelectric composite plate under external voltage excitation, Compos. Struct. 203 (2018) 551–565. [37] M.T. Song, S. Kitipornchai, J. Yang, Free and forced vibrations of functionally graded polymer composite plates reinforced with graphene nanoplatelets, Compos. Struct. 159 (2017) 579–588. [38] X.Y. Guo, P. Jiang, W. Zhang, J. Yang, S. Kitipornchai, L. Sun, Nonlinear dynamic analysis of composite piezoelectric plates with graphene skin, Compos. Struct. 206 (2018) 839–852. [39] X.Y. Guo, P. Jiang, W.H. Yan, S.K. Lai, W. Zhang, Internal resonance responses of rectangular cross-ply composite plates with graphene skins, Int. J. Struct. Stab. Dyn. 19 (2019) 1950057. [40] R.M.R. Reddy, W. Karunasena, W. Lokuge, Free vibration of functionally graded-GPL reinforced composite plates with different boundary conditions, Aerosp. Sci. Technol. 78 (2018) 147–156. [41] H.S. Shen, Y. Xiang, Y. Fan, D. Hui, Nonlinear vibration of functionally graded graphene-reinforced composite laminated cylindrical panels resting on elastic foundations in thermal environments, Composites, Part B, Eng. 136 (2018) 177–186. [42] Y.H. Dong, L.W. He, L. Wang, Y.H. Li, J. Yang, Buckling of spinning functionally graded graphene reinforced porous nanocomposite cylindrical shells: analytical study, Aerosp. Sci. Technol. 82–83 (2018) 466–478.

17 18 19 20 21 22 23 24 25 26 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 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66

JID:AESCTE

20

1 2 3 4 5 6 7 8

AID:105644 /FLA

[m5G; v1.261; Prn:20/12/2019; 11:12] P.20 (1-20)

W. Zhang et al. / Aerospace Science and Technology ••• (••••) ••••••

[43] Y. Niu, W. Zhang, X.Y. Guo, Free vibration of rotating pretwisted functionally graded composite cylindrical panel reinforced with graphene platelets, Eur. J. Mech. A, Solids 77 (2019) 103798. [44] H.S. Shen, H. Wang, Nonlinear vibration of shear deformable FGM cylindrical panels resting on elastic foundations in thermal environments, Composites, Part B, Eng. 60 (2014) 167–177. [45] Y. Kobayashi, A.W. Leissa, Large amplitude free vibration of thick shallow shells supported by shear diaphragms, Int. J. Non-Linear Mech. 30 (1995) 57–66. [46] Y.C. Chern, C.C. Chao, Comparison of natural frequencies of laminates by 3-D theory, Part II: curved panels, J. Sound Vib. 230 (2000) 1009–1030. [47] L. Li, D.G. Zhang, Dynamic analysis of rotating axially FG tapered beams based on a new rigid-flexible coupled dynamic model using the B-spline method, Compos. Struct. 124 (2015) 357–367. [48] L. Li, D.G. Zhang, Free vibration analysis of rotating functionally graded rectangular plates, Compos. Struct. 136 (2016) 493–504.

1 2 3 4 5 6 7 8

9

9

10

10

11

11

12

12

13

13

14

14

15

15

16

16

17

17

18

18

19

19

20

20

21

21

22

22

23

23

24

24

25

25

26

26

27

27

28

28

29

29

30

30

31

31

32

32

33

33

34

34

35

35

36

36

37

37

38

38

39

39

40

40

41

41

42

42

43

43

44

44

45

45

46

46

47

47

48

48

49

49

50

50

51

51

52

52

53

53

54

54

55

55

56

56

57

57

58

58

59

59

60

60

61

61

62

62

63

63

64

64

65

65

66

66