On the stability of classical orbits of the hydrogen ground state in Stochastic Electrodynamics

de la Pe\~na 1980 and Puthoff 1987 show that circular orbits in the hydrogen problem of Stochastic Electrodynamics are stable. Though the Cole-Zou 2003 simulations support the stability, our recent numerics always lead to self-ionisation. Here the de la Pe\~na-Puthoff argument is extended to elliptic orbits. For very eccentric orbits with energy close to zero and angular momentum below some not-small value, there is on the average a net gain in energy for each revolution, which explains the self-ionisation. Next, an $1/r^2$ potential is added, which could stem from a dipolar deformation of the nuclear charge by the electron at its moving position. This shape retains the analytical solvability. When it is enough repulsive, the ground state of this modified hydrogen problem is predicted to be stable. The same conclusions hold for positronium.


Introduction
Stochastic Electrodynamics (SED) is a subquantum theory that considers the quantum vacuum as a true physical vacuum with its zero-point modes being physical electromagnetic modes, see the books [1,2].By construction, the vacuum is filled with classical electromagnetic modes, each carrying an energy 1  2 hω, where Planck's constant is a constant of nature characterising the energy stored in these modes.The task is then to show that SED explains all quantum behaviour of matter at the statistical level.This has been argued to occur on general grounds [1,2].SED is not ruled out by Bell's theorem, since that has an irreparable fatality, the contextually loophole [3].
Testing the working of SED in special cases is needed to get trust in the theory.Since harmonic oscillators perform well, several phenemona are explained: van de Waals forces, the logarithm of the Lamb shift between the hydrogen 1s and 2p states, the Casimir effect, the Unruh effect.A natural next case is the hydrogen ground state.de la Peña [4] and Puthoff [5,6] demonstrate that circular orbits are lead to stable atom: they do not fall onto the centre and they do not self-ionise.The SED theory has been tested numerically on the hydrogen ground state in a 2-d approximation where the orbit remains in its initial plane.Cole and Zou 2003 observed an encouraging agreement with the result from the quantum ground state wave function [7].New simulations have been carried out recently by our team [8] in a 3d approach, benefiting from a decade of progress in computational and programming power and employing a few analytical tricks to make the problem tractable in 3d.However, it was observed that in all runs and all attempts to model the system, there occurred self-ionisation.The latter fact stems with a theoretical prediction [9].Boyer puts forward that the problem may lie in relativistic effects [10], but taking into account relativistic corrections in our computer code did not change much and in particular did not cure the self-ionisation [11].
Faced with the renewed finding of self-ionisation, we search a theoretical underpinning for it and reconsidering the stability argument of de la Peña and Puthoff [4,5], we extend it to include eccentricity.We derive the gain and loss terms, averaged over the stochastic force and averaged over a full period in section 2. In order to study a more general problem with possibly a stable state, we extend in section 3 the problem with an attractive or repulsive 1/r 2 potential.We close with a discussion.In an appendix we consider the case where both masses are finite, such as the positronium.

On the stability of the hydrogen ground state
The nonrelativistic equation of motion of an electron with charge −e around an ion of charge Ze can be written in a dimensionless form where distances are measured in terms of the Bohr radius h/(αZm e c) and times in the Bohr time h/α 2 Z 2 m e c 2 [8], where α = 1/137 the fine structure constant and E(t) the fluctuating electromagnetic field that lies at the heart of SED.The term β 2 ... r arises from damping.We consider a time window during which the orbit remains basically unperturbed and in a plane, which is meaningful because of the small value of β ∼ Z/2000.The correlator of the E field reads in these units, with a factor 3/2 absorbed in E [8], and it is expected that the cutoff τ c ∼ α 2 Z 2 , corresponding in physical units to the Compton time h/m e c 2 , can be taken to 0 at the end.
The unperturbed orbit can be coded as Here L, ε and χ are the angular momentum, eccentricity and angle between long axis and the x-axis, all conserved in the Kepler problem, that is, in the absence of stochastic fields and damping.The energy of this orbit is E = The parameters ε, k and L are related by Let a code the angle φ according to with inverse Here the trace led to inner products of the e's with the connection, following from Eq. (2.5), The orbit may be written in the parametric form The angle, starting with φ = 0 at t = 0, evolves as where arctan takes values between −π/2 and π/2 and is discontinuous at a = (2n + 1)π for integer n; the floor function ⌊x⌋ delivers the largest integer below x.

