2.1. Fluid governing equations

Multi-species hypersonic reacting flows are governed by the Navier–Stokes equations with thermochemical non-equilibrium processes described by a multi-temperature model.

The continuity equation consists of a total mass equation and $(ns - 1)$ species equations:

(2.1)\begin{gather}\frac{{\partial \rho }}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{V}) = 0,\end{gather}

(2.2)\begin{gather}\frac{{\partial \rho {y_s}}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho {y_s}\boldsymbol{V}) = \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho {D_s}\boldsymbol{\nabla }{y_s}) + {\omega _s},\quad s = 1,2, \ldots ,ns - 1,\end{gather}

where the index $ns$ represents the total number of species, $\; t$ is the time, $\rho $ is the mixture density and ${y_s}$ is the mass fraction of species s. Parameter $\; \boldsymbol{V}$ denotes the velocity and ${\omega _s}$ represents the source term of production of species s. Coefficient ${D_s}$ is the effective diffusion coefficient of species s.

The momentum conservation equation is

(2.3)\begin{equation}\frac{{\partial \rho \boldsymbol{V}}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\rho \boldsymbol{VV}) + \boldsymbol{\nabla }p = \boldsymbol{\nabla }\boldsymbol{\cdot }\bar{\boldsymbol{\tau}},\end{equation}

where p is the pressure of the mixture and $\bar{\boldsymbol{\tau}}$ denotes the viscous stress tensor (Garicano-Mena, Lani & Degrez Reference Garicano-Mena, Lani and Degrez2018).

The total energy conservation equation is

(2.4)\begin{align}\frac{{\partial \rho E}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\boldsymbol{V}\rho H) = \boldsymbol{\nabla }\boldsymbol{\cdot }\left( {{k_{tr}}\boldsymbol{\nabla }T + \sum\limits_{m = 1}^{nmol} {{k_{{v_m}}}\boldsymbol{\nabla }{T_{{v_m}}}} + {k_e}\boldsymbol{\nabla }{T_e} + \bar{\boldsymbol{\tau}}\boldsymbol{\cdot }\boldsymbol{V} + \sum\limits_{s = 1}^{ns} {{H_s}\rho {D_s}\boldsymbol{\nabla }{y_s}} } \right),\end{align}

where E and H represent the total energy and total enthalpy per unit mass of mixture, respectively, and ${H_s}$ is enthalpy per unit mass of species s. Index $nmol$ denotes the total number of molecular species. Parameters ${k_{tr}}$, ${k_{{v_m}}}$ and ${k_e}$ represent translational–rotational thermal conductivity of mixture, vibrational thermal conductivity of molecular species m and electron–electronic thermal conductivity of mixture, respectively. It is assumed that the translational temperature is equal to the rotational temperature and T denotes translational–rotational temperature. Here ${T_{{v_m}}}$ is vibrational temperature of molecular species m. In this study, the vibrational temperature of ${\textrm{O}_2}$, ${T_{{v_{{\textrm{O}_2}}}}}$, and the vibrational temperature of ${\textrm{N}_2}$, ${T_{{v_{{\textrm{N}_2}}}}}$, are considered. The electronic–free electron modes are described by electron–electronic excitation temperature ${T_e}$.

The vibrational energy conservation equation for a non-equilibrium molecule of species m is

(2.5)\begin{equation}\frac{{\partial \rho {y_m}{E_{{v_m}}}}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\boldsymbol{V}\rho {y_m}{E_{{v_m}}}) = \boldsymbol{\nabla }\boldsymbol{\cdot }({k_{{v_m}}}\boldsymbol{\nabla }{T_{{v_m}}} + {E_{{v_m}}}\rho {D_m}\boldsymbol{\nabla }{y_m}) + {\omega _{{v_m}}},\end{equation}

where ${E_{{v_m}}}$ is vibrational energy per unit mass of molecular species $m$ and ${\omega _{{v_m}}}$ represents the source term of vibrational energy of species m.

The electron–electronic energy conservation equation is

(2.6)\begin{equation}\frac{{\partial \rho {E_e}}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\boldsymbol{V}\rho {E_e}) = \boldsymbol{\nabla }\boldsymbol{\cdot }\left( {{k_e}\boldsymbol{\nabla }{T_e} + \sum\limits_{s = 1}^{ns} {{H_{{e_s}}}\rho {D_s}\boldsymbol{\nabla }{y_s}} } \right) + {\omega _e},\end{equation}

where ${E_e}$ and ${H_{{e_s}}}$ denote electron–electronic energy per unit mass of mixture and electron–electronic enthalpy per unit mass of species s, respectively, and $\; {\omega _e}$ represents the source term of electron–electronic energy of mixture.

2.2. Thermodynamic relations

The pressure of mixture is the sum of the partial pressures of each species:

(2.7)\begin{equation}p = \sum\limits_{s = 1}^{ns} {{p_s}} = \sum\limits_{s \ne {\textrm{e}^ - }} {\frac{{{\rho _s}{R_u}T}}{{{M_s}}} + \frac{{{\rho _e}{R_u}{T_e}}}{{{M_e}}}} ,\end{equation}

where ${p_s}$ and ${M_s}$ are the partial pressure and molecular weight of species $s$, respectively, ${M_e}$ is the electron molecular weight and ${R_u}$ is the universal gas constant.

The internal energy per unit mass of species $s$, i.e.

(2.8)\begin{equation}{E_s} = {E_{t{r_s}}} + {E_{ro{t_s}}} + {E_{{v_s}}} + {E_{{e_s}}} + \Delta {h_s},\end{equation}

is split between the translational, rotational, vibrational and electron–electronic contributions and the latent chemical energy of species s.

The translational energy and rotational energy per unit mass of species s are given by

