#  Accurate and efficient computation of the Boltzmann equation for Kramer’s
problem

Wei Su11footnotemark: 1222Wei Su and Peng Wang contribute equally Peng
Wang33footnotemark: 344footnotemark: 4 Haihu Liu 55footnotemark: 5 Lei
Wu66footnotemark: 6<EMAIL_ADDRESS>James Weir Fluids Laboratory,
Department of Mechanical and Aerospace Engineering, University of Strathclyde,
Glasgow G1 1XJ, UK  School of Energy and Power Engineering, Xi’an Jiaotong
University, 28 West Xianning Road, Xi’an 710049, China

###### Abstract

The Kramer’s problem is one of the fundamental problems of rarefied gas
dynamics, which has been investigated extensively based on the linearized
Boltzmann equation (LBE) of hard-sphere molecules and simplified kinetic model
equations. However, how the different intermolecular potentials affect the
viscous slip coefficient and the structure of Knudsen layer remains unclear.
Here, a novel synthetic iteration scheme (SIS) is developed for the LBE to
find solutions to Kramer’s problem accurately and efficiently: the velocity
distribution function is first solved by the conventional iterative scheme,
then it is modified such that in each iteration i) the flow velocity is guided
by an ordinary differential equation that is asymptotic-preserving at the
Navier-Stokes limit and ii) the shear stress is equal to the average shear
stress. Based on the Bhatnagar-Gross-Krook model, the SIS is assessed to be
efficient and accurate. Then we investigate the Kramer’s problem for gases
interacting through the inverse power-law, shielded Coulomb, and Lennard-Jones
potentials, subject to diffuse-specular and Cercignani-Lampis gas-surface
boundary conditions. When the tangential momentum accommodation coefficient
(TMAC) is not larger than one, the Knudsen layer function is strongly affected
by the potential, where its value and width increase with the effective
viscosity index of gas molecules. Moreover, the Knudsen layer function
exhibits similarities among different values of TMAC when the intermolecular
potential is fixed. For Cercignani-Lampis boundary condition with TMAC larger
than one, both the viscous slip coefficient and Knudsen layer function are
affected by the intermolecular potential, especially when the “backward”
scattering limit is approached. With the asymptotic theory by Jiang and Luo
(J. Comput. Phys., vol. 316, 2016, pp. 416–434) for the singular behavior of
the velocity gradient in the vicinity of the solid surface, we find that the
whole Knudsen layer function can be well fitted by the power series
$\sum_{n=0}^{2}\sum_{m=0}^{2}c_{n,m}x^{n}(x\ln x)^{m}$, where $x$ is the
distance to the solid surface. Finally, the experimental data of the Knudsen
layer profile are explained by the LBE solution with proper values of the
viscosity index and TMAC.

††journal: Journal of computational physics

##  1 Introduction

The Kramer’s problem is fundamental to solutions of almost all momentum
transfer problems of rarefied gas dynamics [1]. As illustrated in Fig. 1, when
the planar wall at $x_{2}=0$ moves slowly in the horizontal direction, a
nonlinear velocity profile and a finite slip velocity develop near the
surface. This kinetic boundary layer, as known as the Knudsen layer, has a
thickness of several mean free path $\lambda$ of gas molecules. Due to the
infrequent gas-gas interactions, the flow is essentially rarefied, so that the
conventional Navier-Stokes equations only work in the bulk region but break
down in the Knudsen layer. The linearized Boltzmann equation (LBE) can be used
to study this problem. For engineering applications, however, Navier-Stokes
equations are still preferred when the Knudsen number $Kn$ (the ratio of
$\lambda$ to the dimension of flow domain) is small, due to its distinct
computational advantage over the LBE in six-dimensional phase space.

Figure 1: Schematic diagram of the Knudsen layer in the Kramer’s problem. The
velocity defect (Knudsen layer function) $U_{s}$ describes the deviation of
the linearly extrapolated velocity (dash line) in the bulk region from the
true velocity (solid line). The velocity slope in the bulk region is
$(\text{d}U_{1}/\text{d}x_{2})|_{x_{2}\rightarrow\infty}=:k_{1}$, the slip
length is $\zeta$, while the viscous slip coefficient is defined as
$\bar{\zeta}=\zeta/\lambda_{e}$, where $\lambda_{e}$ is the equivalent mean
free path of gas molecules.

Many efforts have been made to predict the rarefied gas flow, through
incorporating the rarefaction effects caused by the presence of solid surface
into hydrodynamic equations [2]. For isothermal flows at small $Kn$, it is
adequate to apply the velocity slip boundary condition (BC) to Navier-Stokes
equations. In this case, the viscous slip coefficient ($\bar{\zeta}$, VSC), as
defined in Fig. 1, is needed. The first estimation
$\bar{\zeta}(\alpha)=(2-\alpha)/{\alpha}$ is proposed by Maxwell using
insightful physical arguments [3]. Here, $\alpha$ is the tangential momentum
accommodation coefficient (TMAC) describing the fraction of diffusely
reflected molecules at the solid surface, while the rest of molecules are
reflected specularly. Almost one hundred years later, using a variational
approach for the LBE and diffuse-specular gas-surface BC, Loyalka obtained the
VSC which is generalized into the following form [4, 5]:

$\bar{\zeta}(\alpha)=\frac{2-\alpha}{\alpha}\left[\bar{\zeta}(1)-0.1211(1-\alpha)\right].$
(1)

Subsequently, lots of investigations were performed to calculate the VSC by
numerically solving the LBE and its simplified model equations. It is found
that VSCs from the LBE and its kinetic model equations have a relative
difference less than $3\%$ when the effective TMAC is fixed [5]. It should be
noted that, although the influence of intermolecular potentials on VSC has
been assessed by Loyalka using the variational results [6] and by Sharipov
comparing the results from different literatures [7], a systematic
investigation of the role of the intermolecular potential on the VSC under
different gas-surface interactions on the basis of the highly accurate
Boltzmann solutions is still absent.

When $Kn$ becomes appreciable, Navier-Stokes equations may be still used, but
in addition to the velocity slip BC the viscosity is modified to be a function
of the distance to solid surfaces. In this case, the structure of the Knudsen
layer provides a critical information to formulate the effective viscosity.
Lockerby et al. first proposed a curve-fitted approximation to the Knudsen
layer function (KLF) as [8]

$U_{s}\left(x\right)\approx\frac{7}{20(1+x)^{2}},$ (2)

where $x$ is the distance to the solid surface normalized by the mean free
path $\lambda$. Although the KLF is fitted from a temperature jump problem
instead of the shear problem, it is found that the Navier-Stokes equations
with the effective viscosity can predict the velocity profiles in Poiseuille
and Couette flows, up to $Kn=0.4$. Later, by fitting the data from the LBE
solution of hard-sphere (HS) gas and the direct simulation Monte Carlo method
for Couette flow [9], Lilley $\&$ Sader obtained a power-law KLF [10, 11]:

$U_{s}\left(x\right)=U_{s}(0)-Cx^{n},$ (3)

where $C$ is a constant and the exponent $n\approx 0.82$. Although the fitting
is carried out in the region $0.1{\lesssim}x\lesssim 1$, they predicted the
power-law divergence of the velocity gradient in the vicinity of the solid
surface, that is, $dU_{s}/dx\rightarrow\infty$ as $x\rightarrow 0$.

The singular behavior of the velocity gradient at the planar surface is
rigorously proved by Takata & Funagane [12], when analyzing the thermal
transpiration based on the LBE of HS molecules. However, instead of the power-
law divergence, they found the logarithmic divergence of the velocity
gradient; that is, the spatial singularity is not stronger than $\ln{}x$ in
the vicinity of the solid surface. This conclusion is confirmed by Jiang & Luo
who, through the asymptotic analysis of the Bhatnagar-Gross-Krook (BGK) model
[13], found that the velocity profile of Couette flow near the solid surface
can be described by the following power series [22]

$U_{s}(x)=\sum_{n=0}^{N}\sum_{m=0}^{M}c_{n,m}x^{n}(x\ln{x})^{m},\quad{}x\rightarrow
0.$ (4)

It should be noted that most contributions to the Kramer’s problem focused
mainly on the diffuse-specular BC and HS molecules (or simplified kinetic
models). How intermolecular potentials (such as the inverse power-law,
shielded Coulomb, and Lennard-Jones potentials) and other gas-kinetic BCs
affect the VSC and KLF remains unclear. This paper is dedicated to addressing
these questions through the numerical simulation of the LBE. We emphasis that,
however, the numerical method to finding the KLF at small values of $Kn$ is
not easy. For instance, the results provided by Takata & Funagane are limited
to $Kn\gtrsim 0.6$, since the computational cost to find the steady-state
solution of the kinetic equations becomes extremely large for small Knudsen
numbers [12]; however, this relative large value of $Kn$ is unfortunately not
small enough to avoid the interference between Knudsen layers. In the present
paper, we first develop an efficient and accurate method to solve the LBE, and
then investigate the role of intermolecular potentials and gas-surface BCs on
the VSC and KLF.

The remainder of the paper is organized as follows. In § 2, the LBE for the
steady Couette flow of a monatomic gas and various kinetic BCs are introduced.
In § 3, a synthetic iteration scheme is developed to boost the convergence in
finding the steady-state solution of the Couette flow in the near-continuum
regime. In § 4, influences of intermolecular potentials and gas-kinetic BCs on
the VSC and KLF as well as the singularity of the velocity gradient near the
solid surface and the similarity of the KLF are investigated. In § 5,
experimental results given by Reynolds et al. [14] are properly explained. The
paper closes with some finial comments in § 6.

##  2 The linearized Boltzmann equation

Consider the steady Couette flow of a monatomic gas between two infinite
parallel plates located at $x_{2}=0$ and $x_{2}=1$. The top plate moves along
the $x_{1}$ direction with the velocity $V_{w}$, while the bottom plate moves
with the opposite velocity. Both plates are maintained at a fixed temperature
$T_{w}$. This Couette flow can be used to study the Kramer’s problem, provided
that the distance between the two plates is large enough so that there is no
interference between the Knudsen layers near the two plates [7].

When $V_{w}$ is far smaller than the most probable speed
(${v_{m}}=\sqrt{2k_{B}T_{w}/m}$, where $k_{B}$ is the Boltzmann constant and
$m$ is the gas molecular mass) of the gas molecules, the velocity distribution
function of gas molecules can be linearized around the global equilibrium
distribution function
$f_{eq}(\textbf{v})={\pi^{-3/2}}{\exp(-|\textbf{v}|^{2})}$ as:

$f(x_{2},\textbf{v})=f_{eq}(\textbf{v})+\frac{V_{w}}{v_{m}}h(x_{2},\textbf{v}),$
(5)

where $\textbf{v}=(v_{1},v_{2},v_{3})$ is the molecular velocity and
$h(x_{2},\textbf{v})V_{w}/v_{m}$ is the small perturbance ($h$ is not
necessary smaller compared to $f_{eq}$). The LBE for $h(x_{2},\textbf{v})$ is:

