Underlying mechanism of the stagnation of positive streamers

Several publications have shown that it is challenging to model the stagnation of positive streamers. They find that as the streamers propagate, the space charge region at the streamer tip reduces in size and that the associated electric field increases towards infinity. In a recent study on streamer propagation in unsteady airflow, we circumvented this problem by relaxing the local density and electric field approximation in the drift-diffusion model, commonly used in the past models. The results, recently published, show that the electric field remains bounded during the streamer propagation. In the present paper, we explore the process of stagnation further with a more rigorous approach. We confirm that the instability in the electric field is an immediate effect of the local density and field approximation and that an extended description of ionization stabilizes the electric field and leads to a decelerating streamer. Finally, we discuss the role of positive ions in the stagnation and we show that the stagnating streamer velocity decreases till it becomes comparable to the ion velocities in the streamer head. This causes a broadening of the streamer head which leads to a sharp decrease in the streamer tip electric field and the streamer stagnation.


Introduction
Streamers are electrically conductive plasma filaments that may propagate in regions where the ambient field is below the breakdown field because of the high space-charge electric field around their tips. Streamers are the building blocks of electric sparks and lightning, and hence a better understanding of their initiation, propagation and termination processes can help us to better understand the dynamics of electric discharges. The stagnation length of streamers is particularly important, both for scientific and industrial purposes, since it defines the length a streamer can reach in a voltage gap, and therefore is a deciding factor for a breakdown event to occur. * Author to whom any correspondence should be addressed.
Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
Several reports have discussed the maximum length or the stagnation length of positive streamers experimentally [1,2] and numerically [3,4]. The numerical studies found that the positive streamers develop a sharp drop in the streamer radius, which in turn leads to an unstable and large electric field in the streamer tip [3][4][5]. In [6], we proposed that this unbounded electric field was a consequence of the local field and density assumptions made for the drift-diffusion equation, which we solved by smoothing the ionization source term using a length scale defined by the local electric field. The smoothing of the ionization front stabilized the electric field and led to a decreasing field in the tip and eventual stagnation of the positive streamer.
In this paper, we investigate further the findings of [6] by using an extended non-local model based on electron density gradients to relax the local density approximation in the calculation of ionization source term [7][8][9]. We show that the more accurate treatment of the ionization source term and its non-local nature allows us to achieve a stagnating positive streamer without the apparently unphysical and indefinite streamer tip electric field. We first simulate a streamer with the local description to illustrate that our numerical model suffers from the same instabilities as observed in [3][4][5], demonstrating that the instability problem is inherent of the local model. We then repeat the simulation with the non-local, extended description and show that the streamer tip electric field is well behaved and fades as the streamer stagnates. Finally, we discuss the physical mechanisms behind streamer stagnation.

