Competition between exchange-driven dimerization and magnetism in diamond(111)

TL;DR

This study investigates the competition between structural dimerization and magnetism in the flat surface states of diamond(111), revealing how electron interactions influence insulating and magnetic phases.

Contribution

It demonstrates that structural dimerization causes insulating behavior, but including spin polarization stabilizes a magnetic, undimerized ground state, highlighting the interplay of structure, electronic, and magnetic properties.

Findings

Dimerization opens an insulating gap in the absence of magnetism.

Including spin polarization stabilizes a magnetic, undimerized ground state.

Calculated band gaps are consistent with experimental optical gaps.

Abstract

Strong electron-electron interaction in ultraflat edge states can be responsible for correlated phases of matter, such as magnetism, charge density wave or superconductivity. Here we consider the diamond(111) surface that, after Pandey reconstruction, presents zig-zag carbon chains, generating a flat surface band. By performing full structural optimization with hybrid functionals and neglecting spin polarization, we find that a substantial dimerization ($0.090$ {\AA} / $0.076$ {\AA} bond disproportionation in the PBE0/HSE06) occurs on the chains; a structural effect absent in calculations based on the LDA/GGA functionals. This dimerization is the primary mechanism for the opening of an insulating gap in the absence of spin polarization. The single-particle direct gap is $1.7$ eV ($1.0$ eV) in the PBE0 (HSE06), comparable with the experimental optical gap of $1.47$ eV, and on the larger(smaller) side of the estimated experimental single particle gap window of 1.57-1.87 eV, after inclusion of excitonic effects. However, by including spin polarization in the calculation, we find that the exchange interaction stabilizes a different ground state, undimerized, with no net magnetization and ferrimagnetic along the Pandey $\pi$-chains with magnetic moments as large as $0.2-0.3~\mu_B$ in the PBE0. The direct single-particle band gap in the equal spin-channel is approximately $2.2$ eV ($1.5$ eV) with the PBE0 (HSE06) functional. Our work is relevant for systems with flat bands in general and wherever the interplay between structural, electronic and magnetic degrees of freedom is crucial, as in twisted bilayer graphene, IVB atoms on IVB(111) surfaces such as Pb/Si(111) or molecular crystals.

Tables

Table 1. Table 1: Bond lengths ( d 1 subscript 𝑑 1 d_{1} and d 2 subscript 𝑑 2 d_{2} ) of the atoms in the Pandey π 𝜋 \pi -chains in the topmost layer, dimerization Δ Δ \Delta , and buckling Δ ​ z Δ 𝑧 \Delta z for different exchange and correlation functionals (XC).

XC	$d_1$(Å)	$d_2$(Å)	$\Delta$	$\Delta z$ (Å)
PBE	1.440	1.440	0.000	0.0048
HSE06	1.476	1.400	0.026	0.0052
PBE0	1.483	1.393	0.031	0.0054
B3LYP	1.482	1.396	0.030	0.0048

Table 2. Table 2: The energy difference between the dimerized and undimerized structures and between magnetic and non-magnetic structures using different hybrid functionals. We use for the undimerized structure the PBE one, as the HSE06 and PBE0 functionals do not have a stable undimerized solution. We then obtain its energy in the HSE06 and PBE0 at fixed atomic positions. The magnetic structure is fully optimized and has no dimerization. The non-magnetic is also completely optimized and has dimerization.

$\Delta$ E (eV/cell)	HSE06	PBE0	B3LYP
dimerized $-$ undimerized	-0.042	-0.062	-0.022
magnetic $-$ non-magnetic	-0.007	-0.008	-0.002

Table 3. Table 3: For each electronic band structure (obtained with XC) calculated @ the atomic structure (relaxed with XC), the fundamental band gap E g subscript 𝐸 𝑔 E_{g} , and the direct band gap at the high-symmetry points of the Brillouin zone, J, K, Γ Γ \Gamma , J’, given in eV.

