Influence of the angular scattering of electrons on the runaway threshold in air

The runaway electron mechanism is of great importance for the understanding of the generation of x- and gamma rays in atmospheric discharges. In 1991, terrestrial gamma-ray flashes (TGFs) were discovered by the Compton Gamma-Ray Observatory. Those emissions are bremsstrahlung from high energy electrons that run away in electric fields associated with thunderstorms. In this paper, we discuss the runaway threshold definition with a particular interest in the influence of the angular scattering for electron energy close to the threshold. In order to understand the mechanism of runaway, we compare the outcome of different Fokker–Planck and Monte Carlo models with increasing complexity in the description of the scattering. The results show that the inclusion of the stochastic nature of collisions smooths the probability to run away around the threshold. Furthermore, we observe that a significant number of electrons diffuse out of the runaway regime when we take into account the diffusion in angle due to the scattering. Those results suggest using a runaway threshold energy based on the Fokker–Planck model assuming the angular equilibrium that is 1.6 to 1.8 times higher than the one proposed by [, ], depending on the magnitude of the ambient electric field. The threshold also is found to be 5 to 26 times higher than the one assuming forward scattering. We give a fitted formula for the threshold field valid over a large range of electric fields. Furthermore, we have shown that the assumption of forward scattering is not valid below 1 MeV where the runaway threshold usually is defined. These results are important for the thermal runaway and the runaway electron avalanche discharge mechanisms suggested to participate in the TGF generation.


Introduction
The interest for the runaway mechanism has increased in the last 25 years due to the discovery of terrestrial gamma-ray flashes (TGFs) by the burst and transient source experiment (BATSE) on board the compton Gamma-ray observatory in 1991 [3]. The emissions are bursts of Γ rays of millisecond duration occurring in association with thunderstorms. The BATSE observation has described millisecond bursts, but the Reuven Ramaty High Energy Solar Spectroscopic Imager revealed later on that they were submillisecond bursts [4]. The TGF emissions are now understood as bremsstrahlung from high energy electrons accelerated in high electric fields associated with thunderstorms (see the review of Dwyer et al [5]), and the idea builds on Wilson's proposal that the electrons accelerated in thunderstorm fields may reach the runaway regime

Plasma Physics and Controlled Fusion
Influence of the angular scattering of electrons on the runaway threshold in air O Chanrion 1 , Z Bonaventura 2 , A Bourdon 3 and T Neubert 1 [6,7]. First, the TGFs were thought to be associated with sprites due to the high absorption of gamma rays in the atmosphere [3]. As the research went on, they appear to originate from cloud altitudes [8,9] where it is estimated that 10 17 runaway electrons are needed to be produced in order to observe a TGF at satellite altitude [8]. From observation, the typical lightning associated with TGFs is a strong positive intracloud [10][11][12]. Today, the two most accepted theories to justify TGF observations are the acceleration of cold electrons in the high field surrounding leaders [13][14][15][16][17][18][19][20] and the relativistic feedback mechanism enhancing the flow of electrons in relativistic runaway electron avalanches (RREAs) associated with the storm electric field [21,22]. It remains a possibility that the two mechanisms are at work concurrently. Modeling associated with TGFs and RREA has started with the works of [23][24][25][26] mostly looking at the generation of RREA from cosmic rays using a Fokker-Planck approach to solve the evolution of relativistic electrons. They addressed the problem analytically and by solving it with finite-difference methods. [27] looked at the RREA beam solving the Fokker-Planck equation with a Monte Carlo method. [28,29] and [8] intensively studied RREA and the feedback process, in particular, they used a Monte Carlo approach, followed by [30,31]. [13-17, 19, 20, 32] looked at the acceleration of cold electrons in discharges to justify TGF emissions. [13] used a 1D Monte Carlo model of a streamer. [14,16] looked at a 2D axisymmetric Monte Carlo model and later on at a hybrid fluid-particle energy-coupled model [20]. [15] used a hybrid fluid-particle space-coupled model. [2,18,19,32] looked at a full Monte Carlo model. In these models, cold electrons have to go through an energy runaway threshold to reach the runaway regime. The definition of the runaway threshold used in [13,14,24] assumes a forward scattering of electrons. [33] discuss the runaway electron mechanism with the nonlocal criterion for runaway electrons in laboratory discharges. In [34], based on a Monte Carlo model in N 2 and He, they have defined a stochastic runaway threshold ε th sto for large electric fields higher than the one corresponding to the maximum of the friction force, and they give the rate of the runaway electron in the absence of ionization. These authors took into account the influence of the angular scattering of electrons. The obtained stochastic runaway threshold is higher than the one with a forward scattering assumption. In the present work, we investigate the influence of the angular scattering of electrons in air on the definition of the runaway energy threshold for a large range of electric fields from 0.3 to 15E k with Fokker-Planck and Monte Carlo models.

