Electronic correlations at paramagnetic $(001)$ and $(110)$ NiO surfaces: Charge-transfer and Mott-Hubbard-type gaps at the surface and subsurface of $(110)$ NiO

We explore the interplay of electron-electron correlations and surface effects in the prototypical correlated insulating material, NiO. In particular, we compute the electronic structure, magnetic properties, and surface energies of the $(001)$ and $(110)$ surfaces of paramagnetic NiO using a fully charge self-consistent DFT+DMFT method. Our results reveal a complex interplay between electronic correlations and surface effects in NiO, with the electronic structure of the $(001)$ and $(110)$ NiO surfaces being significantly different from that in bulk NiO. We obtain a sizeable reduction of the band gap at the surface of NiO, which is most significant for the $(110)$ NiO surface. This suggests a higher catalytic activity of the $(110)$ NiO surface than that of the $(001)$ NiO one. Our results reveal a charge-transfer character of the $(001)$ and $(110)$ surfaces of NiO. Most notably, for the $(110)$ NiO surface we observe a remarkable electronic state characterized by an alternating charge-transfer and Mott-Hubbard character of the band gap in the surface and subsurface NiO layers, respectively. This novel form of electronic order stabilized by strong correlations is not driven by lattice reconstructions but of purely electronic origin. We notice the importance of orbital-differentiation of the Ni $e_g$ states to characterize the Mott-Hubbard insulating state of the $(001)$ and $(110)$ NiO surfaces. The unoccupied Ni $e_g$ surface states are seen to split from the lower edge of the conduction band to form strongly localized states in the fundamental gap of bulk NiO. Our results for the surface energies of the $(001)$ and $(110)$ NiO surfaces show that the $(001)$ facet of NiO has significantly lower energy. This implies that the relative stability of different surfaces, at least from a purely energetic point of view, does not depend on the presence or absence of magnetic order in NiO.


