Molecular Simulations of Cyclic Loading Behavior of Carbon Nanotubes Using the Atomistic Finite Element Method

The potential applications of carbon nanotubes (CNT) in many engineered bionanomaterials and electromechanical devices have imposed an urgent need on the understanding of the fatigue behavior and mechanism of CNT under cyclic loading conditions. To date, however, very little work has been done in this field. This paper presents the results of a theoretical study on the behavior of CNT subject to cyclic tensile and compressive loads using quasi-static molecular simulations. The Atomistic Finite Element Method (AFEM) has been applied in the study. It is shown that CNT exhibited extreme cyclic loading resistance with yielding strain and strength becoming constant within limited number of loading cycles. Viscoelastic behavior including nonlinear elasticity, hysteresis, preconditioning (stress softening), and large strain have been observed. Chiral symmetry was found to have appreciable effects on the cyclic loading behavior of CNT. Mechanisms of the observed behavior have been revealed by close examination of the intrinsic geometric and mechanical features of tube structure. It was shown that the accumulated residual defect-free morphological deformation was the primary mechanism responsible for the cyclic failure of CNT, while the bond rotating and stretching experienced during loading/unloading played a dominant role on the strength, strain and modulus behavior of CNT.


Introduction
Carbon nanotubes (CNTs), originally discovered as a byproduct of fullerene research [1,2], have attracted numerous interests due to the high potential use of these quasi-1D structures as superstrong light-weight materials, novel nanometer-scale mechanical and electronic devices, and so forth.So far, large success has been achieved on the understanding of nanoscale strength, plastic deformation, and failure modes of CNTs under monotonic loading conditions as well as their superior electromagnetic properties [1][2][3][4][5][6][7][8][9][10][11][12][13][14].However, an apparent lack of research efforts exists on the study of fatigue failure and mechanism on CNT subject to cyclic loading conditions.Recently, demand on the knowledge in this topic is increasing rapidly because understanding of fatigue behavior of CNT is essential for their applications in a variety of engineered nanoscale materials and devices where nanotubes may be subject to 10 6 ∼ 10 9 cycles during their life time [15][16][17].For example, when used as electrical contact probes, vertically aligned bundles of CNT will need to be created and required to perform consistently over a large number of touchdown cycles.Another highly promising field is the synthesis of biomaterials using CNT as nanoscale components where soft tissues (e.g., stomach wall or skeletal muscle) will stretch and shrink over a huge number of cycles during their life time.To this end, it will be of great interest to know theoretically whether the superior structural and mechanical properties of CNT will be retained under cyclic loading conditions.
In literature, the limited studies on this topic are focused on the fatigue failure strength and mechanisms of nanotubes in CNT/epoxy composites under cyclic tensile load [18][19][20][21].In these studies, cyclic loading is applied by subjecting CNT to negative and positive voltages at a specified frequency using in situ TEM.A most relevant work recently reported on fatigue resistance of vertically aligned CNT arrays under cyclic compression shows viscoelastic behavior characteristic of soft-tissue materials [22].No fatigue failure is observed under cyclic compressive loading up to half a million cycles at high strain levels.Super-compressible foamlike behavior of CNT films has also been observed under cyclic compression up to 2000 cycles, with the buckling unfolded and recovered to their near original length [23].The high cyclic stress of CNT, as well as its remarkable fatigue resistance, suggests its great potential as building blocks in synthetic fatigue resistant biomaterials and high fracture toughness polymer composites [18,20,22].
In this paper, a theoretical study of the yielding strength and strain of carbon nanotubes subject to cyclic axial compressive and tensile load has been conducted.Classical molecular dynamics simulations employing the Tersoff-Brenner many-body potential [24,25] have been used in this study.This potential accurately represents the lattice constants, binding energies, and elastic constants of both graphite and diamond and has proven capable of reproducing brittle and ductile failure patterns of CNT with very good agreement with results obtained using more accurate ab initiosimulations [4].For cyclic compression, the yielding strength and strain of CNT is governed by structural instability, that is, local and global buckling under axial loads, whose influence on the cyclic behavior of zigzag tubes was investigated.While for cyclic tension, since no thermodynamic effects are considered in the current study, the purely mechanical bond-breaking failure mode governs the yielding behavior of CNT.The cyclic strength and strain was obtained by loading the tubes next to the first bond-breaking, which is irreversible once it occurs.Parametric studies were also made to investigate the effects of chiral symmetry on the cyclic strength and strain of CNT.The justification and consequences of different loading methods employed for tubes with varying chiral symmetries are also discussed.