Local model
We use the drift-diffusion model with the local field and density assumption for the charged particle dynamics commonly used in streamer modeling [10][11][12][13][14]. The ionization source term is calculated based on the local electron density and electric field. We have: ∂n n ∂t = ∇ · (D n ∇n n ) + ∇ · (μ n En n ) , where n e , n p , and n n are the densities of electrons, and positive and negative ions, and E is the electric field vector. μ and D are the mobility and the isotropic diffusion tensor; they are taken from [15,16] for ions and calculated with Bolsig+ [17] for electrons. The recombination loss terms are given by S r e-p and S r n-p with rates coming from [18]. ν i and ν a are the ionization and attachment rates and are calculated using the Bolsig+ software, and depend on the local electric field [17,[19][20][21][22][23]. The attachment rates at low electric fields are important, as they determine the electrical conductivity inside the streamers and consequently affect the maximum streamer length. The attachment rates at electric fields below 40 Td are taken from [18] and for higher electric fields the rates are calculated using Bolsig+ [17].
The photoionization source term S ph in (1.1) and (1.2) has been modeled based on the approximation proposed by [24] and is calculated based on the set of three Helmholtz differential equations described in [25]. A detailed description of the photoionization in our numerical model is provided in [6].

Non-local extended model
In the non-local extended description, the ionization source term is calculated such that it does not depend solely on the local electron density. The extended description is derived according to [7] as follows: The continuity equation for charged particles are derived from the Boltzmann equation: where Γ(r, t) is the charged particle flux and S(r, t) is the source term. The phase-space distribution function of electrons, f (r, v, t), can be expanded in powers of their density gradient as: where the superscript (k) means a tensor of rank k, and ⊗ is an inner dot product of order k. Similarly, S(r, t) and Γ(r, t) can be extended in powers of density gradients [7]. Truncating the expansion of Γ(r, t) and S(r, t), at k = 1, we have [7]: where W b = −μ b E and W f = −μ f E are respectively the bulk velocity vector and the flux velocity vector and D f is the flux diffusion coefficient. Finally and still following [7], we can rewrite the general diffusion equation as: This equation can be interpreted as the drift-diffusion equation (1.1), in which the local ionization source term is replaced by: (1.9) The non-local nature of this description lies in the inclusion of the density gradient term which weakens the locality assumption in the estimation of the ionization source term. Accordingly, the drift-diffusion equation for positive ions is modified: (1.10) A similar model has been developed in [9,26] utilizing the approach proposed by [8], where using perturbation theory they took into account the non-uniformity of the electron density in calculation of the ionization source term. The extended model described in [9,26] is: = ∇ · (D f ∇n e ) + ∇ · (μ f En e ) + n e (ν i − ν a ) ,where k 1 is calculated by [26]: (1.12) where μ f and μ b are respectively the flux mobility and the bulk mobility [26], and α f is the Townsend ionization coefficient [26]. Replacing the expression for k 1 (1.12) in (1.11), we have: which shows that the extended ionization source terms in the two models, (1.8) and (1.13), are identical. The photoionization source term and the transport and rate coefficients in (1.8) are calculated as for the local description. The electron bulk mobilities are calculated using the particle code described in [27], and are consistent with the values reported in [26].
The local and extended descriptions are implemented in the code, previously described in [6], with the same algorithms for adaptive mesh refinement [6] and the implicit time integration scheme with exact linear part integration (implicit-ELP) from [6,28].

Simulations
The simulations are performed in air at standard pressure and temperature. The configuration is axisymmetric and the domain is shown in figure 1. A potential of 7 kV is applied to the anode. The streamer is initiated by placing a Gaussian patch of neutral plasma on the electrode tip. The patch has a maximum density of 2.5 × 10 19 m −3 and an e-folding radius of 0.14 mm. The neutral background ionization is set to 10 9 m −3 for electron and positive ions densities in agreement with the background ionization proposed by [29]. It mimics the natural background ionization due to cosmic rays and radioactivity when attachment is neglected [30]. The photoionization source term included in our model raises the ionization level in front of the streamer head to levels higher than realistic background ionization densities. Other authors have studied the role of background ionization for streamer simulations without photoionization. In [31,32], it is respectively proposed that background ionization levels of 10 15 m −3 and 10 13 m −3 provide results in good agreement with simulations that include the photoionization. These levels are at least four orders of magnitude above our background ionization and hence the photoionization included in our model screens any influence of realistic variations in the background ionization level.
The smallest mesh size for the adaptive mesh refinement is set to 0.25 μm, to accurately determine the density gradient term in (1.8). The largest variable time step size is 5 ps. The implicit-ELP time integration scheme is chosen to remove the time step constraints from the Courant-Friedrich-Levy (CFL) stability condition and from the dielectric relaxation time [6,33]. The upwind spatial discretization scheme has been used. The boundary condition at the electrode surface is completely absorbing.

The local description
The simulation with the local description (1.1) developed the same instability in the tip electric field as described in [3,4]. Figure 2 shows the evolution of the streamer tip electric field as a function of the streamer length (top) and the evolution of the electron density profile along the streamer axis during its propagation (bottom). The electric field at different instances in time is shown on a colour scale in figure 3. Figures 2 and  3 show that the streamer narrows as it propagates and that the electric field at its tip increases. When the streamer reaches 1.65 cm at around 92 ns, the streamer rapidly narrows and the field at its tip gets unbounded, ending the simulation.

The non-local description
The same simulation is repeated but using this time the nonlocal extended model (1.8), and the results are shown in  figures 4 and 5. Firstly, we can note from figure 4 that the propagation velocity of positive streamers are lower when using the extended model than when using the local model as seen on figure 2 despite the fact the propagation velocity of negative streamers have been shown to be faster using the extended model [26,34]. This difference arises from the fact that in (1.8), the sign of the extended term depends on the streamer polarity which lowers the ionization source term in the positive streamer tip, while enhancing the ionization source term in the negative streamer tip and hence a higher electric field and propagation velocity for negative streamers.
Secondly, we see from figures 4 and 5 that the problem of the streamer instability is now solved. While the streamer narrows during its propagation, as in the local description, the electric field remains bounded, allowing the streamer to propagate a longer distance. The ratio of the streamer length to the electrode potential is 3.5 kV cm −1 , which is close to the experimentally evaluated value of 4.4 kV cm −1 reported in [35]. Figure 6 shows the streamer tip electric field (top left) and the streamer velocity (top right) as functions of time, the streamer diameter as a function of the streamer length (bottom left) and the streamer velocity as a function of the streamer diameter (bottom right). We define the streamer diameter as twice the shielding radius [36], which is the distance from streamer axis to the maximum radial electric field.
We see that the streamer velocity increases during the first few nanoseconds corresponding to the avalanche-to-streamer transition phase, and then decreases until the streamer stagnates at around 300 ns. As the streamer propagates, its diameter also decreases in agreement with previous simulations correlating the streamer diameter and its velocity [29,37]. The streamer tip electric field increases as in the local description, but at a lower rate, and reaches a maximum of 230 kV cm −1 when the streamer is 2 cm in length. At this point, the field sharply decreases and the streamer vanishes. This behaviour was also found in our previous work, where we smoothed the ionization source based on a length scale depending on the local electric field [6]. We also note that as the streamer velocity reaches zero the streamer diameter reaches its minimum diameter of 0.12 mm.
The relative importance of the electron density gradient term of the ionization source term (1.9) is illustrated in figure 7. The figure shows a comparison of the local ionization source term, the contribution of the density gradient term, and the total ionization term at different time instants. Note that the gradient term for a positive streamer is negative and that the plots show the absolute value of the gradient term. We see that as time progresses and the streamer becomes narrower, the density gradient in the tip becomes more important. Because the sign of the gradient term is negative, it counteracts the local source term and limits the ionization in high field regions, which stabilizes the electric field. In the model using local density approximation due to the lack of such stabilizing term the narrowing streamer focuses the ionization in an infinitely narrowing region in the tip which leads to an indefinite electric field.

Effects of non-uniform electric field on ionization source term
The non-local extended model described here, takes into account the non-uniformity of the electron density in the calculation of the ionization term. However, there is a rapid spatial variation of the electric field in the streamer tip which can influence the ionization term as well. Here, we investigate the behaviour of the positive streamer when we take into account the effects of non-uniformities in both the electron density and electric field as described in [8]. We have implemented the model described in [8], we have: where N is the gas number density andê is the electric field unit vector. We refer to this model as the fully extended model throughout this section. Similar to [38], we have taken the values of K 2 and K 3 coefficients from [8], while we have used (1.12) to calculate the K 1 coefficient. The ionization source term is: which depends on the spatial variation of the electric field as well as the electron density variations. Accordingly the positive ions equation has been modified to: 3) The term with K 2 coefficient, in our case, is much smaller than the other terms, and hence we have neglected this term. The model parameters and schemes are as described for the previous models. Figure 8 shows the electric field evolution and the electron density profiles along the streamer axis during the propagation of the streamer. The fully extended model (2.1) leads to a slightly faster propagating streamer as compared to the previous extended model (1.8), while they both result in slower streamer propagation compared to the local model. This is due to the positive contribution of the field-gradient source term to the ionization in the fully extended model (2.1). In figure 9, the different terms of the ionization source term (2.2) on the axis of the streamer are shown. Although the source terms due to the density-gradient and the field-gradient are comparable, they are a few micrometers apart and hence they do not cancel out. The field-gradient term leads to a wider ionization source term while the density-gradient term lowers the maximum local ionization. This effect resembles the model described in [6], where the ionization front in the streamer tip was smoothed to account for the non-local nature of the ionization. Moreover, the stabilizing effect of the density-gradient term that limited the local ionization term is not altered which describes the stable propagation of the streamer.
Previously in [39], it was shown for negative streamers that the ionization source term correction from the electron density gradient is more important than the correction from the electric field gradient. Since the results from the fully extended model (2.1) and the extended model (1.8) lead to qualitatively same results, we suggest that the density-gradient term is also more important for positive streamers and hence extended model (1.8) with only the density-gradient term is sufficient for the modelling of positive streamer stagnation. Moreover, as seen in figure 9 there is a sharp variation in the field-gradient term in the streamer tip that needs a very fine mesh to capture otherwise it can lead to numerical instabilities. Therefore, the use of the extended model with only the density-gradient term is numerically recommended.

