15. Macroscopic potential ($\bphi$) and field ($\E$)

https://doi.org/10.48550/arXiv.2403.13981

It is well known that the electric potential is a relative quantity and that, when its value at a point is quoted, this value is always the difference between the potential at the point and a reference potential. In theoretical work the reference point is often taken to be a notional point in the vacuum at infinity, in experimental work it may be a particular electrode, and in engineering it is common to reference potentials to the ‘ground’ or ‘earth’. The difference in meaning between the terms mean inner potential and macroscopic potential is of little relevance to this work. The MIP is the average of $\phi$ over all points in a material; therefore it is a scalar constant. The macroscopic potential is a scalar field, defined at all points in a material, but defined only to a finite precision ${\precphi}$. To within this precision its value at each point in the bulk is equal to the MIP. I mostly refer to $\bphi$ as the macroscopic potential in this section and as the MIP when discussing Bethe’s approximate expression for it in Sec. 16.

It is very important to be able to calculate changes in macroscopic potential. From a conceptual and theoretical viewpoint, our understanding of the relationship between macroscale and microscale electrostatics cannot be considered complete if we do not know how, in principle, to calculate the potential at the macroscale from the charge density or potential at the microscale. From a practical perspective, the MIP, which is believed to be positive, is a key quantity in several areas of experimental and computational science. For example, computer simulations can be used to calculate microscopic charge densities and, up to an unknown constant, microscopic potentials. To improve the designs of devices, such as batteries, fuel cells, chemical sensors, and solar cells, computational scientists simulate their constituent materials independently and, from those simulations, try to calculate the change in the average potential that an electron or ion would experience if it moved from one device component to another [Blumenthal et al., 2017; Hörmann et al., 2019]. The MIP is also used in electron microscopy to analyse and interpret electron diffraction images [Yesibolati et al., 2020]. In both of these contexts, the distinction between the MIP and the time average of the potential felt by the charge carrier is rarely made, but it is the latter that is of interest and the MIP is used as an approximation to it.

In this section I show how to calculate the change ${\Delta\bphi}$ in the macroscale electric potential $\bphi$ between the vacuum above a material’s surface and its bulk from the microscopic charge density $\rho$. The results are easy to generalize to the change in potential across an interface between two materials by treating it as a pair of adjoined surfaces. They also allow the value of $\bphi$ in an isolated material to be calculated relative to a distant point in the vacuum surrounding it.

Some of my results, such as the macroscale potential within a macroscopically-thin film whose surfaces are equally- and oppositely-charged, are well known and serve as a sanity check on my theory and reasoning. However, the main result, which I justified on symmetry grounds in Sec. 6, contradicts most of the literature on this subject over the past century. This result is that ${\bphi=0}$ in the bulk $\bulk$ of an isolated material unless it has charged surfaces or unless its bulk contains charged macroscopic heterogeneities.

It is important to note that ${\Delta\bphi}$ is the difference in potential across a surface at the macroscale. Therefore the derivation that follows is valid in the limit ${a/l\to 0}$, which is the limit in which surfaces become truly planar. Many similar derivations and calculations in existing literature do not assume this limit [Baldereschi et al., 1988; Colombo et al., 1991; Junquera et al., 2007], and this limit is not appropriate for many purposes. For example, it would not be appropriate to assume this limit when calculating variations in the average potential on the nanoscale [Junquera and Ghosez, 2003; Bernardini and Fiorentini, 1998], because it is the limit in which all such variations vanish.

15.1 Change in potential across a surface (${\Delta\bphi}$)

I want to calculate the average potential in a bulk-like region of the material that is much closer to one surface than any other. The plane ${x=x_L}$ is parallel to this surface and in the vacuum just beyond it. Because all other surfaces are further away, it is appropriate to assume that $\abs{x-x_L}$ is much smaller than the material’s lateral dimensions. I will first calculate the microscopic potential due to the charge within a cylindrical region of the material, of radius $R$, whose axis is normal to the surface. I will use the cylindrical coordinates $\rvec=(x,s,\phi)$ or $\rvec=(x,\svec)$, where $\svec\in\realtwo$ is a vector in the plane parallel to the surface, $s=\abs{\svec}$, and $\phi$ is the azimuth.