AFEM Method.
In the current study, the Atomistic Finite Element Method (AFEM), an innovative MD simulation technique developed by Liu et al. [26], is used to optimize the structure subject to external loads until a stable configuration with minimal energy is reached.In AFEM, individual atoms are treated as discrete nodes in conventional FEM, with their equilibrium positions found by minimizing the interatomic potentials.The governing equation of AFEM identical to that of conventional FEM is given by where K, u, and P are stiffness matrix, displacement vector and nonequilibrium force vector, respectively.Specifically, K and P are obtained by making the second and first derivatives of the total interatomic potential U tot of the system with respect to the positions x of all N atoms at the current positions x (0) , respectively, where U tot = N i< j U i j , and F is the external force vector.The most distinct feature of AFEM is its ability to account for the nonlocal, many-body effects of the atomic system.This is achieved by including in a single AFEM element all neighboring atoms that will interact with a certain central atom due to many-body potentials of the system.For example, for the hexagonal lattice structure of CNT, one central carbon atom will interact with both of its three nearest and six second nearest neighbor atoms, resulting in an AFEM element which includes a total number of 10 atoms but only focuses on the central atom.Another desirable feature of AFEM is that it is an N-order method, which means that the computational effort scales with the total number of atoms in the system.This feature makes AFEM much more computationally efficient than the conventional N 2 -order conjugate gradient method [26].
For a nonlinear system, (1) is solved iteratively until P approaches zero.Since the thermodynamic features of the system, such as particle mass, particle velocity, and temperature, are not considered, the time-evolution process of the system change cannot be well represented.Instead, it only gives the final stabilized solution of the system.Specifically, for the current study of CNT under cyclic loading, a stabilized deformed structure will only be found if the structure optimization converges.However, this will satisfy the purpose of the current study since we are only interested in the cyclic strength and strain of CNT before any nonconverging structural failure occurs.

AFEM Method.
The interatomic potential for carbon was proposed by Tersoff [24] and Brenner [25] as where r i j is the distance between atoms i and j; V R and V A are the repulsive and attractive pair terms given by where D (e) , S, β, and R (e) are parameters defined in the modified Morse potential.In particular, D (e) is the equilibrium well depth, and R (e) is the equilibrium interatomic distance.When S = 2, the pair terms reduce to the usual Morse potential.f c (r) is a smooth cutoff function to prescribe the range of the potential and is given by with R (1) = 0.17 nm and R (2) = 0.2 nm.The parameter B i j in (3) represents a multibody coupling between the bond from atom i to atom j and the local environment of atom i and is given by where r ik is the distance between atom i and atom k, θ θi jk is the angle between bonds i-j and i-k, and G is a function given by For atoms i and j having different local environment, Brenner [25] suggested to replace the coefficient B i j in ( 7) by The parameters D (e) , S, β, and R (e) in ( 4) and ( 5), δ in (7), and a 0 , c 0 , and d 0 in (8) have been determined by Brenner [25] by fitting the binding energy and lattice constants of graphite, diamond, simple cubic and facecentered-cubic structures for pure carbon as well as the vacancy formation energy for diamond and graphite.The values of these parameters are D (e) = 6.000eV, S = 122, β = 21 nm −1 and R (e) = 0.1390 nm, δ = 0.50000, a 0 = 0.00020813, c 0 = 330, and d 0 = 3.5.

