1 Introduction

The main motivation of the present paper is to address the issues existing for studying the size-dependent structures at small dimensions within the framework of the nonlocal strain gradient theory proposed by Lim et al.^[1]. These are (ⅰ) what is the higher-order shear resultant expressed in terms of displacement(s), and (ⅱ) whether the results for the same nonlocal parameter and strain gradient parameter can be reduced to the classical solutions for the general boundary value problems.

With a view toward characterizing the softening and stiffening behaviors of structures, an elastic theory considering both the cohesive forces at long range between atoms or molecules and atoms with higher-order deformation mechanism at small dimensions is established. Recently, numerous works have been published with the particular emphasis on the static and dynamic problems of rods^[2], beams^{[1, 3-12]}, plates^[13-16], and shells^[17-18] with their characteristic sizes at small dimensions. In these works, the higher-order boundary conditions are obtained through replacing the classical stress resultants by the nonclassical stress resultants in the balance equations^[19]. Generally speaking, equations of motion of the system can be easily obtained whereas the corresponding boundary conditions may not be physically consistent. The nonclassical shear resultant, Q_n, based on the Timoshenko beam theory is given by^{[3, 20]}

(1)

where f₁ and f₂ are functions of nonlocal parameter l_e, strain gradient parameter l_m, coordinate x, time t, and displacement w. Obviously, the paradox in Eq. (1) is that Q_n becomes an infinite value for l_e = l_m. This paradox exists and it has not attracted deserved attention in the research communities except for the recent work by El-Borgi et al.^[21] who found that the stress resultants are only defined for different values of the nonlocal and strain gradient parameters. Following the procedure for obtaining the higher-order shear resultant, Q_n, similar inconsistence can also be found for the higher-order moment of torque in Ref. [22].

Recently, the integral type of nonlocal elasticity appears to be a useful way to address the inconsistences encountered^[23-24]. Several continuous works have been reported. The closed-form bending deflections of thin beams with various different boundary conditions are presented by Fernández-Sáez et al.^[25], Tuna and Kirca^[26], and Meng et al.^[27]. The buckling and free vibration analyses were studied by Zhu et al.^[28] and Eptaimeros et al.^[29]. Subsequently, the well-posed boundary value problems of nonlocal integral elasticity have been mentioned by numerous researchers^{[24, 30]}. These works reveal that the inconsistences can be overcome by using the nonlocal integral elasticity. Additionally, the nonlocal integral models predict softening effects when they are compared with the corresponding differential ones.

Based on the previous works^[31-32] in addressing the dynamic softening phenomena in nonlocal cantilevers and the buckling behaviors of nonlocal beams undergoing different boundary conditions, the technical method adopted in the present paper is termed as the weighted residual method, which allows one to yield the variationally consistent boundary conditions for the known equations of motion of the system. It will be seen in Subsection 3.2 that the present variationally consistent boundary conditions can avoid the first paradox mentioned above. Additionally, the variational principles have also been utilized by Mousavi and Paavola^[33], Yaghoubi et al.^[34], and Tahaei-Yaghoubi et al.^[35].

The second paradox is that whether the results for the same nonlocal parameter and strain gradient parameter can be reduced to the classical solutions for the general boundary value problems. It is concluded for the same values of the two material length parameters that analytical solutions and numerical results of the nonlocal strain gradient models recover those of classical models^{[4, 20, 22]}. The basis for this conclusion is that, when one uses the classical mode shapes for simply-supported (SS) boundary conditions, their results can be reduced to the classical solutions. However, since the differential order of equations of motion within the framework of the nonlocal strain gradient theory increases, the extra higher-order boundary conditions will simultaneously exist to constitute a well-posed boundary value problem^[36]. In other words, the different selections of the higher-order boundary conditions induced by the strain gradient may vary the numerical results even for structures subject to the given boundary conditions^[37-40]. Since the nonlocal strain gradient theory is the general form of the strain gradient theory, the corresponding boundary value problems may have different numerical results for structures with prescribed constrains. It has been demonstrated by Xu et al.^[32] for studying the bending and buckling of nonlocal strain gradient Euler-Bernoulli beams. In their work, they reformulated the boundary value problems of nonlocal strain gradient beams via the weighted residual method, and studied effects of the two material length parameters on the buckling of Euler-Bernoulli beams. In contrast to the previous works^{[4, 20, 22]}, it is shown that the numerical results of nonlocal strain gradient beams cannot always recover those of the classical models due to the alternative selections of the higher-order boundary conditions.

