Nonlinear Thermal/Mechanical Buckling of Orthotropic Annular/Circular Nanoplate with the Nonlocal Strain Gradient Model

This article presents the nonlinear investigation of the thermal and mechanical buckling of orthotropic annular/circular single-layer/bilayer nanoplate with the Pasternak and Winkler elastic foundations based on the nonlocal strain gradient theory. The stability equations of the graphene plate are derived using higher-order shear deformation theory (HSDT) and first-order shear deformation theory (FSDT) considering nonlinear von Karman strains. Furthermore, this paper analyses the nonlinear thermal and mechanical buckling of the orthotropic bilayer annular/circular nanoplate. HSDT provides an appropriate distribution for shear stress in the thickness direction, removes the limitation of the FSDT, and provides proper precision without using a shear correction coefficient. To solve the stability equations, the differential quadratic method (DQM) is employed. Additionally, for validation, the results are checked with available papers. The effects of strain gradient coefficient, nonlocal parameter, boundary conditions, elastic foundations, and geometric dimensions are studied on the results of the nondimensional buckling loads. Finally, an equation is proposed in which the thermal buckling results can be obtained from mechanical results (or vice versa).


Introduction
Nanoplates including graphene sheets can be named as one category of eminent materials for the next generation of nanodevices.Graphene sheets can be utilized in the manufacture of many nanodevices like sensors and memory devices [1] and there are other applications such as nanosheet resonators [2], mass sensors [3], and gas sensors [4].Because performing highly accurate experiments at the nanoscale is very demanding and high-cost, some techniques such as continuum-based modelling and molecular dynamics (MD) have been employed recently that have the potential to consider the atomic size and length effects.MD modelling of nanoscale structures needs time-consuming procedures and complicated computations.Continuum modelling of nanostructures is computationally less expensive than MD modelling; therefore, continuum theories have been formulated and employed for the investigation of the mechanical properties of nanoscale structures.In other words, continuum modelling of nanostructures can assist in understanding the results of experimental measurements (or MD modelling) more properly, because it can decrease the volume and computational time.In order to make computations easier, continuum mechanics can be employed as an appropriate option for theoretical investigations.The atomic length scales and interatomic force coefficients can be incorporated into the constitutive equations.In the classical continuum model, the influences of size in nanostructures are not regraded [5,6].
Experimental investigations have proved that there are size-dependent effects on the mechanical properties at the nanostructures [7].Therefore, modified continuum theories are required to consider small-scale influences and maintain the ease and computational efficiency of continuum theories.Consequently, new theories have been suggested to take into account nanoscale effects in nanostructures.Various nonclassical continuous models have been presented considering the material length scale coefficient, and multiple modifications to classical elasticity theories including Eringen's nonlocal elasticity theory [8], the strain gradient theory [9], and the nonlocal strain gradient theory [10] have been proposed.
Narendar [11] perused the buckling of rectangular nanoplates on the refined basis of the two-variable plate theory considering nonlocal scale through the Navier technique.He observed that by increasing the size of the square nanoplate, the buckling load ratio decreases significantly.Jamalpoor and his colleagues [12] studied the vibration and buckling of nanoplates exposed to electric as well as magnetic potentials based on the nonlocal plate theory.They noticed that the buckling load increased by increasing the external magnetic potential.Li and Hu [13] used nonlocal continuum theory to study the buckling, vibration, and bending of composite nanobeams.They concluded that the critical buckling load can be increased by reducing the nonlocal coefficient when the nonlocal coefficient is larger than the material characteristic coefficient.Shaban and Alibeigloo [14] investigated the bending and vibration analysis of carbon nanotubes on elastic foundations using the nonlocal three-dimensional theory of elasticity.They observed that the static and vibrational behaviors of carbon nanotubes are significantly affected by elastic foundations.Rokni et al. [15] studied the free vibration of circular/annular plates considering variable thickness with various boundary conditions using three-dimensional elasticity.Shaban and Mazaheri [16] analyzed the electrostatic behavior of microsandwich panels with the functionally graded core using the nonlocal three-dimensional theory of elasticity.They observed that increasing the nonlocal coefficient led to an increase in the displacements.The authors of Ref. [17] perused the free vibration of functionally graded plates on the elastic foundation with elastically restrained edges using FSDT.In their paper, they studied the influences of coefficients including thickness-to-radius, material distribution, various combinations of constraints at edges on the frequency, foundation stiffness coefficients, mode shape, and modal stress.
It is noticed that in the former analyses, the small-scale effect in sort of the stress nonlocality led to a softening-stiffness effect, while the strain gradient size dependency resulted in a hardening-stiffness effect.Consequently, various nonclassical continuum theories of elasticity have been presented to prognosticate various size dependencies in the mechanical features of nanostructures.The theory of nonlocal strain gradients is one of the prominent models for analyzing size-dependent mechanical features.Contrary to classical elasticity theory, the nonlocal strain gradient theory can consider the influences of stiffness-hardening and stiffness-softening by assuming the strain and stress gradient.Therefore, Lim and his colleagues [10] studied a new size-dependent elasticity model named nonlocal strain gradient theory considered both softening and stiffening effects to illustrate the size dependency more accurately.Some analyses have been carried out with respect to the nonlocal strain gradient theory [18].For example, the authors of Ref. [19] investigated the buckling of the nanocrystal shell with the axial loads assuming thermal, magnetic, and electrical conditions with the nonlocal strain gradient theory.They deduced that the nonlocal parameter affects the critical buckling load more notably for nonlocal parameters ranging from 2 to 20 nm.Tanzadeh and Amoushahi [20] perused the buckling of rectangular nanoplates subjected to various uniaxial in-plane loads using the nonlocal strain gradient theory with the higher-order finite strip technique.They concluded that increasing the strain gradient coefficient resulted in an increase in buckling load.Cuong-Le et al. [21] studied the bending, buckling, and vibration of sandwich nanoplates based on the theory of nonlocal strain gradients.They noticed that the nonlocal and strain gradient coefficients result in stiffness reduction and stiffness hardening cases.Therefore, both of these coefficients play a significant effect in the buckling of sandwich nanoplate.Wang and his colleagues [22] perused the buckling of functionally graded nanotubes based on the nonlocal strain gradient theory with high-order theory via the generalized differential quadrature technique.Wang et al. [23] analyzed the buckling of porous functionally graded porous nanobeams considering hybrid effects using the nonlocal strain gradient principle with the aid of the generalized differential quadrature technique.Esen and Özmen [24] studied thermal vibration and buckling functionally graded porous rectangular nanoplates considering magnet electroelastic effects based on the nonlocal strain gradient theory and FDST.They revealed that the temperature increase influences the buckling characteristic of the nanoplate depending on the proportions of the constituents in the nanoplate.Tang and Qing [25] examined the buckling and vibration of functionally graded beams based on the theory of nonlocal strain gradient and using the Laplace transform method.They scrutinized that the effects of length-scale coefficients on the vibration frequencies and buckling loads increase with the rise in vibration and buckling order.Al-Furjan and his colleagues [26] conducted research on the dynamic buckling of carbon nanocones with the nonlocal strain gradient theory and FSDT.They demonstrated that by increasing the strain gradient coefficient, the dynamic instability region occurred at high frequencies.Magorzata Chwa and Muc [27] utilized the nonlocal strain theory and considered higher-order shear deformation theories to analyze buckling and free vibrations of rectangular nanoplates via the Rayleigh-Ritz method.Using nonlocal strain gradient elasticity theory and third-order shear deformation plate theory, Fan et al. [28] studied buckling and postbuckling porous functionally graded nanoplates.They noticed that by increasing the material characteristic gradient index, both nonlocal and strain gradient size influences are noticeable.Sadeghian et al. [29] studied the investigation of large deflection of the annular/circular nanoplate on the basis of nonlocal strain gradient theory.They concluded that both strain gradient and nonlocal coefficients have noticeable effects on decreasing or increasing the deflection of the nanoplate.
One of the important issues in evaluating the stability of various structures is the buckling analysis, which has drawn the attention of many researchers to perform research on it.The term buckling refers to the loss of stability, or in other words, when the structure changes from stable equilibrium to unstable equilibrium, the structure undergoes buckling.In fact, the behavior of the structure under load in which with a small increase in load a disproportionate increase in displacement occurs in the structure is called buckling, and the amount of force for which the structure buckles is called critical buckling load [30].There are many papers that consider buckling analysis.For example, Zghal et al. [31] scrutinized the postbuckling of functionally graded nanotubes considering various mechanical loadings.In their article, they discussed the influences of the carbon-nanotubes (CNT) volume fractions, CNT distributions, gradient indexes, and boundary conditions.In another paper, Zghal and Dammak [32] studied the buckling of functionally graded structures subjected to various compression loads.They noticed that the evenly porous FGM structures have the minimum buckling load while the uneven ones possess the maximum one.Trabelsi et al. [33] studied thermal buckling of functionally graded shells.They demonstrated that the flexural rigidity of the FG cone exposed to thermal loads can be improved with the variation of the powerlaw parameter.Mehr et al. [34] studied functionally graded CNT-reinforced composite shells with thermal buckling.They noticed that increasing aspect ratio and CNT volume fraction cause higher critical buckling temperature.Van Do and Lee [35] investigated the thermal buckling of FGM plates based on the quasi-3D higher-order shear deformation theory.They noticed that the temperature distribution on the plate noticeably influences the thermal buckling responses of FGM plates.
DQM is one of the numerical methods that has been employed by many authors [36,37].The initial purpose of this procedure is to employ Lagrange interpolation polynomials to field constants and to solve the equations at discrete grid nodes.Higher accuracies are gained by considering more grid nodes.With the aid of DQM, Han et al. [38] studied one-electrode microresonators on the basis of a generalized 1DOF principle.Duryodhana et al. [37] employed DQM to study the buckling (and free vibrations) of foamed composites.By means of DQM, Ren et al. [39] investigated the thermal and mechanical buckling behavior of heated rectangular plates considering the characteristics of the temperaturedependent material.
FSDT accounts for the shear deformation influences by the way of linear variation of in-plane displacements through the thickness.Because the FSDT violates the conditions of zero transverse shear stresses on the top and bottom surfaces, a shear correction factor (which can depend on many parameters) is needed to compensate for the error caused by a constant shear strain assumption through the thickness.The HSDTs account for the shear deformation influences and satisfy the zero transverse shear stresses on the top and bottom surfaces; therefore, a shear correction factor is not needed.HSDTs are modelled based on the assumption of higher-order variations of in-plane displacements or both in-plane and transverse displacements through the thickness.Some of HSDTs are computational costs because with each additional power of the thickness coordinate, an additional unknown can be added to the model [40].Also, there is quite the difference between FSDT and HSDT behavior for thick plates, but for thin plates both of the theories predict approximately similar result [41].
Until now, the nonlocal strain gradient theory with HSDT and FSDT has not been carried out for studying the nonlinear thermal as well as mechanical buckling of the orthotropic annular/circular nanoplate via DQM.Furthermore, the buckling analysis of a bilayer annular/circular nanoplate with elastic foundations is also studied.The effects of the nonlocal parameter, strain gradient coefficient, boundary conditions, geometric dimensions, and elastic foundations on the results of the nondimensional buckling loads are investigated.From the results, it can be concluded that the effect of the nonlocal parameter on the buckling of the circular/annular nanoplate is more significant than the types of boundary conditions or elastic foundations.Also, it is noticed that a bilayer nanoplate cannot be assumed as an equivalent single-layer nanoplate (with the same thickness as the bilayer nanoplate) to gain the results of the buckling of the nanoplate.The results of the current article can be a benefit for the design as well as improvement of nanostructured devices, including circular gate transistors, microswitches, and so on.are gained by considering more grid nodes.With the aid of DQM, Han et al. [38] studied one-electrode microresonators on the basis of a generalized 1DOF principle.Duryodhana et al. [37] employed DQM to study the buckling (and free vibrations) of foamed composites.By means of DQM, Ren et al. [39] investigated the thermal and mechanical buckling behavior of heated rectangular plates considering the characteristics of the temperaturedependent material.

