Model of a non-transferred arc cascaded-anode plasma torch: the two-temperature formulation

Under the multi-temperature equilibrium assumption, the thermodynamic state of the system is fully described by the partition functions of each species. They comprise the translational, internal and reactional degrees of freedom of the entities in the plasma-forming gas and can be written as follows:

$$\begin{equation}{Q_i} = Q_i^{{\text{tr}}}\left( {T_i^{\,{\text{tr}}}} \right)Q_i^{{\text{int}}}\left( {T_i^{\,{\text{int}}}} \right){\text{exp}}\left( { - \frac{{e_s^{\,f}}}{{{k_{\text{B}}}T_i^r}}} \right)\end{equation} \tag{ 6 }$$

where $Q_i^{{\text{tr}}}$, $Q_i^{{\text{int}}}$, $e_s^{\,f}$ are respectively the translational, internal partition functions and formation energy of species i. The temperatures $T_i^{\,{\text{tr}}}$, $T_i^{\,{\text{int}}}$ and $T_i^r$ are the translation temperature, temperature linked to the internal degree of freedom of heavy species and reaction temperature, respectively (Capitelli et al 2012).

In the case of a pure argon plasma, mainly the ionization reactions by electron impact have to be considered. In such conditions, it is assumed that the internal degrees of freedom of atoms and ions will be ruled by the electron temperature and ionization reactions as well, so that $T_e^{{\text{tr}}} = {T_e}$, $T_{i \ne e}^{{\text{tr}}} = {T_h}$, $T_i^{\,{\text{int}}} = T_i^r = {T_e}$.

Applying the maximization of the entropy to such a system gives rise to the 2-T Saha equation (Giordano 1998, Capitelli and Giordano 2002), often referred to as the van de Sanden formulation (Van De Sanden et al 1989). It follows that:

$$\begin{align}\begin{array}{l} \frac{{{n_e}{n_{{h^{\left( {z + 1} \right) + }}}}}}{{{n_{{h^{z + }}}}}} = 2{\left( {\frac{{2\pi {m_e}{k_{\text{B}}}{T_e}}}{{h_p^{{\kern 1pt} 2}}}} \right)^{3/2}}\frac{{Q_{{h^{\left( {z + 1} \right) + }}}^{{\text{int}}}\left( {{T_e}} \right)}}{{Q_{{h^{z + }}}^{int}\left( {{T_e}} \right)}}{\text{exp}}\\ \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \times \left( { - \frac{{{E_{{h^{\left( {z + 1} \right) + }}}} - \delta {E_{{h^{\left( {z + 1} \right) + }}}}}}{{{k_{\text{B}}}{T_e}}}} \right) \end{array}\end{align} \tag{ 7 }$$

where ${n_e}$ and $n_{{h^{z + }}}$ are respectively the electron number density and number desnity of heavy species with a charge number $z$, and ${E_{{h^{\left( {z + 1} \right) + }}}} = e_{{h^{\left( {z + 1} \right) + }}}^{\,f} - e_{{h^{{z^ + }}}}^{\,f}$. The same notation is used for the internal partition functions of the heavy species $Q_{{h^{z + }}}^{{\text{int}}}$.

The constants ${m_e}$, ${k_{\text{B}}}$ and $h_p$ are the electron mass, Boltzmann constant and Planck constant, respectively. $\delta {E_{{h^{\left( {z + 1} \right) + }}}} = \left( {z + 1} \right){e^{\,2}}/4\pi {\varepsilon _0}{\lambda _D}$ is the ionization energy lowering where ${\lambda _D}$ is the Debye length and ${\varepsilon _0}$ is the vacuum permittivity (Andre et al 1996).

The internal partition functions are calculated from the NIST data (Kramida et al 2018) and their excited states number taken into account is limited by the ionization energy lowering. The latter is obtained by considering the Debye length involving the electrons only as recommended in (Murphy 2000).

The Dalton's law of partial pressures (8) and electrical neutrality assumption (9) are necessary to obtain the plasma composition:

$$\begin{equation}p = {n_e}{k_{\text{B}}}{T_e} + \sum\limits_{N - 2}^{z = 0} n_{{h^{z + }}}^{ }{k_{\text{B}}}{T_h}\end{equation} \tag{ 8 }$$

$$\begin{equation}{n_e} = \mathop \sum \limits_{z = 0}^{N - 2} z\,n_{{h^{z + }}}\end{equation} \tag{ 9 }$$

where $p$ and ${n_i}$ are the total pressure and the number density of species $i$, respectively.

