2.1. Kinetic description of monatomic gas mixture

We consider an $N$-components mixture of monatomic gases with the velocity distribution functions $f_s(\boldsymbol {x},\boldsymbol {v},t)$ describing their mesoscopic states, where $s$ indicates the species, $t$ is the time, $\boldsymbol {x}\in \mathbb {R}^3$ is the spatial coordinates and $\boldsymbol {v}\in \mathbb {R}^3$ is the molecular velocity. Since all the collisions, either intra-species or inter-species, are binary, the distribution function for species $s$ under external force $\boldsymbol {F}_s$ evolve according to the Boltzmann equation,

(2.1)\begin{equation} \underbrace{\frac{\partial{f_{s}}}{\partial{t}}+ \boldsymbol{v}\boldsymbol{\cdot}\frac{\partial{f_{s}}}{\partial{\boldsymbol{x}}}+ \frac{\boldsymbol{F}_s}{m_s}\boldsymbol{\cdot} \frac{\partial{f_{s}}}{\partial{\boldsymbol{v}}}}_{\mathcal{D}f_{s}} = \sum_{r=1}^N Q_{sr}(\,f_s,f_r), \quad s=1,2,\ldots,N, \end{equation}

where $m_s$ is the molecular mass; the left-hand side is known as the streaming term $\mathcal {D}f_{s}$ and the right-hand side is a sum over all binary Boltzmann collision operators $Q_{sr}$ between molecules of species $s$ and $r$,

(2.2)\begin{equation} Q_{sr}(f_s,f_r) = \int_{\mathbb{R}^3}\int_{4{\rm \pi}} |\boldsymbol{v}_{*}|\sigma_{sr}\left(|\boldsymbol{v}_{*}|, \boldsymbol{\varOmega}\boldsymbol{\cdot}\frac{\boldsymbol{v}_{*}}{|\boldsymbol{v}_{*}|} \right)[f_{s}(\boldsymbol{v}')f_{r}(\boldsymbol{w}')- f_{s}(\boldsymbol{v})f_{r}(\boldsymbol{w})]\, \mathrm{d}\boldsymbol{\varOmega}\,\mathrm{d}\boldsymbol{w}. \end{equation}

Here, $\sigma _{sr}$ is the differential scattering cross-section depending on the intermolecular potential between the two species and $\boldsymbol {\varOmega }$ is the unit vector of the solid angle; $\boldsymbol {v}$ and $\boldsymbol {w}$ are the pre-collision velocities of the two molecules of species $s$ and $r$, respectively. Hence, $\boldsymbol {v}_{*}=\boldsymbol {v}-\boldsymbol {w}$ is the relative velocity, which determines the post-collision velocities $\boldsymbol {v}'$ and $\boldsymbol {w}'$ of the collision pair as

(2.3a,b)\begin{equation} \boldsymbol{v}' = \boldsymbol{v} - \frac{2m_r}{m_s+m_r}(\boldsymbol{v}_{*} \boldsymbol{\cdot}\boldsymbol{\varOmega})\boldsymbol{\varOmega}, \quad \boldsymbol{w}' = \boldsymbol{w} + \frac{2m_s}{m_s+m_r}(\boldsymbol{v}_{*} \boldsymbol{\cdot}\boldsymbol{\varOmega})\boldsymbol{\varOmega}. \end{equation}

The macroscopic variables of each species $s$, namely, the number density $n_s$, mass density $\rho _s$, flow velocity $\boldsymbol {u_s}$, temperatures $T_s$, pressure tensor $\boldsymbol{\mathsf{P}}_s$ and heat flux $\boldsymbol {q}_s$, are obtained by taking the moments of the respective velocity distribution function $f_s$,

(2.4)\begin{equation} \left.\begin{array}{@{}c@{}} n_s = \langle 1,f_s\rangle , \quad \rho_s = \langle m_s,f_s\rangle , \quad \rho_s\boldsymbol{u}_s = \langle m_s\boldsymbol{v},f_s\rangle , \quad \tfrac{3}{2}n_sk_BT_s = \langle \tfrac{1}{2}m_s |\boldsymbol{v}-\boldsymbol{u}_s|^2,f_s\rangle , \\ \boldsymbol{\mathsf{P}}_s = \langle m_s(\boldsymbol{v}-\boldsymbol{u}_s)(\boldsymbol{v}-\boldsymbol{u}_s),f_s\rangle , \quad \boldsymbol{q}_s = \langle \tfrac{1}{2}m_s|\boldsymbol{v}-\boldsymbol{u}_s|^2(\boldsymbol{v}- \boldsymbol{u}_s),f_s\rangle, \end{array}\right\} \end{equation}

where $k_B$ is the Boltzmann constant; and the operator $\langle h,\psi \rangle$ is defined as an integral of $h\psi$ over the velocity space,

(2.5)\begin{equation} \langle h,\psi\rangle \equiv \int_{\mathbb{R}^3}{h\psi}\,\mathrm{d}\boldsymbol{v}. \end{equation}

Then, the corresponding macroscopic quantities for the mixture, the number density $n$, mass $\rho$, flow velocity $\boldsymbol {u}$, temperatures $T$, pressure tensor $\boldsymbol{\mathsf{P}}$ and heat flux $\boldsymbol {q}$, are given by

(2.6) \begin{equation} \left.\begin{array}{@{}c@{}} \displaystyle n = \sum_s\langle 1,f_s\rangle = \sum_{s}n_s, \\ \displaystyle \rho = \sum_s\langle m_s,f_s\rangle = \sum_{s}\rho_s, \\ \displaystyle \rho\boldsymbol{u} = \sum_s\langle m_s\boldsymbol{v},f_s\rangle = \sum_{s}\rho_s\boldsymbol{u}_s, \\ \displaystyle \dfrac{3}{2}nk_BT = \sum_{s}\left\langle \dfrac{1}{2}m_sc^2,f_s \right\rangle =\sum_{s}\dfrac{3}{2}n_sk_BT_s + \dfrac{1}{2}\sum_{s}\rho_s |\boldsymbol{u}_s-\boldsymbol{u}|^2, \\ \displaystyle \boldsymbol{\mathsf{P}} = \sum_{s}\langle m_s \boldsymbol{c}\boldsymbol{c},f_s\rangle = \sum_{s}\boldsymbol{\mathsf{P}}_s + \sum_{s}\rho_s (\boldsymbol{u}_s-\boldsymbol{u})(\boldsymbol{u}_s-\boldsymbol{u}), \\ \begin{aligned} \boldsymbol{q} & = \sum_{s}\left\langle \dfrac{1}{2}m_sc^2\boldsymbol{c},f_s\right\rangle = \sum_{s}\boldsymbol{q}_s + \sum_{s}\dfrac{3}{2}n_sk_BT_s(\boldsymbol{u}_s-\boldsymbol{u})\\ & \quad +\dfrac{1}{2}\sum_{s}\rho_s|\boldsymbol{u}_s-\boldsymbol{u}|^2 (\boldsymbol{u}_s-\boldsymbol{u}) + \sum_{s}\boldsymbol{\mathsf{P}}_s \boldsymbol{\cdot}(\boldsymbol{u}_s-\boldsymbol{u}), \end{aligned} \end{array}\right\} \end{equation}

where $\boldsymbol {c}=\boldsymbol {v}-\boldsymbol {u}$ is the peculiar velocity with respect to the mixture velocity $\boldsymbol {u}$, and therefore, the diffusion velocities $\boldsymbol {u}_s-\boldsymbol {u}$ contribute to the mixture temperature, stress and heat flux. Note that the scalar $c$ denotes the magnitude of the vector $\boldsymbol {c}$, with $c=|\boldsymbol {c}|$, and this notation is consistently used for both the molecular velocity $\boldsymbol {v}$ and the macroscopic velocity $\boldsymbol {u}$ throughout this paper.

2.2. Kinetic model with multi-relaxation collision operators

The proposed kinetic model for gas mixtures adopts relaxation time approximations for each pair of gas components individually to simplify the Boltzmann collision operator (2.2), and thus, the evolution of distribution functions $f_{s}$ can be written as

(2.7)\begin{equation} {\mathcal{D}f_{s}} = {\sum_{r=1}^N\frac{1}{\tau_{sr}}(g_{sr}-f_s)}, \quad s=1,2,\ldots,N, \end{equation}

where $s=r$ indicates a intra-species collision operator, $s\neq r$ are inter-species collisions; $\tau _{sr}$ is the corresponding relaxation time and $g_{sr}$ is the reference distribution function constructed in the form

(2.8) \begin{align} g_{sr}&=\hat{n}_{sr}\left(\frac{m_s}{2{\rm \pi} k_B{T}_{s}}\right)^{3/2} \exp\left(-\frac{m_s(\boldsymbol{v}-\hat{\boldsymbol{u}}_{sr})^2}{2k_B {T}_{s}}\right) \times \left[1+\frac{\hat{T}_{sr}-T_s}{T_s} \left(\frac{m_s(\boldsymbol{v}-\hat{\boldsymbol{u}}_{sr})^2}{2k_B {T}_{s}}-\frac{3}{2}\right)\right.\nonumber\\ &\quad +\left.\frac{2m_s\hat{\boldsymbol{q}}_{sr} \boldsymbol{\cdot}(\boldsymbol{v}-\hat{\boldsymbol{u}}_{sr})}{5\hat{n}_{sr} k_B^2\hat{T}_{s}^2}\left(\frac{m_s(\boldsymbol{v}- \hat{\boldsymbol{u}}_{sr})^2}{2k_B\hat{T}_{s}}- \frac{5}{2}\right)\right], \end{align}

