A Novel Solution to the Klein-Gordon Equation in the Presence of a Strong Rotating Electric Field

The Klein-Gordon equation in the presence of a strong electric field, taking the form of the Mathieu equation, is studied. A novel analytical solution is derived for particles whose asymptotic energy is much lower or much higher than the electromagnetic field amplitude. The condition for which the new solution recovers the familiar Volkov wavefunction naturally follows. When not satisfied, significant deviation from the Volkov wavefunction is demonstrated. The new condition is shown to differ by orders of magnitudes from the commonly used one. As this equation describes (neglecting spin effects) the emission processes and the particle motion in Quantum Electrodynamics (QED) cascades, our results suggest that the standard theoretical approach towards this phenomenon should be revised.

The Klein-Gordon equation in the presence of strong electromagnetic fields with plasma dispersion relation is studied. The stability analysis reveals a bands structure, analogous to the Landau levels, for particles with energy comparable to the electromagnetic field amplitude. Furthermore, a novel analytical solution is derived for particles whose energy is lower than the electromagnetic field amplitude. Substantial deviation from the Volkov wavefucntion is demonstrated. As this equation describes the emission processes and the particle motion in Quantum Electrodynamics (QED) cascades, our results suggest that the standard theoretical approach towards this phenomenon should be revised. Introduction. At present days, several laser infrastructures with expected intensity of 10 24 − 10 25 W/cm 2 are under construction worldwide [1][2][3][4]. The experimental availability of such intense field sources creates exciting opportunities in many research fields [5], such as QED in strong fields [6], Schwinger mechanism [7], Unruh radiation [8,9], Nuclear physics and the search for dark matter candidates [10] as well as novel fast ignition schemes [5,[11][12][13] particles acceleration [14] and high harmonics generation [15].
The fundamental physics underlying all these scientific applications is the interaction of an intense electromagnetic (from now on we shall use the initials EM) field with an electron. The nature of the interaction is determined by the normalized field amplitude ξ ≡ ea/m and the quantum parameter χ ≡ (F µν p ν ) 2 /m 2 . The electron mass and charge are denoted by m, e respectively, the asymptotic momentum of the particle is p µ = (p 0 , p 1 , p 2 , p 3 ) and a = √ −A 2 is the amplitude of the vector potential A µ . The EM field tensor is F µν ≡ ∂ µ A ν − ∂ ν A µ . Natural units are used. If both ξ and χ are larger than one, the electron dynamics is both quantum and non-linear. The appropriate framework is the strong-Field QED [16]. Its basic principle is the inclusion of the term corresponding to the interaction with the classical laser field into the free part of the Lagrangian. As a consequence, the unpertubed states appearing in the cross section calculation are no longer free waves but particles in a classical field, represented by the Volkov wavefunctions [17].
Employing this approach, the properties of QED in the nonperturbative regime were thoroughly investigated through the years [18][19][20][21][22][23][24]. The lowest order strong-field * E-mail address: erez.raicher@mail.huji.ac.il processes are the non linear Compton scattering, where an electron interacts with the laser photons and emits a hard photon and the non-linear Breit-Wheeler process, where a photon decayes into a positron-electron pair in the presence of the EM field [18]. A sequential series of these processes, called "QED cascade" is followed by a rapid formation of a QED plasma whose ingredients are electrons, positrons and gamma photons [25]. Besides its fundamental significance, this phenomena attracts scientific attention as a possible laboratory astrophysics settings [26] and as a potential gamma ray source [27]. Furthermore, it was suggested that spontaneous cascades impose a limit (of about 10 25 W/cm 2 ) on the achievable laser intensity [28].
As the experimental exploration is not yet possible, our knowledge regarding this subject stems from numerical modeling. The kinetic calculation of QED cascades is based on a Monte Carlo technique describing the quantum processes mentioned above, integrated with a Particle-In-Cell code taking into account the collective EM field influence on the classical motion of the electrons [29][30][31][32][33][34][35]. The quantum rates were obtained with the Volkov wavefunctions. Two configurations were particularly investigated -standing waves formed by counterpropagating laser beams and the interaction of laser with a dense plasma. In both cases the vacuum assumption (namely, the vacuum dispersion of the EM wave) lying in the basis of the Volkov solution does not hold.
The main aim of this paper is to question the applicability of the Volkov waveufnction (and hence rates) under the conditions described in the previous paragraph. For this purpose, we shall examine the properties of the exact solutions in the relevant parameters range in comparison with the Volkov solution. We shall argue that these two wavefunctions are substantially different, meaning that the common theoretical treatment towards the QED cascades is fundamentally incorrect and has to be revised.
The Governing Equation and a Stability Analysis. The particles of interest in the context of QED cascades are fermions (electrons and positrons). Therefore, their wavefunction obeys the Dirac equation. However, the spin effect may be neglected in the typical cascades conditions [36]. Hence, for the sake of simplicity we shall treat them as scalars. Thus, the free equation of motion of the particle wavefunction Φ is the familiar Klein-Gordon equation in the presence of a classical EM field.
The center dot stands for Lorentz contraction. The laser field depends upon the spatial coordinate through φ ≡ k · x where k µ = (ω L , 0, 0, k z ) is the wave vector. We assume a circularly polarized field A(φ) = a(φ) e iφ + * e −iφ where the polarization vectors are = (e 1 − ie 2 )/ √ 2 and a(φ) is a slowly-varying amplitude. The dispersion relation of the electromagnetic field is massive-like k 2 = ω 2 L − k 2 ≡ m 2 ph . m ph is the effective mass of the EM wave photons, analogous to the plasma frequency of a classical plasma wave [37]. A selfconsistent calculation of the effective mass may be found elsewhere [38]. Following the standard Volkov derivation [39], we employ the ansatz Φ = e −ipx G(φ).
With the aid of the transformation G = yexp i k·p one can prove the equivalence of (26) to the Mathieu equation [40] y + [λ − 2q cos 2z] y = 0 where the following relations are used and z ≡ [φ + tan −1 (p 2 /p 1 )]/2. It is well known that the solutions space of the Mathieu equation can be divided into stable and unstable regions [41], according to the values of λ, q. For the sake of the stability analysis, we introduce the characteristics values of equation (3) a n (q) and b n (q). For a given q, the equation admits even periodic solution if λ = a n (q) and odd periodic solution if λ = b n (q). The curves a n (q), b n (q) seperate between the stable and unstables regions. An illustration of the characteristic shape of this curves may be found in [41] as well as in the Supplementary Material. The stability criterion determines that if b n (q) < λ < a n (q) the solution is unstable whereas the complementary domain a n (q) < λ < b n+1 (q) is stable. Asymptotic expansions of a n , b n are available for either q/n 2 1 or q/n 2 1.
In the first case we have a n ≈ b n , meaning that the solution is stable for all λ. However, in the opposite limit q/n 2 1 we have b n ≈ a n+1 so that the stable region is limited to narrow bands. The transition region should be studied numerically. This investigation shall be conducted with respect to q/a n (q) instead of q/n 2 for the following reason. The relations between λ, q and our physical quantites (4) clearly show that they are not independant variables. Hence, for given values of p µ , a the relation q/λ is immediate whereas n satisfying a n (q) < λ < b n (q) is unknown and requires numerical calculation. Fig. 1 shows the width of the stable (upper plot) and unstable (lower plot) regions. One can observe that for q/λ < 1/2 the forbidden gap a n −b n is negligible, meaning that the solution is stable. For q/λ exceeding 1/2, however, a n − b n rapidly rises while b n+1 − a n vanishes corresponding to the formation of a bands structure. The calculations were performed for increasing values of n in order to demonstrate that it converges to a universal shape. The quantities are normalized so that widths corresponding to different n values can be compared. Let us evaluate q/λ. For convenience reasons, we shall compute it in the frame of reference where the EM wave is standing (a rotating electric field), namely k µ = (m ph , 0, 0, 0). Therefore, (k · p) 2 /m 2 ph = p 2 0 and |p· | = p 2 1 + p 2 2 / √ 2. We consider a particle with p 3 = 0 so that the value of q/λ will be as high as possible. Consequently, we get p 2 1 + p 2 2 ≈ p 0 (assuming a relativistic particle). Thus the ratio is given by For p 0 ea as well as for the opposite limit p 0 ea we have q λ. The maximal value of the ratio q/λ is 1/ √ 2 and is obtained for ea = p 0 . The domain √ 2 − 1 < p 0 /ea < √ 2 + 1 (6) corresponding to q/λ > 1/2, is characterized by a bands structure (according to the above analysis). The allowed energy values p 0,n obey λ = a n . In the Supplementary Material we use Fig. 1 to deduce the gap between neighboring bands.
As we expect, in the limit m ph → 0 the spacing between the bands vanishes and the whole parameters space is stable. For the expected ELI [1] laser parameters (intensity of ∼ 10 24 W/cm 2 and optical frequency ∼ 1eV ) we have ea ≈ 10 3 m as well as m ph ≈ 10 −6 m leading to ∆ ≈ 0.03m. In order to better understand the above results let us consider m ph < 0, corresponding to a rotating magnetic field or to a wave travelling through a dielectric material. In this case the terms comprising λ (4) have opposite signs. This case was analyzed by [40,42,43] employing the asymptotic expression for q/n 2 1 and was found to exhibit a band structure. This is not surprising, since in the limit ω L → 0 the magnetic field is stationary and the spectrum is known to be quantized (the familiar Landau levels). Our case, however, corresponds to a rotating electric field and hence the emergence of a bands structure is unexpected.
Approximated Solution. Equation (26) admits approximate analytical solutions in several cases. If the asymptotic momentum satisfies p · A = 0, the quantum equation of motion may be solved analytically [44] even with the spin term included (i.e. Dirac equation). For particles with a large axial momentum p 3 (in the laboratory frame) the wavefunction can be found using a WKB approximation [45]. A solution for a discrete set of momentum values was derived in [46,47]. The limit of small m ph (corresponding to low density plasma) was studied as well [48].
In this paper, we construct a solution to (26) in the limit p 0 ea. The detailed derivation appears in the supplemental material. In the following we describe the outline of the proof as well as the final expressions. During the derivation we shall take advantage of the limit m ph → 0. For this reason, our starting point will not be the Mathieu equation, where m ph appears in the denominator and the above-mentioned limit is cumbersome. Instead, we start with (26). It bears an apparent similiarity with the equation of an electron in a crystal, if A µ (φ) is replaced by the periodic potential. We exploit this analogy and employ the Bloch theorem. Accordingly, the solution may be expressed as G = P (φ)e iνφ where P (φ) is a periodic function and ν is called the characteristic exponential. Substitution of this ansatz into (26) leads to − m ph 2 P + µP + δ + Qe iφ + Q * e −iφ P = 0 (8) where we have defined Since P (φ) is periodic it may be expanded in a Fourier series. In terms of the Fourier coefficients c n , the solution of (25) takes the form Recursion relation for c n may be found in the supplementary material. At the moment, the wavefunction spectral width ∆ is unknown. Let us assume that it obeys where In the following we shall obtain an explicit formula for ∆ and thus determine the validity range of our approximation. This assumption results in a symmetric distribution of c n and therefore allows us to obtain analytical expression for ν.
In order to simplify (48) we use the following argument. For m ph = 0, Eq. (26) becomes 1 st order and admits an analytical solution (the familiar Volkov wavefunction).
In the supplementary material we show that the coefficients c n take the same form as c 0n with the simple transformation (k · p) → (k · p)Ω. The analogy with (52) implies that the infinite series in (48) The simplified expression (55) is of great importance for two reasons. The obvious one is that it does not require any further calculation or recursion. The second is that its resembleness to the Volkov wavefunction allows us to conduct cross section calculations with the same analytical technique established in [18]. Furthermore, the analytical form of (55) provides us with an expression for the spectral width ∆ = ea(p · )/ [(k · p)Ω]. Comparing it with (49) we obtain the validity condition, which takes the following form (in the wave framework) In this framework the expression for Ω simplifies to Ω = 1 + (ea/p 0 ) 2 . Suppose the second term under the root is dominant, Ω ≈ ea/p 0 and (57) reduces to p 0 ea. In the opposite case, Ω ≈ 1 so that p 0 ≈ ea and (57) is not satisfied. Hence, our approximation is valid as long as p 0 ea, namely for electrons whose asymptotic-energy is smaller than the EM wave amplitude.
Having completed the solution derivation, its relation to the stability analysis should be addressed. For small asymptotic energies considered here, the characteristic exponential is given by (36). As p 0 increases, ν acquires a correction ∆ν. Beyond the stabilty threshold (6), ∆ν becomes complex and the solution turns unstable.
Results. The final expression (55) demonstrates that the variable that determines the deviation from Volkov solution is Ω, consistently with a previous publication [44]. In particular, the spectral width of (55) is Ω times smaller than the width of a Volkov wavefucntion (50) with the same p µ . For the expected ELI [1] laser parameters mentioned earlier the value of Ω lies in the range 1 − 10 3 . For the sake of demonstration only, the plots appearing below were calculated with different laser parameters. The reason lies in fact that the above parameters result in a huge spectral width making the wavefunction difficult for inspection. Nevertheless, the physical effect remains the same.
The spectral shape of the wavefunction (i.e. the Fourier coefficients c n ) is shown in Fig. 2 for p = (m/5, m/5, 0), ξ = 20, m ph = m/10, corresponding to Ω ≈ 20. The exact value of ν was found numerically and the coefficients were calculated according to the formulas in the supplementary material. The normalization is determined according to The deviation from the Volkov wavefunction is overwhelming. In addition, one can observe that the solid curve is symmetric, in agreement with the approximation employed in the analytical derivation. Fig. 2 depicts the exact wavefunction for p = (5m, 5m, 0), i.e. close to the edge of the stable region. At these conditions, the difference in the spectral width decreases with respect to the previous case. However, the spectral shape is deformed, implying that the symmetry between photon emission and absorption no longer exists. Discussion. As mentioned in the introduction, the emission rates embedded in standard QED cascade calculations relies upon the assumption that the particle wavefunction is approximately similiar to the Volkov wavefunction. In this paper we have challenged the validity this assumption. For an electron whose asymptotic energy is lower than the EM amplitude (p 0 ea) we have obtained an analytical solution (55). Its spectral width is Ω times smaller than the width of a Volkov wavefucntion (50) with the same p µ . It is well established that the spectral width of wavefunction is closely related to the number of laser photons interacting simultaneously with the electron. Consequently, the plot in Fig. 2 implies that the emission of an electron in a plasma will involve much softer photons compared to an electron in vacuum.
Furthermore, we have demonstrated that for an asymptotic energy in the range √ 2 − 1 < p 0 /(ea) < √ 2 + 1, a bands structure is formed, leading to a substantial decrease in the allowed phase space of the particle. This is especially interesting since this is expected to be the active region where the cascades are most probable to occur. For an energy slightly below the threshold corresponding to the bands formation, the wavefunction deviates from Volkov both in the width and in the shape as was shown in Fig. 2. The only regime where the solution of (26) agrees with Volkov is for p 0 ea, corresponding to an energetic electron beam colliding with a laser. The outcome of the above analysis is that the entire theoretical approach towards the laser induced QED cascades should be revised.

Stability Chart
In the main text we have introduced the characteristic values of the Mathieu equation, a n (q) and b n (q). Fig.  1 illustrates these curves for n = 4, 5. As expected, for q = 0 we have a 5 = b 5 = 25 as well as a 4 = b 4 = 16. As q increases, b 5 seperates from a 5 and coalesces with a 4 . The domain restricted by a 4 , b 5 is stable while the one lying between b 4 and a 4 as well as between a 5 and b 5 is unstable. Generally speaking, for q/n 2 1 the solution is stable and for q/n 2 1 the allowed region is limited to narrow bands. For n 1, the transition between both regimes becomes very sharp as demonstrated in the main text.

B. THE ANALYTICAL SOLUTION
In the following we describe the detailed derivation of the analytical solution discussed in the main text. The wavefunction Φ obeys the familiar Klein-Gordon equation in the presence of a classical EM field.
The notation is the same as the one employed in the main text. Following the standard Volkov derivation [39], we employ the ansatz Φ = e −ipx G(φ).
where Q ≡ −2ea(p · ) and * denotes complex conjugate. According to Fluquet theorem, the solution may be expressed as where P (φ) is a periodic function. Substitution of this ansatz into (26) leads to − m ph 2 P + µP + δ + Qe iφ + Q * e −iφ P = 0 (28) where we have defined Since P (φ) is periodic it may be expanded in a Fourier series Substituting (30) into (28) yields with ρ n ≡ m ph 2 n 2 + iµn + δ. Using the definitions of µ, δ we get The equations set represented by (31) may be cast into a matrix form N is a natural number larger than the width ∆ of the spectral distribution. The characteristic exponential is found by the requirement The solution of this equation provides us with the characteristic exponential ν. For Q = 0 it yields where In the vacuum case, ν does not depend on whether Q vanishes or not. We argue that the same reasoning approximately applies to the plasma case. This assumption shall be verified later on. Substituting (36) into (32) we get The next step is to find the coefficients c n . For this purpose, let us look at the first equation of the matrix system (34), i.e.
It yields the following relation From the next equation in the linear system we get The general pattern is apparent the recursion relation we look for is with the definition Analogously, for n > 0 we have as well as During the derivation of the recursion relations (42)(43)(44), the entire matrix system (34) was used besides the following equation Consequently, it can be used as an alternative consistency condition instead of (35). The first term in (46) disappears since ρ 0 (ν) = 0. Writing c 1 , c −1 in terms of c 0 we obtain In terms of the coefficients c n , the solution takes the form At the moment the wavefunction spectral width is unknown. Let us assume that it obeys Therefore, the second term in (38) is negligible and ρ n is anti-symmetric with respect to n. As a result, (47) is automatically satisfied, justifying the expression for ν (36). In the following we shall obtain an explicit formula for ∆ and thus determine the validity range of our approximation. Now we would like to find a compact formula for the infinite series in (48). For this purpose an identity deduced from the Volkov limit shall be exploited. This limit is obained by taking m ph → 0. In this case (26) becomes a 1 st order equation admitting the analytical solution [e(p · A)/(k · p)] dφ (50) where ν 0 ≡ −(e 2 a 2 )/(2(k · p)) (51) The series expansion for the coefficients c 0n corresponding to (50) may be obtained by substituting m ph = 0 in (38), leading to ρ n = −2n(k · p), and utilizing the recursion relation (42)(43)(44). It follows that [e(p · A)/(k · p)] dφ (52) Now we shall take advantage of the fact that the coefficients ρ n (and therefore also c n ) take the same form as in Volkov solution with the simple transformation (k · p) → (k · p)Ω (53) The analogy with (52) implies that the infinite series in (48) converges to the Volkov-like expression The analytical form of (55) provides us with an expression for the spectral width ∆ = ea(p · ) (k · p)Ω (56) Comparing it with (49) we obtain the validity condition, which takes the following form (in the wave framework) In the wave framework the expression for Ω simplifies to Ω = 1 + (ea/p 0 ) 2 (58) Suppose the second term under the root is dominant, Ω ≈ ea/p 0 and (57) reduces to p 0 ea. In the opposite case, Ω ≈ 1 so that p 0 ≈ ea and (57) is inadequate. Hence, our approximation is valid as long as p 0 ea, namely for electrons whose asymptotic-energy is smaller than the EM wave amplitude.
In order to normalize the solution we look at the expression for the electric charge associated with the wavefucntion.
Substituting (55) to (59) yields a constant (up to a negligible term of the order ∆/∆ s ). The normalization factor may be determined by requiring q = 1.