A Two-Dimensional Air Discharge Modified Model Under Unipolar Square Pulse Voltage at Low Temperature and Sub-Atmospheric Pressure

In this paper, a modified COMSOL Multiphysics model with artificial stability is adopted, and the boundary condition developed in the previous paper are applied to the calculation of the photoionization rate. The main processes involved in air discharge were simulated by 39 plasma chemical reactions composed of 12 kinds of particles. The characteristics of low-temperature sub-atmospheric pressure air discharge under unipolar square wave pulse voltage with 13kV applied amplitude, 25 kHz-50 kHz frequency and a duty cycle of 50%-75% are discussed. The results show that: When the duty cycle increases, the average density of electron, N<sub>2</sub><sup>+</sup>, N<sub>4</sub><sup>+</sup>, O<sub>2</sub><sup>+</sup>, O<sub>4</sub><sup>+</sup>, and N<sub>2</sub>O<sub>2</sub><sup>+</sup> show an increasing trend, and the average electron temperature changes little; When the frequency increases, the average electron density, the average density of N<sub>4</sub><sup>+</sup> and O<sub>4</sub><sup>+</sup> show a downward trend, the average density of N<sub>2</sub><sup>+</sup>, O<sub>2</sub><sup>+</sup>, N<sub>2</sub>O<sub>2</sub><sup>+</sup> change little; The average electron temperature shows an increasing trend, but the duration of high electron temperature becomes shorter. When discussing a single period, it is found that the potential difference between the electrode and the plasma determines the direction of the electric field in the space. For the photoionization reaction, even if its influence on the negative streamer is not important, it still needs to be considered in the simulation. The effect of temperature on discharge is also taken as the research emphasis, it includes the following aspects: Firstly, at low temperature, the collision reaction and attachment reaction rates are increased, while the recombination reaction rate is almost unaffected; Secondly, low temperature results in the decrease of secondary electron emission coefficient, which indirectly increases the electric field gradient and the peak value of electric field; Finally, at low temperature, the mobility of carriers, including electrons, positive and negative ions, is increased, but the streamer velocity is decreased.


