Free surfaces recast superconductivity in few-monolayer MgB2: Combined first-principles and ARPES demonstration

Two-dimensional materials are known to harbour properties very different from those of their bulk counterparts. Recent years have seen the rise of atomically thin superconductors, with a caveat that superconductivity is strongly depleted unless enhanced by specific substrates, intercalants or adatoms. Surprisingly, the role in superconductivity of electronic states originating from simple free surfaces of two-dimensional materials has remained elusive to date. Here, based on first-principles calculations, anisotropic Eliashberg theory, and angle-resolved photoemission spectroscopy (ARPES), we show that surface states in few-monolayer MgB2 make a major contribution to the superconducting gap spectrum and density of states, clearly distinct from the widely known, bulk-like σ- and π-gaps. As a proof of principle, we predict and measure the gap opening on the magnesium-based surface band up to a critical temperature as high as ~30 K for merely six monolayers thick MgB2. These findings establish free surfaces as an unavoidable ingredient in understanding and further tailoring of superconductivity in atomically thin materials.

. Surprisingly, the exploitation of the electronic states of the free surfaces of two-dimensional superconductors has remained an almost completely unexplored terrain.
A particularly interesting candidate to examine the effect of thinning on superconductivity is hexagonal magnesium diboride (MgB 2 ), owing to its conveniently layered, graphene-related structure as well as to its high critical temperature of 39 K in bulk form. MgB 2 consists of alternating Mg-and B-planes, with stronger in-plane than out-of-plane electronic bonds. As a result, the electronic structure and superconducting properties of bulk MgB 2 are strongly anisotropic, instigating the formation of two superconducting gaps, each on a different part of the Fermi surface [36][37][38][39][40][41] . A first-principles approach to the Eliashberg formalism has shown that the stronger gap, (0) 7 ∆ ∼ σ meV, in bulk MgB 2 stems from the in-plane σ-bond, while the weaker gap, ∆ ∼ − π (0) 2 3 meV, is connected to the out-of-plane π-bonds [40][41][42] . Superconductivity is fully governed by B-atoms in bulk MgB 2 [36][37][38]41 , since both the σand π-gaps originate from electrons in B-planes that couple to lattice vibrations of the B-atoms. MgB 2 and other multigap superconductors are known to harbour rich new physics [43][44][45][46][47] , as a result of coupling, competition and phase frustration between different band condensates. Owing to the fact that σ-bands lie in-plane and π-bands out-of-plane, the layered, anisotropic crystal structure gives rise to innate multigap superconductivity in MgB 2 , while also enabling growing the material at atomic thickness, layer by layer. It is therefore of particular interest to understand how superconductivity would change in the thinnest limit. Hints of additional superconducting gaps opening in 2D materials have been suggested in recent works, for materials such as two-dimensional Ga and FeSe 48 . They are deduced from fitting of the measured temperature dependence of, e.g., the London penetration depth, which unfortunately does not provide much insight into the origin of the proposed additional gaps. The fitting technique is moreover not fully waterproof, as the measurements can equally result from just anisotropy of the superconducting gap, and no additional gap opening, as shown and discussed in ref. 47 .
Here, we discuss the nucleation of multigap superconductivity in atomically thin MgB 2 , and the crucial role of the electronic states of the free surface therein. The latter remained completely unexplored to date, as the previous studies mentioned above did not consider the role of surface states originating from free surfaces. Our analysis starts from a first-principles study of the structural and electronic properties of monolayer to few-monolayer MgB 2 , that exposes the origin of the emerging surface state. Secondly, we realized experimentally one to eight monolayers (MLs) of MgB 2 using molecular beam epitaxy, and using in situ room-and low-temperature angle-resolved photoemission spectroscopy (ARPES) we validate the theoretical predictions for the surface state. Finally, by comparing ab initio anisotropic Eliashberg theory results and the measured superconducting gap for six MLs MgB 2 , we reveal the surprisingly distinct and influential contribution of the surface states to the superconducting gap spectrum and tunneling properties of few-monolayer MgB 2 . We show furthermore that the emerging superconductivity in the surface state interacts with the contributions due to σand π-bands, leading to novel multigap phenomena in ultrathin MgB 2 , different for each additional monolayer.

