Effect of Confinement and Coulomb Interactions on the Electronic Structure of the (111) LaAlO₃/SrTiO₃ Interface

Abstract

A tight binding supercell approach is used for the calculation of the electronic structure of the (111) LaAlO${}_3$/SrTiO${}_3$ interface. The confinement potential at the interface is evaluated solving a discrete Poisson equation by means of an iterative method. In addition to the effect of the confinement, local Hubbard electron–electron terms are included at the mean-field level within a fully self-consistent procedure. The calculation carefully describes how the two-dimensional electron gas arises from the quantum confinement of electrons near the interface due to the band bending potential. The resulting electronic sub-bands and Fermi surfaces show full agreement with the electronic structure determined by angle-resolved photoelectron spectroscopy experiments. In particular, we analyse how the effect of local Hubbard interactions change the density distribution over the layers from the interface to the bulk. Interestingly, the two-dimensional electron gas at the interface is not depleted by local Hubbard interactions which indeed induce an enhancement of the electron density between the first layers and the bulk.

1. Introduction

Recently, the emergent field of oxide electronics has revealed a rich phenomenology connected to the creation and manipulation of interface electronic states. After the discovery of two-dimensional electron gas (2DEG) at the (001) interface between the perovskite band insulators SrTiO${}_3$ (STO) and LaAlO${}_3$(LAO) [1], which are characterized by high-mobility, much work has been devoted to revealing its properties, like gate-controlled metal-insulator transitions [2], superconductivity [3], including topological ones [4,5], and its possible coexistence with magnetism [6]. Recently, the successful creation of 2DEGs at the (111)-oriented interface of LAO–STO [7,8] has opened the possibility to investigate intriguing phenomena related to topological phase transitions [9], gate tunable anomalous Hall effect [10] and the spin–orbital Edelstein effect [11]. Despite the great number of works, supplementary analysis of the band structure of (111) LAO–STO is still required, especially in relation to the confinement effects and the role of electronic correlations.

A first qualitative understanding of the band structure and Fermi surface has come from a tight-binding (TB) supercell calculation based on an ab initio bulk band structure discussed in [12]. Calculation of the surface electronic structure was performed by introducing a supercell containing 120 Ti atoms stacked along the $\left(111\right)$ direction and using maximally localized Wannier functions with additional on-site potential terms to account for band bending via an electrostatic potential. The TB Hamiltonian was solved self-consistently with Poisson’s equation, incorporating an electric field dependent dielectric constant [13,14] and with only an adjustable parameter, the total magnitude of the band bending at the surface [12]. The derived Fermi surface (FS) consists of three equivalent elliptical sheets oriented along the $\Gamma -M$ direction. The band structure along the $\Gamma -M$ direction shows a single heavy band, corresponding to the long axis of one of the FS ellipses, which is nearly degraded at the band bottom with a more dispersive, doubly deviated band arising from the two remaining FS sheets. The band structure shows three confined 2DEG sub-bands arising from the $t_{2g}$ orbitals and a “ladder” of states with a bulk-like character above the Fermi level due to the finite size of the supercell. The second sub-band was predicted to be just below the Fermi energy, in good agreement with the angle-resolved photoemission spectroscopy (ARPES) experiments. The wave functions of the lowest sub-band at the $\Gamma$ point was predicted to be extended over the ≃15 Ti layers, an order of magnitude more than the lowest bulk sub-band on the (001) STO, due to the lighter effective masses.

In this work we perform a TB supercell calculation for the (111) LAO–STO, and crucially, beyond the effect of the confinement, we also include local Hubbard electron–electron interactions within a fully self-consistent mean-field approach. Furthermore, compared with ref. [12], we also account for the impact of SOC. In particular, the TB supercell Hamiltonian in the (111) direction is obtained by rotating the coordinates and converting the quasi-momentum degree of freedom along the (111) direction to the discrete index numbering the layers of Ti along the (111) axis. Our calculation shows full agreement with the observed electronic structure by ARPES [12] and describes how the 2DEG arises from the quantum confinement of $t_{2g}$ electrons near the surface due to band bending. Moreover, we also demonstrate how the effect of local Hubbard terms changes the density distribution over the layers close to the surface. We show that, contrary to a naif expectation, the 2DEG at the interface is not depleted by local Hubbard interactions which instead induce a modulation of the electron density as a function of the layer number. In fact, we find that local Hubbard terms enhance the electron density between the first layers of the interface and the bulk. The role of the electron potential was previously studied in Ref. [15] for the (001) interface with a different approach.

