Ground state analysis of magnetic nanographene molecules with modified edge

Chemical Physics 415 (2013) 64–68

Contents lists available at SciVerse ScienceDirect

Chemical Physics journal homepage: www.elsevier.com/locate/chemphys

Ground state analysis of magnetic nanographene molecules with modified edge Narjes Gorjizadeh a,⇑, Norio Ota b, Yoshiyuki Kawazoe a a b

Institute for Materials Research, Tohoku University, 2-1-1 Katahira, Aoba-ku, Sendai 980-8577, Japan Pure and Applied Sciences, University of Tsukuba, 1-1-1 Tennodai, Tsukuba 305-8573, Japan

a r t i c l e

i n f o

Article history: Received 13 May 2012 In final form 19 December 2012 Available online 3 January 2013 Keywords: Graphene Molecular magnetism First-principle calculations

a b s t r a c t We study spin states of edge modified nanographene molecules with rectangular and triangular shapes by first principle calculations using density functional theory (DFT) and Hartree–Fock (HF) methods with Møller–Plesset (MP) correlation energy correction at different levels. Anthracene (C14H10) and phenalenyl (C13H9), which contain three benzene rings combined in two different ways, can be considered as fragments of a graphene sheet. Carbon-based ferromagnetic materials are of great interest both in fundamental science and technological potential in organic spintronics devices. We show that non-magnetic rectangular molecules such as C14H10 can become ferromagnetic with high-spin state as the ground state by dihydrogenization of one of the zigzag edges, while triangular molecules such as C13H9 become ferromagnetic with high-spin state by dehydrogenization of one of the zigzag edges. Ó 2012 Elsevier B.V. All rights reserved.

1. Introduction Carbon-based magnetism is a new emerging concept, which attracted the attention of theoretical and experimental scientists. Unlike conventional magnetic materials, whose magnetic properties are due to 3d orbital, in organic materials p orbitals are responsible for magnetization. The advantage of organic materials is their light weight, flexibility, low spin–orbit interaction and low cost compared to the conventional magnets, which make them promising candidates as organic ferromagnets and for organic spintronics. Room-temperature carbon-based magnetism is evidenced in experiments since 1991 [1–9]. Magnetism in p-conjugated systems is considered to be due to the number of unpaired electrons [10]. It was shown that graphite samples with certain defects have spontaneous magnetization [5]. Existence of ferromagnetism at the zigzag edges of graphene nanoribbons and graphene flakes due to localized edge states is also demonstrated in theoretical studies [11]. However, appearance of the edge states in graphene nanoflakes depends on the size of the structure. Depending on the shape of the nanographene, edge states appear at a critical size [12,13]. As the properties of the graphene nanoribbons and nanoflakes are associated with their edge states, modifications of the edge affect their electronic and magnetic properties to a significant value [14,15]. Kusakabe et al. [16,17] have shown that dihydrogenization of one of the edges of zigzag graphene nanoribbons can

⇑ Corresponding author. E-mail address: [email protected] (N. Gorjizadeh). 0301-0104/$ - see front matter Ó 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.chemphys.2012.12.025

make them ferromagnetic. Ota et al. [18,19] have also shown by DFT calculations that edge modified asymmetric nanographenes with one dihydrogenated or dehydrogenated zigzag edge are ferromagnetic structures. The zigzag edges in these structures are large enough to have localized edge states and have singlet ground state before edge modification. In this work we look for strong magnetism in small nanographene molecules by modifications of their zigzag edge. We seek the spin ground state of C14H7,10,13 and C13H7,9,11 using first-principle calculations by density functional theory (DFT) and Hartree–Fock (HF) method with MPn correlation energy corrections. C14H10 and C13H9 contain three benzene rings which are combined in two different ways to yield structures with different shapes and magnetic properties. C14H10 is a rectangular molecule with non-magnetic zigzag edges due to small size of the molecule. While C13H9 is a triangular molecule with one net spin, due to the difference in the number of the two subsets of site A and B, according to Lieb’s rule [20]. Both of these structures can be considered as fragments of a graphene sheet [21,22] and their magnetic properties depend on their shape. We study the stable spin state of these molecules with modified edges by dihydrogenization and dehydrogenization of one of the zigzag edges. Energetically, dihydrogenization of all three carbon atoms of one zigzag edge in C14H13 is not the most stable structure, according to Clar’s sextet rule [23]. The extra hydrogen atoms tend to attach the zigzag and armchair sites of one ring to preserve the sextet. The dihydrogenated structures presented in this study are physically and artificially modified structures which may be induced by artificial proton irradiation discussed in the experimental work by Esquinazi et al. [5].