I. INTRODUCTION
The series of transition metal monoxides MnO, FeO, CoO, and NiO with an electronic configuration ranging from 3d 5 to 3d 8 , respectively, has attracted much attention due to their diverse electronic and magnetic properties 1-5 , allowing for a broad range of applications, e.g., in electronics and spintronics 6-9 , energy storage [10][11][12][13] , and heterogeneous catalysis [14][15][16][17] . At low temperature, these compounds exhibit a correlated Mott-Hubbard or charge-transfer insulating behavior with a large band gap of ∼2-4 eV, associated with a strong localization of the 3d electrons [18][19][20][21][22] . Below the Néel temperature, ranging from T N ∼ 116 to 523 K for MnO to NiO, respectively, these materials display an antiferromagnetic type-II long-range magnetic ordering and undergo a structural phase transition from a rocksalt cubic to a distorted rhombohedral (MnO and NiO) or monoclinic (FeO and CoO) phase 23 .
Over the past decades particular attention has been devoted to understanding the nature of the band gap and excitation spectrum of MnO-NiO [18][19][20][21][22][24][25][26][27][28][29][30] . In fact, due to the strongly correlated nature of electron interactions between the 3d electrons, theoretical computations of the electronic structure of these materials using band-structure methods are particularly difficult. The Coulomb interactions may be modelled by including an on-site Hubbard parameter U to treat the effect of correlations in the partially filled 3d shell, e.g., within the so-called DFT+U method [25][26][27] . While these "beyond standard density functional theory (DFT)" methods often give a reliable description of the electronic properties such as band gaps, magnetic moments, and lattice displacements [30][31][32][33][34][35] , these methods neglect electron dynamics and hence cannot capture correlated electron phenomena related to a Mott transition such as a coherenceincoherence crossover, quasiparticle behavior, and strong quasiparticle mass renormalization 18,19,[36][37][38] . In fact, to determine the electronic properties of these materials one needs to go beyond conventional band-structure methods to include dynamical correlation effects of the 3d electron (e.g., using DFT+U subject to dynamical symmetrybreaking spin and structural effects 28,29 ).
It has been shown that the dynamical correlations and strong localization of the d (and f ) electrons can be described by employing a DFT+dynamical mean-field theory (DFT+DMFT) approach [36][37][38] . DFT+DMFT makes it possible to treat local correlation effects in a self-consistent, numerically exact way, providing a good description of the electronic and magnetic properties of correlated models and materials  . By using DFT+DMFT it becomes possible to compute the material-specific properties of complex correlated materials, e.g., to determine the electronic structure, mag-netic state, and crystal structure of (para-) magnetic materials at finite temperatures, e.g., near the Mott transition 44,46,47,51,53,54,58,59 .
The DFT+DMFT approach has been widely used to study the physical properties of the bulk structure of transition metal monoxides (TMOs), providing a quantitative description of the electronic, magnetic, and structural properties of these materials [43][44][45][46][47][48][49][50][51][52][53][54][55][56][57][58] . While the electronic structure and magnetic state of bulk TMOs are relatively well established nowadays (see however [60][61][62] ), this is not the case for the different surface structures of metal oxides and their interfaces 1-5 . In fact, the electronic structure of a surface or interface can differ significantly from that of bulk, since their electronic states are responsible for a variety of novel (emerging) physical properties. In the case of surfaces, this affects, e.g., molecular adsorption and reactions that determine the catalytic processes at the surface 1-3 . In this respect, NiO is of particular interest among other TMOs due to its importance as catalysts for electrochemical applications, fuel cells, and batteries [10][11][12][13][14][15][16][17] (in addition to being a model system for understanding the Mott-Hubbard metal-insulator transition 18,19 ).
In this paper, we employ a multi-site extension of the DFT+DMFT method 41,42,83,84 implemented within a plane-wave pseudopotential formalism 51,53,54,59 to explore the effects of electronic correlations on the electronic properties of the nonpolar paramagnetic (PM) (001) and (110) surfaces of NiO. We study the electronic structure, magnetic properties, and surface energies of the (001) and (110) NiO with a particular attention given to the effect of structural confinement and its influence on the strength of electronic correlations in the PM (001) and (110) NiO.

II. COMPUTATIONAL DETAILS
We perform a DFT+DMFT study of the electronic structure, magnetic properties, and surface energies of the clean (001) and (110) surfaces of the prototypical correlated insulator NiO. In DFT we employ the generalized gradient Perdew-Burke-Ernzerhof approximation (GGA) 85 within plane-wave pseudopotential formalism 86 . Our DFT+DMFT calculations explicitly include the Ni 3d and O 2p valence states, by constructing a basis set of atomic-centered Wannier functions within the energy window spanned by the p-d band complex 87 . This makes it possible to take into account a charge transfer between the Ni 3d and O 2p states, accompanied by the strong on-site Coulomb correlations in the Ni 3d shell. In order to solve the realistic many-body problem within DMFT we use the continuous-time hybridizationexpansion quantum Monte-Carlo algorithm 88 .
The DFT+DMFT calculations are performed in the paramagnetic state at an electronic temperature T = 390 K. We use the average Hubbard U = 10 eV and Hund's exchange J = 1 eV, in accordance with previous work 25,53,55,89 . The Coulomb interaction U and Hund's J have been treated in the density-density approximation, with no spin-orbit coupling contribution. We employ the fully localized double-counting correction evaluated from the self-consistently determined local occupations, to account for the electronic interactions already described by DFT. We use a fully self-consistent in charge density implementation of DFT+DMFT in order to take into account the effects of charge redistribution caused by electronic correlations and electron-lattice coupling 53,90 . In order to study the effects of structural confinement and relaxations we employ a multi-impurity-site extension of the DFT+DMFT method in order to treat correlations in the 3d bands of the structurally distinct Ni sites 41,42,83,84 . In addition, we perform structural relaxations of the surface and subsurface states of NiO within DFT. The DFT structural relaxation calculations are performed for the (fictitious) ferromagnetic (FM) state of NiO. For the relaxed structures the electronic properties of (001) and (110) PM NiO were evaluated within DFT+DMFT. The spectral functions were computed using the maximum entropy method and the Padé analytical continuation procedure.