with $\hat {n}_{sr}, \hat {T}_{sr}, \hat {\boldsymbol {u}}_{sr}, \hat {\boldsymbol {q}}_{sr}$ being the auxiliary parameters.

Construction of the auxiliary parameters is the most crucial task in building the kinetic model. As the essential constraints, the conservations of mass, momentum and energy have to be guaranteed during any binary collisions. For the intra-species collision operators, the conservations are the same as those for a single-species gas, i.e.

(2.9a–c)\begin{align} \left\langle 1,\frac{1}{\tau_{ss}}(g_{ss}-f_{s})\right\rangle = 0, \quad \left\langle m_s\boldsymbol{v},\frac{1}{\tau_{ss}}(g_{ss}-f_{s})\right\rangle = 0, \quad \left\langle \frac{1}{2}m_sv^2,\frac{1}{\tau_{ss}}(g_{ss}-f_{s})\right\rangle = 0. \end{align}

Therefore, the auxiliary number density, flow velocity and temperature in intra-species collision operators can be simply determined by the macroscopic properties of each species as $\hat {n}_{ss}=n_s$, $\hat {\boldsymbol {u}}_{ss}=\boldsymbol {u}_s$, $\hat {T}_{ss}=T_s$. Thus, the reference distribution function $g_{ss}$ in the intra-species collision term reduces to that in the Shakhov model.

The inter-species collision operators conserve the number density of each inert species, and the total momentum and energy of the collision pairs of species $s$ and $r$:

(2.10) \begin{equation} \left.\begin{array}{@{}c@{}} \displaystyle \left\langle 1,\dfrac{1}{\tau_{sr}}(g_{sr}-f_{s})\right\rangle = 0, \quad \left\langle 1,\dfrac{1}{\tau_{rs}}(g_{rs}-f_{r})\right\rangle = 0, \\ \displaystyle \left\langle m_s\boldsymbol{v},\dfrac{1}{\tau_{sr}}(g_{sr}-f_{s})\right\rangle + \left\langle m_r\boldsymbol{v},\dfrac{1}{\tau_{rs}}(g_{rs}-f_{r})\right\rangle = 0, \\ \displaystyle \left\langle \dfrac{1}{2}m_sv^2,\dfrac{1}{\tau_{sr}}(g_{sr}-f_{s})\right\rangle + \left\langle \dfrac{1}{2}m_rv^2,\dfrac{1}{\tau_{rs}}(g_{rs}-f_{r})\right\rangle = 0. \end{array}\right\} \end{equation}

Then, the auxiliary number density is obtained as $\hat {n}_{sr} = n_s$ and $\hat {n}_{rs} = n_r$, while the auxiliary velocities and temperatures $\hat {\boldsymbol {u}}_{sr}, \hat {\boldsymbol {u}}_{rs}, \hat {T}_{sr}, \hat {T}_{rs}$ cannot be uniquely determined, but yield the constraints

(2.11)\begin{equation} \left.\begin{array}{@{}c@{}} \displaystyle \dfrac{\rho_s}{\tau_{sr}}\hat{\boldsymbol{u}}_{sr} +\dfrac{\rho_r}{\tau_{rs}}\hat{\boldsymbol{u}}_{rs} = \dfrac{\rho_s}{\tau_{sr}}{\boldsymbol{u}}_{s} +\dfrac{\rho_r}{\tau_{rs}}{\boldsymbol{u}}_{r}, \\ \displaystyle\dfrac{1}{\tau_{sr}}\left[\dfrac{3}{2}n_sk_B(\hat{T}_{sr}- T_s)+\dfrac{1}{2}\rho_s(\hat{u}^2_{sr}-u^2_s)\right]\\ \quad + \dfrac{1}{\tau_{rs}}\left[\dfrac{3}{2}n_rk_B(\hat{T}_{rs}-T_r)+ \dfrac{1}{2}\rho_r(\hat{u}^2_{rs}-u^2_r)\right]=0. \end{array}\right\} \end{equation}

Therefore, we further impose the following assumptions to determine the auxiliary velocities and temperatures in inter-species collision operators:

(2.12)\begin{equation} \left.\begin{array}{@{}c@{}} \hat{\boldsymbol{u}}_{sr}-\hat{\boldsymbol{u}}_{rs} = (1-a_{sr})({\boldsymbol{u}}_{s}-{\boldsymbol{u}}_{r}) -b_{sr}(\boldsymbol{\nabla}{\ln T_s}+\boldsymbol{\nabla}{\ln T_r}), \\ \hat{T}_{sr}-\hat{T}_{rs} = (1-c_{sr})(T_{s}-T_{r}) -d_{sr}|{\boldsymbol{u}}_{s}-{\boldsymbol{u}}_{r}|^2. \end{array}\right\} \end{equation}

Here $a_{sr}=a_{rs}, b_{sr}=-b_{rs}, c_{sr}=c_{rs}, d_{sr}=-d_{rs}$ are the adjustable parameters, which describe how rapidly the equilibrium among different gas components can be achieved through inter-species collisions. To be specific, the assumptions imposed to determine the auxiliary velocities and temperatures are designed to recover the correct momentum and energy relaxation, respectively. Since the momentum relaxation is closely related to the diffusion phenomenon, the term with parameters $a_{sr}$ is presented to account for the ordinary diffusion arising from concentration gradients, that is, the Fick effect. The term with parameters $b_{sr}$ is introduced to encapsulate the thermal diffusion induced by temperature gradients, known as the Soret effect. Note that, compared with those in the literature (such as Haack et al. Reference Haack, Hauck and Murillo2017; Bobylev et al. Reference Bobylev, Bisi, Groppi, Spiga and Potapenko2018), (2.12) gives a more general expression for the relations between the auxiliary and macroscopic properties, by adding the temperature gradient to phenomenologically model the thermally induced flow. Then, combined with (2.11), the auxiliary velocities and temperatures are given by

(2.13)\begin{equation} \left.\begin{array}{@{}c@{}} \displaystyle \hat{\boldsymbol{u}}_{sr} ={\boldsymbol{u}}_{s}- \dfrac{\rho_r\tau_{sr}}{\rho_s\tau_{rs}+\rho_r\tau_{sr}}\boldsymbol{{X}}_{sr},\\ \displaystyle \hat{T}_{sr} = T_s-\dfrac{n_r\tau_{sr}}{n_s\tau_{rs}+n_r \tau_{sr}}{Y}_{sr}-\dfrac{\rho_s\rho_r\tau_{sr}\tau_{rs}\boldsymbol{{X}}_{sr} \boldsymbol{\cdot}[\boldsymbol{{X}}_{sr}-2({\boldsymbol{u}}_{s}- {\boldsymbol{u}}_{r})]}{3k_B(n_s\tau_{rs}+ n_r\tau_{sr})(\rho_s\tau_{rs}+\rho_r\tau_{sr})}, \end{array}\right\} \end{equation}

with

(2.14)\begin{equation} \left.\begin{array}{@{}c@{}} \boldsymbol{{X}}_{sr} = a_{sr}({\boldsymbol{u}}_{s}-{\boldsymbol{u}}_{r}) +b_{sr}(\boldsymbol{\nabla}{\ln T_s}+\boldsymbol{\nabla}{\ln T_r}), \\ {Y}_{sr} = c_{sr}(T_{s}-T_{r})+d_{sr} |{\boldsymbol{u}}_{s}-{\boldsymbol{u}}_{r}|^2. \end{array}\right\} \end{equation}

In addition, the auxiliary properties $\hat {\boldsymbol {q}}_{sr}$ are constructed to adjust the relaxation rates of heat fluxes as

(2.15)\begin{equation} \hat{\boldsymbol{q}}_{sr} = (1-{Pr}_{sr}){\boldsymbol{q}}_{s}+\gamma_{sr} ({\boldsymbol{q}}_{sr}-{\boldsymbol{q}}_{s}), \end{equation}

where ${Pr}_{sr}$ is an effective Prandtl number giving the thermal relaxation of species $s$ due to collisions with species $r$, which reduces to ${Pr}_{ss}=2/3$ for the intra-species collisions of a monatomic gas; ${\boldsymbol {q}}_{sr}$ is defined as the heat flux of species $s$ measured relative to auxiliary velocity $\hat {\boldsymbol {u}}_{sr}$,

(2.16)\begin{equation} {\boldsymbol{q}}_{sr} = \langle \tfrac{1}{2}m_s |\boldsymbol{v}-\hat{\boldsymbol{u}}_{sr}|^2(\boldsymbol{v}- \hat{\boldsymbol{u}}_{sr}),f_s\rangle , \end{equation}

and thus, (2.15) yields $\hat {\boldsymbol {q}}_{ss} = (1-{Pr}_{ss}){\boldsymbol {q}}_{s}$ for the intra-species collision term; $\gamma _{sr}=-\gamma _{rs}$ is a dimensionless coefficient taking into account the Dufour effects caused by the diffusive thermal conductivity. Based on the asymptotic analysis of the Boltzmann equation (Chapman & Cowling Reference Chapman and Cowling1970), it is an effect inverse to thermal diffusion and, hence, $\gamma _{sr}$ is not an independent parameter but can be determined by the other adjustable parameters, as shown in the following section.

2.3. Continuum flow limit