Perturbations around the Kepler orbit
The energy changes due to the stochastic field and due to radiation.Due to the small value of α these effects can be calculated separately.In a perturbation expansion to order β 2 the rate of energy change due to the stochastic field follows from r(t) = r (0) (t) + βr (1) as With the leading term vanishing on the average, there remains the statistically averaged energy gain from the field We shall suppose that the stochastic fields and the damping start to act at time t 0 , which is later taken to −∞.The solution for the perturbation r (1) has the form with G the Green's function, to be constructed.Hence the average effect by the field, for times during which the Kepler orbit is not disturbed much, reads in terms of Ġ(t, s) ≡ ∂ t G(t, s), (2.14)

Perturbations of the orbit
We solve the perturbed Kepler dynamics in the rotating frame where the orbit is decomposed on the rotating eigenvectors e 1,2 and the fixed perpendicular vector e 3 .The field then has components Let us search the homogeneous solutions h = (r 1 , r 2 , z 1 ).The first is h (1) (t) ∼ ṙ(t).Another one, h (3) , has r 1 = 0, r 2 ∼ r(t) ∼ ρ a .Two further in-plane solutions h (2,4) are more intricate but can be constructed analytically.h (2) contains oscillations around a term linear in t.For z 1 one finds the solutions r cos φ and r sin φ from components of the orbit on the non-rotating basis.We present these six homogeneous solutions as 1) , ) where τ a = a − ε sin a = k 3 t and ρ a = dτ a /da = 1 − ε cos a.
The solution (2.13) then involves the Greens function G, which on the rotating basis is denoted as Γ, The conditions G(t, t) = 0 and Ġ(t, t) = 1 imply that also Γ(t, t) = 0 and Γ(t, t) = 1.Some effort leads to the elegant exact solution

.18)
Γ 33 comes from the transverse fluctuations, i.e., from h (5,6) , with the result Here φ ab = φ a − φ b with the rotation angle from (2.9).It is easily verified that W(t) = Γ(t, s)| s→t is a matrix valued Wronskian, that is, it is a constant matrix, indeed equal to 1.

Statistical rate of energy change from the stochastic field
The field contribution to the energy change reads from (2.14) (taking from now on t 0 → −∞) where and the dot denotes differentiation to t.The subtracted term in (2.20), cancelling Ġ(t, t) = 1, does not contribute to the integral.
From the result an expansion in b near b = a yields, when expressed as a series in t − s, whereas the in-plane part behaves as As seen in Eq. (2.25) below, the second and third order terms are potentially dangerous for the convergence of the integral at s = t, but they are exactly opposite, so that g(t, s) = O[(t − s) 4 ] for s → t.This behaviour assures that in the limit τ c → 0 the integral behaves well, as it should, since no subtle short-time effects are expected for the energy absorbed from the field.The combined effect takes the compact form with antisymmetric A and symmetric B, We express the energy gain from the field, averaged over one period, as where It is a quite smooth function, ending for circular orbits (ε = 0, κ = 1) at f (1) = 1 2 .The maximal value at ε = 1, κ = 0 is obtained by scaling a → κu, b → κv, which yields for κ → 0 (2.30) To evaluate this, one may analytically perform the v-integral, which leads to several terms odd and one even in u, the latter resulting in a result not far above f (1) = 1 2 : (2.31)

