The Impact of Microstructure Geometry on the Mass Transport in Artificial Pores : A Numerical Approach

Themicrostructure of porous materials used in heterogeneous catalysis determines the mass transport inside networks, which may vary over many length scales. The theoretical prediction of mass transport phenomena in porous materials, however, is incomplete and is still not completely understood. Therefore, experimental data for every specific porous system is needed. One possible experimental technique for characterizing the mass transport in such pore networks is pulse experiments. The general evaluation of experimental outcomes of these techniques follows the solution of Fick’s second law where an integral and effective diffusion coefficient is recognized. However, a detailed local understanding of diffusion and sorption processes remains a challenge. As there is lack of proved models covering different length scales, existing classical concepts need to be evaluated with respect to their ability to reflect local geometries on the nanometer level. In this study, DSMC (Direct Simulation Monte Carlo) models were used to investigate the impact of pore microstructures on the diffusion behaviour of gases. It can be understood as a virtual pulse experiment within a single pore or a combination of different pore geometries.


Introduction
Porous materials are primarily used as catalysts in heterogeneous catalysis.The activity of a catalyst is correlated with its number of active centers, which are located mainly inside a complex pore network.To increase activities accordingly, the inner surface could be increased (for the identical active site density) by reducing the pore diameter and increasing the amount of pores at the same time.The resulting smaller individual pores cause decrease of the rate of mass transfer in the pores.The rate of mass transfer in small pores thus has a great impact on the overall rate of reactions in heterogeneous catalysis.
The understanding of local interactions between diffusive and reactive processes on large inner surfaces, however, is incomplete.This is especially true in highly unstructured porous materials, arising, for example, in support materials prepared by template-assisted.The materials produced in this way exhibit wide geometric complexity.
Regarding complex geometries there are two approaches for describing diffusion phenomena theoretically: first the complex pore system is simplified to a homogeneous material with effective diffusion coefficients (see [1]).Second, local scale limits are defined and used in multiscale models (see [2][3][4][5]).Analysis of the local mass transfer processes occurring at individual geometric characteristics (pore constrictions, pore openings, and transition from macropores to mesopores) is not possible with either of the two approaches.So it is rather difficult to predict the fluid flow around a local edge with only a few nanometers length or its impact on the global mass transfer.
Similar models trying to separate the diffusion in different regimes can be found in the literature.For example, in (1) given by the Pollard and Present approach [6] the effective diffusion coefficient is calculated from diffusion coefficients of various diffusion regimes that can arise depending on the pressure conditions and pore sizes in the porous material: (1) Another approach is given by Crank [7] and provides two equations, which describe the effective diffusion coefficients in combinations of individual layers of a composite transverse and parallel to the flow direction, resulting in series and parallel formulas (see ( 2) and ( 3)): serial: parallel: Here  is the length of the composite,   the length of a layer ,  is the cross-sectional area of the composite, and   the cross-sectional area of a layer .
The porous solid may be divided into thin layers perpendicular to the flow direction and for each of them an effective diffusion coefficient is calculated from the parallel formula.Summing-up all these coefficients with the serial formula results in an effective diffusion coefficient for the entire porous solid.This method is described by Crank [7] as a series-parallel model.The reverse procedure is referred to as a parallel-series model by Crank [7].Here, the solid is first divided into thin layers parallel to the flow direction to get an effective diffusion coefficient for the parallel layers with the serial formula.Afterwards, these layers are combined with the parallel formula to get the total effective diffusion coefficient for the porous solid.
In Bell and Crank [8] weighted balance of the seriesparallel and parallel-series approximation was assumed to predict the ratio of the diffusion coefficient (see ( 4)): In addition to these models, in which the effective diffusion coefficient is calculated from diffusion coefficients of individual regimes and the geometry of individual units, there are further macroscopic approaches.In these models a uniform medium is used to describe the diffusion in heterogeneous structured porous material.For molecules the space for movement is restricted to the pore volume.Therefore, the porosity factor is introduced.All other geometric factors such as effects of the pore connectivity, the extension of the distance traveled by labyrinth effects, the roughness of the pore walls, or variations in pore diameter are expressed by the tortuosity.As such models in [1] describe the "effective medium approximation" [9][10][11], the "parallel pore model" [12], or the "random pore model" [13].
All theoretical models presented so far make use of geometric assumptions and simplifications and are not appropriate to illustrate the internal, local mass transport under the influence of global transport processes in complex structured pore networks.However, these local mass transport processes affect the reaction rates at the individual active surface species.Therefore, the understanding of local geometry, diffusion interactions in solid catalysts, is of central importance.To gain such insights, approaches that are needed are largely geometry and scale independent.The direct numerical simulation of single molecules using Monte Carlo methods (Direct Simulation Monte Carlo (DSMC)) is one such approach and should therefore be the basis for all investigations under this proposed research project.