I. INTRODUCTION
As one of the main forms of gas discharge, streamer discharge has always been a fundamental issue that has attracted much attention in electrical engineering. With the development of science and technology, the links between disciplines have become closer, making gas discharge involve more and more fields. As the author mentioned in the previous paper [1], The associate editor coordinating the review of this manuscript and approving it for publication was Yee Sin Ang .
with the further development of China's aviation industry, under the premise of satisfying safe operation, studying the air discharge characteristics of aircraft in a low temperature sub-atmospheric environment has become a common issue between aviation and high voltage engineering. It is of certain significance to promote the development of China's aviation technology.
No matter in time or space, streamer discharge is a transient process with drastic changes. The time scale is usually on the order of nanoseconds. Although the current experimental technology has met the research of some microscopic parameters of gas discharge to a certain extent [2]- [5], the net charge density, electric field distribution, streamer radius and other parameters cannot be fully recorded. Therefore, physical modeling and numerical simulation have become efficient tools for studying these important parameters. The scholars at home and abroad have been committed to establishing a physical model that can completely describe gas discharge, and by comparing the similarities and differences between numerical simulation results and experimental results, gradually improving the basic theory of streamer formation and development, thereby promoting the development of gas discharge research.
With the rapid development of computer technology, it has gradually become possible to use high-performance computers and optimized numerical algorithms to realize the numerical simulation of the physical process of gas discharge. In the middle of the last century, some scholars performed precise dynamic descriptions of the space charge region in the gas discharge process, and proposed a Monte Carlo particle model [6]- [9]. The model uses the Monte Carlo method to describe the random walk process of particle generation, initial state, space, energy, and direction using random sampling or statistical methods, and records the analysis results to obtain the streamer discharge development characteristics. Although streamers can be simulated more realistically, because electrons and ions are analyzed separately, it is impossible to track the full contribution of particles to the discharge. If we want to consider the motion of various particles in the gas medium, this method is difficult to realize because of its large calculation amount, especially when the density of the gas discharge plasma is high in the atmospheric pressure range. Compared with the Monte Carlo particle model, the fluid model first proposed by Davies has been widely used due to its high computational efficiency [10]- [12]. Passchier and Goedheer [13] used a one-dimensional fluid method to simulate the air discharge process between parallel plates by using the particle migration-diffusion approximation, but the model assumes that the electron energy distribution function is in thermal equilibrium. Lymberopoulos and Economou [14], [15] used a fluid model to study the development process of argon discharge between parallel plates, and the results show that although the number of metastable argon atoms in the discharge process is small, they have a greater impact on the discharge process. Bera et al. [16], [17, [19], Bera [18] studied the discharge process of radio frequency inductively coupled plasma at low pressure using a fluid method. In the radio frequency discharge model, it is assumed that the electron energy is distributed in a Maxwellian mode and the chemical reaction process between particles in the air is not considered. The transport equation, electron energy equation and Poisson equation of electrons and positive ions are adopted in the simulation. The calculated results are consistent with the experimental results. Gordiets et al. [20] used N 2 -O 2 mixture equivalent air to reduce the type of chemical reaction, and studied the role of excited neutral molecules and atoms in the process of gas discharge. Nahorny et al. [21] considered 430 chemical reactions in the research of air discharge, which is currently the research that considers the most types of chemical reactions, and has laid a good foundation for subsequent scientific researchers. Pacheshnyi and Starikovskii [22] selected 10 main particles based on Nahorny, and proposed a two-dimensional physical model of streamer discharge for N 2 :O 2 (9:1) mixed gas at a pressure of 760 Torr. The model is used to analyze the electric field distribution, charged particle concentration, reaction coefficient and other discharge parameters in detail, and it has good consistency with the experimental data. However, in order to simplify the calculation, the above literatures often ignore the average energy transfer equation of electrons, but the average energy of electrons directly determines the reaction rate and energy transfer rate of electrons when they collide and scatter with other particles. Therefore, adding consideration to the average energy of electrons in the mathematical model is helpful to understand the collision and scattering of electrons with other particles during gas discharge, and is of great significance to reveal the microphysical mechanism of gas discharge. Lu et al examined the available results, generalized the physical mechanisms of the guided ionization waves, and described an atmospheric-pressure discharge featuring a plasma bullet behavior as an HSED (high seed electron density) discharge [23]. Throughout the current research results, the external conditions are mainly normal temperature and pressure. Numerical models for studying the characteristics of air discharge at low temperature and sub-atmospheric pressure are still lacking. Therefore, the establishment of this model is of great significance for studying the electron temperature, electron density, streamer radius and other microscopic parameters during air discharge at low temperature and sub-atmospheric pressure.
The selection of chemical reaction rate and transport coefficient is an important issue in air discharge theory, which is often obtained by determining the electron energy distribution. Since air discharge is a plasma discharge at low temperature and non-equilibrium state, its electron energy distribution is different from the Maxwellian distribution and the Druyvesteyn distribution. Therefore, it is necessary to calculate the accurate electron energy distribution function according to the Boltzmann equation to obtain the plasma chemical reaction rate and transport coefficient. BOLSIG+ is an open access electron energy distribution function calculation software [24], which has been widely recognized internationally. In this paper, the electron energy distribution function, as well as the electron mobility, diffusion coefficient and ionization reaction coefficients, are calculated by BOLSIG+. In order to study the characteristics of low-temperature and sub-atmospheric air discharge more accurately, this paper adopts the COMSOL Multiphysics modified model with artificial stabilization items. Due to limited space, the revised model will be introduced in detail in following paper, and its reliability and accuracy will be fully described. VOLUME 9, 2021

II. MODEL FORMULATION A. GOVERNING EQUATION
The electron continuity equation and the average electron energy density equation are [25], [26]: where n e is the electron density (cm −3 ); e is the electron flux (cm −2 s −1 ); S e is the electron yield (cm −3 s −1 ), which mainly involves ionization, excitation, attachment, recombination and photoionization; D e is the electron diffusion coefficient (cm 2 s −1 ); µ e is the electron mobility (cm 2 V −1 s −1 ); ϕ is the electric potential (V); E is the electric field (V/cm);ε is the average electron energy (V); ε e i is the energy loss of a single electron in a collision during the collision reaction (V); ε is the electron energy flux (Vcm −2 s −1 ); e is the elementary charge, e = 1.602176565 × 10 −19 (C); S i is the rate of collision reaction i (cm 3 s −1 ); D ε is the diffusion coefficient of electron energy (cm 2 s −1 ); µ ε is the mobility of electron energy (cm 2 V −1 s −1 ).
The governing equations for ions and neutral particles are [27]: where n + k , n − k , n n is positive ion density, negative ion density and neutral particle density (cm −3 ) respectively; S ion , S ph , S att , S ei rec , S ii rec is the ionizing reaction rate, photoionization reaction rate, attachment reaction rate, ion-electron recombination reaction rate, ionion recombination reaction rate respectively, the unit is cm −3 s −1 .
In this paper, the widely used Helmholtz equation is used to calculate the photoionization rate. The accuracy of this equation has been widely verified [28], the equation is: where pO 2 is the partial pressure of oxygen, and the value is 150 Torr under standard pressure. Parameters such as λ j and A j are consistent with the data in [1]. Equation (9) is the second-order scattering boundary condition developed in [1]. It has been verified that the boundary condition can make the calculation more accurate and efficient. Please refer to the reference [1] for detailed explanation, which will not be repeated here.

n·(∇S
The Poisson equation is [27]: where ε is the dielectric constant of air.

B. BOUNDARY CONDITIONS AND MESH STATISTICS
The model studied in this paper is a two-dimensional axisymmetric configuration based on the COMSOL Multiphysics 5.5 R package. All the calculations included in the model are run on a workstation with dual 2.9 GHz Intel(R) Xeon processors in 64-bit, 256GB RAM, and Windows 7 operating system. The temperature involved in the study is 223K, and the gas pressure is 170 Torr, which is in line with the external environmental conditions of the aircraft (such as an airplane) in the sky. The applied voltage is a unipolar square wave pulse with a peak value of 13kV, a duty cycle of 50%-75%, and a frequency of 25kHz-50kHz. The geometric figures involved in the model are shown in Figure.1. In order to trigger streamer discharge faster, the initial condition set in the calculation area is: where n 0 = 10 14 cm −3 , σ r = 0.1mm, σ z = 0.05mm; n e (r,z), n + i (r,z), n − j (r,z) respectively represents electron density, positive ion density, and negative ion density; z 0 = 0.975cm; the subscript i in n + i (r,z) represents N + 2 , O + 2 , N + 4 , O + 4 , N 2 O + 2 and the subscript j of n − j (r,z) represents O − 2 , O − . At t = 0ns, the electron density was significantly higher than the negative ion density, namely: n e n − j [29]. It can be seen from Figure 1 that B1 is the axis of symmetry, B2 is the anode, B5 is the ground electrode, and B3-B4 are the open boundary. The zero-charge boundary condition [30] is used, namely: The electron density flux on the electrode is expressed as [31]: where n is the boundary normal vector; γ e is the proportion of electrons reflected by the electrode surface; a e is 1 when the electron flux points to the electrode, otherwise it is 0; V e,th is the electron thermal velocity [29]: where T e is the electron temperature; m e is the electron mass; k B is the Boltzmann constant. The electron energy density flux on the electrode is expressed as [32]: The ion flux on the electrode is expressed as [30]: where i is the ion density flux; n i is the ion density; m i is the ion mass; T is the external temperature; γ i is the secondary electron emission coefficient when the ion bombards the electrode. The normal flux of all particles in the open boundary B3-B4 is zero. Figure.2 shows the accuracy of the meshing of this model through some meshes. In numerical calculations, accurate meshing is the guarantee of model convergence.
The mesh parameters are shown in Table 1. The quality of mesh cells is an important consideration when building the model. It can be seen from Table 1 that the minimum cell mass in the needle-plate electrode structure is 0.7821, and the average cell mass is 0.9733. In order to better converge the model, adaptive mesh is added to the basic mesh in this paper. Adaptive mesh is undoubtedly a very economical and effective method to obtain more accurate results for transient processes with rapid spatial changes such as streamer discharge. However, the benefit of this method largely depends on the basic mesh. Since the problem must be solved on the basic mesh for a period, the basic mesh does not need to be too refined under the premise of ensuring the quality of the mesh.

C. PLASMA CHEMICAL REACTION PROCESS
In order to better interpret the microscopic physical process, this paper considers the plasma chemical reaction. There are 12 main particles participating in the plasma chemical reaction, which are electron, Through the analysis of research results at home and abroad, this paper finally selects 39 chemical reactions to describe the air streamer discharge process. The reaction types are divided into the following seven categories: 1) Ionization reaction: 6) Recombination reaction and neutral particle reaction: In the above reaction, the rate of the chemical reaction involving electrons in R1-R11 is a function of the reduced electric field E/N, where N is the density of neutral particles. Based on the reaction cross-section provided by [33], this paper uses BOLSIG+ to obtain the corresponding reaction rate. O 2 (a), O 2 (b), N 2 (A), N 2 (B), N 2 (C) are the simplified symbols of excited state particles O 2 (a 1 g ), O 2 (b 1 + g ), N 2 (A 3 + u ), N 2 (B 3 g ), N 2 (C 3 u ), respectively. k 12 -k 31 are the chemical reaction rates [34], including two-body reaction and three-body reaction, and the units are mol/(s· m 3 ) and mol/(s· m 6 ).

