1 Introduction

Open-cell foams (OCFs) have been used in several industrial applications such as solar energy receivers (Barreto et al. 2020), heat exchangers (Plant and Saghir 2021), heat sinks (Marri and Balaji 2021), thermo-acoustic applications (Caniato et al. 2021), and as monolithic catalyst supports in chemical reactors (Sinn et al. 2021; Frey et al. 2017; Kiewidt and Thöming 2019a). OCFs' complex and porous structure allows for superior performance in the mentioned applications. For instance, OCFs are known as a good choice as monolithic catalyst supports because of low-pressure drop, high thermal conductivity, high catalyst inventory due to large specific surface area, good structural stability, and their excellent heat and mass transfer with fluid (Kiewidt and Thöming 2019b; Ho, et al. 2019). OCFs are generally categorized based on their nominal density of pores per inch (PPI). These pores are connected through channels, also called windows, and form the interconnected porous structure of an OCF (Bracconi et al. 2017). OCF's morphological properties control transport phenomena within them and thus their performance as catalyst support. Therefore, a comprehensive analysis of transport phenomena within OCFs requires detailed information about morphological properties such as open porosity (ε), pore diameter (dp), and window diameter (dw).

CFD simulations have been used for the investigation of transport phenomena within OCF including pressure drop (Inayat et al. 2016; Bracconi et al. 2018a), mass (Meinicke et al. 2017a; Habisreuther et al. 2009), and heat transport (Aghaei et al. 2017; Bracconi et al. 2018b) studies. Microcomputed tomography imaging (μCT) has been widely used to obtain a realistic structure of the OCFs which can then be utilized in pore-scale CFD simulations (Della Torre et al. 2016; Ranut et al. 2015; Sinn et al. 2019).

Mass transport in porous media includes convection and molecular diffusion. Contributing to the latter is concentration gradient-induced diffusion and self-diffusion contribute, both of which are ultimately caused by random molecular movements. While self-diffusion also occurs in fluids without concentration gradients and convection, concentration gradient-induced diffusion can be caused by local variations in flow velocity and direction on the macroscopic scale, resulting in convective dispersion (Dutta 2013; Jha et al. 2010). The ratio of convection and diffusion contribution to mass transport is described by the Péclet number (Pe). In the transition regime with 0.3 < Pe < 5, where the dominant transport mechanism shifts from diffusion to convection, gas flows can be described by a superposition of convection and diffusion resulting in shifted displacement distributions. In the power-law regime (5 < Pe < 300), convection is the dominant transport mechanism, although diffusion should not be ignored (Trunk et al. 2021).

Considering molecular interactions in macroscale computational fluid dynamics (CFD) simulations is currently computationally demanding. For example, simulating molecular diffusion inside the OCF utilized in this work (with a bulk volume of around 2.5 cm3) would require a meshing grid with billions of cells (Della Torre et al. 2014). As an alternative, it was suggested to limit the simulations to a representative elementary volume (REV), i.e., a unit cell of the OCF (Meinicke et al. 2017b), but this neglects the effect of neighboring pores. The impact of such a simplification was recently demonstrated via CFD simulations and magnetic resonance velocimetry (MRV) measurements under REV and full-field conditions which showed that the results of simulations based on the REV of the structure did not match the data of the full-field experiments. On the other hand, under full-field conditions and when the identical OCF is used in both simulations and experimental measurements, cross-validation was achieved for local parameters such as velocity and temperature fields and concentration maps of species (Sadeghi et al. 2020).

Several measurement techniques have been utilized to study gas dispersion. These are, for example, point-wise injection and tracer experiments (Benz et al. 1993; Olsson and Grathwohl 2007), microparticle imaging velocimetry (μPIV) (Winter et al. 2021), and laser-induced fluorescence (LIF) (Schubert et al. 2010). By using a transparent operating fluid such as gas combined with an opaque OCF, such measurements prove to be somewhat difficult and are limited in the number of tools that are suitable for such measurements. Nuclear magnetic resonance (NMR)-based measurement techniques have been developed for macroscale (Sadeghi et al. 2020; Ulpts et al. 2017) and microscopic analysis (Mirdrikvand et al. 2018) of transport phenomena within foam reactors. Tallarek et al. (1998) used pulsed-field gradient (PFG) NMR measurements to investigate dispersion coefficients in packed chromatographic columns in axial and transversal directions. Mirdrikvand et al. (2020) investigated gas dispersion within opaque OCFs via spatially resolved PFG NMR measurements of self-diffusion propagators, i.e., displacement distribution of gas flows. As a result, dispersion coefficients were obtained along axial and radial directions. This was done assuming that propagator measurements can be considered as a first-order approach to transport diffusivities.