Direct Simulation Monte Carlo (DSMC) Method
The DSMC method is based on the Boltzmann equation, a general equation for the description of transport processes in gas dynamics, which describes the velocity distribution of molecules in a fluid (see [14]).
Using the Boltzmann equation, the moving directions and velocities of individual molecules can be predicted, wherein interactions of single molecules with one another and with the system limits are in the form of shocks taken into account.External forces (e.g., gravity) may also be involved.Thus, in gas dynamics for all nonequilibrium and equilibrium states, the Boltzmann equation describes the course of the molecular velocity distribution of a fluid.The molecular velocity distribution can be used for determination of all macroscopic quantities and transport coefficients such as density, pressure, thermal conductivity, and diffusion coefficients.
One possible solution method for the Boltzmann equation is the Monte Carlo method, DSMC [15], which has mainly been developed for rarefied gases but is not limited to this application (see [16]).It can be also extended to gases under nonideal conditions [17].The lower spatial limit is about 10 nm and with adjustments up to 1 nm (of the order of atomic radii).An upper spatial limit does not exist.One drawback of this method is the high computational effort, which can be decreased through the introduction of molecular packets representing a certain number of real molecules.
Pores and pore networks contain only a small number of molecules under normal conditions and therefore the DSMC method can be used for such systems.Then, two additional constraints arise from purely practical circumstances: the insufficient statistics at the small number of molecular packets in a control volume as well as the computation time at too high number of molecular packets.
In literature, the problem of the statistical accuracy of a low molecular packet number has been studied (see, e.g., [18]) and a minimum value of 5 packets per molecule control volume is recommended.The number of molecules that represents a molecular packet must be adjusted accordingly.If the number of control volumes is great, the total number of molecular packets must also be great to meet this requirement.With a small ratio of the mean free path of the gas molecules to the geometrical reference length (i.e., small Knudsen numbers) the molecules density increasesalso leading to higher molecular packet numbers.
All DSMC simulations were done with a modified code based on the DSMC-implementation within OpenFOAM [19] with neon as monatomic gas species with a variable hard sphere model (mass = 33.5E− 27, atom diameter = 2.77E − 10, internal degrees of freedom = 0, and viscosity index = 0.81) that was used.The number of molecules was between 1E + 6 and 2E + 6 represented by 1E + 5 and 2E + 5 particles.All simulations were done under isothermal conditions with a temperature of 323 K and a Maxwellian thermal wall model.

Validation of DSMC Results with Fick's Second Law
For validating the DSMC-code, the characteristic but reduced systems shown in Figure 1 were designed to represent different pore geometries.In the case of the "bisected" system an inner wall is placed to divide the big pore into two half cylinders-introducing an additional inner wall.The pore diameter of this "bisected" pore is therefore the same as in system "big." The simulation procedure starts with the initialization of a small region directly adjacent to the inlet area with a specific number of molecules and then tracks their movements and collisions over the time.The result is a time-dependent molecule distribution inside the pores and the flow rate at the outlet, which is located at the end of the pores.
The pulse response curve, that is, the flow rate of molecules at the outlet, can be described with Fick's second law as shown in (5) for a one-dimensional case with a constant diffusion coefficient.The solution of this equation can either be done analytically as given in [20] or numerically with a suitable software.In this paper, both procedures were applied.The software "Berkeley Madonna" [21] was used for the numerical solution.The fitting procedure was implemented by minimizing the sum of the squares of errors by changing the effective diffusion coefficient: The results of these fitting procedures are shown in Figure 2.
It is obvious that all three systems of geometries shown in Figure 1 can be represented by an effective diffusion coefficient.The pulse response curve of the "bisected" system is situated between the "big" and "small" systems showing the influence of the area/perimeter ratio of the flow cross section and thus reflects the correct physical behavior.
For the very simple case of a cylindrical pore with a given radius, the effective Knudsen diffusion coefficient can be calculated using the analytical expression according to  The analytical (see (6)) and numerical (DSMC) results can then be compared with each other, as shown in Table 1.
Because the pore in the system "bisected" is not of a cylindrical shape, (7) was used to calculate an equivalent pore radius following the definition of the hydraulic diameter in fluid mechanics: With this the deviations between ( 6) and DSMC can be considered to be sufficiently small for all systems.Summarizing these results, it can be concluded that the DSMC model is able to simulate the transport of monatomic gases in small arbitrary geometries with sufficient precision and can be further used to analyze more complex geometries.