Stagnation mechanism
From the previous section, we see that the stagnation process is accompanied by a sharp drop in the streamer tip electric field,  see figure 4 (top), similarly to the behavior described in [6]. In the following, we investigate the mechanism behind this final stage of stagnation which is accompanied by an increase in the streamer diameter and a sharp drop in the streamer tip electric field. Figure 10 shows the evolution of the charge layer thickness along the streamer symmetry axis, streamer diameter and the electric field at the streamer tip. We observe that the streamer diameter gradually decreases during the first 250 ns as the streamer propagates away from the electrode in regions of lower background electric field. This narrowing down of positive streamers has been previously described in [1]. Similarly we observe that the charge layer thickness decreases until 250 ns.
The charge layer thickness is linked with the effective ionization length, δ i , which is defined generally as where v e is the electron velocity and ν i and ν a are respectively the ionization and attachment frequencies [40]. In the Townsend avalanche, the effective ionization length defines the separation of slow positive ions and fast electrons [40]. In a streamer discharge regime, the effective ionization length varies significantly in the streamer front and we estimate its average value in the sheath region in the local description by: where the integration is taken above the breakdown electric field E k . Similarly, we define an equivalent average ionization length for the extended model, where we replace the local ionization source term by the source term from the extended model (1.9) to have: (3.2) Figure 11 shows the evolution ofδ i extended , the charge layer thickness, and the streamer tip electric field as functions of time. As the streamer propagates away from the electrode, δ i extended gradually decreases consistent with the decline in the charge layer thickness, reaching its minimum at 240 ns where it starts to increase. The streamer velocity decreases along with the streamer diameter, see figure 6. Similar relation between positive streamer velocity and diameter has already been shown experimentally in [2] and numerically and analytically in [29,37]. The streamer keeps narrowing down and the charge layer shrinks as it decelerates until the streamer velocity becomes comparable to the positive ions velocity in the streamer tip.   Figure 12 shows the streamer tip charge layer and the ion velocity on the streamer axis. As the propagation of the streamer becomes slow, the propagation speed of the charge layer becomes comparable with the speed of positive ions and the drift motion of the positive ions increases the thickness of the charge layer and the streamer head radius which in turn leads to a drop in the streamer tip electric field.

