Search

Field-scale observations suggest that rock heterogeneities control subsurface fluid flow, and these must be characterised for accurate predictions of fluid migration, such as during \CO2 sequestration. Recent efforts have focused on simulation-based inversion of laboratory observations with X-ray imaging, but models produced in this way have been limited in their predictive ability for heterogeneous rocks. We address the main challenges in this approach through an algorithm that combines: a 3-parameter capillary pressure model, spatial heterogeneity in absolute permeability, the constraint of history match iterations based on marginal error improvement, and image processsing that incorporates more of the experimental data in the calibration. We demonstrate the improvements on five rocks (two sandstones and three carbonates), representing a range of heterogeneous properties, some of which could not be previously modelled. The algorithm results in physically representative models of the rock cores, reducing non-systematic error to a level comparable to the experimental uncertainty.

Abstract Flow of complex fluids in porous structures is pertinent in many biological and industrial processes. For these applications, elastic turbulence, a viscoelastic instability occurring at low Re—arising from a non-trivial coupling of fluid rheology and flow geometry—is a common and relevant effect because of significant over-proportional increase in pressure drop and spatio-temporal distortion of the flow field. Therefore, significant efforts have been made to predict the onset of elastic turbulence in flow geometries with constrictions. The onset of flow perturbations to fluid streamlines is not adequately captured by Deborah and Weissenberg numbers. The introduction of more complex dimensionless numbers such as the M-criterion, which was meant as a simple and pragmatic method to predict the onset of elastic instabilities as an order-of-magnitude estimate, has been successful for simpler geometries. However, for more complex geometries which are encountered in many relevant applications, sometimes discrepancies between experimental observation and M-criteria prediction have been encountered. So far these discrepancies have been mainly attributed to the emergence from disorder. In this experimental study, we employ a single channel with multiple constrictions at varying distance and aspect ratios. We show that adjacent constrictions can interact via non-laminar flow field instabilities caused by a combination of individual geometry and viscoelastic rheology depending (besides other factors) explicitly on the distance between adjacent constrictions. This provides intuitive insight on a more conceptual level why the M-criteria predictions are not more precise. Our findings suggest that coupling of rheological effects and fluid geometry is more complex and implicit and controlled by more length scales than are currently employed. For translating bulk fluid, rheology determined by classical rheometry into the effective behaviour in complex porous geometries requires consideration of more than only one repeat element. Our findings open the path towards more accurate prediction of the onset of elastic turbulence, which many applications will benefit. Article Highlights We demonstrate that adjacent constrictions “interact” via the non-laminar flow fields caused by individual constrictions, implying that the coupling of rheological effects and fluid geometry is more complex and implicit. The concept of characterizing fluid rheology independent of flow geometry and later coupling back to the geometry of interest via dimensionless numbers may fall short of relevant length scales, such as the separation of constrictions which control the overlap of flow fields. By providing direct experimental evidence illustrating the cause of the shortcoming of the status-quo, the expected impact of this work is to challenge and augment existing concepts that will ultimately lead to the correct prediction of the onset of elastic turbulence.

Wettability is one of the key controlling parameters for multiphase flow in porous media, and paramount for various geoscience applications. While a general awareness of the importance of wettability was established decades ago, our fundamental understanding of how wettability influences transport and of how to characterize wettability has improved tremendously in recent years through breakthroughs in imaging technology and modeling techniques. Numerical modeling studies clearly show not only that macroscopic two-phase flow is influenced by the average wettability, but also that the spatial distribution of wetting significantly impacts the macroscopic parameters. Herein, we explore the thermodynamics for porous multiphase systems, and recent breakthroughs in wettability characterization. Our view is that bridging the multiscale characterization of wetting must consider two fundamental perspectives: geometry and energy. Advancing the overall description requires an improved understanding of the operative mechanisms that dominate at various scales, and the development of quantitative approaches to capture these effects. We take a multistage approach, looking at these fundamental perspectives from the sub-pore-to-pore length scales, followed by the pore-to-core length scales using various analytical techniques and numerical simulations. Within this context, there remain many open-ended questions, and we therefore highlight these issues to provide guidance on future research directions. Our overall aim is to provide comprehensive guidance on the multiscale characterization of wettability in porous media, in order to facilitate novel research.