N. Gorjizadeh et al. / Chemical Physics 415 (2013) 64–68

We show that the rectangular molecules can become ferromagnetic by dihydrogenization, while triangular molecules become ferromagnetic by dehydrogenization of one the zigzag edges. Although different methods do not agree in the ground state of the structures between high-spin and low-spin states, we show that the DFT methods give the ground states in agreement with the counting rules for hydrogenated, dihydrogenated and dehydrogenated edges [20,24,25]. HF-based methods even with the higher accuracy of MP4 level do not lead to reliable results due to high spin contamination in at least one of the spin states. Using MPn correlation energy correction of the HF method up to the 4th order of approximation, we show that the C–C bond-length of the molecule can change due to exchange interaction in transition from low-spin to high-spin state, resulting in changing the sign of the exchange interaction and obtaining a different spin ground state, which does not agree with DFT methods. 2. Method The structures are optimized for high spin and low spin states by DFT and HF methods, using Gaussian03 and Gaussian09 [26]. We compare the following methods and basis sets: DFT/UB3LYP/ 6-31+g, DFT/UPBE/6-31+g, UHF/6-31+g, UMP2/6-31+g, UMP2/6311+g(d,p), UMP4/6-311+g(d,p). The spin contamination is zero or negligible for DFT (UPBE) and hybrid DFT (UB3LYP) calculations, but it is large for most of the structures for HF-based methods. All the structures with non-zero net spins are calculated with unrestricted orbitals to allow different spatial components for spin up and spin down orbitals. 3. Results Starting from these molecules, C13H9 and C14H10, we study the ground state of the structures with different modified edges, i.e. the hydrogenated edge carbons, dehydrogenation and dihydrogenation of one of the zigzag edges, as depicted in Fig. 1. After optimizing the geometries for high-spin and low-spin states, the favored spin state and the energy gap between the two states are studied. Geometries are optimized with energy accuracy of 1  e 5 eV. The high spin and low spin states of the structures of Fig. 1(b)–(f) have total magnetic moment of 1 lB and 3 lB, respectively; while the structure of Fig. 1(a) has total magnetic moment of 0 and 2 lB. The two hydrogen atoms of the dihydrogenated carbons of Fig. 1(c) and (f) are out of plane, with similar angle with the

65

plane of molecule. Table 1 shows the ground state and the absolute value of high-spin/low-spin energy gap (H/L gap) predicted by different methods. The spin-squared expectation value, , for different methods is also reported at Table 2, to compare the spin contamination, i.e. the deviation from the exact value of = s(s + 1), which is also shown in the table. UDFT methods have zero or negligible spin contamination, while UHF-based methods have significant spin contamination, which yields error in the ground state wavefunction and the calculated geometry. The problem of spin contamination is common for unrestricted orbitals in the UHF procedure, especially for the systems that the unpaired electrons are delocalized. Increasing electron correlation and size of the basis set may reduce the spin contamination. Hence, we used MPn correlation energy corrections up to the 4th level of approximation to reduce the spin contamination. Yet, The UMP4 method has small or negligible spin contamination for one of the spin states, while still large contamination is observed in the other spin state. All methods find the ground state of C14H10 to be the low-spin state, i.e. non-magnetic, as also predicted by Lieb’s rule [20]. However the energy gap is varied between 1.017 and 3.218 eV. But different methods do not agree with the ground state of C14H7, with three radical carbons. Disregarding UHF-based methods, which have large spin contamination in at least one of the spin states, the other methods, i.e. UB3LYP and UPBE, agree in the ground state and show that the low-spin state is the favored structure energetically. Higher stability of the low-spin state for C14H7 is in agreement with the prediction by the counting rule of radical carbons suggested by Ota et al. [24]. For this structure, the three unpaired r bonds at the edge will have lower exchange energy in the antiparallel configuration. Spin density of the C14H7 for the two spin states calculated by UB3LYP functional is sketched in Fig. 2(a) and (b). The parallel spin configuration in the high-spin state in this figure shows the high unfavorable exchange energy. In case of C14H13, the structure did not converge for UMP2 and UMP4 methods, but the DFT methods predict the high-spin state to be the ground state. However, the structure of a smaller counterpart of the molecule (C10H10), which has two dihydrogenated edge carbon, is optimized with similar methods as Table 1 and the highspin state is the ground state. In low-spin structure of this molecule calculated by UMP2 and UMP4 methods, most of the H atoms are distorted from the plane of the molecule, giving an untidy structure. The counting rule by Maruyama et al. [25] also predicts the high-spin state as the ground state of C10H10 and C14H13. Fig. 2(c) and (d) depicts the spin density of the C14H13 for the