bands	@ structure	$E_g$	J	K	$\Gamma$	J’
PBE	@ PBE	0.000	0.031	0.128	4.392	5.279
HSE06	@ PBE	0.000	0.048	0.157	5.562	6.607
PBE	@ HSE06	0.121	0.560	0.545	4.438	5.322
HSE06	@ HSE06	0.532	1.006	0.994	5.606	6.650
PBE0	@ PBE0	1.194	1.672	1.670	6.339	7.381
B3LYP	@ B3LYP	0.961	1.421	1.407	5.925	7.032

Table 4. Table 4: The magnitude of the spin of each atom in units of the Bohr magneton μ B subscript 𝜇 𝐵 \mu_{B} using different exchange and correlation functionals (XC) using the unrelaxed (dimerized) structure is taken from the non-magnetic calculations, and with the ion-relaxed and cell-relaxed calculations. The labeling of the atoms can be matched to Fig. 3 .

XC	structure	1	2	3	4	5	6	7	8
PBE0	unrelaxed	-0.285(6)	+0.288(1)	+0.022(3)	-0.023(9)	-0.015(2)	+0.014(5)	+0.001(6)	-0.002(5)
PBE0	ion-relaxed	-0.387(4)	+0.390(5)	+0.030(4)	-0.032(5)	-0.020(6)	+0.019(7)	+0.002(2)	-0.003(3)
PBE0	cell-relaxed	-0.374(8)	+0.377(9)	+0.031(4)	-0.033(5)	-0.021(8)	+0.020(8)	+0.002(1)	-0.002(9)
HSE06	unrelaxed	-0.269(1)	+0.271(5)	+0.020(8)	-0.022(4)	-0.014(5)	+0.013(8)	+0.001(5)	-0.002(4)
HSE06	ion-relaxed	-0.348(9)	+0.351(8)	+0.027(1)	-0.029(1)	-0.018(7)	+0.017(9)	+0.001(9)	-0.003(1)

Table 5. Table 5: After inclusion of the magnetism, the bond lengths ( d 1 subscript 𝑑 1 d_{1} and d 2 subscript 𝑑 2 d_{2} ) of the atoms in the Pandey π 𝜋 \pi -chains in the topmost layer, dimerization Δ Δ \Delta , and buckling Δ ​ z Δ 𝑧 \Delta z for different exchange and correlation functionals (XC). Unrelaxed refers to the dimerized structure obtained with the non-magnetic calculations. Changes in the structure is also presented after relaxing the ions only, as well as relaxing the whole cell for the PBE0 functional.

XC	structure	$d_1$(Å)	$d_2$(Å)	$\Delta$	$\Delta z$ (Å)
PBE0	unrelaxed	1.482(9)	1.393(5)	0.031(1)	0.005(4)
PBE0	ion-relaxed	1.439(4)	1.438(9)	0.000(2)	0.006(5)
PBE0	cell-relaxed	1.419(8)	1.419(7)	0.000(1)	0.006(7)
HSE06	unrelaxed	1.475(5)	1.400(1)	0.026(2)	0.005(2)
HSE06	ion-relaxed	1.438(7)	1.438(5)	0.000(1)	0.006(3)

Table 6. Table 6: For the α 𝛼 \alpha (majority) and the β 𝛽 \beta (minority) bands, the fundamental band gap E g subscript 𝐸 𝑔 E_{g} , and the direct band gap at the high-symmetry points of the Brillouin zone, J, K, Γ Γ \Gamma , J’, given in eV, using the unrelaxed (dimerized) structure obtained from the non-magnetic calculations, as well as after ion and cell relaxation.

