# Wall morphology dependence of rare gas Knudsen diffusion in silica and graphite slit nanopore: A molecular dynamics study

Pierre Magnico

► To cite this version:

Pierre Magnico. Wall morphology dependence of rare gas Knudsen diffusion in silica and graphite slit nanopore: A molecular dynamics study. Vacuum, 2025, 242, pp.114756. 10.1016/j.vacuum.2025.114756. hal-05273657

HAL Id: hal-05273657

https://hal.science/hal-05273657v1

Submitted on 23 Sep 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Wall morphology dependence of rare gas Knudsen diffusion in silica and graphite slit nanopore: A molecular dynamics study

Pierre Magnico

<figure><img src="image_2.png" /></figure>

Aix Marseille Univ, CNRS, Centrale Med, M2P2, Marseille, France

## ARTICLE INFO

**Keywords:**

Knudsen diffusion

Accommodation coefficient

Surface diffusion

Molecular dynamics

Anisotropy

Porosity

Stochastic model

## ABSTRACT

Gas/wall collision mechanisms play a key role in Knudsen diffusion process. In particular, the channel wall structure has a major influence in mass transfer. So, we investigate the influence of the wall roughness, anisotropy and porosity on the self-diffusion of helium and neon in nanochannels. Three materials are proposed: graphite and β-cristobalite and amorphous silica. The study makes it possible to analyze, in function of temperature, the correlation between 1/the ballistic/diffusion transition regime of the surface gas transfer, 2/the transition of the bouncing process to a linear increase of the bounce number with time and 3/the shape of the surface residence time distribution characterized by a Fréchet like distribution at short time and an exponential decay at long time. As concerns the amorphous SiO₂, the bounce must be redefined owing to the transfer inside the material which is dominated by a cage effect. The anisotropy effect on collision process and Knudsen diffusion is analyzed by means of a tensorial computation of the tangential momentum accommodation coefficient and of the mean square displacement. Using the Langevin at the channel scale and the Arya model, the ballistic/diffusion transition time of the mean square displacement is related to the collision frequency and the collision number required for the velocity to be uncorrelated. A stochastic model confirms the molecular dynamics results with β-SiO₂ channel: The behavior of the Knudsen diffusion coefficient according to the Arrhenius law and the influence of collision frequency on transition time.

## 1. Introduction

In rarefied gas conditions, gas behavior is not governed by the con-
tinuum flow theories, such as the Navier-Stokes equation, but the
Boltzmann equation. In nanochannels, the mean free path length $\lambda$ is of
the same order of magnitude or larger than the channel width (h), so that
gas/wall collisions must be considered. The degree of rarefaction is
represented by the Knudsen number $Kn = \lambda/h$. In the free molecular
regime, i.e. $Kn > 10$, the molecule trajectory is a series of free flights
interrupted by collisions with channel walls. Therefore, gas/wall mo-
mentum and heat transfer inducing a velocity slip and a temperature
jump wall boundary are the key mechanisms controlling the mass
transfer in nanometric systems implemented in several industrial fields
such as gas storage [1,2], gas sensors [3-5] and advanced material
synthetizing [6,7]. In particular, silica and graphitic materials are often
used or encountered in a wide range of applications including heterog-
enous catalysis [8-12] and shale gas recovery [13-15], fuel cells [16,
17].

Several factors are involved in the gas/wall collision mechanism

such as wall topology characteristics (surface roughness, porosity,
anisotropy), or wall physico-chemical properties (active sites and sur-
face contaminants). The result is a global interaction potential well, the
depth of which depends on position along the wall surface. This
apparent roughness is not only induced by the local wall atom position,
but it is the result of several intricated contributions of different kinds.
Consequently, according to its kinetic energy, an atom of gas may be
adsorbed at the surface for a while during which the atom bounces
several times before desorbing [18,19]. A surface diffusion regime takes
place when the bounce number is high enough [20-22]. In the case of a
porous wall, the gas atom can alternatively move inside the wall and
bounce at the wall surface.

Owing to the complexity of the collision process, numerous empirical
models based on the statistical relation between the incident and the
reflected velocity are proposed like the Yamamoto model or the wash-
board model [23-27]. However, theories such that the Maxwell and
Cercignani-Lampis (CL) models [28,29], necessitating one or two pa-
rameters called accommodation coefficients, are the most popular due to
their simplicity. The first one assumes that a fraction of gas/wall
collision is diffusive, while the other fraction is specular. The second one uses the same principle but separates the transversal from the normal velocity component. As these fractions do not depend on the incident condition, these models have several drawbacks like the velocity component interplay for instance.

Recently, to take into account of wall roughness, authors [30] proposed a new scattering kernel based on the linear combination of the Maxwell model for the hard collision with a static rough wall and the CL model for the weak collision with a smooth wall in thermal equilibrium. The weighting coefficient is determined using the Debye-Waller factor. Considering the intrinsic limits of the Maxwell and the CL kernel, Liang et al. [31] chose the washboard model in which the corrugated surface is composed of vibrating cubes. At the impact point, the tangential component of the molecule's momentum with respect to the local surface orientation is conserved, and the normal component is changed via a hard elastic collision. To solve the problem of reciprocity violation, the authors made improvements which lead to a new model independent of unknown parameters in the case of light gas atoms colliding with heavy solid atoms.

Self-diffusion process in free molecular regime, first considered by Knudsen [32] assuming the diffuse collision only, was better described by Smoluchowski [33] using the Maxwell kernel. The partition into diffuse and specular collisions is controlled by the tangential momentum accommodation coefficient. To provide a more accurate description of the self-diffusion mechanism, Arya et al. [34] took account of the more realistic but simplified Cercignani-Lampis model and Liang et al. [35] improved Washboard model. However, in the two cases the pore walls are rigid, and the adsorption process is not considered when using the hard sphere potential. In addition, the simplification of the CL model lies in the diffuse reflection of the normal velocity component.

In the present work, we propose to analyze the wall heterogeneity effects on self-diffusion and gas/wall collision mechanisms by means of molecular dynamics (MD). Three kinds of heterogeneity are considered: roughness, anisotropy and porosity. In this aim, the silica is used in two forms: the crystalline form such as the β-cristobalite for anisotropy and the amorphous form for porosity. Roughness influence is examined by comparing the MD results with those obtained with the graphitic material. As adsorption has significant effect on mass and heat transfer, all the investigations are carried out over a wide range of temperature. More specifically, the gas atom bounce mechanisms are examined in the case of the porous silica. Owing to the wide characteristic time scale of the processes, three kinds of simulations are carried out: self-diffusion in slit channels, collision dynamics and inside-wall gas transfer.

Owing to the wide temperature range investigated, we have to consider that at small temperature, say $T < 100$ K, the gas/wall collision process (adsorption with trapping, desorption) is dominated by quantum effects [36–38]. In the case of light rare gases, theoretical and experimental studies are mostly carried out at very low temperatures ($T < 50$ K), especially as concerns graphitic substrate [39–41]. For instance, energy accommodation coefficient (EAC) of rare gas is investigated in Ref. [42] for several rare gases (from Ne to Xe). By quantum mechanical forced oscillator calculations, the authors show that the quantum and classical EAC of Ne are close if $T > 50$ K and equal at $T > 100$ K if the solid temperature is null. This is also observed in gas/wall thermal equilibrium with another expression of EAC including the sticking coefficient. However, EAC increases with $T$ as long as $T < 100$ K contrary to the other rare gases. Gorthel et al. [43] performed quantum computations mimicking the time-of-flight experiments of He desorbed from graphite, sapphire and Nichrome surfaces. Desorption mechanism is modeled by a cascade of one-phonon processes. The authors found that the time-of-flight spectra do not follow the Maxwell distribution at $T < 10$ K revealing desorbed atoms with an average energy higher than

that of a particle in thermal equilibrium with the wall. This was confirmed by Brivio et al. [44] by using the distorted-wave Born approximation. However, at larger temperature the angle-resolved distribution does not follow clearly the Knudsen's cosine law at $T > 25$ K with a maximum appearing at $45^\circ$. The authors observed the same trend with Ne/Cu and explained this observation by the increase of population of both closer to zero energy bound states and larger parallel momentum states. More recently, a quantum investigation of trapping/desorption [45] over a long period of time has been conducted with a thermal beam of He (35 K) colliding with a Cu surface ($T_s$ ranging from 100 K to 300 K). Using the reduced density matrix quantum approach [46], the interplay between the diffraction states and the inelastic phonon transition is analyzed. As concerns trapping, it increases with $T_s$ as the probability of phonon excitation increases and decreases with the increasing incident energy as more energy needs to be dissipated into the solid. In particular, trapped particles in DMSA resonance remain in these states for a long time before desorbing and can enhance sticking at low incident energy.

Few investigations of He/graphene interaction and motion in extremely narrow pores at $T \sim 5$ K have been carried out. Strong attractive interaction induces a much higher dense state of He in 0.7 nm slit pores than those in larger slit pores [47] and helium diffuses through a 0.7 nm carbon nanotube less easily than nitrogen molecule due to a larger effective size of He caused by a more pronounced quantum effect at low temperature [48,49].

The paper is organized as follows. Section 2, the simulation set-ups of the channel self-diffusion and of the collision methodology are described. Section 3, the first results of Knudsen diffusion in slit channel are analyzed. In particular, the effect of the anisotropy is discussed as well as the agreement with the Arrhenius law. Section 4, the collision process is investigated. First, anisotropy and roughness effects on the tangential momentum accommodation coefficient are examined in function of temperature. Then the adsorption process is analyzed with the residence time distribution and the mean square displacement. Finally, the gas transfer inside the amorphous silica is studied. Section 5, the transition of the Knudsen mass transfer from the ballistic to the diffusion regime is examined by means of the Langevin model coupled with the Arya model. As the Arya model fails to predict self-diffusion in $SiO_2$ channel, a stochastic model based on the residence time distribution is proposed in section 6. The temperature dependence of the self-diffusion and of the relaxation time of the mean square displacement are analyzed.

## 2. Computational methods

Two kinds of silica materials are used in this study: the β-cristobalite ($\beta$-$SiO_2$) and the amorphous $SiO_2$ (a-$SiO_2$). They are built initially from an ideal β-cristobalite with a volume of $50.2 \AA^3 \times 50.2 \AA^3 \times 50.2 \AA$. This silica is then equilibrated at 300 K during 1 ns using the Amber potential (Eq.SI.E1 of the Supplementary Information) which parameters are proposed by Emami et al. [50] (Table SI.T1-3). The equilibration is carried out with the Nosé-Hoover thermostat in the NVT ensemble and with a temperature damp factor of 1 ps. The time step is 1 fs.