Molecular Simulation Program.
The CNT samples with varying chiral symmetries used in the parametric study are summarized in Table 1.The chiral angle ϕ is defined as the angle between the actual roll-up vector and the zigzag rollup vector.Alternatively, chirality can be specified by a pair of lattice indices (m, n), based on which the chiral angle ϕ, and tube diameter δ can be calculated as [27] The sizes of the samples vary little, with the length and diameter nearly equal to 4 nm and 1 nm, respectively, so that the effects of chiral symmetries on the cyclic behavior can be isolated.Given that the yielding strain of CNT under tension is very high (about 10 times as that under compression), the diameters of the zigzag and armchair tubes used in cyclic tension simulations are selected to be slightly smaller than those used in cyclic compression simulations.This variation, however, has little influence on the simulated cyclic behavior and effects of chiral symmetry.
For all the samples, the tube is first relaxed without any external load to an equilibrium state.Then in each loading cycle, the strain-controlled method is used in the loading process by shifting certain amount of boundary atoms at both ends of the tube along the axis at small displacements.The loading is continued until the structural instability occurs for the first time.The stress-controlled method, on the other hand, is used in the unloading process by releasing in each step a small fraction of the final load exerted on the boundary atoms at the end of loading process.This scheme is implemented to examine exactly how the structural deformation will result upon complete removal of axial loads.The strain-controlled loading will guarantee a stabilized structure configuration which is obtained even after local/global buckling occurs, while stress-controlled unloading will guarantee the monotonic decrease of the loads on each boundary atom, thus avoiding any unwanted structural deformation which otherwise could occur if strain-controlled unloading was employed.
To implement the two loading/unloading methods in AFEM, load steps need to be specified.For the strainedcontrolled loading, a step of 1.39 × 10 −3 nm (1/100 of R (e) ) is selected.In each load step, the prescribed displacement boundary conditions are entered into the displacement vector u properly, and then (1) is solved iteratively until the force vector P becomes sufficiently small.For the stresscontrolled unloading, a step of 1/500 of the final load exerted on the boundary atoms at the end of loading process is selected.These prescribed force boundary conditions are entered into the force vector P properly, followed then by the iterative solving of (1).from cyclic tension and cyclic compression simulations are shown in Figures 1 and 2, respectively.In all the simulations, the cyclic loading is not stopped until the stress path in the final loading cycle starts to repeat that shown in the last loading cycle.In other words, exactly the same cyclic behavior will result in the subsequent loading cycles if the simulation continues.It is interesting to note that the successful implementation of the above criterion in all cyclic loading simulations constitutes in itself a solid proof that CNT with a flawless hexagonal lattice structure exhibits remarkable resistance to cyclic loads.As a matter of fact, except that the zigzag tube in cyclic compression has experienced a progressive loss of its yielding strength and strain during the first nine loading cycles and thereafter reached a steady cyclic behavior, all the other tubes reach their stable cyclic behavior with constant yielding strength and strain in the second loading cycle, as clearly shown in Figures 1 and 2. This result is consistent with the extreme fatigue resistance of vertically aligned CNT arrays that will recover to their near original length in cyclic compressive loading [22,23].

