A study on large amplitude vibration of multilayered graphene sheets

Computational Materials Science 50 (2011) 1043–1051

Contents lists available at ScienceDirect

Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci

A study on large amplitude vibration of multilayered graphene sheets E. Jomehzadeh a,b,⇑, A.R. Saidi b a b

Shahid Bahonar University of Kerman, Young Researchers Society, Kerman, Iran Department of Mechanical Engineering, Shahid Bahonar University of Kerman, Kerman, Iran

a r t i c l e

i n f o

Article history: Received 5 July 2010 Received in revised form 25 October 2010 Accepted 29 October 2010 Available online 26 November 2010 Keywords: Small scale Nonlinear vibration Multilayered graphene sheet Nonlocal continuum

a b s t r a c t In the present article, large amplitude vibration analysis of multilayered graphene sheets is presented and the effect of small length scale is investigated. Using the Hamilton’s principle, the coupled nonlinear partial differential equations of motion are obtained based on the von Karman geometrical model and Eringen theory of nonlocal continuum. The solutions of free nonlinear vibration, based on the harmonic balance method, are found for graphene sheets with three different boundary conditions. For numerical results single, double and triple layered graphene sheets with both armchair and zigzag geometries are considered. The results obtained herein are compared with those available in the literature for linear vibration of multilayered graphene sheets and an excellent agreement is seen. Also, the effects of number of layers, geometric properties and small scale parameter on nonlinear behavior of graphene sheet are discussed in details. Ó 2010 Elsevier B.V. All rights reserved.

1. Introduction Graphene is defined as a flat one-atom-thick of carbon tightly packed into a two-dimensional honeycomb lattice. It is observed to have promising mechanical and electrical properties. Graphene is the thinnest known material in the universe and the strongest ever measured. Its charge carriers exhibit giant intrinsic mobility and can travel for micrometers without scattering at room temperature. Graphene can sustain current densities six orders of magnitude higher than that of copper, shows record thermal conductivity and stiffness, and reconciles such conflicting qualities as brittleness and ductility. Electrons whiz through graphene at rates far beyond those achieved in other materials. These exceptional properties of graphene sheets have led to multiple usages in the field of nanoelectronics, nanodevices, nanosensors, nanooscillators, nanoactuators, nanoresonators, nanocomposites and nanooptomechanical systems [1,2]. Due to its lightness and stiffness, there is also fast-growing interest in graphene as a base material for nano-electro-mechanical systems (NEMS) for sensing applications [3,4]. Since the graphene exhibits the highest stiffness to inertia ratio [4], it can be used in high frequency nanomechanical devices. Thin graphene resonators have recently been demonstrated on small scales and have shown extremely high stiffness to weight ratio [3,5]. Graphene-based resonators offer low inertial masses, ultrahigh frequencies, and in comparison with nanotubes, low-resistance contacts that are essential for matching the impedance of ⇑ Corresponding author. at: Shahid Bahonar University of Kerman, Young Researchers Society, Kerman, Iran. Tel.: +98 341 2111763; fax: +98 341 2120964. E-mail address: [email protected] (E. Jomehzadeh). 0927-0256/$ - see front matter Ó 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.commatsci.2010.10.045

external circuits [1]. All these characteristics have made graphene a playground for researchers including theoretical and high energy physicist looking to lend graphene’s exceptional properties to tomorrow’s ultrasmall technology. Due to the vast computational expenses of nano-structures analyses when using atomic lattice dynamics and molecular dynamic simulations, there is a great interest in applying continuum mechanics for analysis of such structures. Many studies have been carried out for vibration analysis of nano-structures such as nanotubes and graphene sheets or nano-plates. Lima and He [6] developed a von Karman type nonlinear model for ultra-thin, elastically isotropic films with surface effect. Nanoscale vibrational analysis of a multilayered graphene sheet embedded in an elastic medium was investigated by Behfar and Naghdabadi [7] and the corresponding natural frequencies and the associated modes were determined. Kitipornchai et al. [8] used the continuum plate model for mechanical analysis of graphene sheets. He et al. [9] described van der Waals interaction between the layers of multilayered graphene sheets by an explicit formula. Based on the continuum mechanics and a multiple-elastic beam model, Fu et al. [10] investigated the nonlinear free vibration analysis of embedded carbon nanotubes. Liew et al. [11] studied the vibration behavior of multilayered graphene sheets embedded in an elastic matrix based on a continuum-based plate model. Pradhan and Phadikar [12] presented classical and first order shear deformation plate theories for vibration analysis of nano-plates. Their approach was based on the Navier solution and for a nano-plate with all edges simply supported. Ke et al. [13] investigated the nonlinear free vibration of embedded double-walled carbon nanotubes based on the Eringen’s nonlocal elasticity theory and von Karman geometric nonlinearity using differential quadrature method. Dong and Lim [14]

1044

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051 0

point X in the body. The nonlocal stress tensor r at point X can be expressed as

studied the nonlinear free vibrations of a nano-beam with simply supports boundary conditions based on nonlocal elasticity theory. Murmu and Pradhan [15] developed a single elastic beam model for thermo-mechanical vibration of a single-walled carbon nanotube embedded in an elastic medium based on the nonlocal elasticity theory. Duan and Wang [16] characterized the bending and stretching of the circular graphene sheet using von Karman plate theory. Pradhan and Phadikar [17] presented a vibration analysis of multilayered graphene sheets embedded in polymer matrix by employing the nonlocal continuum mechanics. A nonlocal plate model was developed by Ansari et al. [18] to study the vibrational characteristics of multilayered graphene sheets with different boundary conditions embedded in an elastic medium using finite element method. Sadeghi and naghdabadi [19] introduced a hybrid atomistic-structural element for modeling nonlinear behavior of graphene sheets. Pradhan and Kumar [20] investigated the small scale effect on the vibration analysis of orthotropic single layered graphene sheets embedded in an elastic medium. Both Winklertype and Pasternak-type foundation models were employed to simulate the interaction between the graphene sheet and surrounding elastic medium. Nonlinear free vibration of single-walled carbon nanotubes based on the Timoshenko beam model was studied by Yang et al. [21]. An elastic continuum approach for modeling the nonlinear vibration of double-walled carbon nanotubes under harmonic excitation was investigated by Hawwa and Qahtani [22]. Shen et al. [23] studied nonlinear vibration behavior for a simply supported single layer graphene sheet in thermal environments in order to obtain the nonlocal parameter. Recently, Jomehzadeh and Saidi [24] developed a Navier solution for vibration analysis of nano-plates using the three-dimensional nonlocal elasticity theory. Practically, most of the structures in mechanics such as graphene resonators are linear up to displacements on the order of its thickness [3]. If the oscillations of the graphene resonator result in amplitudes which are not very small, then the linear treatment may be too inaccurate for the purpose in view. In such cases, the accuracy can often be improved sufficiently by carrying out further approximations via geometrically nonlinear model. In this paper, a large amplitude vibration analysis of multilayered graphene sheets is presented based on an orthotropic nano-plate model. By considering the small scale effect in constitutive relations and using the von Karman nonlinear model, the governing equations of motion are obtained in form of coupled nonlinear partial differential equations. The nonlinear natural frequencies are presented for orthotropic nano-plates with three different boundary conditions including all edges simply supported, all edges clamped and two edges simply supported two edges clamped. The effects of nonlocal parameter, number of layer and aspect ratio on the nonlinear vibration behavior of multilayered zigzag and armchair graphene sheets are discussed in details.

