Streamer discharge simulations in humid air: uncertainty in input data and sensitivity analysis

In our model, the electron density $n_{\mathrm{e}}$ evolves in time as

where $\mathbf{\Gamma}_{\mathrm{drift}}=-\mu_{\mathrm{e}}\mathbf{E}n_{\mathrm{e}}$ and $\mathbf{\Gamma}_{\mathrm{diff}}=-D_{\mathrm{e}}\nabla n_{\mathrm{e}}$ are the drift and diffusive fluxes of electrons, $\mu_{\mathrm{e}}$ the electron mobility, $D_{\mathrm{e}}$ the diffusion coefficient, and $\mathbf{E}$ the electric field. Furthermore, $S_{\mathrm{e}}$ is an electron source term due to reactions (e.g., ionization or attachment) described in section 2.5, and $S_{\mathrm{ph}}$ is the photo-ionization source term described in section 2.4. Electron transport coefficients are assumed to depend on the local electric field, and they were computed with BOLSIG$-$ [46] using the cross sections described in section 2.3.

We consider time scales up to about $200\,\textrm{ns}$ and found that including ion mobilities made no significant difference in our results. The motion of ions is therefore neglected, so that ions and neutral densities $n_{j}$ (numbered by $j=1,\ldots,n$) evolve in time as

where $S_{j}$ is a source/sink term due to reactions.

The electric field $\mathbf{E}$ is calculated in the electrostatic approximation as $\mathbf{E}=-\nabla\phi$. The electrostatic potential $\phi$ is obtained by solving the Poisson equation

using a geometric multigrid method [47, 48], where $\rho$ is the space charge density and $\epsilon_{0}$ is the vacuum permittivity.

We use a 2D axisymmetric computational domain measuring $200\,\textrm{mm}$ in the $r$ and $z$ directions, which contains a plate-plate geometry with a needle electrode (8 mm in length and 0.8 mm in diameter), as illustrated in figure 1. The bottom plate is grounded, whereas a constant high voltage $V$ is applied on the upper plate and the needle electrode. The axial electric field away from the needle electrode relaxes to the average electric field between the two plate electrodes, which is here defined as the background electric field $E_{\mathrm{bg}}$:

where $d=200$ mm is the distance between two plate electrodes.

For the electrostatic potential, a homogeneous Neumann boundary condition is used on the outer radial boundary. For all species densities, homogeneous Neumann boundary conditions are used on all domain boundaries, including the needle electrode. Secondary electron emission from electrodes due to ions and photons is not included.

As an initial condition, homogeneous background ionization with a density of $10^{11}\,\mathrm{m}^{-3}$ for both electrons and N${}_{2}^{+}$ ions is included. All other ion densities are initially zero. This background ionization provides the first electrons for streamer inception near the electrode tip, where there is significant electric field enhancement.

For computational efficiency the simulations are performed using adaptive mesh refinement (AMR). We use the same refinement criteria as in previous works (see e.g.[44]), namely to refine when $\alpha(E)\Delta x>1.0$ and to de-refine when $\alpha(E)\Delta x<0.125$, where $\alpha(E)$ is the field-dependent ionization coefficient and $\Delta x$ is the grid spacing. The minimum grid spacing in the simulations was $\Delta x=1.5\,\mu\mathrm{m}$. Around the needle electrode, we enforce $\Delta x<50\,\mu\textrm{m}$.

Electron-neutral collisions are the dominant process in streamer discharges, and energy-dependent probabilities of such collisions are described by cross sections. In a fluid model these cross sections are used indirectly through a Boltzmann solver to obtain electron transport coefficients (such as the mobility) and rates for electron-neutral reactions. BOLSIG$-$ [46] with the temporal growth model was used to obtain such data.

For most atoms and molecules, there is still quite some uncertainty in the electron-neutral cross sections. For example, in [49] streamer simulations were compared using different cross section sets for dry air (80% N₂, 20% O₂), and it was found that the streamer velocity could vary by about 10% when compared at the same position. Although the gas mixtures considered here consist mostly of N₂ and O₂, we here want to focus on the uncertainty in cross sections for H₂O. All simulations are therefore performed using Phelps’ cross sections for N₂ and O₂ [50, 51, 52], whereas for H₂O we consider the cross section sets described below.