bands	structure	spin	$E_g$	J	K	$\Gamma$	J’
PBE0	unrelaxed	$\alpha$	1.756	2.183	2.303	6.333	7.387
PBE0	ion-relaxed	$\alpha$	1.789	2.191	2.381	6.332	7.401
PBE0	cell-relaxed	$\alpha$	1.711	2.174	2.363	6.154	7.424
HSE06	unrelaxed	$\alpha$	1.059	1.483	1.589	5.601	6.656
HSE06	ion-relaxed	$\alpha$	1.087	1.493	1.653	5.599	6.665
PBE0	unrelaxed	$\beta$	1.645	2.146	2.083	6.343	7.385
PBE0	ion-relaxed	$\beta$	1.638	2.135	2.088	6.345	7.398
PBE0	cell-relaxed	$\beta$	1.508	2.060	2.022	6.163	7.462
HSE06	unrelaxed	$\beta$	0.943	1.299	1.369	5.611	6.655
HSE06	ion-relaxed	$\beta$	0.939	1.433	1.373	5.612	6.665

Full text

I Introduction

The occurrence of strongly correlated states requires the dominance of the electron-electron interaction over the electronic kinetic energy. In the case of 3$d$ transition metal oxides or high $T_{c}$ superconductors, Mott insulating, magnetic, and superconducting states are stabilized via the strong localization of electrons in 3$d$ orbitals. Recently, it has been shown that a new class of strongly correlated systems can be achieved in ultraflat edge-states or surface bands having small Fermi velocities but not necessarily 3$d$ states, as it happens in twisted bilayer grapheneCao et al. (2018a, b), or multilayer graphene with rhombohedral stacking Lee et al. (2016); Pierucci et al. (2015); Henni et al. (2016); Henck et al. (2018); Pamuk et al. (2017). All these works expand the range of materials hosting strong correlation effects and exotic states of matter and point to the need of understanding exchange and correlation effects in flat edge-states.

One of the simplest and most studied systems hosting a flat edge-state prone to strong exchange-correlation effects is the diamond(111) surface – the structure of which is still under debate more than 100 years after Bragg got the Nobel prize and applied their diffraction technique to determine the structure of bulk diamond BRAGG and BRAGG (1913). The formation of the surface state can be understood by considering that in bulk diamond each carbon atom undergoes $sp^{3}$ hybridization and has four neighbors at distance of $\approx 1.54$ Å  and bond angle at 109.5o. The atoms on the (111) surface have, however, one missing bond and only three nearest neighbors. This dangling bond generates the so-called Pandey Pandey (1982) reconstruction resulting in a $2\times 1$ superstructure forming 1D zig-zag chains (see Fig. 1) and a surface electronic band. The weak, but not negligible, hopping integral in the direction parallel to the surface but perpendicular to the chains is responsible for the $\approx 0.5$ eV energy dispersion of the band. Even if this surface state is not as flat as the one detected in twisted bilayer grapheneCao et al. (2018a, b) or in multilayer graphene with rhombohedral stacking Lee et al. (2016); Pierucci et al. (2015); Henni et al. (2016); Henck et al. (2018); Pamuk et al. (2017), it is substantially more extended in reciprocal space and it holds a larger number of electrons. For this reason diamond(111) should be prone to strong exchange and correlation effects.

Despite its apparent simplicity, the theoretical and experimental description of the structural and ground state properties of the diamond(111) surface has proven to be an exceptionally difficult and yet unsolved problem.

Although it is well accepted (both in theory and experiments) that the surface is a Pandey $\pi$-chain Pandey (1982), the microscopic details, such as dimerization, buckling, deeper layer distortions are still under debate. While X-ray diffraction Huisman et al. (1998a) and ion scattering Huisman et al. (1998b) data suggest buckling of the surface atoms, low-energy electron diffraction measurements Walter et al. (2002) show that the buckling is negligible, but the dimerization is inconclusive within the experimental error. For the electronic structure, angle-resolved photoemission spectroscopy (ARPES) Graupner et al. (1997) measurements find an insulating state with the occurrence of the flat surface state 0.5 eV below the Fermi level. Electron energy loss spectroscopy Pepper (1982) measurements suggest a band gap of $\sim 1$ eV, and reflectance anisotropy spectroscopy Bussetti et al. (2007, 2009) gives larger values of the optical band gap of $\approx 1.47$ eV.

