Chapter 3 Epoch Making Simulation

Large Scale Simulations for Carbon Nanotubes
Project Representative

Syogo Tejima

Research Organization for Information Science & Technology

Authors

Syogo Tejima 1, Satoshi Nakamura 1, Kohji Makino 1, Yutaka Maniwa 3,
Yoshiyuki Miyamoto 4, Yoshikazu Fujisawa 2 and Hisashi Nakamura 1
1 Research Organization for Information Science & Technology (RIST)
2 Honda R&D Co., Ltd
3 Tokyo Metropolitan University
4 Fundamental and Environmental Research Laboratories, NEC Corporation
CARBON NANOTUBE RESARCH GROUPE

Morinobu Endo, Eiji Osawa, Atushi Oshiyama, Yasumasa Kanada, Susumu Saito,
Riichiro Saito, Hisanori Shinohara, David Tomanek, Tsuneo Hirano, Shigeo Maruyama,
Kazuyuki Watanabe, Takahisa Ohno, Yutaka Maniwa, Yoshikazu Fujisawa,
Yoshiyuki Miyamoto, Hisashi Nakamura

Nano carbon materials as nanotubes (CNTs) and fullerenes in nanotechnology have a lot of potential for industrial applications. On the efforts of developing applications, it has been recognized that computational simulations are powerful and efficient tools to find and create new materials from nano scale. Aiming at realistic simulations for nonmaterial, we have developed a large-scale computation technique utilizing tight binding molecular dynamic method, ab initio density functional theory (DFT), and time-dependent DFT method. We have studied various physical properties of nano-carbon and applications e.g.,(1) Growth mechanism of carbon nanotube, (2) Characteristics of Mackay structure, (3) Hole-doped diamond superconductor,(4) Quantum Transport for molecular wire, (5) Irradiation of strong optical field on nano-carbons, (6) One-dimensional magnetic oxygen ordering in CNT, (7) Order-N first principle method for large scale system.
Along these works, we have realized that the Earth Simulator is a very powerful tool for large-scale material simulations.
Keywords: Large scale simulation, TB theory, ab initio theory, DFT, Carbon Nanotube, Fullerenes, GSW, CNT growth,
Hydrogen storage, Hole-doped Diamond superconductor

1. INTRODUCTION

equation linearized with respect to the anomalous Green's
function applied to the prediction of high Tc superconductivity for a hole-doped diamond. We have carried out somesubjects in this work, which are described in the next section. There are three primary objectives with this work: (1)design of innovative nonmaterial with certain desired properties; (2) obtaining fundamental properties in nano-scalematter, and (3) nano-applications. Our purpose is to give theclear explanation of properties and phenomena of nano-scale
events and deduce guiding principle to design new materials
for applications from nanostructures using super-computers.

Carbon materials have been expected to make a breakthrough in material science and nanotechnology. A lot of
potential applications of nanotubes and fullerenes e.g., electronic field emitter and electronic devices have attracted scientific community. In investigating and utilizing their material properties, numerical simulation by means of supercomputer has turned out to be a very efficient tool. A recent
development in nanotechnology has required a more efficient supercomputing capable of a large-scale simulation of
up to ten thousand atoms. Aiming at large-scale simulations
utilizing the Earth Simulator, we have developed computational package based on ab initio DFT theory and parameterized tight-binding method. The TB code we have developed
is shown to be suitable for the very large systems eventhough the lack of symmetrical arrangement. The Eliashberg

2. Physical studies on nano materials
2.1 Growth Mechanism of Carbon nanotube
Controlling growth process of SWCNT in the presence ofa metallic catalyst is an essential requirement to provide high141

Annual Report of the Earth Simulator Center April 2008 - September 2008