Results
Crystal structure. We begin with structural characterizations. Bulk MgB 2 consists of hexagonal Mg layers intercalated with B layers in a honeycomb lattice. It thus adopts the hexagonal space group P6/mmm (No. 191) 49 . Mg occupies Wyckoff position a 1 , i.e., (0,0,0), while B occupies Wyckoff position d 2 , i.e., ( ) . We studied the crystal structure of few-monolayer MgB 2 by performing a structural relaxation within density functional theory (DFT), as implemented in VASP 50 . We found no tendency to break the in-plane symmetry (no buckling, etc.) in one to eight monolayers (MLs) of MgB 2 , neither in freestanding form, nor with a Mg-substrate.
The crystal structure of one ML of MgB 2 on a Mg-substrate is depicted in Fig. 1a. The calculated in-plane lattice parameter of the Mg substrate is slightly larger than that of a freestanding ML of MgB 2 (3.16 Å compared with 3.09 Å), resulting in an advantageously small lattice mismatch (below 3%) in case the ML is grown on this substrate. The calculated interlayer distance evolves from 3.34 Å for one ML to 3.49 Å for eight MLs, steadily converging towards the bulk value, 3.52 Å. Few-monolayer MgB 2 on a substrate can adopt two different forms, namely with B-or Mg-termination. To find out which termination is energetically favoured, we calculate the heat of formation of both possible structures. The heat of formation per atom, ( ) x y f , of the binary structures with stoichiometry Mg x B y , including the Mg-substrate, is calculated as follows: x y x y f t ot Here, the total energies are calculated using DFT. The total energy of elemental metallic Mg 2 was calculated in its hexagonal form (identical to the substrate), while for B 12 the trigonal α-phase was selected. The results are listed in Table 1. Mg-terminated structures consistently have a lower heat of formation compared with B-terminated structures for more than one ML, the difference being in the order of a few tens of meV per atom. For one ML, the difference in heat of formation of both terminations is small. We therefore find that Mg-termination is preferred in few-monolayer MgB 2 . The Mg-atoms at the free surface stand in contrast with the B-dominated fermiology of bulk MgB 2 . This is the physical origin of all findings that will be discussed further in this article.
On the experimental side, we grew MgB 2 films on a Mg(0001) substrate using molecular beam epitaxy, using a technique we developed for high-purity and high-quality films 51,52 . The Mg(0001) substrate was selected because of the minimal lattice mismatch with few-monolayer MgB 2 as mentioned above. The thickness of the samples was monitored in situ by a combination of photoemission and calibrated evaporators, a technique explained in more detail in the Methods section, and in refs 51,52 . To characterize the lattice parameters of the experimental samples, we carried out surface X-ray diffraction (SXRD) experiments. The results are depicted in Fig. 1b. For the thinnest films, the in-plane distance expanded to perfectly match the substrate. For larger thickness, the lattice parameters were found to evolve towards those of bulk MgB 2 . All measurements, including ARPES, were performed in situ, under ultra-high vacuum conditions, to avoid contamination of the samples.
Electronic structure and electron-phonon interaction. The electronic character of the bands is crucial in understanding the contribution of surface states in few-monolayer MgB 2 . The Fermi surface of MgB 2 , calculated using DFT as a function of the number of layers, is displayed in Fig. 2. The σ-band with p x y , character originates from covalent in-plane B-B bonds, while the π-band stems from out-of-plane B-B bonds with p z character. The Fermi level (E F ) is determined by the ionic character of the Mg-B bonds, through charge transfer. Using Bader charge analysis we found that in bulk MgB 2 Mg-atoms donate approx. 0.8 electrons per B atom. In a freestanding ML, Mg donates 0.7 electrons to each B-atom, hence the similarity of its Fermi contour to the bulk Fermi surface. However, in the Fermi contour of few-monolayer MgB 2 , an additional spherical band is present. It is characteristic of the presence of the free Mg-surface, and therefore denoted S, for surface band. Originating directly from the Mg-atoms facing the vacuum, its band character consists predominantly (~90%) of Mg-p states,   as opposed to the B-p character of the other bands. This Mg-p character of band S is shown in Fig. 3, where the orbitals calculated using DFT are projected onto the atomic orbitals. In one ML on the Mg-substrate, the single B-layer is sandwiched between two Mg-layers. As such, each atom in this B-layer receives 1.3 electrons from the two adjacent Mg-layers, in other words, it is overdoped. As a result, the σ-band is pushed down, as can be observed in Fig. 4 for 1 ML. The σ -band is thus eliminated as part of the Fermi surface of one ML of MgB 2 on the Mg-substrate, as seen in Fig. 2. Band S, on the other hand, remains present as a very robust feature of ultrathin MgB 2 films. For two MLs, the σ -band reappears at the Fermi level, and the difference in Fermi contours between freestanding and epitaxial structures vanishes. For four MLs and more, the σ-band starts converging towards its bulk limit. In contrast, the Mg-based surface band persists as a distinct and prominent feature. Owing to its p z -character, the π-band of few-monolayer MgB 2 hybridizes with Mg-p z states of the substrate. This effect, due to Kronig-Penney-like interaction between Mg-layers in the substrate and the films, diminishes with an increasing number of layers. Therefore, hybridization with the substrate becomes weak in six MLs and more, and the π-band re-emerges at the Fermi level.
To verify the appearance of the surface band near E F in the few-monolayer structures, we performed angle-resolved photoemission spectroscopy (ARPES) at room temperature, shown in Fig. 4 together with the calculated σand Sbands. The latter becomes increasingly localised near the surface for thicker films, explaining the increase in intensity of the ARPES signal. The ARPES signal from the σ-bands is significantly weaker than that of the S-band, especially in the one ML case where the σ-bands are not clearly visible in the ARPES data in Fig. 4 (but more clearly so at higher binding energies, such as in Fig. 1(b) of the Suppl. Mat.). The predicted depletion of the σ-bands for 1 ML MgB 2 is nevertheless in line with the lack of intensity at Γ for MgB 2 layers thicker than 1 ML, Figure 4. Room-temperature ARPES measurements of ultrathin MgB 2 . Data are shown for the valence bands of one, two, four, and eight MLs MgB 2 , along Γ-K, with band structure calculations plotted on top. In the latter, thin lines indicate states originating for 75 to 90% from the layers adjacent to the Mg-vacuum interface, namely those states spatially localised within at most 2 nm from the top of the film -corresponding to the depth down to which ARPES can probe -while thick lines indicate states where this contribution exceeds 90%. The ARPES spectra were acquired using ν = .
indicating that the σ-band is shifted towards lower binding energies in the 1 ML case, as pointed out by the calculations. Finally, there is an interesting difference between the electronic states in freestanding MgB 2 films and those attached to the Mg(0001) substrate. We showed in Figs 2-4 that the σ-bands of films of 2 ML thickness and more, deposited on the Mg(001) substrate, lie close together. In freestanding form, on the contrary, one band splits off from the σ-bands, as shown in Fig. 5, for 6 ML thickness (treated further in the next section regarding its superconducting properties). This split-off band, band ′ S , has B-p xy character and originates from a surface state of the free B-surface, in the absence of the substrate.
In order to characterize the nucleation of superconductivity in few-monolayer MgB 2 from electron-phonon interaction, we calculated the phonon spectrum and the electron-phonon coupling for few-monolayer MgB 2 as well as bulk MgB 2 from first-principles. Using density functional perturbation theory (DFPT), as implemented in ABINIT 53,54 , we determined the electron-phonon coupling coefficients ν + g k k q , , where ν denotes the different phonon modes and k and q are wave vectors. The often used isotropic Eliashberg theory does not suffice to study multigap superconductivity and cannot explain the high transition temperatures of MgB 2 41 . Therefore, we employed fully anisotropic Eliashberg theory with the ν + g k k q , as input for the calculation of the superconducting properties, presented in the next section.