From the theoretical point of view, the first density functional theory (DFT) calculations within the standard local density (LDA) and generalized gradient (GGA) approximations gave conflicting results for the surface chains Vanderbilt and Louie (1984); Iarlori et al. (1992); Alfonso et al. (1995). This disagreement is most likely explained by the fact that calculations were very heavy for the time as 12 carbon layers are required for convergence. More recent calculations with the LDA and GGA functionals Kern et al. (1996); Scholze et al. (1996); Bechstedt et al. (2001); Stekolnikov et al. (2002); Marsili et al. (2005, 2007), do not present any buckling or dimerization of the structure. Similarly a calculation with the B3LYP functional also does not show buckling or dimerization Pierre et al. (2014). Consequently, the surface state appears to be metallic within the standard LDA/GGA approximations, in disagreement with all experimental dataGraupner et al. (1997); Scholze et al. (1996); Pepper (1982); Bussetti et al. (2007, 2009). Pioneering works by Marsili et al. Marsili et al. (2005, 2007, 2008) show that quasiparticle GW calculations on top of the GGA crystal structure lead to a gap opening only within the self-consistent $G_{1}W_{1}$ scheme and starting from an artificial band occupation. However, GW calculations were performed with very coarse grids (5 or 9 k-points in the Brillouin zone) that tend to substantially overestimate the band gap (see Appendix A), and in the absence of spin polarization. In these works, the band gap opening is attributed only to an electronic mechanism at fixed ionic coordinates. Finally, the magnitude of excitonic effects has been evaluated to be of the order of $0.1-0.4$ eV Marsili et al. (2008) leading to an experimental single particle gap of the order of $1.57-1.87$ eV, when added to the experimental gap of 1.47 eV.

The main problem of all previous theoretical works is that they rely on the GGA minimized structure, mostly because structural optimization within GW is not possible for solids and non-local exchange calculations in a plane-wave framework are too expensive for such a large system (24 atoms per cell and very dense electronic momentum k-point mesh). Moreover, all calculations neglected magnetism. In this work, we circumvent the difficulties of geometrical optimization with a dense mesh of electronic momentum k-points even in the existence of Hartree-Fock exchange by using a combination of plane wavesGiannozzi et al. (2009, 2017) and Gaussian basis sets Dovesi et al. (2014) that allow for fast structural optimization. We use hybrid functionals with exact exchange and range separation to understand the effects of the exchange interaction on the geometry and electronic structure of diamond(111). Finally, we also explore the occurrence of magnetic solutions.

The structure of the paper is the following: after explaining technical details in sec. II, we present results for non-magnetic (sec. III.1) and magnetic (sec. III.2) calculations.

II Technical details

DFT calculations are performed using the Quantum ESPRESSO Giannozzi et al. (2009, 2017) and CRYSTAL codesDovesi et al. (2014). We use the triple-$\zeta$-polarized Gaussian type basis sets for the C atoms Peintinger et al. (2013), with the PBE Perdew et al. (1996), PBE0 Adamo and Barone (1999), and HSE06 Krukau et al. (2006) functionals. The surface states require an ultradense sampling with an electronic momentum k-point mesh of $90\times 120\times 1$, crucial for an accurate determination of the band gap (as shown in Appendix A) and for the stabilization of magnetism. We used real space integration tolerances of 7-7-7-15-30, and an energy tolerance of $10^{-10}$ Ha for the total energy convergence. Fermi-Dirac smearing for the occupation of the electronic states of 0.001 Ha is used for all of the calculations. In the magnetic case, we further increase the energy tolerance to $10^{-11}$ Ha. We fix the magnetic state in the first iteration of the self-consistent cycle and then we release this constraint.

III Results

III.1 Non-magnetic calculations