Radiative energy loss
The loss terms produce a rate of energy loss averaged over one orbit of period P = 2π/k 3 equal to where, in the averaging over a full period, we could omit a total derivative.Using dt = dφ/ φ = r 2 dφ/L and κ = kL this becomes Combining gain and loss terms leads to the total energy change For spherical orbits (ε = 0, κ = 1) there is the explicit result, first derived qualitatively by de la Peña [4] and Puthoff [5] It exhibits stability: for large negative energy E = − 1 2 k 2 (k ≫ 1) the energy gain is positive on the average, preventing collapse on the nucleus.On the other hand, when E is small negative, its average change is negative, preventing E ≥ 0, i.e., self-ionisation.
But with κ = kL → 0 at fixed L, Eq. (2.34) yields (2.36) Per period in the limit k → 0 at fixed L this implies a fixed amount, (2.37) In words: each revolution produces, on the average, a finite amount of energy, so that self-ionisation happens for orbits that have achieved a small k and an L < f (0).The pericentre of the orbit lies at The critical angular momentum value L c = f (0) = 0.588057 is not particularly small, and neither is the critical pericentre r c = 1 2 f 2 (0) = 0.172921 = 23.7α.Distribution of the dimensionless angular momentum L (in units of h) from the numerical simulation of the SED dynamics reported in ref. [8].Not much weight lies below L = 0.588, confirming that when such a value is reached at near-zero energy, self-ionisation may occur rapidly and the run is ended.Full curve: the distribution of L from the conjecture for the would-be stable ground state distribution, Eq. (3.31) for d, ℓ 0 → 0.
We are faced with the finding that the de la Peña-Puthoff stability argument of circular orbits [4,5] does not hold for all orbits, in particular not for very eccentric ones.This is in full accord with our recent numerical simulations [8].Here for short times orbits with reasonable statistics were observed, sometimes going close to the nucleus or near zero energy, followed by recovery towards less extreme orbits.But at some moment no recovery from a near-zero energy orbit was observed, with self-ionisation as a result.This characteristic is confirmed by the distribution of the dimensionless angular momentum L, obtained from these simulations: In Fig. 1, hardly any weight of orbits occurs lies below L < 0.588.

Adding an 1/r 2 potential
In the hope to find some stable behaviour, we extend the problem without spoiling its analytical solvability.We consider the presence of a "radially-directed dipolar force", stemming from a potential V d (r) = −d/2r 2 , which has the same form as the angular momentum part L 2 /2r 2 of the kinetic energy.Such a term may arise from a dipolar force between electron and nucleus, when the deformation of the nucleus, caused by the electron, is always aligned with the vector between them.

Analysing the orbits
In our classical dynamics approach the unperturbed problem has potential V and Newton force f while the random force and the damping lead to the dynamics For an unperturbed orbit with energy coded by its eccentricity ε ≡ √ 1 − κ 2 .The presence of d in (3.1) now leads to the connections so that now λ = κ.The range of ε and κ is, due to the shape of (3.3), at most the interval (0,1).
For the attractive dipolar force, d > 0, ε and κ indeed range between 0 and 1, while λ takes values between √ δ and √ 1 + δ.However, now there is a second family of orbits having L < √ d, where L 2 eff < 0 and ε > 1 (because κ 2 < 0).While such orbit are hyperbolic (unbound) for d = 0, they are bound for d > 0: these low-angular momentum orbits spiral into and and then out of the nucleus and need to be described within a relativistic framework [12].This property may explain the Darwin term (a δ (3) (r) term) in the relativistic corrections to the H ground state, related to our problem with d = α 2 > 0, as an effect of central spiralling.
The relation φ = L/r 2 now brings for orbits with L eff > 0 One period in a still takes a time P = 2π/k 3 , but it involves µ turns in φ, where µ is in general non-integer.For d > 0 each φ-turn takes the shorter time 2π/µk 3 , so though the solution r(t) remains, r(t) rotates faster.For d < 0 it rotates slower.For µ = 1, sin(φ a − φ b ) and cos(φ a − φ b ) do not reduce to compact forms like (2.7). With , the rate of energy loss from Eq. (2.32) generalises to The result is obtained using the properties (3.3) of the orbit.If d > 0 the orbits with L < √ d the orbits spiral into and out of the centre.They are described by the same solution, after adding absolute values to keep dφ/da and dt/da positive,

Perturbations of the orbit
The equations for perturbations r1 + ∇∇V • r 1 = E now read on the comoving basis The homogeneous solutions are still explicit, where also h (4) picks up a secular term.Indeed, τ a = k 3 t and φ a increase indefinitely with time, with ψ a ≡ φ a /µ − τ a a periodic function in a.

Statistical rate of energy change from the stochastic field
For the field contribution to the energy change we have to inspect the behaviour for s → t.While for d = 0 the problematic terms cancel, the problem is real for d = 0.The terms now add up to To have a well defined τ c → 0 limit, this should be regularised.Technically it can be done by adding at fixed t a contribution to the s-integral, which itself does not contribute to Ė at finite τ c : For the subtraction of the next order O[(t−s) 3 ] term, we point out that it yields a sin a log τ c divergency, which vanishes when integrating over a full orbit.