To synthesize the amorphous material, we follow the methodology proposed by Dufrèche et al. [51]. The initial $\beta$-$SiO_2$ is heated first from 300 K to 6000 K during 1 ns. Then the molten silica is equilibrated at 6000 K during 1 ns. Finally, the a-$SiO_2$ is generated by cooling back to 300 K with a rate of 2.5 K/ps and is equilibrated at 300 K for 1 ns. During all these steps the modified BKS potential [52,53] is used instead of the Amber potential (Eq.SI.E2).

For the collision process and the self-diffusion investigation, the walls are created with a 14.35 Å thick slab cut from the block of $\beta$-$SiO_2$ or
a-SiO₂ along a longitudinal plane xy, the dimension being 50.2 Å in the directions x and y. The oxygens at the opposite surfaces of the β-SiO₂ wall are hydrogenated (Figs. SI.F1A(c/d)). As concerns the a-SiO₂ wall, the two surfaces must be reconstructed to avoid under-coordinates oxygen and silicon atoms. In this aim the methodology of Geske et al. [53] is used. It consists in eliminating the uncoordinated oxygen atoms and silicon atoms which have less than two bonds with oxygen atoms. Then the under-coordinated oxygen and silicon atom are hydrogenated and hydroxylated respectively. With this method, a density of 6.32 silanol/nm is reached. The silanol condensation is not undertaken (Figs. SI.F1B(a/b)).

Comparing with the ideal β-SiO₂, Fig.SI.F1A (c,d) shows that the atoms of the β-SiO₂ at 300 K remain distributed along a direction oriented at 45° to the x,y axes. Owing to the dihedral constraint, anisotropy is statistically preserved. In Fig.SI.F1B (a,b), the a-SiO₂ surface displays no anisotropy. The spatial distribution heterogeneity induces the presence of holes which reveals the presence of porosity inside the wall. In graphite walls, the thermal fluctuation of the carbon atoms position leads to a small roughness because the atoms remain distributed by layer. In the case of silica, the atoms are bonded to each other in the transversal direction and their position at the surface is dependent on the thermal redistribution of atoms in the bulk wall. This also increases the surface roughness. Despite the graphene periodic structure, anisotropy is not considered owing to the very low longitudinal dependence of the gas/wall interaction potential. In particular, this is observed for rare gases [54–56]. This low spatial dependence is amplified by the very high surface flatness.

For the slit channel, the hydrogenated slab is copied and translated in the transversal direction so that the distance (h) between the two inner channel surfaces is around 50 Å (Fig.SI.F1B(c)). It must be specified that the collisional method needs one wall only (Fig.SI.F1C). To impose a temperature of investigation T, the walls, i.e. the slabs, are first equilibrated at 300 K for 200 ps and then they are equilibrated for 800 ps at temperature T, letting the outer surface(s) fixed. The temperature equilibration is carried out with the initial Amber potential and parameters (Table SI.T1-3) in the NVT ensemble.

In the gas transfer investigation as in the temperature equilibrium procedure, the cutoff distance is 12 Å for the 12-6 Lennard Jones (LJ) potential and 8.5 Å for the electrostatic one and the Lorentz-Bertelot mixing rules are used. The LJ parameters of He and Ne (Table SI.T1) are those calculated and used in Refs. [57–59]. To reach the Knudsen condition, the gas/gas interaction is not considered. The O-H bond is constrained with the SHAKE algorithm, and the three-dimensional Ewald summation is used to solve the long-range electrostatics with a tolerance of 10⁻⁴. The time step is set to 1 fs. To minimize the wall/wall electrostatic interaction outside the channel, the dimension of the simulation box in the transversal direction z is around 300 Å for a longitudinal size of 50.16 Å in the x and y direction. Periodic boundary conditions are imposed in the three directions. All the simulations are performed with the large-scale atomic/molecular massively parallel simulator (LAMMPS) package [60].

In the following, the MD results with SiO₂ walls will be compared to those obtained with graphitic walls composed of several graphene layers of surface size 17.04 Å × 17.25 Å. Three layers were used for He [55] and five for Ne in the present work. The Adaptative Intermolecular Reactive Bond Atoms (AIREBO) and the 12-6 LJ potential are used for the interaction between the carbon atoms and for the gas/wall interaction respectively. The carbon LJ potential parameters are as follows:

$$ \epsilon_{CC} = 0.0556 \text{ kcal/mol} \quad \sigma_{CC} = 3.82 \text{ Å} $$

For the Knudsen diffusion in slit channel, the gas phase contains 2000 atoms in the case of SiO₂ walls, 600 and 1200 atoms in the case of

graphite walls with h = 34 Å and 50 Å respectively [55]. The trajectories are recorded during 100 ns every 1 ps after 1 ns gas/wall equilibration at an imposed temperature.

The collision mechanism is investigated with a system composed of one wall of size 50.16 Å × 50.16 Å × 14.35 Å. The gas atoms are thrown to the wall at z > rcut above the wall surface. After the collision with the wall, as soon as an outgoing atom crosses the control plane positioned at rcut above the wall, it is thrown back to the wall with a new imposed velocity (Fig.SI.F1C). The exact transversal position of the control plane is computed by determining each time the position of the highest Hydrogen atom. The gas phase is composed of 49 atoms. As the gas and the solid phase are in thermal equilibrium, the imposed incident velocity distribution is the Maxwell and the Rayleigh distribution for the longitudinal and the transversal components respectively. 2 × 10⁵ collisions are used for the statistics. The gas trajectory and velocity as the control plane position are recorded every 2 fs. The residence time distribution and the surface displacement are examined with the gas trajectory. The accommodation coefficients in Fig. 5 and SI.F.5 are analyzed with the incoming/outgoing velocity statistics over 5 × 10⁴ collisions without trajectory record. During the gas diffusion in the slit channel, the gas gets into the a-SiO₂ walls. Indeed, the record time step of 1 ps used in the self-diffusion simulations is too high to capture all the time scale of the gas transfer inside the wall. Therefore, the gas motion inside the wall is also analyzed with the collision methodology.

### 3. Knudsen self-diffusion

The SiO₂ wall thickness of 14.35 Å is chosen as a compromise between the simulation cpu time and the possibility that gas atoms may go through the walls. In the case of a-SiO₂, the number of He atoms going through the walls increases with temperature from 62 at 200 K to 183 at 300 K and no Ne atoms until 300 K. In the case of β-SiO₂, 14 He atoms cross the walls at 400 K. To examine the diffusional process, the trajectory of these atoms is not considered. For the system Ne/a-SiO₂, the simulations at 75 K and below show too many Ne atoms sticking inside the walls making the results unusable. This is also observed for Ne/β-SiO₂ at 50 K. At this temperature a high number of Ne atoms sticks at the wall surface inducing a non-linearity of mean square displacement. For He/a-SiO₂, the lower temperature bound is 75 K.

Let us first consider the atom transversal profiles in the SiO₂ wall channels. Figs. SI.F2A and SI.F2B display the number atom profiles of all the constituents of the wall channel. The wall profiles are independent of temperature. In these figures, they are represented at T = 300 K and are averaged over time. At the inner channel surfaces, the maximum of the hydroxyl profiles (O₅,H) is located at 16.9 Å and 65.9 Å and the OH layers lie over 2.5 Å. The hydroxyl groups form a superficial dense layer as the overall atom number (dotted black lines) presents a maximum at the same position. Using the maximum of the hydroxyl profile, h is around 49 Å. As concerns the a-SiO₂, the hydroxyl layers have a thickness of around 4.2 Å and their maximum is located at around 14 Å and 69.1 Å. However, the hydroxyl layer is as dense as the bulk wall. Its thickness is much higher than for the β-SiO₂ owing to the disordered structure of the material inducing a higher surface roughness compared to β-SiO₂. This also makes the determination of h value even more inaccurate. h is equal to around 55 Å if the maximum position is used. Taking into account of the layer thickness, h is equal to 51 Å.

As concerns He atoms, Fig.SI.F2C shows that the adsorbed He atoms remain at the β-SiO₂ wall surface. But the time average of the atom number inside the a-SiO₂ walls increases as T decreases. The decrease with T of gas number inside the wall does not contradict the increase of gas number crossing the walls. In fact, as we will see in the next section, the residence time inside the wall decreases with T and therefore the gas
atoms spend less time passing through the walls as T increases. We can also observe that the gas atom number at the surface wall is greater with β-SiO₂ than with a-SiO₂ at the same temperature because the gas atoms cannot enter the β-SiO₂ walls. The Ne profiles in Fig.SI.F2D show that the maximum has a much higher value compared to the He profiles for the same temperature because of the higher interaction magnitude with Ne. However, the number of Ne atoms inside the a-SiO₂ walls seems to be comparable to the He number.

The self-diffusion coefficient is computed with the mean square of the displacement (MSD noted also σ²):

$$D = \lim_{t \to \infty} \frac{\partial}{4\partial t} \langle (\vec{r}(t) - \vec{r}(0))^2 \rangle = \lim_{t \to \infty} \frac{\partial}{4\partial t} \text{MSD}(t) \quad (1)$$

Where $\vec{r}(t)$ and $\vec{r}(0)$ are the longitudinal positions at correlation time t and 0 respectively. This relation assumes that after a high number of wall collision, the diffusion is normal (MSD(t) ∝ t). This is observed after a correlation time of around 10 ns. For times lower than 5 ps, the gas transfer mechanism is mostly ballistic with MSD(t) ~ t² (see Fig. 1). This figure compares the MSD results with the Langevin model [61]. In this model a particle displacing in an effective Newtonian medium is submitted to a random force, so that the particle inertia is taken into account at short time. The equation of motion gives the following bidimensional form of the MSD expression:

$$\text{MSD}(t) = 4D(t - \tau_{\text{Lang}} + \tau e^{-t/\tau_{\text{Lang}}}) \quad (2)$$

Where $τ_{Lang}$ is the relaxation time of the Brownian motion. The model predicts that the intersection between the ballistic interpolation and the diffusion one is located at $2τ_{Lang}$. The typical value of $τ_{Lang}$ is a few dozen picoseconds for silica and hundred picoseconds for graphite. Whatever the silica (see also Fig.SI.F3A for a-SiO₂), the MSD does not follow accurately the Langevin model at the intermediate correlation time contrary to the graphitic channel. This will have an important

consequence for the analysis of $τ_{Lang}$ in section V.

Even though the gas atoms diffuse in the vacuum and in the solid phase, the self-diffusion coefficient follows the Arrhenius law (see Fig. 2):

$$D = D_A \times e^{-(T_A/T)} = D_A \times e^{-(E_A/k_B T)} \quad (3)$$