Viscoelastic
Behavior.Typical viscoelastic behavior characteristic of soft tissues, including hysteresis, preconditioning, nonlinear elasticity, and stress relaxation, has been reported from a study on fatigue resistance of blocks of vertically aligned multiwalled CNTs subject to cyclic compressive loads [22].No signs of fatigue failure have been observed on this material at high strain amplitudes (about 15%) up to half a million cycles.In the current study, similar behavior has been observed in simulations of a single CNT subject to both cyclic compressive and tensile loads.Although a limited number of loading cycles were carried out in the simulations due to high computational loads, well-defined cyclic stressstrain behavior of CNT has already been obtained.However, mechanisms of the observed viscoelastic behavior different from those given in [22] need to be identified.In other words, they should be explained by the intrinsic geometric and mechanical characteristics of the tube structure itself, rather than by the collective behavior (e.g., long-range van der Waals interaction and thermal energy transfer) of nanotube blocks.Further details are given here in after.Nonlinear elastic stress-strain response occurs and holds for all CNT samples loaded to structural yielding in the current study, with the axial strain almost fully recovered at the end of each unloading.For all the cases, the loading and unloading curves are overall parallel to each other but take distinct paths which form closed hysteresis loops.A 3% to 6% residual strain was seen to be left in the tubes upon the complete release of tensile loads for the first time, while for cyclic compression, a 0% to 2% negative (i.e., tensile) residual strain was typical for the tubes at the end of first unloading The residual strain, whether positive or negative, is a result of the adaptation of the tube to the virgin loading, and appears constant during subsequent loading cycles.For the zigzag tube, an apparent enhancement of the yielding strength and strain resulting in an expanded hysteresis loop in the second loading cycle was observed in both cyclic compression and tension.This phenomenon, which however did not occur in other simulations, may arise from the lattice structure reorganization resulting from the specified displacements of selected boundary atoms of the zigzag tube with full symmetry.Hysteresis loops formed by distinct loading and unloading paths were observed for all individual CNTs, as shown in Figures 1 and 2. The strain energy stored during the loading process is greater than the strain energy dissipated during the unloading process, resulting in a net amount of strain energy equal to the area of the hysteresis loop which is permanently accommodated by the residual strain left at the end of unloading.The residual strain alters the intrinsic structure of the tube and acts as an energy cushion of the CNT to external loads.It suggests an alternative state of minimum energy of the tube after its virgin loading cycle In a real world, the residual strain may not be observed since pristine CNTs without any load/deformation history may not exist.When thermodynamic effects are considered, hysteresis could be the result of thermal energy transfer between the tube and its environment.
Preconditioning, which refers to the progressive loss of the yielding strength of soft tissue materials [28,29] during cyclic loading, was only observed in the zigzag tube subject to cyclic compressive loads, as shown in Figure 2(a).The stress-strain hysteresis loop started to shrink from the third loading cycle, and was eventually reduced in the ninth loading cycle to a final-sized loop that was repeated in the later loading cycles.The final yielding strength and strain (including residual tensile strain), as compared to those in the first two loading cycles, has decreased by about 28% and 50%, respectively.Apparently, preconditioning cannot be explained by tube-tube interactions as in [22], which does not exist in the current study.Instead, inherent structural changes caused by cyclic loading need to be examined for this phenomenon.
In its virgin loading, the zigzag tube was loaded to a point where global buckling was pending.To this point, partial postpeak strain softening associated with local buckling  has occurred near the ends of the tube.During the first unloading, although the local buckling deformation was all recovered, the memory of this structural instability was revoked when the tube was reloaded to a certain critical strain level.As a result, failure of CNT in the second and subsequent loading cycles was caused by the localized collapse where local buckling occurred in the virgin loading.
Global buckling as an alternative of structural yielding has been eliminated in this situation.Although the stress path went back to the same point at the end of all unloading processes, the residual pinching deformation at the radial direction started to accumulate upon the completion of second unloading.These accumulated permanent deformations are the major cause of the progressive loss of tube strength, that is, preconditioning of CNT.
A better perspective on this mechanism can be achieved if the cyclic behavior of the zigzag tube is investigated after the global buckling occurs in the virgin loading, as shown in Figure 3(a).The corresponding cyclic stress-strain curve is shown in Figure 3(b).The tube experienced a sudden and large stress loss after the global buckling, which resulted in large irrecoverable pinching deformation at the end of unloading, as shown in Figure 3(c).This permanent structural deformation leads directly to the immediate significant loss of the yielding strength of CNT in the subsequent loading cycles Preconditioning was not observed in other tubes subject to cyclic compressive loads because no plastic morphological change was left at the end of unloading.This could be related to the absence of local buckling deformation associated with the strain-softening behavior in these tubes right before the unloading is started.

