10736 FACTA UNIVERSITATIS Series: Mechanical Engineering Vol. 20, No 3, 2022, pp. 561 - 587 https://doi.org/10.22190/FUME220506029L © 2022 by University of Niš, Serbia | Creative Commons License: CC BY-NC-ND Original scientific paper BENDING, BUCKLING AND FREE VIBRATION ANALYSES OF NANOBEAM-SUBSTRATE MEDIUM SYSTEMS Suchart Limkatanyu1, Worathep Sae-Long2, Jaroon Rungamornrat3, Chinnapat Buachart4, Piti Sukontasukkul5, Suraparb Keawsawasvong6, Prinya Chindaprasirt7,8 1Department of Civil and Environmental Engineering, Faculty of Engineering, Prince of Songkla University, Songkhla, Thailand 2Civil Engineering Program, School of Engineering, University of Phayao, Phayao, Thailand 3Applied Mechanics and Structures Research Unit, Department of Civil Engineering, Faculty of Engineering, Chulalongkorn University, Bangkok, Thailand 4Department of Civil Engineering, Faculty of Engineering, Chiang Mai University, Chiang Mai, Thailand 5Construction and Building Materials Research Center, Department of Civil Engineering, King Mongkut’s University of Technology North Bangkok, Bangkok, Thailand 6Department of Civil Engineering, Faculty of Engineering, Thammasat School of Engineering, Thammasat University, Pathumthani, Thailand 7Sustainable Infrastructure Research and Development Center, Department of Civil Engineering, Faculty of Engineering, Khon Kaen University, Khon Kaen, Thailand 8Academy of Science, Royal Society of Thailand, Dusit, Bangkok, Thailand Abstract. This study presents a newly developed size-dependent beam-substrate medium model for bending, buckling, and free-vibration analyses of nanobeams resting on elastic substrate media. The Euler-Bernoulli beam theory describes the beam-section kinematics and the Winkler-foundation model represents interaction between the beam and its underlying substrate medium. The reformulated strain-gradient elasticity theory possessing three non-classical material constants is employed to address the beam-bulk material small- scale effect. The first and second constants is associated with the strain-gradient and couple- stress effects, respectively while the third constant is related to the velocity-gradient effect. The Gurtin-Murdoch surface elasticity theory is adopted to account for the surface-free energy. To obtain the system governing equation as well as corresponding boundary conditions, Hamilton’s principle is called for. Three numerical simulations are presented to characterize the influences of the material small-scale effect, the surface-energy effect, and the surrounding substrate medium on bending, buckling, and free vibration responses of Received: May 06, 2022 / Accepted June 11, 2022 Corresponding author: Worathep Sae-Long Civil Engineering Program, School of Engineering, University of Phayao, Phayao, 56000, Thailand E-mail: worathep.sa@up.ac.th 562 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. nanobeam-substrate medium systems. The first simulation focuses on the bending response and shows the ability of the proposed model to eliminate the paradoxical characteristic inherent to nanobeam models proposed in the literature. The second and third simulations perform the sensitivity investigation of the system parameters on the buckling load and the natural frequency, respectively. All analytical results reveal that both material small-scale and surface-energy effects consistently stiffen the system response while the velocity-gradient effect weakens the system response. Furthermore, these sized-scale effects are more pronounced when the underlying substrate medium becomes softer. Key words: Reformulated strain gradient elasticity theory, Small-scale effect, Surface- energy effect, Bending analysis, Buckling analysis, Free vibration analysis 1. INTRODUCTION With the discovery of carbon nanotubes (CNTs) by Iijima [1] in 1991, distinct responses of nano-sized structures have drawn high attention from researchers worldwide. As a result, several innovative devices and systems have been possible with superior and unique mechanical characteristics of structures at nanoscale [2-4]. To expedite and support the development of these innovative devices and systems, a thorough understanding of characteristics and behaviours of their nano-sized structural components is of paramount importance. In general, the experimental method is the most straightforward way to characterize mechanical responses of structures. However, there are several limitations on experiments conducted on nano-sized specimens, such as the necessity of a specifically testing apparatus, the requirement of a special testing procedure, and the expensive testing cost, etc. [5-7]. Furthermore, available classical structural models (e.g. bar, beam, shell, and plate models) are inadequate to characterize the structural response at nanoscale since the material small-scale effect and the size dependency inherent to nano-sized structures are not considered [5-13]. Therefore, a rational mathematical model capable of representing the material small- scale effect and the size dependency is deemed necessary and is the main emphasis of the present study. In the research community, both atomistic-based model [14-17] and enhanced structural mechanics-based model [18-20] have been employed to characterize mechanical responses of nano-sized structures for a wide spectrum of applications [21-24]. Although comprehensive details on the mechanical responses of nano-sized structures can be gained with the atomistic- based model [14-17], expensive computational costs and sophisticated modelling efforts are usually required with this numerical model, thus limiting the model access by researchers. Consequently, the enhanced structural mechanics-based model [18-20] have gained increasing popularity due to its good compromise between model efficiency and model accuracy [11]. The enhanced structural mechanics-based model is possible with unification of non-classical elasticity theories [25-33] and classical structural models. When the structural dimension is in the order of nanometer, the material small-scale effect induced by long-range interatomic forces becomes pronounced. A family of higher-order elasticity theories have been proposed in the research community to address this material small-scale effect [34-39]. Among these higher-order elasticity theories, the nonlocal elasticity theory proposed by Eringen [25-26] has been the most popular and was first employed by Peddieson et al. [34] to construct the nano-sized beam model. However, the material small-scale effect was not detected by this beam model in the case of the cantilever beam under an end load [35, 38, 40-42]. This peculiar response was Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 563 defined as a “paradox” in the literature [34, 38, 40-43]. Romano et al. [43] demonstrated that an ill-posed structural-mechanics problem (bounded domain) arises with adoption of the Eringen nonlocal differential model. Furthermore, Koutsoumaris et al. [44] demonstrated that the Eringen nonlocal differential model does not conform to the requirement of the quadratic energy functional form of elasticity. To obtain a more rational small-scale beam model, the strain-gradient elasticity theory has been employed as an alternative [38, 45-49]. The strain-gradient elasticity theory was originally invented by Mindlin [27-28] and possesses three forms, namely: Form I, Form II, and Form III. The difference among these three rudimentary forms is their expressions of the strain energy density function. Form I expresses the strain energy density function in terms of the infinitesimal strain and the second displacement gradient; Form II describes the strain energy density function in terms of the infinitesimal strain and the first strain gradient; and Form III represents the strain energy density function in terms of the infinitesimal strain, the symmetric part of the second displacement gradient; and the rotation gradient (curvature tensor). Several researchers [5, 30, 38, 50-54] have formulated their own strain-gradient elasticity models based on these three forms (Form I, Form II, and Form III). As the simplest variant of the strain-gradient elasticity theory, the simplified strain-gradient elasticity theory of Altan and Aifantis [51] is constructed based on Form II but can account only for the strain gradient along the longitudinal direction. To obtain a more refined variant of the strain-gradient elasticity theory, the modified strain-gradient elasticity theory of Lam et al. [5] is formulated based on either Form I or Form III and contains three material length-scale parameters associated with three strain-gradient measures, namely: the dilatation gradient, the deviatoric stretch gradient, and the symmetric rotation gradient. A thorough study by Papargyri-Beskou et al. [55] has shown the essence of the velocity gradient on wave dispersion in gradient elastic solids and structures. However, the modified strain-gradient elasticity theory of Lam et al. [5] lacks the material length-scale parameter associated with the velocity gradient. Therefore, a more versatile strain-gradient elasticity theory is required to construct a mathematical tool for static and dynamic analyses of nano-sized structures. Recently, Zhang and Guo [30] has proposed the reformulated strain-gradient elasticity theory based on Form I. Strain-gradient measures as well as velocity gradient are both considered in this strain-gradient-type theory. Consequently, the reformulated strain-gradient elasticity theory is preferable to the narrowness of the research gap and is adopted in the present work. In addition to the material small-scale effect, the size dependency invoked by the surface-energy effect is crucial when the structural dimension approaches the range of nanometer [36, 38, 56-59]. To incorporate the surface-energy effect into a conventional structural-mechanics model, the surface elasticity theory proposed by Gurtin and Murdoch [32-33] has been widely employed. In the Gurtin-Murdoch surface elasticity theory, the surface elastic layer is assumed to be a zero-thickness membrane perfectly bonded to its wrapped bulk material. Due to the good compromise between model accuracy and simplicity, several structural-mechanics models have been armed with the ability to account for the size dependency using the Gurtin-Murdoch surface elasticity theory [11, 36, 38, 57-62]. As a result, the present work adopts the Gurtin-Murdoch surface elasticity theory to represent the size dependency induced by the surface-energy effect. The concept of a beam on elastic foundation has found a wide spectrum of applications in nanoengineering and nanoscience [21-24]. Due to its good compromise between model accuracy and simplicity, the Winkler foundation model is the most widely employed foundation model to account for the interactive mechanism between nanobeams and their 564 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. underlying substrate media [9, 38, 59, 63]. During the last decade, several nanobeam- substrate medium models have been proposed in the research community to characterize the nanobeam-substrate medium system [64-68]. For example, Azizi et al. [64] assessed the influence of the surface-free energy on the nonlinear vibration characteristics of the simply-supported nanobeam-substrate medium system; Niknam and Aghdam [65] derived the analytical solution for buckling and vibration analyses of the nonlocal functionally graded (FG) beam on the elastic foundation; Demir [66] employed the differential transform method (DTM) to compute the natural frequencies of simply supported and clamped-clamped nanobeams on elastic foundation; Ponbunyanon et al. [67] extended the Winkler- Pasternak based beam-foundation of Limkatanyu et al. [9] to study the flexural behaviour of the nanobeam-substrate medium system; and Jena et al. [68] performed buckling analyses of single-walled carbon nanotubes on elastic substrate media under both low and high temperature environments. To the best of the authors’ knowledge, the reformulated strain-gradient elasticity theory of Zhang and Guo [30] has not been applied to the problem of nanobeams on substrate media. Therefore, there is still room to add a rational nanobeam-substrate model into the research community due to the merit of the reformulated strain-gradient elasticity theory. The objective of the present work is to propose the rational mathematical model for static and free vibration analyses of nanobeams on substrate media. The presentation of the present work is in the following order. First, the kinematics of the Euler-Bernoulli beam theory is described. Then, the reformulated strain-gradient elasticity theory of Zhang and Guo [30] and the surface elasticity theory of Gurtin and Murdoch [32-33] are briefly discussed. Next, the Winkler-based interactive mechanism between the beam and its underlying substrate medium is presented. Hamilton’s principle is employed to obtain the system governing equation and classical as well as non-classical boundary conditions for bending, buckling, and free vibration analyses. Finally, three numerical simulations are provided to demonstrate behaviours and characteristics of the proposed nanobeam- substrate medium model. The first simulation is employed to demonstrate the ability of the proposed model to overcome the paradoxical response found in the literature [34, 38, 40-42] and to investigate the material small-scale and surface-energy effects as well as the nanobeam-substrate medium interaction on the bending responses. The second and third simulations examine the variation of the system parameters on the critical buckling load and the natural frequency of the simply-supported nanobeam-substrate medium systems. All symbolic calculations in the present work are carried out using the computer software Mathematica [69]. 2. EULER-BERNOULLI BEAM THEORY In the present work, the kinematic assumption of the proposed beam model follows the Euler-Bernoulli beam hypothesis [70] asserting that “Plane section remains plane and still normal to the longitudinal axis”. Based on this kinematical constraint, the displacements of a generic point along the x-, y-, and z- axes are: 0 0 ( ) ( , ) ; ( , ) ( ); and ( , ) 0 x y z v x U x y y U x y v x U x y x  = − = =  (1) Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 565 where Ux(x,y), Uy(x,y), and Uz(x,y) are displacement fields along the x-, y-, and z- axes, respectively; v0(x) is the vertical displacement of the beam section; and y is the vertical distance of the point measured from the reference axis x. 3. REFORMULATED STRAIN-GRADIENT ELASTICITY THEORY The reformulated strain-gradient elasticity theory of Zhang and Gao [30] is modified from Form I of Mindlin’s strain-gradient elasticity theory [27]. The assertion of Form I is that the strain energy density function is expressed as a function of the infinitesimal strain tensor and the second gradient of displacement. For a linear isotropic elastic material, Form I contains five material constants related to strain gradients in addition to two conventional material constants (Lame constants). To simplify Form I of Mindlin’s strain-gradient elasticity theory, Zhang and Gao [30] reduced the number of material length-scale parameters from five to two using the symmetry property of the tensor [71]. The first material length-scale parameter addresses the strain-gradient effect while the second material length-scale parameter considers the couple-stress effect. This simplified strain-gradient elasticity theory is named as the “reformulated” strain-gradient elasticity theory and is employed in this study to address the beam-bulk material small-scale effect. Based on the reformulated strain-gradient elasticity theory proposed by Zhang and Gao [30], the stored strain energy U is defined as [54]: 1 ( ) 2 L ij ij ijk ijk ij ij V U m dV       = + + (2) where L ij  represents the Cauchy stress tensor and is the conjugate-work pair of the strain tensor εij; ijk  represents the symmetric part of the double stress tensor and is the conjugate- work pair of the strain gradient tensor ijk  ; mij defines the symmetric part of the couple-stress tensor and is the conjugate-work pair of the curvature tensor χij; and V is the volume of the elastic body. It is worth mentioning that the stored strain energy U of Eq. (2) is identical to the one associated with the modified couple-stress theory of Yang et al. [29] when the strain- gradient effect is neglected. In the reformulated strain-gradient elasticity theory, the constitutive relations between stress and strain quantities are given by Yin et al. [54] as: 2 L ij kk ij ij    = + (3) 2 2 ijk s ijk l     = (4) 2 2 ij m ij m l = (5) where δij is the Kronecker delta; μ and λ are Lame constants; and ls and lm are the material length-scale parameters associated with the stain-gradient effect and the couple-stress effect, respectively. The strain εij and the strain-gradient measures ijk  , and χij are defined as: , , 1 ( ) 2 ij i j j i u u = + (6) 566 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. , , , 1 ( ) 3 ijk i jk j ki k ij u u u  = + + (7) , , 1 ( ) 2 ij i j j i   = + (8) where ui is the displacement field; and θi is the rotation vector defined as: , 1 2 i ijk k j e u = (9) with eijk being the permutation symbol. 4. SURFACE ELASTICITY THEORY As the structural dimension approaches the range of nanometer, the free energy at the surface of a bulk material can no longer be neglected since it is comparable to the energy stored in the bulk material. This surface-free energy induces the size dependency of nano- sized structures as confirmed by both experiments and numerical simulations [36, 57, 72]. To account for this size-dependent characteristic, the proposed nanobeam-substrate medium model employs the surface elasticity theory of Gurtin and Murdoch [32-33]. Based on the assumptions of this surface elasticity theory [32-33], the beam cross section is assumed to be composed of a bulk core and a wrapping outer surface layer shown in Fig. 1. The bulk core and its wrapping outer surface layer are assumed to be in a perfect bond condition. The wrapping outer surface layer is taken as a mathematically zero-thickness elastic layer. The constitutive relations of the surface layer are: 0 0 , 0 , , , [ ( ) ] ( ) sur sur sur sur sur sur sur sur sur sur res res u u u u               = + + + + − (10) , sur sur sur n res n u  = (11) where sur   and sur n  are, respectively, the in-plane and out-of-plane components of the surface stresses; sur res  is the residual surface stress under unconstrained conditions obtained from the atomistic simulation [72]; α and β are the in-plane Cartesian coordinates of the lateral surface for the nanobeam; usur is the surface-layer deformation; and 0 sur  and 0 sur  are the surface elastic constants. D b hz z y y n t nSurface Layer t Bulk Bulk Fig. 1 Beam cross-sections: Beam bulk and surface layer [32-33] Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 567 5. NANOBEAM-SUBSTRATE MEDIUM INTERACTION In the present work, the interaction between the nanobeam and its underlying substrate medium shown in Fig. 2(a) is simulated based on the concept of the Winkler foundation model [73]. Following this foundation concept, the interactive mechanism of the underlying substrate medium is represented by a set of continuously distributed and non-interconnected springs attached along the nanobeam length shown in Fig. 2(b). The relation between the substrate-medium interactive force Des(x) and the substrate-medium deformation Δes(x) is: ( ) ( ) es es es D x k x=  (12) where kes is the elastic modulus of the substrate medium. x L Elastic Substrate Medium ( )0,y v x ( )yq x Surface Layer Beam-Bulk 4 4 ,P U 1 1 ,P U es k Beam ( )0,y v x x 5 5 ,P U 2 2 ,P U  1 2 3 4 5 6 T P P P P P P=P  1 2 3 4 5 6 T U U U U U U=U Elastic Substrate-Medium Springs 3 3 ,P U 6 6 ,P U (a) (b) Fig. 2 Typical systems: (a) Nanobeam resting on elastic substrate medium; and (b) Beam- substrate medium model 6. REFORMULATED STRAIN-GRADIENT NANOBEAM-SUBSTRATE MEDIUM MODEL 6.1. Compatibility Conditions and Constitutive Relations 6.1.1. Beam bulk In the case of the planar Euler-Bernoulli beam problem, the transverse displacement v0(x,t) serves as a primary variable and is a function of the spatial variable x and the temporal variable t since both static and free vibration responses of nanobeam-substrate medium systems are of interest. Based on the beam kinematics of Eq. (1), non-zero components of the strain εij and the strain-gradient measures ijk  , and χij in Eqs. (6) to (9) are expressed in terms of the transverse displacement v0(x,t) as: 2 0 2 ( , ) ( , ) xx v x t x t y x   = −  (13) 3 0 3 ( , ) ( , ) xxx v x t x t y x    = −  (14) 568 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. 2 0 2 ( , )1 ( , ) 3 xxy v x t x t x    =  (15) 2 0 2 ( , )1 ( , ) ( , ) 2 xz zx v x t x t x t x    = =  (16) Substituting the compatibility relations of Eqs. (13) to (16) into the constitutive relations of Eqs. (3) to (5), the corresponding non-zero stress components in terms of the transverse displacement v0(x,t) are: 2 0 2 ( , ) ( , ) xx xx v x t x t E y x   = −  (17) 3 2 0 3 ( , ) ( , ) 2 xxx s v x t x t l y x     = −  (18) 2 2 0 2 ( , )2 ( , ) 3 xxy s v x t x t l x     =  (19) 2 2 0 2 ( , ) ( , ) ( , ) xz zx m v x t m x t m x t l x   = =  (20) where Exx defines the elastic modulus of the beam bulk. 6.1.2. Outer surface layer For the problem of the planar Euler-Bernoulli beam, the degenerated forms of Eqs. (10) and (11) are given by Guo and Mahmoud [57] as: ( , ) ( , ) sur sur sur sur xx res xx xx x t E x t  − = (21) ( , ) ( , ) sur sur sur nx res nx x t x t  = (22) where ( , ) sur xx x t denotes the in-plane component of the surface stress and is the conjugate- work pair of the in-plane component of the surface strain ( , ) sur xx x t ; ( , ) sur nx x t denotes the out-of-plane component of the surface stress and is the conjugate-work pair of the out-of- plane component of the surface strain ( , ) sur nx x t ; and 0 02 sur sur sur xx E  = + is the surface-layer elastic modulus obtained from two surface elastic constants 0 sur  and 0 sur  [72]. Following the beam kinematics of Eq. (1) and the full compatibility of the composite nanobeam section of Fig. 1, the in-plane ( , ) sur xx x t and out-of-plane surface strain ( , ) sur nx x t components can be expressed in terms of the transverse displacement v0(x,t) as: 2 0 2 ( , ) ( , ) sur xx v x t x t y x   = −  (23) 0 ( , ) ( , ) sur nx y v x t x t n x   =  (24) where ny defines the y component of the unit vector n normal to the lateral surface of beam section as shown in Fig. 1. Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 569 Substituting Eqs. (23) and (24) into Eqs. (21) and (22), the in-plane ( , ) sur xx x t and out- of-plane ( , ) sur nx x t surface stress components are described in terms of the transverse displacement v0(x,t) as: 2 0 2 ( , ) ( , ) ( , ) sur sur sur sur sur xx res xx xx xx v x t x t E x t E y x     − = = −  (25) 0 ( , ) ( , ) ( , ) sur sur sur sur nx res nx res y v x t x t x t n x      = =  (26) 6.1.3. Substrate medium Based on the assumption of the Winkler-foundation model [73], the beam and its supporting Winkler-spring bed is in a full compatibility condition (Δes(x,t) = v0(x,t)). Therefore, the substrate-medium interactive force Des(x) of Eq. (12) can be described in terms of the transverse displacement v0(x,t) as: 0 ( , ) ( , ) ( , ) es es es es D x t k x t k v x t=  = (27) 6.2. Equilibrium Equation: Hamilton’s Principle For an isotropic linear elastic material, the stored strain energy U of the nanobeam- substrate medium system of Fig. 2(a) can be defined as: 2 2 2 2 2 20 0 02 2 0 0 0 Elastic substrate medium contributionLocal stress contribution Couple stress contribution ( , ) ( , )1 1 1 ( ( , )) 2 2 2 1 L L L xx m es v x t v x t U E I dx Al dx k v x t dx x x       = + +         +    In-plane component Out-of-plane component 2 2 23 2 2 20 0 0 3 2 0 0 0 Strain gradient contribution ( , ) ( , ) ( , )1 2 1 2 2 2 3 2 L L L s s v x t v x t v x t Il dx Al dx P xx x          + −               Applied axial load contribution In-plane component Out-of 2 22 0 0 2 0 0 Surface-induced distributed moment contribution ( , ) ( , )1 1 2 2 L L sur sur xx res dx v x t v x t E I dx S dx xx        + +           -plane component Surface-induced shear force contribution (28) where P is the applied axial force; L is the nanobeam length; Г is the beam section perimeter; A A dA=  is the sectional area; 2 A I y dA=  is the second moment of area; 2 I y d   =  is the second moment of perimeter; and 2 y S n d   =  . The first variation of the stored strain energy U of Eq. (28) during the time interval [0, T] is: 570 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. 2 2 2 2 20 0 0 0 2 2 2 2 0 0 0 3 3 2 2 2 20 0 0 0 3 3 2 2 0 0 ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) ( , )2 2 3 T L L xx m L L s s v x t v x t v x t v x t U E I dx Al dx x x x x v x t v x t v x t v x t Il dx Al x x x x                   = +                          + +                   2 2 0 0 0 0 2 2 0 0 0 0 0 0 0 0 ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) ( ( , ))( ( , )) L L sur sur xx res L L es dx v x t v x t v x t v x t E I dx S dx x xx x v x t v x t k v x t v x t dx P dx dt x x                + +                  + −              (29) To remove all differential operators from the virtual transverse displacement δv0(x,t), the integration by parts is applied to Eq. (29), thus leading to the following expression: 6 4 2 0 0 0 0 06 4 2 0 0 5 3 0 0 0 0 5 3 0 0 ( , ) ( , ) ( , ) ( , ) ( ) ( ) ( , ) ( ) ( , ) ( , ) ( , ) ( , ) ( ) ( ) ( ) T L H L sur eff eff es res x L T H L sur eff eff res x v x t v x t v x t U v x t EI EI k v x t S P dxdt x x x v x t v x t v x t v x t EI EI S P dt xx x        =  =     = − + + − −          + − + −      +   4 2 0 0 0 4 2 0 0 2 3 0 0 2 3 0 0 ( , ) ( , ) ( , ) ( ) ( ) ( , ) ( , ) ( ) x L T H L eff eff x x L T H eff x v x t v x t v x t EI EI dt x x x v x t v x t EI dt x x  = = = =    − +         +        (30) where 2( ) 2 H eff s EI Il= denotes the higher-order flexural rigidity associated with the strain- gradient contribution; and 2 22 ( ) 3 L sur eff xx s m xx EI E I Al Al E I   = + + + denotes the lower-order flexural rigidity combining the beam-bulk and surface-layer contributions. Based on the kinetic energy expression given by Papargyri-Beskou et al. [55], the first variation of the kinetic energy K with the inclusion of the velocity gradient during the time interval [0, T] is [30]: 2 4 6 4 2 20 0 0 0 0 2 2 2 4 2 2 2 0 0 3 5 2 20 0 0 2 3 2 0 0 ( , ) ( , ) ( , ) ( , ) ( , ) 2 ( , ) ( , ) ( , ) (2 ) T L v v x L T v v x v x t v x t v x t v x t K v x t A l I l dx dt t x t x t x t v x t v x t v x t Al I Il x t x t         = =           = − − + −                       + − + +             4 20 0 2 2 0 0 ( , ) ( , ) x L T v x dt v x t v x t Il dt x x t   = =     + −          (31) where ρ represents the nanobeam mass density; and lv defines the material length-scale parameter associated with the velocity gradient. It should be noted from the first variation of the kinematic energy in Eq. (31) that the beam configurations are prescribed at t = 0 and t = T, resulting in the vanish of the virtual displacements and its derivatives at t = 0 Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 571 and t = T. Moreover, the mass density ρ over the time interval [0, T] is constant along the cross section of the nanobeam [30]. The first variation of the external work done W by the transverse distributed load qy(x,t) and end forces P shown in Fig. 2(b) during the time interval [0, T] can be defined as: 0 0 0 ( , ) ( , ) T L T y W q x t v x t dx dt     = +      U P (32) where the displacement vector U = [U1 U2 U3 U4 U5 U6] T collects end displacements of the system and the force vector P = [P1 P2 P3 P4 P5 P6] T contains end forces of the system. To establish the governing differential equation of motion and its associated boundary conditions, Hamilton’s principal is called for. As a result, Eqs. (30), (31), and (32) can be written together as: 0 6 4 2 0 0 0 0 6 4 2 0 0 2 4 6 4 2 20 0 0 0 0 2 2 2 4 2 [ ( )] 0 ( , ) ( , ) ( , ) ( , ) ( ) ( ) ( ) ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) 2 T T L H L sur eff eff res y es v v K U W dt v x t v x t v x t v x t EI EI S P x x x v x t v x t v x t v x t q x t k v x t A l I l t x t x t         − − =        = − + −         + − − − − −           2 2 5 3 0 0 0 0 5 3 0 3 5 2 20 0 2 3 2 0 4 0 0 4 ( , ) ( , ) ( , ) ( , ) ( ) ( ) ( ) ( , ) ( , ) (2 ) ( , ) ( , ) ( ) ( ) T H L sur eff eff res x L v v x H L eff eff dxdt x t v x t v x t v x t v x t EI EI S P xx x v x t v x t Al I Il dt x t x t v x t v x t EI EI x x        = =            + − + − −     − + +          + −    2 4 20 0 2 2 2 0 0 2 3 0 0 2 3 0 0 ( , ) ( , ) ( , ) ( , ) ( ) x L T v x x L T H T eff x v x t v x t Il dt x x t v x t v x t EI dt x x    = = = =   −         + − +       U P (33) Considering the arbitrariness of δv0(x,t) during the time interval [0, T], the governing differential equation of motion is unveiled as: 6 4 2 0 0 0 06 4 2 2 4 6 4 2 20 0 0 0 2 2 2 4 2 2 2 ( , ) ( , ) ( , ) ( ) ( ) ( ) ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) 2 H L sur eff eff res es y v v v x t v x t v x t EI EI S P k v x t q x t x x x v x t v x t v x t v x t A l I l t x t x t x t        − + − − +           = − + −              (34) Eq. (34) serves as a fundamental equation for bending, buckling, and free vibration analyses of nanobeam-substrate medium systems. It is worth noting that the governing differential equation of motion for the reformulated strain-gradient nanobeam proposed by Zhang and Gao [30] can be retrieved from Eq. (34) when the surface-energy effect and the underlying substrate medium are neglected. 572 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. Considering the arbitrariness of δv0(x,t) and δU on boundary terms in Eq. (33) provides boundary conditions of the nanobeam-substrate medium system as: Boundary conditions: 5 3 0 0 0 1 5 3 3 5 2 20 0 2 3 2 0 4 2 4 20 0 0 2 4 2 2 2 ( , ) ( , ) ( , ) ( ) ( ) ( ) ( , ) ( , ) (2 ) ; ( , ) ( , ) ( , ) ( ) ( ) H L sur eff eff res v v x H L eff eff v v x t v x t v x t P EI EI S P xx x v x t v x t Al I Il x t x t v x t v x t v x t P EI EI Il x x x t     =     = − − + − −      − + +           = − − −     3 0 3 3 0 0 5 3 0 0 0 4 5 3 3 5 2 20 0 2 3 2 4 0 5 4 ( , ) ; ( ) ; ( , ) ( , ) ( , ) ( ) ( ) ( ) ( , ) ( , ) (2 ) ; ( , ) ( ) ( H eff x x H L sur eff eff res v v x L H eff v x t P EI x v x t v x t v x t P EI EI S P xx x v x t v x t Al I Il x t x t v x t P EI E x   = =  =    = − −           = − + − −      − + +        = −  2 4 3 20 0 0 62 2 2 3 ( , ) ( , ) ( , ) ) ; ( ) L H eff v eff x L x L v x t v x t v x t I Il P EI x x t x  = =       − = −           (35) 6.3. Analytical Solution for Bending Analyses When the applied axial force P and the inertia forces are omitted, the governing equation for bending analyses of nanobeam-substrate medium systems can be obtained from Eq. (34) as: 6 4 2 0 0 0 06 4 2 ( ) ( ) ( ) ( ) ( ) ( ) ( ) 0 H L sur eff eff res es y v x v x v x EI EI S k v x q x x x x      − + − + − =    (36) The general solution to Eq. (36) possesses two parts, namely; the homogeneous solution 0 ( ) h v x and the particular solution 0 ( ) p v x . 0 0 0 ( ) ( ) ( ) h p v x v x v x= + (37) The homogeneous solution 0 ( ) h v x is obtained from Eq. (36) when qy(x) = 0 while the particular solution 0 ( ) p v x depends on the transversely distributed load qy(x). It is noticed that the governing differential equilibrium equation of the beam on Kerr-type foundation and Eq. (36) are in the same form. Consequently, the following homogeneous solution 0 ( ) h v x given by Morfidis [74] and Avramidis and Morfidis [75] can be applied to the proposed nanobeam-substrate medium model. 0 1 1 2 2 3 3 4 4 5 5 6 6 ( ) ( ) ( ) ( ) ( ) ( ) ( ) h v x x C x C x C x C x C x C     = + + + + + (38) where φ1(x), φ2(x), φ3(x), φ4(x), φ5(x) and φ6(x) are the basic displacement functions expressed in Appendix A; and C1, C2, C3, C4, C5 and C6 are constants of integration obtained from imposed boundary conditions. Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 573 6.4. Analytical Solution for Buckling Analyses To determine the buckling load capacity of nanobeam-substrate medium systems, the transversely distributed load qy(x) and the inertia forces present in Eq. (34) are suppressed. Therefore, Eq. (34) becomes: 6 4 2 0 0 0 06 4 2 ( ) ( ) ( ) ( ) ( ) ( ) ( ) 0 H L sur eff eff res es v x v x v x EI EI S P k v x x x x      − + − − + =    (39) For the sake of conciseness, only the simply-supported nanobeam-substrate medium system is considered. The associated classical boundary conditions are defined as: 0 0 ( ) 0 x v x = = and 0 ( ) 0 x L v x = = (40) 4 2 0 0 4 2 0 4 2 0 0 4 2 ( ) ( ) ( ) ( ) ( ) 0 and ( ) ( ) ( ) ( ) ( ) 0 H L eff eff x H L eff eff x L v x v x M x EI EI x x v x v x M x EI EI x x = =    = − + =        = − + =     (41) As suggested by Zhang and Gao [30], the non-classical boundary conditions associated with the simply-supported nanobeam-substrate medium system are: 2 0 2 0 ( ) ( ) 0 x v x x x  =  = =  and 2 0 2 ( ) ( ) 0 x L v x x x  =  = =  (42) The following analytical displacement solution for the critical (lowest) buckling load is provided by Sae-Long et al. [38] as: 0 ( ) sin( )v x A x= (43) where A is the generalized coordinate and ψ = π/L. Substituting Eq. (43) into Eq. (39) and subsequently enforcing the non-trivial condition, the critical buckling load proposed cr P can be gained as: 6 4 2 2 ( ) ( ) H L sur eff eff res esproposed cr EI EI S k P       + + = + (44) It is worth mentioning that Eq. (44) can provide the critical buckling load obtained with either the degenerated strain-gradient model or the modified couple-stress model when the associated material length-scale parameter vanishes (lm = 0 or ls = 0). Furthermore, the classical Euler’s buckling load of the beam-Winkler foundation system can be determined from Eq. (44) when the material length-scale parameters (ls and lm) and the surface-layer parameters ( sur xx E and sur res  ) vanish. 574 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. 6.5. Analytical Solution for Free-Vibration Analyses When the free vibration response of nanobeam-substrate medium systems is of interest, the transversely distributed load qy(x,t) and the applied axial force P present in Eq. (34) are suppressed. Therefore, the governing differential equation of motion becomes: 6 4 2 0 0 0 06 4 2 2 4 6 4 2 20 0 0 0 2 2 2 4 2 2 2 ( , ) ( , ) ( , ) ( ) ( ) ( , ) ( , ) ( , ) ( , ) ( , ) 2 H L sur eff eff res es v v v x t v x t v x t EI EI S k v x t x x x v x t v x t v x t v x t A l I l t x t x t x t        − + −           = − + −              (45) For the sake of simplicity, only the simply-supported nanobeam-substrate medium system is considered in the present study. Classical and non-classical boundary conditions of the simply-supported nanobeam-substrate medium system follow those suggested by Zhang and Gao [30] as: Classical boundary conditions: 0 0 ( , ) 0 x v x t = = and 0 ( , ) 0 x L v x t = = (46) 4 2 4 20 0 0 4 2 2 2 0 4 2 4 20 0 0 4 2 2 2 ( , ) ( , ) ( , ) ( , ) ( ) ( ) 0 and ( , ) ( , ) ( , ) ( , ) ( ) ( ) 0 H L eff eff v x H L eff eff v x L v x t v x t v x t M x t EI EI Il x x x t v x t v x t v x t M x t EI EI Il x x x t   = =     = − − =           = − − =       (47) Non-classical boundary conditions: 2 0 2 0 ( , ) ( , ) 0 x v x t x t x  =  = =  and 2 0 2 ( , ) ( , ) 0 x L v x t x t x  =  = =  (48) To determine the analytical solution to Eq. (45), the method of separation of variables is called for. The following transverse displacement v0(x,t) suggested by Zhang and Gao [30] satisfies all boundary conditions of Eqs. (46) to (48). 0 1 ( , ) sin( ) n i t n n n v x t V x e    = =  (49) where ϕn = nπ/L; Vn is the Fourier coefficient; n is the vibration mode; ωn is the natural frequency associated with the nth vibration mode; and i is an imaginary number. Substituting Eq. (49) into Eq. (45), the natural frequency proposed n  of the nth vibration mode is obtained by enforcing the non-trivial condition (Vn ≠ 0). The resulting explicit form of the nth-vibration mode natural frequency is: 2 4 6 2 2 2 4 2 ( ) ( ) )( ( ) 1 2 sur L H res n eff n eff n v n s v n n eproposed n S EI EI A l I k l            ++ + = + + + (50) Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 575 The natural frequency proposed n  of Eq. (50) becomes identical to the natural frequency classical n  obtained with the local beam-Winkler foundation system when the material length-scale parameters (ls, lm, and lv) and the surface-layer parameters ( sur xx E and sur res  ) vanish. Furthermore, Eq. (50) reduces to the explicit form of the natural frequency obtained with the reformulated strain-gradient beam model of Zhang and Gao [30] when the underlying substrate medium and surface energy are ignored (kes = sur xx E = sur res  = 0). 7. NUMERICAL SIMULATIONS To assess bending, buckling, and free-vibration responses of nanobeam-substrate medium systems, three numerical simulations using the proposed model are presented. The material properties of the lead nanobeam employed in all simulations are obtained from Cuenot et al. [76] and Wan et al. [77] as summarized in Table 1. Table 1 Material properties of the lead nanobeam [76-77]. Material Bulk modulus (GPa) Poisson’s ratio Density (kg/m3) Surface elastic modulus (nN/nm) Residual surface stress (nN/nm) Lead 16 0.4 11,343 8 0.63 7.1. Simulation I: Bending Responses Simulation I focuses on bending responses of the proposed model and considers two analysis cases. The first analysis case demonstrates the ability of the proposed model to obtain a rational response of the cantilever nanobeam under the pure-bending state while the second analysis case assesses the influence of system parameters on the bending response of cantilever nanobeam-substrate medium systems. 7.1.1. Analysis case I Fig. 3 shows a lead cantilever nanobeam of length L = 1,000 nm under an end moment M0 of 3,000 nN-nm, thus inducing the pure-bending state along the whole length. The beam-section shape is square with a dimension h of 100 nm, thus A = 10x103 nm2, I = 8.333x106 nm4, IГ = 6.667x10 5 nm3, and SГ = 200 nm. Only the material small-scale effect is considered in this analysis case. Therefore, the system parameters associated with the underlying substrate medium and the surface energy vanish (kes = sur xx E = sur res  = 0). This cantilever nanobeam is used to demonstrate the capability of the proposed nanobeam model in producing the small-scale dependent response under the pure-bending state. It is worth remarking that the Eringen nonlocal beam model of Limkatanyu et al. [11] and the simplified strain-gradient beam model of Sae-Long et al. [38] fail to represent such a small- scale dependent response under the pure-bending state. Therefore, the cantilever nanobeam of Fig. 3 is also analysed by these two beam models to confirm the ability of the proposed model to suppress the paradoxical behaviour. The material length-scale parameters associated with the strain-gradient effect and the couple-stress effect are set to be identical (ls = lm) and are varied through the following relations: 576 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. s s h l  = and m m h l  = (51) where ξs and ξm are, respectively, the normalized material length-scale parameters associated with the strain-gradient effect and the couple-stress effect and range from 0 to 2. 1, 000L nm= ( )0,y v x x 0 3, 000M nN nm= − z y Beam Section h h 4 3 2 2 ; ; ; 2 12 3 h h A h I I S h   = = = = Fig. 3 Simulation I: Bending analysis case I The transverse-displacement responses associated with all beam models are superimposed in Fig. 4. Obviously, the transverse-displacement response obtained with the local beam model is identical to those obtained with the Eringen nonlocal beam model of Limkatanyu et al. [11] and the simplified strain-gradient beam model of Sae- Long et al. [38] regardless of their values of material length-scale parameters. This observation confirms the inability of these two nanobeam models to represent the material small-scale dependency under the pure-bending state as observed in Limkatanyu et al. [42]. The peculiar transverse-displacement response obtained with the Eringen nonlocal beam model of Limkatanyu et al. [11] stems from the fact that adoption of the Eringen nonlocal differential form leads to an ill-posed structural-mechanics problem [34, 35, 43]. Employment of the simplified strain-gradient model leads to a well-posed structural-mechanics problem [38] and is capable of remedying the well-known paradoxical response of a cantilever beam under an end load [42] but fails to represent the small-scale dependent response under the pure-bending state. This failure relies on the fact that only the axial-strain gradient (∂εxx/∂x) along the beam length is accounted for in the simplified strain-gradient beam model of Sae-Long et al. [38]. Nevertheless, the pure-bending state induces the constant axial-strain variation along the beam length, thus vanishing the axial-strain gradient (∂εxx/∂x = 0). In opposition, the proposed nanobeam model is able to represent the material small-scale dependency under pure bending and the obtained transverse-displacement responses become stiffer with enlarging values of material length-scale parameters. This observation is confirmed with both experimental evidence [5, 78] and analytical results [35, 40-43] available in literature. Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 577 -12 -9 -6 -3 0 0 200 400 600 800 1,000 Distance along the nanobeam (nm) T r a n sv e r se d is p la c e m e n t (n m ) Proposed model ( )0s m = = Proposed model ( )0.5s m = = Proposed model ( )1s m = = Proposed model ( )2s m = = Simplified strain gradient model ( )0For all e a Eringen nonlocal model ( )0For all e a Local model Proposed model ( )0.25s m = = Fig. 4 Transverse displacement versus distance along nanobeam of Simulation I for analysis case I 7.1.2. Analysis case II Fig. 5 shows a lead cantilever nanobeam-substrate medium system under an end moment M0 of 30,000 nN-nm. Geometric and material properties of this nanobeam follow those employed in analysis case I. To assess the influence of underlying substrate media and surface-free energy on the bending response, the nanobeam-substrate medium system of Fig. 5 is analyzed with various values of the substrate-medium stiffness kes and slenderness ratio L/h. Three different beam-bulk constitutive models are employed in this analysis case to address the material small-scale effect, namely: (a) the reformulated strain-gradient model; (b) the degenerated strain-gradient model; and (c) the modified couple-stress model. It is worth mentioning that the reformulated strain-gradient model embraces the degenerated strain-gradient model (ls ≠ 0 and lm = 0) as well as the modified couple-stress model (lm ≠ 0 and ls = 0). z y Beam Section h h 4 3 2 2 ; ; ; 2 12 3 h h A h I I S h   = = = = 1, 000L nm= ( )0,y v x x 0 30, 000M nN nm= − Fig. 5 Simulation I: Bending analysis case II 578 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. Fig. 6(a) compares the variation of the normalized end displacement / Local v End End v v = with the non-dimensional substrate-medium stiffness esK  = kesL 4/(ExxI) obtained from all nanobeam models. The non-dimensional substrate-medium stiffness es K  varies from 1 to 100 while the beam depth h is kept constant at 100 nm. This range of the non- dimensional substrate-medium stiffness esK  follows that employed by Demir et al. [18]. The normalized material-length scale parameters are kept constant at ξs = 1 and ξm = 1. Clearly, the material small-scale effect and the surface-energy effect consistently stiffen the system response (ξv < 1) especially for lower values of non-dimensional substrate- medium stiffness (softer substrate media). Comparison among three nanobeam models indicates that the reformulated strain-gradient beam model leads to the stiffest system response. This superiority of the system stiffness enhancement relies on the fact that both strain-gradient and couple-stress effects are addressed in this beam model and can obviously be observed in Eq. (34). Furthermore, Fig. 6(a) points out that the system stiffness enhancement associated with the couple-stress effect is more pronounced than that associated with the strain-gradient effect. Fig. 6(b) shows the variation of the normalized end displacement / Local v End End v v = with the slenderness ratio L/h obtained from all nanobeam models. Different values of the slenderness ratio L/h are employed to reflect the surface-energy effect as well as the material small-scale effect, ranging from 1 to 100 (very stubby to very slender beams). A specific value of slenderness ratio L/h can be earned by keeping L = 1,000 nm and varying h. As provided by Liew et al. [79], a stiffness coefficient of Kes = 95x10 -3 nN/nm3 is selected to represent the underlying substrate medium as polymer. The material length- scale parameters are kept constant at ls = 100 nm and lm = 100 nm. Fig. 6(b) indicates that with increasing beam slenderness ratio L/h, the system stiffness enhancement induced by the material small-scale effect and the surface-energy effect becomes more pronounced, thus lowering the normalized end displacement / Local v End End v v = . It is worth pointing out that the normalized material-length scale parameters (ξs and ξm) become larger with increasing beam slenderness ratio L/h (decreasing beam depth h), thus magnifying the stiffening characteristic induced by the material small-scale effect. Furthermore, increasing beam slenderness ratio L/h enlarges the beam surface-area/section-area ratio (AS/AB), hence amplifying the stiffness enhancement induced by the surface-energy effect. Among three nanobeam models, the reformulated strain-gradient beam model consistently provides the largest system stiffness. However, the beam depth h dictates the system stiffness obtained with the degenerated strain-gradient model and the modified couple-stress model. As shown in the inset of Fig. 6(b), the response obtained with the degenerated strain-gradient model is stiffer than that obtained with the modified couple- stress model when the beam depth h is larger than 100 nm (L/h < 10). As the beam becomes thinner (h ≤ 100 nm), the system response associated with the modified couple- stress model is stiffer than that associated with the degenerated strain-gradient model. This observation relies on the fact that the higher-order flexural rigidity ( ) H eff EI of the degenerated strain-gradient model provides more contribution to the system stiffness when the beam becomes thicker (h > 100 nm). Zhang et al. [80] also noticed a similar observation for the reformulated strain-gradient Kirchhoff plate model. Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 579 0 0.1 0.2 0.3 0.4 0.5 0.6 0 20 40 60 80 100 Non-Dimensional Substrate-Medium Stiffness es K Higher Substrate-Medium Stiffness es k 100h nm= Reformulated strain-gradient model ( )1s m = = Degenerated strain-gradient model Modified couple-stress model ( )1s = ( )1m = v ξ 0 0.2 0.4 0.6 0.8 1 0 20 40 60 80 100 Slenderness ratio /L h Smaller Beam Depth h 3 3 95 10 / es K x nN nm − = 0 0.2 0.4 0.6 0.8 1 0 20 40 60 80 100 v /v c la ss ic a l Slenderness Ratio Proposed model SGET model MCST model Reformulated strain-gradient model ( )100s ml l nm= = Degenerated strain-gradient model Modified couple-stress model ( )100sl nm= ( )100ml nm= v ξ (a) (b) Fig. 6 Normalized end displacements with variation of: (a) non-dimensional substrate- medium stiffness esK  ; and (b) slenderness ratio L/h 7.2. Simulation II: Buckling Analysis The second simulation assesses the influence of system parameters on the critical buckling load Pcr of the simply-supported lead nanobeam-substrate medium system shown in Fig. 7. Investigated system parameters include the beam depth h, the strain-gradient material parameter ls, and the substrate-medium stiffness kes. The couple-stress material parameter lm is held constant at 100 nm. The normalized beam-depth parameter h/lm is employed to vary the beam depth h and ranges from 0.2 to 4; the normalized material length-scale parameter ls/lm is employed to examine the strain-gradient effect and varies from 0.2 to 1.2; and the non- dimensional substrate-medium stiffness esK = kesL 4/ ( ) L eff EI is employed to vary the substrate- medium stiffness kes and ranges from 1 to 10. 1, 000L nm= ( )0,y v x x es k 4 3 2 2 ; ; ; 2 12 3 h h A h I I S h   = = = = z y Beam Section h hcrP Fig. 7 Simulation II: Buckling analysis The variation of the normalized buckling load / proposed classical cr cr cr P P = with all aforementioned system parameters is presented in Fig. 8. The “classical” buckling load classical cr P is obtained with the beam-substrate medium model in which the material small- scale and surface-energy effects are all neglected. Clearly, decreasing beam depth h and increasing strain-gradient material parameter ls both magnify the normalized buckling load cr  , especially for lower values of non-dimensional substrate-medium stiffness esK (softer substrate media). However, different degrees of magnification are observed for 580 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. these two parameters. The buckling-load magnification associated with decreasing beam depth h is more pronounced than that associated with increasing strain-gradient material parameter ls. It is worth pointing out that decreasing beam depth h results in a larger beam surface-area/section-area ratio (AS/AB) and renders the beam thickness smaller when compared to the couple-stress material parameter lm. The first one amplifies the system stiffening phenomenon induced by the surface-energy effect while the second one magnifies the system stiffening phenomenon associated with the couple-stress effect. 0 0 .2 0 .4 0 .6 0 .8 1 1 .2 0 10 20 30 40 50 60 70 80 0 .2 0 .4 0 .6 0 .8 1 1 .5 2 4 70-80 60-70 50-60 40-50 30-40 20-30 10-20 0-10 1 es K = ξ c r 0 0 .2 0 .4 0 .6 0 .8 1 1 .2 0 10 20 30 40 50 60 70 80 0 .2 0 .4 0 .6 0 .8 1 1 .5 2 4 70-80 60-70 50-60 40-50 30-40 20-30 10-20 0-10 4 es K = ξ c r (a) (b) 0 0 .2 0 .4 0 .6 0 .8 1 1 .2 0 10 20 30 40 50 60 70 80 0 .2 0 .4 0 .6 0 .8 1 1 .5 2 4 70-80 60-70 50-60 40-50 30-40 20-30 10-20 0-10 7 es K = ξ c r 0 0 .2 0 .4 0 .6 0 .8 1 1 .2 0 10 20 30 40 50 60 70 80 0 .2 0 .4 0 .6 0 .8 1 1 .5 2 4 70-80 60-70 50-60 40-50 30-40 20-30 10-20 0-10 10 es K = ξ c r (c) (d) Fig. 8 Normalized buckling load ξcr with variation of normalized parameters ls/lm and h/ lm: (a) esK = 1, (b) esK = 4, (c) esK = 7, and (d) esK = 10 7.3. Simulation III: Free Vibration Analysis The third simulation performs parametric studies on free-vibration responses of the simply-supported lead nanobeam-substrate medium system shown in Fig. 9. The material and geometric properties of this nanobeam follow those employed in the first and second simulations. The mass density of the lead nanobeam is 11,343 kg/m3 as provided in Table 1. Investigated system parameters include the beam depth h, the velocity-gradient material parameter lv, and the substrate-medium stiffness kes. The material length-scale parameters are kept constant at ls = 100 nm and lm = 100 nm. The slenderness ratio L/h is employed to vary the beam depth h and ranges from 5 to 100; the normalized velocity- Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 581 gradient material parameter lv/lm is employed to examine the velocity-gradient effect and varies from 1 to 20; and the non-dimensional substrate-medium stiffness esK = kesL 4/ ( ) L eff EI is employed to vary the substrate-medium stiffness kes and ranges from 1 to 10. A given value of slenderness ratio L/h can be obtained by keeping L = 1,000 nm and varying h. It is worth remarking that the system size-dependency associated with the material small- scale effect and the surface-energy effect is reflected through the variation of the slenderness ratio L/h. 1, 000L nm= ( )0,y v x x es k 4 3 2 2 ; ; ; 2 12 3 h h A h I I S h   = = = = z y Beam Section h h Fig. 9 Simulation III: Free vibration analysis Fig. 10 shows variations of the normalized natural frequency 1 1 1 / proposed classical  = for the first vibration mode (n = 1) with normalized velocity-gradient material parameter lv/lm and slenderness ratio L/h for various values of the non-dimensional substrate- medium stiffness es K . The “classical” natural frequency 1 classical  for the first vibration mode is obtained with the beam-substrate medium model in which material small-scale and surface-energy effects are both neglected. With decreasing beam depth h (increasing slenderness ratio L/h), the normalized natural frequency 1   becomes larger, particularly for lower values of non-dimensional substrate-medium stiffness es K (softer substrate media) as shown in Fig.10. It is worth mentioning that decreasing beam depth h leads to a larger beam surface-area/section-area ratio (AS/AB) and renders the beam thickness smaller compared to the material length-scale parameters (ls = 100 nm and lm = 100 nm). The first one magnifies the system stiffening phenomenon induced by the surface-energy effect while the second one amplifies the system stiffening phenomenon associated with the material small-scale effect, thus resulting in a higher natural frequency 1 proposed  . In opposition, inclusion of the velocity-gradient effect lowers the system stiffness. The normalized natural frequency 1   reduces with increasing normalized velocity-gradient material parameter lv/lm, especially for lower values of non-dimensional substrate- medium stiffness esK (softer substrate media). This observation is consistent with that made by Zhang and Gao [30] and Yin et al. [54]. However, the system stiffening phenomenon associated with the material small-scale effect and the surface-energy effect surpasses the system weakening phenomenon associated with the velocity-gradient effect. Therefore, the combination of the material small-scale effect, the surface-energy effect, and the velocity-gradient effect results in a higher system natural frequency ( 1   >1.0) for specific values of system parameters investigated in the current study. 582 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. 5 2 0 4 0 6 0 8 0 1 0 0 0 1 2 3 4 5 6 7 8 9 10 1 5 1 0 1 5 2 0 9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2 0-1 1 es K = ω ξ 1 5 2 0 4 0 6 0 8 0 1 0 0 0 1 2 3 4 5 6 7 8 9 10 1 5 1 0 1 5 2 0 9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2 0-1 4 es K = ω ξ 1 (a) (b) 5 2 0 4 0 6 0 8 0 1 0 0 0 1 2 3 4 5 6 7 8 9 10 1 5 1 0 1 5 2 0 9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2 0-1 7 es K = ω ξ 1 5 2 0 4 0 6 0 8 0 1 0 0 0 1 2 3 4 5 6 7 8 9 10 1 5 1 0 1 5 2 0 9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2 0-1 10 es K = ω ξ 1 (c) (d) Fig. 10 Normalized natural frequency 1   with variation of normalized parameter lv/lm and slenderness ratio L/h: (a) es K = 1, (b) es K = 4, (c) es K = 7, and (d) es K = 10 8. SUMMARY AND CONCLUSIONS A rational beam-substrate medium model for bending, buckling, and free vibration analyses of nanobeams resting on elastic substrate media is proposed in the present work. The reformulated strain-gradient theory is employed to consider the beam-bulk material small- scale effect. The strain-gradient and couple-stress effects are included in the strain energy density expression while the velocity-gradient effect is included in the kinetic energy expression. The Gurtin-Murdoch surface elasticity theory is used to address the surface- energy effect. The interaction between the nanobeam and its underlying substrate medium is considered through the Winkler-foundation model. To obtain the system governing equation and corresponding boundary conditions, Hamilton’s principle is called for. Three numerical simulations are presented to characterize the influences of the material small-scale effect, the surface-energy effect, and the surrounding substrate medium on bending, buckling, and free vibration responses of nanobeam-substrate medium systems. The first simulation possesses two analysis cases and emphasizes on bending analyses of nanobeam-substrate medium systems. The first one presents the ability of the proposed model to address the size dependency under the pure-bending state and shows the stiffening bending response associated with the material small-scale effect. The second Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 583 one shows that the system-stiffness enhancement associated with the material small-scale and surface-energy effects is more pronounced when the underlying substrate medium becomes softer and the beam becomes more slender. The second simulation assesses the combined impacts of material small-scale, surface energy, and substrate media on the critical buckling load of a simply-supported nanobeam- elastic substrate system. It is indicated that the critical buckling load is enhanced by two distinct stiffening phenomena. The first one is associated with the material small-scale effect while the second one is related to the surface-energy effect. However, these two stiffening phenomena are less pronounced when the underlying substrate medium becomes softer. The third simulation investigates the combined effects of material small-scale, surface energy, and substrate media on the fundamental frequency of a simply-supported nanobeam-elastic substrate system. It is found that the material small-scale and surface- energy effects enhance the system stiffness while the velocity-gradient effect lowers the system stiffness, particularly for softer substrate media. However, the system stiffening phenomenon associated with the material small-scale effect and the surface-energy effect surpasses the system weakening phenomenon associated with the velocity-gradient effect. Therefore, the combination of the material small-scale effect, the surface-energy effect, and the velocity-gradient effect results in a higher system natural frequency ( 1   >1.0) for specific values of system parameters investigated in the current study. The next step in developing the proposed beam model is considering the geometric nonlinearity for analyses of the nanobeam system with large displacement. Moreover, the material nonlinearity is included into the model and applying the model to assess the rational micro/nanobeam structures against the failure scenarios in the future. Acknowledgement: This study was financially supported by Office of the Permanent Secretary, Ministry of Higher Education, Science, Research and Innovation under Research Grant for New Scholar (RGNS 64-134) and by the Thailand Research Fund (TRF) under Grant (RTA6280012). Furthermore, special thanks go to a senior lecturer Mr. Wiwat Sutiwipakorn for reviewing and correcting the English of this paper. REFERENCES 1. Iijima, S., 1991, Helical microtubules of graphitic carbon, Nature, 354, pp. 56-58. 2. Hornyak, G.L., Moore, J.J., Tibbals, H.F., Dutta, J., 2008, Fundamentals of nanotechnology, Taylor & Francis CRC Press, Boca Raton, USA. 3. Abouelregal, A.E., Mohammad-Sedighi, H., Faghidian, S.A., Shirazi, A.H., 2021, Temperature-dependent physical characteristics of the rotating nonlocal nanobeams subject to a varying heat source and a dynamic load, Facta Universitatis Series Mechanical Engineering, 19(4), pp. 633-656. 4. Uzun, B., Civalek, O., Aydogdu, I., 2021, Optimum design of nano-scaled beam using the social spider optimization (SSO) algorithm, Journal of Applied and Computational Mechanics, 7(3), pp. 1348-1361. 5. 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(8), pp. 1477-1508. 6. Jung, W.-Y., Han, S.-C., 2015, Static and eigenvalue problems of Sigmoid Functionally Graded Materials (S-FGM) micro-scale plates using the modified couple stress theory, Applied Mathematical Modelling, 39(12), pp. 3506-3524. 7. Zhu, Z., Wang, Z., Zhou, Y., Chen, Y., Zhou, L., She, A., 2021, Evaluation of the nanostructure of calcium silicate hydrate based on atomic force microscopy-infrared spectroscopy experiments, Nanotechnology Reviews, 10(1), pp. 807-818. 8. Senturia, S.D., 2001, Microsystem design, Springer, Boston, MA, USA. 584 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. 9. Limkatanyu, S., Prachasaree, W., Damrongwiriyanupap, N., Kwon, M., 2014, Exact stiffness matrix for nonlocal bars embedded in elastic foundation media: The virtual-force approach, Journal of Engineering Mathematics, 89, pp. 163-176. 10. Ghodrati, B., Yaghootian, A., Zadeh, A.G., Mohammad-Sedighi, H., 2018, Lamb wave extraction of dispersion curves in micro/nano-plates using couple stress theories, Waves in Random and Complex Media, 28(1), pp. 15-34. 11. Limkatanyu, S., Sae-Long, W., Horpibulsuk, S., Prachasaree, W., Damrongwiriyanupap, N., 2018, Flexural responses of nanobeams with coupled effects of nonlocality and surface energy, ZAMM Journal of Applied Mathematics and Mechanics, 98(10), pp. 1771-1793. 12. Barretta, R., Canadija, M., Marotti de Sciarra, F., 2020, Nonlocal mechanical behavior of layered nanobeams, Symmetry, 12(5), 717. 13. Pourabdy, M., Shishesaz, M., Shahrooi, S., Roknizadeh, S.A.S., 2021, Analysis of axisymmetric vibration of functionally-graded circular nano-plate based on the integral form of the strain gradient model, Journal of Applied and Computational Mechanics, 7(4), pp. 2196-2220. 14. Müller, P., Saúl, A., 2004, Elastic effects on surface physics, Surface Science Reports, 54(5-8), pp. 157-258. 15. Shenoy, V.B., 2005, Atomistic calculations of elastic properties of metallic fcc crystal surfaces, Physical Review B, 71, 094104. 16. Wang, S., Yao, Y., Zhang, B., 2021, Atomistic simulations of the graded residual elastic fields in metallic nanowires, Results in Physics, 25, 104272. 17. Zhang, Y., Cao, P., Deng, B., Huang, L., Shi, Y., 2021, Strain rate-dependent tensile response of glassy silicon nanowires studied by accelerated atomistic simulations, AIP Journal of Applied Physics, 130(8), 085105. 18. Demir, C., Mercan, K., Numanoglu, H.M., Civalek, O., 2018, Bending response of nanobeams resting on elastic foundation, Journal of Applied and Computational Mechanics, 4(2), pp. 105-114. 19. Sedighi, H.M., 2020, Divergence and flutter instability of magneto-thermo-elastic C-BN hetero-nanotubes conveying fluid, Acta Mechanica Sinica, 36, pp. 381-396. 20. Malikan, M., Eremeyev, V.A., Sedighi, H.M., 2020, Buckling analysis of a non-concentric double-walled carbon nanotube, Acta Mechanica, 231, pp. 5007-5020. 21. Pereira, R.S., 2001, Atomic force microscopy as a novel pharmacological tool, Biochemical Pharmacology, 62(8), pp. 975-983. 22. Qi, D., Liu, Z., Leow, W.R., Chen, X., 2017, Elastic substrates for stretchable devices, MRS Bulletin, 42(2), pp. 103-107. 23. Jena, S.K., Chakraverty, S., Malikan, M., Tornabene, F., 2021, Stability analysis of single-walled carbon nanotubes embedded in winkler foundation placed in a thermal environment considering the surface effect using a new refined beam theory, Mechanics Based Design of Structures and Machines, 49(4), pp. 581-595. 24. Qi, D., Zhang, K., Tian, G., Jiang, B., Huang, Y., 2021, Stretchable electronics based on PDMS substrates, Advanced Materials, 33(6), 2003155. 25. Eringen, A.C., 1972, Nonlocal polar elastic continua, International Journal of Engineering Science, 10(1), pp. 1-16. 26. Eringen, A.C., 1983, On differential equations of nonlocal elasticity and solutions of screw dislocation and surface waves, AIP Journal of Applied Physics, 54, pp. 4703-4710. 27. Mindlin, R.D., 1964, Micro-structure in linear elasticity, Archive for Rational Mechanics and Analysis, 16, pp. 51-78. 28. Mindlin, R.D., 1965, On the equations of elastic materials with micro-structure, International Journal of Solids and Structures, 1(1), pp. 73-78. 29. Yang, F., Chong, A.C.M., Lam, D.C.C., Tong, P., 2002, Couple stress based strain gradient theory for elasticity, International Journal of Solids and Structures, 39(10), pp. 2731-2743. 30. Zhang, G.Y., Gao, X.-L., 2020, A new Bernoulli–Euler beam model based on a reformulated strain gradient elasticity theory, Mathematics and Mechanics of Solids, 25(3), pp. 630-643. 31. Cammarata, R.C., 1994, Surface and interface stress effects in thin films, Progress in Surface Science, 46(1), pp. 1-38. 32. Gurtin, M.E., Murdoch, A.I., 1975, A continuum theory of elastic material surface, Archive for Rational Mechanics and Analysis, 57, pp. 291-323. 33. Gurtin, M.E., Murdoch, A.I., 1978, Surface stress in solids, International Journal of Solids and Structures, 14(6), pp. 431-440. 34. Peddieson, J., Buchanan, G.R., McNitt, R.P., 2003, Application of nonlocal continuum models to nanotechnology, International Journal of Engineering Science, 41(3-5), pp. 305-312. 35. Challamel, N., Wang, C.M., 2008, The small length scale effect for a non-local cantilever beam: A paradox solved, Nanotechnology, 19(34), 345703. Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 585 36. Juntarasaid, C., Pulngern, T., Chucheepsakul, S., 2012, Bending and buckling of nanowires including the effects of surface stress and nonlocal elasticity, Physica E: Low-dimensional Systems and Nanostructures, 46, pp. 68-76. 37. Nikam, R.D., Sayyad, A.S., 2020, A unified nonlocal formulation for bending, buckling and free vibration analysis of nanobeams, Mechanics of Advanced Materials and Structures, 27(10), pp. 807-815. 38. Sae-Long, W., Limkatanyu, S., Rungamornrat, J., Prachasaree, W., Sukontasukkul, P., Sedighi, H.M., 2021, A rational beam-elastic substrate model with incorporation of beam-bulk nonlocality and surface-free energy, The European Physical Journal Plus, 136, 80. 39. Barretta, R., Faghidian, S.A., Marotti de Sciarra, F., Pinnola, F.P., 2021, Timoshenko nonlocal strain gradient nanobeams: Variational consistency, exact solutions and carbon nanotube Young moduli, Mechanics of Advanced Materials and Structures, 28(15), pp. 1523-1536. 40. Li, C., Yao, L., Chen, W., Li, S., 2015, Comments on nonlocal effects in nano-cantilever beams, International Journal of Engineering Science, 87, pp. 47-57. 41. Fernández-Sáez, J., Zaera, R., Loya, J.A., Reddy, J.N., 2016, Bending of Euler–Bernoulli beams using Eringen’s integral formulation: A paradox resolved, International Journal of Engineering Science, 99, pp. 107-116. 42. Limkatanyu, S., Sae-Long, W., Sedighi, H.M., Rungamornrat, J., Sukontasukkul, P., Zhang, H., Chindaprasirt, P., 2022, Flexibility-based stress-driven nonlocal frame element: Formulation and applications, Engineering with Computers. 43. Romano, G., Barretta, R., Diaco, M., Marotti de Sciarra, F., 2017, Constitutive boundary conditions and paradoxes in nonlocal elastic nanobeams, International Journal of Mechanical Sciences, 121, pp. 151-156. 44. Koutsoumaris, C.C., Eptaimeros, K.G., Zisis, T., Tsamasphyros, G.J., 2016, A straightforward approach to Eringen’s nonlocal elasticity stress model and applications for nanobeams, AIP Conference Proceedings, 1790(1), 150018. 45. Lu, L., Guo, X., Zhao, J., 2017, A unified nonlocal strain gradient model for nanobeams and the importance of higher order terms, International Journal of Engineering Science, 119, pp. 265-277. 46. Chen, W., Wang, L., Dai, H.-L., 2019, Nonlinear free vibration of nanobeams based on nonlocal strain gradient theory with the consideration of thickness-dependent size effect, Journal of Mechanics of Materials and Structures, 14(1), pp. 119-137. 47. Bian, P.‐L., Qing, H., 2021, On bending consistency of Timoshenko beam using differential and integral nonlocal strain gradient models, ZAMM Journal of Applied Mathematics and Mechanics, 101(8), e202000132. 48. Merzouki, T., Houari, M.S.A., Haboussi, M., Bessaim, A., Ganapathi, M., 2022, Nonlocal strain gradient finite element analysis of nanobeams using two-variable trigonometric shear deformation theory, Engineering with Computers, 38, pp. 647-665. 49. Zhao, X., Zheng, S., Li, Z., 2022, Bending, free vibration and buckling analyses of AFG flexoelectric nanobeams based on the strain gradient theory, Mechanics of Advanced Materials and Structures, 29(4), pp. 548-563. 50. Polizzotto, C., 2017, A hierarchy of simplified constitutive models within isotropic strain gradient elasticity, European Journal of Mechanics - A/Solids, 61, pp. 92-109. 51. Altan, B.S., Aifantis, E.C., 1997, On some aspects in the special theory of gradient elasticity, Journal of the Mechanical Behavior of Materials, 8(3), pp. 231-282. 52. Barretta, R., Marotti de Sciarra, F., 2013, A nonlocal model for carbon nanotubes under axial loads, Advances in Materials Science and Engineering, 360935. 53. Sidhardh, S., Ray, M.C., 2019, Exact solutions for elastic response in micro- and nano-beams considering strain gradient elasticity, Mathematics and Mechanics of Solids, 24(4), pp. 895-918. 54. Yin, S., Xiao, Z., Deng, Y., Zhang, G., Liu, J., Gu, S., 2021, Isogeometric analysis of size-dependent Bernoulli– Euler beam based on a reformulated strain gradient elasticity theory, Computers & Structures, 253, 106577. 55. Papargyri-Beskou, S., Polyzos, D., Beskos, D.E., 2009, Wave dispersion in gradient elastic solids and structures: A unified treatment, International Journal of Solids and Structures, 46(21), pp. 3751-3759. 56. Sapsathiarn, Y., Rajapakse, R.K.N.D., 2018, Mechanistic models for nanobeams with surface stress effects, Journal of Engineering Mechanics, 144(11), 04018098. 57. Gao, X.-L., Mahmoud, F.F., 2014, A new Bernoulli–Euler beam model incorporating microstructure and surface energy effects, Zeitschrift für angewandte Mathematik und Physik, 65, pp. 393-404. 58. Ahouel, M., Houari, M.S.A., Bedia, E.A.A., Tounsi, A., 2016, Size-dependent mechanical behavior of functionally graded trigonometric shear deformable nanobeams including neutral surface position concept, Steel and Composite Structures, 20(5), pp. 963-981. 59. Sae-Long, W., Limkatanyu, S., Prachasaree, W., Rungamornrat, J., Sukontasukkul, P., 2020, A thermodynamics-based nonlocal bar-elastic substrate model with inclusion of surface-energy effect, Journal of Nanomaterials, 8276745. 586 S. LIMKATANYU, W. SAE-LONG, J. RUNGAMORNRAT, ET AL. 60. Foroutan, S., Haghshenas, A., Hashemian, M., Eftekhari, S.A., Toghraie, D., 2018, Spatial buckling analysis of current-carrying nanowires in the presence of a longitudinal magnetic field accounting for both surface and nonlocal effects, Physica E: Low-dimensional Systems and Nanostructures, 97, pp. 191-205. 61. Hosseini-Hashemi, S., Nazemnezhad, R., Rokni, H., 2015, Nonlocal nonlinear free vibration of nanobeams with surface effects, European Journal of Mechanics - A/Solids, 52, pp. 44-53. 62. Yin, F., Chen, C.P., Chen, D.L., 2015, Vibration analysis of nano-beam with consideration of surface effects and damage effects, Nonlinear Engineering, 4(1), pp. 61-66. 63. Limkatanyu, S., Damrongwiriyanupap, N., Prachasaree, W., Sae-Long, W., 2013, Modeling of axially loaded nanowires embedded in elastic substrate media with inclusion of nonlocal and surface effects, Journal of Nanomaterials, 635428. 64. Azizi, S., Safaei, B., Fattahi, A.M., Tekere, M., 2015, Nonlinear vibrational analysis of nanobeams embedded in an elastic medium including surface stress effects, Advances in Materials Science and Engineering, 318539. 65. Niknam, H., Aghdam, M.M., 2015, A semi analytical approach for large amplitude free vibration and buckling of nonlocal FG beams resting on elastic foundation, Composite Structures, 119, pp. 452-462. 66. Demir, C., 2016, Nonlocal vibration analysis for micro/nano beam on Winkler foundation via DTM, International Journal of Engineering & Applied Sciences, 8(4), pp. 108-118. 67. Ponbunyanon, P., Limkatanyu, S., Kaewjuea, W., Prachasaree, W., Chub-Uppakarn, T., 2016, A novel beam- elastic substrate model with inclusion of nonlocal elasticity and surface energy effects, Arabian Journal for Science and Engineering, 41, pp. 4099-4113. 68. Jena, S.K., Chakraverty, S., Tornabene, F., 2019, Dynamical behavior of nanobeam embedded in constant, linear, parabolic, and sinusoidal types of Winkler elastic foundation using first-Order nonlocal strain gradient model, Materials Research Express, 6(8), 0850f2. 69. Wolfram, S., 1992, Mathematica reference guide, Addison-Wesley Publishing Company, Redwood City, California, USA. 70. Bauchau, O.A., Craig, J.I., 2009, Structural analysis, Springer, Dordrecht, Netherlands. 71. Fleck, N.A., Hutchinson, J.W., 1997, Strain gradient plasticity, Advances in Applied Mechanics, 33, pp. 295-361. 72. Miller, R.E., Shenoy, V.B., 2000, Size-dependent elastic properties of nanosized structural elements, Nanotechnology, 11(3), pp. 139-147. 73. Winkler, E., 1867, Die Lehre von der und Festigkeit, Dominicus, Prague, Czechoslovakia. 74. Morfidis, K., 2003, Research and development of methods for the modeling of foundation structural elements and soil, PhD dissertation, Department of Civil Engineering, Aristotle University of Thessaloniki. 75. Avramidis, I.E., Morfidis, K., 2006, Bending of beams on three-parameter elastic foundation, International Journal of Solids and Structures, 43(2), pp. 357-375. 76. Cuenot, S., Frétigny, C., Demoustier-Champagne, S., Nysten, B., 2004, Surface tension effect on the mechanical properties of nanomaterials measured by atomic force microscopy, Physical Review B, 69, 165410. 77. Wan, J., Fan, Y.L., Gong, D.W., Shen, S.G., Fan, X.Q., 1999, Surface relaxation and stress of fcc metals: Cu, Ag, Au, Ni, Pd, Pt, Al and Pb, Modelling and Simulation in Materials Science and Engineering, 7(2), pp. 189-206. 78. Nilsson, S.G., Sarwe, E.-L., Montelius, L., 2003, Fabrication and mechanical characterization of ultrashort nanocantilevers, AIP Applied Physics Letters, 83(5), pp. 990-993. 79. Liew, K.M., He, X.Q., Kitipornchai, S., 2006, Predicting nanovibration of multi-layered graphene sheets embedded in an elastic matrix, Acta Materialia, 54(16), pp. 4229-4236. 80. Zhang, G., Zheng, C., Mi, C., Gao, X.-L., 2021. A microstructure-dependent Kirchhoff plate model based on a reformulated strain gradient elasticity theory, Mechanics of Advanced Materials and Structures. APPENDIX A The general form of the differential equation of Eqs. (36) can be expressed as: 6 4 2 0 0 0 1 2 3 06 4 2 ( ) ( for (0, ) ) ( ) ( ) 0 v x v x v x v x x x x x L      + = + +    (A1) where the parameters λ1, λ2, and λ3 are defined as: Bending, Buckling, and Free Vibration Analyses of Nanobeam-Substrate Medium Systems 587 1 2 3 ; ; ( ) ( ) d ( ) n ( ) a L sur eff res es H H H eff eff eff EI S k EI EI EI     = = =− − (A2) For the sake of conciseness, the parameters λ1, λ2, and λ3 are rewritten in terms of the auxiliary parameters for the expression of the general solution as: 32 3 21 1 2 31 2 (2 / 27) ( / 3)( / 3) ; ; and 3 2          − +− + = =  = + (A3) The general form of homogeneous solution of Eq. (A.1) can be expressed as: 0 1 1 2 2 3 3 4 4 5 5 6 6 ( ) ( ) ( ) ( ) ( ) ( ) ( ) h v x x C x C x C x C x C x C     = + + + + + (A4) As suggested by Morfidis [74] and Avramidis and Morfidis [75], there are two possible homogeneous solutions, which depend on the sign of the parameter Δ as follows: Solution Case I: when Δ = α3 + β2 > 0 1 1 1 2 3 4 5 6 ( ) ( ) ( i ; ; cos[ ]; s n[ ]; cos[ ]; s ] ) ( ) in( ) ( ) [ R x R x Rx Rx Rx Rx e e e Qx e x Qx e x x x x xQx e Qx       − − − = = = = = = (A5) where ( ) ( ) ( ) 2 13 3 1 3 3 1 2 3 2 2 3 ; 2 2 3 ; 2 2 ; 3 2 / 3 R m n m m n m R Q m n         + + + − = =          = − +  + − −  − − +  +    = − = − −  + − + − − −  (A6) Solution Case II: when Δ = α3 + β2 < 0 1 2 1 2 1 2 3 4 5 6 ( ) ; ( ) ; ( ) ; ( ) ; ( ) ; ( ) R x R x R x R xRx Rx x e x e x e x e x e x e      − −− = = = = = = (A7) where 1 1 1 1 31 2 2 2 cos ; 2 cos ; 3 3 3 3 4 2 cos ; cos / 3 3 R R R             − +    = − − = − −        +   = − − = − −      (A8)