Where $k_B$ is the Boltzmann constant. However, given the small number of results with Ne/a-SiO₂, the exponential fit must be considered with caution for this system (E_A = 0.48 kcal/mol i.e. T_A = 220 K, D_A = 724 Å²/ps). The temperature dependance of D, i.e. E_A and D_A, is slightly sensitive to the SiO₂ wall structure. But roughness is responsible for the MSD difference between the graphite and the SiO₂ wall simulations. The Arrhenius law has been also observed with the system He/C in the temperature range 150 K–700 K [55] (see also Fig.SI.F3B). E_A does not depend on h because the system and the gas/wall energy potential are the same. But the ratio of D_A is equal to the wall distance ratio. Comparing the two systems He/SiO₂ and He/C, the value of the E_A ratio is close to the potential ratio $ε_{HeH}/ε_{HeC}$ ($= \sqrt{ε_{HH}/ε_{CC}}$). ε_{HeH} is chosen because H is the most outside atom of the walls. As concerns the two systems He/SiO₂ and Ne/SiO₂, the value of ratio $ε_{HeH}/ε_{NeH}$ is also close to the E_A ratio value. The same observation is made with He and Ne in graphite. However, the ratio agreement is less accurate for the two systems Ne/C and Ne/SiO₂ as the E_A ratio gives a value of 0.38 instead of 0.46 for the LJ energy. Finally, in the case h = 34 Å, D_A does not depend on the gas type. All these observations confirm the consistency of the MD simulations.

At this point of the analysis, Arrhenius law is observed in a temperature range [100, 700] K and [200,700] K as concerns He and Ne in graphitic channel respectively. We will see that Arya's model [34] agrees with the MD results at high temperatures and leads to the same order of magnitude at low temperatures (see Fig.SI.F13). The expression of the Arya's self-diffusion coefficient (D_{Arya}) does not depend explicitly

<figure><img src="image_3.png" /><figcaption>Fig. 1. Means square displacement of He vs. correlation time. Comparison between molecular dynamics simulations and Langevin model. β-SiO₂ wall channel at 300 K. and graphite wall channel at 350 K. Wall distance h = 50 Å.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Channel molecular self-diffusion of He and Ne in silica channel vs. T⁻¹. Comparison between molecular dynamics results (dots) and Arrhenius law D = Dₐe⁻¹/T (solid lines).</figcaption></figure>

on the Arrhenius law. The temperature dependence is due to the
contribution of the mean velocity via T¹/² and of TMAC. In the case of
graphitic channel, Arya's model predicts an agreement with the

Arrhenius law in a temperature range. At high temperatures, TMAC is
constant (see Fig.SI.F5) and D<sub>Arya</sub> increases as T¹/² if T > 750 K. At low
temperatures (T ≤ 100 K for He), TMAC gets closer to 1 and the

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** He/β-SiO₂ system. Channel MSD eigenvalues vs. correlation time. Black lines represent the linear regression of the channel MSD eigenvalues. σ²₁₂ = (σ₁² + σ₂²)/2.</figcaption></figure>
<table><caption>**Table 1**
Self-diffusion coefficient eigenvalue of He in β-SiO₂ channel.</caption><thead><tr><th></th><th>100 K</th><th>200 K</th><th>300 K</th><th>400 K</th></tr></thead><tbody><tr><th>D<sub>1</sub>(Å<sup>2</sup>/ps)</th><td>541</td><td>820</td><td>1180</td><td>1211</td></tr><tr><th>D<sub>2</sub>(Å<sup>2</sup>/ps)</th><td>344</td><td>612</td><td>857</td><td>938</td></tr></tbody></table>

temperature variation of its contribution decreases with T, so that the velocity contribution becomes dominant as T decreases.

The β-SiO₂ surface presents an anisotropy with a preferential orientation at 45° to the axes x and y (see Fig.SI.F1). Therefore, the MSD must be expressed in tensor form:

$$ \sigma^2(t) = \begin{pmatrix} \langle(x(t) - x(0))^2 \rangle & \langle(x(t) - x(0))(y(t) - y(0)) \rangle \\ \langle(y(t) - y(0))(x(t) - x(0)) \rangle & \langle(y(t) - y(0))^2 \rangle \end{pmatrix} \quad (4) $$

$$ = \begin{pmatrix} \sigma_{xx}^2(t) & \sigma_{xy}^2(t) \\ \sigma_{yx}^2(t) & \sigma_{yy}^2(t) \end{pmatrix} $$

The computation of the MSD tensor of He leads to a non-null value of the off-diagonal elements with $\sigma_{xx}^2 = \sigma_{yy}^2$ (see Fig.SI.F4A). The off-diagonal elements $\sigma_{xy}^2$ and $\sigma_{yx}^2$ must be equal. On average over time, $\sigma_{xy}^2$ is equal to 15 % $\sigma_{xx}^2$. In the eigen space, the MSD tensor $\sigma^2(t)$ is of the form:

$$ \sigma^2(t) = \begin{pmatrix} \sigma_1^2(t) & 0 \\ 0 & \sigma_2^2(t) \end{pmatrix} \quad (5) $$

The eigenvalues $\sigma_1^2(t)$ and $\sigma_2^2(t)$, computed with the MSD tensor at each correlation time, is a linear function in time if $T \ge 100$ K. The linear regression of $\sigma_1^2(t)$ and $\sigma_2^2(t)$ leads to the two eigen diffusion coefficients $D_1$ and $D_2$ as shown in Fig. 3. The analysis of the eigenvector shows an

orientation angle of 45° (see Fig.SI.F4B(a)) in accordance with the spatial distribution of Si, O and H shown in Fig.SI.F1A. The two diffusion coefficients $D_1$ and $D_2$ in the eigen space, which values are reported in Table 1, follow also the Arrhenius law (see Fig.SI.F4B(b)). But at low temperature, this orientation is no longer verified. We observe that $\sigma_{xx}^2 \neq \sigma_{yy}^2$ at long time (a deviation between 15 % and 20 % at 50 K) and the covariance value is less than 1 % the value of $\sigma^2$ at 50 K and equal to 7 % at 75 K. The He atoms tend to accumulate and stick at the surface (see Fig.SI.F2C(a)). In this case the MD diffusion method is not accurate enough to capture the surface anisotropy. However, it has no consequence on the self-diffusion averaged over the two longitudinal directions.

As concerns Ne, the orientation angle is not constant whatever T, even if $\sigma_{xy}^2 \approx 13\% \sigma_{xx}^2$ on average over time at high temperature, because $\sigma_{xx}^2 \neq \sigma_{yy}^2$. In the case of amorphous walls, no anisotropy is observed and $\sigma_{xx}^2 = \sigma_{yy}^2$ with $\sigma_{xy}^2 \ll \sigma_{xx}^2$. Therefore, the orientation angle of the eigenvectors is not significant.

## 4. Gas/wall collisional process

### 4.1. Tangential momentum accommodation coefficient

In accordance with numerous collision models which determine the accommodation coefficients with a gas/wall scattering kernel, the momentum accommodation coefficients are based on the incoming and outgoing momentum flux. Based on the Maxwell boundary condition which assumes that a fraction $\alpha_i$ of collision is diffusive and the other is reflexive, the traditional expression of the TMAC is as follows:

$$ \alpha_i = \frac{\langle V_i \rangle - \langle V_i \rangle}{\langle V_i \rangle - \langle V_T \rangle} \quad (6) $$

Where $V_i$ is the incoming flux of the velocity component i, $V_i$ is the outgoing one and $V_T$ is the outgoing flux at the wall temperature. However, this expression is not appropriate for complex surfaces such as

<figure><img src="image_6.png" /><figcaption>**Fig. 4.** TMAC vs. temperature. Dots: Collisional MD method. Lines: fitting function. Four systems are represented: β-SiO₂ wall with He and Ne, a-SiO₂ wall with He and Ne.</figcaption></figure>
anisotropic surfaces for which the velocity components may be cor- related with each other. Furthermore, in the condition of gas/solid ther- mal equilibrium this expression cannot give accurate results the denominator being close to zero. So, to avoid these drawbacks, the method of Spijker et al. [62], based on the least-squares formula for a collision cloud, is chosen and extended to cross correlation:

$$ \alpha_{ij} = 1 - \frac{\mathrm{cov}(V_i, V_j)}{\sigma_{V_i} \sigma_{V_j}} = 1 - \mathrm{cor}_{ij} \quad (7) $$

Where $\sigma_{V_i}$ and $\sigma_{V_j}$ are the standard deviation of $V_i$ and $V_j$ respectively and cor is the correlation coefficient. In Fig. 4, the mean tangential momentum accommodation coefficient (TMAC = $\alpha^t = 0.5(\alpha_{xx} + \alpha_{yy})$) varies from 0.65 to 0.95 as T decreases from 800 K to 75 K, the amplitude of variation depending mainly on the surface type. The high temperature value of TMAC for He/β-SiO$_2$ is in accordance with experimental values [63]. At all temperatures investigated, $\alpha^t$ is much higher than for the system He/C because the surface roughness and the gas/solid interaction are much higher (see Fig.SI.F5).

As regards β-SiO$_2$, $\alpha^t$ is smaller than in the case of a-SiO$_2$ mainly because no gas enters the β-SiO$_2$. It is surprising that $\alpha^t(\text{Ne}) < \alpha^t(\text{He})$ at high temperature. This is also observed with graphite (Fig.SI.F5). Substituting the Ne mass by the He one, simulations carried out at 800 K with Ne/C lead to the same value as $\alpha^t(\text{He})$. So, for Ne/He it is strictly an effect of mass. If we consider the elastic collision between two particles of mass $m_1$ and $m_2$ with $m_1 \gg m_2$, after collision $V_1 \approx V_1$. So, if we take aside the surface roughness and the gas/wall interaction, the longitudinal reflected velocity of Ne should be less modified by the collision than the He one, i.e. the tangential velocity correlation must increase ($\alpha^t$ must decrease) with gas mass. At the same time, $\alpha^t(\text{Ne})$ converges to a constant at higher temperatures than $\alpha^t(\text{He})$ owing to the higher long-range gas/wall interaction.

Authors [64] proposed an exponential expression of the temperature dependence of $\alpha^t$ in the form:

$$ \alpha^t = \alpha_{\infty}^t + \alpha_0^t e^{-T/T_a} \quad (8) $$

As T tends to 0 K, $\alpha^t$ must tend to 1. Therefore $\alpha_0^t \approx 1 - \alpha_{\infty}^t$. Fig. 4 shows that this expression agrees with the DM results. However, the agreement is not as accurate with graphite as concerns the relation between $\alpha_0^t$ and $\alpha_{\infty}^t$ and Eq. (8) underestimates the MD results with Ne at low temperatures (see Fig.SI.F5). The ratio $T_a(\text{Ne})/T_a(\text{He})$ is around 1.8 for the three materials. The higher sensitivity of $\alpha^t(\text{Ne})$ with T is explained mainly by higher gas/wall interaction and mass effect at high temperatures, so that the ratio depends little on the structure of the wall. Owing to the anisotropy of the β-SiO$_2$ surface, TMAC is in fact a tensor. The equality of the diagonal elements is observed with He and Ne as shown in Fig. 5. This condition must be fulfilled because the walls are statistically invariant to a rotation of 90° around the z axis. The high values of $\alpha_{xy}$ mean that the two velocity components are weakly correlated with each other, but the correlation cannot be neglected for the anisotropy analysis because the cross-correlation is clearly different from 0 over most of the temperature range investigated. The eigenvector orientation at 45° is observed because $\alpha_{xx} \approx \alpha_{yy}$. The highest eigenvalue $\alpha_2$ seems not to reach a constant value at high temperatures and is independent of the gas type contrary to $\alpha_1$. $\alpha_1(\text{Ne})$ increases much more when T > 300 K so that $\alpha_1(\text{Ne}) \sim \alpha_2(\text{Ne})$ at T = 75 K. The effect of the surface anisotropy on $\alpha^t$ and the self-diffusion coefficient is consistent as regards He, but the self-diffusion coefficient of Ne is not sensitive to the TMAC anisotropy.