Effects of Chiral Symmetry .
The values of final stationary cyclic strength and strain of all the tubes simulated in this study are listed in Table 2. Their variations with chiral angles are shown in Figure 4.It should be noted that the final yielding strain is measured from the starting point of each enclosed hysteresis loop where the residual strain due to the virgin loading cycle has been properly accounted for (i.e., subtracted if positive or added if negative).This is because the tube structure with the residual strain represents the minimum energy configuration after the initial loading cycle, and could be more close to the real equilibrium structure of CNT under a stress-free condition.
It can be seen from Figure 4(a) that, for both cyclic compression and tension, the final yielding strain increases with increasing chiral angle, with the zigzag type being the weakest and armchair type the strongest.This result is consistent with the observation by Dumitrica et al. [8] and can be explained by the orientation of the C-C bond sustaining the largest strain during axial loading.For a zigzag tube, the bond parallel to the tube axis sustains almost all of the axial strain, while, for an armchair tube, the axial strain is shared by two equivalent bonds oriented 30 • from the axis.This chiral dependency does not change with the number of loading cycles and the axial loading manner, although the yielding strains in cyclic tension are overall much greater than the corresponding values in cyclic compression.
Similar effects of chiral symmetry on the final yielding strength are evident in Figure 4(b), with slight drops occurring at some intermediate ranges of chiral angles.Again, this result is determined by the particular orientation of hexagonal network associated with a certain chiral symmetry.The effect is more pronounced when the chiral angle approaches 30 • , especially in cyclic tension, as evidenced by a sharp increase of the final yielding strength in Figure 4(b).This behavior can be well explained by the stress-strain curves of the two tubes with chiral angles equal to 19.1 • and 30 • (armchair), which are composed of two distinct modulus regions, as shown in Figures 1(c) and 1(d).A lower constant modulus was exhibited by the two tubes up to around 30% strain (a little over 30% for armchair, and a little below 30% for chiral angle of 19.1 • ), followed by a significantly enhanced modulus regime.Similar modulus behavior has been observed by Suhr et al. [22] on CNT blocks and by Nieh [30] on open-cell aluminum foams.In these studies, collapse and densification of individual CNTs or foams were considered as the primary mechanism of the above behavior.This explanation, apparently, is not appropriate for the observed behavior of a single CNT in the current study.Once again, the intrinsic geometric features of the hexagonal lattice structure were carefully studied to give reasonable explanations.It was found that the sharp increase of the modulus was closely related to the bond rotation induced during the cyclic tension.Figure 5(a) shows the variation of average bond rotation during the virgin tension for the armchair and zigzag tube.The rotation directions of bonds and their deformed configurations are also shown in Figure 5(a).It is quite clear that there is a distinct change of rate of bond rotation occurring in the armchair tube at approximately the same strain level where change of modulus behavior was observed, while no such change was found in the zigzag tube.A similar but attenuated rate change of bond rotation was found in the chiral 19.1 • tube, but the effect was reduced to a negligible level for the chiral angle of 10.9 • tube.Besides bond rotation, the variation of average bond stretching with applied strain is also plotted, as shown in Figure 5(b).It is seen that the average bond stretch for a single bond in an armchair tube is much greater than that in a zigzag tube, especially at large strain level.This observation provides further explanation for the excessive yielding strain at the high chiral angle side (Figure 4(b)).Based on the above, it is concluded that the difficulty level of bond rotating and stretching associated with the varying geometries of hexagonal lattice gives full explanation of the effects of chiral symmetry on the stress, strain, and modulus behavior of CNT.

Conclusions
This paper has presented the results of a theoretical study of a single CNT subject to cyclic tensile and compressive loads using molecular simulations under quasistatic conditions.Parametric studies have been conducted to investigate the effects of chiral symmetry on the cyclic loading behavior of CNT.Distinct viscoelastic characteristic of soft tissue materials, including hysteresis, preconditioning, nonlinear elasticity, and large strain, has been observed for an individual CNT during cyclic loading condition.No signs of cyclic failure were observed in the simulations although the tubes were loaded to the maximum strain level they could sustain before yielding in each loading cycle.Extreme cyclic loading resistance of a single CNT is evident.These results agree well with the laboratory observations of cyclic compression of CNT blocks composed of vertically aligned CNT arrays or open-cell foams reported in previous studies.In addition, chiral symmetry was found to have a pronounced effect on the cyclic yielding strain and strength of CNT, with the zigzag type being the weakest and armchair being the strongest.
Mechanisms underlying the above observed behavior of a single CNT, which however differ from those given for the collective behavior of bulk materials, have been reanalyzed and reinterpreted.It was found that the intrinsic geometric and mechanical characteristics of the tube structure itself play the dominant role in the cyclic strength, strain, and modulus behavior.Specifically, accumulated residual morphological deformation due to local/global buckling was the major factor responsible for the preconditioning of the CNT while the amount of bond rotation and stretching experienced during loading/unloading essentially determines the strain and modulus behavior of the CNT with a certain chiral angle.Overall, the current theoretical study has provided clear evidences of the exceptional fatigue resistance of CNT that make this class of material highly potential in the building of synthetic biomaterials, electromagnetic devices, or polymer composites.

