Phase stability and structural distortion of NiO under high pressure
ZHANG Wei-bing(张卫兵), HU Yu-lin(胡玉林), TANG Bi-yu(唐壁玉)
Key Laboratory of Advanced Materials and Ryeological Properties of Ministry of Education,
Department of Physics, Xiangtan University 411105, China
Received 10 April 2006; accepted 25 April 2006
Abstract: The phase stability and structural distortion of NiO under high pressure were investigated using first-principles calculations based on density-functional theory. Different forms of exchange-correlation functional including LDA, GGA and GGA+U were used in the present calculations. All of the three methods predict NiO to be AFM II ordering with the cell slightly compressed along [111] direction and also indicate that there is no structural phase transition of NiO under pressure up to 140 GPa, which are in agreement with the experiment. However, both LDA and GGA incorrectly predict the structural distortion under pressure especially above 60 GPa. Only when strong correlations are included in form of GGA+U, structural distortion under high pressure can qualitatively agree with the experiment. The related mechanism was also analyzed and discussed. These results suggest that the strong electronic correlations still play a very important role in the properties of NiO under high pressure.
Key words: NiO; phase stability; structural distortion; strong correlation; first-principles calculation; high pressure
1 Introduction
As a prototype of Mott insulators and perfect candidate of antiferromagnetic material in FM-AFM exchange bias system [1], NiO has attracted a lot of attention over the past decades. Determining the nature of this classic strongly correlated material is a basic step to understand high Tc superconductivity and exchange bias phenomena.
It’s well known that NiO is a type-II antiferromagnetic insulator with the Néel temperature of 523 K. It crystallizes in the simple cubic rocksalt (B1) structure above TN but becomes rhombohedrally distorted B1 structure (rB1) compressed along [111] direction with rhombohedral angle αrh equal to 60.08? [2] as temperature is below TN. For geophysical interests, equation of state (EOS) has been measured by static compression up to 28 GPa [3] and by shock compression up to 147 GPa [4]. Both experiments indicated that no pressure-induced structure phase transition was observed in experimental pressure range [3, 4]. This is contrasted to other 3d transition-metal oxide (TMO) such as FeO and MnO, in which structural phase transitions such as from B1 to NiAs (B8) structure [5] were observed at high pressures. Recently ETO et al [2] has investigated in detail the pressure effect on the lattice constants expressed in hexagonal lattice up to 141 GPa by static compression using two sets of equipments (named run A and run B), a, c, axial ratio c/a decrease monotonically with increasing pressure without a significant change above 60 GPa which suggested by a early Local Spin Density Function calculation [6].
Because of strong electronic correlations between 3d electrons, the theoretic studies of NiO have long been a challenge to investigators. Using local spin density approximation (LDA), the structural change of NiO with pressure has been investigated theoretically by SASAKI [6]. The calculations suggest that the lattice distortion is enhanced with increasing pressure and becomes extremely large at 60 GPa. A phase transition from rB1 to CsCl (B2) structure at 318 GPa was also predicted. As mentioned previously, a significant change of pressure coefficient of c/a above 60 GPa was not observed in recently experiment. ETO et al [3] suggested this discrepancy between the theory and experiment maybe origin from the fact that LDA was not appropriate to describe the electronic correlations between 3d electrons in NiO.
In order to explain the experimental results, we investigate the structural properties under high pressure in this paper with different exchange correlation functional, and the strong on-site coulomb repulsion of Ni 3d electrons are included in form of GGA+U which are not considered in SASAKI’s calculations.
2 Computational methods
All calculations have been performed with VASP (Vienna Ab initio Simulation Package) [7-11], a first principles plane-wave code based on spin-polarized density functional theory. The interaction between ions and valence electrons was described by the projector augmented-wave (PAW) method [12]. The KOHN- SHAM equations were solved via iterative matrix diagonalization based on the minimization of the norm of the residual vector to each eigenstate and optimized charge- and spin-mixing routines [13-15].
The generalized gradient corrections added in form of PERDEW-WANG functional PW91 [16] were chosen for the exchange correlation function, the spin interpolation of VOSKO et al was also used [17]. To correct the strong electronic correlation, a simple rotationally invariant DFT+U version proposed by DUDAREV et al [18, 19] and implemented in VASP [20] was used. In this method, the parameters U and J did not enter separately, only the difference U-J was meaningful. Parameters U and J represented on-site Coulomb interaction energy and exchange energy respectively. J was kept to 1 eV during all calculations. A detailed description of the DFT+U method can be found in Refs [20, 21]. The change of U with pressure was neglected up to experimental pressure as many authors did in other system like LaMnO3 [22]. A value of U-J=5.3 eV was used in our calculations, U-J= 8 eV suggested by constrained LDA computations [19] and a smaller value U-J=2 eV were also computed for comparison.
All results reported in this paper were carried out on a rhombohedral antiferromagnetic supercell including two NiO formula unit cells. Convergence tests have been checked carefully both for planewave cutoff energy and k points sample, a plane-wave set expanded in energy cutoff 600 eV and k-points sample with a mesh of points 8×8×8 generated by the scheme of MONKHORST and PACK [23] can ensure the total energies difference less than 3 meV/atom. The unit cell shape optimizations for each volume were performed using the conjugate gradient method and a Gaussian-smearing with a width of 0.2 eV. Forces acting on atoms and stress tensors on unit cells were used in the optimization process. For total energy and DOS calculation, the integration over the Brillouin zone was performed using the linear tetrahedron method with bl?chl corrections [24-26]. The calculation total energies as a function of volume were fitted to the MURNAGHAN equation of state [27] to obtain equilibrium bulk properties V0, B0, and its pressure derivatives .
3 Results and discussion
3.1 GGA and LDA calculations
First we give the stability of different magnetic ordering. Different magnetic solutions including AFI, FM, NM and AFII have been taken into account in the GGA calculation and the different magnetic orderings are illustrated in Fig.1. The total energies of per unit cell as a function of the lattice parameters are displayed in Fig.2. The calculations suggest that the full relaxation of the cell shape has a very slight effect on the total energy, which can be neglected within the accuracy of the calculation, so the total energy of different magnetic solution calculated using a cubic B1 structure are shown. The AFII ordering is correctly predicted to be the most stable magnetic configuration followed by FM and then the AFI, which is consistent with experiment and other calculations. It should be noticed that the energy difference between the AFII and the other magnetic orderings is much larger than the energy difference among other three magnetic orderings and the later is about the order of 0.1 eV. The energy differences
Fig.1 Schematic diagram of magnetic structure: the top and bottom panels represent AFI and AFII magnetic ordering, respectively. The gray (big) and black (small) represent Ni and O atoms, respectively. The arrows indicate the spin orientation.
Fig.2 Total energies of double f.u. as function of lattice parameter calculated using GGA. Magnetic solutions are marked in the figure
between different solutions are also very consistent with earlier calculations [28, 29], although our results are slightly larger. The calculated energy difference related to FM state is 0.028 6 eV/double f.u., -0.492 3 eV/double f.u. for AF1 and AFII, separately, compared to 0.011 9 eV/double f.u. and -0.033 1 eV/double f.u. in earlier HF calculations [28]. For AFII magnetic solution, the calculated band gap and magnetic moment are about 0.5 eV and 1.28μB. The calculations also predict other configurations to be metallic which are consistent with other GGA calculations with MnO[30].
The results of ground-state properties including equilibrium volume, bulk modulus and pressure derivative of the bulk properties after full relaxation of the cell shape are listed in Table 1. First we draw our attention to equilibrium volume, which is very important to equation of state. An improper equilibrium volume will lead to a large discrepancy of P-V curve at low pressure. Compared with LDA [6] which give a value of 0.035 22 nm3, 3.3% smaller than experimental one, GGA give a much improved equilibrium volume of 0.036 72 nm3, which is about 0.7% larger than experimental value of 0.036 44 nm3 [2]. Since the B0 is known to depend on the pressure range which is used in the fitting procedure both experimentally and theoretically. First we determine the bulk modulus from the pressure range of 0-60 GPa, as the same pressure range used in the earlier LDA calculations [6]. Compared with the previous LDA results of 236 GPa, our result of 200 GPa is very close to the experimental value of 203(2) GPa [2] in similar pressure volume range. Then we also calculate the bulk modulus in different volume (pressure) range. The bulk modulus in all three pressure ranges are all lower than values found in ETO’s experiment, but agree well with other experiments [3, 4]. In addition, the calculations on cubic rocksalt B1 structure are also performed. Because the lattice distortion is very small, relaxation has very slight effect on the bulk modulus and other properties.
Table 1 Calculated equilibrium volumes (V0), lattice para- meters a, c, axial ratio c/a, bulk modulus B0 and pressure derivative of the bulk modulus using GGA and GGA+U, together with the experimental value and earlier LDA calculation results
Fig.3 gives the volume of a hexagonal unit cell as a function of pressure, together with different experimental results in Refs.[2-4]. The present GGA results are in good agreement with experiments, while previous LDA gave a large discrepancy at low pressure. This may result from the fact that LDA predicts a too small lattice constant mentioned above. Although lattice distortion has a very slight effect on equilibrium properties, a small difference on P-V curves appears in the high pressure ranges.
Fig.3 Pressure dependence of hexagonal unit cell (three f.u.) together with experimental data taken from Refs.[2-4]
Our GGA calculations on AFII ordering predict the ground state to be slightly compressed along [111] direction and the calculated axial ratio c/a is 2.428 75 (rhombohedral angle is 60.38?), overestimating the experimental axial ratio 2.445 7(1) [2] (angle is 60.27 ?) at zero temperature [6] somewhat. Compared with LDA value of 2.42, our result is also more close to the experimental value. The GGA’s overestimated tendency of rhombohedral distortion also is present in other 3d transition-metal monoxides such as MnO [30].
The pressure dependence of the lattice parameters a, c and c/a are shown in Fig.4. For comparison, the experimental data from Ref.[2] and LDA results from Ref.[6] are also included. The lattice parameter c decreases monotonically with increasing pressure, while lattice parameter a exhibits nonmonotic behavior with increasing pressure: it decreases below about 120 GPa and increases above 120 GPa. From Fig.4(c), c/a shows a significant change above 60 GPa, which is all very similar with previous LDA results [6]. In order to check the validity of our results, we also performed LDA calculation and the results are almost identical with SASAKI’s result. On the other side, in ETO’s experiment the lattice parameters a and c decrease monotonically with increasing pressure, the pressure coefficient of c/a is almost constant over the experimental pressure. Moreover, Fig.4 also shows that both GGA and LDA predict too small axial ratio c/a especially above 60GPa and the deviation from experiment becomes more lager at high pressure. However, it’s very clearly seen that GGA gives a more close result of a, c and a/c with experiment at entire experiment pressure range especially at lower pressure. It’s expected that the strong correlation maybe become weak with pressure because of the since the increasing screening [31]. But the discrepancy between theory and experiment even above 100 GPa shows that strong electronic correlations still play a very important role in structural properties under high pressure.
In order to explore the possible structure phase transition from B1 to B2 structure, which has been observed in most of the monoxides with B1 structure, we calculated the total energy of NiO with B2 structure at different volumes. The distorted B1 structure undergoes a phase transition into the B2 structure at volume of 0.029 40 nm3. Fit the total energy curves to a MURNAGHAN equation of state, the obtained phase transition pressure is about 450 GPa, which is much larger than 318 GPa obtained using LDA [6].
Compared with LDA results, using GGA we obtain a much improved structural and bulk properties, but the main discrepancy between experiment and theory still exist. The structural distortion at ambient pressure is still overestimated. The behavior of c/a especially above 60
Fig.4 Pressure dependence of structure parameter a, c and c/a of distorted B1 structure for GGA, together with experimental data taken from Ref.[2] and LDA results from Ref.[6]
GPa doesn’t accord with experiment. These dis- crepancies demonstrate that in order to get a better known of structural distortion under high pressure, the electronic correlation effects must be included.
3.2 GGA+U calculations
The strong electronic correlations are included in form of DFT+U in this paper. ROHRBACH et al [32] have researched the bulk and surface properties in ambient pressure with the same code, who found that DFT+U calculation with U-J=5.3 eV can give reasonable values of band gap, magnetic moment and bulk modulus within experimental uncertainty. Following ROHRBACH’s result, we also check the bulk properties varying with U-J and find that a value of U-J=5.3 eV leads to both improved description of above bulk properties and lattice distortion at ambient pressure. The experimental axial ratio 2.445 7(1) reaches above U-J=4 eV. Meanwhile, as suggested by ROHRBACH, U-J=5.3 eV would lead to a better description of the electronic spectrum with experiment. So, U-J=5.3 eV are used in the present calculation to describe the structural distortion of NiO under high pressure.
The total energies curves with different magnetic solutions calculated using U-J=5.3 eV are shown in Fig.5. Compared with GGA, GGA+U predicts the same order of different magnetic solution and similar energy differences (AFI: 0.028 58 eV, AFII: -0.208 422 eV). And the energy difference between HF calculations still exists. This difference maybe source from the fact that
Fig.5 Total energies of double f.u. as function of lattice parameters calculated using GGA+U (U-J=5.3 eV). Magnetic solutions are marked in the figure
HF [28] exactly considers the exchange energy while the electronic correlations are completely neglected.
As shown in Table 1, the equilibrium volume calculated using GGA+U is slightly larger than GGA result, but the bulk modulus at different pressure are also smaller than GGA and experiment result. The differences seem to source from that two calculations produce different pressure derivative of the bulk modulus. And it’s interesting that the result has very slight effect to p-V relationship.
Fig.6 shows the p-V relationship obtained using U-J=5.3 eV together with two earlier experiment results [2-4]. Our GGA+U results are better consistent with experimental result compared with GGA. The relaxation has a smaller effect under the considering pressure range. The reason is that the p-V relationship is determined by the total energies and the structural distortion is too small to affect the total energy appreciably even at high pressure. Using GGA+U, the calculated bandgap and
Fig.6 Pressure dependence of hexagonal unit cell (three f.u.) calculated using GGA+U together with experimental data taken from Refs.[2, 6]
magnetic moment is 3.1 eV and 1.69μB respectively, which are consistent with other LDA+U calculations in literature.
Now let us focus on the result of structural distortion. The calculated axial ratio under ambient pressure is 2.4493 (rhombohedral angle about 60.01?), which is more close to the experiment. The pressure dependence of the lattice parameters a, c and c/a calculated using GGA+U with U-J=5.3 eV is shown in Fig.7. It can be seen that the GGA+U calculations give a proper qualitative description for the behavior of NiO under pressure. Under experimental pressure range, both a and c decrease monotonically with pressure without an increasing at high pressure as predicted by GGA calculation. The axial ratio c/a decreases with pressure and d(c/a)/dp is almost constant over the entire pressure region which is consistent with experiment. With increasing pressure, it’s also noticed that the nonmonotonic behavior of the extremely large of d(c/a)/dp also appears beyond the experimental pressure.
Fig.7 Pressure dependence of structure parameter a, c and c/a of distorted B1 structure for GGA+U (U-J=5.3 eV) together with experimental data taken from Ref.[2] and LDA results from Ref.[6]
3.3 Theory analysis
The deviation from the B1 structure results from the fact that the magnetic moments on the adjacent [111] planes are coupled antiferromagnetically to each other. With increasing pressure, the strength of the direct nearest-neighbor exchange interaction will increase because the interatomic distance decreases. The increasing exchange interaction will lead to the c/a deviating larger from the cubic B1 structure with pressure, which was observed in experiment. On the other hand, the ordinary DFT inaccurate prediction for structural properties under high pressure has been explained following SASAKI’s model [2, 6], where the rhombohedral distortion was analyzed by expanding the total energy with respect to the shear strain up to the four-order term and determining the coefficient of each term. E(ε,V)≈E0(V)+a(V)ε2+b(V)ε2+c(V)ε3+d(V)ε4+…, the higher-order can be neglected in this equation, where a(V) represents the magnetic distortion, while b(V)ε2 represents the elastic term which functions as the restoring force when b(V) has a positive value. b(V) can be divided into the band structure energy and electrostatic energy, and both increase with pressure. The rhombohedral distortion under pressure is governed mainly by the second-order term b(V) of ε. b(V) decreases with pressure when the sign of b(V) changes, which leads to the large lattice distortion. SASAKI suggested that both LDA and GGA underestimated bbs in entire pressure range. This underestimation will be a decreasing of b(V), i.e., the restoring forces are underestimation. It leads to a larger lattice distortion compared with experiment, as shown Fig.4. GGA seems to give a more proper description for band energy than LDA, but this improvement doesn’t affect the sign of the b(V) which leads to a suddenly change of c/a above 60 GPa. The further underlying mechanisms are being discussed and will be present elsewhere later.
4 Conclusions
1) In the present study, we have investigated the phase stability and structural distortion of NiO as a function of pressure using GGA and GGA+U. For supplement, the LDA are also used in the calculations. It could be shown that although GGA gives much improved structure and bulk properties compared with LDA, the main discrepancy between experiment and theory which was found in early theoretical calculation still exists. The behavior of c/a becomes extremely large under high pressure above 60 GPa, which hasn’t been observed in experiment. This discrepancy is due to the improper description of the electronic correlation between 3d electrons. When the strong electronic correlations are included in form of GGA+U, we obtain overall agreement of structural distortion at entirely pressure range.
2) The increasing distortion with pressure can be induced by increasing direct exchange interactions between Ni atoms on the adjacent [111] planes. Both LDA and GGA underestimate the band-structure energy, which leads to the underestimated restoring forces which results in a large deviation of structure distortion from experiment results. These results show that the strong electronic correlation between 3d electrons still play a very important role in structural and electronic properties under high pressure.
References
[1] NOGU?S J, SORT J, LANGLAIS V, SKUMRYEV V, SURI?ACH S, MU?OZ J S, BAR? M D. Exchange bias in nanostructures [J].Phys Rep, 2005, 422(3):65-117.
[2] ETO T, ENDO S, IMAI M, KATAYAMA Y, KIKEGAWA T. Crystal structure of NiO under high pressure [J]. Phys Rev B, 2000, 61(22):14984-14988.
[3] HUANG E. Compression behavior of NiO in a diamond cell [J]. High Press Res, 1995, 13(6): 307-319.
[4] NOGUCHI Y, UCHINO M, HIKOSAKA H, ATOU T, KUSABA K, FUKUOKA F, MASHIMO T, SYONO Y. Equation of state of NiO studied by shock compression [J]. J Phys Chem Solids, 1999, 60(4): 509-514.
[5] FANG Z, TERAKURA K, SAWADA H, MIYAZAKI T, SOLOVYEV I. Inverse versus normal NiAs structures as high- pressure phases of FeO and MnO[J]. Phys Rev Lett, 1998, 81(5): 1027-1030.
[6] SASAKI T. Lattice distortion of NiO under high pressure [J]. Phys Rev B, 1996, 54(14): R9581-R9584.
[7] KRESSE G, HAFNER J. Ab initio molecular dynamics for liquid metals [J]. Phys Rev B, 1993, 47(1): 558-561.
[8] KRESSE G, HAFNER J. Ab initio molecular-dynamics simulation of the liquid-metal–amorphous-semiconductor transition in germanium [J]. Phys Rev B, 1994, 49(20): 14251-14269.
[9] KRESSE G, FURTHM?LLER J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set [J]. Phys Rev B, 1996, 54(16): 11169-11186.
[10] KRESSE G, FURTHM?LLER J. Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set [J]. Comput Mater Sci, 1996,6(1): 15-50.
[11] KRESSE G, JOUBERT D. From ultrasoft pseudopotentials to the projector augmented-wave method [J]. Phys Rev B, 1999, 59(3): 1758-1775.
[12] BL?CHL P E. Projector augmented-wave method [J]. Phys Rev B, 1994, 50(24): 17953-17979.
[13] WOOD D M, ZUNGER A. A new method for diagonalising large matrices [J]. J Phys A, 1985, 18: 1343-1359.
[14] JOHNSON D D. Modified Broyden’s method for accelerating convergence in self-consistent calculations [J]. Phys Rev B, 1988, 38(18): 12807-12813.
[15] PULAY P. Convergence acceleration of iterative sequences. the case of SCF iteration [J]. Chem Phys Lett, 1980, 73(2): 393-398.
[16] PERDEW J P, CHEVARY J A, VOSKO S H, JACKSON K A, PEDERSON M R, SINGH D J, FIOLHAIS C. Atoms, molecules, solids, and surfaces: Applications of the generalized gradient approximation for exchange and correlation [J].Phys Rev B, 1992, 46(11): 6671-6687.
[17] VOSKO S H, WILK L, NUSAIR M. Accurate spin-dependent electron liquid correlation energies for local spin density calculations: a critical analysis [J]. Can J Phys, 1980, 58(8):1200-1211.
[18] ANISIMOV V I, ZAANEN J, ANDERSEN O K. Band theory and Mott insulators: Hubbard U instead of Stoner I [J]. Phys Rev B, 1991, 44(3): 943-954.
[19] DUDAREV S L, BOTTON G A, SAVRASOV S Y, HUMPHREYS C J, SUTTON A P. Electron-energy-loss spectra and the structural stability of nickel oxide: An LSDA+U study [J]. Phys Rev B, 1998, 57(3): 1505-1509.
[20] BENGONE O, ALOUANI M, BL?CHL P, HUGEL J. Implementation of the projector augmented-wave LDA+U method: Application to the electronic structure of NiO [J]. Phys Rev B, 2000, 62(24): 16392-16401.
[21] ROHRBACH A, HAFNER J, KRESSE G. Electronic correlation effects in transition-metal sulfides [J]. J Phys: Condens Matter, 2003, 15(6): 979-996.
[22] TRIMARCHI G, BINGGELI N. Structural and electronic properties of LaMnO3 under pressure:An ab initio LDA+U study [J]. Phys Rev B, 2005, 71(3): 035101-035109.
[23] MONKHORST H J, PACK J D. Special points for Brillouin-zone integrations [J]. Phys Rev B, 1976, 13(12): 5188-5192.
[24] JEPSEN O, ANDERSON O K. The electronic structure of h.c.p. Ytterbium [J]. Solid State Commun, 1971, 9(20): 1763-1767.
[25] METHFESSEL M, PAXTON A T. High-precision sampling for Brillouin-zone integration in metals [J]. Phys Rev B, 1989, 40(6): 3616-3621.
[26] BL?CHL P E, JEPSEN O, ANDERSEN O K. Improved tetrahedron method for Brillouin-zone integrations [J]. Phys Rev B, 1994, 49(23): 16223-16233.
[27] FU C L, HO K M. First-principles calculation of the equilibrium ground-state properties of transition metals: Applications to Nb and Mo [J]. Phys Rev B, 1983, 28(10): 5480-5486.
[28] TOWLER M D, ALLEN N L, HARRISON N M, SAUNDERS V R, MACKRODT W C, APR? E. Ab initio study of MnO and NiO [J]. Phys Rev B, 1994, 50(8): 5041-5054.
[29] K?DDERITZSCH D, HERGERT W, TEMMERMAN W M, SZOTEK Z, ERNST A, WINTER H. Exchange interactions in NiO and at the NiO(100) surface[J]. Phys Rev B, 2002, 66(6): 064434-064442.
[30] PASK J E, SINGH D J, MAZIN I I, HELLBERG C S, KORTUS J. Structural, electronic, and magnetic properties of MnO[J]. Phys Rev B, 2001, 64(2): 024403-024409.
[31] COHEN R E, MAZIN I I, ISAAK D G. Magnetic collapse in transition metal oxides at high pressure: Implications for the earth [J]. Science, 1997, 275(5300): 654-657.
[32] ROHRBACH A, HAFNER J, KRESSE G. Molecular adsorption on the surface of strongly correlated transition-metal oxides: A case study for CO/NiO (100) [J]. Phys Rev B, 2004, 69(7): 075413-075425.
(Edited by YUAN Sai-qian)
Foundation item: Project(KF0504) supported by Open Program of Key Laboratory of Advanced Materials & Rheological Properties, Ministry of Education
Corresponding author: TANG Bi-yu; Tel: +86-732-2371004; E-mail: tangbiyu@xtu.edu.cn