## 4.2. Gas motion and mass transfer at the β-SiO$_2$ surface

At high temperatures, the gas kinetic energy being much higher than the well depth of the gas/wall potential energy, the gas atoms strike the wall once or two times before leaving the near wall region and $\alpha_{ij}$ depend mainly on surface roughness. At low temperature, the incoming gas atoms bounce several times on the wall, and the outlet velocity is less and less correlated to the inlet velocity ($\alpha_{ij}$ increases) as the bounce

<figure><img src="image_7.png" /><figcaption>Fig. 5. TMAC tensor elements ($\alpha_{xx}$, $\alpha_{yy}$) and the eigenvalues ($\alpha_1$, $\alpha_2$) vs. temperature. Systems: β-SiO$_2$ wall with He and Ne. Dots: collisional simulations.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Fig. 6. Bounce number of He on β-SiO₂ wall vs. residence time at four temperatures. Collisional simulations. Inset: Bounce number over 250 ps. Linear regression slope: 1.3 bounce/ps.</figcaption></figure>

number (BN) increases [19,35,55,65]. This is illustrated in Fig. 6 in the case of He/β-SiO₂. At T = 600 K, the residence time distribution in the interaction region (RTD) lies over 10 ps (see Fig. 7) and the bounce number is lower than 3. As T decreases, the RTD lies over an increasing time range and the BN increases so that at 50 K, the BN is a linear function of the residence time with a rate of 1.3 bounce/ps when the residence time is greater than 20 ps. This rate decreases to 0.4 bounce/ps as T increases to 300 K. The same phenomenon is observed whatever the gas and the walls used in this work as shown Fig.SI.F6 for instance. Obviously, the bounce rate of 0.77 bounce/ps for Ne at low temperature is smaller because of its higher mass. The RTD shows two-time distributions but at 600 K: an exponential decay at long time corresponding to the high number of bounces and its linear increase with time, and a Fréchet distribution at short time where the bounce number is close to one (see also Figs.SI.F.7A-C for comparison). The generalized Fréchet expression is as follows:

$$ Fr(\tau_{Fr}, \beta_{Fr}, t_{Fr}) = \frac{\beta_{Fr}}{\tau_{Fr}} \left( \frac{\tau_{Fr}}{t - t_{Fr}} \right)^{\beta_{Fr}+1} \exp \left[ - \left( \left( \frac{\tau_{Fr}}{t - t_{Fr}} \right)^{\beta_{Fr}} \right) \right] \quad \beta_{Fr} > 0 \quad (9) $$

The reason for the Fréchet distribution has been set out in Ref. [55]. The bounce number at short time is close to one. If T is sufficiently high to neglect the long-range gas/wall interaction and if the wall roughness is small enough, we can consider the gas/wall collision as specular ($V_z \approx -V_z$) and the residence time is around $2r_{cut}/V_z$. As $V_z$ follows the Rayleigh distribution with the scale parameter $\sqrt{2k_B T/m}$, the RTD at short time is close to the Fréchet distribution $Fr(\tau_{Fr}, 2, 0)$.

$\beta_{Fr}$ and $\tau_{Fr}$ are determined at short time by means of a linear regression of the double ln plot of the cumulative distribution function of the RTD vs ln(t) assuming $t_{Fr} = 0$. The time range of agreement increases with T.

The Fréchet of order $\beta_{Fr} = 2$ agrees with the RTD at $T > 100$ K and at $T < 100$ K the generalized Fréchet distribution must be used. For Ne/β-SiO₂ the threshold is 200 K (Fig.SI.7A). In the case of He/C (Fig.SI.F7B) and Ne/C, the RTD at short time follows the Fréchet of order 2 whatever the temperature investigated because the wall has no roughness. All the parameter values of the RTD are gathered in Tables SI.T6 and SI.T7. So, the Fréchet order does not depend on the wall type if $T > 200$ K. This is not the case of $\tau_{Fr}$ which is a decreasing function of T. To confirm the relevance of the collision approach, the position of the distribution maximum in function of $\sqrt{T}$ is reported in Fig. 8. It shows a linear variation for He and Ne whatever the material and a little dependency of the SiO₂ wall type. However, the slope is independent of the SiO₂ wall structure. Whatever the wall type, the He-Ne slope ratio is around 2.3 which is close to the mass ratio square root, i.e. 2.24. The same observation can be made with $\tau_{Fr}$. We will see that gas atoms do not enter the a-SiO₂ wall at short times, which may explain the low sensitivity to SiO₂ wall type.

From Tables SI.T6 and SI.T7, we can notice that the exponential decay is smaller with He and Ne/C than with He and Ne/β-SiO₂. It appears that the He/C momentum and energy exchange leading to desorption is less efficient with graphite walls and needs more bounces. Fig.SI.F8 shows that the characteristic time $\tau_E$ of the decay follows an exponential law of temperature in the form:

$$ \tau_E = \tau'_E e^{T_B/T} \quad (10) $$

In contrast to Fig. 8, which displays a difference between graphite and SiO₂ wall, Fig.SI.F8 shows two groups depending on the gas type. The straight lines are quite parallel in each group. The MD results give the hierarchic inequalities $T_E(C) > T_E(a-SiO₂) > T_E(\beta-SiO₂)$ and $\tau'_E(C) > \tau'_E(a-SiO₂) > \tau'_E(\beta-SiO₂)$ for He and Ne. The inequalities C/SiO₂ are
<figure><img src="image_9.png" /><figcaption>**Fig. 7.** Residence time distribution at four temperatures. System: He/β-SiO₂. Dots: Collision simulations. Lines: fitting functions; (a): Fréchet distribution, (b): exponential distribution.</figcaption></figure>

related to the less efficiency of bounces on graphite walls for the
desorption at long-range residence time. The inequalities a-SiO₂/β-SiO₂
are related to the gas transfer in a-SiO₂ wall.

Fig. 9 (Fig.SI.F9A for Ne/β-SiO2) shows that at T > 100 K the main transfer mechanism is ballistic. Whatever T, the ballistic regime lies over a time range lower than 1 ps. At low temperatures, i.e. 50 K/75 K for He and 200 K/100 K for Ne, the diffusion processus appears at t > 10 ps. The correlation time t corresponds to a time interval inside a residence time t<sub>R</sub> = t<sub>exit</sub> - t<sub>entry</sub> > t. So, the MSD(t) is the result of the contribution of all the near wall trajectories during this time interval. This means that the MSD is more sensitive to the long residence time, i.e. when the bounce number increases linearly with t<sub>R</sub>, because the number of intervals increases with t<sub>R</sub>. Therefore, as concerns the diffusion regime, the time threshold of 10 ps means a bounce number greater than 13 (10 × 1.3). For the ballistic regime the bounce number is lower than 1.3.

These BN values are much smaller than in the case of He/C at 50 K
[55]. Indeed, considering the time threshold and the bounce rate, the
diffusion and the ballistic regime are established if the bounce number is
greater than 95 and lower than 20 respectively. With Ne/C at 150 K, the
BN must be greater than 45 and lower than 9 respectively. These find-
ings confirm the previous observation about the low bounce efficiency
with graphite wall.

At low temperatures, the surface anisotropy of β-SiO₂ is also captured
by means of the collision-MSD tensor (see Figs.SI.F.9B,C). As σ²ₓₓ ≈ σ²ᵧᵧ,
the orientation angle of the eigenvectors is 45°. The eigen diffusion
coefficients shown in Table 2 are far much smaller than in the case He/C
[55] where the surface diffusion coefficient is around 300 Å²/ps at 50 K
for instance. As a last comment, the MD simulations reveal that the
Langevin model agrees with the collision-MSD in the case of the graphite
wall but not in the case of SiO₂ walls, which is consistent with the
<figure><img src="image_10.png" /><figcaption>Fig. 8. Time of the RTD maximum vs. T⁻¹/². Systems: He and Ne/β-SiO₂/a-SiO₂ and/C. Dots: collision simulations. Lines: linear regression.</figcaption></figure>

<table><caption>**Table 2**
Collision-diffusion coefficient eigenvalue of He and Ne on β-SiO₂ wall.</caption><thead><tr><th></th><th>50 K (He)</th><th>75 K (He)</th><th>100 K (Ne)</th></tr></thead><tbody><tr><td>D₁(Å²/ps)</td><td>5.07</td><td>15.9</td><td>5.5</td></tr><tr><td>D₂(Å²/ps)</td><td>2.58</td><td>16.1</td><td>3.2</td></tr></tbody></table>

channel simulations (not shown).

### 4.3. Gas motion and mass transfer at the a-SiO₂ surface

For the a-SiO₂, the bounce must be defined in a different way because the gas atoms can enter the wall. A bounce is characterized by a sign change of the transversal velocity $V_z$ but inside the wall, the gas velocity fluctuates owing to the contact with the silicon and oxygen atoms without leaving the wall. Two constraints can be used to avoid inner bounces. The first one consists in imposing a minimum transversal position ($z_{min}^B$) below which a sign change of the transversal velocity $V_z$ cannot be considered as a bounce process. The second one consists in imposing a condition on the maximum transversal position ($z_{max}^B$) of a trajectory. A sign change is considered as a bounce if the maximum of transversal position reached by an atom immediately after the sign change is greater than ($z_{max}^B$). We have already observed the presence of a superficial hydroxyl layer located at 14 Å and 3.5 Å thick. So, a value of 14 Å and 16 Å for $z_{min}^B$ and $z_{max}^B$ is chosen respectively. For the first constraint (C1), the transversal bounce position is inside this layer. For the second constraint (C2), the maximum trajectory after bouncing is at least at the upper side of this layer.

The effects of these constraints on the mean transversal position of the bounces are displayed in Fig. 10(a) for He at 50 K. The average is carried out over all the bounces which happened during a collision

process and averaged over all the collisions having the same residence time. The minimum transversal position ($z_{min}^C$) reached by a gas atom at each collision process and averaged over all the collisions having the same residence time is also reported in this figure (black curve). The $z_{min}^C$ curve shows that as the gas atoms explore the wall more deeply the residence time increases and shows that, on average, He atoms cross the superficial hydroxyl layer at long residence time. On contrary, with the constraints and in average, only the bounces at the wall surface are considered. At short time, the constraints have no influence because of very quick bounces contribution only.