Dispersion coefficients can also be obtained from analytical calculations (Salles et al. 1993; Edwards et al. 1991), Lagrangian and Eulerian random-walk particle tracking methods (Trunk et al. 2021; Salamon et al. 2006), or Navier–Stokes CFD simulation data used for tracking the fluid (Fan et al. 2016). Fluid dynamic dispersion of flow within OCFs was investigated in several numerical studies (Saber et al. 2012; Voltolina et al. 2017). As discussed before, the trilemma of grid resolution, size of the investigated geometry, and computational cost forced these conventional CFD simulations to compromise between the size of the examined geometry and the grid. For instance, Fan et al. (2016) used stationary Navier–Stokes finite-element CFD simulations for evaluation of the axial and radial fluid distribution within a 30 PPI OCF. They considered 100 massless particles on the streamline that passed through a 2-mm circle at the inlet of the foam reactor during 30 ms. Since diffusion was not considered in their simulations, the reported axial distribution graphs in their work did not show the movement of gas molecules in directions opposite to the flow caused by molecular interaction. The reported distribution graphs were also not smooth because the number of particles in the simulations was small. In addition, instead of performing tracer simulations over the converged stationary solution, they tracked 100 massless particles over the stationary streamlines within OCF during certain observation time. Such conventional stationary CFD simulations cannot predict the realistic distribution propagators since the movement of particle is limited to the convection. Chandra et al. (2019) performed direct numerical simulations (DNS) to investigate the axial dispersion coefficients of flow within representative elementary volumes (REV) of OCFs based on Kelvin’s ideal unit cell. Also, they compared simulation results with calculations based on the Taylor–Aris dispersion equation (Salles et al. 1993; Edwards et al. 1991). They studied the effect of flow rate, porosity, and flow direction on dispersion coefficients. Parthasarathy et al. (2013) investigated air flow longitudinal dispersion coefficients within OCFs via pore-scale Navier–Stokes standard tracer simulations. They implemented the stationary solutions as the initial condition and performed tracer simulations. Their simulation results were compared with the experimental data from the literature on conventional packed beds. Pore-scale CFD simulations of transport phenomena within OCFs usually employ a mesh grid with millions of meso-sized cells (~ 100–200 µm).

Conventional tracer simulations which implement the stationary solution of Navier–Stokes simulations as initial condition cover diffusion driven by concentration gradients. However, without concentration gradients being present, mass transport by thermal motion is spatially random and not covered by Fick’s laws. This self-diffusion of gas atoms or molecules is, therefore, not directly accounted for in conventional Navier–Stokes-based tracer equations.

We hypothesize that it is this missing consideration of self-diffusion that can cause an underestimation of diffusion's contribution to dispersion coefficients obtained from conventional Navier–Stokes-based tracer simulations. This can be expected in resting fluids and in near-wall regions where transport is dominated by diffusion. Therefore, due to this incapability of resolving such microscopic molecular interactions, conventional Navier–Stokes-based tracer simulations alone are not adequate for predicting flow transport in flow regimes with a pronounced contribution of diffusion (i.e., for Pe < 16). This is especially true for flows within porous structures such as OCFs which encounter high solid–fluid interface density and near-wall regions.

The combination of standard pore-scale CFD simulations with particle tracking (PT) solvers can overcome this issue. In such coupled simulations, pore-scale velocity fields are provided to the particle tracking solver which then performs accurate microscopic tracking of molecules during a certain observation time. This allows us to consider microscopic molecular movement of gas in dispersion analysis without the implementation of a micro-sized computational network. For instance, Manz et al. (1999) investigated the longitudinal and transversal displacement distributions of water flow via random-walk particle tracking in conjunction with lattice Boltzmann (LB) simulations and NMR measurements. They used three different packings of spheres and analyzed the influence of varying the observation time from Δ = 20, 50, and 100 ms to Δ = 200, 500, and 1000 ms by direct comparison of experimental and simulated displacement distributions. Although the simulations and experiments agreed well, the authors restricted the study to the Darcy flow regime with very low Re numbers (0.4 < Re < 0.77) and only to one flow rate per sample. Hussain et al. conducted a cross-validation study comparing the results of NMR experiments with lattice Boltzmann and pore network simulations for water transport in Bentheimer sandstone (Hussain et al. 2013). Freund et al. (2005) studied the pore-scale mass transport of isothermal flows within a packed bed reactor by performing LB simulation in conjunction with a particle tracking method. Firstly, they compared the obtained averaged velocity profile and porosity from LB simulations with experimental NMR measurement data. Subsequently, they coupled the LB simulation of flow within a simple cube packing with a particle tracking module to calculate the axial and radial dispersion coefficients. The comparison of calculated dispersion coefficients from simulations and their corresponding experimental tracer measurements showed good agreement.