We emphasize that the choice of H₂O cross section sets included in our comparison was rather arbitrary, and that we probably could have included other data as well (from e.g., [53, 54, 55]). However, our goal here is simply to illustrate how much variation there will be in streamer simulations using different H₂O cross section sets. For a more extensive discussion of H₂O cross sections and their accuracy, we refer to [56, 53, 57, 55].

Photoionization in humid air resembles that in dry air: excited N₂ molecules can emit photons with enough energy to non-locally ionize O₂ molecules, see e.g. [67]. We consider two approaches for modeling photoionization in humid air, which are described below. For comparison, we also briefly describe the dry air case.

In addition to electron-neutral cross sections, chemical reactions are another important input data for streamer fluid simulations. For many reactions, there are significant uncertainties in both the reaction form and the reaction rate coefficient, see e.g. [76, 77]. Furthermore, previous authors have used different sets of chemical reactions for humid air, and it is often unclear which set is the most suitable and which reactions are most significant.

In this work, we use the three chemistry sets that are described below. We have made small modifications to these sets so that they are compatible with the electron-neutral cross sections that we use. We have also included reactions to obtain the emitted light from the second positive system ($\rm N_{2}(C)\to N_{2}(B)+{\it h\nu}$), which is responsible for most of the optical emission in air at atmospheric pressure [78]. Furthermore, we have made corrections to reaction rates when necessary. All these changes are documented in B.

We remark that the choice of chemistry sets was also arbitrary, with the sets obtained from recent papers for streamer simulations in dry/humid air. There are other datasets, e.g., from [39, 41, 42]. Our purpose here is simply to demonstrate how sensitive streamer simulations are to the use of different datasets.

We use so-called sensitivity tests [83] to better understand which chemical reactions are important. These tests are performed by multiplying individual reaction rate coefficients with a factor $c_{i}$ for $i=1,\dots,N_{c}$, using the following values: $(0.0,0.25,0.5,0.75)$. With $N_{r}$ reactions, the total number of additional simulations is thus $N_{c}N_{r}$. We then pick a “quantity of interest”, here called $Q$, that we compare at some particular time to the case with unmodified rate coeffients $Q_{0}$. For reaction $n$ there are then results $Q_{n,i}$ for $i=1,\dots,N_{c}$, from which we compute normalized derivatives as

where $\Delta c=c_{i}-1$. Finally, we determine the mean of the normalized derivatives

and the sample standard deviation of the $Q^{\prime}_{n,i}$

Note that there are more sophisticated ways to do sensitivity tests in which multiple coefficients are varied at the same time, see e.g. [77, 84, 85]. However, our purpose here is simply to understand which reactions are important, and not to quantify the uncertainty in results due to uncertainties in the rate coefficients.

As described in section 2.3, five cross section sets for H₂O were considered, namely Kawaguchi, Triniti, Morgan, Phelps, and Itikawa. We remind the reader that rotational cross sections were not included in the Morgan and Phelps data, and that we left them out for the Itikawa data. In this section, we investigate the effect of these H₂O cross sections on positive streamers in humid air with 3% and 10% H₂O. The simulations were performed in a background field of 10 kV/cm, using the ‘dry air’ Zheleznyak photoionization model and the chemistry set of Malagón.

We first illustrate a typical time evolution in these simulations on the left of figure 4, which shows the electric field and electron density profiles for the streamer in air with 3% H₂O using the Kawaguchi H₂O cross sections. The streamer starts from the needle electrode and accelerates downwards, with the velocity and the optical diameter increasing approximately linearly with time, from about $0.8\times 10^{6}$ m/s to $2.3\times 10^{6}$ m/s and from about 1.2 mm to 3.7 mm, respectively. However, the maximal electric field at the streamer head remains approximately constant at about 110 kV/cm.

Simulations with different H₂O cross sections are compared on the right of figure 4, which shows the electric field and electron density profiles for air with 3% and 10% H₂O, at the moment the streamers reach $z=50$ mm. In addition, the maximal electric field, streamer velocity, and optical diameter are shown versus the streamer head position in figure 5.

