Nonlocal Strain Gradient Anisotropic Elastic Shell Model for Vibration Analysis of Single‒walled Carbon Nanotubes

1. Introduction

The use of classical continuum mechanics models in the study of carbon nanotubes (CNTs) dynamics can lead to inaccurate results. This is due to the actual discrete structure of CNTs and to their reduced dimensions. Therefore, in order to accurately investigate vibrations and stability of carbon nanotubes, they should be considered non‒classical continuum mechanics models based on anisotropic and size‒dependent theories.

Various anisotropic elastic shell theories were developed by the researchers. A very effective one was the theory first proposed by Chang [1,2], where the prediction of chirality and size‒dependent elastic properties of single‒walled carbon nanotubes was obtained via molecular mechanics model. The most important result obtained by Chang was that, for CNTs, the classical relationship from the isotropic elastic theory of continuum mechanics between Young’s modulus and shear modulus is not retained, and a more refined relationship taking into account the effect of tube diameter and chiral angle was proposed.

Starting from this theory, Ghavanloo and Fazelzadeh [3] proposed an anisotropic elastic shell model including chirality effect to investigate the vibration characteristics of SWCNTs. Considering Flügge shell theory and using complex variable method, they studied the effect of tube chirality and diameter on the natural frequencies of SWCNTs, together with the influence of external loads. On the basis of the anisotropic elastic shell model [3], the same Authors of the present paper carried out in Ref. [4] a comparison of shell theories for the vibration analysis of SWCNTs, specifically Donnell, Sanders and Flugge shell theories. Assuming as reference molecular dynamics results available in literature, they obtained that Flügge shell theory is the most accurate but also the most computationally expensive; on the other hand, they found that Donnell shell theory is not accurate while Sanders shell theory is very accurate in the vibration modelling of SWCNTs for all geometries and wavenumbers.

In addition to anisotropic models, several size‒dependent theories were introduced in literature, where the first one was the nonlocal elasticity theory developed by Eringen [5,6]. In the nonlocal differential constitutive relations of Eringen, the stress tensor at a reference point of a body is written as a function not only of the strain tensor at that point but also of the strain tensor at all other points of the body. To this aim, in the nonlocal elasticity equations of Eringen, it is inserted a nonlocal parameter, which is a small length scale constant appropriate to each material, whose value must be obtained by means of comparisons with the results of molecular dynamics simulations.

Starting from Eringen nonlocal theory, and taking into account the anisotropic model [3], Fazelzadeh and Ghavanloo [7] proposed a nonlocal anisotropic elastic shell model to study the linear vibrations of CNTs with arbitrary chirality. They investigated the effect of the nonlocal parameter on the natural frequencies of zigzag, armchair and chiral SWCNTs with different geometries and wavenumbers.

Another important size‒dependent theory was the strain gradient theory developed by Mindlin [8,9],

which represents an extension of the classical elasticity theory by considering additional higher‒order strain gradient terms with respect to stress tensor. Specifically, Mindlin strain gradient theory is based on the assumption that the materials cannot be simply modelled as collections of points, but they have to be considered as atoms with higher‒order deformation mechanisms at small (micro/nano) scale.

Combining nonlocal elasticity and strain gradient theories, Lim [10] proposed a new refined nonlocal strain gradient theory. Two different small length scale parameters, i.e., nonlocal and strain gradient (material) parameters, are adopted to account for the size‒dependent characteristics of nanomaterials. Dispersion relations based on the nonlocal strain gradient model with different values of nonlocal and material parameters with respect to wave propagation in the case of Euler–Bernoulli and Timoshenko nanobeams are analysed.

The nonlocal strain gradient theory developed by Lim was adopted to investigate the linear dynamics of beams in the framework of an isotropic elastic beam model in Refs. [11,12,13,14]. The effect of nonlocal and material parameters on vibrations and stability of beams was analysed, where nonlocal parameter was introduced to consider the influence of nonlocal elasticity, and material parameter was introduced to consider the significance of strain gradient. It was found that, when the material parameter is lower than the nonlocal parameter, the beam provides a stiffness‒softening effect on critical buckling force and natural frequencies, while, when the material parameter is higher than the nonlocal parameter, the beam exerts a stiffness‒hardening effect on critical buckling force and natural frequencies. Some

interesting reviews of size‒dependent continuum mechanics models for the linear vibration analysis of nanostructures can be found in Refs. [15,16,17].

Mehralian et al. [18] developed a nonlocal strain gradient isotropic elastic shell model to analyse the effect of nonlocal and material parameters on the linear vibrations of SWCNTs. The model reported in this paper is similar to the one reported in the present work. However, there is a relevant difference, that is the use of an isotropic (instead of an anisotropic) elastic shell model, where the adoption of an isotropic model for simulating SWCNT vibrations is proven to be non‒accurate, due to the anisotropic characteristics that are inherent of nanostructures [19,20].

Adopting classical continuum mechanics theories to model the dynamic behaviour of nanostructures needs a careful choice of equivalent parameters. Specifically, in the case of SWCNTs, their effective discrete structure can be modelled by means of a continuous cylindrical shell if equivalent parameters, i.e., Young’s modulus, Poisson’s ratio and thickness, are properly selected. To this aim, Yakobson et al. [21] obtained the values of the equivalent parameters comparing them with the values of the strain energy of discrete SWCNTs derived through molecular dynamics simulations. Interesting results on molecular dynamics simulations for SWCNT vibration analysis can be found in Refs. [22,23,24,25,26].

Readers interested in deepening shell theories are invited to refer to the fundamental books [27,28,29,30,31,32]. In particular, Leissa [27] studied the linear vibrations of cylindrical shells under different geometries, boundary conditions and wavenumbers. In addition, readers that are interested in nonlinear vibrations and energy exchanges in CNTs are invited to refer to the related papers [33,34,35,36,37,38,39,40,41], where also the effect of resonance interactions between different vibration modes, e.g., radial breathing and circumferential flexural modes, is evaluated, together with the influence of anisotropy and nonlocality. Finally, static models of CNTs considering the effect of nonlocal elasticity together with pull-in instability problems can be found in Refs. [42,43].