Theory and Formulation
FSDT accounts for the shear deformation influences by the way of linear variation of in-plane displacements through the thickness.Because the FSDT violates the conditions of zero transverse shear stresses on the top and bottom surfaces, a shear correction factor (which can depend on many parameters) is needed to compensate for the error caused by a constant shear strain assumption through the thickness.The HSDTs account for the shear deformation influences and satisfy the zero transverse shear stresses on the top and bottom surfaces; therefore, a shear correction factor is not needed.HSDTs are modelled based on the assumption of higher-order variations of in-plane displacements or both inplane and transverse displacements through the thickness.Some of HSDTs are computational costs because with each additional power of the thickness coordinate, an additional unknown can be added to the model [40].Also, there is quite the difference between FSDT and HSDT behavior for thick plates, but for thin plates both of the theories predict approximately similar result [41].
Until now, the nonlocal strain gradient theory with HSDT and FSDT has not been carried out for studying the nonlinear thermal as well as mechanical buckling of the orthotropic annular/circular nanoplate via DQM.Furthermore, the buckling analysis of a bilayer annular/circular nanoplate with elastic foundations is also studied.The effects of the nonlocal parameter, strain gradient coefficient, boundary conditions, geometric dimensions, and elastic foundations on the results of the nondimensional buckling loads are investigated.From the results, it can be concluded that the effect of the nonlocal parameter on the buckling of the circular/annular nanoplate is more significant than the types of boundary conditions or elastic foundations.Also, it is noticed that a bilayer nanoplate cannot be assumed as an equivalent single-layer nanoplate (with the same thickness as the bilayer nanoplate) to gain the results of the buckling of the nanoplate.The results of the current article can be a benefit for the design as well as improvement of nanostructured devices, including circular gate transistors, microswitches, and so on.Taking into account the HSDT, the displacement field can be written in the r, θ, and z axes clarified by U, V, and W, respectively.

