Comparative study of the Boltzmann and McCormack equations for Couette and Fourier flows of binary gaseous mixtures

We evaluate the accuracy of the McCormack model by comparing its solutions for Couette and Fourier ﬂows of binary gaseous mixtures with results from the linearized Boltzmann equation. Numerical simulations of Ne–Ar and He–Xe gas mixtures are carried out from slip to near free-molecular ﬂow regimes for different values of the molar concentration. Our numerical results show that while there are only small differences in the shear stress in Couette ﬂow and the heat ﬂux in Fourier ﬂow, calculated from the two kinetic equations, differences in other macroscopic quantities can be very large, especially in free-molecular ﬂow regime. Moreover, the difference between results from the two models increases with the molecular mass ratio and the molar concentration of the heavier species. Finally, the applicability of the McCormack model, which was derived for linearized ﬂows only, is investigated by comparing its solutions with those from the Boltzmann equation for Fourier ﬂow with large wall-temperature ratios.


Introduction
In practical applications like vacuum technology, porous media, and the chemical industry, information about the heat/mass transfer in rarefied gaseous mixtures is indispensable. Benchmark test cases are of great importance as they can validate new numerical models developed to describe gas flows or test the validity of existing approaches under various physical conditions. In this paper, Couette and Fourier flows between two parallel plates are chosen as benchmark test cases as they are classical problems in fluid mechanics. Although solutions for single-species gases have been well studied, few papers have investigated gaseous mixtures.
The plane Couette flow of a binary gaseous mixture was first studied [1][2][3][4][5] using kinetic models for the Boltzmann equation (BE), such as the Hamel model [6] for Maxwellian molecules and the McCormack model for general intermolecular potentials [7]. Notably, following the McCormack model, the influence of intermolecular interactions on the velocity and shear stress in three mixtures (Ne-Ar, He-Ar, and He-Xe) [4] and the influence of the gas-surface interaction on the flow properties were investigated [5]. The linearized Boltzmann equation (LBE) for hard-sphere (HS) molecules has been solved by an analytical version of the discrete-ordinates (ADO) method [8], and the accuracy of the McCormack model has been assessed for a He-Ar mixture: the McCormack model produces accurate shear stress for each species, and the velocity of the heavier species [5,8]; however, the velocity of the lighter species and especially the heat flux significantly deviate from the LBE results (by over 100% for some case).
Very few papers have tackled the heat transfer through a binary gaseous mixture. Plane Fourier flow was first simulated by solving the BE using the numerical kernel method [9]. Later, the heat transfer between two plates with a small temperature difference was studied using the McCormack model [10] and the LBE [11]. Surprisingly, the normalized heat flux for Ne-Ar and He-Xe mixtures obtained from the linearized equations were found to agree with results from the BE, with the maximum relative deviation being about 4%. However, there were large differences in the density and temperature between the McCormack model and the LBE: for density, up to 15% difference in Ne-Ar mixture and 51% in He-Xe mixture were observed, while for temperature the maximum differences were 12% and 20% for Ne-Ar and He-Xe mixtures, respectively. The influence of intermolecular potentials on the heat flux between two parallel plates, for three binary mixtures of noble gases (Ne-Ar, He-Ar, He-Xe), has been studied using the McCormack model [12]: the heat flux is sensitive to the intermolecular potential, and the difference between the HS and realistic potential [13] reached 15% near the hydrodynamic regime.
To summarize, only two papers have compared the McCormack model and the LBE for mixtures of Ne-Ar, He-Xe [11], and He-Ar [8]. Therefore, systematic new comparisons between the McCormack model and the LBE will be useful for further development of numerical tools for the simulation of gaseous mixture flows. In this paper, Couette and Fourier flows are considered with two types of binary gaseous mixture composition. Different species accommodation coefficients are considered, and their influence on the flow parameters is analyzed. A large temperature difference between the two plates in Fourier flow is also investigated, where the results from the McCormack model are compared with those from the BE in order to establish the limits of the linearized approach.