(a) C14H10

(b) C14H7

(c) C14H13

(d) C13H9

(e) C13H7

(f) C13H11

Fig. 1. Schematic picture of graphene molecules, showing connectivities and hydrogen counts on each carbon: (a) C14H10 , with hydrogenated edge carbons (b) C14H7 , with three dehydrogenated edge carbons (c) C14H13 , with three dihydrogenated edge carbons (d) C13H9 , with hydrogenated edge carbons (e) C13H7 , with two dehydrogenated edge carbons (f) C13H11 , with two dihydrogenated edge carbons. The two H atoms of the dihydrogenated edges in (c) and (f) are out of plane with the same angle with the plane of molecule.

66

N. Gorjizadeh et al. / Chemical Physics 415 (2013) 64–68

Fig. 2. Spin density isosurface of C14H7 with magnetization of 1 lB (a), and 3 lB (b), and C14H13 with magnetization of 1 lB (c), and 3 lB (d).

Table 1 Ground state and High-spin/Low-spin energy gap of the structures optimized with different levels of accuracy. Structure

Method

Ground state (magnetic moment)

|High-spin/Low-spin energy gap| (eV)

C14H10

B3LYP/6-31+g PBEPBE/6-31+g HF/6-31+g MP2/6-31+g MP2/6-311+g(d,p) MP4/6-311+g(d,p) UB3LYP/6-31+g UPBEPBE/6-31+g UHF/6-31+g UMP2/6-31+g UMP2/6-311+g(d,p) UMP4/6-311+g(d,p) UB3LYP/6-31+g UPBEPBE/6-31+g UHF/6-31+g UMP2/6-31+g UMP2/6-311+g(d,p) UMP4/6-311+g(d,p) UB3LYP/6-31+g UPBEPBE/6-31+g UHF/6-31+g UMP2/6-31+g UMP2/6-311+g(d,p) UMP4/6-311+g(d,p) UB3LYP/6-31+g UPBEPBE/6-31+g UHF/6-31+g UMP2/6-31+g UMP2/6-311+g(d,p) UMP4/6-311+g(d,p) UB3LYP/6-31+g UPBEPBE/6-31+g UHF/6-31+g UMP2/6-31+g UMP2/6-311+g(d,p) UMP4/6-311+g(d,p)

Low spin (0 lB) Low spin (0 lB) Low spin (0 lB) Low spin (0 lB) Low spin (0 lB) Low spin (0 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) High spin (3 lB) High spin (3 lB) Low spin (1 lB) High spin (3 lB) High spin (3 lB) High spin (3 lB) –a –a –a Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) High spin (3 lB) High spin (3 lB) High spin (3 lB) High spin (3 lB) High spin (3 lB) High spin (3 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB) Low spin (1 lB)

1.838 1.645 1.017 3.218 2.797 2.426 1.232 1.271 1.269 0.247 0.048 0.060 0.501 0.501 1.451 –a –a –a 3.879 3.881 6.074 2.402 2.303 22.683 0.320 0.594 0.738 0.188 0.066 0.583 2.758 2.689 2.331 2.718 2.874 2.735

C14H7

C14H13