Fig. 1 shows the diamond structure along the [111] direction in a $2\times 2\times 1$ cell. We consider 12 carbon layers with the bottom layer saturated by hydrogen and the top unhydrogenated. We choose an in-plane lattice parameter of $a=4.369$ Å, and $b=2.522$ Å, as derived from the experimental lattice constant of bulk diamond $a_{0}=3.567$ Å  Huisman et al. (1998a), but we also perform full structural optimization (cell and internal coordinates) although the results are weakly affected. A vacuum of 50 Å  is placed between the periodic images along the $z$-direction. We first optimize the structure within the PBE using both plane waves and Gaussian basis sets, finding practically indistinguishable results. We then use the PBE optimized structure as a starting guess for geometrical optimization with hybrid functionals. In the Gaussian basis set calculations we do not use any symmetry so that no a priori guess on the crystal structure is retained.

The Pandey $\pi$-chain surface atomic structure can be parametrized by the dimerization $\Delta$ and the buckling $\Delta z$ (see Fig. 1),where $d_{1}$ and $d_{2}$ label the two distinct bond lengths of the atoms that make up the zigzag chain on the topmost layer. The dimerization is then defined as: $\Delta=|d_{1}-d_{2}|/(d_{1}+d_{2})$. The buckling of the atoms, $\Delta z$, on this layer is simply the difference in their position along the $z$-direction (see Fig. 1).

Table 1 shows the calculated values for the topmost layer for the optimized atomic structures. With all three functionals, the buckling of the top layer is estimated to be very small, with $\Delta z\sim 0.005$ Å, which agrees well with the low-energy electron diffraction measurements of about 0.01 Å  Walter et al. (2002). With the PBE functional, the two bond lengths are equal, $d_{1}=d_{2}$, hence there is no dimerization. The inclusion of unscreened exchange with the PBE0 functional, or of screened exchange via the HSE06 functional, gives a significant imbalance between the two bond lengths, predicting a dimerization of the surface structure with $\Delta\sim 0.026$ in the HSE06 case. The dimerization is slightly larger with the PBE0 functional than with the HSE06 functional. In general, we find that the larger the amount of Hartree-Fock exchange included in the calculation and the more unscreended the exchange, the larger the dimerization. The energy gain induced by the dimerization is substantial, as shown in Table 2.

In a previous calculation using similar settings and the B3LYP functional, the author of Ref. Pierre et al., 2014 found no dimerization in the carbon chains. We repeated this calculation starting (i) from the undimerized PBE structure and (ii) from the dimerized HSE06 structure. In the first case, the simulation remains in the undimerized structure as in Ref. Pierre et al., 2014. However, in the second case, the structural optimization with the B3LYP functional converges to a dimerized structure that is lower in energy than the undimerized one, on the same line of what has been obtained with the other hybrid functionals.

Having demonstrated the crucial effect of the exchange interaction on the atomic structure, we now investigate its effect on the electronic spectrum. We first consider the PBE approximation on top of the PBE optimized structure (labeled PBE@PBE). We find, in agreement with all previous calculations Kern et al. (1996); Scholze et al. (1996); Bechstedt et al. (2001); Stekolnikov et al. (2002); Marsili et al. (2005, 2007), a metallic solution with no gap and practically very small direct gaps, as shown in Fig. 2. Interestingly, the use of the HSE06 on top of the PBE crystal structure (HSE06@PBE) still leads to a metallic solution with no indirect gap and tiny direct gaps at J and K. On the contrary, if the PBE functional is used on top of the HSE06 geometry (PBE@HSE06), a gap opens and the electronic structure is insulating with an indirect gap of 0.121 eV and direct gaps at J and K of 0.560 eV and 0.545 eV, respectively. The fact that HSE06 on top of the PBE structure leads to a metallic insulating solution, while even the PBE functional on top of the dimerized HSE06 geometry is successful in inducing an insulating state, demonstrates unambiguously that, in the absence of spin polarization, gap opening is mostly driven by the dimerization of the Pandey $\pi$-chains. Contrary to all previous works that tried to stabilize an insulating solution at fixed atomic coordinatesMarsili et al. (2005, 2007, 2008), our work underlines the crucial importance of the atomic distortion. The complete HSE06 calculation (i.e. HSE06 on top of the HSE06 structure, HSE06@HSE06) leads to a larger direct gap of $\approx 1$ eV (see Table 3) both at K and J and to a fundamental indirect gap of $\approx 0.532$ eV. The electronic bands of the full Brillouin zone as well as a comparison with the ARPES data are given in Appendix B. A larger direct gap of $1.7$ eV can be obtained using unscreened functionals such as PBE0, as shown Fig. 2. The direct band gaps at different high-symmetry points and for all the used approximations are reported in Table 3.