(2.9a,b)\begin{equation}{E_{t{r_s}}} = \left\{ {\begin{array}{@{}ll} {\dfrac{{3{R_u}T}}{{2{M_s}}}}&{s \ne {\textrm{e}^ - }}\\ 0&{s = {\textrm{e}^ - }} \end{array}} \right.\quad \textrm{and}\quad {E_{ro{t_s}}} = \left\{ {\begin{array}{@{}ll} {\dfrac{{{R_u}T}}{{{M_s}}}}&{s = dia}\\ 0&{s \ne dia,} \end{array}} \right.\end{equation}

where index $dia$ represents diatomic molecule.

The vibrational energy and electron–electronic energy per unit mass of species s are

(2.10a,b)\begin{align}{E_{{v_s}}} = \left\{ {\begin{array}{@{}ll} {\dfrac{{{R_u}{\theta_{{v_s}}}/{M_s}}}{{{\textrm{e}^{{\theta_{{v_s}}}/{T_{{v_s}}}}} - 1}}}&{s = dia}\\ 0&{s \ne dia} \end{array}} \right.\quad \textrm{and}\quad {E_{{e_s}}} = \left\{ {\begin{array}{@{}ll} {\dfrac{{\left( {\dfrac{{{g_{{1_s}}}}}{{{g_{{0_s}}}}}} \right)\,{\textrm{e}^{ - {\theta_{{e_s}}}/{T_e}}}}}{{1 + \left( {\dfrac{{{g_{{1_s}}}}}{{{g_{{0_s}}}}}} \right)\,{\textrm{e}^{ - {\theta_{{e_s}}}/{T_e}}}}}\dfrac{{{R_u}{\theta_{{e_s}}}}}{{{M_s}}}}&{s \ne {\textrm{e}^ - }}\\ {\dfrac{{3{R_u}{T_e}}}{{2{M_e}}}}&{s = {\textrm{e}^ - },} \end{array}} \right.\end{align}

where ${\theta _{{v_s}}}$ and ${\theta _{{e_s}}}$ denote the characteristic vibrational temperature and the characteristic electronic temperature of species s, respectively, and ${g_{{0_s}}}$ and ${g_{{1_s}}}$ are the degeneracy of the ground and the first electronic level of species s, respectively.

2.3. Chemistry model

This study selects the Gupta chemical reaction model (Gupta et al. Reference Gupta, Yos, Thompson and Lee1990) and the Park model (Park Reference Park1993) which are applied to 11 air species $\textrm{(}{\textrm{O}_2},{\textrm{N}_2},\,\textrm{N},\,\textrm{O},\textrm{ NO},\, {\textrm{N}^ + },\, {\textrm{O}^ + }, \textrm{N}_2^ + ,\; \; \textrm{O}_2^ + ,\textrm{ N}{\textrm{O}^ + },\, {\textrm{e}^ - })$. The mass production rate ${\omega _s}$ is calculated by

(2.11)\begin{equation}{\omega _s} = {M_s}\sum\limits_{z = 1}^{nr} {(\nu _{z,s}^b - \nu _{z,s}^f)R_z^{fb}} = {M_s}\sum\limits_{z = 1}^{nr} {(\nu _{z,s}^b - \nu _{z,s}^f)(R_z^f - R_z^b)} ,\end{equation}

where $\nu _{z,s}^f$ and $\nu _{z,s}^b$ are the stoichiometric coefficients of species s for the reactants and products in reaction z, respectively, and $R_z^{fb}$ denotes the total reaction rate of reaction z. The forward and backward reaction rates can be written as

(2.12a,b)\begin{equation}R_z^f = {\varGamma _z}k_z^f\prod\limits_{l = 1}^{ns} {{{\left( {\frac{{{\rho_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^f}}} ,\quad R_z^b = {\varGamma _z}k_z^b\prod\limits_{l = 1}^{ns} {{{\left( {\frac{{{\rho_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^b}}} ,\end{equation}

where ${\varGamma _z}$ is given as follows:

(2.13)\begin{equation}{\varGamma _z} = {\left[ {\sum\limits_{l = 1}^{ns} {{\alpha_{z,l}}} \left( {\frac{{{\rho_l}}}{{{M_l}}}} \right)} \right]^{{L_z}}},\end{equation}

where ${\alpha _{z,l}}$ denotes the third-body efficiency for species l in reaction z. Here ${L_z}$ equals 1 when a third body M is needed in reaction z. Otherwise, ${L_z}$ equals 0.

When the Gupta model is employed, the forward and backward reaction rate coefficients, $k_z^f$ and $k_z^b$, are expressed in modified Arrhenius form as

(2.14a,b)\begin{equation}k_z^f = {A_{f,z}}T_f^{{B_{f,z}}}\,{\textrm{e}^{ - {C_{f,z}}/{T_f}}},\quad k_z^b = {A_{b,z}}T_b^{{B_{b,z}}}\,{\textrm{e}^{ - {C_{b,z}}/{T_b}}},\end{equation}

where ${T_f}$ and ${T_b}$ represent the controlling temperature of the forward and backward reactions, respectively. When the Park model is applied, $k_z^f$ is still calculated by (2.14a,b) and $k_z^b$ is obtained by the equilibrium constant (Park Reference Park1985):

(2.15)\begin{equation}k_z^b = \frac{{k_z^f({T_b})}}{{{K_{eq,z}}({T_b})}},\end{equation}

where the equilibrium constant ${K_{eq,z}}$ is given as the following curve-fitting expression (Park Reference Park1990):

(2.16)\begin{equation}{K_{eq,z}}({T_b}) = exp \left( {\frac{{A_1^z}}{S} + A_2^z + A_3^z\ln S + A_4^zS + A_5^z{S^2}} \right),\end{equation}

where $S = 10\;000/{T_b}$ and the constants $A_1^z,A_2^z,A_3^z,A_4^z,A_5^z$ are given by Park (Reference Park1990). Those five constants for ${\textrm{N}_2} + {\textrm{N}^ + }\rightleftharpoons \textrm{N} + \textrm{N}_2^ + $ are not included in Park (Reference Park1990), so the equilibrium constant ${K_{eq,z}}$ for ${\textrm{N}_2} + {\textrm{N}^ + }\rightleftharpoons \textrm{N} + \textrm{N}_2^ + $ is calculated using the similar method proposed by Park (Reference Park1985) and Ghezali, Haoui & Chpoun (Reference Ghezali, Haoui and Chpoun2019). Radiative recombination reactions in the Park model are not considered, because radiation is ignored in this paper.

2.4. Energy exchange model

The vibrational energy source term of molecular species m can be decomposed as

(2.17)\begin{equation}{\omega _{{v_m}}} = {Q_{t - {v_m}}} + {Q_{v - {v_m}}} + {Q_{e - {v_m}}} + {Q_{r - {v_m}}},\end{equation}

where ${Q_{t - {v_m}}}$ denotes the energy exchange between translation and vibration energies. Energy exchange ${Q_{t - {v_m}}}$ is evaluated using Landau–Teller theory (Landau & Teller Reference Landau and Teller1936) as

(2.18)\begin{equation}{Q_{t - {v_m}}} = \rho {y_m}\frac{{{E_{{v_m}}}(T) - {E_{{v_m}}}({T_{{v_m}}})}}{{{\tau _{vt,m}}}},\end{equation}

and the relaxation time of vibration–translation exchanges ${\tau _{vt,m}}$ modelled by the Millikan–White empirical formula (Millikan & White Reference Millikan and White1963; Miller et al. Reference Miller, Tannehill, Lawrence and Edwards1995) is expressed as

(2.19)\begin{gather}\tau _{vt,m}^{MW} = {{\sum\limits_{s \ne {\textrm{e}^ - }} {{X_s}} } / {\sum\limits_{s \ne {\textrm{e}^ - }} {\frac{{{X_s}}}{{{\tau _{vt,m,s}}}}} }},\end{gather}

(2.20)\begin{gather}{\tau _{vt,m,s}} = \frac{{101\;325}}{p}exp [{a_{m,s}}({T^{ - 1/3}} - {b_{m,s}}) - 18.42],\end{gather}

where ${X_s}$ is the mole fraction of heavy particle s. Here ${a_{m,s}}$ and ${b_{m,s}}$ are given by

(2.21)\begin{gather}{a_{m,s}} = 0.00116\mu _{m,s}^{1/2}\theta _{{v_m}}^{4/3},\end{gather}

(2.22)\begin{gather}{b_{m,s}} = 0.015\mu _{m,s}^{1/4},\end{gather}

where the reduced mass ${\mu _{m,s}}$ is defined as

(2.23)\begin{equation}{\mu _{m,s}} = \frac{{1000{M_m}{M_s}}}{{{M_m} + {M_s}}}.\end{equation}

When the temperature is higher than 8000 K, Park (Reference Park1990) suggests that the correction term should be added on the basis of $\tau _{vt,m}^{MW}$ and the total relaxation time ${\tau _{vt,m}}$ can be expressed as

(2.24)\begin{equation}{\tau _{vt,m}} = \tau _{vt,m}^{MW} + {({\sigma _m}{\bar{C}_m}{n_m})^{ - 1}},\end{equation}

where the average molecular velocity ${\bar{C}_m}$, the effective cross-section ${\sigma _m}$ and number density ${n_m}$ of species m are given by

(2.25)\begin{gather}{\bar{C}_m} = {\left( {\frac{{8{R_u}T}}{{{\rm \pi} {M_m}}}} \right)^{1/2}},\end{gather}

(2.26)\begin{gather}{\sigma _m} = {10^{ - 21}}{\left( {\frac{{50\;000}}{T}} \right)^2},\end{gather}

(2.27)\begin{gather}{n_m} = \frac{{\rho {y_m}{N_0}}}{{{M_m}}},\end{gather}

where ${N_0}$ is Avogadro's constant. Time ${\tau _{vt,m}}$ can be expressed as

(2.28)\begin{equation}{\tau _{vt,m}} = {f_{{\tau _{vt,m}}}}(\rho ,T,{y_1}, \ldots ,{y_{ns}}) = \frac{1}{\rho }{f^{\prime}_{{\tau _{vt,m}}}}(T,{T_e},{y_1}, \ldots ,{y_{ns}}),\end{equation}

where ${f_{{\tau _{vt,m}}}}$ stand for the functions by which ${\tau _{vt,m}}$ can be calculated. The term ${f^{\prime}_{{\tau _{vt,m}}}}$ does not involve the density $\rho $. Thus, the term ${Q_{t - {v_m}}}$ can be written as

(2.29)\begin{equation}{Q_{t - {v_m}}} = {\rho ^2}{f^{\prime}_{{Q_{t - {v_m}}}}}(T,{T_{{v_m}}},{T_e},{y_1}, \ldots ,{y_{ns}}),\end{equation}

where ${f^{\prime}_{{Q_{t - {v_m}}}}}$ is independent of $\rho $.

The vibrational energy exchange between the different molecules ${Q_{v - {v_m}}}$ is described according to Roberts (Reference Roberts1994):

(2.30)\begin{align}{Q_{v - {v_m}}} &= \frac{{\rho {y_m}}}{{{M_m}}}\mathop \sum \limits_{s = 1,s \ne m}^{mol} \left\{ \frac{{{X_s}}}{{{\tau_{vv,m,s}}}}\left[ {\frac{{1 - exp \left( { - \frac{{{\theta_{{v_m}}}}}{T}} \right)}}{{1 - exp \left( { - \frac{{{\theta_{{v_s}}}}}{T}} \right)}}} \right]\frac{{{e_{{v_s}}}({T_{{v_s}}})}}{{{e_{{v_s}}}(T)}}[{e_{{v_m}}}(T) - {e_{{v_m}}}({T_{{v_m}}})]\right.\nonumber\\ &\quad - \left.\frac{{{X_s}}}{{{\tau_{vv,m,s}}}}\frac{{{e_{{v_m}}}({T_{{v_m}}})}}{{{e_{{v_s}}}(T)}}[{e_{{v_s}}}(T) - {e_{{v_s}}}({T_{{v_s}}})] \vphantom{\left\{ \frac{{{X_s}}}{{{\tau_{vv,m,s}}}}\left[ {\frac{{1 - exp \left( { - \frac{{{\theta_{{v_m}}}}}{T}} \right)}}{{1 - exp \left( { - \frac{{{\theta_{{v_s}}}}}{T}} \right)}}} \right]\frac{{{e_{{v_s}}}({T_{{v_s}}})}}{{{e_{{v_s}}}(T)}}[{e_{{v_m}}}(T) - {e_{{v_m}}}({T_{{v_m}}})]\right.}\right\},\end{align}

where ${e_{{v_m}}}$ denotes the vibrational energy per mole of species m. Note that species m can only be ${\textrm{N}_2}$ in (2.30). For other species, the transfer of vibrational energy can be obtained according to the assumption of conservation of vibrational energy (Miller et al. Reference Miller, Tannehill, Lawrence and Edwards1995). The relaxation time of vibration–vibration exchanges ${\tau _{vv,m,s}}$ can be written as

(2.31)\begin{equation}{\tau _{vv,m,s}} = \frac{1}{{{P_{m,s}}{{(8{\rm \pi} KT{N_0}/{\mu _{m,s}})}^{1/2}}\sigma _{m,s}^2{n_s}}},\end{equation}

where K and ${\sigma _{m,s}}$ represent Boltzmann's constant and the collision cross-sections (Miller et al. Reference Miller, Tannehill, Lawrence and Edwards1995), respectively. The probability of energy transfer during each collision ${P_{m,s}}$ is provided by Park (Reference Park1993) as

(2.32)\begin{gather}{P_{{\textrm{N}_2},\textrm{NO}}} = 5.5 \times {10^{ - 5}} {\left( {\frac{T}{{1000}}} \right)^{2.32}},\end{gather}

(2.33)\begin{gather}{P_{{\textrm{N}_2},{\textrm{O}_2}}} = 3.0 \times {10^{ - 6}} {\left( {\frac{T}{{1000}}} \right)^{2.87}}.\end{gather}

In this paper, only collision between ${\textrm{N}_2}$ and ${\textrm{O}_2}$ is taken into consideration. Like ${Q_{t - {v_m}}}$, the term ${Q_{v - {v_m}}}$ can be expressed as

(2.34)\begin{equation}{Q_{v - {v_m}}} = {\rho ^2}{f^{\prime}_{{Q_{v - {v_m}}}}}(T,{T_{{v_m}}},{y_m},{y_{s1}}, \ldots ,{y_{smol}}),\end{equation}

where ${y_{s1}}, \ldots ,{y_{smol}}$ represent the molecular species except species m.

The energy exchange between vibrational and electronic modes ${Q_{e - {v_m}}}$ is modelled as (Candler & Maccormack Reference Candler and Maccormack1991)

(2.35)\begin{equation}{Q_{e - {v_m}}} = \rho {y_m}\frac{{{E_{{v_m}}}({T_e}) - {E_{{v_m}}}({T_{{v_m}}})}}{{{\tau _{e,m}}}},\end{equation}

where ${\tau _{e,m}}$ denotes the relaxation time of electronic–vibrational exchanges for species m. Parameter $\; {\tau _{e,{\textrm{N}_2}}}$ is given based on research of Lee (Reference Lee1993) and Candler (Reference Candler1988) by

(2.36)\begin{gather}{\log _{10}}\left( {\frac{{{p_e}{\tau_{e,{\textrm{N}_2}}}}}{{101\;325}}} \right) = 3.91\log _{10}^2({T_e}) - 30.36{\log _{10}}({T_e}) + 48.90,\quad 1000\;\textrm{K} \le {T_e} \le 7000\;\textrm{K},\end{gather}

(2.37)\begin{gather}{\log _{10}}\left( {\frac{{{p_e}{\tau_{e,{\textrm{N}_2}}}}}{{101\;325}}} \right) = 1.30\log _{10}^2({T_e}) - 9.09{\log _{10}}({T_e}) + 5.58,\quad 7000 < {T_e} \le 50\;000\;\textrm{K},\end{gather}

where ${p_e}$ is the electron partial pressure.

Parameter ${\tau _{e,{\textrm{O}_2}}}$ is obtained from the formula of Park & Lee (Reference Park and Lee1995) based on ${\tau _{e,{\textrm{N}_2}}}$:

(2.38)\begin{equation}{\tau _{e,{\textrm{O}_2}}} = 300 \times {\tau _{e,{\textrm{N}_2}}}(1.492{T_e}).\end{equation}

The term ${Q_{e - {v_m}}}$ can also be written as

(2.39)\begin{equation}{Q_{e - {v_m}}} = {\rho ^2}{f^{\prime}_{{Q_{e - {v_m}}}}}(T,{T_e},{y_m},{y_{{\textrm{e}^ - }}}),\end{equation}

where ${f^{\prime}_{{Q_{e - {v_m}}}}}$ does not depend on $\rho $.

The change in vibrational energy due to chemical reactions ${Q_{r - {v_m}}}$ can be expressed as (Miller et al. Reference Miller, Tannehill, Lawrence and Edwards1995)

(2.40)\begin{equation}{Q_{r - {v_m}}} = {\omega _m}{E_{{v_m}}}.\end{equation}

The electron–electronic energy source term can be decomposed as

(2.41)\begin{equation}{\omega _e} = {Q_{t - e}} - \sum\limits_{m = mol} {{Q_{e - {v_m}}}} + {Q_{r - e}} + {Q_{i - e}} - {p_e}\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{V}.\end{equation}

The energy exchange between translation and electron modes ${Q_{t - e}}$ can be written as (Miller et al. Reference Miller, Tannehill, Lawrence and Edwards1995; Roberts Reference Roberts1996)

(2.42)\begin{equation}{Q_{t - e}} = 3{\rho _e}{R_u}(T - {T_e})\sum\limits_{s \ne {\textrm{e}^ - }} {\frac{{{v_{e,s}}}}{{{M_s}}}} ,\end{equation}

where ${v_{e,s}}$ denotes the collision frequency between electrons and species s. For ions and electrons, it follows Coulomb collisions and can be expressed as (Appleton & Bray Reference Appleton and Bray1964)

(2.43)\begin{equation}{\nu _{e,s}} = \frac{8}{3}{\left( {\frac{{\rm \pi} }{{{m_e}}}} \right)^{1/2}}{n_s}e_c^4\frac{1}{{{{(2K{T_e})}^{3/2}}}}\ln \left( {\frac{{{K^3}T_e^3}}{{{\rm \pi} {n_e}e_c^6}}} \right),\end{equation}

where ${m_e}$ is the electron mass, ${n_e}$ denotes number density of electrons and ${e_c}$ is the magnitude of the electronic charge. For collisions between neutral species and electrons, the collision frequency can be expressed as (Appleton & Bray Reference Appleton and Bray1964)

(2.44)\begin{equation}{\nu _{e,s}} = {n_s}{\sigma _{e,s}}{\left( {\frac{{8K{T_e}}}{{{\rm \pi} {m_e}}}} \right)^{1/2}},\end{equation}

where the effective cross-section ${\sigma _{e,s}}$ is evaluated by curve fits (Gnoffo, Gupta & Shinn Reference Gnoffo, Gupta and Shinn1989) as

(2.45)\begin{equation}{\sigma _{e,s}} = {a_s} + {b_s}{T_e} + {c_s}T_e^2.\end{equation}

The collision frequency ${v_{e,s}}$ and the effective cross-section ${\sigma _{e,s}}$ can also be calculated by similar formulae obtained by Petschek & Byron (Reference Petschek and Byron1957) or Matsuzaki (Reference Matsuzaki1988) and the simulation in § 4 uses Matsuzaki's model. Regardless of the model mentioned, the term ${Q_{t - e}}$ can be written as

(2.46)\begin{align}\begin{aligned}{Q_{t - e}} & = {Q_{t - e\_neu}} + {Q_{t - e\_ions}}\\ & = {\rho ^2}\textrm{[}f_{t - e}^{neu}(T,{T_e},{y_{{\textrm{e}^ - }}},{y_{neu1}}, \ldots ,{y_{neun}}) + f_{t - e}^{ion}(\rho ,T,{T_e},{y_{{\textrm{e}^ - }}},{y_{ion1}}, \ldots ,{y_{ionn}})\textrm{] ,} \end{aligned}\end{align}

where ${Q_{t - e\_neu}}$ denotes the translation–electron energy exchange from collisions between electrons and neutral particles and ${Q_{t - e\_ions}}$ denotes that from collisions between electrons and ions. The term $f_{t - e}^{neu}$ represents the contribution of the electron–neutral particle collisions and is not correlated with $\rho $; $f_{t - e}^{ion}$ stands for the contribution of the electron–ion collisions and depends on $\rho $.

The added or removed electron energy due to chemical reactions ${Q_{r - e}}$ can be specified as (Hao, Wang & Lee Reference Hao, Wang and Lee2016)

(2.47)\begin{equation}{Q_{r - e}} = \sum\limits_{s = 1}^{ns} {{\omega _s}{E_{{e_s}}}} .\end{equation}

The term ${Q_{i - e}}$ represents the energy loss of free electrons by the impact ionization reactions (Scalabrin Reference Scalabrin2007), which gives

(2.48)\begin{equation}{Q_{i - e}} ={-} {\textstyle{1 \over 3}}({I_\textrm{O}}{\omega _{{\textrm{O}^ + }}} + {I_\textrm{N}}{\omega _{{\textrm{N}^ + }}}),\end{equation}

where ${I_\textrm{O}}$ and ${I_\textrm{N}}$ are the first energy of ionization of $\textrm{O}$ and $\textrm{N}$; ${\omega _{{\textrm{O}^ + }}}$ and ${\omega _{{\textrm{N}^ + }}}$ are the formation rate of ${\textrm{O}^ + }$ and ${\textrm{N}^ + }$ in the electron-impact ionization reactions.

The term $- {p_e}\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{V}$ denotes the contribution of electron partial pressure in electron–electronic energy change (Gnoffo et al. Reference Gnoffo, Gupta and Shinn1989; Lee Reference Lee1984).

2.5. Transport properties

The dynamic viscosity coefficient of each species is obtained by curve fits taken from Gupta et al. (Reference Gupta, Yos, Thompson and Lee1990) in the temperature range (1000–30 000 K) and Lennard-Jones potential (Anderson Reference Anderson2006) in the other temperature range. The translational–rotational, vibrational and electron–electronic thermal conductivities of each species are given by Eucken's relation (Anderson Reference Anderson2006). The dynamic viscosity $\mu $ and thermal conductivities of gas mixture are derived from Wilke's rule (Wilke Reference Wilke1950). So $\mu $, ${k_{tr}}$, ${k_{{v_m}}}$ and ${k_e}$ can be expressed as

(2.49)\begin{gather}\mu = {f_\mu }(T,{T_e},{y_1}, \ldots ,{y_{ns}}),\end{gather}

(2.50)\begin{gather}{k_{tr}} = {f_{{k_{tr}}}}(T,{y_1}, \ldots ,{y_{ns}}),\end{gather}

(2.51)\begin{gather}{k_{{v_m}}} = {f_{{k_{{v_m}}}}}(T,{T_{{v_m}}},{y_1}, \ldots ,{y_{ns}}),\end{gather}

(2.52)\begin{gather}{k_e} = {f_{{k_e}}}(T,{T_e},{y_1}, \ldots ,{y_{ns}}),\end{gather}

where ${f_\mu }$, ${f_{{k_{tr}}}}$, ${f_{{k_{{v_m}}}}}$ and ${f_{{k_e}}}$ stand for the functions by which $\mu $, ${k_{tr}}$, ${k_{{v_m}}}$ and ${k_e}$ can be calculated, respectively. One can observe that $\mu $, ${k_{tr}}$, ${k_{{v_m}}}$ and ${k_e}$ depend on temperatures and mass fraction of each species.

It is assumed that the binary diffusion coefficients between each pair of components ${D_{sj}}$ are the same and the electric current of the flow field is zero. The multi-component diffusion coefficient of each species ${D_s}$ is given as follows:

(2.53)\begin{align}{D_s} = \left\{ {\begin{array}{@{}ll} {\dfrac{{(1 - {y_s})\mu }}{{\rho (1 - {X_s})Sc}}}&{s \ne {\textrm{e}^ - },}\\ {\dfrac{{\sum\limits_{s = 1}^{nion} {{y_s}{D_s}} }}{{\sum\limits_{s = 1}^{nion} {{y_s}} }}}&{s = {\textrm{e}^ - },} \end{array}} \right.\end{align}

where $Sc$ is the Schmidt number and index $nion$ denotes the total number of ion species. For neutral particles, $Sc = 0.5$. For ions, $Sc = 0.25$. In a similar way, ${D_s}$ can also be expressed as

(2.54)\begin{equation}{D_s} = {f_{{D_s}}}(T,{T_e},\rho ,{y_1}, \ldots ,{y_{ns}}) = {f^{\prime}_{{D_s}}}(T,{T_e},{y_1}, \ldots ,{y_{ns}})/\rho ,\end{equation}

where the functions ${f_{{D_s}}}$ and ${f^{\prime}_{{D_s}}}$ show that ${D_s}$ depends on density, temperature and mass fraction of each species.

3.1. Similarity analysis for continuity equation

The analysis based on the dimensionless continuity equation (Zhang Reference Zhang1990) is reviewed first:

(3.1)\begin{equation}\frac{{\partial \bar{\rho }{y_s}}}{{\partial \bar{t}}} + \frac{{\partial \bar{\rho }{y_s}{{\bar{u}}_j}}}{{\partial {{\bar{x}}_j}}} = \frac{1}{{R{e_\infty }}}\frac{{\partial \left( {\frac{{\bar{\rho }{{\bar{D}}_s}\partial {y_s}}}{{\partial {{\bar{x}}_j}}}} \right)}}{{\partial {{\bar{x}}_j}}} + \frac{L}{{{\rho _\infty }{U_\infty }}}{\omega _s},\end{equation}

where the overbar is used to denote non-dimensional variables. Parameter $\; {\bar{D}_s}$ can be expressed in the following form (see (2.49), (2.54) and (A4), (A5) in Appendix A):

(3.2)\begin{align}{\bar{D}_s} = \frac{{{{f^{\prime}}_{{D_s}}}(\bar{T} {T_\infty },{{\bar{T}}_e} {T_\infty },{y_1}, \ldots ,{y_{ns}})/\rho }}{{{f_\mu }( {T_\infty },{T_{e,\infty }},{y_{1,\infty }}, \ldots ,{y_{ns,\infty }})/{\rho _\infty }}} = \frac{{{{f^{\prime}}_{{D_s},{T_\infty }}}(\bar{T},{{\bar{T}}_e},{y_1}, \ldots ,{y_{ns}})}}{{{f_{\mu ,{T_\infty }}}(1,{T_{e,\infty }}/{T_\infty },{y_{1,\infty }}, \ldots ,{y_{ns,\infty }}) \bar{\rho }}},\end{align}

where the index ${T_\infty }$ represents that the functions ${f_{\mu ,{T_\infty }}}$ and ${f^{\prime}_{{D_s},{T_\infty }}}$ are relevant to the temperature of free stream ${T_\infty }$. Formula (3.2) shows that ${\bar{D}_s}$ only depends on dimensional parameter ${T_\infty }$ in addition to non-dimensional parameters.

The second item on the right-hand side of (3.1) can be written as

(3.3)\begin{align}\begin{aligned} rh{s_{c2}} & = \dfrac{L}{{{\rho _\infty }{U_\infty }}}{M_s}\sum\limits_{z = 1}^{nr} {(\nu _{z,s}^b - \nu _{z,s}^f)} {\left[ {\sum\limits_{l = 1}^{ns} {{\alpha_{z,l}}} \left( {\dfrac{{{\rho_l}}}{{{M_l}}}} \right)} \right]^{{L_z}}}\left[ {k_z^f\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{\rho_l}}}{{{M_l}}}} \right)}^{\nu_{z,l}^f}}} - k_z^b\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{\rho_l}}}{{{M_l}}}} \right)}^{\nu_{z,l}^b}}} } \right]\\ & = \dfrac{{{\rho _\infty }L}}{{{U_\infty }}}{M_s}\sum\limits_{z = 1}^{nr} {(\nu _{z,s}^b - \nu _{z,s}^f)} {\left[ {\sum\limits_{l = 1}^{ns} {{\alpha_{z,l}}} \left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)} \right]^{{L_z}}}\left[ {k_z^f{{\bar{\rho }}^{{L_z} + \sum\nolimits_{l = 1}^{ns} {\nu_{z,l}^f} }} \rho_\infty^{{L_z} + \sum\nolimits_{l = 1}^{ns} {\nu_{z,l}^f - 2} }\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)}^{\nu_{z,l}^f}}} } \right.\\ & \left. \quad -\, k_z^b {{\bar{\rho }}^{{L_z} + \sum\nolimits_{l = 1}^{ns} {\nu_{z,l}^b} }} \rho_\infty^{{L_z} + \sum\nolimits_{l = 1}^{ns} {\nu_{z,l}^f - 2} }\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)}^{\nu_{z,l}^b}}} \right]\\ &= \dfrac{{{\rho _\infty }L}}{{{U_\infty }}}{M_s}\sum\limits_{z = 1}^{nr} {(\nu _{z,s}^b - \nu _{z,s}^f)} {\left[ {\sum\limits_{l = 1}^{ns} {{\alpha_{z,l}}} \left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)} \right]^{{L_z}}}{{R^{\prime}}_z}, \end{aligned}\end{align}