A. DISCHARGE CHARACTERISTICS UNDER UNIPOLAR SQUARE PULSE VOLTAGE
In this section, the characteristics of low temperature sub atmospheric pressure air discharge under unipolar square wave pulse voltage are studied. The voltage amplitude is 13kv, the frequency is 25kHz-50kHz, and the duty cycle is 50%-75%. By changing the duty cycle and frequency, the effects on the average electron density, average electron temperature and average ion density in the whole calculation domain were studied.
The first is the relationship between electron density and electric potential. When the voltage is very low, the electron density in the plasma region is at a low level. When the voltage increases gradually, the electron density also increases gradually. As shown in Figure 3, with the time evolution of square wave voltage, the average electron density changes in the form of triangle wave. When the frequency remains unchanged and the duty cycle increases from 50% to 75%, the electron density increases. The reason for the increase in the electron density is that the increase in the duty cycle leads to an increase in the duration of the high level, which promotes the ionization reactions R1 and R2, thereby increasing the electron density. In addition to the duty cycle, the frequency of the pulse voltage also affects the change of electron density. When the duty cycle remains unchanged and the frequency increases from 25kHz to 50kHz, the electron density decreases. The main reason is that when the number of periods is the same, the duration of high-level decreases with the increase of frequency, which weakens the ionization reactions R1 and R2, and then decreases the electron density. In addition, when the number of periods is the same, the increase in frequency will also shorten the duration of the low level, which facilitates the accumulation of charge, thus increasing the pre-ionization level. With the increase of pre-ionization level, the breakdown field decreases gradually, and the discharge phenomenon is more likely to occur. Similarly, with the increase of pre-ionization level, the space charge density increases and the space charge distribution becomes more uniform, which reduces the probability of bifurcation phenomenon.
The second is the relationship between electron temperature and duty cycle and frequency. As shown in Figure 4, the electron temperature is different from the electron density. When the frequency remains unchanged and the duty cycle increases from 50% to 75%, the peak average electron temperature remains almost unchanged at 70eV. When the duty cycle remains unchanged and the frequency increases from 25kHz to 50kHz, the average electron temperature shows an increasing trend, but the duration of the high electron temperature becomes shorter. The main reason is that as the frequency increases, the duration of the high-level decreases and the electron density decreases, resulting in a thicker plasma sheath and an increase in the average electron temperature. However, due to the short duration of high level, the duration of high electron temperature is also shortened.
Finally, the relationship between positive ion density and duty cycle and frequency. As shown in Figure 5, the positive ion density is basically the same as the electron density. When the frequency remains unchanged and the duty cycle increases from 50% to 75%, the density changes of the main positive ions such as N are basically the same, which will not be repeated in this paper. The main reason for the increased density of N + 2 and O + 2 is the enhancement of ionization reaction R1 and R2. When the voltage is 0V, the O + 2 density in R11 is low, which leads to the weakening of R11 reaction, and then the density of N 2 O + 2 density is low. As the voltage gradually increases, the ionization reaction R2 gradually strengthens, leading to the strengthening of the R11 and the increase of the density of N 2 O + 2 . However, since the chemical reactions related to N 2 O + 2 are only R10 and R11, the density of N 2 O + 2 is five orders of magnitude lower than that of other positive ions. When duty ratio remains unchanged and frequency increases from 25kHz to 50kHz, the density relationship between N + 2 and N + 4 is different. According to the reaction rate of R12-R15, the inconsistency between the production rate and the consumption rate leads to N + 2 density being an order of magnitude higher than N + 4 . As mentioned above, the increase of frequency is conducive to the accumulation of charge. Besides the positive ions mentioned above, there is O − 2 , and  the increase of charge number indirectly affects the plasma chemical reaction rate. With the time evolution of voltage, N 2 O + 2 has a square wave pattern, but the high-level duration becomes shorter due to the increase in frequency, and the duration of high density N 2 O + 2 is also shortened. In order to further analyze the discharge characteristics under unipolar square pulse voltage with a frequency of 25kHz and a duty cycle of 50%, this paper selects a period for discussion and divides it into five time points: a, b, c, d, and e. As shown in Figure 6, the five time points are divided into four stages, which are rising stage, high-level stage, falling stage, and low-level stage.
Firstly, in the voltage rising stage, the sharp rise of the voltage on the electrode leads to the instantaneous increase of the electric field strength, which is prone to discharge. If the y-axis is positive, the electric field in the space is negative, as shown in the lower part of Figure 7. After the transition to the high-level stage, the electric potential remains basically unchanged, and the electric field distribution in the space VOLUME 9, 2021 FIGURE 7. Space electric field distribution: the upper part shows that the plasma potential is greater than the electrode potential; the lower part shows that the plasma potential is less than the electrode potential. is like that in the DC stage. It is worth mentioning that the plasma in space has a certain potential. When the voltage on the electrode falls sharply, the electric field changes dramatically, which is also a stage prone to discharge. With the decrease of the electrode voltage, when the electrode potential is less than the plasma potential, the electric field in the space changes from negative to positive, as shown in the upper part of Figure 7. In the last low-level stage, the charged particles gradually disappear through diffusion motion and recombination reaction, and the insulation performance of air gradually recovers.
For the high frequency discharge in this paper, when the time goes to millisecond level, the air has been in a conducting state. In Figure 8, the spatial structure of streamer in space is shown by five plasma parameters: electric field, potential, electron density, N + 2 density and O − 2 density.