Superconducting properties.
To study superconductivity stemming from the surface band in few-monolayer MgB 2 , we performed low-temperature ARPES measurements, down to 20 K. The results are shown in Fig. 6, together with the spectrum of the Fermi level region of the polycrystalline Ta foil, measured at 10 K, which was used as a reference for the position of E F . The spectra are acquired in the point where band S crosses the Fermi level along the Γ -K direction. Given the very strong ARPES signal on the S-band (cf. Fig. 4), the superconducting gap opening on this band can be determined very accurately. Therefore, our experiment focusses on superconductivity in this band, unique to ultrathin MgB 2 . The position of the spectral edge in the valence band spectra of the one to four ML thick films remained the same at all temperatures and the saddle point along the slope coincided with the E F position of the Ta reference. In contrast, a change in the position of the spectral edge with temperature is observed in Fig. 6a for the six and eight ML thick films. Moreover, this modification occurred only at low temperatures, as can be witnessed in the figure. Such shift towards higher binding energy of the spectral leading edge from the position of the reference E F indicates a vanishing density of states at the Fermi level and therefore the opening of a superconducting gap 55 , ∆, in samples of six and eight ML thick. For samples thinner than six MLs no superconducting gap is observed in the temperature range attainable with our experimental set-up. Also, the proximity effect of the substrate is expected to play a role in suppressing superconductivity in the thinnest limit. Effectively, this comprises two kinds of proximity effects due to the substrate. Firstly, there is the purely electronic proximity effect for one ML MgB 2 , shown in Fig. 2, where charge transfer from the substrate to the film eliminates the σ-bands from the Fermi surface. The second kind of proximity effect stems from a transfer of Cooper pairs from the superconducting film to the metallic substrate, which also depletes superconductivity.  h 9 0 eV incident radiation. The data were acquired at the k-points where the S band crosses E F . One notices the superconducting gap for six and eight MLs at the lowest measured temperature, 20 K. The spectrum acquired on a Ta-foil is shown as a reference for the position of E F . (b) The evolution of the gap on the S -band with temperature, as measured using ARPES. A fit of the profile yields ∆ = . Since purely from theory there is no reason why there would be no superconductivity in the thinner films 56 it is very likely that this second kind of proximity effect inhibits superconductivity in films thinner than six MLs if grown on a metallic substrate.
For six MLs, we performed measurements at several temperatures in the range 20-30 K, shown in Fig. 6b. The gap width, ∆, was obtained by fitting the valence band spectra with a BCS spectral function multiplied by a Fermi-Dirac distribution and convoluted with the experimental Gaussian broadening, the width of which was determined from the Ta spectra. The measurement enables us to fit the superconducting gap as a function of temperature according to ∆ = ∆ − T TT ( ) (0)(1 ( / ) ) p c 1/2 , from which we obtain = .
p 2 4, (0) 3 6 ∆ = . meV, and the critical temperature of this gap, = T 31 S c, K. This superconductivity originating from the surface state cannot be explained by previously proposed models of few-monolayer MgB 2 based on the tight-binding formalism, in which surface states (electronic as well as vibrational) were entirely omitted 57 . So, to unravel the origin of the measured superconductivity in the surface band, and its relation to the other gaps, we employed fully anisotropic Eliashberg theory on a freestanding, six ML MgB 2 film, using the electron-phonon coupling coefficients obtained from first-principles as input (more detailed information on which can be found in Methods). The calculated distribution of the gap opening specifically on band S as a function of temperature is shown in Fig. 6a (dashed line). It is observed that this S-gap closes at T 33 0 13 ⁎ to model the electron-electron repulsion that counteracts Cooper pair formation 58 . From the close agreement between our theoretical predictions and experimental data it can be concluded that the screening modelled by µ ⁎ does not change in the ultrathin limit in the case of MgB 2 , as opposed to, e.g., ultrathin TaS 2 where it was proposed that a renormalization of this repulsion lies at the base of the increase in critical temperature in thinner samples 59 .
The gap spectrum displayed in Fig. 7a is very rich, and the S-gap is certainly not the only contribution. In order to understand the physical origin of this rich spectrum, we plot the calculated distribution of the gap on the Fermi surface, ( ) ρ ∆ , at 1 K, in Fig. 7b. It consists of three domes of which two (with values ~3-3.5 meV and ~6-7 meV) are comparable to the πand σ-gaps of bulk MgB 2 respectively. However, the free surfaces provide additional contributions that phonons couple to, akin to, e.g., the interlayer state crucial for superconductivity in functionalized graphene 17,19,20 . The Mg-based band S, resulting from the Mg-vacuum interface, contributes to the first dome by hybridizing with the π-gap (denoted Sπ-dome). The gap on band S′ is in the range 4-5 meV and is also mixed with the π-gap, forming the central dome of ρ ∆ ( ) denoted π ′ S . In the experiment, this B-surface is bound to the substrate, and as a result band S′ approaches the σ-bands.
As bands S and ′ S stem from the free surfaces of opposite sides of the film, what happens to band S′ does not influence the gap opening on the S-band of the Mg-surface. Therefore, separate control of the contributions of these two surface states to the gap spectrum can be achieved by means of chemical functionalization of one or either of these free surfaces. Specifically, we envisage that the surface bands can be eliminated from the Fermi surface by means of adatoms on the Mg-side (for the S-band) and/or the B-side (for the ′ S -band). With elements with a low atomic number, such as hydrogen, the charge transfer between the adatoms and the film can be limited, such that the σand π-states would be largely unaltered under the influence of this functionalization. This provides for a unique (possibly local) control of the superconducting gap spectrum through easily accessible surfaces.
Furthermore, the different domes in the superconducting gap spectrum with their distinct contributions in six MLs MgB 2 , Sπ, S π ′ and σ, are robust under the influence of temperature, persisting even at temperatures close to the overall critical temperature of six MLs of MgB 2 , = T 35 c K. This is concluded from the temperature dependence of the full gap distribution, ( ) ρ ∆ , displayed in Fig. 7a, together with the weighted averages of the different domes. The effect of the surface contributions on the density of states (DOS) in the superconducting state, N S , calculated from the superconducting gap spectrum, is also very pronounced, as displayed in Fig. 7c. N S exhibits three sharply resolved peaks, corresponding to the three domes in ρ ∆ ( ), whereas in bulk MgB 2 the DOS in the superconducting state consists of only two peaks, corresponding to the σand π-gaps (as shown both in theory 42 and in scanning-tunneling microscopy experiments 60 ). As such, the tunneling properties undergo equally radical change going from bulk to ultrathin MgB 2 .