r¼

KðjX 0  Xj; sÞr0 ðX 0 Þ dX

0

ð1Þ

V

where r0 is the classical stress tensor and K(jX0  Xj) is the Kernel function represents the nonlocal modulus. While the constitutive equations of classical elasticity is an algebraic relation between stress and strain tensors, that of nonlocal elasticity involves spatial integrals which represent weighted averages of contributions of the strain of all points in the body to the stress at the given point. Eringen [25] showed that it is possible to represent the integral constitutive relation in an equivalent differential form as

ð1  lr2 Þr ¼ r0

ð2Þ

2

where l = (e0 a0) is the nonlocal parameter, a0 an internal charac@2 @2 teristic length and e0 a constant. Also, r2 ¼ @x 2 þ @y2 is the Laplacian operator. 3. Formulation Consider a thin orthotropic nano-plate of total thickness h with dimension a  b for modeling the single layer graphene sheet. The origin of the Cartesian coordinate system is located in the middle of the plate (Fig. 1). Since the graphene sheet is assumed to have large amplitude motion, the von Karman type strain–displacement relations are used as

ex ¼ ex0 þ zjx ey ¼ ey0 þ zjy cxy ¼ cxy0 þ zjxy

ð3Þ

where

8 > <

9

8

  1 @w 2

9

8

9

8

@2 w

9

> > > @x þ 2 @x > ex0 > >  @x2 > > > >j > > = < = >  2 = < x = < @2 w ey0 ¼ @@yv þ 12 @w  j ¼ ; y 2 @y @y > > > > > > :c > : > ; > > > 2 : ; > jxy ; > : @u @ v @w @w > ; xy0 2 @ w þ þ @u

@y

@x

ð4Þ

@x@y

@x @y

Here the time dependent variables u, v and w are the mid-plane displacement components in the x, y and z directions, respectively. As it can be seen, the strain–displacement relations are nonlinear with respect to transverse displacement. According to Hamilton’s principle, the equations of motion of the orthotropic nano-plate can be given by

Z

t2

d

L dt ¼ 0

ð5Þ

t1

where L is the Lagrangian and t is the time variable. Expressing the Lagrangian parameter based on the von Karman theory and considering the small scale effect from Eq. (2) in constitutive relations, the nonlinear equations of motion for a nano-plate can be obtained as follows:

2. Constitutive relations of nonlocal continuum According to nonlocal elasticity theory, the stress at a reference point X is considered to be a function of the strain field at every

(a)

Z

(b)

y

b

x a

Fig. 1. Geometry of a single layer graphene sheet. (a) discrete model (graphene sheet) and (b) continuum model (nano-plate).

1045

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

  € @Nx @Nxy @w €  I1 þ ¼ ð1  lr2 Þ I0 u @x @x @y   € @Nxy @Ny @w 2 þ ¼ ð1  lr Þ I0 v€  I1 @x @y @y

  @ 2 Mx @ 2 Mxy @ 2 My @w @w 2 @ þ 2 þ ð1  l r Þ þ N þ N x xy @x2 @x@y @y2 @x @x @y   @ @w @w N xy þ Ny þ ð1  lr2 ÞP þ ð1  lr2 Þ @y @x @y !!   € @2w € € @ v€ @u @2w 2 € 1 þ þ 2 ¼ ð1  lr Þ I0 wþI  I2 @x2 @y @x @y

Nx

3

2

A11

A12

A13

B11

B12

7 6A 7 6 12 7 6 7 6 A13 7¼6 7 6B 7 6 11 7 6 5 4 B12 Mxy B13

A22

A23

B12

B22

A23

A33

B13

B23

B12

B13

D11

D12

B22 B23

B23 B33

D12 D13

D22 D23

6 Ny 6 6 6 Nxy 6 6M 6 x 6 4 My

ex0 3 7 6 B23 7 76 ey0 7 76 7 B33 76 cxy0 7 76 7 7 6 D13 7 76 jx 7 76 7 D23 54 jy 5 jxy D33 B13

Z

h=2

Substituting these relations into Eq. (6c) and expressing the resultant moments in terms of transverse displacement yields 4

4   @4u @4u ð1Þ ð1Þ ð1Þ @ u þ 2 A12 þ 2A33 þ A22 4 2 2 @x @x @y @y4 !2 @2w @2w @2w ¼  2 @x@y @x @y2

ð1Þ

A11

3

x

v¼

Z

y 0

 2 ! 2 2 @2u 1 @w ð1Þ @ u ð1Þ @ u þ A12 þ A13  dx @x2 @y2 @x@y 2 @x !  2 2 2 2 1 @w ð1Þ @ u ð1Þ @ u ð1Þ @ u dy A12 þ A þ A  22 23 @x2 @y2 @x@y 2 @y

ð1Þ A11

ð14bÞ

4. Modeling of a multilayered graphene sheet (MLGS)

ð9Þ

Let us consider a multilayered graphene sheet with N layers (Fig. 2). The pressure at a point between any two adjacent layers depends on the difference of their deflections at that point. Thus, the van der Waals interaction can be expressed as

Pi ¼

N X

cij ðwi  wj Þ

j¼1

h=2

qð1; z; z2 Þ dz

ð14aÞ

Thus, by determining the transverse displacement and stress function, the in-plane displacements can also be defined.

where s = sin h and c = cos h and h is 0° and 90° for armchair and zigzag graphene sheets, respectively. Also, the inertia parameters (I0, I1, I2) in Eq. (6) can be expressed in term of density of the graphene sheet, q, as follows:

Z

Z 0

E1 =ð2 þ 2m12 Þ