The models
In this section, we introduce the four models that will be used to investigate the influence of the scattering on the runaway threshold. The first two models are based on the Fokker-Planck equation describing the movement of an electron drifting in the electric field. The collisional operator is simplified by averaging over many collisions, introducing an effective friction force that simply slows the electrons down, and a diffusion coefficient that affects the electron velocity angular distribution. For the first Fokker-Planck model, we impose forward scattering at all electron energies, thus, neglecting the diffusion coefficient. This reduced model gives the definition of the runaway energy threshold in which all electrons have forward velocities, i.e., antiparallel to the electric field. For the second Fokker-Planck model, we keep the diffusion coefficient, but we assume the angular equilibrium [1]. This allows us to reduce the problem and to define an energy threshold for runaways. Being interested in the runaway threshold, we drop the source terms in both Fokker-Planck models. For the sake of correctness, we use the same set of collision cross sections to derive the Fokker-Planck equation coefficients as in the Monte Carlo model. The Monte Carlo model uses the method of [35] to implement the scattering in the velocity angle knowing the ratio of the momentum transfer to the total cross section of individual collisions. To cancel the source terms in the Monte Carlo models, we remove the attachment process, and we do not create new electrons when ionization occurs. Instead, we simply take into account the energy loss and the deviation in angle.

Fokker-Planck models
Following [24,27], the Fokker-Planck equation reads where e, p, and μ are the electron charge, the magnitude of the electron momentum, and the cosine of the angle θ, respectively, between the electric field and the momentum vector. E is the magnitude of the applied field. F D is the friction force that represents the change in electron momentum due to collisions, and D represents the diffusion in angle due to collisions. The friction force F D is reconstructed as in [13,14] from a set of cross sections, obtained thanks to the work of [36,37] with the exception that we use the relativistic binary encounter Bethe (RBEB) model for ionization as in [2,38]. The diffusion coefficient D is approximated by one-fourth the rate of change in the electron mean square angle [27,39] where N g is the neutral gas density, v is the electron velocity, and χ is the scattering angle. The average square deviation angle ⟨ ⟩ χ 2 is calculated by an independant Monte Carlo scheme using the set of cross sections used in the present paper.

Fokker-Planck model with forward scattering.
Under the assumption of forward scattering, the diffusion coefficient D is null. Moreover, following [24], (1) can be rewritten in the form Particle orbits describing the movement of an electron through the phase space can be obtained by dropping the right hand side of (3). The electron orbits then are given by the Langevin equation, With a further assumption of electron drifting in the direction of the electric force, i.e., µ = −1, these equations become simply This equation gives the definition of the runaway threshold with the forward scattering assumption ε th , corresponding to the electron momentum or energy for which F D ( p ) = eE in the energy range from 100 eV to 1 MeV [13,14,24]. The plot of the friction force is given in figure 1, similar to [13,14]. For comparison, we plot on top, the friction force used in Lehtinen et al [1] (red curve) and the one given by the National Insitute of Standards and Technology (NIST) [40] (dashed black curve). We observe good agreement with both. The definition of the runaway threshold ε th is obtained when the friction force F D ( p ) (the blue curve in figure 1) intersects the magnitude of the Lorentz force eE( p ) (horizontal green curves in figure 1) for energies above 100 eV. A list of runaway threshold energies is given in table 1 for various electric fields expressed as a function of E k , the conventional breakdown field in air (E k = 32 kV cm −1 at ground pressure).