Discussion
Our combined theoretical and experimental study of few-monolayer MgB 2 unambiguously reveals that surface states drastically change superconductivity in atomically thin superconductors. The surface states in question stem from the outer layers of the material itself, i.e., the free surfaces facing vacuum, without need for any specific adatoms, intercalants or substrates that were previously shown to stimulate superconductivity in the two-dimensional limit. Moreover, the impact of free surfaces on superconductivity is inherent to atomically thin materials, but still susceptible to further chemical functionalization and nanoscale engineering. Our measured and calculated enriched gap spectrum of six MLs of MgB 2 proves that these surface states can be exploited to qualitatively change the multigap nature of superconductivity in the atomically thin limit. We showed furthermore that the condensate originating from the surface state hybridizes with the π-condensate in few-monolayer MgB 2 . This is completely different from bulk-sized MgB 2 where ARPES experiments indicate near-degeneracy of a surface state gap with the σ-gap, without any clear influence on the measured two-gap superconductivity 39 . This comparison implies a strong evolution of the superconducting gap spectrum, and strong changes in the interplay and hybridization between condensates when thinning MgB 2 layer by layer. Therefore, as an outlook of this work, measurements below 20 K are foreseen to explore possible superconductivity in the thinnest structures of MgB 2 -in line with theoretical predictions presented in ref. 56 . In that limit, our results indicate a considerable role of the substrate since, e.g., the Mg-substrate eliminates several electronic bands in case of a single ML of MgB 2 (cf. Fig. 2a). Non-metallic substrates would exclude the possibility of Cooper-pairs escaping and would thus make the superconducting state even more robust. On a related note, recent experimental advances enable non-epitaxial fabrication of not only weakly bound layered materials, e.g., NbSe 2 14 , but also of ultrathin compounds with ionic interlayer bonds, like MgB 2 61 . It is thus envisioned that the full potential of superconductivity from free surfaces can be explored using the available technology and can be further incorporated in, e.g., atomically thin devices aimed at cryogenic computing.