The manuscript is organized as follows. In Section 2 we introduce the Hamiltonian, the TB supercell approach and present the results for the band structure, Fermi surface, and self-consistent band bending potential. We analyse the effects of both the absence and presence of local Hubbard electron–electron interactions. In Section 3 we discuss our results and give a comparative discussion of previous studies.

2. Methods and Results

In this section we present the model and results obtained within a TB supercell approach, both in the absence and presence of local Hubbard electron-electron interactions.

2.1. Model

STO has a cubic perovskite structure as shown in Figure 1.

The conductance bands form out of the $t_{2g}=\{d_{yz},d_{zx},d_{xy}\}$ orbitals of the Ti atoms in the bulk structure. Therefore, we focus only on the Ti lattice, which has a simple cubic structure at room temperature [17] with a lattice constant $a_0=0.3905$ nm. The presence of a thin film of LAO over STO leads to the formation of a 2DEG and a thin positive charge density at its top. This is mostly ascribed to oxygen vacancies in STO [12] leading to an electronic reconfiguration which neutralizes this positive charge. As a consequence, the conduction band is partially filled, so that the electronic properties are determined by the low-energy region of such bands. In the present work we will take the positive charge at the interface as a free parameter of the model. In this sense we manipulate the number of oxygen vacancies.

In order to reconstruct the self-consistent electronic band structure, we adopt a TB Hamiltonian framed in the basis of atomic orbitals, using hopping parameters which fit the available ARPES data for the lowest bands [12]. The bulk Hamiltonian for the conductance bands in STO, expressed in the quasi-momentum $(K_x,K_y,K_z)$ directed along the cubic axes, is

$$H_{\mathrm{TB}}^{\mathrm{Bulk}}=\sum _{\overrightarrow{K}}\sum _{i\ne j\ne k}\sum _{\sigma }\left(-t_{D}cos\left(K_i{a}_0\right)-t_{D}cos\left(K_j{a}_0\right)-t_{I}cos\left(K_k{a}_0\right))\right)d_{ij,\sigma ,\overrightarrow{K}}^{†}{d}_{ij,\sigma ,\overrightarrow{K}}$$

(1)

where we truncated the nearest-neighbour hopping. Here $\{i,j,k\}$ runs over $\{x,y,z\}$, $d_{ij,\sigma ,\overrightarrow{K}}$ is the annihilation operator of the electron characterized by the $d_{ij}$ orbital, spin $\sigma$ and quasi-momentum $\overrightarrow{K}$. $t_D$ and $t_I$ are the direct and indirect hopping parameters, which we choose to be $t_D=0.25$ eV and $t_I=0.02$ eV [9] in agreement with ARPES data. Since the electric field produced from the interfacial charge breaks the translational invariance along the (111) direction, the Hamiltonian expressed in terms of the quasi-momentum component along the (111) is not the optimal choice for the description of the two-dimensional gas. Therefore, from this Hamiltonian, we construct a TB supercell Hamiltonian in the (111) direction by a rotation of coordinates and converting the quasi-momentum degree of freedom along the (111) direction to the discrete index numbering the layers of Ti along the (111) axis. By this procedure, we convert the $6\times 6$ bulk Hamiltonian (considering the spin degree of freedom) to a $6N\times 6N$ Hamiltonian, for which N represents the number of layers considered (in the bulk system $N\to \infty$). We include two other local terms in real coordinates, which are therefore independent of $\overrightarrow{K}$: the spin–orbit coupling (SOC) $H_{\mathrm{SOC}}$, and a trigonal crystal field along the (111) direction $H_{\mathrm{TRI}}$ [9,10,11,18]. The matrix for the TB supercell Hamiltonian has the form

$$H=\left(\begin{array}{cccccc}{H}_0& H_t& 0& 0& 0& \dots \\ H_t^{†}& H_0& H_t& 0& 0& \dots \\ 0& H_t^{†}& H_0& H_t& 0& \dots \\ 0& 0& H_t^{†}& H_0& H_t& \dots \\ \dots & \dots & \dots & \dots & \dots & \dots \end{array}\right),$$

