Nonlinear resonant dynamics of geometrically imperfect higher-order shear deformable functionally graded carbon-nanotube reinforced composite beams

Composite Structures 174 (2017) 45–58

Contents lists available at ScienceDirect

Composite Structures journal homepage: www.elsevier.com/locate/compstruct

Nonlinear resonant dynamics of geometrically imperfect higher-order shear deformable functionally graded carbon-nanotube reinforced composite beams Raheb Gholami a,⇑, Reza Ansari b,⇑, Yousef Gholami b a b

Department of Mechanical Engineering, Lahijan Branch, Islamic Azad University, P.O. Box 1616, Lahijan, Iran Department of Mechanical Engineering, University of Guilan, P.O. Box 3756, Rasht, Iran

a r t i c l e

i n f o

Article history: Received 4 March 2017 Revised 10 April 2017 Accepted 17 April 2017 Available online 20 April 2017 Keywords: Geometrically imperfect shear deformable FG-CNTRC beams Nonlinear resonant dynamics Unified higher-order shear deformable beam model Initial imperfection

a b s t r a c t This study aims at numerically analyzing the nonlinear resonant dynamics of geometrically imperfect higher-order shear deformable functionally graded carbon nanotube-reinforced composite (FG-CNTRC) beams with various end conditions subjected to a harmonic transverse load. Introducing a generalized displacement field including various beam theories, employing Hamilton’s principle and taking into account geometrical nonlinearity and initial imperfection, three nonlinear coupled equations and associated boundary expressions are obtained for geometrically imperfect FG-CNTRC beams. These equations formulate the longitudinal, transverse and rotational motions of FG-CNTRC beams. An efficient multistep numerical solution approach based on the generalized differential quadrature (GDQ) method, a numerical Galerkin-based scheme and time periodic discretization is employed to convert the time-dependent nonlinear partial differential equations (PDEs) into a Duffing-type nonlinear set of ordinary differential equations (ODEs) which can be solved via the pseudo arc-length continuation technique. Nonlinear resonant dynamics characteristics are illustrated in the form of frequency-response and force-response curves; highlighting the influences of initial geometrical imperfection, geometrical parameters, excitation frequency and boundary conditions. Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction Remarkable mechanical, thermal and electrical properties of carbon nanotubes (CNTs) including their high stiffness and strength, light weight and high aspect ratio [1–4] enable them to extensively utilized in reinforcing polymer [5,6], ceramic [7] and metal [8] matrixes instead of traditional fibers. Accordingly, a new class of advanced materials known as functionally gradedcarbon nanotube reinforced (FG-CNTRC) materials have received a lot of attention of the researchers and engineers due to their potential applications in the novel nano- and micro-electromechanical systems (NEMs and MEMS) and electronic devices as well as numerous industrial fields such as aerospace, automotive, sports and electronics [9,10]. In the last few years, a wide range of examinations are carried out by scientists and researchers to experimentally and theoretically distinguish various aspects of materials properties as well as mechanical behaviors of CNTRC ⇑ Corresponding authors. E-mail addresses: [email protected] (R. Gholami), [email protected] (R. Ansari). http://dx.doi.org/10.1016/j.compstruct.2017.04.042 0263-8223/Ó 2017 Elsevier Ltd. All rights reserved.

materials [11–17]. Especially, the dynamics of mechanical structures made of CNTRC materials has been investigated expensively in the literature; however, a thorough literature review on this subject will not be undertaken in this work. Most of these examinations are limited to linear mechanical characteristics; few investigations have been carried out on the geometrically nonlinear dynamics of CNTRC structures. The linear studies generally concentrated on investigation of the linear free vibration characteristics, linear static and dynamic stability, determining the critical buckling loads and natural fundamental frequencies and linear bending [18–22]. For instance, using the finite element approach, the buckling and free vibration of FG-CNT reinforced polymer composite beams subjected to nonuniform temperature fields was examined by George and Murigendrappa [23]. Ghorbani Shenas et al. [24] analyze the free vibration behavior of the pre-twisted FG-CNTRC beams with different sets of boundary conditions using the Chebyshev–Ritz method. It was illustrated that in addition to the type of boundary condition, the pre-twist angle on the fundamental natural frequencies depend on the vibration mode number. Utilizing a unified formulation of finite prism method and a variational

46

R. Gholami et al. / Composite Structures 174 (2017) 45–58

approach, the three-dimensional free vibration of FG-CNTRC plates as well as laminated fiber-reinforced composite plates was analyzed by Wu and Li [25]. Recently, Ansari et al. [26] tried to use an efficient and novel numerical variational approach to examine the buckling of axially-loaded FG-CNTRC conical panels. A further few investigations on this subject incorporated the geometric nonlinearities in the mathematical formulation as well as nonlinear mechanical behaviors [27–31]. For example, Ansari et al. [32] developed the nonlinear Timoshenko beam model in order to numerically analyze the nonlinear forced vibration of shear deformable FG-CNTRC beams with different edge supports. Based on the first-order shear deformable beam model and employing the Ritz method and an iteration approach, the effect of various possible types of geometric imperfection on the nonlinear free vibration of FG-CNTRC beams was investigated by Wu et al. [33]. On the basis of the first-order shear deformation plate theory and a numerical solution strategy, the large amplitude vibration analysis of FG-CNTRC rectangular plates subjected to a harmonic excitation was carried out by Ansari et al. [34]. The imperfection sensitivity on the nonlinear postbuckling equilibrium path of FG-CNTRC beams was analyzed by Wu et al. [35]. Furthermore, a third-order shear deformable plate model was proposed by Ansari and Gholami [14] to investigate the effect of important design parameters such as CNT distribution, CNT volume fraction, geometry, boundary condition and excitation frequency on the frequency-response and force-response curves of FG-CNTRC rectangular plates. An improved perturbation technique was utilized by Shen and Xiang [36] to describe the large amplitude free vibration characteristics of nanocomposite cylindrical shells reinforced by single-walled carbon nanotubes (SWCNTs) in thermal environments. Ansari et al. [37] presented an analytical solution to study the nonlinear buckling and postbuckling of FG-CNTRC circular cylindrical shells with piezoelectric layers subjected to the combined electro-thermal loadings, axial compression and lateral forces. Regarding the state of the art of mechanical analyses performed on the FG-CNTRC composites, the large amplitude resonant dynamics of geometrically imperfect higher-order shear deformable FG-CNTRC beams with various edge conditions subjected to a harmonic transverse excitation is still unexplored. Therefore, the objective of the present paper is to investigate the effect of the initial geometrical imperfection on the nonlinear resonant dynamics of FG-CNTRC higher-order shear deformable beams with various edge conditions subjected to the distributed harmonic transverse load. By defining a generalized displacement field as well as applying the von Kármán hypotheses and taking into account the initial imperfection, three coupled nonlinear governing equations associated with the longitudinal, transverse and rotational motions of geometrically imperfect FG-CNTRC beams can be achieved by means of Hamilton’s principle. These equations include the influences of transverse shear deformation, rotary inertia, initial imperfection and geometrical nonlinearity. Furthermore, by selecting an appropriate shape function, the unified developed relations can be reduced to simpler beam models based on the existing beam theories such as the Euler-Bernoulli and Timoshenko beam theories as well as third-order, parabolic, trigonometric, hyperbolic and exponential shear deformation beam theories. Afterwards, an efficient multistep numerical solution approach based on the generalized differential quadrature (GDQ) method, numerical-Galerkin-based scheme and time periodic discretization is utilized to transform the time-dependent nonlinear continuous partial differential equations (PDEs) into a Duffing-type nonlinear set of ordinary differential equations (ODEs). Then, the pseudo arc-length continuation technique is utilized to the nonlinear resonant dynamics responses of FGCNTRC by plotting the frequency-response and force-response