B. THE INFLUENCE OF PHOTOIONIZATION ON STREAMER DISCHARGE CHARACTERISTICS
In the second section, the calculation of photoionization is briefly introduced. In this section, the role of photoionization  is discussed. As mentioned in [1], although the reaction rate of photoionization is very small compared with the ionization rate, photoionization provides seed electrons for discharge and promotes the development of discharge. As shown in Figure 9, in the simulation of double-headed streamer discharge in [1], it is found that photoionization is mainly concentrated in the region of streamer head, while the photoionization rate of positive streamer head is slightly higher than that of negative streamer head. According to [1], it can see that the development of positive streamer depends on photoionization. When photoionization cannot produce enough secondary electrons, positive streamer cannot develop. For the negative streamer, as shown in Figure 10, the planar front velocity of the negative streamer has little to do with whether photoionization is considered, that is, photoionization reaction is not important for the negative streamer, which is consistent with the conclusion of reference [35]. But even so, the photoionization reaction cannot be ignored in the calculation of air streamer discharge.
The random occurrence of photoionization leads to the random distribution of photoelectrons, and then the secondary electron avalanche is also random, so the bifurcation phenomenon is easy to occur. In unipolar square pulse discharge, metastable particles tend to remain in space because of their long life. As mentioned above, with the increase of frequency, the number of charges in space increases, the background charge becomes more uniform, and the probability of bifurcation phenomenon decreases. When the pre-ionization level of background charge reaches 10 7 cm −3 , the photoionization effect is no longer obvious.