ðI0 ; I1 ; I2 Þ ¼

ð13Þ

where Aij represent the components of the inverse of stretching stiffness matrix. Eqs. (12) and (13) are the starting equations for analyzing the nonlinear vibration analysis of graphene sheets. As it can be seen, these equations are two nonlinear partial differential equations with a total degree of eight. The nonlinear behavior of graphene sheets can be considered as two kinds, one from the inplane forces due to the nature of edge restrains and the other due to the interaction of displacement components affected in the inplane differential equations as well as compatibility conditions. Also, the in-plane displacement components can be expressed in terms of the transverse deflection and stress function by help of Eqs. (4), (7) and (11) by the following relations:

ð8Þ

e 11 c4 2c2 s2 s4 4c2 s2 Q 7 6e 7 6 2 2 4 4 2 2 7 6 Q 12 7 6 c s c þs c s 4c2 s2 7 6 7 6 2 2 4 2 2 7 6 e 7 6 s4 2c s c 4c s Q 7 6 22 7 6 7¼6 3 7 6 3 3 3 2 2 e 7 6 c s cs  c s cs 2csðc  s Þ 7 6Q 7 6 13 7 6 6 e 7 6 cs3 c3 s  cs3 c3 s 2csðc2  s2 Þ 7 5 4 Q 23 5 4 e 33 c 2 s2 ðc2  s2 Þ2 c2 s2 2c2 s2 Q 3 2 E1 =ð1  m12 m21 Þ 6 E =ð1  m m Þ 7 12 21 7 6 2 6 7 4 m21 E1 =ð1  m12 m21 Þ 5

4

Eq. (12) is a nonlinear fourth order partial differential equation in terms of transverse displacement and stress function, and it needs to be augmented with a compatibility equation. Eliminating u and v from Eq. (4) and expressing the in-plane strain components in terms of stress function, one can obtained the compatibility equation as

u¼

2

4

ð11Þ

ð1Þ

ð7Þ

In which the transformed coefficients are expressed in terms of material properties (Young modulus and Poison’s ratio) as [23]

3

Nxy ¼ 

ð12Þ ð6cÞ

h=2

2

Ny ¼

2 2 w @ w € € D11 @@xw4 þ 2ðD12 þ 2D33 Þ @x@2 @y 2 þ D22 @y4 þ ðI 0  I 2 Þð1  lr Þðw þ r wÞ 2 2  2 2 2 2 @ w @ u þ @@yw2 @@xu2 ð1  lr2 ÞP ¼ ð1  lr2 Þ @@xw2 @@yu2  2 @x@y @x@y

32

e ij ð1; z2 Þ dz Q

@2u @x@y

ð6bÞ

where Aij and Dij (i, j = 1, 2, 3) are called the stretching and bending stiffness, respectively and Bij the stretching–bending coupling stiffness which is zero for the symmetric graphene sheet. For a generally single layered anisotropic graphene sheet, the stiffness coefficients are expressed in terms of the reduced transformed coefe ij as ficients Q

ðAij ; Dij Þ ¼

@2u ; @x2

Nx ¼

where a dot denotes differentiation with respect to time and P is the external load acting on the graphene sheet in z direction. The modeling of the graphene sheet is based on the fact that the elastic moduli in two perpendicular orientations are different due to armchair or zigzag configurations, so the nano-plate is considered to be anisotropic. For such a case, the resultant forces (Nx, Ny, Nxy) and resultant moments (Mx, My, Mxy) can be defined in terms of strains as

2

@2u ; @y2

ð6aÞ

ð10Þ

h=2

As it can be concluded, the parameter I1 is identically zero for a symmetric graphene sheet with respect to z-axis. Raju et al. [26] showed that the effect of longitudinal or in-plane inertia on large amplitude vibration of thin-walled structures is negligible. Vanishing the in-plane inertia terms, it can be seen that two first Eq. (6) will be exactly satisfied if a stress function u is defined such as

Fig. 2. A sample of the multilayered graphene sheet.

ð15Þ

1046

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

where Pi is the external load in the form of van der Waals pressure between the ith layer and the jth layer and cij is the coefficient of the van der Waals interaction between the ith layer and the jth layer. This coefficient can be obtained using the van der Waals interaction potential energy, as a function of the inter-atomic spacing using the Lennard–Jones model as

 12  6 r r U LJ ¼ 4e  d d

ð16Þ

where ULJ is the Lennard–Jones potential, r and e are parameters changing for different types of carbon–carbon van der Waals interaction and d is distance between interacting atoms that have the Lennard–Jones interaction. The van der Waals force can be determined by taking the derivative of the potential energy. Differentiating the potential energy (16) with respect to distance d and expressing the Taylor expansion of the results around the equilibrium position, the coefficient of the van der Waals interaction can be defined as [9]

4 @ 4 u1 @4u ð1Þ ð1Þ ð1Þ @ u1 þ 2ðA12 þ 2A33 Þ 2 12 þ A22 4 @x @x @y @y4 !2 @ 2 w1 @ 2 w1 @ 2 w1 ¼  @x@y @x2 @y2

ð1Þ

A11

ð19bÞ

@ 4 w2 @ 4 w2 @ 4 w2 þ 2ðD12 þ 2D33 Þ 2 2 þ D22 4 @x @x @y @y4 2 2 € 2 Þ  ð1  lr2 Þc12 ðw2  w1 Þ €2 þ r w þ ðI0  I2 Þð1  lr Þðw ! @ 2 w2 @ 2 u2 @ 2 w2 @ 2 u2 @ 2 w2 @ 2 u2 2 ð19cÞ 2 ¼ ð1  lr Þ þ @x2 @y2 @x@y @x@y @y2 @x2

D11

4 @ 4 u2 @4u ð1Þ ð1Þ ð1Þ @ u2 þ 2ðA12 þ 2A33 Þ 2 22 þ A22 4 @x @x @y @y4 !2 @ 2 w2 @ 2 w2 @ 2 w2 ¼  @x@y @x2 @y2

ð1Þ

A11

ð19dÞ

( pffiffiffi!2   6  8 5 4 3 24e r 3003p X ð1Þk 5 r 1 cij ¼  12 2 9acc r acc 256 2k þ 1 a k cc hij k¼0 )   2 35p X ð1Þk 2 1 ð17Þ  8 k¼0 2k þ 1 k h6ij

where w1 and w2 are the transverse deflections of bottom and top layers of the graphene sheet, respectively. Two different boundary conditions are considered for the DLGSs as simply supported and clamped which are both movable in-plane edges.