(2)

where

$$H_0=H_{\mathrm{SOC}}+H_{\mathrm{TRI}}$$

(3)

and $H_t$ is the tunnelling Hamiltonian describing the hopping between two neighbouring layers for a given state of defined quasi-momentum parallel to the interface. In the Appendix A we give all the details of the calculation and the explicit forms of the in-plane contributions and out-of-plane hoppings. In order to model a slab of the material, we cut the block matrix to a finite size, which in this paper is fixed to 51 layers.

On top of this matrix, we introduce a potential $\phi$, which includes both a contribution from the interfacial charge and a screening contribution from the electrons themselves which populate the interface. Therefore, this component has to be determined self-consistently. In order to do this, we fix the positive charge density $\rho$ at the beginning of the slab, and solve the classical equations of the electrodynamics along the $\left(111\right)=\widehat{Z}$ direction

$${\partial }_Z\left(\epsilon \left(F\right){\partial }_Z\phi \right)=-\frac{\rho }{{\epsilon }_0}$$

(4)

which for a discrete system of infinite charged planes becomes

$$\left\{\begin{array}{c}{\phi }_l=-\frac{a_0}{\sqrt{3}}\sum _{l^{\prime }=1}^l{F}_{l^{\prime }}\hfill \\ {\epsilon }_0\epsilon \left(F_l\right)F_l=D_l\hfill \\ D_l=\left|e\right|(n_{2D}-n_l)\hfill \end{array}\right.$$

(5)

where $D_l$ is the electric displacement, $F_l$ the electric field, ${\epsilon }_0$ is the absolute dielectric constant value, $\epsilon \left(F\right)$ is the relative dielectric constant, $n_{2D}$ is the total positive density charge at the interface divided by the in-plane elementary unit cell surface $a_0^2$, while $n_l$ is the 2D density charge on the layer l. Equation (5) involve $\epsilon$, which leads to solutions which are sensitive to the choice of dielectric constant model. We choose $\epsilon$ at zero temperature as indicated in ref. [19].

$$\epsilon \left(F\right)=1+\frac{{\chi }_0}{{(1+{\left(\frac{F}{F_0}\right)}^2)}^{1/3}},$$

(6)

where ${\chi }_0=21,000$ and $F_0=80,000$ V/m, which for $F=0$ tends to the standard order of magnitude in STO at low temperatures [20]. The choice of a $F^{-2/3}$ dependence represents the STO ferroelectric behaviour at low temperatures [21]. This particular behaviour is motivated by the Barret formula [22]. Other parametrization in the literature are adopted in ref. [20].

We adopt the following procedure to reach self-consistency: we fix a value of $n_{2D}$ and a trial potential ${\phi }_0$, and include it in the Hamiltonian and diagonalize it. We find the chemical potential at which the total electron density is $n_{2D}$ and compute the electron density on each layer. We use this density to solve the system (5) and obtain the potential ${\tilde{\phi }}_0$. At this point the input potential in the Hamiltonian is ${\phi }_1=\alpha {\tilde{\phi }}_0+(1-\alpha ){\phi }_0$, where $\alpha$ is chosen to guarantee a stable convergence. We repeat the procedure until ${\phi }_{i+1}\approx {\phi }_i$. Typical values of $\alpha$ are $0.8$, $0.9$ and $0.95$, and the stopping criterion is $|{\phi }_{i+1}-{\phi }_i|<{10}^{-2}$ eV. In the following we choose three benchmark choices of $n_{2D}$ ($n_{2D}=1\times {10}^{14}$ cm${}^{-2}$, $n_{2D}=2\times {10}^{14}$ cm${}^{-2}$, and $n_{2D}=3\times {10}^{14}$ cm${}^{-2}$) in order to study changes of the electronic confinement induced by increasing values of electron density. This analysis is rather relevant because the self-consistent densities of 2DEG are in agreement with the expectations, of the order of $1.5\times {10}^{14}$ cm${}^{-2}$ [12]).

2.2. Results in the Absence of Local Hubbard Interaction Terms