C. EFFECT OF TEMPERATURE ON AIR DISCHARGE AT LOW TEMPERATURE SUB ATMOSPHERIC PRESSURE
Gas pressure and temperature are the two main factors that affect the discharge of air at low temperature and subatmospheric pressure. In the previous paper [1], the corresponding conclusions about the influence of gas pressure on discharge have been given. This paper will discuss the influence of temperature and gradually improve the basic theory of air discharge at low temperature and sub-atmospheric pressure.
First, analyze the effect of temperature on the rate of 39 plasma chemical reactions. In order to better explain the influence of temperature on plasma reaction rate, this paper only takes 1cm needle-plate discharge under DC voltage as an example. The interaction between various plasma chemical reactions and the sequence of their occurrence are reflected only by the change of the density of charged species containing O element with time.
It can be seen from Figure 11 that the earliest occurrence is the electron collision reaction, which includes the ionization reaction of R1 and R2, the elastic collision reaction of R3 and R4, and the excitation reaction of R5-R9. In this paper, f is used to represent the number of collisions experienced by the particle per unit time, that is, the collision frequency. The collision frequency is proportional to the collision reaction rate. The more the number of collisions, the greater the collision frequency and the greater the collision reaction rate.
where f i , f ex , f el are ionization frequency, excitation frequency, and elastic collision frequency, respectively. The probability of interaction between groups in the gas phase depends on the cross-section of the interaction, and each interaction needs to be characterized by a specific crosssection, such as σ i , σ ex , σ el denote ionization cross-section, excitation cross-section, and elastic collision cross-section respectively. N is the density of neutral particles.ῡ is the average value of the random movement speed of the molecule. According to the ideal gas equation, formula (20) is derived as: The neutral particle density N is proportional to the gas pressure p and inversely proportional to the temperature T. Under the premise of constant gas pressure, when the temperature decreases, the density of neutral particles increases. The collision cross section itself is a quantity independent of the neutral particle density. Therefore, as the temperature decreases, the ionization frequency, excitation frequency, and elastic collision frequency all increase, and then the ionization rate, excitation rate, and elastic collision rate increase.
What happens next is the charge transfer reaction between particles. The reaction rates k 12 -k 20 are constants, and temperature has almost no effect on this type of reaction rate. When the time exceeds 10 −8 s, the electric field inside the plasma channel is very low due to the shielding effect of the space charge [1], and the particle movement speed is slow. Currently, the recombination reaction between electrons and positive ions including R21, R22, R30, and R31 is likely to occur. Although the recombination reaction rate between electrons and positive ions has little relationship with temperature, the decrease in temperature reduces the conductivity in the plasma channel. Analysis shows that ''cluster ions'' such as N + 4 and O + 4 play a major role. In the recombination reaction, the recombination rate of ''cluster ions'' is much higher than that of normal N + 2 and O + 2 . When the temperature decreases, the ''cluster ions'' are difficult to decompose, which leads to the increase of electron loss and the decrease of conductivity. When the time reaches 10 −7 s, due to the presence of electronegative gas oxygen in the air, the electron attachment reactions R10 and R11 occur successively. In this paper, R10 is named ''three-body attachment reaction'' and R11 is named ''two-body attachment reaction''. Like the collision frequency, the attachment frequency increases with the decrease of temperature. Both three-body attachment and two-body attachment reaction rates increase with the increase of the attachment frequency. The difference is that the three-body attachment reaction often occurs in the main streamer, while the two-body attachment reaction often occurs in the subsequent secondary streamer. According to the reaction mechanism, the three-body attachment reaction rate is proportional to the square of the neutral particle density N 2 , and the two-body attachment reaction rate is proportional to the neutral particle density N. When the temperature decreases, the neutral particle density N increases, and the three-body attachment reaction rate and the two-body attachment reaction rate both increases, but the increase in the three-body attachment reaction rate is significantly higher than the two-body attachment reaction rate. With the increase of attachment reaction rate, the density of negative ion O − 2 gradually increases, and the recombination reaction between positive ion and negative ion including R23-R27 takes place successively. However, temperature has little influence on this kind of recombination reaction rate.
Secondly, in addition to the photoionization reaction during the air discharge, the secondary electron emission caused by the positive ion bombarding the cathode metal also provides some seed electrons for the discharge. Different metal materials have different secondary electron emission coefficients γ , and the same metal material has different secondary electron emission coefficients at different temperatures.
As shown in Figure 12, the secondary electron emission coefficient increases exponentially with the increase of temperature [36]. At low temperature, the secondary electron emission coefficient is small, and the number of electrons emitted by positive ions bombarding cathode metal decreases, which is consistent with the conclusion that the conductivity in the plasma channel is low at low temperature.
The change of the secondary electron emission coefficient also affects the distribution of current and axial electric field in the plasma channel. Taking the needle-plate streamer discharge under 10kV DC voltage as an example, when the secondary electron emission coefficient is smaller, it means that the positive ion bombards the cathode to emit fewer electrons, so it needs a longer propagation time to generate enough electrons. As shown in Figure 13, a longer propagation time will increase the electric field gradient and peak electric field. In Figure 14, the larger the secondary electron emission coefficient, the shorter the propagation time, so the  electrons required for streamer propagation will accumulate faster. Temperature indirectly affects the distribution of current and axial electric field in the plasma channel by changing the secondary electron emission coefficient.
Finally, according to the conclusions obtained above, the conductivity in the plasma channel decreases at low temperatures, and the conduction current also decreases accordingly. The formation of conduction current is due to the directional movement of carriers, which mainly include electrons, positive ions, and negative ions. Temperature also affects electron mobility and ion mobility. This paper uses COMSOL Multiphysics to construct a transient model based on the migration tube method to measure ion mobility, calculates the change curve of positive and negative ions in the temperature range of −50 • −50 • (223K-323K), and fits the corresponding formula. As shown in Figure 15, with the increase of temperature, the mobility of positive and negative ions declines, while that of positive ions decline slowly, while that of negative ions decline relatively quickly. When the temperature is between 223K and 318K, the mobility of negative ions is higher than that of positive ions, and after 318K, the mobility of negative ions is gradually lower than that of positive ions. The reason may be that with the change of temperature, the type and proportion of the product ions in different types of plasma reactions change. According to the curve in Figure 15, the equation (21) where µ + and µ − represent positive and negative ion mobility (cm 2 V −1 s −1 ) respectively, and T is the temperature ( • C). It can also be seen from the fitting formula that the ion mobility is negatively correlated with temperature. In order to facilitate analysis, this paper divides ions into light ions and heavy ions. Compared with heavy ions, light ions are small in weight and size, so the mobility of light ions is obviously greater than that of heavy particles. The result of simulation calculation is the weighted average of light and heavy ions.
where µ L and µ H represent the mobility of light ions and heavy ions, and α is the proportion of heavy ions in the air. According to the reference [37], when in an ideal state, the heavy ion mobility of positive and negative polarities in the air is 1.55 × 10 −2 cm 2 V −1 s −1 and 2.0 × 10 −2 cm 2 V −1 s −1 ; The mobilities of light ions are 1.55 cm 2 V −1 s −1 and 2.0 cm 2 V −1 s −1 respectively. If the effect of temperature on the ion mobility of a single component is ignored, the value of ion mobility in Figure 15 is substituted into equation (23)  mentioned above belong to heavy ions. According to the measurement principle of the migration tube method, when the temperature increases, the light ions migrate under the action of the electric field and diffuse around. Because the light ions are very small in mass and volume, the diffusion caused by thermal motion is obvious, so the number of light ions that reach the measuring point is less. Due to the large mass and volume of heavy ions, the diffusion caused by thermal motion is relatively weak. The number of heavy ions near the measurement point does not decrease significantly, so the proportion α of heavy ions near the measurement point increases, resulting in ion mobility µ decreases. The main reason why the negative ion mobility is greater than the positive ion mobility is that the proportion of nitrogen in the air is 78%, and the proportion of N + 2 and N + 4 generated by nitrogen participating in the ionization reaction is also large, resulting in the volume and mass of positive ions larger than that of negative ions, and then the mobility becomes smaller.
The influence of temperature on the electron mobility mainly reflects the change of the reduced electric field E/N caused by the change of temperature, where the unit of the reduced electric field strength E/N is Td (1Td = 10 −17 Vcm 2 ). According to the formula obtained in the reference [36], when 39.9 Td ≤E/N ≤206.1 Td, the calculation formula of electron mobility is: When 206.1 Td ≤ E/N ≤ 1000 Td, the calculation formula of electron mobility is: According to the curve fitted in Figure 16, with the increase of reduced electric field E/N, the electron mobility decreases. According to equation (20), it can also be found that the increase of temperature causes the neutral particle density to decrease, and the reduced field E/N increases. Temperature and reduced field strength are positively correlated. When in low temperature conditions, the reduced field strength E/N is small, and the electron mobility is large currently. According to the above analysis, when the temperature is low, both the electron mobility and the ion mobility are at a high value, but the conductivity in the channel is low. The main reason for this phenomenon is that although the carrier migration speed is faster at low temperature, the rate of electron attachment and positive and negative ion recombination reaction also increases, resulting in the decrease of carrier number, conduction current and conductivity.
In addition to affecting the mobility of carriers, the streamer propagation velocity is also affected by temperature. This paper compares the electron density distribution and the electric field distribution under the temperature of 223K and 293K to illustrate the influence of temperature on the streamer propagation velocity. Figure 17 shows the electron density distribution at the same time under different temperatures. Figure 18 shows the electric field distribution at the same time under different temperatures. When the temperature is 223K, the streamer propagation speed becomes slower, and the propagation distance is shorter than 293K in the same time. The main reason is that when the temperature is low, the reduced electric field E/N is smaller, and the  streamer head field also decreases, which leads to the slow ionization reaction of streamer head and the difficulty of streamer propagation.
In order to further analyze the difference of temperature and gas pressure in affecting the velocity of streamer, this paper discusses them separately. One group keeps the gas pressure constant and changes the temperature; the other  group keeps the temperature constant and changes the gas pressure, and observe their respective changes. As shown in Figure 19, when the same time is selected, the higher the temperature, the greater the deviation of the streamer head position, that is, the faster the streamer spreads, which is consistent with the conclusion that the streamer propagates slowly at low temperatures. Compared with p = 170Torr, T = 323K in Figure 19 and p = 144Torr, T = 273K in  Figure 20, it is found that when the reduced field E/N is the same, the higher the temperature, the faster the streamer propagation speed. The main reason is the effect of temperature on the conductivity in the plasma channel. When the temperature is higher, the conductivity in the streamer channel is higher, the voltage drop in the plasma channel is smaller, and the electric field at the streamer head is larger, which promotes the ionization reaction and accelerates the streamer propagation.

