Energy exchange and localization of low-frequency oscillations in single-walled carbon nanotubes
Valeri V. Smirnov, Leonid I. Manevitch, Matteo Strozzi, Francesco Pelicano
EEnergy exchange and localization of low-frequency oscillations in single-walled carbonnanotubes
V.V. Smirnov , ∗ L.I. Manevitch , M.Strozzi , and F.Pelicano Institute of Chemical Physics, RAS, 4 Kosygin str.,119991 Moscow, Russia Department of Engineering ”Enzo Ferrari”,University of Modena and Reggio Emilia, Strada Vignolese 905, 41125 Modena, Italy (Dated: October 10, 2018)We present the results of analytical study and Molecular Dynamics simulation of low energynonlinear non-stationary dynamics of single-walled carbon nanotubes (CNTs). New phenomena ofintensive energy exchange between different parts of CNT and weak energy localization in a partCNT are analytically predicted in the framework of the continuum shell theory. These phenomenatake place for CNTs of finite length with medium aspect ratio under different boundary conditions.Their origin is clarified by means of the concept of Limiting Phase Trajectory, and the analyticalresults are confirmed by the MD simulation of simply supported CNTs.
PACS numbers: 61.48.De, 63.22.Gh, 63.20.D-,05.45.-a
I. INTRODUCTION
From a modern point of view, carbon nanotubes areattractive subjects for two reasons. In the first place,they are associated with great hopes for the creation ofsuper-small and ultra-fast electronic and electromechani-cal devices with unique physical properties [1–4]. On theother hand, they are quasi-one-dimensional objects, thatallows to check out some of the fundamentals of the mod-ern solid-state physics. In particular, variuos computa-tional and in-situ measurements of thermoconductivityof CNT [5–9] are directly related with the problem offiniteness of thermoconductivity of one-dimensional an-harmonic lattices. This problem has been formulatedmore than fifty years ago by Fermi, Pasta and Ulam [10]and it is not completely resolved till now. The stationarydynamics of CNTs or nonstationary, but non-resonancedynamics of CNT, can be treated in terms of linear ornonlinear normal modes. Using their combinations, onecan describe the CNT oscillations under arbitrary initialconditions. However, the situation drastically changesif we deal with non-stationary resonance processes suchas energy transfer. In the framework of the linear the-ory, the energy transfer requires the formation of a wavepacket, the time evolution of which depends strongly onthe dispersion properties of the system. The dispersionleads to the wave packet spreading that strongly affectsthe energy transfer efficiency. In the nonlinear systems,the dispersive spreading, can be compensated by non-linearity. As a result, a soliton (breather) mechanismof energy transfer in the infinite quasi-one-dimensionalnonlinear lattices arises. However, it was recently shown[11–13] that the resonant interaction of nonlinear normalmodes in the finite lattices leads to the existence of signif-icant non-stationary phenomena, which disappear in theinfinite case: i) the intensive energy exchange between ∗ [email protected] different parts of the system, which can be observed at asmall enough excitation level, ii) the existence of an in-stability threshold for the zone-boundering mode, iii) thetransition to weak energy localization in some part of thesystem. All these phenomena can be understood and effi-ciently described by a unified viewpoint in the frameworkof Limiting Phase Trajectory (LPT) concept. The LPTscorrespond to strongly non-stationary processes, whichare characterized by the maximum possible (under givenconditions) energy exchange between different parts ofthe system. In the present paper we show that the pro-cesses of the intensive energy exchange and the transitionto the energy capture in the some part of the CNT canbe explained by the transformation of the LPT with thegrowth of the excitation level. The problem of nonlinearresonant interaction of low-frequency vibrational modesof CNT together with the description of phenomena men-tioned above are the main subjects of this paper. A briefpreliminary discussion of the revealed phenomena waspresented in [14]. II. THE MODELA. Sanders-Koiter thin shell theory and itsmodification
The dynamics of carbon nanotubes (CNT) is one ofthe few areas of solid state physics, in which the classicaltheory of thin elastic shells (TTES) can be legitimatelyapplied. It is noteworthy that, in contrast to macro-scopic mechanics, where the fundamental limits of TTESare restricted by the possibility of plastic deformation,this theory can also be used for large displacements ofCNTs, even in the analysis of their collapse. The onlycomplicating factor is the uncertainty of the parametercharacterizing the thickness of the CNT [15]. The ap-plicability of a well-designed TTES allows us to obtainan effective description of the vibrational spectrum inthe framework of the linear approximation. This can be a r X i v : . [ c ond - m a t . m e s - h a ll ] D ec easily performed for the simplest of the boundary con-ditions, when a CNT of finite size can be considered asa part of an infinite CNT. However, the modified theorypresented below admits efficient study of both linear andnonlinear dynamics of CNTs under arbitrary boundaryconditions.We have found analytically the corresponding part ofthe oscillation spectrum and estimated the effect of theboundary conditions taking into account the existence ofboundary layers [16]. As for the non-stationary nonlinearproblems for CNT, in particular, resonance intermodalinteraction, they can be, in principle, studied by numer-ical methods [17–19]. However, such approach is unsuf-ficient when dealing with prediction of new phenomenaand the range of their manifestation. On the contrary, weshow that the analytical approach to nonlinear dynamicsof the CNT turns out to be efficient from this viewpoint.In the considered case an interaction of lowest (by eigenfrequency) CFMs, which is essential due to the effectivecrowding of the eigenvalues in the left range of the spec-trum, leads to existence of resonance non-stationary pro-cess. This analytical investigation is based on the re-duced nonlinear Sanders-Koiter thin shell theory, and itis a far going extension of our recent study relating tothe origin of the oscillation localization in the nonlinearlattices of various types [11–13].Two optical-type vibration branches in the CNT spec-trum (fig. 1) are of interest from the viewpoint of theenergy exchange processes mentioned above: the one ofthem is the well-known Radial Breathing Mode (RBM),which is associated with the circumferential wave num-ber n = 0 and corresponds to uniform radial extension-compression. The lowest optical mode in the CNT spec-trum (Circumferential Flexure Mode - CFM) [18] is spec-ified by n = 2, and the main deformation is a deviationof the CNT cross-section from the initial circular one[20]. To the best of our knowledge, the nonlinear dy-namic processes on CNTs were analytically studied onlyon the basis of a simplest modal analysis (RBM and itsparametric instability) [17–19].The CFM oscillations are characterized by relativesmallness of the ring and shear deformations (in particu-lar, the contour length of the lateral section is not essen-tially changed during deformation). So, we can assumethat only the bending, torsion and longitudinal deforma-tion contribute to the potential energy.The energy of elastic deformation of CNT in the di-mensionless units is written as follows: E el = 12 (cid:90) π (cid:90) (cid:18) ε ξ + ε ϕ + 2 νε ξ ε ϕ + 1 − ν ε ξϕ (cid:19) dξdϕ ++ β (cid:90) π (cid:90) (cid:18) κ ξ + κ ϕ + 2 νκ ξ κ ϕ + 1 − ν κ ξϕ (cid:19) dξdϕ, (1)where ε ξ , ε ϕ and ε ξϕ are the longitudinal, circumfer- FIG. 1. The CNT spectrum according to the exact Sanders-Koiter thin shell theory: solid curves correspond to circum-ferential wave number n = 0, dashed ones - to n = 1 anddot-dashed one - to n = 2. The insert shows the small wavenumber part of the CFM branch. All the frequencies ω aremeasured in dimensionless units and k - denotes the numberof longitudinal half-waves along the CNT. ential and shear deformations, and κ ξ , κ ϕ , κ ξϕ are thelongitudinal and circumferential curvatures, and torsion,respectively.Considering the small-amplitude oscillations of theCNT in the limiting case of a large aspect ratios, onecan write the following equation of motion in terms ofthe radial component of the displacement (see AppendixA for details) ∂ W∂τ + W − ε µω ∂ W∂ξ − ε γ ∂ W∂ξ ∂τ + ε κω ∂ W∂ξ + a W (cid:32)(cid:18) ∂W∂τ (cid:19) + W ∂ W∂τ (cid:33) + ε a ω (cid:18) ∂W∂ξ (cid:19) ∂ W∂ξ = 0 , (2)where W characterizes the radial displacement of theshell, α is the inverse aspect ratio of the CNT (i.e. theratio of the CNT radius R to its length L ), β is theratio of the thickness of the CNT wall h to its radius and ω is the gap frequency. Other parameters depend onthe circumferential wave number n and the Poisson ratio ν . ξ and τ = ω τ are the dimensionless coordinatesalong the CNT axes (0 ≤ ξ ≤
1) and the dimensionlesstime reduced to the gap frequency ω , respectively. Oneshould note that, taking into account the small value of ω ∼ β (cid:28) α (cid:28)
1, we clearlyassign the order of smallness of the different terms inequation (2). At the same time we formally consider thatthe coefficients of equation (2) to be of the zero order bysmall parameter.Equation (2) is a useful tool to analyze the effect ofvarious boundary conditions on the spectrum of natu-ral oscillations of the CNT [21] (see brief discussion inthe Appendix C). The frequency spectrum in the case ofsimply supported edges is determined by the followingexpression: ω = ω + µ π k + κπ k γπ k , (3)where k is a longitudinal wave number correspondingto the number of half-waves along the CNT axis.It is convenient to rewrite the equation (2) using com-plex variables:Ψ = 1 √ (cid:18) ∂W∂τ + iW (cid:19) W = − i √ − Ψ ∗ ) ∂W∂τ = 1 √ ∗ ) , (4)where the asterisk denotes the complex conjugation.Performing the multiscale expansion procedure (seeAppendix B) one can get the equation for the amplitudeof the main order in the ”slow” time τ = ε τ : i ∂χ ∂τ − µ − ω γ ω ∂ χ ∂ξ + κ ω ∂ χ ∂ξ − a | χ | χ = 0 , (5)where the main order value χ is coupled with the com-plex function Ψ = εχ exp ( − iτ ), and the small param-eter ε ∼ α .First of all, equation (5) admits the plane-wave solu-tion χ = A exp ( − i ( ωτ − kξ )) with the dispersion ratio ω = ( µ − ω γ ) k + κk ω − a A , (6)where A is the amplitude. As it can be seen, thisdispersion relation is in accordance with the relation (3).Equation (5) is the modified Nonlinear Schr¨odingerEquation (NLSE), in the the standard version of whichthe fourth derivative is absent. As it is well known, thestandard NLSE admits the localized solution - the en-velope soliton or the breather. The presence of fourth (a)(b)FIG. 2. (Color online) (a) The ”soliton” solutionof Eq. (5) at the variuos values of the ”frequency” ω ( − . , − . , − . , − . ω . derivative complicates the problem, but using Pade ap-proximation [22] one can obtain the following localizedsolution χ = X e − iωτ sech( λξ ) , (7)where λ = (cid:118)(cid:117)(cid:117)(cid:116) µ − γω − (cid:113) ( µ − γω ) + 8 κωω κX = 2 √ ω (cid:115) − λ ( µ − γω ) − ω ωa (8)Solution (7) describes a set of soliton-like excitations,which are parametrized by the ”frequency” parameter ω .The permissible values of ω are determined by the con-ditions of the reality of the magnitude (1 /λ ) (the solitonwidth) as well as of the amplitude X . Therefore, thesevalues have to be negative. It is a natural requirementbecause the localized solutions can exist in the gap of thevibration spectrum.Figures (2) show the solution (7) and its parameters λ and X at the various values of the frequency ω . Equa-tion (5) describes the nonlinear dynamics in the asymp-totic limit of the infinitely long CNT. Commonly speak-ing, it is the only limit when the localized soliton-likeexcitations occur in the nonlinear one-dimensional sys-tems, keeping in mind the boundary conditions in theinfinity. The modification of equation (5) is needed tocompare the analytical results and the data of numericalor physical experiments. One should note that the sourceof the nonlinearity in equation (5) is the inertial part ofthe energy. However, the finiteness of the CNT requirestaking into account the nonlinear part of the elastic defor-mation energy also. As it is mentioned in the AppendixA the most essential contribution arises from the nonlin-ear term ( ∂W/∂ξ ) ∂ W/∂ξ . Taking into account thiscontribution one can modify the equation (5): i ∂χ ∂τ − µ − ω γ ω ∂ χ ∂ξ + κ ω ∂ χ ∂ξ − a | χ | χ + a ∂∂ξ (cid:32)(cid:12)(cid:12)(cid:12)(cid:12) ∂χ ∂ξ (cid:12)(cid:12)(cid:12)(cid:12) ∂χ ∂ξ (cid:33) = 0 , (9)Equation (9) like equation (5) admits the plane-waveas well as the soliton-like solution (7). Only the effectiveamplitude of soliton X is modified as: X = 4 ω (cid:115) κ (9 λ ( µ − γω ) + 20 ω ω )3 a λ ( µ − γω ) − a κ − a ω ω ) (10)One should note that the low-frequency limit ω → λ = ω (cid:115) ωγω − µX = 4 (cid:114) − ωa (11)Nonlinear equation (9) can be used for the analysisof nonlinear normal modes interaction and, in particu-lar, for finding the transition between two regimes - theintensive energy exchange and energy localization [14].To perform this, one should take into account that thevibration spectrum for any CNT with a finite length isdiscrete, i.e. the longitudinal wave numbers are integers.To consider the intermodal interaction let us use the sumof the plane waves with the wave numbers k and k . χ = χ ( τ ) sin ( πk ξ ) + χ ( τ ) sin ( πk ξ ) (12)Substituting solution (12) into equation (9) one shoulduse the Galerkin procedure to obtain the equations forcomplex amplitudes χ and χ : i ∂χ ∂τ + δω χ − σ | χ | χ − σ (cid:16) | χ | χ + χ χ ∗ (cid:17) = 0 i ∂χ ∂τ + δω χ − σ | χ | χ − σ (cid:16) | χ | χ + χ χ ∗ (cid:17) = 0 , (13) where δω i = µ − ω γ ω π k i + κ ω π k i , i = 1 , ω of the considered brunch and σ ij = (cid:0) a + π a k i k j (cid:1) ( i, j = 1 , . (One can estimate that the frequency shift betweenthe lowest modes ( k = 1 , k = 2) is approximately twicesmaller than that for the next pair of modes ( k = 2 , k =3).)It is easy to see that the nonlinear terms in equations(13) are separated into two groups: the terms | χ j | χ i ( i, j = 1 ,
2) determine the nonlinear frequency shift,while the terms χ i χ ∗ j ( i (cid:54) = j ) describe the nonlinear in-teraction between modes. The Hamiltonian correspond-ing to equations (13) can be written as H = δω | χ | + δω | χ | − (cid:16) σ | χ | + σ | χ | (cid:17) − σ (cid:16) | χ | | χ | + (cid:0) χ ∗ χ + χ ∗ χ (cid:1)(cid:17) (14)Equations (13), besides the obvious energy integral(14), possess another integral X = | χ | + | χ | , (15)which characterises the excitation level of the system,and it is an analogue of the occupation number integralin quantum-mechanical terminology. B. LPT and the localilization of the CNTvibrations
As it was shown in [11, 14], the modal analysis be-comes inadequate at the resonance conditions. Thereforewe introduce new variables as the linear combinations ofresonating modes with preservation the integral X : φ = 1 √ χ + χ ); φ = 1 √ χ − χ ) . (16)The new variables describe the dynamics of some partsof the CNT [14] (similary to some groups (clusters) of theparticles in the effective discrete one-dimensional chain[11–13]). Considering the distribution of energy alongthe nanotube one can see that such a linear combinationof NNMs describes a predominant energy concentrationin certain region of the CNT, while the other part ofCNT has a lower energy. Because of small difference be-tween frequencies of the modes, the selected parts of CNTdemonstrate a coherent behavior similar to beating in thesystem of two weakly coupled oscillators. Therefore wecan consider these regions as new large-scale elementaryblocks, which can be identified as unique elements of thesystem - the ”effective particles” [11]. The existence ofintegral of motion (15) allows to reduce the dimension ofthe phase space up to 2 variables - θ and ∆, which char- acterize the relationship between the amplitudes and thephase shift between the effective particles, respectively: φ = √ X cos θe − i ∆ / ; φ = √ X sin θe i ∆ / . (17)Substituting these expressions into equations (13), theequations of motion in the terms of ”angular” variables( θ, ∆) can be obtained:sin 2 θ ( ∂θ∂τ −
12 ( δω − δω ) sin ∆ − X (3 σ (sin 2∆ sin 2 θ + 2 sin ∆)+3 σ (sin 2∆ sin 2 θ − − σ sin 2∆ sin 2 θ )) = 0sin 2 θ ∂ ∆ ∂τ + ( δω − δω ) cos ∆ cos 2 θ − X cos 2 θ (cos ∆(3 σ ( cos ∆ sin 2 θ + 1)+3 σ (cos ∆ sin 2 θ − − σ (cos 2∆ + 5) sin 2 θ ) = 0 (18)First of all, it is easy to show that equations (18) havetwo stationary points with coordinates ( θ = π/ , ∆ = 0)and ( θ = π/ , ∆ = π ). Taking into account the relations(17, 16), one can observe that these points correspondto the steady states χ and χ , respectively. All tra-jectories surrounding the stationary points describe theevolution of ”mixed” states with different contributionsof χ and χ . In particular, the lines θ = 0 and θ = π/ θ = 0 and θ = π/ π/ ± mπ ), where m = 0 , , . . . . On the otherside these points lie in the trajectories, which separatethe the normal modes attraction domains. Such a tra-jectory is the most distant from the stationary pointsand is named as the Limiting Phase Trajectory (LPT).As it was mentioned above the last notes the extremelynonuniform distribution of the energy (from a possibleones).The numerical solutions of Eq. (18) with the initialconditions corresponding to the immovable point ( θ = 0,∆ = π/
2) for the various values of the excitation X areshown in the Fig. 3 (a-f).Figs. 3(a, b) show the evolution of θ ( τ ) and ∆( τ ) forsmall value of X , when the system is close to the linearone. In this case one can see the non-smooth behavior ofthe relative amplitudes as well as of the phase shift of the”effective particles” φ and φ . Such a behavior corre-lates with that the any states belonging to the lines θ = 0or θ = π/ (cid:54) = ( π/ ± mπ ), in fact, are some ”vir-tual” ones, and they must be passed in the infinitesimal a Τ Θ b Τ- - - D c Τ Θ d Τ- - - D e Τ Θ f Τ- - - - - D FIG. 3. Time evolution of the variables θ, ∆ for the differentlevels of the CNT excitation: (a-b) X (cid:28) X loc ; (c-d) X =0 . X loc ; (e-f) X = 1 . X loc time.Figs. 3(c,d) demonstrate the behavior of the functions θ and ∆ if the excitation level X is large enough and isextremely close, but smaller than some threshold value X loc , the origin of which we will consider below. Thesefigures show that qualitetively this behavior qualitativelydoes not different from that in fig. 3(a,b).However, Figs.3(e,f) exhibit the drastic changes in theevolution of the functions θ and ∆ in spite of excitationof the system was increased on 0 .
5% only. First of all theinteval of variation of the function θ becomes twice less.It is the most important fact, which shows, that the statewith θ = 0 is inaccessible, if the initial conditions corre-spond to θ = π/ θ and ∆ and consider the topologyof the phase space. H = X δω (cos ∆ sin 2 θ + 1)+ δω (1 − cos ∆ sin 2 θ ) − X
64 ( σ (cos ∆ sin 2 θ + 1) + σ (1 − cos ∆ sin 2 θ ) + σ (cos 2∆ cos 4 θ + 4 sin ∆ + 10 cos θ ))) (19)Fig. 4 shows the phase portraits for various valuesof the parameter X . The initial structure of the phasespace for the small X , when the system is close to thelinear, is clearly seen in the left panel of Fig. 4. Therepresentative domains of the phase space are boundedby the intervals 0 ≤ θ ≤ π/ − π/ ≤ ∆ ≤ π/ χ and χ . The trajectory, which separates theattraction domains of different steady states and roundsthe normal mode χ , contains two lines θ = 0 and θ = π/
2, and two fragments, which connect the pairsof points: (( θ = 0 , ∆ = − π/ θ = π/ , ∆ = − π/ θ = 0 , ∆ = π/ θ = π/ , ∆ = π/ χ . The motionalong these trajectories leads to the non-smooth behavioras it is shown in the Fig. 3However, the steady state χ becomes unstable if theparameter X exceeds some threshold. Its value X ins canbe calculated from the condition of the instability: ∂ H∂θ |{ ∆=0 ,θ = π/ } = 0 ,X ins = 16 ( δω − δω ) σ − σ . (20)Two new stationary points arise after loosing the mode χ its stability. They are new nonlinear normal modes.The distance between them grows while the parameter X increases. These new stationary points correspond tosome non-uniform distribution of the energy along theCNT, however, this non-uniformity is weak. The mainfeatures of these states consists in that no trajectory sur-rounding them cannot attain the separatrix, which passesthrough the unstable stationary state χ . Therefore, thenon-uniformity of energy distribution remains for the in-finite time. Nevertheless, any trajectories, which are sit-uated in the gap between the separatrix and the LPT,preserve the possibility to pass from the vicinity of φ state ( θ = 0) into the vicinity of φ state ( θ = π/
2) (seeFig. 3(c, d). These process is accompanied with the slowenergy transfer from one part of the CNT to another one. However, the behavior of the solution of equations (18)is changed drastically if the value of X overcomes nextthreshold X loc . The existence of this threshold resultsfrom that the new stationary states move away from theunstable state and the separatrix grows while the LPTmoves to the unstable state in the vicinity of θ = π/ θ = π/ , ∆ = 0). At this moment the gap betweenthe LPT and the separatrix disappears and the only tra-jectory passed from θ = 0 to θ = π/ X leads to that newseparatrix, which is passed through the unstable station-ary points ( θ = π/ , ∆ = 0) and ( θ = π/ , ∆ = 2 π ),arises (see Fig. 4(c)). It separates the phase space of thesystem into uncoupled parts and any trajectories, whichstart near the θ = 0, cannot attain the value θ = π/ χ , φ and φ ,i.e. H ( θ = 0 , δ = 0) = H ( θ = π/ , δ = 0)= H ( θ = π/ , δ = 0) . (21)So, the value of the localization threshold turns out tobe X loc = 64 ( δω − δω )3 σ − σ − σ (22)Fig. 6 demonstrates the dependence of the localiza-tion threshold in the terms of the radial displacement w from the inverse aspect ratio of the CNT. The solid curveshows the threshold value in the accordance with equa-tion (22). To compare with it, the threshold, which hasbeen estimated on the base of equation (5), is drawn bythe dashed line. One can see that the asymptotic valuesfor the long CNTs (at α →
0) of both thresholds are thesame, but the difference between them becomes essentialfor the finite-length CNT.
C. Comparison of two-mode approximation withother numerical methods
The two-mode approximation in the framework of non-linear Sanders-Koiter theory of thin elastic shells allowsus to predict the bifurcation of dynamical behavior oflow-frequency CNT vibrations as well as to estimate thethreshold values of oscillation amplitude. However, the (a) (b) (c)FIG. 4. (Color online) Phase portraits of the system with Hamiltonian (19) for different values of the excitation: (a) X (cid:28) X loc ,(b) X = 0 . X loc , (c) X = 1 . X loc .(a) (b) (c)FIG. 5. (Color online) (a) The energy distribution along the CNT during the MD-simulation with the small excitation level X ,that corresponded to phase portrait in the Fig. 4(a). (b - c) The same as in the panel (a) with initial excitation level X (cid:46) X loc and X > X loc , respectively. Α X loc FIG. 6. The threshold of vibration localization for the CNTwith radius 0.79 nm vs inverse aspect ratio α : solid black andthe dashed red curves are the pedictions, which based on theEqs. (9) and (5), respectively. The dotted blue curve showsthe threshold value estimated by the numerical method (seesection II C for detail). influence of the other part of the spectrum is very impor-tant for the estimation of the reliability of the obtainedresults. Therefore, they should be verified by the in-dependent numerical methods. One of the approachesconsists in the direct numerical integration of the modalnonlinear equations of the Sanders-Koiter thin shell the-ory.In order to carry out the numerical analysis of the CNTdynamics, a two-step procedure was used: i) the displace-ment field was expanded by using a double mixed series,then the Rayleigh-Ritz method was applied to the lin-earized formulation of the problem, in order to obtain an approximation of the eigenfunctions; ii) the displacementfields are re-expanded by using the linear approximatedeigenfunctions, the Lagrange equations were then con-sidered in conjunction with the nonlinear elastic strainenergy to obtain a set of nonlinear ordinary differentialequations of motion.So, to satisfy the boundary conditions the displace-ment field was expanded into series r ( ξ, ϕ, t ) = (cid:34) M u (cid:88) m =0 N (cid:88) n =0 R m,n T ∗ m ( ξ ) cos nϕ (cid:35) f ( t ) (23)where function r ( ξ, ϕ, t ) substitutes the displacements u , v or w .In equations (23) T ∗ m ( ξ ) = T m (2 ξ −
1) are the Cheby-shev orthogonal polynomials of the m − th order, n is thenumber of nodal diameters, and f ( t ) describes the timeevolution of the CNT vibrations.The maximum number of variables needed for describ-ing a general vibration mode with n = 2 nodal diame-ters (Circumferential Flexural Mode) is obtained by therelation ( N max = M u + M v + M w + 3 − p ), where ( M u = M v = M w ) denote the maximum degree of theChebyshev polynomials and p describes the number ofequations for the boundary conditions to be respected.A specific convergence analysis was carried out to se-lect the degree of the Chebyshev polynomials: degree 11was found suitably accurate, ( M u = M v = M w = 11).In the cases of a SWCNT with simply supported orclamped edges ( p = 8), the maximum number of degreesof freedom of the system with is equal to ( N max = 33 +3 − p = 0), the maximum number of degrees of freedom ofthe system is equal to ( N max = 33 + 3 = 36).The equations (23) are inserted into the expressionsof the potential energy E el and kinetic energy T to com-pute the Rayleigh quotient R ( q ) = E el, max /T ∗ , where E el,max = max( E el ) is the maximum of the poten-tial energy during a modal vibration, T ∗ = T max /ω , T max = max( T ) is the maximum of the kinetic en-ergy during a modal vibration, ω is the circular fre-quency of the synchronous harmonic motion and q =[ ..., U m, , V m, , W m, , ... ] T represents a vector containingall the unknown variables.After imposing the stationarity to the Rayleigh quo-tient, one obtains the eigenvalue problem( − ω M + K ) q = 0 (24)which gives approximate natural frequencies (eigenval-ues) and modes of vibrations (eigenvectors). The re-sults of performed calculation show that the eigenspec-trum values are in the good accordance with the estima-tions made in the framework of reduced Sanders-Koitertheory discussed above. The specific difference betweeneigenvalues amounts to the values 2 −
4% for the long-wave modes and reachs up to 20% while the longitudinalwavenumber grows [16].In the nonlinear analysis, the full expression of thedimensionless potential energy E el containing terms upto the fourth order (cubic nonlinearity), is considered.In the cases of simply supported and clamped bound-ary conditions, the two low-frequency optical-type cir-cumferential flexure modes ( m = 1 , n = 2) and ( m =2 , n = 2) are considered.Using the Lagrange equations ddτ (cid:18) ∂T∂q (cid:48) i (cid:19) + ∂E el ∂q i = 0 , (25)a set of nonlinear ordinary differential equations is ob-tained; these equations must be completed with suitableinitial conditions on displacements and velocities. Thissystem of nonlinear equations of motion was finally solvedby using the implicit Runge-Kutta numerical methodswith suitable accuracy, precision and number of steps.The solution of nonlinear equations with initial condi-tions in the vicinity of the bifurcation threshold showsthe coincidence of the threshold values in the analyticalmodel and the numerical one for the wide interval of as-pect ratios (see Fig. 6). The procedure and results willbe discussed in the nearest future. FIG. 7. (Color online) Snapshot of typical energy distributionalong the deformed CNT during the vibration associated withCFM spectrum branch.
D. MD simulation
To verify the results of analytical model the simulationof the low-frequency vibrations of CNTs was performedby molecular dynamics (MD) techniques using realisticinter-atomic potential functions. Classical molecular dy-namics technique which uses predefined potential func-tions (force fields) was applied for the calculation of thetotal potential energy of the system. The typical MD ex-periment consisted of several stages. At the first stage theCNT was kept at high temperature ( (cid:39) K ) for struc-tural relaxation. Then the termostat temperature wasdecreased with a constant rate down to approximately1 K with a subsequent low-temperature relaxation. Thethird stage dealt with CNT deformation according to an-alytical solution with subsequent relaxation. The secondversion of initial conditions were given by initial veloci-ties of atoms at zero initial displacements. After that theexternal field was turned off, and the free natural oscil-lations of CNT with the fixed boundary conditions wererealized. In accordance with analytical description, theatoms at the edges of CNT were fixed by the force fieldagainst any radial displacements ( W (0 , t ) = W (1 , t ) = 0)that corresponds to the boundary conditions similar tosimply supported shell. The typical snapshot of distri-bution of the CNT deformation energy during the MDsimulation is shown in fig.(7).The consequent analysis of MD simulation data in-cluded the control of natural frequencies and energy dis-tribution along the CNT axis via variation of the oscilla-tion amplitude. The 3D pictures of energy distributionalong the CNT axis measured during the MD simulationshave been discussed in the section II B (see the figures 5(a-c)).Fig. 8 shows the variation of the CNT vibration spec-trum with changing of the initial excitation level. Thedot, black and red curves correspond to the excitation X < X loc .One should note, that the large narrow peak near thefrequency ω (cid:39) .
05 corresponds to the lowest eigenvalueof the system under consideration. Therefore no vibra-tion with the smaller frequency exists in the gap for theany initial excitation, if they do not exceed the local-ization threshold X loc . However, the overcoming thethreshold of localization changes the spectrum drasti- FIG. 8. (Color online) CNT vibration spectra at the variousexcitation level: X = (0 . − . X loc . cally. From one side, one can see in Fig. 5(c), that theshape of the initial excitation is deformed essentially dur-ing the MD simulation and the local temperature reachesa great value ( ∼ K ). Therefore, the onset of high-frequency modes is naturally sufficient. On the otherside, the intensive oscillations in the gap of the spec-trum can not be explained by the increasing of the localtemperature. To fill this gap, the excitation of anotherlow-frequency modes is required (they may be acoustictype modes like the bending or the longitudinal stretch-ing modes - see Fig. 1). The another possibility is form-ing of the localized excitation, because the Fourier spec-trum is wide enough. Unfortunately, the two-mode ap-proximation used in our analysis can not answer in thisquestion. The accurate study of this problem needs inthe consideration of the inter-brunch mode interaction.This problem will be formulated in our future studies. III. CONCLUSION
Instability and bifurcation of the edge-spectrum opti-cal NNM at the value X = X ins leads to appearance ofthe localized NNMs with stationary (in slow time scale)energy localization in some part of the CNT. In contrastto a breather, this is not a strong localization becausethe two-mode approximation, which is valid for consid-ered aspect ratio, can reveal a weak localization only.We demonstrate that instability of edge-spectrum opti-cal modes of CNT vibrations is the preliminary conditionof non-stationary (in slow time scale) energy localizationin the some domain of CNT. The energy capture in oneof the CNT parts can be achieved, if the excitation levelexceeds the specific threshold X = X loc , which corre-sponds to merging two trajectories, which are the LPTand the separatrix appeared at X = X ins . When this threshold is exceeded the phase portrait of the systemunder consideration changes drastically: the separatrixpassing through the unstable stationary point ( θ = π/ θ = π/
4, ∆ = π ) and prevents full energy ex-change between effective particles ψ and ψ . Simulta-neously a set of transit-time trajectories, which involveany values of phase difference ∆, is created. It meansthat initial conditions corresponding to identical veloci-ties or displacements of both modes lead to the energycapture by the effective particles. Then only a partialenergy exchange becomes possible along the trajectories,surrounding the stable stationary point and situated in-side the separatrix. One should keep in mind that theprocess of energy capture does not suggest the creationof strongly localized solutions whose formation requires aparticipation of more components of the spectrum. Thiscan be achieved for CNT with larger aspect ratio.It should be noted once more that the developmentand the use of the analytical framework based on theLPT concept is motivated by the fact that resonant non-stationary processes occurring in a broad variety of fi-nite dimensional physical models are beyond the well-known paradigm of nonlinear normal modes (NNMs),fully justified only for quasi-stationary processes and non-stationary processes in non-resonant case. While theNNMs approach has been proved to be an effective toolfor the analysis of instability and bifurcations of station-ary processes (see, e.g., [23]), the use of the LPTs conceptprovides the adequate procedures for studying stronglymodulated regimes as well as the transitions to energylocalization and chaotic behavior [11]. Such an approachclarifies also the physical nature of the breathers forma-tion in infinite discrete or continuum systems.As a conclusion we would like to note that the phe-nomenon of energy localization considered above has uni-versal character and it is the common peculiarity of thesystems possessing the optical-type branches of vibra-tional spectrum. However, as it was studied early [11–13],the occurence of the localization depends on the typesof nonlinearity as well as on the relations between coef-ficients σ i,j in the Hamiltonian (14). If some ratios be-tween coefficients σ i,j are sa tisfied, an additional integralof motion arises that leads to the effective linearizationof the equations of motion and, as consequence, to theabsence of the localization processes [24]. So, in spiteof that the interaction of resonating nonlinear modes de-scribes by the Hamiltonian (14) for a wide class of thenonlinear systems, the results of this interaction may varyconsiderably. In any case, the analysis of the Hamilto-nian (14) in combination with the LPT concept gives ususefull tool for the study of nonlinear systems. ACKNOWLEDGMENTS
The work was supported by Russia Basic ResearchFoundation (grant 08-03-00420a) and Russia Science0Foundation (grant 14-17-00255)
Appendix A: The reduced Sanders- Koiter thin shelltheory
It is convenient to use the dimensionless variableswhich determine the elastic deformation of circular thinshell. In such a case all components of the displacementfield ( u - longitudinal along the CNT axis, v - tangentialand w - radial displacement, respectively) are measuredin the units of CNT radius R . The displacements and re-spective deformations refer to the middle surface of theshell. The coordinate along the CNT axis ξ = x/L ismeasured via the length of nanotube and variates from 0up to 1, and θ is the azimuthal angle.One can define the dimensionless energy and time vari-ables, which are measured in the units E = Y RLh/ (1 − ν ) and t = 1 / (cid:112) Y /ρR (1 − ν ), respectively. Here Y isthe Young modulus of graphene sheet, ρ - its mass den-sity, ν - the Poisson ratio of CNT, and h is the effectivethickness of CNT wall. There are two dimensionless geo-metric parameters which characterize CNT: the first oneis inverse aspect ratio α = R/L and the second - effectivethickness shell β = h/R .The energy of elastic deformation of CNT in the di-mensionless units is written as follows: E el = 12 (cid:90) π (cid:90) ( N ξ ε ξ + N ϕ ε θ + N ξϕ ε ξϕ + M ξ κ ξ + M ϕ κ ϕ + M ξϕ κ ξϕ ) dξdϕ (A1)where ε ξ , ε ϕ and ε ξϕ are the longitudinal, circumferen-tial and shear deformations, and κ ξ , κ ϕ and κ ξϕ are thelongitudinal and circumferential curvatures, and torsion,respectively. The respective forces and momenta may bewritten in the physically linear approximation: N ξ = ε ξ + νε ϕ , N ϕ = ε ϕ + νε ξ ,N ξϕ = 1 − ν ε ξϕ M ξ = β
12 ( κ ξ + νκ ϕ ) , M ϕ = β
12 ( κ ϕ + νκ ξ ) ,M ξϕ = β
24 (1 − ν ) κ ξϕ (A2)One should note that both curvatures and torsion arethe dimensionless variables in accordance with our defi-nition of dispacement field ( u, v, w ).The Sanders-Koiter approximation of defectless thinshell allows to write the nonlinear deformations ( ε ) andcurvatures ( κ ) in the following form ε ξ = α ∂u∂ξ + α ∂w∂ξ ) + 18 ( α ∂v∂ξ − ∂u∂ϕ ) ε ϕ = ∂v∂ϕ + w + 12 ( ∂w∂ϕ − v ) + 18 ( ∂u∂ϕ − α ∂v∂ξ ) ε ξϕ = ∂u∂ϕ + α ∂v∂ξ + α ∂w∂ξ ( ∂w∂ϕ − v ) (A3) κ ξ = − α β ∂ w∂ξ ,κ ϕ = β (cid:18) ∂v∂ϕ − ∂ w∂ϕ (cid:19) ,κ ξϕ = β (cid:18) − α ∂ w∂ξ∂ϕ + 3 α ∂v∂ξ − ∂u∂ϕ (cid:19) . (A4)One should make some physically grounded relation-ships between the displacement components to simplifythe description of the CNT nonlinear dynamics. We con-sider the low-frequency optical-type vibrations which arespecified by circumferential wave number n = 2. Thisbranch is characterized by relatively small circumfer-ential and shear deformations, while the displacementsthemselves may not be small. In such a case we canwrite: ε ϕ = 0; ε ξϕ = 0 (A5)The components of displacement field are u ( ξ, ϕ, τ ) = U ( ξ, τ ) + U ( ξ, τ ) cos( nϕ ) v ( ξ, ϕ, τ ) = V ( ξ, τ ) sin( nϕ ) w ( ξ, ϕ, τ ) = W ( ξ, τ ) + W ( ξ, τ ) cos( nϕ ) (A6)These relations allow us to express the longitudinal andtangential components, and axially symmetric part of theradial displacement via the radial one. Correspondingrelationships can be written as folows: V ( ξ, τ ) = − n W ( ξ, τ ); U ( ξ, τ ) = − αn ∂W ( ξ, τ ) ∂ξW ( ξ, τ ) = − n (( n − W ( ξ, τ ) + α ( ∂W ( ξ, τ ) ∂ξ ) ); ∂U ( ξ, τ ) ∂ξ = − n + 14 n α ( ∂W ( ξ, τ ) ∂ξ ) (A7)Because the kinetic energy contains the inertial termscorresponding to all components of the deformation field E kin = 12 (cid:90) π (cid:90) (cid:18) ( ∂u∂τ ) + ( ∂v∂τ ) + ( ∂w∂τ ) (cid:19) dξdϕ (A8)1we need in taking into account the relations (A7) also.Omitting the calculation details one can write the fi- nal equation of motion in terms of radial displacement W ( ξ, t ): ∂ W∂τ + ω W − µ ∂ W∂ξ − γ ∂ W∂ξ ∂τ + κ ∂ W∂ξ + a W (cid:32)(cid:18) ∂W∂τ (cid:19) + W ∂ W∂τ (cid:33) + a (cid:18) ∂W∂ξ (cid:19) ∂ W∂ξ + a (cid:32) ∂W∂τ ∂W∂ξ ∂ W∂ξ∂τ − W (cid:18) ∂ W∂ξ∂τ (cid:19) + (cid:18) ∂W∂τ (cid:19) ∂ W∂ξ (cid:33) + a (cid:20) ∂ W∂ξ∂τ (cid:18) ∂ W∂ξ∂τ ∂ W∂ξ + 2 ∂W∂ξ ∂ W∂ξ ∂τ (cid:19) + ∂W∂ξ (cid:18) ∂W∂ξ ∂ W∂ξ ∂τ + 2 ∂ W∂ξ∂τ ∂ W∂ξ (cid:19)(cid:21) = 0 (A9)where ω = β n ( n − n + 1) , µ = α β ( n − n − ν )6( n + 1) ; γ = α n ( n + 1) , κ = α (12 + n β )12 n ( n + 1) ∼ α n ( n + 1) ; a = ( n − n ( n + 1) , a = 2 α ( n − n ( n + 1) ; a = α ( n − n ( n + 1) , a = α n ( n + 1)(A10)The estimation of the different terms of equation (A9)shows that the essential contribution get the first and thesecond nonlinear terms only. In further we skip the lastnonlinear terms.Eq. (A9) allows us to calculate the eigenfrequencies inthe linear approximation as well as to estimate the effectof nonlinearity on these frequencies at different boundaryconditions. The detail analysis shows that two first non-linear terms in the equation (A9) give the dominant con-tribution in the low-frequency dynamics of CNTs. Oneshould note, that the parametrs α and β are small enoughin most cases. Therefore, the parameter κ is very closeto the γ . Moreover, the own frequency of the gap ( ω issmall due to the smalness of the effective thickness of theCNT wall). In such a case it is convenient to introducethe ’new’ time scaled by the gap frequency ω : τ = ω τ .Taking into account that the coefficients µ, γ, a , a areof the order of unity, one can rewrite the equation (A9)keeping the dominant terms with the order of small pa-rameters, which do not exceed two: ∂ W∂τ + W − µω ∂ W∂ξ − γ ∂ W∂ξ ∂τ + κω ∂ W∂ξ + a W (cid:32)(cid:18) ∂W∂τ (cid:19) + W ∂ W∂τ (cid:33) − a ω (cid:18) ∂W∂ξ (cid:19) ∂ W∂ξ = 0 , (A11) Appendix B: The multiscale expansion
Because we consider the small-amplitude oscillations,one can represent the complex amplitude Ψ as a series ofsmall parameter ε :Ψ = ε (cid:0) ψ + εψ + ε ψ + . . . (cid:1) (B1)Next we should introduce the ’time’ series: τ = ετ , τ = ε τ , . . . and the respective time derivatives: ∂∂τ = ∂∂τ + ε ∂∂τ + ε ∂∂τ + . . . (B2)Substituting the expansion (B1) into (2) and takinginto account the hierarchy of the times, we get the equa-tions in the different orders by small parameter ε . ε : i ∂ψ ∂τ − ψ = 0So, we get ψ = χ ( ξ, τ , τ ) e − iτ ε : i ∂ψ ∂τ − ∂ψ ∂τ − ψ = 0Then we get: ψ ( ξ, τ , τ , τ ) = χ ( ξ, τ , τ ) e − iτ ψ ( ξ, τ , τ , τ ) = χ ( ξ, τ ) e − iτ . (B3)2 ε : i ∂ψ ∂τ − ψ + i ∂ψ ∂τ + i ∂ψ ∂τ − µ ω ∂ ψ ∂ξ − γ ∂ ψ ∂ξ − i γ ∂ ψ ∂τ ∂ξ + κ ω ∂ ψ ∂ξ + µ ω ∂ ψ ∗ ∂ξ − γ ∂ ψ ∗ ∂ξ − κ ω ∂ ψ ∗ ∂ξ + a (cid:18) | ψ | ψ − i ∂ψ ∂τ (cid:0) ψ + ψ ∗ − | ψ | (cid:1)(cid:19) = 0 , (B4)Taking into account the relations (B3), one can inte-grate the equation (B4) with respect to τ and τ . Thenthe condition of the secular terms absence give us theequation for the main approximation amplitude χ (5). Appendix C: Influence of boundary conditions
The presence of boundary conditions different fromthe simple supporting affects on the NNM and theirfrequency. In this Appendix we consider the effectivemethod for the solution of boundary problem and demon-strate the procedure of the normal mode construction onthe example of CNT with free edges.It is intuitively clear that the strong boundary condi-tions like the clamping lead to frequency growth whilethe more ”soft” conditions can decrease the frequencies.To estimate the variation of normal modes we used thelinear approximation of RSKTST (reduced Sanders Koi-ter thin shell theory).Let us assume that the solution of linearized equationof the the CNT vibrations ∂ W∂τ + ω W − µ ∂ W∂ξ − γ ∂ W∂ξ ∂τ + κ ∂ W∂ξ = 0 (C1)is represented as the periodic process W ( ξ, τ ) ∼ f ( ξ ) cos( ωτ ) . (C2)Taking into account expression (C2) explicitly, one canrewrite equation (C1) with the help of the product of twodifferential operators: κ (cid:18) d dξ + k (cid:19) (cid:18) d dξ − λ (cid:19) f = 0 , (C3) where the parameters µ , γ , and ω are linked by therelationships κλ k = ω − ω κ (cid:0) λ − k (cid:1) = µ − γω . (C4)Because the operators (cid:0) d /dξ + k (cid:1) and (cid:0) d /dξ − λ (cid:1) are commutative ones, any function f ( ξ ), which satisfies one of the equations (cid:18) d dξ + k (cid:19) f ( ξ ) = 0 (C5) (cid:18) d dξ − λ (cid:19) f ( ξ ) = 0 , (C6)is a solution of equation (C3).So a general solution of equation (C3) is a linear com-bination of the solutions of equations (C5), (C6): f ( ξ ) = ( C sin k ( ξ − ξ )+ C e λ ( ξ − ξ ) + C e − λ ( ξ − ξ ) ) , (C7)where k , λ , C j , j = 1 , , ξ j , j = 0 , W ( ξ, τ ): α ∂ W ( ξ, τ ) ∂ξ − ν (cid:0) n − (cid:1) W ( ξ, τ ) = 0 α ∂ W ( ξ, τ ) ∂ξ − (2 − ν ) (cid:0) n − (cid:1) ∂W ( ξ, τ ) ∂ξ = 0 ξ = 0 , [1] C. Li and T.-W. Chou, Phys. Rev. B , 073405 (2003).[2] V. Sazonova, Y. Yaish, H. stnel, D. Roundy, T. A. Arias, and P. L. McEuen, Nature , 284 (2004).[3] H. Peng, C. W. Chang, S. Aloni, T. D. Yuzvinsky, and FIG. 9. The comparison of CNT vibration spectra for thedifferent boundary conditions. CNT parameters - L = 10 nmand R = 0 .
79 nm.A. Zettl, Phys Rev Lett , 087203 (2006).[4] M. Anantram and F.Leonard, Rep. Prog. Phys. , 507(2006).[5] S. Berber, Y.-K. Kwon, and D. Tomanek, Phys. Rev.Lett. , 4613 (2000).[6] J. Wang and J.-S. Wang, Appl. Phys. Lett. , 111909(2006).[7] N. Mingo and D. A. Broido, NanoLetters , 1221 (2005).[8] B.Li, J.Wang, L.Wang, and G.Zhang, Chaos , 015121(2005).[9] A.V.Savin, B.Hu, and Yu.S.Kivshar, Phys. Rev. B ,195423 (2009).[10] E.Fermi, J.Pasta, and S.Ulam, in Collected Papers ofEnrico Fermi , Vol. 2, edited by E. Segre (University ofChicago Press, Chicago, 1965) p. 978.[11] L.I.Manevitch and V.V.Smirnov, Phys. Rev. E ,036602 (2010).[12] L.I.Manevitch and V.V.Smirnov, Doklady Physics ,324 (2010). [13] V.V.Smirnov and L.I.Manevitch, Acoustical Physics ,271 (2011).[14] V. V. Smirnov, D. S. Shepelev, and L. I. Manevitch,Phys. Rev. Lett. , 135502 (2014).[15] Y. Huang, J. Wu, and K. C. Hwang, Phys. Rev. B ,245413 (2006).[16] M. Strozzi, L. I. Manevitch, F. Pellicano, V. V. Smirnov,and D. S.Shepelev, Journal of Sound and Vibration ,2936 2957 (2014).[17] M. X. Shi, Q. M. Li, and Y. Huang, Proc. R. Soc. A , 1941 (2008).[18] M. X. Shi, Q. M. Li, and Y. Huang, Proc. R. Soc. A , 03069 (2009).[19] P. Soltani, J. Saberian, R. Bahramian, and A. Farshidi-anfar, IJFPS , 47 (2011).[20] A.M.Rao, E. Richter, S. Bandow, B. Chase, P. C. Eklund,K. A. Williams, S. Fang, K. R. Subbaswamy, M. Menon,A. Thess, R. E. Smalley, G. Dresselhaus, and M. S.Dresselhaus, Science , 187 (1997).[21] V.V.Smirnov, D.S.Shepelev, and L.I.Manevitch, in XIVconference on Polymers and Composite Materials , editedby S.A.Dubrovski, Institute of Chemical Physics RAS(ICP RAS, Moscow, Russia, 2013) p. 42, in Russian.[22] J. George A.Baker and P. Graves-Morris,
Pade Approx-imants , edited by G.-C. Rota, Encyclopedia of Mathe-matics and Its Application, Vol. 13, 14 (Addison-WesleyPublishing Co, Reading, Massachusetts, 1981).[23] A. F. Vakakis, L. I. Manevitch, Y. V. Mikhlin, V. N.Pilipchuk, and A. A. Zevin,
Normal Modes and Local-ization in Nonlinear Systems (Wiley, New York, 1996) p.552.[24] L. I. Manevitch and V. V. Smirnov, “Resonant energy ex-change in nonlinear oscillatory chains and limiting phasetrajectories: from small to large system.” in