Many properties of complex porous media such as reservoir rocks are strongly affected by heterogeneity at different scales. Complex depositional and diagenetic processes have a strong control on the pore structures, leading to systems with a wide range of pore sizes covering many orders of magnitude in length scales. This poses a significant challenge for digital rock analysis since a single resolution image and associated simulation model cannot capture all the relevant length scales in sufficient detail due to limitations in computer memory and speed. The scale-transgressive effects of heterogeneity must therefore be accounted for through a multiscale digital rock workflow. We introduce a generalized multiscale imaging and pore-scale modelling workflow to derive transport properties of complex rocks having broad pore size distributions. A dry/wet micro-CT imaging sequence is used to spatially map the porosity and the connectivity of resolved and unresolved porous regions. The unresolved porosity regions are classified into different porosity classes or rock types. The resulting 3D rock-type map and the porosity map are combined and transformed into a multiscale pore network model. Resolved pores are treated in a conventional pore network manner while unresolved network elements are treated as a continuum Darcy-type porous medium. Similar to conventional continuum models, each Darcy pore is populated with single and multiphase flow properties. These properties are derived from high-resolution rock-type models constructed from backscatter SEM images and/or high-resolution micro-CT images of sub-samples. The multiscale digital rock workflow is applied to two heterogeneous rock samples: a mixed wet thinly laminated reservoir sandstone and an oil wet reservoir carbonate. Experimentally measured mercury–air primary drainage and oil–water imbibition capillary pressure curves (after ageing to restore wettability) are used to anchor the multiscale pore network model. Waterflood relative permeability is calculated in a blind test and compared with high-quality experimental data. A very encouraging agreement between computed and measured properties is found.

The long-standing question on the adequate description of multiphase flow in porous media may be ultimately decided based on the ability to estimate model parameters with sufficient accuracy that make the models distinguishable. Since the most-common Darcy scale multiphase flow models all use a somewhat phenomenological relative permeability or resistance factor formulation, the key question is what the associated uncertainty really is when derived from flow experiments by inverse modeling. In this work, a recently developed workflow for systematic assessment of uncertainty was used to analyze the impact of the choice of relative permeability models and associated uncertainty. In an exemplary fashion, the Corey and LET relative permeability parameterizations were compared. The choice of Corey and LET models in the inverse modeling workflows showed differences in the derived relative permeability relations. The Corey parameterization is found to be more restrictive and imposed additional constraints on parameters. For example, varying the connate water saturation and residual oil saturation did not improve the match with experimental data. The pressure drop, saturation profiles and the capillary pressure–saturation relationship constrained the solution and imposing additional constraints on, e.g., residual oil saturation has very little impact on the result. In contrast, the LET function provided more degrees of freedom in order to accommodate the shape of the relative permeability curves. The findings also suggested that both Corey and LET models may not necessarily provide optimum parameterizations of the experimental data. The cross-correlations of fit parameters and non-Gaussian residuals indicated that we were still dealing with a phenomenological parameterization that is not yet the fully adequate description of the data. This may be the starting point for a comparison of different flow models beyond the uncertainty imposed by the choice of model parameterizations. Future work is aimed at assessing whether better choices in interpretation workflows and optimized experimental workflows can minimize these issues.

Multiphase flow in porous media is strongly influenced by the pore-scale arrangement of fluids. Reservoir-scale constitutive relationships capture these effects in a phenomenological way, relying only on fluid saturation to characterize the macroscopic behavior. Working toward a more rigorous framework, we make use of the fact that the momentary state of such a system is uniquely characterized by the geometry of the pore-scale fluid distribution. We consider how fluids evolve as they undergo topological changes induced by pore-scale displacement events. Changes to the topology of an object are fundamentally discrete events. We describe how discontinuities arise, characterize the possible topological transformations and analyze the associated source terms based on geometric evolution equations. Geometric evolution is shown to be hierarchical in nature, with a topological source term that constrains how a structure can evolve with time. The challenge associated with predicting topological changes is addressed by constructing a universal geometric state function that predicts the possible states based on a non-dimensional relationship with two degrees of freedom. The approach is validated using fluid configurations from both capillary and viscous regimes in ten different porous media with porosity between 0.10 and 0.38. We show that the non-dimensional relationship is independent of both the material type and flow regime. We demonstrate that the state function can be used to predict history-dependent behavior associated with the evolution of the Euler characteristic during two-fluid flow.

Hypothesis: The wetting behaviour is a key property of a porous medium that controls hydraulic conductivity in multiphase flow. While many porous materials, such as hydrocarbon reservoir rocks, are initially wetted by the aqueous phase, surface active components within the non-wetting phase can alter the wetting state of the solid. Close to the saturation endpoints wetting phase fluid films of nanometre thickness impact the wetting alteration process. The properties of these films depend on the chemical characteristics of the system. Here we demonstrate that surface texture can be equally important and introduce a novel workflow to characterize the wetting state of a porous medium. Experiments: We investigated the formation of fluid films along a rock surface imaged with atomic force microscopy using zeta-potential measurements and a computational model for drainage. The results were compared to spontaneous imbibition test to link sub-pore-scale and core-scale wetting characteristics of the rock. Findings: The results show a dependency between surface coverage by oil, which controls the wetting alteration, and the macroscopic wetting response. The surface-area coverage is dependent on the capillary pressure applied during primary drainage. Close to the saturation endpoint, where the change in saturation was minor, the oil-solid contact changed more than 80%. (C) 2019 The Authors. Published by Elsevier Inc.

