Nanoscale Smoothing and the Analysis of Interfacial Charge and Dipolar Densities

The interface properties of interest in multilayers include interfacial charge densities, dipole densities, band offsets, and screening-lengths, among others. Most such properties are inaccesible to direct measurements, but are key to understanding the physics of the multilayers. They are contained within first-principles electronic structure computations but are buried within the vast amount of quantitative information those computations generate. Thus far, they have been extracted from the numerical data by heuristic nanosmoothing procedures which do not necessarily provide results independent of the smoothing process. In the present paper we develop the theory of nanosmoothing, establishing procedures for both unpolarized and polarized systems which yield interfacial charge and dipole densities and band offsets invariant to the details of the smoothing procedures when the criteria we have established are met. We show also that dipolar charge densities, i. e. the densities of charge transferred across the interface, and screening lengths are not invariant. We illustrate our procedure with a toy model in which real, transversely averaged charge densities are replaced by sums of Gaussians.


I. INTRODUCTION
One of the central problems of physics from the mid nineteenth century on has been how to make the transition from a microscopic to a macroscopic theory of matter. This problem has two aspects. The first is the task of deriving the equations governing the macroscopic behaviour of matter from the underlying microscopic equations, exemplified by the derivation of Maxwell's macroscopic equations from the microscopic theory of charges and fields in vacuum. This problem is elegantly solved by a coarse-graining procedure which takes an average over "physically infinitesimal" regions. Such an elementary region is chosen to be small enough to let the average of a quantity follow all the changes that are observable at the macroscopic level, but large enough compared with characteristic atomic dimensions for it to contain so many particles that the behaviour of an individual particle has a negligible effect on the average quantity. This coarse-graining procedure smooths over the atomic-scale fluctuations in physical quantities, leaving only the slow spatial variation of their macroscopic components. A beautifully clear derivation of the macroscopic Maxwell's equations, first derived by H. A. Lorentz in 1902, 1 can be found in Rosenfeld's Theory of Electrons, a now forgotten classic. 2 In the coarse-graining procedure there are three distance scales: λ 1 , the scale on which the macroscopic quantities vary; λ 2 , the scale on which the smoothing is carried out; and λ 3 , the microscopic scale, that is, the atomic scale. For the procedure to work, λ 1 must be sufficiently larger than λ 3 that the pair of inequalities λ 1 >> λ 2 >> λ 3 can both be satisfied. Indeed, this criterion makes clear the distinction between macroscopic and microscopic.
The second aspect is the task of deriving macroscopic constitutive equations from microscopic properties. An early example is the derivation of the Claussius-Mossotti equation 3,4 which provides the link between the microscopic polarizability (response of the atoms or molecules to the local electric field) and the macroscopic dielectric constant.
The systems of interest in the present paper are those with interfaces between quite different materials. A planar capacitor comprised of an insulating layer sandwiched between metallic electrodes is a good example of such a system. Further examples are heterojunctions between different semiconductors and Schottky barriers at semiconductor-metal interfaces. In all such systems, there are multiple causes of charge transfer across or to the interfaces. These can include charge transfer to establish spatial uniformity of the chemical potential (the Fermi level), charge transfer in response to the local change in chemical composition across the interface and to the interface-induced atomic relaxation of the structure, charge accumulation to screen the interface charge density associated with the termination of bulk polarization at the interface, and charge accumulation attendant to charging or shorting 5,6 of a capacitor. Interface dipole densities arise from such transfer of charges and are responsible for offsets of the average electrostatic potentials across the interfaces, a dominant factor in determining Schottky barriers and valence and conduction band offsets in semiconductor heterojunctions. 7,8 In all such cases, the charge density of each material is perturbed, and the perturbation is localized near the interfaces, usually within a few interatomic distances.
The case of thin ferroelectric films between metallic electrodes provides a good illustrative example of this general class of systems. Since the early seventies a phenomenological model has been developed 9,10,11,12 to explain the modification of the polar phases (substantial reduction of the spontaneous polarization for small thicknesses, or even the complete suppression of ferroelectricity below a certain critical thickness) and of their thermodynamic properties (depression of the transition temperature with respect to that of the bulk material). The model, mainly due to Batra and coworkers, relies on three basic assumptions: (i) the polarization charge lies in a sheet right at the interface, (ii) the surface polarization charge density equals the magnitude of the polarization inside the thin film, and (iii) the free compensation charge spreads out at least over a finite distance λ within the electrode, decaying exponentially towards its interior as in the Thomas-Fermi approximation. In this model, the screening length λ is dependent only on intrinsic properties of the electrode, such as the density of free carriers or the dielectric constant. All effects are neglected which might come from a particular choice of the electrode/ferroelectric interface, such as the different chemical bondings formed at the junction or the interpenetration of the electrode and dielectric/ferroelectric wave functions that might screen the polarization charge in part within the insulator, reducing therefore the magnitude of the interface dipole density. Atomic level charge fluctuations are neglected as well, implying a smoothing on a microscopic scale distinct from the smoothing on the coarse-graining scale in the derivation of macroscopic equations and properties.
The need to go beyond such simple models is well illustrated by this interpenetration of wave functions across the interface. As pointed out first by Heine, 13 and later by Tejedor et al. 14 and Tersoff,15 the bulk Bloch-states of the metal with energies below the Fermi level of the metal and within the semiconductor band-gap and its valence band decay exponentially inside the semiconductor (and, indeed, might have a significant amplitude for a few layers from the interface), creating a continuum of gap states [the so-called metal-induced gap states (MIGS)]. Achieving deep understanding of interface properties with quantitative predictive power and free of adjustable parameters when they are determined by such effects occurring at the atomic scale requires first-principles computations. In recent years it has become possible to carry out first-principles calculations for systems of the complexity of those under discussion here. These simulations provide a wealth of information at the atomic level about the structural and electronic properties of materials and their responses to various external perturbations. 16 Some quantities, such as the microscopic charge density distribution ρ ( r) or the corresponding electrostatic potentials, are routinely available from first-principles calculations.
The question becomes how to extract from the immense detail provided by the first-principles computations reliable values of the physical quantities of interestinterface charge and dipole densities, screening lengths, etc. -which enter the pseudo-macroscopic models currently used. Two major difficulties arise. First, coarse graining is inapplicable because the relevant distance scale for interface properties is the atomic scale, i. e. λ 2 ∼ λ 3 . Second, the chargedensity changes associated with maintaining the constancy of the Fermi levels can be orders of magnitude smaller than the unperturbed bulk charge densities, themselves very rapidly varying functions of position, reflecting the underlying atomic structure. The relative magnitudes of these changes in the microscopic charge densities are illustrated in Fig. 1. Moreover, the polarization-induced charge densities and their screening charge densities can be smaller by additional orders of magnitude than those arising from imposing Fermi-level constancy. Therefore, all the interface-related dipole densities are overwhelmed by much larger variations of the total microscopic charge density.
It is still possible to carry out smoothing of the computed charge density at the nanoscale where the conditions for coarse-graining are not met, providing sufficient care is taken. A heuristic smoothing procedure has been introduced 8,17,18 to extract the quantities of interest from the results of first-principles charge-density computations as described in more detail below. However, that procedure has not yet been systematically analyzed to establish the conditions under which it accurately extracts the quantities of interest: surface charge densities, surface dipole densities, etc. In the present paper, we introduce such an analysis. In addition, in order to focus only on the perturbations introduced by the nanosmoothing procedure, avoiding other sources of numerical noise coming from the first-principles simulations, we illustrate the analysis with a toy model whose accuracy can be arbitrarily improved. The application of this theory to first-principles computations on realistic ferroelectric capacitors is the subject of a forthcoming paper. Despite this focus on the ferroelectric capacitor, our analysis is of general utility for the extraction of interface properties for all multilayer systems.
The rest of the paper is organized as follows. We set the grounds of our discussion of simulations of interfaces from first-principles and define the microscopic behaviour of the different quantities that are the targets of our study in Sec. II. In Sec. III we define the interface quantities of interest for an unpolarized interface. In Sec. IV we describe the difficulties encountered in defining precisely the location of a reference interface. We develope the theory of nanosmoothing of unpolarized systems in Sec. V with particular attention to questions of the sensitivity of quantities of interest to the smoothing procedure. We describe in Sec. VI a toy model used to illustrate all of the essential elements of the smoothing theory. In Sec. VII, we present the results of the toy model computations for the unpolarized case. We generalize the nanosmoothing theory of Sec. V for polarized systems in Sec. VIII   criteria smoothing functions must satisfy to yield interfacial properties insensitive to their parameters and specifying which of our results are new.