curves; illustrating the effects of initial geometric imperfection, geometric parameters, excitation frequency and boundary conditions. The paper is organized as follows. Section 2 deals with developing a unified coupled nonlinear governing equations of motion of geometrically imperfect higher-order shear deformable FGCNTRC beams. Section 3 sets up an efficient multistep solution procedure used in solving the nonlinear vibration problems. Section 4 provides the numerical results corresponding to the nonlinear resonant dynamics behaviors of imperfect FG-CNTRC beams with various edge conditions. The study is completed by a summary of the main findings in Section 5. 2. Mathematical formulation A schematic representation of shear deformable FG-CNTRC beams with the rectangular cross-section, length L, thickness h and considered Cartesian coordinate system (i.e., 0 6 x 6 L and
h=2 6 z 6 h=2) are illustrated in Fig. 1. The motion of any point of considered beam in the longitudinal and transverse directions are denoted by ux ðt; x; zÞ and uz ðt; x; zÞ. Also, the displacements of any point located on the middle-axis in the longitudinal and transverse directions are represented by uðt; xÞ and wðt; xÞ, respectively; wx ðt; xÞ signifies the rotation of the transverse normal. Moreover, It is assumed that the FG-CNTRC beam subjects an initial geometric imperfection in the positive out-of-plane direction indicated by w ¼ w ðxÞ. Furthermore, a harmonic uniformly-distributed load

t is applied on the FG-CNTRC per unit length Fðt; xÞ ¼
f 0 sin X
and t are forcing beams in the transverse direction, in which
f 0 , X amplitude, excitation frequency and time, respectively.

Fig. 1. Schematic view of a shear deformable FG-CNTRC beams: geometry, kinematic parameters and Cartesian coordinate system.

R. Gholami et al. / Composite Structures 174 (2017) 45–58

2.1. Kinematics of deformations of shear deformable FG-CNTRC beam In order to derive a unified higher-order shear deformable beam model for FG-CNTRC beams, the components of displacement vector are considered as follows

z @wðt;xÞ @x

ux ðt; x; zÞ ¼ uðt; xÞ


  þ !ðzÞ wx ðt; xÞ þ @wðt;xÞ ; @x

uz ðt; x; zÞ ¼ wðt; xÞ:

ð1Þ

where !ðzÞ is generalized shape function; indicating the transverse shear deformation and stress distribution through the thickness of the FG-CNTRC beam. By considering !ðzÞ as following form

Euler—Bernoulli beam theory ðEBTÞ : !ðzÞ ¼ 0;

the functionally graded distribution can be categorized into three main types, namely FGA, FGO and FGX. The volume fraction of CNTs corresponding to any type of distributions in the thickness direction can be written as follows

UD : V cnt ¼ V cnt ; FGA : V cnt

  2z  ¼ 1
V cnt ; h

  2jzj  FGO : V cnt ¼ 2 1
V cnt ; h

Timoshenko beam theory ðTBTÞ : !ðzÞ ¼ z;   4z2 Reddy beam theory ðRBTÞ : !ðzÞ ¼ z 1
3h 2 ;

FGX : V cnt ¼

Parabolic shear deformable beam theory ðPSDBTÞ :   5z2 !ðzÞ ¼ z 54
3h 2 ;

V cnt ¼ ð2Þ

Hyperbolic shear deformable beam theory ðHSDBTÞ :

 !ðzÞ ¼ h sinh hz
z cosh 12 ;

the unified displacement field defined in Eq. (1) can be converted to the displacement filed corresponding to various types of beam theories. 2.2. Strain-displacement relations

ð7cÞ

ð7dÞ

ð3Þ

where 2 @w x e0xx ¼ @u þ 12 @w þ dw ; e1xx ¼ @@xw2 ; e2xx ¼ @w ; @x @x dx @x @x 0 @w cxz ¼ wx þ @x ; !1 ¼ !ðzÞ
z; !2 ¼ !ðzÞ; !3 ¼ d!dzðzÞ :

Kcnt   cnt  :
qqm Kcnt

qcnt qm

ð8Þ

Also, Kcnt stands for the mass fraction of CNT, and qcnt and qm respectively indicates the mass densities of CNT and matrix. Furthermore, it is clear that the volume fraction of matrix phase can be determined as

m E11 ¼ g1 V cnt Ecnt 11 þ V m E ;

g2 E22

g3

Utilizing the displacement field expressed in Eq. (1) and considering the von Kármán geometric nonlinearity and initial geometric imperfection, the non-zero components of strains in terms of displacements can be computed as



ð7bÞ

ð9Þ

On the basis of the rule of mixture, the effective material properties of CNTRC are approximately determined as [38]

2

exx ¼ e0xx þ !1 e1xx þ !2 e2xx ; cxz ¼ !3 c0xz :

Kcnt þ



V m ¼ 1
V cnt :

Exponential shear deformable beam theory ðESDBTÞ :




4jzj  V : h cnt

ð7aÞ

in which

Trigonometric shear deformable beam theory ðTSDBTÞ :
 !ðzÞ ¼ ph sin phz ;

!ðzÞ ¼ ze
2ðz=hÞ

47

2

ð4Þ

G12

ð10aÞ

¼

V cnt V m þ m; E Ecnt 22

ð10bÞ

¼

V cnt V m þ m: G Gcnt 12

ð10cÞ

in which the subscripts/superscript m and cnt represent the matrix phase and carbon nanotubes, respectively; the Young modulus and shear modulus are symbolized by E and G, respectively and the coefficients gj ðj ¼ 1; 2; 3Þ stand for the CNT efficiency parameters describing the scale-dependent material properties. It is worthy remarked that gj ðj ¼ 1; 2; 3Þ must be evaluated by matching the elastic modulus of CNTRCs computed via the molecular dynamics simulations with ones estimated by means of the rule of mixture. Moreover, the effective Poisson’s ratio m and mass density q can be determined as

2.3. Constitutive relations

m ¼ V cnt mcnt þ V m mm :

ð11aÞ

In a continuum made of a linear elastic material, the non-zero components of the stress tensor are expressed through the kinematic parameters as follows

q ¼ V cnt qcnt þ V m qm :

ð11bÞ






rxx ¼ Q 11 ðzÞ e þ ! e þ ! e ; rxz ¼ ks Q 55 ðzÞ! c 0 xx

1 1 xx

2 2 xx

0 3 xz :

ð5Þ

in which ks appeared in Eq. (5) denotes the shear correction factor and equals to 5=6 for Timoshenko beam theory and 1 for higherorder shear deformation beam theories and

Q 11 ðzÞ ¼

E11 ðzÞ 1
mðzÞ2

; Q 55 ðzÞ ¼ G12 ðzÞ:

ð6Þ

In Eq. (6), E11 , G12 and m are the effective material properties of FG-CNTRC beams which should be determined. The considered composite beam is composed of the mixture of the SWCNTs and isotropic matrix. Furthermore, the carbon nanotubes can be dispersed as a uniform distribution (UD) or functionally graded distribution in the thickness direction of considered beam. Also,

where mcnt and mm respectively indicate the Poisson’s ratios of carbon nanotube and matrix. 2.4. Coupled nonlinear governing equations of motion In this subsection, a variational approach on the basis of the Hamilton’s principle is utilized to derive the unified nonlinear governing equations as well as corresponding boundary conditions for FG-CNTRC beams. Based upon the Hamilton’s principle, the variation of total energy of system can be expressed as

Z

t2

t1

ðdPT
dPs þ dPw Þdt ¼ 0:

ð12Þ

where dPT ; dPs and dPw are the variation of kinetic energy, strain energy and work done by external loads, respectively.

48

R. Gholami et al. / Composite Structures 174 (2017) 45–58

By considering the components of strain and stress tensors obtained in Eqs. (3) and (5), the strain energy of higher-order shear deformable FG-CNTRC beam can be expressed as

dPs ¼

Z

8

Z

0

 Nxx de0xx þ Mxx de1xx þ Pxx de2xx þ Q x dc0xz dx:

ð13Þ