The multi-relaxation kinetic model equations have been constructed for a monatomic gas mixture, which contains several adjustable parameters: relaxation rates $\tau _{sr}$, effective Prandtl number ${Pr}_{sr}$ and the set of coefficients $a_{sr}$, $b_{sr}$, $c_{sr}$, $d_{sr}$, $\gamma _{sr}$ used for calculating auxiliary properties. Here, we perform the Chapman–Enskog analysis to get the asymptotic limit of the proposed model equation and determine the relevant parameters to recover the transport coefficients, so that the kinetic model and macroscopic fluid dynamics are consistent in the continuum flow regime, which is the basic requirement to achieve accurate kinetic modelling.

Without loss of generality, we consider a binary mixture of monatomic gases. The extension to a multi-species mixture is expected to be relatively easy, due to the fact that our kinetic model is built on each binary collision pair. Therefore, although the expressions of the transport coefficients of the whole mixture become more complex for mixtures with more than two components, only the properties of binary mixtures are required to determine all the model parameters in each possible binary collision term.

In the continuum limit, when all the relaxation times are considerably smaller than the characteristic time of the gas flow, only the species number density $n_s$, total momentum $\rho \boldsymbol {u}$ and energy $\frac {3}{2}nk_BT+\frac {1}{2}\rho u^2$ are the collisional invariants. Then, the set of equations for the conserved macroscopic variables $n_s, \boldsymbol {u}, T$ can be obtained by taking momentums of the kinetic equations (2.7) and summing over the species ($s=1,2$) for momentum and energy equations,

(2.17) \begin{equation} \left.\begin{array}{@{}c@{}} \displaystyle \dfrac{\partial{n_s}}{\partial{t}} + \boldsymbol{\nabla}\boldsymbol{\cdot} (n_s\boldsymbol{u}) + \boldsymbol{\nabla}\boldsymbol{\cdot} (n_s(\boldsymbol{u}_s-\boldsymbol{u})) = 0, \\ \displaystyle \dfrac{\partial}{\partial{t}}(\rho\boldsymbol{u}) + \boldsymbol{\nabla}\boldsymbol{\cdot}(\rho\boldsymbol{u}\boldsymbol{u}) + \boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{\mathsf{P}} = n_1\boldsymbol{F}_1+n_2\boldsymbol{F}_2, \\ \begin{aligned} & \dfrac{\partial}{\partial{t}} \left(\dfrac{3}{2}nk_BT+\dfrac{1}{2}\rho u^2\right) +\boldsymbol{\nabla}\boldsymbol{\cdot}\left[\left(\dfrac{3}{2}nk_BT+\dfrac{1}{2}\rho u^2\right)\boldsymbol{u}\right] \\ & \quad +\boldsymbol{\nabla}\boldsymbol{\cdot}(\boldsymbol{\mathsf{P}}\boldsymbol{\cdot}\boldsymbol{u}) + \boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{q} = n_1\boldsymbol{F}_1\boldsymbol{\cdot}\boldsymbol{u}_1+n_2\boldsymbol{F}_2\boldsymbol{\cdot}\boldsymbol{u}_2, \end{aligned} \end{array}\right\} \end{equation}

where $\boldsymbol {u}_s$ is the species macroscopic velocity defined in (2.4), $\boldsymbol{\mathsf{P}}$ and $\boldsymbol {q}$ are the mixture pressure tensor and heat flux, respectively, given in (2.6). To close the above set of equations, it is necessary to find the approximations to distribution functions $f_s$, and hence, the properties $\boldsymbol {u}_s, \boldsymbol{\mathsf{P}}, \boldsymbol {q}$ can be expressed as functions of the macroscopic variables $n_s, \boldsymbol {u}, T$, giving the constitutive relations and transport properties of the gas mixture (see the details in Appendix A).

The Navier–Stokes-type equation for a gas mixture is then obtained by the second approximation of $f_s$, and the species velocities $\boldsymbol {u}_s^{{NS}}$, stress tensor $\boldsymbol{\mathsf{P}}^{{NS}}$ and heat flux $\boldsymbol {q}^{{NS}}$ are given by

(2.18) \begin{align} \left.\begin{array}{@{}c@{}} \displaystyle \boldsymbol{u}_1^{{NS}} =\boldsymbol{u} -\dfrac{\rho_1\tau_{21}+\rho_2\tau_{12}}{a_{12}} \dfrac{p}{\rho\rho_1}\boldsymbol{d}_{12}- \dfrac{2b_{12}\rho_2}{a_{12}\rho}\boldsymbol{\nabla}{\ln T}, \\ \displaystyle \boldsymbol{u}_2^{{NS}} = \boldsymbol{u}+ \dfrac{\rho_1\tau_{21}+\rho_2\tau_{12}}{a_{12}} \dfrac{p}{\rho\rho_2}\boldsymbol{d}_{12}+\dfrac{2b_{12}\rho_1}{a_{12} \rho}\boldsymbol{\nabla}{\ln T}, \\ \displaystyle \boldsymbol{\mathsf{P}}^{{NS}} = nk_BT\boldsymbol{\mathsf{I}} - k_BT\left(\dfrac{\tau_{11}\tau_{12}}{\tau_{11}+\tau_{12}}n_1+ \dfrac{\tau_{22}\tau_{21}}{\tau_{22}+\tau_{21}}n_2\right) \left(\boldsymbol{\nabla}\boldsymbol{u}+\boldsymbol{\nabla}\boldsymbol{u}^{\mathrm{T}}- \dfrac{2}{3}\boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{u}\boldsymbol{\mathsf{I}}\right) , \\ \begin{aligned} \boldsymbol{q}^{{NS}} & = \dfrac{5}{2}k_BT(n_1\boldsymbol{u}_1^{{NS}} +n_2\boldsymbol{u}_2^{{NS}}-n\boldsymbol{u}) + \dfrac{2b_{12}\rho_1\rho_2}{\rho_1\tau_{21}+\rho_2\tau_{12}} ({\boldsymbol{u}}_{1}^{{NS}}-{\boldsymbol{u}}_{2}^{{NS}}) \\ & \quad -\left[\left(\dfrac{n_1}{m_1}\dfrac{\tau_{11}\tau_{12}}{{Pr}_{12} \tau_{11}+{Pr}_{11}\tau_{12}}+\dfrac{n_2}{m_2}\dfrac{\tau_{22} \tau_{21}}{{Pr}_{21}\tau_{22}+{Pr}_{22}\tau_{21}}\right) \dfrac{5k_B^2T}{2} \right.\\ & \quad - \left.\dfrac{4b_{12}^2\rho_1\rho_2}{a_{12} (\rho_1\tau_{21}+\rho_2\tau_{12})T}\right]\boldsymbol{\nabla} T, \end{aligned} \end{array}\right\} \end{align}

where $p=nk_BT$ is the pressure, $\boldsymbol{\mathsf{I}}$ is a $3\times 3$ identity matrix and $\boldsymbol {d}_{12}$ is given as

(2.19)\begin{equation} \boldsymbol{d}_{12} = \boldsymbol{\nabla} \left(\frac{n_1}{n}\right) +\frac{n_1n_2(m_2-m_1)}{n\rho}\boldsymbol{\nabla}{\ln p} -\frac{\rho_1\rho_2}{\rho p}\left(\frac{\boldsymbol{F}_1}{m_1}-\frac{\boldsymbol{F}_2}{m_2}\right), \end{equation}

which yields $\boldsymbol {d}_{12}=-\boldsymbol {d}_{21}$, since $\boldsymbol {\nabla } ({n_1}/{n}) =-\boldsymbol {\nabla } ({n_2}/{n})$.

2.3.1. Viscosity

The shear viscosity $\mu$ of a binary mixture can be obtained from the off-diagonal components of the non-equilibrium stress tensor $\boldsymbol{\mathsf{P}}$ in (2.18), which leads to

(2.20)\begin{equation} \mu = k_BT\left(\frac{\tau_{11}\tau_{12}}{\tau_{11}+ \tau_{12}}n_1+\frac{\tau_{22}\tau_{21}}{\tau_{22}+\tau_{21}}n_2\right). \end{equation}

Clearly, the mixture shear stress depends on all the relaxation times. It is known that the dependence of species shear viscosity $\mu _s$ on mean collision time $\tau _{ss}$ of intra-species collisions is $\mu _s=n_sk_BT\tau _{ss}$. However, the relaxation time $\tau _{sr} (s\neq r)$ in the kinetic model is no longer a mean molecular collision time of inter-species interactions, but measures the time scale approaching the reference states $g_{sr}$ for the component $s$ due to the collisions with component $r$. In other words, $n_r\tau _{sr}$ should be smaller than $n_s\tau _{rs}$ when $m_s< m_r$, because the collisions between molecules with disparate mass have a more significant influence on the light one than on the heavy one. For example, by considering a mixture consisting of ions and electrons with mass $m_i$ and $m_e$, respectively, the approximated relation between the inter-species relaxation times given in the literature is $n_i\tau _{ei}/n_e\tau _{ie}=m_e/m_i\ll 1$ (Bellan Reference Bellan2006).

Therefore, here we define

(2.21)\begin{equation} \phi_{sr}=\frac{n_s\tau_{ss}}{n_r\tau_{sr}}, \end{equation}

to quantify the ratio between the relaxation times of each intra- and inter-species collision for species $s$, where $\phi _{sr}$ is only determined by the intermolecular interactions and temperature, but independent of the concentrations of the components. Then the mixture viscosity (2.20) can be rewritten in the form

(2.22)\begin{equation} \mu = \frac{\mu_1}{1+\dfrac{n_2}{n_1}\phi_{12}} + \frac{\mu_2}{1+\dfrac{n_1}{n_2}\phi_{21}}, \end{equation}