Understanding the interaction between competing fluids in the pore space of rocks is key for predicting subsurface flow and trapping, such as with CO2 in a saline aquifer. These processes occur over a large span of timescales (from seconds to thousands of years), and length scales (from microns to kilometres). Understanding the link between these temporal and spatial scales will enable us to interpolate between observations made at different resolutions. In this work we explore the temporal scales present during macroscopically steady-state multiphase flow in a porous rock using differential pressure measurements. We observe a cascade of timescales in the pressure differential i.e. a continuous range of frequencies, with lower frequencies having greater amplitudes. We demonstrate a scaling of the spectral density with frequency of S ∼ 1/f^2, or red noise, to describe the dynamics. This scaling is independent of the flow rate of the fluids or the fraction of the flow taken by water. This red, or Brownian, noise indicates a stochastic process where pressure fluctuations are seen throughout the pore space, resulting in intermittent filling of pores over a wide range of time-scales, from seconds to minutes in these experiments. This observation will aid future modelling of subsurface flow as it suggests self-organised critically of the system with no characteristic time or length scale.

Multiphase flow in porous media systems is a critical element of many processes in the energy industry. The characteristics of the simultaneous flow of the immiscible phases can be quantified using relative permeability relations. In geoscience applications, these relations are determined in coreflooding studies that often comprise coreflood tests of oil–water mixtures performed on centimetre-scale rock samples. The outcomes of these are subject to uncertainty, which ultimately influences how accurately the parameters from small-scale tests translate to the upscaled estimations. To assess this uncertainty, Shell researchers have developed an inverse modelling workflow for the uncertainty analysis of relative permeability functions derived from coreflood tests. The results suggest that, even at a small scale, the uncertainty can be significant.

HypothesisHigh molecular weight polymers are widely used in oilfield applications, such as in chemical enhanced oil recovery (cEOR) technique for hydrocarbon recovery. However, during flow in a porous rock, polymer retention is usually a major challenge, as it may result in the decrease of polymer concentration or lead to plugging of pores with significant permeability reduction and injectivity loss. Hence, an understanding of the retention mechanisms will have a profound effect in optimizing the process of polymer flooding, in particular, for carbonate rocks, which hold more than half of the world's oil reserves. Therefore, in this study, the retention of hydrolysed polyacrylamide (HPAM) polymer, a commonly used chemical for EOR, is investigated during flow in Estaillades carbonate rock.ExperimentsA novel approach of investigating HPAM retention in Estaillades carbonate rock was carried out using Atomic force microscopy (AFM). Since Estaillades carbonate rock is ∼98% calcite, HPAM retention was first characterised on a cleaved flat calcite mineral surface after immersing in HPAM solution. Afterwards, HPAM was flooded in Estaillades carbonate to observe the effect of flow dynamics on the retention mechanisms.FindingsWe find that the dominant mechanism for retention of HPAM on calcite after fluid immersion is polymer adsorption, which we believe is driven by the electrostatic interaction between the calcite surface and the solution. The thickness of the adsorbed layer on calcite is beyond 3 ​nm suggesting it is not adsorbed only flat on the surface. Different types of adsorbed layers were formed representing trains, and the more extended loops or tails with the largest polymer layer thickness about 35 ​nm, representing the longer loops or tails. Layers of this thickness will begin to impair the permeability of the rock. However, in Estaillades, thicker adsorbed layers are observed in different regions of the rock surface ranging between 50 and 350 ​nm. We suggest that this is due to either mechanical entrapment and/or polymer entanglement during flow in Estaillades carbonate rock, which will cause the major permeability impairment in porous rocks.

Hypothesis: The macroscopic movement of subsurface fluids involved in CO2 storage, groundwater, and petroleum engineering applications is controlled by interfacial forces in the pores of rocks. Recent advances in modelling these systems has arisen from approaches simulating flow through a digital representation of the complex pore structure. However, further progress is limited by difficulties in characterising the spatial distribution of the wetting state within the pore structure. In this work, we show how observations of the fluid coverage of mineral surfaces within the pores of rocks can be used as the basis for a quantitative 3D characterisation of heterogeneous wetting states throughout rock pore structures. Experiments: We demonstrate the approach with water–oil fluid pairs on rocks with distinct lithologies (sandstone and carbonate) and wetting states (hydrophilic, intermediate wetting, and heterogeneously wetting). Findings: Fluid surface coverage the within rock pores is a robust signal of the wetting state across varying rock types and wetting states. The wetting state can be quantified and the resulting 3D maps can be used as a deterministic input to pore scale models. These may be applied to multiphase flow problems in porous media ranging from soil science to fuel cells.

Hypothesis Emerging energy-related technologies deal with multiscale hierarchical structures, intricate surface morphology, non-axisymmetric interfaces, and complex contact lines where wetting is difficult to quantify with classical methods. We hypothesise that a universal description of wetting on multiscale surfaces can be developed by using integral geometry coupled to thermodynamic laws. The proposed approach separates the different hierarchy levels of physical description from the thermodynamic description, allowing for a universal description of wetting on multiscale surfaces. Theory and Simulations The theoretical framework is presented followed by application to limiting cases of wetting on multiscale surfaces. Limiting cases include those considered in the Wenzel, Cassie-Baxter, and wicking state models. Wetting characterisation of multiscale surfaces is explored by conducting simulations of a fluid droplet on a structurally rough surface and a chemically heterogeneous surface. Findings The underlying origin of the classical wetting models is shown to be rooted within the proposed theoretical framework. Integral geometry provides a topological-based wetting metric that is not contingent on any type of wetting state. The wetting metric is demonstrated to account for multiscale features along the common line in a scale consistent way; providing a universal description of wetting for multiscale surfaces.