where N xx , M xx and Pxx signify the resultant in-plane load, bending moment and higher-order moment, respectively. Also, the shear resultant is denoted by Q x . These parameters are defined as follows

Z

Z

ð1; !1 ; !2 Þrxx dA; Q x ¼ ks

ðNxx ; Mxx ; Pxx Þ ¼ A

rxz !3 dA

ð14Þ

A

Z

fA11 ; B11 ; C 11 ; D11 ; F 11 ; H11 g ¼ Z A55 ¼ Q 55 !23 dA;

A

Q 11 f1; !1 ; !2 ; !21 ; !1 !2 ; !22 gdA; ð15Þ

A

the resultants defined in Eq. (14) can be expressed as follows



@u @x



þ 12


@w2 @x



þ dw dx



@w @x



x ; þ B11 @@xw2 þ C11 @w @x 2


2 dw @w 2 x M xx ¼ B11 @u þ 12 @w þ dx @x þ D11 @@xw2 þ F11 @w ; @x @x @x  
  2 2 @w x Pxx ¼ C11 @u þ 12 @w þ dw ; þ F11 @@xw2 þ H11 @w @x @x dx @x @x
 @w Q x ¼ ks A55 wx þ @x : Z L ( " 2  2 # @u @w @u @ 2 w @u @wx þ 2I1 I0 þ þ 2I2 @t @t @t @t@x @t @t 0 9 !2  2 @w @ 2 w @2w @wx = dx þ 2I3 x þ I4 þ I5 ; @t@x @t @t@x @t Z A

n

o

q 1; !1 ; !2 ; !1 !2 ; !21 ; !22 dA;

ð17Þ

ð18Þ

Finally, the variation of work done due to the external transverse load Fðt; xÞ can be obtained as

dPw ¼

Z

L

ð19Þ

Fðt; xÞdwdx 0

Now, by substituting Eqs. (13), (17) and (19) into (12) and consequently utilizing the fundamental lemma of the calculus of variations, the unified geometrically nonlinear governing equations for imperfect higher-order shear deformable FG-CNTRC beams are achieved in terms of the resultants as follows 2

3

2

@Nxx @ u @ w @ w ¼ I 0 2 þ I 1 2 þ I 2 2x ; @x @t @t @x @t

ð20aÞ

@Pxx @2u @2w @3w
Q x ¼ I2 2 þ I5 2x þ I3 2 ; @x @t @t @x @t

A11

! 2  @ 2 u @w @ 2 w d w @w dw @ 2 w @3w @ 2 wx þ þ þ þ B11 3 þ C11 2 2 2 2 @x @x @x @x dx @x @x2 dx @x

@2u @3w @2w þ I1 2 þ I 2 2x ; ð22aÞ 2 @t @t @x @t   2 2 x x þ @@xw2
D11 @@xw2
F11 @w ks A55 @w @x @x     2 3  2  2   3 2 2 2 @3 w @3 w
B11 @@xu3 þ @@xw2 þ @w þ ddxw3 @w þ ddxw2 @@xw2 þ ddxw2 @@xw2 þ dw @x @x3 @x dx @x3 n  o 2 
2 dw @w 2  2 @ w x þ 12 @w þ dx @x þ B11 @@xw2 þ C11 @w þ ddxw2 þ A11 @u @x @x @x @x2 n 2  o
2 2    3 @2 w @2 w þ ddxw2 @w þ dw þ dw þ A11 @@x2u þ @w þ B11 @@xw3 þ C11 @@xw2x @w @x @x2 @x @x dx @x2 dx 3

u w þFðt; xÞ ¼ I0 @@tw2
I1 @t@2 @x
I3 @t@ 2w@xx
I4 @t@2 @x 2 ; 3

3

ð22bÞ

It is noted that using the developed unified shear deformable beam model, extracting various beam models on the basis of Euler-Bernoulli, Timoshenko and any existing higher-order shear deformable beam theory becomes readily. Using the natural and essential boundary conditions given in Eq. (21), the mathematical expression of various types of end supports can be defined as a. FG-CNTRC beam with simply supported- simply supported end conditions (SS-SS)

u ¼ w ¼ Mxx ¼ Pxx ¼ 0 at edges x ¼ 0; L

ð23Þ

b. FG-CNTRC beam with clamped–clamped end conditions (C–C)

u ¼ w ¼ wx ¼

@w ¼ 0 at edges x ¼ 0; L @x

ð24Þ

c. FG-CNTRC beam with simply supported-clamped end conditions (SS-C)

u ¼ w ¼ Mxx ¼ Pxx ¼ 0 at edges x ¼ 0

ð25Þ

3. Multistep numerical solution approach

ð20bÞ

ð20cÞ

Also, the associated essential and natural boundary conditions can be obtained as

du ¼ 0 or Nxx ¼ 0

Now, inserting the resultants obtained in Eq. (16) into Eq. (20) results in the nonlinear governing equations of imperfect FGCNTRC beams in terms of displacements as follows

u ¼ w ¼ wx ¼ @w ¼ 0 at edges x ¼ L @x

    @Q x @ Mxx @ @w @ dw þ N þ N
þ Fðt; xÞ xx xx @x @x @x @x @x2 dx 2

@2w @3u @3w @3w ¼ I 0 2
I1 2
I3 2 x
I 4 2 2 ; @t @t @x @t @x @t @x

ð21dÞ

! 2  @ 2 u @w @ 2 w d w @w dw @ 2 w @3w @ 2 wx þ þ þ þ H C11 þ F 11 11 2 @x2 @x @x2 @x3 dx @x2 @x2 dx @x   2 2 3 @w @ u @ w @ w
ks A55 wx þ ¼ I 2 2 þ I 5 2x þ I 3 2 ; ð22cÞ @x @t @t @x @t

1 PT ¼ 2

fI0 ; I1 ; I2 ; I3 ; I4 ; I5 ; I6 g ¼

ð21cÞ

@w ¼ 0 or M xx ¼ 0 @x

2

ð16Þ

Moreover, the kinetic energy is calculated as follows

where

d

ð21bÞ

¼ I0

Introducing the following quantities

Nxx ¼ A11

  @w dw @Mxx þ
¼0 @x dx @x

dwx ¼ 0 or Pxx ¼ 0

ðrxx dexx þ rxz dcxz Þd8 L


¼

dw ¼ 0 or Q x þ Nxx

ð21aÞ

In the present study, to investigate the nonlinear resonant dynamics of imperfect FG-CNTRC beams at primary resonance, a
t (
f 0 signifies uniform periodic transverse force Fðt; xÞ ¼
f 0 cos X
indicates the external excitation) is the forcing amplitude, X applied to the FG-CNTRC beam in the lateral direction. Also, it is assumed that the system subjected to a viscous damping and consequently the energy dissipation of system includes into formulation by considering the non-conservative damping load of

49

R. Gholami et al. / Composite Structures 174 (2017) 45–58

Table 1 Convergence study of the present solution for nondimensional natural frequencies of perfect and imperfect clamped FGX-CNTRC beams with different nanotube volume fractions ðL=h ¼ 25; c0 ¼ 0:1Þ. V cnt

Numbers of grid points

Refs. [33]

5

7

9

11

13

15

17

21

0.12

Perfect Imperfect

1.20173 1.21718

1.23050 1.24323

1.20544 1.21736

1.19482 1.20665

1.18782 1.19970

1.18280 1.19447

1.17927 1.18995

1.17509 1.18664

1.1733 1.1845

0.17

Perfect Imperfect

1.49488 1.51280

1.52858 1.54344

1.49748 1.51131

1.48476 1.49848

1.47655 1.49032

1.47463 1.48998

1.47441 1.48897

1.47135 1.48776

1.4750 1.4878

0.28

Perfect Imperfect

1.72066 1.74628

1.76702 1.78786

1.73081 1.75062

1.71413 1.73375

1.70276 1.72249

1.69479 1.71412

1.68934 1.70868

1.68382 1.70614

1.6887 1.7065

viscous type CX_ in which C indicates the damping matrix and can
L ÞK [39] (
c is the damping coefficient). be assumed as as C ¼ ð2
c=x Furthermore, in the course of the numerical calculations, the initial imperfection is considered in the form of the first vibration mode of perfect FG-CNTRC beam associated with each boundary condition with given value of deflection c0 . A multistep numerical solution approach employed in this section help is in predicting the nonlinear frequency-response and force-response curves of geometrically imperfect higherorder shear deformable FG-CNTRC beams with various edge conditions. At this step, first of all, by considering the energy dissipation force and applying the GDQ method on the nonlinear governing equations of motion appeared in Eq. (20), the discretized equations of motions as functions of generalized coordinates un , wn and wn can be derived as follows


t ¼ 0 € þ KX þ CX_ þ Knl ðXÞ þ F cos X MX

ð26Þ

 T where X ¼ uTn ; wTn ; wTn is the displacement vector including the 3n elements of time-dependent generalized coordinates and n is the number of grid points along the x-axis. Also, K and M are the stiffness and mass matrices with 3n  3n elements. Also, F and Knl ðXÞ are the force and nonlinear stiffness vectors, respectively.
t, an eigen_ Knl ðXÞ and F cos X Now, by neglecting the terms CX, value problem corresponding the free vibration of FG-CNTRC beam is obtained in which can be solved to give the mode shapes associated with various mode numbers for the FG-CNTRC beams with different edge supports. By employing the obtained mode shapes and applying the numerical Galerkin approach on Eq. (26), one can obtain the following time-varying Duffing-type set of ordinary differential equations

~ q_ þ K
t ~q ~ þC ~ nl ðqÞ ¼ F ~ cos X € þ Kq M

geometrically imperfect FG-CNTRC beams as the frequencyresponse and force-response curves. It is assumed that the matrix phase is composed of Poly methyl methacrylate (PMMA) with Em ¼ 2:5 GPa; qm ¼ 1150 kg=m3 and mm ¼ 0:34 [38]. Furthermore, the material properties of singlewalled (10, 10) armchair CNTs as the reinforcement are considered cnt cnt as [38]: Ecnt 11 ¼ 5:6466 TPa; E22 ¼ 7:0800 TPa; G12 ¼ 1:9445 TPa; cnt 3 cnt q ¼ 1400 kg=m ; m ¼ 0:175. Moreover, for the three selected carbon nanotube volume fractions, on the basis of the molecular dynamics simulations and extended rule of mixture, the carbon nanotube efficiency parameters appeared in Eq. (10) are considered as follows [5,46]

V cnt ¼ 0:12 : g1 ¼ 0:137; g2 ¼ 1:022; g3 ¼ 0:715 ; V cnt ¼ 0:17 : g1 ¼ 0:142; g2 ¼ 1:626; g3 ¼ 1:138; V cnt ¼ 0:28 : g1 ¼ 0:141; g2 ¼ 1:585; g3 ¼ 1:109 It is remarked that for the subsequent results presented in this section, the non-dimensional (i.e., Non. Dim.) frequencies is pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
;x
ÞL I00 =A110 in which A110 and I00 stands denoted by ðX; xÞ ¼ ðX for the values of A11 and I0 corresponding to the homogeneous matrix beam. Also, the nondimensional applied forcing and maximum (i.e., Max.) vibration amplitudes are defined as f ¼
f 0 L2 =hA110 and W max ¼ wmax =h, respectively. Furthermore, the pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi nondimensional damping coefficient is obtained as c ¼
cL= A110 I00 . First of all, a convergence test is carried out. In this regard, the nondimensional natural frequencies of perfect and imperfect clamped FGX-CNTRC beams with various nanotube volume

ð27Þ

More discussion about the numerical-based Galerkin approach can be found in works done by Ansari et al. [40–42]. The new set of equations can be solved by means of: (1) the time periodic discretization method [43,44] is employed to discretize Eq. (27) on the time domain; (2) By utilizing the linear frequency and the linear solution as the initial guess, the pseudo-arc length continuation technique [45] is utilized for the analysis of geometrically nonlinear resonant responses and specifying the stability and bifurcations. These multistep numerical approach helps us in constructing frequency–responses and force-responses of imperfect shear deformable FG-CNTRC beams with different end conditions. 4. Results and discussions In this section, the developed mathematical formulation and numerical-based multistep solution procedure presented in the previous sections are complemented by some numerical results via plotting the nonlinear resonant dynamics behaviors of

Fig. 2. Nonlinear frequency-response curve of a perfect FGX-CNTRC beam with
 clamped edge condition L=h ¼ 10; V cnt ¼ 0:12; c ¼ 0:01; f ¼ 0:1 .

50

R. Gholami et al. / Composite Structures 174 (2017) 45–58

fractions are tabulated in Table 1 for various numbers of grid points. It can be found that n ¼ 17 is enough to obtained the converged results. Furthermore, the converged frequencies are compared with those given by Wu et al. [33]. It can be seen that the

Fig. 3. Nonlinear frequency-response curve of an imperfect FGX-CNTRC beam with
 clamped edge conditions L=h ¼ 10; V cnt ¼ 0:12; c ¼ 0:01; f ¼ 0:1; c0 ¼ 0:6 .

Fig. 4. Nonlinear force-response curve of an imperfect FGX-CNTRC beam with
 simply-supported boundary conditions L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:01; X=xL ¼ 0:8 .

Fig. 5. Nonlinear force-response curve of an imperfect FGX-CNTRC beam with
 simply-supported boundary conditions L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:01; X=xL ¼ 1:01 .

Fig. 6. Geometric imperfection effect on the nonlinear frequency-response curve of
FGX-CNTRC beam with different edge conditions L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:01; f ¼ 0:08Þ.

R. Gholami et al. / Composite Structures 174 (2017) 45–58

present results are in excellent agreement with those provided by of Wu et al. [33]; illustrating the accuracy and validity of developed model, utilized solution procedure and provided results. Figs. 2 and 3 illustrate the frequency-response curves of the prefect and imperfect FG-CNTRC beams with clamped edge conditions

Fig. 7. Initial geometric imperfection effect on the force-response curve of FGX
 CNTRC beam with different edge conditions L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:01 .

51

at the primary resonance, respectively. The limit point bifurcations as well as the stable and unstable solutions are indicated in these figures. It is remarked that the nonlinear portion of system includes the quadratic and cubic terms and, depending on the dominating term, the system may exhibit the hardening or softening behavior. It can be seen that, due to the cubic nonlinear term which comes into formulation because of the middle-plane stretching, the frequency-response curve of the perfect FG-CNTRC beam exhibits a strong nonlinear hardening spring-type behavior including two limit point bifurcations. However, as displayed in Fig. 3, the frequency-response curve of imperfect FG-CNTRC beam initially illustrates a softening behavior and then follows a hardening-type response. Also, it includes four limit point bifurcations as well as two unstable solutions. It indicates that for an imperfect FG-CNTRC beam with low vibration amplitude, the effective nonlinearity may be dominated by the quadratic nonlinear terms and subsequently the beam exhibits the softening behavior. On the other hand, at the high vibration amplitudes, the cubic nonlinearities dominate the overall geometric nonlinearities; resulting in hardening behavior. Furthermore, it is noted that these hardening and softening behaviors strongly depend on various parameters including the amount of imperfection, geometries and type of boundary condition.

Fig. 8. Effect of length-to-thickness ratio on the nonlinear frequency-response
curve of FGX-CNTRC beam at the primary resonance V cnt ¼ 0:12; c ¼ 0:005; f ¼ 0:1; c0 ¼ 0:6Þ.

52

R. Gholami et al. / Composite Structures 174 (2017) 45–58

Displayed in Figs. 4 and 5 are the force-response curves of an imperfect FG-CNTRC beam with simply-supported edge conditions corresponding to X=xL 6 1 and X=xL P 1, respectively. As illustrated in Fig. 4, for X=xL 6 1, four limit point bifu1rcations with two unstable branches between points A and B as well as points C and D can be found. It is noted that unlike the imperfect systems, the force-response curves of perfect systems for X=xL 6 1 includes no bifurcation point and unstable solution, as reported by Ansari and Gholami [14]. Furthermore, over considered range of forcing, depending on the type of forcing sweep, two jumps from higher to lower or lower to higher values of vibration amplitude exist. While, as represented in Fig. 5, for the case of X=xL P 1, the force-response curve includes two limit point bifurcations with one unstable solution between points A and B. It is noteworthy that in addition to the excitation frequency ratio, the accordance and number of the limit point bifurcations strongly depend on different factors such as the imperfection, geometry and boundary conditions. This issue will be examined in next figures. Illustrated in Fig. 6 is the influence of geometric imperfection on the nonlinear frequency-response curve of higher-order shear deformable FGX-CNTRC beams with various edge supports at the primary resonance. Also, the fundamental natural frequencies as well as the frequency-response curves corresponding to the nonlinear free vibration are provided in Fig. 6. The fundamental fre-

Fig. 9. Nonlinear free vibration response of FGX-CNTRC beam with various length
 to-thickness ratio V cnt ¼ 0:12; c ¼ 0:005; f ¼ 0:1; c0 ¼ 0:6 .

quency of system with various end conditions increases as the magnitude of initial imperfection increases. Also, it is evident from this figure that the increase in the value of magnitude of initial geometric imperfection results in weakening the hardening-type nonlinearity and transforming the geometric nonlinearity in the FG-CNTRC beams from hardening type to softening type. Additionally, depending on the value of the magnitude of the initial imperfection, the frequency-response curve of system indicates either two or four limit point bifurcations as well as one or two unstable branches. Furthermore, it can be found that the initial geometric imperfection leads to more intensifying the softening behavior of FG-CNTRC beams with SS–SS boundary conditions compared with those with C–C and C–SS edge supports. Also, Fig. 7 displays the force-response curves of shear deformable FGX-CNTRC beams with various boundary conditions corresponding to various magnitudes of initial geometric imperfections. Ii is noted that the curves are provided in the case of X=xL 6 1. It is illustrated that for an perfect FG-CNTRC beam ðc0 ¼ 0Þ in the case of X=xL 6 1, no limit bifurcation point, unstable solution and jump phenomena can be found in the force-response curve of system. Therefore, the maximum vibration amplitude of FG-CNTRC beam gradually increases with increasing the forcing amplitude. But, for an imperfect FG-CNTRC

Fig. 10. Effect of length-to-thickness ratio on the force-response curve of
FGX-CNTRC beam at the primary resonance V cnt ¼ 0:12; c ¼ 0:01; X=xL ¼ 0:94; c ¼ 0:01; c0 ¼ 0:6Þ.

R. Gholami et al. / Composite Structures 174 (2017) 45–58

Fig. 11. The frequency-response curve of CNTRC beams associated with various
 CNT distribution L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:005; f ¼ 0:1; c0 ¼ 0:3 .

53

Fig. 12. The force-response curve of CNTRC beams associated with various CNT
 distribution L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:005; X=xL ¼ 0:9; c0 ¼ 0:6 .

54

R. Gholami et al. / Composite Structures 174 (2017) 45–58

Fig. 13. Nonlinear force-response curve of imperfect FGX-CNTRC beams associated
 with various excitation frequency ratio ðX=xL Þ L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:01 .

Fig. 14. Comparisons of the frequency-response curve of FGX-CNTRC beams
predicted by various beam theories L=h ¼ 14; V cnt ¼ 0:12; c ¼ 0:004; f ¼ 0:1; c0 ¼ 0:6Þ.

55

R. Gholami et al. / Composite Structures 174 (2017) 45–58

beam, depending on the value of imperfection, the force-response curve may be displays either two or four limit point bifurcations as well as one or two unstable solutions. However, for the case of small initial imperfections, the force-response curves include no limit bifurcation point or unstable solution. Figs. 8 and 9 show the effect of length-to-thickness ratio ðL=hÞ on the frequency-response curves corresponding the geometrically nonlinear forced and free vibrations of the imperfect shear deformable FGX-CNTRC beams. It can be seen that decreasing the lengthto-thickness ratio strengthens the nonlinear initial softening at the small vibration amplitude and the later softening behavior at the higher vibration amplitudes. Moreover, for the FGX-CNTRC beams, at the primary resonance, due to decreasing the total stiffness of the system, the peak of vibration amplitude increases as the length-to-thickness ratio increases. Furthermore, Fig. 10 highlights the effect of length-to-thickness ratio on the force–response curves of the imperfect FGX-CNTRC beams. This figure illustrates that increasing the length-to-thickness ratio results in occurring the limit-point bifurcations at lower forcing amplitudes and larger vibration amplitudes. Illustrated in Fig. 11 is the frequency-response curve of imperfect CNTRC beams with various types of CNT distribution. It is evident from this figure that in the case of distribution of CNTs near the top and bottom surfaces, due to more intensifying the effective stiffness of beams, the system displays a weaker nonlinear behavior. For this reason, it can be seen that CNTRC beam with FGX-type distribution possesses highest strength and lowest nonlinear hardening-type behavior compared with those of with UD, FGA and FGO distributions. Furthermore, the influence of CNT distribution on the force-response curve of imperfect CNTRC beam with various boundary conditions are plotted in Fig. 12. It can be seen that depending on the type of boundary condition and CNT distribution, the force-response curves illustrate either two or four limit point bifurcations as well as one or two unstable branch/branches. Also, it can be found that for the CNT distribution with higher effective material properties, the bifurcation limit points occur at higher forcing amplitudes. The influence of excitation frequency ratio ðX=xL Þ on the forceresponse curve of imperfect FGX-CNTRC beams with various end conditions is addressed in Fig. 13. It is seen in the figure that increasing the excitation frequency ratio results in occurrence of

the limit point bifurcations at higher forcing amplitudes. However, the nonlinear force response curves of imperfect FG-CNTRC beams in the cases X=xL P 1 and X=xL < 1 may be completely different. In the case of X=xL P 1, regardless the type of boundary condition and excitation frequency ratio, the system includes two limit point bifurcations as well as a unstable branch. While, in the case of X=xL < 1, depending on the geometry, boundary condition and excitation frequency ratio, the force response curve may be consists of either two or four limit point bifurcations as well as either one or two unstable branch/branches. Finally, the comparisons of the frequency-response curve of imperfect FGX-CNTRC beams corresponding to various beam theories are provided in Fig. 14. It can be seen that the nondimensional natural frequencies predicted by TBT is higher than those given by higher-order shear deformation beam theories. Moreover, the difference between the frequency-response curve of TBT and shear deformation beam theories especially in the case of C–C edge condition is significant. Furthermore, the curves given by all higherorder beam theories are almost same with some negligible differences. Furthermore, the comparisons of nondimensional natural frequencies of imperfect FGX-CNTRC beams with various edge conditions and length-to-thickness ratios predicted by various beam theories are provided in Table 2. It can be found that in the thick and moderately thick beams a considerable difference can be found in the natural frequencies. However, increasing the lengthto-thickness ratio of beams results in decreasing these differences and consequently the simple theories such as TBT or EBT can be used to predict the vibrational behavior of considered system. 5. Concluding remarks The geometrically nonlinear resonant dynamics of an initiallyimperfect higher-order shear deformable FG-CNTRC beams with various boundary conditions subjected to a transverse harmonic excitation was numerically analyzed. In this regard, a unified higher-order shear deformable beam model was developed through use of Hamilton’s principle along with utilizing the von Kármán kinematic nonlinearity and taking into account the initial geometric imperfection. It was considered that the initial imperfection is in the form of the first vibration mode of perfect FG-CNTRC beam associated with each boundary condition. To numerically

Table 2 Comparisons of the nondimensional natural frequencies of imperfect FGX-CNTRC beams corresponding to different length-to-thickness ratios, types of theories and boundary
 conditions V cnt ¼ 0:12; c0 ¼ 0:2 . Boundary condition

L=h

TBT

RBT

PSDBT

TSDBT

HSDBT

ESDBT

C–C

10 12 14 16 20 30

1.60752 1.52073 1.45327 1.39588 1.29619 1.08775

1.74679 1.62919 1.53976 1.46603 1.34406 1.10900

1.74472 1.62689 1.53728 1.46341 1.34130 1.10639

1.76416 1.64272 1.55059 1.47486 1.35016 1.11178

1.56703 1.51107 1.56234 1.47711 1.35547 1.11244

1.78420 1.65857 1.56349 1.48555 1.35776 1.11545

SS-C

10 12 14 16 20 30

1.50187 1.39244 1.30516 1.23115 1.10657 0.87230

1.57433 1.44744 1.34807 1.26504 1.12831 0.88063

1.57352 1.44646 1.34696 1.26384 1.12705 0.87953

1.58386 1.45477 1.35384 1.26964 1.13134 0.88185

1.59144 1.46474 1.36042 1.27186 1.13059 0.88544

1.59526 1.46372 1.36102 1.27547 1.13530 0.88354

SS-SS

10 12 14 16 20 30

1.42713 1.28586 1.17040 1.07280 0.91577 0.66081

1.43376 1.28989 1.17296 1.07449 0.91658 0.66102

1.43376 1.28989 1.17296 1.07449 0.91658 0.66102

1.43617 1.29158 1.17420 1.07541 0.91713 0.66121

1.43818 1.29021 1.17403 1.07144 0.91655 0.66101

1.43978 1.29424 1.17620 1.07695 0.91808 0.66157

Type of theory

56

R. Gholami et al. / Composite Structures 174 (2017) 45–58

analyze the nonlinear resonant dynamics characteristics of FGCNTRC beams, an efficient multistep numerical solution approach including the GDQ method, numerical-Galerkin-based scheme, time periodic discretization and pseudo-arc length continuation technique was employed to obtain the nonlinear frequencyresponse and force-response curves of considered system at the primary resonance; illustrating the effects of initial geometric imperfection, geometric parameters, CNT distribution, excitation frequency and boundary condition. Investigating the frequency–response curves of the FG-CNTRC beams indicated that in the case of perfect and imperfect beams with small magnitude of initial geometric imperfection, the considered beams exhibit only a nonlinear hardening-type behavior and the frequency-response curve includes two limit point bifurcations and an unstable branch. However, for higher values of initial imperfection, both softening and hardening behaviors can be observed in the frequency response curve. Moreover, the frequency-response curve includes four limit point bifurcations and two unstable branches. Furthermore, increasing the magnitude of the initial imperfection results in intensifying the initial softening behavior of the system, especially for the beams with SS-SS end conditions. Moreover, analyzing the force-response curves of the FG-CNTRC beams with various edge conditions revealed that in the case of X=xL < 1, depending on the magnitude of initial geometric imperfection and boundary conditions, the system includes either two or four limit point bifurcations. But, in the case of X=xL P 1, regardless the value of initial imperfection, the system shows two limit-point bifurcations and one unstable solution.


t ¼ 0 € þ KX þ CX_ þ Knl ðXÞ þ F cos X MX where

2

kuu

6 K ¼ 4 kwu kwu

The details about the multistep numerical solution procedure of the nonlinear resonant dynamics of imperfect FG-CNTRC beams at primary resonance are given herein. In this regard, first of all, using the shifted Chebyshev–Gauss–Lobatto grid points, the grid points in x-direction are generated as follows

xi ¼

  L i
1 1
cos p ; i ¼ 1; 2; . . . ; n: 2 n
1

ðA-1Þ

Furthermore, assuming the values of displacement variables at the defined grid points as ui ¼ uðxi Þ; wi ¼ wðxi Þ and wi ¼ wx ðxi Þ, the column vectors un ; wn and wn indicating the displacement vectors corresponding to the axial motion, transverse motion and rotation are considered as follows

un ¼ ½ u1 wn ¼ ½ w1

u2 w2

Moreover,

   un T ; wn ¼ ½ w1

w2

   wn T ;

T

   wn  : the

imperfection

ðA-2Þ vector

is

defined

as

kuw

kuw

3

kww

7 kww 5;

kwu

kwu

kuu ¼ A11 Dð2Þ x þ A11



    ð1Þ ð1Þ  ð2Þ  Dð2Þ ; x wn }Dx þ Dx wn }Dx

