Nonlinear thermal stability and snap-through buckling of temperature-dependent geometrically imperfect graded nanobeams on nonlinear elastic foundation

In this research, thermal postbuckling and nonlinear thermal bending of size-dependent functionally graded (FG) perfect/imperfect nanobeams in-contact with a nonlinear elastic foundation and exposed to thermal loading are first scrutinized. The second objective of this research is to investigate snap-through instability in a thermally postbuckled FG nanobeam due to uniform lateral load. The geometrical imperfection of the nanobeam is taken into account. Thermo-mechanical material properties are temperature-dependent and are assumed to vary continuously throughout the nanobeam thickness on the basis of the power-law model. The nonlinear equilibrium equations are derived according to the Euler-Bernoulli beam theory in conjunction with the nonlinear von-Karman assumptions based on the nonlocal elasticity theory. Using Chebyshev polynomial of the first kind, Ritz method is utilized into the principle of virtual displacement to form nonlocal governing equations. Three different methodologies, including direct displacement control approach, Newton-Raphson iterative method, and cylindrical arch-length scheme are utilized to investigate the nonlinear thermal stability curves and snap-through phenomenon through limit points of a thermally postbuckled FG nanobeam. The primary purpose of this research is to explore the influences of the nonlocal parameter, nonlinear elastic foundation coefficients, power-law index, imperfection amplitude as well as the different types of boundary conditions on the nonlinear thermal stability and snap-through phenomenon of the FG nanobeam.


