Changing and predicting the frequency of double wall carbon nanotubes oscillator

Multi-wall carbon nanotubes have been prepared and intensively studied, as promising candidates for nanoscale molecular motors, switches, bearings, springs and oscillators.^1–6 In 2000, Cumings and Zettl¹ reported that, when a inner tube was pulled out from an outer tube and released, it would quickly turn back, driven by the Van der Waals Force. Then proposed by Q.S. Zheng, Q. Jiang et al. that, a double wall carbon nanotube (DWCNT) with both open ends could work as oscillators with ultra-high frequency up to several GHz.^7–9 Since then, many studies were established to explore the preparation, the mechanics and dynamics properties of DWCNT oscillators.^3–5,9–20 For example, W.L. Guo et al. studied the energy dissipation during the oscillation,¹⁰ L.H. Wong et al. found the oscillation energy leakage can be reduced by placing a single defect on the outer tube,¹³ and R.F. Zhang et al. detected the superlubricity²¹ between the ultra-long outer and inner tubes.²² 

Although many experimental and simulation works have been performed in this area, most of them just focus on one or a few kinds of nanotubes, which have fixed size and chirality.^{4,9–12,15–18} How the oscillation frequency changes exactly with the DWCNT’s shape has not been well discussed. On the other hand, theoretical work has given the expression of oscillation period based on continuum model or discrete model,^7,8,23 but the formula usually contains complex integral, which can’t be easily calculated.

In this article, we perform molecular dynamics (MD) simulation for a series of DWCNT oscillators with different tube length, different radius and different chirality. The tendency for the frequency change is carefully studied, and the mechanism behind is revealed.

MD simulation has provided many details for the oscillation process of the DWCNTs.^{3,9–12,15–18} In our MD simulation, carbon nanotubes with different shape are built, and the COMPASS forcefield²⁴ is selected for describing the inter-atomic interactions. It’s a universal forcefield provided by Materials Studio and LAMMPS. As shown by C.Y. Guo et al., it can best reproduce the phonon frequencies in the room temperature range, comparing with Tersoff, CHARMM, CVFF and PCFF forcefields.²⁵ It can also reproduce the surface energy of graphene, which consists well with the earlier experimental studies.¹⁹ In this forcefield, bond stretch, bond angle bend, inversion and tensional rotation terms are considered, and a Lennard-Jones type potential is used for determining the non-bond Van der Waals interactions.²⁴ 

In MD simulation, the selection of cut-off distance is very important, which may change the oscillating behaviors. Therefore, a sensitivity test is needed. Five identical (5,0),(8,8) DWCNTs are built. The lengths of the outer and inner tubes are 49.2Å and 25.6Å, respectively. Before the sensitivity test start, the positions of the outer tubes are fixed. Then the inner tubes are released with the initial amplitudes, 47.9Å. Fig. 1(a) shows the oscillating curves of these DWCNTs with different cut-off distances, ranging from 5.5Å to 13.5Å. From Fig. 1(a) we can see the cyan (13.5Å), blue (11.5Å) and green (9.5Å) curves overlap with each other, while the red (7.5Å) and black (5.5Å) curves move far away. By measuring the distance between the valleys, we find the oscillation periods are very similar, if the cut-off distance is larger than 9.5Å. The periodic error is smaller than 1.5%. Similar things happen when the length of the outer tube increases to 147.5Å, as shown in Fig. 1(b). Considering the calculation results above, we set the cut-off distance to 9.5Å, in order to balance the calculation efficiency and accuracy. This value consists with the former studies.¹⁹ 

The initial temperature (T₀) is another important factor which needs to be carefully chosen. Here we also carry out a test for six identical (9,0),(18,0) oscillators working at different T₀. T₀ ranges from 1K to 500K. The total simulation time is 500ps. As shown in Fig. 2, The oscillating curves look the same at the first few time steps. But as time goes on, the oscillators working at higher T₀ couldn’t keep their initial amplitudes. Especially, when T₀ reaches 500K, the amplitude decreases by 25.5Å after 500ps simulation (see the blue curve in Fig. 2). And the oscillation frequencies become different as well. This phenomenon can be attributed to the faster energy dissipation at higher T₀. As mentioned by the former studies, the higher the T₀ is set, the faster the energy dissipation will become.¹⁵ As a result, the amplitude decreases. In this article, our aim is to study the oscillating behavior with a given shape and amplitude. So T₀ is set to 1K, which may avoid the influence from the amplitude change.