The role of positive ions motion
To understand better the role of the positive ion motion on the streamer stagnation, we have investigated the stagnation with and without the ion motion. Accordingly, we removed the drift term, ∇ · (μ p En p ), from the equation of motion of the positive ions (1.10), and we ran the simulation again. Figure 13 shows the evolution of the streamer diameter and of the charge layer thickness from 150 to 300 ns. We can see that without considering the ion motion, the streamer keeps narrowing down. With considering the ion motion, the streamer head broadens from around 250 ns. This result further emphasizes the role of positive ion motion on the behaviour of the positive streamer stagnation. Streamer models commonly assume that ions are stationary, we show that it is valid when the streamer velocity is large as compared to the drift velocity of ions. For stagnating streamer this is no more valid.
We should note that uncertainties in the ion mobility parameter does not alter the mechanism of the stagnation process or the final streamer length. The latter comes from the fact that the ion velocity becomes important only when the streamer is already too slow to propagate significantly before the full stagnation.
If we go back to the local model description for the study of streamer stagnation, we observe that the average ionization length and the streamer diameter keep narrowing down until the streamer tip electric field becomes unstable, see figure 14. The streamer velocity decreases until the end of the simulation, however it never reaches the velocities that are comparable to the ion drift velocities. Therefore, the local description does not capture the role of ions in the streamer stagnation and inevitably experiences an instability in the streamer tip electric field.

Conclusion
In this paper, we investigated the stagnation of positive streamers using two drift-diffusion models with different descriptions of the ionization source term. The first model uses the local density approximation in the ionization source term which result in an instability in the streamer head, where the electric field grows towards infinity, similar to [3][4][5][6]. The second model is a non-local fluid model [7,8,26] where the field remains limited and fades as the streamer stagnates. This result suggests that the instability of the streamer originates in the non-local nature of ionization in the high electric field regions which the local model fails to capture. Therefore, the use of the extended model or models involving higher moments of the Boltzmann equation [41] is suggested, which provide a more accurate estimation of ionization in the streamer tip. We also showed that the results from the fully extended model including both the electron density and electric field gradients agree well with the extended model including only the density gradient term. Moreover, we found that the motion of positive ions plays a major role in the stagnation of positive streamers. It suggests that streamer models with the assumption of immobile ions are not suitable for the stagnation of positive streamers.
This study led to the following observations on the stagnation process: • The streamer propagation is accompanied by a decrease in equivalent ionization length, charge layer thickness, and streamer diameter along with an increase of streamer tip electric field. The streamer velocity continuously decreases. • When the streamer velocity becomes comparable to the positive ions velocity in the streamer tip, the streamer charge layer thickness and streamer diameter start to increase leading to a drop in the streamer tip electric field eventually leading to the streamer stagnation.
To summarize the above-mentioned points, in this paper we have presented for the first time a stable modelling of positive streamer stagnation, and have described the underlying mechanism of stagnation. The model allows the calculation of the maximum length a streamer can reach in a laboratory discharge gap. One can note that the results presented here are obtained in the framework of non-local extended fluid models. Further investigation of the complex problem of positive streamer stagnation with higher order fluid or particle models might be required.