Introduction
Known as a novel class and microscopically inhomogeneous composite materials, functionally graded materials (FGMs) are composed of a controlled mixture of at least two constituents such that the volume fraction of each alters continuously in one or more desired direction(s). FGMs have decreased stress concentrations and thermal stresses and have the capability of resisting intense temperature gradients [1,2]. For this reason, extensive studies on thermal stability investigation of FGM structures have been performed. Hence, within the framework of the conventional elasticity, a survey of the research works on dynamic and static features of FGM structures under thermo-mechanical loads is presented hereinafter.
Ma and Lee [3] analysed the thermal postbuckling and vibration analysis of FGM Timoshenko beam exposed to in-plane thermal loading using the shooting method. They concluded that the bifurcation buckling does not happen for simply-supported FG beams by reason of the bending-stretching coupling effect. Kanani et al [4] investigated the influence of nonlinear elastic foundation and various types of boundary conditions on the free and forced vibration of graded beam with large amplitude assumption. In this study, governing equations of motion are according to Euler-Bernoulli theory and are solved using the variational iteration method (VIM). It is observed that the nonlinearity in the time and frequency responses of forced vibration enhances by increasing the nonlinear coefficient of elastic foundation. Fallah and Aghdam [5] performed the postbuckling analysis and nonlinear free vibration behavior of graded beams placed on the nonlinear elastic foundation. In this work, using the approach of He's variational iteration, nonlinear natural frequency and postbuckling load-deflection curves of FG beams subjected to axial force are analysed. By adopting the Timoshenko beam theory and generalized differential quadrature method (GDQM), Esfahani et al [6] examined an analysis on the nonlinear thermal stability of FG beams embedded on nonlinear elastic foundations with hardening nonlinearity behavior. Temperature dependency of material properties of the graded beam is also considered. It is shown that the behavior of the FG beam may be of the bifurcation type of buckling or unique and stable equilibrium curve depending on the assumed type of boundary conditions and the thermal loading. Komijani et al [7] presented the vibration response of thermo-electrically postbuckled functionally graded piezoelectric material (FGPM) beam with a rectangular cross-section and various types of boundary conditions employing the Ritz finite element method (FEM). It is observed that near the bifurcation points of the FG structures, the numerical values of the fundamental frequencies reach to zero. She et al [8] explored the thermal buckling and postbuckling of FGM beams with or without piezoelectric layers by adopting the higher-order shear deformation beam theory and the concept of physical neutral surface. Recently, Dehrouyeh-Semnani [9] developed an analytical solution on the nonlinear behavior of shear-deformable beams made of FGMs with initial curvature subjected to thermal loading. Through vibration and thermal stability investigation of the beams, it is confirmed that temperature-dependent material properties and nonlinear elastic foundation have a notable effect on the obtained numerical results.
Due to the absence of a material length scale parameter, modeling of nano/micro-structural components, like nano/micro-scale beams and plates via conventional elasticity theory is inadequate to describe a phenomenon often named size effect. In general, experimental approaches, molecular dynamics (MD) simulations, and continuum mechanics methods are often applied for examining the mechanical properties of nano/micro-scaled structural systems. Experimental observation and MD simulations have shown notable small-scaled influences on the mechanical properties of nano/micro-scaled systems. There is, however, often some difficulty in MD simulations and experimental investigations (often the MD simulations are complicated and time-consuming and micro/nano-scaled controlled experiments are difficult to implement). In addition, the capability of this approach is much limited by the astronomical amount of data produced in the calculation.
Furthermore, the structural stiffness softening is observed at nanoscale levels, whereas the behavior of microstructures is usually directed by the stiffness hardening. Accordingly, size-dependent theories, including the strain gradient and couple stress models, are often employed to analyse the mechanical response of microstructures while the nonlocal elasticity theory is utilized to nanoscale structures, including nanobeams and nanoplates [10]. To clarify this deficiency, several models incorporating size-dependency have been developed including modified couple stress theory (MCST) [11], modified strain gradient theory (MSGT) [12] and nonlocal elasticity theory [13,14]. Also, another class of size-dependent theory which is called nonlocal strain gradient theory [15] has been recently developed based on a combination of the nonlocal elasticity theory and the strain gradient theory. According to Eringen's nonlocal theory, at a given point in a continuum, the stress depends on the strain at all other adjacent points of the continuum, and so the size effect is achieved through constitutive equations utilizing a nonlocal parameter. Moreover, to obtain high sensitivity and desired efficiency, FGMs have found potential utilization in nano/micro-electro mechanical systems (NEMS/MEMS) [16,17] such as atomic force microscopes (AFMs) [18], shape memory alloy films [19,20], micro sensors, micro piezo actuator, and nano-motors [21]. Therefore, from the point of view of theoretical analysis, it is essential to establish the nano/micro-mechanical models of the system made of FGMs and to examine the mechanical properties of the system. Also, the basic elements such as the beams, plates and shells in the nano/micro-scale are commonly utilized as the components in the NEMS/MEMS. FG nanostructures in these devices could be manufactured based on the multilayer process which combines both chemical vapour deposition and highgrowth rate plasma-enhanced chemical vapour deposition bulk layer [16].
In recent years, some investigators have applied modified couple stress and modified strain gradient theories to analyse FG microstructures and to examine the influence of size-dependency on their mechanical behavior. For instance, Reddy [22] employed the modified couple stress for bending, vibration, and buckling analysis of microstructure-dependent simply-supported FG beams. It is concluded that the natural frequency anticipated by the classical continuum model is lower than that by the microstructure-dependent model. By adopting the Timoshenko theory and differential quadrature method (DQM), a numerical solution is presented by Ke et al  [23] to obtain the nonlinear vibration behavior of FG beams using the modified couple stress theory. By means of the modified strain gradient theory, Ansari et al [24] explored free vibration responses of Timoshenko FG microbeams according to the Mori-Tanaka scheme and Navier solution. Li et al [25] reported static bending and free vibration analysis of simply-supported FGP microbeams employing the modified strain gradient theory. It is revealed that natural frequency obtained by the MSGT is larger than the classical model. Karami et al [26] examined the influences of different homogenization models on the free vibration of FG curved micro-beam based on the first-order shear deformation beam theory in conjunction with the MSGT.
Simsek [27] investigated the nonlinear vibration behavior of microbeams placed on a nonlinear elastic foundation according to the modified couple stress theory and Euler-Bernoulli beam theory. Komijani et al [28] analysed the thermal postbuckling and small amplitude vibration in the pre/postbuckling configurations of FG microbeams embedded on a nonlinear elastic foundation and under in-plane thermal loads using modified couple stress theory. It is demonstrated that the nonlinear structural behavior of the FG beams may be affected by consideration of microstructure-dependency. Recently, Shenas et al [29] conducted a study on thermal postbuckling and free vibration analysis of postbuckled rotating FG pre-twisted microbeams with a mid-plane symmetric distribution of material properties using MCST. As a conclusion, they have reported that the raise of the angular speed reduces the mid-point deflection, leading to the increase in the critical buckling temperature.
Similar to microbeams, the literature on the problem of FG circular and rectangular microplates is presented. Reddy et al [30] directed a research study on an axisymmetric bending of circular graded microplates applying the modified couple stress theory with regard to the finite element method and Newton-Raphson iterative procedure. Ashoori and Sadough Vanini [31] anticipated the size-dependent free vibration of circular FG piezoelectric plates in thermal environment. In this study, it is owing to the lateral loading of thermally preloaded microplates the buckled states are considered to be emerging from either thermo-electrical bifurcation type of buckling or limit load buckling. Moreover, Kim et al [32] presented bending, free vibration, and buckling response of simply-supported porous FG microplates employing the first-order shear deformation and classical plate theories. It is concluded that the porosity distribution in the FGMs has an important effect on the bending-stretching coupling coefficients.
Lately, some outstanding researches have been performed in which the nonlocal elasticity theory has put to employ to capture size dependency in mechanical responses of nanostructures. Shahsavari et al [33] performed an investigation on dynamic deflections of viscoelastic orthotropic nanoplates subjected to the moving load embedded within visco-Pasternak substrate and hygrothermal environment. It is concluded that at high values of the nonlocal parameter, the visco-Pasternak foundation is more efficient than Pasternak and Winkler foundations. Static characteristics of FGM plates made of hexagonal beryllium crystals as an anisotropic material are analysed by Karami et al [34] based on the nonlocal strain gradient theory. It is verified that the nonlocality and strain gradient size-dependency increase leads to decrement and increment the stiffness of the square plate, respectively. Also, Karami et al [35] used a closed-form solution to obtain the resonance deflection of FG rectangular nanoplate where the size dependency of the proposed model is based on the bi-Helmholtz nonlocal strain gradient theory. Moreover, based on the quasi-3D bi-Helmholtz nonlocal strain gradient model, in which the proposed model includes three-fine scale parameters in addition to the elastic constants of the anisotropic material, the elastic bulk waves of of FG triclinic nanoplates are studied by Karami et al [36]. Niknam and Aghdam [37] introduced a closed-form solution for nonlocal large amplitude free vibration and buckling analysis of FG nanobeams on nonlinear elastic foundation. As a result, it is shown that the fully clamped nanobeams are more affected by the nonlocal parameters than other types of boundary conditions. Based on the Euler-Bernoulli and Timoshenko beam theories, Ebrahimi and Salari [38][39][40] researched the vibration and thermal buckling analysis of graded nanobeams. In these studies, linear thermal buckling and free vibration responses of the FG beams with temperature-dependent material properties in thermal environment are studied using the nonlocal theory and differential transform method (DTM). Gao et al [41] examined the nonlinear bending and thermal postbuckling behavior of FG piezoelectric nanobeams on the basis of the two-step perturbation method and nonlocal theory of elasticity. In this research, only immovable clamped types of boundary condition are examined, and the uniform transverse load is also included in the formulation. Karami et al [42] presented the numerical results of the nonlocal porous FG nanoplate model for the wave propagation analysis of rectangular graded plates with fully clamped boundary conditions. It is shown that, increasing the nonlocal parameter will decrease the wave frequencies as well as phase velocities. Karami et al [43] also used the nonlocal elasticity theory to analyse the propagation of guided waves in a functionally graded carbon nanotubereinforced composite (FG-CNTRC) plate based on the first-order shear deformation plate theory. Solution method of this research is based on the analytical solution which satisfies fully clamped boundary conditions and the material properties is supposed via the Mori-Tanaka approach. In another study, Karami et al [44] applied the second-order shear deformation theory to obtain the size-dependent static, stability and dynamic analysis of FG-CNTRC plates resting on Winkler-Pasternak foundation. It is observed that the frequency and buckling obtained by the nonlocal model are smaller than those by the classical model regard to effect of nonlocal parameter. In the category of FG nanoshells, the vibrational behavior of porous doubly-curved nanoshells with simply-supported boundary conditions and based on the nonlocal strain gradient theory is investigated by Karami et al [45]. Also, a variational approach is developed by Karami et al [46] to obtain the wave dispersion in anisotropic doubly-curved nanoshells based on the nonlocal strain gradient higher-order shell theory. Jouneghani et al [47] investigated the bending response of simply-supported porous FG nanobeams under a hygro-thermo-mechanical loading based on the Navier method and Eringen's nonlocal theory. Yang et al [48] performed research on the buckling, vibration, and nonlinear bending of bi-directional Euler-Bernoulli FG nanobeams subjected to axial load by means of the DQ method. Karami and Janghorban [49] analysed the free vibration of nonlocal strain gradient FG nanobeams resting on elastic foundation using a new shear strain shape function and Navier method. It is deduced that the length-to-thickness ratio rising leads to an increment in the natural frequencies. Recently, Aria et al [50] proposed the nonlocal finite element formulation for the linear thermal buckling and vibration analysis of FG nanobeams with porosities. It is seen that the natural frequency increases by the presence of porosity distribution in the graded nanobeam.
In general, the complexity of the snap-through type of instability depends on the geometric and loading parameters. Hence, there are only a few available papers on the topic of the snap-through phenomenon in the thermally postbuckled FGM beams and plates due to uniform lateral load. For instance, Fallah and Nosier [51] presented the snap-through analysis of thermally postbuckled FG circular plates under asymmetric transverse loading based on the conventional elasticity theory. In this study, for modeling the nonlinear stability response of circular FG plates with different types of boundary conditions, a two-parameter perturbation procedure along with Fourier series method is implemented. Also, Ashoori and Sadough Vanini [52] utilized the modified couple stress theory within the framework of the classical circular FG piezoelectric microplates to explore the thermal postbuckling and snap-through behavior of thermally preloaded microplates under concentrated and uniform lateral loads. It is observed that FGP microplates reveal the higher value of upper and lower limit loads and also more enhanced snap-through feature by greater levels of deformed state due to thermal preloading.
For curved tubes and beams without thermally postbuckled configuration and subjected to lateral mechanical loading, Babaei et al [53] studied snap-through behavior of FGM curved tubes subjected to the uniform lateral pressure by means of the conventional elasticity theory and the higher-order shear deformation tube theory. She et al [54,55] selected the nonlocal strain gradient theory to examine snap-through buckling response of porous FG curved nanotubes and curved nanobeams exposed to a uniform transverse load. In the mentioned works, governing equations that are based on the higher-order shear deformation theory are solved using the two-step perturbation method. They concluded that small size parameters have a considerable influence on the nonlinear bending analysis of curved structures. In another study, She et al [56] studied the influence of the thermal environment on the snap-buckling phenomena of temperature-dependent FG curved nanobeams subjected to uniform transverse load and based on the nonlocal strain gradient theory. It is observed that the power-law index has important influences on the snap-buckling response. However, the snap-through instability of thermally postbuckled FG nanobeams on nonlinear elastic foundations has not been scrutinized. Snap-through type of instability may occur in antipathetic case where the direction of the applied lateral load and lateral deflection due to thermal postbuckling is opposite. In sympathetic case, where the direction of applied load and thermal deflection are the same, snap-through phenomenon may not occur. Hence, the problem of the present paper is a two-step solution. In the first step, thermal postbuckling analysis should be carried out, which is a nonlinear eigenvalue problem. Afterwards, the displacement vector and the associated temperature are obtained, and the response of the thermally postbuckled nanobeam under the action of a uniform lateral load of antipathetic type will be investigated. In the second step, a nonlinear bending problem should be solved, which consists the possibility of the limit loads and snap-through phenomenon.
Despite that various endeavors to investigate the static and dynamic analysis of FG nanostructures in the literature, there is no research reported on the nonlinear temperature-dependent thermal stability and the snapthrough behavior of thermally postbuckled FG nanobeams resting on nonlinear elastic foundations. To this end, the current research paper aims to analyse an examination of nonlinear thermal stability of perfect/imperfect temperature-dependent FG nanobeams under in-plane thermal load. Three-parameter elastic foundation with cubic nonlinearity, temperature dependency of properties, and initial geometrical imperfections are also taken into account. Thermo-mechanical properties of the FG nanobeam are assumed to vary in the nanobeam thickness and are determined with regard to the power-law function. The nonlinear nonlocal governing equations and related boundary conditions of graded size-dependent nanobeams are derived utilizing the Euler-Bernoulli beam theory in conjunction with the geometrically nonlinear von-Karman strain field by means of the nonlocal elasticity theory and virtual displacement principle. The Ritz method on the basis of Chebyshev polynomial of the first kind is implemented to form the matrix representation of the highly coupled nonlocal nonlinear governing equations associated with different types of boundary conditions. In the cases that bifurcation type buckling may occur, the thermal postbuckling path is traced using a direct displacement control strategy. Afterwards, the Newton-Raphson iterative procedure as a standard load control method is used to extract the nonlinear thermal bending curves. Furthermore, it should be noted that tracing the nonlinear equilibrium paths through limit points is not achievable using standard load control procedure and so pathfollowing techniques must be applied to examine snap-through instability of FG nanobeams due to lateral loading opposite to the direction of deflection. Hence, to follow the nonlinear equilibrium path of a thermally postbuckled nanobeam, the cylindrical arch-length method is used. After comparing the numerical outcomes of this work with the available results in the open literature, a comprehensive parametric study is conducted to investigate influence of the size-dependent nonlocal parameter, nonlinear elastic foundation coefficients, power-law index, imperfection amplitude, and different types of boundary conditions on the nonlinear thermal stability and snap-through response of the thermally preloaded FG nanobeam.