There are other factors which need to be determined before MD simulation start. Since it’s a nonequilibrium system, there are no reliable simulation tools for the temperature control until now.¹⁰ We just simply use the NVE ensemble, to take a first step for understanding and predicting the oscillation behavior with these assumptions. The time step is set to 1fs, and the total simulation time is 50ps.

Firstly, we want to investigate the effect when the length of the outer tube changes separately. A zigzag (5,0) inner tube with six repeated units (RU) is built, whose length is 25.6Å (For zigzag tube, 1RU= 4.26Å, while for armchair tube, 1RU = 2.46Å), as shown in Fig. 3(a). Then five (8,8) outer tubes are built with five different lengths, ranging from 20 RU (49.2Å) to 60 RU (147.5Å). The radius difference between the inner and outer tube is 3.47Å. According to former studies, this value (∼3.4Å) is appropriate for providing stable oscillation.^9,23

After structure optimization, we combine the outer and inner tubes together along z axis, and fix the position of the outer tubes.⁹ Then the inner tubes are pulled out, until their mass centers reach the edge of the outer tube. So the initial extrusion is half the length of the inner tube, as shown in Fig. 3(b).

After MD simulation, the trajectory of the inner tubes along the z-axis can be tracked. In Fig. 3(c) we plot the displacement for each inner tube. Different curves refer to different outer tube length. From Fig. 3(c) we can clearly see that all these curves arrange regularly, implying a stable oscillation. The amplitudes are different because of the increasing outer tube length. The red curve stands for the inner tube in the shortest outer tube (20 RU). As seen, it moves forward and turns back for three complete cycles during the 50.0 ps simulation time. By checking the time difference between the peaks, the oscillation period can be determined to be 16.0 ps. Correspondingly, the oscillation frequency is 62.5 GHz. This frequency is much faster than the CPU dominant frequency nowadays. When the outer tube becomes longer, the oscillation frequency slows down. As shown in Fig. 3(d), the oscillation period changes from 16.0 ps to 38.0 ps, and the frequency changes from 62.5 GHz to 26.3 GHz, respectively. In order to show the tendency in the period change, we plot the best fitting line through the data points (red line in Fig. 3(d)). We find that, all the data points are very close to the red curve, implying a good linear feature.

To better understand the linear behavior, firstly, we need to deduce the oscillation period expression using the maximum mean force approximation(effective force approximation), as suggested by the former researchers.²³ In this approximation, the Van der Waals force F act on the inner tube can be divided into two parts. When the inner tube moves out from the outer tube, the force F reaches its maximum value, F_m; when the inner tube completely moves into the outer tube, F suddenly decreases to zero. This approximation is based on the fact that the Van der Waals interaction is mainly cut off upon 9.5Å. As shown in Fig. 4 (a), when the inner tube immerges deeply into the outer tube, the Van der Waals range remains to be constant. So the inner tube can not feel any difference during its movement. The resultant force F is zero, and the inner tube’s velocity v₀ remains unchanged (shown in Fig. 4(a)). This period can be considered as the uniform motion period. On the other hand, when the inner tube moves out from the outer tube, the area of the outer tube providing the effective pull-in force is fixed (shown in Fig. 4(b)). So the resultant force F=F_m is also constant, where |F_m| is only determined by the radius of the outer and inner tubes.¹ This period can be considered as the uniformly accelerated motion period.

We need to strengthen that, this maximum main force approximation above neglect the fluctuation of the real force acting on the inner tube, which caused by the discrete distribution of the mass. The approximation also neglect the force change when the inner tube is about to move completely into the inner tube, during which period the force shift from its maximum value to zero gradually due to the edge energy barrier. However, our calculation result shows that, this simplified theoretical model consists well with our MD simulation results.