However, the bounce number is much more sensitive to the constraint type. At long residence time, the constraints break the linearity, and Fig. 10(b) shows that most bounces (sign changes of $V_z$ without constraint) occur without the gas atoms leaving the wall. With C2, the bounce number converges to around 60. So, over a time threshold, 50 ps in the case of He, the bounce number being quite constant during one collision with the wall, the gas atoms remain more times inside the wall. In other words, the bounces are more and more spaced over time because their average transversal bounce position is around 15.8 Å and the gas atoms explore more and more deeply the amorphous wall between two bounces. The threshold of 50 ps corresponds to the threshold value of the exponential decay of the RTD (see Fig.SI.F7C), which reveals the relevance of the constraint C2. All these results are also observed with Ne (see Fig.SI.F.10).

The collision-MSD with β-SiO₂ (Fig. 9) and a-SiO₂ (Fig.SI.F11) is qualitatively identical However, the surface diffusion is 0.6 Å²/ps with a-SiO₂ instead of 3.8 Å²/ps with β-SiO₂. So, it suggests that the sub-diffusional mechanism inside wall has an influence on the collision-MSD. As concerns Ne, the diffusion regime at 100 K is not clearly reached at long time.
<figure><img src="image_11.png" /><figcaption>Fig. 9. Collision-MSD vs. correlation time. System: He/β-SiO₂ at five temperatures.</figcaption></figure>

4.4. *Mass transfer and residence time distribution inside the a-SiO₂ wall*

To investigate the gas transfer inside the channel wall during the Knudsen diffusion; the MSD inside the wall, noted wall-MSD, is also computed with the collision methodology. To determine the wall-MSD, only trajectories which cross a first control plane located at z = 13.75 Å are considered and the parts of these gas trajectories under a second control plane located at z = 13.5 Å are used for the wall-MSD calculation. These plane positions are located under the OH layer. The first plane is to ensure that gas atoms leave the bulk wall region when crossing the second wall. We recall that the collision-MSD is calculated all along a trajectory contrary to the wall-MSD.

Fig. 11 shows the wall-MSD of He from 50 K to 600 K (see Fig.SI.F12B for Ne). These curves are compared to results obtained by self-diffusion simulations of 1000 gas atoms inside the initial amorphous silica block (50 Å × 50 Å × 50 Å) in thermal equilibrium (internal MD method). These simulations are carried out for more than 10 ns and the self-diffusion MSD (noted int-MSD) is averaged over the three directions. As concerns the collision simulations, the wall-MSD is averaged over the two longitudinal directions.

The internal self-diffusion methodology allows us to reach long correlation time and to capture the three transfer regimes: the ballistic process at short time with MSD $\propto t^2$, the sub-diffusion process with MSD $\propto t^\alpha$ and $\alpha < 1$ at intermediate time range and the normal diffusion process (MSD $\propto t$) at long time. The sub-diffusion process results from a cage effect, i.e. the first neighbors (Si,O) of the gas atoms prevent them from spreading [66–69]. In the ballistic range, MD simulations with He and Ne agree with the linear dependance of MSD on the thermal velocity square $k_B T/m$. The internal self-diffusion simulations reveal a somewhat diffusion transfer ($\alpha = 0.8$) at 600 K for He when $t > 5$ ns. Of course, the transition happens at higher correlation time for Ne owing to its higher mass and gas/wall interaction. In the sub-diffusion regime, $\alpha$ varies from 0.14 to 0.41 as T increases. The ballistic/sub-diffusion transition happens at the same time range for the two gas types. Owing to the smaller time range, the collision methodology cannot capture the diffusion

regime. It agrees with the internal methodology in the ballistic regime but in the sub-diffusion regime, wall-MSD is almost constant in time ($\alpha \sim 0$) and the transition to the sub-diffusion regime is more pronounced in the same time range. So, the collision methodology reveals a higher contribution of the ballistic transfer and the gas stops spreading as the time correlation is greater than 1 ps.

The RTD inside the wall is also investigated with the same collision method. In particular, the trajectory parts located below $z = 13.5$ Å are only considered. Figs.SI.F12A,C show typical He and Ne wall-RTDs with the same shape as the collision RTD while one would expect an exponential distribution. The residence time at short time is related to the ballistic gas transfer. But the decay follows a power law derived from the sub-diffusion regime. The wall-RTD leads to characteristic times 20 times smaller than the collision-RTD ones. This does not mean that the time spent inside the wall is negligible during the collision process. As observed above, the rate number of bounce outside the wall, computed with the constraint C2 decreases with the near wall residence time as this residence time is greater than 50 ps. The time spent in the wall ($z < 13.5$ Å) between two bounces increases with the wall residence time inducing an increase of its contribution to the collision-RTD at long time.

Considering the collision-MSD time range lower than 200 ps and the wall-MSD time range lower than 50 ps in accordance with the wall-RDT, we can conclude that the diffusion regime inside the a-SiO₂ cannot be reached with the collision and the channel diffusion methodology.

5. 2D analysis of the channel diffusion

In this section the Langevin relaxation time is analyzed (see §3.). Knudsen diffusion coefficient is computed from the MSD at long time where MSD $\propto t$. At short time normal diffusion process does not occur because the collision number with the walls is too small to reach the asymptotic stochastic behavior. At the early stage of gas transfer, the ballistic process is dominant as MSD $\propto t^2$. The MD simulations with graphite channel show that the transition from the ballistic to the diffusional regime can be modeled by means of the Langevin model. In
<figure><img src="image_12.png" /><figcaption>Fig. 10. Transversal bounce position (a) and bounce number (b) of He vs. residence time. Wall: a-SiO₂. T = 50 K. Comparison with several constraints. (a): black curve (Z<sub>min</sub><sup>C</sup>(t)): minimum transversal position reached during a collision process. (b): Inset bounce number over 250 ps.</figcaption></figure>

the case of SiO₂ channel, the agreement is not fulfilled even if the transition is clearly established (see Fig. 1 and SI.F3A). However, we will see that the Langevin model can give some insights with silica wall.

This model links the diffusion coefficient, the relaxation time and the velocity correlation with the particle mass, the effective medium viscosity (η) and the temperature. In bidimensional systems, the Langevin model gives the following relations:

$$D = \frac{k_B T}{\eta} \quad \text{and} \quad \tau_{\text{Lang}} = \frac{m}{\eta} \qquad (11)$$

$$\langle V(t)V(0) \rangle = \langle V(0)^2 \rangle e^{-\frac{t}{\tau_{\text{Lang}}} \qquad (12)}$$

$$\text{With } \langle V(0)^2 \rangle = 2D/\tau_{\text{Lang}} = 2k_B T/m = V_T^2 \qquad (13)$$
<figure><img src="image_13.png" /><figcaption>Fig. 11. Mean square displacement of He in a-SiO₂ vs. correlation time at four temperatures. Comparison between the collision method (wall-MSD) and the internal method (int-MSD).</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Fig. 12. Velocity vs. temperature. Channel molecular simulations. Comparison between the Langevin Model (lines) and the kinetic theory of gases (Dots).</figcaption></figure>

Approximately, we can reason that the friction and random force comes from the gas/wall interaction, and it is interesting to verify if the same behavior can be found with the MD simulations and atomistic wall. In Fig. 1 the ballistic and the diffusion processes are represented with a

straight line of slope 2 and 1 respectively which intersects at $t_i = 2\tau_{\text{Lang}} = 335$ ps.

Using the MD results for He/C and Ne/C, the variation of $V_{\text{Lang}} = \sqrt{2D_{\text{MD}}/\tau_{\text{Lang}}}$ with temperature is displayed in Fig. 12 for two wall
distances. It is clear that $V_{Lang}$ is equivalent to $V_T$ in this case. This observation and the above comments prove the coherence between the MD, the Langevin model, the gas kinetic theory, and the isotropic gas/wall collision process at all time scales even at very low temperatures. At this stage, a parallel can be drawn with the rarefied gas diffusion. In an infinite medium without walls, the self-diffusion coefficient $D_G$ is related to the mean free path $l_G$ and the 3D mean velocity:

$$D_G = \frac{l_G}{3} \sqrt{\frac{8k_B T}{\pi m}} \qquad (14)$$

In the present study, the 3D velocity is replaced by the 2D root mean-square velocity (or thermal velocity $V_T$) and the mean free path by a relaxation length ($l_{Lang} = V_{Lang}\tau_{Lang}$):

$$D_{Lang} = \frac{l_{Lang}}{2} \sqrt{\frac{2k_B T}{m}} \qquad (15)$$

The factor 2 in the denominator is related to the dimensionality of the process.

In the case of $SiO_2$ walls, $\tau_{Lang}$ is computed in the same way. Despite the overestimation of the Langevin MSD, $V_{Lang}$ and $V_T$ have the same order of magnitude for He and Ne for the two $SiO_2$ walls. The anisotropy of the $\beta$-SiO$_2$ wall has no effect.

Since $V_{Lang}$ is close to the longitudinal thermal velocity, an analytic expression of relaxation time can be derived from the Arya model [34]:

$$\tau_{Arya} = 2 \frac{D_{Arya}}{V_T^2} = -h \sqrt{\frac{\pi m}{2k_B T} \left[ \frac{\ln \epsilon_C}{2\pi} - \frac{1-\alpha_t}{\alpha_t} \right]} \qquad (16)$$

Where $\epsilon_C = mV_C^2/(2k_B T)$. The authors used a stochastic model coupled with the Cercignani-Lampis collision kernel [29], in the case of a slit pore delimited by two rigid planes. The gas/wall interaction is modeled by the hard sphere potential, and the gas atoms do not interact. The analytic self-diffusion expression $D_{Arya}$ requires a velocity cutoff owing to the planar wall hypothesis. In the following, a cutoff time $\tau_C = h/V_C$ is

used instead of $V_C$.

Let us consider the case of the graphitic channel. As concerns the diffusion coefficient, Fig.SI.F13 compares the MD results with the Arya model. The model is computed with the TMAC values obtained by MD. At $T > 800$ K the TMAC value at $T = 700$ K is used because it remains constant at high temperatures (Fig.SI.F5). The figure shows the agreement between the model and the MD simulations from 800 K to 100 K with $\tau_C = 10$ ns for He and down to 250 K with $\tau_C = 25$ ns for Ne. We can notice that the value of $\tau_C(Ne)/\tau_C(He)$ is close to that of $\sqrt{m_{Ne}/m_{Ne}}$. Below these temperatures, the Arya model is no longer predictive because the gas/wall long range interaction is not considered.