quality CNT for electrical and thermal transport wire inapplications. There has been large progress in the experimental mass production of SWCNT. Nevertheless, thegrowth mechanism of high quality SWCNT is still poorly
understood. On the basis of such experimental results, we
presented a scenario for an acceptable SWCNT growth
process, which are constructed from four stages. In this scenario, the metal cluster plays two important roles in thegrowth process. The first is to act as a hydrocarbon decomposition catalyst to produce C atoms (stage 1). The second isto act as a solvent for the C atoms (stage 2, 3, 4).
Using a first principle molecular dynamics calculation,
we examined a decomposition of Methane molecules bycatalysis. The reason we used quantum rather than classicalmethods is that it is difficult for classical method to consider strong quantum effects of H atoms having small mass.
Stage 1: From a simulation at 1100 K, methane molecule isdissociated into C and H atoms on the surface of the catalysis as shown in Fig. 1. The isolated flying H atoms tend to
become hydrogen molecules by binding. On the contrary, C
atoms tend to remain on the surface of catalyst during thesimulation. We speculate that due to the lack of simulation
time, we could not observe the diffusion of C atoms entering the catalyst.
After decomposition of hydrocarbon atoms, C atoms for
the CNT growth are placed at the surface of the catalyst.
According to our assumption, H molecules do not play animportant role in growth of SWCNT; therefore we consider
a behavior of C atoms only. On the basis of such a reason,
we can introduce an empirical classical molecular dynamic

simulation to deal with C and Fe atoms except for H atoms.
From the classical molecular dynamics simulations basedon empirical potential, we obtained very interesting behaviorof C atoms in the presence of Fe cluster composed of 55
atoms. Stage 2: until containing about 25 atoms of C in cluster, C atom remain inside the cluster with diffusion. Stage 3:
Over about 25 atoms of C, C atom saturates inside the cluster
and precipitate on the surface of Fe cluster. Stage 4: The precipitated carbon atom creates a carbon cloud near the surface.
At last, caped carbon structure lifts up from Fe cluster. Suchcarbon structure is certainly tubular structure not carbon nanotube. The tubular structure possesses sp3 hybridization andtube hollow are filled with C atoms as shown in Fig. 2.
Simulation results depend on conditions like a hydrocarbon flow constraint, system temperature, and size of Fe cluster. From the detailed research for the effect of these conditions on the SWCNT growth process, we will find the best
condition which helps understanding of the growth of CNT.
2.2 Characteristics of the Mackay structure [1]We have focused on the Mackay structure that consists of
only carbon and expects useful characteristics (stable, very
rigid, magnetic body, and so on). Three production processes
from connecting carbon nanotubes to the Mackay structure
were found out by using only sequential GSW rearrangements in last year. However, only small size structure's energy was obtained, because the amount of calculation is veryhuge. Then, the each energy of larger size structures is calculated newly. Figure 3 shows the energy of each structure.
Although there is energy barrier in the production process,Fig. 1 Fore methane molecules are placed on the Fe cluster of (001) and (111) surface in a initial structure (left figure). After 54ps using first principle MD simulation, methane molecules are decomposed into H and C , and one hydrogen molecule appears (right figure).

142

Chapter 3 Epoch Making Simulation

the barrier is not so high. These results confirm that there is
the possibility of the realization of the Mackay structure.
Furthermore, we clarified that the Mackay structure can be
distinguished into two types: zigzag and armchair type. Thisdifference is the direction of six-membered ring (Fig. 4).
Zigzag type has been known, but armchair type has not.
Some researchers have analyzed smaller zigzag types withcomputer simulation. In this research, the characteristics oflarge size and/or armchair type Mackay structure and are
examined. However, in order to investigate characteristics

Fig. 2 After 10 ns using an empirical classical simulation, carbon tubular structure is created on the Fe catalyst. Green circle shows Fe
atom red isolated carbon atom, blue bonding carbon atom.

Fig. 3 Energy of each structure.

Fig. 4 Zigzag and armchair type of the Mackay structures.

143

Annual Report of the Earth Simulator Center April 2008 - September 2008of large size structure, a large scale computer is necessary.
Therefore, large structures have not been dealt with. First,
bulk modulus is investigated. The results are shown in
Fig. 5. The smallest structure is almost as rigid as diamondand is more rigid than large one. Moreover, there was no difference according to the types. Next, it has been guessed thatthe carbon atoms in the Mackay structure are arranged along
the faces of Schwarz's minimal p-surface (Fig. 4:upper-left).
However, the carbons are arranged along the faces of truncated octahedron (Fig. 6). The results are useful to analyze thecharacteristics of the Mackay structure. Finally, the electric is

calculated with PWscf code. E-dos of 3 Mackay structures
are illustrated in Fig. 7. As a consequence, the Mackay structures have semiconducting property.
2.3 Hole-doped diamond superconductor
We have been evaluating the critical transition temperature (Tc) of hole-doped diamond on the basis of the charge
fluctuation mechanism. The result of our calculation lastyear showed that pure diamond with 0.03 holes per atom
was possibly the superconductor with 75 (K). The Tc's
observed in experiments, however, are lower than those esti-