which is a linear combination of the species viscosities of the mixture components and shares the exact same form as that given by Wilke's mixture rule (Wilke Reference Wilke1950). Nevertheless, the more accurate values of $\phi _{sr}$ can be obtained by fitting the mixture viscosity measured experimentally.

2.3.2. Diffusion

In a gas mixture the ordinary diffusion coefficient $D_{sr}$ is usually measured when the gas mixture is uniform in temperature and pressure and without external forces acting on the gas molecules; while the thermal diffusion with coefficient $D_{T,sr}$ also contributes to the diffusion velocity when there is a temperature gradient present. In general, the mass flux $J_s$ of the species $s$ in a binary mixture caused by gradients of concentration, pressure and temperature is given by (Hirschfelder, Curtiss & Bird Reference Hirschfelder, Curtiss and Bird1954)

(2.23)\begin{equation} J_s =\rho_s(\boldsymbol{u}_s-\boldsymbol{u}) ={-}\frac{n^2}{\rho}m_sm_r[D_{sr}\boldsymbol{d}_{sr}-D_{T,sr}\boldsymbol{\nabla}{\ln T}]. \end{equation}

Meanwhile, the diffusion velocity in the proposed kinetic model is obtained as (2.18) from the Chapman–Enskog method,

(2.24)\begin{equation} \boldsymbol{u}_s-\boldsymbol{u} ={-}\frac{\rho_1\tau_{21}+ \rho_2\tau_{12}}{a_{sr}}\frac{p}{\rho_1\rho_2}\boldsymbol{d}_{sr}- \frac{2b_{sr}}{a_{sr}}\boldsymbol{\nabla}{\ln T}. \end{equation}

Therefore, the binary diffusion coefficient $D_{12}$ is

(2.25)\begin{equation} D_{12} = \frac{k_BT}{m_1m_2n}\frac{\rho_1\tau_{21}+\rho_2\tau_{12}}{a_{12}}, \end{equation}

which determines the parameter $a_{12}=a_{21}$ when the relaxation times are known from viscosity. Also, the thermal-diffusion coefficient $D_{T,12}$, as well as the thermal-diffusion ratio $k_{T,12}=D_{T,12}/D_{12}$ are obtained:

(2.26a,b)\begin{equation} D_{T,12} = \frac{2b_{12}n_1n_2}{a_{12}n^2}, \quad k_{T,12} =\frac{2b_{12}\rho_1\rho_2}{p(\rho_1 \tau_{21}+\rho_2\tau_{12})}. \end{equation}

Thus, the parameter $b_{12}=-b_{21}$ is adjusted to match the thermal-diffusion ratio and recover the Soret effect in the continuum limit. The Soret effect describes mass separation due to a temperature gradient and it may be positive or negative depending on the mass difference and intermolecular potentials (Chapman & Cowling Reference Chapman and Cowling1970). In general, thermal diffusion becomes stronger in molecules with a larger mass difference and, thus, is crucial in the modelling of gas mixtures with significant mass disparity. Meanwhile, although the absolute value of the thermal-diffusion ratio $k_{T,12}$ usually has an order of magnitude less than $10^{-1}$ for a neutral gas mixture, it can be greatly increased in ionized gases (Chapman Reference Chapman1958).

2.3.3. Thermal conductivity

In addition to the direct transfer of kinetic energy during collisions, the diffusional migration of molecules also carries thermal energy and contributes to the heat transport in a gas mixture, which is known as the Dufour effect. Based on the Chapman–Enskog method, the mixture heat flux $\boldsymbol {q}$ is obtained as

(2.27)\begin{align} \boldsymbol{q} &= \frac{5}{2}k_BT[n_1(\boldsymbol{u}_1- \boldsymbol{u})+n_2(\boldsymbol{u}_2-\boldsymbol{u})] + \frac{5}{2}k_BTA({\boldsymbol{u}}_{1}-{\boldsymbol{u}}_{2})\nonumber\\ &\quad -\frac{5}{2}k_BT\left(\frac{n_1}{m_1} \frac{\tau_{11}\tau_{12}}{{Pr}_{12}\tau_{11}+{Pr}_{11}\tau_{12}}+ \frac{n_2}{m_2}\frac{\tau_{22}\tau_{21}}{{Pr}_{21}\tau_{22}+{Pr}_{22}\tau_{21}} - A\frac{2b_{12}}{a_{12}k_BT}\right)k_B\boldsymbol{\nabla} T. \end{align}

with

(2.28)\begin{equation} A = \gamma_{12} \left(\frac{1}{m_1}\frac{\tau_{11} \tau_{12}}{{Pr}_{12}\tau_{11}+{Pr}_{11}\tau_{12}}+ \frac{1}{m_2}\frac{\tau_{22}\tau_{21}}{{Pr}_{21}\tau_{22} +{Pr}_{22}\tau_{21}}\right)\frac{a_{12}\rho_1\rho_2}{\rho_1 \tau_{21}+\rho_2\tau_{12}}. \end{equation}

It can be seen that there are three factors contributing to the total heat flux. (i) The first term occurs since the heat flux is measured relative to the mixture flow velocity, instead of the species flow velocity, and thus, represents energy carried by the molecular flux in the presence of diffusion $\boldsymbol {u}_s-\boldsymbol {u}$. (ii) The second term with binary diffusion velocity $\boldsymbol {u}_1-\boldsymbol {u}_2$ arises as an inverse process to thermal diffusion, and has a coefficient $\frac {5}{2}k_BTA=nk_BTk_{T,12}$ based on the asymptotic analysis of the original Boltzmann equation (Chapman & Cowling Reference Chapman and Cowling1970), leading to

(2.29)\begin{equation} \gamma_{12} =\left(\frac{1}{m_1}\frac{\tau_{11}\tau_{12}}{{Pr}_{12} \tau_{11}+{Pr}_{11}\tau_{12}}+\frac{1}{m_2}\frac{\tau_{22} \tau_{21}}{{Pr}_{21}\tau_{22}+{Pr}_{22}\tau_{21}}\right)^{{-}1} \frac{4b_{12}}{5a_{12}k_BT}. \end{equation}

(iii) The third term $-\kappa \boldsymbol {\nabla } T$ is generated by a temperature gradient, where $\kappa$ is the ordinary thermal conductivity of the mixture that is usually measured experimentally in the absence of any diffusion velocity:

(2.30)\begin{align} \kappa = \left(\frac{n_1}{m_1}\frac{\tau_{11}\tau_{12}}{{Pr}_{12} \tau_{11}+{Pr}_{11}\tau_{12}}+\frac{n_2}{m_2}\frac{\tau_{22} \tau_{21}}{{Pr}_{21}\tau_{22}+{Pr}_{22}\tau_{21}}\right) \frac{5k_B^2T}{2} - \frac{4b_{12}^2\rho_1\rho_2}{a_{12} (\rho_1\tau_{21}+\rho_2\tau_{12})T}. \end{align}

In analogy with the parameter $\phi _{sr}$ used to measure the ratio between the relaxation times of intra- and inter-species collisions, we also define

(2.31)\begin{equation} \varphi_{sr}=\frac{{Pr}_{sr}}{{Pr}_{ss}}, \end{equation}

to represent the ratio of the thermal relaxation rates. Given the species thermal conductivity $\kappa _s ={5n_sk_B^2T\tau _{ss}}/{2m_s{Pr}_{ss}}$, the mixture thermal conductivity (2.30) can be rewritten in the form

(2.32)\begin{equation} \kappa =\frac{\kappa_1}{1+\dfrac{n_2}{n_1}\phi_{12}\varphi_{12}} +\frac{\kappa_2}{1+\dfrac{n_1}{n_2}\phi_{21}\varphi_{21}} - \frac{D_{12}k_T^2n^3k_B}{n_1n_2}, \end{equation}

where the last term indicates the effect of thermal diffusion on the thermal conductance of a gas mixture. Similar to the viscosity, the values of $\varphi _{sr}$ can be determined by matching the mixture thermal conductivity measured at different proportions of gas components.

2.4. Inter-species energy relaxation

We have shown that in the continuum flow limit the temperatures of different components stay the same up to the Navier–Stokes approximation of the proposed model, when all the relaxation times are considerably smaller than the characteristic time of gas flow. Thus, the transport coefficients in the continuum flow limit are not affected by the parameters $c_{sr}$ and $d_{sr}$. However, these parameters determine the auxiliary temperatures and, hence, the energy relaxation rates between different species, which may have a significant impact in strong non-equilibrium situations. In the kinetic modelling of a single-species gas, we have found that having all the transport coefficients is not enough to exactly describe the underlying relaxation processes (Li et al. Reference Li, Zeng, Su and Wu2021, Reference Li, Zeng, Huang and Wu2023; Zeng, Li & Wu Reference Zeng, Li and Wu2022). The same is true for the energy relaxation during inter-species collisions of a gas mixture.

Therefore, we determine the parameters $c_{sr}$ and $d_{sr}$ for calculating auxiliary temperature by imposing that the energy exchange rates of the inter-species collision operator coincide with that of the Boltzmann collision operator,

(2.33)\begin{equation} \left\langle\frac{1}{2}m_sv^2,\frac{1}{\tau_{sr}}(g_{sr}-f_{s})\right\rangle = \left\langle \frac{1}{2}m_sv^2,Q_{sr}\right\rangle . \end{equation}

The exchange rates of the kinetic model (left-hand side) can be calculated straightforwardly, while that of the Boltzmann collision operator may only be explicitly evaluated for the Maxwellian intermolecular potential,