Here acc is the carbon–carbon bond length and hij is the distance between two layers. Let us assume that all layers have the same thickness and effective material constants. Applying Eqs. (12) and (13) to each layer of the multilayered graphene sheet by considering the van der Waals interaction, it can be obtained

Let us consider a double layered graphene sheet with all edges movable simply supported with the following boundary conditions:

D11

@ 4 wi @ 4 wi @ 4 wi þ 2ðD12 þ 2D33 Þ 2 2 þ D22 4 @x @x @y @y4

€ i þ r2 w € i Þ  ð1  lr2 Þ þ I0  I2 Þð1  lr2 Þðw

N X

cij ðwi  wj Þ

j¼1

¼

@ 2 wi @ 2 ui @ 2 wi @ 2 ui @ 2 wi @ 2 ui 1  lr Þ 2 þ @x2 @y2 @x@y @x@y @y2 @x2 2

4   @4u @ 4 ui ð1Þ ð1Þ ð1Þ @ ui i þ 2 A12 þ 2A33 þ A22 @x4 @x2 @y2 @y4 !2 @ 2 wi @ 2 wi @ 2 wi ¼  @x@y @x2 @y2

! ð18aÞ

5.1. A simply supported DLGS (2-SSSS)

2

R b=2

2

R a=2

ui wi ¼ ðM x Þi ¼ @@x@y ¼ ui wi ¼ ðM y Þi ¼ @@x@y ¼

ð1Þ

ð18bÞ

where i = 1, 2, . . . , N. Therefore, the large amplitude vibration of a Nlayered graphene sheet can be described by 2N coupled nonlinear Eq. (18).

dy ¼ 0 at x ¼  2a

@ 2 ui a=2 @x2

dx ¼ 0 at y ¼  2b

ð20Þ

where in this case the index i is 1 or 2. The movable in-plane edge is the edge which is kept straight by a distribution of normal stresses and therefore the resultant stresses on the edge is zero. The nonlinear free vibration response of double layered simply supported graphene sheets can be obtained by introducing the following admissible function for transverse deflection of each layer:

wi ¼ hW i ðtÞ cos

A11

@ 2 ui b=2 @y2

npx a

cos

mpy b

ð21Þ

where n and m are the numbers of half cosine waves in the x and y directions, respectively and Wi(t) is a function of time only. It is obvious that Eq. (21) satisfies the first two boundary conditions in Eq. (20) of each edge. Substituting Eq. (21) into the right side of Eqs. (19b) and (19d), the general solutions for the stress function ui can be obtained as 2

u1 ¼

Eh

2 32n2 m2 a2 b

  4  m4 a4 cosð2npx=aÞ þ n4 b cosð2mpy=bÞ W 1 ðtÞ2 ð22aÞ

5. Free vibration analysis of double layered graphene sheets (DLGSs) Let us consider the large amplitude free vibration analysis of double layered graphene sheets. In this case, the following four nonlinear governing equations should be considered:

@ 4 w1 @ 4 w1 @ 4 w1 þ 2ðD12 þ 2D33 Þ 2 2 þ D22 4 @x @x @y @y4 2 2 €1 þ r w € 1 Þ  ð1  lr2 Þc12 ðw1  w2 Þ þ ðI0  I2 Þð1  lr Þðw ! @ 2 w1 @ 2 u1 @ 2 w1 @ 2 u1 @ 2 w1 @ 2 u1 ð19aÞ  2 ¼ ð1  lr2 Þ þ @x2 @y2 @x@y @x@y @y2 @x2

D11

2

u2 ¼

Eh

2 32n2 m2 a2 b

  4  m4 a4 cosð2npx=aÞ þ n4 b cosð2mpy=bÞ W 2 ðtÞ2 ð22bÞ

These functions exactly satisfy two last conditions of Eq. (20) (in-plane boundary conditions). Introducing the transverse displacement wi and the stress function ui from Eqs. (21) and (22) into the basic Eqs. (19a) and (19c) and then using the Galerkin method, the modal equations of the problem can be obtained. In according to the Galerkin procedure, the integral

Z A

Ci ðw; uÞWi dA

ð23Þ

1047

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

can be computed over the area of the graphene sheet. Ci (w, u) is the nonlinear Eqs. (19a) or (19c) and Wi is the corresponding spatial part of admissible function (21). A simple but lengthy calculation of the above integral leads to modal time differential equations of Duffing’s type, for both layers, which can be written as

5.2. A two edges clamped two edges simply supported DLGS (2-CCSS) It is assumed that the graphene sheet has clamped edges in x direction and simply supported edges in y direction. The boundary conditions for such a nano-plate may be identified as 2

2

d W 1 ðtÞ 2

þ aW 1 ðtÞ þ bW 1 ðtÞ3 þ cW 2 ðtÞ ¼ 0

ð24aÞ

2

þ aW 2 ðtÞ þ bW 2 ðtÞ3 þ cW 1 ðtÞ ¼ 0

ð24bÞ

dt 2 d W 2 ðtÞ dt

where the parameters a, b and c are defined in terms of graphene sheet properties as

b¼

2

2 ð12a2 b

þ

2 a2 h

2

2 2

p2 þ b h p

2 Þða2 b2

2 e 11 Q e2 Þ e 22 þ b4 Q e 11 Þð Q e 22  Q 3p4 h ða4 Q 12 2 2 2 2 e 11 Q e 22 a2 b ð12a2 b þ a2 h p2 þ b h2 p2 Þ 4q Q

þ a2 lp2 þ b

2

12c12 a2 b

q

2 hð12a2 b

þ

2 a2 h

wi ¼ hW i ðtÞ cos2

4

ð25bÞ

For linear vibration of the graphene sheet in which the term bWi(t)3 can be neglected, the corresponding linear natural frepffiffiffiffiffiffiffiffiffiffiffiffi quency is given by xL ¼ a  c. In order to find the nonlinear natural frequencies of DLGSs, the harmonic balance (HB) method is employed. The HB method is an analytical approach for solving nonlinear oscillators, in which the initial conditions are generally simplified by setting velocity or displacement to be zero. Here, the periodic solutions with frequency of x are considered for Eq. (24). Introducing a new variable as s = xt, Eq. (24) can be transformed into

XW 001 þ aW 1 þ bW 31 þ cW 2 ðtÞ ¼ 0

ð26aÞ

XW 002 þ aW 2 þ bW 32 þ cW 2 ðtÞ ¼ 0

ð26bÞ