Relative permeability is commonly used to model immiscible fluid flow through porous materials. In this work, we derive the relative permeability relationship from conservation of energy, assuming that the system to be non-ergodic at large length scales and relying on averaging in both space and time to homogenize the behavior. Explicit criteria are obtained to define stationary conditions: (1) there can be no net change for extensive measures of the system state over the time averaging interval; (2) the net energy inputs into the system are zero, meaning that the net rate of work done on the system must balance with the heat removed; and (3) there is no net work performed due to the contribution of internal energy fluctuations. Results are then evaluated based on direct numerical simulation. Dynamic connectivity is observed during steady-state flow, which is quantitatively assessed based the Euler characteristic. We show that even during steady-state flow at low capillary number ([Formula: see text]), typical flow processes will explore multiple connectivity states. The residence time for each connectivity state is captured based on the time-and-space average. The distribution for energy fluctuations is shown to be multi-modal and non-Gaussian when terms are considered independently. However, we demonstrate that their sum is zero. Given an appropriate choice of the thermodynamic driving force, we show that the conventional relative permeability relationship is sufficient to model the energy dissipation in systems with complex pore-scale dynamics that routinely alter the structure of fluid connected pathways.

We describe the underlying mathematics, validation, and applications of a novel Helmholtz free-energy—minimizing phase-field model solved within the framework of the lattice Boltzmann method (LBM) for efficiently simulating two-phase pore-scale flow directly on large 3D images of real rocks obtained from micro-computed tomography (micro-CT) scanning. The code implementation of the technique, coined as the eLBM (energy-based LBM), is performed in CUDA programming language to take maximum advantage of accelerated computing by use of multinode general-purpose graphics processing units (GPGPUs). eLBM’s momentum-balance solver is based on the multiple-relaxation-time (MRT) model. The Boltzmann equation is discretized in space, velocity (momentum), and time coordinates using a 3D 19-velocity grid (D3Q19 scheme), which provides the best compromise between accuracy and computational efficiency. The benefits of the MRT model over the conventional single-relaxation-time Bhatnagar-Gross-Krook (BGK) model are (I) enhanced numerical stability, (II) independent bulk and shear viscosities, and (III) viscosity-independent, nonslip boundary conditions. The drawback of the MRT model is that it is slightly more computationally demanding compared to the BGK model. This minor hurdle is easily overcome through a GPGPU implementation of the MRT model for eLBM. eLBM is, to our knowledge, the first industrial grade–distributed parallel implementation of an energy-based LBM taking advantage of multiple GPGPU nodes. The Cahn-Hilliard equation that governs the order-parameter distribution is fully integrated into the LBM framework that accelerates the pore-scale simulation on real systems significantly. While individual components of the eLBM simulator can be separately found in various references, our novel contributions are (1) integrating all computational and high-performance computing components together into a unified implementation and (2) providing comprehensive and definitive quantitative validation results with eLBM in terms of robustness and accuracy for a variety of flow domains including various types of real rock images. We successfully validate and apply the eLBM on several transient two-phase flow problems of gradually increasing complexity. Investigated problems include the following: (1) snap-off in constricted capillary tubes; (2) Haines jumps on a micromodel (during drainage), Ketton limestone image, and Fontainebleau and Castlegate sandstone images (during drainage and subsequent imbibition); and (3) capillary desaturation simulations on a Berea sandstone image including a comparison of numerically computed residual non-wetting-phase saturations (as a function of the capillary number) to data reported in the literature. Extensive physical validation tests and applications on large 3D rock images demonstrate the reliability, robustness, and efficacy of the eLBM as a direct visco-capillary pore-scale two-phase flow simulator for digital rock physics workflows.

Pore-to-Darcy scale upscaling of multiphase flow is one of the major unresolved problems in many fields of porous media research. While this problem involves very fundamental aspects, there are many practical and application-driven challenges as well, such as the accurate prediction of Darcy-scale multiphase effective properties, e.g., relative permeability by pore-scale flow simulation on the basis of the imaged pore geometry, e.g., via X-ray computed micro-tomography. Validation of pore-scale modeling methods against experimental data by comparison of measured against simulated relative permeability curves has proven to be insufficient. Comparison of the fluid topology, in particular, the non-wetting phase topology, is a much more reliable criteria as relative permeability shows a very strong correlation with connectivity. While percolation-based quasi-static modeling approaches operating in the capillary limit have proven moderately successful in drainage, they largely fail to predict fluid connectivity in imbibition. We show that a fully coupled visco-capillary simulation, using a free-energy based lattice Boltzmann approach, correctly predicts fluid connectivity in imbibition, which was not predicted correctly by a quasi-static approach in a previous study on the same dataset. The respective relative permeability data shows a close match with results from Darcy-scale core flooding experiments while there is a major mismatch from quasi-static approach. The results once more highlight the close connection between relative permeability and fluid topology and suggest that topological measures are much more meaningful validation criteria for pore-scale simulation. The results also show that the direct simulation approach described in this paper, using appropriate initial and boundary conditions, correctly predicts Darcy-scale effective properties.