ð2Þ kuw ¼ B11 Dð3Þ x ; kuw ¼ C 11 Dx ;      ð1Þ ð1Þ  ð2Þ  kwu ¼ B11 Dð3Þ Dð2Þ ; x þ A11 x wn }Dx þ Dx wn }Dx      ð2Þ ð2Þ ð3Þ  ð1Þ ð2Þ  kww ¼ ks A55 Dx
D11 Dx
B11 Dx wn }Dx þ Dx wn }Dð2Þ x           ð1Þ ð1Þ  ð1Þ  ð2Þ  þA11 2 Dxð1Þ wn  Dð2Þ ; x wn }Dx þ Dx wn  Dx wn }Dx      ð1Þ ð1Þ ð2Þ  ð1Þ ð1Þ  ð2Þ kww ¼ ks A55 Dx
F11 Dx þ C11 Dx wn }Dx þ C 11 Dx wn }Dx ;      ð1Þ ð1Þ  ð2Þ  kwu ¼ C 11 Dð2Þ Dð2Þ x ; kww ¼ C11 x wn }Dx þ Dx wn }Dx ð1Þ þF11 Dð3Þ x
ks A55 Dx

kww ¼ H11 Dð2Þ x
ks A55 Ix : ðA-4Þ

2

I 0 Ix 6 ð1Þ M ¼
6
I 4 1 Dx I2 I x F ¼ ½ 0 Fw