The aim of the present work is to develop an advanced elastic shell model for the vibration analysis of SWCNTs considering all three fundamental effects inherent to nanostructures previously reported, i.e., anisotropy, nonlocal elasticity and strain gradient. This is very important since a comprehensive hybrid anisotropic and size‒dependent model can lead to more realistic and therefore accurate results.

Sanders‒Koiter shell theory is used to obtain the strain‒displacement relationships. Eringen nonlocal elasticity theory and Mindlin strain gradient theory are used to get the constitutive equations, where the anisotropic elastic constants are expressed by means of Chang molecular mechanics model. The complex variable method is considered to analytically solve the dynamic equations of motion and to obtain the natural frequencies of SWCNTs with simply supported boundary conditions. The present model is first validated in an anisotropic form (i.e., without size‒dependent effects) via comparisons with the results of molecular dynamics simulations from the literature. Then, a parametric analysis is performed on the complete developed model to analyse the influence of size‒dependent effects (i.e., nonlocal and material parameters) on the natural frequencies of SWCNTs with different geometries and wavenumbers.

2. Sanders‒Koiter Shell Theory for SWCNTs

In the present paper, the actual discrete SWCNT is modelled by means of an equivalent continuous cylindrical shell, see Figure 1(a, b), with radius $R$, length $L$ and thickness $h$. A cylindrical coordinate system $(O,x,\theta ,z)$ is used, where the origin $O$ of the reference system is located at the centre of one end of the shell. Three displacements are measured: longitudinal $u(x,\theta ,t)$, circumferential $v(x,\theta ,t)$ and radial $w(x,\theta ,t)$, where the radial displacement $w$ is assumed as positive outward, $(x,\theta )$ are the longitudinal and circumferential coordinates of an arbitrary point on the middle surface of the shell, $z$ is the radial coordinate along the thickness $h$ of the shell and $t$ is the time.

In this paper, Sanders‒Koiter linear shell theory is considered to model SWCNT dynamics. The linear relationships between strains and displacements in Sanders‒Koiter theory are based on “Kirchhoff‒Love’s assumptions”, see Ref. [27] for the details.

The consequences of these geometric assumptions are that, in presence of thin cylindrical shell, the transverse shear deformations may be neglected $({\gamma }_{xz}={\gamma }_{\theta z}=0)$ in the expression of the constitutive equations, and the rotary inertia of the shell can be neglected in the expression of the kinetic energy. By considering the previous assumptions, in the Sanders‒Koiter linear shell theory the middle surface strains $({\epsilon }_{x,0},{\epsilon }_{\theta ,0},{\gamma }_{x\theta ,0})$ of the cylindrical shell are written as a function of the displacements $(u,v,w)$ in the following form [27]:

$${\epsilon }_{x,0}=\frac{\partial u}{\partial x}\hspace{1em}\hspace{1em}\hspace{1em}{\epsilon }_{\theta ,0}=\frac{1}{R}\frac{\partial v}{\partial \theta }+\frac{w}{R}\hspace{1em}\hspace{1em}\hspace{1em}{\gamma }_{x\theta ,0}=\frac{1}{R}\frac{\partial u}{\partial \theta }+\frac{\partial v}{\partial x}$$

(1)

Again, by considering the previous assumptions, in the Sanders‒Koiter linear shell theory the middle surface changes in curvature and torsion $(k_x,k_{\theta },k_{x\theta })$ of the cylindrical shell are expressed as [27]:

$$k_x=-\frac{{\partial }^{2}w}{\partial x^2}\hspace{1em}\hspace{1em}\hspace{1em}{k}_{\theta }=\frac{1}{R^2}\frac{\partial v}{\partial \theta }-\frac{1}{R^2}\frac{{\partial }^{2}w}{\partial {\theta }^2}\hspace{1em}\hspace{1em}\hspace{1em}{k}_{x\theta }=-\frac{2}{R}\frac{{\partial }^{2}w}{\partial x\partial \theta }+\frac{1}{2R}\left(3\frac{\partial v}{\partial x}-\frac{1}{R}\frac{\partial u}{\partial \theta }\right)$$

(2)

According to the Sanders‒Koiter shell theory, the strain components $({\epsilon }_x,{\epsilon }_{\theta },{\gamma }_{x\theta })$ at an arbitrary point of the surface of the cylindrical shell are related to the middle surface strains and changes in curvature and torsion via the radial coordinate $z$ by means of the following relationships [27]:

$${\epsilon }_x={\epsilon }_{x,0}+z{k}_x\hspace{1em}\hspace{1em}\hspace{1em}{\epsilon }_{\theta }={\epsilon }_{\theta ,0}+z{k}_{\theta }\hspace{1em}\hspace{1em}\hspace{1em}{\gamma }_{x\theta }={\gamma }_{x\theta ,0}+z{k}_{x\theta }$$

(3)

The adoption of Sanders‒Koiter shell theory to model the SWCNT dynamics in the present work is justified on the basis of the results obtained by the same Authors in a previous paper, see Ref. [4]. In that paper, a comparison of shell theories for vibration analysis of SWCNTs based on an anisotropic elastic shell model, specifically Donnell, Sanders and Flugge shell theories, was achieved, where the results of molecular dynamics simulations available in the literature were considered as references to check the accuracy of the three different shell theories. Flügge shell theory was found to be the most accurate and, actually, this theory was adopted in several papers investigating the linear vibrations of SWCNTs based on anisotropic elastic shell models, see Refs. [3,7]. On the other hand, it was proven that the additional terms present in the expressions of force and moment resultants, which give Flügge shell theory greater accuracy than Sanders and Donnell ones, lead to very high computational effort in the numerical simulations of SWCNT dynamic behaviour. Moreover, on the basis of the parametric analyses performed, it was found that Donnell shell theory is not accurate for several geometries and wavenumbers, while Sanders shell theory is very accurate for all geometries and wavenumbers. This is the reason why, in the present paper, Sanders shell theory is adopted instead of the more accurate but also more complex Flügge shell theory for the vibration modelling of SWCNTs.

In the framework of Sanders‒Koiter linear shell theory, in the following Section, a novel advanced anisotropic elastic shell model will be proposed, taking into account both nonlocal elasticity and strain gradient, i.e., two relevant small length scale effects characterising the SWCNT dynamic behaviour.

3. Nonlocal Strain Gradient Anisotropic Elastic Shell Model

