J. Phys. Chem. C 2007, 111, 13087-13091
13087
Diffusion Dynamics of the Li Ion on C60: A Direct Molecular Orbital-Molecular Dynamics Study Hiroto Tachikawa† DiVision of Materials Chemistry, Graduate School of Engineering, Hokkaido UniVersity, Sapporo 060-8628, Japan ReceiVed: March 28, 2007; In Final Form: June 29, 2007
The diffusion dynamics of the Li+ ion on fullerene (C60) have been investigated by means of the direct molecular orbital-molecular dynamics (MO-MD) method. The total energy and energy gradient on the full dimensional potential energy surface of the Li+C60 system were calculated at each time step in the trajectory calculation. The optimized structure, where the Li+ ion is located at the hexagonal site of C60, was used as an initial structure at time zero. Simulation temperatures were chosen in the range of 10-300 K. The dynamics calculations showed that the Li+ ion vibrates around the initial equilibrium point below 40 K, while the ion can move above 50 K. At low temperature (below 300 K), the diffusion coefficients for the Li+ ion on the C60 surface are larger than those of the graphite surface. The diffusion coefficients on both C60 and the graphite surface were almost equivalent at medium temperatures around 300 K. At higher temperatures (T > 300 K), the coefficients for the graphite surface were significantly larger than those of C60. On the basis of theoretical results, we designed an ion-switching molecular device composed of C60 and graphite sheet.
1. Introduction The chemistry of fullerenes consists of wide fields from basic science to technological application.1-5 Also, it has potential utilities from electronic materials to medical pharmacy.5 The variety of utilities is originated in the specific properties of C60 which can interact easily with several atoms, ions, and molecules inside and outside of C60. For example, the fullerene turns into superconducting material on doping by alkali atoms (K, Rb, and Cs) with the transition temperature (Tc) exceeding 30 K. Also, alkali atoms (Li, Na, and Mg) doped C60 play as semiconductors.6-9 To elucidate the electronic properties of doped fullerenes, a lot of experiments and theoretical works have been carried out. Grpmov et al. measured UV-vis, IR, and Raman spectra of LiC60 produced by low-energy ion implantation of C60.10 They found that new peaks appeared at 618 and 1100 nm after interaction of C60 with Li. These transitions were assigned to the electronic excitations within C60-. The IR spectrum of LiC60 was very close to that of C60-, indicating that the LiC60 complex exists as a charge-separated state (Li)+(C60)-. Zhao et al. calculated electronic states of Li+C60 using the density functional theory (DFT) method. They showed that the Li+ ion can bind in the hexagonal site of C60 with a binding energy of 41.7 kcal/ mol.11 Santos et al. suggested that the Li+ ion is bound at about 2.0 Å from the surface of C60.12 Hannongbua calculated potential energy curve for the C60-Li+ system. The Li+ ion is bound at 2.5 Å from the surface of C60. Pitarch-Ruiz et al. showed that Na is bound to the hexagonal site of C60 at the 2A1 state.13 For a dynamics feature, collisions of Li to C60 have been investigated by several groups. Ohno et al. calculated the collision process of Li+ to C60- with high collision energy regions (collision energy was ca. 5 eV). The Li+ ion is bound on the surface of C60- after the collision to C60.14 Thus, †
E-mail:
[email protected]. Fax: +81-11-706-7897.
information on the static properties of M-C60 (where M means alkali metal) and collision reaction of Li to C60 at higher energy regions have been accumulated. However, there is no information on the diffusion process of M inside and outside of C60. In the present study, diffusion dynamics of Li+ outside of C60 have been investigated by means of the direct molecular orbital-molecular dynamics (MO-MD) method. In this method, full dimensional potential energy surface is calculated quantum mechanically, and then trajectory is run by including all degrees of freedom. The electronic states are calculated and updated in each time step, so this method gives the information that cannot be obtained by the usual MD calculation. In previous papers,15,16 we investigated the diffusion dynamics of lithium ion and atom (Li+ and Li) on the graphite sheet using the direct MO-MD method developed by us. In this method, the potential energy surface and energy gradient of all atoms of the system are explicitly calculated quantum chemically. Therefore, different characteristics from the usual classical MD calculation are obtained. It was found that the Li+ ion diffuses preferentially along the node of the highest occupied molecular orbital (HOMO) of the graphite sheet. This implies that modification of the HOMO leads to control of the diffusion route of the lithium ion. Also, it was found that the diffusion of lithium atom is much different from that of Li+ ion. The Li atom moves together with spin density distribution generated on the graphite surface. Hence, the diffusion rate for Li atom is slower than that of the Li+ ion. The diffusion coefficients calculated were in good agreement with the experiments. In this work, the similar technique has been applied to the diffusion of Li+ on C60. 2. Methods of Calculation First, the fullerene (C60) was fully optimized at the AM1 level. For the interaction system (LiC60+), one Li+ ion was put on the hexagonal site of C60, and then the structures of LiC60+ were
10.1021/jp072451c CCC: $37.00 © 2007 American Chemical Society Published on Web 08/15/2007
13088 J. Phys. Chem. C, Vol. 111, No. 35, 2007
Tachikawa
fully optimized. It should be noted that the AM1 calculation represents reasonably the structure and electronic states of the lithium-graphite system: charge and activation energy for the diffusion are in good agreement with those of B3LYP/ LANL2MB as described in previous studies.15,16 Diffusion processes of the Li+ ion on the C60 surface were investigated by means of the direct MO-MD method.17-19 The total energy and energy gradient on the multidimensional potential energy surface were calculated at each time step at the AM1-MO level of theory, and then the classical equation of motion is full-dimensionally solved. Therefore, charges and electronic states of the Li ion and all carbon atoms are exactly treated within the level of theory at each time step. This point is much different from the usual classical MD calculation where the charges of all atoms and ion are constant during the diffusion. Hence, one can obtain details of the diffusion processes of lithium ion on amorphous carbon using the direct MO-MD method. We carried out the direct MO-MD calculations under constant temperature condition. Mean temperature of the system is defined by
T)
1 3kN
〈∑ 〉 miVi2
(1)
i
where N is number of atoms, Vi and mi are the velocity and mass of the ith atom, and k is Boltzmann’s constant. We choose temperatures in the range of 100-1100 K. The velocities of atoms at the starting point were adjusted to the selected temperature. In order to keep a constant temperature of the system, bath relaxation time (τ) was introduced in the calculation.20 We have chosen τ ) 0.01 ps in all trajectory calculations. The equations of motion for n atoms in the system are given by
dQj ∂Η ) dt ∂Pj ∂Pj ∂U ∂H ))∂t ∂Qj ∂Qj
(2)
where j ) 1 - 3N and H is the classical Hamiltonian. Qj, Pj, and U are Cartesian coordinates of the jth mode, conjugated momentum, and potential energy of the system. These equations were numerically solved by the Runge-Kutta method. No symmetry restriction was applied to the calculation of the gradients. The time step size was chosen by 0.2 fs, and a total of 10 000 steps were calculated for each dynamics calculation. The drift of the total energy is confirmed to be less than 1 × 10-3% throughout at all steps in the trajectory. The momentum of the center of mass and the angular momentum are assumed to zero. Static DFT calculation, B3LYP/LANL2DZ, B3LYP/6-31G(d), and semiempirical AM1 calculations for the potential energy curves (PECs) along a diffusion path on carbon cluster were carried out using the Gaussian 03 program package.21 3. Results A. Structures of C60-Li+. First, the structures of Li+C60 are fully optimized at the B3LYP/6-31G(d), B3LYP/LANL2DZ, and AM1 levels of theory. The optimized structure of the Li+ ion adsorbed on a hexagonal site of C60 (denoted by h-site) is illustrated in Figure 1, while the optimized geometrical parameters are given in Table 1. The bond distances between Li+ ion
Figure 1. Optimized structure of LiC60+ (hexagonal site) obtained at the AM1 calculation and geometrical parameters. The notations R1 and R2 mean the distances of the Li+ ion from carbons C1 and C2, respectively.
and carbon atoms (C1 and C2 atoms) are calculated to be R1 ) 2.463 and R2 ) 2.466 Å, respectively, at the B3LYP/LANL2DZ level, which are slightly different from each other. This difference is caused by a slight difference of the electronic states of the C1 and C2 atoms of C60. These carbons are not equivalent to each other. The ion is located at 1.999 Å above the hexagonal carbon plane. In the pentagonal site, the Li+ ion is slightly higher in the position than that of the hexagonal site. The Li+ ion can bind also to the pentagonal site (denoted by p-site). In this site, the Li+ ion is located at 2.054 Å above the pentagonal carbon plane of C60. The energy of the Li+ ion in the p-site is 0.97 kcal/mol lower in energy than that of the h-site. The charges of Li+ ion are calculated to be +0.69 (h-site) and +0.72 (p-site), indicating that about 30% of the electron is transferred from C60 to the Li+ ion. The optimized geometrical parameters and relative energies calculated at the AM1 level are given in Table 1. The AM1 calculation gives the similar results of the B3LYP/LANL2DZ level of theory. Potential energy curves along the diffusion path of Li+ between hexagonal and pentagonal sites are plotted in Figure 2. The calculations are carried out at the B3LYP/LANL2DZ and AM1 levels of theory. The optimized structures of Li+C60 in the hexagonal sites are used as the geometries of C60 and height of Li+ ion from the surface. The height of Li+ ion is fixed to the optimized distance in both calculations, and then the angle (φ) is changed from 8.0 to 73.0. Two minimum points corresponding to p- and h-sites and one maximum point corresponding to the saddle point between p- and h-sites are found along the diffusion path. The barrier heights at the saddle point relative to the hexagonal site are calculated to be 9.3 (B3LYP) and 9.8 kcal/mol (AM1). The saddle point is located above the center of C-C bond between the p- and h-sites. To elucidate the structure of the saddle point in more detail, the search of the transition state (TS) between the p- and h-sites is carried out at the AM1 level. The structure of the TS is slightly deviated from the center of the C-C bond. The imaginary frequencies corresponding to a translational mode of Li+ from p- to h-sites are calculated to be 208i cm-1. B. Diffusion of Li+ Ion on C60. First, the trajectory calculations were carried out at lower temperatures below 50 K. However, the Li+ ion does not move on the surface of C60 at 10-40 K. The diffusion of Li+ ion takes place from 50 K. Snapshots of the C60-Li+ system at 300 K are illustrated in Figure 3. At time zero, the Li+ ion is located in the hexagonal site, and it moves gradually on the surface. The ion reaches the
Diffusion of Li+ on C60
J. Phys. Chem. C, Vol. 111, No. 35, 2007 13089
Figure 2. Potential energy curves along the diffusion path of Li+ between pentagonal and hexagonal sites. The calculations are carried out at the B3LYP/LANL2DZ and AM1 methods. The heights of the Li+ ion are fixed to those of the optimized structure. Inset figures show the structures of LiC60+ at the pentagonal and hexagonal sites and saddle point.
backside of the starting point on C60 at 0.60 ps. The diffusion path in this trajectory at 300 K is approximately given by
hexagonal site (h) f C-C bond f h f CC bond f h f p f C-C bond f p f h f ...
Figure 3. Snapshots of the LiC60+ system at 300 K obtained by direct MO-MD calculation. The Li+ ion is started from the hexagonal site of C60.
where C-C bond means that the Li+ moves along the C-C bond of the hexagonal or pentagonal site. After 1.20 ps, the ion returns again near the starting point. Time propagations of the position of the Li+ ion (R1 and R2) are plotted in Figure 4 as a function of time. At time zero, the distances R1 and R2 are almost equivalent to each other (2.624 and 2.623 Å, respectively). Both distances are gradually increased with increasing time. The distance R1 reaches a maximum value (R1 ) 9.18 Å) at 0.60 ps (backside of the starting point). The distance R1 decreases and reaches a minimum value at 1.30 ps. The profile of distance R2 is very similar to that of R1. These results indicate that the Li+ ion rotates around C60 with a time period of about 1.30 ps at 300 K. From the calculation, it is found that the diffusion coefficient at 300 K is calculated to be 1.67 × 10-10 m2/s for this sample trajectory.
Diffusion dynamics of lithium ion on C60 are essentially similar to that of graphite. The Li+ ion moves preferentially along the node of the HOMO in both cases. The difference is activation energies along the diffusion path. The Li+ ion can diffuse more easily on C60 than on the graphite surface. This difference is important to develop a new material as described in the Discussion section. C. Check of Level of Theory. To check the level of theory used in the direct MO-MD calculation, the potential energies of the system are calculated at both AM1 and B3LYP/ LANL2DZ levels of theory. The results for 300 K are potted in Figure 5 as a function of time. As clearly seen in this figure, the result of the AM1 calculation is in good agreement with that of B3LYP/LANL2DZ. This agreement strongly indicates that the direct MO-MD calculation at the AM1 level of theory gives a reasonable value for the potential energy of the system.
TABLE 1: Optimized Parameters of Li+ Adsorbed on C60 Molecules and Relative Energies between Hexagonal and Pentagonal Sites (∆E in kcal/mol)a B3LYP/LANL2DZ
a
B3LYP/6-31G(d)
AM1
parameter
hexagonal
pentagonal
hexagonal
pentagonal
hexagonal
pentagonal
R1 R2 θ hb ∆E
2.4626 2.4657 54.3 1.9988 0.0
2.3726 2.3739 60.0 2.0539 -0.97
2.4486 2.4351 54.3 1.978 0.0
2.3190 2.3203 57.4 1.9528 +0.18
2.6209 2.6261 56.7 2.195 0.0
2.5737 2.5737 62.5 2.283 -0.87
Bond lengths and angles are in angstroms and degrees, respectively. b Height of Li+ ion from the C60 surface.
13090 J. Phys. Chem. C, Vol. 111, No. 35, 2007
Tachikawa
Figure 4. Time propagation of positions of the Li+ ion (R1 and R2) at 300 K plotted as a function of time.
Figure 5. Comparison of the potential energies plotted as a function of time calculated by AM1 and with those calculated by the B3LYP/ LANL2DZ methods. The direct MO-MD calculation was carried out at 300 K.
Figure 7. Ion-switching cell material designed theoretically on the basis of the direct MO-MD calculation. The arrow means a trajectory of Li+ superimposed on the material. At low temperature (below 300 K), the Li+ ion moves preferentially along the C60 one-dimensional nanowire. On the other hand, the Li+ ion moves on both the C60 wire and graphite sheets at medium temperature near 300 K. At higher temperatures (above 300 K), the Li+ ion diffuses preferentially on the graphite sheets. The diffusion path of the Li+ ion is selectively changed by thermal activation of the system.
(50-300 K). The activation energy is estimated to be 0.2 kcal/ mol. For comparison, the diffusion coefficients for the Li+ ion on the graphite surface15 are plotted in Figure 5 as an Arrhenius plot. The results can be summarized as follows. (1) At low temperature (T < 300 K), the diffusion coefficients for the Li+ ion on the C60 surface are larger than those of the graphite. (2) The diffusion coefficients of the Li+ ion on both C60 and graphite surfaces are almost equivalent at medium temperatures around 300 K. (3) At higher temperatures (T > 300 K), the coefficients for graphite are significantly larger than those of C60. 4. Discussion Li+
Figure 6. Arrhenius plots of the diffusion coefficients of the ion on C60 in the temperature range of 80-300 K. For comparison, the diffusion coefficients for the Li+ ion on the graphite surface are plotted. The data for the graphite surface are cited from ref 15.
D. Arrhenius Plots of Diffusion Coefficients. The dynamics calculations are performed at 50, 100, 150, 200, 250, and 300 K. The Arrhenius plot of the diffusion coefficients is given in Figure 6. The plot is almost linear in the temperature range
The present results showed that the thermal behavior of the Li+ ion on C60 is much different from that on the graphite surface. On the basis of these evidences, we design a highperformance molecular device composed of C60 and graphite sheet. The designed molecular device is illustrated in Figure 7. This device is composed of one-dimensional C60 nanowire and graphite sheets. The present study suggested that the diffusion coefficient of Li+ on C60 is larger than that of Li+ on the graphite surface at low temperature below 300 K. Therefore, the Li+ ion diffuses preferentially on the C60 nanowire at lower temperatures. At intermediate temperature around 300 K, the
Diffusion of Li+ on C60 Li+ ion can move with the same diffusion rate on both C60 wire and graphite sheets. Therefore, two-dimensional diffusion of Li+ ion is possible at the intermediate temperature. At higher temperatures (T > 300 K), on the other hand, the diffusion on the graphite sheets is faster than that on C60. Hence, the Li+ ion moves preferentially on the graphite surface. Thus, the diffusion path of the Li+ ion on the designed molecular device is drastically changed by the control of temperature. Acknowledgment. The author is indebted to the Computer Center at the Institute for Molecular Science (IMS) for the use of the computing facilities. H.T. also acknowledges a partial support from a Grant-in-Aid from the Ministry of Education, Science, Sports and Culture of Japan. Supporting Information Available: Figure of temperature and potential energy of the Li+C60 system plotted as a function of simulation time (Figure S1). This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Ellison, E. H. J. Phys. Chem. B 2006, 110, 11406. (2) Amusia, M. Y.; Liverts, E. Z.; Mandelzweig, V. B. Phys. ReV. 2006, 74, No. 042712. (3) Wu, J. W.; Gan, Z. Z. Phys. ReV. B 2005, 72, No. 014528. (4) Kato, H.; Takemura, S.; Iwasaki, K.; Watanabe, Y.; Nanba, N.; Hiramatsu, T.; Nishikawa, O.; Taniguchi, M. J. Vac. Sci. Technol. 2006, 24, 1500. (5) Radomski, A.; Jurasz, P.; Alonso-Escolano, D.; Drews, M.; Morandi, M.; Malinski, T.; Radomski, M. W. Br. J. Pharmacol. 2005, 146, 882. (6) Lavrentiev, V.; Naramoto, H.; Narumi, K.; Sakai, S.; Avramov, P. Chem. Phys. Lett. 2006, 423, 366. (7) Chen, Y.; Stepniak, F.; Weaver, J. H.; Chibante, L. P. F.; Smalley, R. E. Phys. ReV. B 1992, 45, 8845.
J. Phys. Chem. C, Vol. 111, No. 35, 2007 13091 (8) Gu, C.; Stepniak, F.; Poirier, D. M.; Jost, M. B.; Benning, P. J.; Chen, Y.; Ohno, T. R.; Martins, J. L.; Weaver, J. H.; Fure, J.; Smalley, R. E. Phys. ReV. B 1992, 45, 6348. (9) Poirier, D. M.; Ohno, T. R.; Kroll, G. H.; Benning, P. J.; Stepniak, F.; Weaver, J. H.; Chibante, L. P. F.; Smalley, R. E. Phys. ReV. B 1993, 47, 9870. (10) Grpmov, A.; Krawez, N.; Lassesson, A.; Ostrovskii, D. I.; Campbell, E. E. B. Curr. Appl. Phys. 2002, 2, 51. (11) Zhao, Y. L.; Pan, X. M.; Zhou, D. F.; Su, Z. M.; Wang, R. S. Synth. Met. 2003, 227, 135. (12) Santos, J. D.; Bulhoes, L. O. S.; Longo, E.; Varela, J. A. J. Mol. Struct. (THEOCHEM) 1995, 335, 149. (13) Pitarch-Ruiz, J.; Evangelisti, S.; Maynau, D. J. Chem. Theory Comput. 2005, 1, 1079. (14) Ohno, K.; Maruyama, Y.; Esfarjani, K.; Kawazoe, Y. Phys. ReV. Lett. 1996, 76, 3590. (15) Tachikawa, H.; Shimizu, A. J. Phys. Chem. B 2005, 109, 13255. (16) Tachikawa, H.; Shimizu, A. J. Phys. Chem. B 2006, 110, 20445. (17) Tachikawa, H. J. Chem. Phys. 2006, 125, 133119. (18) Tachikawa, H. J. Chem. Phys. 2006, 125, 144307. (19) Tachikawa, H. J. Phys. Chem. A 2006, 110, 153. (20) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; di Nola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (21) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A., Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision B.04; Gaussian, Inc.: Pittsburgh, PA, 2003.