$v_{2}\frac{\partial{h}}{\partial{x_{2}}}=L(h,f_{eq}),$ (6)

where the linearized Boltzmann collision operator is [15]:

$L=\underbrace{\iint
B(\theta,|\textbf{u}|)[f_{eq}(\textbf{v}^{\prime})h({\textbf{v}}^{\prime}_{\ast})+f_{eq}(\textbf{v}^{\prime}_{\ast})h({\textbf{v}}^{\prime})-f_{eq}(\textbf{v})h({\textbf{v}}_{\ast})]d\Omega
d{\textbf{v}}_{\ast}}_{L^{+}}-\nu_{eq}(\textbf{v})h(\textbf{v}),$ (7)

and the equilibrium collision frequency is

$\nu_{eq}(\textbf{v})=\iint
B(|\textbf{v}-\textbf{v}_{\ast}|,\theta)f_{eq}(\textbf{v}_{\ast})d\Omega{}d{\textbf{v}}_{\ast}.$
(8)

Note that in the above LBE, the coordinate $x_{2}$ has been normalized by the
distance between the two plates $H$, the molecular velocity v has been
normalized by the most probable speed $v_{m}$, and velocity distribution
functions $f_{eq}$ and $h$ have been normalized by $n_{0}/v_{m}^{3}$, where
$n_{0}$ is the average number density of the gas molecules between the
parallel plates. The relative velocity of the two molecules before binary
collision is $\textbf{u}=\textbf{v}-\textbf{v}_{\ast}$, and $\Omega$ is a unit
vector along the relative post-collision velocity
$\textbf{v}^{\prime}-\textbf{v}^{\prime}_{\ast}$. The deflection angle
$\theta$ between the pre- and post-collision relative velocities satisfies
$\cos\theta=\Omega\cdot{\textbf{u}}/|\textbf{u}|$, $0\leq\theta\leq\pi$.
Finally, $B(\theta,|\textbf{u}|)=|\textbf{u}|\sigma$ is the collision kernel,
with $\sigma$ being the differential cross-section that is determined by the
intermolecular potential.

In the present paper, we consider the inverse power-law potentials, where the
collision kernels are modeled as [16, 15]

$B(|{\textbf{u}}|,\theta)=\frac{|{\textbf{u}}|^{2(1-\omega)}}{K}{\sin^{\frac{1}{2}-\omega}\left(\frac{\theta}{2}\right)\cos^{\frac{1}{2}-\omega}\left(\frac{\theta}{2}\right)},$
(9)

with $\omega$ being the viscosity index (i.e. the shear viscosity $\mu$ of the
gas is proportional to $T^{\omega}$) and $K$ some normalization constants
[15]. HS and Maxwell molecules have $\omega=0.5$ and 1, respectively. Note
that this type of collision kernel cannot describe the charged molecules
interacting through the Coulomb potential with $\omega=2.5$ [16]. As discussed
in the Chapter 10 of Ref. [17], in reality, however, charged molecules
interact through the shielded Coulomb potential:

$U^{\prime}(\rho^{\prime})=\epsilon\frac{\lambda_{d}}{\rho^{\prime}}\exp\left(-\frac{\rho^{\prime}}{\lambda_{d}}\right),$
(10)

where $\epsilon$ is related to the strength of the potential, $\rho^{\prime}$
is the intermolecular distance, and $\lambda_{d}$ is the Debye shielding
length. For simplicity, we only consider the single-species charged molecules
interacting through the repulsive force. We also consider noble gases
interacting through the following Lennard-Jones potentials:

$U^{\prime}(\rho^{\prime})=4\epsilon\left[\left(\frac{d}{\rho^{\prime}}\right)^{12}-\left(\frac{d}{\rho^{\prime}}\right)^{6}\right],$
(11)

where $d$ is the distance at which the potential is zero. The details of
implementation of the Lennard-Jones potentials in FSM can be found in Ref.
[18].

The differential cross-section for the above shielded Coulomb and Lennard-
Jones potentials can be calculated according to Sharipov & Bertoldo [19]. Then
the linearized Boltzmann collision operator Eq. (7) can be solved by the fast
spectral method developed by the authors [18].

To fully determine the gas dynamics in spatially-inhomogeneous problems, the
gas-surface BC is needed. The general form of the BC, which specifies the
relation between the velocity distribution function $f(\textbf{v})$ of the
reflected and incident gas molecules at the solid surface, is given below:

$v_{n}f(\textbf{v})=\int_{v_{n}^{\prime}<0}|v_{n}^{\prime}|R(\textbf{v}^{\prime}\rightarrow\textbf{v})f(\textbf{v}^{\prime})d\textbf{v}^{\prime},\quad
v_{n}>0,$ (12)

where $\textbf{v}^{\prime}$ and v are velocities of the incident and reflected
molecules, respectively, $v_{n}$ is the normal component of the molecular
velocity v directed into the gas, and
$R(\textbf{v}^{\prime}\rightarrow\textbf{v})$ is the non-negative scattering
kernel.

The most popular gas-surface BC is the diffuse-specular one, with the
scattering kernel reading:

$R_{M}(\textbf{v}^{\prime}\rightarrow\textbf{v})=\alpha_{M}\frac{m^{2}{v}_{n}}{2\pi(kT_{w})^{2}}\exp\left(-\frac{m\textbf{v}^{2}}{2kT_{w}}\right)+\left(1-\alpha_{M}\right)\delta\left(\textbf{v}^{\prime}-\textbf{v}+2\textbf{n}{v}_{n}\right),$
(13)

where the constant $\alpha_{M}$ is the TMAC, with a value in the range of
$0\leq\alpha_{M}\leq 1$, and $\delta$ is the Dirac delta function. Purely
diffuse reflection has $\alpha_{M}=1$. The BC proposed by Cercignani & Lampis
[20] has also been widely used, which reads:

$\displaystyle R_{CL}\left(\textbf{v}^{\prime}\rightarrow\textbf{v}\right)=$
$\displaystyle\frac{m^{2}v_{n}}{2\pi\alpha_{n}\alpha_{t}\left(2-\alpha_{t}\right)\left(kT_{w}\right)^{2}}I_{0}\left(\frac{\sqrt{1-\alpha_{n}}mv_{n}v_{n}^{\prime}}{\alpha_{n}kT_{w}}\right)$
(14)
$\displaystyle\times\exp\left\\{-\frac{m\left[v_{n}^{2}+(1-\alpha_{n})v_{n}^{\prime}\right]^{2}}{2kT_{w}\alpha_{n}}-\frac{m\left[\textbf{v}_{t}-(1-\alpha_{t})\textbf{v}_{t}^{\prime}\right]^{2}}{2kT_{w}\alpha_{t}(2-\alpha_{t})}\right\\},$

where $\textbf{v}_{t}$ is the tangential velocity,
$I_{0}(x)=\int_{0}^{2\pi}\exp\left(x\cos\phi\right)\text{d}\phi/2\pi$, and
$\alpha_{n}\in[0,1]$ and $\alpha_{t}\in[0,2]$ are the energy and momentum
accommodation coefficients, respectively. When $\alpha_{n}=\alpha_{t}=1$ or
$\alpha_{n}=\alpha_{t}=0$, the fully diffuse or specular BCs are recovered,
respectively, while for $\alpha_{n}=0$ and $\alpha_{t}=2$, the Cercignani-
Lampis scattering kernel descries “backward” scattering. Other types of BCs
have also been proposed and discussed [21], but for the Kramer’s problem, as
will be shown below, the two BCs are adequate to explain the experimental data
of Reynolds et al. [14].

The macroscopic quantities of interest are the flow velocity normalized by
$V_{w}$ and the shear stress normalized by $n_{0}k_{B}T_{w}V_{w}/v_{m}$, which
can be calculated as

$\displaystyle U_{1}=$ $\displaystyle\int{v_{1}h}d\textbf{v},$ (15)
$\displaystyle P_{12}=$ $\displaystyle\int{2v_{1}v_{2}h}d\textbf{v}.$

##  3 Numerical method: the synthetic iteration scheme

To resolve the singular behavior of the velocity gradient that occurs in the
vicinity of the solid surface [12, 22], high spatial resolution is required.
This means that it is better to solve kinetic equations by time-implicit
deterministic numerical method, otherwise the restriction on the Courant-
Friedrichs-Lewy condition will render the time step extremely small and hence
the computational cost enormous; also, the direct simulation Monte Carlo
method will be expensive to resolve the velocity profile in very small cells
near the solid surface.

To avoid the interference of the two Knudsen layers, the mean free path of gas
molecules should be sufficiently smaller than the distance between two
parallel plates; or equivalently, the rarefaction parameter

$\delta=\frac{H}{\lambda_{e}},~{}\lambda_{e}=\frac{\mu{(T_{w})}v_{m}}{n_{0}k_{B}T_{w}},$
(16)

should be sufficiently large. Here the equivalent mean free path
$\lambda_{e}=(2/\sqrt{\pi})\lambda$, where $\lambda$ is the mean free path of
the gas molecules. Therefore, the Knudsen number is $Kn=\sqrt{\pi}/2\delta$.

The integro-differential system Eq. (6) is usually solved by the conventional
iteration scheme. Given the value of $h^{(k)}(x_{2},\textbf{v})$ at the $k$-th
iteration step, the velocity distribution function at the next iteration step
is calculated by solving the following equation [9, 16, 18]:

$\nu_{eq}h^{(k+1)}+v_{2}\frac{\partial{h}^{(k+1)}}{\partial{x_{2}}}=L^{+}(h^{(k)},f_{eq}),$
(17)

where the derivative with respect to $x_{2}$ is usually approximated by a
second-order upwind finite difference, and the collision operator in Eq. (7)
can be calculated by the fast spectral method [15, 18] based on the velocity
distribution function at the $k$-th iteration step. The process is repeated
until relative differences between successive estimates of macroscopic
quantities are less than a convergence criterion $\epsilon$.

The conventional iteration scheme is efficient for highly rarefied gas flows
(when $\delta$ is very small), where converged solutions can be quickly found
after several iterations. However, the number of iteration increases
significantly when $Kn$ decreases (or $\delta$ increases), especially when the
gas flow is in the near-continuum regime [23, 24]. These behaviors are in fact
a result of the competition between the molecular collision and streaming. In
the free-molecular flow regime, gas molecules move in straight way (except the
collision with solid surfaces) so that any disturbance at one point can be
quickly felt by all other spatial points, so the exchange of “information” and
hence the convergence is fast. However, for near-continuum flows, binary
collisions dominate so that the exchange of information through streaming
becomes very inefficient: the perturbance decays rapidly due to frequent
binary collisions and takes a long time to be felt by other points. The
implicit unified gas-kinetic scheme may be used to achieve fast convergence
[25], however, currently there is no version developed for the LBE.