According to the nonlocal strain gradient theory developed by Lim et al. [12], the general constitutive equation for size‒dependent structures is expressed as:

$$\left(1-{\mu }^2{\nabla }^2\right)\mathit{t}=\mathbf{C}:\mathit{\epsilon }-l^2\nabla \mathbf{C}:\nabla \mathit{\epsilon }$$

(4)

where $\mathit{t}$ is the stress tensor, $\mathbf{C}$ is the fourth‒order elasticity tensor, $\mathit{\epsilon }$ is the strain tensor, $\nabla \mathbf{C}$ is the elasticity gradient tensor, $\nabla \mathit{\epsilon }$ is the strain gradient tensor, ${\nabla }^2$ is the Laplace operator, $\mu$ is the nonlocal parameter, which is introduced to investigate the effect of the nonlocal elasticity, and $l$ is the material parameter, which is introduced to analyse the influence of the strain gradient.

For a shell‒type structure, the size‒dependent behaviour must be taken into consideration in the axial, circumferential and radial direction. Therefore, starting from the general equation (4), the constitutive equation of the nonlocal strain gradient anisotropic shell theory is given by:

$$\left(1-{\mu }^2{\nabla }^2\right)\mathit{t}=\frac{1}{h}\left(1-l^2{\nabla }^2\right)\mathbf{Y}\mathit{\epsilon }$$

(5)

where $\mathit{t}$ and $\mathit{\epsilon }$ are the stress and strain vectors, respectively, which for an elastic shell‒type structure under plane stress hypothesis are expressed as:

$$\mathit{t}={\left[{\sigma }_x,{\sigma }_{\theta },{\tau }_{x\theta }\right]}^{\mathrm{T}}\hspace{1em}\hspace{1em}\hspace{1em}\mathit{\epsilon }={\left[{\epsilon }_x,{\epsilon }_{\theta },{\gamma }_{x\theta }\right]}^{\mathrm{T}}$$

(6)

$h$ is the thickness of the shell, $\mathbf{Y}$ is the corresponding anisotropic elastic matrix, and:

$${\nabla }^2=\frac{{\partial }^2}{\partial x^2}+\frac{1}{R^2}\frac{{\partial }^2}{\partial {\theta }^2}$$

(7)

is the Laplace operator in the polar coordinate system.

The constitutive equation (5) can be projected in the $(x,\theta )$ plane in the form:

$$\left(1-{\mu }^2{\nabla }^2\right){\sigma }_x=\frac{1}{h}\left(1-l^2{\nabla }^2\right)\left(Y_{11}{\epsilon }_x+Y_{12}{\epsilon }_{\theta }+Y_{13}{\gamma }_{x\theta }\right)$$

$$\left(1-{\mu }^2{\nabla }^2\right){\sigma }_{\theta }=\frac{1}{h}\left(1-l^2{\nabla }^2\right)\left(Y_{21}{\epsilon }_x+Y_{22}{\epsilon }_{\theta }+Y_{23}{\gamma }_{x\theta }\right)$$

(8)

$$\left(1-{\mu }^2{\nabla }^2\right){\tau }_{x\theta }=\frac{1}{h}\left(1-l^2{\nabla }^2\right)\left(Y_{31}{\epsilon }_x+Y_{32}{\epsilon }_{\theta }+Y_{33}{\gamma }_{x\theta }\right)$$

The surface elastic constants $Y_{ij}$ (8), as elements of the anisotropic elastic matrix $\mathbf{Y}$, are given by [1]:

$$Y_{ij}=\frac{2}{3\sqrt{3}}\left(K_{\rho }{G}_{li}{G}_{lj}+\frac{{2K}_{\theta }}{a^2}{H}_{li}{H}_{lj}\right) i,j,l=\mathrm{1,2},3 \left(\mathrm{s}\mathrm{u}\mathrm{m} \mathrm{o}\mathrm{v}\mathrm{e}\mathrm{r} l\right)$$

(9)

where $a$ is the carbon‒carbon bond length, $(K_{\rho },K_{\theta })$ are force constants associated with the stretching and angular distortion of the carbon‒carbon bond, respectively, see Ref. [1] for the details.

The corresponding matrices $\mathbf{G}$ and $\mathbf{H}$ can be calculated as follows [3]:

$$\mathbf{G}={\mathbf{B}}^{-1}\left(\mathbf{I}-\mathbf{D}\mathbf{F}\right), \mathbf{H}=\mathbf{Q}\mathbf{F}$$

(10)

where $\mathbf{I}$ is the identity matrix, matrix $\mathbf{F}$ is given by [3]:

$$\mathbf{F}={\left[\mathbf{U}{\mathbf{B}}^{-1}\mathbf{D}-\left(\frac{2K_{\theta }}{K_{\rho }{a}^2}\mathbf{V}\mathbf{A}+\mathbf{W}\right)\right]}^{-1}\mathbf{U}{\mathbf{B}}^{-1}$$

(11)

and matrices $\left(\mathbf{A},\mathbf{B},\mathbf{D},\mathbf{U},\mathbf{V},\mathbf{W},\mathbf{Q}\right)$ are given by [2]:

$$\mathbf{A}=\left\{A_{ij}\right\}=\left\{-\cos {\omega }_{ik}\cos {\omega }_{jk}\right\} i,j,k=\mathrm{1,2},3 \left(\mathrm{s}\mathrm{u}\mathrm{m} \mathrm{o}\mathrm{v}\mathrm{e}\mathrm{r} k\right)$$

(12)

$$\mathbf{B}=\frac{1}{3\sqrt{n^2+nm+m^2}}\left(\begin{array}{ccc}(2n+m)\cos {\varphi }_1& -(n-m)\cos {\varphi }_2& -(n+2m)\cos {\varphi }_3\\ \sqrt{3} m\sin {\varphi }_1& -\sqrt{3}(n+m)\sin {\varphi }_2& \sqrt{3} n\sin {\varphi }_3\\ (2n+m)\sin {\varphi }_1& -(n-m)\sin {\varphi }_2& -(n+2m)\sin {\varphi }_3\end{array}\right)$$

(13)