(2.34)\begin{align} &\frac{n_s}{\tau_{sr}}\left[\frac{3}{2}k_B(\hat{T}_{sr}-T_s)+ \frac{m_s}{2}(\hat{\boldsymbol{u}}_{sr}-\boldsymbol{u}_s)^2\right]\nonumber\\ &\quad =\lambda_{sr}\frac{m_sm_r}{(m_s+m_r)^2}n_sn_r [3k_B(T_r-T_s)+m_r(\boldsymbol{u}_r-\boldsymbol{u}_s)^2], \end{align}

where $\lambda _{sr}=\lambda _{rs}$ are constants related to the integral of collision cross-sections. By substituting the auxiliary velocity and temperature (2.13) into (2.34), we immediately get

(2.35)\begin{equation} \left.\begin{array}{@{}c@{}} \displaystyle c_{sr} =\dfrac{2\lambda_{sr}m_sm_r}{(m_s+m_r)^2}(n_s\tau_{rs}+n_r\tau_{sr}),\\ \displaystyle d_{sr} =\dfrac{\lambda_{sr}m_sm_r}{3k_B(m_s+m_r)^2}[{\lambda_{sr} (n_r\rho_r\tau_{sr}^2-n_s\rho_s\tau_{rs}^2)} -{2(\rho_r\tau_{sr}-\rho_s\tau_{rs})}]. \end{array}\right\} \end{equation}

Note that the same form of (2.34) can be obtained for non-Maxwellian intermolecular potentials with non-constant $\lambda _{sr}$, when it is subject to the restriction that the distribution functions are Maxwellian at different temperatures but with small diffusion velocities (Morse Reference Morse1964). Therefore, we calculate the parameters $c_{sr}$ and $d_{sr}$ using (2.35) for any type of intermolecular potential.

Furthermore, variable $\lambda _{sr}$ can be approximated by matching the momentum exchange rates of the inter-species collision operator with that of the Boltzmann collision operator for the Maxwellian intermolecular potential,

(2.36)\begin{equation} \left\langle m_sv,\frac{1}{\tau_{sr}}(g_{sr}-f_{s})\right\rangle = \langle m_sv,Q_{sr}\rangle , \end{equation}

which leads to

(2.37)\begin{equation} \frac{n_sm_s}{\tau_{sr}}(\hat{\boldsymbol{u}}_{sr}-{\boldsymbol{u}}_{s}) =\lambda_{sr}\frac{m_sm_r}{m_s+m_r}n_sn_r({\boldsymbol{u}}_{r}- {\boldsymbol{u}}_{s}). \end{equation}

Therefore, $\lambda _{sr}$ can be evaluated by parameter $a_{sr}$ as

(2.38)\begin{equation} \lambda_{sr} = \frac{a_{sr}(m_s+m_r)}{\rho_s\tau_{rs}+\rho_r\tau_{sr}}, \end{equation}

which is thus determined by the binary diffusion coefficient as given by (2.25).

2.5. Indifferentiability principle

The indifferentiability principle states that, when all the molecules are mechanically identical (e.g. they have the same mass and scattering cross-section), the model equation reduces to a single one by adding the distribution functions (Garzó, Santos & Brey Reference Garzó, Santos and Brey1989), although in reality one does not have this situation unless one uses two Boltzmann equations to simulate single-species flow. This property holds for the Boltzmann equation because of the bilinearity of its operators. It is however non-trivial for a model equation to inherit since the operators constructed are usually highly nonlinear. Historically, several kinetic models using a single-relaxation collision operator have been proven to fulfil this principle. It should be noted that some of them require a condition that the diffusion velocities vanish for the indifferentiable molecules, when the models contain parameters to recover the Fick law (Brull, Pavan & Schneider Reference Brull, Pavan and Schneider2012; Todorova & Steijl Reference Todorova and Steijl2019).

Therefore, we also adopt the assumption that $\boldsymbol {u}_s=\boldsymbol {u}$ for all the indifferentiable species $s$, to demonstrate that our multi-relaxation model with the linearised collision operator complies with the indifferentiability principle. Consider a gas mixture system that slightly deviates from an equilibrium state with flow velocity $\boldsymbol {u}$, temperature $T$ and number density $n_s$ of each component, the reference distribution in the linearised collision operator is given as

(2.39) \begin{align} g_{sr}^{linear} &= f_{s}^{eq}\left[1 + \frac{m_s (\hat{\boldsymbol{u}}_{sr}-\boldsymbol{u}) \boldsymbol{\cdot}\boldsymbol{c}}{k_BT} + \frac{\hat T_{sr}-T}{T} \left(\frac{m_sc^2}{2k_B{T}}-\frac{3}{2}\right)\right.\nonumber\\ &\quad +\left. \frac{2m_s\hat{\boldsymbol{q}}_{sr}\boldsymbol{\cdot}\boldsymbol{c}}{5{n}_{s} k_B^2{T}^2}\left(\frac{m_sc^2}{2k_B{T}}-\frac{5}{2}\right)\right], \end{align}

with

(2.40)\begin{equation} f_s^{eq} = n_s\left(\frac{m_s}{2{\rm \pi} k_BT}\right)^{3/2} \exp{\left(-\frac{m_s(\boldsymbol{v}-{\boldsymbol{u}})^2}{2k_BT}\right)}. \end{equation}

For the indifferentiable molecules, we have (i) $m_s=m$ and $n_r\tau _{sr}$ is constant for any $s$ and $r$ due to the identical mass and scattering cross-sections, respectively; (ii) $\hat {\boldsymbol {u}}_{sr}=\boldsymbol {u}$ based on (2.13) with the thermal-diffusion coefficient vanishing; (iii) $\hat T_{sr}+\hat T_{rs}=T_s+T_r$ and $nT=\sum n_sT_s$ from total energy conservation (2.11) and calculation of the mixture temperature (2.6), respectively; (iv) $\hat {\boldsymbol {q}}_{sr}=(1-{Pr})\boldsymbol {q}_{s}$ from (2.15). Then, the sum of the kinetic equations over all species yields

(2.41)\begin{align} {\mathcal{D}\left(\sum_{s=1}^Nf_{s}\right)} &= \sum_{s=1}^N{\sum_{r=1}^N\frac{1}{\tau_{sr}} (g_{sr}^{linear}-f_s)} \nonumber\\ &= \frac{1}{\tau}\left\{\sum_sf_s^{eq}\left[1+\frac{2m(1 -{Pr}){\left(\sum_s\boldsymbol{q}_{s}\right)}\boldsymbol{\cdot} \boldsymbol{c}}{5{n}k_B^2{T}^2}\left(\frac{mc^2}{2k_B{T}}- \frac{5}{2}\right)\right]-\sum_{s=1}^Nf_{s}\right\}, \end{align}

where $\tau =n_r\tau _{sr}/n$ is the overall relaxation time. Clearly, the kinetic model equation reduces to the Shakhov model of single-species monatomic gas with the distribution function $f=\sum f_{s}$.

2.6. Dimensionless forms

Let $L_0, T_0, n_0, m_0$ be the reference length, temperature, number density and mass, respectively, then the most probable speed is

(2.42)\begin{equation} v_m=\sqrt{\frac{2k_BT_0}{m_0}}, \end{equation}

and the reference pressure is $p_0=n_0k_BT_0$. The dimensionless variables are introduced as

(2.43)\begin{equation} \left.\begin{array}{@{}c@{}} \tilde{\boldsymbol{x}}=\boldsymbol{x}/L_0, \quad \tilde{n}=n/n_0, \quad \tilde{m}=m/m_0, \quad\tilde{T}=T/T_0, \\ \tilde{\boldsymbol{v}}=\boldsymbol{v}/v_m, \quad \tilde{\boldsymbol{c}}=\boldsymbol{c}/v_m, \quad \tilde{t}=v_mt/L_0, \quad \tilde{\tau}=v_m\tau/L_0, \\ \tilde{p}=p/p_0, \quad \tilde{\boldsymbol{q}}=\boldsymbol{q}/(p_0v_m), \quad \tilde{f}_s=v_m^{3}f_s/n_0. \end{array}\right\} \end{equation}

The Knudsen numbers ${Kn}_{s}$ of each species $s$ is defined as

(2.44)\begin{equation} {Kn}_{s}=\frac{\mu_s(T_0)}{n_0L_0}\sqrt{\frac{\rm \pi}{2m_sk_BT_0}}. \end{equation}

It is noted that the species-specific Knudsen numbers are correlated as ${Kn}_{r}={Kn}_{s}\beta ^{\mu }_{rs}/\sqrt {\beta ^m_{rs}}$, with $\beta ^{\mu }_{rs}=\mu _r(T_0)/\mu _s(T_0)$ being the viscosity ratio at the reference temperature and $\beta ^m_{rs}=m_r/m_s$ the mass ratio. Therefore, the dimensionless relaxation times can be written in terms of the Knudsen numbers,

(2.45a,b)\begin{equation} \tilde{\tau}_{ss} = 2{Kn}_{s}\sqrt{\frac{\tilde{m}_s}{\rm \pi}} \frac{\tilde{T}_s^{\omega_s-1}}{\tilde{n}_s}, \quad \tilde{\tau}_{sr} = \tilde{\tau}_{ss}\phi_{sr}^{{-}1} \frac{\tilde{n}_s}{\tilde{n}_r}, \end{equation}

where $\omega _s$ is the viscosity index of species $s$ in