Theoretical formulations
Assume a FG nanobeam with rectangular cross-section, height h, length L, and width b in the Cartesian coordinates system (x, z), as illustrated in figure 1. It is assumed that the FG nanobeam during deformation is supported by a nonlinear elastic foundation and is exposed to uniform lateral and thermal loads.

Material properties descriptions
The thermo-mechanical material properties are estimated based on the power-law model distribution across the nanobeam thickness for the volume fractions of the constituents as [57] where V c and V m are the volume fractions of ceramic and metal, respectively. The power-law index ζ is a positive constant z ¥   0 ( ) which demonstrates the dispersion profile of material properties. The effective material properties considered for the FG nanobeam, including Young's modulus E and thermal expansion coefficient α, are obtained according to the Voigt rule. Based on this rule, each nonhomogeneous properties P(z, T) may be functionally graded in the following form in which the corresponding metal and ceramic properties are denoted by subscripts m and c, respectively. Moreover, the nonlinear function of temperature according to the Touloukian formulation is implemented to considering the temperature-dependent material properties as follows [57] = + + + + --

P T P P T P T P T P T
Here, for each of the constituents, T is temperature calculated in Kelvin and P i ʼs represent temperaturedependence coefficients.

Nonlocal elasticity theory in the presence of thermal loading
Within the framework of the Eringen's nonlocal elasticity theory [13,14], at a point x in an elastic continuum body the stress tensor is not only calculated based on the strain at its point but also based on the strains at all other points. Therefore, the nonlocal stress tensor σ ij at point x is written as follows where ¢ x refers to all neighbors of a point, is the kernel function defining the nonlocal modulus, and C ijkl are the fourth-order elasticity tensor with 81 coefficients. Moreover, ε kl is the strain tensor, and τ is the material constant, which is specified as τ=e 0 a/l;a/l being a dimensionless ratio between the internal characteristic length (e.g. lattice parameter and granular distance) and external characteristic length (e.g. crack length) of the nanostructures, and by conforming the dispersion paths of plane waves with those of atomic lattice dynamics, the material constant e 0 is obtained experimentally or approximated by atomistic and molecular dynamic simulation. By extending the nonlocal elasticity theory in the presence of the thermal effect, equation (4) takes the form as where α ij and T 0 are the thermal expansion coefficient and reference temperature, respectively. Furthermore, equation (5) is simplified to make the integral equation become solvable by utilizing a linear differential operator m = -  1 2 stated by Eringen [13], where ∇ 2 is the Laplacian operator and m = e a 0 2 ( ) is called the nonlocal parameter exhibiting the size effect on the nanostructure response.
The value of the nonlocal parameter depends on boundary condition, chirality, mode shapes, number of walls, and the nature of motions. Huang et al [58] calibrated the value of e 0 for the static bending analysis of single-layered graphene sheet using molecular dynamic simulations. The nonlocal parameter is experimentally obtained for various materials; for instance, a conservative estimate of μ<4 nm 2 for a single-walled carbon nanotube is proposed by Wang [59]. There is no rigorous study made on estimating the value of nonlocal parameter to simulate mechanical behavior of functionally graded micro/nanobeams. Hence all researchers who worked on the size-dependent mechanical behavior of FG nanostructures based on the nonlocal elasticity theory investigated the size effect on the mechanical behavior of FG nanostructures by changing the value of the nonlocal parameter. In the present research, a conservative estimate of the nonlocal parameter is considered to be in the range of μ=0-4 nm 2  [37,40].
Then, equation (5) can be described equivalently in differential form as Based on the symmetry of the stress and strain tensors [60], the elasticity tensor must have the following However, relations 7 reduce the number of independent components to 36 constants. Also, the concept of strain energy [60] leads to the relation C ijkl =C klij or equivalently C ij =C ji , which provides further reduction to 21 independent elastic components for general anisotropic elastic materials.
In order to continue further, we must address the issues of material isotropy. If we assume isotropic behavior, it can be shown that the most general form that satisfies this isotropy condition is given by [61] g d d g d d where γ 1 , γ 2 , and γ 3 are arbitrary constants. Using the general form 8 in the stress-strain relation 6, and letting g = n n n + - where ν is the Poissonʼs ratio. Relation 9 is called the generalized Hookeʼs law for linear isotropic elastic solids. Notice the significant simplicity of the isotropic form when compared to the general stress-strain law originally given by 6. It should be noted that only two independent elastic constants are needed to describe the behavior of isotropic materials.
In the expanded form, equation 9 can be written as In the present research, the governing equations and related boundary conditions of graded nanobeam are derived utilizing the Euler-Bernoulli beam theory. Hence, the Euler-Bernoulli beam theory assumptions are as follows • The cross-section is infinitely rigid in its own plane.
• The cross-section of a beam remains plane after deformation.
• The cross-section remains normal to the deformed axis of the beam.
The above assumptions state that the shear and transverse components of stress tensor are assumed to vanish, σ yy =0 and σ zz =0. Substituting these results into equation (10) determines ε yy and ε zz , and by substitution of obtained expressions into σ xx component of stress-strain relation gives the simplified constitutive relation along with the nonlocal theory and linear thermoelasticity for an isotropic linear elastic material where σ xx and ε xx are the axial components of the stress and strain tensors, respectively, and ΔT=T−T 0 represents the temperature rise of the graded nanobeam.