II. SIMULATION OF INTERFACES FROM FIRST PRINCIPLES.
First-principles calculations of interfaces between two materials, where there is no periodicity in at least one direction, are almost universally done by means of the supercell approximation. 19 Within this approach a basic unit cell that contains a suitable number of multiatom layers of the two materials is periodically repeated over all space (Fig.  2). For the interfaces within a nanostructured multilayer to be well defined, with properties distinct from those of the bulk-like regions between them, the widths of the layers of each material introduced in the construction of the basic unit cell must be large enough to avoid the interaction between adjacent interfaces through the bulk materials, so that the calculation accurately represents an isolated interface.
Throughout this work, we shall assume that the interface is oriented along the z axis, and each material is periodic in the plane parallel to the interface, referred to as the (x, y) plane.
Precisely the same methodology is used to treat nanostructured multilayer materials, and in this paper we do not distinguish between the cases, specializing to multilayers in which the individual material thicknesses are large enough for the interfaces to be noninteracting.
The microscopic charge densities ρ ( r) [see Fig. 1(c)] and electrostatic potentials V ( r) provided by the first-principles computations for the previously described supercells are continuous functions periodically repeated in space with the periodicity of the supercell (a sc in Fig. 2). A few interatomic distances away from the interfaces, the microscopic quantities recover their bulk features. In other words, if the layer widths are large enough, it is possible to identify "bulk-like" regions in the middle of each of the layers that constitute the superlattice, with large variations of the microscopic quantities with the same periodicity as in the bulk, reflecting the underlying atomic structure.
Many of the interface-related quantities that we shall define below depend only on what happens in the direction perpendicular to the interface, where the actual discontinuities of the physical structure occur. Since we assume inplane periodicity, we can trivially eliminate the in-plane dependence by taking a planar average of the corresponding  0000  0000  0000 0000  0000  0000  0000  0000 0000  0000  0000  0000  0000  0000 0000  0000  0000  0000  0000 0000  0000  0000  0000  0000 0000  0000  0000 0000  0000  0000  0000  0000 0000   1111  1111  1111 1111  1111  1111  1111  1111 1111  1111  1111  1111  1111  1111 1111  1111  1111  1111  1111 1111  1111  1111  1111  1111 1111  1111  1111 1111  1111  1111  1111  1111 1111   z  x   a   a   a   sc   2   1 PSfrag replacements FIG. 2: Schematic view of a supercell used in first-principles simulations of interfaces. The z axis is taken as normal to the interface, whereas the plane parallel to the discontinuity is taken as the (x, y) plane. n1 and n2 multiatomic layers of two materials, whose bulk lattice constants in the out-of-plane direction are respectively a1 and a2, are stacked to build a basic unit cell (hatched rectangle) that is periodically repeated in space. asc is the length of the supercell in the out-of-plane direction, with asc = n1a1 + n2a2. In the figure, n1 = 2, and n2 = 3.
microscopic quantity, e. g., for the electron density, where S is the area of the interface unit cell. As is proven in Appendix A, this transverse averaging has no effect on the Poisson equation, which reads after averaging

III. INTERFACE QUANTITIES OF INTEREST; THE UNPOLARIZED CASE.
One of the most important physical properties of heterojunction devices is the band offset or Schottky barrier at an interface, that is, the relative positions of the energy levels on both sides of the interface. In the case of semiconductor heterojunctions, the valence-band offset (VBO) [conduction-band offset (CBO)] is defined as the difference between the positions of the tops of the valence bands [the bottoms of the conduction bands] of the two materials. In the case of metal-semiconductor contacts, we can define the p-type (n-type) Schottky barrier as the difference between the Fermi level of the metal and the top of the valence band (bottom of the conduction band) of the semiconductor. These differences in the band positions determine the effective barrier for electron or hole transport across the junction.
The computation of such effects from first-principles cannot be achieved by a direct comparison of the corresponding single-particle energies (tops of the valence bands, bottoms of the conduction bands and/or Fermi level of the metal) in the two compounds as obtained from two independent bulk band-structure calculations. The reason is the lack of an intrinsic energy origin to which to refer all the energies: in a first-principles simulation, the Hamiltonian eigenvalues are referred to an average of the electrostatic potential that is ill defined for an infinite system 20 where, due to the long-range nature of the Coulomb interaction, it is defined only to within an arbitrary constant. Consequently, together with the eigenvalue difference, we must consider both the shift of this average between the two materials and the redefinition of the averaging process so that it is appropiate to the multilayer system under consideration. As was mentioned in Sec. I, the coarse-graining procedure used historically fails when applied to the results of first-principles computations. We designate the averaged potential as V u and its shift as ∆ V u . The brackets indicate that the averaging process is not yet defined. We shall label all quantities of physical interest of the unpolarized system by a subscript u. This potential shift depends on the dipole induced by the electronic charge transferred from one side of the interface to the other after interfacial hybridization. As the charge transfer depends not only on the materials that constitute the interface, but also on intrinsic interface effects such as the chemical composition (termination of each material at the interface), on particular orientation and on other structural details, the shift can only be obtained from a self-consistent calculation on a supercell including both materials. This ensures that the averaged electrostatic potentials of both materials in the "bulk-like" regions defined in the previous section, where all the physical quantities recover the bulk features, are expressed with respect the same reference and allows a direct extraction of the shift.
This shift of the averaged electrostatic potential should be directly related to an averaged interface dipole-moment density p u , We prove in Appendix C that Eq. (3) holds for the specific definition of the average procedure introduced in Sec. V A below and discussed in this context in Sec. V B. From a fundamental point of view, 21 the charge Q u , and electric dipole moment at the interface p u are defined respectively as the zero and the first moment of the total microscopic charge density, For the previous equations to be meaningful and truly represent an interface quantity, the thicknesses of the adjacent layers must be wide enough so that they contain regions within which the charge density is essentially unaffected by the presence of interfaces within which z 1 and z 2 could be located. However, from a practical point of view, such a definition poses serious questions. Indeed, both Q u and p u are ill-defined due to the large and rapid oscillations of the microscopic charge density [see Fig. 1(c)]. In fact, different choices of the integration limits z 1 and z 2 yield widely different values of the charge and interface dipole moment. Only in the extreme case of a Clausius-Mossotti model, in which the total charge is unambiguously decomposed into an assembly of localized and neutral charge distributions, so that a unit cell can be chosen with no charge at the surface, the dipole moment of a periodic charge distribution would be well-defined as the integral of the first moment of the charge density. However, any Claussius-Mossotti approach does not correspond to reality, particularly in materials where delocalized covalent charge is present. 22 This ambiguity in the definition of p u with respect to the boundaries of the region within which the dipole moment is computed is closely connected to the problem of defining the polarization of a periodic system from the charge density. 23