(2.46)\begin{equation} \mu_s(T)=\mu_s(T_0)\left(\frac{T}{T_0}\right)^\omega_s. \end{equation}

Then, the kinetic model equations are non-dimensionalised as

(2.47)\begin{equation} \frac{\partial{\tilde{f}_{s}}}{\partial{\tilde{t}}}+\tilde{\boldsymbol{v}} \boldsymbol{\cdot} \frac{\partial{\tilde{f}_{s}}}{\partial{\tilde{\boldsymbol{x}}}}+ \frac{\tilde{\boldsymbol{F}}_s}{\tilde{m}_s}\boldsymbol{\cdot} \frac{\partial{\tilde{f}_{s}}}{\partial{\tilde{\boldsymbol{v}}}} =\sqrt{\frac{\rm \pi}{\tilde{m}_s}}\frac{\tilde{T}_s^{1-\omega_s}}{2{Kn}_{s}} \left[\tilde{n}_s(\tilde{g}_{ss}-\tilde{f}_s)+\sum_{r \neq s}\tilde{n}_r\phi_{sr}(\tilde{g}_{sr}-\tilde{f}_s)\right], \end{equation}

with the dimensionless reference velocity distribution function

(2.48) \begin{align} \tilde{g}_{sr}&=\widetilde{{n}}_{sr}\left(\frac{\tilde{m}_s}{{\rm \pi} \tilde{T}_{s}}\right)^{3/2}\exp\left(-\frac{\tilde{m}_s(\tilde{\boldsymbol{v}}- \tilde{\hat{\boldsymbol{u}}}_{sr})^2}{\tilde{T}_{s}}\right) \times \left[1+\frac{\tilde{\hat{T}}_{sr}-\widetilde{T_s}}{\widetilde{T_s}} \left(\frac{\tilde{m}_s(\tilde{\boldsymbol{v}}- \tilde{\hat{\boldsymbol{u}}}_{sr})^2}{\tilde{T}_{s}}- \frac{3}{2}\right)\right.\nonumber\\ &\quad +\left.\vphantom{1+\frac{\tilde{\hat{T}}_{sr}-\widetilde{T_s}}{\widetilde{T_s}} \left(\frac{\tilde{m}_s(\tilde{\boldsymbol{v}}- \tilde{\hat{\boldsymbol{u}}}_{sr})^2}{\tilde{T}_{s}}- \frac{3}{2}\right)}\frac{4\tilde{m}_s\tilde{\hat{\boldsymbol{q}}}_{sr} \boldsymbol{\cdot}(\tilde{\boldsymbol{v}}-\tilde{\hat{\boldsymbol{u}}}_{sr})}{5 \tilde{\hat{n}}_{sr}\tilde{T}_{s}^2}\left(\frac{\tilde{m}_s (\tilde{\boldsymbol{v}}-\tilde{\hat{\boldsymbol{u}}}_{sr})^2}{\tilde{T}_{s}}- \frac{5}{2}\right)\right]. \end{align}

It clearly shows that the strengths of intra- and inter-species collisions are indicated by the magnitudes of $\tilde {n}_s$ and $\tilde {n}_r\phi _{sr}$, respectively. Therefore, in a gas mixture with a large disparity in concentration or mass, the intra-species collisions become dominant for the component $s$ with a major proportion of number density ($\tilde {n}_s \gg \tilde {n}_r$) or significantly heavier mass ($\phi _{sr} \ll 1$ when $m_r\ll m_s$).

4.1. Normal shock wave

We investigate the structure of normal shock waves under various free-stream conditions. The conditions are defined by the species mole fraction $\chi _s$ in the upstream, as well as the Mach number (${Ma}$) calculated based on the speed of sound ${v}_{mix}=\sqrt {5k_BT_u/3m_{mix}}$ in the upstream, where $T_u$ is the upstream temperature and $m_{mix}=\sum {m_s\chi _s}$ is the averaged mass of the mixture. Given the upstream conditions, the macroscopic quantities at the downstream end are determined by the classical Rankine–Hugoniot relation. The mass of the lighter species (denoted as species 1), the mixture number density $n_u$ and the temperature $T_u$ of the upstream flow are taken as reference values, namely, $m_0=m_1$, $n_0=n_u$, $T_0=T_u$. The characteristic length $L_0$ is set equal to the mean free path of the lighter species in the upstream flow, thus leading to ${Kn}_1=1$.

The numerical investigations of normal shock wave structures in gas mixtures based on the Boltzmann equation and the DSMC method have been conducted in the literature, taking into account various intermolecular interactions, including hard-sphere and ab initio potentials. We validate the proposed kinetic model by comparing its results with existing published data (Kosuge, Aoki & Takata Reference Kosuge, Aoki and Takata2001; Raines Reference Raines2002; Sharipov & Dias Reference Sharipov and Dias2018), as given in Appendix B. This comparison not only demonstrates the accuracy of our model but also highlights its applicability in accommodating diverse intermolecular potentials. Furthermore, we apply both the proposed kinetic model and the model of Kosuge (Reference Kosuge2009) to simulate a weak and a strong shock wave with ${Ma}=2$ and 10, respectively. The comparative analysis with DSMC simulations indicates that both models can accurately capture the structure of the weak shock wave, while our model outperforms the Kosuge model in modelling strong non-equilibrium shock waves, as detailed in Appendix C.

We then simulate the binary mixtures with large mass ratios (up to 1000) to determine the limits of our kinetic model. The simulation domain $[-L_x,L_x]$ is selected to ensure that the upstream and downstream are in equilibrium states, and then the specific values of $L_x$ are chosen as $40L_0,500L_0,30L_0$ for mixtures 1, 2 and 3, respectively, due to the significant differences in the properties of these gas mixtures and the upstream conditions. Numerical results of both the kinetic model and the DSMC method are compared in figures 2–4 for mixtures 1, 2 and 3, respectively. We present the normalized values for number density, flow velocity, temperature and dimensionless heat flux for each species and also the mixture, under the specified conditions with different concentrations $\chi _1=0.1,0.5,0.9$ and ${Ma}=3,4,5$. For mixtures 1 and 2 consisting of Maxwell gases, excellent agreements between the results of our kinetic model and the DSMC method are achieved (figures 2 and 3), even when the mass ratio is as high as 1000. For mixture 3 consisting of hard-sphere molecules, the number density and velocity of the lighter species predicted by the kinetic model deviate from DSMC results, when the mole fraction of the lighter species is small ($\chi _1=0.1$; see figures 4a and 4d). However, the average properties of the mixture remain highly accurate, since the lighter species, due to its low concentration, has a minimal impact on the overall behaviour of the mixture. Generally, three reasons may contribute to the possible deviation predicted by the kinetic model for mixtures consisting of non-Maxwell gases: (i) the parameters $c_{sr}$ and $d_{sr}$ accounting for the energy relaxation are derived approximately for non-Maxwell molecules; (ii) the term with parameter $b_{sr}$ is designed to capture the correct thermal-diffusion phenomena in the continuum limit, and thus, may have a discrepancy in strong non-equilibrium cases; (iii) the velocity-dependent collision frequency for non-Maxwell molecules is not recovered in the BGK-type operators (Yuan & Wu Reference Yuan and Wu2022). Nevertheless, as shown in figure 4, the kinetic model gives good overall agreement with DSMC simulations for hard-sphere gas mixtures. Therefore, the agreement suggests that the kinetic model can predict accurate results for realistic molecular models, whose behaviour usually lies between that of hard-sphere and Maxwell molecules.

Figure 2.

Comparisons of the normalized (a–c) number density, (d–f) flow velocity, (g–i) temperature and (j–l) dimensionless heat flux of the gas mixture between the kinetic model (lines) and the DSMC method (symbols) for the normal shock wave at ${Ma}=5$. The binary mixture consists of Maxwell molecules with a mass ratio $m_2/m_1=10$, a diameter ratio $d_2/d_1=1$ and the mole fraction of light species $\chi _1=0.1,0.5,0.9$.

Figure 3.

Comparisons of the normalized (a–c) number density, (d–f) flow velocity, (g–i) temperature and (j–l) dimensionless heat flux of the gas mixture between the kinetic model (lines) and the DSMC method (symbols) for the normal shock wave at ${Ma}=3$. The binary mixture consists of Maxwell molecules with a mass ratio $m_2/m_1=1000$, a diameter ratio $d_2/d_1=1$ and the mole fraction of light species $\chi _1=0.1,0.5,0.9$.

Figure 4.

Comparisons of the normalized (a–c) number density, (d–f) flow velocity, (g–i) temperature and (j–l) dimensionless heat flux of the gas mixture between the kinetic model (lines) and the DSMC method (symbols) for the normal shock wave at ${Ma}=4$. The binary mixture consists of hard-sphere molecules with a mass ratio $m_2/m_1=100$, a diameter ratio $d_2/d_1=2$ and the mole fraction of light species $\chi _1=0.1,0.5,0.9$.

Although the shock waves in mixtures with large mass ratios (e.g. 32.8 for Helium–Xenon) have been studied in the literature, and have shown unique characteristics absent in those composed of similar gas molecules. The shock structures can be significantly altered by substantial mass disparity. By comparing the shock structures in a wide range of mass ratios and species concentrations, the following features can be observed and correctly captured by our kinetic model.

1. (i) The shock wave thickness of a gas mixture is markedly thicker than that of a pure gas, especially when the mass ratio is large, as the mixture viscosity and diffusivity become stronger and the relaxation between components gets slower. For example, as shown in figure 3, the Maxwell gas mixture with a mass ratio of 1000 and diameter ratio of 1 form a significantly large transition zone from the upstream to the downstream, which spans several hundreds of the molecular mean free path.