To have a convergence-accelerated scheme for the LBE, synthetic equations for
the evolution of macroscopic flow variables that are asymptotic preserving the
Navier-Stokes limit provides an alternative way to enhance the information
exchange across the whole computational domain [24]. For the Navier-Stokes
equations which can be derived from the Boltzmann equation through the
Chapman-Enskog expansion to the first-order of the Knudsen number, the
governing equation for the flow velocity is

$\frac{\partial{U_{1}}}{\partial{}x_{2}}=-\delta{P_{12}},$ (18)

where $P_{12}$ is a constant across the whole domain. For the LBE, the shear
stress remains a constant (this can be easily proven by multiplying Eq. (6)
with $v_{1}$ and then integrating with respect to v), but the equation for
$U_{1}$ contains high-order terms beyond the Navier-Stokes level. That is, the
governing equation is in general can be expressed as

$\frac{\partial{}U_{1}}{\partial{}x_{2}}=-\delta{P_{12}}+\text{High-order
terms}.$ (19)

To obtain the synthetic equation Eq. (19) that will facilitate the fast
convergence to the steady-state, we first rewrite Eq. (7) as
$L=(L-L_{BGK})+L_{BGK}$, where

$L_{BGK}=\delta[2U_{1}v_{1}f_{eq}-h]$ (20)

is the linearized collision operator of the BGK equation for Couette flow
between two parallel plates [24]. Then multiplying Eq. (6) by $2v_{1}v_{2}$
and integrating the resulting equation with respect to the molecular velocity
v, we obtain

$\displaystyle\frac{\partial U_{1}}{\partial
x_{2}}=-\delta{}P_{12}+\underbrace{\int{}2v_{1}v_{2}(L-L_{BGK})d\textbf{v}-\frac{\partial}{\partial
x_{2}}\int(2v_{2}^{2}-1)v_{1}fd\textbf{v}}_{\text{High-order terms}}.$ (21)

It is obvious that, in the near-continuum regime where
$\delta\rightarrow\infty$, the high-order terms in the right-hand side of Eq.
(21) are negligible compared to $\delta{}P_{12}$, so that the derived
synthetic equation is asymptotic preserving the Navier-Stokes limit. With this
macroscopic equation to update the flow velocity, we devise the following new
iteration scheme to find the steady-state solution of the LBE Eq. (6) quickly:

  * 1.

Due to the symmetry condition
$h(x_{2},v_{1},v_{2},v_{3})=-h(1-x_{2},v_{1},-v_{2},v_{3})$, the computational
domain will be limited to $0\leq{x_{2}}\leq{1/2}$. When $h^{(k)}$ and
$U_{1}^{(k)}$ are known at the $k$-th iteration, we calculate one of the high-
order terms $H_{1}(x_{2})=\int 2v_{1}v_{2}(L-L_{BGK})d\textbf{v}$. We also
calculate the velocity distribution function $h^{(k+1/2)}$ according to the
conventional iteration scheme Eq. (17), that is, we solve the following
equation:

${\nu_{eq}}h^{(k+1/2)}+v_{2}\frac{\partial{h}^{(k+1/2)}}{\partial{x_{2}}}=L^{+}(h^{(k)},f_{eq}),$
(22)

by a second-order upwind finite difference in the bulk and a first-order
upwind scheme at the solid surface [9].

  * 2.

From $h^{(k+1/2)}$, we calculate the flow velocity $U_{1}^{(k+1/2)}(x_{2})$,
the shear stress $P_{12}^{(k+1/2)}(x_{2})$, and one of the high-order terms
$H_{2}(x_{2})=\int(2v_{2}^{2}-1)v_{1}fd\textbf{v}$. We also calculate the
average shear stress as

$\bar{P}=2\int_{0}^{1/2}{}P_{12}^{(k+1/2)}dx_{2}.$ (23)

  * 3.

We obtain the flow velocity $U_{1}^{(k+1)}$ by solving Eq. (21) with the
symmetrical boundary condition $U_{1}(1/2)=0$, where $P_{12}$ is replaced by
$\bar{P}$. That is,

$\displaystyle\frac{\partial U_{1}^{(k+1)}(x_{2})}{\partial
x_{2}}=-\delta{}\bar{P}+H_{1}(x_{2})-\frac{\partial{H_{2}(x_{2})}}{x_{2}}.$
(24)

  * 4.

The velocity distribution function $h(x_{2},\textbf{v})$ is modified to
incorporate the change of the macroscopic flow velocity. Meanwhile, the shear
stress is adjusted to its mean value $\bar{P}$, for all spatial points. That
is,

$h^{(k+1)}=h^{(k+1/2)}+2\left(U_{1}^{(k+1)}-U_{1}^{(k+1/2)}\right)v_{1}f_{eq}+2\left(\bar{P}-P_{12}^{(k+1/2)}\right)v_{1}v_{2}f_{eq}.$
(25)

  * 5.

The above steps are repeated until convergence.

Since the gas kinetic equation is solved together with the macroscopic
equation Eq. (24) for flow velocity, the above scheme is called the synthetic
iterative scheme (SIS). Note that although the SIS has been widely applied to
the radiation transport processes [26] and rarefied gas flows driven by local
pressure, temperature, and concentration gradients to overcome the slow
convergence in the near-continuum flow regime [23, 27, 28], it is the first
time that the SIS is developed for the linearized Couette flow.

###  3.1 Numerical tests of efficiency and accuracy

Figure 2: Profiles of the flow velocity at different iteration steps obtained
from the conventional iteration scheme (a) and SIS (b), when $\delta=100$.
Circles show the converged solution. (c) The total iteration number needed to
obtain the converged solution as a function of the rarefaction parameter
$\delta$, where circles and squares are the results from the conventional
iteration scheme and SIS, respectively. The initial condition is
$h(x_{2},\textbf{v})=0$. The iteration is terminated when the maximum relative
difference in the flow velocity between two consecutive iterations is less
than $10^{-5}$.

Numerical simulations are carried out to assess the efficiency and accuracy of
the SIS. We consider the simple BGK kinetic model with the diffuse boundary
condition as it has recently been solved with high accuracy [22, 29, 30].

We first test the efficiency of the SIS. We choose the rarefaction parameter
$\delta=100$ and discretize the half spatial space into 50 even-spaced points.
The molecular velocity space $v_{1}$ and $v_{3}$ are discretized by the roots
of the physicists’ version of the fourth-order Hermite polynomial, while the
molecular velocity $v_{2}$ is truncated to $[-6,6]$ and approximated by the
non-uniform points [15, 31]:

$v_{2}=\frac{6}{(N_{v}-1)^{\imath}}(-N_{v}+1,-N_{v}+3,\cdots,{N_{v}-1})^{\imath},$
(26)

which is useful to capture the discontinuity in the velocity distribution
function near $v_{2}\sim 0$. In this test we take $\imath=3$.

Figure 2 compares the convergence history and speed (in terms of the number of
iterative steps to reach the converged solution) of the SIS to the
conventional iteration scheme. Starting from the initial guess
$h(x_{2},\textbf{v})=0$, the perturbance from the solid surface quickly
adjusts the flow velocity near the solid surface towards the surface velocity
in the conventional iteration scheme: from Fig. 2(a) we find that $U_{1}(0)$
is already very close to the final converged solution after 100 iterations.
However, due to the frequent binary collision, such a perturbance slowly
penetrates the bulk regime. Since the symmetry condition always guarantees
$U_{1}(1/2)=0$, a large number of iterations are needed to alter the velocity
profile between the solid surface and the center of the channel to be nearly
linear. This situation is completely changed in the SIS, where the flow
velocity is corrected to be nearly linear at each iteration according to the
synthetic equation Eq. (24), which can be approximated by
$\partial{U_{1}}/\partial{}x_{2}=-\delta\bar{P}$ when $\delta$ is large. Such
a macroscopic governing equation allows the efficient exchange of information,
and therefore fast convergence is realized in the whole computational domain,
see Fig. 2(b).

As far as the convergence speed is concerned, we see from Fig. 2(c) that, when
$\delta$ is small, i.e. in the free-molecular flow regime, the conventional
iteration scheme and SIS are as efficient as each other, where the converged
solutions are obtained within 5 iterations. As $\delta$ increases so that the
flow enters the transition and near-continuum regimes, the iteration number of
the conventional iteration scheme increases rapidly, while that of the SIS
quickly reaches the saturation number of about 20. At $\delta=200$, the SIS is
about 500 times more efficient than the conventional iteration scheme. The
gain of using SIS becomes larger and larger as $\delta$ further increases.

Table 1: Comparisons of the velocity at the solid surface $x_{2}=1$, the velocity derivative at the channel center $x_{2}=1/2$, and the shear stress between the results of Jiang & Luo [22] and SIS. The linearized BGK equation is used. At each value of $\delta$, the data of Jiang & Luo [22] and SIS are shown in the first and second rows, respectively.  $1/\delta$ | $U_{1}(1)/2$ | $dU_{1}(1/2)/2dx_{2}$ | $-P_{12}/4$
---|---|---|---
$0.003$ | 0.497891535 | 0.993939801 | 1.490909702$\times 10^{-3}$
| 0.497891548 | 0.993939827 | 1.490909741$\times 10^{-3}$
$0.01$ | 0.493069780 | 0.980081002 | 4.900405010$\times 10^{-3}$
| 0.493069792 | 0.980081024 | 4.900405118$\times 10^{-3}$
$0.1$ | 0.441224641 | 0.835285766 | 4.155607783$\times 10^{-2}$
| 0.441224646 | 0.835285536 | 4.155607809$\times 10^{-2}$
$1$ | 0.251861340 | 0.444228470 | 1.694625753$\times 10^{-1}$
| 0.251861372 | 0.444228442 | 1.694625700$\times 10^{-1}$
$10$ | 0.072922113 | 0.132195579 | 2.611624603$\times 10^{-1}$
| 0.072922127 | 0.132195588 | 2.611624596$\times 10^{-1}$
$100$ | 0.013430729 | 0.025200983 | 2.796682147$\times 10^{-1}$
| 0.013430736 | 0.025200817 | 2.796682147$\times 10^{-1}$

We then assess the accuracy of the SIS by comparing the solution of the
integral equation derived from the linearized BGK equation, which has the
accuracy of at least 12 significant digits [22]. In order to capture the
Knudsen layer near the solid surface, the spatial domain $0\leq{}x_{2}\leq
1/2$ is divided into $N_{s}$ nonuniform sections, with most of the discrete
points placed near the wall:

$x_{2}=(10-15s+6s^{2})s^{3},$ (27)

where $s=(0,1,\cdots,N_{s})/2N_{s}$. The size of the smallest section is
$1.2\times 10^{-9}$ when $N_{s}=1000$. The iterations terminate when the
maximum relative error in the flow velocity between two consecutive iterations

$\epsilon=\text{max}\left|\frac{U_{1}^{(k+1)}(x_{2})}{U_{1}^{(k)}(x_{2})}-1\right|$
(28)

is less than $10^{-10}$; the point $U_{1}(1/2)$ is excluded since the velocity
is always zero.