Svidrytski et al. (2019) utilized the LB method and a random-walk particle tracking module to analyze longitudinal dispersion coefficients in particulate beds of different heterogeneity and morphology. Brosten et al. (2011) used NMR measurements in conjunction with LB simulations with random-walk PT to study the dispersion of water flow within polymer OCFs. However, they used different OCFs in simulations (50 PPI) and NMR measurements (80 and 110 PPI, respectively) and thus did not present a direct comparison between experiments and measurements. A direct comparison of NMR-measured displacement propagators of gas flow within OCF with the CFD-PT method is still missing in the open literature and therefore in the scope of our current contribution.

In this study, we present an experimental cross-validation of full-field CFD-PT simulations of methane flow within two different ceramic OCFs. CFD-PT simulations are performed based on µCT measurements of OCFs used in NMR measurements published earlier (Mirdrikvand et al. 2020) to analyze gas flow microscopic displacement distribution. Particularly, we hypothesize that the concept of random placement and stepwise particle tracking provides a better estimation of diffusion-dominated transport in near-wall regions with fluid-wall boundary conditions such as rejection and reflection. Using such wall boundary conditions enables CFD-PT solvers to predict gas molecules displacement in the opposite direction of flow and better estimation of diffusion in near-wall regions.

Four different flow rates in the Darcy (Rep < 1), Darcy–Forchheimer (1 < Rep < 10), and Forchheimer (10 < Rep < 150) regime are investigated to cross-validate simulations with PFG NMR measurements over a wide range of operating conditions.

2 Materials and Methods

For coupling pore-scale CFD simulations based on the Navier–Stokes equations with a random-walk particle tracking (PT) module within the OpenFOAM framework, a novel toolkit named disTrackFoam was recently introduced (Report et al. 2016). This represents a numerical tool to investigate microscopic displacement distribution and dispersion of gas flows within porous media, such as ceramic sponges or periodic open cellular structures (POCS) (Trunk et al. 2021). Here, we perform CFD-PT simulations using geometries based on µCT measurements of two OCFs. The µCT-based CFD simulations predict the pore-scale flow pattern, while the particle tracking module simulates the microscopic movement of the gas molecules during observation time. The velocity fields obtained from µCT-based CFD simulations were imported into the PT module. The PT module then considers a suitable number of massless particles moving through the flow fields. It tracks them during the respective observation time. In this work, the results of the coupled CFD and particle tracking (CFD-PT) simulation are compared with recently published PFG NMR data (Mirdrikvand et al. 2020) (Fig. 1).

Fig. 1
figure 1

Flowchart of the numerical and experimental methods used in this work to compare the displacement distribution of gas flows within open-cell foams

Full size image