High-resolution X-ray imaging combined with a steady-state flow experiment is used to demonstrate how pore-scale displacement affects macroscopic properties in an altered-wettability microporous carbonate, where porosity and fluid saturation can be directly obtained from the grey-scale micro-CT images. The resolvable macro pores are largely oil-wet with an average thermodynamic contact angle of 120°. The pore-by-pore analysis shows locally either oil or brine almost fully occupied the macro pores, with some oil displacement in the micro-porosity. We observed a typical oil-wet behaviour consistent with the contact angle measurement. The brine tended to occupy the larger macro pores, leading to a higher brine relative permeability, lower residual oil saturation, than under water-wet conditions and in a mixed-wet sandstone. The capillary pressure was negative and seven times larger in the carbonate than the sandstone, despite having a similar average pore size. These different displacement patterns are principally determined by the difference in wettability.

Summary A key step in de-risking chemical enhanced oil recovery (cEOR) projects is to assess the incremental recovery for the field of interest in customised laboratory experiments that mimic conditions within target reservoirs. Any deviation from these conditions, as is oftentimes the case, leads to discrepancies which call the reliability of laboratory results into question, thereby affecting the economics of the cEOR projects. Concerning iron-bearing formations, one approach is to treat samples with a reducing fluid in order to mimic native reservoir redox conditions. In this study, investigations into the effect of a solution of the reducing agent, Sodium Dithionite, in brine on surface wettability were performed using Atomic Force Microscopy (AFM) to quantify interactions between model crude oil components and an iron-bearing sandstone under varying redox conditions. Results show that the adhesion of the oil components to the sandstone surface decreased in the order -NH2 (∼70%) > -COOH (∼36%) > -CH3 (∼3%) on introduction of the reducing fluid, potentially providing a basis for deployment in core floods to ascertain the suitability of cEOR procedures.

When oil fields fall during their lifetime below the bubble point gas comes out of solution. The key questions are at which saturation the gas becomes mobile (“critical gas saturation”) and what the gas mobility is, because mobile gas reduces the production of oil significantly. The traditional view is that the gas phase becomes mobile once gas bubbles grow or expand to a size where they connect and form a percolating path. For typical 3D porous media the saturation corresponding to this percolation limit is on the order of 20%. However, significant literature report on gas mobility below lower limits of percolation thresholdsi.e.below 0.1%. A direct experimental insight for that is lacking because laboratory measurements are notoriously difficult since the formation of gas bubbles below the bubble point includes thermodynamic and kinetic aspects, and the pressure decline rates achievable in laboratory experiments are orders of magnitude higher than the decline rates in the field. Here we study the nucleation and transport of gas coming out of solutionin-situin 3D rock using X-ray computed micro tomography which allows direct visualization of the nucleation kinetics and connectivity of gas. We use either propane or a propane–decane mixture as model system and conduct pressure depletion in absence of flow finding that – consistent with the literature – observation of the bubble point in the porous medium is decreased and becomes pressure decline rate dependent because of the bubble nucleation kinetics. That occurs in single-component systems and in hydrocarbon mixtures. Pressure depletion in absence of flow results in critical gas saturations between 20 and 30% which is consistent with typical percolation thresholds in 3D porous structures. That does not explain experimentally observed critical gas saturations significantly below 20%. Also, the respective pore level fluid occupancy where pores are filled with either gas or liquid phase but not partially with both as in normal 2-phase immiscible systems rather diminishes connectivity of gas and liquid phases. This observation indicates that likely other mechanisms play a role in establishing gas mobility at saturations significantly below 20%. Experiments under flow conditions, where gas is injected near the bubble point suggest that diffusion may significantly contribute to the transport of gas and may even be the dominant transport mechanism at field relevant flow rates. The consequence of diffusive transport are compositional gradients where locally the composition is such gas nucleation may occur. That would lead to a disconnected but mobile gas distribution ahead of the convective front. Furthermore, diffusive exchange leads to ripening and anti-ripening effects which influences the distribution for which we see evidence in pressure depletion experiments but not so much at low rate gas injection. Respective relative permeability computed from the imaged fluid distributions using a lattice Boltzmann approach show distinctly different behavior between pressure depletion and flowing conditions. These findings suggest that capillarity in a gas–liquid hydrocarbon mixture is far more complex than in a 2-phase immiscible system. Capillarity is coupled to phase behavior thermodynamics and kinetics on a fast time scale and diffusion-dominated mechanisms such as ripening and anti-ripening effects at a slow time scale. While the consequences for the current experimental and field modelling approaches are not yet fully clear, this shows that more research is needed to fully understand these effects and their implications.