III. RESULTS AND DISCUSSION
a. Electronic structure of (001) NiO We begin with the electronic structure and equilibrium lattice volume calculations of bulk NiO in the paramagnetic state using DFT+DMFT 51,53 . In agreement with experiment, our calculations yield a correlated insulator with a large Ni d-d energy gap of about 2.9 eV and equilibrium lattice constant a = 4.233Å (lattice volume of 128 a.u. 3 and bulk modulus 187 GPa). The Ni 2+ d 8 ions are in a high-spin (S = 1) state with an instantaneous local magnetic moment of 1.83 µ B , corresponding to the fluctuating moment of 1.7 µ B . This result is in good agreement with the experimental estimate of 1.7-1.9 µ B 91 . We note that the top of the valence band of NiO has a mixed Ni e g and O 2p character, with a resonant peak in the filled Ni e g band located at about 1 eV below the Fermi level. The latter can be ascribed to the formation of a Zhang-Rice bound state 92 . Our results suggest a mixture of a Mott-Hubbard Ni d-d and charge transfer type of the band gap which is caused by the effect of Coulomb correlations in the Ni 3d shell 20,21,24 . The bottom of the conduction band is mainly of the Ni 3d character with a broad paraboliclike Ni 4s band at the Brillouin zone Γ point.
Next, we use the equilibrium lattice constant obtained by DFT+DMFT for bulk NiO to construct the (001) symmetric slab consisting of 7-NiO-monolayers (ML), with a thickness of about 12.59Å and vacuum of ∼21 A. Each ML in the slab contains two Ni sites. The top two NiO MLs are considered as a surface and subsurface, respectively, with 3-ML-thick NiO as a quasi-bulk. In our calculations we set the kinetic energy cutoff for the plane wave basis to 65 Ryd for the wave function and 650 Ryd for the charge density. We employ the three-impuritysite DFT+DMFT to treat structurally distinct the surface, subsurface, and quasi-bulk Ni 3d ions 93 . We first perform structural optimization of the surface and subsurface (001) NiO layers within FM GGA. It gives a small inward relaxation of the top ML of δ 12 −2.5% (0.054 A) and a weak outward relaxation of the subsurface ML δ 23 ∼ 0.06%, with nearly absent surface buckling. Here, surface relaxations are characterized as a percent change of the spacing between layers i and j versus the equilibrium inter-layer spacing d 0 , δ ij = (z Ni i − z Ni j )/d 0 . We also find a small inward inter-plane relaxation between the first and third MLs, δ 13 ∼ −1.0%. Interestingly, the nomagnetic GGA yields -2.7% and -2.4% relaxation for the top and subsurface NiO-MLs, respectively, with a small surface buckling of ∼0.11Å.
In Fig. 1 we show the orbitally-resolved spectral functions of (001) NiO obtained by DFT+DMFT for the relaxed within FM GGA surface of NiO. Our results for the unrelaxed (bulk-terminated) (001) NiO are shown in Supplementary Fig. S1. In accordance with experiment, we find a Mott-Hubbard insulator with a large (fundamental) energy gap of about 2.1 eV. The calculated local magnetic moments are 1.84, 1.82, and 1.82 µ B for the surface, subsurface, and quasi-bulk Ni ions, respectively. The corresponding fluctuating moments are 1.78, 1.74, and 1.75 µ B . Our results reveal a weak charge redistribution between the surface and bulk Ni sites. The total Wannier Ni 3d occupancies are nearly the same for all the Ni sites and differ by less than 0.03. In addition, we find no sizeable distinction between the Ni 3z 2 − r 2 and x 2 − y 2 orbital occupancies, which differ by less than 0.03 for all the Ni sites. At the same time, we observe a substantial splitting of the occupied O 2p and the empty Ni e g bands at the surface caused by the effect of surface symmetry breaking (see the top panel of Fig. 1). It leads to a sizeable reduction of the band gap from about 2.6 eV in the quasi-bulk to 2.1 eV at the surface. Moreover, surface relaxation (strain) is seen to result in a small increase of the band gap by less than ∼5%, from 2.0 to 2.1 eV. This result underlines the importance of the crystalfield splitting caused by the effect of surface truncation and relaxations for establishing the Mott insulating state of (001) NiO.
In the surface layer the top of the (001) NiO valence band shows a mixed Ni x 2 −y 2 and O 2p z character, with a large contribution from the O 2p states, caused by the surface symmetry breaking and strain effects. Moreover, in the subsurface ML the top of the valence band has a mixed Ni 3z 2 − r 2 and O 2p x /p y character. Our results therefore suggest a charge-transfer character of the (001) surface. Moreover, we also notice the importance of orbital-differentiation of the Ni e g states (for the surface and subsurface NiO-MLs) to characterize the correlated insulating state of (001) NiO. In addition to this, the unoccupied Ni 3z 2 − r 2 surface states are seen to split from the lower edge of the conduction band and form strongly localized (and dispersionless) band at about 2 eV above the E F (see our results for the k-resolved spectra of (001) NiO in Fig. 2). In contrast to this, at the top of the valence band, all bands merge with the continuum of bulk states.
We also calculated the local (dynamical) susceptibility χ(τ ) = m z (τ )m z (0) for the surface, subsurface, and quasi-bulk Ni 2+ ions of PM (001) NiO. In Fig. 3 we display our results for the N e g contributions in χ(τ ). τ denotes imaginary times. The value of this quantity around τ = β 2 indicates the long time limit, and deviations from this value for small and large β indicate strong magnetic fluctuations. In the present case, all the e g contributions are seen to be almost independent of τ , suggesting that the 3d electrons are localized to form fluctuating local moments. In fact, χ(τ ) is seen to be nearly constant and close to its maximal value S = 1 for the Ni e g states. However, in the surface layer χ(τ ) is seen to be remarkably larger at all τ , implying higher localization of the surface e g states. This leads to a higher charge-transfer character of the (001) NiO surface, while for the bulk our results suggest a mixture of a Mott-Hubbard d-d and charge transfer character of the band gap 20,21,24 .
For the (001) NiO we compute the surface energy within DFT+DMFT as γ (E N slab − N E bulk )/(2S) 94 . Here, E N slab is the total energy of an N -formula-unit slab of NiO. E bulk is the total energy of bulk NiO per formula unit. S is the surface area and the factor two accounts for the two surfaces in the slab unit-cell. In order to minimize numerical differences between the bulk and the slab DFT+DMFT calculations we adopt the same setup parameters within the DFT+DMFT total energy calculations (such as kinetic energy cutoffs, k-point sampling, etc.). In particular, we compute the 4-formula-unitsupercell to estimate E bulk for the (001) NiO. Our result 59 meV/Å 2 is in good agreement with previous estimates in the antiferromagnetically ordered phase 74,76,79 . We note that the nonmagnetic GGA yields 44 meV/Å 2 .