The potential at ${\rvec_b\equiv(x_b,s_b,\phi_b)=(x_b,\svec_b)}$ due to the charge density within a cylinder bounded by the surfaces ${\abs{\svec-\svec_b}=R}$, ${x=x_L}$, and ${x=x_b}$ is

The average volumetric charge density on a disc of radius $R$, which is parallel to the surface and centered at $(x,\svec)$, is

and similarly for $\Phi_r$. It is easy to see that the surface excess of ${\Delta \rho(\uvec;x,\svec_b,R)}$ vanishes when averaged over the plane parallel to the surface. The planar averages of ${\Phi^{[\Delta\rho]}_L}$ and ${\Phi^{[\Delta\rho]}_r}$ also vanish in the large $R$ limit, i.e.,

15.1.1 Mesoscale average of ${\Phi_L^{[\rhobar]}}$

Integrating over $u$ in the expression for ${\Phi_L^{[\rhobar]}}$, choosing ${R>\depth}$, and using a Taylor expansion gives

I assume that there exists a well-defined macroscopic average of the volumetric charge density on every plane parallel to the surface. By this I mean that, although $\bar{\rho}(x,\svec_b;R)$ may exhibit microscopic fluctuations as $R$ increases, it converges to a well-defined value rather than systematically growing or shrinking. Furthermore, as $\svec_b$ is varied at fixed $R$, ${\bar{\rho}(x,\svec_b;R)}$ fluctuates microscopically about the value to which it converges in the large $R$ limit. If I also assume that the bulk of the material is charge-neutral, only the first term of the series expansion in Eq. (79) can survive the large $R$ limit. In anticipation of this limit, and with the understanding that ‘$=$’ means ‘$\approx$’ until the limit is taken, I write

The mesoscale average over $x_b$ of $\int_{x_L}^{x_b}x\,\bar{\rho}(x\,;R)\dd{x}$ is $X_{\bsigma}^L \,\bsigma_L(R)$, where

15.1.2 Mesoscale average of $\Phi_r$

15.1.3 Test case I: Thin film, ${\abs{x_L-x_R}\ll R}$

Defining ${\bsigma_R(R)}$ as the average areal charge density of the surface at $\mxr$, and ${X_\bsigma^R\in\mxr}$ as the center of the microscale distribution of excess charge at the right-hand surface, allows me to express the mesoscale average of Eq. (82), when ${x_r=x_R}$, as

15.1.4 Test case II: Macroscopic sample, $\abs{x_L-x_R}\not\ll R$

If the bulk of the material is charge neutral, there is no contribution to the potential in the vicinity of $x_b$ from points $x_r>x_b$ sufficiently far from it. However, the potential at $x_b$ will depend on the precise choice of $x_r$ because the integrated charge between $x_b$ and $x_r$ depends sensitively on its value. Therefore, as before, I will take the mesoscale average over $x_r$ of the microscopic potential, $\Phi_r(x_b,x_r;R)$, at $x_b$. Eq. (99) and Eq. (101), which are derived in Sec. 9, can be used in Eq. (82) to write the following expression for the mesoscale average of $\Phi_r(x_b,x_r;R)$ over $x_r$.

Now let us assume that ${\bsigma_L=0}$. Eq. (87) becomes ${\bphi=0}$, which implies that the macroscopic potential, which is the mesoscale average of the microscopic potential, is zero in the bulk of any material that is mesoscopically charge-neutral in the bulk and which does not have charged surfaces or contain any charged macroscale heterogeneities. There are two very important points to note about this result.

The first is that it contradicts the prevailing view that the MIP is finite and positive [Bethe, 1928; Yesibolati et al., 2020]. It also contradicts a view commonly expressed or implied in textbooks on electromagnetism [Jackson, 1998] and solid state physics [Ashcroft and Mermin, 1976; Kittel, 2004], namely, that it is possible for the symmetry of a crystalline microstructure to endow a material with a macroscale electric field.