The structure of this paper is as follows. Section 2 presents a short review of the nonlocal strain gradient theory. Section 3 gives the boundary value problems of nonlocal strain gradient Timoshenko beams in two different ways. The first way is the usual formulation by numerous researchers that the first paradox can be observed. To solve the first paradox, the second way derives the variationally consistent boundary conditions by using the weighted residual method. It is found that the differential order of equations of motion of the system and their variationally boundary conditions will not be varied for the same values of the two nonzero material length parameters, indicating that the numerical results for this case will not be, in general, reduced to those of classical models. After that, the analytical solutions for buckling of Timoshenko beams subject to three common boundary conditions are determined in Section 4. The effects of aspect ratios, two material length parameters, and boundary conditions on the buckling loads of nonlocal strain gradient Timoshenko beams are investigated in Subsection 5.1. Subsection 5.2 proposes three semi-empirical formulae for capturing the buckling loads of nonlocal strain gradient Timoshenko beams. Finally, the main results of the present paper are summarized in Section 6.

2 Nonlocal strain gradient theory

Based on the nonlocal strain gradient theory established by Lim et al.^[1], the total strain energy of an elastic body occupying domain V may be given by

(2)

where and are the strains at the points x and y, C_ijkl is the conventional elastic tensor, l_e0 and l_e1 are the nonlocal parameters, l_m is the material length scale parameter, is the attenuation kernel function capturing the nonlocal effect, and all the other symbols follow the usual notation unless otherwise specified.

Then, by using the variational principle with respect to Eq. (2), Lim et al.^[1] derived the following constitutive relation:

(3)

with t_ij the total stress tensor of the nonlocal strain gradient theory.

The stress tensor σ_ij and the higher-order stress tensor are defined by

(4)

(5)

In general, great challenges will be encountered when one solves the above integral constitutive equations. Therefore, by assuming , Lim et al.^[1] rewrote Eq. (3) into a simple differential constitutive equation,

(6)

For the static one-dimensional problems, Eq. (6) further degenerates to

(7)

where E is the Young's modulus of the body.

It is noted that Eq. (7) is the combined form of the nonlocal constitutive equation^[41] and the strain gradient constitutive equation^[42-44]. As a result, it has the spectacular advantages of characterizing the materials that both exhibit the softening effect and stiffening effect.

3 Boundary value problems of Timoshenko beams 3.1 Review of boundary value problems of Timoshenko beams

In this section, we briefly review the main procedures given by Li and Hu^[3] for obtaining the boundary value problems of Timoshenko beams based on the nonlocal strain gradient theory developed by Lim et al.^[1].

Consider a straight, homogeneous, isotropic, thick circular beam of the length L subject to the lateral load q(x) and the axial load f(x), as shown in Fig. 1. The x-axis is taken along the geometric centroid of the beam cross section, and the z-axis is taken positive downward. Based on the conventional beam theory, the displacement field () along the xz-plane can be expressed by

(8)

where u and w are the axial and transverse displacements of a point on the mid-plane of the beam, respectively, ϕ is the rotation angle, and the prime denotes the differentiation with respect to x.

The von Kármán nonlinear strains for a beam under large displacements can be given by Reddy^[45],

(9)

where and are the normal and shear strains, respectively.