Appendix A

ðA-3Þ

0 T

I1 Dxð1Þ I0 Ix
I4 Dxð2Þ I3 Dð1Þ x

I2 I x

3

7 7;
I3 Dð1Þ x 5

ðA-5Þ

I 5 Ix

Fw ¼ ½. . . ;
f 0 ; . . .1n ;

h iT Knl ðXÞ ¼ KTu ðXÞ KTw ðXÞ KTw ðXÞ :

ðA-6Þ ðA-7Þ

Also, the components of Knl ðXÞ are expressed as

    ð2Þ Ku ðXÞ ¼ A11 Dð1Þ x wn  Dx wn ;         ð2Þ ð1Þ ð3Þ Kw ðXÞ ¼
B11 Dð2Þ x wn  Dx wn þ Dx wn  Dx wn    1  ð1Þ   ð1Þ   ð1Þ    ð1Þ  þ A11 Dð1Þ D u þ w w w w  D þ D  D n n n n x x x x n 2 x   A     11 ð1Þ ð1Þ   Dð2Þ þ B11 Dð2Þ Dð2Þ x wn þ C11 Dx wn x wn þ x wn  Dx wn 2         ð1Þ  ð1Þ ð2Þ  Dð1Þ x wn þ A11 Dx wn  Dx wn  Dx wn ; n          ð1Þ ð2Þ ð1Þ ð2Þ  þ A11 Dð2Þ x un þ Dx wn  Dx wn þ Dx wn  Dx wn     o   ð2Þ ð2Þ ð1Þ  þ Dð1Þ þ B11 Dð3Þ x wn  Dx wn x wn þ C11 Dx wn  Dx wn ;     ð2Þ Kw ðXÞ ¼ C11 Dð1Þ ðA-8Þ x wn  Dx wn :