Problem statement
Consider a binary mixture of monatomic gases, where the mass of a molecule of the first (second) species is m 1 (m 2 ), and the corresponding number density is n 1 (n 2 ). Without loss of generality, we assume m 1 < m 2 . The gaseous mixture is confined between two parallel plates situated at y 0 ¼ AEH=2, see Fig. 1. In Couette flow, the two plates with temperature T 0 move in opposite directions with a speed U=2. In Fourier flow, both plates are stationary, but the plate at y 0 ¼ ÀH=2 has a temperature of T C ¼ T 0 À DT=2, while the other one has T H ¼ T 0 þ DT=2. We assume that, in the case of Couette flow, the relative speed of the two plates U is much smaller than the most probable molecular velocity v 0 of the mixture. We also assume that the temperature difference DT, in the case of Fourier flow, is much smaller than the equilibrium gas temperature T 0 , so the gaseous mixture deviates only slightly from thermodynamic equilibrium and the McCormack model can be applied.
The most probable molecular velocity of the mixture is: where k is the Boltzmann constant and m ¼ C 0 m 1 þ 1 À C 0 ð Þ m 2 is the mean molecular mass of the mixture. Here, C 0 is the equilibrium molar concentration of the lighter species, and n 0a is the equilibrium number density of species a (a ¼ 1; 2Þ.
The equilibrium number density of the mixture is n 0 ¼ n 01 þ n 02 . In Couette flow, we are interested in the profiles of the species velocity u 0 ax and shear stress P 0 axy , while in Fourier flow, we are interested in the deviated density n 0 a , concentration C 0 ¼ ðn 01 þ n 0 1 Þ=ðn 01 þ n 02 þ n 0 1 þ n 0 2 Þ, deviated temperature T 0 a , and heat flux q 0 ay .

Kinetic equations
To describe the gas dynamics at various conditions, the gas kinetic theory is necessary. In this section, the BE for the binary gaseous mixture is first introduced. Then, the McCormack model is used to simplify the Boltzmann collision operator for linearized Couette and Fourier flows. Finally, the gas-wall boundary condition is specified and numerical techniques to solve the kinetic equations are briefly described.

The Boltzmann equation
Let f a ðt; x; vÞ be the distribution function of specie a with molecular velocity v at spatial location x and time t. In the absence of external forces, the following BE describes the evolution of f 1 and f 2 : where the Boltzmann collision operator Q ab ðf a ; f b Þ is In the above equations, v and v Ã are the pre-collision velocities of molecules of species a and b, respectively, while 0 v ab ; 0 v ab Ã are the corresponding post-collision velocities. Conservation of momentum and energy yield the following relations: is the relative pre-collision velocity, X is the unit vector in the sphere S 2 with the same direction as the relative post-collision velocity, and h is the deflection angle between the two relative velocities, i.e. cos h ¼ X Á v r =jv r j. Finally, for hardsphere molecules, the collision kernels C ab are given by where d a is the molecular diameter of species a.
For systems that only slightly deviate from equilibrium, the BE can be linearized. Introducing the equilibrium Maxwellian distribution function: with the linearized Boltzmann collision operator L ab ðf When the perturbation function f d a is known, the deviated (from equilibrium values) macroscopic flow quantities of species a, such as the density, velocity, shear stress, temperature, and heat flux, are calculated as follows:

The McCormack model
The McCormack model was proposed for linearized flows of multicomponent monatomic mixture [7], where the linearized collision operator is obtained by requiring that its first three velocity moments be the same as the corresponding moments of the lin- where the linearized collision operator L ab h is given in A. Note that in writing Eq. (8), we have used the following dimensionless quantities: where c a is the molecular velocity and p 0a ¼ n 0a kT 0 is the partial pressure of species a.
When the perturbation functions h a are known, the velocity and a dc a , respectively, while the deviated density, deviated temperature, and the heat flux in Fourier flow are The dimensionless macroscopic quantities of the binary gaseous mixture in Couette flow are defined as follows: while those in Fourier flow are defined as follows: where p 0 ¼ n 0 kT 0 is the equilibrium pressure of the mixture. Expressions for the macroscopic quantities of the binary gaseous mixture in the LBE can be given similarly.

Boundary conditions
The Maxwell diffuse-specular boundary condition is adopted to describe the gas-wall interaction. When molecules hit the plate, an a a portion of them (0 < a a 1) are diffusely reflected, while the rest are specularly reflected. There is complete accommodation when the accommodation coefficient a a is unitary. In this paper, if not otherwise specified, complete accommodation is assumed.
In the linearized Couette flow the boundary condition for the LBE has the following form while that for the McCormack model is [4] h þ a ðy ¼ AE1=2Þ where the superscripts þ and À of the perturbation functions in Eqs. (11) and (12) refer to the outgoing and incoming molecules with respect to the plates' surfaces, respectively. In the linearized Fourier flow the boundary condition for the LBE reads while that for the McCormack model is [12] h þ a ðy ¼ AE1=2Þ

Numerical techniques
The Boltzmann collision operators in the BE (3) and LBE (7) are solved by the fast spectral method (FSM) [14]. The main idea of the FSM is to expand the distribution function and collision operator into Fourier series, and handle the binary collision in a corresponding frequency space. The method can deal with highly rarefied gas flows, where the distribution function has large discontinuities. The number of discretized velocities can be large to capture the discontinuities, but the number of frequency components is relatively small, resulting in high computational accuracy and efficiency [15]. We discretize the three-dimensional molecular velocity space by 32 Â 64 Â 32 points (the points are uniformly distributed in the v x and v z directions, while the discretization in the v y direction is non-uniform, with most of them concentrating on v y $ 0 to capture the discontinuities in the molecular distribution function), while the corresponding frequency space is uniformly discretized by 32 Â 32 Â 32 points.
When the Boltzmann collision operator is obtained, the LBE (7) is solved by an iterative scheme [14], where the spatial derivative is approximated by a second-order upwind finite-difference. The physical half-space À1=2 6 y 0 is discretized by 50 points nonuniformly, with most points located near the plate to capture the velocity slip and temperature jump. The iterations terminate when the maximum relative difference (excluding the point y ¼ 0) in macroscopic quantities (such as the velocity and shear stress in Couette flow, and the density, temperature, and heat flux in Fourier flow) between two consecutive steps is less than 10 À5 .
The discrete velocity method [4]  Á dc az . The entire physical space is uniformly discretized into 400 points, and a first-order finite-difference is used to approximate the spatial derivative. The two-dimensional reduced molecular velocity space is discretized by 50 Â 50 points according to Gaussian-Hermit quadrature. Grid-independence is checked with R 1=2 À1=2 jA lþ1 =A l À 1j < 10 À10 is employed, where l is the iteration step, and A ¼ P xy in the Couette flow and A ¼ q y in the Fourier flow.

Numerical results
The binary gaseous mixtures Ne-Ar and He-Xe are chosen in this comparative study in order to investigate the influence of the molecular mass ratio. Three values of the molar concentration C 0 ¼ 0:1; 0:5, and 0.9 are considered at three values of the rarefaction parameter d 0 ¼ 0:1; 1, and 10. For Ne-Ar and He-Xe mixtures at an equilibrium temperature T 0 ¼ 300K, the molecular diameter ratios are d 2 =d 1 ¼ 1:406 and 2:226, respectively, while the corresponding molecular mass ratios m 2 =m 1 are 1:979 and 32:8. Both Couette and Fourier flows are characterized by the following rarefaction parameter [17]: where l ¼ l 1 þ l 2 is the viscosity of the mixture, see Appendix A.

Couette flow
The velocity in the Ne-Ar mixture with molar concentration C 0 ¼ 0:5 at d 0 ¼ 0:1; 1, and 10 is shown in the left column of Fig. 2. When d 0 ¼ 0:1, the Ar gas velocity is different by 9% between the McCormack model and the LBE, while the difference in Ne is less than 1%. When d ¼ 1, the difference in the Ar velocity quickly decreases, while that in the Ne velocity slightly increases. By d ¼ 10, there are practically no differences in all the respective velocities. The influence of the molecular mass ratio on the velocity profile is observed by comparing the left (Ne-Ar) and right (He-Xe) columns in Fig. 2 The influence of the molar concentration C 0 on the gas velocity at the plate is summarized in Table 1. As C 0 varies from 0.1 to 0.9, the difference in Ar velocity between the McCormack model and the LBE increases from 6% to 13% when d 0 ¼ 0:1. For d 0 ¼ 10, this difference is less than 1% for all considered values of C 0 . On the other hand, the relative difference in the Ne velocity between the two kinetic equations is less than 4% for all values of C 0 and d 0 : the maximum 4% is when C 0 ¼ 0:9 and d 0 ¼ 0:1.
The shear stresses P xy in the Ne-Ar mixture are listed in Table 2 for different values of the molar concentration and the rarefaction parameter. Theoretically, the mixture shear stress has to be constant, but numerical results vary slightly across the channel. Therefore, the average shear stress P av xy ¼ R 1=2 À1=2 P xy dy is presented. The maximum variation in the shear stress across the channel, max i jðP xy ðiÞ À P av xy Þ=P av xy j, is less than 0.4%, which indicates the good numerical accuracy of both the FSM and the discrete velocity method. Good agreement between the McCormack model and the LBE is observed: the relative difference increases from 0.1% when d 0 ¼ 0:1 up to 2% when d 0 ¼ 10, for all molar concentrations. We have also calculated the shear stress of the He-Xe mixture with C 0 ¼ 0:5. They are, respectively, 0.2163, 0.1482, and 0.0400 when d 0 ¼ 0:1; 1, and 10 from the McCormack model, and 0.2150, 0.1491, and 0.0411 from the LBE. The relative difference in the shear stress in the He-Xe mixture between the two kinetic equations is slightly higher than that in the Ne-Ar mixture.

Fourier flow
The number density in the Ne-Ar mixture with C 0 ¼ 0:5 when d 0 ¼ 0:1; 1, and 10 is shown in the left column of Fig. 3. The difference in the Ne number density predicted by the McCormack model and the LBE reaches 31% when d 0 ¼ 0:1. This difference decreases as the rarefaction parameter increases, with only 5% difference by d 0 ¼ 10. Fig. 3 also shows that for mixtures with disparate molecular masses (i.e. He-Xe) the disagreement in the number density between the McCormack model and the LBE is larger than for the Ne-Ar mixture, especially at small values of the rarefaction parameter. The influence of the molar concentration on the number density is shown in Table 3. Large differences, i.e. 38% and 24%, are found in Ne and Ar, respectively, when d 0 ¼ 0:1 and C 0 ¼ 0:1. This difference is reduced when the molar concentration of Ne is increased: at C 0 ¼ 0:9, the maximum differences in Ne and Ar are 22% and 3%, respectively.
The molar concentration C in the Ne-Ar and He-Xe mixtures, which is defined as the deviation of the concentration of the lighter species from its initial equilibrium state, is shown in Fig. 4. Positive values of C mean that the concentration of the lighter species increases near the hotter plate due to thermodiffusion. Although the McCormack model, in contrast to several other kinetic models such as the Hamel model [6], can describe the thermodiffusion phenomenon in gaseous mixtures, comparisons in Fig. 4 show that the McCormack model underestimates the thermodiffusion significantly when compared to the LBE. For instance, when d 0 ¼ 0:1, the McCormack model underpredicts the molar concentration at the hotter plate by 40% and 67% for Ne-Ar and He-Xe mixtures, respectively. Even when d 0 ¼ 10, the relative differences between the McCormack and LBE results are as high as 14% and 20% for Ne-Ar and He-Xe mixtures, respectively.
The temperature variation in the Ne-Ar mixture with C 0 ¼ 0:5 when d 0 ¼ 0:1; 1, and 10 is shown in the left column of Fig. 5. Good agreement between the McCormack model and the LBE is found when d 0 ¼ 10. However, when d 0 ¼ 0:1, the difference in the temperature results for Ne reaches 13%, while the difference for Ar is negligible. The influence of the molecular mass ratio on the temperature is seen by comparing the left and right columns of The influence of the molar concentration on the gas temperature at the hotter plate is shown in Table 4. When d 0 ¼ 0:1, the maximum difference between results from the two kinetic equations. (18%) is found for Ne when C 0 ¼ 0:1, and this reduces to 7% when C 0 ¼ 0:9. So the difference between the two kinetic equations decreases as C 0 increases.
The heat flux in the Ne-Ar mixture is given in Table 5. As with the shear stress in Couette flow, here the average heat flux across the channel is reported. From Table 5 we see that the heat flux has its maximum value for the molar concentration C 0 ¼ 0:5 for all rarefactions. The maximum difference between the McCormack model and the LBE results is about 2%. We have also calculated the heat flux in the 3622. This means that the relative difference in the He-Xe mixture between the two models is greater than in the Ne-Ar mixture.

Effect of the incomplete accommodation
We studied the influence of the gas-wall interaction by simulating the incomplete accommodation. In the following, simulations of a Ne-Ar mixture with C 0 ¼ 0:5 are carried out for both Couette and Fourier flows, where the accommodation coefficients for Ne and Ar are a Ne ¼ 0:6 and a Ar ¼ 0:8, respectively.

Couette flow with incomplete accommodation
The influence of the incomplete accommodation on the velocity can be seen by comparing data in Table 6 with those in Table 1. Clearly, the species velocities decrease with a decrease in the accommodation coefficient. Profiles of the shear stress are shown in Fig. 6. We find the shear stress of each species is quite different,  Table 2, it is clear that the shear stress decreases with the accommodation coefficient. The relative difference (about 2%) between the results from the two kinetic equations remains the same as in the complete accommodation case.

Fourier flow with incomplete accommodation
The number density and temperature of each species and the mixture at the hot plate are presented in Table 6. The influence of incomplete accommodation can be seen by comparing the corresponding data in Tables 3 and 4. We find that the agreement between the McCormack model and the LBE is essentially improved when there is incomplete accommodation. For instance, with incomplete accommodation, the differences in the density of Ne and Ar when d 0 ¼ 0:1 are reduced from 31% to 21%, and from 14% to 5%, respectively. For the temperature, the difference in Ne when d 0 ¼ 0:1 is reduced from 13% to 7%, but that in Ar increases from 2% to 7%.
The influence of incomplete accommodation on the concentration is shown in Fig. 7. The absolute value of the concentration increases when the accommodation coefficient decreases, especially in the near free-molecular regime (d 0 ¼ 0:1), and the agreement between the McCormack model and the LBE becomes better for incomplete accommodation.
The absolute value of the heat flux is shown in Fig. 8.   Table 5 when C 0 ¼ 0:5 and there is complete accommodation, the heat flux decreases with a decrease in the accommodation coefficient, but the difference between the McCormack model and the LBE results remains unchanged.

The limit of the linearized approach
In the above sections, numerical results obtained from the McCormack model are compared with results from the LBE because the McCormack model requires the deviations of the gas parameters from their equilibrium values to be small. Therefore, it is useful to determine the applicability limit of this linearized approach, especially for the case of the heat flux through a gaseous mixture with large wall-temperature ratios. We carried out simulations on the heat transfer between two parallel plates with the temperature differences DT=T 0 ¼ 0:6 and 1, which leads to the large temperature ratios T H =T C ¼ 1:8 and 3, respectively. Under these conditions, the assumption of a small deviation of the plate's temperature from its equilibrium value (i.e. DT ( T 0 ) is not fulfilled. In this case the full BE is solved by the FSM [14] for the equal-mole Ne-Ar gas mixture with complete accommodation at three values of the rarefaction parameter. Good agreement in the number density is observed in Fig. 9 when d 0 ¼ 0:1 and 10, except in a small region in the vicinity of the wall (especially in the slip-flow regime and near the colder plate).
For the comparison of temperature profiles between the McCormack model and the BE, we note that near the free-molecular regime the gas temperature tends to a constant value However, in the linearized approach we calculate the temperature deviation from its reference value T 0 ¼ ðT C þ T H Þ=2. Therefore, to compare the temperature obtained from the McCormack model with the results from the BE when d 0 ¼ 0:1, we take T FM as the reference temperature. From Fig. 10b it is clear that the use of T FM as the reference temperature makes the comparison meaningful, and the temperature profiles from the McCormack model display good agreement with these from the BE. This is not the case when the mean temperature T 0 is chosen as the reference value, and large shifts between the two kinetic model results are observed, see Fig. 10a. In the slip-flow regime (d 0 ¼ 10), however, choosing the mean temperature T 0 as the reference value (Fig. 10c), instead of the free-molecular value T FM (Fig. 10d) enables a meaningful comparison, and perfect agreement between the two kinetic equations is found.
The dimensionless heat fluxes obtained from the McCormack model, the LBE, and the BE are compared in Table 7. Good agreement of the McCormack results with the solution of the BE, of the order of 4%, is found at DT=T 0 ¼ 0:6. However, for the temperature difference, DT=T 0 ¼ 1, the deviation between the McCormack model and the BE increases, to 9%.

Conclusion
Plane Couette and Fourier flows of binary gaseous mixture have been simulated using the McCormack kinetic model and the linearized Boltzmann equation, over a wide range of the molar concentration and rarefaction parameter. Two gaseous mixtures, one with similar molecular masses (Ne-Ar) and the other with disparate molecular masses (He-Xe), have been considered. Our numerical results showed that when only the shear stress in Couette flow and the heat flux in Fourier flow are required, the McCormack model can be used, as the differences in the results from the two kinetic models are within 2%. However, difference in other Fig. 6. The shear stress in the linearized Couette flow of an equal-mole Ne-Ar gas mixture with complete (aNe ¼ aAr ¼ 1:0) and incomplete (aNe ¼ 0:6; aAr ¼ 0:8) accommodation.

Table 6
Velocity ux in Couette flow, and number density n with temperature T in Fourier flow, at y ¼ 0:5 in the Ne-Ar gas mixture with C0 ¼ 0:5. The accommodation coefficients are aNe ¼ 0:6 and aAr ¼ 0:8.   Table 5 Total heat flux q y in the linearized Fourier flow of a Ne-Ar gas mixture. On the other hand, if the molecular mass ratio and the rarefaction parameter are fixed, the difference increases with an increase in the molar concentration of the heavier species. The fact that the  McCormack model is not accurate for highly rarefied gas is because the collision term of McCormack model only recovers the relaxation rate of low order momentum (such as pressure tensor, heat flux), while at small value of d 0 the relaxation of higher order moments is important.
Finally, we investigated the applicability of the McCormack model in nonlinear Fourier flow. We found that the McCormack model can provide temperature profiles in relatively good agreement with those from the Boltzmann equation, even for large differences in the plates' temperatures. However, the reference  temperature has to be adapted to the rarefaction parameter: in the slip-flow and hydrodynamic regimes the arithmetic mean of the two plates' temperatures should be chosen as the reference value, while near the free-molecular regime the geometric average (complete accommodation) should be chosen. In the intermediate region, either the arithmetic or geometric temperature can be used. As far as the heat flux is concerned, the use of the McCormack model is acceptable, with a maximum error of about 9% even when the temperature ratio of the two plates is 3. Table 7 The heat flux q y in a Ne-Ar gas mixture with C0 ¼ 0:5 for different wall-temperature ratios.