Shell scientists are making the most of advancing imaging technology to reveal what happens in a 3D porous medium during emulsification.

At the pore scale, slow invasion of a wetting fluid in porous materials is often modelled with quasi-static approximations which only consider capillary forces in the form of simple pore filling rules. The appropriateness of this approximation, often applied in pore network models, is contested in literature, reflecting the difficulty of predicting imbibition relative permeability with these models. However, validation by sole comparison to continuum-scale experiments is prone to induce model overfitting. It has therefore remained unclear whether difficulties generalizing the model performance are due to network extraction, the pore filling rules, or whether a quasi-static description is useful at all. Here, we address this by examining whether such a model can predict the pore-scale fluid distributions underlying the behaviour at the continuum scale. To this end, we compare the fluid arrangement evolution measured in fast synchrotron micro-CT experiments on two rock types to quasi-static simulations which implement capillary-dominated pore filling and snap-off, including a sophisticated model for cooperative pore filling. The results indicate that pore network models with appropriate pore filling rules can in principle obtain a good first-order prediction of the upscaled flow properties of strongly-wetted rocks at low capillary numbers.

The wetting state is an important control on flow in subsurface multi fluid phase systems, e.g., carbon storage and oil production. Advances in X-ray imaging allow us to characterise the wetting state using imagery of fluid arrangement within the pores of rocks. We derived a model from equilibrium thermodynamics relating fluid coverage of rock surfaces to wettability and fluid saturation. The model reproduces the behaviour measured in a water-wet, nearly all-quartz, Bentheimer sandstone imaged during steady-state imbibition. A shift in fluid surface coverage is observed when the rock is altered to a new wetting state with crude oil. In two multi-mineralogical (Berea) samples, one water-wet and the other altered with crude oil, the analysis of fluid surface coverage after imbibition revealed mineral specific wetting preferences only in the altered system. Clays and calcite preferentially alter to an oil wet state, leading to mixed wettability in the rock.

Low salinity waterflooding is a low-cost method of enhancing oil recovery although, no consistent concept has been established explaining why some oil-fields show an increase in oil production when the salinity of the injected brine is reduced, while others do not. Various studies were conducted investigating the underlying mechanisms of the ‘low salinity effect’ using different crude oil, brine and rock compositions. Core floods of sandstone rock and analyses of molecular interactions using model systems indicate that clay content may play a dominant role. However, the spatial configuration of the sheet-like clay particles, which may vary from rock to rock, complicate comparisons of these model scenarios with reality. In the present study, we report the development of a pre-screening method using Atomic Force Microscopy (AFM) to assess rock-fluid interactions, which has previously only been used either on artificial model systems or minerals from crushed rock, by exploring the capability to operate in-situ in complex rock without crushing. The orientation of clay particles within a pore of an outcrop sandstone, Bandera Brown, was investigated with AFM and these particles were further assessed for changes in adhesion in brines of differing salinity. The results show a decrease in adhesions between CH3-functionalised AFM tips and the rock surface in low salinity brine, predominantly at the clay edges. This demonstrates that the edges of the clay particles, which may pin the oil phase after wettability alteration and therewith prevent oil from getting produced, lose this capacity when exposed to low salinity brine.

Capillary energy barriers have important consequences for immiscible fluid flow in porous media. We derive time-and-space averaging theory to account for non-equilibrium behavior and understand the role of athermal capillary fluctuations in the context of their relationship to larger scale phenomenological equations. The formulation resolves several key challenges associated with two-fluid flow in porous media: (1) geometric and thermodynamic quantities are constructed as smooth functions of time based on time-and space averages; (2) averaged thermodynamics are developed for films; (3) multi-scale fluctuation terms are identified, which account for transient behaviours of interfaces and films that occur due to pore-scale events; (4) geometric constraints are derived and imposed on the averaged thermodynamics; (5) a new constitutive model is proposed for capillary pressure dynamics that includes contributions from films; and (6) a time-and-space criterion for representative elementary volume (REV) is established based on capillary fluctuations. Capillary fluctuations are assessed quantitatively based on pore-scale simulations and experimental core-flooding data.

Hypothesis Understanding wetting behavior is of great importance for natural systems and technological applications. The traditional concept of contact angle, a purely geometrical measure related to curvature, is often used for characterizing the wetting state of a system. It can be determined from Young's equation by applying equilibrium thermodynamics. However, whether contact angle is a representative measure of wetting for systems with significant complexity is unclear. Herein, we hypothesize that topological principles based on the Gauss-Bonnet theorem could yield a robust measure to characterize wetting. Theory and Experiments We introduce a macroscopic contact angle based on the deficit curvature of the fluid interfaces that are imposed by contacts with other immiscible phases. We perform sessile droplet simulations followed by multiphase experiments for porous sintered glass and Bentheimer sandstone to assess the sensitivity and robustness of the topological approach and compare the results to other traditional approaches. Findings We show that the presented topological principle is consistent with thermodynamics under the simplest conditions through a variational analysis. Furthermore, we elucidate that at sufficiently high image resolution the proposed topological approach and local contact angle measurements are comparable. While at lower resolutions, the proposed approach provides more accurate results being robust to resolution-based effects. Overall, the presented concepts open new pathways to characterize the wetting state of complex systems and theoretical developments to study multiphase systems., Comment: 11 pages, 9 figures, 1 table