As hybrid functionals give only a slight underestimation of the experimental band gap in bulk diamondHummer et al. (2009); Garza and Scuseria (2016), the $32\%$ ($0.47$ eV) underestimation in HSE06 of the optical gap with respect to the experiments Bussetti et al. (2007, 2009) is fairly surprising. The situation is much better in the PBE0 leading to a somewhat larger gap than the optical direct gap measured in experiments ($0.1$ eV larger). However, this value is in better agreement with experiments; if we add the excitonic effects of $0.1-0.4$ eV Marsili et al. (2008) to the experimental gap of 1.47 eV, this leads to an experimental single particle gap of $1.57-1.87$ eV. The PBE0 value is thus on the lower side of the window ($1.7$ eV). However, all calculations presented up to now have been carried out neglecting spin-polarization. The occurrence of a flat band could also lead to magnetic solutionsCao et al. (2018a, b), as it happens in multilayer graphene with rhombohedral stackingHenck et al. (2018); Pamuk et al. (2017), a very similar system. For this reason, we investigate below the occurrence of magnetism using the HSE06 and PBE0 functionals.

III.2 Magnetic calculations

We first perform magnetic calculations at atomic coordinates fixed at the dimerized solution, referred to as “unrelaxed” in the remainder of the figures and tables. We choose as initial condition of the simulation a fully ferromagnetic (FM) and antiferromagnetic (AFM) configuration on the two surface atoms of the topmost layer along the zig-zag chain. As expected, the PBE functional does not stabilize a magnetic state, therefore we focus on the PBE0 and HSE06 functionals. In both cases, we find that the most stable solution has global zero magnetic moment, and is ferrimagnetic within each layer of atoms with large atomic magnetic moments. With the HSE06, the magnetic moments on the two atoms in the chain are $+0.271~{}\mu_{B}$ and $-0.269~{}\mu_{B}$. A similar solution is obtained with the PBE0, namely we obtain an atomic spin of $+0.288~{}\mu_{B}$ and $-0.285~{}\mu_{B}$ left on the atoms. The small imbalance of $0.003~{}\mu_{B}$ between the majority and minority spin electrons is partially linked to the dimerization in the atomic structure and in great part to the structure of deeper layers making the two atoms in the chain inequivalent with respect to deeper layers, as shown in Fig. 3. The magnetic moments decrease of about one order of magnitude between every two layers, going towards bulk diamond, as shown in Table 4.

As at the end of the magnetic simulation one of the atoms in the chain experience a restoring force towards the non-dimerized solution, we perform structural optimization in the presence of the magnetic solution. We find that the spin-polarized structural optimization restores the non-dimerized solution with a negligible dimerization remaning on the surface atoms, as shown in Table 5. Similarly, the dimerization of the atoms on the second layer, i.e. atoms 4 and 5, also becomes negligible, while the bond lengths remain unchanged in the deeper layers. We have also checked the effect of the cell relaxation with the PBE0 functional. The in-plane lattice parameter decreases to $a=4.289$ Å  from the experimental value of $a=4.369$ Å. However, the conclusions after the ionic relaxation remain, that the dimerization of the surface atoms is still negligible.

Our calculations show that even in the absence of dimerization, within hybrid functionals the ground state is a zero magnetization state, weakly ferrimagnetic on the two surface atoms, as shown in Table 4. After the ionic relaxation, the magnetic moment of the atoms changes by $\sim 0.1~{}\mu_{B}$, and the small imbalance on the magnetic moments of the surface atoms remains to be $\sim 0.002-0.003~{}\mu_{B}$. We have further checked the effect of the cell relaxation on the magnetic moments using the PBE0 functional. While the magnetic moment of the surface atoms decreased by $\sim 0.01~{}\mu_{B}$ with respect to the ion-relaxed calculation, the conclusion that the ferrimagnetism of the top surface atoms is stabilized remains.