The boundary value problems of nonlocal strain gradient Timoshenko beams were presented by Li and Hu^[3]. For completeness of the context, we summarize their main results below. It was assumed in their paper that the derivation with respect to the thickness direction was omitted for simplicity.

The virtual work of the Timoshenko beams within the framework of the nonlocal strain gradient theory is given by

(10)

which, by using the integration by parts, becomes

(11)

where V and A are the volume and the cross-sectional area of the beam, respectively. The various stress resultants defined in Eq. (11) are given by

(12)

The virtual work done by the external forces is expressed by

(13)

where we have assumed that there is no work done by other higher-order forces acting on the boundaries.

The principle of virtual displacements states, for an elastic body, that the total virtual work of the system must be zero,

(14)

In view of Eqs. (11), (13), and (14), after integrating by parts with respect to x, we have the following expressions:

(15)

and the boundary conditions at the ends of the beam as

(16a)

(16b)

(16c)

(16d)

(16e)

(16f)

It is noted that the boundary conditions (16) do not coincide with those given by Li and Hu^[3]. More precisely, the stress results N_c and N_n given by Eqs. (16c) and (16d) are not included in their boundary conditions.

By taking Eq. (9) into account, the constitutive equations within the framework of the nonlocal strain gradient theory can be expressed in terms of displacements as

(17)

where G is the shear modulus of the material.

It immediately follows from Eq. (12) that

(18)

(19)

(20)

where I is the second moment of inertia of the beam about the neutral axis, and K is the shear correction factor of the Timoshenko beam.

In view of Eq. (15) and Eqs. (18)-(20), we have

(21)

(22)

(23)

By substituting Eqs. (21)-(23) into Eq. (15), we obtain the equations of motion in terms of displacements as

(24)

(25)

(26)

Equations (24)-(26) are the general governing equations of Timoshenko beam models based on the nonlocal strain gradient theory. They can be reduced to the governing equations of nonlocal Timoshenko beam models (e.g., ) and those of strain gradient Timoshenko beam models (e.g., ). Recently, Li and Hu^[3] derived the displacement forms of higher-order shear stress resultant Q_n using straightforwardly algebraic manipulations. Their formulation, however, appears to fail when the two material length parameters are the same. More precisely, it is found from Eq. (1) that the higher-order shear stress resultant for . In other words, they have not presented a well-posed boundary value problem of nonlocal strain gradient beams. To clarify this issue, we are encouraged to provide the variationally consistent boundary value problems, as detailed in Subsection 3.2.

3.2 Variationally consistent boundary conditions

The weighted residual method allows us to formulate the variationally consistent boundary conditions provided that the governing equations of motion of the system are known in advance. This motivation of this subsection is to derive the corresponding variationally consistent boundary conditions of nonlocal strain gradient Timoshenko beams.

By making use of the following classical stress resultants:

(27)

(28)

(29)

and Eqs. (24)-(26), we can easily rewrite equations of motion in terms of the classical stress resultants as

(30)

(31)

(32)

Next, we will derive the variationally consistent boundary conditions of the equations of motion for nonlocal strain gradient Timoshenko beams.

Firstly, we can easily obtain the following weak form of Eqs. (30)-(32):

(33)

Then, by integrating by parts with respect to Eq. (33), we obtain

(34)

The virtual strain energy expressed in Eq. (34) is defined by

(35)

The virtual work done by the external forces is defined by

(36)

The variationally consistent boundary conditions are

(37)

where N and N_h are the force resultant and the higher-order force resultant, respectively, Q, M, Q_h, and M_h are the shear force, moment, higher-order moment, and higher-order shear force resultants, respectively.

By using Eqs. (21) and (27)-(29), the boundary conditions given above may be rewritten in the following forms:

(38)