where X = x2 and the superscript (0 ) denotes differentiation with respect to s. Considering the initial velocity of the graphene sheet equal to zero, the (K + 1)th HB solutions can be described as [27]

W Kþ1 ¼ i

Kþ1 X

Ski cos½ð2k  1Þs

R a=2

dy ¼ 0 at x ¼  2a

@ 2 ui a=2 @x2

ð30Þ

dx ¼ 0 at y ¼  2b

npx mpy cos a b

ð31Þ

3

ð25aÞ

lp2 Þqh

ð25cÞ

p2 þ b2 h2 p2 Þ

@ 2 ui b=2 @y2

Assume that a solution of transverse displacement can be expressed in the following form:

2

2

c¼

2

R b=2

ui wi ¼ ðM y Þi ¼ @@x@y ¼

e 12 a2 b þ 4 Q e 33 a2 b þ Q e 22 a4 þ Q e 11 b Þh 12a2 b ð1 þ p2 la2 þ p2 lb Þc12  p4 ð2 Q 2

a¼

ui i wi ¼ @w ¼ @@x@y ¼ @x

ð27Þ

k¼1

It is easy to show that the admissible function (31) exactly satisfies the two first boundary conditions (30). Upon substituting Eq. (31) into the right side of Eqs. (19b) and (19d), the general solution for the stress function can be obtained. Doing the same procedure as the previous section, a Duffing’s equations with the different coefficients can be obtained. Its coefficients are presented in Appendix A. 5.3. A clamped DLGS (2-CCCC) The boundary conditions for a clamped edges nano-plate with movable in-plane edges may be identified as 2

R b=2

2

R a=2

ui i wi ¼ @w ¼ @@x@y ¼ @x ui i wi ¼ @w ¼ @@x@y ¼ @y

@ 2 ui b=2 @y2 2

@ ui a=2 @x2

dy ¼ 0 at x ¼  2a

ð32Þ

dx ¼ 0 at y ¼  2b

Let us consider the following relation for the transverse displacement:

wi ¼ hW i ðtÞ cos2

npx mpy cos2 a b

ð33Þ

Above admissible function exactly satisfies the two first boundary conditions of clamped edges. Substituting Eq. (33) into the right side of Eqs. (19b) and (19d), the general solution for the stress function is obtained. Doing the same procedure as the previous cases, the Duffing’s equations can be obtained as Eq. (24) with different coefficients.

Substituting the above equation into Eq. (26), expanding results as Fourier series about s and letting the coefficients of {cos s,cos 3s, . . . , cos (2N + 1)s} be zero, it follows that:

6. Free vibration analysis of triple layered graphene sheets (TLGSs)

CSk k ¼ 1; 2; . . . ; K þ 1 i ðSi ; xi Þ ¼ 0

For large amplitude free vibration analysis of graphene sheets with three layers, there are six nonlinear governing equations as follows:

K þ1 X

Si2k1 ¼ Ai

ð28Þ

k¼1

where Ai is the initial displacement of the ith layer. From the above equations, the frequency and parameter S can be obtained for each step. For example, first step of the HB method have the following solution:

pffiffiffiffiffi First step W i ¼ Ai cos xt

3 4

x2 ¼ a þ bA2i  c

@ 4 w1 @ 4 w1 @ 4 w1 þ 2ðD12 þ 2D33 Þ 2 2 þ D22 4 @x @x @y @y4 2 2 €1 þ r w € 1Þ þ ðI0  I2 Þð1  lr Þðw

D11

 ð1  lr2 Þ½c12 ðw1  w2 Þ  c13 ðw1  w3 Þ @ 2 w1 @ 2 u1 @ 2 w1 @ 2 u1 @ 2 w1 @ 2 u1  2 ¼ ð1  lr Þ þ @x2 @y2 @x@y @x@y @y2 @x2 2

ð29Þ

! ð34aÞ

1048

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

4   @4u @ 4 u1 ð1Þ ð1Þ ð1Þ @ u1 1 þ 2 A12 þ 2A33 þ A22 4 2 2 @x @x @y @y4 !2 @ 2 w1 @ 2 w1 @ 2 w1 ¼  @x@y @x2 @y2

2

d W 3 ðtÞ

ð1Þ

A11

2

2

D11

þ a3 W 3 ðtÞ þ b2 W 3 ðtÞ3 þ c3 W 1 ðtÞ þ k3 W 2 ðtÞ ¼ 0

2

2

4

3

2

2

2

ð12a2 b þ a2 h

2

4

3

e 12 a2 b þ 4 Q e 33 a2 b þ Q e 22 a4 þ Q e 11 b Þh 12a2 b ð1 þ p2 la2 þ p2 lb Þðc23 þ c13 Þ  p4 ð2 Q 2

2

2

2

ð12a2 b þ a2 h

2

4

3

@ 2 w2 @ 2 u2 @ 2 w2 @ 2 u2 @ 2 w2 @ 2 u2 ¼ ð1  lr Þ  2 þ @x2 @y2 @x@y @x@y @y2 @x2

b1 ¼ b2 ¼ b3 ¼

!

2

ð34cÞ

4   @4u @ 4 u2 ð1Þ ð1Þ ð1Þ @ u2 2 þ 2 A þ 2A þ A 12 33 22 @x4 @x2 @y2 @y4 !2 2 2 2 @ w2 @ w2 @ w2 ¼  @x@y @x2 @y2

2 e 11 Q e2 Þ e 22 þ b4 Q e 11 Þð Q e 22  Q 3p4 h ða4 Q 12 2 2 2 2 e 11 Q e 22 a2 b ð12a2 b þ a2 h p2 þ b h2 p2 Þ 4q Q

c1 ¼ c2 ¼

12c12 a2 b 2 hð12a2 b

q

ð1Þ

c3 ¼ k1 ¼

@ 2 w3 @ 2 u3 @ 2 w3 @ 2 u3 @ 2 w3 @ 2 u3 2 þ 2 2 @x @y @x@y @x@y @y2 @x2

k2 ¼ k3 ¼

!

4   @4u @ 4 u3 ð1Þ ð1Þ ð1Þ @ u3 3 þ 2 A12 þ 2A33 þ A22 4 2 2 @x @x @y @y4 !2 @ 2 w3 @ 2 w3 @ 2 w3 ¼  @x@y @x2 @y2

ð34eÞ

ð34fÞ