A comparison between the SIS and accurate results of [22] is tabulated in
Table 1 for the linearized Couette flow. The molecular velocity $v_{2}$ is
discretized according to Eq. (26) with $N_{v}=64$ and $\imath=5$, while in the
spatial discretization Eq. (27) we choose $N_{s}=500$. Clearly our SIS has an
accuracy of at least 6 significant digits. The accuracy can be further
increased when more refined velocity and spatial grids are used.

##  4 Numerical results of the linearized Boltzmann equation

Using the accurate and efficient SIS, the LBE is solved for different
intermolecular potentials, under different gas-surface BCs. In the numerical
simulation, we set the rarefaction parameter to be $\delta=100$, so that the
distance between two plates is about 100 times as large as the mean free path
of gas molecules; thus, the interference between the Knudsen layers near each
plate is avoided. The molecular velocity $v_{2}$ is discretized according to
Eq. (26) with $N_{v}=128$ and $\imath=5$, while $v_{1}$ and $v_{3}$ are
discretized by $32\times 32$ uniform grids in the range of $[-6,6]$; in the
spatial discretization we choose $N_{s}=500$ in Eq. (27). In the fast spectral
approximation of the linearized Boltzmann collision operator Eq. (7), the
integral with respect to the solid angle $\Omega$ is calculated by the Gauss-
Legendre quadrature with $M=8$, see equation (39) in Ref. [16]. All these
measures enable our results holding an accuracy of at least 6 significant
digits.

When the steady-state solution is obtained, the velocity profile in the bulk
region is linearly fitted by $U_{NS}=k_{1}x_{2}+k_{0}$ in the dimensionless
form, where $k_{0}$ and $k_{1}$ are coefficients from the least square
fitting. Then the KLF is calculated according to the following equation:

$U_{s}\left(\frac{x_{2}}{Kn}\right)=\frac{U_{NS}(x_{2})-U_{1}(x_{2})}{k_{1}Kn},$
(29)

and the VSC is calculated as

$\bar{\zeta}=-\frac{1+k_{0}}{\bar{P}}.$ (30)

In the numerical simulation, we find that $\delta=100$ is accurate enough to
recover the KLF and VSC, when compared to the solution of $\delta=1000$.
However, when $\delta=10$, that is, the distance between two plates is roughly
10 times of the mean free path, two Knudsen layers interact with each other,
which leads to an inaccurate KLF by using Eq. (29).

###  4.1  The viscous slip coefficient

Although a large number of VSCs have been computed from kinetic model
equations [7], very few data are available based on the LBE for various
intermolecular potentials, in particular of the highly accurate solutions. In
this section, we study how the intermolecular potentials (including the
inverse power-law, shielded Coulomb, and Lennard-Jones potentials) and gas-
kinetic BCs (including the diffuse-specular and Cercignani-Lampis BCs) affect
the Kramer’s problem.

Table 2:  The VSCs $\bar{\zeta}$ for HS, VHS with $\omega=0.81$, and Maxwell molecules under the diffuse-specular BC with different TMACs. LBE results of the present paper, Siewert [32], and Wakabayashi et al. [33] are denoted by a, b, and c, respectively.  $\alpha_{M}$ |  | HS | VHS | Max. | $\alpha_{M}$ |  | HS | VHS | Max.
---|---|---|---|---|---|---|---|---|---
$0.1$ | a | 17.04836 | 17.09319 | 17.12847 | $0.6$ | a | 2.215672 | 2.245067 | 2.267942
| b | 17.04780 | – | – |  | c | 2.214780 | – | –
| c | 17.00580 | – | – |  | b | 2.209300 | – | –
$0.2$ | a | 8.173130 | 8.214580 | 8.247129 | $0.7$ | a | 1.781936 | 1.808636 | 1.829366
| b | 8.172480 | – | – |  | b | 1.780980 | – | –
| c | 8.152400 | – | – |  | c | 1.776600 | – | –
$0.3$ | a | 5.206345 | 5.244574 | 5.274524 | $0.8$ | a | 1.453926 | 1.478044 | 1.496725
| b | 5.205630 | – | – |  | b | 1.452920 | – | –
| c | 5.192800 | – | – |  | c | 1.449400 | – | –
$0.4$ | a | 3.716862 | 3.752014 | 3.779490 | $0.9$ | a | 1.196466 | 1.218108 | 1.234829
| b | 3.716090 | – | – |  | b | 1.195400 | – | –
| c | 3.706900 | – | – |  | c | 1.192500 | – | –
$0.5$ | a | 2.818444 | 2.850653 | 2.875774 | $1.0$ | a | 0.988451 | 1.007717 | 1.022560
| b | 2.817610 | – | – |  | b | 0.987328 | – | –
| c | 2.810700 | – | – |  | c | 0.984900 | – | –

####  4.1.1 The influences of intermolecular potential and gas-kinetic BC

Table 2 tabulates the VSCs obtained from the LBE for HS, variable hard-sphere
(VHS) with $\omega=0.81$, and Maxwell molecules, when the diffuse-specular BC
of different TMACs is used. Results of Wakabayashi et al. [33] using a
discrete velocity method and Siewert [32] using a polynomial expansion
technique to solve the LBE for HS molecules are also listed for comparison. It
is noticed that the three groups of data agree well with each other,
especially the relative difference between our results and those of Siewert
[32] is less than $10^{-4}$. As expected, the VSC increases as the TMAC
decreases. Also, the VSC is insensitive to the intermolecular potential, which
only slightly increases with the viscosity index $\omega$, where the relative
difference between HS and Maxwell molecules is less than $4\%$. This results
confirm the statement in previous studies [5, 7, 34].

Table 3:  The VSC for HS, VHS with $\omega=0.81$, and Maxwell molecules under Cercignani-Lampis BC with different effective TMAC $\alpha_{t}$ and energy accommodation coefficient $\alpha_{n}$. $\alpha_{t}$ | $\omega$ | $\alpha_{n}=0.25$ | $0.5$ | $0.75$ | $1$
---|---|---|---|---|---
$0.25$ | 0.5 | 6.365427 | 6.343336 | 6.324267 | 6.307321
| 0.81 | 6.400178 | 6.372316 | 6.347638 | 6.325202
| 1 | 6.426786 | 6.394845 | 6.366237 | 6.339971
$0.5$ | 0.5 | 2.799516 | 2.785158 | 2.772602 | 2.761338
| 0.81 | 2.829277 | 2.811279 | 2.795092 | 2.780207
| 1 | 2.851688 | 2.831142 | 2.812430 | 2.795028
$0.75$ | 0.5 | 1.598122 | 1.591127 | 1.584932 | 1.579323
| 0.81 | 1.622629 | 1.613906 | 1.605945 | 1.598540
| 1 | 1.641138 | 1.631215 | 1.622031 | 1.613380
$1.0$ | 0.5 | 0.988451 | 0.988451 | 0.988451 | 0.988451
| 0.81 | 1.007717 | 1.007717 | 1.007717 | 1.007717
| 1 | 1.022560 | 1.022560 | 1.022560 | 1.022560
$1.25$ | 0.5 | 0.615670 | 0.622315 | 0.628343 | 0.633906
| 0.81 | 0.629985 | 0.638188 | 0.645891 | 0.653221
| 1 | 0.641382 | 0.650657 | 0.659514 | 0.668067
$1.5$ | 0.5 | 0.361248 | 0.374217 | 0.386121 | 0.397213
| 0.81 | 0.371198 | 0.387115 | 0.402275 | 0.416866
| 1 | 0.379386 | 0.397328 | 0.414721 | 0.431729
$1.75$ | 0.5 | 0.174178 | 0.193187 | 0.210840 | 0.227456
| 0.81 | 0.180631 | 0.203809 | 0.226193 | 0.247988
| 1.0 | 0.185886 | 0.211919 | 0.237535 | 0.262897
$2$ | 0.5 | 0.028851 | 0.053665 | 0.076984 | 0.099153
| 0.81 | 0.032881 | 0.062901 | 0.092298 | 0.121255
| 1.0 | 0.035527 | 0.069105 | 0.102637 | 0.136246

In order to study the Kramer’s problem with a more sophisticated gas-surface
interaction, the LBE is then solved with the Cercignani-Lampis BC. Results are
summarized in Table 3, for the effective TMAC $\alpha_{t}\in[0.25,2]$ and the
energy accommodation coefficient $\alpha_{n}\in[0.25,1]$. When the value of
$\alpha_{n}$ and the intermolecular potential are fixed, the VSC increases
rapidly when $\alpha_{t}$ decreases, which is consistent with that in the
diffuse-specular BC. The additional free parameter $\alpha_{n}$ in Cercignani-
Lampis BC introduces new interesting results. When $\alpha_{t}<1$, for a fixed
$\alpha_{t}$ and intermolecular potential, the VSC decreases slightly as
$\alpha_{n}$ increases, where the maximum drop in the VSC is less than $2\%$.
When $\alpha_{t}=1$, the Cercignani-Lampis BC is reduced to the fully diffuse
one in this problem, and the VSC does not vary with $\alpha_{n}$. When
$\alpha_{t}>1$, the variation of VSC on $\alpha_{n}$ reverses when compared to
that of $\alpha_{t}<1$; and it is strongly influenced by $\alpha_{n}$,
especially when $\alpha_{t}$ is large. For instance, for HS molecules at
$\alpha_{t}=2$, the VSC is increased by more than three times when
$\alpha_{n}$ changes from $0.25$ to $1$. For fixed $\alpha_{n}$ and
$\alpha_{t}$, the change in the VSC is insensitive to the intermolecular
potentials when $\alpha_{t}\lesssim 1.75$. However, when $\alpha_{t}$ is close
to two (i.e. the “backward” scattering), the influence of the intermolecular
potential becomes considerable. For example, when $\alpha_{t}=2$ and
$\alpha_{n}=1$, Maxwell molecules have a VSC that is about 37% higher than
that for HS molecules.

####  4.1.2 The viscous slip coefficient as a function of the effective TMAC

The variation of the VSC with respect to the effective TMAC $\alpha$ (for
diffuse-specular and Cercignani-Lampis BCs, $\alpha=\alpha_{M}$ and
$\alpha_{t}$, respectively) could be generalized to some simple expressions.
By considering the VSCs at the two limit ends of $\alpha=1$ and
$\alpha\rightarrow 0$, the VSC is fitted by Eq. (1) for the diffuse-specular
BC. For Cercignani-Lampis BC, Sharipov proposed a similar equation which is a
linear combination of the VSCs at $\alpha_{t}=1$ and $\alpha_{t}=2$ [35].
However, the estimation shows a large error when $\alpha_{t}\rightarrow 0$.
Here we construct a more accurate expression for the VSC with respect to the
effective TMAC $\alpha$.

We find from Tables 2 and 3 that the VSC can be fitted by a general function
as the one used by Lilley & Sader for both diffuse-specular and Cercignani-
Lampis BCs [11]:

$\bar{\zeta}(\alpha)=\frac{a}{\alpha}-b\alpha-c,$ (31)