Fig. 5 The relationship between primitive cell size and bulk modulus.

Fig. 6 Stable shape of the Mackay structure.

144

Chapter 3 Epoch Making Simulation

Fig. 7 E-DOS of 48, 144 and 192 atoms of the Mackay structures.

mated by simulation. At the present, the maximum value of
Tc is 7 (K). This discrepancy seems to be attributed mainly
to the method of incorporating the bare Coulomb interaction
within the tight-binding approximation. In the simulationachieved last year, the spatial spread of atomic orbitals has
not been taken into account at all.
According to the consideration of this spatial spread, thebare Coulomb interaction V has to have the following two
limits:nents. As a function that interpolates between these two limits, we have adopted.

We show in Fig. 8 the value of Tc's against hole-concentration. The figure represents that for relatively high holeconcentration, Tc's are lower than those simulated last year.
In this way, the discrepancy between the simulation andexperimental results has been improved.
2.4 Quantum Transport for molecular wire
Application of carbon nanotubes as electronic devices had
been expected due to their unique properties. Devices suchas random access memory and radio receiver have alreadybeen developed up to a considerable stage. In order to fully
achieve control on the properties of carbon nanotube-based
devices, it is important to perform analyses of the carbon
nanotubes of which sizes are the same as those observed inexperiments.
To this aim, we have developed a extended tight-binding
molecular dynamics method based on conventional tightbinding molecular dynamics method. A conventional tightbinding molecular dynamics method neglects the following

where

and U at is intra-atomic bare Coulomb interaction. For thesake of simplicity, the value of U at was set at 8 (eV) for the
diagonal components and 0 (eV) for the off-diagonal compo-

Fig. 8 Tc's against hole-concentration. (a) Result simulated last year. (b) Result simulated this year.

145

Annual Report of the Earth Simulator Center April 2008 - September 2008

Fig. 9 Electrical current distribution in the nanotube attached to one-dimensional electrodes. The diameter of the arrow
indicates the magnitude of the current in log scale. (a) The line connecting two electrode-junctions is parallel tothe axis of the nanotube. (b) The line connecting two electrode-junctions is not parallel to the axis of the nanotube.

effects.: A change in long range Coulomb interaction due toa change in the electronic density distribution and polarization within carbon atoms. We have formulated, therefore, atight-binding molecular dynamics method incorporatingthese effects. Furthermore, to treat carbon atoms in the electro-magnetic field, this method was combined with themicroscopic Maxwell's equations.
As a test calculation, we show in Fig. 9 the electrical current distribution in the nanotube attached to one-dimensional
electrodes. A detailed investigation is left for future work.

time, the field varies with respect to time aswhere, w and t respectively correspond to 800 nm of lightwavelength and 50 fs of pulse width. The maximum strengthof field E0 has been set as 1.6 V/Å as arbitrary test and we
found with this condition the surface structure of graphite
shows structural change shown in Fig.10 as below.
We expect the current condition will be realized by femtosecond laser technique and also expect other structuralchange, like graphite-diamond transition, will also be
observed by modifying the pulse shape and wavelength oflaser shot.
The other application of TDDFT with time-varying external field is penetration of optical field inside semiconducting
carbon nanotubes. It is well known that semiconducting nanotubes shows fluorescence depending on their chiralities.
The dependence is owing to their unique high density ofstate of one-dimensionality. Recently, precise experiments
can realize optical transition with polarization of opticalfield perpendicular to tube axis which opposed to previous
understanding that nanotubes completely depolarize the optical field in such polarization.
Motivated by the recent experiments, we investigate how
depolarization is strong by monitoring induced dielectricfield at the center of semiconducting nanotubes under external dielectric field. We have tested two semiconducting nanotubes, (8,0) and (14,0) tubes. In both cases, we found thatsignificant enhancement of field tentatively occurs when frequencies of the applied dielectric field are close to the resonance of excitation energies of these nanotubes. In contrary,when the frequency is out-of-resonance, the dielectric field
inside nanotubes is significantly attenuated. Figure 11,
shows the current results. From top to bottom, blue lines
denote time-evolution of applied field while the red lines
denote the allied field plus induced field by nanotubes(called as total field). Note that only the top of Fig.11 shows
off-resonance condition. This finding raises further interestfor optical excitation of molecules being encapsulated inside