It is emphasized from the boundary conditions (38) that either stress resultants Q_h, M_h) or its work-conjugate displacements should be prescribed. As a result, one may encounter different selections of higher-order boundary conditions for beams with typical constrained ends. As the motivation of the present work, effects of different selections of higher-order boundary conditions on the critical buckling load are studied in Section 5. In addition, the above variationally consistent boundary conditions contain, as simplified cases, various Timoshenko beam models. For instance, when the strain gradient parameter l_m vanishes, the boundary conditions above can be reduced to those of nonlocal Timoshenko beam models^{[31, 46]}. On the other hand, the boundary conditions above can be reduced to those of strain gradient Timoshenko beam models if the nonlocal parameter is zero^[47-48]. For the comparison purpose, we present the main differences of the boundary stress resultants for the nonlocal strain gradient Timoshenko beams in Table 1. As expected, the present results can be reduced to those of nonlocal models and those of strain gradient models.

4 Buckling solutions

In this section, we derive the closed-form solutions for buckling of nonlocal strain gradient Timoshenko beams. For the buckling case, we take and . The governing equations (25) and (26) can be rewritten as

(39)

Similarly, it follows from Eq. (38) that the variationally consistent boundary conditions are

(40)

It can be inferred from Eq. (40) that the general solutions of Eq. (39) may contain 8 unknown integration constants for boundary value problems of beams, as shown in Subsections 4.1 and 4.2.

4.1 Nondimensional boundary value problems

For simplicity of the illustration, we define the following nondimensional parameters:

(41)

Then, Eq. (39) may be converted into the following nondimensional equations of motion based on the Timoshenko beam theory:

(42a)

(42b)

In this paper, spatial derivatives of nondimensional parameters denote results with respect to the nondimensional coordinate X.

Analogously, it follows from Eq. (40) that the variationally consistent nondimensional boundary conditions are

(43a)

(43b)

(43c)

(43d)

It is noted that Eqs. (42a) and (42b) are coupled with the displacements W and ϕ, and as a result, we should decouple Eqs. (42a) and (42b). To this end, we can do the the following steps: (ⅰ) we solve from Eq. (42a); (ⅱ) we differentiate Eq. (42b) with respect to X; (ⅲ) we substitute and from Step (ⅰ) into the resulting expression of Step (ⅱ); (ⅳ) we rearrange the expression of Step (ⅲ) and yield the following equations:

(44a)

(44b)

where

(45)

By taking , it is interestingly found that Eq. (44) coincides with Eq. (47) of Xu et al.^[31].

4.2 Analytical solutions

For the decoupled set of the ordinary differential equations (44a) and (44b), the general solutions of the nondimensional displacements W(X) and ϕ(X) depend on 15 integration parameters (C_i, (i = 1, 2, · · ·, 8) and D_j (j = 1, 2, · · ·, 7)),

(46)

where

(47)

However, not all the above 15 constants are independent of each other. In fact, their relationships can be found by substituting Eq. (46) into Eq. (42b), namely,

(48)

where

(49)

It is worth mentioning that the last relation in Eq. (48) can only be obtained by inserting Eq. (46) into Eq. (42a). It is noted that Eq. (48) greatly simplifies the boundary value problems for Timoshenko beams, as shown below.

4.2.1 CC boundary conditions

For the CC boundary conditions, the lower-order boundary conditions require that the displacement W(X) and the rotation angle ϕ(X)should be zero at the fixed ends of the beam,

(50)

The higher-order boundary conditions may have 4 different selections.

(Ⅰ) The first selection (BC1) we consider is

(51)

Equation (51) means that the displacement and rotation slopes (e.g., W ′ and φ) are prescribed. In other words, the higher-order shear stress resultant and the higher-order moment resultant may not necessarily vanish.

By substituting the general solutions (46) into the boundary conditions in Eqs. (50) and (51), and remembering Eq. (48), we obtain the following set of linear equations in terms of the coefficients C_i (i = 1, 2, · · ·, 8):

(52)

For non-trivial solutions, the determinant of Eq. (52) must be zero. Then, we have