Fokker-Planck model without forward scatter-
ing. Without the assumption of forward scattering, the diffusion term in (1) cannot be neglected. Following the strategy of [27], we assume that the equilibrium in angles is achieved. This should hold, particularly, around the runaway threshold, where p varies slowly and the value of D is relatively high [1]. This allows setting the µ ∂ ∂ / term in (1) to zero and, furthermore, leads to the distribution function for the angular equilibrium in the form By integrating over all angles, equation (1) recasts to is the reduced distribution function and is the averaged cosine μ. From the expression for M( p ), we can see that the beam of electrons at angular equilibrium is forward when ( ) ∼ − M p 1, which requires that ( ) eE pD p . The comparison of the diffusion force pD( p )/e and the Lorentz force for different electric fields is presented in figure 2. For comparison, we present the force pD L /e with the diffusion coefficient based on the screened Rutherford model of cross section used in Lehtinen et al [1].
It can be seen from figure 2 that electrons drifting in high fields ( E 6 k ) must have an energy above 1 MeV and electrons drifting in weak fields (0.5E k ) must have an energy even above 100 MeV in order to satisfy the condition ( ) eE pD p for the forward beam. Clearly, these conditions cannot be fulfilled for models that intend to model the generation of runaway electrons from thermal electrons nor for models involving MeV electrons in a RREA.
Furthermore, the forward beam assumption is not valid close to the runaway threshold energy ε th , which is always below 1 MeV. The assumption that electrons drift forward as well as the forward scattering assumption is implicitly included in the definition of ε th . In the following, we, therefore, adopt the definition of runaway threshold energy [1,2] that does not rely on the forward scattering nor on the forward drift of electrons. gives the runaway threshold energy, th ε for a given electric field E. E k is the conventional breakdown field in air. For comparison, we present the friction force F DL used in Lehtinen et al [1] (red curve) and the one given by the NIST [40] (dashed black curve).
Looking back to equation (6), we now see that the distribution f p (t, p) follows the 1D transport equation inside a force field given by F D ( p ) + eEM( p ). Proceeding as in the previous section to obtain the electron orbits, we can write (6) as The sign of this force gives the definition of a runaway threshold in a similar way as in the case with the forward scattering assumption. Here, it means that electrons at angular equilibrium and above a certain threshold ε th eq will accelerate, and those below the threshold ε th eq will slow down. The comparison of the different forces present in this expression are given in figure 3.
The threshold definition ε th eq is obtained when the friction force F D ( p ) intersects the effective driving force −eEM( p ) for energies above 100 eV as shown in figure 3. The comparison of both the threshold ε th and the threshold ε th eq is presented in table 1 for various values of the electric field E. For the threshold ε th eq , we give both our value and the one obtained by Lehtinen et al [1] noted ε L th eq, . Note that the difference in thresholds is relatively significant. First, the threshold ε th eq is 5 to 26 times higher than the one with the forward scattering assumption, which emphasizes the role of diffusion in the threshold definition. Second, our threshold is 1.6 to 1.8 times higher than the one proposed by Lehtinen et al [1] and Celestin et al 19, depending on the magnitude of the ambient electric field. The latter difference comes mostly from the choice of cross sections from which the friction force and the diffusion coefficient are derived. It is important to recall the good agreement of our friction force observed in figures 1 and 3.
On the other hand, the comparison of our diffusion coefficient with the one used by Lehtinen et al [1] presented in figure 3 shows that our diffusion coefficient is 2.5 to 3.8 times higher. As an example, we plot the effective driving force −eEM L ( p ) calculated with the diffusion coefficient D L in figure 3 for a field of 0.5E k . It shows clearly that the combination of friction force and effective driving force from Lehtinen et al [1] gives a lower threshold than with our data and that the small deviations observed in the friction force are responsible for a small part of the difference observed in the threshold, the major part of which comes from the deviation in the diffusion coefficient.
The cross-sectional data we used were obtained from the work of Biagi [37], which gives a data set of total and momentum transfer cross sections validated experimentally by means of transport coefficients. The anisotropy function comes from the work of [41] calculating the cross section by using the Dirac partial wave method with screened potentials obtained from Dirac-Hartree-Fock atomic electron density. The anisotropy model is expected to be more accurate and more diffusive than the Rutherford model used in Lehtinen et al [1].
In the following, we give a fitting formula for runaway threshold energy ε th eq valid from 0.3 E k to 15 E k : keV.
k th eq 1.3292 (8) Since the bremsstrahlung process strongly affects the friction force above 1 MeV, a correction of the formula for fields less than 0.3 E k would have to be considered.   Note: A formula approximating ε th eq is given by equation (8).

Monte Carlo models
In the following, we briefly introduce the Monte Carlo model we use to compare with the preceding models. The model is described in detail in [14] and corresponds to a way to approximate the Boltzmann equation, where m e is the rest mass of an electron, γ = −v c 1/ 1 / 2 2 , v is the electron speed, c is the speed of light, and Q e is the Boltzmann collision term accounting for both elastic and inelastic collisions.
The Monte Carlo method allows for simulating the trajectory of electrons drifting in a constant electric field while modeling individual collisions by taking into account their stochastic nature, accounting also for the fact that each different collisional process, e.g., elastic collisions, excitation, or ionization, absorbs a different amount of energy of the electron. The model uses the set of cross sections of [37] with the exception that we use the model of RBEB for ionization as in [2,38]. The scattering angle after the collision is obtained following the scheme of [35]. Comparing with the Fokker-Planck model, the Monte Carlo model corresponds to a description of collisions at a timescale short enough not to be able to average over all collisions. Being interested in the runaway threshold and the electron orbits, we cancel the source terms in the model by removing the attachment process and by not creating new electrons when ionization occurs. We simply take into account the energy loss and the deviation in angle for this process. The model will be used in the following: (a) with an assumption of forward scattering, forcing particle velocity not changing direction in a collision; and (b) without this assumption as the classical Monte Carlo model.