2. (ii) In mixture 1 with a mass ratio of 10 (moderate mass difference), a pronounced temperature overshoot of the heavy species (higher than the downstream temperature) is observed when the heavier gas has only a small proportion (figure 2i), which has been shown in the literature from kinetic modelling (Bird Reference Bird1968; Kosuge et al. Reference Kosuge, Aoki and Takata2001; Sharipov & Dias Reference Sharipov and Dias2018) and hydrodynamic equations (Schmidt, Seiler & Wörner Reference Schmidt, Seiler and Wörner1984). However, when the mass ratio increases to 1000, the temperature overshoot gradually vanishes when $\chi _1=0.9$. As illustrated in figure 3, a comparable or even larger proportion of the heavier gas is required for the temperature overshoot to occur.

3. (iii) Flow velocity undershoot of the lighter gas happens in a mixture with a large mass ratio (100 and 1000) and a small proportion of lighter species ($\chi _1=0.1$) (figures 3d and 4d), which confirms the phenomena predicted by the multi-temperature hydrodynamic equations on weak shock (Goldman & Sirovich Reference Goldman and Sirovich1969). In these situations, the lighter gas decelerates to a velocity lower than that of the downstream flow, and thus, is compressed to a density above the downstream one.

4. (iv) A two-stage shock structure with distinct gradients of the mixture properties can be observed. In mixture 1 with a mass ratio of 10 and a large concentration of light species ($\chi _1=0.9$), as shown in the third column of figure 2, the properties of the shock exhibit a very steep change on the upstream side, while followed by a sudden change in the form of a long tail downstream. This phenomenon happened in a mixture with a moderate mass difference in experiments of Helium–Xenon (Gmurczyk, Tarczynski & Walenta Reference Gmurczyk, Tarczynski and Walenta1979) and found by hydrodynamic equations (Ruyev, Fedorov & Fomin Reference Ruyev, Fedorov and Fomin2002). However, when the mass ratio becomes significantly large in mixtures 2 and 3 (figures 3 and 4), this phenomenon disappears. On the contrary, the two-stage shock structure with an opposite trend occurs when there is only a small proportion of light species in the mixture ($\chi _1=0.1$). This structure consists of a smooth change in gas properties on the upstream side followed by a sudden and dramatic change on the downstream side.

4.2. Planar Fourier flow

A stationary rarefied mixture confined between two infinite parallel plates located at $x=0$ and $x=L_0$ is considered. The surfaces of the plates have different temperatures $T_w=T_0 \pm \Delta T$, and reflect the gas molecules in a fully diffuse way. The Knudsen number is determined in terms of the averaged number density of the mixture $n_0$, the temperature $T_0$ and the distance between the two plates $L_0$. A variety of cases were considered for the mixtures with different mole fractions $\chi _s$, temperature difference $\Delta T$ and Knudsen numbers, while a selection of the representative results is shown in figures 5 and 6. It can be seen that the solutions given by the kinetic model agree well with the DSMC results.

Figure 5.

Comparisons of the normalized number density (first column), dimensionless temperature (second column) and heat flux (third column) of the gas mixture between the kinetic model (lines) and the DSMC method (symbols) for the planar Fourier flow with $\Delta T=0.2T_0$. The binary mixture consists of Maxwell molecules with a diameter ratio $d_2/d_1=1$ and the mole fraction of light species $\chi _1=0.5$. Each row belongs to a specific flow condition: (a–c) mixture 1 ($m_2/m_1=10$), ${Kn}_1=0.1$; (d–f) mixture 1, ${Kn}_1=1$; (g–i) mixture 2 ($m_2/m_1=1000$), ${Kn}_1=0.1$; (j–l) mixture 2, ${Kn}_1=1$.

Figure 6.

Comparisons of the normalized number density (a,d), dimensionless temperature (b,e) and heat flux (c, f) of the gas mixture between the kinetic model (lines) and the DSMC method (symbols) for the Fourier flow, when ${Kn}_1=0.1$ and $\Delta T=0.8T_0$. The binary mixture consists of hard-sphere molecules with a mass ratio $m_2/m_1=100$, a diameter ratio $d_2/d_1=2$ and the mole fraction of light species $\chi _1=0.1$ (a–c) and 0.9 (d–f). The black lines are the kinetic model solutions without thermal-diffusion effects ($b_{12}=0$).

The temperatures of the components in mixture 1 ($m_2/m_1=10$) stay close when ${Kn}$ is up to 1 (figures 5a–5f), while pronounced temperature separation and concentration variation can be observed for mixture 2 ($m_2/m_1=1000$) when ${Kn}=0.1$ (figures 5g–5l). It is noteworthy that all these Maxwell molecules under consideration have the same size of mean free path, due to their identical diameter and interaction potential. However, the inter-species relaxation is much slower for the gases with disparate mass, although the spatial Knudsen number is the same. In other words, the mixtures with larger mass ratios may exhibit significant non-equilibrium phenomena even at small ${Kn}$, and hence, shrink the applicable range of the hydrodynamic description of the mixtures. The same observation has been found in a previous work solving the linearised Boltzmann equation and McCormack model (Ho et al. Reference Ho, Wu, Graur, Zhang and Reese2016).

Unlike Maxwell gas mixtures, where the thermal-diffusion effect is absent, mixtures of hard-sphere molecules have a significant thermal-diffusion effect. As shown in figure 6 (${Kn}_1=0.1$ and $\Delta T=0.8T_0$), although the temperature of the two species remains the same, the concentration ratio between components varies across the domain due to the temperature gradient, and the heavy gas molecules tend to concentrate in the cold region. Meanwhile, the highly nonlinear feature arising from the pronounced temperature difference between the two plates is accurately captured. We also solve the kinetic model without the thermal-diffusion effect by setting parameter $b_{12}=0$ (other transport properties and relaxation rates remain unchanged), corresponding results are shown in figure 6 by black lines. Clearly, the species separation phenomenon cannot be reproduced, thus leading to an incorrect prediction of concentration and, hence, the heat flux of each species.

Additionally, we also compare our kinetic model with the DSMC results reported by Strapasson & Sharipov (Reference Strapasson and Sharipov2014), where the heat transfer in a He–Ar mixture confined between two parallel plates are investigated across a wide range of Knudsen numbers. Figure 7 presents the dimensionless heat flux $q^*$ for a He–Ar mixture with hard-sphere molecules as functions of the rarefaction parameter $\delta$ ranging from 0.01 to 40, for various mole fractions of helium $\chi _{\textrm {He}}=0.25,0.5,0.75$, and temperature differences between two plates $\Delta T/T_0=0.1,0.75$. The maximum relative error, observed in the moderate rarefaction regime, is consistently below 2.8 %. Note that the rarefaction parameter $\delta$ is related to the Knudsen number through the equation

(4.3)\begin{equation} {Kn}_s = \frac{\sqrt{\rm \pi}\mu_s}{2\mu_{mix}}\sqrt{\frac{m_{mix}}{m_s}}\frac{1}{\delta}, \end{equation}

where $\mu _{mix}$ is the viscosity of the mixture. The dimensionless heat flux $q^*$ is defined as

(4.4)\begin{equation} q^* = \frac{q}{p_0}\sqrt{\frac{m_{mix}}{2k_BT_0}}\frac{T_0}{\Delta T}. \end{equation}

Figure 7.

Comparison of the dimensionless heat flux between the kinetic model (lines) and the DSMC method (symbols) (Strapasson & Sharipov Reference Strapasson and Sharipov2014) for He–Ar mixture with hard-sphere molecules at rarefaction parameter $\delta$ ranging from 0.01 to 40, with mole fractions $\chi _{\textrm {He}}=0.25,0.5,0.75$ and temperature differences between the two plates $\Delta T/T_0=0.1,0.75$. The corresponding parameters in the kinetic model fitted by matching the mixture viscosity and thermal conductivity are $\phi _{12}=2.149$, $\phi _{21}=0.3775$, $\varphi _{12}=1.373$, $\varphi _{21}=1.642$.

4.3. Planar Couette flow

The Couette flow shares the same configuration as that of the Fourier flow, but the temperatures of both plates are kept the same at $T_0$, and the lower and upper plates move along the $y$ direction with velocity $u_{1,y}=-v_w$ and $u_{2,y}=v_w$, respectively, where $v_w=\sqrt {2k_BT_0/m_{mix}}$. Due to the symmetry, only half of the domain $(L_0/2\le {}x_2\le {}L_0)$ is used in the simulations.

The simulation results of our kinetic model and the DSMC method are shown in figure 8 for the three mixtures with $\chi _1=0.5$, when ${Kn}_1=1$ for mixture 1 and ${Kn}_1=0.1$ for the others. Good agreements are found for all the macroscopic properties. We analyse the mixture shear stress predicted by two methods across various mole fractions $\chi _1$ and Knudsen numbers ${Kn}_1$ (table 2). Note that the actual plate's velocities vary significantly for mixtures with different mass ratios and mole fractions, and the boundary velocity is supersonic or even hypersonic for the heavier species, but subsonic for lighter gas in most cases (exclude the one for mixture 1 with $\chi _1=0.9$). Nevertheless, the maximum relative error of shear stress in all the considered cases remains below $3.2\,\%$. It should be noted that the moderate rarefaction level extends to a wider range of Knudsen numbers for a mixture with a larger mass difference, because of the presence of multiscale species ${Kn}$ and relaxation times of inter-species collisions. Importantly, the data demonstrate a consistent level of accuracy for our kinetic model, regardless of variations in the mass ratio of the gas species.

