Tracer dispersion in two-dimensional rough fractures

Prediction of the moments in advection-diffusion lattice Boltzmann method. I. Truncation dispersion, skewness, and kurtosis

The effect of the heterogeneity in the soil structure or the nonuniformity of the velocity field on the modeled resident time distribution (RTD) and breakthrough curves is quantified by their moments. While the first moment provides the effective velocity, the second moment is related to the longitudinal dispersion coefficient (k_{T}) in the developed Taylor regime; the third and fourth moments are characterized by their normalized values skewness (Sk) and kurtosis (Ku), respectively. The purpose of this investigation is to examine the role of the truncation corrections of the numerical scheme in k_{T}, Sk, and Ku because of their interference with the second moment, in the form of the numerical dispersion, and in the higher-order moments, by their definition. Our symbolic procedure is based on the recently proposed extended method of moments (EMM). Originally, the EMM restores any-order physical moments of the RTD or averaged distributions assuming that the solute concentration obeys the advection-diffusion equation in multidimensional steady-state velocity field, in streamwise-periodic heterogeneous structure. In our work, the EMM is generalized to the fourth-order-accurate apparent mass-conservation equation in two- and three-dimensional duct flows. The method looks for the solution of the transport equation as the product of a long harmonic wave and a spatially periodic oscillating component; the moments of the given numerical scheme are derived from a chain of the steady-state fourth-order equations at a single cell. This mathematical technique is exemplified for the truncation terms of the two-relaxation-time lattice Boltzmann scheme, using plug and parabolic flow in straight channel and cylindrical capillary with the d2Q9 and d3Q15 discrete velocity sets as simple but illustrative examples. The derived symbolic dependencies can be readily extended for advection by another, Newtonian or non-Newtonian, flow profile in any-shape open-tabular conduits. It is established that the truncation errors in the three transport coefficients k_{T}, Sk, and Ku decay with the second-order accuracy. While the physical values of the three transport coefficients are set by Péclet number, their truncation corrections additionally depend on the two adjustable relaxation rates and the two adjustable equilibrium weight families which independently determine the convective and diffusion discretization stencils. We identify flow- and dimension-independent optimal strategies for adjustable parameters and confront them to stability requirements. Through specific choices of two relaxation rates and weights, we expect our results be directly applicable to forward-time central differences and leap-frog central-convective Du Fort-Frankel-diffusion schemes. In straight channel, a quasi-exact validation of the truncation predictions through the numerical moments becomes possible thanks to the specular-forward no-flux boundary rule. In the staircase description of a cylindrical capillary, we account for the spurious boundary-layer diffusion and dispersion because of the tangential constraint of the bounce-back no-flux boundary rule.

Prediction of the moments in advection-diffusion lattice Boltzmann method. II. Attenuation of the boundary layers via double-Λ bounce-back flux scheme

Impact of the unphysical tangential advective-diffusion constraint of the bounce-back (BB) reflection on the impermeable solid surface is examined for the first four moments of concentration. Despite the number of recent improvements for the Neumann condition in the lattice Boltzmann method-advection-diffusion equation, the BB rule remains the only known local mass-conserving no-flux condition suitable for staircase porous geometry. We examine the closure relation of the BB rule in straight channel and cylindrical capillary analytically, and show that it excites the Knudsen-type boundary layers in the nonequilibrium solution for full-weight equilibrium stencil. Although the d2Q5 and d3Q7 coordinate schemes are sufficient for the modeling of isotropic diffusion, the full-weight stencils are appealing for their advanced stability, isotropy, anisotropy and anti-numerical-diffusion ability. The boundary layers are not covered by the Chapman-Enskog expansion around the expected equilibrium, but they accommodate the Chapman-Enskog expansion in the bulk with the closure relation of the bounce-back rule. We show that the induced boundary layers introduce first-order errors in two primary transport properties, namely, mean velocity (first moment) and molecular diffusion coefficient (second moment). As a side effect, the Taylor-dispersion coefficient (second moment), skewness (third moment), and kurtosis (fourth moment) deviate from their physical values and predictions of the fourth-order Chapman-Enskog analysis, even though the kurtosis error in pure diffusion does not depend on grid resolution. In two- and three-dimensional grid-aligned channels and open-tubular conduits, the errors of velocity and diffusion are proportional to the diagonal weight values of the corresponding equilibrium terms. The d2Q5 and d3Q7 schemes do not suffer from this deficiency in grid-aligned geometries but they cannot avoid it if the boundaries are not parallel to the coordinate lines. In order to vanish or attenuate the disparity of the modeled transport coefficients with the equilibrium weights without any modification of the BB rule, we propose to use the two-relaxation-times collision operator with free-tunable product of two eigenfunctions Λ. Two different values Λ_{v} and Λ_{b} are assigned for bulk and boundary nodes, respectively. The rationale behind this is that Λ_{v} is adjustable for stability, accuracy, or other purposes, while the corresponding Λ_{b}(Λ_{v}) controls the primary accommodation effects. Two distinguished but similar functional relations Λ_{b}(Λ_{v}) are constructed analytically: they preserve advection velocity in parabolic profile, exactly in the two-dimensional channel and very accurately in a three-dimensional cylindrical capillary. For any velocity-weight stencil, the (local) double-Λ BB scheme produces quasi-identical solutions with the (nonlocal) specular-forward reflection for first four moments in a channel. In a capillary, this strategy allows for the accurate modeling of the Taylor-dispersion and non-Gaussian effects. As illustrative example, it is shown that in the flow around a circular obstacle, the double-Λ scheme may also vanish the dependency of mean velocity on the velocity weight; the required value for Λ_{b}(Λ_{v}) can be identified in a few bisection iterations in given geometry. A positive solution for Λ_{b}(Λ_{v}) may not exist in pure diffusion, but a sufficiently small value of Λ_{b} significantly reduces the disparity in diffusion coefficient with the mass weight in ducts and in the presence of rectangular obstacles. Although Λ_{b} also controls the effective position of straight or curved boundaries, the double-Λ scheme deals with the lower-order effects. Its idea and construction may help understanding and amelioration of the anomalous, zero- and first-order behavior of the macroscopic solution in the presence of the bulk and boundary or interface discontinuities, commonly found in multiphase flow and heterogeneous transport.