where ${\rho _\infty }L$ and ${U_\infty }$ are each required to be the same between two different flows according to the binary scaling law. Terms $k_z^f$ and $k_z^b$ can be expressed as a function of dimensionless parameters like ${\bar{D}_s}$:

(3.4)\begin{equation}k_z^{f/b} = {f_{k_z^{f/b}}}(T,{T_e},{T_{{v_1}}}, \ldots ,{T_{{v_{mol}}}}) = {f_{k_z^{f/b},{T_\infty }}}(\bar{T},{\bar{T}_e},{\bar{T}_{{v_1}}}, \ldots ,{\bar{T}_{{v_{mol}}}}).\end{equation}

As far as elementary reactions of high-temperature air are concerned (Gupta et al. Reference Gupta, Yos, Thompson and Lee1990; Park Reference Park1993), the forward reactions can be only second-order reactions and the backward reactions can be either second-order or third-order reactions. When the reverse of reaction z is a second-order reaction, ${L_z} + \sum\nolimits_{l = 1}^{ns} {\nu _{z,l}^f} = 2$ and ${L_z} + \sum\nolimits_{l = 1}^{ns} {\nu _{z,l}^b} = 2$ in (3.3). Thus, ${R^{\prime}_z}$ in (3.3) can be written as (refer to (3.4))