The second important point is that, because $\Phi[\rho]$ is a linear functional of $\rho$, this result is the only result that could emerge from an internally-consistent theory of structure homogenization. Mathematically, the spatial average ${\Rho}$ of $\rho$ is simply the weighted sum (integral) of the infinite set of charge densities ${\{\rho_u: \rho_u(x+u)\equiv \rho(x), \forall u\in\realone \;\;\text{and}\; \forall x \in \realone\}}$. It follows from linearity that ${\bphi=\expval{\Phi[\rho]}_\prectheo = \Phi[\expval{\rho}_\prectheo] = \Phi[\Rho]}$. Therefore, if ${\Rho=0}$ everywhere, as is the case for an isolated uniform material whose surfaces are uncharged, then ${\bphi=0}$ everywhere.

15.2 Lorentz's fallacy: the macroscopic local field

Nevertheless, most textbooks posit the existence of a non-vanishing $\E$ emanating from the bulks of crystals that lack inversion symmetry at the microscale. Furthermore, it is commonly believed that the net field acting at each point $x$ in a material’s bulk is ${\E + \pp/3\epsilon_0 + \me(x)}$, where $\me$ is a microscopic field and ${\pp/3\epsilon_0}$ is another macroscopic (${\bk=0}$) contribution. The purpose of this section is to critically examine the reasoning used to infer the existence of the contribution ${\pp/3\epsilon_0}$. Almost all derivations of this term are based on a construction and line of reasoning first presented by Hendrik A. Lorentz in a series of lectures given at Columbia University in 1906, which were subsequently published in book form [Lorentz, 1916] (p137). His construction, which is illustrated in Fig. 15, is sometimes known as the Lorentz cavity. This construction has been used in many textbooks [Griffiths, 1999; Jackson, 1998; Born and Huang, 1954; Kirkwood, 1936; van Vleck, 1937; Kirkwood, 1940], but as I will now explain, both Lorentz’s original argument and all of its descendents that I am aware of are fatally flawed.

Lorentz set out to calculate the average force or electric field acting on a microscopic particle $A$ in the bulk of a dielectric in which there exists a macroscopic electric field $\E$. The particle could be a molecule, an atom, or an electron. He expressed the electric field at $A$ as ${\me_A = \EAfar + \EAnear}$, where $\EAfar$ is the field emanating from all material beyond a mesoscopic spherical region of radius $R$ centered at $A$ and $\EAnear$ is the field emanating from all charges within the region, except $A$ itself. As the boldface notation suggests, $\EAfar$ is calculated by treating the material as a continuum; therefore it is a macroscopic quantity. This is reasonable because the length scale on which fluctuations of the microscopic charge density occur is much smaller than the distance ($>R \sim l$) to $A$. On the other hand $\EAnear$ can be expressed as ${\EAnear = \EAnearmacro + \Delta\EAnearmicro}$, where ${\EAnearmacro}$ and ${\Delta\EAnearmicro}$ are macroscopic and microscopic contributions, respectively. It is a microscopic quantity.

In some presentations of this approach, and as illustrated in Fig. 15, these two separate contributions are imagined in different and separated material systems. To calculate $\EAfar$ a continuous material with a cavity in its bulk is imagined, with $A$ at the center of the cavity. To calculate $\EAnear$ a microscopically-varying spherical charge distribution is imagined, with $A$ at its center. This is the charge that was evacuated to form the cavity and it is frozen in the arrangement it had prior to being evacuated.

From here, different authors have derived the term ${\pp/3\epsilon_0}$ in different ways, but it tends to be thought of as arising either from the charge on the cavity’s surfaces or from the dipole moment of the material evacuated from it. It is interesting that the version of the argument that appears in the first edition of Jackson’s book [Jackson, 1962] differs substantially from the one appearing in its 1975 second edition [Jackson, 1975] and that the latter is very similar to the one in Ashcroft and Mermin’s 1976 book [Ashcroft and Mermin, 1976]. However, it is not necessary to go in detail into these differences because we have already introduced the fatal flaw in Lorentz’s reasoning and, to my knowledge, all variants of his derivation suffer from it.