In order to clarify the effect of screening, we self-consistently obtain the eigenstates of the Hamiltonian (2) with the potential of Equation (5) for the three different benchmark parameters defined above. The band structures for low fillings, expressed in terms of the dimensionless quasi-momentum $\overrightarrow{k}=\overrightarrow{K}{a}_0\sqrt{\frac{2}{3}}$, are shown in Figure 2.

The original $6N$ degrees of freedom are structured in N subsets of 6 bands each, which are progressively less confined at the interface. By increasing the positive charge at the interface, the splitting between the subsets of bands increases. For $n_{2D}=1\times {10}^{14}$ cm${}^{-2}$, two subsets of bands intersect with one another, while for $n_{2D}=2\times {10}^{14}$ cm${}^{-2}$, and $n_{2D}=3\times {10}^{14}$ cm${}^{-2}$, the first subset of bands is separated from the higher bands. The identification of the subset to which a band belongs is more easily performed by looking at the number of nodes of the wave function for each band evaluated at $\overrightarrow{k}=0$. We show this result in Appendix B. The splitting, and the confinement in turn, is proportional to the slope of the potential close to the first layer. Figure 3 shows the potential for all benchmark choices and the 2D charge density for each layer.

Independently of the filling density $n_{2D}$, the chemical potential always lies well above the maximum value of the potential well. Therefore, a bulk contribution is always present, as also visible in the right panel of Figure 3. However, the higher the 2D density, the smaller the confinement region of the quasi-2DEG, which shrinks to almost 10 layers for $n_{2D}=3\times {10}^{14}$ cm${}^{-2}$. In this case, the 2D charge density of the 2DEG reaches the value of ${\sum }_{l=1}^{10}{n}_l\approx 2.5\times {10}^{14}$ cm${}^{-2}$, the typical order of magnitude for 2DEG densities [12]. We point out that the way we extract the density of the 2DEG is different from the one which is used from the typical (001) interface, since in the (111) interface the first band is not really 2D due to the fact that the electron hopping happens between the different planes. Therefore, the estimation of the 2DEG density computed as the area of the FS of the external band is naive. An interesting feature which appears above the chemical potential is the band splitting at $\overrightarrow{k}=0$ for the higher bands. This splitting disappears by removing the atomic SOC and thus resembles a huge linear Rashba splitting induced by the combined effect of the electric potential, which naturally breaks the inversion symmetry, and the atomic SOC. Even if the Fermi energy is at the filling energy, and lower than splitting, the chemical potential can be changed by using an external gate voltage, in order to investigate these kind of bands. The application of such a gate does not change the splitting between the bands, since this is mostly influenced by the positive charge at the interface.

In Figure 4 we show the FSs of the band for the benchmark density $n_{2D}=3\times {10}^{14}$ cm${}^{-2}$. The contour shows a six-fold symmetry of the energy spectra and, in agreement with the analysis of ref. [10,11,12], every ellipse-shaped band possesses a strong d-orbital character away from $\overrightarrow{k}\sim 0$, while nearby the $\Gamma$ point the $a_g$ or $e_g^{\pi }$ character of the band is restored (i.e., the spherical symmetry of the inner bands), due to the trigonal crystal field.

2.3. Effect of Local Hubbard Interaction Terms

In this section we explore the effect of local Hubbard electron–electron terms on the results we presented above. These interactions can be expressed in a mean-field approximation as

$$H_{\mathrm{C}}=\sum _{\overrightarrow{K}}\sum _{\alpha ,l}{n}_{\alpha \uparrow ,l,\overrightarrow{K}}\left(U\langle n_{\alpha \downarrow ,l}\rangle +\sum _{\beta \ne \alpha }{U}^{\prime }\langle n_{\beta \uparrow ,l}\rangle +U^{\prime }\langle n_{\beta \downarrow ,l}\rangle \right)+(\uparrow \leftrightarrow \downarrow ),$$

(7)

