There’s an identity in electromagnetism which has been bugging me since college.

Gauss’s law says that the divergence of the electric field is equivalent to the charge distribution: $\nabla \cdot \mathbf{E}=\rho$. But in order to use this for a point charge—which is the most basic example in the subject!—we already don’t have the mathematical objects we need to calculate the divergence on the left or to represent the charge distribution on the right.

After all, the field of a point charge has to be $\mathbf{E}=q\hat{\mathbf{r}}/4\pi r^2$, and since its charge should be concentrated at a point it has to be a delta function: $\nabla \cdot (q\hat{\mathbf{r}}/4\pi r^2)=q\delta (\mathbf{x})$. In your multivariable-calculus-based E&M class you might mention this briefly, at best. Yet it is… kinda weird? And important? It feels like it should make a basic fact that lives inside of a larger intuitive framework of divergences and delta functions and everything else.

Why, in the first place, are we using this divergence operator that we didn’t know how to actually calculate—are we missing something? Are there other divergences that we don’t know how to calculate? Does it work the same way in other dimensions? What about other powers of $\frac{1}{r}$? Are there other derivative operators we don’t know about that do similar tricks? Is there an equivalent version for the curl and by extension the magnetic field? Is there an equivalent version for dipoles, or multipoles? Etc. (The answer to all of these questions is ‘yes’, by the way.)

Not only is it unsatisfying, it’s also hard to learn about. For years I’ve been referring back to this one rather confusing physicsforum.com post, and I’m pretty tired of reading that. It’s not even good! Griffith’s and other E&M textbooks also mention it but they’re obscured by pedagogy and most of the interesting parts are left as exercises… and even then they don’t have much to say. Meanwhile venerable Wikipedia’s treatment is very slim and spread out over many hard-to-navigate articles; the best one is probably here but it’s still not great.

So today’s the day: I’m going to figure this out in all the generalization I want and write myself the reference I have wanted so I never have to visit that forum post, or that one page of Griffiths, ever again.

---

1. The Basic Argument

The first thing we learn in electrostatics is that the electric field of a point particle is

$$\mathbf{E}=\frac{q\hat{\mathbf{r}}}{4\pi r^2}$$

That is, the field points radially out in every direction from the ‘infinitely concentrated’ point charge, and the magnitude falls off proportional to $4\pi r^2$. Non-coincidentally, $4\pi r^2$ is the formula for the surface area of a sphere of radius $r$. Evidently electric flux lines get weaker exactly in proportion to how much they “spread out”. It is as though you had a pipe whose input has to be equal to its output, except the input is at the origin and the output is “every direction at once”. Put differently, an electric charge is the source of a flux and then that flux fluxes around in exactly the way a flux has to flux around, which is: conservatively. A source of nonzero electric flux is what a charge is.

Which means that you can detect the presence of charges by measuring the flux around a volume. This is Gauss’s Law: that summing the electric flux through any closed surface measures the total charge contained within it.

$$\iint \mathbf{E}\cdot d\mathbf{A}=q_{\text{enclosed}}$$

The divergence theorem turns Gauss’s Law into

$${\iiint }_V\nabla \cdot \mathbf{E}dV=q_{\text{enclosed}}$$

We also learn the differential form of Gauss’s Law, which says that the divergence $\nabla \cdot \mathbf{E}$ equals the charge distribution $\rho (\mathbf{x})$. For a point particle the integral’s value is entirely concentrated at the origin, so $\rho (\mathbf{x})$ has to be a delta function:

$$\rho (\mathbf{x})=q\delta (\mathbf{x})$$

But we also know the functional form of $\mathbf{E}$ for a point charge: it’s $q\hat{\mathbf{r}}/4\pi r^2$. Hence at least in ${\mathbb{R}}^3$ it must be true that:

$$\nabla \cdot \frac{\hat{\mathbf{r}}}{r^2}=4\pi {\delta }^3(\mathbf{x})$$

Equivalently:¹

¹Sorry but I am stubbornly opposed to the Laplacian symbol $\Delta ={\nabla }^2$ ↩

$$-{\nabla }^2\frac{1}{r}=4\pi {\delta }^3(\mathbf{x})$$

We can also write this delta function in terms of $r$:²

²By the way (because I definitely didn’t know this off the top of my head) you can’t just replace $\delta (\mathbf{x})$ with $\delta (r)$. Translating $\delta (\mathbf{x})$ to spherical coordinates requires some extra care because it has to be true that ${\int }_V\delta (x,y,z)d^3\mathbf{x}={\int }_V\delta (r)(r^2\sin \theta )drd\theta d\varphi$. The two angular integrals integrate to $4\pi$, so $\delta (r)$ has to have a $\frac{1}{4\pi r^2}$ to cancel everything out. ↩