IV. THE DIFFICULTY OF DEFINING A REFERENCE INTERFACE.
In order to get rid of bulk effects and extract interface-related features, some authors 24,25,26,27 have defined an ideal interface by stacking alternate slabs, each of them made from slicing the planar average of the bulk charge density of the corresponding material perpendicular to a particular direction. Let us define ρ where a (s) bulk is the lattice constant of the bulk unit cell of each material in the z direction, and s refers to the side of the interface considered, 1 or 2. The bulk unit cell boundaries are chosen so as to preserve inversion symmetry.
Replicating the bulk charge densities up to the as yet unspecified interface from each side, we could define the reference charge density ρ 0 for all z as where z int is the coordinate assigned to the position of the interface. ρ 0 (z) is discontinous at the interfacial plane z int by its very definition. We now define the interface-induced deformation of the charge density ∆ρ u (z) as If the thicknesses of the two layers are wide enough, ∆ρ u (z) becomes negligibly small over the ranges R 1 in the left material and R 2 in the right material [see Fig. 6(b)]. The interface region can thus be identified as comprised of those ranges where ∆ρ u (z) differs significantly from zero, in other words where the microscopic charge density differs from the relevant bulk values. The interface charge and dipole density associated with ∆ρ u (z) are defined as, The advantage of this approach is that both ∆Q u and ∆p u are well defined quantities with respect the location of the integration limits z 1 and z 2 in Eqs. (8a)-(8b), provided that z 1 lies in R 1 and z 2 lies in R 2 . However, this approach has pitfalls. In particular (i) the position of the interface, z int in Eq. (6), is not yet specified. Some recipes have been given for how to cut the bulk slabs, but they have limited applicability. One case is for common anion heterostructures with non-relaxed interfaces along high-symmetry planes, such as the (001) 24 , (110) 25 , or (111) 26 interfaces of GaAs/AlAs superlattices. As soon as an interface-induced rippling of the atomic layers is introduced, for instance after an atomic relaxation of the interface geometry, the problem of defining the position of the interface worsens. (ii) Therefore, the interface charge and dipole densities are not unique, since they depend critically on where the mathematical surface representing the interface is chosen, and no objective criterion for locating it has been established. In particular ∆Q u and ∆p u [Eqs. (8a)-(8b)] equal Q u and p u [Eqs. (4a)-(4b)] if and only if z int − z 1 contains an integer number of unit cells of the left material and z 2 − z int contains an integer number of unit cells of the right material. The location of the interface determines where on one side, the bulk charge density of the left material is subtracted, and on the other side that of the right material. A different choice of the mathematical interface can produce very different charge and dipole densities. In addition, comparison between different interface orientations makes little sense with this definition. 7 (iii) Since the interface dipole moment is dependent on the reference charge density, the corresponding potential drop at the interface (∆ dipole in the notation of Refs. 24,25,26,27) must be too. However, it is important to note that the potential drop generated by ∆ρ u (z) is only part of the total potential shift. The total charge density of the interface is given by thus the potential shift associated with ρ 0 (z) must be included as well. Although the existence of a potential shift generated by ρ 0 (z) is general, we shall explain its origin only in the particular case where z int − z 1 contains an integer number of unit cells of the left material, and z 2 − z int contains an integer number of unit cells of the right material. In this particular situation both slabs used to construct the reference charge density in Eq. (6) have neither a charge nor a dipole moment. Under these circumstances, the potential shift is the difference in the locally averaged potentials produced by the ρ (s) 0 of each material in the regions R 1 and R 2 respectively. For this shift to be an interface property, the layer width of each material must be large enough that the local average V (s) slab has approached the unperturbed value within the center of an slab. The planar averaged potential within any point of the slab will be given by 20 so the local average can be computed as Since V (s) slab depends on the charge density distribution, it differs for the left and right slabs used in the construction of the reference charge density and produces an additional shift. Therefore, the total potential drop at the interface is the sum of the potential drop generated by ∆ρ u (z), ∆ dipole , and the difference of the average potential of the two reference slabs ∆ ref . Only the sum ∆ ref + ∆ dipole is independent of the reference charge density chosen and is a physically measurable property of the interface. Since each term in the sum is sensitive to the arbitrary location of the interface, each must be computed accurately enough for the sensitivity to disappear from the sum. As the charge density shifts of interest are so small, this is an unnecessary burden, removed by the use of a proper nanosmoothing 17 procedure as shown in Sec. V.
A procedure to eliminate charge fluctuations in the regions of the material which do not contribute to the interfacial hybridization, thereby localizing the physically relevant charge densities to the interface, consists of filtering out the periodic oscillations of microscopic quantities, which typically follow the underlying atomic structure, preserving only those features that emerge in the vicinity of the interface.
To obtain this smoothed charge density, we have followed the recipe given by Baldereschi et al. in Ref. 17 and generalized by Colombo and coworkers for lattice mismatched heterostructures in Ref. 18. Starting from the planaraveraged charge density ρ u (z), we construct the smoothed density ρ u (z) by convoluting it with a smoothing function f (z), which has the following properties: 28 so that In addition, f (z) should be monotonic in |z| and sufficiently smooth itself. L should be chosen on the scale of the unit cell length or larger, but smaller than the widths of the left and right layers. A sharper criterion for L is introduced below. We define this as the averaging procedure left unspecified above in Sec. III. The particular smoothing function we have used, following Refs. 8, 17 and 18, is the convolution of two square-wave filter functions: where Giustino and coworkers 29,30 propose convolution with a Gaussian kernel that can be an approximation to the asymptotic limit of a convolution of a large number of square-wave filter functions. This method is best suited for superlattices where crystal deviates from perfect periodicity far away from the interface so that it is not possible to define regions where the interface-induced charge density vanishes, or in disordered three-dimensional systems with short-range order. Even more general functions can be used, providing the criteria established above are met.
The explicit dependence of f (z), defined in Eq. (15), on z is: 0 as shown in Fig. 3. l > is the greater of l 1 and l 2 . After nanosmoothing the planar averaged charge density ρ u (z), the resulting charge density ρ u (z) is a continuous function that joins smoothly at the interface. Even though we have not subtracted any reference charge, the smoothed charge density becomes negligibly small over ranges R Providing the filter function f (z) satisfies the following additional conditions, The Poisson equation remains invariant after the nanosmoothing (see Appendix B) and transforms into The nanosmoothing function f (z) defined in Eq. (16a) and Eq. (16b) violates condition Eq. (18a) at its end points z = ± l 2 , where f (z) is discontinuous so that condition (18c) cannot be unambiguously applied. Similarly, the f (z) of Eq. (17) violates condition (18a) at the four points z = ± l1±l2 2 , and df /dz is discontinuous at the points z = ± l1±l2 2 so that condition (18b) cannot be unambiguously applied. Nevertheless, one can think of the f (z) as a distribution, a family of smooth functions all of which meet conditions (18a) -(18c) and which approach the f (z) of Eq. (17) as their limit. In practice, because the smoothing operation is a convolution, Eq. (12), one carries out smoothing via fast Fourier transformations. The family of the finite Fourier series involved is thus a distribution which converges to f (z) in the limit, meeting conditions (18a) -(18c) along the way. Eq. (19) holds in general for both the unpolarized and the polarized cases, so subscripts have been omitted. Otherwise the subscript u is used because we are presently treating the unpolarized case.
We prove in Appendix C that the full electrostatic potential shift ∆V u is given by the nanosmoothed dipole density p u , where C. Insensitivity of the dipole-moment density and potential shift to the smoothing function.
The microscopic charge density ρ u must return to the bulk microscopic charge densities ρ (1) 0 and ρ (2) 0 in the regions R 1 and R 2 for it to be possible to ascribe physical properties specifically to individual interfaces. However rapidly it approaches those values, the approach must be complete within R 1 and R 2 , as summarized by Eq. (22), with g (z) such that ρ u (z) and all its derivatives are continuous (the cusps at the nuclei are washed out by lateral averaging). Consequently, locally within R 1 and R 2 , the microscopic charge density ρ u can be represented by the Fourier transforms of ρ (1) 0 and ρ In Eq. (23) A (s) 0 vanishes because of charge neutrality. The smoothing function f (z − z ′ ), Eq. (15), is a convolution of two square wave functions ω ls (z − z ′ ), Eqs. (16a)-(16b), s = 1, 2. The order of the ω ls within the convolution is immaterial, so the ω ls can be applied to the nanosmoothing of ρ Similarly the electrostatic potential V s (z) can be expressed within R s as a comparable Fourier series, with the B  2), which is influenced by the charge density outside of R ′ s . Upon nanosmoothing, all contributions to V u (z) for z within R ′ s vanish except that for n = 0, which is invariant to the smoothing process. The potential shift ∆V u , defined in Eq. (C2) as is thus invariant to the smoothing procedure, Moreover, according to Eq. (20), the nanosmoothed dipole-moment density p u is invariant as well, D. The transferred charge density and the dipolar density.
It is of considerable physical interest to establish the value of the charge transferred across the interface, a difficult task. A criterion for establishing the position of the interface is needed, the difficulty of which is discussed in Sec. IV. The rapid, large oscillations of ρ u (z) and the relative smallness of the pertinent component of ρ u (z) make using it impractical. On the other hand, if one uses a criterion based on ρ u (z), z int can be sensitive to the smoothing function. Nevertheless, we shall attack the problem using ρ u (z) and attempt to overcome the resulting sensitivity to the smoothing function of the position of z int and the amount of charge transferred.
We start by defining two cumulative charge densities For the unpolarized case now under consideration, Q u = 0. Thus, as holds for ∀z ∈ (z 1 , z 2 ), and it is sufficient to consider either one or the other. Define q as the magnitude of the charge transferred per unit area of the interface, the transferred charge density. We estimate q as and estimate z int as Now both q and z int are sensitive to the choices of l 1 and l 2 in the smoothing function f . As l 1 and l 2 increase, ω l1 can reach across the interface from material 1 to material 2 bringing contributions from ρ(z ′ ), z ′ within material 2, to ρ(z), z within material 1, and vice versa, thus returning part of the transferred charge back to its origin and reducing the value of q. Moreover, since ρ (z) contains components which oscillate strongly with z, Q (z) could develope multiple suprema or maxima as l 1 and l 2 increase in multiples of the lattice constants of materials 1 and 2, respectively. This would vitiate the utility of the definitions (33) and (34) of the transferred charge q and the interface location z int , respectively, should it happen. We have found that it does happen in the toy model described in Sec. VI and studied in Sec. VII, cf. Fig. 8 below, in the case where the interatomic distance remains unchanged across the entire superlattice. Accordingly, as a precaution, the smallest acceptable values of l 1 and l 2 should be used for f (z), a single lattice constant of each material, an important additional condition on l 1 and l 2 .
If the transferred charge density were concentrated equally on two surfaces at either side of the interface, separated by a distance λ, a dipole moment density of magnitude qλ would be created. Setting qλ equal to the actual dipolar density p u allows us to define λ u as the dipolar length where we have restored the subscript u to q as we are dealing with the unpolarized case.
E. Loss of invariance of physical magnitudes of interest with nanosmoothing.
Now we can ask whether the physical magnitudes of interest, such as the interfacial charge [Eq. (4a)] or dipole moment densities [Eq. (4b)] remain unchanged if the microscopic charge density is replaced with the nanosmoothed charge density. In other words, if we define the interfacial charge Q u and dipolar densities p u computed from ρ u (z) as then the question is whether Q u = Q u and p u = p u for given integration limits z 1 and z 2 . Replacing Eq. (14) into Eq. (36a),  The region of integration in Eq. (37) within the z, z ′ plane is shaded in Fig. 4(a). Then, the integral of Eq. (37) can be decomposed into the integral on the central square plus the integrals on the cross hatched triangles On the other hand, Q is defined as where the region of integration in the z, z ′ plane is now shaded in Fig. 4(b). Decomposing the domain of integration as for Q u above gives Therefore, the equality Q u = Q u is verified if and only if the integrals over the shaded regions on either side of the line z = z ′ are the same, that is if the following pair of equations hold If in Eq. (41a) we apply the following change of variables z = z 2 − u and z ′ = z 2 + v in the left hand side and z = z 2 + u and z ′ = z 2 − v in the right hand side, then Eq. (41a) transforms into Due to the parity conditions of the filter function, we know that if and only if ρ is even about z 2 . Completely analogous reasoning can be applied to Eq. (41b) and is omitted here. Thus, for Q u to equal Q u , it must be possible to find a z 1 and a z 2 about which ρ(z) is symmetric for |z − z 1,2 | ≤ L. We now show for the unpolarized case that this symmetry condition can be satisfied. For a multilattice consisting on alternating layers of two different materials s, s = 1 or 2, there are two distinct interfaces within the supercell bounded by z 1 and z 1 + a sc , interface 1,2 and interface 2,1. We define total charge Q sc and Q sc which are the sums of the charges associated with each individual interface, Proceeding in analogy with Eq. (38) and Eq. (40), we obtain The fact that ρ u (z) is periodic in z with period a sc allows as to rewrite Eqs. (44) and (45) so as to establish the equality of Q sc and Q sc , Since the supercell is electrically neutral, so must the nanosmoothed supercell be, Consequently, from Eq. (43b) it follows that implying that there would be a smooth electrostatic field within each layer if Q 12 is nonzero. The existence of such a field would polarize the system in contradiction to the initial condition that the system is unpolarized. We conclude that Thus, for the interface charge to be invariant to nanosmoothing, that is, for Q 12 = Q 12 to hold, z 1 and z 2 must be positioned in R 1 and R 2 so that Q 12 vanishes in the unpolarized case. To do this, one could make an arbitrary choice of z 1 in R 1 , say, and then integrate ρ u (z) from z 1 up to some z 2 in R 2 at which the integral vanishes. There is no need to do this, as it is the nanosmoothing quantities themselves which are of interest.
Repeating the reasoning for the dipole moment density, we arrive at the conclusion that for p u = p u , the following condition must be satisfied Applying the same change of variables as before, that is z = z 2 − u and z ′ = z 2 + v in the left hand side and z = z 2 + u and z ′ = z 2 − v in the right hand side, then Eq. (50) transforms into Even in the case of a function ρ u symmetric around z 1 and z 2 , the previous condition does not hold in general, and the difference between p u and p u amounts to plus a similar term that comes from the difference in the lower triangles in Fig. 4(a) and Fig. 4(b). This has to be evaluated for each particular case.
In conclusion, the interface dipole-moment density is not invariant to nanosmoothing and the interface charge density can be made so only by exquisite case in the choice of z 1 and z 2 . This is of no concern, as it is the nanosmoothed quantites which are of physical interest.

