Jtam.dvi JOURNAL OF THEORETICAL AND APPLIED MECHANICS 49, 3, pp. 621-639, Warsaw 2011 VARIATIONAL PRINCIPLES AND NATURAL BOUNDARY CONDITIONS FOR MULTILAYERED ORTHOTROPIC GRAPHENE SHEETS UNDERGOING VIBRATIONS AND BASED ON NONLOCAL ELASTIC THEORY Sarp Adali School of Mechanical Engineering, University of KwaZulu-Natal, Durban, South Africa e-mail: adali@ukzn.ac.za Variational principles are derived for multilayered orthotropic graphene sheetsundergoing transversevibrationsbasedonthenonlocal elastic the- ory of orthotropic plates which provide a continuummodel for graphene sheets. The variational formulation allows the derivation of natural bo- undary conditions which are expressed in the form of a set of coupled equations formultilayeredsheets as opposed touncoupledboundarycon- ditions applicable to simply supported and clamped boundaries and also in the case of a formulation based on the local (classical) elasticity the- ory. For the free vibrations case, the Rayleigh quotient is derived. The methods for the variational formulation use techniques of calculus of va- riations and the semi-inverse method for deriving variational integrals. Variational formulations provide the basis for a number of approximate and numerical methods of solutions and improve the understanding of the physical phenomena. Key words: variational formulation, multilayered graphene sheets, non- local theory, vibrations, semi-inversemethod 1. Introduction Graphene is a two-dimensional carbon nanostructure with many applications in several fields. The covalent bond of carbon atoms makes a graphene sheet one of the stiffest and strongestmaterials with aYoung’smodulus in the range of 1TPaandhigher as supportedby the results givenbyPoot andvanderZant (2008), Sakhaee-Pour (2009), Gao and Hao (2009) and Shokrieh and Rafiee (2010) on themechanical properties of graphene sheets. Its superior properties have already been put into use in a number of applicationswhich include their 622 S. Adali use as sensing devices (Arsat et al., 2009; Wu et al., 2010), in lithium-ion batteries (Lian, 2010), for desalination of seawater (Mishra andRamaprabhu, 2011), in electrochemical capacitors (Yuan et al., 2011) for electrooxidation (Choiet al., 2011) aswell as for sensors for thedetection of cancer cells (Yang et al., 2010;Feng et al., 2011).Theyarealsousedas reinforcements in composites, and a review of the graphene based polymer composites is given by Kuilla et al. (2010). Further applications and potential applications of graphene in the information technology and in other fields are discussed in the review articles by Soldano et al. (2010) and Terrones et al. (2010). Experimental study of nano-scale structures has been a difficult field due to the size of the phenomenon. Similarly, molecular dynamics approach has its drawbacks in the form of extensive computer time and memory required to investigate even relatively small nano structures using nano time scales. This situation led to the development of continuum models for nano-sized components, e.g., carbon nanotubes, and in particular, graphene sheets to investigate their mechanical behaviour (He et al., 2004; Kitipornchai et al., 2005; Hemmasizadeh et al., 2008), and these models were used extensively to investigate the mechanical behaviour of graphene sheets. However, the nano- scale thickness of the sheets leads to inaccurate results when the models are based on classical elastic constitutive relations. Classical elasticity is a scale free theory and as such neglects the size effects which become prominent at atomistic scale. Size effects have been observed in experimental andmolecular dynamic simulations of carbon nanotubes due to the influences of interatomic and intermolecular interaction forces (ChangandGao, 2003; SunZhang, 2003, Ni et al., 2010). The most often used continuum theory to analyze nano-scale structures is nonlocal elasticity developed in 70s to take small scale effects into account by formulating a constitutive relation with the stress at a point expressed as a function of the strains at all points of the domain instead of the strain at the same point as in the case with the classical elasticity theory (Edelen and Laws, 1971; Eringen, 1972, 1983). The recent book by Eringen (2002) provi- des a detailed account of the nonlocal theory. Continuum models were also implemented to study themechanical behaviour of grahene sheets, and in par- ticular, the buckling of single-layered graphene sheets byPradhan andMurmu (2009), Sakhaee-Pour (2009) andPradhan (2009) where nonlocal theories we- re employed. Vibrational behaviour of graphene sheets has been the subject of several studies due to its importance in many applications. Vibrations of single-layered graphene sheets using nonlocal models were studied byMurmu and Pradhan (2009), Shen et al. (2010) and Narendar and Gopalakrishnan Variational principles and natural boundary conditions... 623 (2010). Vibrations of multilayered graphene sheets were investigated byHe et al. (2005), Behfar and Naghdabadi (2005), Liew et al. (2006) and Jomehza- deh and Saidi (2011) based on the classical elasticity theory. More recently nonlocal continuum models were used in the study of vibrations of multilay- ered graphene sheets by Pradhan and Phadikar (2009), Pradhan and Kumar (2010), Ansari et al. (2010), and Pradhan and Kumar (2011). The objective of the present study is to derive the variational principles and the applicable boundary conditions involving the transverse vibrations of multilayered orthotropic graphene sheets using the nonlocal theory of ela- sticity discussed above. Previous studies on variational principles involving nano-structures include multi-walled nanotubes under buckling loads (Adali, 2008), and undergoing linear and nonlinear vibrations (Adali, 2009a,b) which are based on nonlocal theory of Euler-Bernoulli beams. The corresponding results based on nonlocal Timoshenko theory are given in Adali (2011) for nanotubes under buckling loads and in Kucuk et al. (2010) for nanotubes un- dergoing vibrations. In the present study, these results are extended to the case of multilayered graphene sheets undergoing transverse vibrations, and natural boundary conditions are derived which are fairly involved due to co- upling between the sheets and small size effects. Moreover Rayleigh quotient for freely vibrating graphene sheets is obtained. The governing equations of the vibrating multilayered graphene sheets constitute a system of partial dif- ferential equations and the variational formulation for this system is obtained by the semi-inverse method developed by He (1997, 2004). This method was applied to several problems of mathematical physics governed by a system of differential equations some examples of which can be found inHe (2005, 2006, 2007), Liu (2005), Zhou (2006). The variational formulations given in Adali (2008, 2009a,b, 2011) and in Kucuk et al. (2010) were also obtained by the semi-inverse method. 2. Governing equations A continuum model of multilayered graphene sheets is shown in Fig.1a with the van der Waals interaction between the adjacent layers depicted as elastic springs. For an n-layered graphene, the top layer is numbered as i = 1 and the bottom layer as i= n. Top view of a graphene sheet is shown in Fig.1b where a and b are the dimensions of the sheets in the x and y directions, respectively. Bending stiffnesses of the orthotropic sheets are given by D11, D12,D22 and D66 which are defined as (Pradhan and Phadikar, 2009) 624 S. Adali Fig. 1. Multileyered graphene sheets, (a) side view, (b) top view D11 = E1h 3 12(1−ν12ν21) D12 = ν12E2h 3 12(1−ν12ν21) D22 = E2h 3 12(1−ν12ν21) D66 = G12h 3 12 (2.1) where h is the thickness of the graphene sheet, E1 and E2 areYoung’smodu- li in the x and y directions, respectively, G12 is the shear modulus, and ν12 and ν21 are Poisson’s ratios. Let wi(x,y,t) indicate the transverse deflection of the i-th layer and η the small scale parameter of the nonlocal elastic the- ory as defined in Pradhan and Phadikar (2009), Pradhan and Kumar (2010, 2011). Then the differential equations governing the transverse vibrations of multilayered graphene sheets in the time interval t1 ¬ t ¬ t2 and based on the nonlocal theory of elasticity (Pradhan and Phadikar, 2009) are given as D1(w1,w2)=L(w1)−η 2N(w1)− c12∆w1+η 2c12∇ 2(∆w1) = f(x,y,t)−η2∇2f(x,y,t) Variational principles and natural boundary conditions... 625 Di(wi−1,wi,wi+1)=L(wi)−η 2N(wi)+ c(i−1)i∆wi−1− ci(i+1)∆wi (2.2) −η2c(i−1)i∇ 2(∆wi−1)+η 2ci(i+1)∇ 2(∆wi)= 0 for i=2,3, . . . ,n−1 Dn(wn−1,wn)=L(wn)−η 2N(wn)+c(n−1)n∆wn−1 −η2c(n−1)n∇ 2(∆wn−1)= 0 where f(x,y,t) is a transverse load acting on the topmost layer (i=1)which can also be taken as acting at the bottommost layer (i = n) due to the symmetry of the structure, the symbol ∆wi is defined as ∆wi ≡wi+1−wi (2.3) and L(wi) and N(wi) are differential operators given by L(wi)=D11 ∂4wi ∂x4 +2(D12+2D66) ∂4wi ∂x2∂y2 +D22 ∂4wi ∂y4 +m0 ∂2wi ∂t2 −m2 ( ∂4wi ∂x2∂t2 + ∂4wi ∂y2∂t2 ) (2.4) N(wi)=∇ 2 [ m0 ∂2wi ∂t2 −m2 ( ∂4wi ∂x2∂t2 + ∂4wi ∂y2∂t2 )] with ∇2 = ∂ 2 ∂x2 + ∂ 2 ∂y2 . In Eqs. (2.4), m0 = ρh and m2 = ρh 3/12. The coeffi- cient c(i−1)i is the interaction coefficient of van derWaals forces between the (i−1)-th and i-th layerswith i=2, . . . ,n. The constant η= e0α is amaterial parameter defining the small scale effect in the nonlocal elastic theory where e0 is an experimentally determined constant andhas to be determined for each material independently (Eringen, 1983). α is an internal characteristic length such as lattice parameter, size of grain, granular distance, etc. (Narendar and Gopalakrishnan, 2010). 3. Variational functional In the present section, the semi-inverse method (He, 1997, 2004) will be em- ployed in order to derive the variational formulation of the problem. For this purpose, we first define a trial variational functional V (w1,w2, . . . ,wn) given by V (w1,w2, . . . ,wn)=V1(w1,w2)+V2(w1,w2,w3)+ . . . (3.1) +Vn−1(wn−2,wn−1,wn)+Vn(wn−1,wn) 626 S. Adali where V1(w1,w2)=U(w1)−Ta(w1)−Tb(w1) + t2 ∫ t1 b ∫ 0 a ∫ 0 [(−f+η2∇2f)w1+F1(w1,w2)]dxdydt Vi(wi−1,wi,wi+1)=U(wi)−Ta(wi)−Tb(wi) (3.2) + t2 ∫ t1 b ∫ 0 a ∫ 0 Fi(wi−1,wi,wi+1)dxdydt for i=2,3, . . . ,n−1 Vn(wn−1,wn)=U(wn)−Ta(wn)−Tb(wn)+ t2 ∫ t1 b ∫ 0 a ∫ 0 Fn(wn−1,wn)dxdydt with the functionals U(wi), Ta(wi) and Tb(wi) defined as U(wi)= 1 2 t2 ∫ t1 b ∫ 0 a ∫ 0 [ D11 (∂2wi ∂x2 )2 +2D12 ∂2wi ∂x2 ∂2wi ∂y2 +D22 (∂2wi ∂y2 )2 +4D66 (∂2wi ∂x∂y )2] dxdydt Ta(wi)= 1 2 t2 ∫ t1 b ∫ 0 a ∫ 0 { m0 (∂wi ∂t )2 +m2 [(∂2wi ∂x∂t )2 + (∂2wi ∂y∂t )2]} dxdydt (3.3) Tb(wi)= η2 2 t2 ∫ t1 b ∫ 0 a ∫ 0 { 2m0 (∂2wi ∂x2 ∂2wi ∂t2 + ∂2wi ∂y2 ∂2wi ∂t2 ) +m2 [( ∂3wi ∂x2∂t )2 +2 ( ∂3wi ∂x∂y∂t )2 + ( ∂3wi ∂y2∂t )2]} dxdydt where i = 1,2, . . . ,n. It is observed that U(wi) represents the strain energy and Ta(wi) the kinetic energy of the i-th layer of the multilayered graphene sheet. The functional Tb(wi) arises due to small scale effects, i.e. the nonlo- cal theory used in modeling of the graphene sheet. Similarly, the expression ∫ t2 t1 ∫ b 0 ∫a 0 Fi(wi−1,wi,wi+1)dxdydt in equation (3.2)1 represents the potential energy due to van der Waals forces between the layers. Similarly the term ∫ t2 t1 ∫ b 0 ∫a 0 [(−f+η 2∇2f)w1]dxdydt represents thework done by external forces where the second term arises due to small scale effects. In equations (3.2), Fi(wi−1,wi,wi+1) denotes the unknown functions of wi−1,wi and wi+1, and Variational principles and natural boundary conditions... 627 their derivatives and should be determined such that the Euler-Lagrange equ- ations of variational functional (3.2)1 correspondtodifferential equations (2.2). These equations are given by L(w1)−η 2N(w1)+ 2 ∑ j=1 δFj δw1 =L(w1)−η 2N(w1) + 2 ∑ j=1 ∂Fj ∂w1 − 2 ∑ j=1 ∂ ∂x ( ∂Fj ∂w1x ) − 2 ∑ j=1 ∂ ∂y ( ∂Fj ∂w1y ) =0 L(wi)−η 2N(wi)+ i+1 ∑ j=i−1 δFj δwi =L(wi)−η 2N(wi) (3.4) + i+1 ∑ j=i−1 ∂Fj ∂wi − i+1 ∑ j=i−1 ∂ ∂x ( ∂Fj ∂wix ) − i+1 ∑ j=i−1 ∂ ∂y ( ∂Fj ∂wiy ) =0 L(wn)−η 2N(wn)+ n ∑ j=n−1 δFj δwn =L(wn)−η 2N(wn) + n ∑ j=n−1 ∂Fj ∂wn − n ∑ j=n−1 ∂ ∂x ( ∂Fj ∂wnx ) − n ∑ j=n−1 ∂ ∂y ( ∂Fj ∂wny ) =0 where i=2,3, . . . ,n−1 and the subscripts x, y and t denote differentiations with respect to that variable, and δFi/δwi is the variational derivative defined as δFi δwi = ∂Fi ∂wi − k=3 ∑ k=1 ∂ ∂ξk ( ∂Fi ∂wiξk ) + k=3 ∑ k=1 j=3 ∑ j=k ∂2 ∂ξk∂ξj ( ∂Fi ∂wiξkξj ) + . . . (3.5) where ξ1 = x, ξ2 = y and ξ3 = t. It is noted that the variational derivative δFi/δwi of Fi(wi−1,wi,wi+1) follows from the Euler-Lagrange equations of the functional ∫ t2 t1 ∫ b 0 ∫a 0 Fi(wi−1,wi,wi+1)dxdydt. Comparing equations (3.4) with equations (2.2), we observe that the following equations have to be satisfied for Euler-Lagrange equations (3.4) to represent the governing equ- ations (2.2) 2 ∑ j=1 δFj δw1 =−c12∆w1+η 2c12 (∂2∆w1 ∂x2 + ∂2∆w1 ∂y2 ) 628 S. Adali i+1 ∑ j=i−1 δFj δwi = c(i−1)i∆wi−1− ci(i+1)∆wi−η 2c(i−1)i (∂2∆wi−1 ∂x2 + ∂2∆wi−1 ∂y2 ) (3.6) +η2ci(i+1) (∂2∆wi ∂x2 + ∂2∆wi ∂y2 ) n ∑ j=n−1 δFj δwn = c(n−1)n∆wn−1−η 2c(n−1)n (∂2∆wn−1 ∂x2 + ∂2∆wn−1 ∂y2 ) From equations (3.6), it follows that F1(w1,w2)= c12 4 (∆w1) 2+ c12 4 η2 [(∂∆w1 ∂x )2 + (∂∆w1 ∂y )2] Fi(wi−1,wi,wi+1)= c(i−1)i 4 (∆wi−1) 2+ ci(i+1) 4 (∆wi) 2 + η2c(i−1)i 4 [(∂∆wi−1 ∂x )2 + (∂∆wi−1 ∂y )2] (3.7) + η2ci(i+1) 4 [(∂∆wi ∂x )2 + (∂∆wi ∂y )2] for i=2,3, . . . ,n−1 Fn(wn−1,wn)= c(n−1)n 4 (∆wn−1) 2+ η2c(n−1)n 4 [(∂∆wn−1 ∂x )2 + (∂∆wn−1 ∂y )2] With Fi given by equations (3.7), we observe that equations (3.4) are equiva- lent to equations (2.2), viz. D1(w1,w2)=L(w1)−η 2N(w1)+ 2 ∑ j=1 δFj δw1 = f−η2∇2f Di(wi−1,wi,wi+1)=L(wi)−η 2N(wi)+ i+1 ∑ j=i−1 δFj δwi =0 (3.8) Dn(wn−1,wn)=L(wn)−η 2N(wn)+ n ∑ j=n−1 δFj δwn =0 4. Free vibrations In the present Section, the variational principle and the Rayleigh quotient are given for the case of freely vibrating graphene sheets. Let the harmonicmotion of the i-th layer be expressed as wi(x,y,t) =Wi(x,y)e √ −1ωt (4.1) Variational principles and natural boundary conditions... 629 where ω is the vibration frequency and Wi(x,y) is the deflection amplitude. The equations governing the free vibrations are obtained by substituting equ- ation (4.1) into equations (2.2) with f(x,y,t)= 0 and replacing the deflection wi(x,y,t) by Wi(x,y).Theoperators L(wi) and N(wi) nowbecome LFV (Wi) and N(Wi) given by LFV (Wi)=D11 ∂4Wi ∂x4 +2(D12+2D66) ∂4Wi ∂x2∂y2 +D22 ∂4Wi ∂y4 −m0ω 2Wi +m2ω 2 (∂2Wi ∂x2 + ∂2Wi ∂y2 ) (4.2) N(Wi)=∇ 2 [ −m0ω 2Wi+m2ω 2 (∂2Wi ∂x2 + ∂2Wi ∂y2 )] The variational principle for the case of free vibrations is the same as the one given by equations (3.1) and (3.2) with the deflection wi(x,y,t) replaced by Wi(x,y), the triple integrals replaced by the double integrals with respect to x and y, i.e., ∫ b 0 ∫a 0 Fi(Wi−1,Wi,Wi+1)dxdy, and U(wi),Ta(wi) and Tb(wi) replaced by UFV (Wi), TFVa(Wi) and TFVb(Wi) given by UFV (Wi)= 1 2 b ∫ 0 a ∫ 0 [ D11 (∂2Wi ∂x2 )2 +2D12 ∂2Wi ∂x2 ∂2Wi ∂y2 +D22 (∂2Wi ∂y2 )2 +4D66 (∂2Wi ∂x∂y )2] dxdy TFVa(Wi)= 1 2 b ∫ 0 a ∫ 0 { m0W 2 i +m2 [(∂Wi ∂x )2 + (∂Wi ∂y )2]} dxdy (4.3) TFVb(Wi)=− η2 2 b ∫ 0 a ∫ 0 { m0 [(∂Wi ∂x )2 + (∂Wi ∂y )2] +m2 [(∂2Wi ∂x2 )2 +2 (∂2Wi ∂x∂y )2 + (∂2Wi ∂y2 )2]} dxdy The functions Fi(Wi−1,Wi,Wi+1) are of the same form as given by equations (3.7) since the functions Fi(Wi−1,Wi,Wi+1) are independentof time.Next the Rayleigh quotient is obtained for the vibration frequency ω from equations (3.1), (3.2) and (4.3) as ω2 =min Wi ∑n i=1UFVi(Wi)+ ∑n i=1 b ∫ 0 a ∫ 0 Fidxdy ∑n i=1[TFVa(Wi)+TFVb(Wi)] (4.4) 630 S. Adali where Fi (i=1,2, . . . ,n) are given by equations (3.7)with wi(x,y,t) replaced by Wi(x,y). 5. Boundary conditions After substituting equations (3.2) into functional (3.1), we take its first varia- tion with respect to wi in order to derive the natural and geometric boundary conditions. The first variations of V (w1,w2, . . . ,wn) with respect to wi, de- noted by δwiV , are given by δw1V (w1,w2, . . . ,wn)= δw1V1(w1,w2)+ δw1V2(w1,w2,w3) δwiV (w1,w2, . . . ,wn)= δwiVi−1(wi−2,wi−1,wi)+ δwiVi(wi−1,wi,wi+1) (5.1) +δwiVi+1(wi,wi+1,wi+2) for i=2,3, . . . ,n−1 δwnV (w1,w2, . . . ,wn)= δwnVn−1(wn−2,wn−1,wn)+ δwnVn(wn−1,wn) The first variation of Vi(wi−1,wi,wi+1) with respect to wi is given by δwiVi(wi−1,wi,wi+1)= δwiU(wi)− δwiTa(wi)− δwiTb(wi) (5.2) +δwi t2 ∫ t1 b ∫ 0 a ∫ 0 [ Fi(wi−1,wi,wi+1) ] dxdydt for i = 1,2, . . . ,n− 1,n. Let δwi denote the variation of wi satisfying the boundary conditions δwix(x,0, t) = 0 δwix(x,b,t) = 0 δwiy(0,y,t) = 0 δwiy(a,y,t)= 0 (5.3) where the following notation was used δ(∂wi/∂x) = δwix, δ(∂wi/∂y) = δwiy. Moreover, the deflections wi(x,y,t) and their space derivatives vanish at the end points t = t1 and t = t2, i.e., δwi(x,y,t1) = 0, δwi(x,y,t2) = 0, δwix(x,y,t1)= 0, δwix(x,y,t2)=0, etc. Next using the subscript notation for differentiation, i.e., wix = ∂wi/∂x, wiy = ∂wi/∂y etc., we derive the first variations δwiU(wi), δwiTa(wi), δwiTb(wi) and δwi ∫ t2 t1 ∫ b 0 ∫a 0 [Fi(wi−1,wi,wi+1)]dxdydt by integration by parts to obtain Variational principles and natural boundary conditions... 631 δwiU(wi)= t2 ∫ t1 b ∫ 0 a ∫ 0 (D11wixxδwixx+D12wixxδwiyy +D12wiyyδwixx +D22wiyyδwiyy +4D66wixyδwixy ) dxdydt = t2 ∫ t1 b ∫ 0 a ∫ 0 [D11wixxxx+2(D12+2D66)wixxyy +D22wiyyyy]δwidxdydt +B1(wi,δwi) δwiTa(wi)= t2 ∫ t1 b ∫ 0 a ∫ 0 [m0witδwit+m2(wixtδwixt+wiytδwiyt)]dxdy dt =− t2 ∫ t1 b ∫ 0 a ∫ 0 [m0witt−m2(wixxtt+wiyytt)]δwidxdydt+B2(wi,δwi) (5.4) δwiTb(wi)= η 2 t2 ∫ t1 b ∫ 0 a ∫ 0 m0(wittδwixx+wixxδwitt+wittδwiyy +wiyyδwitt)dxdydt. . . +η2 t2 ∫ t1 b ∫ 0 a ∫ 0 m2(wixxtδwixxt+2wixytδwixyt+wiyytδwiyyt)dxdydt = η2 t2 ∫ t1 b ∫ 0 a ∫ 0 ∇2[m0witt−m2(wixxtt+wiyytt)]δwidxdydt+B3(wi,δwi) t2 ∫ t1 b ∫ 0 a ∫ 0 δwi[Fi(wi−1,wi,wi+1)]dxdydt = 1 2 t2 ∫ t1 b ∫ 0 a ∫ 0 (c(i−1)i∆wi−1− ci(i+1)∆wi)δwidxdydt. . . + η2 2 t2 ∫ t1 b ∫ 0 a ∫ 0 [−c(i−1)i∇ 2(∆wi−1)+ ci(i+1)∇ 2(∆wi)]δwidxdydt +B4(wi,δwi) where Bk(wi,δwi), k=1, . . . ,4 are the boundary terms. B1(wi,δwi) is given by 632 S. Adali B1(wi,δwi)= k=3 ∑ k=1 B1k(wi,δwi) (5.5) where B11(wi,δwi)= t2 ∫ t1 b ∫ 0 [D11(wixxδwix−wixxxδwi) +D12(wiyyδwix−wixyyδwi)] ∣ ∣ ∣ x=0 x=0 dydt B12(wi,δwi)= t2 ∫ t1 a ∫ 0 [D12(wixxδwiy −wixxyδwi) (5.6) +D22(wiyyδwiy −wiyyyδwi)] ∣ ∣ ∣ y=b y=0 dxdt B13(wi,δwi)=−2D66 t2 ∫ t1 b ∫ 0 wixyyδwi ∣ ∣ ∣ x=0 x=0 dydt−2D66 t2 ∫ t1 a ∫ 0 wixxyδwi ∣ ∣ ∣ y=b y=0 dxdt Similarly, B2(wi,δwi) and B3(wi,δwi) are given by B2(wi,δwi)=−m2 t2 ∫ t1 b ∫ 0 wixttδwi ∣ ∣ ∣ x=0 x=0 dydt−m2 t2 ∫ t1 a ∫ 0 wiyttδwi ∣ ∣ ∣ y=b y=0 dxdt (5.7) B3(wi,δwi)= k=3 ∑ k=1 B3k(wi,δwi) where B31(wi,δwi)= η 2m0 t2 ∫ t1 b ∫ 0 (wittδwix−wixttδwi) ∣ ∣ ∣ x=0 x=0 dydt +η2m0 t2 ∫ t1 a ∫ 0 (wittδwiy −wiyttδwi) ∣ ∣ ∣ y=b y=0 dxdt (5.8) B32(wi,δwi)= η 2m2 t2 ∫ t1 b ∫ 0 [−wixxttδwix+(wixxxtt+wixxytt)δwi] ∣ ∣ ∣ x=0 x=0 dydt B33(wi,δwi)= η 2m2 t2 ∫ t1 a ∫ 0 [−wiyyttδwiy +(wiyyytt+wixyytt)δwi] ∣ ∣ ∣ y=b y=0 dxdt Variational principles and natural boundary conditions... 633 Finally, we have B4(w1,δw1)= η2 2 ( t2 ∫ t1 b ∫ 0 c12∆w1xδw1 ∣ ∣ ∣ x=0 x=0 dydt+ t2 ∫ t1 a ∫ 0 c12∆w1yδw1 ∣ ∣ ∣ y=b y=0 dxdt ) B4(wi,δwi)= η2 2 t2 ∫ t1 b ∫ 0 (c(i−1)i∆w(i−1)x+ ci(i+1)∆wix)δwi ∣ ∣ ∣ x=0 x=0 dydt (5.9) + η2 2 t2 ∫ t1 a ∫ 0 (c(i−1)i∆w(i−1)y + ci(i+1)∆wiy)δwi ∣ ∣ ∣ y=b y=0 dxdt for i=2, . . . ,n−1 B4(wn,δwn)= η2 2 ( t2 ∫ t1 b ∫ 0 c(n−1)n∆w(n−1)xδwn ∣ ∣ ∣ x=0 x=0 dydt + t2 ∫ t1 a ∫ 0 c(n−1)n∆w(n−1)yδwn ∣ ∣ ∣ y=b y=0 dxdt ) Using the fundamental lemma of calculus of variations, the boundary condi- tions at x = 0,a and y = 0,b are obtained from equations (5.5)-(5.9) for i=2, . . . ,n−1. The boundary conditions at x=0,a are given by D11wixx+D12wiyy +η 2(−m0witt+m2wixxtt)= 0 or wix =0 −D11w1xxx−D12w1xyy −2D66w1xyy +m2w1xtt +η2[m0w1xtt−m2(w1xxxtt+w1xxytt+ c12∆w1x)] = 0 or w1 =0 −D11wixxx−D12wixyy −2D66wixyy +m2wixtt +η2[m0wixtt−m2(wixxxtt+wixxytt)] (5.10) +η2(c(i−1)i∆w(i−1)x+ ci(i+1)∆wix)= 0 o wi =0 for i=2, . . . ,n−1 −D11wnxxx−D12wnxyy −2D66wnxyy+m2wnxtt +η2[m0wnxtt−m2(wnxxxtt+wnxxytt+ c(n−1)n∆w(n−1)nx)] = 0 or wn =0 and at y=0,b by D12wixx+D22wiyy +η 2(−m0witt+m2wiyytt)=0 or wiy =0 −D12w1xxy −D22w1yyy −2D66w1xxy+m2w1ytt +η2[m0w1ytt−m2(w1yyytt+w1xyytt)]+η 2c12∆w1y =0 or w1 =0 634 S. Adali −D12wixxy−D22wiyyy −2D66wixxy+m2wiytt +η2[m0wiytt−m2(wiyyytt+wixyytt)] (5.11) +η2(c(i−1)i∆w(i−1)y + ci(i+1)∆wiy)= 0 or wi =0 for i=2, . . . ,n−1 −D12wnxxy−D22wnyyy −2D66wnxxy+m2wnytt +η2[m0wnytt−m2(wnyyytt+wnxyytt+ c(n−1)n∆w(n−1)ny)] = 0 or wn =0 It is observed that when η 6=0, the natural boundary conditions are coupled, that is, the nonlocal formulation of the problem leads to natural boundary conditionswhich containderivatives of wi−1 and wi+1 in the expression for wi, e.g. see the first equations of (5.10)3 and (5.11)4. 6. Conclusions The variational formulations for the free and forced vibrations of multilayered graphene sheets were derived using a continuum formulation based on the nonlocal orthotropic plate theory. The nonlocal theory used in the formula- tion allows the inclusion of small size effects and as such improves the accuracy of the model. A semi-inverse approach was employed in the derivation of the variational principles and the Rayleigh quotient for free vibrations was obta- ined. The formulation was used to obtain the natural boundary conditions. The variational principles presented here may form the basis of approxima- te and numerical methods of solution such as the Rayleigh-Ritz and finite elementmethods based on the energy functional of the problem andmay faci- litate the implementation of complicated boundaryconditions. Itwas observed that the nonlocal theory leads to coupled boundary conditions as opposed to uncoupled natural boundary conditions in the case of local theory of graphene sheets. Acknowledgement The research reported in the present paper was supported by a grant from the National ResearchFoundation (NRF) of SouthAfrica. This support is gratefully ack- nowledged. Variational principles and natural boundary conditions... 635 References 1. Adali S., 2008, Variational principles for multi-walled carbon nanotubes un- dergoing buckling based on nonlocal elasticity theory,Physics Letters A, 372, 5701-5705 2. Adali S., 2009a, Variational principles for transversely vibratingmulti-walled carbonnanotubesbasedonnonlocalEuler-Bernoullibeammodel,NanoLetters, 9, 5, 1737-1741 3. Adali S., 2009b, Variational principles formulti-walled carbon nanotubes un- dergoingnonlinear vibrations by semi-inversemethod,Micro andNanoLetters, 4, 198-203 4. Adali S., 2011, Variational formulation for buckling of multi-walled carbon nanotubes modelled as nonlocal Timoshenko beams, Journal of Theoretical an Applied Mechanics, to appear 5. Ansari R., Rajabiehfard R., Arash B., 2010, Nonlocal finite elementmo- del for vibrations of embedded multi-layered graphene sheets, Comp. Mater. Sci., 49, 831-838 6. Arsat R., Breedon M., Shafiei M., Spizziri P.G., Gilje S., Kaner R.B., Kalantar-Zadeh K., Wlodarski W., 2009, Graphene-like nano- sheets for surface acoustic wave gas sensor applications,Chemical Physics Let- ters, 467, 4/6, 344-347 7. BehfarK., Naghdabadi R., 2005,Nanoscale vibrational analysis of amulti- layered graphene sheet embedded in an elastic medium, Composites Science and Technology, 65, 1159-1164 8. Chang T., Gao H., 2003, Size-dependent elastic properties of a single-walled carbon nanotube via a molecular mechanics model, Journal of Mechanics and Physics of Solids, 51, 1059-1074 9. Choi S.M., Seo M.H., Kim H.J., Kim W.B., 2011, Synthesis of surface- functionalizedgraphenenanosheetswithhighPt-loadingsandtheir applications to methanol electrooxidation,Carbon, 49, 3, 904-909 10. Edelen D.G.B., Laws N., 1971, On the thermodynamics of systems with nonlocality,Archive for Rational Mechanics and Analysis, 43, 24-35 11. Eringen A.C., 1972, Linear theory of nonlocal elasticity and dispersion of plane waves, International Journal of Engineering Science, 10, 425-435 12. Eringen A.C., 1983, On differential of non-local elasticity and solutions of screw dislocation and surfacewaves, Journal of Applied Physics, 54, 4703-4710 13. EringenA.C., 2002,Nonlocal ContinuumField Theories, Springer,NewYork 636 S. Adali 14. Feng L., Chen Y., Ren J., Qu X., 2011, A graphene functionalized electro- chemical apta-sensor for selective label-free detection of cancer cells,Biomate- rials, In press, Available online 22 January 2011 15. Gao Y., Hao P., 2009, Mechanical properties of monolayer graphene under tensile and compressive loading,Physica E, 41, 1561-1566 16. He J.-H., 1997, Semi-inverse method of establishing generalized variational principles for fluidmechanics with emphasis on turbomachinery aerodynamics, International Journal of Turbo Jet-Engines, 14, 23-28 17. He J.-H., 2004, Variational principles for some nonlinear partial differential equations with variable coefficients,Chaos, Solitons and Fractals, 19, 847-851 18. He J.-H., 2005,Variational approach to (2+1)-dimensional dispersive longwa- ter equations,Phys. Lett. A, 335, 182-184 19. He J.-H., 2006, Variational theory for one-dimensional longitudinal beam dy- namics,Phys. Lett. A, 352, 276-277 20. He J.-H., 2007, Variational principle for two-dimensional incompressible invi- scid flow,Phys. Lett. A, 371, 39-40 21. He L.H., Lim C.W., Wu B.S., 2004, A continuummodel for size-dependent deformation of elastic films of nano-scale thickness, Int. J. Solids Struct., 41, 847-857 22. He X.Q., Kitipornchai S., Liew K.M., 2005, Resonance analysis of multi- layeredgraphenesheetsusedasnanoscale resonators,Nanotechnology,16, 2086- 2091 23. Hemmasizadeh A., MahzoonM., Hadi E., Khandan R., 2008, Amethod for developing the equivalent continuummodel of a single layer graphene sheet, Thin Solid Films, 516, 7636-7640 24. Jomehzadeh E., Saidi A.R., 2011, A study on large amplitude vibration of multilayered graphene sheets,Comp. Mater. Sci., 50, 1043-1051 25. Kitipornchai S., He X.Q., Liew K.M., 2005, Continuum model for the vibration of multilayered graphene sheets,Phys. Rev. B, 72, 075443 26. Kucuk I., Sadek I.S.,Adali S., 2010,Variational principles formulti-walled carbon nanotubes undergoing vibrations based on nonlocal Timoshenko beam theory, Journal of Nanomaterials, V. 2010, 1-7 27. Kuilla T., Bhadra S., Yao D., Kim N.H., Bose S., Lee J.H., 2010, Recent advances in graphene based polymer composites, Progress in Polymer Science, 35, 11, 1350-1375 28. LianP., ZhuX., Liang S., Li Z., YangW.,WangH., 2010, Large reversi- ble capacityofhighqualitygraphene sheets as ananodematerial for lithium-ion batteries,Electrochimica Acta, 55, 12, 3909-3914 Variational principles and natural boundary conditions... 637 29. Liew K.M., He X.Q., Kitipornchai S., 2006, Predicting nano vibration of multi-layered graphene sheets embedded in an elastic matrix,ActaMater., 54, 4229-4236 30. Liu H.-M., 2005, Generalized variational principles for ion acoustic plasma waves by He’s semi-inversemethod, Chaos, Solitons and Fractals, 23, 573-576 31. Mishra A.K., Ramaprabhu S., 2011, Functionalized graphene sheets for arsenic removal and desalination of seawater,Desalination, In press, Available online 11 February 2011 32. Murmu T., Pradhan S.C., 2009, Vibration analysis of nano-single-layered graphene sheets embedded in elastic medium based on nonlocal elasticity the- ory, J. Appl. Phys., 105, 064319 33. Narendar S., Gopalakrishnan S., 2010, Strong nonlocalization induced by small scale parameter on terahertz flexural wave dispersion characteristics of a monolayer graphene,Physica E, 43, 423-430 34. Ni Z., BuH., ZouM.,YiH., BiK., ChenY., 2010,Anisotropicmechanical properties of graphene sheets frommolecular dynamics,Physica B: Condensed Matter, 405, 5, 1301-1306 35. Poot M., Van Der Zant H.S.J., 2008, Nanomechanical properties of few- layer graphenemembranes,Appl. Phys. Lett., 92, 063111 36. PradhanS.C., 2009,Buckling of single layer graphene sheet basedonnonlocal elasticity and higher order shear deformation theory,Phys. Lett. A, 373, 4182- 4188 37. Pradhan S.C., Kumar A., 2010, Vibration analysis of orthotropic graphene sheets embedded in Pasternak elastic medium using nonlocal elasticity theory and differential quadraturemethod,Comp. Mater. Sci., 50, 239-245 38. Pradhan S.C., Kumar A., 2011, Vibration analysis of orthotropic graphe- ne sheets using nonlocal elasticity theory and differential quadrature method, Composite Structures, 93, 774-779 39. Pradhan S.C.,MurmuT., 2009, Small scale effect on the buckling of single- layered graphene sheets under biaxial compression via nonlocal continuumme- chanics,Comput. Mater. Sci., 47, 268-274 40. Pradhan S.C., Phadikar J.K., 2009, Small scale effect on vibration of em- bedded multilayered graphene sheets based on nonlocal continuum models, Phys. Lett. A, 373, 1062-1069 41. Sakhaee-Pour A., 2009, Elastic properties of single-layered graphene sheet, Solid State Communications, 149, 1/2, 91-95 42. Sakhaee-Pour A., 2009, Elastic buckling of single-layered graphene sheet, Comput. Mater. Sci., 45, 266-270 638 S. Adali 43. Shen L., ShenH.-S., ZhangC.-L., 2010,Nonlocal platemodel for nonlinear vibrationof single layergraphenesheets in thermalenvironments,Comp.Mater. Sci., 48, 680-685 44. ShokriehM.M.,RafieeR., 2010,PredictionofYoung’smodulus of graphene sheets and carbon nanotubes using nanoscale continuum mechanics approach, Materials and Design, 31, 2, 790-795 45. SoldanoC., MahmoodA., Dujardin E., 2010, Production, properties and potential of graphene,Carbon, 48, 8, 2127-2150 46. Sun C.T., Zhang H.T., 2003, Size-dependent elastic moduli of platelike na- nomaterials, Journal of Applied Physics, 93, 1212-1218 47. Terrones M., Botello-Méndez A.R., Campos-Delgado J., López- Uŕıas F., Vega-Cantú Y.I., Rodŕıguez-Maćıas F.J., Eĺıas A.L., MUñoz-Sandoval E., Cano-Márquez A.G., Charlier J.-C., Terro- nes H., 2010, Graphene and graphite nanoribbons: Morphology, properties, synthesis, defects and applications,Nano Today, 5, 4, 351-372 48. WuW., Liu Z., Jauregui L.A.,YuQ., PillaiR., CaoH., Bao J., Chen Y.P., Pei S.-S., 2010, Wafer-scale synthesis of graphene by chemical vapor deposition and its application in hydrogen sensing, Sensors and Actuators B: Chemical, 150, 296-300 49. Yang M., Javadi A., Gong S., 2010, Sensitive electrochemical immuno- sensor for the detection of cancer biomarker using quantum dot functionali- zed graphene sheets as labels, Sensors and Actuators B: Chemical, In Press, Available online 2 December 2010 50. Yuan C., Hou L., Yang L., Fan C., Li D., Li J., Shen L., Zhang F., Zhang X., 2011, Interface-hydrothermal synthesis of Sn3S4/graphene sheet composites and their application in electrochemical capacitors, Mater. Lett., 65, 2, 374-377 51. ZhouW.X., 2006,Variational approach to the Broer-Kaup-Kupershmidt equ- ation,Phys. Lett. A, 363, 108-109 Zasady wariacyjne i naturalne warunki brzegowe dla wielowarstwowych ortotropowych paneli grafenowych poddanych drganiom, sformułowane w ramach nielokalnej teorii sprężystości Streszczenie W pracy zajęto się problemem drgań poprzecznych ortotropowych paneli grafe- nowych, dla których sformułowano zasady wariacyjne na podstawie nielokalnej teo- rii sprężystości, co pozwoliło na budowę ciągłego modelu takich struktur. Formuła Variational principles and natural boundary conditions... 639 wariacyjna umożliwiła konstrukcję naturalnych warunków brzegowych wyrażonych zbiorem sprzężonych równań opisujących grafenowe panele wielowarstwowew odróż- nieniu od rozprzężonych warunków brzegowych stosowanych jedynie do zamocowań typu swobodne podparcie lub zamurowanie, jednocześnie przy zastosowaniu lokalnej (klasycznej) teorii sprężystości. Dla przypadku drgań swobodnych wyznaczono iloraz Rayleigha układu z grafenu. W prezentowanym sformułowaniu użyto odpowiednich technik obliczania funkcjonałów i półodwrotnej metody wyznaczania całek. Wyka- zano, że postać wariacyjna stanowi podstawę dla numerycznych metod poszukiwa- nia przybliżonych rozwiązań i pogłębia zrozumienie zachodzących zjawisk fizycznych w takich układach. Manuscript received April 18 2011; accepted for print June 1, 2011