Investigation of Linear Pore Combinations
To investigate combinations of pore geometries, two combinations ("big-small" and "small-big" shown in Figure 3) were modelled.The resulting pulse response curves simulated with DSMC and the curves obtained by the fitting procedure with 2nd Fick's law are illustrated in Figure 4. Additionally, the results of "small" and "big" (see Figures 2 and 3) are shown for comparison.

Big-small
Small-big It is obvious (see Figure 5) that DSMC in case of "bigsmall" and "small-big" differs significantly from the analytical result using Fick's 2nd law.The pulse responses deviate substantially from each other.This observation can be discussed as follows.Geometries representing pore combinations with a sudden diameter change cannot be represented by one characteristic effective diffusion coefficient as it would result from Fick's 2nd law.Any calculation of such geometry using the analytical equation fails to reproduce the real behaviour and it is not possible to describe geometry by a characteristic effective diffusion coefficient.Nevertheless, the curves of the DSMC results can be analyzed further.Comparing the diffusion rates, it can be seen that the diffusion rate in the system "big-small" is slower than in all of the other systems-even "small." This indicates that molecules are moving slower in geometries with a transition from a big to a small diameter compared to geometries with a constant small diameter.Also, the system "small-big" reaches diffusion rates comparable with the system "big, " which is also an interesting observation.However, the general trend would have been also reflected by Fick's 2nd law, thus giving tendencies.Figure 5: Molecule distributions at 0.00084 s for cases "small" (a), "big" (b), "small-big" (c), and "big-small" (d).All geometries were scaled by a factor of 100 in the flow direction.DSMC makes it possible to visualize distributions (see Figure 6).Thus, a detailed view with respect to different time steps is possible and may help to understand the deviations, which have been found.

Bisected-big
At comparable time steps, all densities differ and show significant different profiles.At this specific time step, the densities of molecules for the systems "big" and "small" are comparable and system "small-big" results in lower densities, whereas system "big-small" ends at much higher densities.Apparently, the sudden diameter reduction in "big-small" system represents some kind of obstacle, which influences the mass transport significantly.
To further investigate the influence of sudden diameter changes, two pore widths have been combined.These pore combinations are based on the system "bisected, " shown in Figure 6.Both systems consist of a constant identical diameter as system "big" but with an inner wall as system "bisected" in one half length of the pore.Therefore, a diameter  change was achieved without creating an obstacle in the geometry.
Figure 7 shows the pulse response curves in comparison to the results of the systems "big" and "bisected" (see Figure 3).All pore geometries possess similar properties, for example, identical diameters and length.Only the inner wall starts at different positions.As the inner wall is no obstacle which reduces the diameter in a significant way, Fick's 2nd law represents the transport adequately.
The two new systems (see Figure 7) result in response curves proceeding between those of the two original geometries.However, both new systems behave differently, thus indicating the influence of the position of the inner wall.
Effective diffusion coefficients have been evaluated for the systems "big-bisected" and "bisected-big" and are compared with the results of the serial formula of Crank [7] (see (2)-combining the effective diffusion coefficients of the systems "big" and "bisected").Table 2 summarizes the values and their deviations.The serial equation does not reflect sequences and so the effective diffusion coefficients for both combinations are equal.The effective diffusion coefficients obtained by evaluating the DSMC curves and by using the serial formula of Crank [7] show only small differences with opposite signs of the deviation.The DSMC results show an influence of sequence of the pore geometries and two different effective diffusion coefficients could be derived.It can be concluded that the serial formula from Crank [7] delivers reasonable results; however, it is not able to consider effects of sequence-contribution in pore geometries.Therefore an over-or underestimation of effective diffusion coefficients might occur, which would influence the local behaviour.As all microkinetic steps are highly nonlinear, the error in estimating would be nonpredictable.