where the fitting coefficients $a,~{}b$, and $c$ are shown in Table 4 for
typical inverse power-law intermolecular potentials. Fig. 3 shows that the
fitted curve (constructed from the data when $\alpha\geq 0.2$) can accurately
predict the VSC even in the limit $\alpha\rightarrow 0$. For instance, when
the TMAC is 0.05, relative differences between the fitted VSC and the LBE
solutions are less than $0.1\%$ for both diffuse-specular and Cercignani-
Lampis BCs.

Table 4: Fitting coefficients in Eq. (31) for HS, VHS with $\omega=0.81$, and Maxwell molecules under the diffuse-specular and Cercignani-Lampis BCs.  BC | $\alpha_{n}$ | $\omega$ | $a$ | $b\cdot 10$ | $c$
---|---|---|---|---|---
Diffuse-specular | n/a | 0.5 | 1.773 | 1.1660 | 0.6687
| n/a | 0.81 | 1.773 | 1.4270 | 0.6238
| n/a | 1 | 1.773 | 1.6370 | 0.5885
Cercignani-Lampis | $0.25$ | 0.5 | 1.774 | 0.7266 | 0.7127
|  | 0.81 | 1.775 | 0.8889 | 0.6781
|  | 1 | 1.776 | 1.0130 | 0.6516
Cercignani-Lampis | $0.5$ | 0.5 | 1.773 | 0.4772 | 0.7367
|  | 0.81 | 1.774 | 0.5867 | 0.7069
|  | 1 | 1.774 | 0.5742 | 0.6837
Cercignani-Lampis | $0.75$ | 0.5 | 1.772 | 0.2434 | 0.7597
|  | 0.81 | 1.773 | 0.2915 | 0.7358
|  | 1 | 1.773 | 0.3373 | 0.7164
Cercignani-Lampis | $1.0$ | 0.5 | 1.772 | 0.0218 | 0.7818
|  | 0.81 | 1.766 | 0.0543 | 0.7370
|  | 1 | 1.772 | 0.0000 | 0.7499
Figure 3:  The VSC as a function of the effective TMAC $\alpha$ for HS
molecules, when the (a) diffuse-specular BC and (b) Cercignani-Lampis BC with
$\alpha_{n}=0.25$ are used. Solid lines: the numerical fitting of Eq. (31)
using the data from the LBE solutions (circles). Squares: the data from the
LBE but no used for fitting. Note that other values of $\alpha_{n}$ and other
types of intermolecular potentials show a similar behavior.

###  4.2 The Knudsen layer function

####  4.2.1 The influence of the intermolecular potential

Figure 4:  KLFs from the LBE solutions for the inverse power-law potentials,
the Lennard-Jones potentials of Helium and Xenon when the gas temperature is
$300$K, and the shielded Coulomb potential of charged molecules when
$\epsilon=k_{B}T_{w}$. The result from BGK equation is also included for
comparison. The diffuse BC is used.

Figure 4 illustrates the KLFs obtained from the LBE with the diffuse BC, for
the inverse power-law potential with various values of viscosity index
$\omega$, the Lennard-Jones potential of helium and xenon, and the shielded
Coulomb potential. It is found that, for the inverse power-law potential, the
KLF increases with the viscosity index in the whole Knudsen layer. For the
Lennard-Jones potential, the KLF of xenon molecules is larger than that of
helium, but the results of both helium and xenon lie between those of HS and
Maxwell molecules. This is comprehensible because the effective viscosity
indexes of helium and xenon at a temperature of $300$K are 0.66 and 0.85 [36],
respectively. The KLF predicted by the BGK is even larger than that from the
Maxwell molecules, but is smaller than that of $\omega=1.5$ where the gas
molecules interact with soft potentials. The shielded Coulomb potential has
the largest KLF, since its effective viscosity is close to 2.5 [17].

Thus, contrary to the VSC whose value is insensitive to the intermolecular
potential, the KLF is strongly affected by the intermolecular potential. That
is, when the effective viscosity index increases, (i) the value of the KLF
increases, and (ii) the KLF decays more slowly, or equivalently, the Knudsen
layer becomes wider. For example, at the solid surface, the relative
difference between KLFs of Maxwell and HS molecules is approximate $20\%$, and
that between the shielded Coulomb and HS potentials reaches 60%. Relative
differences at distances one to two mean free path away from the solid surface
are even larger, say, when $x_{2}/\lambda=2$ the value of the KLF of the
shielded Coulomb potential is about 4 times of that of HS potential. On the
other hand, when $U_{s}$ is decreased to 0.01 of its value at the solid
surface, the corresponding distances to the solid surface for the HS and
Maxwell molecules are about $2.7\lambda$ and $3.5\lambda$, respectively.

####  4.2.2 The influence of the gas-kinetic BC

Figure 5: The KLF for HS molecules. (a) The diffuse-specular BC. Along the
arrow, $\alpha_{M}=0.2,0.4,0.6,0.8,1$. (b) The Cercignani-Lampis BC with
$\alpha_{n}=0.25$. Along the arrow, $\alpha_{t}=0.5,~{}1,~{}1.5$ and 2. (c)
The Cercignani-Lampis BC with $\alpha_{t}=2$. Along the arrow,
$\alpha_{n}=0.25,~{}0.5,~{}0.75$ and 1. Dots: LBE solutions. Solid lines:
fitted curves using Eq. (4) with $M,N=2$. Insets: zoomed regions in the
vicinity of the solid surface.

For the sake of clarity, we only focus on HS molecules. The diffuse-specular
BC is first considered in Fig. 5(a). It is found that the KLF increases as
$\alpha_{M}$ decreases. For a fixed $\alpha_{M}$, the KLF decreases rapidly
when moving away from the solid surface, say, its value decays by roughly
$85\%$ of the value on the solid surface when $x_{2}$ is about one mean free
path away from the solid surface.

Typical KLFs under Cercignani-Lampis BC are included in Figs. 5(b) and (c).
When $\alpha_{n}$ is fixed, for example, for $\alpha_{n}=0.25$, the KLF
decreases as $\alpha_{t}$ increases. The relative reduction in $U_{s}(0)$ is
about $40\%$ when $\alpha_{t}$ rises from 0.25 to 1. However, the variation of
the KLF with respect to $\alpha_{t}$ becomes weaken as $\alpha_{n}$ increases,
such that at $\alpha_{n}=1$ the reduction in $U_{s}(0)$ with $\alpha_{t}$
falls below $3\%$ (this is not visualized here but can be deducted from Table
6 below). When $\alpha_{t}$ ($\neq 1$) is fixed, the influence of $\alpha_{n}$
on the KLF becomes larger when $\alpha_{n}$ increases. And the greater the
TMAC $\alpha_{t}$ exceeds 1, the more pronounced the change in the KLF with
$\alpha_{n}$. As an example, when $\alpha_{t}=2$, the KLF is increased by
three times, as the $\alpha_{n}$ varies from 0.25 to 1, see Fig. 5(c).