Three types of boundary conditions are considered for the triple layered graphene sheets as all edges simply supported, two edges simply supported two edges clamped and all edges clamped which are all movable in-plane edges. Similar to the double layered graphene sheets, the admissible functions for transverse deflections of each layer are assumed as Eqs. (21), (31) and (33) with i = 1, 2, 3 for TLGSs. Doing the same procedure as the previous case, the three Duffing equations can be obtained for this case. For example, the modal equations of triple layered simply supported graphene sheets are 2

dt 2 d W 2 ðtÞ 2

þ a1 W 1 ðtÞ þ b1 W 1 ðtÞ3 þ c1 W 2 ðtÞ þ k1 W 3 ðtÞ ¼ 0 3

ð36eÞ

p2 þ b2 h2 p2 Þ

12c13 a2 b

2

ð36fÞ

qhð12a2 b2 þ a2 h2 p2 þ b2 h2 p2 Þ

þ a2 W 2 ðtÞ þ b2 W 2 ðtÞ þ c2 W 1 ðtÞ þ k2 W 3 ðtÞ ¼ 0

12c23 a2 b 2 hð12a2 b

q

þ

2 a2 h

2

ð36gÞ

p2 þ b2 h2 p2 Þ

Considering the harmonic vibration, the mode shapes and frequencies of this case can also be obtained using HB method.

7. Numerical results

ð1Þ

A11

d W 1 ðtÞ

þ

2 a2 h

ð34dÞ

@ 4 w3 @ 4 w3 @ 4 w3 þ 2ðD þ 2D Þ þ D D11 12 33 22 @x4 @x2 @y2 @y4 2 2 € € þ ðI0  I2 Þð1  lr Þðw3 þ r w3 Þ  ð1  lr2 Þ½c13 ðw3  w1 Þ  c23 ðw3  w2 Þ

ð36dÞ

2

A11

¼ ð1  lr2 Þ

ð36cÞ

p2 þ b2 h2 p2 Þða2 b2 þ a2 lp2 þ b2 lp2 Þqh

 ð1  lr2 Þ½c12 ðw2  w1 Þ  c23 ðw2  w3 Þ

dt

ð36bÞ

p2 þ b2 h2 p2 Þða2 b2 þ a2 lp2 þ b2 lp2 Þqh

@ 4 w2 @ 4 w2 @ 4 w2 þ 2ðD12 þ 2D33 Þ 2 2 þ D22 4 @x @x @y @y4 2 2 €2 þ r w € 2Þ þ ðI0  I2 Þð1  lr Þðw

2

ð36aÞ

p2 þ b2 h2 p2 Þða2 b2 þ a2 lp2 þ b2 lp2 Þqh

2

ð35cÞ

where the coefficients ai, bi and ci are defined in terms of graphene sheet properties as

e 12 a2 b þ 4 Q e 33 a2 b þ Q e 22 a4 þ Q e 11 b Þh 12a2 b ð1 þ p2 la2 þ p2 lb Þðc12 þ c23 Þ  p4 ð2 Q

2

a3 ¼ 

2

ð12a2 b þ a2 h 2

a2 ¼ 

ð34bÞ

2

e 12 a2 b þ 4 Q e 33 a2 b þ Q e 22 a4 þ Q e 11 b Þh 12a2 b ð1 þ p2 la2 þ p2 lb Þðc12 þ c13 Þ  p4 ð2 Q 2

a1 ¼ 

dt

In order to verify the accuracy of the present formulations, the linear natural frequencies are compared with the available results in literature for a simply supported classical double and triple isotropic nano-plates in which the nonlocal effect is neglected (l = 0). In Table 1, linear natural frequencies are presented and compared with those reported in Ref. [9]. It can be observed that the results of the present study are sufficiently accurate. For numerical results, six different graphene sheets with the following properties are considered [23]: Table 1 Comparison of the linear natural frequencies of a simply supported nano-plate. Number of layers

Linear natural frequency (THz)

Ref. [9]

Present

Ref. [9]

Present

2

xL1 xL2

0.069 2.683

0.0667 2.6807

0.173 2.688

0.1659 2.6812

3

xL1 xL2 xL3

0.069 1.865 3.286

0.0667 1.8632 3.2828

0.173 1.872 3.29

0.1659 1.8667 3.2817

ð35aÞ ð35bÞ

m = 1, n = 1

m = 1, n = 2

1049

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

Armchair1: a ¼ 9:519 nm; b ¼ 4:844 nm; h ¼ 0:129 nm; q ¼ 6316 kg=m3 ;

Also, the van der Waals coefficients are considered as follows [9]:

e0 a0 ¼ 0:67 nm

m12 ¼ 0:197; E1 ¼ 2:434 TPa; E2 ¼ 2:473 TPa;G12 ¼ 1:039 TPa ð37aÞ Armchair2: a ¼ 6:995 nm; b ¼ 4:847 nm; h ¼ 0:143 nm; q ¼ 5727 kg=m3 ; e0 a0 ¼ 0:47 nm

m12 ¼ 0:202; E1 ¼ 2:154 TPa; E2 ¼ 2:168 TPa; G12 ¼ 0:923 TPa ð37bÞ Armchair3: a ¼ 4:888 nm; b ¼ 4:855 nm; h ¼ 0:156 nm; q ¼ 5295 kg=m3 ; e0 a0 ¼ 0:27 nm

m12 ¼ 0:201; E1 ¼ 1:949 TPa; E2 ¼ 1:962TPa; G12 ¼ 0:846 TPa

c12 ¼ 108:651 GPa=nm; c13 ¼ 1:87206 GPa=nm

ð38Þ

The negative sign of van der Waals coefficient represents an attraction between two layers and the positive sign represents repulsion. The fundamental linear natural frequencies for different types of graphene sheets with three different boundary conditions are presented in Table 2. The results are obtained for single, double and triple layered graphene sheets. Based on the results in this table, it can be concluded that the first natural frequency of the graphene sheets is independent of the number of layers. The results indicate that the increase of the number of the layers results in a higher linear natural frequencies as all the layers of the graphene sheet vibrate as a whole. It can also be seen that the clamped graphene sheets have the higher linear

ð37cÞ Zigzag1:

3 3

a ¼ 9:496 nm; b ¼ 4:877 nm; h ¼ 0:145 nm; q ¼ 5624 kg=m ;

2.75

e0 a0 ¼ 0:47 nm

2.5

m12 ¼ 0:223; E1 ¼ 2:145 TPa; E2 ¼ 2:097 TPa; G12 ¼ 0:938 TPa ð37dÞ

a ¼ 7:065 nm; b ¼ 4:877 nm; h ¼ 0:149 nm; q ¼ 5482 kg=m3 ;

2.25

ω / ωL

Zigzag2:

zigzag1-CCSS zigzag1-SSSS armchair1-CCSS armchair1-SSSS zigzag1-CCCC armchair1-CCCC

2 1.75

e0 a0 ¼ 0:32 nm

m12 ¼ 0:204; E1 ¼ 2:067 TPa; E2 ¼ 2:054 TPa; G12 ¼ 0:913 TPa

1.5

ð37eÞ Zigzag3:

1.25 3

a ¼ 1:987 nm; b ¼ 1:974 nm; h ¼ 0:154 nm; q ¼ 5363 kg=m ;

1

e0 a0 ¼ 0:22 nm

0

1

m12 ¼ 0:205; E1 ¼ 1:987 TPa; E2 ¼ 1:974 TPa; G12 ¼ 0:857 TPa

2

Wmax

ð37fÞ

Fig. 3. The backbone curves for single layer graphene sheets.

Table 2 Linear natural frequencies (THz) of graphene sheets.

1.012

Number of layer

CCSS

CCCC

2

3

Armchair1 Armchair2 Armchair3

0.0569 0.0763 0.1143

0.0569 0.0763 0.1143

2.5988 2.5923 2.5818

0.0569 0.0763 0.1143

1.8062 1.8020 1.7957

3.1827 3.1744 3.1610

Zigzag1 Zigzag2 Zigzag3

0.0453 0.0426 0.0491

0.0453 0.0426 0.0491

2.5972 2.5948 2.5799

0.0453 0.0426 0.0491

1.8049 1.8032 1.7929

3.1808 3.1778 3.1596

Armchair1 Armchair2 Armchair3

0.0642 0.0958 0.1695

0.0642 0.0958 0.1695

2.5999 2.5927 2.5843

0.0642 0.0958 0.1695

1.8063 1.8028 1.7998

3.1827 3.1747 3.1629

Zigzag1 Zigzag2 Zigzag3

0.0505 0.0641 0.1225

0.0505 0.0641 0.1225

2.5972 2.5950 2.5820

0.0505 0.0641 0.1225

1.8050 1.8037 1.7961

3.1807 3.1779 3.1610

Armchair1 Armchair2 Armchair3

0.1126 0.1450 0.2127

0.1126 0.1450 0.2127

2.5999 2.5944 2.5871

0.1126 0.1450 0.2127

1.8083 1.8057 1.8040

3.1833 3.1758 3.1649

Zigzag1 Zigzag2 Zigzag3

0.1166 0.1258 0.1634

0.1166 0.1258 0.1634

2.5988 2.5968 2.5838

0.1166 0.1258 0.1634

1.8076 1.8066 1.7991

3.1818 3.1792 3.1624

1.01

3-CCCC 3-CCSS 2-CCCC 3-SSSS 2-CCSS 2-SSSS

1.008

ω / ωL

SSSS

1

1.006

1.004

1.002

1

0

1

2

3

4

Wmax Fig. 4. The backbone curves for double and triple layer armchair graphene sheets.

1050

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

1.012

2 1.9

1.01

3-CCCC 3-CCSS 2-CCCC 3-SSSS 2-CCSS 2-SSSS

c12 =0 c12 =10 16

1.7

c12 =10 17

1.6

ω / ωL

ω / ωL

1.008

1.8

1.006

1.5 1.4

1.004

1.3 1.2

1.002

1.1 1

0

1

2

3

4

1

0

Wmax

2.6

-9 2

μ=(0.75×10 ) -9 2

μ=(0.5×10 )

-9 2

μ=(0.25×10 ) μ=0

ω / ωL

2 1.8 1.6 1.4 1.2 1

0

1

2

3

4

Fig. 7. The backbone curves of a simply supported double layer graphene sheets for some van der Waals coefficients.

nonlinear ratio increases. To study the effect of small length scale on nonlinear behavior of graphene sheets, the backbone curves of a simply supported single layer armchair graphene sheet (type 1) are depicted in Fig. 6 for various values of nonlocal parameters. It can be seen that the nonlocal parameter has a considerable effect on the nonlinear vibration behavior of SLGSs. At a given vibration amplitude, the nonlinear frequency ratio increases with the increase of the nonlocal parameter. During calculations, it has been found that the nonlocal parameter does not have a significant effect on nonlinear frequency ratio of multilayered graphene sheets. In order to study the effect of van der Waals interaction on nonlinear behavior of graphene sheets, the backbone curves of a simply supported double layer graphene sheet are shown in Fig. 7 for some values of van der Waals coefficient. It can be seen that the van der Waals interaction decreases the nonlinear behavior of the graphene sheets.

-9 2

μ=(1×10 )

2.2

2

Wmax

Fig. 5. The backbone curves for double and triple layer zigzag graphene sheets.

2.4

1

3

4

Wmax Fig. 6. The backbone curves of a simply supported single layer graphene sheets for different nonlocal parameters.

natural frequencies than the other boundary conditions in all cases since the boundary support is the stronger for the clamped graphene sheets. To study the effects of physical properties of graphene sheets, the variation of nonlinear frequency ratio (nonlinear to linear frequencies) versus maximum amplitude of the bottom surface is investigated. Backbone curves (nonlinear to linear frequency– amplitude curve) are plotted for single, double and triple layered graphene sheets in Figs. 3–5 with all edges clamped, all edges simply supported and two simply supported two clamped edges boundary conditions. It can be found that the geometry nonlinearity makes graphene sheets have hardening behavior. i.e., the nonlinear frequency ratio increases as the vibration amplitude increases. Also, it can be seen from Fig. 3 that zigzag graphene sheets have higher values of nonlinear to linear frequencies than armchair graphene sheets. In addition, it can be concluded from Figs. 4 and 5 that as the number of layer in MLGSs increases, the

8. Conclusion The nonlinear free vibration analysis of multilayered graphene sheets has been studied based on the nonlocal elasticity theory. Considering the Eringen nonlocal theory and von Karman hypothesis, the nonlinear equations of motion have been obtained using the Hamilton’s principle. For illustrative examples, the large amplitude vibration solutions of single, double and triple layered graphene sheets with simply supported or clamped boundary conditions have been found using the harmonic balance method. The results in the form of linear and nonlinear frequencies have been presented for different types of armchair and zigzag graphene sheets. It has been seen that the small length scale increases the nonlinear behavior of the graphene sheets but the van der Waals interaction decreases it. Also, it can be seen the van der Waals interaction decreases the nonlinear behavior of DLGSs but increases the nonlinear behavior of TLGSs.