Circular orbits
At ε = 0 the cubic term of (3.13) vanishes, so after accounting for (3.14) the problem is regular.For d < 0 the maximal k is k m = 1/ |d|.The total energy loss rate is which prevents to pass below the lowest energy coded by k = k m and λ = 0.For small k and hence λ ≈ 1, the loss term dominates, exhibiting Puthoff-stability. For At small k this always shows the stability Ėtot < 0. For large-k the stability condition Ėtot > 0 demands that d < 1/4, a condition met curiously already in the quantum approach.So the dipole force does not essentially modify the stability of circular orbits.
In conclusion, circular orbits retain their de la Peña-Puthoff stability in the presence of the −d/2r 2 potential.

Very eccentric orbits
The d = 0 case (pure hydrogen problem) taught us that the remaining interest lies in the limit of very eccentric orbits with energy close to zero, where the possibility of self-ionisation looms.In the limit ε → 1, κ → 0 and small E = − 1 2 k 2 , we set k = kκ, λ = µκ, with the relation µ 2 = 1 + d k2 from Eqs. (3.3), (3.5).Expressed in angular momentum L, Eq. (3.5) As usual, the easy part is the loss term, which can be taken directly from (3.6).It scales as k 3 k5 , which implies in the limit k → 0 a finite contribution per period P = 2π/k 3 , as before, The energy gain from the field during one period is We scale a → κx, b → κy.Similar to Eq. (2.37) of the case d = 0, there results per period a finite limit of the energy gain, Even in this scaling limit this involves a lengthy expression, with φ(x, y) = 2µ(arctan x − arctan y) = 2µ arctan x−y 1+xy and numerators and the common denominator For the regulator needed to cancel the (t − s) 3 term in (3.13) we may take Indeed, when integrating it over y, a logarithmically divergent term appears, but the result vanishes over a full period, i.e., upon integration over x.Both gain and loss terms being proportional to k 3 implies that per period π ≤ a ≤ π there appears, averaged over disorder and over a period, an energy change Inserting k = (µ 2 − 1)/d, this becomes a function of µ alone, Self-ionisation is likely prevented when H(µ) < |d| for all orbits, that is, for all relevant µ.In the repulsive case d < 0 the range for µ is 0 < µ < 1.Since H(0) = 5.99 and H(1) = 0, our statistical argument suggest a stable bound state for d < d c = −35.8.This d c is finite, though rather large.
For d > 0 the physical domain is µ 2 ≥ 1 and µ 2 < 0. For µ 2 ≥ 1 we find that H(µ) is an increasing function.Its asymptotic behaviour can be analyzed.The limit µ → ∞ describes the orbits with lowest possible angular momentum L eff = 0, L = √ d.We can scale y − x → w(1 + x 2 )/2µ.(The leading and subleading terms in 1/µ can just be evaluated; for the second order correction a regularisation is needed, a subtraction of total derivatives of the form d[w + (w 3 + w) cos w + (w 2 + 1) sin w]/dw, with coefficients that depend on u.)This ends up with With this shape of H(µ), Eq. (3.27) predicts that these orbits remain stable only for the extremal value d c = 1 4 , but then (3.16) predicts that spherical orbits sink to the centre, so no stable cases are found for d > 0. The role of orbits spiralling into and out of the centre (the regime µ 2 < 0) is left as an open question.

Quantum mechanics
In a quantum approach one would introduce the angular momentum operator L2 and an effective one, L2 eff = L2 − d, the latter taking the eigenvalues l(l + 1) − d ≡ ℓ l (ℓ l + 1) for l = 0, 1, 2, • • • , so that ℓ l = − For d = α 2 Z 2 it has an analogy with the Schrödinger equation for the spinless relativistic electron in an H-atom [13].Indeed, it produces the Dirac square-root formula for the eigenenergies, be it that the total angular momentum operator J = L + S for the spinning electron is reduced to L and hence its eigenvalues j → l.The ground state is It is normalised according to 4π ∞ 0 dr r 2 ψ 2 0 (r) = 1.For the relativistic "Klein Gordon" H atom, d = α 2 Z 2 affects only the relativistic corrections; here we shall keep d as a parameter of order unity, positive or negative.