Table 5:  Fitted coefficients corresponding to Eq. (4) with $M,N=2$ for the KLF obtained from the LBE with HS, VHS with $\omega=0.81$, and Maxwell molecules, when the diffuse-specular BC is used. $\omega$ | $\alpha_{M}$ | $c_{0,0}$ | $c_{0,1}$ | $c_{0,2}$ | $c_{1,0}$ | $c_{1,1}$ | $c_{1,2}\cdot 10$ | $c_{2,0}$ | $c_{2,1}$ | $c_{2,2}\cdot 10^{2}$
---|---|---|---|---|---|---|---|---|---|---
$0.5$ | 0.1 | 0.6502 | 1.2720 | -0.3624 | 1.4640 | 0.9950 | -0.6074 | -2.0090 | 0.1593 | 0.2246
| 0.2 | 0.6084 | 1.1760 | -0.3257 | 1.3260 | 0.9093 | -0.5358 | -1.8360 | 0.1410 | 0.1967
| 0.3 | 0.5675 | 1.0830 | -0.2915 | 1.1970 | 0.8283 | -0.4699 | -1.6720 | 0.1241 | 0.1712
| 0.4 | 0.5277 | 0.9941 | -0.2598 | 1.0760 | 0.7517 | -0.4095 | -1.5170 | 0.1085 | 0.1479
| 0.5 | 0.4889 | 0.9091 | -0.2302 | 0.9630 | 0.6792 | -0.3543 | -1.3710 | 0.0943 | 0.1267
| 0.6 | 0.4509 | 0.8277 | -0.2029 | 0.8570 | 0.6108 | -0.3040 | -1.2330 | 0.0813 | 0.1075
| 0.7 | 0.4139 | 0.7498 | -0.1776 | 0.7582 | 0.5463 | -0.2585 | -1.1030 | 0.0695 | 0.0902
| 0.8 | 0.3778 | 0.6752 | -0.1543 | 0.6661 | 0.4856 | -0.2173 | -0.9805 | 0.0588 | 0.0748
| 0.9 | 0.3424 | 0.6038 | -0.1328 | 0.5805 | 0.4284 | -0.1805 | -0.8652 | 0.0492 | 0.0610
| 1.0 | 0.3079 | 0.5356 | -0.1132 | 0.5012 | 0.3746 | -0.1480 | -0.7568 | 0.0406 | 0.0489
$0.81$ | 0.1 | 0.7436 | 1.6500 | -0.6657 | 2.3940 | 1.4840 | -1.3430 | -3.0070 | 0.3430 | 0.5216
| 0.2 | 0.6935 | 1.5130 | -0.5949 | 2.1550 | 1.3470 | -1.1840 | -2.7270 | 0.3030 | 0.4571
| 0.3 | 0.6450 | 1.3840 | -0.5294 | 1.9330 | 1.2180 | -1.0380 | -2.4640 | 0.2664 | 0.3983
| 0.4 | 0.5980 | 1.2610 | -0.4691 | 1.7270 | 1.0980 | -0.9046 | -2.2190 | 0.2329 | 0.3449
| 0.5 | 0.5523 | 1.1450 | -0.4135 | 1.5360 | 0.9858 | -0.7836 | -1.9890 | 0.2024 | 0.2966
| 0.6 | 0.5081 | 1.0360 | -0.3624 | 1.3590 | 0.8806 | -0.6739 | -1.7750 | 0.1747 | 0.2531
| 0.7 | 0.4651 | 0.9315 | -0.3155 | 1.1940 | 0.7823 | -0.5747 | -1.5750 | 0.1495 | 0.2139
| 0.8 | 0.4233 | 0.8331 | -0.2727 | 1.0430 | 0.6907 | -0.4855 | -1.3890 | 0.1269 | 0.1789
| 0.9 | 0.3827 | 0.7400 | -0.2336 | 0.9035 | 0.6052 | -0.4056 | -1.2160 | 0.1065 | 0.1478
| 1.0 | 0.3433 | 0.6519 | -0.1981 | 0.7753 | 0.5257 | -0.3345 | -1.0550 | 0.0883 | 0.1203
$1$ | 0.1 | 0.8123 | 1.9930 | -1.0180 | 3.3540 | 1.9460 | -2.3230 | -4.0160 | 0.5808 | 0.9405
| 0.2 | 0.7559 | 1.8180 | -0.9055 | 3.0050 | 1.7580 | -2.0430 | -3.6190 | 0.5118 | 0.8235
| 0.3 | 0.7016 | 1.6530 | -0.8023 | 2.6820 | 1.5820 | -1.7880 | -3.2510 | 0.4490 | 0.7174
| 0.4 | 0.6491 | 1.4980 | -0.7076 | 2.3840 | 1.4190 | -1.5570 | -2.9100 | 0.3918 | 0.6214
| 0.5 | 0.5983 | 1.3530 | -0.6210 | 2.1090 | 1.2670 | -1.3470 | -2.5930 | 0.3400 | 0.5349
| 0.6 | 0.5493 | 1.2170 | -0.5419 | 1.8570 | 1.1260 | -1.1580 | -2.3000 | 0.2931 | 0.4571
| 0.7 | 0.5019 | 1.0880 | -0.4699 | 1.6240 | 0.9955 | -0.9880 | -2.0290 | 0.2508 | 0.3874
| 0.8 | 0.4560 | 0.9682 | -0.4044 | 1.4120 | 0.8745 | -0.8354 | -1.7790 | 0.2127 | 0.3252
| 0.9 | 0.4116 | 0.8553 | -0.3450 | 1.2170 | 0.7625 | -0.6991 | -1.5480 | 0.1787 | 0.2700
| 1.0 | 0.3685 | 0.7495 | -0.2914 | 1.0390 | 0.6590 | -0.9010 | -1.3350 | 0.1483 | 0.2212
Table 6: Fitting coefficients corresponding to Eq. (4) with $M,N=2$ for the LBE solutions of the KLF, when HS molecules and Cercignani-Lampis BC are used. Since the KLF is independent of $\alpha_{n}$ when $\alpha_{t}=1$, only the fitting coefficients for $\alpha_{t}=1$ and $\alpha_{n}=0.25$ are tabulated.  $\alpha_{n}$ | $\alpha_{t}$ | $c_{0,0}$ | $c_{0,1}$ | $c_{0,2}$ | $c_{1,0}$ | $c_{1,1}$ | $c_{1,2}\times 10$ | $c_{2,0}$ | $c_{2,1}\times 10$ | $c_{2,2}\times 10^{3}$
---|---|---|---|---|---|---|---|---|---|---
$0.25$ | 0.25 | 0.4756 | 0.7342 | -0.1390 | 0.5628 | 0.4537 | -0.2119 | -0.9529 | 0.5551 | $0.7804$
| 0.5 | 0.4174 | 0.6648 | -0.1299 | 0.5416 | 0.4265 | -0.1888 | -0.8847 | 0.5012 | $0.6781$
| 0.75 | 0.3616 | 0.5988 | -0.1214 | 0.5210 | 0.4001 | -0.1676 | -0.8195 | 0.4520 | $0.5814$
| 1 | 0.3079 | 0.5356 | -0.1132 | 0.5012 | 0.3746 | -0.1476 | -0.7568 | 0.4056 | $0.4893$
| 1.25 | 0.2559 | 0.4746 | -0.1053 | 0.4823 | 0.3501 | -0.1283 | -0.6965 | 0.3607 | $0.4004$
| 1.5 | 0.2051 | 0.4151 | -0.0975 | 0.4643 | 0.3265 | -0.1089 | -0.6379 | 0.3156 | $0.3134$
| 1.75 | 0.1551 | 0.3563 | -0.0897 | 0.4473 | 0.3037 | -0.0891 | -0.5806 | 0.2687 | $0.2270$
| 2 | 0.1056 | 0.2980 | -0.0816 | 0.4316 | 0.2818 | -0.0683 | -0.5245 | 0.2188 | $0.1405$
$0.5$ | 0.25 | 0.4080 | 0.6424 | -0.1285 | 0.5242 | 0.4152 | -0.1848 | -0.8585 | 0.4961 | $0.6048$
| 0.5 | 0.3735 | 0.6048 | -0.1231 | 0.5159 | 0.4011 | -0.1714 | -0.8228 | 0.4635 | $0.5872$
| 0.75 | 0.3403 | 0.5694 | -0.1180 | 0.5084 | 0.3877 | -0.1592 | -0.7891 | 0.4338 | $0.5372$
| 1.25 | 0.2762 | 0.5023 | -0.1085 | 0.4942 | 0.3619 | -0.1363 | -0.7252 | 0.3781 | $0.4524$
| 1.5 | 0.2446 | 0.4689 | -0.1037 | 0.4871 | 0.3491 | -0.1247 | -0.6934 | 0.3498 | $0.3951$
| 1.75 | 0.2130 | 0.4345 | -0.0987 | 0.4797 | 0.3363 | -0.1124 | -0.6608 | 0.3198 | $0.3462$
| 2 | 0.1810 | 0.3987 | -0.0935 | 0.4719 | 0.3232 | -0.0992 | -0.6271 | 0.2873 | $0.2953$
$0.75$ | 0.25 | 0.3545 | 0.5836 | -0.1203 | 0.5099 | 0.3926 | -0.1647 | -0.8019 | 0.4479 | $0.5554$
| 0.5 | 0.3383 | 0.5658 | -0.1177 | 0.5061 | 0.3861 | -0.1582 | -0.7852 | 0.4319 | $0.5304$
| 0.75 | 0.3229 | 0.5502 | -0.1154 | 0.5034 | 0.3802 | -0.1526 | -0.7705 | 0.4182 | $0.5090$
| 1.25 | 0.2931 | 0.5211 | -0.1111 | 0.499 | 0.3691 | -0.1426 | -0.7433 | 0.3932 | $0.4698$
| 1.5 | 0.2780 | 0.5058 | -0.1088 | 0.4964 | 0.3634 | -0.1372 | -0.7290 | 0.3798 | $0.4489$
| 1.75 | 0.2625 | 0.4888 | -0.1063 | 0.4929 | 0.3572 | -0.1310 | -0.7130 | 0.3645 | $0.4251$
| 2 | 0.2462 | 0.4696 | -0.1035 | 0.4884 | 0.3504 | -0.1238 | -0.6951 | 0.3466 | $0.3979$
$1$ | 0.25 | 0.3093 | 0.5405 | -0.1138 | 0.504 | 0.3762 | -0.1499 | -0.7617 | 0.4109 | 0.4980
| 0.5 | 0.3083 | 0.5370 | -0.1134 | 0.5020 | 0.3751 | -0.1483 | -0.7582 | 0.4071 | 0.4918
| 0.75 | 0.3080 | 0.5357 | -0.1132 | 0.5013 | 0.3747 | -0.1477 | -0.7570 | 0.4058 | 0.4896
| 1.25 | 0.3079 | 0.5354 | -0.1132 | 0.5011 | 0.3746 | -0.1475 | -0.7567 | 0.4055 | 0.4890
| 1.5 | 0.3075 | 0.5342 | -0.1130 | 0.5004 | 0.3742 | -0.1470 | -0.7555 | 0.4042 | 0.4869
| 1.75 | 0.3066 | 0.5311 | -0.1127 | 0.4987 | 0.3732 | -0.1455 | -0.7524 | 0.4008 | 0.4813
| 2 | 0.3050 | 0.5255 | -0.1120 | 0.4955 | 0.3715 | -0.1430 | -0.7469 | 0.3949 | 0.4714

####  4.2.3 Fitting the Knudsen layer function and the singularity of velocity
gradient

The KLF is essential not only in determining the nonlinear constitution in the
Knudsen layer [8], but also in defining the singularity of the velocity
gradient near the solid surface. In a recent work based on the BGK model,
Jiang & Luo [22] have rigorously shown that the velocity near the solid
surface can be described by Eq. (4), whose gradient possesses a logarithmic
divergence. However, in their work, it was only numerically demonstrated that
the first four leading terms of Eq. (4) can capture the velocity profile in an
extremely small interval $0\leq{}x_{2}\leq 1.5\times 10^{-7}$.

In this work, based on the highly accurate results of the LBE, surprisingly,
we find that the entire KLF can be described by Eq. (4), provided that more
high-order terms are included. The associated fitting coefficients in Eq. (4)
with $M,N=2$ for the diffuse-specular and Cercignani-Lampis BCs are tabulated
in Tables 5 and 6, respectively. We note that when $\alpha_{n}$ is fixed, the
absolute value of the fitting coefficient decreases as $\alpha_{t}$ increases.
Meanwhile, the dependency of each fitting coefficient on $\alpha_{t}$ becomes
weaker and weaker as $\alpha_{n}$ increases. For instance, when $\alpha_{n}$
is increased from 0.25 to 1, the maximum relative difference in $c_{0,0}$ for
different $\alpha_{t}$ is reduced from $350\%$ to $2\%$. From the insets in
Fig. 5, we observe that the fitted curves agree quite well with the numerical
results. Note that Eq. (4) with $M,N=2$ can also describe the KLF very well
when the distance to the solid surface reaches $10\lambda$.

Figure 6: The absolute value of the velocity gradient $|\lambda
dU_{s}/dx_{2}|$ for HS (solid line), Maxwell (dash-dot line), and soft-
potential with $\omega=1.5$ (dash line) molecules, when the diffuse boundary
is used. Insets are the zoomed velocity gradient and the equilibrium collision
frequency $\nu_{eq}(0,v_{2},0)$ normalized by the rarefaction parameter
$\delta$.

Next, the singularity of the velocity gradient in the vicinity of the solid
surface is investigated through the deviation of Eq. (4) with respect to
$x_{2}$. This singularity is dominated by the term with $n=0$ and $m=1$ in Eq.
(4), that is, the velocity gradient near the solid surface is
$c_{0,1}\ln{}x_{2}$ [12, 22]. From Table 5, it is found that for a fixed TMAC,
$c_{0,1}$ increases with the viscosity index, indicating that the rate of
divergence is faster for the gas molecules with a larger value of the
viscosity index, see Fig. 6. However, this trend reverse at $x_{2}\approx
0.015\lambda$. This behavior is somehow related to the variation of the
equilibrium collision frequency $\nu_{eq}$. From the left inset in Fig. 6 we
see that, when the rarefaction parameter $\delta$ is fixed, $\nu_{eq}(0,0,0)$
increases with the viscosity index $\omega$, which means that the collision
frequency is larger for larger values of $\omega$, so that the gas approaches
to the equilibrium quicker and hence the velocity defect decreases faster.
Similarly, the velocity gradient at $x_{2}>0.015\lambda$ seems to be
proportional to $\nu_{eq}(0,v_{2}>3,0)$. It should be noted that, however,
this explanation is phenomenological; one may resort to the rigorous
mathematical analysis to have a deep understanding [12, 37]. When the
intermolecular potential is fixed, a smaller effective TMAC will produce a
larger velocity gradient near the solid surface, see Tables 5 and 6.