where $n_{\alpha \uparrow ,l,\overrightarrow{K}}=d_{\alpha \uparrow ,l,\overrightarrow{K}}^{†}{d}_{\alpha \uparrow ,l,\overrightarrow{K}}$, we exploited the spatial homogeneity in the interfacial plane so that mean densities are independent of positions, $\alpha$ runs over the orbital degree of freedom, and U and $U^{\prime }$ parametrize the strengths of the interaction. In the absence of any term breaking the $C_{3v}$ symmetry and the time-reversal invariance, for each layer the electron density is equal for every orbital and spin. In such a regime, the $U^{\prime }$ terms provide only a renormalization of the U term; therefore, we can neglect them reducing at the same time the computational effort for the simultaneous self-consistent calculation of the local Hubbard potential $U\langle n_l\rangle$ and the potential ${\phi }_l$. Equation (7) shows that the local Coulomb interaction introduces an effective potential varying over each layer, proportional to the local particle density at that layer. We expect that this leads to a broadening of the electron density over the whole slab of material, since it energetically favours the lowest occupied layers.

We choose a benchmark value of $U=4$ eV, as chosen in ref. [8], and compute the bands shown in Figure 5 for the same benchmark values of $n_{2D}$, in order to clearly discriminate the effects originating from local Hubbard interactions on the band structure. For $n_{2D}=1\times {10}^{14}$ cm${}^{-2}$ and $n_{2D}=2\times {10}^{14}$ cm${}^{-2}$, the band structures change only slightly from the situations without Hubbard terms. However, at $n_{2D}=3\times {10}^{14}$ cm${}^{-2}$ the local Hubbard terms significantly enhance the separation between the confined and free bands. Moreover, due to the effect of the Hubbard terms, an additional sub-band crosses the chemical potential close to the $\Gamma$ point. As a result, as shown in Figure 4, the FS shows a reconfiguration in the same region of the Brillouin zone. We can clarify interesting features of the electronic band structure by comparing the self-consistent potential with and without the effects of the Hubbard terms, in the former case also accounting for the effective potential $U\langle n_l\rangle =U{n}_l/6$. We show these results in Figure 6. For $n_{2D}=3\times {10}^{14}$ cm${}^{-2}$, the effective potential exhibits a small peak which enhances the separation between the high-slope region to the plateau.

Even in this case, however, the chemical potential is above the highest value of the potential, which results in a bulk component. The electron density for the case of $U=4$ eV is reported in Figure 7. Actually, it shows an intermediate region between the layers 10 and 30 where the local Hubbard terms induce an enhancement of the local density. This is a very interesting result since the Hubbard terms do not deplete the 2DEG but adds a modulation of the electron density in the intermediate region. This is mainly due to the fact that the interfacial electron charge concentration repels electrons in the intermediate region, creating a slight inflexion in the effective potential felt by each electron, mildly favouring occupation of the intermediate region.

3. Discussion and Conclusions

In this manuscript we have systematically discussed a TB supercell method to describe the band structure of the (111) LAO–STO interface taking into account the effect of atomic SOC, trigonal strain, the electronic confinement and local Hubbard electron–electron interactions within a fully self-consistent procedure. The sub-band energy and the FS show full agreement with the observed electronic structure by ARPES [12] and describes how the 2DEG arises from the quantum confinement of $t_{2g}$ electrons near the surface due to the band bending potential, for different values of the positive charge density at the interface, i.e., by varying the concentration of the oxygen vacancies. In particular, we have shown how the effect of local Hubbard electron–electron interactions change the density distribution over the layers close to the surface. Indeed, the 2DEG at the interface is not depleted by local Hubbard interactions which, instead, induce a modulation of the electron density as a function of the layers number. The net effect of the Hubbard interactions is to enhance the electron density in the intermediate spatial region between the first layers at the interface and bulk. This effect was previously analyzed in Ref. [23] for the (001) interface for which the Coulomb interactions have the effect of helping the formation of the 2DEG. The main difference compared to the (001) is the orbital structure of the electrons at the interface which causes a different impact of the Coulomb interactions on the results. However, we do not find any qualitative changes in the behaviour of the eigenstates in the confined region due to the effect of the Hubbard interaction, as found in ref. [8]. This is possibly due to different numerical strategies used in the band structure computation. Moreover, despite in our model a Rashba coupling or any odd hopping parameters in the quasi-momentum was not inserted explicitly, an interesting huge linear splitting for the higher bands above the Fermi level appears, probably due to the combined effect of SOC and the electric potential. Its origin can be the subject of further studies.