C13H9

C13H7

C13H11

a

Structure did not converge in relaxation.

two spin states calculated by UB3LYP functional. The hybridization of the dihydrogenated carbons changes to sp3 which results in moving the p orbitals towards the other carbons, which have sp2

hybridization. The stable spin state of the zigzag edge with p orbital is the parallel spin configuration on the carbons of the same subset (A or B), as achieved in the high-spin state in Fig. 2(d). The major disagreement of UDFT and UMP4 methods is in finding the ground state of C13H9. Both functionals in the UDFT method favor the low-spin state energetically, while UMP4 method predicts the high-spin state as the ground state, with the H/L gap of 22.683 eV. While the energy difference of the two states is around 3.88 eV for UDFT methods. The optimized structure of C13H9 by UMP4 indicates that the exchange interaction stretches the C–C bonds of the high-spin state compared to the low-spin. The three C–C bonds in the middle of the molecule changes from 1.41 to 1.50 Å, in transition from low-spin to high-spin state, so that the molecule is stretched in its plane. Fig. 3 shows the C–C bond lengths for the two spin states of C13H9 optimized by UMP4. UDFT method, on the other hand, distorts one of the hydrogen atoms to go out of plane to converge C13H9 in the high-spin state. This distortion of the geometry raises the energy of the structure, as can be noticed from the spin densities sketched in Fig. 4(b), and makes the high-spin state less favorable energetically. The low-spin state on the other hand, has a planar structure (Fig. 4(a)). Using UDFT calculation with UB3LYP/6-31G(d), Radovic et al. [27] show that C13H8 with one radical carbon at the edge, has the stable edge. The ground state of C13H8 is shown to be the high-spin state with H/L gap of 0.533 eV [27]. These radical carbons at the edges are produced in the heat treatment process of Highly Ordered Pyrolytic Graphite (HOPG) production [28]. However, there is not an agreement in the stability or reactivity of radical edge carbons with rdangling bonds [27–31]. While UDFT method predicts that C13H9 does not have stable edge with 9 hydrogen atoms, UMP4 method finds this molecule stable with a planar structure. The middle C– C bonds of the optimized structure of low-spin (high-spin) states using UMP4 is stretched up to 2% (7%) compared to the edge C–C bonds. On the other hand, when the structure goes from low-spin to high-spin state, the middle (edge) C–C bonds are stretched for 3% (6%) due to the larger exchange interaction (Fig. 3). The spin contamination in the high-spin state, which is the ground state in UMP4 method, is zero, while the low-spin state has a large spin contamination. Hence spin densities are not calculated properly for the low spin state and comparison of the energy of the two spin states does not lead to a reliable result. The large H/L gap is due

67

N. Gorjizadeh et al. / Chemical Physics 415 (2013) 64–68 Table 2 Spin-squared expectation value, , for different methods, and the exact value of = s(s + 1). Structure

Magnetic moment (lB)

s(s + 1)

C14H10

0 2 1 3 1 3 1 3 1 3 1 3

0 2 0.75 3.75 0.75 3.75 0.75 3.75 0.75 3.75 0.75 3.75

C14H7 C14H13 C13H9 C13H7 C13H11 a

B3LYP/6-31+g

PBEPBE/6-31+g

HF/6-31+g

MP2/6-31+g

MP2/6-311+g(d,p)

MP4/6-311+g(d,p)

0 2.00 0.82 3.75 0.90 3.75 0.75 3.75 0.85 3.75 0.75 3.75

0 2.00 0.80 3.75 0.78 3.75 0.75 3.75 0.75 3.75 0.75 3.75

0 2.81 11.77 4.95 3.62 4.29 3.56 3.75 7.89 5.59 3.20 4.15

0 2.09 11.21 8.02 –a 4.28 3.50 3.75 7.76 5.54 2.56 4.15

0 2.02 9.46 6.87 –a 4.13 2.76 3.75 6.77 5.03 1.83 4.036

0 2.81 1.48 3.78 –a 4.13 2.96 3.75 0.75 5.04 1.83 4.04

Structure did not converge in relaxation.

(b)

(a)

1.378 1.410 1.410

1.418