Governing equations
The displacement field of the geometrically imperfect FG nanobeam is written with regard to the Euler-Bernoulli theory as where in equation (12) u and w are displacements through the axial and transverse directions, respectively. Also, u 0 and w 0 are the mid-plane displacement components of the nanobeam.
Postbuckling has a nonlinear behavior, and consequently, a geometrically nonlinear analysis should be conducted. According to the von-Karman geometrical nonlinearity consistent with moderate rotations and small strains, the nonlinear strain-displacement relation can be written as [62] in which w * is the initial geometrical imperfection function through the nanobeam which denotes a deviation concerning to the flat position.
To obtain the governing nonlinear equations of the FG nanobeam involving the region of volume  under thermal loading, the static form of the generalized Hamilton's principle identified as virtual displacement principle may be utilized. Regarding the potential energy of the external uniform lateral load, the virtual total potential energy of the graded nanobeam is equivalent to the sum of the virtual strain energy of the nanobeam and the virtual energy of elastic foundation. Thus, in an equilibrium state, one may write where d is virtual work done by uniform lateral load q on nanobeam. Also, the virtual strain energy δU s and the virtual energy of nonlinear elastic foundation δU f of the nanobeam per unit width are equal to ò ò ò d sd e d d d where K w , K g , and K nl indicate the Winkler, Pasternak, and nonlinear coefficients of elastic foundation, respectively. Substitution of equations (13) and (15) in which the stress resultants may be generated upon the integration of stress field through the thickness of nanobeam. Stress resultant components, within the framework of the Euler-Bernoulli beam theory, can be defined as By substituting equation (11) into equation (17), the corresponding stress resultants of FG nanobeam in terms of the displacements are obtained as following forms In equation (18), the stiffness components A 11 , B 11 , and D 11 represent the stretching, coupling bendingstretching, and bending stiffnesses, respectively. Furthermore, N T and M T denote the axial thermal force and thermal moment resultants of the graded nanobeam calculated by T T The equilibrium equations of FG nanobeams are achieved within the standard scheme, involving integrating by parts in equation (16) to relieve the virtual primary variables and then setting the coefficients of δu 0 and δw 0 equal to zero results in the next equations The complete set of associated boundary conditions, arising from the integration procedure at x=−L/2 and x=L/2, are revealed as