Methods
Sample Preparation and characterization. The Mg(0001) substrate was prepared using cycles of sputtering and annealing at 493 K. During molecular beam epitaxy (MBE) of Mg and B from pure sources, with atomic flux ratio Mg:B = 3:2, the pressure was kept below − 10 9 mbar, while the substrate was held at 475 ± 15 K. For the thinnest films (<4 ML), deposition of B was sufficient, whereas for thicker layers co-deposition of Mg and B was needed. We monitored the film thickness through the evaporation rate, using the calibrated attenuation of the substrate core level photoemission peaks when evaporating B on Mg(0001) and Mg on a copper plate. Surface X-ray diffraction (SXRD) was carried out at the ALOISA beamline of the Elettra Synchrotron, Trieste, Italy 62 .
Angle Resolved Photoemission Spectroscopy. Angle-resolved photoelectron spectroscopy (ARPES) was performed at the BaDElPh beamline of the Elettra Synchrotron 63 to map the variations in the electronic band structure of MgB 2 as a function of film thickness. The spectra were collected both at room temperature and at low temperatures with a SPECS Phoibos-150 hemispherical electron analyzer, equipped with a 2D CCD sensor. The nominal energy resolution of the detector was 5 meV and its nominal angular resolution was 0.15°. In the low-temperature ARPES experiments, the samples were cooled using a liquid He-cryostat, and the temperature was measured by means of silicon diodes. The energy resolution and stability of the photon beam energy throughout the course of the experiment was monitored using photoemission spectra around the Fermi level of the clean Mg(0001) at the Γ point of the Brillouin zone and that of a polycrystalline Ta-foil, in electric contact with the sample. More detailed information on stability and resolution in the ARPES measurements can be found in Supplementary Information, section 1.