1.378 1.410

1.414 1.414 1.414

1 .3 7 8

1.410 1 .3 7 8

1.396 1.396

1.378

1.396 1.496 1.496 1.496

1 .4 1 8

1.410

1.418

1.396 1 .4 1 8

1.396 1.410

1.378

1.418

1.396

1.418

Fig. 3. C–C bond lengths (Å) of C13H9 for the low-spin state (a) and high-spin state (b) optimized by UMP4 method. The H atoms are not shown for more clarity of the bond lengths.

Fig. 4. Spin density isosurface of C13H9 with magnetization of 1 lB (a), and 3 lB (b), and C13H11 with magnetization of 1 lB (c), and 3 lB (d).

to this error in the energy of the low spin state. Prediction of the low-spin state as the ground state by DFT method, is in agreement with the counting rule of Maruyama et al. [25].

The two methods (UDFT and UMP4) disagree on the ground state of C13H7 also. In this molecule, high-spin state is favorable by UDFT, while the low-spin state has lower energy in UMP4. In

68

N. Gorjizadeh et al. / Chemical Physics 415 (2013) 64–68

this case, UMP4 method lowers the energy of the low-spin state by breaking the symmetry of the molecule. However, the spin contamination for the high-spin state in UMP4 method result in raising its total energy and making the low-spin as the ground state. There is a spin contamination of 6%–9% for the low-spin states calculated by UDFT methods, but it is fairly small. The Counting rule of Ota et al. [24] also predicts the high-spin state as the ground state of C13H7, in agreement with UDFT methods. C13H11 on the other hand, has lower energy for the low-spin state as predicted by all methods, and the counting rule by Maruyama et al. [25], with energy difference of 2.331–2.874 eV. As depicted in Fig. 4(c) and (d), the up–up spin pair in the high-spin state results in a larger exchange energy and increases its total energy. The sp3-like carbons in this structure also have less p character, with a very small spin density. Comparison of the UHF, UMP2 and UMP4 methods indicate that electron correlation plays a significant role in obtaining the electronic structure of the magnetic molecules. Correction of the correlation function up to UMP4 level helps to decrease the spin contamination, but still it remains high for one of the spin states. Increasing the size of the basis set also affects the electronic structure calculation, but the problem with the spin contamination still remains. This deficiency of the UHF-based methods in calculating the spin states of was also discussed by Suhai [32]. 4. Conclusion UDFT and UHF-based calculations have been performed for C13H9 and C14H10 molecules with modified edges to obtain the most favorable spin state. UHF-based methods even with the higher accuracy of the 4th level of perturbation do not lead to reliable results due to large spin contamination in at least one of the spin states. UDFT based calculations, on the other hand, are in agreement with the ground states predicted by the counting rules in references [20,24,25]. Based on our DFT results, we show that the non-magnetic rectangular-shape graphene molecules, such as C14H10, can become ferromagnetic with high spin state ground state by dihydrogenization of one of the edges. Triangular-shape graphene molecules, such as C13H9, on the other hand, which are low-spin states at their ground state, will become ferromagnetic structures with high-spin ground state by dehydrogenization of one of the edges. These light-weight ferromagnetic nanographene molecules are applicable as organic magnets and also in organic spintronic devices. Acknowledgements Narjes Gorjizadeh would like to thank the crew of the Center for Materials Research of Tohoku University for their support of the