b. (110) NiO surface
We now turn to our results for the PM (110) NiO. To model (110) NiO we design the (110) 11-ML-thick slab with symmetry-equivalent surfaces. It is of 14.59 A thick, with a vacuum spacer of about 21Å. In our DFT+DMFT calculations, we treat the top three NiO MLs as surface, subsurface, and sub-subsurface layers, respectively. The remaining 5-ML-thick layer in the middle is taken as quasi-bulk NiO. We use the four-impuritysite extension of DFT+DMFT to treat the effect of correlations in the 3d shell of the structurally distinct Ni ions. In DFT the kinetic energy cutoff for the plane wave basis was set to 50 Ryd for the wave function and 500 Ryd for the charge density. All the calculations are performed in the local basis set determined by diagonalization of the corresponding Ni 3d occupation matrices. Structural optimization of the surface and two subsurface (110) NiO layers within FM GGA suggests that surface relaxations in (110) NiO are long-range. It gives a large inward relaxation of the top layer of about -11.0% (-0.16 A) accompanied by a sizable outward relaxation of the subsurface layer by 3.4%. The inter-plane relaxation between the first and third layer δ 13 is about -3.8%, while δ 15 ∼ −2.9%.
In Fig. 4 we summarize our results for the Ni 3d and O 2p spectral functions of PM (110) NiO. Our results for the k-resolved spectra are shown in Fig. 5. We obtain a Mott-Hubbard insulating solution with a relatively large band gap of 0.9 eV. Similarly to (001) NiO, we observe a weak variation of the local moments among the different Ni layers. In particular, the calculated local magnetic moments are 1.83, 1.79, 1.82, and 1.82 µ B for the surface, two subsurfaces, and quasi-bulk Ni-ions, respectively. Our estimate of the fluctuating moments is 1.77, 1.70, 1.73, and 1.75 µ B . The Ni 3d occupancies are nearly the same for all the Ni sites, differing by less than 0.04. In addition, we observe no sizeable difference between the Ni 3x 2 − r 2 and y 2 − z 2 orbital occupancies in (110) NiO. The deviations are less than 0.03. Our analysis of the layer-dependent Ni 3d and O 2p spectral functions reveals a remarkable reduction of the band gap at the surface layer. In fact, it drops by about three times, from 2.9 eV in the quasi-bulk to 0.9 eV in the surface layer, in the relaxed (110) NiO. Most importantly, our result for the (fundamental) gap value of (110) NiO is more than twice smaller than that in the (001) NiO (2.1 eV). In close similarity to (001) NiO, structural relaxations of the surface result in an increase of the band gap from 0.8 to 0.9 eV.
Our results for the Ni 3d and O 2p spectral functions note however that the situation significantly differs from that in the (001) NiO or in bulk NiO. In fact, in contrast to the bulk and (001) NiO we observe that the top of the valence bands of (110) NiO has almost pure O 2p character at the surface and in the sub-subsurface layer, implying a pure charge-transfer type of the band gap of the (110) surface of NiO. In contrast to this, at the subsurface layer, it is nearly purely of Ni 3d (y 2 − z 2 ) character, suggestive of a Mott-Hubbard type of the band gap in the subsurface layer of (110) NiO. Our results therefore document a novel electronic state characterized by an alternating pure charge-transfer and Mott-Hubbard character of the band gap in the surface and subsurface (110) NiO, respectively. This behavior suggests the absence of prominent Zhang-Rice physics at the (110) surface of NiO that reappear in bulk NiO 43,45,53 . We notice that this remarkable alternating charge-transfer and Mott-Hubbard state is robust with respect to the lattice reconstructions of the (110) surface as seen in the spectral function of the bulk-terminated (110) NiO (see Fig. S2). This suggests that it is of purely electronic origin. Moreover, we notice a remarkable orbital-differentiation of the Ni e g states at the (sub-) surface which seems to be important to characterize the transport properties of (110) NiO. In Fig. 6 we display our results for the layer-dependent orbitally-resolved local susceptibility χ(τ ). We observe that all the e g contributions are seen to be almost independent of τ , suggesting a strongly localized nature of the Ni e g states in (110) NiO. Similarly to (001) NiO the surface layer χ(τ ) is somewhat larger for all τ , implying a higher degree of localization of the surface e g states. On the other hand, counter-intuitive to this is the reduction of the band gap at the surface layer seen in the Ni 3d and O 2p spectral function for the (110) NiO surface, which is caused by charge-transfer effects.
We also calculate the surface energy of the relaxed (110) NiO using DFT+DMFT. In our calculations we use the 8-formula-unit supercell NiO as a reference to estimate E bulk . Our result for the surface energy of (110) NiO 125 meV/Å 2 is substantially larger than that for the (001) NiO (59 meV/Å 2 ), implying that (001) NiO has significantly lower energy facet. The nonmagnetic GGA calculation gives 85 meV/Å 2 for (110) NiO.

