Time-dependent one-dimensional simulation of atmospheric dielectric barrier discharge in N2/O2/H2O using COMSOL Multiphysics

The results of time-dependent one-dimensional modelling of a dielectric barrier discharge (DBD) in a nitrogen–oxygen–water vapor mixture at atmospheric pressure are presented. The voltage–current characteristics curves and the production of active species are studied. The discharge is driven by a sinusoidal alternating high voltage–power supply at 30 kV with frequency of 27 kHz. The electrodes and the dielectric are assumed to be copper and quartz, respectively. The current discharge consists of an electrical breakdown that occurs in each half-period. A detailed description of the electron attachment and detachment processes, surface charge accumulation, charged species recombination, conversion of negative and positive ions, ion production and losses, excitations and dissociations of molecules are taken into account. Time-dependent one-dimensional electron density, electric field, electric potential, electron temperature, densities of reactive oxygen species (ROS) and reactive nitrogen species (RNS) such as: O, O−, O+, $$ {\text{O}}_{2}^{ - } $$O2-, $$ {\text{O}}_{2}^{ + } $$O2+, O3, $$ {\text{N}}, {\text{N}}_{2}^{ + } $$N,N2+, N2s and $$ {\text{N}}_{2}^{ - } $$N2- are simulated versus time across the gas gap. The results of this work could be used in plasma-based pollutant degradation devices.


Introduction
Identifying and studying the plasma and active species characteristics have attracted much attention because of successful experimental studies on atmospheric pressure plasma applications such as processing technology, engineering [1], sterilization and surface treatment [2][3][4][5][6]. Many works have been devoted to simulation of the atmospheric pressure plasmas. For example, Gadkari et al. modeled a co-axial DBD plasma reactor in pure helium using a 2-D fluid model in COMSOL Multiphysics [7]. They investigated the influence of partial packing on the discharge characteristics of the dielectric barrier discharge in helium. Pan et al. used the fluid model to carry out numerically the evolution features of the atmospheric-pressure CF 4 plasma in a dielectric barrier discharge [8]. Abidat et al. studied a one-dimensional model of atmospheric pressure helium gas dielectric barrier discharge using COMSOL Multiphysics software [9]. Golubovskii et al, investigated spatio-temporal characteristics of the homogeneous barrier discharge in helium, numerically. They studied the influence of the elementary processes on the discharge by rate constants using a one-dimensional fluid model [10]. Various types of configurations for the derivation of electrical discharge at different working pressures have been reported. Dielectric barrier discharge (DBD) is one of the most famous working setup for a surface non-thermal electrical discharge. To prevent fast transition from atmospheric pressure electrical discharge to arc mode, one or two electrodes are covered by dielectric layers, therefore, they are kept away from the discharge contact [6]. Different processes such as electron impact ionization, photoionization, ionization by excited atoms and so on lead to stability of an electrical discharge [11]. Atmospheric pressure DBD as a process to generate active species, ozone, free radicals have been used for biological and medical applications such as sterilization and decontamination of environmental pollutions [12,13].
In this work, electrical characteristics of DBD in N 2 /O 2 / H 2 O (0.78/0.21/0.01) together with the plasma species were investigated by COMSOL Multiphysics 5.2 package. We performed time-dependent one-dimensional simulation of the DBD using reaction rates up to 99 reactions including effective processes such as electron attachment and detachment, surface charge accumulation, charged species recombination, conservation of negative and positive ions, ion production and losses, excitations and dissociations of molecules. Time-dependent one-dimensional electric field, electron density and temperature, densities of reactive oxygen and nitrogen species such as: N 2s and N À 2 , were evaluated versus time across the gas gap. It is well-known that the initiated reactions by reactive oxygen species (ROS) and reactive nitrogen species (RNS) play the key role in plasma disinfection and surface processing. Therefore, the results of this work could be used in plasma-based pollutant degradation devices.

Geometry description
In this paper, atmospheric pressure DBD is simulated. The plasma is driven by a sinusoidal alternating high voltagepower supply at 30 kV with a frequency of 27 kHz. The initial value for the electron density assumed to be 10 13 m -3 , the gas pressure 1.01 9 10 5 Pa and the gas temperature 300 K. The reduced electron mobility and corresponding energy were obtained by executing BOLSIG? [14,15]. The results were 4 9 10 24 [1/v.s.m] and 5 [eV], respectively. The mass fraction of the discharge gas was N 2 /O 2 /H 2 O (0.78/0.21/0.01). As shown in Fig. 1, the gas gap between the two copper electrodes is 0.2 mm. The radius of round electrodes is assumed to be 5 cm, thicknesses of the electrode and quartz glass were assumed to be 0.2 and 0.1 mm, respectively.