The argon plasma composition is calculated up to 35,000 K at atmospheric pressure considering the electron, Ar, Ar⁺, Ar²⁺, Ar³⁺ and Ar⁴⁺ for different values of the disequilibrium degree $\theta = {T_e}/{T_h}$. For computational purposes, the data were generated for $\theta$ from 1 up to 25 with a step of 0.1 for the reasons given in the section about the implementation of the 2T model.

Figure 2 displays the number densities as function of the electron temperature for different values of $\theta$. The results are in good agreement with those obtained in Colombo et al (2008).

Zoom In Zoom Out Reset image size

It has to be underlined that the widely used representation of figure 2, i.e. the electron temperature dependence of all number densities, is confusing because it is seen that the number densities of the major species (e, Ar, Ar+) increase as $\theta$ increases while the total pressure is kept constant. The number densities of the heavy species are T_h -dependent: in a T_h-representation their number densities are shifted to lower T_h temperatures as $\theta$ increases. Note that the 2-T equilibrium constants of ionization reactions (right-hand side of equation (7)) are almost independent of $\theta$ even though the translation partition function of the heavy species $Q_h^{{\text{tr}}}$ are T_h -dependent. Consequently, if the reference temperature is T_e , the partial pressure of electrons increases as $\theta$ increases at the expense of the partial pressure of the heavy species according to the total pressure conservation. At very low temperature, ${n_{{\text{Ar}}}}\sim \theta p/{k_{\text{B}}}{T_e}$. At high T_e and $\theta \gg 1$ the electron number density tends to ${n_e}\sim p/{k_{\text{B}}}{T_e}$ and becomes $\theta$-independent as seen in figure 2 for $\theta = 10$.

Figure 3 depicts the electron enthalpy with the second formulation of enthalpy, referred to as 2T(II) in the previous section. The insert displays the low temperature values in a log-scale that is very useful for computational purposes. Only the 2T(II) enthalpy is shown because the other formulation (2T(I)) will lead to inconsistent results as shown in the next sections. As already known, the main contribution is due to the ionization terms highlighted by the changes in the slopes. For a given T_e , the electron enthalpy decreases with an increase of $\theta$ due to the increase of the mass density conditioned by the increase of the number density of ions seen in figure 2.

Zoom In Zoom Out Reset image size

The last thermodynamic properties useful for the 2-T model is the specific heat derived from the specific enthalpy. Figure 4(a) displays a comparison of the specific heat with Colombo et al (2008) calculated as ${c_p} = {\text{d}}\left( {{h_e} + {h_h}} \right)/{\text{d}}{T_e}$; a very good agreement is found between both calculations. However, we have adopted a different definition, i.e. ${c_{pe}} = {\text{d}}{h_e}/{\text{d}}{T_e}$ and ${c_{ph}} = {\text{d}}{h_h}/{\text{d}}{T_h}$. Figure 4(b) shows the electron temperature dependence of ${c_{pe}}$, derived from the 2T(II) enthalpy formulation, for different values of $\theta$. The specific heat of heavy species in the 2T(II) formulation is almost constant, i.e. ${c_{ph}}\sim 5{k_{\text{B}}}/2{m_{{\text{Ar}}}}$.

Zoom In Zoom Out Reset image size

This section briefly displays the 2-T argon transport coefficients calculated for the model. The details of formulae used are not shown here but can be found for example in Colombo et al (2008) and references therein. A comparison of the equilibrium and 2-T calculations with previous published results is also carried out. The collision integrals ${\bar Q_{ij}}$ for species collisions in the argon plasma are calculated at the T_e temperature when a collision involves an electron and at T_h otherwise. The temperature dependence of the collision integrals (e-Ar, Ar-Ar, Ar-Ar^z+) are taken from Aubreton and Fauchais (1983). The collision integrals for charged species are calculated from the approach given by Mason et al (1967), Devoto (1973) considering the screened Coulomb potential involving the Debye length involving the electrons only.