GPU accelerated Monte Carlo simulation of pulsed-field gradient NMR experiments

The simulation of diffusion by Monte Carlo methods is often essential to describing NMR measurements of diffusion in porous media. However, simulation timescales must often span hundreds of milliseconds, with large numbers of trajectories required to ensure statistical convergence. Here we demonstrate that by parallelising code to run on graphics processing units (GPUs), these calculations may be accelerated by over three orders of magnitude, opening new frontiers in experimental design and analysis. As such cards are commonly installed on most desktop computers, we expect that this will prove useful in many cases where simple analytical descriptions are not available or appropriate, e.g. in complex geometries or where short gradient pulse approximations do not hold, or for the analysis of diffusion-weighted MRI in complex tissues such as the lungs and brain.

Three-dimensional modeling of the brain's ECS by minimum configurational energy packing of fluid vesicles

The extracellular space of the brain is the heterogeneous porous medium formed by the spaces between the brain cells. Diffusion in this interstitial space is the mechanism by which glucose and oxygen are delivered to the brain cells from the vascular system. It is also a medium for the transport of certain informational substances between the cells (called volume transmission), and for drug delivery. This work involves three-dimensional modeling of the extracellular space as void space in close-packed arrays of fluid membrane vesicles. These packings are generated by minimizing the configurational energy using a Monte Carlo procedure. Both regular and random packs of vesicles are considered. A random walk algorithm is then used to compute the geometric tortuosities, and the results are compared with published experimental data. For the random packings, it is found that although the absolute values for the tortuosities differ, the dependence of the tortuosity on pore volume fraction is very similar to that observed in experiment. The tortuosities we measure are larger than those computed in previous studies of packings of convex polytopes, and modeling improvements, which require higher resolution studies and an improved modeling of brain cell shapes and mechanical properties, could help resolve remaining discrepancies between model simulations and experiment. It is also shown that the specular reflection scheme is the appropriate technique for implementing zero-flux boundary conditions in random walk simulations commonly encountered in diffusion problems.

Lattice Boltzmann method for non-Newtonian (power-law) fluids

We study an ad hoc extension of the lattice Boltzmann method that allows the simulation of non-Newtonian fluids described by generalized Newtonian models. We extensively test the accuracy of the method for the case of shear-thinning and shear-thickening truncated power-law fluids in the parallel plate geometry, and show that the relative error compared to analytical solutions decays approximately linear with the lattice resolution. Finally, we also tested the method in the reentrant-flow geometry, in which the shear rate is no longer a scalar and the presence of two singular points requires high accuracy in order to obtain satisfactory resolution in the local stress near these points. In this geometry, we also found excellent agreement with the solutions obtained by standard finite-element methods, and the agreement improves with higher lattice resolution.

Self-affine fronts in self-affine fractures: large and small-scale structure

The evolution and spatial structure of displacement fronts in fractures with self-affine rough walls are studied by numerical simulations. The fractures are open and the two faces are identical but shifted along their mean plane, either parallel or perpendicular to the flow. An initially flat front advected by the flow is progressively distorted into a self-affine front with a Hurst exponent equal to that of the fracture walls. The lower cutoff of the self-affine regime depends only on the aperture, while the upper cutoff grows with the lateral shift and linearly with the width of the front.

Boundary conditions for stochastic solutions of the convection-diffusion equation

Stochastic methods offer an attractively simple solution to complex transport-controlled problems, and have a wide range of physical, chemical, and biological applications. Stochastic methods do not suffer from the numerical diffusion that plagues grid-based methods, but they typically lose accuracy in the vicinity of interfacial boundaries. In this work we introduce some ideas and algorithms that can be used to implement boundary conditions in stochastic simulations of the convection-diffusion equation with accuracies comparable to the bulk phase. The algorithms have been tested in two-dimensional channel flows over a range of Peclet numbers, and compared with independent finite-difference calculations.

Transport in rough self-affine fractures

Transport properties of three-dimensional self-affine rough fractures are studied by means of an effective-medium analysis and numerical simulations using the Lattice-Boltzmann method. The numerical results show that the effective-medium approximation predicts the right scaling behavior of the permeability and of the velocity fluctuations, in terms of the aperture of the fracture, the roughness exponent, and the characteristic length of the fracture surfaces, in the limit of small separation between surfaces. The permeability of the fractures is also investigated as a function of the normal and lateral relative displacements between surfaces, and it is shown that it can be bounded by the permeability of two-dimensional fractures. The development of channel-like structures in the velocity field is also numerically investigated for different relative displacements between surfaces. Finally, the dispersion of tracer particles in the velocity field of the fractures is investigated by analytic and numerical methods. The asymptotic dominant role of the geometric dispersion, due to velocity fluctuations and their spatial correlations, is shown in the limit of very small separation between fracture surfaces.