Hellmann-Feynman Forces within the DFT+U in Wannier functions basis
aa r X i v : . [ c ond - m a t . s t r- e l ] J a n Hellmann–Feynman Forces within the DFT+U inWannier functions basis
D. Novoselov , Dm.M. Korotin , V.I. Anisimov , Institute of Metal Physics, S.Kovalevskoy St. 18, 620137 Yekaterinburg, Russia Ural Federal University, Mira St. 19, 620002 Yekaterinburg, RussiaE-mail: [email protected]
Abstract.
The most general way to describe localized atomic-like electronic statesin strongly correlated compounds is to use Wannier functions. In the present paperwe continue the development of widely-used DFT+U method onto Wannier functionbasis set and propose a technique to calculate a Hubbard contribution to an atomicforces. The technique was implemented as a part of plane-waves pseudopotential codeQuantum-ESPRESSO and tested on a charge transfer insulator NiO.PACS numbers: 71.15.-m
1. Introduction
Quantitative description of micro- and macroscopic properties for materials with strongelectron-electron interactions is a challenge for condensed matter physics during few lastdecades. The Density Functional Theory extensions such as DFT+U and DFT+DMFT,that allow to take into account electronic correlations, are under intensive developmentat the moment. Core of the DFT+correlations methods is a merging of the densityfunctional and a model approaches. Hubbard model for correlated states is build on aresult of ab-initio
DFT calculation. Some localized atomic wavefunctions are chosen todescribe correlated electrons. Then various Hubbard corrections for potential, chargedensity, total energy, atomic forces, etc. are added to the DFT calculation results. Animportant point in this approach is a choice of the localized atomic wavefunctions.Many researchers used Wannier function basis for strongly correlated materialscalculations in recent years, for example, see Refs. [1, 2, 3, 4]. This choice is convenientto describe a correlation effects in compounds with mixed character partially filled energybands due to the fact that it is localized in the direct space and is a superposition ofatomic orbitals of neighboring atoms. In our previous works [5, 6] we have formulated theDFT+U approach with Wannier functions (WF) for correlated states description. In thepresent paper we extend the approach with an atomic forces calculation technique. Thiswill allow to evaluate a phonon spectra and perform a molecular dynamics simulationsfor a correlated materials with complex energy bands structure near a Fermi level. ellmann–Feynman forces within Wannier functions basis
2. Method
The Hubbard contribution to the total energy could be expressed [8] as: E U [ n n Iσmn o ] = U X Imσ (cid:16) n Iσmm − X n n Iσmn n Iσnm (cid:17) , (1)where U is the on-site Coulomb interaction parameter and n Iσmn is the correlated statesoccupation matrix for site I and spin σ . By using Hellmann-Feynman theorem one canwrite a Hubbard contribution to an atomic forces as [9]: F Uαi = − ∂E U ∂τ αi = − X Imnσ ∂E U ∂n Iσmn ∂n Iσmn ∂τ αi , (2)where δτ αi is the displacement of the atom α of the unit cell in the i -th direction.Therefore the Hubbard contribution to the atomic forces is: F Uαi = − U X Imnσ ( δ mn − n Iσnm ) ∂n Iσmn ∂τ αi . (3)The key point in the atomic forces calculation is an accurate computation of theoccupation matrix derivative. The occupation number operator for k -point in reciprocalspace is: ˆ n σ k = X ν θ ( ε ν − E f ) | ψ σ k ν ih ψ σ k ν | , (4)where ǫ ν is the ν -th band energy, | ψ σ k ν i is ν -th eigenvector of the Hamiltonian matrixand E f is the Fermi energy. The occupation matrix element in WFs basis is ( k -pointweight is omitted here and below for simplicity): n Iσmn = X k h W Im k | ˆ n σ k | W In k i = X k ν θ ( ε ν − E f ) h W Im k | ˆ S | ψ σ k ν ih ψ σ k ν | ˆ S | W In k i . (5)WFs used here are defined as a projections of Bloch sums of the atomic orbitals | φ Im k i onto a subspace of the Bloch functions (the detailed description of WFs constructionprocedure within pseudopotential method is given in Ref. [10]): | ˜ W Im k i = N X µ = N | ψ σ k µ ih ψ σ k µ | ˆ S | φ Im k i , (6)where ˆ S is an overlap operator of the ultrasoft pseudopotential formalism [11] and ithas the following form:ˆ S = ˆ1 + X stI | β Is i q Ist h β It | , (7)where | β Is i is a projector function that satisfies the condition h β Is | ˜ φ It i r 3. Results and discussion The proposed approach was tested on one of typical objects for calculation schemesverification for correlated materials - nickel oxide. NiO is a charge-transfer insulatorwherein the partially filled bands are formed by Ni-d and O-p orbitals of neighboringatoms O. It is known that spin-polarized DFT+U calculation is able to successfullyreproduce bands structure of the compound [14]. Therefore we calculated an atomicforces acting on atoms at the end of self-consistent cycle within spin-polarized DFT+Uapproach in WF basis [5, 6].For the density-functional calculations, we used the Perdew–Burke–Ernzerhof GGAexchange–correlation functional together with Vanderbilt ultrasound pseudopotentials.We used a kinetic energy cutoff of 45 Ryd (360 Ryd) for the plane-wave expansion ofthe electronic states (core-augmentation charge). The self-consistent calculations wereperformed with a (6, 6, 6) Monkhorst–Pack k-point grid. Calculations were performedfor a cell containing two formula units. One-site effective Coulomb interaction parameterU = 8.0 eV for Ni 3 d states was chosen [14]. Constructed Wannier functions had asymmetry of Ni-d and O-p states. The total and partial densities of states for Ni ionsare shown on figure 1. Obtained energy gap values equals 4 eV and one is in agreementwith previous works [14] and photoelectron and XAS measurements [15],[16]. The upperHubbard band above the Fermi level consists of Ni-d states and the lower Hubbard bandis a mixture of Ni-d and O-p states due to a strong hybridization.To test the forces calculation procedure we displaced one nickel atom from itsequilibrium position in the x -direction (the displacement direction within a cell is shownas a blue-color vector on figure 2). Then the Hubbard contribution to the total forcewere computed using equations (3-13). ellmann–Feynman forces within Wannier functions basis -5-4-3-2-1 0 1 2 3 4 5 -5 -4 -3 -2 -1 0 1 2 3 D O S ( s t a t e s / e V / c e ll ) Energy (eV) total spin-uptotal spin-downNi1-d spin-upNi1-d spin-downNi2-d spin-upNi2-d spin-down Figure 1. (color online) NiO density of states. Upper and lower panels correspond tothe of the up and down projection of the spin moment. Shaded area corresponds tothe total density of states. Zero energy corresponds to Fermi level. The resulting total force, obtained within DFT+U in WF basis calculation, actingon the Ni atom is shown on figure 3 with blue-color dashed line. The force dependencewas obtained for 10 various displacement and then interpolated with a straight line. The”analytical” force is compared with a numerical one (red solid line on figure 3) obtainedas a numerical derivative on total energy dependence shown on the inset of the figure. Figure 2. (color online) Schematic NiO crystal structure view with displacement ofthe one Ni atom (dark-gray sphere) in the x -direction represented by blue-color vector.Light-gray and dark-gray spheres denoted Ni atoms with different spin. Red spherescorresponding an oxygen atoms. One can see that the atomic force obtained analytically are in a good agreementwith the numerical derivative. A slight difference in the lines slope could be explainedas a result of neglect of some terms in equation (13).The most clear and straightforward method for phonon frequencies calculation isthe ”frozen phonon” approach. The phonon frequency equals to the second derivativeof the total energy over atom displacement. On the other hand it could be computed asthe first derivative of the total force acting on atom. Since in both cases the derivative ellmann–Feynman forces within Wannier functions basis -0.04-0.03-0.02-0.01 0 0.01 0.02 0.03 0.04-0.2 -0.15 -0.1 -0.05 0 0.05 0.1 0.15 0.2 T o t a l f o r c e ( R y / B oh r) Atom displacement (Bohr) AnalyticalNumerical -234.902-234.899-234.896 -0.2-0.1 0 0.1 0.2 T o t a l ene r g y ( R y ) Atom displacement (Bohr) Figure 3. (color online) The dependence of total force and total energy (on inset) ondisplacement of one Ni atom in the x -direction. is computed numerically, the second approach is more accurate and it is used nowadaysmore intensively than the first one. From the data presented on the figure 3 we computedthe second derivative of the toal energy and the first derivative of the total force overatom displacement. The obtained values 0.21 Ry/(Bohr) and 0.2 Ry/(Bohr) are in agood agreement. Therefore the presented technique could be used not only for forcescalculation but additionally for vibrational properties computation within the ”frozenphonons” approach. 4. Conclusion In the present work the approach to calculate the Hubbard term contribution to anatomic forces within Wannier functions basis into DFT+U framework with ultrasoftpseudopotential formalism is proposed. We have performed a calculation of an atomicforce acting on slightly displaced nickel atom in NiO. The good agreement betweenatomic force evaluated analytically and numerically is obtained that confirms theapplicability and reliability of the proposed method. Acknowledgments The present work was supported by the grant of the Russian Scientific Foundation(project no. 14-22-00004). References [1] Leonov I, Binggeli N, Korotin D, Anisimov V I, Stoji´c N and Voll-hardt D 2008 Physical Review Letters http://link.aps.org/doi/10.1103/PhysRevLett.101.096405 [2] Pavarini E, Biermann S, Poteryaev A, Lichtenstein A I, Georges A and An-dersen O K 2004 Physical Review Letters http://link.aps.org/doi/10.1103/PhysRevLett.92.176403 ellmann–Feynman forces within Wannier functions basis [3] Yin W G, Volja D and Ku W 2006 Physical Review Letters http://link.aps.org/doi/10.1103/PhysRevLett.96.116405 [4] Lechermann F, Georges A, Poteryaev A, Biermann S, Posternak M, Yamasaki Aand Andersen O K 2006 Physical Review B http://link.aps.org/doi/10.1103/PhysRevB.74.125120 [5] Korotin D, Kukolev V, Kozhevnikov A V, Novoselov D and Anisimov V I 2012 Journal of Physics:Condensed Matter Preprint ) URL http://arxiv.org/abs/1201.6154 [6] Korotin D M, Novoselov D and Anisimov V I 2014 Journal of Physics: Condensed Matter Preprint )[7] Giannozzi P, Baroni S, Bonini N, Calandra M, Car R, Cavazzoni C, Ceresoli D, Chiarotti G L,Cococcioni M, Dabo I, Corso A D, Gironcoli S D, Fabris S, Fratesi G, Gebauer R, GerstmannU, Gougoussis C, Kokalj A, Lazzeri M, Martin-samos L, Marzari N, Mauri F, Mazzarello R,Paolini S, Pasquarello A, Paulatto L, Sbraccia C, Smogunov A and Umari P 2009 Journal ofPhysics: Condensed Matter Phys. Rev. B (3) 035105 URL http://link.aps.org/doi/10.1103/PhysRevB.71.035105 [9] Feynman R P 1939 Phys. Rev. (4) 340–343 URL http://link.aps.org/doi/10.1103/PhysRev.56.340 [10] Korotin D, Kozhevnikov A, Skornyakov S, Leonov I, Binggeli N, Anisimov V andTrimarchi G 2008 The European Physical Journal B http://dx.doi.org/10.1140/epjb/e2008-00326-3 [11] Vanderbilt D 1986 Applied optics [12] Martin R M 2004 Electronic Structure Basic Theory and Practical Methods (Cambridge UniversityPress) ISBN 9780521782852[13] Himmetoglu B, Floris A, de Gironcoli S and Cococcioni M 2014 International Journal of QuantumChemistry http://doi.wiley.com/10.1002/qua.24521 [14] Anisimov V I, Zaanen J and Andersen O K 1991 Physical Review B Phys. Rev. Lett. (20) 1951–1953 URL http://link.aps.org/doi/10.1103/PhysRevLett.53.1951 [16] Sawatzky G A and Allen J W 1984 Phys. Rev. Lett. (24) 2339–2342 URL(24) 2339–2342 URL