Why are S-state solutions of Dirac equation for hydrogen atom allowed to be unbounded?

Question

In this Wikipedia page there's a statement about 1S orbital as solved from Dirac equation:

Note that $\gamma$ is a little less than $1$, so the top function is similar to an exponentially decreasing function of $r$ except that at very small $r$ it theoretically goes to infinity (but this behavior appears at a value of $r$ smaller than the radius of a proton!).

For comparison, when solving radial part of Schrödinger's equation for hydrogen, we explicitly strive to avoid solution which is singular at the singular point of potential.

In Dirac equation it seems that the only justification of allowing unbounded solution is a physical one, not mathematical. It looks very inconsistent with what I've seen when solving various problems with Schrödinger's equation.

So why is Dirac equation's solution allowed to be unbounded? What are the boundary and smoothness conditions on the wavefunction?

I suspect that Dirac equation simply doesn't admit S states with bounded wavefunction in Coulomb potential. Can it then be shown that if we start from a spherically-symmetric wave packet with a zero in its center, then its evolution in Coulomb field will inevitably lead to formation of similar $r^{\gamma-1}$ singularity?

Answer 1

The whole idea of always looking for non-singular solution of a PDE when solving physical problem is flawed. There are legitimate cases when a singular solution does solve the equation and satisfies all the boundary conditions originally imposed on it — and — solves the physical problem as well.

Similarly, the condition of square integrability of wavefunction, mentioned in akhmeteli's answer, is not enough to get e.g. discrete energy spectrum of a bound particle. Consider e.g. spherical/cylindrical Neumann function as a solution of radial equation for a particle in spherical/cylindrical box with constant potential inside. Such a solution is square integrable with corresponding Jacobian, but clearly results in continuous energy spectrum if we take it as a solution of such a problem.

What we should look for is the solution which solves the PDE as a distribution. Namely, it must take into account all distributional singularities like Dirac delta arising when applying the differential operator to it. As an example, consider spherical Neumann function as a candidate solution for the problem mentioned above for S state. It has the following series expansion:

$$y_0(r)=-\frac1r+\frac r2-\frac {r^3}{24}+\cdots.$$

It's easy to see that it's normalizable when integrated in the box of radius $R$, since the Jacobian $r^2$ cancels the $\frac1{r^2}$ term in $y_0(r)^2$. But this function has a serious problem: applying the kinetic energy operator to it yields a Dirac delta, which isn't compensated by anything in the equation. It becomes obvious if you consider that $\frac1r$ is (up to a constant factor) a solution of Poisson's equation for point charge. Another way to see it is to verify that

$$\lim_{\varepsilon\to0}\int_{\mathbb R^3}\nabla^2\left(\frac1{\sqrt{r^2+\varepsilon^2}}\right)\,\mathrm dV\ne0,$$

recalling the definition of Dirac delta distribution.

To summarize: the premise of the question is false: we don't try to avoid singular solution to Schrödinger's or other similar equations — we try to make sure that our solution solves our equation in a distributional sense, not leaving any stray distributional singularities. Thus it's completely OK to allow solutions of Dirac equation to be unbounded, provided that they do solve it exactly as distributions.

Answer 2

I would say the bound state solution should have integrable probability density near the center, which condition does hold (a factor of $r^2$ appears when you integrate in the spherical coordinates), although the wave function (and the probability density) is indeed singular. This condition is used at page 20 of the article @LewisMiller refers to in a comment ("Schrödinger and Dirac equations for the hydrogen atom, and Laguerre polynomials").