$$4\pi {\delta }^3(\mathbf{x})=-{\nabla }^2\frac{1}{r}=\nabla \frac{\hat{\mathbf{r}}}{r^2}=\frac{\delta (r)}{r^2}$$

Which is neat, and also rather suspicious-looking. Seems like the more interesting identity here is that ${\delta }^3(\mathbf{x})=\delta (r)/4\pi r^2$ where the numerator is the surface area of a $2$-sphere.

It’s pleasing (since it’s pleasing when any integral is easy) that you can simply plug that into the equation for the electric field of an arbitrary charge distribution and recover Gauss’s law:

$$\begin{array}{rl}\mathbf{E}(\mathbf{x})& =\frac{1}{4\pi }\int \frac{\rho (x^{\prime })}{|\mathbf{x}-{\mathbf{x}}^{\prime }{|}^2}d{x}^{\prime }\\ \nabla \cdot \mathbf{E}(\mathbf{x})& =\frac{1}{4\pi }\int [\nabla \cdot \frac{1}{|\mathbf{x}-{\mathbf{x}}^{\prime }{|}^2}]\rho (x^{\prime })d{\mathbf{x}}^{\prime }\\ & =\frac{1}{4\pi }\int [4\pi \delta (\mathbf{x}-{\mathbf{x}}^{\prime })]\rho (x^{\prime })d{\mathbf{x}}^{\prime }\\ & =\int \delta (\mathbf{x}-{\mathbf{x}}^{\prime })\rho (x^{\prime })d{\mathbf{x}}^{\prime }\\ & =\rho (\mathbf{x})\end{array}$$

---

2. The Other Definition of Divergence

Producing this result by working backwards from physics is good enough for most purposes, but it’s a bit perplexing. Maybe there’s a cleaner derivation?

I’ve looked around and there are some formal-ish ways to do it, by a procedure they call ‘regularizing’ $\frac{\mathbf{r}}{r^2}$ as a limit of a more complicated function like $\mathbf{r}/(r^2+a^2{)}^{\frac{3}{2}}$, which is a way of producing distributions as a limit of non-distributions. I guess it’s rigorous, but I don’t want to do it. It doesn’t teach me anything new about divergences or delta functions at all. Plus it just feels unnecessary.³

³In general there are two ways of defining delta functions: either you make them out of a limit of functions you know how to do analysis on and prove that the limit is well-defined, which is this regularization procedure… or you define them to have certain properties by fiat, and then show that they exist. The latter, IMO, is the “right” way. I think the approximations are only to satisfy people who are unnecessarily fixated on classical functions that have definite values at points. (There’s a rather nice book called “Theory of Distributions: A Non-Technical Introduction” by Richards & Youn which I like because of how much it commits to the better approach.) ↩

Others claim that the divergence of $\mathbf{r}/r^2$ “is” undefined according to the usual definition, and that we’re just assigning a value to make the divergence theorem work. They’re obviously wrong: we’re not inventing a value; we’re discovering the actual value and it just requires delta functions to express. For the purposes of physics we don’t care at all about confining the space of objects we consider to just the standard-issue smooth functions. Evidently multivariable calculus wants distributions to get involved; we may as well let it happen.

The most satisfying explanation, in my opinion, is based on a different definition of divergence which isn’t used as much:

Recall that in multivariable calculus class we initially define divergence as a sum of partial derivatives ${\partial }_x\hat{\mathbf{x}}+{\partial }_y\hat{\mathbf{y}}+{\partial }_z\hat{\mathbf{z}}$ (or whatever it becomes in other coordinate systems). But there’s another definition which is really a more direct extension of the one-variable derivative⁴

⁴The sense in which this is a 1d derivative: $1/V$ can be written as $\int dV$, with the integral is over a ball of radius $ϵ$. In the numerator it’s over the boundary of that ball. So divergence is $\underset{e\to 0}{\lim }[{\int }_{\partial B_{ϵ}}fdA]/({\int }_{B_{ϵ}}dV).$ In one dimension a ball of radius $ϵ$ is just a line segment, so this is literally the same as the 1d derivative ${\lim }_{ϵ\to 0}\frac{f(x+ϵ)-f(x)}{ϵ}$. ↩

. It looks like this:

$$\nabla \cdot \mathbf{F}=\underset{V\to 0}{\lim }\frac{1}{|V|}{\oint }_{\partial V}\mathbf{F}\cdot d\mathbf{n}$$

That is, it’s a ratio of the flux through a volume surrounding the point divided by the volume itself, as the volume goes to zero. It’s actually a standard definition and is at the top of the Wikipedia page on divergence, but for whatever reason it doesn’t come up as often. To use it, you compute the volume in the denominator as a sphere or a cube or whatever you want. For instance if $\mathbf{F}=x\hat{\mathbf{x}}+y\hat{\mathbf{y}}+z\hat{\mathbf{z}}=r\hat{\mathbf{r}}$ and we integrate over a sphere, then