IV. CONCLUSIONS
In conclusion, we have calculated the electronic structure, magnetic properties, and surfaces energies of the (001) and (110) surfaces of the prototypical correlated insulating material, PM NiO, using a multi-site extension of the DFT+DMFT method. Our results reveal a complex interplay between electronic correlations and surface effects in NiO. We obtain that the electronic structure of the (001) and (110) NiO differs significantly from that of bulk NiO. In both (001) and (110) NiO we observe a sizable reduction of the band gap at the surface. The latter is most significant for the (110) NiO, from 2.9 in the quasi-bulk to 0.9 eV at the surface. Surface relaxations (surface strain) are seen to result in a remarkable increase of the band gap in comparison with the bulk-terminated NiO.
We observe a substantial splitting of the O 2p and Ni e g bands at the surface caused by the effect of surface symmetry breaking. Our results suggest a charge-transfer character of the (001) and (110) surfaces of NiO. Most importantly, for the (110) NiO surface we observe a remarkable electronic state characterized by an alternating pure charge-transfer and Mott-Hubbard character of the band gap in the surface and subsurface NiO layers, respectively. This novel form of electronic order stabilized by strong correlations is not driven by lattice reconstructions but of purely electronic origin. Moreover, we notice the importance of orbital-differentiation of the Ni e g states to characterize the Mott insulating state of the (001) and (110) NiO. The unoccupied Ni e g surface states are seen to split from the lower edge of the conduction band, to form strongly localized states in the fundamental gap of bulk NiO.
Our DFT+DMFT calculations reveal a remarkable difference between the electronic structure of the (001) and (110) NiO surfaces. For example, we obtain a large difference in their (fundamental) band gap value (by about a factor of two) due to the surface effects, ∼2.1 eV in the (001) and 0.9 eV for the (110) NiO surfaces. This suggests a higher catalytic activity of the (110) NiO surface than that of the (001) NiO surface. In agreement with previous estimates in the antiferromagnetically ordered phase, our DFT+DMFT calculations yield significantly different surface energies for the (001) and (110) NiO surfaces, of 59 and 125 meV/Å 2 , respectively. Thus, the (001) NiO surface is found to have a significantly lower energy facet. We find that the effect of electron-electron correlations in NiO results in a sizable enhancement of the NiO surface energies by 34-48% with respect to the nonmagnetic GGA. Overall, our results for the electronic structure, magnetic properties, and surface energies of (001) and (110) NiO agree well with experimental data measured in the antiferromagnetic phase. Our results further reveal an intriguing independence of properties such as the energetical hierarchy of different surfaces of the presence or absence of antiferromagnetic order.
In our study we demonstrate the ability of multi-site DFT+DMFT to compute the properties of transition metal oxides surfaces. The DFT+DMFT computations can be used to promote our understanding of the physics of surfaces, which are promising for applications. Further studies could for example deal with NiO deposited on SrTiO 3 , elucidating the increased reactivity with water of such systems 95 or the design of highly active catalysts, e.g., for developing hybrid nickel-zinc and zinc-air batteries 13 . NiO is also being actively investigated for its resistive switching properties [96][97][98] for the description of which extensions of DFT+DMFT might also prove useful.

Supplementary Material
In this supplementary material, we display additional information concerning the spectral properties of the investigated NiO surfaces.