Fig. 13 compares $\tau_{Lang}$ calculated by MD simulations in graphite channel with the correlation time $\tau_{Arya}$. $\tau_{Arya}$ is computed in two ways: with Eq. (8) ($\alpha^t(\alpha_\infty, \alpha_0, T_\alpha)$) and with the MD values of TMAC ($\alpha^t(MD)$). The overall agreement between the Arya model and the MD results is observed throughout the temperature range for the two gases. This is surprising at small temperature because the adsorption is not actually considered by the Arya model. However, as concerns Ne the Arya model with the two TMAC disagrees with the MD results at around 150 K because the model underestimates the channel MSD at this temperature (see Fig.SI.F13).

Two regions are shown with a transition temperature of 270 K for He and 500 K for Ne, located inside the Arrhenius law region. The variation of $\tau_{Arya}$ results of the increase of $(1-\alpha^t)/\alpha^t$ and the decrease of $T^{-1/2}$ as T increases. $T^{-1/2}$ represents the collision frequency ($\approx V_T/\hbar$) and the TMAC ratio is related to the number of collisions ($N_C$) for the tangential velocity component to be uncorrelated.

If we consider gas transfer as a series of free flights interrupted by wall collisions, $N_C$ is the ratio of the correlation time $\tau_{Arya}$ to the mean time between two collisions ($\tau_F$). As the transversal velocity is described by the Rayleigh distribution, the expression $\langle \tau_F \rangle$ is computed as follows:

<figure><img src="image_15.png" /><figcaption>Fig. 13. Channel MSD relaxation time vs. temperature. (a): He and Ne with graphite wall. Comparison between molecular dynamics results and Arya model. O: TMAC determined by collision methodology. -: TMAC computed with Eq. (8).</figcaption></figure>
$$ \langle \tau_F \rangle = h \sqrt{\frac{m}{k_B T}} \times \frac{\int_0^\infty e^{-x}/\sqrt{2x}dx}{\int_0^\infty e^{-x}dx} = h \sqrt{\frac{\pi m}{2k_B T}} \quad (17) $$

With $x = mV_z^2/k_B T$. Therefore:

$$ N_C = - \left[ \frac{\ln \epsilon_C}{2\pi} - \frac{1 - \alpha^t}{\alpha^t} \right] \quad (18) $$

$N_C$ can be also interpreted as a kind of Peclet number because it is related to the ratio of the self-diffusion to the thermal velocity and the channel width:

$$ Pe = \frac{D}{hV_T} = \frac{\sqrt{\pi}}{2} N_C \quad (19) $$

The ratio $D/V_T \sim N_C h$ represents the diffusion length or correlation length. Using the estimated values of $\epsilon_C$ and $\alpha^t$ (MD), $N_C$ increases from 4 to 35 and from 3 to 40 as T increases from 75 K to 700 K for He and Ne respectively.

At high temperature, the wall collision frequency per gas atom increases with T. Therefore, the time needed for a particle to have a random velocity orientation decreases with T because $\alpha^t$ is constant. If the surface is perfectly smooth, the collisions are only specular, and the tangential velocity of a gas atom always remains unchanged so that the mass transfer regime is ballistic for an infinite time whatever the collision frequency. The value of $\tau_{Arya}$ is related to the surface roughness.

At small temperature, the correlation time $\tau_{lang}$ increases with T because even the adsorption probability and the adsorption time decreases, i.e. $\alpha^t$ decreases, the number of gas/wall collision for the reflected velocity to reach a random orientation must increase. In addition, the term $(1 - \alpha^t)/\alpha^t$ is dominant and the increase of the collision frequency is not enough to balance the decrease of the adsorption contribution. At the transition temperature, where the collision frequency is balanced by the collision number $N_c$, a regime transition takes place from an adsorbed atom to a hard sphere. The transition temperature range of $\tau_{Arya}$ covers the Arrhenius one of the Knudsen diffusion.

As concerns He and SiO₂ walls, the Arya model agrees with the MD results if the cutoff time is 475 fs in the case of β-SiO₂ at T ≥ 100 K, but no agreement is found with a-SiO₂ owing to gas absorption. For the β-SiO₂ wall, the cutoff parameter $\epsilon_C$ is no longer negligible at all temperatures. So, the value of the cutoff time is questionable. The same observations are made for Ne/β-SiO₂ with a cutoff value of 500 fs. In Fig. SI.F14, the maximum of $\tau_{lang}$ seems to be located at 75 K with a value around 20 ps for β-SiO₂. The Arya model predicts the decrease of the correlation time with T at moderate and high temperatures but not the time threshold at low temperatures. In other words, the model prediction is a cutoff effect. However, the decrease with temperature is induced by the $T^{-1/2}$ contribution of the frequency term without any effects of adsorption time and incident/reflected velocity correlation. Results with a cutoff value of 360 fs are reported for the amorphous wall which leads to a Knudsen diffusion not too far from the MD results. Using the MD values of $\tau_{lang}$ and the expression of $\langle \tau_F \rangle$, the collision number increases from 1 to 2 in the temperature range investigated. This means that the diffuse collision as wall boundary condition is consistent from the perspective of the transition regime. As concerns Ne, the MD

results are not accurate enough to show a clear dependance of the relaxation time with temperature for the two SiO₂ walls.

## 6. Stochastic approach

In the following, a stochastic model of self-diffusion is proposed for the systems (He and Ne)/β-SiO₂. Gas particles at temperature T move in a 2D channel with a longitudinal and a transversal velocity following the Maxwell and the Rayleigh distribution respectively. The channel, delimited by two flat planes, has a width noted $h$. The particle displacement is computed at each time step $dt$. In order to consider the near wall RTD computed by collisional simulations, the particles move freely over a transversal distance of $h - 2 \times r_{cut}$. As soon as a particle reaches the distance $r_{cut}$ from a wall, a residence time ($t_{RTD}$) is chosen randomly in accordance with the RTD. Collisions are diffusive for the transversal velocity component. Therefore, when the particle leaves the near wall region, a value of transversal velocity is chosen randomly from the Rayleigh distribution. Owing to the surface diffusion at low temperature, two kinds of treatment are carried out for the longitudinal displacement in the near wall region.

At high temperature ($T \ge 100$ K for He and $T > 100$ K for Ne), the wall MSD does not reveal a well-established diffusion process. So, the ballistic process is taken into account only, and we assume that the RTD follows a Fréchet distribution of order 2 ($Fr(\tau_{Fr}, 2, 0)$ in Eq. (9)) as explained in §4.2. In the near wall region, we consider that the particles move in the longitudinal direction with the mean velocity $\frac{1}{2}(V_x + V_x)$ where $V_x$ is the incoming velocity and $V_x$ the outgoing velocity following the Cercignani-Lampis boundary condition:

$$ V_x = (1 - \alpha^t) \times V_x + \sqrt{\alpha^t(2 - \alpha^t)} \times V_x^M \quad (20) $$

Where $V_x^M$ is the Maxwell distribution. Therefore, the ballistic displacement $dx_{bal}$ during $t_{RTD}$ is $\frac{1}{2}(V_x + V_x) \times t_{RTD}$. So, during the time step $dt$

$$ dx_{bal} = \frac{1}{2} (V_x + V_x) \times dt \quad (21) $$

At low temperature ($T < 100$ K for He and $T = 100$ K for Ne), a fraction of the gas/wall collision is diffusive. Diffusive collision happens at long time, i.e when the RTD follows the exponential distribution characterized by the decay time $\tau_E$. RTD is modeled by the Fréchet $Fr(\tau_{Fr}, 2, 0)$ and the exponential distributions at $t_{RTD} < \tau_{sep}$ and $t_{RTD} > \tau_{sep}$ respectively. Considering the wall-MSD, $\tau_{sep}$ is the transition time from the ballistic process to the diffusion one.

So, two cases must be considered:

1. If $t_{RTD} < \tau_{sep}$, the ballistic displacement $dx_{bal}$ during $dt$ is considered only and is computed with Eqs.20 and 21.

2. If $t_{RTD} > \tau_{sep}$, the gas/wall collisions are mainly diffusive. The particle displacement is the sum of the ballistic ($dx_{bal}$) and of the diffusive ($dx_{diff}$) displacement. Owing to the diffusive gas/wall collision, we impose $V_x^M = V_x$ for the ballistic process which occurs during $\tau_{sep}$. The diffusion process occurs during $t_{RTD} - \tau_{sep}$. It induces a displacement

<caption>**Table 3**</caption>

Stochastic computation. Self-diffusion coefficient and relaxation time of He and Ne in β-SiO₂ channel. In parenthesis: results with the high temperature method. $\tau_{Lang}$ is the MD relaxation time.

<table><thead><tr><th></th><th>50 K</th><th>75 K</th><th>100 K</th><th>200 K</th><th>300 K</th><th>400 K</th><th>600 K</th></tr></thead><tbody><tr><td><strong>He/β-SiO<sub>2</sub></strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>D<sub>Sto</sub>(Å<sup>2</sup>/ps)</td><td>177 (232)</td><td>270 (390)</td><td>462</td><td>660</td><td>837</td><td>942</td><td>2515</td></tr><tr><td>τ<sub>Sto</sub> (ps)</td><td>33 (29)</td><td>26 (33)</td><td>29</td><td>20</td><td>15</td><td>14</td><td>11</td></tr><tr><td>τ<sub>Lang</sub> (ps)</td><td>18</td><td>21</td><td>20</td><td>16</td><td>15</td><td>13</td><td></td></tr><tr><td><strong>Ne/β-SiO<sub>2</sub></strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>D<sub>Sto</sub>(Å<sup>2</sup>/ps)</td><td></td><td></td><td>108 (174)</td><td>301</td><td>385</td><td></td><td>690</td></tr><tr><td>τ<sub>Sto</sub> (ps)</td><td></td><td></td><td>59 (57)</td><td>46.4</td><td>39</td><td></td><td>29.5</td></tr><tr><td>τ<sub>Lang</sub> (ps)</td><td></td><td></td><td>31</td><td>29</td><td>31</td><td></td><td></td></tr></tbody></table>
<figure><img src="image_16.png" /><figcaption>Fig. 14. Self-diffusion of He and Ne in β-SiO₂ channel vs. T⁻¹. Comparison between the stochastic model and the MD computation. Lines: Arrhenius law.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. 15. Relaxation time of He and Ne in β-SiO₂ channel vs. T. Comparison between the stochastic model (τSto) and the channel MD computation (τLang). Lines: Fit with the function T⁻⁰.⁵.</figcaption></figure>

N(0, √2D_S(t_RTD - τsep)) where N is the gaussian distribution of zero mean value and D_S is the surface diffusion coefficient. The longitudinal displacement dxsurf during dt has the following expression:

$$dx_{surf} = \frac{1}{2}(V_x + V_x) \times dt + V_{diff} \times dt \quad (22)$$

Where Vdiff = -N(0, √2D_S(t_RTD - τsep)). In the expression of Vdiff, the diffusion displacement is assumed to be linear in time.

The distribution parameters, TMAC and D_S are reported in Table SI. T.8. τFr is computed by fitting the RTD with the Fréchet distribution by means of the maximum likelihood estimation and τsep is determined
from the wall MSD curve (Fig. 9).