Many subsurface fluid flows, including the storage of CO2 underground or the production of oil, are transient processes incorporating multiple fluid phases. The fluids are not in equilibrium meaning macroscopic properties such as fluid saturation and pressure vary in space and time. However, these flows are traditionally modeled with equilibrium (or steady-state) flow properties, under the assumption that the pore-scale fluid dynamics are equivalent. In this work, we used fast synchrotron X-ray tomography with 1 s time resolution to image the pore-scale fluid dynamics as the macroscopic flow transitioned to steady state. For nitrogen or decane, and brine injected simultaneously into a porous rock, we observed distinct pore-scale fluid dynamics during transient flow. Transient flow was found to be characterized by intermittent fluid occupancy, whereby flow pathways through the pore space were constantly rearranging. The intermittent fluid occupancy was largest and most frequent when a fluid initially invaded the rock. But as the fluids established an equilibrium the dynamics decreased to either static interfaces between the fluids or small-scale intermittent flow pathways, depending on the capillary number and viscosity ratio. If the fluids were perturbed after an equilibrium was established, by changing the flow rate, the transition to a new equilibrium was quicker than the initial transition. Our observations suggest that transient flows require separate modeling parameters. The time scales required to achieve equilibrium suggest that several meters of an invading plume front will have flow properties controlled by transient pore-scale fluid dynamics.

We develop nonequilibrium theory by using averages in time and space as a generalized way to upscale thermodynamics in nonergodic systems. The approach offers a classical perspective on the energy dynamics in fluctuating systems. The rate of entropy production is shown to be explicitly scale dependent when considered in this context. We show that while any stationary process can be represented as having zero entropy production, second law constraints due to the Clausius theorem are preserved due to the fact that heat and work are related based on conservation of energy. As a demonstration, we consider the energy dynamics for the Carnot cycle and for Maxwell's demon. We then consider nonstationary processes, applying time-and-space averages to characterize nonergodic effects in heterogeneous systems where energy barriers such as compositional gradients are present. We show that the derived theory can be used to understand the origins of anomalous diffusion phenomena in systems where Fick's law applies at small length scales, but not at large length scales. We further characterize fluctuations in capillary-dominated systems, which are nonstationary due to the irreversibility of cooperative events.

The wetting state is an important control on flow in subsurface multi fluid phase systems, e.g., carbon storage and oil production. Advances in X-ray imaging allow us to characterise the wetting state using imagery of fluid arrangement within the pores of rocks. We derived a model from equilibrium thermodynamics relating fluid coverage of rock surfaces to wettability and fluid saturation. The model reproduces the behaviour measured in a water-wet, nearly all-quartz, Bentheimer sandstone imaged during steady-state imbibition. A shift in fluid surface coverage is observed when the rock is altered to a new wetting state with crude oil. In two multi-mineralogical (Berea) samples, one water-wet and the other altered with crude oil, the analysis of fluid surface coverage after imbibition revealed mineral specific wetting preferences only in the altered system. Clays and calcite preferentially alter to an oil wet state, leading to mixed wettability in the rock.

SummaryAdvances in pore-scale imaging, increasing availability of computational resources, and developments in numerical algorithms have started rendering direct pore-scale numerical simulations of multiphase flow on pore structures feasible. In this paper, we describe a two-phase-flow simulator that solves mass- and momentum-balance equations valid at the pore scale (i.e., at scales where the Darcy velocity homogenization starts to break down). The simulator is one of the key components of a molecule-to-reservoir truly multiscale modeling work flow.A Helmholtz free-energy-driven, thermodynamically based diffuse-interface/phase-field method is used for the effective simulation of numerous advecting interfaces, while honoring the interfacial tension (IFT). The advective Cahn-Hilliard (CH) (mass-balance, energy dissipation) and Navier-Stokes (NS) (momentum-balance, incompressibility) equations are coupled to each other within the phase-field framework. Wettability on rock/fluid interfaces is accounted for by means of an energy-penalty-based wetting (contact-angle) boundary condition. Individual balance equations are discretized by use of a flexible discontinuous Galerkin (DG) method. The discretization of the mass-balance equation is semi-implicit in time using a convex/concave splitting of the energy term. The momentum-balance equation is split from the incompressibility constraint by a projection method and linearized with a Picard splitting. Mass- and momentum-balance equations are coupled to each other by means of operator splitting, and are solved sequentially.We discuss the mathematical model and its DG discretization, and briefly introduce nonlinear and linear solution strategies. Numerical-validation tests show optimal convergence rates for the DG discretization, indicating the correctness of the numerical scheme and its implementation. Physical-validation tests demonstrate the consistency of the phase distribution and velocity fields simulated within our framework. Finally, two-phase-flow simulations on two real pore-scale images demonstrate the usefulness of the pore-scale simulator. The direct pore-scale numerical-simulation methodology rigorously considers the flow physics by directly acting on pore-scale images of rocks without remeshing. The proposed method is accurate, numerically robust, and exhibits the potential for tackling realistic problems.