(53)

which is a transcendental equation, and the infinite numerical solutions may be found. The smallest positive value of Eq. (53) is termed in the present work as the critical buckling load.

(Ⅱ) The second selection of the higher-order boundary conditions (BC2) is

(54)

Similarly, the following condition is obtained for non-trivial solutions:

(55)

where

(56)

(Ⅲ) The third selection of the higher-order boundary conditions (BC3) is

(57)

By combining Eqs. (57) and (43c), we get the following equivalent equation:

(58)

The boundary value problem for this case is determined by the following characteristic equation:

(59)

(Ⅳ) The fourth selection of the higher-order boundary conditions (BC4) is

(60)

Then, we arrive at the following characteristic equation:

(61)

4.2.2 SS boundary conditions

For nonlocal strain gradient Timoshenko beams with the SS boundary conditions, there also exist 4 possible boundary value problems. The lower-order boundary conditions are given by

(62)

which indicates that the displacement and moment should vanish at the constrained ends of the beam.

(Ⅰ) The first selection of the higher-order boundary conditions (BC1) is

(63)

In view of the first equation of Eq. (63), the second equation of Eq. (63) can be simplified as

(64)

Then, by combining Eqs. (62), (64), and the first equation of Eq. (63), we have the following simplified boundary conditions:

(65)

The critical buckling load can be determined by solving the following characteristic equation:

(66)

The above expression is equivalent to

(67)

which, together with the second equation of Eq. (47), yields the analytical solution for the buckling load,

(68)

Note that Eq. (68) contains various buckling results. For nonlocal Timoshenko beams, it can be reduced to the results of Xu et al.^[31] and Zhang et al.^[46] in the absence of the strain gradient parameter . On the other hand, the above equation is equivalent to that of Wang et al.^[49]. Moreover, it can be concluded from Eq. (68) that the inclusion of the strain gradient parameter has the ability to increase the buckling load. However, for the nonlocal parameter and the shear deformation parameter s, it has the opposite tendency. Here and in the following, we choose the physical boundary conditions based on the concept that the buckling loads increase as s increases.

(Ⅱ) The second selection of the higher-order boundary conditions (BC2) is

(69)

Following the similar procedure as that of BC1, we can obtain the following characteristic equation:

(70)

where some parameters are defined in Eq. (56), and the undefined parameters are given by

(71)

(Ⅲ) The third selection of the higher-order boundary conditions (BC3) is given by

(72)

The corresponding boundary value problem enables one to obtain the following characteristic equation:

(73)

(Ⅳ) The last selection of the higher-order boundary conditions (BC4) is expressed by

(74)

Then, we can arrive at the following expressions for buckling loads of beams:

(75)

4.2.3 CF boundary conditions

For cantilever beams with CF boundary conditions, we assume that the left end is clamped at X = 0 and the right end is free at X = 1. The lower-order boundary conditions of Timoshenko beams are

(76)

There may exist 16 (= 4 × 4) possible higher-order boundary conditions for beams with CF boundary conditions. As a result, choosing the physically reasonable boundary conditions that agree with the intuitional results is of great significance.

(Ⅰ) The first selection of the higher-order boundary conditions (BC1) is

(77)

In view of Eqs. (46), (76), and (77), we have the following characteristic equation:

(78)

where

(79)

It is emphasized that in simplifying the fifth row of Eq. (78), we have used Eq. (49) and the following relations:

(80)

(Ⅱ) The second selection of the higher-order boundary conditions (BC2) is

(81)

The critical buckling load can be determined from

(82)

(Ⅲ) The third selection of the higher-order boundary conditions (BC3) is

(83)

which can be simplified by considering Eqs. (43a)–(43d) as

(84)

Inserting Eq. (46) into Eqs. (76) and (84) allows one to have the following characteristic equation:

(85)

(Ⅳ) For the fourth selection of the higher-order boundary conditions (BC4), we have