$$\begin{array}{rl}\nabla \cdot \mathbf{F}& =\frac{1}{4/3\pi r^3}\oint (r\hat{\mathbf{r}})\cdot \hat{\mathbf{r}}(r^{2}d\theta d\varphi )\\ & =\frac{4\pi r^3}{4/3\pi r^3}\\ & =3\\ & =({\partial }_x\hat{\mathbf{x}}+{\partial }_y\hat{\mathbf{y}}+{\partial }_z\hat{\mathbf{z}})\cdot (x\hat{\mathbf{x}}+y\hat{\mathbf{y}}+z\hat{\mathbf{z}})\end{array}$$

In many ways this is more intuitive! On the other hand I have no idea how to prove that it’s equivalent to ${\partial }_x\hat{\mathbf{x}}+{\partial }_y\hat{\mathbf{y}}+{\partial }_z\hat{\mathbf{z}}$ in general, and it’s hard to google for because you just get results about proving the divergence theorem. Sigh. But it makes some sense. $\nabla \cdot \mathbf{F}=({\partial }_x,{\partial }_y,{\partial }_z)\cdot \mathbf{F}$ acts like the same formula but implemented on a cube instead of a sphere.

Using this definition we can derive the weird equation from E&M as follows. The flux of $\hat{\mathbf{r}}/r^2$ through a sphere of radius $ϵ$ is $4\pi {ϵ}^2/{ϵ}^2=4\pi$ as long as the volume contains the origin and therefore the limit is $\nabla \cdot \hat{\mathbf{r}}/r^2={\lim }_{V\to 0}4\pi 1_{\mathcal{O}\in V}/V$ which, if you integrate that against test functions, acts like $4\pi \delta (\mathbf{x})$. Something like that.

---

3. Dealing with $\delta (r)$ in other dimensions

One nice thing about the integral definition is that it makes generalizations of the delta function divergence to other dimensions very natural: just integrate over different types of objects. In each case the coefficient is given by the surface area of an $(n-1)$-sphere of radius $R=1$, which you can look up.

$$\begin{array}{rl}\nabla \cdot \frac{\hat{\mathbf{r}}}{r^{n-1}}& =S_{n-1}\delta (\mathbf{x})\\ & =\frac{\delta (r)}{r^{n-1}}\end{array}$$

e.g. in ${\mathbb{R}}^2$ with polar coordinates (so ${\mathbf{r}}_{xy}=x\hat{\mathbf{x}}+y\hat{\mathbf{y}}$):

$$\nabla \cdot \frac{{\hat{\mathbf{r}}}_{xy}}{r_{xy}}=2\pi \delta (x,y)=\frac{\delta (r_{xy})}{r_{xy}}$$

Note that you can totally compute a 2-divergence in a plane in ${\mathbb{R}}^3$, or a 3-divergence in ${\mathbb{R}}^4$, etc. I guess we could write it as ${\nabla }_{xy}\cdot \mathbf{F}$. I do vaguely recall seeing that object in formulas in my life but can’t remember where.

In fact this construction works in ${\mathbb{R}}^1$ also, but it’s kinda weird: the 1d version of $\hat{\mathbf{r}}/r^2$ in ${\mathbb{R}}^3$ and ${\hat{\mathbf{r}}}_{xy}/r_{xy}$ in ${\mathbb{R}}^2$ is ${\hat{\mathbf{r}}}_x$, the “one dimensional radius function”, also written less strangely as $\text{sgn}(x)\hat{\mathbf{x}}$. That is, it’s a unit vector pointing in the $+\mathbf{x}$ direction in the positive numbers and the $-\mathbf{x}$ direction in the negative numbers. Then:

$$\begin{array}{rl}\nabla \cdot ({\hat{\mathbf{r}}}_x)& =({\partial }_x\hat{\mathbf{x}})\cdot (\text{sgn}(x)\hat{\mathbf{x}})\\ & ={\partial }_x\text{sgn}(x)\\ & =2\delta (x)\\ & =\delta (r_x)\end{array}$$

The factor of $2$ can be regarded as the “surface area” of a $0$-sphere, that is, of a line segment. Admittedly it’s kinda weird to write $2\delta (x)=\delta (r_x)$. One way of thinking about it is that ${\hat{\mathbf{r}}}_x$ “acts like” $\delta (x)\theta (x)$, covering only half the displacement at the origin with a step function, whereas $\hat{\mathbf{x}}(\text{sgn}(x))$ acts like $\delta (x)\text{sgn}(x)$ and covers the full displacement. Hence the factor of $2$.