Pore-scale CFD simulations were conducted to determine the flow pattern within OCFs. For this purpose, the OCF structure was reconstructed from μCT images using image processing software (ImageJ version 1.52a, https://imagej.nih.gov/ij/). The surface-triangulated geometry of the OCFs was prepared using suitable meshing software (MeshLab version 1.3.2, http://www.meshlab.net/). The morphological parameters of the OCFs such as open porosity (ε), pore diameter (dp), window diameter (dw), strut diameter (ds), and specific surface area (Sv) are reported in Table 1. The reconstructed triangulated structures were scaled and inserted into the tubular reactor via CAD software (FreeCAD version 0.18.4, https://www.freecadweb.org/). To prepare the computational grid, the meshing software CfMesh (version 1.1.2, https://cfmesh.com/) was utilized. Four different flow rates were investigated, ranging from 0.1 SL/min to 1.5 SL/min. Here, we calculated the Péclet number by using the pore diameter as the characteristic length (\({\text{Pe}}_{p} = ud_{{\text{p}}} /D_{{{\text{eff}}}}\)). For effective self-diffusion, the diffusion coefficient is considered a macroscopic property that reflects the transport behavior of particles without convection and concentration gradients present in the system. It may be influenced by microstructural factors such as pore size, tortuosity, and surface interactions. The effective self-diffusion coefficient (Deff) of methane, measured in the OCFs without flow, was reported as 1.73·10–5 and 1.53·10–5 m2/s for the 10 and 20 PPI OCFs, respectively (Mirdrikvand et al. 2018). The reported effective diffusion coefficients were utilized for calculating the Péclet numbers of flow within OCFs. The corresponding Péclet numbers were calculated with the OCF pore diameter. The Knudsen numbers of the system were calculated from Re numbers (\(Kn= \sqrt{\frac{\gamma \Pi }{2} }\frac{{\text{Ma}}}{{\text{Re}}}\)) as 1.34·10−5 and 2.08·10−5, respectively. Here γ is the ratio of heat capacities, and Ma is Mach number, i.e., the ratio of velocity to sound speed (Agrawal 2011).

Table 1 Morphological properties of the utilized open-cell foams of 25-mm length (Mirdrikvand et al. 2020)
Full size table

At the highest inlet velocity (\({u}_{{\text{max}}}=0.51\) mm/s), the Mach number is still lower than 0.2, so that the flow was considered incompressible. The flow regime was regarded laminar, since the highest pore Reynolds number (\({\text{Re}}_{{\text{p}}} = \rho u_{z} d_{{{\text{in}}}} /\mu _{f}\)) was 76.49. Pore-scale flow patterns were calculated by numerically solving the incompressible Navier–Stokes equation using the well-known SIMPLE algorithm within the OpenFOAM 4.1 CFD toolbox (simpleFoam). Grid independence tests were conducted, and the optimum cell numbers were determined for each OCF CFD case (Table 2). For this purpose, the pressure loss caused by foam structures at an inlet velocity of 0.051 m/s (1.5 SL/min) was determined for several grids with increasing cell numbers. The optimum number of cells was found for each OCF at which the pressure drop starts to become independent of the number of cells used (Figure S1). For the calculation meshes used, the resolution (i.e., mesh size) was about 130–190 µm.

Table 2 Information on CFD and particle tracking simulations
Full size table

A homogeneously distributed constant inlet velocity was considered as an inlet boundary condition. Also, a no-slip boundary condition (u = 0) was implemented for all solid surfaces, including the OCFs and the foam reactor wall. An ambient atmospheric pressure (p = 0 Pa gauge pressure) was considered at the outlet of the foam reactor to match the conditions during the NMR measurements (Jha et al. 2010), which served as a reference.

The actual particle displacement calculation was performed in a second simulation step with the OpenFOAM disTrackFoam toolkit (Trunk et al. 2021) based on the precalculated stationary flow field. The displacement of the non-interacting particles is described by a stochastic differential equation which contains a deterministic term for advection and a random diffusion term (modeling Brownian motion) (Kloeden and Platen 1992),

$$\vec{x}\left( {t + \Delta t} \right) = \vec{x}\left( t \right) + \underbrace {{\vec{u}\left( {\vec{x}\left( t \right),t} \right)\Delta t}}_{{{\text{Convection}}}} + \underbrace {{\vec{\zeta }\left( t \right)\sqrt {2D_{{{\text{mol}}}} \Delta t} }}_{{{\text{Diffusion}}}}$$
(1)

The calculation of the particle trajectory in the Eulerian approximation in the form of Eq. (1) for each time step Δt is dependent on a convective part caused by the particle velocity \(\overrightarrow{u}\) at the location \(\overrightarrow{x}\) and time \(t\), and the additional diffusive part described by a random walk, depending on the molecular diffusion coefficient \({D}_{{\text{mol}}}\) and a stochastic randomly oriented vector \(\overrightarrow{\zeta }\left(t\right)\) with a fixed length (\(\left|\overrightarrow{\zeta }\left(t\right)\right|=\sqrt{3}\)), zero mean, and unit variance (Kloeden and Platen 1992). The random vector is generated with a random number uniformly distributed within the interval [0, 1] for each of the three vector components via the OpenFOAM® random number generation method osRandomDouble. Afterwards, the three components of the vector are projected into the interval [ − 1, 1]. In the last step, the random vector is normalized to a length of \(\sqrt{3}\) (Freund et al. 2005; Heidig et al. 2014; Maier et al. 2000). To approximate a physically correct treatment of a solid wall contact, the stepwise particle tracking (SPT) of the disTrackFoam toolkit was selected to calculate diffusive particle movement in close vicinity to walls where the fluid velocity and, therefore, the convective particle velocity are zero due to wall adhesion. Using stepwise wall boundary conditions enabled the CFD-PT solver to predict gas transport from the wall surface caused by diffusion. Since flow within OCFs encounters considerable solid interfaces, consideration of such diffusive transports is necessary for accurate gas displacement prediction.

As illustrated in Fig. 2, the particle hitting a solid wall receives a recalculated random vector leading away from the actual surface, so that the particle leaves the surface with no flow velocity. This boundary condition offers the most accurate representation of reality.

Fig. 2
figure 2

Schematic description of the generation of particle trajectories for solid–wall interactions

Full size image

The number of particles as well as the molecular diffusion coefficient Dmol at ambient temperature and atmospheric pressure used in particle tracking simulations are listed in Table 2.

For the distribution method of the initial particle positions, a homogeneous and constant particle density within the ROI for both cases was targeted. The generation was conducted via the dedicated particle placement utility included within the disTrackFoam toolkit (Trunk et al. 2021) that considers the investigated simulation mesh for enabling various convenient placement strategies even in complex porous domains. The inherent tracking of the particles enables the disTrackFoam toolkit to conveniently provide the positions of all particles together with an—for each particle distinct—identification label at desired time steps as data output. Therefore, it is easily possible to determine the total displacement of each particle individually as the magnitude of the difference between the investigated start and end time positions.

To finally determine the full width at half maximum (FWHM) of the displacement distributions, the two known neighboring points of the curve were first determined at the height of the half maximum, and the curve was interpolated between them. The width of the curve could then be determined from the resulting coordinates.

3 Results and Discussion

To describe the mass transport characteristics of OCFs by performing CFD-PT simulations, we first investigated the axial displacement propagators. Such propagators can be understood as a description of particles that were initially placed at an origin position and distributed later longitudinally via diffusion. This leads to a symmetrical curve of the displacement probability, which can be represented by a normal distribution. Such a normal distribution can be observed when self-diffusion occurs or when identical concentration gradients are present starting from the original position in both directions. Using conventional stationary CFD simulations alone cannot present diffusional contribution in displacement distributions. To illustrate this incapability, we considered 6300 massless particles within ROI for flow rates of 0.5 and 1 SL/min in 10 PPI OCF. The particles were tracked over the streamline and during certain observation time equal to 10 ms. The comparison of the obtained axial displacement distribution of the tracked particles from the stationary CFD simulations with NMR data revealed their incapability and low degree of agreements because of neglecting the diffusional contribution (see Figure S2). The left tails of the obtained displacement propagators from stationary CFD simulations considerably presented a shorter displacement in opposite direction of flow compared to NMR data. Such considerable deviation confirmed the large contribution of diffusion in displacement in opposite direction of flow. By increasing the flow rate, the deviation of the right tail of the obtained displacement propagators from stationary CFD simulations from corresponding NMR data decreased due to increasing the contribution of convection in transport of particles. Overall, the deviation between the results of stationary CFD simulations and NMR data necessitated using CFD-PT simulation for a better quality of agreement.

As the first step, the steady-state flow field within 10 and 20 PPI OCFs was obtained by performing μCT-based pore-scale CFD simulations. In the following step, the motion of massless particles was tracked via PT solver over steady-state flow fields and within the determined region of interests (ROIs) during certain observation periods. Each tracked particle had different velocities due to their position in the flow field. Figure 3 illustrates the method used to observe particle evolution within OCF by CFD-PT simulations.

Fig. 3
figure 3

The graphical illustration of the evolution of particles by CFD-PT simulations during different time steps (Δ)

Full size image

The data on the particle velocities were experimentally obtained by PFG NMR measurements in an earlier work (Mirdrikvand et al. 2020). Here, the CFD-PT simulations were performed with identical parameters (flow rates, ROI, and observation time) and identical OCF structures as used in the experiments. This allowed us to compare simulated displacement distributions with experimental ones.

In a case close to restricted self-diffusion, a superimposed small convection was tested, and the data obtained at a flowrate of 0.1 SL/min are still almost normally distributed, as illustrated in the first curve in the upper left panel of Fig. 4. They reveal the dominating contribution of restricted self-diffusion to the displacement height and width (Fig. 5). As the flow velocity increases, the curve peaks widen, indicating an increasing share of convection in the dispersion.

Fig. 4
figure 4

Comparison between the obtained displacement distributions of gas flow within a 10 PPI open-cell foam at different flow rates of 0.1, 0.5, 1, and 1.5 SL/min from PFG NMR experimental measurements and CFD-PT simulations. Blue symbols represent CFD-PT data, and orange symbols represent NMR data. Observation time is 10 ms

Full size image
Fig. 5
figure 5

Effect of flow rate on the obtained displacement distribution of gas flows within a 10 PPI open-cell foam from a. PFG NMR measurements and b. CFD-PT simulations. Observation time is 10 ms

Full size image

A comparison of displacement distributions derived from PFG NMR and CFD-PT for a 10 PPI OCF and 10-ms observation time (Fig. 4) displays good agreement at all investigated flow rates. However, the right tail of the curves from the CFD-PT data tends to shift further to the right with increasing flow rate. In all these instances, the distributions are mainly caused by effective diffusion. This can be most clearly seen in the lowest convection tested (0.1 SL/min), which shows mainly effective self-diffusion in the pores with negligible convectional dispersion. Characteristics of this case are that the peak of the displacement probability curve occurs at zero displacement.

Here, the effective molecular diffusion is the molecular diffusion reduced by the factor of the pores' tortuosity.

The comparison of simulations and measurements revealed a deviation that increases with flow rate. This also applies to a similar comparison using a 20 PPI OCF (Figure S3). This is in agreement with earlier findings that NMR measurements underrepresent convection to a small extent which led to lower velocities (Sadeghi et al. 2020). This underrepresentation could be due to the NMR measurement method, where convection can cause labeled protons to leave the field of view too early before they can be measured.

In both CFD-PT and NMR data, for each flow rate, the illustrated displacement distributions of 10 PPI OCF were broader than their corresponding ones within 20 PPI OCF. Comparing the displacement distributions of 20 PPI OCF from simulations with 10 PPI, the right tails of the distribution graphs were less shifted to the flow direction by an increasing flow rate. Also, the peak locations where maximum probably occurred were shifted toward the flow direction. As a result of their higher porosity (smaller shift of the right tail) and smaller window sizes (larger shift of peak position), these two differences may be explained by more solid interfaces and higher velocity within 20 PPI OCF, respectively. NMR and CFD displacement distributions shift to the right and become broader with rising flow rates. This is shown for a 10 PPI OCF during 10 ms of observation time in Fig. 5. The maximum displacement probability also decreases with the flow rate. By increasing the flow rate, convection's contribution to dispersion increases and becomes dominant compared to diffusion. This decreases the number of stagnant molecules with zero displacement and increases the number of molecules with higher displacement in the flow direction.

To investigate the effect of observation time Δ on the quality of agreement, the displacement distributions for the flow rate of 1.5 SL/min within the 10 PPI OCF were compared to the corresponding experimental ones at different observation times of 5 and 10 ms. By reducing the observation time from 10 ms (Fig. 6b) to 5 ms (Fig. 6a), the deviation between simulation and measurements decreases from 0.21 to 0.089 mm. In PFG NMR measurements, molecules are labeled and tracked over a certain observation time within a region of interest (ROI). At shorter observation times, less labeled molecules leave the ROI; decreasing the observation time increases the accuracy. As NMR measurements do not track molecules leaving the ROI, this leads to higher agreement between simulations and measurements.

Fig. 6
figure 6

The effect of the observation time ∆ on the quality of agreement between the normalized displacement probabilities of the CFD-PT simulations (blue circles) and the PFG NMR measurements (orange squares) at a flow rate of 1.5 SL/min within a 10 PPI open-cell foam. ∆ = 5 ms (a) and ∆ = 10 ms (b)

Full size image

Even at the fastest velocities, the methane signal should be fully polarized, as shown by the following estimate of the residence time: The displacement curves (Fig. 4) show a maximum displacement of 3.5 mm at an observation time of Δ = 10 ms, corresponding to velocity v = 0.35 m/s and leading to the residence time tR = 0.04 m/(0.35 m/s) = 100 ms for fast-flowing components in the FOV of max. 40 mm. Thus, the time available for the spins to see the 7 Tesla and built-up magnetization is significantly longer than the T1 relaxation of 20–25 ms for methane at ambient temperature and normal pressure. In addition, a high magnetic field is also present outside the FOV, so the spins enter the 7T region in a pre-polarized state.

In the CFD-PT method, this is not an issue as particles are tracked regardless of their final position. To increase the accuracy of the comparison, the CFD-PT method was modified so that particles leaving the defined ROI are ignored. This shifts the displacement distribution to the left, but less than we expected (Fig. 7). While it increases the quality of agreement slightly it does not eliminate general deviation. During 10 ms of observation time, 16.7% of particles left the ROI. This was 7.85% in positive axial direction, 0.65% in the negative direction, and 8.16% in the radial direction. The considerable contribution of lateral flow causes only a small difference between the two CFD distributions.

Fig. 7
figure 7

The effect of considering all particles (black line) versus neglecting the particles which leave the ROI during 10 ms of observation time (blue line) at a flow rate of 1.5 SL/min using a 10 PPI OCF. The right panel shows the magnified view of the shaded region

Full size image

We further investigated whether a different size or position of the considered ROI in the CFD-PT simulation would influence the quality of agreement (Figs. S3 and S4). Changing the size of the ROI did not have a considerable influence on the quality of agreement. Only with the largest ROI (15 × 15 × 14 mm), the relative deviation of the propagators improved slightly in comparison with the original size which is 12 × 12 × 12 mm (Figure S4). The negligible effect of changing ROI size was also reported during experimental measurements (Mirdrikvand et al. 2020). Shifting the center of the ROI ± 1 mm in each direction did not show any improvement compared to the original ROI position (Figure S5). Finally, the displacement distribution of gas with a flow rate of 0.5 SL/min within an empty tube was compared with experimental PFG NMR measurements (Figure S6).

For the following evaluation, it is helpful to consider the different flow regimes investigated. They range from a Darcy regime (i.e., Rep < 1) with dominating viscous (diffusional) transport forces to a Darcy–Forchheimer regime with an increasing contribution of inertial (convective) forces (i.e., 1 < Rep < 10) and with predominant inertial forces (i.e., Rep > 10). For instance, the calculated Rep numbers for 10 and 20 PPI OCFs (Table 3) show that the flow regime is dominated by diffusion at a flow rate of 0.1 SL/min.

Table 3 Rep and Peeff values of different flow rates within 10 and 20 PPI OCFs for the three regimes, the superposition regime (diffusion), the transition to predominant mechanical dispersion (transition), and predominant mechanical dispersion (dispersion)
Full size table

The full width at half maximum (FWHM), a representative value also used for determining dispersion coefficients, is extracted from Figs. 4 and S3 distributions. It shows similar trends for PFG NMR measurements and CFD-PT simulations and increases with flow rate (Fig. 8, top row). The deviation of the FWHM obtained from both methods increases with flow rate (Fig. 8c). While the deviation of the FWHM obtained from the 20 PPI OCF shows only a slight increase from 0.19 to 0.20 mm when the flow rate rises from 0.1 to 1 SL/min, the deviation of the FWHM of 10 PPI OCF increases more strongly from 0.07 to 0.21 mm. This might be caused by a difference in flow velocity determined by both methods. The PFG NMR measurements usually show lower velocities than CFD calculations (not shown here) (Sadeghi et al. 2020). Either PFG NMR measurements underestimate convection contribution or the CFD workflow overpronounces it. At low flow rates (0.1 and 0.5 SL/min), the system operates in the transition range from the Darcy (Rep < 1) to the Darcy–Forchheimer flow regime (1 < Rep < 10). In this area, the flow is dominated by viscous forces rather than inertial forces. The 10 PPI OCF has a lower porosity and a smaller specific surface area than the 20 PPI sponge. Therefore, at 10 PPI and at low flow rates, the contribution of convection to mass transport is small, and thus, the deviation of the obtained FWHM between CFD and NMR is small (Fig. 6c). From 1 to 1.5 SL/min, the FWHM deviation plateaus. In this range of higher flow rates, the 20 PPI OCF shows lower deviation than the 10 PPI OCF. Here, the deviation can be explained by the higher averaged flow velocity within 10 PPI compared to 20 PPI. Thus, in the 10 PPI OCF, convection is more influential than in the 20 PPI OCF. This causes the higher deviations between simulation and NMR data at 10 PPI and high flow rates. This is particularly because NMR measurements underestimate convection, as mentioned above (Sadeghi et al. 2020).

Fig. 8
figure 8

Dispersion analysis: Full width at half maxima (FWHM) obtained from CFD-PT simulations and PFG NMR measurement at the different flow rates within a 10 PPI (a) and 20 PPI (b) open-cell foam, deviations between the FWHM obtained from simulation and measurement (c), and dependence of normalized dispersion coefficients Dax/Deff on Pe numbers for all cases (d). The observation time ∆ = 10 ms; axial dispersion coefficients (Dax) are calculated from Eq. (2)

Full size image

From the FWHM, it is possible to calculate the axial dispersion coefficient Dax (Mirdrikvand et al. 2020),

$$D_{{{\text{ax}}}} = \frac{{({{FWHM}})^{2} }}{16 \Delta \ln \left( 2 \right)},$$
(2)

where Δ is the observation time. Like the FWHM, the normalized dispersion coefficients (Dax/Deff) from CFD simulations are larger than those from PFG NMR measurements (Fig. 8d). When plotting Dax/Deff against Peeff, the curves do not collapse, and an offset between NMR and CFD data remains. All four curves, however, show a similar monotonic increase in Dax/Deff with Peeff. Therefore, this comparison clearly shows that both methods give the same dependence of axial dispersion on convection. This is, both convection and axial dispersion increase with flow rate but with an offset between each method and between the two OCFs. Nevertheless, at two high flow rates of 1 and 1.5 SL/min, both methods exhibit stronger convection in 10 PPI OCF. This results in higher axial dispersion coefficients Dax at the same velocity as 20 PPI OCF. This can be explained by a higher contribution of convection to dispersion in the case of the 10 PPI OCF. In addition, the two OCFs may also differ in the way they form convective dispersion due to their different pore or web geometries. For example, the ratio between the averaged window size and pore size is 0.6 and 0.7 for OCFs with 10 and 20 PPI, respectively. The smaller windows of the 10 PPI OCF resulted in a further increase in convection contribution compared to the profiles in the 20 PPI OCF.

The offset, on the other hand, might be a consequence of the PFG method used. This method shows increasing deviations with increasing fluid velocity, since it was developed for diffusion measurements and underestimates convection.

Fluctuations of FWHM values observed in the literature for 0.68 < Rep < 4.23 (Della Torre et al. 2014) might point toward measurement uncertainties, but there is also the option that the evaluation method used in the reference (fitting the non-normal distributed data to a Gaussian function) caused a certain scattering of the FWHM values.

In order to evaluate the influence of the flow rate on the discrepancy between CFD-PT and measured data, the differences between the integral values of the simulated and measured normalized displacement distributions were calculated from Figs. 4 and S2 (Fig. 9). It can be seen that in the case of 10 PPI, the difference between simulated and measured integral values increased significantly with increasing flow rate. However, in the case of 20 PPI, the difference remained almost constant. These slight changes are consistent with the trend in FWHM deviation observed for 10 and 20 PPI OCFs (see Fig. 8c). However, the integral differences remain below 0.04% and can be considered marginal.

Fig. 9
figure 9

Comparison of the difference between the obtained integral value of the normalized displacement distributions graphs from Figs. 4 and S3

Full size image

As stated in Table 2, the measured self-diffusion coefficient of methane (Dm = 2.24 × 10–5 m2/s) was implemented for performing CFD-PT simulations. This value was measured in a closed tube filled solely with methane by NMR (Mirdrikvand et al. 2018). In addition, to investigate also the effect of utilized diffusion coefficients, we performed CFD-PT simulations by implementing the diffusion coefficient of methane in air equal to Dm = 2.21 × 10–5 m2/s (Tang et al. 2015). For the 20 PPI OCF, a comparison revealed an increase in the deviation between CFD-PT and NMR propagators with a Dm rising from 2.21 × 10–5 to 2.24 × 10–5 m2/s, while the effect for 10 PPI was negligible (Figure S7). This can be explained by a more pronounced contribution of diffusion for the 20 PPI OCF as indicated by smaller Pe numbers (Table 3).

4 Summary and Conclusion

Coupling μCT-based pore-scale CFD simulations in OpenFOAM with the random-walk particle tracking solver disTrackFoam resulted in a numerical tool that enables microscopic analysis of gas flows within OCFs at much lower expense than experimental measurements or full-field molecular simulations such as Monte Carlo analysis. When using OCFs in CFD simulations identical to those employed in full-field PFG NMR measurements, we could obtain similar displacement distributions of gas flows at flow rates of 0.1, 0.5, 1, and 1.5 SL/min that show diffusive contributions to mass transport also in the opposite direction of the flow and allow for quantifying dispersion at low Pe numbers. However, match quality declined slightly with flow rate, which could be traced back to an underestimation of convection by PFG NMR measurements. By reducing the observation time from 10 to 5 ms, the discrepancy was almost eliminated. This could be explained by a reduced loss of magnetically labeled molecules that leave the region of interest before being detected. By further reducing this observation time, it should be possible to improve the matching and strive for even higher Péclet numbers.

Even though the integral difference of the displacement data could be shown to be smaller than 0.04%, it was sufficient to impact the full-width half maxima and dispersion coefficients obtained by CFD and PFG NMR. This demonstrates the sensitivity of these parameters to slight changes in the determined propagators.

Nevertheless, in the Darcy flow regime, which features a small contribution of inertia forces, the deviation between simulations and measurements is still low. At higher velocities, the discrepancy was more pronounced. In addition, the match was better for the 20 PPI OCF than for the 10 PPI OCF. This can be explained by the higher averaged mean velocities inside the 10 PPI OCF and thus an increased convection contribution.

Also, a lateral flow was determined by CFD-PT when particles were observed leaving the region of interest. A higher percentage of particles left in the radial direction than in the axial flow direction. Despite this considerable contribution of lateral flow, its impact on the deviations of CFD-PT simulations and PFG NMR measurements was small.

The cross-validation of CFD-PT and PFG NMR measurements in terms of local parameters such as displacement distribution of gas verified the capability of both complementary methods for microscopic analysis of gas flows within OCFs and demonstrated the advantage of the CFD-PT approach at larger flow velocities.

Similar to PFG NMR, where gradient pulses can be applied in all of the three directions so that both longitudinal and radial displacement distributions can be measured (Mirdrikvand et al. 2020), CFD-PT simulations have the potential to study displacement distributions in a directional manner as well. In this way, spatially specified dispersion coefficients in porous media such as open-cell foams can also be determined using this method.