Just as the charges at the surfaces of the materials depicted in Fig. 11 and Fig. 12 depend sensitively, in magnitude and sign, on how the surfaces are terminated, so too does the charge on the surface of the cavity and the dipole moment of the evacuated material. They both depend sensitively and microscopically on the cavity radius $R$ and vanish when averaged over a continuous mesoscopic range of radii. Therefore the true value of the macroscopic field at $A$ is the sum of only two contributions: the applied field and the field from charge at the material’s surfaces.

15.3 LO-TO splitting

However, it is well known that the frequency of a ${\bk\to 0}$ longitudinal optical (LO) phonon is increased by the electric field that is intrinsic to it, and which opposes its motion [Royo and Stengel, 2021; Pick et al., 1970; Cochran and Cowley, 1962; Born and Huang, 1954; Ashcroft and Mermin, 1976; Jones and March, 1973; Coiana et al., 2024; Lyddane et al., 1941]. Were it not for this field, the frequencies of some crystals’ LO and TO phonons would be equal, by symmetry, in the long wavelength limit (${\bk\to 0}$). The breaking of this degeneracy by the LO phonon’s intrinsic field is commonly referred to as LO-TO splitting [Born and Huang, 1954; Ashcroft and Mermin, 1976; Jones and March, 1973].

Therefore I am claiming that, in the ${\bk\to 0}$ limit, an LO phonon of wavevector $\bk$ creates an electric field of wavevector $\bk$, but that a ${\bk=0}$ phonon does not create an electric field in a crystal whose surfaces are earthed.

A ${\bk=0}$ phonon does create a uniform (i.e., ${\bk=0}$) electric field if the surfaces are not earthed, because the polarization current that flows during the rigid relative motion of sublattices changes the areal charge densities on parallel opposing surfaces at equal and opposite rates. If this charge accumulates, a field emanates from it.

15.3.1 Why finite-wavevector LO phonons create electric fields

Suppose that the crystal’s microstructure ($\rho$) is modulated along $\hat{x}$ by an LO phonon, of finite wavelength $\lambda$, propagating along $\hat{x}$. At each instant, this modulation creates regions of excess positive charge and regions of excess negative charge, which alternate along $\hat{x}$ with a wavelength of $\lambda$. This excess charge density wave, in turn, creates an electric field of the same wavelength, $\lambda$, which opposes the LO mode’s motion.

To help understand why excesses of charge are created, it is instructive to consider the perfect crystal, without any perturbation, and to calculate the excess charge, ${\bsigma(\tx)}$, on plane ${\plane(\tx)}$, which is perpendicular to ${\hat{x}}$ at ${\tx\in\bulk}$. To do so, we can treat the plane as a pair of adjoined surfaces, use Eq. (60) to calculate the excess charge on each one, and add them to give ${\bsigma(\tx)}$. Because ${\tx}$ is in the bulk, and the bulk of an unperturbed crystal is uniform, we can choose ${x_L=x_b=\tilde{x}}$. We find that the excess charge at $\tx$ is

Now let us consider what happens when an LO phonon breaks periodicity by modulating the structure along the $x$ axis. When this happens, ${\bsigma_+(\tx)}$ and ${\bsigma_-(\tx)}$ are no longer exactly equal in magnitude, in general, which means that ${\bsigma(\tx)}$ does not vanish. Calculating its value is more complicated than in the periodic case because the crystal is no longer uniform. Therefore it is no longer valid to regard the point $\tx$ as both defining the position of our imaginary surfaces and as points in the ‘bulk’ beneath them. However, as an illustration, let us calculate the net charge in the interval ${\interval(\tx+u,a)}$ averaged over all ${u}$ between ${-a/2}$ and ${a/2}$.