(86)

Following the similar procedure, we have

(87)

from which the buckling loads of the beam can be numerically solved.

The remaining 12 selections of the higher-order boundary conditions are summarized in Table 2. Analogously, these 12 boundary value problems result in 12 characteristic equations that are not presented in this paper for brevity.

5 Numerical results and discussion

In this section, the critical buckling loads are numerically solved based on the expressions derived in Section 4. Since the present formulations for Timoshenko beams are based on the nonlocal strain gradient theory which, as special cases, can degenerate to either the corresponding nonlocal ones or the strain gradient ones, the results are presented in the three theories as follows: (ⅰ) the nonlocal theory by taking ^[31]; (ⅱ) the strain gradient theory by setting ; (ⅲ) the nonlocal strain gradient theory for . In the numerical simulations, we consider a circular cross-sectional beam with

(88)

where is Poisson's ratio, and d is the diameter of the beam.

5.1 Numerical results and discussion

As is well-known, the aspect ratio has an important effect on the buckling behaviors of Timoshenko beams. Therefore, we make the first attempt to show how the critical buckling loads vary with respect to the aspect ratios for beams with all the possible boundary conditions, as presented in Fig. 2. In Fig. 2(a), numerical results for all the other 12 boundary value problems show the similar tendency as BC1 or have no numerical values. For BC4 in Fig. 2(b), the numerical results cannot be obtained. However, for BC4 in Fig. 2(c), the determinant is uniformly zero. Consequently, these mentioned cases are not depicted in Fig. 2. Interestingly, it is found from Fig. 2 that three distinct behaviors can be clearly observed. For small values of aspect ratios, it is seen from Fig. 2 that the critical buckling loads increase or decrease rapidly with increasing aspect ratios. For large values of aspect ratios, the curves for CF and SS boundary conditions approach the same constant values, but the curves for CC boundary conditions converge to different constant values. Moreover, our numerical results show that different selections of the higher-order boundary conditions may yield identical numerical results (e.g., BC3 or BC4 in Fig. 2(a)). As benchmark results, Table 3 shows the critical buckling loads versus aspect ratios of the Timoshenko beams.

It is accepted from the mechanics of structures that buckling loads increase as the aspect ratio of the beam increases, and they finally converge to those of Euler-Bernoulli beam models^[45]. We therefore follow the notation to identify the physically acceptable higher-order boundary conditions. It is observed from Fig. 2 that the physically acceptable higher-order boundary conditions we identify are BC3, BC1, and BC1-BC3 for beams with CF boundary conditions, SS boundary conditions, and CC boundary conditions, respectively.

Recently, it has been demonstrated that, for the same material length parameters, the buckling^{[14, 50]} and vibration^{[4, 20]} behaviors of structures may recover the classical cases. However, the present analytical solutions for buckling problems of nonlocal strain gradient Timoshenko beams formulated in Section 4 do not agree with the previous works^{[4, 14, 20, 50]}. Therefore, it is necessary to study the effects of the selected higher-order boundary conditions on the critical buckling loads of Timoshenko beams.

Firstly, we make our first attempt by comparing the present results with those reported for Euler-Bernoulli beam models with a view toward demonstrating the differences between the Timoshenko beam models and Euler-Bernoulli beams (see Table 4). It should be mentioned that the results of Akgöz and Civalek^[51] are obtained based on the strain gradient Euler-Bernoulli beam models. As observed from this table, the following conclusions can be drawn.

(ⅰ) The buckling loads calculated by the present work are lower than those predicted by the corresponding Euler-Bernoulli beam models, as expected.

(ⅱ) The differences between the present work and those reported in the literature are prominent for higher buckling modes.

Having demonstrated the differences of different beam models, we next focus on the critical buckling of Timoshenko beams subject to three typical boundary conditions. Unless otherwise stated, we take the aspect ratio L/d = 20.