In air with 3% H₂O, all five different cross sections do not lead to obvious differences in streamer properties, with variations below 10% when compared at the same streamer length. This is not that surprising, as differences in the electron transport coefficients are also small with 3% H₂O, as mentioned in section 2.3.5. Differences become larger with 10% H₂O in air, and then three groups can be distinguished: ① the Kawaguchi, Triniti and Phelps cases, which show good agreement; ② the Morgan case; and ③ the Itikawa case. Compared to the cases with 3% H₂O, streamers in group ① become about 40% slower and thinner, due to a much faster decay of conductivity in the streamer channel. The Morgan case shows even slower propagation at about $1.0\times 10^{6}$ m/s and a smaller diameter of about 1.4 mm, which is attributed to a higher attachment coefficient, as demonstrated in figure 3. An interesting phenomenon occurs in the Morgan case where a new positive streamer is generated from the needle electrode behind the previous one. This phenomenon, recently observed in strongly attaching gases in [15] as well, will be further illustrated in A. In contrast, the Itikawa case propagates much faster with a larger diameter, similar to those with 3% H₂O, but with a slightly lower $E_{\mathrm{max}}$. Although the attachment coefficient in the Itikawa data is nearly identical to that in the Morgan data, the ionization coefficient is about 5% higher. This suggests that ionization reactions play a more significant role than attachment, as discussed in A.

In the rest of the paper, we will use the Kawaguchi H₂O cross sections, as they can reproduce experimental swarm measurements and because they can be used with an isotropic scattering model, see section 2.3.1.

We now compare different photoionization models to understand how important changes in photoionization due to humidity are. Three approaches are considered: the Naidis and Aints models for humid air discussed in section 2.4, and for comparison also the standard Zheleznyak model for dry air described in section 2.4.1. We performed simulations with these models using the Malagón chemistry and the Kawaguchi H₂O cross sections, in a background field of 10 kV/cm. Figure 6 shows the electric field and electron density profiles when the streamers have reached $z=50\,\textrm{mm}$, and corresponding axial profiles are shown in figure 7. Furthermore, figure 8 shows streamer properties as a function of the streamer head position $z$.

There is a clear effect of photoionization on $E_{\mathrm{max}}$ and other streamer properties. If the standard Zheleznyak model is used, $E_{\mathrm{max}}$ is almost the same for air with 0%, 3% and 10% H₂O. However, other streamer properties such as velocity and diameter vary significantly, with a higher humidity leading to a reduced velocity and diameter. With the Naidis and Aints models there is less photoionization, which leads to an increase in $E_{\mathrm{max}}$ of up to about 15%. Due to the higher field at the streamer tip, the electron density in the streamer channel is also higher. Furthermore, the streamer velocity is substantially higher, for example about 30% higher with the Naidis model and 10% H₂O than with the Zheleznyak model. The effect on the streamer diameter is smaller, with differences of up to about 15% between the Naidis and Zheleznyak cases for 10% H₂O.

In addition, the reduction in photoionization with the Naidis and Aints models can cause the axisymmetric streamers to become unstable and eventually branch, which happens with the Aints model and 10% H₂O. Since we cannot simulate branching with the axisymmetric fluid model used here, this case was not included in the figures. Performing realistic simulations of the branching dynamics of streamers in humid air would require fully 3D simulations, as well as a stochastic treatment of photoionization [70, 86, 69], which we leave for future work.

That a reduction in photoionization can lead to faster streamers was also found in previous work [49]. Perhaps this is most easily understood from the scenario in which there is a very large amount of photoionization. In such cases, the streamer will have less field enhancement, and thus become less conductive, which means that it will propagate more slowly.

In this section, we study the effect of chemical reactions on positive streamers in humid air. We used the chemistry sets described in section 2.5 with varying humidity levels to perform a total of ten streamer simulations, for which the conditions are summarized in table 3. Figure 9 shows the electric field and electron density profiles when the streamers have reached $z=50$ mm, and figure 10 shows streamer properties as a function of the streamer head position $z$.

The chemistry set also has an important effect on streamer properties. We first look at the maximal electric field at the streamer head. For each humidity considered here, there are almost no differences in $E_{\mathrm{max}}$ with different chemistry sets. However, an increase in air humidity leads to a higher $E_{\mathrm{max}}$.

We then focus on the streamer velocity and diameter. As a first test, we compare results in dry air (0% H₂O) with the four considered chemistry sets, since these chemistries also differ in the included reactions for oxygen and nitrogen species. Streamers with the Komuro and Malagón chemistries agree very well in this case, whereas the Starikovskiy and Guo chemistries lead to slower streamers. At $z=50$ mm, the Starikovskiy case is about 8% slower, whereas the Guo case is about 12% slower. In streamer diameter these differences are about 10% and 14%, respectively.