$$\mathbf{D}=\frac{1}{3\sqrt{n^2+nm+m^2}}\left(\begin{array}{ccc}-(2n+m)\sin {\varphi }_1& (n-m)\sin {\varphi }_2& (n+2m)\sin {\varphi }_3\\ \sqrt{3} m\cos {\varphi }_1& -\sqrt{3}(n+m)\cos {\varphi }_2& \sqrt{3} n\cos {\varphi }_3\\ (2n+m)\cos {\varphi }_1& -(n-m)\cos {\varphi }_2& -(n+2m)\cos {\varphi }_3\end{array}\right)$$

(14)

$$\mathbf{U}=\left(\begin{array}{ccc}\sin {\varphi }_1& \sin {\varphi }_2& \sin {\varphi }_3\\ \cos {\varphi }_1& \cos {\varphi }_2& \cos {\varphi }_3\\ m\cos {\varphi }_1& -(n+m)\cos {\varphi }_2& n\cos {\varphi }_3\end{array}\right)$$

(15)

$$\mathbf{V}=\left(\begin{array}{ccc}-\cos {\varphi }_1& -\cos {\varphi }_2& -\cos {\varphi }_3\\ \sin {\varphi }_1& \sin {\varphi }_2& \sin {\varphi }_3\\ 0& 0& 0\end{array}\right)$$

(16)

$$\mathbf{W}=\left(\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ -m\sin {\varphi }_1& (n+m)\sin {\varphi }_2& -n\sin {\varphi }_3\end{array}\right)$$

(17)

$$\mathbf{Q}=\left\{Q_{ij}\right\}=\left\{-\cos {\omega }_{ji}\right\} i,j=\mathrm{1,2},3$$

(18)

where [2]:

$$\cos {\omega }_{ij}=\left\{\begin{array}{c}(\cos {\varphi }_i\sin {\varphi }_k\cos {\phi }_j-\sin {\varphi }_i\cos {\varphi }_k)/\sin {\theta }_j\\ 0\end{array}\right. \begin{array}{c}i\ne j\ne k\\ i=j\end{array}$$

(19)

and $(n,m)$ are the chirality indices of the SWCNT, which define its radius via the relation [1]:

$$R=\frac{\sqrt{3} a}{2\pi } \sqrt{n^2+nm+m^2}$$

(20)

The structural parameters of the SWCNT, i.e., chiral angles $({\varphi }_1,{\varphi }_2,{\varphi }_3)$, torsion angles $({\phi }_1,{\phi }_2,{\phi }_3)$ and bond angles $({\theta }_{1,i},{\theta }_{2,i},{\theta }_{3,i})$, can be calculated by means of the equations [1]:

$${\varphi }_1=\arccos \frac{2n+m}{2\sqrt{n^2+nm+m^2}}\hspace{1em}\hspace{1em}\hspace{1em}{\varphi }_2=\frac{4\pi }{3}+{\varphi }_1\hspace{1em}\hspace{1em}\hspace{1em}{\varphi }_3=\frac{2\pi }{3}+{\varphi }_1$$

(21)

$$\begin{gathered} {\phi }_1=\frac{\pi }{\sqrt{n^2+nm+m^2}}\cos {\varphi }_1\hspace{1em}\hspace{1em}\hspace{1em}{\phi }_2=\frac{\pi }{\sqrt{n^2+nm+m^2}}\cos \left(\frac{\pi }{3}+{\varphi }_1\right) \\ {\phi }_3=\frac{\pi }{\sqrt{n^2+nm+m^2}}\cos \left(\frac{\pi }{3}-{\varphi }_1\right) \end{gathered}$$

(22)

$$\cos {\theta }_i=\sin {\varphi }_j\sin {\varphi }_k\cos {\phi }_i+\cos {\varphi }_j\cos {\varphi }_k i,j,k=\mathrm{1,2},3 i\ne j\ne k$$

(23)

4. Force and Moment Resultants

In the present nonlocal strain gradient anisotropic elastic shell model, the force $(N_x,N_{\theta },N_{x\theta })$ and moment $(M_x,M_{\theta },M_{x\theta })$ resultants per unit length are derived by integrating the stress components of the constitutive equations (8) and considering the thin shell assumption $\left(z/R\ll 1\right)$ of Sanders‒Koiter shell theory as follows:

$$\begin{gathered} \left(1-{\mu }^2{\nabla }^2\right)N_x=\left(1-l^2{\nabla }^2\right)·\left[Y_{11}\frac{\partial u}{\partial x}+\right. \\ \left.+\frac{Y_{12}}{R}\left(\frac{\partial v}{\partial \theta }+w\right)+Y_{13}\left(\frac{\partial v}{\partial x}+\frac{1}{R}\frac{\partial u}{\partial \theta }\right)\right] \end{gathered}$$

(24)

$$\begin{gathered} \left(1-{\mu }^2{\nabla }^2\right)N_{\theta }=\left(1-l^2{\nabla }^2\right)·\left[Y_{21}\frac{\partial u}{\partial x}+\right. \\ \left.+\frac{Y_{22}}{R}\left(\frac{\partial v}{\partial \theta }+w\right)+Y_{23}\left(\frac{\partial v}{\partial x}+\frac{1}{R}\frac{\partial u}{\partial \theta }\right)\right] \end{gathered}$$

(25)

$$\begin{gathered} \left(1-{\mu }^2{\nabla }^2\right)N_{x\theta }=\left(1-l^2{\nabla }^2\right)·\left[Y_{31}\frac{\partial u}{\partial x}+\right. \\ \left.+\frac{Y_{32}}{R}\left(\frac{\partial v}{\partial \theta }+w\right)+Y_{33}\left(\frac{\partial v}{\partial x}+\frac{1}{R}\frac{\partial u}{\partial \theta }\right)\right] \end{gathered}$$

(26)

$$\begin{gathered} \left(1-{\mu }^2{\nabla }^2\right)M_x=\left(1-l^2{\nabla }^2\right)·\left[-X_{11}\frac{{\partial }^{2}w}{{\partial x}^2}+\frac{X_{12}}{R^2}\right.· \\ \left.·\left(\frac{\partial v}{\partial \theta }-\frac{{\partial }^{2}w}{{\partial \theta }^2}\right)+X_{13}\left(-\frac{2}{R}\frac{{\partial }^{2}w}{\partial x\partial \theta }-\frac{1}{2R^2}\frac{\partial u}{\partial \theta }+\frac{3}{2R}\frac{\partial v}{\partial x}\right)\right] \end{gathered}$$

(27)

$$\begin{gathered} \left(1-{\mu }^2{\nabla }^2\right)M_{\theta }=\left(1-l^2{\nabla }^2\right)·\left[-X_{21}\frac{{\partial }^{2}w}{{\partial x}^2}+\frac{X_{22}}{R^2}·\right. \\ ·\left.\left(\frac{\partial v}{\partial \theta }-\frac{{\partial }^{2}w}{{\partial \theta }^2}\right)+X_{23}\left(-\frac{2}{R}\frac{{\partial }^{2}w}{\partial x\partial \theta }-\frac{1}{2R^2}\frac{\partial u}{\partial \theta }+\frac{3}{2R}\frac{\partial v}{\partial x}\right)\right] \end{gathered}$$

(28)

$$\begin{gathered} \left(1-{\mu }^2{\nabla }^2\right)M_{x\theta }=\left(1-l^2{\nabla }^2\right)·\left[-X_{31}\frac{{\partial }^{2}w}{{\partial x}^2}+\frac{X_{32}}{R^2}·\right. \\ ·\left(\frac{\partial v}{\partial \theta }-\frac{{\partial }^{2}w}{{\partial \theta }^2}\right)\left.+X_{33}\left(-\frac{2}{R}\frac{{\partial }^{2}w}{\partial x\partial \theta }-\frac{1}{2R^2}\frac{\partial u}{\partial \theta }+\frac{3}{2R}\frac{\partial v}{\partial x}\right)\right] \end{gathered}$$

(29)

where the elements of the bending stiffness matrix $\mathbf{X}$ can be defined as [3]:

$$X_{ij}=\frac{Y_{ij}{h}^2}{12} i,j=\mathrm{1,2},3$$

(30)

It must be stressed that the force $(N_x,N_{\theta },N_{x\theta })$ and moment $(M_x,M_{\theta },M_{x\theta })$ resultants per unit length written by considering Sanders‒Koiter shell theory are different from the ones of Flügge shell theory since, in the last theory, the ratio $z/R$ is not neglected, i.e., the thin shell assumption is not taken into consideration, see Ref. [27] for the details.

5. Equations of Motion

The classical dynamic equilibrium equations in terms of force and moment resultants are written as (external forces and moments are neglected) [3]:

$$\frac{\partial N_x}{\partial x}+\frac{1}{R}\frac{\partial N_{x\theta }}{\partial \theta }-\frac{1}{2R^2}\frac{\partial M_{x\theta }}{\partial \theta }-\rho h\frac{{\partial }^{2}u}{{\partial t}^2}=0$$

(31)

$$\frac{1}{R}\frac{\partial N_{\theta }}{\partial \theta }+\frac{\partial N_{x\theta }}{\partial x}+\frac{3}{2R}\frac{\partial M_{x\theta }}{\partial x}+\frac{1}{R^2}\frac{\partial M_{\theta }}{\partial \theta }-\rho h\frac{{\partial }^{2}v}{{\partial t}^2}=0$$

(32)

$$\frac{{\partial }^2{M}_x}{{\partial x}^2}+\frac{2}{R}\frac{{\partial }^2{M}_{x\theta }}{\partial x\partial \theta }+\frac{1}{R^2}\frac{{\partial }^2{M}_{\theta }}{{\partial \theta }^2}-\frac{N_{\theta }}{R}-\rho h\frac{{\partial }^{2}w}{{\partial t}^2}=0$$

(33)

where $\rho h$ is the mass density per unit area (i.e., surface density) of the SWCNT.

Applying the nonlocal elasticity operator $(1-{\mu }^2{\nabla }^2)$ on the dynamic equilibrium equations (31‒33), and then substituting equations (24‒29) into the modified equations (31‒33), we obtain:

$$\begin{gathered} \left(1-l^2{\nabla }^2\right)\left\{\left[Y_{11}\frac{{\partial }^2}{{\partial x}^2}+\frac{2Y_{13}}{R}\frac{{\partial }^2}{\partial x\partial \theta }+\left(\frac{Y_{33}}{R^2}+\frac{X_{33}}{{4R}^4}\right)\frac{{\partial }^2}{{\partial \theta }^2}\right]u\right. \\ +\left[Y_{13}\frac{{\partial }^2}{{\partial x}^2}+\left(\frac{Y_{12}+Y_{33}}{R}-\frac{{3X}_{33}}{4R^3}\right)\frac{{\partial }^2}{\partial x\partial \theta }+\left(\frac{Y_{23}}{R^2}-\frac{X_{23}}{{2R}^4}\right)\frac{{\partial }^2}{{\partial \theta }^2}\right]v \\ \left.+\left[\frac{Y_{12}}{R}\frac{\partial }{\partial x}+\frac{Y_{23}}{R^2}\frac{\partial }{\partial \theta }+\frac{X_{13}}{{2R}^2}\frac{{\partial }^3}{{\partial x}^2\partial \theta }+\frac{X_{33}}{R^3}\frac{{\partial }^3}{\partial x{\partial \theta }^2}+\frac{X_{23}}{{2R}^4}\frac{{\partial }^3}{{\partial \theta }^3}\right]w\right\} \\ =\rho h\left(\frac{{\partial }^{2}u}{{\partial t}^2}-{\mu }^2\frac{{\partial }^{4}u}{{\partial x}^2{\partial t}^2}-{\mu }^2\frac{1}{R^2}\frac{{\partial }^{4}u}{{\partial \theta }^2{\partial t}^2}\right) \end{gathered}$$

(34)

$$\left(1-l^2{\nabla }^2\right)\left\{\left[Y_{13}\frac{{\partial }^2}{{\partial x}^2}+\left(\frac{Y_{12}+Y_{33}}{R}-\frac{{3X}_{33}}{{4R}^3}\right)\frac{{\partial }^2}{\partial x\partial \theta }+\left(\frac{Y_{23}}{R^2}-\frac{X_{23}}{{2R}^4}\right)\frac{{\partial }^2}{{\partial \theta }^2}\right]u\right.$$

$$\begin{gathered} +\left[\left(Y_{33}+\frac{{9X}_{33}}{{4R}^2}\right)\frac{{\partial }^2}{{\partial x}^2}+\left(\frac{{2Y}_{23}}{R}+\frac{{3X}_{23}}{R^3}\right)\frac{{\partial }^2}{\partial x\partial \theta }+\left(\frac{Y_{22}}{R^2}+\frac{X_{22}}{R^4}\right)\frac{{\partial }^2}{{\partial \theta }^2}\right]v \\ +\left[\frac{Y_{23}}{R}\frac{\partial }{\partial x}+\frac{Y_{22}}{R^2}\frac{\partial }{\partial \theta }-\frac{{3X}_{13}}{2R}\frac{{\partial }^3}{{\partial x}^3}-\left(\frac{X_{12}+{3X}_{33}}{R^2}\right)\frac{{\partial }^3}{{\partial x}^2\partial \theta }-\frac{{7X}_{23}}{2R^3}\frac{{\partial }^3}{\partial x{\partial \theta }^2}\right. \end{gathered}$$

(35)

$$\left.\left.-\frac{X_{22}}{R^4}\frac{{\partial }^3}{{\partial \theta }^3}\right]w\right\}=\rho h\left(\frac{{\partial }^{2}v}{{\partial t}^2}-{\mu }^2\frac{{\partial }^{4}v}{{\partial x}^2{\partial t}^2}-{\mu }^2\frac{1}{R^2}\frac{{\partial }^{4}v}{{\partial \theta }^2{\partial t}^2}\right)$$

$$\begin{gathered} \left(1-l^2{\nabla }^2\right)\left\{\left[-\frac{Y_{12}}{R}\frac{\partial }{\partial x}-\frac{Y_{23}}{R^2}\frac{\partial }{\partial \theta }-\frac{X_{13}}{{2R}^2}\frac{{\partial }^3}{{\partial x}^2\partial \theta }-\frac{X_{33}}{R^3}\frac{{\partial }^3}{\partial x\partial {\theta }^2}-\frac{X_{23}}{2R^4}\frac{{\partial }^3}{{\partial \theta }^3}\right]u\right. \\ +\left[-\frac{Y_{23}}{R}\frac{\partial }{\partial x}-\frac{Y_{22}}{R^2}\frac{\partial }{\partial \theta }+\frac{{3X}_{13}}{2R}\frac{{\partial }^3}{{\partial x}^3}+\left(\frac{X_{12}+{3X}_{33}}{R^2}\right)\frac{{\partial }^3}{{\partial x}^2\partial \theta }+\frac{7X_{23}}{2R^3}\frac{{\partial }^3}{\partial x{\partial \theta }^2}\right. \\ \left.+\frac{X_{22}}{R^4}\frac{{\partial }^3}{{\partial \theta }^3}\right]v+\left[-\frac{Y_{22}}{R^2}-X_{11}\frac{{\partial }^4}{{\partial x}^4}-\frac{{4X}_{13}}{R}\frac{{\partial }^4}{{\partial x}^3\partial \theta }-\left(\frac{{{2X}_{12}+4X}_{33}}{R^2}\right)\frac{{\partial }^4}{\partial x^2{\partial \theta }^2}\right. \\ \left.\left.-\frac{{4X}_{23}}{R^3}\frac{{\partial }^4}{\partial x{\partial \theta }^3}-\frac{X_{22}}{R^4}\frac{{\partial }^4}{{\partial \theta }^4}\right]w\right\}=\rho h\left(\frac{{\partial }^{2}w}{{\partial t}^2}-{\mu }^2\frac{{\partial }^{4}w}{{\partial x}^2{\partial t}^2}-{\mu }^2\frac{1}{R^2}\frac{{\partial }^{4}w}{{\partial \theta }^2{\partial t}^2}\right) \end{gathered}$$

(36)

which represent the equations of motion for an arbitrary chiral SWCNT in terms of longitudinal $u$, circumferential $v$ and radial $w$ displacements of the middle surface of the SWCNT.

From the equations of motion (34‒36), it can be observed that the strain gradient operator $(1-l^2{\nabla }^2)$ is applied on the strains and changes in curvature and torsion, while the nonlocal elasticity operator $(1-{\mu }^2{\nabla }^2)$ is applied on the accelerations.

6. Solution Method

In this paper, the complex variable method is considered to analytically solve the partial differential equations of motion (34‒36) by setting the real and imaginary parts equal to zero and therefore to get the natural frequencies of the SWCNT.

In the present work, simply supported boundary conditions are adopted. These boundary conditions, for the complex variable method, impose the relations $\mathrm{R}\mathrm{e}\left(v\right)=\mathrm{R}\mathrm{e}\left(w\right)=0$ at both ends $x=(0, L)$ of the SWCNT. The displacement field that satisfies these boundary conditions can be written as [7]:

$$\begin{gathered} u(x,\theta ,t)=\overline{U}\exp \left(i{\lambda }_{q}x\right)\cos \left(s\theta \right)\cos \left(\omega t\right) \\ v(x,\theta ,t)=-i\overline{V}\exp \left(i{\lambda }_{q}x\right)\sin \left(s\theta \right)\cos \left(\omega t\right) \\ w(x,\theta ,t)=-i\overline{W}\exp \left(i{\lambda }_{q}x\right)\cos \left(s\theta \right)\cos \left(\omega t\right) \end{gathered}$$

(37)

where $\left(\overline{U},\overline{V},\overline{W}\right)$ denote the displacement amplitudes along the longitudinal $u$, circumferential $v$ and radial $w$ directions, respectively, $i$ is the imaginary unit, ${\lambda }_q$ is the wavenumber along the longitudinal direction, with ${\lambda }_q=q\pi /L$, where $q$ is the number of longitudinal half‒waves and $L$ is the length of the SWCNT, $s$ is the number of circumferential waves and $\omega$ is the circular frequency.

Substituting equation (37) into equations (34‒36), a set of algebraic equations for the displacement amplitudes $\left(\overline{U},\overline{V},\overline{W}\right)$ is obtained, which can be written in the matrix form [7]:

$$\mathbf{E}{({\lambda }_q,s,\omega )}_{3\times 3}\left[\begin{array}{c}\overline{U}\\ \overline{V}\\ \overline{W}\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]$$

(38)

where $\mathbf{E}$ is a non‒symmetric matrix, whose elements are reported in Appendix.

For a non‒trivial solution $\left(\overline{U}=\overline{V}=\overline{W}=0\right)$, the determinant of matrix $\mathbf{E}$ (38) must be equal to zero:

$$\det {\mathbf{E}({\lambda }_q,s,\omega )}_{3\times 3}=0$$

(39)

Solving equation (39) we get a third‒degree algebraic equation in ${\omega }^2$; this last equation provides three different eigenfrequencies for each number of waves $(q,s)$ that give three different vibration modes (i.e., longitudinal, torsional and radial modes). Since the highest natural frequency corresponds to the radial vibration mode, then in the numeric results only the radial natural frequencies will be computed.

7. Numeric Results

In this paper, the effect of nonlocal and material parameters on the natural frequencies of SWCNTs is considered. Sanders-Koiter shell theory is used to obtain the strain‒displacement relationships. An anisotropic elastic shell model is adopted to take into account the intrinsic chirality effects of CNTs. Simply supported boundary conditions are imposed. Vibration modes with different number of waves along the longitudinal and circumferential directions are analysed. SWCNTs with different chiralities and geometries are investigated.

In Table 1, the values of the carbon‒carbon bond parameters $(a,k_{\rho },k_{\theta })$ and the equivalent parameters $(h,\rho )$ retrieved from the pertinent literature are reported. In particular, parameters $k_{\rho }$ and $k_{\theta }$, which are force constants related to the variance of carbon‒carbon bond length $a$ and angle $\theta$, respectively, are considered to obtain the surface elastic constants $Y_{ij}$ (9) of the SWCNT by means of the molecular mechanics model developed by Chang [1,2].

Moreover, in order to study the dynamics of the actual discrete SWCNT via an equivalent continuous cylindrical shell, an equivalent thickness $h$, which is obtained from MD simulations of CNT energy, and an equivalent mass density $\rho$, resulting from the surface density of graphite, are considered, see Ref. [21] for the details.

Finally, in the parametric analyses, a SWCNT with radius R = 1.34 nm (i.e., thickness ratio R/h = 20) will be considered. It should be stressed that this value of thickness ratio respects the hypothesis that allows a thin shell theory (in this paper Sanders‒Koiter shell theory) to be applied, which is R/h > 10, see Ref. [27] for the details.

7.1. Validation of the Anisotropic Elastic Shell Model

The first step of the present work is the validation of the anisotropic elastic model based on Sanders–Koiter shell theory that will be adopted in the following to analyse the effect of nonlocal and material parameters on the natural frequencies of SWCNTs.

This validation is carried out by comparing the results of the present anisotropic elastic shell model with the ones of molecular dynamics simulations available in literature [22]. The natural frequencies of the radial breathing mode (i.e., the undeformed vibration mode characteristic of CNTs presenting no longitudinal and circumferential waves) under different chirality indices are considered.

From the comparisons it can be observed that the percentage difference is relatively low (maximum value ≈ 2.7%, medium value ≈ 1.6%) for all the considered chirality indices, see Table 2, and therefore the present anisotropic elastic shell model can be considered as accurate.

On the other hand, it must be underlined that the development and implementation of the anisotropic elastic shell model present high analytical complexity and computational effort. Therefore, it is useful to check if the corresponding isotropic elastic shell model, which presents low analytical complexity and computational effort, can provide similar results.

From the comparisons between the results of the isotropic elastic shell model and molecular dynamics simulations [22] it can be observed that the percentage difference is higher than the anisotropic elastic shell model (maximum value ≈ 3.6%, medium value ≈ 2.9%) for all the considered chirality indices, see Table 3.

Since the anisotropic elastic shell model was demonstrated to be significantly more accurate than the corresponding isotropic one, then it will be adopted in the following parametric analyses.

7.2. Effect of Nonlocal and Material Parameters

In this Section, the effect of nonlocal and material parameters on the natural frequencies of the simply supported SWCNT of Table 1 is studied. Different chiralities and geometries are analysed. Vibration modes with different wavenumber along the longitudinal and circumferential directions are evaluated.

The first goal is to analyse the effect of nonlocal and material parameters on the natural frequencies of the simply supported SWCNT of Table 1 for a generic vibration mode.

In Figure 2, the natural frequencies of the vibration mode (q = 1, s = 2) with one longitudinal half‒wave and two circumferential waves are presented. The chirality indices (n = 34, m = 0) (i.e., zigzag SWCNT) are considered. Thickness ratio R/h = 20 and aspect ratio L/R = 10 are adopted. From Figure 2 it is observed that, for a fixed value of material parameter $l$, the natural frequencies decrease as the nonlocal parameter $\mu$ increases. On the contrary, for a fixed value of nonlocal parameter $\mu$, the natural frequencies increase as the material parameter $l$ increases. Therefore, an opposite effect between the two small length scale parameters on the natural frequencies is found.

The second goal is to study the effect of nonlocal and material parameters on the natural frequencies of the simply supported SWCNT of Table 1with different chiralities.

In Figure 3, the natural frequencies of the same vibration mode (q = 1, s = 2) of Figure 2 are shown. A SWCNT with the same geometry but with different chirality indices (n = 20, m = 20) (i.e., armchair SWCNT) is considered. From Figure 3 it can be observed that the effect of nonlocal $\mu$ and material $l$ parameters on the natural frequencies of the vibration mode (q = 1, s = 2) is the same of that of Figure 2. Therefore, it can be deduced that the effect of the two small length scale parameters on the natural frequencies is independent of SWCNT chirality. Starting from this result, in the following simulations the chirality indices (n = 34, m = 0) (i.e., zigzag SWCNT) will be considered.

The third goal is to analyse the effect of nonlocal and material parameters on the natural frequencies of the simply supported SWCNT of Table 1 with different geometries.

In Figure 4, the natural frequencies of the same vibration mode (q = 1, s = 2) of Figure 2 are shown. A SWCNT with the same aspect ratio L/R = 10 but different thickness ratios R/h is considered. From Figure 4, it can be noted that the natural frequencies decrease as the value of the thickness ratio R/h increases. For relatively low thickness ratios 20 < R/h < 50, the decrease is exponential as the value of material parameter l increases, while it is linear as the value of nonlocal parameter $\mu$ increases. For relatively high thickness ratios R/h > 80, the decrease is linear and the natural frequencies are similar for all values of nonlocal and material parameters. Therefore, the effect of the two small length scale parameters on the natural frequencies is strongly dependent of SWCNT radius.

In Figure 5, the natural frequencies of the same vibration mode (q = 1, s = 2) of Figure 2 are shown. A SWCNT with the same thickness ratio R/h = 20 but different aspect ratios L/R is considered. From Figure 5 it can be observed that, for every value of nonlocal $\mu$ and material $l$ parameters, the natural frequencies are constant as the value of aspect ratio L/R increases, and therefore the effect of the two small length scale parameters on the natural frequencies is independent of SWCNT length.

The last goal is to analyse the effect of nonlocal and material parameters on the natural frequencies of the simply supported SWCNT of Table 1 for vibration modes with different wavenumbers. First it is evaluated the effect of the number of longitudinal half‒waves q.

In Figure 6, the natural frequencies of the axisymmetric modes (s = 0) of the SWCNT of Table 1 with thickness ratio R/h = 20 and aspect ratio L/R = 10 are considered. The effect of nonlocal μ and material l parameters for a different number of longitudinal half‒waves q is investigated. From Figure 6, first of all, it is derived that the natural frequency of the undeformed mode (q = 0, s = 0) is independent of the two small length scale parameters (i.e., it is constant). Moreover, increasing the value of nonlocal parameter μ, it is found a linear decrease of the natural frequencies within the range q = (0 ‒ 3) and a subsequent linear increase by further increasing the number of longitudinal half‒waves. On the other hand, increasing the value of material parameter $l$, it is observed an exponential increase of the natural frequencies as the number of the longitudinal half‒waves increases within the range q = (0 ‒ 5).

In Figure 7, the natural frequencies of the beam‒like modes (s = 1) of the same SWCNT of Figure 6 are analysed. The effect of nonlocal μ and material l parameters for a different number of longitudinal half‒waves q is investigated. From Figure 7, as the value of nonlocal parameter μ increases, it is found a little decrease of the natural frequencies within the range q = (0 ‒ 3) and a subsequent little increase by further increasing the number of longitudinal half‒waves. On the other hand, increasing the value of material parameter $l$, it is obtained an exponential increase of the natural frequencies as the number of the longitudinal half‒waves increases within the range q = (0 ‒ 5), where this exponential increase is lower than the one of Figure 6 for the axisymmetric modes (s = 0).

In Figure 8, the natural frequencies of the shell‒like modes (s = 2) of the same SWCNT of Figure 6 are analysed. The effect of nonlocal μ and material l parameters for a different number of longitudinal half‒waves q is studied. From Figure 8 it is obtained that, increasing the number of longitudinal half‒waves within the range q = (0 ‒ 5), the natural frequencies remain quasi‒constant increasing the value of nonlocal parameter, while they increase exponentially increasing the value of material parameter, where this exponential increase is lower than the one of Figure 7 for the beam‒like modes (s = 1).

Therefore, by comparing Figure 6, Figure 7 and Figure 8, it is obtained that, as the number of longitudinal half‒waves q increases, the natural frequencies first linearly decrease and then linearly increase with increasing the nonlocal parameter μ, while they exponentially increase with increasing the material parameter l. The magnitude of these opposite behaviours reduces with increasing the number of circumferential waves s (i.e., the linear first decrease and then increase of natural frequencies with increasing the nonlocal parameter μ becomes quasi‒constant, the exponential increase of natural frequencies with increasing the material parameter l becomes quasi‒linear).

Then it is investigated the effect of the number of circumferential waves s.

In Figure 9, the natural frequencies of the modes with zero longitudinal half‒waves (q = 0) of the SWCNT of Table 1 with thickness ratio R/h = 20 and aspect ratio L/R = 10 are considered. The effect of nonlocal μ and material l parameters for a different number of circumferential waves s is evaluated. From Figure 9, first of all, it is confirmed that the natural frequency of the undeformed mode (q = 0, s = 0) is independent of the two small length scale parameters. Additionally, with increasing the value of nonlocal parameter μ, it is observed a little increase of the natural frequencies as the number of the circumferential waves increases within the range s = (0 ‒ 5). On the other hand, increasing the value of material parameter $l$, it can be observed a strongly exponential increase of the natural frequencies as the number of the circumferential waves increases within the range s = (0 ‒ 5).

In Figure 10 and Figure 11, the natural frequencies of the modes with respectively one (q = 1) and two (q = 2) longitudinal half‒waves of the same SWCNT of Figure 9 are analysed. The effect of nonlocal μ and material l parameters for a different number of circumferential waves s is investigated. From Figure 10 and Figure 11, it is noted exactly the same behaviour of the natural frequencies of that in Figure 9.

Therefore, by comparing Figure 9, Figure 10 and Figure 11, it is observed that, as the number of circumferential waves s increases, the natural frequencies have a little increase with increasing the nonlocal parameter μ, while they have a strongly exponential increase with increasing the material parameter l, and this behaviour is independent of the number of longitudinal half‒waves q.

8. Conclusions

The main contribution of this paper is the development of a novel nonlocal strain gradient anisotropic elastic shell model to analyse the linear vibrations of SWCNTs. Based on Eringen nonlocal elasticity and Mindlin strain gradient theories, and taking into account the inherent anisotropic elastic behaviour of nanostructures by adopting Chang molecular mechanics model, for the first time, to the Authors’ knowledge, the combined effect of nonlocal and material (strain gradient) parameters on the natural frequencies of SWCNTs with different chirality, geometry and wavenumber is investigated.

The most relevant findings of the present work are summarised as follows.

• The natural frequency of the undeformed vibration mode is independent of both nonlocal and material parameters.
• For a generic linear vibration mode, the natural frequencies decrease as the nonlocal parameter increases, while they increase as the material parameter increases.
• The decrease of the natural frequencies with increasing SWCNT radius is exponential as the material parameter increases, while it is linear as the nonlocal parameter increases.
• The effect of nonlocal and material parameters on the natural frequencies is independent of SWCNT chirality and length.
• As the number of longitudinal half‒waves increases, the natural frequencies linearly decrease with increasing the nonlocal parameter, while they exponentially increase with increasing the material parameter.
• As the number of circumferential waves increases, the natural frequencies little increase with increasing the nonlocal parameter, while they strongly exponentially increase with increasing the material parameter.