Figure 7: The rescaled KLF $U_{s}/U_{s}(x_{2}=0)$ for $\alpha_{M}=0.2$ and 1, when the HS molecules and the diffuse-specular BC are used. For clarity, results at other values of $\alpha_{M}$ are not shown. Inset: the relative difference ($R\%$) of $U_{s}/U_{s}(x_{2}=0)$ for various $\alpha_{M}$ compared to that of $\alpha_{M}=0.2$.  Figure 8: The rescaled KLF $U_{s}/U_{s}(x_{2}=0)$ of HS molecules when (a) $\alpha_{n}=0.25$, (b) $\alpha_{n}=0.5$, and (c) $\alpha_{n}=1$, when the Cercignani-Lampis BC with $\alpha_{t}=0.25,0.5,1.0,1.5$, and 2.0 are used. Inset: the relative difference ($R\%$) of $U_{s}/U_{s}(x_{2}=0)$ at various TMAC, when compared to that of $\alpha_{t}=0.25$.  Table 7:  Fitting coefficients of the rescaled KLF $U_{s}/U_{s}(0)$ for inverse power-law potentials with different values of the viscosity index $\omega$, when the diffuse BC is used.  $\omega$ | $c_{0,0}$ | $c_{0,1}$ | $c_{0,2}$ | $c_{1,0}$ | $c_{1,1}$ | $c_{1,2}\cdot 10$ | $c_{2,0}$ | $c_{2,1}$ | $c_{2,2}\cdot 10^{2}$
---|---|---|---|---|---|---|---|---|---
$0.5$ | 1.0000 | 1.739 | -0.3677 | 1.628 | 1.217 | -0.4794 | -2.458 | 0.1317 | 0.1589
$0.75$ | 1.0000 | 1.864 | -0.5256 | 2.113 | 1.462 | -0.8429 | -2.932 | 0.2245 | 0.2976
$1$ | 1.0000 | 2.034 | -0.7907 | 2.820 | 1.788 | -1.5680 | -3.623 | 0.4024 | 0.6003
$1.25$ | 1.0000 | 2.275 | -1.2390 | 3.864 | 2.217 | -2.9800 | -4.648 | 0.7380 | 1.2300
$1.5$ | 1.0000 | 2.807 | -2.5280 | 6.238 | 2.821 | -8.4950 | -6.999 | 1.9310 | 4.4710

####  4.2.4 The similarity of the Knudsen layer function

In the above section, the details of the KLFs under several specific TMACs
have been presented, which can serve as benchmark solution. In this section we
investigate the similarity in the structure of the Knudsen layer.

We first study the KLF normalized by its value on the solid surface $x_{2}=0$,
when the HS molecules and the diffuse-specular BC are used. Results of other
types of molecules are similar. Fig. 7 shows the rescaled KLF
$U_{s}/U_{s}(x_{2}=0)$ and their relative difference at different TMAC, when
compared with that at $\alpha_{M}=0.2$. We notice that the rescaled KLF for
$\alpha_{M}=0.2$ and 1 almost overlap; as shown in the inset of figure 7, the
maximum relative difference among all TMACs is less than $7\%$. Thus, the KLF
for diffuse-specular BC possesses a good similarity between different values
of TMAC.

For Cercignani-Lampis BC, as can be seen from Fig. 8, when $\alpha_{n}=1$, the
maximum relative discrepancy for all $\alpha_{t}$ is less than $10\%$. When
$\alpha_{n}$ decreases, however, the deviation of the rescaled KLF between
different $\alpha_{t}$ increases. For instance, when $\alpha_{n}=0.25$, the
maximum relative difference for $\alpha_{t}=2$ is about $30\%$, as compared
with $\alpha_{t}=0.25$, see the inset in Fig. 8. Nevertheless, it should be
noted that, for all $\alpha_{n}$ with $\alpha_{t}\leq 1$, the relative
difference of the rescaled KLF is less than $7\%$.

Approximately, the KLFs is defined to have similarity if the relative
difference of the rescaled KLF for different TMAC is less than $10\%$.
Therefore, as shown in Figs. 7 and 8, the KLF has the similarity when the
diffuse-specular BC and the Cercignani-Lampis BC with $\alpha_{n}=1$ are
considered, in the full range of the effective TMAC; for Cercignani-Lampis BC
with other values of $\alpha_{n}$, the similarity is preserved when
$\alpha_{t}\leq 1$.

Under the diffuse-specular BC, the rescaled KLF can be fitted using Eq. (4),
with the fitting coefficients for different intermolecular potentials
tabulated in Table 7. Furthermore, the corresponding KLF on the solid surface
$x_{2}=0$ can be fitted using an exponential function of the effective TMAC
$\alpha$ as

$U_{s}(x_{2}=0)=c_{1}\exp(-c_{2}\alpha)+c_{3},$ (32)

where $c_{1},c_{2}$ and $c_{3}$ are the fitting coefficients tabulated in
Table 8 for different intermolecular potentials. As a consequence, the KLF at
arbitrary TMAC can be roughly estimated by multiplying Eq. (32) and the
rescaled KLF, with the maximum relative error being smaller than $10\%$. The
KLF for the Cercignani-Lampis BC can also be rescaled according to the data in
Table 6.

Table 8:  Fitting coefficients of the velocity defect $U_{s}(0)$ at the solid surface by Eq. (31) for inverse power-law potentials with different values of the viscosity index $\omega$, when the diffuse-specular BC is used.  | $\omega$ | $c_{1}$ | $c_{2}$ | $c_{3}$
---|---|---|---|---
| 0.5 | 1.798 | 0.2410 | -1.1050
| 0.75 | 1.776 | 0.2825 | -1.0030
| 1 | 1.750 | 0.3378 | -0.8800
| 1.25 | 1.820 | 0.3822 | -0.8353
| 1.5 | 1.895 | 0.4368 | -0.7721

##  5 Comparison with the experiment

Figure 9:  Comparisons of the KLFs between the experiments and the LBE
solution with $\alpha_{t}=0.88$, when the air molecules of $\omega=0.75$ and
various $a_{n}$ are applied. Exp. A1 and Exp. A2 were measured by Reynolds et
al. [14]. The diffuse-specular BC with $\alpha_{M}=0.88$ is used in the BGK
model.

Reynolds et al. measured the VSC and KLF for air passing along the surface of
a highly polished aluminum plate [14] . They found that the KLF is different
from the results predicted by the BGK model. Loyalka pointed out that such
discrepancy is due to the deficiency of the BGK model [38], where the
collision frequency does not depend on the molecular velocity; by using a
kinetic model with a variable collision frequency, a reasonable agreement of
the velocity profile with the experimental data was observed. Given the
apparent deficiency of the model equation, results from the LBE of HS
molecules were also compared with the experimental data [9]. However, all the
previous works were based on the HS gas with an viscosity index of
$\omega=0.5$, while air has an effective viscosity index of 0.75 at the room
temperature. Moreover, the TMAC used in the numerical simulations was one,
which results in a VSC of about one, while that measured by Reynolds et al.
has an average value of $\bar{\zeta}_{Exp}=1.1$ (which has been corrected by
multiplying a factor of $\sqrt{\pi}/2$) [14] .

In this section, we try to explain the experimental data using the LBE
solutions for the inverse power-law potential with $\omega=0.75$. Although air
is a mixture of oxygen and nitrogen, we treat it as a single-species monatomic
gas, since (i) the molecular masses of oxygen and nitrogen are close to each
other and (ii) for isothermal flow the mass flow rate (and hence the VSC) is
insensitive to the rotational degrees of freedom [39, 40].

Figure 9 shows the KLF obtained from the LBE with $\alpha_{t}=0.95$, and
$\alpha_{n}=0.1$ and 1 under the Cercignani-Lampis BC, as well as the
experimental data. The result from the BGK equation is also included for
comparison. We use the value of TMAC $\alpha=0.95$, as our numerical
calculation in the previous section suggests that the predicted VSC from the
LBE agrees well with the experimental value of 1.1 [14]. It is found that the
KLF changes slightly under different $\alpha_{n}$, and the results of
$\alpha_{n}=1$ seems better than the others in the agreement with the
experimental data, while the solution of the BGK equation has a visible
deviation from the experimental results. Note that when using the diffuse-
specular BC, similar results can also be obtained for $\alpha_{M}=0.95$.

We note that the KLF from the experiments are scattered, which is inconsistent
with the theoretical analysis that the normalized velocity near the solid
surface should be independent of the mean free path and shear gradient.
Reynolds et al. argued that the most possible reason was the inaccurate
determination of the mean free path [14]. Therefore, intuitively, in order to
interpret the experimental results, one should take this factor into account.
To this end, we first assume the actual TMAC for the interaction of air with
the polished aluminum plate is $\alpha_{M}$ in the diffuse-specular BC. Then
we calculate the VSC $\bar{\zeta}(\alpha_{M})$ from the LBE. If
$\bar{\zeta}_{Exp}<\bar{\zeta}(\alpha_{M})$, the mean free path in the
experiment has been overestimated due to the inaccuracy in measuring the gas
pressure. Therefore, the value of the KLF from the experimental should be
multiplied by $1/\sigma=\bar{\zeta}(\alpha_{M})/\bar{\zeta}_{Exp}$, while the
width of the KLF should be stretched by a factor of $1/\sigma$. In the
numerical simulation, various values of $\alpha_{M}$ are attempted, until good
agreement between the results of experiment and numerical simulation are
achieved.

Figure 10:  Comparisons of the KLF between the experiments and the LBE
solution with $\alpha_{M}=0.8,~{}0.9$ and 1, when the air molecules with
$\omega=0.75$ are used. The LBE results are scaled by a factor
$\sigma=\bar{\zeta}_{Exp}/\bar{\zeta}(\alpha_{M})$, where
$\bar{\zeta}_{Exp}=1.1$ is the average VSC from the experiments. Exp. A1 and
Exp. A2 are measured by [14].

To show all the results in one figure, however, the KLF $U_{s}(x_{2})$
obtained from the numerical simulation of the LBE has been rescaled to
$\sigma{U_{s}}(\sigma{}x_{2})$. Comparisons between the numerical and
experimental results are depicted in Fig. 10. It is seen that, when the TMAC
varies from 0.8 to 1, the results of LBE can cover almost all the experimental
data. In other words, the TMAC of the aluminum plate used in the air
experiments is most likely $0.9\pm 0.1$. If the TMAC is 0.9, we have
$\sigma=0.9$, this means that the mean free path in the experiment is
overestimated by $10\%$, which seems reasonable due to the accuracy of the
micro-manometers at that time.

##  6 Conclusions

In summary, we have proposed a synthetic iteration scheme to expedite the
convergence of finding the steady-state solution of the linearized Boltzmann
equation for the Couette flow between two parallel plates. In the free
molecular flow regime, both the conventional and synthetic schemes lead to the
same converged solution after several iterations. However, the synthetic
iteration scheme converges significantly faster than the conventional one in
the transition and near-continuum gas flow regimes, which is about two to
three orders of magnitude faster than the conventional iterative scheme. Based
on the Bhatnagar-Gross-Krook kinetic model, the synthetic iteration scheme is
assessed to be accurate at least with six significant digits.

