Simulating flow and reactive transport in porous media

In this work, we developed and applied computational methods for simulating flow, transport and reactive transport in porous media. This comprised four main components: single-phase flow calculation; chemical transport calculation; the coupling to reaction kinetics at mineral surfaces and resulting structural changes; and the use of parallel and GPU computing to make the calculation practicable on realistic rock geometries. Single phase flow was calculated using the Lattice Boltzmann (LB) method. We used the multi-relaxation-time (MRT) operator for its superior stability and viscosity-independence. A sparse memory approach was employed which improves the efficiency of calculations performed on low-permeability rock pore-space images. We also extended an idea proposed by Skordos in which the lattice Boltzmann densities were transformed to increase the number of floating point bits retained in the calculation. We showed that this enhances the numerical precision of the calculation considerably, where the original paper found no appreciable benefit. We showed how this now permits the 4-byte datatype to be used reliably in slow flowing, heterogeneous domains. The LB algorithm was implemented for the use of parallel GPUs (Graphics Processors) using the MPI (Message Passing Interface) and shown to give strong scaling on a cluster of 24 Tesla K40 GPUs. A study of single phase permeability on micro-CT images of sandstone and carbonate rock pore structures of varying degrees of heterogeneity was carried out. Good agreement with experiment was found for the simpler pore spaces, while discrepancies in the micro-porous samples was attributed to two causes: 1) the exclusion of flow through unresolved micro-porosity and 2) unrepresentative sample sizes used in the simulation. The effect of image resolution and segmentation was studied by comparing single phase permeability computed in 1) scans of the same volume obtained at different voxel sizes, individually segmented and 2) numerically coarsened images from a high resolution segmented image. Numerical coarsening from a high resolution segmented image was found to be much more consistent than 1) and was shown to preserve porosity and permeability down to lower voxel size images unlike the images scanned and segmented at different voxel sizes. Finally, representative elementary volume (REV) was investigated for the rock samples. A statistical method was used in which porosity and permeability were obtained from sub-volumes sampled from the domain. The convergence of these parameters with sub-volume size was used to obtain characteristic length scales and measures of heterogeneity. The image sizes used were found to be unrepresentative for the complex microporous carbonates. Transport curves (propagators) were computed in three different porous media samples of increasing heterogeneity (a bead-pack; sandstone; and carbonate) and found to agree with experiment. Questions about the origins of stagnant transport zones in the microporous carbonate were pursued by investigating the effects of image segmentation. The effects of the image segmentation techniques, in which grey-scale micro-porosity in a scanned pore image is binarised into fluid or mineral, were quantified by computing the fraction of trapped solute (stagnant zones) for segmentations of varying porosity. Physical differences between experiment and calculation were clarified, and we suggest alternative approaches for the treatment of micro-porous rocks. A pore-scale reactive flow model was put together by coupling flow calculation and solute transport methods with changes in pore-structure through chemical kinetics. Convection and diffusion in this model was solved using a finite-volume approach: a second order transport model with a flux limiter function made the model suitable for high Peclet number transport calculations. We also proposed a method for counteracting errors associated with the staircase representation of diagonal surfaces in the Cartesian grid in which exposed grid surfaces are associated with a rescaling factor. First order reaction kinetics were included at mineral surfaces and the dissolution of a sphere was shown to give different dissolution profiles with different dimensionless transport and reaction parameters. The dissolution model was applied to the reaction between HCl acid and calcite mineral under the assumption that products of the reaction could be neglected. An experimental system in which HCl acid was injected through a flow cell containing a calcite block was simulated and the normalised volume of undissolved calcite was compared with the experimental data, as well as resulting morphologies obtained by micro-CT scanning. Good agreement with the experimental dissolution rate was obtained, however some differences in the resulting morphologies were found. This was attributed to neglecting the influence of product ions on the diffusion behaviour of the reactant and was discussed. By obtaining the concentration of H+ reactant on the surface of calcite block, the process could be concluded to be strongly transport-controlled. This enabled the definition of a new effective Damkohler number in terms of the reactant surface concentration which no longer required approximating length scales or separating convection or diffusion rates. Finally, the dissolution of a Ketton carbonate sample was computed. The injection process mirrored that of a strong acid flowing through the pore-space at a given flow rate, and having an intrinsic surface reaction rate with the rock mineral. It was found that the flow rate strongly affected the resulting dissolution pattern, in line with experimental observation. This lead to drastically altered flow properties, including single-phase permeability which was quantified.