Ritz method
For the practical problems with general geometry, loading, and boundary conditions, exploring their analytical solutions is impossible because of the mathematical complexity of the size-dependent models compared to the classical ones. Consequently, numerical strategies such as finite element method, differential quadrature method, Ritz method, Galerkin method, etc become the most proper ones for solving such problems. Among different numerical techniques, the Ritz method is a powerful tool and commonly employed for the analysis of structures. The Ritz method needs the expansion of unknown functions in infinite series, and by assuming the adequate number of functions, it is possible to achieve a high accuracy solution to the problem examined. Also, it is well known that the weak form of equations includes both the governing differential equation and the natural boundary conditions of the problem, and it puts less stringent continuity requirements on the approximate solution than the original differential equation. Therefore, the Ritz method with Chebyshev basis polynomials of the first kind is applied to solve the nonlocal nonlinear governing equations of the graded nanobeam on nonlinear elastic foundation in which the creation of a weak form of equations (20) At this step, to achieve the spatial approximation, each of the essential variables of this research, i.e., u 0 and w 0 may be approximated as where M is the total number of shape functions to ensure the convergence of the series expansion. Furthermore, N m u (x) and N x m w ( ) represent the approximation functions of Ritz method which should be determined with regard to the essential type of boundary conditions [62].
Upon the substitution of the series expansion 23 for essential variables into 22, results in the following matrix representation of the nonlocal nonlinear equilibrium equations In a more compact form, equation (24) may be represented as where T m (x) is the m-th order Chebyshev polynomial function that is nonzero at x=−L/2 and x=L/2. In this regard, table 1 provides mentioned shape functions associated with three types of boundary conditions, namely; clamped-clamped (C-C), clamped-simply supported (C-S) and simply supported-simply supported (S-S). Furthermore, there arise three circumstances needing different approaches to handle in designating equilibrium curve of size-dependent FG nanobeams, which are provided next.

Thermal postbuckling analysis
In order to analyse thermal postbuckling of the graded nanobeam, the governing equation is achieved for the thermal postbuckling behavior when the force vector is omitted from equation (25). In such conditions, the nonlocal governing equation for thermal postbuckling analysis is written as To solve the above system of equations, which is in the form of a nonlinear eigenvalue problem, a direct displacement control method is utilized to trace the thermal postbuckling path. At first, to achieve the buckled configuration and critical buckling temperature of the nanobeam a linear analysis should be done. To this end, the critical temperature difference is extracted for measured thermo-mechanical properties at reference temperature. This procedure which corresponds to the buckled state of nanobeam proceeds on to get the converged critical temperature difference and also the related eigenvector. After implementing the linear analysis, the next process is applied to obtain the thermal postbuckling path. Accordingly, it is assumed that an element of the displacement vector is known. Mid-point deflection of the nanobeam is an appropriate choice for this end. The mid-point deflection of nanobeam is incremented at the first iteration of each stage. Next, determining the eigenvalue problem in which the incremental displacement vector updates the nonlinear stiffness matrix supplies new estimations to the both of eigenvalue and eigenvector. After making suitable iterations, the convergence is obtained. For the next step, the mid-point deflection of nanobeam is now incremented again. The procedure is continued until the desired load-deflection path is reached.