(3.5)\begin{align}{R^{\prime}_z} &= k_z^f {\bar{\rho }^2}\prod\limits_{l = 1}^{ns} {{{\left( {\frac{{{y_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^f}}} - k_z^b {\bar{\rho }^2}\prod\limits_{l = 1}^{ns} {{{\left( {\frac{{{y_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^b}}}\nonumber\\ &\quad = {f_{{{R^{\prime}}_z},{T_\infty }}}(\bar{\rho },\bar{T},{\bar{T}_e},{\bar{T}_{{v_1}}}, \ldots ,{\bar{T}_{{v_{mol}}}},{y_1}, \ldots ,{y_{ns}}),\end{align}

where ${f_{{{R^{\prime}}_z},{T_\infty }}}$ only depends on the dimensional parameter ${T_\infty }$ in addition to dimensionless parameters.

When the reverse of reaction z is a third-order reaction, ${L_z} + \sum\nolimits_{l = 1}^{ns} {\nu _{z,l}^f} = 2$ and ${L_z} + \sum\nolimits_{l = 1}^{ns} {\nu _{z,l}^b} = 3$. Therefore, ${R^{\prime}_z}$ can be expressed as

(3.6)\begin{align}{R^{\prime}_z} &= k_z^f {\bar{\rho }^2}\prod\limits_{l = 1}^{ns} {{{\left( {\frac{{{y_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^f}}} - k_z^b {\bar{\rho }^3}{\rho _\infty }\prod\limits_{l = 1}^{ns} {{{\left( {\frac{{{y_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^b}}}\nonumber\\ &\quad = {f_{{{R^{\prime}}_z},{T_\infty },{\rho _\infty }}}(\bar{\rho },\bar{T},{\bar{T}_e},{\bar{T}_{{v_1}}}, \ldots ,{\bar{T}_{{v_{mol}}}},{y_1}, \ldots ,{y_{ns}}),\end{align}

where ${f_{{{R^{\prime}}_z},{T_\infty },{\rho _\infty }}}$ depends on both ${T_\infty }$ and ${\rho _\infty }$ in addition to dimensionless parameters. In this case, the similarity between two different hypersonic airflows not only needs the same ${\rho _\infty }L$ but also needs the same ${\rho _\infty }$. In other words, the scale of the flow field L needs to be equal everywhere and the two flows are actually the same.

Generally, the density $\rho $ of the flow field around hypersonic vehicles is relatively low because of the low density of the free stream ${\rho _\infty }$ at high altitude. In this case, the backward third-order reaction rate is far less than the forward second-order reaction rate for reaction z in most areas of the flow field. Thus, for high-altitude flow, formula (3.6) can be written as follows:

(3.7)\begin{align}\begin{aligned} {{R^{\prime}}_z} & = k_z^f {{\bar{\rho }}^2}\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^f}}} - \rho \left[ {k_z^b {{\bar{\rho }}^2}\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)}^{\nu_{z,l}^b}}} } \right] \approx k_z^f {{\bar{\rho }}^2}\prod\limits_{l = 1}^{ns} {{{\left( {\dfrac{{{y_l}}}{{{M_l}}}} \right)}^{\nu _{z,l}^f}}} \\ & = {{f^{\prime}}_{{{R^{\prime}}_z},{T_\infty }}}(\bar{\rho },\bar{T},{{\bar{T}}_e},{{\bar{T}}_{{v_1}}}, \ldots ,{{\bar{T}}_{{v_{mol}}}},{y_1}, \ldots ,{y_{ns}}), \end{aligned}\end{align}

where ${f^{\prime}_{{{R^{\prime}}_z},{T_\infty }}}$ is unrelated to ${\rho _\infty }$. In this case, should ${\rho _\infty }L$, ${T_\infty }$, $\; {U_\infty }$, $\; {y_{1,\infty }}, \ldots ,{y_{ns,\infty }}$ be identical, the airflows over geometrically similar objects will be similar and the binary scaling law can be applied to the continuity equation.

3.2. Similarity analysis for vibrational energy conservation equation

The dimensionless vibrational energy conservation equation is written as

(3.8)\begin{equation}\frac{{\partial \bar{\rho }{y_m}{{\bar{E}}_{{v_m}}}}}{{\partial \bar{t}}} + \frac{{\partial \bar{\rho }{y_m}{{\bar{E}}_{{v_m}}}{{\bar{u}}_j}}}{{\partial {{\bar{x}}_j}}} = \frac{1}{{R{e_\infty }}}\frac{{\partial \left( {{{\bar{k}}_{{v_m}}}\frac{{\partial {{\bar{T}}_{{v_m}}}}}{{\partial {{\bar{x}}_j}}} + \bar{\rho }{{\bar{D}}_m}{{\bar{E}}_{{v_m}}}\frac{{\partial {y_m}}}{{\partial {{\bar{x}}_j}}}} \right)}}{{\partial {{\bar{x}}_j}}} + \frac{L}{{{\rho _\infty }U_\infty ^3}}{\omega _{{v_m}}},\end{equation}

where ${\bar{k}_{{v_m}}}$ and ${\bar{E}_{{v_m}}}$ can be expressed as follows (see (2.49), (2.51) and (A3), (A6) in Appendix A):

(3.9)\begin{gather}{\bar{k}_{{v_m}}} = \frac{{{f_{{k_{{v_m}}}}}(T,{T_{{v_m}}},{y_1}, \ldots ,{y_{ns}})}}{{{\mu _\infty }U_\infty ^2/{T_\infty }}} = \frac{{{f_{{k_{{v_m}}},{T_\infty }}}(\bar{T},{{\bar{T}}_{{v_m}}},{y_1}, \ldots ,{y_{ns}})}}{{{f_{\mu ,{T_\infty }}}(1,{T_{e,\infty }}/{T_\infty },{y_{1,\infty }}, \ldots ,{y_{ns,\infty }})U_\infty ^2/{T_\infty }}},\end{gather}

(3.10)\begin{gather}{\bar{E}_{{v_m}}} = \frac{{{f_{{E_{{v_m}}},{T_\infty }}}({{\bar{T}}_{{v_m}}})}}{{U_\infty ^2}}.\end{gather}

Formulae (3.9) and (3.10) indicate that ${T_\infty }$ and ${U_\infty }$ are the only dimensional parameters on which ${\bar{k}_{{v_m}}}$ and ${\bar{E}_{{v_m}}}$ depend. The second item on the right-hand side of (3.8) can be written as (refer to (2.29), (2.34), (2.39) and (2.40))

(3.11)\begin{align}rh{s_{v2}} &= \frac{{{\rho _\infty }L}}{{\rho _\infty ^2U_\infty ^3}}({Q_{t - {v_m}}} + {Q_{v - {v_m}}} + {Q_{e - {v_m}}} + {Q_{r - {v_m}}}) = \frac{{{\rho _\infty }L}}{{U_\infty ^3}}({\bar{\rho }^2}{f^{\prime}_{{Q_{t - {v_m}}}}} + {\bar{\rho }^2}{f^{\prime}_{{Q_{v - {v_m}}}}}\nonumber\\ &\quad + {\bar{\rho }^2}{f^{\prime}_{{Q_{e - {v_m}}}}}) + \frac{L}{{{\rho _\infty }U_\infty ^3}}{\omega _m}{E_{{v_m}}},\end{align}

where ${f^{\prime}_{{Q_{t - {v_m}}}}}$ can be written as a function of dimensionless parameters by (2.29):

(3.12)\begin{equation}{f^{\prime}_{{Q_{t - {v_m}}}}} = {f^{\prime}_{{Q_{t - {v_m},{T_\infty }}}}}(\bar{T},{\bar{T}_{{v_m}}},{\bar{T}_e},{y_1}, \cdots ,{y_{ns}}),\end{equation}

where ${f^{\prime}_{{Q_{t - {v_m},{T_\infty }}}}}$ only depends on dimensional parameter ${T_\infty }$ besides dimensionless parameters. Similarly, ${f^{\prime}_{{Q_{v - {v_m}}}}}$ and ${f^{\prime}_{{Q_{e - {v_m}}}}}$ can also be expressed as functions ${f^{\prime}_{{Q_{v - {v_m},{T_\infty }}}}}(\bar{T},{\bar{T}_{{v_m}}},{y_m},{y_{s1}}, \ldots ,{y_{smol}})$ and ${f^{\prime}_{{Q_{e - {v_m},{T_\infty }}}}}(\bar{T},{\bar{T}_e},{y_m},{y_{{\textrm{e}^ - }}})$, which are only related to dimensional parameter ${T_\infty }$. These three terms can ensure the similarity of the vibrational energy conservation equation when the preconditions of the binary scaling law are satisfied. The term $(L/{\rho _\infty }U_\infty ^3){\omega _m}{E_{{v_m}}}$ can be analogous to the term $(L/{\rho _\infty }{U_\infty }){\omega _s}$ in continuity equation (3.1) and depends on both ${T_\infty }$ and ${\rho _\infty }$. Thus, the term $(L/{\rho _\infty }U_\infty ^3){\omega _m}{E_{{v_m}}}$ can cause dissimilarity in two hypersonic airflows which have different scales but meet the conditions of the binary scaling law. The fundamental cause of this dissimilarity also lies in the three-body reactions of the reverse reaction. Just as in the continuity equation (3.1), the binary scaling law can be applied to the vibrational energy conservation equation as long as the incoming flow density ${\rho _\infty }$ is low enough.

Thus, only ${Q_{r - {v_m}}}$ in the vibrational energy conservation equation can theoretically result in the failure of the binary scaling law and the root cause is still three-body reactions.

3.3. Similarity analysis for electron–electronic energy conservation equation

The dimensionless electron–electronic energy conservation equation is written as

(3.13)\begin{equation}\frac{{\partial \bar{\rho }{{\bar{E}}_e}}}{{\partial \bar{t}}} + \frac{{\partial \bar{\rho }{{\bar{E}}_e}{{\bar{u}}_j}}}{{\partial {{\bar{x}}_j}}} = \frac{1}{{R{e_\infty }}}\frac{{\partial \left( {{{\bar{k}}_e}\frac{{\partial {{\bar{T}}_e}}}{{\partial {{\bar{x}}_j}}} + \bar{\rho }\sum\limits_{s = 1}^{ns} {{{\bar{D}}_s}{{\bar{H}}_{{e_s}}}} \frac{{\partial {y_s}}}{{\partial {{\bar{x}}_j}}}} \right)}}{{\partial {{\bar{x}}_j}}} + \frac{L}{{{\rho _\infty }U_\infty ^3}}{\omega _e},\end{equation}

where ${\bar{k}_e}$, ${\bar{E}_e}$ and ${\bar{H}_{{e_s}}}$ can be expressed as follows (see (2.49), (2.52) and (A3), (A6) in Appendix A):

(3.14)\begin{gather}{\bar{k}_e} = \frac{{{f_{{k_e}}}(T,{T_e},{y_1}, \ldots ,{y_{ns}})}}{{{\mu _\infty }U_\infty ^2/{T_\infty }}} = \frac{{{f_{{k_e},{T_\infty }}}(\bar{T},{{\bar{T}}_e},{y_1}, \ldots ,{y_{ns}})}}{{{f_{\mu ,{T_\infty }}}(1,{T_{e,\infty }}/{T_\infty },{y_{1,\infty }}, \ldots ,{y_{ns,\infty }})U_\infty ^2/{T_\infty }}},\end{gather}

(3.15)\begin{gather}{\bar{E}_e} = \frac{{{f_{{E_e},{T_\infty }}}({{\bar{T}}_e},{y_1}, \ldots ,{y_{ns}})}}{{U_\infty ^2}},\end{gather}

(3.16)\begin{gather}{\bar{H}_{{e_s}}} = \frac{{{f_{{H_{{e_s}}},{T_\infty }}}({{\bar{T}}_e})}}{{U_\infty ^2}}.\end{gather}

Similar to ${\bar{k}_{{v_m}}}$ and ${\bar{E}_{{v_m}}}$ in § 3.2, the dimensional parameters on which ${\bar{k}_e}$, ${\bar{E}_e}$ and ${\bar{H}_{{e_s}}}$ depend are only ${T_\infty }$ and ${U_\infty }$. The second item on the right-hand side of (3.13) can be written as (refer to (2.39), (2.46), (2.47) and (2.48))

(3.17)\begin{align} rh{s_{e2}} & = \dfrac{{{\rho _\infty }L}}{{\rho _\infty ^2U_\infty ^3}}\left( {{Q_{t - e}} - \mathop \sum \limits_{m = mol} {Q_{e - {v_m}}} + {Q_{r - e}} + {Q_{i - e}} - {p_e}\dfrac{{\partial {u_j}}}{{\partial {x_j}}}} \right)\nonumber\\ &\quad = \dfrac{{{\rho _\infty }L}}{{U_\infty ^3}}\left( {{{\bar{\rho }}^2}f_{t - e}^{ion} + {{\bar{\rho }}^2}f_{t - e}^{neu} - {{\bar{\rho }}^2}\mathop \sum \limits_{m = mol} {{f^{\prime}}_{{Q_{e - {v_m}}}}}} \right)\nonumber\\ & \quad + \dfrac{L}{{{\rho _\infty }U_\infty ^3}}\left( {\mathop \sum \limits_{s = 1}^{ns} {\omega_s}{E_{{e_s}}} - \dfrac{1}{3}({I_\textrm{O}}{\omega_{{\textrm{O}^ + }}} + {I_\textrm{N}}{\omega_{{\textrm{N}^ + }}})} \right) - \dfrac{{L{p_e}}}{{{\rho _\infty }U_\infty ^3}}\dfrac{{\partial {u_j}}}{{\partial {x_j}}}, \end{align}

where $f_{t - e}^{ion}$ and $f_{t - e}^{neu}$ can be written as functions $f_{t - e,{T_{\infty ,}}{\rho _\infty }}^{ion}(\bar{\rho },\bar{T},{\bar{T}_e},{y_{{\textrm{e}^ - }}},{y_{ion1}}, \ldots ,{y_{ionn}})$ and $f_{t - e,{T_\infty }}^{neu}(\bar{T},{\bar{T}_e},{y_{{\textrm{e}^ - }}},{y_{neu1}}, \ldots ,{y_{neun}})$, respectively, based on (2.46), similar to $\; {f^{\prime}_{{Q_{e - {v_m},{T_\infty }}}}}(\bar{T},{\bar{T}_e},{y_m},{y_{{\textrm{e}^ - }}})$ in § 3.2. One can notice that $f_{t - e,{T_\infty }}^{neu}$ only depends on dimensional variable ${T_\infty }$ like ${f^{\prime}_{{Q_{e - {v_m}}},{T_\infty }}}$, whereas $f_{t - e,{T_{\infty ,}}{\rho _\infty }}^{ion}$ not only depends on ${T_\infty }$ but also is dependent on ${\rho _\infty }$. Therefore, it is observed that the energy exchange between translation and electron mode caused by electron–ion collisions, ${Q_{t - e\_ions}}$, can invalidate the binary scaling law. If the energy exchange caused by electron–ion collisions is small enough compared with electron–electronic energy change caused by other means, the binary scaling law will be approximately applicable.

Due to the preceding three-body reactions in the reverse reaction, the electron–electronic energy change resulting from chemical reactions (including $\sum\nolimits_{s = 1}^{ns} {{\omega _s}{E_{{e_s}}}} $, ${I_\textrm{O}}{\omega _{{\textrm{O}^ + }}}$ and ${I_\textrm{N}}{\omega _{{\textrm{N}^ + }}}$) can make the binary scaling law invalid. Similar to the analysis for term $(L/{\rho _\infty }U_\infty ^3){\omega _m}{E_{{v_m}}}$ in (3.11), the binary scaling law can be applied to the electron–electronic energy conservation equation if the incoming flow density ${\rho _\infty }$ is sufficiently low.

The term $- (L{p_e}/{\rho _\infty }U_\infty ^3)(\partial {u_j}/\partial {x_j})$ can be expressed as

(3.18)\begin{equation}- \frac{{L{p_e}}}{{{\rho _\infty }U_\infty ^3}}\frac{{\partial {u_j}}}{{\partial {x_j}}} ={-} \frac{{{T_\infty }{R_u}{y_{{e^ - }}}}}{{U_\infty ^2{M_{{e^ - }}}}}\frac{{\bar{\rho }{{\bar{T}}_e}\partial {{\bar{u}}_j}}}{{\partial {{\bar{x}}_j}}}.\end{equation}

One can observe from (3.18) that the term $- (L{p_e}/{\rho _\infty }U_\infty ^3)(\partial {u_j}/\partial {x_j})$ only depends on dimensional parameters ${T_\infty }$, ${U_\infty }$ and does not affect the feasibility of the binary scaling law.

Thus, ${Q_{t - e}}$, ${Q_{r - e}}$ and ${Q_{i - e}}$ in the electron–electronic energy conservation equation are the terms that could potentially invalidate the binary scaling law. For ${Q_{r - e}}$ and ${Q_{i - e}}$, the primary cause is still three-body reactions. For ${Q_{t - e}}$, the primary cause is the translation–electron energy exchange resulting from electron–ion collisions, ${Q_{t - e\_ions}}$.

4.1. Code validation

The flight experimental data of RAM-C II vehicle (Jones & Cross Reference Jones and Cross1972) are used for code validation. The geometry of RAM-C II vehicle and grid of the computation field are illustrated in figure 1. The detailed computational conditions are listed in table 1.

Figure 1. Geometry and grid of RAM-C II vehicle (top: grid; bottom: geometry of vehicle).

Table 1. Computational conditions for RAM-C II vehicle.

A comparison of peak electron number density along the RAM-C II vehicle surface at 71 km between CFD predicted result and flight data is plotted in figure 2. It is observed that the distribution of maximum electron number densities calculated by the Gupta model presents a good agreement with the experimental data. Although the distribution of peak electron number density predicted by the Park model results in a poorer agreement with the flight data compared with that predicted by the Gupta model, the result obtained by the Park model is no more than an order of magnitude different from the experimental data on the whole and has the same trend as the experimental data.

Figure 2. Peak electron number density along the surface of RAM-C II vehicle at 71 km.

4.2. The failure of the binary scaling law for electron distribution

The applicability of the binary scaling similarity law for thermochemical non-equilibrium air was examined using the classical case of flow around cylinders. The original cylinder which had a radius of 0.15 m and its 10-times-smaller model were used for numerical simulation. The original cylinder was subjected to a free-stream condition of 70 km and Mach 27 while the scaled-down cylinder had its free-stream density magnified 10 times to ensure that the condition of the binary scaling law $({\rho _\infty }R = \textrm{const}\textrm{.)}$ is satisfied. The detailed computational conditions are listed in table 2.

Table 2. Computational conditions for original (R = 0.15 m) and 10-times-scaled (R = 0.015 m) cylinders.

Grid independence studies of the original and scaled-down cylinders were implemented with three different grids: A, B and C. The grid quality increases from A to C by increasing both the circumferential and radial grid numbers. Detailed information of the grids is given in table 3. Grids A, B and C are constructed with the same normal size of the first cell next to surfaces. All the meshes are refined near the bow shock and Grid B with medium quality for two cylinders is illustrated in figure 3. The distributions of electron number density along the stagnation line for both the original and the scaled-down cylinders using Grids A, B and C are shown in figure 4. The results with different grid resolutions apparently demonstrate good agreement. According to the results with Grids B and C, the flow parameters no longer rely on the quality of the grid. Therefore, only the results calculated with Grid B are used for the later analysis.

Figure 3. Grid B of cylinders at 70 km and $M{a_\infty } = 27$. (a) Original cylinder of R = 0.15 m. (b) Scaled cylinder of R = 0.015 m.

Figure 4. Comparisons of distributions of electron number density along the stagnation line using different grids at $M{a_\infty } = 27$ with the Gupta model. (a) Original cylinder. (b) Scaled cylinder.

Table 3. Mesh information for cylinders at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$.

Figure 5 shows the Mach number and ${\textrm{e}^ - }$ mass fraction contours of the original and scaled cylinders at $M{a_\infty } = 27$ using the Gupta model. Although there is a good similarity in the bow-shock shape between the two cylinders, the distributions of free electrons are significantly different between the two cylinders. Figure 6 shows the distributions of T, ${T_e}$, ${T_{{v_{{\textrm{N}_2}}}}}$ and ${T_{{v_{{\textrm{O}_2}}}}}$ along the stagnation streamline at $M{a_\infty } = 27$ using the Gupta model. It is observed that the temperatures of the two flows maintain a good similarity on the whole despite the relatively small differences in local regions. The mass fractions of free electrons and ions along the stagnation line are plotted in figure 7. As shown figure 7(a,b), the gaps of maximum mass fraction between original and scaled cylinders for ${\textrm{N}^ + }$, $\textrm{N}{\textrm{O}^ + }$ and $\textrm{N}_2^ + $ are all no more than 6 %. Nevertheless, figure 7(c,d) indicates that there are tremendous differences between the original and scaled cylinders for mass fractions of ${\textrm{e}^ - }$ and ${\textrm{O}^ + }$. The difference in peak mass fraction on the stagnation line between the two cylinders reaches 50.4 % for ${\textrm{e}^ - }$ and 59.0 % for ${\textrm{O}^ + }$. As such, the binary scaling law fails completely for distribution of ${\textrm{e}^ - }$ at $M{a_\infty } = 27$ with the Gupta model.

Figure 5. Contours of (a) Mach number and (b) mass fraction of ${\textrm{e}^ - }$ at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Gupta model (top: original cylinder; bottom: scaled cylinder).

Figure 6. Distributions of temperatures along the stagnation line at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Gupta model. (a) Translational–rotational and electron–electronic temperatures and (b) ${\textrm{N}_2}$ and ${\textrm{O}_2}$ vibrational temperatures.

Figure 7. Distributions of electrons and ions along the stagnation line at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Gupta model. Mass fractions of (a) $\textrm{N}{\textrm{O}^ + }$ and $\textrm{O}_2^ + $, (b) $\textrm{N}_2^ + $ and ${\textrm{N}^ + }$, (c) ${\textrm{e}^ - }$ and (d) ${\textrm{O}^ + }$.

Figure 8 shows the distributions of ${\textrm{e}^ - }$, ${\textrm{O}^ + }$ and temperatures along the stagnation line at $M{a_\infty } = 27$ using the Park model. It is observed that the binary scaling law holds well for the distributions of ${\textrm{e}^ - }$ and ${\textrm{O}^ + }$ when the Park model is employed. Furthermore, the similarity in temperatures between the two cylinders using the Park model is greater than that using the Gupta model. However, the binary scaling law was found to be invalid for modelling electron distribution using the Park model when $M{a_\infty }$ increases to 39 (other boundary conditions remain unchanged). The distributions of ${\textrm{e}^ - }$ and ${\textrm{N}^ + }$ at $M{a_\infty } = 39$ using the Park model are illustrated in figure 9 (the detailed information of grids is listed in table 4). It is observed that the peak mass fraction of ${\textrm{e}^ - }$ and ${\textrm{N}^ + }$ in the original cylinder is obviously higher that in the scaled-down cylinder at $M{a_\infty } = 39$ using the Park model. Zeng et al. (Reference Zeng, Lin, Liu and Qu2009) also found that, at ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}$ and ${U_\infty } = 7650\;\textrm{m}\;{\textrm{s}^{ - 1}}$, the binary scaling law failed for distribution of ${\textrm{O}^ + }$ by numerical simulation of RAM-C II and its scaled model using the Gupta model.

Figure 8. Distributions of ${\textrm{e}^ - }$, ${\textrm{O}^ + }$ and temperatures along the stagnation line at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Park model. Mass fractions of (a) ${\textrm{e}^ - }$ and (b) ${\textrm{O}^ + }$. (c) Translational–rotational and electron–electronic temperatures and (d) ${\textrm{N}_2}$ and ${\textrm{O}_2}$ vibrational temperatures.

Figure 9. Distributions of ${\textrm{e}^ - }$ and ${\textrm{N}^ + }$ at $M{a_\infty } = 39$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Park model. (a) Contours of mass fraction of ${\textrm{e}^ - }$, (b) contours of mass fraction of N⁺, (c) mass fraction of ${\textrm{e}^ - }$ along the stagnation line and (d) mass fraction of ${\textrm{N}^ + }$ along the stagnation line.

Table 4. Mesh information for cylinders at $M{a_\infty } = 39$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$.

4.3. Effect of the electron-impact ionization reactions

According to the analysis in § 3.1, the three-body reaction is one of the main reasons for the invalidation of the binary scaling law. We try to find out which three-body reactions caused the failure of the binary scaling law for the electron distribution at high altitude with extremely high Mach number. The formation of free electrons is governed by seven elementary reactions in the Gupta air-chemistry model: (i) associative ionization reactions ($\textrm{N} + \textrm{O}\rightleftharpoons \textrm{N}{\textrm{O}^ + } + {\textrm{e}^ - }$, $\textrm{O} + \textrm{O}\rightleftharpoons \textrm{O}_2^ +{+} {\textrm{e}^ - }$ and $\textrm{N} + \textrm{N}\rightleftharpoons \textrm{N}_2^ +{+} {\textrm{e}^ - }$), (ii) electron-impact ionization reactions ($\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ and $\textrm{N} + {\textrm{e}^ - }\rightleftharpoons {\textrm{N}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$) and (iii) other reactions (${\textrm{O}_2} + {\textrm{N}_2}\rightleftharpoons \textrm{N}{\textrm{O}^ + } + \textrm{NO} + {\textrm{e}^ - }$ and $\textrm{NO} + \textrm{M}\rightleftharpoons \textrm{N}{\textrm{O}^ + } + {\textrm{e}^ - } + \textrm{M}$). Noticing that the shape of the distribution curve of ${\textrm{e}^ - }$ is very similar to that of ${\textrm{O}^ + }$ from figure 7(c,d), it is speculated that the dissimilarity of electron distribution might be related to $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ in the Gupta model. The dimensionless reaction rates can be defined as follows:

(4.2)\begin{equation}\overline {R_z^{fb}} ,\overline {R_z^f} ,\overline {R_z^b} = \frac{{R_z^{fb},R_z^f,R_z^b}}{{{\rho _\infty }{U_\infty }/(R{M_{mix,\infty }})}},\end{equation}

where ${M_{mix,\infty }}$ is the average molar mass of free stream and R is the radius of the cylinder. The dimensionless total reaction rates of the foregoing seven elementary reactions are plotted in figure 10. It is obvious that the reaction $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ plays the most dominant role in formation of free electrons. Figure 11 shows the dimensionless reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ on the stagnation line. As shown in figure 11(a), the dimensionless total reaction rates of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ for the original and scaled models are dissimilar and the difference of maximum reaction rates between the two reaches 26.3 %. As described in § 3.1, the precondition for the binary scaling law is that the two-body forward reaction rate is far larger than the three-body backward reaction rate, whereas both the original and scaled cylinders demonstrate a backward reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ of the same order of magnitude as the forward reaction rate, as illustrated in figure 11(b).

Figure 10. The dimensionless total rates of reactions related to generation of ${\textrm{e}^ - }$ along the stagnation line at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Gupta model.

Figure 11. The dimensionless reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ along the stagnation line at $M{a_\infty } = 27$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Gupta model. (a) Total reaction rate. (b) Forward and backward reaction rate.

In light of the fact that the binary scaling law did not fail at $M{a_\infty } = 27$ using the Park model, reaction rate coefficients of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ in the Gupta model were replaced with those in the Park model to verify the preceding assumption. Leaving the computational conditions listed in table 2 unchanged, the Gupta model was still adopted for coefficients of other reactions and the results are illustrated in figure 12. It can be observed that the binary scaling law works well for distributions of electrons, ions and temperatures. As a consequence, it can be inferred that the electron-impact ionization reaction $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ is the main reason causing the dissimilarity in flow field using the Gupta model.

Figure 12. Distributions of electrons, ions and temperatures along the stagnation line at $M{a_\infty } = 27$ with only the reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ calculated using the Park model (the Gupta model is employed for other reactions). (a) Mass fraction of ${\textrm{e}^ - }$, (b) mass fraction of ${\textrm{O}^ + }$, (c) mass fraction of $\textrm{N}_2^ + $ and ${\textrm{N}^ + }$, (d) mass fraction of $\textrm{N}{\textrm{O}^ + }$ and $\textrm{O}_2^ + $, (e) translational–rotational and electron–electronic temperatures and (f) ${\textrm{N}_2}$ and ${\textrm{O}_2}$ vibrational temperatures.

Figure 13 compares the dimensionless reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ along the stagnation streamline between the two cylinders with only reaction rate coefficients of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ calculated using the Park model. Although there are obvious differences between the two cylinders for the three-body reverse reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ from figure 13(a), the two-body forward reaction rates of that between the two cylinders have good similarity and the three-body reverse reaction rate can be negligible compared to the forward reaction rate from figure 13(b). Thus, $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ satisfies the necessary conditions of the binary scaling law in this case.

Figure 13. The dimensionless reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ along the stagnation line at $M{a_\infty } = 27$ with only reaction rate coefficients of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ calculated using the Park model (the Gupta model is employed for other reactions). (a) Backward reaction rate. (b) Forward and backward reaction rate.

The ratio of the forward reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ to the reverse reaction rate ${\gamma _{\textrm{O} + {\textrm{e}^ - }}}$ is defined as follows:

(4.3)\begin{equation}{\gamma _{\textrm{O} + {\textrm{e}^ - }}} = \frac{{\frac{{{\rho _\textrm{O}}}}{{{M_\textrm{O}}}}\frac{{{\rho _{{\textrm{e}^ - }}}}}{{{M_{{\textrm{e}^ - }}}}}k_{\textrm{O} + {\textrm{e}^ - }}^f({T_e})}}{{\frac{{{\rho _{{\textrm{O}^ + }}}}}{{{M_{{\textrm{O}^ + }}}}}\frac{{{\rho _{{\textrm{e}^ - }}}}}{{{M_{{\textrm{e}^ - }}}}}\frac{{{\rho _{{\textrm{e}^ - }}}}}{{{M_{{\textrm{e}^ - }}}}}k_{\textrm{O} + {\textrm{e}^ - }}^b({T_e})}} = \frac{{{M_{{\textrm{O}^ + }}}{M_{{\textrm{e}^ - }}}}}{{{M_\textrm{O}}\rho }}\frac{{{y_\textrm{O}}}}{{{y_{{\textrm{O}^ + }}}{y_{{\textrm{e}^ - }}}}}{K_{eq,\textrm{O} + {\textrm{e}^ - }}}({T_e}).\end{equation}

The reverse reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ calculated by the Gupta model is of the same order of magnitude as the forward reaction rate $({\gamma _{\textrm{O} + {\textrm{e}^ - }}}\sim 1)$ and the precondition for the binary scaling law $({\gamma _{\textrm{O} + {\textrm{e}^ - }}} \gg 1)$ is no longer satisfied. Here ${T_e}$ is the controlling temperature of the forward and backward reaction rates for the electron-impact ionization reactions. At $M{a_\infty } = 27$, the maximum ${T_e}$ in cases using the Gupta model or the Park model for $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ are close from figures 6(a) and 12(e). One can notice that the equilibrium constant of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ in the Gupta model is close to that in the Park model from figure 14(a). Combined with (4.3), figures 7(c,d) and 12(a,b), it can be seen that higher mass fraction of ${\textrm{O}^ + }$ and ${\textrm{e}^ - }$ is the reason that ${\gamma _{\textrm{O} + {\textrm{e}^ - }}}$ with the Gupta model is far less than that with the Park model. Figure 14(b) compares the forward reaction rate coefficient of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ from Gupta and Park models as a function of the controlling temperature according to (2.14a,b). It is clear that the forward reaction rate coefficient of the Gupta model is more than one order of magnitude larger than that of the Park model at the corresponding peak ${T_e}$. Higher forward reaction rate with the Gupta model caused the reaction to move in the positive direction and produce more ${\textrm{O}^ + }$ and ${\textrm{e}^ - }$, which increased the backward reaction rate and reduced ${\gamma _{\textrm{O} + {\textrm{e}^ - }}}$ to around 1.

Figure 14. The equilibrium constants and forward reaction rate coefficients of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ in the Gupta and Park models. (a) Equilibrium constants (particle number density equals 1.0 × 10¹⁵ cm³ in the Park model). (b) Forward reaction rate coefficients (the arrows represent the maximum electron–electronic excitation temperatures in figures 6a and 12e).

The formation of free electrons is governed by five elementary reactions in the Park model: (i) associative ionization reactions ($\textrm{N} + \textrm{O}\rightleftharpoons \textrm{N}{\textrm{O}^ + } + {\textrm{e}^ - }$, $\textrm{O} + \textrm{O}\rightleftharpoons \textrm{O}_2^ +{+} {\textrm{e}^ - }$ and $\textrm{N} + \textrm{N}\rightleftharpoons \textrm{N}_2^ +{+} {\textrm{e}^ - }$) and (ii) electron-impact ionization reactions ($\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ and $\textrm{N} + {\textrm{e}^ - }\rightleftharpoons {\textrm{N}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$). The dimensionless total reaction rates of these five elementary reactions along the stagnation line at $M{a_\infty } = 39$ using the Park model are depicted in figure 15. It can be seen that the electron-impact ionization reaction $\textrm{N} + {\textrm{e}^ - }\rightleftharpoons {\textrm{N}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ is dominant in the formation of ${\textrm{e}^ - }$. The other electron-impact ionization reaction $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ also plays an important role in electron formation near the location where a significant difference in electron mass fraction began to appear (x/R = −0.09; see figure 9c). Figure 16 shows the dimensionless reaction rates of the two electron-impact ionization reactions along the stagnation streamline at $M{a_\infty } = 39$ using the Park model. It is noticed that the backward reaction rate is of the same order of magnitude as the forward reaction rate for both reactions, which leads to the significant differences in total reaction rate between the original and scaled cylinders in figure 16(b,d).

Figure 15. The dimensionless total rates of reactions related to generation of ${\textrm{e}^ - }$ along the stagnation line at $M{a_\infty } = 39$ using the Park model.

Figure 16. The dimensionless reaction rate of two electron-impact ionization reactions along the stagnation line at $M{a_\infty } = 39$ using the Park model. (a) Forward and backward reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$. (b) Total reaction rate of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$. (c) Forward and backward reaction rate of $\textrm{N} + {\textrm{e}^ - }\rightleftharpoons {\textrm{N}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$. (d) Total reaction rate of $\textrm{N} + {\textrm{e}^ - }\rightleftharpoons {\textrm{N}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$.

Figure 17 compares the forward reaction rate coefficient of two electron-impact ionization reactions in the Gupta and Park models. Considering the foregoing results at $M{a_\infty } = 27$ and 39 with the Gupta and Park models, one can see that the binary scaling law is invalid for electron distribution when the forward reaction rate coefficient ${k^f}$ of electron-impact ionization reactions rises to a certain degree with an increase of $M{a_\infty }$. In order to further verify this, $M{a_\infty }$ is reduced to 17 (other boundary conditions remain unchanged and the detailed information of grids is listed in table 5). Figure 18 shows the distributions of ${\textrm{e}^ - }$ and ${\textrm{O}^ + }$ along the stagnation line at $M{a_\infty } = 17$ using the Gupta model. One can observe that there is good agreement in the distributions of ${\textrm{e}^ - }$ and ${\textrm{O}^ + }$ between the original and scaled cylinders (see figure 18) with a decrease in ${k^f}$ of electron-impact ionization reactions at $M{a_\infty } = 17$ using the Gupta model (see figure 17).

Figure 17. Forward reaction rate coefficient of two electron-impact ionization reactions in the Gupta and Park models (the red arrow represents the maximum electron temperature of the flow field around the cylinder under the corresponding $M{a_\infty }$ using the Gupta model and the green arrow represents that using the Park model).

Figure 18. Distributions of ${\textrm{e}^ - }$ and ${\textrm{O}^ + }$ along the stagnation line at $M{a_\infty } = 17$ with ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ using the Gupta model. Mass fractions of (a) ${\textrm{e}^ - }$ and (b) ${\textrm{O}^ + }$.

Table 5. Mesh information for cylinders at ${\rho _\infty }R = 1.3131 \times {10^{ - 5}}\;\textrm{kg}\;{\textrm{m}^{ - 2}}$ with $M{a_\infty } = 17$.

Due to effects of the electron avalanche, it is hard to measure the forward reaction rates of electron-impact ionization reactions (Park Reference Park1990; Hao et al. Reference Hao, Wang and Lee2016). However, the reverse reactions of electron-impact ionization are very nearly in quasi-steady-state distribution and the forward reaction rates can be calculated using the reverse reaction rates and equilibrium constants (Park Reference Park1990). Park measured the ionic recombination rate of $\textrm{N}$ atom ${\textrm{N}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - } \to \textrm{N} + {\textrm{e}^ - }$ by observing the rate of decrease of electron density in an expanding nozzle (Park Reference Park1990, Reference Park1968, Reference Park1973). Although the rates of ${\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - } \to \textrm{O} + {\textrm{e}^ - }$ have never been measured, the ionic recombination rate coefficients of various atoms (such as $\textrm{C}$, $\textrm{Ar}$, $\textrm{He}$) are fairly similar (Park Reference Park1990). Thus, Park assumed the recombination rates of ${\textrm{O}^ + }$ to be the same as that of ${\textrm{N}^ + }$ and obtained the rate coefficient of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$. This assumption is also embodied in the Gupta model. The rate coefficients of electron-impact ionization reactions in the Gupta model are taken from the model of Dunn and Kang (Dunn & Kang Reference Dunn and Kang1973; Gupta et al. Reference Gupta, Yos, Thompson and Lee1990). Only the three-body ionic recombination rate coefficient for ${\textrm{C}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - } \to \textrm{C} + {\textrm{e}^ - }$ was measured by Dunn (Reference Dunn1971). Hao et al. (Reference Hao, Wang and Lee2016) found an abnormal non-equilibrium relaxation process in the FIRE II case at $M{a_\infty } = 39$ simulated using the Gupta model, which proved to be a consequence of electron-impact ionization reaction rates calculated using the Gupta model. The problem was solved using modified coefficients of electron-impact ionization reactions based on the Park model. Hao considered that parameters for electron-impact ionization reactions in the Park model might be more accurate because Park's data were taken from more recent experiments. In addition, a more recent chemical kinetic model was developed by Ozawa, Zhong & Levin (Reference Ozawa, Zhong and Levin2008) and the forward rate coefficient of $\textrm{O} + {\textrm{e}^ - }\rightleftharpoons {\textrm{O}^ + } + {\textrm{e}^ - } + {\textrm{e}^ - }$ at high temperature calculated using the Ozawa model is closer to that calculated using the Park model (Niu et al. Reference Niu, Yuan, Dong and Tan2018). Although the Gupta and Park models have their own advantages in predictions of species formations (Niu et al. Reference Niu, Yuan, Dong and Tan2018), based on the points discussed above, the result at $M{a_\infty } = 27$ obtained using the Park model might be more reliable. The rates of electron-impact ionization reactions require to be further studied experimentally and theoretically because they strongly affect the applicable range of the binary scaling law for modelling electron distribution. Although the differences in forward electron-impact ionization reaction rates between the Gupta and Park models can reach as much as an order of magnitude or more (Niu et al. Reference Niu, Yuan, Dong and Tan2018), the electron-impact ionization reactions can lead to the failure of the binary scaling law in electron distribution at extremely high Mach number no matter what chemical model is employed, which reflects that it is probably a general conclusion in nature.