Hitachi SR11000 (model K2) supercomputer, and Global COE program ‘‘Materials Integration (International Center of Education and research), Tohoku university’’ MEXT, Japan, for financial support. References [1] M. Tamura, Y. Nakazawa, D. Shiomi, K. Nozawa, Y. Hosokoshi, M. Ishikawa, M. Takashi, M. Jnoshita, Chem. Phys. Lett. 186 (1991) 401. [2] J.M. Manriquez, G.T. Yee, R.S. McLean, A.J. Epstein, J.S. Miller, Science 252 (1991) 1415. [3] M. Allemand, K.C. Khemani, A.F. Wudl, K. Holczer, S. Donovan, G. Gruner, J.D. Thompson, Science 253 (1991) 301. [4] Y. Nakazawa, D. Shiomi, K. Nozawa, Y. Hosokoshi, M. Ishikawa, M. Takahashi, M. Kinoshita, Chem. Phys. Lett. 186 (1991) 401. [5] P. Esquinazi, D. Spemann, R. Höhne, A. Setzer, K. Han, T. Butz, Phys. Rev. Lett. 91 (2003) 227201. [6] K. Kamishima, T. Noda, F. Kadonome, K. Kakizaki, N. Hiratsuka, J. Magn. Magn. Mater. 310 (2007) e346. [7] Y. Wang, Y. Huang, Y. Song, X. Zhang, Y. Ma, J. Liang, Y. Chen, Nano Lett. 9 (2009) 220. [8] X. Yang, H. Xia, X. Qin, W. Li, Y. Dai, X. Liu, M. Zhao, Y. Xia, S. Yan, B. Wang, Carbon 47 (2009) 1399. [9] T. Saito, D. Nishio-hamane, S. Yoshii, T. Nojima, Appl. Phys. Lett. 98 (2011) 052506. [10] P.P. Letinen, A.S. Foster, Y. Ma, A.V. Krasheninnikov, R. Nieminen, Phys. Rev. Lett. 93 (2004) 187202. [11] M. Fujita, K. Wakabayashi, K. Nakada, K. Kusakabe, J. Phys. Soc. Jpn. 65 (1996) 1920. [12] D.E. Jiang, B.G. Sumpter, S. Dai, J. Chem. Phys. 127 (2007) 124703. [13] D.E. Jiang, S. Dai, Chem. Phys. Lett. 466 (2008) 72. [14] N. Gorjizadeh, A.A. Farajian, K. Esfarjani, Y. Kawazoe, Phys. Rev. B 78 (2008) 155427. [15] N. Gorjizadeh, Y. Kawazoe, J. Nanomater. 2010 (2010) 513501. [16] K. Kusakabe, M. Maruyama, Phys. Rev. B 67 (2003) 092406. [17] M. Maruyama, K. Kusakabe, J. Phys. Soc. Jpn. 73 (2004) 656. [18] N. Ota, N. Gorjizadeh, Y. Kawazoe, J. Magn. Soc. Jpn. 34 (2010) 573. [19] N. Ota, N. Gorjizadeh, Y. Kawazoe, J. Magn. Soc. Jpn. 36 (2012) 36. [20] E.H. Lieb, Phys. Rev. Lett. 62 (1989) 1201. [21] Y. Morita, S. Suzuki, K. Sato, T. Takui, Nat. Chem. 3 (2011) 197. [22] H. Nagai, M. Nakano, K. Yoneda, R. Kishi, H. Takahashi, A. Shimizu, T. Kubo, K. Kamada, K. Ohta, E. Botek, B. Champagne, Chem. Phys. Lett. 489 (2010) 212. [23] E. Clar, Polycyclic Hydrocarbon, Academic, London, 1964. [24] N. Ota, N. Gorjizadeh, Y. Kawazoe, J. Magn. Soc. Jpn. 35 (2011) 414. [25] M. Maruyama, K. Kusakabe, S. Tsuneyuki, K. Akagi, Y. Yoshimoto, J. Yamauchi, J. Phys. Chem. Solids 65 (2004) 119. [26] M.J. Frisch, G.W. Trucks, H.B. Schlegel et al., GAUSSIAN 09, Revision B.01, Gaussian, Inc., Wallingford, CT, 2009. [27] L.R. Radovic, B. Bockrath, J. Am. Chem. Soc. 127 (2005) 5917. [28] D. Ingram, D. Austen, in: Proceedings of Conferences on Industrial Carbon Graphite, London, 1958. [29] M. Shibyama, H. Sato, T. Enoki, X.X. Bi, M.S. Dresselhaus, M. Endo, J. Phys. Soc. Jpn. 69 (2000) 754. [30] T. Kawai, Y. Miyamoto, O. Sugino, Y. Koga, Phys. Rev. B 62 (2000) R16346. [31] S. Mrozowski, Phys. Rev. 85 (1952) 609. [32] S. Suhai, Phys. Rev. B 51 (1995) 16553.