Classical phase space density
If a stationary ground state exists in the SED problem, it should be expressable as a function of the conserved quantities E and L, more precisely, as functions of R ≡ −1/E and L eff ≡ √ L 2 − d.Since this task was worked out by us for the ground state and excited states of the relativistic hydrogen atom [13], we can adjust the approach here.For (3.30) we coin the shape Let us verify this and fix the normalisation C; the result will be given in (3.36).The value of the Hamiltonian allows to denote the radial velocity p r and the effective angular momentum The fact that L eff rather than L itself enters here will imply its presence in (3.31).(The angle µ is this section should not be confused with the short hand µ = λ/κ in the remainder of the text.)As this implies dp r dL eff = dR dµ (r/R 2 ), the volume element in momentum space, with 0 ≤ ν ≤ 2π the azimuthal angle, reads Using L 2 eff R = 2r(R − r) sin 2 µ, we have Comparing with ψ 2 0 (r) from (3.30) we see that it has the proper shape, with Another quantity of interest is the distribution function of the conserved quantities E and L or L eff .As worked out in Eqs. ( 47) -(49) of Ref. [8], the relation between P EL eff and P pr = f (E, L eff ) is (3.37) With κ = L eff /L max eff = kL eff , this becomes which is normalised to unity, when taking all 0 < κ < 1 and 0 < R < ∞.Despite the setback for the SED program for d = 0 and all d > 0, it would be interesting to test this distribution for the regime d < d c = −35.8where a stable ground state of the problem should occur.

Discussion
It was put forward by de la Peña 1980 and Puthoff 1987 that circular orbits lead to stability in the hydrogen problem of Stochastic Electrodynamics (SED) [4][5][6].This gave hope for stability of the full problem, supported by the 2003 Cole-Zou numerical simulations of the dynamics [7].Our own, recent simulations improved on these, benefitting from new computer power, new architectures and analytical tricks.However, it was found that self-ionisation always occurs [8].This is explained by the present analytical approach.
We have developed the de la Peña-Puthoff statistical theory, where the average energy gain from the field per period of the orbit is evaluated and compared with the average energy lost by radiation.
Since both parts are defined by the unperturbed problem, the derivation is elegant and prone for study by students.
Our approach predicts indeed that with the electron and the nucleus modelled as point charges, self-ionisation takes place in the SED description of the hydrogen atom.Technically, the problem arises from orbits with nearly vanishing energy and moderately small angular momentum, below L = 0.5880h (we restore physical units).The perihelion then lies at distance r min = 1 2 (L/h) 2 a 0 = 0.173a 0 from the nucleus.The speed takes here its maximum of 0.024c, so the problem likely is insensitive to relativistic corrections, as we verified numerically [11].Indeed, when one looks at Fig. 1, the problem is not only the self-ionisation for L < 0.588, but moreover that the whole distribution is different from what one would expect from a conjecture based on the shape of the quantum ground state wavefunction.Hence if the SED program can be saved, either the gain or the loss term, or both, need to have a different shape.
Next we have added an −d/2r 2 potential, which for d = α 2 Z 2 has connections with the relativistic H problem for a spinless electron.The problem remains exactly solvable.In the repulsive case d < 0, |d| > d c = 28.6 we predict stability, which would be an interesting test for our numerical approach.In the attractive situation d > 0 we always find instability, even without accounting for effects from orbits which spiral into and out of the origin in a finite time [12].Those orbits nevertheless offer in principle a connection with the Darwin term, a relativistic delta-function correction to the Hamiltonian for the hydrogen problem, so that for some scholars the SED theory will retain a magic spell.
Our results support the conclusions reached by a number of authors regarding the failure of the "old" approach in SED to account for atomic stability: atomic orbits are not simply classical orbits perturbed by the stochastic field [9].For a proper functioning of SED, an intricate equilibrium state of both the electron and the stochastic field seems to be needed [1,2,14,15].

Figure 2 . 1 .
Figure 2.1.Distribution of the dimensionless angular momentum L (in units of h) from the numerical simulation of the SED dynamics reported in ref.[8].Not much weight lies below L = 0.588, confirming that when such a value is reached at near-zero energy, self-ionisation