Lattice screening of the polar catastrophe and hidden in-plane polarization in KNbO$_{3}$/BaTiO$_{3}$ interfaces

We have carried out first-principles simulations, based on density functional theory, to obtain the atomic and electronic structure of (001) BaTiO$_{3}$/KNbO$_{3}$ interfaces in an isolated slab geometry. We tried different types of structures including symmetric and asymmetric configurations, and variations in the thickness of the constituent materials. The spontaneous polarization of the layer-by-layer non neutral material (KNbO$_{3}$) in these interfaces cancels out almost exactly the"built-in"polarization responsible for the electronic reconstruction. As a consequence, the so-called polar catastrophe is quenched and all the simulated interfaces are insulating. A model, based on the modern theory of polarization and basic electrostatics, allows an estimation of the critical thickness for the formation of the two-dimensional electron gas between 42 and 44 KNbO$_{3}$ unit cells. We also demonstrate the presence of an unexpected in-plane polarization in BaTiO$_{3}$ localized at the $p$-type TiO$_2$/KO interface, even under in-plane compressive strains. We expect this in-plane polarization to remain hidden due to angular averaging during quantum fluctuations unless the symmetry is broken with small electric fields.


I. INTRODUCTION
The surprising discovery by Ohtomo and Hwang 1 of a metallic state at the interface between two good band insulating oxides, LaAlO 3 and SrTiO 3 , has triggered a large amount of new studies on polar oxide interfaces. 2 Indeed, the quasi two-dimensional electron gas (2DEG) that forms when LaAlO 3 is grown on top of a TiO 2 terminated (001)-surface of SrTiO 3 , i.e. when the interface between the two materials is LaO/TiO 2 , displays very different properties from those generated at interfaces between standard III-V semiconductors (such as GaAs and Al x Ga 1−x As). Among them we find conducting carrier densities and electron effective masses orders of magnitude larger than those found at semiconductor interfaces. 3 It is also fascinating how, depending on growth conditions, magnetic 4 and superconducting 5 ground states have been experimentally identified at this interface between non-magnetic insulating oxides. Very recently, two independent groups have proven how both magnetic and superconducting states might even coexist on the same sample, 6,7 a very unexpected result since magnetic order is usually considered detrimental to superconductivity. As a consequence of all these phenomena, interfaces in polar oxides can open the door to novel implementations of field effect transistors and to a new era of oxide electronics. 3,8 Despite this recent activity, many fundamental questions regarding the origin and confinement of the 2DEG remain highly debated. Different models have been proposed to explain the experimental results. The pioneering one invokes the so called "polar catastrophe" 9 that arises from the polarization discontinuity 10  and Ti +4 O −2 2 0 layers of the perovskite structure of SrTiO 3 are charge neutral. But, aside this first rationalization, other explanations can be found in the literature for the origin of the 2DEG, among them: (i) the interlayer mixing between LaAlO 3 and SrTiO 3 and non-abruptness of the interfaces (with the formation of a few monolayers of metallic La 1−x Sr x TiO 3 11 ), (ii) doping due to oxygen vacancies 12,13 (including those produced in surface redox reactions 14 ), and/or (iii) the presence of charged defects and adsorbates. 15,16 All these models highlight the importance of the growth conditions of these structures for the appearance and the behavior of the functional properties of the 2DEG.
One of these properties, that is well reproduced by different experimental groups on many samples grown with a variety of techniques, is the existence of a critical thickness, t c , in the number of layers of LaAlO 3 for the formation of the 2DEG. Thiel and coworkers 17 have demonstrated that, for the interfaces to be conducting, the number of layers of LaAlO 3 has to be larger than four unit cells. The thickness of the polar layer increases up to five unit cells in n-type LaVO 3 /SrTiO 3 interfaces. 18 These observations are consistent with the fact that the conductivity of SrTiO 3 -LaAlO 3 -SrTiO 3 heterostructures with dissimilar interfaces is reduced if their p-type (AlO 2 /SrO) and n-type (LaO/TiO 2 ) interfaces are spaced by less than six unit cells. 19,20 Remarkably, even below the critical thickness, a metal-insulator transition can be driven by an external electric field. 15,17 These studies suggested the possibility of the design of new polar interfaces where the appearance of the 2DEG could be switched on and off by the action of an external perturbation.
Along this line, the replacement of one (or the two) materials at the interface by ferroelectric perovskites is particularly attractive. The spontaneous polarization present in these materials is very sensitive to electric fields and could be used to create a bound charge at the interface that could reinforce/deplete the 2DEG. Previous theoretical works have been focused on I-V/II-IV interfaces [NaNbO 3 /SrTiO 3 , 21,22 and KNbO 3 /ATiO 3 interfaces, where A = Sr, Ba or Pb [22][23][24] ]. From the formal ionic charge point of view, I-V ferroelectric perovskite oxides, such as NaNbO 3 and KNbO 3 , are made of alternating positive (B +5 O −2 2 ) +1 and negative (A +1 O −2 ) −1 charged layers along the [001] direction (essentially as LaAlO 3 , although now the AO layers of the perovskite ABO 3 structure are negative, while the BO 2 layers are positive). Therefore, the layer-by-layer electrostatic of the previous I-V/II-IV interfaces is analogous to that in the LaAlO 3 /SrTiO 3 interface. Density functional theory simulations on non-stoichiometric (i.e. with a noninteger number of unit cells of the layer-by-layer non neutral perovskite), symmetric superlattices indeed suggested the existence of a 2DEG in KNbO 3 /ATiO 3 interfaces, switchable between two conducting states by the ferroelectric polarization orientation of the titanate layer. 23,24 However, with the simulation boxes used in the previous works only the n-type or p-type interfaces are present. It can be proved (see Sec. 4 of the Supplemental Material of Ref. 25) that, within this configuration, the local interface properties exactly reproduce those of the infinite isolated slab geometries, obviously beyond the critical thickness for the formation of the 2DEG. In other words, the calculations of Refs. 23 and Ref. 24 show the charge distribution and properties after the electronic reconstruction has taken place, but nothing is said about the magnitude of the critical thickness for the formation of the 2DEG. 26 In this work we carry out first principles calculations of BaTiO 3 /KNbO 3 interfaces where we explicitly avoid the issue of lack of stoichiometry in the simulation box. We have found that, due to the large lattice screening provided by the KNbO 3 layer, the critical thickness for the formation of the 2DEG is one order of magnitude larger than in LaAlO 3 /SrTiO 3 interfaces. An unexpected result of our simulations is that an in-plane polarization develops on the BaTiO 3 side of a (TiO 2 /KO) p-type BaTiO 3 /KNbO 3 interface even when BaTiO 3 is subject to in-plane compressive strains. We explain this effect using basic electrostatic arguments.
The rest of the paper is organized as follows. After summarizing the basic theory behind the electronic reconstruction in Sec. II, we present the computational details used in our simulations in Sec. III. The firstprinciples results, together with the relevant comparisons to the model, can be found in Sec. IV.