Yes, that sounds weird and made up. I’m happy with it mostly I realized that it gives a satisfying result in in ${\mathbb{R}}^3$ as well: recall that in spherical coordinates the radial term of the divergence looks like $\nabla \cdot \mathbf{f}=\frac{1}{r^2}{\partial }_r(r^2{f}_r)$. Well, suppose $f=\hat{\mathbf{r}}\theta (r)/r^2$ where once again we imagine that we need the $\theta (r)$ in there to deal with how $r$ switches signs at the origin. Then $\nabla \cdot \mathbf{f}=\frac{1}{r^2}{\partial }_r[r^2\frac{\theta (r)}{r^2}]=\frac{1}{r^2}{\partial }_r\theta (r)=\delta (r)/r^2$ is the right value. Not bad, eh?

By the way, there is some information about all of this on the Wikipedia article for Newtonian potential. They call the function which is the fundamental solution to ${\nabla }^{2}f=\delta$ in ${\mathbb{R}}^d$ the “Newtonian Kernel” $\Gamma$, and write

$$\Gamma (x)=\{\begin{array}{ll}2\pi \log r& d=2\\ \frac{1}{d(2-d)V_d}{r}^{2-d}& d\ne 2\end{array}$$

Where $V_d$ is the volume of the $d$-sphere. That’s a bit confusing. It’s easier to follow with the identity $V_d=\frac{S_d}{d}$ where $S_d$ is the surface area of the $d$-sphere. Then this is really

$$\Gamma (x)=\{\begin{array}{ll}\frac{1}{2\pi }\log r& d=2\\ \frac{1}{(2-d)S_d}{r}^{2-d}& d\ne 2\end{array}$$

And its gradient is given by the same formula in all dimensions:

$$\nabla \Gamma (x)=\frac{1}{S_d}\frac{\hat{\mathbf{r}}}{r^{d-1}}$$

This agrees with what we wrote above, and even works in $d=1$ if you consider the “surface area of the 1-sphere” to be $S_1=2$.

---

4. Other Shapes and Multipoles

In fact every example charge distribution from elementary E&M has an expression as a delta function. It’s just that we’re not very… good… at using delta functions so we don’t normally write them that way.

Here are two classic examples from intro E&M, plus the electric fields that you get from applying Gauss’s law to their symmetries:

• an infinite line of charge in the $z$-direction with linear charge density $\mu$ has electric field $\mathbf{E}(\mathbf{x})=\mu {\hat{\mathbf{r}}}_{xy}/(2\pi r_{xy})$.
• an infinite plane of charge in the $xy$ plane with area charge density $\sigma$ has constant electric field $\mathbf{E}(\mathbf{x})=\sigma \text{sgn}(z)\hat{\mathbf{z}}/2$. The $\text{sgn}(z)$ makes this valid on both sides of the plane.

In each case it should be that $\nabla \cdot \mathbf{E}=\rho (\mathbf{x})$. Evidently:

$$\begin{array}{rl}{\rho }_{\text{line}}(\mathbf{x})& =\mu \delta (y,z)\\ & =\mu \frac{\delta (r_{xy})}{2\pi r_{xy}}\\ {\rho }_{\text{plane}}(\mathbf{x})& =\sigma \delta (z)\\ & =\sigma \frac{\delta (r_z)}{2}\end{array}$$

Here are the forms of $V$, $\mathbf{E}$, and $\rho$ side-by-side:

$$\begin{array}{rl}{\rho }_{\text{line}}(\mathbf{x})& =\mu \delta (y,z)\\ {\mathbf{E}}_{\text{line}}(\mathbf{x})& =\frac{\mu }{2\pi }\frac{{\hat{\mathbf{r}}}_{xy}}{r_{xy}}\\ V_{\text{line}}(\mathbf{x})& =\frac{\mu }{2\pi }\ln r_{xy}\\ & \\ {\rho }_{\text{plane}}(\mathbf{x})& =\sigma \delta (z)\\ {\mathbf{E}}_{\text{plane}}(\mathbf{x})& =\frac{\sigma }{2}\text{sgn}(z)\hat{\mathbf{z}}\\ V_{\text{plane}}(\mathbf{x})& =\frac{\sigma }{2}|z|\hat{\mathbf{z}}\end{array}$$

How about some other interesting charge distributions?

A perfect electric dipole is the limiting case of a positive and negative charge next to each other, so that their net charge is zero but there is a nonzero dipole moment $\mathbf{p}$ along a certain axis axis. The potential, electric field, and charge distributions of dipoles are given by the limit as we press two point charges together while keeping the product $qd$ fixed. But in fact this limit is just a directional derivative:

$${\rho }_{\text{dipole}}(\mathbf{x})=\underset{d\to 0}{\lim }[(+q)\delta (\mathbf{x}-\mathbf{d}/2)+(-q)\delta (\mathbf{x}+\mathbf{d}/2)]=-\mathbf{p}\cdot \partial [\delta (\mathbf{x})]=-{\partial }_{\mathbf{p}}\delta (\mathbf{x})$$

So the charge distribution of a dipole is the gradient of a delta-function. That makes sense: the net charge is zero, but there’s two infinite spikes at the origin infinitesimally close to each other, which is what ${\delta }^{\prime }$ looks like.

We can immediately write down the electric field and potential also:

$$\begin{array}{rl}{V}_{\text{dipole}}(\mathbf{x})& =-{\partial }_{\mathbf{p}}[\frac{1}{4\pi r}]=\mathbf{p}\cdot [\frac{\hat{\mathbf{r}}}{4\pi r^2}]\\ {\mathbf{E}}_{\text{dipole}}(\mathbf{x})=-\nabla V_{\text{dipole}}(\mathbf{x})& =-{\partial }_{\mathbf{p}}[\frac{\hat{\mathbf{r}}}{4\pi r^2}]=\frac{3(\mathbf{p}\cdot \hat{\mathbf{r}})\hat{\mathbf{r}}-\mathbf{p}}{4\pi r^3}\\ {\rho }_{\text{dipole}}(\mathbf{x})=\nabla {\mathbf{E}}_{\text{dipole}}(\mathbf{x})& =-{\partial }_{\mathbf{p}}[\delta (\mathbf{x})]\end{array}$$

(Although see the next aside for a correction to $\mathbf{E}$: there’s apparently supposed to be a delta function term there also.)

Another way of looking at dipoles is to consider manually placing a bunch of charges at a distance $h$ apart and then taking $h\to 0$. Write ${\Delta }_{\mathbf{v}}$ for a finite difference at a distance $h$ along the $\mathbf{v}$. For instance ${\Delta }_{\mathbf{v}}=f(\mathbf{x}+\mathbf{v}h)-f(\mathbf{x})$. Note that ${\partial }_{\mathbf{v}}f(\mathbf{x})={\lim }_{h\to 0}{\Delta }_{\mathbf{v}}f(\mathbf{x})$. Also, we can write $T_{\mathbf{z}}f\equiv f(\mathbf{x}+\mathbf{z}h)$, such that $D_{\mathbf{z}}=T_{\mathbf{z}}-1$ and $D_{\mathbf{z}}f=(T_{\mathbf{z}}-1)f=T_{\mathbf{z}}f-f$.

Then a “physical” dipole (where the charges are a small but finite distance apart) is proportional to

$${\Delta }_{\mathbf{z}}\delta (\mathbf{x})=\delta (\mathbf{x}+\mathbf{z}h)-\delta (\mathbf{x})$$

Then the infinitesimal dipole charge distribution is given by $\rho (\mathbf{x})=-\frac{q}{h}{\Delta }_{\mathbf{z}}\delta (\mathbf{x})$, which in the limit where $h\to 0$ with $hq=p$ fixed gives

$${\rho }_{\text{dipole}}(\mathbf{x})=(-p){\partial }_{\mathbf{z}}\delta (\mathbf{x})=(-p\mathbf{z})\cdot \partial \delta (\mathbf{x})$$

A physical quadrupole is given by the “second finite difference” (so, second derivative). We can consider the case along a single axis:

$$\begin{array}{rl}{\Delta }_{\mathbf{z}}^2\delta & =(T_{\mathbf{z}}-1{)}^2\delta \\ & =T_{\mathbf{z}}^2\delta -2T_{\mathbf{z}}\delta +\delta \\ & \equiv \delta (\mathbf{x}+2h\mathbf{z})-2\delta (\mathbf{x}+h\mathbf{z})+\delta (\mathbf{x})\end{array}$$

In the limit we take $h\to 0$ while holding $q^{2}h=Q$ and get ${\rho }_{\mathbf{zz}\text{-quadrupole}}=Q{\partial }_{\mathbf{z}}^2\delta =\hat{Q}\cdot {\partial }^2\delta$ (where $\hat{Q}$ is a quadrupole tensor which only has an $\mathbf{zz}$ component). Or, we can do an $\mathbf{y}$-$\mathbf{z}$ quadrupole:

$$\begin{array}{rl}{\Delta }_{\mathbf{y}}{\Delta }_{\mathbf{z}}\delta (\mathbf{x})& =(T_{\mathbf{y}}-1)(T_{\mathbf{z}}-1)\delta \\ & =T_{\mathbf{y}}{T}_{\mathbf{z}}\delta -T_{\mathbf{y}}\delta -T_{\mathbf{z}}\delta +\delta \\ & \equiv \delta (\mathbf{x}+h(\mathbf{z}+\mathbf{y}))-\delta (\mathbf{x}+h\mathbf{xz})-\delta (\mathbf{x}+h\mathbf{z})+\delta (\mathbf{x})\end{array}$$

The limit with $q^{2}h=Q$ is ${\rho }_{yz\text{-quadrupole}}=Q{\partial }_y{\partial }_z\delta (\mathbf{x})=\hat{Q}\cdot {\partial }^2\delta$ (where $\hat{Q}$ is now a quadrupole tensor which only has an $yz$ component).

This construction is nicely easy to generalize, for instance to any charge distribution that’s a mix of points and multipoles at any separation from each other.

We can also make lines and planes and other shapes out of multipoles. For instance a “line of dipoles” looks like a positively-charged line infinitesimally close to a negatively-charged line. The result is just that we take an additional $-\partial$ of every term, which is equivalent to forcing the two charged surfaces to be $d\to 0$ apart with opposite signs while holding the ratio $qd$ constant. For instance a line of charge on the $z$ axis had charge distribution ${\rho }_{\text{line}}=\mu \delta (x,y)$. A dipole along the $x$-axis made out of two of thsee has charge distribution

$${\rho }_{\text{line of dipoles}}=\underset{h\to 0,h\mu =p}{\lim }[{\rho }_{\text{line}}(\mathbf{x}+h\hat{x})-{\rho }_{\text{line}}(\mathbf{x})]=p{\partial }_x\delta (x,y)$$

More generally, the multipole distribution of a potential gives an Taylor expansion of $V$ away from the charges, in terms of the increasingly higher-order moments $(q,\mathbf{p},\hat{Q},\dots )$ of the underlying charge distribtion.⁵

⁵$\hat{Q}$ here is the rank-2 quadrupole tensor. Equations using it and higher-order multipoles are best unpacked in index notation: ${\partial }_{\hat{Q}}^2\delta (\mathbf{x})=Q^{ij}{\partial }_i{\partial }_j\delta (\mathbf{x})$. By the way, I haven’t learned a ton about $\hat{Q}$, and I’m a bit confused about when it ought to have a factor of $1/2$ or not. It might be a convention. Definitely double-check before using this. ↩

$$\begin{array}{rl}V(\mathbf{x}{)}_{\text{multipole}}& =\frac{1}{4\pi }[\frac{q}{r}+\mathbf{p}\cdot \frac{\hat{\mathbf{r}}}{r^2}+\frac{1}{2}\hat{Q}\cdot \frac{{\hat{\mathbf{r}}}^{\otimes 2}}{r^3}+\dots ]\\ & =\frac{1}{4\pi r}[q+\mathbf{p}\cdot (-\partial )+\frac{1}{2}\hat{Q}\cdot (-\partial {)}^2+\dots ]\end{array}$$

Each of these terms is going to have a $-{\nabla }^2$ that gives a delta-function derivative of some order.

$$\begin{array}{rl}-\nabla V=\mathbf{E}(\mathbf{x}{)}_{\text{multipole}}& =[q+\mathbf{p}\cdot (-\partial )+\frac{1}{2}\hat{Q}\cdot (-\partial {)}^2+\dots ]\frac{\hat{\mathbf{r}}}{4\pi r^2}\\ -{\nabla }^{2}V=\rho (\mathbf{x}{)}_{\text{multipole}}& =[q+\mathbf{p}\cdot (-\partial )+\frac{1}{2}\hat{Q}\cdot (-\partial {)}^2+\dots ]{\delta }^3(\mathbf{x})\end{array}$$

---

5. Other powers of $r$

The multipole examples imply that in general, there are lots of objects in ${\mathbb{R}}^3$ that have delta function divergences and it’s not just $\hat{\mathbf{r}}/r^2$, but the results are going to involve derivatives of delta functions instead… which are even harder to detect with the usual implementations of divergence.

For instance we can can compute $\nabla \cdot \hat{\mathbf{r}}/r^3$ in two ways. Everywhere except the origin, we can use $\nabla \cdot f(r)=\frac{1}{r^2}{\partial }_r(r^2{f}_r)$ to get

$$\nabla \cdot \frac{\hat{\mathbf{r}}}{r^3}=-\frac{1}{r^4}$$

And around the origin we use the integral definition of divergence:

$$\begin{array}{rl}\nabla \cdot \frac{\hat{\mathbf{r}}}{r^3}& =\underset{R\to 0}{\lim }{\oint }_R\frac{1}{r}\frac{d\Omega }{r^2}\\ & =\frac{4\pi }{r}{\delta }^3(\mathbf{x})\\ & =\frac{\delta (r)}{r^3}\end{array}$$

It seems like we should be able to make the delta function into a derivative, similar to what showed up in the multipole distribution. But it’s a little weird. Normally we can replace $\delta /x^n$ with $\frac{(-1{)}^n}{n!}{\delta }^{(x)}$. But it seems like the identity is probably a little bit different in radial coordinates, since after all we expect this to be true:

$$\frac{4\pi }{r}{\delta }^3(\mathbf{x})=\frac{1}{r}\frac{\delta (r)}{r^2}=\frac{-{\delta }^{\prime }(r)}{r^2}$$

That is, $\frac{1}{r^3}$ should give a first derivative, not a third derivative:

$$\frac{4\pi }{r}{\delta }^3(\mathbf{x})=\frac{1}{r}\frac{\delta (r)}{r^2}\stackrel{!}{\ne }\frac{-{\delta }^{(3)}(r)}{3!}$$

The problem, I presume, is basically that ${\partial }_r\delta (r)$ is a weird object, because in an integral against a test function $⟨-{\partial }_r\delta (r),f⟩$, the normal integration-by-parts that lets us move the derivative across doesn’t work: $⟨-{\partial }_r\delta (r),f⟩\ne ⟨\delta ,{\partial }_{r}f⟩$: since the radial integral has bounds $(0,r)$, we can’t ignore the boundary. This integration by parts is what justifies $\delta (x)/x=-{\partial }_x\delta (x)$ normally, since $⟨\delta (x)/x,f⟩=⟨\delta ,f/x⟩=-f^{\prime }(0)$ (in a principal-value sense?). Therefore it is probably best to leave $\delta (r)/r^3$ as-is instead of trying to turn it into a radial derivative.

Nevertheless I’m pretty sure there are ways to do it, but it’s a lot more than I want to figure out right now. Roughly speaking, though, we can expect that a term like $\delta (r)/r^k$ is going to turn into a delta-function derivative that is comparable to $\frac{1}{r^2}\frac{(-1{)}^{k-2}}{(k-2)!}{\delta }^{(k-2)}$, that is, it will act like a delta-derivative but it acts like a factor of $\frac{1}{r^{(k-2)}}$ instead. But I hope to figure out the actual details in a future article.

---

6. Curl and Magnetic Fields

One last question. How does this work for magnetism and curl?

The equivalent Maxwell equation is Ampère’s law, which establishes that the curl of the magnetic field is proportional to the current density (in units with ${\mu }_0=1$):

$$\nabla \times \mathbf{B}=\mathbf{J}$$

The integral form is:

$${\oint }_{\partial A}\mathbf{B}\cdot d\ell ={\iint }_A\mathbf{J}\cdot d\mathbf{A}$$

Like divergence, there’s an integral form for the curl, which is basically the same idea except that it is computed each plane instead of in the overall volume. The coordinates on the plane with normal $\mathbf{u}$ are given by:

$$(\nabla \times \mathbf{F})\cdot \hat{\mathbf{u}}=\underset{A\to 0}{\lim }\frac{1}{|A|}{\oint }_{\partial A}\mathbf{F}\cdot d\ell$$

We can use this to justify delta-function formulas for currents on various surfaces, although I’m going to skip most of the steps. The equivalent identity is going to be for a current which is entirely concentrated in a line, which we’ll assume is at the origin of the $(x,y)$ plane and directed up the $z$ axis.

$$\mathbf{J}=j\hat{\mathbf{z}}\delta (x,y)=j\hat{\mathbf{z}}\frac{\delta ({\mathbf{r}}_{xy})}{r_{xy}}$$

(That’s in $(r_{xy},\theta ,z)$ cylindrical coordinates; the $\frac{1}{r_{xy}}$ factor is the same as the one for the 2d divergence up above.)

Of course the magnetic field due to a infinitely thin wire is a basic textbook example, so we know immediately what function has this as its curl:

$$\begin{array}{rl}\mathbf{B}& =j\frac{\hat{\theta }}{2\pi r_{xy}}\\ \mathbf{J}=\nabla \times \mathbf{B}& =j\hat{\mathbf{z}}\delta (r_{xy})\end{array}$$

The $\hat{\theta }/r_{xy}$ vector field is the classic ‘twist around the origin’ vector field that points along $\hat{\theta }$ and at a right angle to ${\hat{\mathbf{r}}}_{xy}$ everywhere. It might look more familiar in cartesian coordinates:

$$\hat{\theta }=\frac{x\hat{\mathbf{y}}-y\hat{\mathbf{x}}}{\sqrt{x^2+y^2}}=\frac{x\hat{\mathbf{y}}-y\hat{\mathbf{x}}}{r_{xy}}$$

Then its curl is:

$$\nabla \times \frac{\hat{\theta }}{r_{xy}}=2\pi \hat{\mathbf{z}}\delta (x,y)=\hat{\mathbf{z}}\frac{\delta (r_{xy})}{r_{xy}}$$

That’s neat, I guess.

We can also do the magnetic field due to a single magnetic dipole (an infinitesimal magnetic or loop of current or classically-interpreted particle with spin) with magnetic dipole moment $\mathbf{m}$. We’ll use Gsponer’s notation of including a $\text{sgn}(r)$ term and see if it gives us some good delta functions terms. The vector potential is:

$$\mathbf{A}=\frac{1}{4\pi }\frac{\mathbf{m}\times \mathbf{r}}{r^3}\text{sgn}(r)$$

The magnetic field is (yes, definitely had to look up some identities for this):

$$\begin{array}{rl}4\pi \mathbf{B}=4\pi \nabla \times \mathbf{A}& =\nabla \times [(\mathbf{m}\times \frac{\mathbf{r}}{r^3})\text{sgn}(r)]\\ & =[(-\mathbf{m}\cdot \nabla )\frac{\mathbf{r}}{r^3}+\mathbf{m}(\cancel{\nabla \cdot \frac{\mathbf{r}}{r^3}})]\text{sgn}(r)-\frac{\mathbf{m}\times \mathbf{r}}{r^3}\times \nabla \text{sgn}(r)\\ & =[\mathbf{m}\cdot (-\frac{1}{r^3}+\frac{3{\mathbf{r}}^{\otimes 2}}{r^5})]\text{sgn}(r)-\frac{\mathbf{m}\times \mathbf{r}}{r^3}\times \hat{\mathbf{r}}\delta (r)\\ \mathbf{B}& =\frac{3(\mathbf{m}\cdot \mathbf{r})\mathbf{r}-r^2\mathbf{m}}{4\pi r^5}\text{sgn}(r)+(\hat{\mathbf{r}}\times \mathbf{m}\times \hat{\mathbf{r}})\frac{\delta (r)}{4\pi r^2}\end{array}$$

(Since we included the $\text{sgn}$ term that should track the delta functions for us, it seemed like $\nabla \cdot \frac{\mathbf{r}}{r^3}=\nabla \cdot \frac{\hat{\mathbf{r}}}{r^2}$ could be ignored now.) The latter term is the delta-function correction to the magnetic dipole field. The internet tells me that its integral over space is $+\frac{8\pi }{3}\mathbf{m}$, compared to $-\frac{4\pi }{3}\mathbf{p}$ for the scalar dipole field, and, as mentioned earlier, Jackson says this is responsible for the specific wavelength in the hyperfine splitting of hydrogen. Weird.

---

7. Summary

Most of these equations are the same toy examples from an intro electromagnetism course, written in a different way. But it is satisfying to see them written “explicitly”, which is what the delta functions let us do, instead of “working around” the delta function formulation by computing with e.g. Gauss’s law. I think things would have been easier to learn, back then, if the delta-function forms of these objects were made explicit from the start.

For posterity here’s a summary of the talked we’ve talked about:

Delta Functions in Radial Coordinates

$$\begin{array}{rl}{\delta }^n(\mathbf{x})& =\delta (r_n)/S_{n-1}\\ {\delta }^3(\mathbf{x})& =\delta (r_3)/4\pi r^2\\ {\delta }^2(\mathbf{x})& =\delta (r_2)/2\pi r\\ \delta (x)& =\delta (r_1)/2\end{array}$$

For the $1d$ case, recall that the “0-sphere” is a line segment whose “surface” area is usefully understood to be $S_0=2$, giving $\delta (x)=\frac{1}{2}\delta (r)$. That does seem a bit weird—why are we chopping our delta function in two?—but it does seem to work.

Integral Form of Divergence and Curl

Divergence in general is given by

$$\nabla \cdot \mathbf{F}=\underset{V\to 0}{\lim }\frac{1}{|V|}{\oint }_{\partial V}\mathbf{F}\cdot d\mathbf{n}$$

Curl is given by

$$(\nabla \times \mathbf{F})\cdot \hat{\mathbf{u}}=\underset{A\to 0}{\lim }\frac{1}{|A|}{\oint }_{\partial A}\mathbf{F}\cdot d\ell$$

for any plane with normal $\mathbf{u}$; chose $\mathbf{u}=\{\mathbf{x},\mathbf{y},\mathbf{z}\}$ to get the usual vector projections.

Functions whose divergence/curl/exterior derivative are delta functions

In ${\mathbb{R}}^n$ the divergence’s integrand includes a factor of $r^{n-1}$ from the coordinates, while all the angular coordinates naturally integrate to $S_{n-1}$. Therefore it’s $1/r^{n-1}$ that cancels the radial part out and produces a delta function at the origin:

$$\nabla \cdot \frac{\hat{\mathbf{r}}}{r^{n-1}}=S_{n-1}{\delta }^n(\mathbf{x})=\frac{\delta (r)}{r^{n-1}}$$

Which in ${\mathbb{R}}^3$ is

$$\nabla \cdot \frac{\hat{\mathbf{r}}}{r^2}=4\pi {\delta }^3(\mathbf{x})=\frac{\delta (r)}{r^2}$$

Meanwhile if curl is integrated around a loop in e.g. the $\mathbf{xy}$ plane, then the integrand includes a factor of the radius $\rho$ in that plane and is therefore canceled out by ${\rho }^{-1}$.

$$\nabla \times \frac{\hat{\theta }}{{\rho }^{-1}}=(0,0,2\pi \delta (x,y))=(0,0,\frac{\delta (\rho )}{\rho })$$

The analogs on other dimensions of $n$-spheres can be used to generalize these to higher dimensions, or to charge or current distributions that take lower-dimensional forms like lines or planes of charge.