Equations for performing DBD simulations
In this section, we introduce the governing equations of fluid dynamics that were used in this simulation. Surface chemistry was used to take into account reactions of different species and rate of production and losses at the electrodes surfaces [16]. The electron density and the electron mean energy were calculated by solving the pair of the propulsion and propagation equations. Equations (1) and (2) show the electron continuity and flux equation, respectively.
on e ot þ r Á C e ¼ R e À ðũ Á rÞn e ð1Þ where n e is the electron density, D e electron diffusion coefficient, and C e is the electron flux, ũ average fluid velocity of the species, R e the rate of electron production. The electron flux consists of two terms: a term caused by the electric field and the other by the density gradient. The electron energy density equation is: where the term Ẽ Á C e indicates the amount of energy gained by the electron from the electric field. R e the energy rate stemmed from inelastic collisions which is obtained from the following equation: S en is the power dissipation, Q gen is the main source of heat and q is the electron charge. D e is the electron diffusion coefficient, l e is the energy mobility, and D e is energy diffusion coefficient: We used the Townsend coefficients with the assumption that the electron source can be calculated from the following equation: where M is the number of reactions, x j the molar fraction of the target species for the reaction j, a j the Townsend coefficient for the jth reaction and N n the total number of neutral particles. Assuming that p is the number of nonelastic collisions of an electron, then we have: where De j is the energy dissipation from the jth reaction. In non-electron species, the following equations are solved for each mass fraction: where w k is the density of ions,r Á j k the energy flux of the ions. The electrostatic field is calculated by the following equation: where e 0 is the permittivity of vacuum, e r is the relative dielectric constant. According to the boundary conditions for electron flux and energy flux, we have the following relations: The second term on the right side of Eq. (11) indicates that the electron is produced according to the secondary emission and c is the secondary emission coefficient. On the surface of the electrodes, ions and excited species become neutralized by surface reactions. Surface interactions on the electrodes are simulated using the coefficient b j , which it determines the probability of operating j species. For the ion species in the discharge, the conjugation equation is as follows: u and S i show the electrostatic potential and rate of change of electron density, respectively.

Simulation results
The evolutions of active species density, electron temperature and electric potential in the DBD gap are investigated in this work. Figure 2 shows the V-I characteristics of the electrical discharge across the gap in one period. In this figure, the change in current depicts breakdown of the gas gap. That is, in each cycle, there are two breakdown events.
As shown in Fig. 3, due to aggregation of positive charges on dielectric covering the grounded electrode, electric potential between two electrodes is changed.
The evolutions of the applied voltage across the gap are indicated in Fig. 4 in different times.
The spatio-temporal evolutions of the electric field and potential between the two electrodes are indicated in Fig. 5. As shown in this figure, in the presence of electrical discharge, electric field in the gap is reduced because of the plasma diamagnetic property.
Time evolutions of the electron temperature are appeared in Fig. 6 for center of the gap. It is seen that the electron temperature changes twice in each period due to the two discharge events in a cycle as indicated in Fig. 2. Spatio-temporal evolutions of the electron density are depicted in Fig. 7 versus the gap spacing. It is seen that the discharge occurs twice in each cycle, one in the positive and the other in the negative half cycle of the voltage. The  first discharge event initiates near the powered electrode corresponds to the right hand side of the picture. The maximum electron density is about 1.2 9 10 18 m -3 . It is seen that the electric field and the potential changes periodically versus position. The root mean square of the electric field can be used to investigate the sheath region [17,18].
The evolutions for the electron density and temperature distribution across the gap in different times are shown in Figs. 8 and 9, respectively. Figure 10 shows the two-dimensional simulation image of the electron temperature between the two electrodes. The electron temperature grows and diminishes corresponding to each discharge event.
Variations of N 2 ? and N 2 s (excited nitrogen) densities are indicated in Figs. 11. Figure 12 shows the distribution of atomic nitrogen and oxygen at the time of T.  The reactions due to electron impact inside the gas gap are given in Tables 1 and 2. The reactions of electron impact with active species of nitrogen including: metastable nitrogen N 2s , nitrogen molecule N 2 and singly ionized nitrogen molecule N 2 ? are depicted in Table 1 [14,19]. The reactions of electron impact with oxygen and water vapor molecules are shown in Table 2. In Table 2, the reaction rates R 1 till R 14 are for oxygen molecules while the rates R 15 till R 22 are for the water vapor molecules [14,19]. Two and three body reaction rates among the atoms and molecules are depicted in Table 3. The surface reactions are also given in Table 4.
Spatio-temporal simulation results of N 2s , N À 2 and N 2 ? densities across the discharge gap are shown in Figs. 15, 16 and 17, respectively. The variation of ozone concentration as an important radical species in process of sterilization and pollution decontamination is shown in Fig. 18. It is seen that the ozone concentration increases versus time.
The accumulation of the ozone inside the discharge gap is due to its long enough life time. The discharge event duration in the DBD is in the order of micro-seconds, while the ozone life time is in the order of several minutes [20].

Conclusions
Atmospheric DBD discharge has drawn great interest as a device for sterilization and decontamination of pollutions.
To identify and study the DBD characteristics and active species densities, we performed a simulation work by COMSOL Multiphysics package. To realize the density evolutions more precise we took into account 99 rate equations. The rates include elastic scattering of the electrons by the gas mixture molecules, electron impact excitation and ionization of the nitrogen, oxygen and water vapor molecules, two and three body collisions among the atoms, molecules and ions. It was revealed by examining the densities of the produced ROS and RNS that they are responsible for plasma disinfection and surface treatment processing based on their substantial number densities. Electron density and temperature evolutions were examined across the gap versus time. Spatio-temporal evolutions of the electric field and potential showed significant changes inside the discharge gap due to presence of the charged plasma species. One of the important results was the growing density of the ozone. It is produced in each discharge event and accumulated inside the gap. Its density