VI. DESCRIPTION OF THE TOY MODEL.
As highlighted in the introduction, the components of the density which give rise to the interface-related dipole densities are nearly obscured by the much larger variations of the total microscopic charge density. This atomicscale charge density, routinely provided by any density-functional-based first-principles code, is affected by numerical noise and convergence problems inherent in some of the standard approximations in the practical implementations of density functional theory (DFT). Therefore, the accuracy of the computations required for extracting the actual charge transferred from one side of the interface to the other must be high enough so that the numerical noise of the calculations is orders of magnitude smaller than the relevant interface-related charge densities. In this paper, in order to illustrate all of the essential elements of the theory of smoothing while avoiding these practical problems, we shall define a toy interface model that resembles closely a realistic multilayer material but whose computational accuracy can be systematically improved.
The requirements that such a toy model should meet are: (i) its electron density must be a continuous function; (ii) far away from the interfaces, where the interface-induced perturbation of the charge density becomes negligible, the toy electron density must tend to two distinct periodic functions on the left and on the right of each interface, mimicking the differing behaviour at the bulk level of the materials that constitute the multilayer system containing the interfaces; and (iii) the interlayer spacing at the interfaces should be distorted with respect those at bulk so as to simulate better the interface induced relaxations that happen in real interfaces.
In the toy model we propose here, we represent only the laterally averaged density, a one-dimensional function. As before, the direction perpendicular to the interface is referred to as the z axis.
We shall consider atomic-like charge densities g (s) i (again as before, the superindex s = {1, 2} refers to the side of the interface, left or right, where a given "atom" i is located) of the sum of two gaussians, centered at each atomic site z The positive (negative) gaussian, whose standard deviation is denoted by σ n (σ e ), mimics the nuclear charge density (electronic charge density) of atom i centered on position z i . By imposing σ n < σ e we make the "nuclear" charge more confined than the "electronic" charge. Since the gaussians are normalized, a net charge per atomic site can be simulated by making A n,s i = A e,s i . The parameter δ (s) i allows us to displace the electronic clouds with respect to the nuclei and thereby produce a net dipole moment on a particular atom.
The atomic-like charge densities are arranged in bulk unit cells that might contain one single atom or a more complicated polyatomic basis. In this Section we shall assume that the bulk unit cells of each of the materials that form the superlattice contains a single atom that does not carry any charge (A n,s For the polarized case, we refer the reader to Sec. IX. Thus, for the unpolarized interface each material can be considered as a one-dimensional monoatomic chain, where consecutive "atoms" are separated by a distance a (s) . The interatomic distance at the interface, a int , is taken as a int = a (1) +a (2) 2 . Then, consecutive interatomic distances from the interface evolve smoothly towards the bulk value as a function of the distance to the interface. In our simulations, when we move from the interface towads material s, s = {1, 2} , the second interatomic distance is set up to 1 4 is the lattice constant of the other material. The bulk value is recovered only at the third interatomic distance from the surface.
A supercell is then built, as described in Sec. II. The basic unit cell, periodically repeated in space, contains a suitable number N 1 and N 2 of bulk unit cells of the two materials. The microscopic charge density, ρ (z), is defined as the superposition of all the atomic-like charge densities The resulting model is illustrated for the non polar case (δ  PSfrag replacements In Fig. 6(a) we illustrate the specific toy model we shall analyze in detail for the unpolarized case. The parameters of the model are A n,1 = A e,1 = 1, σ n,1 = 0.5, σ e,1 = 2.0, A n,2 = A e,2 = 2, σ n,2 = 0.7, σ e,2 = 4.5, δ (1) = δ (2) = 0 for all the "atoms", a (1) = 6, and a (2) = 8. Atomic units are used throughout the paper. First, we consider a reference density ρ 0 (z) defined as in Eq. (6), locating the interface at z int half way between the rightmost atomic layer of material 1 and the leftmost atomic layer of material 2. ∆ρ u (z), defined as in Eq. (7), is plotted in Fig. 6(b). Next, we construct ρ u (z) from ρ u (z) directly using a nanosmoothing function f (z − z ′ ) defined as in Eq. (15), in which l 1 = a (1) = 6, and l 2 = a (2) = 8. The smoothed charge density ρ u (z) is shown in Fig. 6(c). As we have already discussed in Sec. V C, the regions R ′ s within which ρ u vanishes lie within the regions R s within which ∆ρ u vanishes, because smoothing ρ u within R s brings into ρ u values of ρ u for z outside R s . Due to the symmetry of the microscopic charge density at the center of each layer, the interfacial charge density between material 1 and material 2 is the mirror image of the interfacial charge density between material 2 and material 1 [see the center and the edges of the Fig. 6(b) and Fig. 6(c)].
A closer look at the interface region is shown in Fig. 7. The discontinuity of ∆ρ u at the interface plane is clearly observed. Also, we can see how ∆ρ u displays large fluctuations at the atomic scale in the neighborhood of the interface due to the interface-induced relaxations of the atomic layers. The positions of the atoms in the layers close to the interface do not coincide with the positions of the atoms after the cleaving of the bulk to define ρ 0 (see Sec. IV). Therefore, in the computation of ∆ρ u we are subtracting charge densities centered on different positions. The nanosmoothing procedure eliminates not only the contributions associated with ρ 0 in the bulk regions R 1 and R 2 , but it filters out the oscillations due to the interface induced relaxations while producing a continuous charge distribution. Nanosmoothing ρ u is clearly superior to forming ∆ρ u .
In   Table I). This striking sensitivity of the nanosmoothed charge density to the smoothing function impedes detailed physical interpretation of its features. Even though we are dealing with an unpolarized interface made of the juxtaposition of neutral and non-polar atoms, the interface charge density Q u , Eq. (4a), vanishes if and only if the integration limits z 1 and z 2 are taken midway between atoms inside each material layer, meeting the requirements for symmetry described in Sec. V E (see Table  I). ∆Q u , Eq. (8a), on the other hand, does not vanish since, with our criterion for locating the interface plane z int , we do not conform to the requirement that z int − z 1 contains an integer number of unit cells of the left material, and z 2 − z int contains an integer number of unit cells of the right material. Consequently, the integral of the reference charge density ρ 0 between z 1 and z 2 carries a net charge. On the other hand, as it should, Q u vanishes and is independent of the integration limits provided that z 1 lies in R ′ 1 and z 2 lies in R ′ 2 . We have also calculated the corresponding interface dipole densities p u , Eq. (4b); ∆p u , Eq. (8b); and p u , Eq. (36b), with the results displayed in Table I. By definition, p u and ∆p u are independent of the nanosmoothing procedure. As long as l 1 and l 2 equal an integer number of lattice constants of the bulk unit cell of the material along z, p u is insensitive to the shape and range of the smoothing function. Note that, as expected after the discussion in Sec. V E, the nanosmoothed dipole density p u differs from both ∆p u , computed from ∆ρ u (z), and p u computed from ρ u (z). The latter two also differs because ρ 0 (z) defined in Eq. (6) generates an extra dipole-moment density when integrated within our integration limits. Thus the electrostatic potential shifts ∆V u and ∆V u differ correspondingly. The shift ∆V u is the physically meaningful one because what enters in the band offsets are the local averages of the electrostatic potentials in each material, which are independent of the location of z 1 and z 2 within R shows that ∆V is the difference of the two potentials at specific points within regions in which V (z) varies rapidly. In Fig. 9 we show V (z) and V (z) and indicate the positions z 1 and z 2 . It is clear that ∆V is about 20% larger Interface charge densities  than ∆V , explaining the relation between p u and p u in Table I, the difference arising from the way V is sampled by nanosmoothing and by the selection of z 1 and z 2 .
In Table I we also show two cases where the nanosmoothing does affect the value p u . First, the value of p u departs from the correct value 2.157 when the range of the smoothing function, L in Eqs. (13a) and (13b), is of the same magnitude as the width of the layer of one of the materials, as for l 1 = 5a (1) = 30 and l 2 = 5a (2) = 40, so L = l 1 + l 2 = 70, slightly larger than the width of material 1, made of 10 layers with an interlayer distance a (1) = 6. Under these circumstances, it is not possible to define regions R ′ 1 and R ′ 2 where ρ u (z) vanishes, impeding a proper location of the integration limits z 1 and z 2 . Second, the same occurs when the range of every filter function entering Eq. (15) does not equal an integer number of lattice constants of the bulk unit cell along z, as for l 1 = a (1) 2 = 3, and l 2 = a (2) 2 = 4. In such a case, the charge density after nanosmoothing, ρ u (z), still shows large and rapid oscillations. In Table II we report the amount of charge transferred from one material to the other, the screening length and the interface position computed with the method summarized in Sec. V D for the nanosmoothed charge densities plotted in Fig. 8. The negative cumulative charges are displayed in Fig. 10. As expected the value of q is sensitive to the  Fig. 6(a). The interface charge Qu and the dipole-moment pu densities are defined respectively in Eqs. (4a) and (4b) and the integration limits z1 and z2 taken midway between the atoms at the center of each material layer. When a reference charge density of an ideal interface is subtracted, the corresponding charge ∆Qu and dipole moment ∆pu are defined as in Eqs. (8a) and (8b) with the same integration limits as before. Since the nanosmoothed charge density does not play any role in the definition of Qu, pu, ∆Qu, and ∆pu, these magnitudes are insensitive to nanosmoothing. Q u and p u are defined in Eqs. (36a) and (36b). l1 and l2 are the lengths of the square-wave functions entering in the definition of the smoothing function, Eq. (15).  smoothing function, decreasing with the increasing of its range. Accordingly, the screening length is also sensitive to the range of the smoothing function L, taking a value that is roughly the greater of l 1 and l 2 . Nevertheless, the position of the interface seems to be insensitive to the filtering function.
The results displayed in Fig. 8 for ρ(z) demonstrate that the larger the width of the smoothing function, the more complex the spatial dependence of ρ(z) and the less it resembles a simple interface charge density. The introduction of multiple extrema is artificial and caused by excessive transfer of charge back across the interface by smoothing. This back transfer of charge is manifested clearly in the values of q in Table II. Accordingly the smallest allowable widths l s = a s , should be used for nanosmoothing. The resulting values of q and λ are the best estimates which can be extracted by nanosmoothing. Now, as a final test of the insensitivity of the dipole-moment density and potential shift to the shape of the smoothing function, we shall simplify the toy model to one with a common atomic separation in both materials. This opens FIG. 10: Negative cumulative charges defined in Eq. (30a) for the nanosmoothed charge densities shown in Fig. 8. We can estimate the magnitude of the charge transferred per unit area of the interface q from the supremum of this cumulative charge and the interface position from the position where this supremum appears. Numerical values are given in Table II.  three simple options for the filter function. We can set f equal to a single ω when the filter function adopts a square shape [see Fig. 3(a)], we can set f equal to the convolution of two ω's with the same l, a triangle function, or we can set f equal to the convolution of two ω's with different l, a trapezoidal function, as in Eq. (17) and in Fig.  3(b). Comparing the results of the smoothing with the square wave, the triangle, and the trapezoidal function for various values of l allows us to establish which physical properties are insensitive both to the range and to the shape of the smoothing function. Again, as soon as l 1 and/or l 2 equal an integer number of lattice constants of the bulk unit cell of the relevant material along z, p u is insensitive to the shape (triangular, square or trapezoidal) and range of the smoothing function.

VIII. POLARIZED SYSTEMS.
Up to now we have dealt with materials unpolarized except at the interfaces. In such systems it was possible to identify bulk-like regions in the middle of each layer where the microscopic charge density is unaffected by the presence of interfaces and returns to the bulk microscopic charge density, which is centrosymmetric within a unit cell. We now turn our attention to polarized interfaces, where at least one of the materials has a non-vanishing polarization P .
Basic electrostatic arguments that can be found in any textbook 21,32 show that a non uniform polarization in a dielectric generates a volume charge density, the polarization charge ρ pol ( r), whose value at any point of space is given by III: Interfacial charge and dipole moment densities for a microscopic interface-like charge density obtained with the following parameters: A n,1 = A e,1 = 1, σn,1 = 0.5, σe,1 = 2.0, A n,2 = A e,2 = 2, σn,2 = 0.7, σe,2 = 4.5, δ (1) = δ (2) = 0 for all the "atoms", a (1) = 10, and a (2) = 10. The meaning of the symbols and integration limits is as in Table I.
Even in the case of uniformly polarized material, where P is constant and therefore its divergence vanishes inside the dielectric, the discontinuity of the polarization at the surface or interface gives rise to a net surface or interface charge density σ pol given by the familiar form wheren is a unit vector normal to the surface or interface pointing outwards.
In the case of an interface, the other material at the interface responds to this perturbation in order to minimize the electrostatic energy cost associated with the build up of the polarization charge ∇ · P at the interfaces. If the second material is a dielectric, a uniform polarization is induced within it as well. 33,34 If the second material is a metal, a screening charge is induced that spreads over a finite distance (the screening length) in the electrode, producing additional dipole layers at each polar dielectric/metal interface. If the screening of the polarization charge is not perfect, a residual depolarizing field appears inside the dielectric. 5, 6 We now face the same question as for the unpolarized case: how to extract from the results of the first-principles calculations values of the physical quantities of interest (polarization and screening charge densities, screening lengths, depolarizing fields, etc).
In parallel to what was done for the unpolarized case, we can define the net charge density Q p and dipole moment density p p associated with the microscopic charge density of the polar interface ρ p as The computation of these two quantities present the same difficulties as for the unpolarized case, summarized in Sec. III. The approach of determining a reference interface charge density ρ 0 (z) as in Eq. (6) and defining an interfaceinduced deformation of the charge density ∆ρ p (z) by subtracting ρ 0 (z) from the microscopic charge density, discussed in Sec. IV for the unpolarized case, still allows for the identification of two regions R 1 and R 2 where ∆ρ p (z) vanishes [see Fig. 11(b)]. Therefore we can define both the interface charge and dipolar densities invariant with respect the position of the integration limits z 1 and z 2 , as long as z 1 lies in R 1 and z 2 lies in R 2 . However, on top of all the drawbacks to that approach presented in Sec. IV, there is an additional pitfall now: an extra dependence on the position of the integration limits z 1 and z 2 appears in the computation of the interface dipole moment from the definition of the reference interface position z int . Indeed, for at least one of the materials the dipole moment density does not vanish within the bulk unit cell, Therefore, the dipole moment density ∆p p associated with ∆ρ p (z) does not equal p p defined in Eq. (59b), not even in the case where z int − z 1 contains an integer number M 1 of unit cells of the left material and z 2 − z int contains an integer number M 2 of unit cells of the right material, Nevertheless, the nanosmoothing procedure developed in Sec. V A remains a useful tool, since the methodology for constructing the nanosmoothed charge density, Eq. In Fig. 11(a) we illustrate the first specific toy model we shall analyze in detail for the polarized case. The parameters of the model are the same as in the caption of Fig. 5 with the exception of the parameter δ (1) i , that now takes a non-zero value δ (1) i = 1.0 where i runs over all the "atoms" of material 1. This toy model simulates a multilayer in which a net polarization has been induced in material 1 by displacing the "electronic clouds" rigidly 1.0 length units towards the right. The neutral but polarized atom of this toy model can be thought of as representing the overall neutral atomic planes of a ferroelectric within which a dipole-moment density is generated by buckling or puckering. In principle, the charge density of material 2 would be modified as a response to the presence of a polarization in material 1. This polarization-induced response, which should be computed self-consistently, is not considered in the present simple toy model. The displacement of the negative charges translates into the asymmetry of the charge density inside material 1 represented in Fig. 11(a). The interface-induced deformation of the charge density ∆ρ p (z), defined as in Eq. (60), is shown in Fig. 11(b) where, as for the unpolarized case, the position of the interface z int has been located half way between the rightmost left atomic layer and the leftmost right atomic layer. The nanosmoothed charge density ρ p (z) is displayed in Fig. 11(c), obtained with the same smoothing function as the the one used in the construction of ρ u (z) in Fig. 6(c). In the last two panels, the two regions inside materials 1 and 2 where ∆ρ p (z) and ρ p (z) vanish are clearly shown.
Some of the features already discussed in Sec. VII for the unpolarized interface still remain valid for the polarized case. In particular: (i) the regions R ′ s within which ρ p vanishes lie within R s within which ∆ρ p vanishes, because smoothing ρ p within R s brings into ρ p values of ρ p for z outside R s ; (ii) the discontinuity of ∆ρ p at the interface and the large fluctuations at the atomic scale in the neighborhood of the interface; (iii) the smooth and continuous profile [note the change of scale in the charge density in Fig. 11(b) and Fig. 11(c)] of ρ p , highlighting its advantages with respect ∆ρ p . However, some other issues are new. Now, as a consequence of the asymmetry of ρ p (z), the charge density at the interface regions, defined as the ranges where ∆ρ p (z) and ρ p (z) differ significantly from zero, are different in the adjacent interfaces contained in our simulated supercell [see the center and the edge of Fig. 11(b) and 11(c)]. They are not simply reflections of one another as in Fig. 6(b) and 6(c).
In Table IV we report the values of the interface charges computed from ρ p (z) (Q p ), ∆ρ p (z) (∆Q p ), and ρ p (z) (Q p ). The integration limits z 1 and z 2 are taken midway between the atoms at the center of each material layer. In contrast to the unpolarized interface, Q p does not vanish. The surplus of negative charge is due to the charge density entering region z 2 − z int from the left from the region between z int − z 1 which is not compensated by the departure of any charge. Since ρ p (z) is not even about z 1 , Q p = Q p . As happened in the unpolarized case, as long as l 1 and l 2  equal an integer number of lattice constants of the bulk unit cell of the material along z, Q p and p p are insensitive to the shape and range of the smoothing function. Only when the range of every filter function entering in Eq. (15) does not equal an integer number of lattice constants of the bulk unit cell, or when the range of the smoothing function L is of the same order as the width of the layer of one of the materials, do Q p and p p differ from the correct value.
Subtracting from ρ p (z) the related ρ u (z), a microscopic charge density constructed by nanosmoothing with filter functions defined with the same parameters with the exception of the displacement of the electronic charge, we obtain the profile of the nanosmoothed polarization-induced charge density, Fig. 12, whose integral between z 1 and z 2 gives σ pol as the accumulation of charge at the interface. Note that since the resulting Q u vanishes, Q p is identical to σ pol . Although, as listed in Table IV, the amount of charge accumulated at the interface is independent of the shape and range of the filter function used (provided that l 1 and l 2 are integer numbers of lattice constants of the bulk unit cell of the material along z), the profile of the polarization-induced charge density is not. The wider the nanosmoothing function, the larger the range of the polarization charge density. Again, this sensitivity of the shape of the smoothed charge density to the nanosmoothing function spoils a direct physical interpretation of its features. The smallest values of l 1 and l 2 , a single lattice constant of each material, yields, as for the unpolarized case, the simplest and physically most relevant profile for ρ p .
Solving the nanosmoothed Poisson equation, Eq. (19), for the nanosmoothed charge density of Fig. 11(a) we obtain the nanosmoothed potential V p (z) displayed in Fig. 13. This electrostatic potential has two distinct features: (i) a jump at each interface due to the interface dipole moment already present in the unpolarized case (see Fig. 9 and Table I); (ii) a field generated by the polarization charge σ pol at the interface. We can isolate this field by subtracting the nanosmoothed potential for the polarized (Fig. 13) and unpolarized (Fig. 9) systems. The result is shown in Fig. 14. The electric field can be computed from the slope of the nanosmoothed potential.
The magnitude of the electric field is related to the magnitude of the polarization charge σ pol . The periodic boundary conditions enforces that E (1) N 1 a (1) bulk + E (2) N 2 a (2) bulk = 0, where N 1 and N 2 are the number of unit cells of materials 1 and 2 stacked to build the superlattice (see Sec. VI). Since each interface carries a charge of magnitude σ pol , the electric field changes its value at the interface as E (2) − E (1) = 4πσ pol , yielding to values for the fields  Potential (2), with the charge density ρ p (z) shown in Fig. 11(a), V p (z), and by solving Eq. (19) with the charge density ρ p (z) shown in Fig. 11(c), V p (z). The nanosmoothing function is defined as in Eq. (15) with l1 = a (1) = 6, and l2 = a (2) = 8, single interplanar distances along z for each material. Atomic units are used.
. Taking σ pol from Table IV, we infer a value of E (1) = 1.196, and E (2) = −0.897, in excellent agreement with the slopes computed in Fig. 14.
From the knowledge of the change of the macroscopic electric field across the interface we can determine the difference in the zero-field polarization of the two materials. Assuming that in our toy model the dielectric constant of the two materials are the same and equal to 1, then − 4π P (2) − P (1) = ǫ E (2) − E (1) .
Since E (2) − E (1) = 4πσ pol , FIG. 14: Difference of the nanosmoothed potential for the polarized (dashed line in Fig. 13) and unpolarized (dashed line in Fig. 9) cases. The magnitude of the electric field E generated by a periodic array of polarization charge distributions σ pol can be computed from the slope of the profile of the potential within material 1 and material 2. Atomic units are used.  Fig. 11(a). The meaning of the different symbols is as in the caption of Table I, with the subscript u replaced by p for the polar case.
But our material 2 is unpolarized, so P (2) = 0, and we arrive to the conclusion that This conclusion can be checked analytically in our particular toy model, where the charge density is the juxtaposition of "atomic-like " charge densities. Integrating the first-moment of the microscopic charge density for the slab of material 1 used to build the superlattice, defined in Eq. 54), we arrive to the conclusion that P (1) = − δ (1) a (1) bulk , whose numerical value in our numerical example equals the polarization charge shown in Table IV. A generalization for the dielectrically mismatched interface can be found in Sec. III. F of Ref. 35.
A second toy model we shall analyze for the polarized case is shown in Fig. 15. It is made of a one dimensional chain in which the unit cell of material 1 has two atoms per unit cell. Inside the unit cell, one of the atoms is charged negatively and the other positively so that the overall charge in the unit cell is neutral. Thus, besides the electronic polarization discussed in the previous case, polarization can be induced by ionic displacements. The parameters of the model are as in Fig. 11, with the exception of A e,1 = 1.5 for the odd "atoms" (net charge, -0.5 per site), A e,1 = 0.5 for the even "atoms" (net charge, +0.5 per site), and δ   1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2   The main conclusions that can be drawn from Fig. 15 are the same as in Fig. 11. As we made for the electronic polarized case, we can infer the value of the bulk zero field polarization from the polarization charge at the interface, that in this case amounts to +0.25.