For the purpose of this illustration, let us suppose that each unit cell contains a single anion-cation pair and that the distance between the pair in interval ${[\tx-a,\tx]}$ is smaller by ${\delta x}$ than the distance between the pair in interval ${[\tx,\tx+a]}$, such that the difference between the dipole moments of these unit cells is ${\Delta d = q\delta x}$, where ${q}$ is the cation’s charge. Then, the average over ${u\in(-a/2,a/2)}$ of the net charge in interval ${\interval(\tx+u,a)}$ is ${q\delta x/ a =\Delta d/a}$.

Now suppose that we have calculated the same quantity for every pair of adjacent unit cells in an interval ${\interval(\tx,\ell)}$, where ${a\ll \ell\ll\lambda}$, and then repeated this calculation for a continuous range of values of $\tx$. Let us denote ${1/\ell}$ times the sum of all net charges in interval ${\interval(x,\ell)}$ by ${\expval{\rho}^*(x)}$ and ${1/\ell}$ times the sum of the cells’ dipole moments by ${\expval{\mpp}^*(x)}$. Then, it can be shown that ${\expval{\rho}^*(x)=-\div \expval{\mpp}^*(x)}$. The similarity of this expression to the relation ${\rho=-\div\pp}$ is not coincidental: Maxwell used a similar line of reasoning to deduce it, albeit with displacements of charges replaced by displacements of the ether.

This example illustrates that ${\bsigma(\tx)}$ does not vanish in the presence of an LO perturbation because the symmetry reason for it vanishing no longer exists. Furthermore, because a more realistic charge density $\rho$ would be a smooth function of position, ${\bsigma(\tx)}$ would be a smooth function of $\tx$, with the same periodicity $\lambda$ as the LO perturbation that created it. Therefore, there would be an excess charge density wave of periodicity $\lambda$, from which would emanate an electric field of periodicity $\lambda$.

15.3.2 Zero-wavevector LO phonons

The fact that ${\melo(\bk,u)}$ is discontinuous at ${\bk=0}$ is well known when expressed in a different way: Squared phonon frequencies are eigenvalues of a crystal’s dynamical matrix. Therefore, if ${\melo(\bk,u)}$ vanishes suddenly at ${\bk=0}$, causing LO phonon frequencies at ${\bk=0}$ to be smaller than their values in the ${\bk\to 0}$ limit, the dynamical matrix must be discontinuous at ${\bk=0}$. It is very well known that it is discontinuous, and non-analytic corrections are commonly applied to the ${\bk=0}$ dynamical matrix to calculate the ${\bk\to 0}$ dynamical matrix [Pick et al., 1970; Cochran and Cowley, 1962; Cochran, 1960; Giannozzi et al., 1991; Gonze and Lee, 1997; Baroni et al., 2001; Born and Huang, 1954; Jones and March, 1973]. The ${\bk\to 0}$ dynamical matrix is not one matrix, in general, because both it, and the LO phonon frequency, depend on the direction in reciprocal space from which the point ${\bk=0}$ is approached.

To understand why the ${\bk\to 0}$ limit is singular, it is easier to think about the LO phonon’s wavelength in the ${\lambda\to\infty}$ limit than its wavevector in the ${\bk\to 0}$ limit: Imagine a microscopic or mesoscopic neighbourhood of a point in the bulk of an arbitrarily large perfect crystal, and then imagine that the crystal is perturbed by displacing it from equilibrium along the eigenvector of an LO phonon of wavelength $\lambda$. Now imagine increasing $\lambda$.

As ${\lambda}$ becomes much larger than the size of the neighbourhood, and continues to increase, the microstructure within the neighbourhood looks more and more like it would look if the crystal’s sublattices had been displaced rigidly relative to one another. Therefore it looks more and more like a crystal that has been perturbed by displacing it along the eigenvector of a ${\bk=0}$ phonon. Nevertheless, no matter how large $\lambda$ becomes, if one moves a distance ${\lambda/2}$ in the direction of ${\bk}$, the relative displacements of the atoms in the direction of $\bk$ are reversed.

In other words, microscopically, increasing $\lambda$ brings the structure closer to the ${\bk=0}$ structure, but macroscopically it does not; and it is the macroscopic structure that determines whether or not there is a macroscopic $\E$ field.