The stochastic self-diffusion ($D_{Sto}$) and the stochastic relaxation time ($\tau_{Sto}$) are gathered in Table 3. Results obtained with the high temperature method at low temperatures are inserted in parenthesis. Fig.SI.F15, the stochastic MSD shows an underestimation in the intermediate regime and confirms the MD results (see Fig. 1). The underestimation is observed over all the temperature range investigated. The mean time between two wall collisions (Eq. (17)) is 8.3 ps and 4.8 ps at 75 K and 300 K respectively. These times are close to the stochastic relaxation time (Table .3).

Fig. 14, $D_{Sto}$ (stochastic model) is compared to the MD self-diffusion coefficient. At high temperatures, the stochastic method tends to underestimate the self-diffusion coefficient, but they agree the MD results in the case of He at low temperatures. If the high temperature method is used at low temperatures, a clear overestimation is observed. For He, the stochastic and the MD activation energies are rather close. So, the stochastic approach allows us to obtain the same order of magnitude as the MD computation.

The relaxation time computed by the stochastic method $\tau_{Sto}$ is a decreasing function of temperature (see Fig. 15). The results with the RTD at T = 600 K are inserted in the Figure. It decreases roughly as $T^{-1/2}$ and confirms that $\tau_{Sto}$ is linked to the mean time between two wall collisions. However, owing to the small number of results, this confirmation must be considered with caution as concerns Ne. Comparing the MD results ($\tau_{Lang}$) with He at low temperatures, the stochastic computation overestimates the relaxation time at $T \le 200$ K. For Ne, $\tau_{Lang}$ seems to be constant in the small temperature range owing to its fluctuation coming from the MD uncertainties. It was not possible to capture the wall collision frequency dependence of $\tau_{Lang}$.

In the case of a-SiO$_2$, the same trends are observed. However, $D_{Sto}$ is not estimated with sufficient accuracy at low temperatures.

## 7. Conclusion

We have presented the investigation of the Knudsen diffusion of He and Ne in nano-slit channels through the gas/wall collision process. Influence of surface roughness and anisotropy as well as wall porosity is analyzed with several materials such as graphite, β-cristobalite and amorphous SiO$_2$. The gas/wall interaction process depending on temperature, the investigation has been carried out from 50 K to 800 K.

The molecular dynamics results show that the Knudsen diffusion follows the Arrhenius law in a temperature range depending on the wall material. The anisotropy of β-cristobalite has been captured with He. The Knudsen diffusion eigenvalues also obey the Arrhenius law, and the two eigenvector directions are parallel and perpendicular to the hydrogen surface distribution. This was also observed for the TMAC. Determining the characteristic transition time according to Langevin model, a 2D approach of the mass transfer is proposed by replacing the effective viscosity by a kind of gas/wall surface friction. Combined with the Arya model, the Langevin model reveals the temperature dependence of the transition time on the gas/wall collision frequency and the collision number needed for the tangential velocity to be uncorrelated. Considering the RTD, a stochastic model confirms the channel self-diffusion analysis in the case of β-SiO$_2$ wall. In particular, the Arrhenius law is verified and the temperature variation of velocity correlation time.

At a smaller scale in space and time, the collision mechanism has been investigated with the analysis of the residence time distribution in the vicinity of the wall surface, the bouncing process and the gas transfer inside the channel wall. The RTD is composed of two distributions: the Fréchet distribution at short time and an exponential decay at long time. The Fréchet distribution is the results of the ballistic transfer and collision process composed of one to three surface bounces. The exponential decay is related to a high bounce number and to the surface diffusion at low temperatures. In the case of a-SiO$_2$, a new accurate bounce criterium

shows that the gas atoms spend a time between two bounces inside the amorphous wall which increases with the residence time. The study of gas transfer within the wall showed short-term ballistic transfer and a long-term cage effect. These regimes are consistent with the inside RTD shape. The diffusion regime at very long time (higher than 5 ns) is expected to be observable at high temperature only.

The relevance of this work lies in considering the Knudsen diffusion as an interplay of processes at multiple length and time scales, taking into account both RTD in the region near the wall and inside it, as well as the mass transfer mechanism at the surface and inside the wall. In addition, with a two-step stochastic model, we were able to predict the Knudsen diffusion and its ballistic/diffusion transition time.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgments

Centre de Calcul Intensif d'Aix-Marseille is acknowledged for granting access to its high-performance computing resources.

## Appendix A. Supplementary data

Supplementary data to this article can be found online at https://doi.org/10.1016/j.vacuum.2025.114756.

## Data availability

Data will be made available on request.

## References

[1] E. Rivard, M. Trudeau, K. Zaghib, Hydrogen storage for mobility: a review, Materials 12 (12) (2019) 1973, https://doi.org/10.3390/ma12121973.

[2] H. Pan, C. Yu, X. Suo, L. Yang, X. Cui, H. Xing, Emerging porous materials for carbon dioxide adsorptive capture: progress and challenges, Mater. Chem. Front. 7 (2023) 6463, https://doi.org/10.1039/d3qm00705g.

[3] A. Staerz, U. Weimar, N. Barsan, Current state of knowledge on the metal oxide based gas sensing mechanism, Sensor. Actuator. B Chem. 358 (2022) 131531, https://doi.org/10.1016/j.snb.2022.131531.

[4] R. Alrammouz, J. Podlecki, P. Abboud, B. Sorli, R. Habchi, A review on flexible gas sensors: from materials to devices, Sensor. Actuator. B 284 (2018) 209–231, https://doi.org/10.1016/j.sna.2018.10.036.

[5] J.H. Choi, J. Lee, M. Byeon, T.E. Hong, H. Park, C.Y. Lee, Graphene based gas sensors with high sensitivity and minimal sensor-to-sensor variation, ACS Appl. Nano Mater. 3 (3) (2020) 2257–2265, https://doi.org/10.1021/acsanm.9b02378.

[6] Z. Cai, B. Liu, X. Zou, H.M. Cheng, Chemical vapor deposition growth and applications of two-dimensional materials and their heterostructures, Chem. Rev. 118 (13) (2018) 6091–6133, https://doi.org/10.1021/acs.chemrev.7b00536.

[7] S. Bhowmik, A.G. Rajan, Chemical vapor deposition of 2D materials: a review of modeling, simulation, and machine learning studies, iScience 25 (3) (2022) 103832, https://doi.org/10.1016/j.isci.2022.103832.

[8] B. Biesemans, J. De Clercq, C.V. Stevens, J.W. Thybaut, J. Lauwaert, Recent advances in amine catalyzed aldol condensations, Catal. Rev. 66 (2) (2022) 401–483, https://doi.org/10.1080/01614940.2022.2048570.

[9] A. Garcia, I.I. Slowing, J.W. Evans, Pore diameter dependence of catalytic activity: nitrobenzaldehyde conversion to an aldol product in amine-functionalized mesoporous silica, J. Chem. Phys. 149 (2018) 024101, https://doi.org/10.1063/1.5037618.

[10] I.C. Gerber, P. Serp, A theory/experience description of support effects in carbon-supported catalysts, Chem. Rev. 120 (2) (2019) 1250–1349, https://doi.org/10.1021/acs.chemrev.9b00209.

[11] N. Talukder, Y. Wang, B.B. Nunna, E.S. Lee, Nitrogen-doped graphene nanomaterials for electrochemical catalysis/reactions: a review on chemical structures and stability, Carbon 185 (2021) 198–214, https://doi.org/10.1016/j.carbon.2021.09.025.

[12] B.R. Goldsmith, B. Peters, J.K. Johnson, B.C. Gates, S.L. Scott, Beyond ordered materials: understanding catalytic sites on amorphous solids, ACS Catal. 7 (2017) 7543–7557, https://doi.org/10.1021/acscatb.7b01767.

[13] S. Wang, F. Javadpour, Q. Feng, Molecular dynamics simulations of oil transport through inorganic nanopores in shale, Fuel 171 (2016) 74–86, https://doi.org/10.1016/j.fuel.2015.12.071.
[14] P. Asai, J. Jin, M. Deo, J.D. Miller, D. Butt, Non-equilibrium molecular dynamics simulation to evaluate the effect of confinement on fluid flow in silica nanopores, Fuel 317 (2022) 123373, https://doi.org/10.1016/j.fuel.2022.123373.

[15] D. Feng, Z. Chen, K. Wu, J. Li, X. Dong, Y. Peng, X. Jia, X. Li, D. Wang, A comprehensive review on the flow behaviour in shale gas reservoirs: multi-scale, multi-phase, and multi-physics, Can. J. Chem. Eng. 100 (11) (2022) 3084–3122, https://doi.org/10.1002/cjce.24439.

[16] H. Su, Y.H. Hu, Recent advances in graphene-based materials for fuel cell applications, Energy Sci. Eng. 9 (7) (2021) 958–983, https://doi.org/10.1002/ese3.833.

[17] F. Yang, J. Gu, L. Ye, Z. Zhang, G. Rao, Y. Liang, K. Wen, J. Zhao, J.B. Goodenough, W. He, Justifying the significance of Knudsen diffusion in solid oxide fuel cells, Energy 95 (2016) 242–246, https://doi.org/10.1016/j.energy.2015.12.022.

[18] A. Filinov, M. Bonitz, D. Loffhagen, Microscopic modeling of gas-surface scattering. I. a combined molecular dynamics-rate equation approach, plasma sour, Sci. Technol. 27 (6) (2018) 064003, https://doi.org/10.1088/1361-6595/aac61e.

[19] A. Filinov, M. Bonitz, D. Loffhagen, Microscopic modeling of gas-surface scattering: II. application to argon atom adsorption on a platinum (111) surface, Plasma Sources Sci. Technol. 27 (6) (2018) 064002, https://doi.org/10.1088/1361-6595/aac620.

[20] J.V. Barth, Transport of adsorbates at metal surfaces: from thermal migration to hot precursors, Surf. Sci. Rep. 40 (2000) 75–149, https://doi.org/10.1016/S0167-5729(00)00002-9.

[21] T. Ala-Nissilayz, R. Ferrando, S.C. Ying, Collective and single particle diffusion on surfaces, Adv. Phys. 51 (3) (2002) 949–1078, https://doi.org/10.1080/00018730110107902.

[22] G. Antczaka, G. Ehrlich, Jump processes in surface diffusion, Surf. Sci. Rep. 62 (2007) 39–61, https://doi.org/10.1016/j.surfrep.2006.12.001.

[23] M. Logan, R.E. Stickney, Simple classical model for the scattering of gas atoms from a solid surface, J. Chem. Phys. 44 (1) (1966) 195, https://doi.org/10.1063/1.1726446.

[24] M. Logan, James C. Keck, Classical theory for the interaction of gas atoms with solid surfaces, J. Chem. Phys. 49 (2) (1968) 860, https://doi.org/10.1063/1.1670153.