The weak ferrimagnetism of the two surface atoms is visible also in the electronic structure as it splits the degeneracy of the spin bands in $\alpha$ (majority) and $\beta$ (minority) spin bands. Table 6 shows the band gap values at the high-symmetry points of the Brillouin zone with the magnetic calculations. The optical direct band gap between equal spin states is now $\approx 2.174(2.060)$ eV for majority(minority) spins at the J-point with the PBE0 functional. Finally, if we add the excitonic effects of $0.1-0.4$ eV Marsili et al. (2008) to the experimental gap of 1.47 eV, this leads to an experimental single particle gap of $1.57-1.87$ eV. The PBE0 gap value with the magnetism is thus on the higher side of the experimental window. With the HSE06 functional, the direct band gap at the J-point is lower 1.493(1.433) eV for majority (minority) spin electrons, hence on the lower side of the experimental window.

As the magnetic solution is a mean field solution, it is in principle possible that the full many-body multi-determinant ground state is actually non-magnetic, particularly in low dimensional materials. A many-body multi-determinant calculation with inclusion of structural optimization is not possible for this system. Thus, in order to test the reliability of hybrid functionals in predicting the competition between dimerization and magnetism, we consider carbyne, the linear carbon chain that is a prototype of dimerization in one dimension. This system is very pathological for what concerns structural and magnetic instabilities and is known to be non-magnetic, but dimerized. Our calculations show that the HSE06 functional favors the polyyne (dimerized carbyne) structure rather than cumulene (non-dimerized carbyne) structure – unlike the standard LDA/PBE functionals. This is in agreement with the literature Liu et al. (2013). Furthermore, we have started with an initial AFM configuration on the two atoms and found that the final state is non-magnetic. Therefore, we show that in a similar carbon-based system, magnetism is not stabilized even when we use a hybrid functional. Hence the magnetism in diamond(111) surface can be a physical effect and our calculations predict that the magnetic solution is the most stable with a slight energy difference of a few meV/cell from the non-magnetic solution as shown in Table 2. Further experiments are needed to verify this hypothesis.

Our work demonstrates that, within hybrid functionals, the ground state of the diamond(111) surface is then insulating with zero net magnetization and ferrimagnetic order along the top surface atoms of the Pandey $\pi$-chains; a very surprising result given that diamond is non-magnetic and the atomic orbitals forming the surface state are of $p$ character and thus, at the atomic level, not as localized as 3$d$ orbitals.

IV Conclusion

We have studied the diamond(111) surface by using hybrid-functionals with different degrees of screened exchange. Contrary to all previous theoretical works, we include the exchange interaction at all levels in the calculation, both in the structural optimization and in the calculation of electronic and properties. Moreover, we allowed for magnetism in calculations.

In the absence of spin polarization, the primary effect responsible for the gap opening is the dimerization of the Pandey $\pi$-chains, that is enough to lead to an insulating state. This is at odds with all previous spinless calculations that were either finding no gap Kern et al. (1996); Scholze et al. (1996); Bechstedt et al. (2001); Stekolnikov et al. (2002), or claimed that gap opening was purely an electronic mechanism Marsili et al. (2005, 2007, 2008). The PBE0 band gap of 1.672 eV is on the higher side, while the HSE06 band gap of 1.006 eV is on the lower side of the experimental window for single particle gap of 1.57-1.87 eV obtained by summing excitonic effects to the experimental gap. Thus, in the absence of spin polarization the system could be classified as a Peierls-Slater insulator.