Theory and Formulation
where u 0 and w 0 define the displacements of the midplane on the r and z axes, respectively.Furthermore, φ(r) is the rotation component around the θ axis.Also, the function g(z) can be mentioned as: f (z) and y* can be assumed as various functions utilized in different papers which are listed in Table 1 (for instance, the Ambartsumian [42] model can be considered as − 1 6 z 3 Taking into account the HSDT, the displacement field can be written in the r, θ, and z axes clarified by U, V, and W, respectively. where 0 u and 0 w define the displacements of the midplane on the r and z axes, respec- tively.Furthermore, ( ) r φ is the rotation component around the θ axis.Also, the function g(z) can be mentioned as: f (z) and y* can be assumed as various functions utilized in different papers which are listed in Table 1 (for instance, the Ambartsumian [42] model can be considered as The nonlinear strains, considering von Karman's presumptions, can be written as Equations ( 3)-( 7): ( ) ( ) rz dg z dW dU dr dz dz It is noted since the symmetric assumption is considered, Equations ( 5) and ( 7) are equal to zero.Taking into account the HSDT, the displacement field can be written in the r, θ, and z axes clarified by U, V, and W, respectively.

Model
where 0 u and 0 w define the displacements of the midplane on the r and z axes, respectively.Furthermore, ( ) r φ is the rotation component around the θ axis.Also, the function g(z) can be mentioned as: f (z) and y* can be assumed as various functions utilized in different papers which are listed in Table 1 (for instance, the Ambartsumian [42] model can be considered as The nonlinear strains, considering von Karman's presumptions, can be written as Equations ( 3)-( 7): ( ) ( ) rz dg z dW dU dr dz dz It is noted since the symmetric assumption is considered, Equations ( 5) and ( 7) are equal to zero.

Model
Table 1.Recommended functions used in HSDT by various authors.
The nonlocal form (NL) of the force and momentum resultants can be clarified as: The potential energy of the system is the sum of the strain energy created due to the work of internal forces as well as the potential energy created by external forces.
Π, U, and Ω are described as the potential energy of the entire system, the strain energy of the structure, and the potential energy of external forces, respectively.According to the concept of minimum potential energy, for a structure in equilibrium, the variation in the potential energy is equal to zero.
By zeroing δΠ, the coefficients of δu 0 , δw 0 , and δφ must be equal to zero.The Euler-Lagrange equations (which are in non-local form, so they are marked using superscript NL) are obtained as Also, Equation ( 16) can be rewritten as The nonlocal strain gradient theory (assumed as a combination of both the strain gradient theory and the nonlocal stress field) was suggested by Lim et al. [10]: in which C ijkl , l, and µ define elastic, strain gradient (internal material length scale), and nonlocal parameters, respectively.In addition, the constitutive stress-strain equation in the very small scale is defined as [49] (1 It should be mentioned in Equation ( 19), E 1 and E 2 indicate the Young modulus along the 1 and 2 directions.Additionally, v 12 and v 21 are Poisson ratios in those directions.Also, G 13 signifies the shear modulus.
The nonlocal form can be mentioned as The force and moment resultants in the local form can be written as Moreover, the resultants in terms of displacements are gained as follows: The equilibrium equations of a monolayer axisymmetric circular/annular nanoplate considering Pasternak and Winkler elastic foundations can be locally written as To derive the stability equations and determine the critical buckling force, the adjacent equilibrium method is used.Therefore, displacements are replaced by the following equations: Therefore, for the resultants of force and moment By substituting Equation (37) into equilibrium Equations ( 33)- (35), Equations ( 38)-( 40) are obtained: In Equations ( 38)- (40), by omitting the terms that have only zero subscript and setting terms related to the preload state equal to zero and considering the uniform compressive load, we have Therefore, the stability equations in terms of local stresses are written in the form of Equations ( 42)-( 44):

Derivation of the Governing Equations for the Bilayer Annular/Circular Nanoplate
Numerous graphene layers can be employed to fix their weak buckling strengths.Therefore, graphene plates are set on top of each other by means of van der Waals bonds, generating graphene layers [42].Van der Waals force is a general term used to explain the attraction of intermolecular forces among molecules (for more information please refer to Refs.[50,51]).Figure 2 shows a bilayer nanoplate (considering van der Waals forces between layers) on the Pasternak and Winkler elastic foundations.

Derivation of the Governing Equations for the Bilayer Annular/Circular Nanoplate
Numerous graphene layers can be employed to fix their weak buckling strengths Therefore, graphene plates are set on top of each other by means of van der Waals bonds generating graphene layers [42].Van der Waals force is a general term used to explain the attraction of intermolecular forces among molecules (for more information please refer to Refs.[50,51]).Figure 2 shows a bilayer nanoplate (considering van der Waals forces between layers) on the Pasternak and Winkler elastic foundations.Equilibrium equations considering the bilayer annular/circular nanoplate with the Pasternak-Winkler elastic foundation can be procured similar to the procedure used for the single-layer nanoplate.Furthermore, the annular/circular displacement fields of the axisymmetric bilayer nanoplate can be written as (i = 1 shows the top layer and i = 2 illustrates the bottom layer) Equilibrium equations considering the bilayer annular/circular nanoplate with the Pasternak-Winkler elastic foundation can be procured similar to the procedure used for the single-layer nanoplate.Furthermore, the annular/circular displacement fields of the axisymmetric bilayer nanoplate can be written as (i = 1 shows the top layer and i = 2 illustrates the bottom layer) In addition, the strain equations are the same as those procured in the single-layer annular/circular nanoplate.It should be noted that considering the minimum potential energy method to achieve the equilibrium equations (as well as boundary conditions) of the top and bottom layers, the following equations are considered: k o is the van der Waals stiffness constant between layers.Furthermore, the equilibrium equations in terms of local stresses for both layers can be mentioned as follows: In addition, the stability equations can be presented for the thermal buckling of circular and annular plates.Moreover, the only difference between the equations of thermal buckling compared to mechanical buckling is in the strains.Consequently, the stress and moment resultants change.The nonlinear strains, taking into account von Karman's assumptions, can be written as Equations ( 56)-(58): The local stress resultants (considering thermal effects) according to displacements are as follows: In order to derive the stability equations from the equilibrium equations and determine the critical buckling temperature using the adjacent equilibrium technique, the equilibrium equations are obtained in the form of Equations ( 38)- (40).In those equations, by ignoring the terms that only have a zero subscript, zeroing the terms related to the preload state, and knowing that in thermal buckling Stability equations in terms of local stresses are written in the form of Equations ( 62)-(64): 2

Equations of Thermal Stability for a Bilayer Annular/Circular Plate
The stability equations in terms of local stresses (considering the first layer (i = 1) and the second layer (i = 2)) can be written as follows:

Boundary Conditions
The boundary conditions for the annular/circular plate can be assumed as Clamped (C): Simply supported (S): Free (F):

Nondimensional Assumptions
Because of the very small values that exist on the nanoscale, and for the ease of calculations, Equation ( 72) is introduced to convert the equations of the nanoplate nondimensional:

Numerical Procedure
DQM is one of the most efficient and elegant techniques for the numerical solution of partial differential equations used by many authors (which are briefly explained in the Introduction).This method is obtained from the quadratic integration technique, which expresses the integral at one point in the direction of the domain and relies on all points along that direction.Moreover, the value of dependency can be obtained using weight constants: where, w 1 , w 2 , . . ., w n as well as f 1 , f 2 , . . ., f n can be defined as weight constants and function values at discrete points, respectively.Belman et al. [52] (regarding quadratic integration) proposed that the derivative at one node of the function domain can be based on the function values at the entire nodes of the domain via weight constants.
in which A ij and N are the weight constant and the total number of nodes in the direction of r, respectively.For the first-order derivative, the weighting constants can be gained as For higher-order derivatives, the following equation can be defined as For the second-and higher-order derivatives, the weighting constants can be defined as follows: In this study, the grid points are distributed according to the Chebyshev-Gauss-Lubato distribution as follows: where r i and r o are nodes at the start and end of the function.
It should be noted that for each node (considering a single-layer plate) there are three equations and four unknowns, whose unknowns are displacement components of the plate (w 0 , ϕ, u 0 ) and critical buckling load in the case of mechanical/thermal axial symmetric analysis.In other words, there are 27 equations and 28 unknowns for 9 nodes.The desired unknown (mechanical buckling critical load or thermal buckling critical term) is obtained after applying the DQM.Therefore, by inserting Equations ( 25)- (32) (which are the resultants in terms of displacements) in the stability equations, the stability equations in terms of displacements are obtained.Then, by using Equation (74) (for first-order derivatives) or Equation (78) (for higher-order derivatives) and by using appropriate weighting constants (Equations (75)-(77) for the first-order derivative or Equations (79)-(80) for higher derivatives) the discretized form of equations will be obtained in terms of displacements.The matrix form of these equations can be written as follows: [Ce f ] is the matrix form of coefficients of equations.For N nodes, the matrix form of equations includes 3*N algebraic equations and 3*N + 1 unknowns.Since solving these algebraic equations and unknowns are long and time-consuming, Maple software 2023 was used to solve it.So, in the computer program, the determinant of the matrix of coefficients should be equal to zero (|Ce f | = 0) to obtain the dimensionless critical load of mechanical buckling or the dimensionless term of thermal buckling.Then, an equation will be obtained that contains a Polynomial in terms of buckling load.By solving it, and in order to find the critical buckling load, the minimum value should be chosen.Also, for a better understanding of applying DQM, an example of the discrete form of the equations can be seen in Appendix A.

Results and Discussion
In this part, various factors are studied to identify their effect on the nondimensional buckling load of the annular/circular nanoplate.For this purpose, HSDT, FSDT, and the nonlocal strain gradient theory are utilized with the aid of DQM.Moreover, to validate the procedure, the current results are compared with the results of references.Also, the following assumptions are considered (if not mentioned in the text): Figure 3 reveals the influence of the number of nodes employed in the DQM to achieve the results of this paper.As can be seen, after nine nodes the appropriate convergence is gained.So, the number of nine nodes is considered to obtain the numerical results.

Results and Discussion
In this part, various factors are studied to identify their effect on the nondimensional buckling load of the annular/circular nanoplate.For this purpose, HSDT, FSDT, and the nonlocal strain gradient theory are utilized with the aid of DQM.Moreover, to validate the procedure, the current results are compared with the results of references.Also, the following assumptions are considered (if not mentioned in the text): Figure 3 reveals the influence of the number of nodes employed in the DQM to achieve the results of this paper.As can be seen, after nine nodes the appropriate convergence is gained.So, the number of nine nodes is considered to obtain the numerical results.To verify accuracy and validity, the present results (gained using both FSDT and HSDT) in the clamped boundary condition are compared with the reference, which shows good accuracy (Table 2).In addition, it can be concluded that by increasing the radius of the plate, the nondimensional buckling loads are decreased.To verify accuracy and validity, the present results (gained using both FSDT and HSDT) in the clamped boundary condition are compared with the reference, which shows good accuracy (Table 2).In addition, it can be concluded that by increasing the radius of the plate, the nondimensional buckling loads are decreased.Table 3 compares the nondimensional buckling load gained by HSDT by considering various shape functions for a circular nanoplate at clamped as well as simply supported boundary conditions.Therefore, different functions are used to distribute the shear stress along the thickness, which can be summarized as: From Table 3, it can be deduced that using these shape functions results in obtaining almost the same values of dimensionless buckling loads for both clamped and simply supported boundary conditions.

Number of nodes
Figure 4 illustrates the changes in the nondimensional buckling loads versus the nonlocal constant for the annular nanoplate at various boundary conditions.It can be seen that by increasing the nonlocal parameter the non-dimensional buckling load decreases.This may be persuaded by increasing the nonlocal parameter; the stiffness of the nanoplate is increased and, therefore, the nondimensional buckling of the plate is reduced.Also, by increasing the rigidity of the boundary condition, increasing the nonlocal parameter results in a significant decrease in the nondimensional buckling load.Moreover, as the nonlocal coefficient increases the effect of the boundary condition on the nondimensional load decreases.
From Figures 5 and 6, it is distinguished that generally increasing the elastic foundation results in an increase in the nondimensional buckling load, which is due to the fact that assuming an elastic medium results in a stiffer structure.Also, it is noticed that in local analysis, the nondimensional buckling loads with the Winkler foundation are higher than values with the Pasternak foundation, but by increasing the nonlocal parameter the Pasternak foundation has more values of nondimensional buckling loads.Moreover, In the local analysis, by increasing the flexibility of the boundary condition, the nondimensional buckling load decreases, but it can be noticed in Figures 5 and 6, by increasing the rigidity of the boundary condition and increasing the nonlocal coefficient, the nondimensional buckling load decreases.This can highlight that in nonlocal analyses, the influence of the nonlocal coefficient on the buckling analysis of nanoplate is more significant than the type of boundary conditions or elastic foundations.
This may be persuaded by increasing the nonlocal parameter; the stiffness of the nanoplate is increased and, therefore, the nondimensional buckling of the plate is reduced.Also, by increasing the rigidity of the boundary condition, increasing the nonlocal parameter results in a significant decrease in the nondimensional buckling load.Moreover, as the nonlocal coefficient increases the effect of the boundary condition on the nondimensional load decreases.From Figures 5 and 6, it is distinguished that generally increasing the elastic foundation results in an increase in the nondimensional buckling load, which is due to the fact that assuming an elastic medium results in a stiffer structure.Also, it is noticed that in local analysis, the nondimensional buckling loads with the Winkler foundation are higher than values with the Pasternak foundation, but by increasing the nonlocal parameter the Pasternak foundation has more values of nondimensional buckling loads.Moreover, In the local analysis, by increasing the flexibility of the boundary condition, the nondimensional buckling load decreases, but it can be noticed in Figures 5 and 6, by increasing the rigidity of the boundary condition and increasing the nonlocal coefficient, the nondimensional buckling load decreases.This can highlight that in nonlocal analyses, the influence of the nonlocal coefficient on the buckling analysis of nanoplate is more significant than the type of boundary conditions or elastic foundations.Figures 7 and 8 illustrate the nondimensional buckling load versus the nonlocal parameter obtained from HSDT and FSDT for two thickness-to-radius ratios including h/r = 0.05 and h/r = 0.06.As can be deduced by increasing the nonlocal parameter, the nondimensional buckling load obtained from HSDT and FSDT differ significantly.Also, by comparing these two figures, it can be seen that in the higher thickness-to-radius ratios the difference between FSDT and HSDT increases more significantly.Figures 7 and 8 illustrate the nondimensional buckling load versus the nonlocal parameter obtained from HSDT and FSDT for two thickness-to-radius ratios including h/r = 0.05 and h/r = 0.06.As can be deduced by increasing the nonlocal parameter, the non-dimensional buckling load obtained from HSDT and FSDT differ significantly.Also, by comparing these two figures, it can be seen that in the higher thickness-to-radius ratios the difference between FSDT and HSDT increases more significantly.Figure 9 shows the nondimensional buckling load versus the radius ratios (inner radius-to-outer radius ratios) for different boundary conditions and considering different values of strain gradients.It can be concluded that by increasing radius ratios the buckling loads decrease.Also, increasing the strain gradient results in a decrease in buckling loads which specifies that inclusion of the strain gradient coefficient makes the plate stiffer than that of the classical plate [54].Figure 9 shows the nondimensional buckling load versus the radius ratios (inner radius-to-outer radius ratios) for different boundary conditions and considering different values of strain gradients.It can be concluded that by increasing radius ratios the buckling loads decrease.Also, increasing the strain gradient results in a decrease in buckling loads which specifies that inclusion of the strain gradient coefficient makes the plate stiffer than that of the classical plate [54].

Thermal Analysis
In this section, to validate the research, the current results of thermal buckling of circular/annular plates are compared with the results of valid reference, and then the results of the relevant diagrams are presented.Figure 10

Thermal Analysis
In this section, to validate the research, the current results of thermal buckling of circular/annular plates are compared with the results of valid reference, and then the results of the relevant diagrams are presented.Figure 10 reveals the effect of the number of nodes employed in the DQM to achieve the results of thermal analysis.It is noted that T * is the critical buckling temperature, which is considered as follows: of the relevant diagrams are presented.Figure 10 reveals the effect of the number of n employed in the DQM to achieve the results of thermal analysis.It is noted that T the critical buckling temperature, which considered as follows: T Δ is the difference between the critical buckling temperature and the refer temperature.The coefficient of thermal expansion is considered as ∆T is the difference between the critical buckling temperature and the reference temperature.The coefficient of thermal expansion is considered as In order to validate the numerical results of the present article and the solution technique and compare them with the results of Ref. [55], the dimensionless thermal buckling parameter (λ T ) is defined as follows: In this comparison, the dimensionless parameter of the thermal buckling of the isotropic graphene plate for the annular nanoplate assuming various boundary conditions, radius ratios, and thickness-to-radius ratios is given in Table 4.According to this table, it is clear that in smaller thickness-to-radius ratios, the results agree well with the results of Ref. [43], and as the h/r o and r i /r o ratios increase, a small amount of difference is observed.Moreover, the difference between the results is caused by the different types of theories.
To better observe the differences between the present study with the reference [55], the percentage errors are gathered in Table 5, which can be defined as follows: where E 1 (%) and E 2 (%) represent the HSDT and FSDT errors compared to the reference's values, respectively.Also, E HSDT , E FSDT , and E Ref represent values of the dimensionless thermal buckling parameter (λ T ) gained by the present study HSDT, FSDT, and the reference (obtained from Table 4), respectively.Figure 11 shows the critical buckling temperature of the annular nanoplate versus the nonlocal parameter for different boundary conditions.As can be distinguished, the critical buckling temperature decreases with the rise of the nonlocal parameter (which has the same trend as the mechanical analysis).
Finally, by comparing the equations and checking the results for both thermal and mechanical analyses, it is found that with the following equation, thermal and mechanical buckling analyses can be related to each other: It is noted that in Equation (93), T * is critical buckling temperature and N * is mechanical load.In other words, the results of the thermal buckling can be gained using Equation ( 93 Figure 11 shows the critical buckling temperature of the annular nanoplate versus the nonlocal parameter for different boundary conditions.As can be distinguished, the critical buckling temperature decreases with the rise of the nonlocal parameter (which has the same trend as the mechanical analysis).Finally, by comparing the equations and checking the results for both thermal and mechanical analyses, it is found that with the following equation, thermal and mechanical buckling analyses can be related to each other:

Bilayer Analysis
Figures 12-15 demonstrate the variation in the critical buckling temperature in terms of nonlocal constants for the annular bilayer nanoplate with a thickness of h (0.34 nm) and the single-layer nanoplate considering a thickness of 2*h (0.68 nm) for different boundary conditions.In fact, this part studies the probability of replacing a single-layer nanoplate (considering double thickness) with a bilayer nanoplate.It can be concluded that a bilayer nanoplate cannot be considered a single-layer nanoplate with the same thickness (as the bilayer nanoplate) for obtaining the exact results of the buckling analysis.Furthermore, the results obtained from the equivalent single-layer nanoplate with double thickness overestimate the real results gained using a bilayer nanoplate.This can be justified because of the van der Waals interaction that exists between the layers of bilayer nanoplate, which is an approximately weak force.If the buckling loads are strong enough, then molecules break free of the van der Waals forces that hold them together.Thus, the results of a single-layer plate are not similar to those of a bilayer plate.

Bilayer Analysis
Figures 12-15 demonstrate the variation in the critical buckling temperature in terms of nonlocal constants for the annular bilayer nanoplate with a thickness of h (0.34 nm) and the single-layer nanoplate considering a thickness of 2*h (0.68 nm) for different boundary conditions.In fact, this part studies the probability of replacing a single-layer nanoplate (considering double thickness) with a bilayer nanoplate.It can be concluded that a bilayer nanoplate cannot be considered a single-layer nanoplate with the same thickness (as the bilayer nanoplate) for obtaining the exact results of the buckling analysis.Furthermore, the results obtained from the equivalent single-layer nanoplate with double thickness overestimate the real results gained using a bilayer nanoplate.This can be justified because of the van der Waals interaction that exists between the layers of bilayer nanoplate, which is an approximately weak force.If the buckling loads are strong enough, then molecules break free of the van der Waals forces that hold them together.Thus, the results of a single-layer plate are not similar to those of a bilayer plate.

Conclusions
In this paper, the nonlinear mechanical and thermal buckling of the circular and annular (single layer/bilayer) nanoplate using the nonlocal strain gradient theory consider-

Conclusions
In this paper, the nonlinear mechanical and thermal buckling of the circular and annular (single layer/bilayer) nanoplate using the nonlocal strain gradient theory considering both HSDT and FSDT is carried out.DQM is employed to solve the governing equations of the annular/circular nanoplate.Also, the results were compared with other references that showed good agreement.Some of the results of this article can be summarized as follows:

•
In the higher thickness-to-radius ratios, the difference between FSDT and HSDT increases more significantly.

•
The influence of the nonlocal coefficient on the buckling analysis of circular/annular nanoplate is more significant than the types of boundary conditions or elastic foundations.

•
A bilayer nanoplate cannot be considered an equivalent single-layer nanoplate (with the same thickness as the bilayer nanoplate) to obtain the results of the buckling analysis of nanoplate.Furthermore, the results obtained from the equivalent singlelayer nanoplate with double thickness overestimate the real results gained using a bilayer nanoplate.
Author Contributions: M.S.: Conceptualization, methodology, software, validation, formal analysis, investigation, and writing: preparation of original draughts.A.P.: Supervision, conceptualization, methodology, and formal analysis.G.J.: validation investigation.All authors have read and agreed to the published version of the manuscript.
Funding: This research received no external funding.
Data Availability Statement: Data sharing is unavailable due to privacy or ethical restrictions.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A
As an example, in this part, the discrete forms of the dimensionless equations of stability for circular nanoplate on the Winkler elastic foundation considering FSDT using DQM are considered.The obtained FSDT stability equations in terms of displacements are as follows:

Figure
Figure 1a,b demonstrate an annular graphene plate with inner radius r i , outer radius r o , and constant thickness h on the Winkler-Pasternak elastic foundation (in which k p and k w define the Pasternak and Winkler elastic foundations) and under uniform extended buckling loads.

Figure
Figure 1a,b demonstrate an annular graphene plate with inner radius i r , outer ra-

Figure 2 .
Figure 2. A schematic of the bilayer nanoplate considering the Pasternak and Winkler elastic foundations.

Figure 2 .
Figure 2. A schematic of the bilayer nanoplate considering the Pasternak and Winkler elastic foundations.

Figure 3 .
Figure 3.The effect of the number of nodes on the dimensionless mechanical buckling load.

Figure 3 .
Figure 3.The effect of the number of nodes on the dimensionless mechanical buckling load.

Figure 4 .
Figure 4.The nondimensional buckling load of the annular nanoplate with respect to the nonlocal parameter in diverse boundary conditions.

Figure 4 . 29 Figure 5 .
Figure 4.The nondimensional buckling load of the annular nanoplate with respect to the nonlocal parameter in diverse boundary conditions.Micromachines 2023, 14, x FOR PEER REVIEW 18 of 29

Figure 5 .
Figure 5. Nondimensional buckling loads of the annular nanoplate with respect to the nonlocal parameter for different values of elastic foundations at various boundary conditions (including C-C and S-S).

Figure 5 .
Figure 5. Nondimensional buckling loads of the annular nanoplate with respect to the nonlocal parameter for different values of elastic foundations at various boundary conditions (including C-C and S-S).

Figure 6 .
Figure 6.The nondimensional buckling load of the annular nanoplate with respect to the nonloca parameter for different values of elastic foundations at various boundary conditions (including C-S and S-C).

Figure 6 .
Figure 6.The nondimensional buckling load of the annular nanoplate with respect to the nonlocal parameter for different values of elastic foundations at various boundary conditions (including C-S and S-C).

Micromachines 2023 , 29 Figure 7 .Figure 7 .
Figure 7.The nondimensional buckling load of the annular nanoplate with respect to the nonlocal parameter for FSDT and HSDT considering h/r = 0.05 for the simply supported boundary condition.

Figure 7 .
Figure 7.The nondimensional buckling load of the annular nanoplate with respect to the nonlocal parameter for FSDT and HSDT considering h/r = 0.05 for the simply supported boundary condition.

Figure 8 .
Figure 8.The nondimensional buckling load of the annular nanoplate with respect to the nonlocal parameter for FSDT and HSDT considering h/r = 0.06 for the simply supported boundary condition.

Figure 8 .
Figure 8.The nondimensional buckling of the annular nanoplate with respect to the nonlocal parameter for FSDT and HSDT considering h/r = 0.06 for the simply supported boundary condition.

Figure 9 .
Figure 9.The nondimensional buckling load of the annular nanoplate versus the radius ratios.
reveals the effect of the number of nodes employed in the DQM to achieve the results of thermal analysis.It is noted that * T is the critical buckling temperature, which is considered as follows:

Figure 9 .
Figure 9.The nondimensional buckling load of the annular nanoplate versus the radius ratios.

1 1Figure 10 .Figure 10 .
Figure 10.The critical buckling temperature for the annular nanoplate with respect to the nu of nodes.

Figure 11 .
Figure 11.The critical buckling temperature for the annular nanoplate versus the nonlocal parameter.

Figure 11 .
Figure 11.The critical buckling temperature for the annular nanoplate versus the nonlocal parameter.

Figure 12 .Figure 12 .
Figure 12.Critical buckling temperature for single-layer and bilayer annular nanoplate with respect to the nonlocal parameter for the C-C condition.

Figure 13 .
Figure 13.Critical buckling temperature for single-layer and bilayer annular nanoplate versus the nonlocal parameter for S-S condition.

Figure 14 .Figure 13 . 29 Figure 13 .
Figure 14.Critical buckling temperature for single-layer and bilayer annular nanoplate versus the nonlocal parameter for S-C condition.

Figure 14 .Figure 14 . 29 Figure 15 .
Figure 14.Critical buckling temperature for single-layer and bilayer annular nanoplate versus the nonlocal parameter for S-C condition.

Figure 15 .
Figure 15.Critical buckling temperature for single-layer and bilayer annular nanoplate versus the nonlocal parameter for C-S condition.

Table 1 .
Recommended functions used in HSDT by various authors. 2

Table 1 .
Recommended functions used in HSDT by various authors.

Table 2 .
Comparison of the nondimensional buckling loads of the circular nanoplate gained by the present study with a reference.

Table 2 .
Comparison of the nondimensional buckling loads of the circular nanoplate gained by the present study with a reference.

Table 3 .
Comparison of the dimensionless buckling load gained by considering various functions of HSDT for a circular nanoplate under both simply supported and clamped boundary conditions.

Table 4 .
Comparison of the dimensionless thermal buckling parameter (λ T ) gained by the present study (HSDT and FSDT) with the reference for the annular nanoplate under clamped boundary conditions.

Table 5 .
[55]differences in percentage errors between the present study with the Ref.[55].
) which can save the time of computations.