T

wn ¼ ½ w1 w2    wn  . Applying the GDQ method on Eqs. (22a)–(22c) and considering the uniform periodic transverse force
t and CX, _ and damping load of viscous type as Fðt; xÞ ¼
f 0 cos X respectively, result in the following equation

h

DðrÞ x

i

ðrÞ

ij

¼ Wij

where  and } identify the Hadamard and SJT product [47]. Furthermore, the weighting coefficients associated with the first derivative and higher-order derivatives can be calculated via the following relation

8 where Iij is a n  n identity matrix r ¼ 0 Iij > > > Pðx Þ > i > ; i; j ¼ 1; . . . ; n and i – j and r ¼ 1 > > < ðxi
xj ÞPðxj Þ

 n ðr
1Þ X ¼ Wij ð1Þ ðr
1Þ ðrÞ > Wik ; i ¼ j i; j ¼ 1; . . . ; n and r P 2 ; i–j and
W
r W > ij ii xi
xj > > > > k¼1 > : k–i

ðA-9Þ

R. Gholami et al. / Composite Structures 174 (2017) 45–58

where Pðxi Þ ¼

n Q

ðxi
xk Þ. The mathematical expressions of

k¼1;i–k

boundary conditions can be discretized similarly. Now, the numerical-based Galerkin method can be used to reduce Eq. (A-3) to a time-dependent ordinary differential equa_ Knl ðXÞ and tion of Duffing type. In this regard, first, the terms CX;
F cos Xt appeared in Eq. (A-3) are neglected. Then, considering a