Forward scattering models
3.1.1. Fokker-Planck with forward scattering. In this section, we present the orbits of the Fokker-Planck equation (3) with the forward scattering assumption, described by the Langevin equation (4). The solution of (4) is performed numerically as an initial value problem by lsode solver from ODEPACK [42]. The sample orbits are started by setting the initial condition for (4) on a circle with energy of 100 keV. Individual orbits are then followed as long as the electron energy stays in the range in between 100 eV and 100 keV. Electron orbits in phase space are shown in figure 4 for electric fields of 0.5, 2, 4, and 6E k . In each panel of figure 4, we observe that some electrons are trapped and slow down (thin blue orbits), while others run away (golden orbits). A curve made of points ( ) µ p, that realize the equilibrium between the friction force F D ( p ) and the Lorentz force µ eE is shown with a thick blue line on top of each panel of figure 4. It is clear that some electrons slow down until they cross this curve and then accelerate again to finally run away. Note that electrons drifting strictly against the field such that µ = −1 hold the definition of the runaway threshold with the forward scattering assumption (see table 1). The delimitation of golden orbits and thin blue ones can be used to define a probability to run away as a function of the initial angle θ 0 . The probability is simply one for the golden orbits (these will run away) and zero for the thin blue orbits (these are trapped), and the result for the electric field of 2E k is plotted in red in figure 5 as a function of the initial angle θ 0 for an initial electron energy of 20 keV. The steplike simplicity of the probability function comes from the fact that we have omitted random changes in the velocity angle and we consider the friction force instead of individual collisional events. Both effects normally cause single electrons to "jump" randomly from one orbit to another. As discussed in section 2.1, the assumption of forward scattering often is valid when an electron has a really high energy, i.e., high above 1 MeV, but it is questionable for electrons with energies around keV as emphasized below.

Monte Carlo model with forward scattering.
In the following, we present results from two simulations. In the first, electrons are injected into the homogeneous background electric field of 2 E k with an initial angle of θ 0 and an energy of 20 keV. Figure 5 shows the converged probability of running away as a function of the initial angle θ 0 . The probability is derived by counting the number of runaway electrons and dividing it by the number of electrons injected. Electrons are classified as runaway if their energy is above the threshold energy ε th . In figure 5, results from the Monte Carlo model are compared with the results from the Fokker-Planck model with the forward assumption described in the preceding section. Note that the inclusion of the stochastic nature in the treatment of collisions into the Monte Carlo model allows electrons to jump from one orbit to another and, thus, smooths the sharp probability obtained from the Fokker-Planck model as we can observe in figure 5.
In the second simulation, we initiated several monoenergetic beams of electrons injected forward with different initial energies. Figure 6 shows the probability to run away as a function of the initial energy of electrons drifting in homogeneous electric fields of 0.5, 2, 4, and 6E k . For 0.5, 2, and 4E k , the Monte Carlo results converged in less than 200 ps, and we show results obtained at convergence. For 6E k , as shown in figure 1, the driving force is very close to the maximum of the friction force. With the forward scattering assumption, low energy electrons are gradually lost and run away. To illustrate this, in figure 6, we show results every 30 ps from 0 to 330 ps, showing the probability increasing in time. Results from the Fokker-Planck model appear as vertical dashed lines, and these are simply the threshold energies for runaway ε th , delimiting the probability to run away, that is, by definition, equal to 0 when the electron is below the energy threshold ε th and 1 when above. Once again, we observe good agreement on the probability to run away between the two models. Again, introducing the stochastic nature of collisions by the Monte Carlo model simply smooths the probability curve without introducing too many other differences. This indicates that, for the case of forward scattering, the runaway threshold ε th is valid even for the Monte Carlo models. (6) is, in fact, the continuity equation for quantity p f p 2 and can be solved numerically by the monotonic upstream-centered scheme for conservation laws [43] using the finite volume method. The population of electrons was initiated by a set of  sharp Gaussian peaks centered at energies around the runaway threshold ε th eq (given in table 1) and separated by 0.1 keV. Looking at the trajectory of the electron population, we observe similar behavior for different electric fields. As expected, the electrons accelerate if their energy is above the threshold ε th eq and decelerate if their energy is below the threshold. Naturally, the force field F D ( p ) + eEM( p ) sign and magnitude in (6) drive the acceleration or deceleration of the electron population. As an example, in figure 7, we present the time evolution of p f p 2 for an applied electric field of 4E k .