Digital rock technology and pore-scale physics have become increasingly relevant topics in a wide range of porous media with important applications in subsurface engineering. This technology relies heavily on images of pore space and pore-level fluid distribution determined by X-ray microcomputed tomography (micro-CT). Digital images of pore space (or pore-scale fluid distribution) are typically obtained as gray-level images that first need to be processed and segmented to obtain the binary images that uniquely represent rock and pore (including fluid phases). This processing step is not trivial. Rock complexity, image quality, noise, and other artifacts prohibit the use of a standard processing workflow. Instead, an array of strategies of increasing sophistication has been developed. Typical processing pipelines consist of filtering, segmentation, and postprocessing steps. For each step, various choices and different options exist. This makes selection and validation of an optimum processing pipeline difficult. Using Darcy-scale quantities as a benchmark is not a good option because of rock heterogeneity and different scales of observation. Here, we present a conceptual workflow where noisy images are derived from a ground truth by systematically including typical image artifacts and noise. Artifacts and noise are not simply added to the images. Instead, tomographic forward projection and reconstruction steps are used to incorporate the artifacts in a physically correct way. A proof of concept of this workflow is demonstrated by comparing seven different image-segmentation pipelines ranging from absolute thresholding to a machine-learning approach (Trainable Weka Segmentation). The Trainable Weka Segmentation showed the best performance of the tested methods.

Digital rock physics (DRP) is a rapidly evolving technology targeting fast turnaround times for repeatable core analysis and multi-physics simulation of rock properties. We develop and validate a rapid and scalable distributed-parallel single-phase pore-scale flow simulator for permeability estimation on real 3D pore-scale micro-CT images using a novel variant of the lattice Boltzmann method (LBM). The LBM code implementation is designed to take maximum advantage of distributed computing on multiple general-purpose graphics processing units (GPGPUs). We describe and extensively test the distributed parallel implementation of an innovative LBM algorithm for simulating flow in pore-scale media based on the multiple-relaxation-time (MRT) model that utilizes a precise treatment of body force. While the individual components of the resulting simulator can be separately found in various references, our novel contributions are (1) the integration of all of the mathematical and high-performance computing components together with a highly optimized code implementation and (2) the delivery of quantitative results with the simulator in terms of robustness, accuracy, and computational efficiency for a variety of flow geometries including various types of real rock images. We report on extensive validations of the simulator in terms of accuracy and provide near-ideal distributed parallel scalability results on large pore-scale image volumes that were largely computationally inaccessible prior to our implementation. We validate the accuracy of the MRT-LBM simulator on model geometries with analytical solutions. Permeability estimation results are then provided on large 3D binary microstructures including a sphere pack and rocks from various sandstone and carbonate formations. We quantify the scalability behavior of the distributed parallel implementation of MRT-LBM as a function of model type/size and the number of utilized GPGPUs for a panoply of permeability estimation problems.

We generate a novel reference dataset to quantify the impact of numerical solvers, boundary conditions, and simulation platforms. We consider a variety of microstructures ranging from idealized pipes to digital rocks. Pore throats of the digital rocks considered are large enough to be well resolved with state-of-the-art micro-computerized tomography technology. Permeability is computed using multiple numerical engines, 12 in total, including, Lattice-Boltzmann, computational fluid dynamics, voxel based, fast semi-analytical, and known empirical models. Thus, we provide a measure of uncertainty associated with flow computations of digital media. Moreover, the reference and standards dataset generated is the first of its kind and can be used to test and improve new fluid flow algorithms. We find that there is an overall good agreement between solvers for idealized cross-section shape pipes. As expected, the disagreement increases with increase in complexity of the pore space. Numerical solutions for pipes with sinusoidal variation of cross section show larger variability compared to pipes of constant cross-section shapes. We notice relatively larger variability in computed permeability of digital rocks with coefficient of variation (of up to 25%) in computed values between various solvers. Still, these differences are small given other subsurface uncertainties. The observed differences between solvers can be attributed to several causes including, differences in boundary conditions, numerical convergence criteria, and parameterization of fundamental physics equations. Solvers that perform additional meshing of irregular pore shapes require an additional step in practical workflows which involves skill and can introduce further uncertainty. Computation times for digital rocks vary from minutes to several days depending on the algorithm and available computational resources. We find that more stringent convergence criteria can improve solver accuracy but at the expense of longer computation time.