With this efficient and accurate method, the influences of the intermolecular
potentials (i.e. the inverse power-law, Lennard-Jones, and shielded Coulomb
potentials) and the kinetic boundary BCs on the Knudsen layer have been
investigated based on the linearized Boltzmann equation, where the Boltzmann
collision operator for general intermolecular potentials is solved by the fast
spectral method. Both the diffuse-specular and Cercignani-Lampis boundary
conditions are considered. It has been found that, although different
intermolecular potentials lead to roughly the same value of the viscous slip
coefficient, the KLF is strongly affected by the potential, whose value and
width increase with the effective viscosity index of the gas.

The highly accurate VSC and its general relation to the TMAC are presented for
different intermolecular potentials and gas-surface boundary conditions. In
addition, the KLF is found to be perfectly fitted by the series
$\sum_{n=0}^{2}\sum_{m=0}^{2}c_{n,m}x^{n}(x\ln{}x)^{m}$, where $x$ is the
distance to the solid surface. Correspondingly, based on the obtained KLF, the
macroscopic flow velocity gradient exhibits a logarithmic divergence on the
boundary. The strength of this divergence depends on the coefficient
$c_{0,1}$, whose value also increases with the viscosity index. Furthermore,
the similarity of the KLF has been established by rescaling the KLF by the
defect velocity at the solid surface. Consequently, the KLF at arbitrary TMAC
can be predicted by multiplying the rescaled KLF and the defect velocity at
the solid surface which is accurately fitted by an exponential function of the
TMAC. These results are useful to formulate the effective shear viscosity [8]
and slip boundary condition to be used in the framework of Navier-Stokes
equations [2].

The experimental data of the viscous slip coefficient and KLF measured by [14]
has been interpreted fairly well by the linearized Boltzmann equation with a
realistic viscosity index. We concluded that the TMAC for the interaction of
air with the polished aluminum is most likely $0.9\pm 0.1$, instead of $1.0$
as used in previous studies for a comparison with the experiment. This result
suggests that mean free path in the experiment has been overestimated by about
$10\%$.

Finally, it should be noted that the accurate and efficient synthetic
iterative scheme developed in this paper are readily to be extended to multi-
species gas mixtures [24]. The influence of intermolecular potentials and gas
kinetic boundary conditions on the Kramer’s problem of gas mixtures are
subject to future studies.

## Acknowledgments

This work is founded by joint project from the Royal Society of Edinburgh and
National Natural Science Foundation of China under Grant No. 51711530130, the
Carnegie Research Incentive Grant for the Universities in Scotland, and the
Engineering and Physical Sciences Research Council (EPSRC) in the UK under
grant EP/R041938/1.

## References

## References

  * [1] H. A. Kramers, On the behaviour of a gas near a wall, Il Nuovo Cimento (1943-1954) 6 (1949) 297–304.
  * [2] D. A. Lockerby, J. M. Reese, On the modelling of isothermal gas flows at the microscale, J. Fluid Mech. 604 (2008) 235–261.
  * [3] J. C. Maxwell, VII. On stresses in rarified gases arising from inequalities of temperature, Proc. Royal Soc. Lond. 170 (1879) 231–256.
  * [4] S. K. Loyalka, Momentum and temperature-slip coefficients with arbitrary accommodation at the surface, J. Chem. Phys. 48 (1968) 5432–5436.
  * [5] F. Sharipov, V. Seleznev, Data on internal rarefied gas flows, J. Phys. Chem. Ref. Data 27 (1998) 657–706.
  * [6] S. Loyalka, Slip and jump coefficients for rarefied gas flows: variational results for Lennard-Jones and n (r)-6 potentials, Physica A: Statistical Mechanics and its Applications 163 (3) (1990) 813–821.
  * [7] F. Sharipov, Data on the velocity slip and temperature jump on a gas-solid interface, Journal of Physical and Chemical Reference Data 40 (2) (2011) 023101\.
  * [8] D. A. Lockerby, J. M. Reese, M. A. Gallis, Capturing the Knudsen layer in continuum-fluid models of nonequilibrium gas flows, AIAA journal 43 (2005) 1391–1393.
  * [9] T. Ohwada, Y. Sone, K. Aoki, Numerical analysis of the shear and thermal creep flows of a rarefied gas over a plane wall on the basis of the linearized Boltzmann equation for hard-sphere molecules, Phys. Fluids A 1 (1989) 1588–1599.
  * [10] C. R. Lilley, J. E. Sader, Velocity gradient singularity and structure of the velocity profile in the Knudsen layer according to the Boltzmann equation, Phys. Rev. E 76 (2007) 026315.
  * [11] C. R. Lilley, J. E. Sader, Velocity profile in the Knudsen layer according to the Boltzmann equation, Proc. Royal Soc. A 464 (2008) 2015–2035.
  * [12] S. Takata, H. Funagane, Singular behaviour of a rarefied gas on a planar boundary, J. Fluid Mech. 717 (2013) 30–47.
  * [13] P. L. Bhatnagar, E. P. Gross, M. Krook, A model for collision processes in gases. I. Small amplitude processes in charged and neutral one-component systems, Phys. Rev. 94 (1954) 511.
  * [14] M. A. Reynolds, J. J. Smolderen, J. F. Wendt, Velocity profile measurements in the Knudsen layer for the Kramers problem., Tech. rep., In Proceedings of the Ninth International Symposium of Rarefied Gas Dynamics, M. Becker and M. Fiebig, eds., Dfvlr-Press (1974).
  * [15] L. Wu, J. M. Reese, Y. H. Zhang, Solving the Boltzmann equation by the fast spectral method: application to microflows, J. Fluid Mech. 746 (2014) 53–84.
  * [16] L. Wu, C. White, T. J. Scanlon, J. M. Reese, Y. H. Zhang, Deterministic numerical solutions of the Boltzmann equation using the fast spectral method, J. Comput. Phys. 250 (2013) 27–52.
  * [17] S. Chapman, T. G. Cowling, The Mathematical Theory of Non-uniform Gases, Cambridge University Press, 1970.
  * [18] L. Wu, H. H. Liu, Y. H. Zhang, J. M. Reese, Influence of intermolecular potentials on rarefied gas flows: Fast spectral solutions of the Boltzmann equation, Phys. Fluids 27 (2015) 082002.
  * [19] F. Sharipov, G. Bertoldo, Numerical solution of the linearized Boltzmann equation for an arbitrary intermolecular potential, J. Comput. Phys. 228 (2009) 3345–3357.
  * [20] C. Cercignani, M. Lampis, Kinetic models for gas-surface interactions, Transp. Theory Stat. Phys. 1 (1971) 101–114.
  * [21] L. Wu, H. Struchtrup, Assessment and development of the gas kinetic boundary condition for the Boltzmann equation, J. Fluid Mech. 823 (2017) 511–537.
  * [22] S. Jiang, L. S. Luo, Analysis and accurate numerical solutions of the integral equation derived from the linearized BGKW equation for the steady Couette flow, J. Comput. Phys. 316 (2016) 416–434.
  * [23] D. Valougeorgis, S. Naris, Acceleration schemes of the discrete velocity method: Gaseous flows in rectangular microchannels, SIAM J. Sci. Comput. 25 (2003) 534–552.
  * [24] L. Wu, J. Zhang, H. H. Liu, Y. H. Zhang, J. M. Reese, A fast iterative scheme for the linearized Boltzmann equation, J. Comput. Phys. 338 (2017) 431–451.
  * [25] Y. Zhu, C. Zhong, K. Xu, Implicit unified gas-kinetic scheme for steady state solutions in all flow regimes, J. Comput. Phys. 315 (2016) 16–38.
  * [26] M. L. Adams, E. W. Larsen, Fast iterative methods for discrete-ordinates particle transport calculations, Progr. Nucl. Energy 40 (2002) 3–159.
  * [27] S. Naris, D. Valougeorgis, D. Kalempa, F. Sharipov, Flow of gaseous mixtures through rectangular microchannels driven by pressure, temperature, and concentration gradients, Phys. Fluids 17 (2005) 100607.
  * [28] L. Szalmás, Accelerated discrete velocity method for axial-symmetric flows of gaseous mixtures as defined by the McCormack kinetic model, Comput. Phys. Commun. 184 (2013) 2430–2437.
  * [29] Y. W. Yap, J. E. Sader, High accuracy numerical solutions of the Boltzmann Bhatnagar-Gross-Krook equation for steady and oscillatory Couette flows, Phys. Fluids 24 (2012) 032004.
  * [30] W. Li, L.-S. Luo, J. Shen, Accurate solution and approximations of the linearized BGK equation for steady Couette flow, Computers & Fluids 111 (2015) 18–32.
  * [31] W. Su, S. Lindsay, H. H. Liu, L. Wu, Comparative study of the discrete velocity and lattice Boltzmann methods for rarefied gas flows through irregular channels, Phys. Rev. E 96 (2017) 023309.
  * [32] C. E. Siewert, The linearized Boltzmann equation: concise and accurate solutions to basic flow problems, Z. Angew. Math. Phys. 54 (2003) 273–303.
  * [33] M. Wakabayashi, T. Ohwada, F. Golse, Numerical analysis of the shear and thermal creep flows of a rarefied gas over the plane wall of a Maxwell-type boundary on the basis of the linearized Boltzmann equation for hard-sphere molecules, Eur. J. Mech. B-Fluids 15 (1996) 175–201.
  * [34] S. K. Loyalka, N. Petrellis, T. S. Storvick, Some numerical results for the BGK model: Thermal creep and viscous slip problems with arbitrary accomodation at the surface, Phys. Fluids 18 (9) (1975) 1094–1099.
  * [35] F. Sharipov, Application of the Cercignani–Lampis scattering kernel to calculations of rarefied gas flows. II. Slip and jump coefficients, Eur. J. Mech. B-Fluids 22 (2003) 133–143.
  * [36] G. A. Bird, Molecular Gas Dynamics and the Direct Simulation of Gas Flows, Oxford Science Publications, Oxford University Press Inc, New York, 1994.
  * [37] S. Takata, H. Funagane, Poiseuille and thermal transpiration flows of a highly rarefied gas: over-concentration in the velocity distribution function, J. Fluid Mech. 669 (2011) 242–259.
  * [38] S. K. Loyalka, Velocity profile in the Knudsen layer for the Kramer’s problem, Phys. Fluids 18 (1975) 1666–1669.
  * [39] L. Wu, C. White, T. J. Scanlon, J. M. Reese, Y. H. Zhang, A kinetic model of the Boltzmann equation for non-vibrating polyatomic gases, J. Fluid Mech. 763 (2015) 24–50.
  * [40] S. K. Loyalka, T. S. Storvick, Kinetic theory of thermal transpiration and mechanocaloric effect. iii. flow of a polyatomic gas between parallel plates, J. Chem. Phys. 71 (1979) 339–350.