The critical buckling loads of nonlocal strain gradient Timoshenko beams versus the nonlocal parameter are depicted in Fig. 3 with several different selections of the higher-order boundary conditions. The numerical results show that the critical buckling loads decrease as the nonlocal parameter increases. It means that the nonlocal parameter has the softening effect. It is also found from Fig. 3 that the higher-order boundary conditions for the beams with the CC boundary conditions have a more profound effect than the other two boundary conditions.

Figure 4 displays the influences of the selected higher-order boundary conditions and strain gradient parameters on the critical buckling loads of nonlocal strain gradient Timoshenko beams. It is seen that the critical buckling loads increase with increasing material length parameters, demonstrating that the present model can capture the stiffening effect.

5.2 Formulae

Subsection 5.1 shows that solving the complex boundary value problems of beams within the framework of the nonlocal strain gradient theory is a challenging task. Therefore, it is necessary to develop the approximately analytical formulae. Since the present size-dependent models may be reduced to those of the conventional ones when the material length parameters are not included in the constitutive equations, the corresponding buckling loads will also converge to the conventional ones. In view of the simplified expression (68) for beams with the SS boundary conditions, we take

(89)

as the semi-empirical formula for beams with the CF boundary conditions, where B is the curve fitting parameter.

The numerical value for B in Eq. (89) is 2.467, 4 by curve fitting with data plotted in Figs. 2(a), 3(a), and 4(a). As can be seen, the fitted curves and the numerical results are all in excellent agreement, demonstrating the efficiency of the proposed formulae. Moreover, the fitting parameter is found to be independent of material length parameters and the shear deformation parameter s

For beams with the CC boundary conditions, two independent formulae may be possible. One is given by

(90)

The other one is given by

(91)

In the above formulae, A and B are the curve fitting parameters to be determined.

Figures 2(c), 3(c), and 4(c) display the fitting curves. Good agreement between the numerical results and the proposed formulae can be found. For benchmarks, Table 5 tabulates the fitting parameters depicted in Figs. 2–4.

In order to verify the proposed formulae, we depict the critical buckling parameter versus the aspect ratio, L/d, for the (5, 5) and (7, 7) armchair carbon nanotubes in Fig. 5. It should be stated that the critical strains reported in Zhang et al.^[52] are now converted into the critical buckling parameter, , based on the following relations^[31]:

(92)

In view of Eq. (90) and the third equation of Eq. (41), we can obtain the following generalized equation:

(93)

It can be seen that Eq. (93) has an asymptotic value, which corresponds to the result of the Euler-Bernoulli beam model.

Numerical values for the fitting parameters A and B in Eq. (93) are obtained by fitting the molecular dynamics (MD) results given by Zhang et al.^[52]. Good agreement between the present fitting curves and the MD results is observed. In conclusion, the present expression of Eq. (93) can be used to model the size-dependent beam structures.

6 Concluding remarks

The nonlocal strain gradient theory, featuring a combined nonlocal theory and strain gradient theory, can capture both the softening and stiffening effects of materials/structures observed in the experiments or atomic simulations. However, when one uses this theory in analyzing the static and dynamic behaviors of materials/structures, the variationally consistent boundary conditions have not been correctly used. With a view toward exploring the variationally consistent boundary conditions, this paper reformulates the static boundary value problems of nonlocal strain gradient Timoshenko beams (see Eqs. (39) and (40)) with the help of the weighted residual method. Then, the analytical solutions for buckling loads of nonlocal strain gradient Timoshenko beams subject to three typical boundary conditions are obtained. The numerical results demonstrate clearly the effects of the two material length parameters and different selections of the higher-order boundary conditions on the buckling loads of Timoshenko beams. It is also found that the different selections of the higher-order boundary conditions of the beams with the CC boundary conditions play an important role in the buckling loads. In addition, the present analytical solutions may be useful for the choices of appropriate shape functions with the Galerkin method or the finite element method.