Figure 5 presents the electrical conductivity for different $\theta$ values. A good agreement is found at equilibrium with Murphy (2000), Colombo et al (2008) and also for non-equilibrium conditions with Colombo et al (2008). Below roughly 12,000 K, it decreases when $\theta$ increases as shown in the insert. Actually, the plasma is dominated by the collisions between heavy species so that the collision integrals ${\bar Q_{ij}}\left( {{T_h} = {T_e}/\theta } \right)$ turn out to be higher as $\theta$ increases, hence the decrease of electrical conductivity. Above 12,000 K, the electrical conductivity follows the electron number density $\theta$-dependence of figure 2. Figure 6 depicts the thermal conductivities, including the translational thermal conductivity of electrons and heavy species, ${\lambda _e}$ and ${\lambda _h}$ respectively, and the reactive thermal conductivity ${\lambda _r}$. The viscosity is shown in figure 7. Again, a good agreement is found with (Murphy and Arundell 1994, Murphy 2000, Colombo et al 2008). The dependence of the coefficients ${\lambda _e}$ on $\theta$ is consistent with the electron number density. The coefficient ${\lambda _h}$ decreases with $\theta$ because of the collision integrals ${\bar Q_{ij}}\left( {{T_h} = {T_e}/\theta } \right)$. The viscosity has the same $\theta$- dependence as ${\lambda _h}$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

A good agreement is also found for ${\lambda _r}$ at equilibrium (Murphy and Arundell 1994). It has to be underlined that different approaches were published and discussed in the 2-T assumption, as for example in (Capitelli et al 2012) or (Santos et al 2019). It is beyond the scope of this paper to debate the effect of the different methods, which would require a close examination together with the methods of plasma composition calculation. Instead, ${\lambda _r}$ is calculated according to the method described in Bonnefoi and Fauchais (1983) and exhibits a decreasing evolution with $\theta$ as shown in figure 6(b).

The exchange coefficient is calculated with the equation (10) as in Trelles et al (2007), Freton et al (2012):

$$\begin{equation}{K_{{\text{exch}}}} = \sum\limits_{i \ne e}^{{ }N} 3{k_B}\frac{{{m_e}}}{{{m_i}}}{n_i}{v _{ei}}.\end{equation} \tag{ 10 }$$

${v _{ei}}$ is the elastic collision frequency between the electron and heavy species i and ${v _{ei}} = {n_i}{\bar v_e}{\bar Q_{ei}}$. ${\bar v_e}$ and ${\bar Q_{ei}}$ are the electron mean speed and average elastic cross section, respectively, calculated according to Freton et al (2012).

Figure 8 shows the exchange term ${E_{{\text{exch}}}} = {K_{{\text{exch}}}}\left( {{T_e} - {T_h}} \right)$ for different $\theta$ values and different scales. It is seen that ${E_{{\text{exch}}}}$ increases with $\theta$ like the electron number density.

Zoom In Zoom Out Reset image size

The energy conservation equations require the calculation of the divergence of the radiative flux which represents the net balance between the emission and absorption radiative processes. In the commonly-used approach of the net emission coefficient (NEC), the absorption coefficient, which encompasses all the continuum and lines contributions and includes the integration over all the spectrum wavelengths assuming the radiation isotropically emerges from an isothermal sphere of radius ${R_p}$ (Lowke 1974, Liebermann and Lowke 1976). For reasons linked to the computational cost (Trelles and Modirkhazeni 2014) and availability of ready-to-use NEC data as a function of temperature, the NEC approach is used in many models (Colombo et al 2011, Boselli et al 2013, Liang and Groll 2018, Modirkhazeni and Trelles 2018, Liang 2019, Liang and Trelles 2019). It has to be noted that the NEC coefficient cannot consider the absorption in the fringes of the arc or in the cold gas surrounding the arc, which would bring about a negative divergence of the radiative flux in these low temperature regions (Kovitya and Lowke 1985, Randrianandraina et al 2011). In other models (Trelles et al 2007, Baeva et al 2016, 2019, Baeva 2017), the radiation power is estimated by means of the Beulens' et al formula (1991) that gives a dependence on the electron number density.

The NEC given by Erraki (1999), Freton et al (2012) that considers the contributions of continuum and lines radiation of argon where the self-absorption of lines is taken into account in this paper. Figure 9 displays the NEC of argon at atmospheric pressure and at LTE. It shows that the absorption has to be considered when comparing ${R_p} = 0$ (Menart and Malik 2002) and ${R_p} = 5\,\;{\text{mm}}$ (Erraki 1999). Figure 9 also shows the continuum NEC from (Erraki 1999). The lines emission contribution (including the self-absorption lines), NEC (total)—NEC (continuum), is important especially at high temperatures.

Zoom In Zoom Out Reset image size

Note that Beulens et al (1991) give a valid approximation only up to 15,000 K for the total and continuum NEC as shown in figure 9 by the extrapolation of the Beulens' et al formula obtained with the plasma composition calculated up to 35,000 K. The approximation does not consider any absorption which does not hold at high temperature where radiation losses are the highest. However, this approximation is easy to handle for different $\theta$ since it depends on the electron number densities. However, it has to be underlined that the radiation losses are dominant for high current arcs leading to high temperatures and also to thermal equilibrium that is why one can assume the radiation losses are independent of $\theta$ at least in the NEC approximation.