Classical Monte Carlo model.
In this last comparison, we present the results of similar test cases as above but performed with the full Monte Carlo model. We recall that this model contains not only the stochastic nature of collisions, but also considers the scattering in the angle. In the first simulation, we, again, run the code initiated with 20 keV electrons injected with different angles into the homogeneous background field of 2 E k . The green curve in figure 5 shows the converged probability of running away as a function of the initial angle θ 0 . The results of the full Monte Carlo model are compared with the results of the Fokker-Planck model with the forward scattering assumption and with the Monte Carlo model with the forward assumption as described in preceding sections. We observe a radical change between the full Monte Carlo model that includes the complete description of the scattering processes and the other models where the assumption of forward scattering is employed. Even when shot forward, i.e., for θ π = 0 and with an initial energy of more than ten times the threshold ε th , 15% of electrons are lost from the runaway population by diffusion. On the other hand, about 15% of electrons shot antiforward will turn back inside the runaway regime due to diffusion.
If we integrate the results presented in figure 5 over all incident angles, we find that 31%, 33%, and 49.8% of electrons are lost from the runaway population for the Fokker-Planck model with the forward scattering assumption, for the Monte Carlo model with the forward scattering assumption, and for the full Monte Carlo model, respectively. This indicates that a significant number of electrons will abandon the runaway regime due to the angular diffusion from scattering.
In the second simulation, we initiated several beams of monoenergetic electrons injected forward with different initial energies around the threshold ε th eq , see table 1. Figure 8 shows the converged probability to run away as a function of the initial energy for electrons drifting in different electric fields. This time, the Monte Carlo results are compared with those from the Fokker-Planck model without the forward scattering assumption but with the angular equilibrium assumption.
Once again, we observe good agreement for the probability to run away between the two models. Again, for the Fokker-Planck model, the probability to run away is, by definition, 0 when the electron is below the energy threshold ε th eq and 1 above. Clearly, introducing the stochastic nature of collisions in the Monte Carlo model simply smooths the probability curve, indicating that the runaway threshold ε th eq is adequate for both models.
It is interesting to see that the probabilities obtained from the models with and without forward scattering are extremely different. This was expected according to the differences in runaway threshold ε th and ε th eq presented in table 1. For instance, it shows that less than 10% of electrons initiated at 10 keV (i.e., five times above the runaway threshold ε th ), really will run away in a field of 2E k (i.e., meaning that 90% of electrons will abandon the runaway regime due to the angular diffusion from the scattering).

Conclusion
In this paper, we have used several models describing the motion of electrons in a velocity space drifting in a constant  background electric field. The first class of models is constructed around the Fokker-Planck equation describing scattering with a friction force and a diffusion coefficient. After reduction, these models allow for defining the runaway threshold energy in an intuitive manner. The second class of models is constructed around the classical Monte Carlo method including the stochastic nature of the scattering. The latter models do not allow for defining the runaway threshold simply but are used to test the thresholds obtained by the Fokker-Planck model. Plotting the probability to run away for test cases assuming a forward scattering and not, we have found good agreement in the definition of the runaway threshold for both model classes. The assumption of forward scattering leads to the simple definition of the threshold.
Without this assumption, the Fokker-Planck model also allows for defining a threshold to run away if we assume an angular equilibrium. The threshold is found to be 5 to 26 times higher than the one with the forward scattering assumption, depending on the magnitude of the ambient electric field. The runaway threshold values are 1.6 to 1.8 times higher than those derived in Lehtinen et al [1] and Celestin et al [19]. For high electric fields, the runaway rates in air obtained with a Monte Carlo model by following the method of Bakhov et al [34] that deviates significantly from those derived in pure N by [34]. Both deviations are attributed to differences in the cross-sectional data between models and indicate that our model is more diffusive than the others. Furthermore, we have shown that the assumption of forward scattering should not be used below 1 MeV where the runaway threshold usually is obtained.
Those results suggest that the definition for the runaway threshold energy ε th eq , based on the Fokker-Planck model assuming angular equilibrium, is more appropriate than one assuming forward scattering. The outcome is particularly important for the acceleration of thermal electrons to the runaway regime and to the runaway electron avalanche. Those processes are suggested to participate as a source of runaway electrons that are vital for theories explaining observations of recently discovered terrestrial Γ-ray flashes.