Assuming the length difference between the outer tube L_out and the inner tube L_in is ΔL, the oscillation period T can be calculated by the following formula:

In the formula above, T consists with two terms: the uniformly accelerated motion term, and the uniform motion term, as mentioned before. In Eq. (1), L_ex corresponds to the initial extrusion of the inner tube (shown in Fig. 4 (c)), and a corresponds to the inner tube acceleration. a can be calculated by:

where F_m is the maximum main pull-in force, and m_in is the total mass of the inner tube. Meanwhile, v₀ in the second term of Eq. (1) can be calculated by:

where v₀ corresponds to the uniformly motion speed inside the outer tube.

After deduce the oscillation period expression, we now try to explain why the oscillation period T increases linearly with the outer tube length L_out. According to Eq. (1), T is determined by a, L_ex, v₀, and ΔL. We would like to discuss them separately. As shown in Fig. 3(b), when L_out becomes longer, the radius difference and the inner tube length L_in remain unchanged. So F_m and m_in in Eq. (2) are constant, and the acceleration a is constant, too. On the other hand, since L_ex never change, the v₀ in Eq. (3) keeps constant. As a result, the oscillation period T depends only upon ΔL. In our case, L_out changes linearly, while L_in doesn’t change. Therefore ΔL increases, which leads to the linearly increasing of the oscillation period, as shown in Fig. 3(d).

In this section, we want to study the case that the inner tube and the outer tube have the same length (L_in=L_out). Five (5,0)@(8,8) DWCNTs with different length are chosen. The length of the armchair-type outer tubes ranges from 24.6 Å(10RU) to 54.12 Å(22RU). The initial condition for each oscillator is similar to Section III A, that the initial extrusion L_ex is half of the inner tube length L_in.

After MD simulation, we plot the trajectory of each inner tube in Fig. 5(a). The red, green, blue, cyan and purple curves correspond to the tube length of 10RU, 13RU, 16RU, 19RU and 22RU. The corresponding oscillation period T is shown in Fig. 5(b). We can see T shifts from 10.6ps to 23.2ps. The red line in Fig. 5(b) is the best fitting line, and the fitting expression is written bottom-right. In this expression, L_s is in RU, and T is in ps. In most cases, the length of DWCNT oscillator is much longer than 15Å, so the first term -0.093 in the expression can be negligible. Approximately, T is proportional to L_in.

This phenomenon above can also be interpreted by our effective force model. In this case, L_in=L_out, ΔL=0, and L_ex=L_in/2. Substitute ΔL, L_ex to Eq. (1), we obtain:

In the formula above, a = F_m/m_in. Similar to Section III A, F_m is constant because of the unchanged radius. However, m_in is no longer constant because L_in changes. m_in is proportional to L_in and can be written as:

where σ is the areal density of the nanotube, and r_in is the radius of the inner tube. Substitute m_in into a= F_m/m_i, and substitute a into Eq. (4), we get the final expression:

In the equation above, F_m, r_in and σ are constant. So the terms inside the radical sign remain constant, too. As a result, the phenomenon that the oscillation period T increases proportionally to the inner tube length L_in can be well explained. Our simulation result corresponds well with the earlier studies.¹⁴ 

In this section, we consider the effect induced by the radius change. We build five inner tubes with different radius, ranging from 4.1 Å to 9.5 Å. Then we choose five outer tubes, whose radius is larger than the inner tubes by 3.39 Å. Similarly to Section III B, all the outer and inner tubes have the same length, 24.60 Å (10RU for armchair tubes). The initial extrusion L_ex is half of L_in, and then MD simulation starts.

Fig. 6(a) shows the displacement of these inner tubes. As seen, the amplitudes of these tubes are similar, while the oscillation period T is different. The larger the radius r_in becomes, the longer the period T becomes. T changes with r_in is shown in Fig. 6(b). We can see T increase from 10.5 ps to 12.2 ps, which is not as obvious as in Section III A and Section III B. Moreover, the period change is neither linear nor proportional versus the inner tube radius.