II. BACKGROUND ON THE "POLAR DISCONTINUITY" MODEL
In order to establish the nomenclature and the basic theory that will be used later, we review the most important points of the "polar discontinuity" model. Although this model has been invoked since the discovery of the 2DEG at polar oxide interfaces, 9 only recently it has been rigorously rationalized with explanations firmly rooted on the modern theory of polarization (for a recent review, see Ref. 27 and references therein). This has been developed by Stengel and Vanderbilt in Ref. 10 for insulating interfaces, and later generalized by Stengel for the case of a non-zero surface density of "free" charge in Ref. 28, and to the case of surfaces in Ref. 25. The theory presented in these works is absolutely general, and we strongly point the interested reader to those milestone papers. Here, we particularize it to the conditions considered in this work, and estimate the critical thickness for the formation of 2DEG in the case where any of the two materials forming the interface is ferroelectric.
The standard nomenclature used in the literature of polar oxide interfaces denote the material that is nonneutral layer-by-layer (i.e. LaAlO 3 ) as polar, and the material that is neutral layer-by-layer (i.e. SrTiO 3 ) as non-polar. Rigorously, this notation does not apply here since the two materials that constitute our interfaces are ferroelectric and, therefore, might undergo polar phase transitions with the appearance of a non-vanishing spontaneous polarization. Nevertheless, for the sake of consistency with previous works we will maintain the convention and refer to the ferroelectric non-neutral layer by layer material (i.e. KNbO 3 ) as the polar material and the ferroelectric neutral layer-by-layer material (i.e. BaTiO 3 ) as the non-polar one.
During the development of the model we assume a ntype interface, simulated within an isolated slab geometry, with the non-polar material at the left and the polar material (with a formal ionic charge of ±e alternating from layer to layer, where e is the magnitude of the electronic charge) at the right [see Fig. 1(a)]. The generalization for other configurations is straightforward, changing the appropriate signs when required.
Within the modern theory of polarization, we can compute the "formal" bulk polarization from the positions of the atomic nuclei and the center of localized Wannier functions. This decomposition of the charge (nuclear and electronic) into localized contributions allows for a simple classical interpretation of the bulk polarization in terms of a point charge model, and rescue the Clausius-Mossotti formulation.
In the perfectly ideal structure without rumpling, where an atomically sharp junction in the absence of defects is supposed, all the atoms at a given layer lie at the same plane [ Fig. 1(a)]. Then, we can always choose unit cells that tile the crystal under appropriate primitive translations, and that leave the left-over interface region charge neutral. 10 (It is important to note that, for the moment, we are assuming that the thickness of the polar layer is below the critical thickness for the formation of the 2DEG.) Then, the magnitude of the dipole of an individual (AO)-(BO 2 ) unit in this material is d = ea/2, where a is the out-of-plane lattice constant. The sign of the dipole is always directed from the negative to the positive layer, so in the considered configuration the sign is negative [see Fig. 1(a)]. This dipole corresponds to a "built-in formal" polarization (calculated for the previous choice of unit cell for the primitive basis of atoms and Wannier functions) of where Ω is the volume of the unit cell of the polar material and S is the cell surface. Analogously, since the layers are formally charged neutral in the non-polar material, P 0 non−polar = 0. Now, we can wonder what would happen if the materials that constitute the interface are ferroelectric, with a non-vanishing ferroelectric contribution to the polarization P S (this might be also the case when a compressive epitaxial strain is applied to the SrTiO 3 /LaAlO 3 interface 29 ). In this situation, the ferroelectric polarization must be added to the "built-in" polarization of the polar layer. 25 The difference in the polarization at the interface between the two materials produces a surface density of bound charge, Classical electromagnetic theory (Gauss's theorem) teaches us that this sheet of interfacial charge gives rise to a change in the macroscopic electric field in the two materials of where 0 is the dielectric permittivity of vacuum [ Fig. 1(b)]. The exact magnitude of the fields depends on the electrostatic boundary conditions, extremely linked with the geometry of the simulation boxes used in the computations (see Supplemental Material of Ref. 28 for a complete review). The field might induce strong structural changes in the materials and polarizes the electronic Wannier functions. Both facts translate into the development of a field-induced polarization ∆P [ Fig. 1(c)] that tends to screen the discontinuity of the total polarization P , and, consequently, of the macroscopic electric field. At the end, a self-consistent solution of Eqs. (2)-(4) is achieved where the polarizations in the two materials are in equilibrium with the corresponding macroscopic electric fields, and the total energy of the system is minimized.
For a sufficiently thin polar layer thickness, before the electronic reconstruction takes place, the absence of free charge at the interface requires the normal component of the electric displacement field D to be preserved, The finite electric displacement is an input parameter of the model. In a first-principles simulations its value can be set by hand using the virtual crystal approximation to introduce external fractional charges in the surface atoms layers, while constraining the macroscopic electric field to be strictly zero in the vacuum region. 28 Other authors fix the atomic positions on the surface unit cell to some specific values. 29 From the definition of the electric displacement field, and assuming that the polar material behaves as a linear dielectric of susceptibility χ polar around the spontaneous polarization structure, then and Eq. (6) transforms into where polar = 1 + χ polar is the dielectric constant of the polar material. From Eqs. (5) and (8), This electric field tilts the electronic bands of the polar layer [ Fig. 1(d)]. At a given critical thickness, t c , the top of the valence band of the polar material reaches the level of the bottom of the conduction bands. 30 Beyond t c , a Zener breakdown takes place, with the concomitant transfer of charge from the surface of the polar material to the interface. The magnitude of t c can be easily computed from Eq. (9) as where ∆ is the interfacial potential step. ∆ will depend on the type of band alignment and on the particular interface (p or n). Here, according to the band alignment of Fig. 1(d) for the n interface is given by 31 where E polar gap is the band-gap of the polar material. Note that the band alignment in our interfaces, where the band gaps of BaTiO 3 and KNbO 3 are very similar, is of type II, different from that in the prototypical LaAlO 3 /SrTiO 3 case (type I). For schematic views of the types of interfaces according to the band offset see Ref. 32.
The first conclusion that can be drawn from Eq. (10) is that the polarizability of the polar layer is essential to determine the critical thickness for the formation of the 2DEG (it is directly proportional to polar ). The role played by the polar distortions to avoid the polar catastrophe in LaAlO 3 /SrTiO 3 interfaces has been confirmed by first-principles simulations. 33 The importance of the extra screening due to lattice relaxations has also been discussed in the strongly related LaTiO 3 /SrTiO 3 interfaces. [34][35][36] Indeed, the success for explaining the critical thickness for the formation of the 2DEG, 17,20,26 together with the electrostrictive effect on the polar LaAlO 3 films, 37 are between the most important achievements of the polar discontinuity model.
The second conclusion is related with the relationship between D and the ferroelectric polarization in the nonpolar layer. In many works, the simulations try to reproduce the behavior of a thin polar layer on top of a thick non-polar substrate. In such cases, the macroscopic field in the non-polar materials is forced to be zero, either by symmetry or by using a dipole correction in vacuum. As a consequence, D = P non−polar . If the ferroelectric contribution to the polarization in the nonpolar layer points in the same direction as the "built-in" polarization in the polar one, then it contributes to increase the critical thickness. This has been proven in Ref. 29 for the case of a SrTiO 3 /LaAlO 3 interface subject to epitaxial strain [Eq. (10) is equivalent to Eq. (4) in Ref. 29 taking into account that, in this particular system P S polar = P S LaAlO3 = 0 and D = P S SrTiO3 ]. Finally, the third conclusion is that there is also a strong influence of an eventual spontaneous polarization of the polar material in the value of t c . In particular, if P S polar is close in magnitude to P 0 polar and points in the opposite direction, so the "built-in" polarization can be almost compensated by the ferroelectric contribution to the polarization, a large cancellation of the term in parentheses in the denominator of Eq. (10) is produced, with the concomitant increase in the critical thickness. In a previous work, Murray and Vanderbilt 38 have estimated how in SrTiO 3 /KNbO 3 superlattices the system would not become metallic until the number of layers of KNbO 3 is larger than 32.
To further validate Eq. (10) in an isolated slab geometry, we have carried out simulations on BaTiO 3 /KNbO 3 interfaces. The motivation for this choice is four fold: (i) under appropriate compressive in-plane strains, KNbO 3 is a ferroelectric polar material, with the spontaneous polarization pointing along the [001] direction, and with a "built-in" polarization of P 0 polar = P 0 KNbO3 = e/2S ≈ 53 µC/cm 2 (value computed at the theoretical in-plane lattice constant of an hypothetical SrTiO 3 substrate, a = 3.874Å), (ii) under this mechanical boundary condition, both KNbO 3 and BaTiO 3 can be stabilized with the same tetragonal P 4mm symmetry, (iii) the theoretical spontaneous polarization of our tetragonal KNbO 3 in bulk is P S polar = P S KNbO3 = 48 µC/cm 2 , close to the built-in polarization, and (iv) we can also test to which extent the ferroelectric contribution to the polarization in the BaTiO 3 layer is dominated by the imposed value of the displacement field in the simulations. First-principles results and comparison with the previous model will be presented in Sec. IV.

III. COMPUTATIONAL DETAILS
We have carried out density functional first-principles simulations based on a numerical atomic orbital method as implemented in the Siesta code. 39 All the calculations have been carried out within the local density approximation (LDA), using the Perdew and Zunger 40 parametrization of the Ceperley and Alder functional 41 to simulate the electronic exchange and correlation. This choice avoids the systematic overestimation of the ferroelectric character 42 of perovskite oxides found in other commonly used functionals 43 based on the generalized gradient approximation. This is an important point in this study, since the dielectric properties of the oxides at the bulk level might determine the behavior of the interfaces, in particular a tendency for "overscreening" of the 2DEG when the ferroelectric properties are favored. 28 Core electrons were replaced by ab-initio norm conserving pseudopotentials, generated using the Troullier-Martins scheme, 44 in the Kleinman-Bylander fully nonlocal separable representation. 45 Due to the large overlap between the semicore and valence states, the semicore 3s and 3p electrons of Ti, 3s and 3p electrons of K, 4s and 4p electrons of Nb, and 5s and 5p electrons of Ba were considered as valence electrons and explicitly included in the simulations. K, Ti, Nb and Ba pseudopotentials were generated scalar relativistically. The reference configuration and cutoff radii for each angular momentum shell for the pseudopotentials used in this work can be found in Ref. 46 for Ba, Ti, and O, and in Table I for K and Nb.
The one-electron Kohn-Sham eigenstates were expanded in a basis of strictly localized numerical atomic orbitals. 47, 48 We used a single-ζ basis set for the semicore states of K, Ti, Nb, and Ba and double-ζ plus polarization for the valence states of all the atoms. For K (Ba), an extra shell of 3d (5d) orbitals was added. All the parameters that define the shape and range of the basis functions were obtained by a variational optimization of the energy 49 in bulk cubic BaTiO 3 (for Ba, Ti, and O), and of the entalphy 50 (with a pressure P = 0.03 GPa) in bulk cubic KNbO 3 (for K and Nb, the basis set of O was frozen to that obtained in BaTiO 3 ).
The electronic density, Hartree, and exchange correlation potentials, as well as the corresponding matrix elements between the basis orbitals, were calculated in a uniform real space grid. An equivalent plane wave cut-off of 1200 Ry was used to represent the charge density. For the Brillouin zone integrations we use a Monkhorst-Pack sampling 51 equivalent to 12 × 12 × 12 in a five atom perovskite unit cell.
To avoid the problem of artificially charging the interface with the use of non-stoichiometric superlattices with periodic boundary conditions, we follow the proposal of Lee and Demkov. 26 Within this approach, the calculations were performed on vacuum-terminated (KNbO 3 ) m /(BaTiO 3 ) l /(KNbO 3 ) m slabs, where the number of KNbO 3 cells, m, is always an integer number to guarantee the electrical neutrality of the interface. In particular, we have focused on two kinds of systems: symmetric slabs where both interfaces between BaTiO 3 and KNbO 3 are of the same kind (either p-type TiO 2 /KO, or n-type BaO/NbO 2 interfaces), and asymmetric interfaces where one interface is TiO 2 /KO and the other BaO/NbO 2 . In both cases we have relaxed structures containing a different number of KNbO 3 and BaTiO 3 unit cells, given by the subscripts m and l, respectively. In the asymmetric interface a dipole slab correction is used to guarantee that the electric field in vacuum vanishes.
The in-plane lattice constant was fixed to the theoretical one of an hypothetical SrTiO 3 substrate (a = 3.874 A). Under this constraint, both BaTiO 3 and KNbO 3 are under conditions of compressive epitaxial strains. This will be relevant for the discussion of Sec. IV B.
Starting from ideal reference structures, built by piling up the corresponding unit cells of bulk strained materials without rumpling, the atomic coordinates are relaxed until the maximum component of the force on any atom was smaller than 0.01 eV/Å. In the symmetric slabs, the minimization is performed in a two step process: First, mirror symmetry planes are imposed on the central layer of BaTiO 3 to avoid the development of a polarization in any direction. Then, symmetry is broken displacing coherently the cations by hand, and a second relaxation is carried out without any imposed symmetry. In the asymmetric slabs, the relaxations are carried out in a single step, since the symmetry is broken directly in the initial reference structure.
To stablish the notation, we will call the plane parallel to the interface the (x, y) plane, whereas the perpendicular direction will be referred to as the z-axis.

A. Out-of-plane lattice relaxations and screening
In order to characterize the atomic displacements induced by the relaxation, we define the "out-ofplane" rumpling parameter along z of layer i as η z i = [z(M i ) − z(O i )] /2, where z(M i ) and z(O i ) are, respectively, the z coordinates of the cations and the oxygens at a given layer i [ Fig. 2(a)]. We also define the "in-plane" rumplings (η x , η y ), in an equivalent way, as represented TABLE I. Reference configuration and cutoff radii of the pseudopotential used in our study. Because of the inclusion of the semicore states in valence within the Troullier-Martin scheme, K and Nb pseudopotentials must be generated for ionic configurations (ionic charge of +1). However, these are more suitable than the neutral ones, given the oxidation numbers of these atoms in the perovskites. Units in Bohr.  Fig. 2(b). In all the studied systems a large out-of-plane rumpling is observed within the KNbO 3 layers, with a magnitude that is essentially independent of their thickness and the kind of interface: symmetric p [ Fig. 3(a)], symmetric n [ Fig. 3(b)] or asymmetric [ Fig. 3(c)]. Two oxide layers away from the interface, the layer-by-layer rumpling converges to a rather uniform sawtooth pattern, with values slightly larger than those observed in bulk KNbO 3 under the same epitaxial conditions. This fact is consistent with P 0 KNbO3 > P S KNbO3 , with ∆P KNbO3 tending to compensate for the difference in order to screen the polarization discontinuity at the interface. Only in the neighborhood of the free surface, a small deviation from this trend is obtained due to the larger relaxations on the surface atoms.
The dipole slab correction ensures that the displacement field in vacuum vanishes, D = 0. Due to the absence of surface external charges in the simulations, this value is preserved at the interfaces. This implies that a polarization in the BaTiO 3 layer induces a depolarizing field that is responsible for a large electrostatic energy. Therefore, under this electrical boundary conditions, no polarization is expected on BaTiO 3 . Both facts are well reproduced in our simulations, where we observe that the out-of-plane polarization vanishes within the BaTiO 3 layer in all cases. This happens even when the mirror symmetries are lifted by hand (in the symmetric interfaces), or spontaneously (in the asymmetric slabs).
According to the model developed in Sec. II, the polarization in the polar layer tends to screen the discontinuity of the polarization at the interface and, therefore, its sign opposes that of the "formal" polarization. The lattice screening avoids the development of an electric field in the polar material that would result in the tilting of its bands and, for sufficiently large thicknesses, to a Zener breakdown and accumulation of charge at the interfaces. The lattice screening provided by KNbO 3 is much stronger than the one anticipated in LaAlO 3 , where the structural distortion sustains the insulating behavior up to only 5 overlayers of LaAlO 3 on SrTiO 3 . 33 As a consequence, in our simulations all the computed structures are insulating. The reason behind this is that the LaAlO 3 is a wide band gap insulator with a low dielectric constant ( r = 25) and no ferroelectric instability (P S LaAlO3 =0). It costs some energy to polarize it. On the contrary, KNbO 3 is a ferroelectric oxide that polarizes spontaneously and contributes to reduce the field [making smaller the numerator of Eq. (9)], and increase the critical thickness for the formation of the 2DEG [making larger the denominator of Eq. (10)]. The reduction of the internal field within KNbO 3 can be directly checked from the nanosmoothed 52,53 electrostatic potential in the slabs. Independently of the geometry, the magnitude of the field amounts to 0.024 eV/Å [see red dashed lines in Fig. 4] to be compared with the roughly constant electric field of 0.24 eV/Å (one order of magnitude larger) found in LaAlO 3 /SrTiO 3 interfaces. 26,37 This is in very good agreement with the prediction of the electrostatic model developed in Sec. II: using a dielectric constant of bulk KNbO 3 around the ferroelectric structure under the same mechanical boundary condition of 25.0, 54 Eq. (9) yields a value of 0.023 eV/Å.
Since the interfacial potential steps between KNbO 3 and BaTiO 3 are 1.08 eV for the TiO 2 /KO interface and

B. In-plane polarization
Even though the ground state of bulk KNbO 3 and BaTiO 3 at zero temperature is rhombohedral, where the polarization displays both in-plane (x, y) and outof-plane (z) components, when a compressive in-plane strain is applied the polarization in the (x, y)-directions is strongly reduced or even suppressed. 56 In particular, when KNbO 3 or BaTiO 3 thin films are grown on top of a SrTiO 3 substrate, the in-plane polarization of these materials vanishes and the tetragonal c-phase is stabi-  lized. 57,58 However, contrary to current thought, in our calculations of the BaTiO 3 /KNbO 3 interfaces, we can observe in Fig. 3 the appearance of a moderate in-plane polarization in these systems along the [110] direction. The layer-by-layer in-plane rumpling profile plotted in Fig. 3 reveals that the effect is highly localized at the ptype TiO 2 /KO interface, quickly decaying upon moving into BaTiO 3 , and it is completely absent at the n-type BaO/NbO 2 interfaces.
The origin of the hidden interfacial in-plane polarization can be easily understood with a simple electrostatic model based on formally charged ions. In Fig. 5  interfaces. The cleavage of bulk BaTiO 3 and KNbO 3 to form the p-type TiO 2 /KO interface is accompanied by a change in the local electrostatic potential felt by the Ti cations. At the interface, some of the Ba cations (nominal charge +2) in the first neighbour atomic layer are replaced by K cations (nominal charge +1). This implies that at that interface the in-plane Ti cation movement is less constrained than in bulk, due to the reduced repulsions with K ions with respect to Ba ones. Analogously, in the n-type interface the electrostatic potential felt by the Nb atoms is also altered. In this case the K +1 cations are replaced by Ba +2 ions, leading to an enhanced repulsion that hinders the appearance of any in-plane polarization. As this effect is directly related to the nominal charges of the ions in both interfacial materials we would expect it to be present in other similar systems and could be used to induce in-plane polarizations in other ferroelectric nanostructures. To further characterize the development of the in-plane polarization, we take the relaxed coordinates of the asymmetric slab with m = 5 and l = 8 and remove by hand the in-plane displacements. We use this new structure as a reference configuration. Then, we compute the in-plane distortions required to go from the reference configuration to the relaxed structure and decompose into their x and y components. Finally, given fractions of these distortions are frozen in on top of the reference structure. The energy surface obtained for different fractions of the in-plane rumplings is represented in Fig. 6. From this energy landscape we can extract (i) the stabilization energy (the energy difference between the polar state and the reference one), that amounts to 43.5 meV/slab, and (ii) the energy barriers that prevents the system to rotate the in-plane polarization from [110] direction to any of the other symmetry-equivalent positions, passing through the transition state at [100] positions. This is only 16 meV, slightly larger than those found in artificial Ruddlesden-Popper-type superlattices. 59 To gauge the magnitude of these barriers, we compare them with the rotation zero-point-energy (ZPE) in this slab. This value is obtained through the nuclear motion Hamilto-nianĤ where V (η x , η y ) is the energy represented in Fig. 6 and M * can be shown to be by writing the kinetic energy operator expressed in terms of atomic coordinates and masses (M i ) using the linear trasformation relating the effective modes η x and η y with the atomic displacements through the coefficients c i . To solve the Schrödinger equation associated to Eq. (12) and find the vibrational levels associated to Fig. 6 we follow the recipe given in Ref. 60. Then, we compute the ZPE by comparing the position of the first level to the minimum of the energy surface. The resulting value for the ZPE, 13.6 meV, is much smaller than the global stabilization energy but comparable to the height of the rotational barrier, indicating that quantum fluctuations will be important in these nanostructures. These fluctuations prevent the observation of any in-plane spontaneous polarization, since the system will be delocalized over the four equivalent minima in [110] directions in a similar way to what happens in a dynamic Jahn-Teller problem. 60 This case can be compared with what happens in an incipient ferroelectric. In the latter case quantum fluctuations make the ZPE larger than the double well stabilization energy making the maximum of the probability density associated to the distortion (and the polarization) to be localized at the centrosymmetric state (origin in Fig. 6). On the other hand, in the present case, the maximum of the probability density corresponds with a nonzero polarization region around the origin (see dashed line in Fig. 6) but directional averaging results in a null net polarization. However, the coherent dynamics between the wells could be disrupted by small electric fields along the plane, that would induce large changes in the directionality of the in-plane polarization. The signature of such potential energy surface for polarization rotation would be a high dielectric constant. 59

V. CONCLUSIONS
In summary, using accurate first-principles simulations we have studied the influence of the ferroelectric polarization of the polar layer in the formation of 2DEG at BaTiO 3 /KNbO 3 interfaces. The most important conclusions that can be drawn are: (i) the spontaneous polarization of the polar KNbO 3 layer cancels out almost exactly the "built-in" polarization discontinuity at the interface; (ii) as a consequence of this compensation, the critical thickness for the formation of the 2DEG is estimated to be between 42 and 44 unit cells of KNbO 3 , one order of magnitude larger than in SrTiO 3 /LaAlO 3 interfaces; (iii) this behavior can be easily explained in terms of a simple model based on the modern theory of polarization and basic electrostatics; and (iv) surprisingly, BaTiO 3 displays an in-plane component of the polarization at the p-type TiO 2 /KO interface, even when BaTiO 3 is under in-plane compressive strains. However we do not expect this in-plane polarization to be experimentally observable as the barriers for its rotation are very small and quantum fluctuations will prevent it from being localized in a particular direction. This situation could be easily modified by the application of small electric fields that will break the symmetry and reveal the hidden polarization.
We thank Massimiliano Stengel for useful discussions and a critical reading of the manuscript. This work was supported by the Spanish Ministery of Science and Innovation through the MICINN Grant FIS2009-12721-C04-02, by the Spanish Ministry of Education through the FPU fellowship AP2006-02958 (PAP), and by the European Union through the project EC-FP7, Grant No. CP-FP 228989-2 "OxIDes". The authors thankfully acknowledge the computer resources, technical expertise and assistance provided by the Red Española de Supercomputación. Calculations were also performed at the ATC group of the University of Cantabria.