Density Functional (Perturbation) Theory.
Our density functional theory (DFT) calculations make use of the Perdew-Burke-Ernzerhof (PBE) functional implemented within a planewave basis in the VASP code. Electron-ion interactions are treated using projector augmented wave (PAW) potentials, taking into account Mg-2p 6 3s 2 and B-2s 2 2p 1 as valence electrons. An energy cutoff of 450 eV for the planewave basis was used, to achieve convergence of the total energy below 1 meV per atom. To obtain a very accurate description of the Fermi surface, a very dense × × 35 35 31 Γ -centered Monkhorst-Pack k-point grid is used for bulk MgB 2 , and a × × 35 35 1 grid for few-monolayer structures. The lattice parameters were obtained using a conjugate-gradient algorithm so that forces on each atom were minimized below 1 meV/Å. Subsequently, density functional perturbation theory (DFPT) calculations were carried out within the framework of ABINIT, keeping the same valence electrons, and also using the PBE functional. The crystal structure was optimized again in ABINIT, with no significant differences with the VASP results. The total number of perturbations due to atomic displacements to be treated (in other words, the number of phonon branches) amounts to ⋅ = N 3 54 atoms for six MLs. In this way, the phonon spectrum and electron-phonon coupling coefficients, matrix elements of the perturbative part of the Hamiltonian 54 are obtained. We carried out the DFPT calculations on a 22 22 1 × × electronic k-point grid and a × × 11 11 1 q-point grid, a subgrid of the k-point grid, for the phonons. The DFPT calculations are performed on freestanding six MLs, thus the interface state between film and substrate, as well as interface phonons are deliberately omitted. The phonons and electron-phonon coupling of six ML and bulk MgB 2 are discussed in Supplementary Information, section 2.
Anisotropic Eliashberg Theory. The anisotropic Eliashberg equations were solved self-consistently in Matsubara space using the Uppsala superconductivity (UppSC) code, starting from the electron and phonon band structures and electron-phonon coupling obtained with DFPT, after which the converged solutions were analytically continued to real frequencies (see ref. 42 for further details). In this scheme, we iterated until convergence better than 10 −6 on the relative gap values was reached. In all calculations, we employed standard µ = . Data availability statement. The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.