Next, we compare results in humid air with 3% and 10% H₂O using the Komuro, Malagón and Starikovskiy chemistries. The differences in streamer velocity then become larger, with the Komuro case leading to the fastest streamers and the Starikovskiy case leading to the slowest ones. With 3% H₂O, the Malagón case is about 15% slower than the Komuro case at $z=50$ mm, and the Starikovskiy case is about 35% slower. With 10% H₂O, these differences are about 15% and 50%, respectively. The differences in streamer diameter are approximately the same as those in the velocity, with the slower streamer having a smaller diameter. Note that the Starikovskiy case with 10% H₂O is approximately a steady positive streamer, whose properties do not vary in time, as discussed in [14].

Regardless of the particular chemistry, an increase in the percentage of H₂O leads to slower streamers with a smaller diameter, as was also found in previous work [33, 34]. In the results presented here, streamers have propagated for about a hundred ns. On such short timescales, the main effect of a chemistry on streamer propagation is through reactions that affect the electron density inside the streamer channel, such as electron attachment, detachment and electron-ion recombination reactions. Below, in section 3.4, we will present a sensitivity analysis to identify which reactions in each of the chemistry sets are most important in this respect.

Using the method described in section 2.6, we have studied the sensitivity of streamer simulations to individual chemical reactions. The total number of electrons $N_{\mathrm{e}}$ at $t=100\,\textrm{ns}$ was here used as the “quantity of interest”. We exclude electron impact ionization reactions from this analysis; their effect is discussed in A. All simulations were performed in air containing 3% H₂O, using a background field of 10 kV/cm, the Aints photoionization model and the Kawaguchi cross sections for H₂O.

In figures 11–13, the ten most important reactions from the Komuro, Malagón and Starikovskiy chemistry sets are shown. ‘Most important’ is here defined as those reactions that lead to largest $Q^{\prime}_{n}$ (see equation (14)), which is an indication of the relative change in the quantity of interest with respect to a change in the reaction rate coefficient.

Almost all of the important reactions decrease $N_{\mathrm{e}}$, as indicated by the red color in the bar plots. The two three-body attachment reactions given by equations (5–6) are about equally important for all considered chemistries. However, there is a major difference in the importance of water clusters around positive ions. In the Komuro chemistry, the formation of single clusters through $\rm O_{2}^{+}+H_{2}O+M\to O_{2}^{+}(H_{2}O)+M$ and $\rm O_{4}^{+}+H_{2}O\to O_{2}^{+}(H_{2}O)+O_{2}$ is included, but these reactions have no strong effect on the later discharge evolution. In contrast, water clusters play a very important role in the Malagón and Starikovskiy chemistries. The Malagón chemistry includes reactions up to the formation of $\rm H_{3}O^{+}(H_{2}O)_{3}$, whereas the Starikovskiy chemistry includes reactions up to $\rm H_{3}O^{+}(H_{2}O)_{6}$. For both chemistries, the recombination of electrons with such clusters are the reactions with the highest sensitivity.

The standard deviations $\sigma$ of the ‘sensitivities’, as defined by equation (15), are given by the numbers next to the bars in figures 11–13. Note that there are several reactions related to water clusters that have a large standard deviation. This is caused by the following phenomenon: if there is a chain of reactions (as there is for water cluster formation), the simulation might not be sensitive to changes in the rate of one reaction in the chain. But if one reaction in the chain is fully disabled, the whole chain is stopped, which has a large effect on the result. In our tests, we include such a ‘disabling’ of reactions (namely by multiplying the rate coefficient by zero), which therefore leads to the large $\sigma$ values.

Sensitivity analysis can be a useful tool for reducing the number of reactions in a set. If we are only interested in the total number of electrons $N_{\mathrm{e}}$ at $100\,\textrm{ns}$, our results show that the effect of humidity in the Komuro chemistry can be described by two three-body attachment reactions. For the Malagón and Starikovskiy chemistries, the formation of water clusters around positive ions and the consecutive recombination of electrons with these clusters is also important. Since electrons predominantly recombine with the largest clusters, it might be possible to simplify the chemistry by leaving out the intermediate states.