IV. CONCLUSON
This paper is based on the modified COMSOL Multiphysics model with artificial stability, and applies the boundary condition developed in the previous [1] to the calculation of photoionization reaction rate. The main processes involved in air discharge were simulated by 39 plasma chemical reactions composed of 12 kinds of particles, and an improved model of air discharge under needle-plate electrode structure was constructed. The characteristics of low-temperature sub-atmospheric pressure air discharge under unipolar square wave pulse voltage with external application amplitude of 13kV, frequency of 25kHz-50kHz and duty ratio of 50%-75% are discussed. The main conclusions are as follows: 1) The average electron density during the discharge process is positively correlated with the duty cycle of the pulse voltage, and negatively correlated with the pulse voltage frequency; the average electron temperature is positively correlated with the pulse voltage frequency, but the increase in frequency will shorten the duration of the high electron temperature. The average electron temperature is hardly affected by the pulse voltage duty cycle; the average density of N + 2 , N + 4 , O + 2 , O + 4 , N 2 O + 2 are positively correlated with the pulse voltage duty cycle. When the duty cycle remains unchanged and the frequency increases, the average density of N + 2 and O + 2 differs from N + 4 and O + 4 by an order of magnitude. With the time evolution of the square wave voltage, N 2 O + 2 changes in a square wave form, but the increase in frequency also shortens the duration of the N 2 O + 2 high-density area. The increase of the pulse voltage frequency will also increase the number of space charges, make the background charge more uniform, and prevent bifurcation. 2) When the electrode potential is greater than the plasma potential, the space field is in the negative y-axis direction, otherwise it is in the positive y-axis direction.
It is easy to discharge in the two stages of voltage rise and fall. The high-level stage is similar to the DC discharge, while the air gradually recovers its insulation performance in the low-level stage.
3) The photoionization reaction is mainly concentrated in the head of the streamer, and the photoionization reaction is not important for the negative streamer. Even so, the photoionization reaction still needs to be considered in the streamer discharge calculation. When the spatial pre-ionization level reaches 10 7 cm −3 , the photoionization effect is no longer obvious. 4) The influence of low temperature on air discharge is mainly manifested in the following aspects: the decrease of temperature leads to the increase of collision reaction rate including elastic collision, ionization, and excitation. The rate of attachment reaction also increased, including ''three-body attachment'' and ''two-body attachment''. The results show that the recombination reaction rate is almost not affected by temperature; the low temperature leads to the decrease of the secondary electron emission coefficient, which means that the discharge propagation time is prolonged, and the longer the propagation time will enhance the electric field gradient and the peak value of electric field; the decrease of temperature leads to the increase of carrier migration rate including electron, positive and negative ions, and the mobility of negative ions is significantly higher than that of positive ions. The decrease in temperature will also decrease the streamer propagation velocity. When the reduced field strength E/N is the same, the lower the temperature, the slower the streamer propagation speed. Due to limited space, this paper only discusses the discharge characteristics under unipolar square wave pulse voltage. The research on the change rule under bipolar square wave pulse voltage and other electrode structures will be reflected in following paper. ZHONGHUA ZHANG was born in Henan, China, in 1998. He received the B.S. degree in electrical engineering from the Harbin University of Science and Technology, Harbin, China, in 2020, where he is currently pursuing the Ph.D. degree in high voltage and insulation. His research interest includes epoxy impregnated paper insulation in HVDC bushing.
KAILUN YANG was born in Jilin, China, in 1996. He received the B.S. degree in electrical engineering from the Changsha University of Science and Technology, Changsha, China, in 2018. He is currently pursuing the Ph.D. degree in high voltage and insulation with the Harbin University of Science and Technology. His research interests include high thermal conductivity insulating materials and dielectric failure mechanisms of polymer insulating materials. VOLUME 9, 2021