Moreover, the radiation power is split between the energy equations of electron and heavy species assuming that the free-free and free-bound transitions leading to radiative emission correspond to electron energy loss mechanisms. The bound-bound transitions produce radiation losses from the heavy species (Freton et al 2012). The last statement can be questionable and suggests that the electron kinetic energy is weakly affected by inelastic collisions populating the excited states of heavy species. We suppose that only the tail of the electron distribution function may have departure from the Maxwellian distribution which weakly affect the mean electron kinetic energy electrons for ${T_e}$ temperatures below 35,000 K. Consequently, as also supposed by Freton et al (2012), the radiation losses due to the continuum transitions are ascribed to the electron energy equation while the line transitions are attributed to the heavy species equations.

The effect of the radiative properties (Beulens et al 1991, Erraki 1999) on modeling results will be shown in the following.

For the sake of clarity, the thermal boundary conditions are expressed in terms of temperature. In fact, in the model, the values imposed in the Dirichlet boundary conditions are the values of enthalpies, taken from the lookup tables of the plasma properties.

The only boundary where the electron temperature is explicitly imposed (Dirichlet boundary condition) is the inlet. The electron temperature at the inlet is assumed to be 1950 K. Below this temperature, the electron number density is too low resulting in much less than one electron per fluid cell, which makes the concepts of Maxwellian distribution and electron temperature inapplicable. An artificial increase in the electron number density at low temperatures, which is used in some 2T arc simulations (Trelles et al 2007, Freton et al 2012) to equalize the electron and heavy species temperature in the cold periphery, does not change the distributions of the arc plasma properties in our model while it produces unwanted fluctuations in the exchange term and predicted electron temperature. The heavy species temperature at the inlet is assumed to be 300 K. Thus, the disequilibrium degree $\theta = {T_e}/{T_h}$ at the inlet is 6.5.

The electron temperature on the other boundaries has a zero-flux boundary condition. Zero flux of electron temperature on the neutrode wall is assumed due to the low sensitivity of the model to the imposed boundary condition. Even ${T_e}$ of 300 K imposed on the neutrode wall results in noticeable changes of the electron temperature only in the boundary cells. This is because the simulated electron enthalpy balance in the periphery of the arc is dominated mostly by the Joule power, exchange term and, to a much lesser extent, by convection and dissipation. The zero-flux boundary condition for electron temperature on the electrode surface in this case is a simplification while in general the electron temperature on the electrode-plasma interface should be defined by a sheath model, which is omitted in the present study. Due to the lack of cathode sheath, the electron temperature in plasma is not coupled with the cathode temperature and does not affect the thermionic emission on the cathode surface.

The heavy species temperature on the surface of the neutrodes is assumed to be 500 K between the cathode tip and anode due to the heating of the neutrodes by the plasma radiation and 300 K upstream from the cathode tip where the injected gas is still cold. The model has little sensitivity to the exact value of the thermal boundary condition on the neutrode surface between the electrodes. An increase of the assumed neutrode temperature by 100 K in the model results in barely noticeable changes in other model indicators. Since the model does not include the copper neutrodes, the neutrode surface temperature of 500 K is merely an assumption needed to take into consideration the neutrode heating by plasma radiation. In the real SinplexPro plasma torch, a water-cooling circuit ensures a constant wall temperature.

The plasma-forming gas is injected through 24 small spots behind the cathode as shown in figure 10. It is intended to imitate the gas injection through the 24 holes of the gas injector ring of the actual torch. The gas injection angle is set at 25°, which is controlled through the momentum source term at the boundary cells. The injection through small holes increases the angular momentum and, thus, intensifies the flow swirling. The velocity at the boundary condition is controlled in order to maintain the gas flow rate, which in this study was set at 60 normal liters per minute (NLPM).

The voltage ${V_{{\text{recomputed}}}}$ applied on the external surface of the cathode as boundary condition is recomputed after each time step in the model in order to maintain the imposed electric current intensity through the plasma (Électricité de France R&D 2017, Zhukovskii et al 2020b). As a part of the voltage scaling algorithm the total Joule power in the model from the previous time step is computed and compared with the product of the imposed current intensity (I) and voltage (V) imposed at the previous time step. If the total Joule power is higher than the product I·V, the voltage is decreased proportionally, and vice versa. In order to avoid detrimental fluctuations at the beginning of the simulation, when all the parameters change a lot from one time step to the other, the increment of the recomputed voltage is limited by 10% per time step. The time-evolution of the arc voltage computed this way is shown in figure 11. A similar method with the recalculation of the imposed voltage was employed in Baeva and Uhrlandt (2011).