2.5 Application of time-dependent density functional theory
for irradiation of strong optical field on nano-carbons [2]In the last year, we developed our time-dependent density
functional theory (TDDFT) code being available for computing material response under time-varying strong external
field. We develop formalisms to compute 'work' done by
the external field and confirm that increase of total energy(including kinetic energy of ions) is the same as the work.
This means our numerical simulation goes under stablecondition.
This year, we apply the code for two subjects. The one isstructural change of graphite surface irradiated with pulse
laser. Recently, mono-atomic layer of graphene sheet has
attracts high attention since mono-layer graphene has been
expected to show ballistic transport property, mechanicaland thermal robustness. However, the massive production of
mono-layer graphene is still challenging and usual processes(scotch tape, CVD, solvent treatment, etc.) suffer from contaminant. We were searching a way to make structural transition of graphite surface into diamond phase by strong laser
pulse, but just by accident, we discover a condition of exfoliation of single-layer graphene sheet.
In this calculation, we use periodic boundary condition
consisting of 10-layer AB-stacked 2 × 2 unit cell of graphitewith vacuum region of 10Å. The optical field by laser pulse
was mimicked by fictitious dipole charge layers uniformly
expand in parallel direction to graphite layer at the middle ofthe vacuum region making sawtooth type potential, which isa standard way of investigating the field effect. But this146

Chapter 3 Epoch Making Simulation

Fig.10 Time evolution of graphite surface after irradiation of laser shot.

Fig.12 Amplitude of the magnetic moment of oxygen molecule inside(m,0) and (n,n) CNTs.

multi-helix structures with the increasing CNT diameter.
Then, the first-principal calculations using pwscf code wereconducted on the narrow CNTs with one-dimensional arrayof oxygen molecules. It was suggested that amplitude of themagnetic moment of oxygen molecule decreases withdecreasing the CNT diameter below ~7Å as shown in
Fig.12. The results will give important suggestion for the
low-dimensional magnetic material design.

Fig.11 Field enhancement in semiconducting nanotubes. Blue and redlines are time-evolution of applied and total electric field.

semiconducting nanotubes. By tuning the chiralty of nanotubes, nanotube excitation energy can be adjusted to that ofmolecule. If this is the case, we obtained higher quantum
efficiency of optical transition because of the enhancementof optical field inside nanotubes.

2.7 Large scale calculation for application using a sophisticated first principle order-N calculation method [3].
We succeeded in structural optimization by a first principle self-consistent calculation for the nano structure whichcontained about five thousands atoms. By this code the stability of three dimension pyramid structure is obtained
shown in Fig. 13. This sample contains about 23,000 atoms.
It is understood to be able to optimize stability and the structure to such a hugeness system.

2.6 One-dimensional magnetic Oxygen ordering in CNT
Oxygen molecule is an interesting functional moleculewith spin S = 1. Therefore, it would be possible to develop
novel nanoscale magnetic materials by encapsulation ofoxygen molecules inside the carbon nanotubes (CNTs).
Here, we examined the structures and electronic states ofoxygen-encapsulating CNTs using the computer simulation
technique. First, classical molecular dynamics calculations
clarified that oxygen molecules form one-dimensional to147

Annual Report of the Earth Simulator Center April 2008 - September 2008

3. SUMMARY
The large-scale simulations on nonmaterial have been carried out by ab initio density functional method and the parameterized tight-binding calculations. The optimized codes
showed that the computation on the Earth Simulator could
give an exceptional performance and enables us morechance for large-scale and realistic simulations. Our large
scale simulations can provide the nanotechnology industries
valuable information on novel nano material properties andon nano electrical designs for application.

Fig.13 Optimized nano structure in about 23000 atom system.

References[1] K. Makino, et al, RIST NEWS No.45.
[2] Y. Miyamoto, and J. Zhang, "Testing the numerical stability of time-dependent density functional simulations
using the Suzuki-Trotter formula", Phys Rev B77,165123 (2008).
[3] T. Miyazaki, et-all, "The energetics of hut-cluster selfassembly in Ge/Si(001) from linear-scaling DFT calculations", J. Phys. Soc. Jpn. in press.

148

Chapter 3 Epoch Making Simulation

1

1

1

3

4

2

1

1
2
3
4 NEC

CNT

1.
CNT

2.
H14
H16 H17 H18 H19

H15
1

2

3
TC
4
2
5
TDDFT
1

6
1

CNT

CNT
N

7
DNA

N
1

149