Acknowledgment The support for this work, provided by the Iran Nanotechnology Initiative Council, is gratefully acknowledged.

E. Jomehzadeh, A.R. Saidi / Computational Materials Science 50 (2011) 1043–1051

e 12 a2 b þ 16 Q e 33 a2 b þ 3 Q e 22 a4 þ 16 Q e 11 b Þh a2 b ð36a2 b þ 36p2 la2 þ 48p2 lb Þc12  p4 ð8 Q 2

a¼

2

2

4 b4

½ðð3p4 a4 þ 16p

þ 8p

4 a2 b2 Þh2

þ 36p

2 a4 b2

þ 48p

2

2

2 a2 b4 Þ

2 a4 b2

l þ ð3p

4

2 a2 b4 Þh2

þ 4p

þ

1051

3

4 36a4 b 

qh

ðA:1Þ

2 e 2 4 6 2 2 8 4 8 e e e2 b ¼ f3p4 h ð Q 11 Q 22  Q 12 Þ½ðð16p a b l þ 64p a b l þ 16a b Þ Q 22 þ 10 8 10 e 2 8 2 2 4 4 8 4 e2 e2 ð64p2 b l þ 16p2 a2 b l þ 16a2 b Þ Q 66 Þ Q 11 þ ðð17p a b l þ 36p a b l þ 17a b Þ Q 22 þ 4 6 6 6 8 e 22  ð16p2 a4 b l þ 64p2 a2 b l þ 16a4 b8 Þ Q e 66 Q e2  ð41p2 a6 b l þ 116p2 a4 b l þ 41a6 b Þ Q 12

l þ 128p2 a2 b8 l þ 32a4 b8 Þ Qe 66 Qe 12 Þ Qe 11 þ ð17p2 a10 l þ 36p2 a8 b2 l þ 17a10 b2 Þ Qe 66 Qe 222  2 4 4 e2 2 8 2 2 6 4 8 4 e e e ðð17p2 a8 b l þ 36p2 a6 b l þ 17a8 b Þ Q 12  ð34p a b l þ 72p a b l þ 34a b Þ Q 66 Q 12 Þ Q 22 g= 4 2 2 2 2 4 2 2e 2 2 2 4 2 4 4 e 22 þ b Q e 11 Q e 66  a b Q e  2a b Q e 12 Q e 66 þ a Q e 66 Q e 22 Þð3p a b h þ 36a b þ fða b Q 11 Q 12 4 2 4 4 2 2 4 2 4 2 2 2 2 2 4 2 2 2 4 3p a h l þ 36p a b l þ 8p a b h l þ 4p a b h þ 48p a b l þ 16p4 b h lÞg

ð32p2 a4 b

6

2

c¼

2

12c12 a2 b ð3a2 b þ 3p2 a2 l þ 43p2 b ½ðð3p

4 a4

þ 16p

4 b4

4 a2 b2 Þh2

þ 8p

þ 36p

2 a4 b2

þ 48p

2 a2 b4 Þ

2

l þ ð3p

Appendix A The Duffing’s coefficients for a double layered CCSS graphene sheet References [1] A.K. Geim, Science 324 (2009) 1530–1534. [2] A.K. Geim, P. Kim, Scientific American 298 (2008) 68–75. [3] J.S. Bunch, A.M. van der Zande, S.S. Verbridge, I.W. Frank, D.M. Tanenbaum, J.M. Parpia, H.G. Craighead, P.L. McEuen, Science 315 (2007) 490–493. [4] J.T. Robinson, M. Zalalutdinov, J. W Baldwin, E.S. Snow, Z. Wei, P. Sheehan, B.R. Houston, Nano Letters 8 (2008) 3441–3445. [5] I.W. Frank, D.M. Tanenbaum, A.M. van der Zande, P.L. McEuen, Journal of Vacuum Science and Technology B 25 (2007) 2558–2561. [6] C.W. Lima, L.H. He, International Journal of Mechanical Sciences 46 (2004) 1715–1726. [7] K. Behfar, R. Naghdabadi, Composites Science and Technology 65 (2005) 1159– 1164. [8] S. Kitipornchai, X.Q. He, K.M. Liew, Physical Review B 72 (2005) 075443. 6pp.. [9] X.Q. He, S. Kitipornchai, K.M. Liew, Nanotechnology 16 (2005) 2086–2091. [10] Y.M. Fu, J.W. Hong, X.Q. Wang, Journal of Sound and Vibration 296 (2006) 746– 756.

lÞ

2 a4 b2

4

2

4

þ 4p2 a2 b Þh þ 36a4 b qh

ðA:2Þ

ðA:3Þ

[11] K.M. Liew, X.Q. He, S. Kitipornchai, Acta Materialia 54 (2006) 4229–4236. [12] S.C. Pradhan, J.K. Phadikar, Journal of Sound and Vibration 325 (2009) 206– 223. [13] L.L. Ke, Y. Xiang, J. Yang, S. Kitipornchai, Computational Materials Science 47 (2009) 409–417. [14] Y.X. Dong, C.W. LIM, Science in China Series E: Technological Sciences 52 (2009) 617–621. [15] T. Murmu, S.C. Pradhan, Computational Materials Science 46 (2009) 854–859. [16] W.H. Duan, C.M. Wang, Nanotechnology 20 (2009) 075702. 7pp. [17] S.C. Pradhan, J.K. Phadikar, Physics Letters A 373 (2009) 1062–1069. [18] R. Ansari, R. Rajabiehfard, B. Arash, Computational Materials Science 49 (2010) 831–838. [19] M. Sadeghi, R. Naghdabadi, Nanotechnology 21 (2010) 105705. 10pp. [20] S.C. Pradhan, A. Kumar, Computational Materials Science 50 (2010) 239–245. [21] J. Yang, L.L. Ke, S. Kitipornchai, Physica E 42 (2010) 1727–1735. [22] M.A. Hawwa, H.M. Al-Qahtani, Computational Materials Science 48 (2010) 140–143. [23] L. Shen, H.S. Shen, C.L. Zhang, Computational Materials Science 48 (2010) 680– 685. [24] E. Jomehzadeh, A.R. Saidi, Composite Structures 93 (2011) 1015–1020. [25] A.C. Eringen, Journal of Applied Physics 54 (1983) 4703–4710. [26] I.S. Raju, G.V. Rao, K. Raju, Journal of Sound and Vibration 49 (1976) 415–422. [27] S.K. Lai, C.W. Lim, B.S. Wub, C. Wang, Q.C. Zeng, X.F. He, Applied Mathematical Modelling 33 (2009) 852–866.