By including spin polarization, we find that the flatness of the diamond(111) edge-state stabilizes an insulating state with zero net magnetic moment and with ferrimagnetic ordering along the top surface atoms of the chain with sizable magnetic moments of the order of $0.2-0.3~{}\mu_{B}$. As the magnetic moment depends weakly on the underlying crystal structure, the electronic structure depends weakly on the amount of dimerization. Interestingly, structural optimization in the presence of magnetism converges to a ground state with a negligible dimerization on the surface atoms. We find that the PBE0 gap of $\approx 2.174(2.060)$ eV is on the higher side, while the HSE06 gap of 1.493(1.433) eV for majority (minority) spin electrons is on the lower side of the experimental window. Thus, within a hybrid functional approach the ground state is essentially antiferromagnetic with negligible dimerization, i.e. a Slater insulator.

As diamond(111) can be seen as formed from buckled graphene layers with rhombohedral (ABC) stacking (see Fig. 1), it is instructive to compare our magnetic state with the one recently detected in multilayer graphene with ABC stacking Lee et al. (2016); Pierucci et al. (2015); Henni et al. (2016); Henck et al. (2018); Pamuk et al. (2017). In ABC graphene multilayers the state is globally antiferromagnetic, but weakly ferrimagnetic on the outer layers, exactly as in the present case. However, the magnetic moment per carbon atom is much smaller than in diamond(111). It is, however, important to recall that the flat edge states in ABC graphene extends in an extremely small part of the Brillouin zone and hosts less electrons than the flat band of diamond(111). The similarity of these two states suggests the occurrence of a magnetic state even in diamond(111).

The fundamental point underlined by our work is that there are two competing mechanisms for opening of a gap in diamond(111), namely magnetization or dimerization. Experimentally, it would be possible to detect the occurrence of magnetism via spin-resolved scanning tunneling spectroscopy by using magnetic tips.

Finally, our work demonstrates that in order to describe the correlated states in flat edge bands, it is necessary to include the electron-electron interaction at all levels in the calculations, both in the structural and electronic properties. This is relevant far beyond the case of diamond(111), and it is most likely also crucial to describe the phase diagram of twisted bilayer grapheneCao et al. (2018a, b), or other low dimensional system presenting edge states such as IVB atoms on top of IVB(111) surfaces such as Pb/Si(111)Tresca et al. (2018) or one dimensional polyenesGu et al. (2017).

Acknowledgements.

This work is supported by the Graphene Flagship, by RhomboG grant (ANR-17-CE24-0030) from Agence Nationale de la Recherche and by PRACE. Calculations were performed at IDRIS, CINES, CEA, BSC TGCC, and CAC(Tardis). B. P. acknowledges National Science Foundation (Platform for the Accelerated Realization, Analysis, and Discovery of Interface Materials (PARADIM)) under Cooperative Agreement No. DMR-1539918.

Appendix A Band gap convergence with respect to k-points

We test the convergence of the band gap with respect to the electronic momentum k-point mesh shown in Fig. 4. For this purpose, we perform electronic structure calculations with the PBE0 functional on top of the undimerized structure obtained with the PBE functional (PBE0@PBE). The band gap can be similarly compared with the second panel of Fig. 2 (HSE06@PBE) of the manuscript. Our tests show that the band gap strongly depends on the k-point grid. A coarse k-point mesh of $4\times 6\times 1$ clearly overestimates the band gap with respect to the denser k-point meshes. This mesh is comparable with that used in previous GW calculationsMarsili et al. (2005, 2007).

Appendix B Comparison with ARPES data

We calculate the electronic band structure using the HSE06 and PBE0 functionals. The electronic structure along the high-symmetry points of the Brillouin zone is given in Fig. 5 and Fig. 6, respectively. The valence band maximum is at the K-point and the conduction band minimum is at the J-point. The experimental data are taken from the ARPES measurements Graupner et al. (1997). We calculate the full path both for the non-magnetic and magnetic calculations with ionic relaxations.

The slope of the calculated path from the K-point towards the $\Gamma$-point depends on the amount of exchange and on the range of the interaction. For the PBE0 functional the best agreement is found with inclusion of magnetism.