Call for Papers
The idea of magnetoelectricity has been around for decades, and several materials and concepts have been studied; however, recent discoveries of compounds in which the magnetoelectric (ME) coupling is strong have sparked a renaissance.ME symmetry allows magnetization (resp., polarization) to be induced by the application of an electric (resp., magnetic) field.The initial impetus driving ME materials was the prospect of employing a single field excitation, for example, voltage gradient, to tune both polarization and magnetization characteristics.Significant progress in ME nanostructures has been made in recent years, mainly driven by materialist use in new class of multifunctional device applications, for example, electric field-controlled magnetic data storage, antenna miniaturization using passive substrates, and field control for spintronics applications.Manipulation of electric polarizations by external magnetic fields (and vice versa) has been demonstrated in some of single-phase multiferroic materials, for example, BiFeO 3 .In this case, the intrinsic ME effect arises from the coupling between electric and magnetic order parameters.Intrinsic ME coupling occurs in compounds in which time-reversal and space-inversion symmetries are absent.A significant ME effect is also observed in the case of two-phase (thin films) composite multiferroics made by combining ferroelectric and ferromagnetic substances together, for example, BaTiO 3 /CoFe 2 O 4 , and lead-zirconate-titanate (PZT)/ferrite (or Terfenol-D).In these previous studies, it is commonly believed that elastic coupling amongst the ferroelectric and magnetic order parameters is at the origin of the ME effect, that is, an electric field E induces strain in the ferroelectric; this strain is passed on to the ferromagnet, where it causes a net magnetization M in the sample.A basic understanding of the distributions of electric, magnetic, and internal stress fields within the material allows one to manipulate the relative stability of domain structures.One would expect that the ME effect is large if the coupling at the interface is large; thus, composites with large surface area (such as multilayer systems and nanogranular mixtures) and strongly ferroelastic components are especially effective.The theoretical understanding of the microscopic mechanism that causes the strong coupling between the electric and magnetic degrees of freedom is quite limited.We therefore anticipate that ME nanostructures may open up a whole new dimension for a wide range of materials having cost effective, small size planar, and broadband multifunctional device applications.
Papers are solicited in, but not limited to, the following areas: • Physical mechanisms that control the ME property in layered and granular nanostructures

Figure 3 :
Figure 3: Behavior of the zigzag tube under cyclic compressive loads after global buckling occurs: (a) global buckling; (b) stress-strain curve; (c) structure morphology at the end of the first unloading.

Figure 4 :Figure 5 :
Figure 4: Effects of chiral symmetry on the final yielding strength and strain: (a) final yielding strain versus chiral angle; (b) final yielding strength versus chiral angle.

Table 1 :
CNT samples used in the parametric study of MD simulations.

Table 2 :
Steady state cyclic strength and strain of all CNTs obtained from simulations.
• Multiferroic composites • Multifunctional applications of nanostructures • Giant ME effects • Magnetoelectric interactions at electromechanical resonance • Magnetization and polarization tuning by external fields • Coexistence of magnetic and ferroelectric orders • Static and dynamic magnetoelasticity • Crystallographic engineering by assembling ME building blocks • ME effect in superconductors Before submission authors should carefully read over the journal's Author Guidelines, which are located at http://www .hindawi.com/journals/jnm/guidelines.html.Prospective authors should submit an electronic copy of their complete manuscript through the journal Manuscript Tracking System at http://mts.hindawi.com/according to the following timetable: Department of Physics, Université de Bretagne Occidentale, Université Européenne de Bretagne, Brest, France; brosseau@univ-brest.fr