Nonlinear thermal bending
In the circumstances that simply-supported boundary conditions are selected or initial geometrical imperfection for nanobeam is assumed, there is a unique and stable equilibrium curve for thermally loaded FG nanobeams which is directed by equation (25). Owing to the appearance of nonlinear terms in this equation, an iterative procedure has to be done for each load increment. In this paper, the Newton-Raphson iterative procedure that is a standard load control method is employed. The left-hand side of equation (24) is expanded by Taylor's series about the known solution, after assembling the element equations, results in where T X T , [ ( )]is called the tangent stiffness matrix. Based on the Newton-Raphson method, the converged solution of each load increment is employed as the new guess vector for the first iteration of the next step in which the ultimate applied load level is divided into a number of small load increments. The procedure is continued until the desired load level is obtained. The elements of this matrix in conjunction with the nonlinear stiffness matrix in equation (24) are defined in appendix.

Snap-through phenomenon
In solution of nonlinear thermal bending and thermal postbuckling problem, temperature and the associated displacement vector are obtained at the end of each solution procedure. Since the nanobeam is curved due to the thermally postbuckled configuration and edge is of immovable type, there is the possibility of limit load type of instability when nanobeam is subjected to uniform lateral load of antipathetic type. It should be pointed out that tracing the nonlinear equilibrium paths within limit points is not obtainable using standard load control strategies such as Newton-Raphson scheme or Picard method. Hence, path-following methods are employed in tracing the equilibrium trajectories of a structure with limit load points. In the present work, the cylindrical arclength method developed by Crisfield [63] is adopted to trace the load-deflection path through the limit loads. The governing equation to trace the equilibrium path in this case can be expressed as The details of the process according to this procedure are omitted here for the sake of brevity, and the reader is referred to [63].

Results and discussions
The solution procedure mentioned in the previous section is applied herein to examine the nonlinear thermal stability behavior and snap-through features of size-dependent FG nanobeams embedded on nonlinear elastic foundation. To this end, the upper surface of the FG nanobeam is Silicon Nitride (Si 3 N 4 ) (as the ceramic-rich), and the lower surface of the nanobeam is Stainless Steel (SUS304) (as the metal-rich). Following equation (3), thermo-mechanical material properties of the aforementioned constituents are assumed to be temperaturedependent, and the coefficients P i for each property of mentioned constituents are given in table 2. Accordingly, for temperature-dependent (TD) case, each property should be estimated at the current temperature, while for temperature-independent (TID) material case, properties are obtained at reference temperature, i.e. T 0 =300 K, using equation (3). It is worth mentioning that the nonlocal effects disappear after a certain length. For instance, scale effects predicted by the nonlocal elasticity on the vibration and buckling of FG nanobeams gradually disappear for L>50 nm, and when L is greater than 50 nm local and nonlocal results are the same [37]. Unless otherwise stated, in all of the numerical results, h=2 nm [42,64] and L/h=25 (L=50 nm) are assumed in the rest. Additionally, the following parameters are introduced and employed to generalize the subsequent results where E c ref stands for the elasticity modulus at reference temperature for ceramic constituent. Furthermore, the initial geometrical imperfection through the nanobeam is assumed by the following function in which η is known as the imperfection amplitude parameter.

Comparison studies
Firstly, two comparison studies are conducted to demonstrate the accuracy and proficiency of the proposed model. To the best of the authors' knowledge, nonlinear thermal stability and snap-through results of thermally postbuckled size-dependent FG nanobeams are not available in open literature. However, when the nonlocal parameter is assumed to be zero, FGM beams may be considered as a special case of size-dependent beams. For the first comparison study, in figure 2, for both TID and TD cases of material, thermal postbuckling path of C-C FGM beams based on the conventional elasticity theory is compared with those of Esfahani et al [6] employing the Timoshenko beam theory and the generalized differential quadrature method. The considered parameters such as length to thickness ratio and the power-law index values of the contact-less beam are assumed as L/h=25 and ζ=1, respectively. It can be revealed from this figure that results obtained by the present research are in excellent agreement with those of Esfahani et al [6].
In another comparison study, the critical buckling temperature difference of homogeneous nanobeams is evaluated from the present research and compared with those given by Ebrahimi and Salari [40]. Comparison study of S-S nanobeam is presented in table 3 for different values of nonlocal parameters. It is seen that, results are in excellent agreement with those provided by Ebrahimi and Salari [40] which accepts the accuracy and correctness of the present formulation. The small divergences are due to the different beam theories.

Parametric studies
Eventually, a parametric and comprehensive study is performed to investigate the influence of the power-law index, nonlocal parameter, imperfection amplitude, nonlinear elastic foundation coefficients, and different Table 2. Temperature dependent coefficients for Si 3 N 4 and SUS304 [57].