Zoom In Zoom Out Reset image size

For this study, the open-source software code_saturne was adapted for the 2-T modeling of electric arc. For that an extra user scalar was added with its own diffusivity and source terms. The user-added scalar was solved by the same solver as the default thermal variable. The first step in development of the 2-T model was similar to the one described in Freton et al (2012): default and user-added thermal equations are duplicated as LTE, in order to verify that these two equations are numerically identical. Both thermal scalars, default and user-defined, are the same in the way their initialization, boundary conditions and source terms can be treated in the user subroutines.

The thermal variable in code_saturne by default is the enthalpy because it is convenient for the computation of the energy balance: the thermal source terms and thermal flux, which are usually expressed in terms of energy, can be used as they are. In addition, the use of enthalpy as the thermal variable allows to simulate the phase transition in the cells of the electrodes included into the computational domain. During the phase transition from solid to liquid the energy accumulated in the heated cells of the electrode increases, while the temperature remains the same. Therefore, the temperature cannot be used to capture this process. In the 2-T model, the thermal equations solved inside the electrodes are identical for the default thermal variable and the used-added one.

However, when enthalpy is used as the thermal variable, a problem arises: all the plasma properties are defined as functions of temperature and not of enthalpy. These properties are usually stored in a lookup table, which associates the precalculated values of plasma properties to each value of some certain set of temperature values. Therefore, in order to calculate the plasma properties for a given cell, the temperature is first needed. In the LTE model, the translation of the enthalpy to temperature is done using a simple lookup table with one column for temperature and one column for enthalpy. Meanwhile, the 2-T model solves the enthalpy of electrons ${h_e}$ and enthalpy of heavy species ${h_h}$, which must be translated into the temperature of electrons and temperature of heavy species. The solution of this problem is not as straightforward as in the LTE case, since both enthalpies are functions of both temperatures. In order to solve this issue a set of disequilibrium degrees $\theta = {T_e}/{T_h}$ is selected prior the running of simulation. The set of $\left\{ {{\theta _i}} \right\}$ must be fine enough to guarantee the smooth transitions of all the properties from one ${\theta _i}$ to another. For each ${\theta _i}$ a lookup table with plasma properties as functions of electron temperature is built. Thus, the number of elements in the set $\left\{ {{\theta _i}} \right\}$ is equal to the number of lookup tables used in the model. The interpolation between the values of $\left\{ {{\theta _i}} \right\}$ is not employed in order to optimize the computational cost of the plasma properties recalculation at each time step for each cell.

In order to obtain the electron and heavy species temperatures of a given cell for calculation of the plasma properties the following procedure is performed. The enthalpies ${h_e}$ and ${h_h}$ are taken either from the previous time step or from the initialization values for the first time step at the start of the model run. For each value ${\theta _i}$ from the set $\left\{ {{\theta _i}} \right\}$ a simple translation ${h_e} \to T_e^{\,i}$ and ${h_h} \to T_h^{\,i}$ is done. Following that, the initial set $\left\{ {{\theta _i}} \right\}$ is compared with the obtained set of ratios $\left\{ {T_e^{\,i}/T_h^{\,i}} \right\}$. Some obtained ratios of temperatures are higher than the initial ${\theta _i}$, some are lower. The closest obtained ratio $T_e^{\,i}/T_h^{\,i}$ to the initial ${\theta _i}$ is selected as the best fit. The best fit ${\theta _{{\text{selected}}}}$ is assigned to the given cell and, using the obtained $T_e^{\,{\text{selected}}}$, all the other plasma properties are calculated for the given cell.

Each simulation instance starts with the assignment of the initial values of all the solved variables. After the variable initialization, the following key stages are repeated during each time step in the developed 2-T model in code_caturne:

• (a)  
  Translation of enthalpies to temperatures: $\left\{ {{h_e},{h_h}} \right\} \to \left\{ {{\theta _{{\text{selected}}}},T_e^{\,{\text{selected}}},T_h^{\,{\text{selected}}}} \right\}$.
• (b)  
  Update of the plasma properties as functions of ${\theta _{{\text{selected}}}}$ and $T_e^{\,{\text{selected}}}$ for each plasma cell (density, thermal and electrical conductivities, viscosity, exchange term, radiation heat loss, specific heat).