~ jx
t and harmonic solution for the considered system as X ¼ Xe assuming the initial imperfection as the first vibration mode of perfect FG-CNTRC beam associated with each boundary condition and given value of deflection, the following eigenvalue problem can be achieved



~ ¼
x ~ X ~¼ u
2 MX; ~T KX

~T w

~T w

T

ðA-10Þ

:

It is remarked that in Eq. (A-10), the discretized boundary conditions components associated with boundary grid points are substituted with the identical components in the earlier matrices. Solving Eq. (A-10) gives the linear natural frequency and corresponding mode shapes (eigenvectors). Selecting the first m eigenvectors, the solution of system can be expressed as

X ¼ Uq:

ðA-11Þ

where q and U are the reduced generalized coordinates and base functions used in the numerical-based Galerkin technique and can be written as

h iT ðmÞ ð1Þ ð2Þ ðmÞ ð1Þ ð2Þ ðmÞ q3m1 ¼ quð1Þ qð2Þ :    qu qw qw    qw qw qw    qw u ðA-12Þ 2

Uu

6 U¼4 0 0

0 Uw 0

0

3

7 0 5 : Uw 3n3m

ðA-13Þ

with

~ ð1Þ Uu ¼ ½ u n1

~ ðmÞ ~ ð1Þ  ; Uw ¼ ½ w  u n1 nm n1

~ ð1Þ Uw ¼ ½ w n1

~ ðmÞ  :  w n1 nm

~ ðmÞ  ;  w n1 nm

Substituting Eq. (A-11) into (A-3) results in the following residual vector

ðA-15Þ

Based upon the numerical-based Galerkin approach, by defining a matrix operator as following form which enable us to simultaneously multiply each equation by associated eigenvector and integrate over the space domain

Gm3n ¼ UT diagðSÞ; S ¼ ½ Sx

Sx

Sx 13n

ðA-16Þ

where Sn indicates the integral operator [41], one can obtain the following relation

~ q_ þ K
t ~q ~ þC ~ nl ðqÞ ¼ F ~ cos X € þ Kq M

57

ðA-19Þ

where Q is the discretized form of q in Eq. (A-17) on the time domain, DðÞ t symbolizes the time differentiation matrix operator. Furthermore, v ecðYÞ signifies the vectorization of matrix Y and  indicates the Kronecker product. Also, It represent an identity matrix. Eq. (A-19) can be expressed as the following set of nonlinear parameterized equations


Þ ¼ 0: H : R3mNt þ1 ! R3mNt ; HðvecðQ Þ; X

ðA-20Þ

where N t denotes the number of discrete points in the time domain. The pseudo-arc length continuation technique [45] can be used to construct the frequency–responses and force-responses of imperfect shear deformable FG-CNTRC beams with different end conditions. It is remarked that the pseudo-arc length continuation technique is a appreciate approach to approximately obtain the solution set of a system of nonlinear equations such as FðX; pÞ ¼ 0, for various values of the system parameter p. To proceed with the algorithm, an initial point on the solution path should be chosen. Each point on this procedure contains two separate steps known as the prediction step and the correction step. At the first step, utilizing the unit length tangent vector X_ i at Xi , the next point on the path, Xpiþ1 is predicted (i.e., Xpiþ1 ¼ Xi þ DsX_ i , where Ds > 0 signifies the step size). Then, using the Newton iteration method, the prediction made at the previous stage is improved. In this step, in addition to improving the predicted point Xpiþ1 , one can found the next solution point Xiþ1 . Furthermore, it should be noted that the convergence of pseudo-arc length continuation technique is strongly dependent on the initial guess which can be assess using various schemes including the random guess, globalized Newton approach and homotopy method. References

ðA-14Þ


t € þ KUq þ CUq_ þ Knl ðUqÞ þ F cos X R ¼ MUq

        
2 ð2Þ X ~ þ It  K ~ þ X
Dð1Þ  C ~ v ecðQ Þ D  M t t 2p 2p     ~ nl ðQ Þ
It  F ~ AT ¼ 0: þv ec K

ðA-17Þ

where

~ ¼ GCU; K ~ ¼ GMU; K ~ ¼ GKU; C ~ nl ðqÞ ¼ GKnl ðUqÞ; F ~ ¼ GF: M ðA-18Þ It is remarked that using the numerical-based Galerkin approach enables us to reduce the general coordinates form 3n to 3 m. Now, by means of the time periodic discretization method [43,44], Eq. (A-17) can be discretized on the time domain and can be expressed as

[1] Sumfleth J, Prehn K, Wichmann MH, Wedekind S, Schulte K. A comparative study of the electrical and mechanical properties of epoxy nanocomposites reinforced by CVD-and arc-grown multi-wall carbon nanotubes. Compos Sci Technol 2010;70:173–80. [2] Manchado ML, Valentini L, Biagiotti J, Kenny J. Thermal and mechanical properties of single-walled carbon nanotubes–polypropylene composites prepared by melt processing. Carbon 2005;43:1499–505. [3] Ho Y, Chang C, Shyu F, Chen R, Chen S, Lin M-F. Electronic and optical properties of double-walled armchair carbon nanotubes. Carbon 2004;42:3159–67. [4] Yakobson BI, Avouris P. Mechanical properties of carbon nanotubes. In: Carbon nanotubes. Springer; 2001. p. 287–327. [5] Han Y, Elliott J. Molecular dynamics simulations of the elastic properties of polymer/carbon nanotube composites. Comput Mater Sci 2007;39:315–23. [6] Zare Y. Effects of interphase on tensile strength of polymer/CNT nanocomposites by Kelly-Tyson theory. Mech Mater 2015;85:1–6. [7] Aguilar-Elguézabal A, Bocanegra-Bernal M. Fracture behaviour of a-Al2O3 ceramics reinforced with a mixture of single-wall and multi-wall carbon nanotubes. Compos B Eng 2014;60:463–70. [8] Morsi K, Esawi A, Borah P, Lanka S, Sayed A. Characterization and spark plasma sintering of mechanically milled aluminum-carbon nanotube (CNT) composite powders. J Compos Mater 2010. [9] Ashrafi B, Hubert P, Vengallatore S. Carbon nanotube-reinforced composites as structural materials for microactuators in microelectromechanical systems. Nanotechnology 2006;17:4895. [10] Lau K-T, Gu C, Gao G-H, Ling H-Y, Reid SR. Stretching process of single-and multi-walled carbon nanotubes for nanocomposite applications. Carbon 2004;42:426–8. [11] Kashyap K, Koppad PG, Puneeth K, Ram HA, Mallikarjuna H. Elastic modulus of multiwalled carbon nanotubes reinforced aluminium matrix nanocomposite – a theoretical approach. Comput Mater Sci 2011;50:2493–5. [12] Ci L, Ryu Z, Jin-Phillipp NY, Rühle M. Investigation of the interfacial reaction between multi-walled carbon nanotubes and aluminum. Acta Mater 2006;54:5367–75.

