Jtam-A4.dvi JOURNAL OF THEORETICAL AND APPLIED MECHANICS 55, 3, pp. 823-838, Warsaw 2017 DOI: 10.15632/jtam-pl.55.3.823 LINEAR FREE VIBRATION OF GRAPHENE SHEETS WITH NANOPORE VIA AIFANTIS THEORY AND RITZ METHOD Sima Ziaee Yasouj University, Mechanical Engineering Department, Yasouj, Iran e-mail: ziaee@yu.ac.ir This article aims to study the natural frequency of defective graphene sheets since the existence of cut-outs in plates may be essential on the basis of their desired functionality. A combination of the Aifantis theory and Kirchhoff thin plate hypothesis is used to derive governing equations ofmotion. TheRitzmethod is employed to derive discrete equations of motion.Themolecular structuralmechanicsmethod is also employed to specify the effective length scale parameter. In the ‘numerical results’ Section, the effects of different parameters suchasboundaryconditionsanddiameter of thehole-to-side length ratio on the fundamental frequency of graphene sheets are studied. Keywords: free vibration, defective graphene sheet, Aifantis theory, molecular structural mechanics 1. Introduction Graphene sheets, the two-dimensional carbon nanostructure, amazingly exhibit the electrical, mechanical, thermal and optical properties (Jomehzadeh et al., 2015) that have been considered as a promisingmaterial for a wide range of applications (Lebedeva et al., 2012) such as compo- sites, chemical sensors, ultra capacitors, transparent electrodes, photovoltaic cells, bio-devices (Zandiatashbar et al., 2014) and the gigahertz oscillator suggested based on the telescopic oscil- lation of graphene layers (Lebedeva et al., 2012). It is of great importance to simulatemechanical behavior of single/muti-layered graphene sheets accurately due to this wide scope of applica- tions. One important topic which has recently received great attention is the study of free and forced vibration of graphene sheets used in nano-scale devices which may experience vibration. Among different theoretical modelings used to simulate carbon nano-structures, continuum mechanics modeling andmolecular structural mechanics (MSM)methods are the most popular ones employed to study themechanical behavior of carbon nanostructures because their compu- tational tasks are much more time-effective than different atomistic modeling such as classical molecular dynamics, tight-binding molecular dynamics and density functional theory (Li and Chou, 2003; Hu et al., 2007). The linear vibration of zigzag and armchair single-layered gra- phene sheets (SLGS) was studied by Sakhaee-Pour et al. (2008) via the MSM approach. They obtained fundamental frequencies and correspondingmode shapes for different boundary condi- tions and showed that the natural frequencies of SLGS are independent of chirality and aspect ratio (Sakhaee-pour et al., 2008), while the Hashemnia et al. study (2009) indicated that the fundamental frequency of SLGS with fixed-fixed or fixed-free end conditions decreased as the aspect ratio increased. They proposed the use of graphene sheets with a lower aspect ratio to prevent resonance and dynamic damage (Hashemnia et al., 2009). Gupta and Batra (2010) im- plemented MSM simulation to investigate the effect of pretention on the natural frequencies of SLGS. The comparison between the results obtained by the MSMmethod and those predicted by an equivalent linear elastic isotropic continuum model showed that there was a noticeable 824 S. Ziaee difference between themode shapes corresponding to the several lowest frequencies of the SLGS and those of an equivalent linear elastic isotropic continuum model (Gupta and Batra, 2010). Sadeghi and Naghdabadi (2010) introduced a hybrid atomistic structural element to model the nonlinear behavior of SLGS.Their results verifiedwith experimental observations demonstrated that natural frequencies predicted by nonlinear analysis were far higher than those obtained by the linear approach. Wang et al. (2013) proposed a pseudo beam model with 3-node beam ele- ment to simulate themodal behavior of thewrinkledSLGSbyusing theMSMmethod.Theynot only obtained the features of formation and evolution of wrinkles but also offered a prediction model to estimate the natural frequency of wrinkle SLGS (Wang et al., 2013). The importance of incorporating the size effect into continuum mechanics, in order to inve- stigate the mechanical behavior of micro- or nano-scale devices, is well known and higher-order continuum theories containing additional material constants have been developed to this end (Askes and Aifantis, 2011). Mindlin (1964) developed a theory of elasticity withmicrostructure in which strain energy was considered as a function of macroscopic strain, the difference betwe- en macroscopic and microscopic deformation and the gradient of the microscopic deformation. In this theory, the strain energy contained 16 constitutive coefficients in addition to Lame’s constants.Mindlinwrote the kinetic energy density in terms of quantities atmicroscale andma- croscale as well (Askes and Aifantis, 2011). Mindlin also simplified this theory and introduced three new versions which differed in the assumed relation between themicroscopic deformation gradient and themacroscopic displacement (Mindlin, 1964; Askes andAifantis, 2011). TheMin- dlin simplified theory in which strain energy is only a function of first-order gradient of strain tensor contains five new constants as well as Lame’s constants for an isotropic linear elastic ma- terial (Mindlin, 1964; Askes and Aifantis, 2011). It can be shown that these five new constants can be grouped in two new constants (Askes andAifantis, 2011). In this way, the number of new material constants is reduced from 5 down to 2 (Askes and Aifantis, 2011). A simplified kinetic energy whose constitutive equation consists of velocity and velocity gradient was also sugge- sted by Mindlin (1964). Therefore, Mindlin’s simplified theory contains one additional inertia parameter as well as two additional elastic parameters (Askes and Aifantis, 2011). Lam et al. (2003) proposed a modified strain gradient theory in which three new length scale parameters were introduced into isotropic linear elastic materials. In this theory, the total deformation energy density is independent of the anti-symmetric rotation gradient tensor and it is only a function of the symmetric strain tensor, the dilatation gradient vector, the deviatoric stretch gradient tensor and the symmetric rotation gradient tensor (Lam et al., 2003). The other well-known higher-order continuum theories are classical couple stress theory (Mindlin and Tiersten, 1962) with two material length scale parameters for an isotropic elastic material, and themodified couple stress theory (Yang et al., 2002) whose modified constitutive equation contains only a new length scale parameter. It can be shown that the modified couple stress theory is a special case of themodified strain gradient elasticity theory if two of the three material length scale parameters of themodified strain gradient theory are taken to be equal to zero (Ashoori Movassagh andMahmoodi, 2013). On the basis of Aifantis’ studies in plasticity and nonlinear elasticity, Aifantis and his cowor- kers proposed another gradient elasticity theory, the constitutive equation of which is a function of strain and Laplacian of strain with one internal length (Aifantis, 1992; Altan and Aifantis, 1997). It can be shown that the associated equilibrium equations obtained by simplifiedMidlin’s theory can be converted to those obtained by the Aifantis gradient elasticity theory if two new elastic parameters of simplifiedMindlin’s theory are taken to be equal to each other (Askes and Aifantis, 2011). Another powerful gradient elasticity theory proposed to simulate dynamical behavior of micro/nanostructures is achievedbycombining stable straingradientswithacceleration (inertia) gradients (Askes and Aifantis, 2011) such as Mindlin’s simplify theory. Metrikine and Askes Linear free vibration of graphene sheets... 825 (2002, 2006) derived such a theory from a discrete lattice. It is worth mentioning that these theories knownas dynamically consistentmodels (Askes andAifantis, 2011), incorporate at least two length scales: one of them related to strain gradients and another related to acceleration (inertia) gradients. The nonlocal strain gradient theory is another Laplacian based gradient elasticity whose constitutive equation is a function of strain, stress, Laplacian of strain and Laplacian of stress (Askes and Aifantis, 2011; Aifantis, 2011). Meany researchers employed different gradient elasticity theories mentioned above in con- junction with different classical continuummechanic theories to study the mechanical behavior of size-dependentmicro/nanostructures so far (Ansari et al., 2011, 2013; Reddy, 2011; Rahmani and Pedram, 2014; Akgöz and Civalek, 2012, 2014, 2015; Gholami et al., 2016; Binglei et al. 2011, 2016; Askes and Aifantis, 2009; Şimşek, 2016; Ebrahimi et al., 2016; Li and Hu, 2016). On the basis of the strain gradient Timoshenko beam theory, the free vibration characteri- stics of functionally graded microbeams were investigated by Ansari et al. (2011). Ansari et al. (2013) also combined the most general strain gradient elasticity theory containing five additio- nal material length scale parameters with the classical Timoshenko beam theory to investigate bending and buckling of functionally gradedmicrobeams. Akgöz andCivalek (2014) proposed a new size-dependent trigonometric beammodel based on the strain gradient theory.They (Akgöz andCivalek, 2014) employed their newmodel to estimate load-bearing capacity of microbeams. A new size-dependent sinusoidal platemodel to predictmechanical behavior of thin,moderately thick and thick microplate was suggested by Akgöz and Civalek (2015) as well. They employed a modified strain gradient elasticity hypothesis to incorporate the size effect into the classical continuum plate theory. In an earlier work, Akgöz and Civalek (2012) estimated vibratory be- havior of single-layered graphene sheets embedded in an elastic matrix. To this end, they used the Krichhoff plate hypothesis in conjunction with themodified couple stress theory to develop governing equations of motion (Akgöz and Civalek, 2012). Gholami et al. (2016) combinedMindlin’s simplified theory hypothesis with first-order shear deformation shell theory to investigate free vibration and axial buckling of circular cylindrical micro-/nano-shells. They ignored the effects of inertia gradients in their model and derived the kinetic energy in terms ofmacroscopic kinematic quantities. They compared three different size- dependent shell models on the basis of strain gradient theory, modified strain gradient theory, andmodified couple stress theory in predicting the natural frequency and load-bearing capacity ofmicroshells aswell. They showed that the values of natural frequency and the critical buckling force predicted by the strain gradient theory are higher than those predicted by the modified strain gradient theory because of inducing a higher stiffness. Binglei et al. (2011) simulated static bending, static instability and vibratory behavior of a simply-supported micro-plate via the modified strain gradient theory and the modified couple stress theory. Their results clearly show that the strain gradient theory induces a higher stiffness than the modified couple stress theory.They (Binglei et al., 2016) re-derived the governing equation of the nonclassicalKirchhoff micro-plate as well as the general boundary conditions based on the strain gradient elasticity to simulate mechanical behavior of a micro-plate with any reasonable boundary conditions. Askes and Aifantis (2009) employed a dynamically consistent model with two length scales to simulate flexural wave dispersion in carbon nanotubes. They verified validity of the proposed modelbycomparing theirfindingswith thoseobtainedviamolecular dynamics simulationand/or Eringen’s nonlocal elasticity theory. Using nonlocal strain gradient elasticity, Li and Hu (2016) developed a model to analyze wave propagation in fluid-conveying carbon nanotubes. They investigated the effects of different parameters such as nonlocal parameter, small scale parameter, damping coefficient and flow velocity on wave propagation properties of fluid-conveying carbon nanotubes. The otherwell-known gradient theory isEringen’s nonlocal elasticity theory inwhich there is 826 S. Ziaee only one length scale parameter and its constitutive equation is expressed based on the nonlocal stress tensor and Laplacian stress tensor (Askes and Aifantis, 2011). Eringen’s nonlocal elasti- city is widely employed to simulate mechanical behavior of carbon nano-structures because of being a good agreement between obtained results via Eringen’s theory andmolecular dynamics simulation. Murmu and Pradhana (2009) employed a nonlocal elasticity theory to investigate the vibrational behavior of SLGS. Numerical results clearly revealed that the small scale value could affect the fundamental frequencies of SLGS significantly (Murmu and Pradhana, 2009). Shen et al. (2010) employed a nonlocal orthotropic plate model which contained small scale effects to simulate nonlinear vibration behavior of SLGS in thermal environmental. They used temperature-dependentmaterial properties obtained bymolecular dynamic simulation. To esti- mate the value of the small-scale parameter, they reconciled the natural frequencies of graphene sheets obtained via molecular dynamics (MD) simulation with those predicted by the nonlocal plate model (Shen et al., 2010). Their results clearly showed the significant impact of vibration amplitude and temperature change on the nonlinear vibration response of both armchair and zigzag graphene sheets (Shen et al., 2010). Ansari et al. (2010) employed the nonclassical Min- dlin plate theory to study the vibrational behavior of SLGS. To incorporate the size effect into the continuum platemodel, they used Eringen’s nonlocal elasticity theory.Matching the results obtained from MD simulation with those taken from the nonlocal model, they determined the proper value of the nonlocal parameter. They showed that the nonlocal parameter depends on chairality and boundary conditions. Using nonlocal elasticity, Farajpour et al. (2011) estimated the load-bearing capacity of circular graphene sheet subjected to uniform radial compression. To study the vibratory behavior of annular and circular graphene sheet in thermal environment, Mohammadi et al. (2014) employed Eringen’s theory. They investigated the effect of different parameters such as elastic medium, boundary conditions and small scale parameter on the na- tural frequencies of circular and annular graphene sheets. Farajpour et al. (2012) estimated load-bearing capacity of a rectangular single-layered graphene sheet by combining the nonlocal elasticity theory with Kirchhoff’s plate hypothesis as well as orthotropic material stress-strain relations. Although there have been numerous studies carried out on the simulation of vibrational behavior of perfect graphene sheets, the notable studies shown that the effects of nanopores on linear and/or nonlinear vibration responses of graphene sheets do not exist, while the existence of nanopores in the graphene lattice can be essential based upon the desired functionality of graphene sheets in different nano-devices such as bio-devices and the DNA-decorated graphene (Zandiatashbar et al., 2014). On the other hand, the production process used or environmental and operating conditions under which the graphene device operates can lead to the appearance of vacancy defects in the graphene lattice (Zandiatashbar et al., 2014). Therefore, simulation of the dynamic response of graphene sheets with nanopores to determine how the size, the position and the shape of cut-outs influence the natural frequencies is important. Accordingly, the main purpose of this study is to model the dynamic response of SLGS with cut-outs. Different me- thods, such as Ritz, finite difference and finite element can be employed to this end. Rajamani and Prabhakaran (1977) used Lagrange’s equations to obtain discrete equations of motion of a composite plate with a central cut-out. For this purpose, they employed potential energy and ki- netic energy of a uniformplate and considered the cut-out as a displacement-dependent external loading on the plate (Rajamani and Prabhakaran, 1977). The finite differencemethodwas used byAksu and Ali (1976) to determine vibration characteristics of rectangular plates with one or two central cut-outs. To develop a simple method to study the effects of the rectangular cut- out on natural frequencies of rectangular plates, Ali and Atwal (1980) used Rayleigh’s method and employed trigonometric functions to estimate the lateral deflection. They showed that this method can predict the fundamental frequency with a quite good accuracy although the accu- racy could be increased if one added a correction function to the estimated lateral displacement Linear free vibration of graphene sheets... 827 function (Ali andAtwal, 1980). Lam et al. (1989) presented amodification of the Rayleigh-Ritz method to investigate the vibrational behavior of rectangular plates with one or two cut-outs. They used an orthogonal polynomial function generated based on theGram-Schmidt process to estimate the lateral displacement of plate segments (Lam et al., 1989). Similar methods were used by Liew et al. (2003) to investigate the effect of different boundary conditions on vibratio- nal behavior of a rectangular plate with the central rectangular cut-out. The negative stiffness method was applied to analyze free vibration and buckling of plates with cut-outs by Tham et al. (1986). One of the popular variational methods employed to derive discrete equations of motion is the Ritz method because the displacement field is approximated by a linear combination of shape functions which are only satisfied by essential boundary conditions. On the other hand, it is shown that one can reduce the complexity of the mechanical behavior modelling of structures with cut-outs by combining the construction of energy functional with the Ritz method (Malekzadeh et al., 2013).Therefore, in this study, theRitzmethod is employed toderive discrete equations of motion which yield the eigenvalue problem. Because Eringen’s nonlocal elasticity theorydoes not allow the construction of energy functional (Reddy, 2011), the dynamic consistentmodel is used to simulate vibrational behavior of the graphene sheetwith a nanopore. To this end, Aifantis’ theory, as well as inertia gradients is combined withKirchhoff’s thin plate hypothesis to construct the energy functional and incorporate the size effect into the classical continuum plate model. Also, the length scales related to strain gradients and inertia gradients are set equal to each other. Themolecular structuralmechanics method is also used to simulate vibrational behavior of defective SLGS to estimate the length scale parameter which must be used in the proposed nonclassical plate model. 2. Governing equation 2.1. Nonclassical thin plate theory According toHamilton’s principle, one can formulate discretized free vibration equations of a nano-plate with free edges cut-out via theRitzmethod.TheAifantis theory in conjunctionwith inertia gradients is also used to incorporate the size effect into the equation governing natural frequencies of the nano-plate with the cut-out. For this reason, the variation of strain energy and kinetic energy can be written as δPE = ∫ VI (σijδε L ij) dV − ∫ VII (σijδε L ij) dV δKE = ∫ VI [ρu̇δu̇+ρv̇δv̇+ρẇδẇ+ρl21(u̇,iδu̇,i+ v̇,iδv̇,i+ ẇ,iδẇ,i)] dV − ∫ VII [ρu̇δu̇+ρv̇δv̇+ρẇδẇ+ρl21(u̇,iδu̇,i+ v̇,iδv̇,i+ ẇ,iδẇ,i)] dV (2.1) where σij, ε L ij and l1 are components of the stress tensor, linear strain tensor and length scale related to inertia gradients, respectively. Vi (i = I,II) is the volume of i-th segment of the plate. The displacement field is shown by u, v andw. ρ is the specific mass density of the plate material and the time derivative of displacements are shown by overhead dot. According to the classical plate theory hypothesis, if the xy-plane of the Cartesian coordinate system (x,y,z) coincides with the geometrical mid-plane of the undeformed rectangular micro-/nano-plate, the displacement field can be expressed as (Reddy, 1999) 828 S. Ziaee u(x,y,z,t) =u0(x,y,t)−z ∂w(x,y,t) ∂x v(x,y,z,t) = v0(x,y,t)−z ∂w(x,y,t) ∂y w(x,y,z,t) =w(x,y,t) (2.2) where u0 and v0 are the x- and y-components of the displacement vector of a point in the mid-plane of the plate at time t, respectively. According to the classical plate theory hypothesis, the linear components of the strain tensor can be found on the basis of Eqs. (2.2) (Reddy, 1999) εLxx = ∂u0 ∂x −z ∂2w ∂x2 εLyy = ∂v0 ∂y −z ∂2w ∂y2 εLxy = ε L yx = 1 2 (∂v0 ∂x + ∂u0 ∂y −2z ∂2w ∂x∂y ) (2.3) According to Aifantis’ theory, the stress tensor components are related to the linear strain tensor components as (more details can be found in Aifantis (1992, 2011), Askes and Aifantis (2011), Gitman et al. (2005)) σij =Cijkl(εkl− l 2εkl,mm) (2.4) which can be obtained on the basis of the implicit gradient elasticity model given by tr(α1ε+α2σ)1+α3ε+α4σ+∇ 2[tr(α5ε+α6σ)1+α7ε+α8σ] = 0 through a proper choice of the constants α1, . . . ,α8 (Aifantis, 2011).C can be defined based on the classical plate theory hypothesis as C=       E 1−ν2 Eν 1−ν2 0 Eν 1−ν2 E 1−ν2 0 0 0 G       (2.5) whereE,G and ν are themodulus of elasticity, modulus of rigidity and Poisson’s ratio, respec- tively. l is the length scale parameter related to strain gradients. By substitutingEqs. (2.2) to (2.5) into Eqs. (2.1), one can find the variation of strain energy andkinetic energyon thebasis of thedisplacementfield.Byusing theappropriate approximation of dependent unknowns (u0,v0,w) satisfying the essential boundary conditions and employing Hamilton’s principle, one can find the following system of ordinary differential equations gover- ning the lateral vibration of the plate. It is completely independent of the other two equations governing in-planemotion of the plate (MImnpq −M II mnpq)Ẅpq+(K I mnpq−K II mnpq)Wpq =0 (2.6) in which MI,IImnpq = ∫ VI,VII [ −1 12 ρh3 (∂φp ∂x ∂φm ∂x ϕqϕn+ ∂ϕq ∂y ∂ϕn ∂y φpφm ) −ρhφpϕqφmϕn ] dxdy − l21 ∫ VI,VII 2ρh (∂φp ∂x ∂φm ∂x ϕqϕn+ ∂ϕq ∂y ∂ϕn ∂y φpφm ) dxdy − l21 ∫ VI,VII 1 12 ρh3 (∂2φp ∂x2 ∂2φm ∂x2 ϕqϕn+ ∂2ϕq ∂y2 ∂2ϕn ∂y2 φpφm+2 ∂φp ∂x ∂φm ∂x ∂ϕq ∂y ∂ϕn ∂y ) dxdy (2.7) Linear free vibration of graphene sheets... 829 and KI,IImnpq = ∫ VI,VII [ −D (∂2φp ∂x2 ϕq +ν ∂2ϕq ∂y2 φp )∂2φm ∂x2 ϕn −D (∂2ϕq ∂y2 φp+ν ∂2φp ∂x2 ϕq )∂2ϕn ∂y2 φm− [2D(1−ν)] ∂2(φpϕq) ∂x∂y ∂2(φmϕn) ∂x∂y ] dxdy − l2 ∫ VI,VII [ −D (∂4φp ∂x4 ϕq +(1+ν) ∂2ϕq ∂y2 ∂2φp ∂x2 +ν ∂4ϕq ∂y4 φp )∂2φm ∂x2 ϕn −D (∂4ϕq ∂y4 φp+(1+ν) ∂2φp ∂x2 ∂2ϕq ∂y2 +ν ∂4φp ∂x4 ϕq )∂2ϕn ∂y2 φm −[2D(1−ν)] (∂4(φpϕq) ∂x3∂y + ∂4(φpϕq) ∂x∂y3 )∂2(φmϕn) ∂x∂y ] dxdy (2.8) whereD=Eh3/[12(1−ν2)] and h is thickness of the plate. TheGram-Schmidt process is used to generate the polynomial functions ϕn(y) (ϕq(y)) and φm(x) (φp(x)) as well. Equation (2.6) can be solved on the basis of the eigen-value problem to find the natural frequencies corresponding to lateral vibration of the nano-plate. 2.2. Molecular structural mechanics Themolecular structuralmechanics approach has been built on the basis of some similarities between themolecular model of carbon nanostructures (nanotube and graphene sheet) and the structure of a space frame building (Li and Chou, 2003; Hu et al., 2007; Sakhaee-pour et al., 2008; Hashemnia et al., 2009; Gupta and Batra, 2010, Sadeghi and Naghdabadi, 2010; Wang et al., 2013). In this method, the atomic lattice of carbon nanostructures is rebuilt by beams and point masses replacing C-C covalent bonds and carbon atoms, respectively (Li and Chou, 2003; Hu et al., 2007; Sakhaee-pour et al., 2008; Hashemnia et al., 2009; Gupta andBatra, 2010; Sadeghi and Naghdabadi, 2010; Wang et al., 2013). To determine stiffness parameters of the equivalent beam, i.e., tensile stiffness (EA), bendingstiffness (EI) and torsional rigidity (GJ), one can use the concepts of energy equivalence between the total potential energies in computational chemistry (Eq. (2.9)) and the total elemental strain energies in structuralmechanics (Eq. (2.10)) (Li and Chou, 2003; Hu et al., 2007; Sakhaee-pour et al., 2008; Hashemnia et al., 2009; Gupta and Batra, 2010; Sadeghi and Naghdabadi, 2010; Wang et al., 2013). The total potential function among carbon atoms due to bonded andnonbonded interactions are (Li andChou, 2003;Hu et al., 2007; Sakhaee-pour et al., 2008;Hashemnia et al., 2009;Gupta and Batra, 2010; Sadeghi and Naghdabadi, 2010; Wang et al., 2013) UP = ∑ Ur + ∑ Uθ+ ∑ (Uϕ+Uω)+ ∑ Uvdw (2.9) in which Ur is the bond stretching, Uθ is the bond angle bending, Uϕ is the dihedral angle torsion, Uω is the out of plane torsion, Uvdw is the nonbonded van der Waals interaction and Ue is the nonbonded electrostatic interactions. On the other hand, the total strain energy for a beam element is (Li and Chou, 2003; Hu et al., 2007; Sakhaee-pour et al., 2008; Hashemnia et al., 2009; Gupta and Batra, 2010; Sadeghi and Naghdabadi, 2010; Wang et al., 2013) USE = ∑ UA+ ∑ UM + ∑ UT + ∑ UV (2.10) whereUA is the strain energyof axial tension,UM is the strain energyof bending,UT is the strain energy of torsion and Uv is the strain energy of the shear force. One can obtain a relationship between themolecular mechanics force field constants and the structural mechanics parameters 830 S. Ziaee by using the energy equivalence between Eq. (2.9) andEq. (2.10) (Li andChou, 2003; Hu et al., 2007; Sakhaee-pour et al., 2008; Hashemnia et al., 2009; Gupta and Batra, 2010; Sadeghi and Naghdabadi, 2010; Wang et al., 2013) as follows EA L = kr EI L = kθ GJ L = kτ (2.11) where kr, kθ and kτ are force field constants in molecular mechanics. On the basis of the previous research (Li and Chou, 2003; Hu et al., 2007; Sakhaee-pour et al., 2008; Hashemnia et al., 2009; Gupta andBatra, 2010; Sadeghi andNaghdabadi, 2010;Wang et al., 2013), the following numerical values of the force field constants used in this research are selected kr =65.2nNÅ −1 kθ =2.78nNÅrad −2 kτ =2.78nNÅrad −2 (2.12) Then, the diameter-to-length ratio of the beam element estimated by Eq. (2.11) shows that this element is too thick andmust bemodeled based on shear deformation beam theories. Wang et al. (2013) used the study of Scarpa et al. (2009) who considered the C-C bond as a Timoshenko beam, to obtain the mechanical properties and optimized cross-section diameter of the equivalent beam element via the following relationships (Wang et al., 2013) d< 2 √ 6 kθ kr E = k2rL 4πkθ G= kτ 2kθ (2.13) where kθ is expressed as kθ = krd 2 16 4C1+C2 C1+C2 (2.14) and (Wang et al., 2013) C1 =112L 2kτ +192L 2kτν+64L 2kτν 2 C2 =9krd 4+18krd 4ν+9krd 4v2 (2.15) and L is length of the C-C bond. They obtained the diameter-to-length ratio of the equivalent Timoshenko beam around 0.704 (Wang et al., 2013), which is used in this study. Based on the finite element method and using modal analysis concepts, one can find the natural frequencies of this space frame-like structure constructed byTimoshenko beam elements and point masses. 3. Verification 3.1. Non-classical thin plate theory Based on the best knowledge of the author, it is the first attempt to analyze the vibrational behavior of graphene sheets with nanopores via Aifantis’ theory combined with acceleration (inertia) gradients and the Ritzmethod. Then, to verify the correctness of the presentedmodel, the author first compares the obtained natural frequencies of classical thin plates with central holes with those obtained based on the finite element method (Table 1). To this end, the value of length scale parameters (l and l1) is taken equal to zero. According toTable 1, althoughnatural frequencies of perfect plates canbe accurately predic- ted by the presentmodel, those of the defective plate are overestimated. However, the calculated percentage error (lower than 4% and 3% for the simply supported and clamped plate, respecti- vely) shows that there is an acceptable agreement between the results. Linear free vibration of graphene sheets... 831 Table 1. Comparisons of the first dimensionless natural frequencies Ω = (ωa2/π2) √ ρh/D for the fully simply supported (SSSS)and fully clamped (CCCC)plate (a/b=1,a/h=29.4, d– the diameter of the central hole) d/a SSSS CCCC Present FEM Present FEM 0 1.9986 1.9986 3.6445 3.6445 0.1 2.0021 1.9913 3.6907 3.6533 0.2 2.0304 1.9950 3.8482 3.7644 0.3 2.0918 2.0454 4.1673 4.0898 0.4 2.1979 2.1727 4.7690 4.7286 0.5 2.4988 2.4106 5.8615 5.8571 Next, after setting the length scale parameters to be equal to each other (i.e. l = l1), the effects of an increase in the length scale parameter on the variation of the first natural frequency of perfect graphene sheets by increasing the side-length of nano-plate is compared with the data available in Ansari et al. (2010), showing the influence of the small scale parameter on the natural frequency of perfect graphene sheets by using Eringen’s nonlocal elasticity theory. Figure 1 shows that, regardless of the boundary conditions, the first natural frequency of perfect graphene sheets decreases with an increase in the length scale parameter, although a rise in the side-length of the nano-plate reduces the importance of the length scale value, so that the non-classical natural frequencies tend to the classical ones. The observed behavior is completely in agreement with that reported in Ansari et al. (2010). Fig. 1. The effects of the length scale value on the natural frequency of perfect graphene sheets; (a) fully simply supported nano-plates, (b) fully clamped nano-plates (l= l1) 3.2. Molecular structural mechanics method To verify the accuracy of predicted natural frequencies of graphene sheets obtained via molecular structural mechanics (MSM), the estimated first natural frequency of fully simply supported and fully clamped single-layered graphene sheets are compared with those obtained by the molecular dynamics method and are available in Ansari et al. (2010), see Table 2. As it is seen, there is a quite good agreement between the results. 4. Numerical results In this Section, first, the effects of different parameters such as length scale and diameter of the central hole to the side length ratio, which is not more than the half of side length of nano- plate, on the dimensionless natural frequency of graphene nano-sheet with a central pore are 832 S. Ziaee Table 2.Comparison of the first natural frequency of a square single-layered graphene sheet Side-length SSSS CCCC a [nm] MSM MD MSM MD 10 0.0576347 0.0587725 0.1148990 0.1146223 15 0.0292556 0.0273881 0.0523147 0.0517078 20 0.0167249 0.0157524 0.0310192 0.0306219 25 0.0118420 0.0099840 0.0182260 0.0179975 investigated (seeFig. 2).Then, thevariation of dimensionless natural frequencies of the graphene nano-sheet with the pore diameter to side length ratio equal to 0.5 versus the length scale and side length is studied. Finally, it is tried to estimate the effective length scale of the defective graphene sheet viamatching the results obtained from the non-classical plate theory with those found on the basis of molecular structural mechanics. It should be mentioned that the length scales related to strain gradients and inertia gradients are set equal to each other (i.e. l= l1) in this Section. Fig. 2. Defective monolayer graphene sheet modeled with a plate with a circular hole InFig. 3, the effects of the length scale anddiameter of the central hole to the side length ratio on the dimensionless natural frequency of the graphene nano-sheet with a central pore for two different boundary conditions are investigated. Regardless of boundary conditions, for a small value of the length scale parameter, the first natural frequency of defective graphene sheets rises as diameter of the central hole to the side length ratio increases, while with a rise in the length scale, the rate of increase of the first natural frequencywithporediameter decreases. For a larger value of the length scale, a decrease in the first natural frequency may be seen by increasing the pore diameter (Figs. 3a and 3d). Generally, however, the sensitivity of simply supported graphene sheets to the existence of the central hole is less than in the fully clamped nano-sheets. The sensitivity of the second and third natural frequency of a fully clamped graphene sheet to an increase in pore diameter is less than the first one, although the length scale value can change the increasing or decreasing trend of the curves. The impact of side length on frequencies can be seen in Fig. 4 for fully clamped and simply supported graphene sheets. It is assumed that graphene has a central hole the diameter of which is half of the side length of the sheet. It is clearly seen that the importance of the length scale role in decreasing the natural frequencies reduceswith an increase in the side length of the defective sheet. It can be concluded that the influence of the length scale on higher order frequencies is more than the first natural frequency as well. As expected, the natural frequencies of simply-supported sheets are less than the fully clamped ones. In order to show if the existence of vacancy defects makes a difference to the effective length scale used in the non-classical plate theory, the results obtained from the molecular structural Linear free vibration of graphene sheets... 833 Fig. 3. The effects of the length scale and diameter of the hole to side length of plate ratio on the first three dimensionless natural frequencies of fully clamped graphene sheets (a), (b) and (c) and fully simply-supported graphene sheets (d), (e) and (f) (l= l1) mechanics are used. First, the effective length scale for perfect square graphene nano-sheets the side length of which is 10nm is obtained by matching the first natural frequencies predicted by the molecular structural mechanics with the non-classical plate model. They are around 1.15nmand 0.5nm for simply supported and fully clamped nano-sheets, respectively. According to these effective length scales, the first three natural frequencies of nano-sheets are 0.0578181, 0.133181 and 0.193658THz, 0.1150115, 0.2230178 and 0.3152479THz for simply supported and fully clamped nano-plates, respectively. The comparison of these results with those determined bymolecular structuralmechanics (0.05763, 0.12857 and0.18252THz for simply supported sheet and 0.11490, 0.23047 and 0.34411 for clamped sheet) shows that the percentage error in higher order frequencies is more than that of the first natural frequencies. Similar results can be seen in Table 3 and 4 in which the samemethod is used to estimate the effective length scale of the defective sheet. It can be concluded that the effective length scale may vary with an increase 834 S. Ziaee Fig. 4. The effect of the length scale and side length on the first three dimensionless natural frequencies of the fully clamped graphene sheet (a), (b) and (c) and fully simply-supported graphene sheet (d), (e) and (f) with a central hole the diameter to side length ratio of which is 0.5 (l= l1) in diameter of the vacancy defect as well. It seems that the dependency of the effective length scale on diameter of the vacancy defect changes with boundary conditions. The effective length scale of a clamped defective graphene is more sensitive to diameter of the vacancy than that of the simply supported one. 5. Conclusion This article attempts to study the natural frequency of defective graphene sheets because the existence of cut-outs in platesmay be essential on the basis of their desired functionality. In this study, Aifantis’ theory in conjunction with inertia gradients is combined with Kirchhoff’s thin plate hypothesis to incorporate the size effect into the classical continuum plate theory because Eringen’s nonlocal elasticity theory does not allow the construction of an energy functional. Linear free vibration of graphene sheets... 835 Table 3. The estimated length-scale parameter of simply supported defective graphene sheets on the basis of matching natural frequencies obtained from the present method against MSM results No. of Mode No. MSM [THz] Present Length scale lost method parameter atoms [THz] l= l1 [nm] ω11 0.057446 0.05738 6 ω12 0.12840 0.14056 1.2 ω22 0.18224 0.19349 ω11 0.057170 0.05711 24 ω12 0.12797 0.13944 1.145 ω22 0.18157 0.19140 ω11 0.057008 0.05706 54 ω12 0.12718 0.13934 1.155 ω22 0.18049 0.19096 Table 4. The estimated length-scale parameter of fully clamped defective graphene sheets on thebasis ofmatchingnatural frequencies obtained fromthepresentmethodagainstMSMresults No. of Mode No. MSM [THz] Present Length scale lost method parameter atoms [THz] l= l1 [nm] ω11 0.11444 0.11439 6 ω12 0.23712 0.21931 0.55 ω13 0.34343 0.30796 ω11 0.11381 0.11383 24 ω12 0.23650 0.21525 0.61 ω13 0.34196 0.29900 ω11 0.11356 0.11355 54 ω12 0.23650 0.21236 0.71 ω13 0.34196 0.29157 The Ritz method is employed to derive discrete equations of motion which yield the eigenvalue problem. The molecular structural mechanics method is also employed to specify the effective length scale parameter. In the ‘numerical results’ Section, the effects of different boundary conditions, length scale, diameter of hole to side length ratio and side length of nano-sheets on the fundamental frequency of graphene sheets are studied. The results demonstrate that: • Regardless of boundary conditions, for a small value of the length scale parameter, the first natural frequency of defective graphene sheets rises as diameter of the central hole to the side length ratio increases, while with a rise in the length scale, the rate of increase of the first natural frequency with pore diameter decreases. • For larger values of the length scale, by increasing pore diameter, a decrease in the first natural frequencymay be seen. • The sensitivity of simply supported graphene sheets to the existence of the central hole is less than that of the fully clamped nano-sheets. • The importance of the length scale role in decreasing the natural frequencies reduces with an increase in the side length of the defective sheet. 836 S. Ziaee • The influence of the length scale on higher order frequencies is stronger than that on the first natural frequencies. • Although the dependency of the effective length scale on diameter of the vacancy defect changeswith boundary conditions, the effective length scale does not varywith an increase in diameter of the vacancy defect significantly. Acknowledgements Theauthorgratefullyacknowledgesthe supportofYasoujUniversityunderGrantNo.Gryu-89111109. References 1. Aifantis E.C., 1992, On the role of gradients in the localization of deformation and fracture, International Journal of Engineering Science, 30, 1279-1299 2. AifantisE.C., 2011,On thegradientapproach-relationtoEringen’snonlocal theory, International Journal of Engineering Science, 49, 1367-1377 3. Akgöz B., Civalek Ö., 2012, Free vibration analysis for single-layered graphene sheets in an elastic matrix via modified couple stress theory,Materials and Design, 42, 167-171 4. Akgöz B., Civalek Ö., 2014, A new trigonometric beam model for buckling of strain gradient microbeams, International Journal of Mechanical Sciences, 81, 88-94 5. Akgöz B., Civalek Ö., 2015, A microstructure-dependent sinusoidal plate model based on the strain gradient elasticity theory,Acta Mechanica, 226, 2277-2294 6. Aksu G., Ali R., 1976, Determination of dynamic characteristics of rectangular plates with cu- touts using a finite difference formulation, Journal of Sound and Vibration, 44, 147-158 7. Ali R., Atwal S.J., 1980, Prediction of natural frequencies of vibration of rectangular plates with rectangular cutouts,Composite Structures, 12, 819-823 8. Altan B., Aifantis E.C., 1997, On some aspects in the special theory of gradient elasticity, Journal of Mechanical Behavior of Material, 8, 231-282 9. Ansari R., Gholami R., Faghih Shojaei M., Mohammadi V., Sahmani S., 2013, Size- -dependent bending, buckling and free vibration of functionally graded Timoshenko microbeams based on the most general strain gradient theory,Composite Structures, 100, 385-397 10. AnsariR.,GholamiR., SahmaniS., 2011,Freevibrationanalysis of size-dependent functionally gradedmicrobeams based on the strain gradient Timoshenko beam theory,Composite Structures, 94, 221-228 11. AnsariR., Sahmani S.,ArashB., 2010,Nonlocal platemodel for free vibrationof single-layered graphene sheets,Physics Letters A, 375, 53-62 12. Ashoori Movassagh A., Mahmoodi A.M.J., 2013, A micro-scale modeling of Kirchhoff plate based on modified strain-gradient elasticity theory,European Journal of Mechanics A/Solids, 40, 50-59 13. Askes H., Aifantis E.C., 2009,Gradient elasticity and flexural wave dispersion in carbon nano- tubes,Physical Review, B80, 195412 14. Askes H., Aifantis E.C., 2011, Gradient elasticity in statics and dynamics: An overview of formulations, length scale identificationprocedures,finite element implementationsandnewresults, International Journal of Solids and Structures, 48, 1962-1990 15. Binglei W., Shenghua H., Junfeng Z., Shenjie Z., 2016, Reconsiderations on boundary conditions of Kirchhoff micro-plate model based on a strain gradient elasticity theory, Applied Mathematical Modelling, 40, 7303-7317 Linear free vibration of graphene sheets... 837 16. Binglei W., Shenjie Z., Junfeng Z., Xi C., 2011, A size-dependent Kirchhoff micro-plate model based on strain gradient elasticity theory, European Journal of Mechanics A/Solids, 30, 517-524 17. Ebrahimi F., Barati M.R., Dabbagh A., 2016, A nonlocal strain gradient theory for wave propagation analysis in temperature-dependent inhomogeneous nanoplates, International Journal of Engineering Science, 107, 169-182 18. Farajpour A., Mohammadi M., Shahidi A.R., Mahzoon M., 2011, Axisymmetric buckling of the circular graphene sheets with the nonlocal continuumplatemodel,Physica E, 43, 1820-1825 19. Farajpour A., Shahidi A.R., Mohammadi M., Mahzoon M., 2012, Buckling of orthotropic micro/nanoscale plates under linearly varying in-plane load via nonlocal continuum mechanics, Composite Structures, 94, 1605-1615 20. Gholami R., Darvizeh A., Ansari R., Sadeghi F., 2016, Vibration and buckling of first- -order shear deformable circular cylindrical micro-/nano-shells based on Mindlin’s strain gradient elasticity theory,European Journal of Mechanics A/Solids, 58, 76-88 21. Gitman I.M., Askes H., Aifantis E.C., 2005, The representative volume size in static and dynamic micro-macro transitions, International Journal of Fracture, 135, L3-L9 22. Gupta S.S., Batra R.C., 2010, Elastic properties and frequencies of free vibrations of single- -layered graphene sheets, Journal of Computational and Theoretical Nanoscience, 10, 2151-2164 23. Hashemnia K., Farid M., Vatankhah R., 2009, Vibrational analysis of carbon nanotubes and graphene sheets usingmolecular structuralmechanics approach,Computational Materials Science, 47, 79-85 24. Hu N., Nunoya K., Pan D., Okabe T., Fukunaga H., 2007, Prediction of buckling characte- ristics of carbon nanotubes, International Journal of Solids and Structures, 44, 6535-6550 25. Jomehzadeh E., Saidi A.R., Jomehzadeh Z., Bonaccorso F., Palermo V., Galiotis C., Pugno N.M., 2015,Nonlinear subharmonic oscillation of orthotropic graphene-matrix composite, Computational Materials Science, 99, 164-172 26. Lam D.C.C., Yang F., Chong A.C.M., Wang J., Tong P., 2003, Experiments and theory in strain gradient elasticity, Journal of the Mechanics and Physics of Solids, 51, 1477-1508 27. Lam K.Y., Hung K.C., Chow S.T., 1989, Vibration analysis of plates with cutouts by the modified Rayleigh-Ritz method,Applied Acoustics, 28, 49-60 28. Lebedeva I.V., Knizhnik A.A., Popov A.M., Lozovik Yu. E., Potapkin B.V., 2012,Mo- deling of graphene-basedNEMS,Physica E, 44, 949-954 29. Li C., ChouT.-W., 2003,A structuralmechanics approach for the analysis of carbon nanotubes, International Journal of Solids and Structures, 40, 2487-2499 30. Li L., Hu Y., 2016,Wave propagation in fluid-conveying viscoelastic carbon nanotubes based on nonlocal strain gradient theory,Computational Materials Science, 112, 282-288 31. Liew K.M., Kitipornchai S., Leung A.Y.T., Lim C.W., 2003, Analysis of the free vibration of rectangular plates with central cut-outs using the discrete Ritz method, International Journal of Mechanical Sciences, 45, 941-959 32. Malekzadeh P., Bahranifard F., Ziaee S., 2013, Three-dimensional free vibration analysis of functionally graded cylindrical panels with cut-out using Chebyshev-Ritz method, Composite Structures, 105, 1-13 33. Metrikine A.V., Askes H., 2002, One-dimensional dynamically consistent gradient elasticity models derived from a discrete microstructure. Part 1: Generic formulation, European Journal of Mechanics A/Solids, 21, 555-572 34. Metrikine A.V., Askes H., 2006, An isotropic dynamically consistent gradient elasticitymodel derived from a 2D lattice,Philosophical Magazine, 86, 3259-3286 838 S. Ziaee 35. Mindlin R., 1964, Micro-structure in linear elasticity, Archive for Rational Mechanics and Ana- lysis, 16, 52-78 36. Mindlin R.D., Tiersten H.F., 1962, Effects of couple-stresses in linear elasticity, Archive for Rational Mechanics and Analysis, 11, 415-448 37. Mohammadi M., Farajpour A., Goodarzi, Dinari F., 2014, Thermo-mechanical vibration analysis of annular and circular graphene sheet embedded in an elastic medium, Latin American Journal of Solids and Structures, 11, 659-682 38. Murmu T., Pradhana S.C., 2009, Vibration analysis of nano-single-layered graphene sheets embedded in elastic medium based on nonlocal elasticity theory, Journal of Applied Physics, 105, 064319 39. Rahmani O., Pedram O., 2014, Analysis and modeling the size effect on vibration of functio- nally graded nano-beams based on nonlocal Timoshenko beam theory, International Journal of Engineering Science, 77, 55-70 40. Rajamani A., Prabhakaran R., 1977, Dynamic response of composite plates with cut-outs. Part 1: Simply-supported plates, Journal of Sound and Vibration, 54, 549-564 41. Reddy J.N., 1999,Theory and Analysis of Elastic Plates, Taylor and Francis, USA 42. Reddy J.N., 2011,Microstructure-dependent couple stress theories of functionally graded beams, Journal of Mechanics and Physics of Solids, 59, 2382-2399 43. Sadeghi M., Naghdabadi R., 2010, Nonlinear vibrational analysis of single-layered graphene sheets,Nanotechnology, 21, 105705-105710 44. Sakhaee-pour A., Ahmadian M.T., Naghdabadi R., 2008, Vibrational analysis of single- layered graphene sheets,Nanotechnology, 19, 85702-85707 45. ScarpaF., Adhikari S., PhaniA.S., 2009,Effective elasticmechanical properties of single layer graphene sheets,Nanothecnology, 20, 065709 46. Shen L., Shen H.-S., Zhang C.-L., 2010, Nonlocal platemodel for nonlinear vibration of single layer graphene sheets in thermal environments,Computational Materials Science, 48, 680-685 47. ŞimşekM., 2016,Nonlinear free vibration of a functionally gradednanobeamusingnonlocal strain gradient theory and a novel Hamiltonian approach, International Journal of Engineering Science, 105, 12-27 48. Tham L.G., Chan A.H.C., Cheung Y.K., 1986, Free vibration and buckling analysis of plates by negative stiffness method,Composite Structures, 22, 687-692 49. WangC.G., Lan L., LiuY.P., TanH.F., HeX.D., 2013,Vibration characteristics of wrinkled single-layered graphene sheets, International Journal of Solids and Structures, 50, 1812-1823 50. YangF., ChongA.C., LamD.C.C., TongP., 2002,Couple stress based strain gradient theory for elasticity, International Journal of Solids and Structurs, 39, 2731-2743 51. Zandiatashbar A., Lee G.-H., An S.J., Lee S., Mathew N., Terrones M., Hayashi T., PicuC.R., Hone J., KoratkarN., 2014,Effect of defects on the intrinsic strength and stiffness of graphene,Nature Communications, 5, 3186, DOI: 10.1038/ncomms4186 Manuscript received May 15, 2016; accepted for print February 6, 2017