• (c)  
  Assignment of the boundary conditions for all the solved variables. The boundary condition for the magnetic vector potential is recomputed once per 2 time steps during the initial 30 time steps and once per 100 time steps afterwards.
• (d)  
  Assignment of the source terms for velocity: user source terms and Lorentz force $\vec J \wedge \vec B$.
• (e)  
  Solution of the equations of mass and momentum conservation, computation of pressure.
• (f)  
  Assignment of the source terms for the electron enthalpy: exchange term, Joule power, heat transport by electric current, continuum radiation.
• (g)  
  Solution of the equation of the electron enthalpy.
• (h)  
  Solution of the equation of electric potential
• (i)  
  * Recalculation of the electric field $\vec E$, electric current density $\vec J$ and Joule power $\vec J \cdot \vec E$
• (j)  
  * Voltage scaling algorithm.
• (k)  
  Solution of the equation of the magnetic vector potential
• (l)  
  * Recalculation of the magnetic field $\vec B$ and Lorentz force $\vec J \wedge \vec B$.
• (m)  
  Assignment of the source terms for the heavy species enthalpy: exchange term, line radiation.
• (n)  
  Solution of the equation of the heavy species enthalpy.
• (o)  
  Output of the computed values and start of the following iteration.

The computational time cost of the 2-T model is not significantly higher than the time of the LTE model used in this laboratory: 2 h 45 min per 1000 time steps for the 2-T model against 2 h per 1000 time steps for the LTE model. The additional computational time is mostly associated to the extra solved scalar and translation of enthalpies to temperatures. On the other hand, due to considered thermal non-equilibrium and hence lower gradient of electron temperature compared to the gradient of LTE temperature, the electrical conductivity of plasma along the entire path of electric current from the cathode to the anode is much smoother. Thus, the average number of iterations necessary for the solution of the equation for the electric potential decreases from 12 in the LTE model down to 2 in the 2-T model. After the right initialization parameters are found and all the local instabilities are eliminated, the 2-T model becomes more numerically stable than the LTE model.

The velocity inside the electrodes is zeroed through the momentum source terms as suggested in Patankar (1980).

The anode and cathode sheaths are beyond the scope of this study. Thus, the cathode and anode voltage drops are assumed to be constant in time and uniform over the whole surface of each electrode. It is important to mention that, when the cathode sheath is considered, the predicted cathode arc attachment tends to be more realistic (Liang 2019). The cathode and anode sheath voltage drops are assumed to be 10 V and 0 V respectively.

The electrode heating by the electric current is implemented as described in Zhukovskii et al (2021). For that purpose, the heavy species enthalpy field is coupled with the enthalpy of the electrode by the continuity of corresponding temperatures. The continuous transition of the plasma heavy species temperature to the electrode temperature in this case is a simplification, since the model neglects the anode and cathode sheath. Thus, the enthalpy balance at the electrode surface is rather approximate and requires further development.

The line for this enthalpy balance analysis extends from the cathode tip center to the outlet center. Figures 21 and 22 depict the enthalpy balances in electron and heavy species along the line 1 for the 2T(I) and 2T(II) models, respectively. Since the line includes the torch axis and the predicted plasma arc is axisymmetric, the maxima of all the radial temperature profiles lie on this line. Three intervals can be highlighted on the torch axis:

• (a)  
  Area of the highest Joule power from the cathode tip until around 10 mm. Due to the highest constriction of the arc in this area, the Joule power is the highest resulting in the highest plasma temperature. Therefore, the electron enthalpy dissipation is high and mostly negative except in the very vicinity of the cathode tip. The electron enthalpy dissipation stays positive only in the ∼200 μm layer in the cathode tip vicinity, where electrons are heated up not only by Joule power but also by enthalpy dissipation. Another result of the highest plasma temperature is that the continuum radiation of electrons and line radiation of heavy species have the highest values in this area too.
• (b)  
  The heavy species enthalpy dissipations of both formulations are surprisingly very close despite the higher heavy species enthalpy of the 2T(I) formulation than that of the 2T(II) formulation. Actually, the higher heavy species enthalpy gradient in the 2T(I) formulation is mitigated by a low thermal diffusivity of heavy species which is defined as $\left[ {\left( {{\lambda _h} + {\lambda _r}} \right)/C_p^h} \right]$. The specific heat of heavy species with the 2T(I) formulation contains the ionization contribution leading to higher $C_p^h$ value than with the 2T(II) formulation. It results in a lower heavy species enthalpy diffusivity with the 2T(I) formulation compared to the 2T(II) formulation.