58

R. Gholami et al. / Composite Structures 174 (2017) 45–58

[13] Alfonso I, Navarro O, Vargas J, Beltrán A, Aguilar C, González G, et al. FEA evaluation of the Al4C3 formation effect on the Young’s modulus of carbon nanotube reinforced aluminum matrix composites. Compos Struct 2015;127:420–5. [14] Ansari R, Gholami R. Nonlinear primary resonance of third-order shear deformable functionally graded nanocomposite rectangular plates reinforced by carbon nanotubes. Compos Struct 2016;154:707–23. [15] Zhu P, Lei Z, Liew KM. Static and free vibration analyses of carbon nanotubereinforced composite plates using finite element method with first order shear deformation plate theory. Compos Struct 2012;94:1450–60. [16] Rafiee M, Yang J, Kitipornchai S. Large amplitude vibration of carbon nanotube reinforced functionally graded composite beams with piezoelectric layers. Compos Struct 2013;96:716–25. [17] Kamarian S, Salim M, Dimitri R, Tornabene F. Free vibration analysis of conical shells reinforced with agglomerated carbon nanotubes. Int J Mech Sci 2016;108:157–65. [18] Huang B, Guo Y, Wang J, Du J, Qian Z, Ma T, et al. Bending and free vibration analyses of antisymmetrically laminated carbon nanotube-reinforced functionally graded plates. J Compos Mater 2016. p. 0021998316685165. [19] Song Z, Zhang L, Liew K. Active vibration control of CNT-reinforced composite cylindrical shells via piezoelectric patches. Compos Struct 2016;158:92–100. [20] Ansari R, Torabi J, Shakouri AH. Vibration analysis of functionally graded carbon nanotube-reinforced composite elliptical plates using a numerical strategy. Aerosp Sci Technol 2017;60:152–61. [21] Lei Z, Liew KM, Yu J. Buckling analysis of functionally graded carbon nanotubereinforced composite plates using the element-free kp-Ritz method. Compos Struct 2013;98:160–8. [22] Mehrabadi SJ, Aragh BS, Khoshkhahesh V, Taherpour A. Mechanical buckling of nanocomposite rectangular plate reinforced by aligned and straight singlewalled carbon nanotubes. Compos B Eng 2012;43:2031–40. [23] George N, Jeyaraj P, Murigendrappa S. Buckling and free vibration of nonuniformly heated functionally graded carbon nanotube reinforced polymer composite plate. Int J Struct Stab Dyn 2016. p. 1750064. [24] Shenas AG, Malekzadeh P, Ziaee S. Vibration analysis of pre-twisted functionally graded carbon nanotube reinforced composite beams in thermal environment. Compos Struct 2017;162:325–40. [25] Wu C-P, Li H-Y. Three-dimensional free vibration analysis of functionally graded carbon nanotube-reinforced composite plates with various boundary conditions. J Vib Control 2016;22:89–107. [26] Ansari R, Torabi J, Shojaei MF, Hasrati E. Buckling analysis of axially-loaded functionally graded carbon nanotube-reinforced composite conical panels using a novel numerical variational method. Compos Struct 2016;157:398–411. [27] Mirzaei M, Kiani Y. Nonlinear free vibration of temperature-dependent sandwich beams with carbon nanotube-reinforced face sheets. Acta Mech 2016;227:1869–84. [28] Wu H, Kitipornchai S, Yang J. Thermal buckling and postbuckling analysis of functionally graded carbon nanotube-reinforced composite beams. Appl Mech Mater 2016;846. [29] Guo X, Zhang W. Nonlinear vibrations of a reinforced composite plate with carbon nanotubes. Compos Struct 2016;135:96–108. [30] Chaudhari VK, Lal A. Nonlinear free vibration analysis of elastically supported nanotube-reinforced composite beam in thermal environment. Proc Eng 2016;144:928–35.

[31] Zhang L, Lei Z, Liew K, Yu J. Large deflection geometrically nonlinear analysis of carbon nanotube-reinforced functionally graded cylindrical panels. Comput Methods Appl Mech Eng 2014;273:1–18. [32] Ansari R, Shojaei MF, Mohammadi V, Gholami R, Sadeghi F. Nonlinear forced vibration analysis of functionally graded carbon nanotube-reinforced composite Timoshenko beams. Compos Struct 2014;113:316–27. [33] Wu H, Yang J, Kitipornchai S. Nonlinear vibration of functionally graded carbon nanotube-reinforced composite beams with geometric imperfections. Compos B Eng 2016;90:86–96. [34] Ansari R, Torabi J, Shakouri AH. Vibration analysis of functionally graded carbon nanotube-reinforced composite elliptical plates using a numerical strategy. Aerosp Sci Technol 2016. [35] Wu H, Kitipornchai S, Yang J. Imperfection sensitivity of thermal post-buckling behaviour of functionally graded carbon nanotube-reinforced composite beams. Appl Math Model 2017;42:735–52. [36] Shen H-S, Xiang Y. Nonlinear vibration of nanotube-reinforced composite cylindrical shells in thermal environments. Comput Methods Appl Mech Eng 2012;213:196–205. [37] Ansari R, Pourashraf T, Gholami R, Shahabodini A. Analytical solution for nonlinear postbuckling of functionally graded carbon nanotube-reinforced composite shells with piezoelectric layers. Compos B Eng 2016;90:267–77. [38] Shen H-S. Nonlinear bending of functionally graded carbon nanotubereinforced composite plates in thermal environments. Compos Struct 2009;91:9–19. [39] Ibrahim S, Patel B, Nath Y. Modified shooting approach to the non-linear periodic forced response of isotropic/composite curved beams. Int J NonLinear Mech 2009;44:1073–84. [40] Ansari R, Gholami R. Size-dependent nonlinear vibrations of first-order shear deformable magneto-electro-thermo elastic nanoplates based on the nonlocal elasticity theory. Int J Appl Mech 2016;8:1650053. [41] Ansari R, Gholami R. Surface effect on the large amplitude periodic forced vibration of first-order shear deformable rectangular nanoplates with various edge supports. Acta Astronaut 2016;118:72–89. [42] Gholami R, Ansari R. A most general strain gradient plate formulation for sizedependent geometrically nonlinear free vibration analysis of functionally graded shear deformable rectangular microplates. Nonlinear Dyn 2016;84:2403–22. [43] Ansari R, Gholami R, Shahabodini A. Size-dependent geometrically nonlinear forced vibration analysis of functionally graded first-order shear deformable microplates. J Mech 2016;32:539–54. [44] Ansari R, Mohammadi V, Shojaei MF, Gholami R, Sahmani S. On the forced vibration analysis of Timoshenko nanobeams based on the surface stress elasticity theory. Compos B Eng 2014;60:158–66. [45] Keller HB. Numerical solution of bifurcation and nonlinear eigenvalue problems, applications of bifurcation theory. In: Proc. advanced sem., Univ. Wisconsin, Madison, Wis.. New York: Academic Press; 1977. p. 359–84. [46] Griebel M, Hamaekers J. Molecular dynamics simulations of the elastic moduli of polymer–carbon nanotube composites. Comput Methods Appl Mech Eng 2004;193:1773–88. [47] Ansari R, Gholami R, Shojaei MF, Mohammadi V, Darabi M. Coupled longitudinal-transverse-rotational free vibration of post-buckled functionally graded first-order shear deformable micro-and nano-beams based on the Mindlin’s strain gradient theory. Appl Math Model 2016;40:9872–91.