[25] T. Yan, W.L. Hase, J.C. Tully, A washboard with moment of inertia model of gas-surface scattering, J. Chem. Phys. 120 (2004) 1031–1043, https://doi.org/10.1063/1.1628674.

[26] N. Mateljevic, J. Kerwin, S. Roy, J.R. Schmidt, J.C. Tully, Accommodation of gases at rough surfaces, J. Phys. Chem. C 113 (2009) 2360–2367, https://doi.org/10.1021/jp8077634.

[27] K. Yamamoto, H. Takeuchi, T. Hyakutake, Scattering properties and scattering kernel based on the molecular dynamics analysis of gas-wall interaction, Phys. Fluids 19 (2007) 087102, https://doi.org/10.1063/1.2770513.

[28] J.C. Maxwell, On stresses in rarefied gases arising from inequalities of temperature, Phil. Trans. Roy. Soc. Lond. 170 (1879) 231–256.

[29] C. Cercignani, M. Lampis, Kinetic models for gas-surface interactions, Transport Theor. Stat. Phys. 1 (1971) 101–114, https://doi.org/10.1080/00411457108231440.

[30] Y. Chen, L. Gibelli, M.K. Borg, Impact of random nanoscale roughness on gas-scattering dynamics, Phys. Rev. E 109 (2024), https://doi.org/10.1103/PhysRevE.109.065308.

[31] T. Liang, Q. Li, W. Ye, A physical-based gas-surface interaction model for rarefied gas flow simulation, J. Comput. Phys. 352 (2018) 105–122, https://doi.org/10.1016/j.jcp.2017.08.061.

[32] M. Smoluchowski, Zur kinetischen theorie der transpiration und diffusion verdißneter Gase, Ann. Phys. (Leipzig) 33 (1910) 1559–1570.

[33] M. Knudsen, Die Gesetze der Molekularstromung und der inneren Reibungsstromung der Gase durch Rohren, Ann. Phys. (Leipzig) 28 (1909) 75–130.

[34] G. Arya, H.C. Chang, E.J. Maginn, Knudsen diffusivity of a hard sphere in a rough slit pore, Phys. Rev. Lett. 91 (2) (2003) 026102, https://doi.org/10.1103/PhysRevLett.91.026102.

[35] T. Liang, Q. Li, Accurate modeling of Knudsen diffusion in nanopores using a physical-based boundary model, J. Appl. Phys. 126 (8) (2019) 084304, https://doi.org/10.1063/1.5110197.

[36] G.P. Brivio, T.B. Grimley, Dynamics of adsorption/desorption at solid surfaces, Surf. Sci. Rep. 17 (1993) 1–84, https://doi.org/10.1016/0167-5729(93)90023-I.

[37] J.A. Barker, D.J. Auerbach, Gas-surface interactions and dynamics, Surf. Sci. Rep. 4 (1984) 1–99, https://doi.org/10.1016/0167-5729(84)90005-0.

[38] V. Buch, Sticking probability of light atoms on the graphite (0001) surface, J. Chem. Phys. 91 (1989) 4974–4984, https://doi.org/10.1063/1.456738.

[39] P. Taborek, L.J. Senator, Helium on graphite: low-temperature desorption kinetics and sticking coefficient, Phys. Rev. Lett. 56 (1986) 628, https://doi.org/10.1103/PhysRevLett.56.628.

[40] J. Böheim, Sticking and accommodation of He on metal surfaces at low temperatures, Surf. Sci. 148 (1984) 463–477, https://doi.org/10.1016/0039-6028 (84)90593-4.

[41] M. Sivani, P. Taborek, D. Goodstein, Direct and thermal desorption of ⁴He films, Phys. Lett. 95A (1983) 59, https://doi.org/10.1016/0375-9601(83)90782-X.

[42] T. Brunner, W. Brenig, Elastic scattering, sticking and accommodation of noble gas atoms at a cold surface, Surf. Sci. 291 (1993) 192–206, https://doi.org/10.1016/0039-6028(93)91491-7.

[43] Z.W. Gortel, H.J. Kreuzer, M. Schaff, G. Wedler, Time of flight spectra in thermal desorption of helium, Surf. Sci. 134 (1983) 577–600, https://doi.org/10.1016/0039-6028(83)90060-2.

[44] G.P. Brivio, T.B. Grimley, V. Bortolani, G. Santoro, Theory of thermal desorption of light noble gases from metals: application to Ne-Cu(111), Surf. Sci. 245 (1991) 17–28, https://doi.org/10.1016/0039-6028(91)90464-4.

[45] Z. Medina, B. Jackson, Quantum studies of light particle trapping, sticking, and desorption on metal and graphite surfaces, J. Chem. Phys. 128 (2008) 114704, https://doi.org/10.1063/1.2890043.

[46] B. Jackson, Reduced density matrix approach for gas-solid scattering, Chem. Phys. Lett. 270 (1997) 484–490, https://doi.org/10.1016/S0009-2614(97)00404-1.

[47] T. Ohba, Excess adsorption of helium in extremely narrow slit pores, J. Low Temp. Phys. 177 (2014) 274–282, https://doi.org/10.1007/s10909-014-1214-5.

[48] T. Ohba, Limited quantum helium transportation through nano-channels by quantum fluctuation, Sci. Rep. 6 (2016) 28992, https://doi.org/10.1038/srep28992.

[49] A.W. Hauser, A.O. Mitrushchenkov, M.P. de Lara-Castells, Quantum nuclear motion of helium and molecular nitrogen clusters in carbon nanotubes, J. Phys. Chem. C 121 (2017) 3807–3821, https://doi.org/10.1021/acs.jpcc.6b12959.

[50] F.S. Emami, V. Puddu, R.J. Berry, V. Varshney, S.V. Patwardhan, C.C. Perry, H. Heinz, Force field and a surface model database for silica to simulate interfacial properties in atomic resolution, Chem. Mater. 26 (2014) 2647–2658, https://doi.org/10.1021/cm500365c.

[51] M.F. Döpke, J. Lützenkirchen, O.A. Moultos, B. Siboulet, J.F. Dufreche, J. T. Padding, R. Hartkamp, Preferential adsorption in mixed electrolytes confined by charged amorphous silica, J. Phys. Chem. C 123 (27) (2019) 16711–16720, https://doi.org/10.1021/acs.jpcc.9b02975.

[52] B.W.H. Van Beest, G.J. Kramer, R.A. Van Santen, Force fields for silicas and aluminophosphates based on Ab initio calculations, Phys. Rev. Lett. 64 (1990) 1955–1958, https://doi.org/10.1103/PhysRevLett.64.1955.

[53] J. Geske M. Vogel, Creating realistic silica nanopores for molecular dynamics simulations, Mol. Simul. 43 (1) (2017) 13–18, https://doi.org/10.1080/08927022.2016.1221072.

[54] M. Bartolomei, E. Carmona-Novillo, M.I. Hernández, J. Campos-Martínez, F. Pirani, Global potentials for the interaction between rare gases and graphene-based surfaces: an atom–bond pairwise additive representation, J. Phys. Chem. C 117 (2013) 10512–10522, https://doi.org/10.1021/jp401635t.

[55] P. Magnico, Q.D. To, Collisions, adsorption and self diffusion of gas in nanometric channels by molecular dynamics and stochastic simulation and the case of helium gas in graphitic slit pore, Int. J. Heat Mass Tran. 214 (2023) 124371, https://doi.org/10.1016/j.ijheatmasstransfer.2023.124371.

[56] L. Zhang, Z. Song, B. Zhao, E. Villarreal, H. Ban, Fast atom effect on helium gas/graphite interfacial energy transfer, Carbon 161 (2020) 206–218, https://doi.org/10.1016/j.carbon.2020.01.058.

[57] O. Talu, A.L. Myers, Reference potentials for adsorption of helium, argon, methane, and krypton in high-silica zeolites, Colloids Surf. A Physicochem. Eng. Asp. 187–188 (2001) 83–93, https://doi.org/10.1016/S0927-7757(01)00628-6.

[58] W.A. Steele, The physical interaction of gases with crystalline solids I. Gas-solid energies and properties of isolated adsorbed atoms, Surf. Sci. 36 (1973) 317–352, https://doi.org/10.1016/0039-6028(73)90264-1.

[59] H. Sha, S. Zhang, R. Faller, Molecular investigation of gas adsorption, separation, and transport on carbon nanoscrols: a combined grand canonical monte carlo and molecular dynamics study, Carbon 132 (2018) 401–410, https://doi.org/10.1016/j.carbon.2018.02.078.

[60] S. Plimpton, Fast parallel algorithms for short-range molecular dynamics, J. Comput. Phys. 117 (1995) 1–19, https://doi.org/10.1006/jcph.1995.1039.

[61] P. Langevin, Sur la théorie du mouvement Brownien, Comptes rendus de l'académie des sciences de Paris (1908) 530–533.

[62] P. Spijker, A.J. Markvoort, S.V. Nedea, P.A.J. Hilbers, Computation of accommodation coefficients and the use of velocity correlation profiles in molecular dynamics simulations, Phys. Rev. E 81 (2010) 011203, https://doi.org/10.1103/PhysRevE.81.011203.

[63] T. Ewart, P. Perrier, I. Graur, J.G. Méolans, Tangential momentum accommodation in microtube, Microfluid. Nanofluidics 3 (2007) 689–695, https://doi.org/10.1007/s10404-007-0158-3.

[64] B.-Y. Cao, M. Chen, Z.-Y. Guo, Temperature dependence of the tangential momentum accommodation coefficient for gases, Appl. Phys. Lett. 86 (2005) 091905, https://doi.org/10.1063/1.1871363.

[65] Y. Chen, J. Li, S. Datta, S.Y. Docherty, L. Gibelli, M.K. Borg, Methane scattering on porous kerosen surfaces and its impact on mesopore transport in shale, Fuel 316 (2022) 123259, https://doi.org/10.1016/j.fuel.2022.123259.

[66] C. Tang, G. Sun, Y. Liu, Is hydrogen diffusion in amorphous metals non-Arrhenian? Int. J. Hydrogen Energy 47 (2022) 9627–9634, https://doi.org/10.1016/j.ijhydene.2022.01.023.

[67] C. Donati, S.C. Glotzer, P.H. Poole, W. Kob, S.J. Plimpton, Spatial correlations of mobility and immobility in a glass-forming Lennard-Jones liquid, Phys. Rev. E 60 (1999) 3107, https://doi.org/10.1103/PhysRevE.60.3107.

[68] M. Micoulaut, Y. Guissani, B. Guillot, Simulated structural and thermal properties of glassy and liquid germania, Phys. Rev. E 73 (2006) 031504, https://doi.org/10.1103/PhysRevE.73.031504.

[69] Y. Lu, X.Y. Liu, G.H. Hu, Double-spring model for nanoparticle diffusion in a polymer network, Macromolecules 55 (11) (2022) 4548–4556, https://doi.org/10.1021/acs.macromol.2c00257.