• (c)  
  Near the cathode tip, the convection and heat transport by electrons change their signs from negative (0–3 mm) to positive (beyond 3 mm). The plasma in the cathode vicinity being the coldest on the whole torch axis, the motion of electrons and heavy species from the cathode tip downstream until 3 mm has a cooling effect where $\vec J \cdot \vec \nabla {h_e} < 0$. At around 3 mm from the cathode tip the plasma reaches its highest temperature, thus both convection and heat transport by electrons have a positive contribution to the enthalpy balance beyond 3 mm.
• (d)  
  In both 2T(I) and 2T(II) formulations the area of thermal non-equilibrium with $\theta$ from 1.1 to 3.0 is limited to the ∼200 μm layer in the cathode tip vicinity, while the cell size in this area is around 17 μm.
• (e)  
  It has to be underlined that even in the high temperature plasma core where thermal equilibrium is reached, θ has an insignificant excess over unity (values around 1.01) which still makes the exchange term a significant loss of electron enthalpy close to the cathode due to the high electron temperatures with both 2T(I) and 2T(II) formulations.
• (f)  
  The main difference between the 2T(I) and 2T(II) formulations in the area of high Joule power is the fluctuations of the exchange term in the 2T(I) formulation, which are not found in the 2T(II) formulation. These fluctuations are conditioned by the fluctuations of the heavy species temperature due to unstable translation of enthalpies to temperatures, and therefore fluctuations of $\theta$ from 1.001 to 1.1. In addition, the exchange term goes to zero at around 8 mm from the cathode tip with the 2T(I) formulation as a result of the issue of $\theta < 1$.
• (g)  
  Note that the electron enthalpy gain through the Joule power is too high for the 2T(I) formulation due to the low enthalpy storage in electrons with this formulation, hence resulting in a higher exchange term.
• (h)  
  Area of the medium Joule power from around 10 mm from the cathode tip until the end of the arc core at the anode arc attachment at around 35 mm.
• (i)  
  In this area, the Joule power in both 2T(I) and 2T(II) formulations is significant enough to create a temperature decrease towards the channel wall and, thus, high electron enthalpy dissipation losses in both formulations. The electron enthalpy dissipation decreases continuously from the cathode tip to the outlet because of the continuous widening and flattening of the radial profile of the electron temperature as seen in figures 14, 17 and 18.
• (j)  
  The striking differences between figures 21 and 22 are the changes in the exchange term and convection. Since the ionization energy in the 2T(I) formulation is assigned to the heavy species, all the energy losses have a higher impact on electrons and lower impact on the heavy species. It makes the cooling of the heavy species on the way towards the outside domain slower with the 2T(I) formulation. In addition, the convection source term has higher values for the heavy species due to the higher storage of enthalpy in the heavy species with the 2T(I) formulation as the convection heats up the heavy species too much with the 2T(I) formulation. Hence, the exchange term gets either lower than in the 2T(II) formulation or even equal to zero due to$\,\theta < 1$.
• (k)  
  On the contrary, the heavy species with the 2T(II) formulation gain energy mostly from the exchange term which is mainly balanced by the line radiation.
• (l)  
  Area of the plasma jet without Joule power beyond 35 mm from the cathode tip. The main source of electron enthalpy in this area is the convection, while this energy is lost to dissipation, continuum radiation and for the 2T(II) formulation to the exchange with the heavy species. The main source of heavy species enthalpy with the 2T(I) formulation is convection, while with the 2T(II) formulation it is the exchange term. This difference is still explained by the characteristic behavior of both formulations discussed above.
• (m)  
  In both formulations, the heavy species lose their energy to dissipation and line radiation. However, the radiation is defined as a function of the electron temperature which decreases more rapidly over the course of the plasma jet with the 2T(I) formulation. Hence, the line radiation decreases more rapidly with the 2T(I) formulation. A decrease of the electron temperature below 18,000 K beyond 25 mm from the cathode tip, as seen in figure 15, results in an increase of the heavy species translational and reactive thermal conductivities as can be seen in figure 6. Hence, the heavy species enthalpy dissipation keeps growing along the torch axis towards the outlet, as the electron temperature keeps decreasing. Due to the rapid decrease of radiation with the 2T(I) formulation, the dissipation becomes the main loss of heavy species enthalpy around 44 mm from the cathode tip as seen in figure 21. Meanwhile, the electron temperature does not decrease so fast with the 2T(II) formulation. Therefore, even the decreasing line radiation stays the highest loss with the 2T(II) formulation all the way to the outlet of the domain 74 mm downstream of the cathode tip.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In order to draw a parallel between the other studies of 2-T electric arc models and the present study, the 2T(II) formulation can be compared with the transferred arc modeling results from (Freton et al 2012) for the configuration denoted as 'case 5' with the radiation distributed between electrons and heavy species. The following similarities can be found:

• (a)  
  Joule power is the highest electron enthalpy gain near the cathode, but is surpassed by convection in the arc part closer to the anode.
• (b)  
  Convection is a loss near the cathode and gain in the rest of the arc.
• (c)  
  Line radiation is the main loss for the heavy species.

The main difference with the 'case 5' from Freton et al (2012) is that the electron enthalpy dissipation is significantly higher than the other losses in 'case 5'. In this study, the continuum radiation and exchange term are surpassed by the electron enthalpy dissipation only near the cathode, while in the rest these three terms are very close. The assignation of the reactive thermal conductivity to the heavy species explains this difference.

In the previous section, it has been shown that the ionization energy should be assigned to the electrons which otherwise leads to inconsistent results regarding the exchange term and $\theta$. Moreover, the enthalpy balance of heavy species is also dependent on the radiation loss term which may have different data sources as explained above. Consequently, in order to test the influence of the radiation heat loss data on the predictions of the 2-T model, the data from Beulens et al (1991) was embedded to the 2T(II) configuration of the model similarly to Trelles et al (2007), Baeva et al (2016), (2019), Baeva (2017). Figure 23 depicts the electron and heavy species temperatures calculated in the same conditions as figures 14, 19, and 20 (500 A) but using the radiation data from Beulens et al (1991). The thermal disequilibrium degree $\theta = {T_e}/{T_h}$ is also shown in figure 23 where the blue areas correspond to $\theta < 1$.

Zoom In Zoom Out Reset image size

In comparison with figures 14 and 19 where the radiation data from Erraki (1999) are used, significantly higher temperatures are observed in the whole domain due to lower radiation heat loss as seen in figure 24. Moreover, $\theta < 1$ in a small region at the nozzle exit occurs but with less extension than in figure 20 (2T(I) formulation).

Zoom In Zoom Out Reset image size

The radial profile of the electron temperature at nozzle exit displayed in figure 25 confirms the overestimation of the electron temperature obtained with the Beulens data close to the torch axis. Actually, as shown in figure 9 and also in figure 25, the Beulens radiation data are valid only up to 15,000 K.

Zoom In Zoom Out Reset image size

The enthalpy balances shown in figure 24 highlight lower radiation in the hot region compared to figure 22. The heat dissipation has a much higher role in the enthalpy balance than with the data from Erraki (1999) for both electrons and heavy species due to the higher gradients of both temperatures.

Beyond 30 mm from the cathode tip, where the Joule power decreases to zero, the enthalpy balance of electrons is mainly made between convection and enthalpy dissipation. Meanwhile, the exchange term collapses due to $\theta < 1$ because of too low line radiation losses of heavy species in the hot region upstream of the anode. In other words, the radiation losses of heavy species are not sufficient to establish a proper enthalpy balance in the model.

Attention is hereafter paid to the influence of the arc current on the disequilibrium degree keeping the same other calculation conditions as in figure 23. Figure 26 shows the electron and heavy species temperatures calculated with I = 100 A. As expected, the maximum temperature is significantly lower than with I = 500 A and remains below 15,000 K. Therefore, the radiation heat loss from Beulens et al (1991) in the obtained temperature range is close to the data from Erraki (1999) and the thermal disequilibrium degree is never below unity at such a low arc current.

Zoom In Zoom Out Reset image size

As seen in figure 26, some axial asymmetry of the electron and heavy species temperature is observed near the anode. Such asymmetry was not found in the simulations with 500 A with both formulations of enthalpy and both sources of radiation loss data (figures 14, 19 and 23). The asymmetry seems to be conditioned by the low electron temperature in the vicinity of the anode surface and hence too low electrical conductivity, which causes fluctuations of the anode arc attachment. In order to illustrate these fluctuations, the electric current density over several cross-sections perpendicular to the torch axis is shown in figure 27. It is seen that the electric current density is axisymmetric upstream of the anode, while at the anode some deviation from axial symmetry is found. The electric current is still spread over the whole circumference of the anode, but it is much less uniform than with 500 A. It is worth noting that, in reality, the SinplexPro plasma torch operated at the low current of 100 A is in fact not stable. The measured voltage does not oscillate but drifts almost randomly with the current of 100 A.

Zoom In Zoom Out Reset image size