Figure 8.

Comparisons of the (a–c) normalized number density and dimensionless (d–f) flow velocity, (g–i) temperature and (j–l) heat flux between the kinetic model (lines) and the DSMC method (symbols) for the Couette flow with the mole fraction of light species $\chi _1=0.5$, when ${Kn}_1=1$ (mixture 1 in the first column) and 0.1 (mixture 2 in the second column and mixture 3 in the third column).

Table 2.

The dimensionless shear stress $P_{xy}/n_0k_BT_0$ of the mixture calculated by the kinetic model equation for the Couette flow, with the mole fraction $\chi _1=0.1,0.5,0.9$ and ${Kn}_1=0.1,1,10$. The values in parentheses are the relative errors between the results given by the kinetic model and the DSMC method.

5.1. Supersonic flow around a circular cylinder

We consider the supersonic gas flow with density $n_0$ at ${Ma}_{\infty }=3$ passing a cylinder with diameter $L_0$. The temperatures of both the free stream and isothermal surfaces of the cylinder are maintained at $T_0$. The Knudsen number of the lighter species in the free stream is ${Kn}_{1}=0.5$. Based on the diameter ratio of the mixtures’ components, the Knudsen number of the heavier species in mixtures 1 and 2 is ${Kn}_{2}=0.5$ as well, while that in mixture 3 of hard-sphere molecules is ${Kn}_{2}=0.125$. Besides, the mole fraction considered in this problem is $\chi _1=0.5$ for all mixtures. The simulations are conducted only in the upper half-domain $[-L,L]\times [0,L]$ due to symmetry, with $L/L_0=6,12,10$ for mixtures 1, 2 and 3, respectively.

The detailed flow fields about species number density and flow velocity of the surrounding gas are presented in figure 9 for the gas mixture with a mass ratio of 1000, where the result given by our kinetic model matches the DSMC data well. Figure 10 compares kinetic model results of the windward side number density, flow velocity, temperature and heat flux along the stagnation line with those solved by the DSMC method, and the overall agreement is very good.

Figure 9.

Comparisons of the dimensionless (a) number density and (b) flow velocity between the results of the kinetic model (black lines) and the DSMC method (background contours) for a supersonic mixture flow passing a cylinder, when the mole fraction of the free stream is $\chi _1=0.5$, ${Ma}_{\infty }=3$ and ${Kn}_{1}=0.5$. The binary mixture consists of Maxwell molecules with a mass ratio $m_2/m_1=1000$ and a diameter ratio $d_2/d_1=1$.

Figure 10.

Comparisons of the dimensionless (a–c) number density, (d–f) flow velocity in the $x$ direction, (g–i) temperature and (j–l) heat flux in the $x$ direction, along the windward side stagnation line between the results of the kinetic model and the DSMC method for a supersonic mixture flow passing a cylinder, when the mole fraction of the free stream is $\chi _1=0.5$, ${Ma}_{\infty }=3$ and ${Kn}_{1}=0.5$. The first, second and third columns correspond to mixture 1, 2 and 3, respectively.

Despite having the same free-stream flow velocity that exceeds the sound speed of the mixture, individual species in a gas mixture with disparate mass experience vastly different flow characteristics. This is due to their distinct species Mach numbers ${Ma}_{\infty,s}$, which are defined based on their own individual sound speeds $\sqrt {5k_BT_0/3m_{s}}$. For instance, ${Ma}_{\infty,1}=0.13$ and ${Ma}_{\infty,2}=4.24$ for the mixture shown in figure 9. Therefore, the lighter component forms a subsonic flow field and exhibits significantly less compression compared with the heavier species. As shown by the density profiles in figure 10 for the mixture with a mass ratio of 1000, the number density of the heavier gas gets nearly 10 times that of the lighter one near the stagnation point.

The overall properties of the shock are mainly determined by the heavier species, because of its higher number density and molecular mass, as long as the components’ mole fractions in the free stream are comparable. However, it is found that, as the mass ratio changes from 10 to 1000 for mixtures 1 and 2, the thickness of the shock wave and the peak values of the mixture temperature increase only slightly (10 %). Comparatively, mixture 3 has a thinner shock structure since its effective Knudsen number is lower (${Kn}_{2}={Kn}_{1}/4$). Note that the actual free-stream velocities of the mixtures vary significantly due to the distinct average mixture mass, though ${Ma}_{\infty }$ keeps constant.

Figure 11 shows the pressure and heat flux along the surface of the cylinder. The pressure predicted by the kinetic model matches the DSMC results very well, while the heat fluxes given by the two methods have a 6 % relative difference around the windward side stagnation region for mixture 3 consisting of hard-sphere molecules. Interestingly, the forces acting on the object are found to be very close in the three types of mixtures, despite the disparity in average mixture mass, and even the intermolecular potential. Particularly, for the flow of mixtures 1 and 2, the aerodynamic forces are nearly the same and, hence, roughly independent of the mass ratio, which is the only different dimensionless variable in the two flows. Also, the force is found to be insensitive to the intermolecular potential and molecular diameter ratio (mixture 3). On the other hand, the values of heat flux on the surface are not only inversely scaled by the square root of the mass ratio, but also significantly influenced by the intermolecular potential and molecular diameter ratio.

Figure 11.

Comparisons of the dimensionless (a–c) pressure and (d–f) heat flux along the surface of the cylinder between the results of the kinetic model and the DSMC method for the mixtures, when the mole fraction of the free stream is $\chi _1=0.5$, ${Ma}_{\infty }=3$ and ${Kn}_{1}=0.5$. Here $p_x$ and $p_y$ denote the pressures exerted on the surface along the $x$ and $y$ directions, respectively. The angle is measured from the windward to the leeward side. The first, second and third columns correspond to mixture 1, 2 and 3, respectively.

5.2. Nozzle flow

Our kinetic model is applied to simulate a two-dimensional rarefied gas mixture flowing through a nozzle into the vacuum. The structure of the nozzle is shown in figure 12, which has a straight channel with width $L_0$, a converging section shrinking the width to $L_0/2$ at the throat and a diverging section. At the inlet of the nozzle ($x=0$), the flow was assumed to be considerably subsonic ${Ma}_{in}=0.05$ with Maxwellian velocity distributions at temperature $T_0$. The gas molecules are reflected on the cold walls of the nozzle ($T_w=T_0/2$) with complete thermal accommodation, and then go through the outlet ($x=3L_0$) into the vacuum. The Knudsen number, defined in terms of the gas properties at the inlet, is ${Kn}_1=0.1$ for lighter species for all the mixtures.

Figure 12.

The dimensionless (a,c) number density and (b,d) local Mach number solved by the kinetic model (black lines) and the DSMC method (background contours) for a gas mixture flowing through a nozzle, when ${Kn}_{1}=0.1$ at the inlet. The inlet mole fraction is $\chi _1=0.05$ for mixture 2 (a,b) and $\chi _1=0.5$ for mixture 3 (c,d).

We consider a very small proportion of lighter molecules ($\chi _1=0.05$) mixed with heavier ones for mixture 2 and a half-half mix ($\chi _1=0.5$) of hard-sphere molecules for mixture 3. Figure 12 shows the species number density and local Mach number distribution in the nozzle solved by the kinetic model and the DSMC method. Also, figure 13 shows the density, temperature and heat flux along the centreline of the nozzle. The good accuracy of our model equation is demonstrated.

Figure 13.

Comparisons of the dimensionless (a,d) flow velocity, (b,e) temperature and (c, f) heat flux along the centreline between the results of the kinetic model and the DSMC method for a gas mixture flowing through a nozzle, when ${Kn}_{1}=0.1$ at the inlet. The inlet mole fraction is $\chi _1=0.05$ for mixture 2 (a–c) and $\chi _1=0.5$ for mixture 3 (d–f).

It is found that, when a gas mixture with disparate mass flows through a nozzle, the density and temperatures of each component do not experience remarkable separation across the nozzle. However, the species velocities at the outlet become noticeably different. As illustrated in figures 13(a) and 13(d), mixture 2, having a mass ratio of 1000, exhibits an outlet velocity for the lighter species that is 20 times higher than that of the heavier one; similarly, mixture 3 shows a velocity ratio of 3.7 between its lighter and heavier components. Consequently, the significant diffusion velocities make a noticeable contribution to the total heat flux of the mixture, particularly in the diverging section of the nozzle, where the heat flux due to conductance becomes negligible. On the other hand, the two components achieve close values of species Mach number ${Ma}_{s,local}$ calculated using their respective local sound speeds, as plotted in figures 12(b) and 12(d), which also indicates similar degrees of compression of the disparate species flowing through the nozzle. Meanwhile, the relatively higher ${Ma}_{s,local}$ of the heavier components primarily arises from inter-species collisions with lighter molecules, which accelerate the heavier ones.

We also calculate the flow rate of individual components passing through the nozzle outlet, as plotted in figure 14. The molecular number flux is found to be comparable for the mixture components, while the heavier species will hence dominate the mass flow rate due to its much higher molecular mass. Compared with DSMC results, the mass flow rate predicted by the kinetic model yields less than 1.7 % relative difference for all the species in mixture 2 and 3.2 % relative difference for those in mixture 3.

Figure 14.

Comparisons of the dimensionless molecular number flux of each species passing through the outlet of the nozzle between the results of the kinetic model (lines) and the DSMC method (symbols). The results of mixtures 2 and 3 are given in (a) and (b), respectively.