X. SUMMARY AND DISCUSSION.
The problem of how to extract from the immense detail provided by the first-principles calculations, with resolution at the atomic scale, reliable values of physical quantities of interest which enter into nanoscale electrostatic analysis has been reviewed.
This problem is particularly challenging in the case of interfaces between quite different materials, since all the relevant magnitudes (interface charge and dipole densities, screening lengths, etc.) are overwhelmed by the large and rapid oscillations of the microscopic charge density.
The different procedures to filter out the periodic oscillations of the microscopic quantities, which typically follow the underlying atomic structure, preserving only those features that change in the vicinity of a surface or interface are critically analyzed, and the criteria under which they accurately extract the quantities of interest are discussed.
The approach of defining a reference interface charge density from the bulk unperturbed charge densities of each material is spoilt by the fact that the interface position is undetermined. The profile of the charge density at the interface is discontinuous and displays large fluctuations in the neighborhood of the interface due to the interfaceinduced atomic relaxations. Moreover, both the interface charge and dipole densities (and so the potential drop at the interface) are not unique in this approach, since they depend critically on the position of the interface.
A clearly superior method consists in nanosmoothing the microscopic charge density by taking its convolution with a filter function. This procedure eliminates the contributions coming from the bulk together with the oscillations due to the interface-induced relaxations while producing a continuous charge distribution.
In this work, we prove rigorously that the interface charge and dipolar densities are independent of the nanosmooth-ing function used provided the following conditions are met: (i) it is positive definite, even and normalized within the support region of space within which it is defined [Eqs. (13a)-(13d)]; (ii) it must be sufficiently smooth itself for the Poisson equation to remain invariant after nanosmoothing. By "sufficiently smooth" we understand that the second derivative must exist, and both the function and the first derivative must vanish at the end points. [Eqs. (18a)-(18c)]; (iii) in practice the convolution is done using fast Fourier transforms. Therefore, although the nanosmothing function might violate some of the previous conditions, the family of the finite Fourier series involved in the convolutions is a distribution which converges to the smoothing function in the limit, meeting all the requirements along the way.
(iv) if square-wave filter functions (or convolutions of them) are chosen as filtering functions, the width of each filter function l 1 and l 2 must be integral numbers of the lattice constant of each material in order to have charge and dipole moment densities (and therefore potential shifts across the interface) insensitive to the smoothing function; (v) the smallest acceptable value for l 1 and l 2 should be used in order to avoid the presence of multiple maxima and minima in the smoothed charge density; (vi) the total width of the filter function L should be chosen on the scale of the unit cell length or larger, but significantly smaller than the regions R s of each material where the microscopic charge density is unaffected by the presence of the interface; (vii) after nanosmoothing, the charge density displays two regions R ′ s , within the R s , where the smoothed charge density ρ vanishes. The integration limits z 1 and z 2 used to compute the interface charge and dipole densities must lie in these regions; (viii) once the integration limits are chosen, the charge density at the interface computed from the microscopic and the nanosmoothed charge density are equal if and only if the microscopic charge density is symmetric for |z − z 1,2 | ≤ L. The interface dipole density computed from the microscopic and the nanosmoothed charge density are never equal; (ix) nanosmoothing is only valid for computing charge and dipole moment densities. Nothing can be said about the shape of the nanosmoothed charge density at the interface since it depends critically on the filter function; (x) by using the smallest acceptable values of l 1 and l 2 , generally one unit cell length, reasonable estimates of the density of the charge transferred across the interface and of dipole-layer width can be made.
The nanosmoothing procedure is a powerful technique that allows us to extract the information relevant for computing the change in the average potentials and charge densities from one side of the interface to the other, opening the door to the calculations of band offsets, 7,8,17,18 polarization 29 and dielectric permittivity 30 profiles, effective charges, 36 and force constants 37 in semiconductor-semiconductor interfaces, and depolarizing electric fields and screening lengths in real ferroelectric capacitors. 5,6 The present work sets the basis for a rigorous selection of the minimum size of the supercell required to obtain accurate and smoothing-independent results.
The first two integrals on the right-hand-side can be performed trivially, where boundary x refers to the two intersections of the boundary of the unit cell with the x axis. Since the potential and all its derivatives are periodic in the plane, the previous integral vanishes. The same holds for the second integral in the right-hand-side of Eq. (A3).
Regarding the third integral in Eq. (A3), we can take the second derivative with respect z out of the integral, Gathering together the results of Eqs. (A2) through (A5) we can conclude that that is, the transverse average of the Poisson equation yields a one-dimensional Poisson equation for the transverse average of the potential.

APPENDIX B: NANOSMOOTHING OF THE POISSON EQUATION.
Applying the nanosmoothing procedure to the electrostatic potential V ( r) results in where V (z) is the planar average of the of V ( r) in planes parallel to the interface, and f (z − z ′ ) vanishes for z ′ ≥ z + L and z ′ ≤ z − L. Taking the second derivative of the nanosmoothed potential and requiring that d 2 f (z−z ′ ) dz 2 exists yields Requiring that makes the first term in the right hand side of Eq. (B3) vanish. Integrating the second term in the right hand side of that equation yields so that we arrive at the conclusion that The formal solution of the Poisson equation, Eq. (19), is The electrostatic potential shift across the interface is, therefore, The location of the points at ±∞ are thus far unspecified. They can both be chosen at an image of z 1 or of z 2 within the unit cell without loss of generality since all the relevant computational procedures require periodicity. In that case, all four integrals within the square brackets of Eq. (C2) vanish because of overall electrical neutrality and of the choice of z 1 and z 2 at the microscopic unit cell boundaries, with Eq. (20) following.