Material
Properties 12.330e-6 0 8.086e-4 0 0 types of boundary conditions on the thermal postbuckling, nonlinear thermal bending and snap-through response of nonlocal graded nanobeams. Table 4 provides the influences of temperature dependency of materials, various power-law indices, nonlocal parameters, and contact conditions on the critical temperature difference ΔT cr of uniformly heated C-C FG nanobeams. From this table, it is evident that the predicted critical buckling temperature under TID material properties is higher than the predicted one following TD case, where in TID case the buckling temperature is overestimated. Moreover, it is obvious that, as the nonlocal parameter increases, the critical buckling temperature of the nanobeam decreases. In other words, it designates that the nonlocality causes a reduction in the stiffness of an FG nanobeam. Different contact conditions are considered, and the nonlinear parameter of the elastic foundation is put equal to zero. This is due to the fact that the nonlinear elastic foundation coefficient does not affect the buckling temperature of FG nanobeam. As observed, the critical temperature of the nanobeam increases with increasing the Winkler constant of the elastic foundation, which is due to the increase in the total stiffness of the structure foundation interplay. Further, with the introduction of a shear layer on Pasternak elastic foundation, results in an increase in buckling temperature and stability enhancement. Through the presented results, it is also concluded that, as the power-law index increases, the buckling temperature difference of the nanobeam decreases. This is anticipated because of less bending rigidity results from more metal-rich FG nanobeams. Figure 3 is presented to study the influence of the power-law index and temperature dependency on the thermal postbuckling equilibrium paths of the contact-less C-C FG nanobeams in which the bifurcation type buckling appears when power-law index is an arbitrary value. As a consequence, due to stretching-bending coupling and nonsymmetric distribution of material properties concerning the mid-plane of the FG nanobeams, thermal moment resultant is generated through the FG nanobeam. In such condition, for perfect nanobeams, when both edges are clamped, the supports handle the induced thermal moments, and FG nanobeam undergoes primary-secondary equilibrium path. It can also be revealed from this figure that the postbuckling paths for the FG nanobeam with TD case are lower than that of the FG nanobeam with TID material properties and continuing temperature dependency of the considered material properties remarkably Figure 2. A comparison on thermal postbuckling path of FGM beams between the results of the present study and those reported by Esfahani et al [6] for C-C Timoshenko beams. enhance the accuracy of solutions. For this reason, in the consequent results of the present research, only TD material is taken into the consideration. Figure 4 illustrates the nonlocal nonlinear thermal stability of contact-less FG nanobeams corresponding to various values of power-law indices for C-S and S-S types of boundary conditions. It is evident that for each type of boundary condition, only the case of homogeneous nanobeams exhibits the bifurcation type behavior while for nonhomogeneous ones, the response of the nanobeam is the unique and stable nonlinear bending type. Therefore, FG nanobeam initially deflects downward since the thermal expansion coefficient of the ceramic constituent is lower than the metal one. It is worthwhile to notice that, simply-supported edge is unable to apply the additional bending moments at the support and hence, the edge support does not handle the induced thermal moments. It also is verified that the buckling temperature of a homogeneous nanobeam with C-S type of boundary conditions is higher than a homogeneous nanobeam with S-S type of boundary conditions. Figure 5 is exhibited to examine the effect of nonlocal parameter on the nonlinear responses of perfect and imperfect FG nanobeams with both edges clamped. In the results of this figure, two cases of linear contact condition, i.e. (k w , k g )=(0, 0), (100, 10) are considered. It can be inferred from this figure that for nanobeams with initial geometrical imperfection, a primary-secondary equilibrium path is not observed for the considered problem. As the nonlocal parameter increases, results in a decrease of critical buckling temperature and  postbuckling strength of the nanobeam, and the nonlocal parameter has a softening influence on the thermal stability behavior of FG nanobeams. Further observations from this figure indicate that due to the improvement of the elastic stiffness of the nanobeam with the existence of Winkler and Pasternak constants of the elastic foundation, the critical temperature difference of the nanobeam increases and thermal postbuckling deflection decreases. Figure 6 is depicted to elaborate on the effect of the nonlinear elastic foundation parameters on the thermal postbuckling path of linearly graded nanobeams. It is demonstrated that thermal postbuckling curve of incontact FG nanobeam is higher than its contact-less condition. Again it is confirmed that by increasing each coefficient of elastic foundation, the elastic stiffness of nanostructure heightens. Moreover, the nonlinear coefficient of elastic foundation does not affect critical temperature difference while this constant notably influences the postbuckling strength of the nanobeam. Especially, this influence is more prominent in deep postbuckling configuration.
The influence of nonlinear elastic foundation constants on nonlinear thermal bending curves of linearly graded nanobeams is indicated in figure 7. This examination is prepared for two types of boundary conditions, i.e., C-S and S-S. Based on the previous descriptions, when nanobeam edges are simply-supported, the additional thermal moments are not handled by the nanobeam edges. In this case, the obtained result is not of the bifurcation type of instability and nanobeam undergoes a stable and unique nonlinear bending path for C-S and S-S types of boundary conditions. Results of this figure reveal that nonlinear elastic foundation constants raise the nonlinear bending curve and result in less lateral deflection by increasing the stiffness of the nanostructure. Unlike the C-C type, it is also observed that the effect of nonlinear constant of elastic foundation in nonlinear thermal bending response emerges at the beginning of thermal loading.
To investigate the thermal stability behavior of the contact-less C-C FG nanobeam, the load-deflection curves are depicted in figure 8, and the effects of initial imperfection parameter and length to thickness ratio are     studied. Six values of initial geometrical imperfection parameter for both positive and negative magnitudes are considered. Evaluation of the obtained results reveals the fact that the behavior of perfect FG nanobeam is of bifurcation type of buckling while the nanobeams with initial imperfection experience nonlinear bending behavior as unique and stable path. However, the load-deflection curves are much closer to the bifurcation type of instability response by considering smaller imperfection amplitude parameters. In addition, it can be concluded that with increasing the L/h ratio from 25 to 35, the critical temperature of the graded nanobeam decreases and also thermal postbuckling equilibrium curve decreases. The physical reason for this observation is that when this parameter increases, the bending rigidity of the nanobeam decreases significantly. Figure 9 aims to compare the importance of the amplitude of initial imperfection on the nonlinear thermal bending of the contact-less linearly graded nanobeams for two types of boundary conditions, namely, C-S and S-S. Based on the prior discussions in figure 4, perfect FG nanobeams subjected to thermal loading with C-S and S-S types of boundary conditions deflect downward. However, the FG nanobeams may deflect upward with positive magnitudes of initial imperfection (which indicates that the nanobeam is originally curved upward before inducing a thermal load) when the amplitude of this parameter reaches a critical value. Hence, it is shown that FG nanobeams with C-S types of boundary conditions deflects downward with η=0, 0.01, and deflects upward when the critical value of imperfection parameter is η=0.03. Similarly, for FG nanobeams with both edges simply-supported, the critical value of initial geometrical imperfection is η=0.05.
In figures 10 and 11, the influences of the various nonlinear elastic foundation conditions on end-shortening force and moment at mid-point of linearly graded nanobeams for C-S and S-S types of boundary conditions are depicted, respectively. Only TD case of material properties is considered. It is observed that end-shortening force changes linearly with respect to temperature and there is no bifurcation point within the end-shortening force and moment-temperature paths of the nanobeam which confirms the uniqueness of a stable load-deflection curve. As revealed, by increasing each constant of the elastic foundation for a prescribed amount of temperature,  the end-shortening force of the FG nanobeam increases and the moment of the nanobeam decreases. These anticipations are the same for the C-S and S-S boundary conditions. The next parametric study is presented to scrutinize the effect of the deformed configurations and thermal preloading on the snap-through instability of contact-less FG nanobeams. To this end, the snap-through phenomenon of the linearly graded nanobeam is depicted in figure 12. The left figure exhibits snap-through response of C-C type, while the right one presents C-S type. Four different thermal preloading associated with the maximum mid-point deflections and values of W/h=±0.4, ±0.6, ±0.8 and ±1 are considered. It is revealed that when the maximum deflection due to the employed temperature is less than W/h=0.4, the snapthrough phenomenon does not take place due to small deflections caused by the associated temperature rise parameters. However, for more thermal preloading where the maximum mid-point deflection is greater than W/h=0.4, the induced deflection is large enough to existence of the snap-through type of instability. Moreover, it can be concluded from this figure that the intensity of the snap-through instability and also the upper and lower limit loads increase as the maximum deflection due to thermal preloading increases. Figure 13 investigates the snap-through buckling of linearly graded nanobeam subjected to the uniform lateral load with S-S type of boundary condition. Two different cases of nonlinear contact conditions, i.e. (k w , k g , k nl )=(0, 0, 0), (100, 10, 40) as well as four deformed configurations due to thermal preloading with the negative values of maximum deflections W/h=0.4, 0.6, 0.8 and 1 are considered. It is verified that for the maximum deflection higher than W/h=0.4, the snap-through feature is observant. Further consideration of this figure indicates that as the maximum deflection due to the thermal preloading enhances, upper and lower limit loads both increases and causes in more intensified snap-through instability. It is also found that as the constants of nonlinear elastic foundation increase, both of the upper and lower limit loads enhance, particularly for W/h=−1.   Figure 14 is exhibited to analyse the effect of the power-law index on snap-through behavior of contact-less FG nanobeams. It is supposed that FG nanobeam which is thermally preloaded reaches a postbuckled configuration equal to W/h=1. It is observed that as the value of power-law index increases, results in a decrease of upper limit loads of the nanobeam, and the power-law index has an important effect on the snapthrough phenomenon of the graded nanobeams. In addition, all trajectories cross through the point (0, −1) as an equilibrium point in which this point indicates that the bifurcation buckling is symmetric.
In figure 15, the effects of various power-law indices and thermal preloading on snap-through features of contact-less FG nanobeams are presented. In this example, four various values of power-law index such as 0.2, 0.5, 2 and 5 are considered and the nonlocal parameter is set equal to μ=4 nm 2 . The left figure represents the snap-through behavior of C-S type, which corresponds to W/h=−0.8, while the right one concerns S-S type with W/h=−1. Results of this figure reveal that the snap-through buckling behaviors of the thermally preloaded FG nanobeams are affected by this parameter. Therefore, it is shown that higher values of power-law index causes a reduction of upper limit loads.