Noticing L_in=L_out=constant, L_ex=L_in/2, the period change can only be attributed to the change of F_m, induced by the tube radius change. In Section III A and Section III B, we have shown that the effective force approximation is a good mathematical model for understanding the oscillation behavior of the DWCNTs. Here we also use this model, to estimate the effective force F_m acting on each inner tube with different radius r_in. From Eq. (6), F_m can be written as:

Assuming $σ$ = 7.606 × 10⁻⁷ kg/m³, and substituting L_in, r_in and T into Eq. (7), the value of F_m can be obtained. Then the F_m change versus r_in is plotted in Fig. 6(c). The red line is the best fitting line. Surprisingly, we find that F_m increases linearly with the radius r_in. The expression can be written as:

where F_m is in nN, and r_in is in nm. This linear feature is interesting, implying that we can estimate the oscillation period and frequency of different DWCNT oscillators using a simple formula. The formula will be deduced in Section III E.

From Eq. (8) we can also obtain the surface energy density γ. Referring to the former research,¹⁹ γ can be estimate by γ =F_m/4πr_in. If r_in is large enough, γ will converged at 0.18 N/m². This value corresponds well with the former studies for nanotubes and graphene, which give a value ranging from 0.12 to 0.18 N/m².^19,22,26

In this section, we want to study the chirality effects on the oscillation frequency. Two outer tubes (14,0) (8,8) with the similar radius and length (∼25.6Å) are built, while two inner tubes (3,3) (5,0) are built as well. In Table I we plot the frequency of these oscillators. As shown, all these frequencies are very similar. The frequency difference between these oscillators is smaller than 0.7%. That means, the chirality takes tiny influence to the frequency of the DWCNT oscillator, similarly to Xiao et al.’s studies.¹⁵ However, as discussed by W.L. Guo. et al.,¹⁰ the similar chirality between the outer and inner tubes will lead to faster energy dissipation. As a result, we’d better choose different chirality for the inner and the outer tubes.

From the calculation above, we can conclude that the effective force approximation can provide good description for the behavior of different kinds of DWCNT oscillators. Also we find the effective force estimate by our model seems to be proportional to the radius of the inner tube. By substituting Eq. (8), Eq. (5) into Eq. (2), and substituting Eq. (2), Eq. (3) into Eq. (1), we can finally get the formula below:

where the inner tube radius r_in is in nm, T is in s, and L_in, L_ex, ΔL is in m. This formula can be used for predicting the oscillation period of different DWCNT oscillators, if the distance between the inner and the outer tube is ∼ 3.4Å.

In order to check the precision of the period formula obtained above, we choose different kinds of DWCNTs with different length, different amplitude different radius and different chirality for comparison. As shown in Table II, all the simulation results fit well with our theoretical predictions, with an error smaller than 5%. That means, our simplified theoretical model has good accuracy.

However, we need to mention that, in order to avoid the error induced by the Van der Waals range, L_out, L_in should not be shorter than 20 Å, and L_ex should not be shorter than 10 Å. Fortunately, nearly all the DWCNT oscillators studied by the former researchers meet this requirement.

In this article, the oscillation frequency changes with the DWCNT’s shape is carefully studied using MD simulation. We find that, the oscillation period can be easily influenced by the tube length, the radius, and the extrusion distance. By employing the maximum main force approximation, the linear behavior of the period change can be well understood. Moreover, we find the effective maximum force increases linearly with the inner tube radius. After fitting the data point, a universal formula is deduced, which can predict the oscillation period with an error smaller than 5%. This work may provide theoretical guidance for changing and predicting the frequency of future ultra-fast DWCNT oscillators.

We would like to thank Prof. X.H. Yan, Prof. D.N. Shi, Dr. R. Ma, Dr. K. Zhou and Dr. K.X. Zhang for helpful discussions. This work was supported by National Natural Science Foundation of China (NSFC11574155 and NSFC11605091), Pre-research Project of NUIST (2014x034) and the Startup Foundation for Introducing Talent of NUIST (2012x062).