Local Effects at Pore Diameter Changes
To further investigate the influence of diameter changes or the transition shapes on the response curves, pore geometries illustrated in Figure 8 were designed and simulated with DSMC.In the system "big-obstacle1-big" an infinitely thin wall with a centric hole is placed at the middle of the pore.An increased thickness of the obstacle is used in the system "big-obstacle2-big" to get kind of an orifice.
Figure 9 compares the resulting response curves of the systems "big-cone-small" and "small-cone-big" with the ones from systems "big-small" and "small-big." It can be seen that the curves of systems with cone-shapes are very similar to the curves of the corresponding systems with the sudden change of the diameter, that is, without a cone-shape.This indicates that cone-shaped diameter changes do not influence the results.Thus, the angle of the wall, separating the two diameters from each other, seems to have no influence on the investigated geometries.
Figure 10 illustrates the effect of the length of the obstacle in the pore on the simulated pulse response curve.The result of the infinitely thin wall (system "big-obstacle1-big") is identical to the response of the system "big." The curve of the system "big-obstacle2-big" is slightly shifted towards the Normalized intensity Time (s)  response of the system "small." These results indicate that an infinitesimal thin obstacle does not influence the diffusion rate and with increasing length of the obstacle the mass transport is decreased.
It can be concluded that the decrease of the diffusion coefficient in case of reduced pore diameters is not like a bulging effect, but more of a plugging effect.If the diameter gets smaller, the local concentration is increased, which hinders following molecules and even more induces a backward movement of molecules.The driving force of diffusion in the forward direction is lowered.
Conversely, if the pore diameter is suddenly increased, the local concentration gets smaller, which boosts the driving force of diffusion and makes the movement of molecules in the forward direction much higher.

Summary and Conclusion
Within this work the diffusion in pore geometries and pore geometry combinations were investigated by numerical experiments with the DSMC method and compared with Fick's 2nd law.The model allows simulating mass transport in small pores with arbitrary size and transition between different pore diameters-only limited by computing power.
It has to be noted that the knowledge of the pore diameter alone is not sufficient to predict mass transport in heterogeneous pore networks but also the knowledge of the connectivity and local pore diameter distributions is also necessary.Furthermore, the kind and shape of the transition between different pore geometries have crucial impact on the diffusion processes and must be regarded.
Ongoing simulations consider geometries where the cone-shape is extended to the whole length of the pore.Also increasing the length of the obstacle in the system "big-obstacle2-big" will be regarded in future work.The question for this is as follows: can the attenuation of diffusion at diameter-reduction be equalized by a diameter-enhancement downstream?
The obtained results are valuable for discussion of flow phenomena at molecular dynamics level and allow a better insight to diffusion in small pores.

Figure 2 :
Figure 2: Simulated (DSMC) pulse response curves in comparisons to Fick's 2nd law for three model geometries.

Figure 7 :
Figure 7: Pulse response of half-cylindrical shaped pore combinations and the individual pores.

Figure 9 :Figure 10 :
Figure 9: Pulse response curves for the investigation of the effect of a cone-shaped diameter change.

Nomenclature𝐶:
Concentration in mol/m 3 : Pressure in Pa : Ideal gas constant 8.3145 J/(mol⋅K) : Temperature in K Kn: Knudsen number Λ: Mean free path length in m : Characteristic length of geometry in m   : Boltzmann-constant 1.3406488 − 23 J/K : Diffusion coefficient in m 2 /s : Molar mass in kg/kmol.Greek Letters  : Atom diameter in m.

Table 1 :
Effective diffusion coefficients from DSMC and theoretical calculations.