Concluding remarks
In the present investigation, nonlinear thermal bending and thermal postbuckling of perfect/imperfect FG nanobeam in-contact with a nonlinear elastic foundation are first analysed. Afterwards, the snap-through  phenomenon of thermally deformed FG nanobeams due to uniform lateral load is presented. Thermomechanical material properties of the nanobeam are estimated according to a power-law form, which are considered to be temperature-dependent and also graded in the nanobeam thickness. Considering the nonlocal elasticity theory and the von-Karman assumptions, the nonlocal nonlinear governing equations and the related boundary conditions of the Euler-Bernoulli nanobeam are obtained. The Ritz method on the basis of Chebyshev polynomial is applied into the principle of virtual displacement to form the matrix representation of nonlocal governing equations. The solution for the nonlinear thermal stability and snapthrough behavior is divided into three different methodologies, including direct displacement control method, Newton-Raphson iterative procedure, and cylindrical arch-length scheme. After validating the proposed model, the parametric investigations are performed to examine the effects of different parameters on the nonlinear thermal stability and snap-through features of the FG nanobeams. In summary, the following general conclusions can be obtained: • The behavior of nonlinear thermal stability indicates that the equilibrium path of homogeneous /nonhomogeneous perfect nanobeam with C-C type of boundary condition is of the bifurcation type buckling. For C-S and S-S types of boundary conditions, only the case of homogeneous nanobeam exhibits the bifurcation type behavior while for nonhomogeneous ones, the unique and stable nonlinear thermal bending curve is observed.
• It is found that the stiffness and the critical temperature difference of the FG nanobeam decrease as the nonlocal parameter increases.
• Increasing the nonlinear elastic foundation constants enhances the elastic stiffness of the FG nanobeam and results in higher thermal postbuckling path and nonlinear bending curve. Also, it is revealed that as the constants of the nonlinear elastic foundation increase, the end-shortening force of the FG nanobeam increases, and the moment at mid-point of the nanobeam decreases.
• It is shown that C-C FG nanobeams with initial imperfection experience a unique and stable nonlinear thermal bending behavior. However, the load-deflection curves are much closer to the bifurcation type of buckling by considering smaller imperfection amplitude parameters.
• It is concluded that the deformed configurations due to thermal preloading play an essential role in the existence of the snap-through phenomenon of FG nanobeam. Increasing the maximum mid-point deflection due to thermal preloading increases the upper and lower limit loads, and also the intensity of the snap-through instability enhances.
• The snap-through behavior of the thermally preloaded FG nanobeam is affected by power-law index parameter. As the value of power-law index increases, the upper limit load of nanobeam decreases.
The components of the tangent stiffness matrices included formerly in equation (29)