Structure exploiting methods for fast uncertainty quantification in multiphase flow through heterogeneous media

The multiphase fluid flow in porous media is one of the most fundamental phenomena in various physical processes, such as oil/gas flow in reservoir, subsurface contamination dispersion, chemical separation, etc. Due to its importance, the efficient and accurate solution and prediction of multiphase flow in porous media is highly required in engineering applications and mechanism studies, which has been a research hot spot with increasing interest in recent years. However, the strong nonlinearity implicated in the multiphase flow model has brought great challenges for the computation and analysis. In addition, the permeability in Darcy-type pressure equation is always represented as a high contrast coefficient, which further complicates the simulation difficulty especially for high fidelity prediction and decision-making process. In this study, we attempt to apply one popular global model reduction method, namely Proper Orthogonal Decomposition (POD), to accelerate the simulation of incompressible multiphase flow in porous media. The cornerstone of POD is processing information from a sequence of snapshots (or instantaneous solutions) to identify a set of basis functions to construct a low-dimensional space, the most energetic part of these basis functions is then selected to represent the solution of the original system. In our work, one key point is how to establish the accurate POD reduced order model (ROM) for incompressible mul-tiphase porous flow problems. According to Li et al. [1], the traditional Galerkin projection method (directly project the original governing equation into the low-dimensional space spanned by basis functions) is incapable of introducing the fluid saturation and permeability located on the cell face in POD-ROM for incompressible multiphase flow whether in homogeneous or heterogeneous porous media, thus we adopt the matrix operation method (project the discrete governing equation instead of the original governing equation into the low-dimensional space spanned by basis functions) to build the POD-ROM in this study. Based on the proposed POD-ROM framework, we use the implicit method to solve the (non-) wetting-phase pressure and saturation equations. However, in the implicit solution of the pressure and saturation equations, the standard POD-ROM alone sometimes cannot yield an expected acceleration owing to two main reasons: (1) the nonlinearity between pressure and saturation is sensitive, thus the saturation would change largely even with slight pressure change. And the small error of pressure results would lead to a large deviation of saturation, which in turn exerts negative impacts on the computation of pressure; (2) the complexity for computing a projected nonlinear term in POD-ROM still depends on the dimension of the original full-order system, it is unavoidable in conventional POD-ROM. Therefore, more POD modes are needed to be used in order to achieve the prescribed numerical accuracy, which could worsen the speed up efficiency of proposed POD-ROM. To further reduce the computational complexity of POD-ROM to speed up the implicit solution of pressure and saturation equations, another model reduction method called Discrete Empirical Interpolation Method (DEIM) [2] is applied to treat the projected nonlinear term in the proposed incompressible multiphase porous flow POD-ROM. Different from previous studies, we develop an adaptive DEIM (ADEIM) in this work to achieve a more flexible selection of the interpolation points. Combined the POD-ROM with ADEIM, a POD-ADEIM-ROM for incompressible multiphase flow in porous media is established. The computational efficiency and numerical accuracy of the proposed model is validated through several representative numerical examples. Numerical results indicate that an outstanding acceleration (2~3 orders of magnitude faster) without sacrificing numerical accuracy significantly is obtained from the proposed model when comparing to the traditional solution method that without any acceleration technique. In addition, the effects of POD mode number, the distribution of permeability field, and the mesh density on the overall performances of the proposed model are investigated in detail.

Multiphase flow in porous media is involved in various natural and industrial applications, including water infiltration into soils, carbon geosequestration, and underground hydrogen storage. Understanding the invasion morphology at the pore scale is critical for better prediction of flow properties at the continuum scale in partially saturated permeable media. The deep learning method, as a promising technique to estimate the flow transport processes in porous media, has gained significant attention. However, existing works have mainly focused on single-phase flow, whereas the capability of data-driven techniques has yet to be applied to the pore-scale modeling of fluid–fluid displacement in porous media. Here, the conditional generative adversarial network is applied for pore-scale modeling of multiphase flow in two-dimensional porous media. The network is trained based on a data set of porous media generated using a particle-deposition method, with the corresponding invasion morphologies after the displacement processes calculated using a recently developed interface tracking algorithm. The results demonstrate the capability of data-driven techniques in predicting both fluid saturation and spatial distribution. It is also shown that the method can be generalized to estimate fluid distribution under different wetting conditions and particle shapes. This work represents the first effort at the application of the deep learning method for pore-scale modeling of immiscible fluid displacement and highlights the strength of data-driven techniques for surrogate modeling of multiphase flow in porous media.

The focus of the underlying research work is on the macroscopic modeling of unstable multiphase fluid flow in deformable porous media, where a lower‐viscous fluid is displaced by a more viscous fluid. This process leads to the formation of channel‐like networks, called viscous fingering. The instability effect is involved in a wide range of different fields in engineering. Some of the most common applications include carbon sequestration to store carbon dioxide (CO2) in underground reservoirs, contaminant transport in geostructures, in industrial processes, such as filtration, catalytic reactions, and in the operation of polymer electrolyte membrane fuel cells with multiphase flow in the gas diffusion layers. In this work, ideal miscible water–glycerin fluids are considered. It is assumed that the interacting fluids in the deformable porous media are incompressible. In addition, a dispersion–diffusion law is applied to capture the fluid–fluid interactions. The presence of a deformable porous material adds additional effects to the problem of the multiphase flow in porous media. The stresses in the porous solid matrix lead to changes in the porosity and influence the flow velocity of the fluids. To couple the deformable porous media with the multiphase flow, a macroscopic approach is used that relies on the theory of porous media. For the porous solid matrix, a linear elastic material model within small strains assumption is applied. The influence of the deformation‐dependent porosity on the instability is studied for 2D simulations, such as a multilayered geometry with different elastic parameters. The presented coupled nonlinear system of differential equations is simulated with the finite element method. Furthermore, a stabilization technique based on the quasi‐compressibility method is used.

A phase-field lattice Boltzmann model for simulating multiphase flows in porous media: Application and comparison to experiments of CO2 sequestration at pore scale

A 6M digital twin for modeling and simulation in subsurface reservoirs

Multiphase flow in porous media is relevant to amount of engineering processes, such as hydrocarbon extraction from reservoir rock, water contamination, CO2 geological storage and sequestration. Pore scale modeling, as an alternative approach to lab measurement, firstly serves as an effective bridge to link the pore scale properties (pore geometry and wettability) and displacement mechanisms to continuous scale multiphase flow in porous media; and secondly allows us to determine essential flow functions, such as capillary pressure and relative permeability curves, which are required for continuous scale modeling. In the literature, three methodologies, Bundle of Capillary Tube Modeling (BCTM), Direct Pore Scale Modeling (DPSM) and Pore Network Modeling (PNM), have appeared to be mostly widely adopted in the investigation of the pore-scale mechanics of fluid-fluid and fluid-solid interactions in porous media by numerical simulation. In this review article, a comprehensive review is provided to show their strengths and weaknesses and to highlight challenges that are faced in modelling of multiphase flow, key challenges include: are contact angle characterization, validation and upscale pore scale findings to core, or even field scale.

In this research study, we conducted a coupled thermalstress-fluid flow numerical modeling on the UTF Phase A SAGD project in order to investigate the fundamental geomechanical behavior involved in the SAGD process and to gain insight into the reservoir response to temperature and pore pressure changes. The numerical simulation is carried out using a self-developed coupled finite element model which incorporates our proposed strain-induced permeability model. The obtained simulation results were compared with the measured data. Introduction Thermal recovery processes involve coupling between heat transfer, multiphase flow and stress/deformation, which has become an increasingly important subject in the petroleum field(1). Particularly, the coupling is crucial in problems such as borehole stability, hydraulic fracturing and injection/production induced deformation of the ground surface during the thermal recovery processes in the heavy oil or oil sand reservoirs. Numerical modeling of the coupled processes is historically carried out in the areas of geomechanics modeling and the reservoir simulation. The former is to compute the stress-strain behavior, therefore the deformation; the latter is to essentially model the multiphase flow and heat transfer in porous media. Each of these disciplines simplifies the part of the problem that is not of primary interest. These approaches are unacceptable in situations where the coupling is strong and the changes of porosity and permeability cannot be accounted by rock compressibility alone. Gutierrez and Lewis(2) extend Biot's theory to multiphase fluid flow in deformable porous media. Based on their formulation, they conclude that the coupling between the geomechanics and the multiphase flow occurs simultaneously. Thus, fully coupled system equations of deformations, multiphase flow and heat transfer should be solved simultaneously. Development of such kinds of fully coupled geomechanics-multiphase flow-heat transfer simulators needs tremendous effort, since the existing FEM geomechanics codes and the FDM reservoir simulators cannot be used. Settari and Mourits(3) present an approach to couple the stress-strain behavior to multiphase flow, heat transfer using porosity as a coupling parameter. The geomechanics module and the thermal reservoir simulator are used in a staggered manner. Pore pressure and temperature changes are calculated from the thermal reservoir simulator and transferred to the geomechanics module. The stress and the displacement changes are then calculated in the geomechanics simulation. An iterative algorithm is used to ensure that the porosity calculated from the geomechanics module is the same as that from the thermal reservoir simulator. The staggered technique employed to solve the coupled system equations allows for the use of the existing geomechanics codes in conjunction with a standard reservoir simulator. Currently, most of the commercial coupled geomechanics-multiphase flow-heat transfer simulators are developed in this way. The disadvantage of these kinds of coupled simulators is that the thermal reservoir module, usually developed using finite difference method (FDM), cannot accommodate the full permeability tensor, since they adopt the standard discretization scheme such as 5-spot for 2-D problems and 7-spot for 3-D problems. In this paper, the development of a coupled geomechanicsmultiphase flow-heat transfer simulator using finite element method (FEM) is described with the use of the Galerkin's least squares (GLS) technique(4) to stabilized the saturation equation.

In this research study, we conducted a coupled thermal-stress-fluid flow numerical model on the UTF Phase A SAGD project in order to investigate the fundamental geomechanical behaviour involved in the SAGD process, and to gain insight into the reservoir response to temperature and pore pressure changes. The numerical simulation is carried out by using a self-developed coupled finite element model which incorporates our proposed strain-induced permeability model. The obtained simulation results were compared with the measured data. Introduction Thermal recovery processes involve coupling between heat transfer, multiphase flow and stress/deformation, which has become an increasingly important subject in the petroleum field(1). Particularly, the coupling is crucial in problems such as borehole stability, hydraulic fracturing and injection/production induced deformation of the ground surface during the thermal recovery process in heavy oil or oil sand reservoirs. Numerical modelling of the coupled processes is historically carried out in the areas of geomechanics modelling and reservoir simulation. The former is to compute the stress-strain behaviour; therefore, the deformation. The latter is to essentially model the multiphase flow and heat transfer in porous media. Each of these disciplines simplifies the part of the problem that is not of primary interest. These approaches are unacceptable in situations where the coupling is strong and the changes of porosity and permeability cannot be accounted for by rock compressibility alone. Gutierrez and Lewis(2) extend Biot's theory to multiphase fluid flow in deformable porous media. Based on their formulation, they conclude that the coupling between the geomechanics and the multiphase flow occurs simultaneously. Thus, fully coupled system equations of deformations, multiphase flow and heat transfer should be solved simultaneously. Development of such kinds of fully coupled geomechanics-multiphase flow-heat transfer simulators needs tremendous effort, since the existing FEM geomechanics codes and FDM reservoir simulators cannot be used. Settari and Mourits(3) present an approach to couple the stress-strain behaviour to multiphase flow-heat transfer using porosity as a coupling parameter. The geomechanics module and the thermal reservoir simulator are used in a staggered manner. Pore pressure and temperature changes are calculated from the thermal reservoir simulator and transferred to the geomechanics module. The stress and the displacement changes are then calculated in the geomechanics simulation. An iterative algorithm is used to ensure that the porosity calculated from the geomechanics module is the same as that from the thermal reservoir simulator. The staggered technique employed to solve the coupled system equations allows for the use of the existing geomechanics codes in conjunction with a standard reservoir simulator. Currently, most of the commercial coupled geomechanics-multiphase flow-heat transfer simulators are developed in this way. The disadvantage of these kinds of coupled simulators is that the thermal reservoir module, usually developed using the finite difference method (FDM), cannot accommodate the full permeability tensor, since they adopt the standard discretization scheme, such as 5-spot for 2D problems and 7-spot for 3D problems. In this paper, the development of a coupled geomechanics-multiphase flow-heat transfer simulator using the finite element method (FEM) is described with the use of Galerkin's least squares (GLS) technique(4) to stabilize the saturation equation.

Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries

This is the first year in the proposed three-year effort to develop high-resolution numerical methods for multi-phase flow in hierarchical porous media. The issues being addressed in this research are: Computational efficiency: Field-scale simulation of enhanced oil recovery, whether for energy production or aquifer remediation, is typically highly under-resolved. This is because rock transport properties vary on many scales, and because current numerical methods have low resolution. Effective media properties: Since porous media are formed through complex geologic processes, they involve significant uncertainty and scale-dependence. Given this uncertainty, knowledge of ensemble averages of flow in porous media can be preferable to knowledge of flow in specific realizations of the reservoir. However, current models of effective properties do not represent the observed behavior very well. Relative permeability models present a good example of this problem. In practice, these models seldom provide realistic representations of hysteresis, interfacial tension effects or three-phase flow; there are no models that represent well all three effects simultaneously. Wave propagation: It is common in the petroleum industry to assume that the models have the same well-posedness properties as the physical system. An example of this fallacy is given by the three-phase relative permeability models; they were widely assumed by the petroleum community to produce hyperbolic systems for the Buckley-Leverett equations, but later the mathematics community proved that these models inherently produce local elliptic regions. Since numerical methods must use the models for computations, oscillations that develop could erroneously be attributed to numerical error rather than modeling difficulties. During this year, we have made significant progress on several tasks aimed at addressing these issues.

Modeling multiphase flow in porous structures remains a challenge due to the complexity of handling multiple interfaces. This paper presents a one-domain pore-resolved simulation approach for immiscible two-phase flows in porous media, using a monolithic fluid–solid coupling framework to implicitly consider the existence of solid objects, with the fluid–fluid interfaces captured through solving an algebraic volume of fluid equation. Fluid interfacial tension is considered by integrating a continuum surface force, and the wall wettability condition is imposed by modifying the contact angle of the fluid interface at the embedded solid surface. The resulting equations are simple and stable, as there are no empirical models or parameters involved for the interface representation. This approach has been validated through performing a series of test-case simulations, including capillary-dominated flow, capillary rise with gravity, Taylor film formation, and finally two-phase flow in a heterogeneous porous structure. The numerical approach is demonstrated to be well suited for investigating pore-scale two-phase flows in realistic porous media.

We present a distributed active subspace method for training surrogate models of complex physical processes with high-dimensional inputs and function valued outputs. Specifically, we represent the model output with a truncated Karhunen-Lo\`eve (KL) expansion, screen the structure of the input space with respect to each KL mode via the active subspace method, and finally form an overall surrogate model of the output by combining surrogates of individual output KL modes. To ensure scalable computation of the gradients of the output KL modes, needed in active subspace discovery, we rely on adjoint-based gradient computation. The proposed method combines benefits of active subspace methods for input dimension reduction and KL expansions used for spectral representation of the output field. We provide a mathematical framework for the proposed method and conduct an error analysis of the mixed KL active subspace approach. Specifically, we provide an error estimate that quantifies errors due to active subspace projection and truncated KL expansion of the output. We demonstrate the numerical performance of the surrogate modeling approach with an application example from biotransport.

Fostered by advances in computing, simulation and data science have emerged as the third and fourth paradigms enabling scientific discovery. Many mathematical challenges arise when integrating these paradigms to produce new knowledge. One of the most important challenges is the inverse or calibration problem: find the inputs of a simulation such that its outputs agree with a given set of measurements. As both the cost of the complex simulations and the volume of measurement data sets increase, intuition-driven trial-and-error methods for calibration quickly become untenable; modern scientists need more computationally efficient and mathematically rigorous methods. When measurements and simulated predictions contain uncertainty, the Bayesian perspective on statistical inference provides a popular mathematical framework for developing computational methods to solve inverse problems. Markov chain Monte Carlo (MCMC) methods are a mainstay of this framework; they characterize a probability distribution on the space of model inputs that is conditioned on the given measurements. However, despite recent advances, MCMC methods remain impractical for expensive simulations with more than a handful of inputs. At its heart, MCMC searches the space of inputs for corresponding outputs that match measurements, and so it is afflicted with the curse of dimensionality. Additionally, MCMC methods become less efficient as the volume of the measurement data increases, particularly if the measurements are redundant. To make MCMC—and Bayesian methods, in general—amenable to today’s computer models with millions of inputs and terabytes of measurement data, one must devise methods to reduce the dimension of both the simulation inputs and the measurements. For many models used in practice, the outputs vary primarily on a low-dimensional manifold of the nominally high-dimensional input space. Recently developed active subspace methods discover this low-dimensional manifold and exploit it to enable analyses typically suited for low-dimensional settings (e.g., statistical inversion of complex systems). We hypothesize that active subspace methods can provide the necessary dimension reduction, enabling existing Bayesian methods to find solutions for a large class of otherwise intractable statistical inverse problems. We will test the hypothesis with a series of mathematical analyses and computational experiments that will proceed along three research veins: (i) integrating active subspace methods into methods for statistical inversion, (ii) addressing current challenges in active subspace methods, and (iii) incorporating methods for nonlinear dimensionality reduction into active subspace methods. We will apply the methods we develop to real inverse problems in chemical kinetics and turbulent flame modeling. Our team combines expertise in methods for inverse problems, dimension reduction, and reduced order modeling. With a history of active collaboration and thoughtful engagement, we bring enthusiasm and fresh perspectives to a well studied subject. We are excited by the potential opportunity to work together to extend and apply active subspace methods to data-intensive inverse problems.

This study aims to bridge length scales in immiscible multiphase flow simulation by connecting two published governing equations at the pore-scale and continuum-scale through a novel validation framework. We employ Niessner and Hassnaizadeh's [“A model for two-phase flow in porous media including fluid-fluid interfacial area,” Water Resour. Res. 44(8), W08439 (2008)] continuum-scale model for multiphase flow in porous media, combined with the geometric equation of state of McClure et al. [“Modeling geometric state for fluids in porous media: Evolution of the Euler characteristic,” Transp. Porous Med. 133(2), 229–250 (2020)]. Pore-scale fluid configurations simulated with the lattice-Boltzmann method are used to validate the continuum-scale results. We propose a mapping from the continuum-scale to pore-scale utilizing a generalized additive model to predict non-wetting phase Euler characteristics during imbibition, effectively bridging the continuum-to-pore length scale gap. Continuum-scale simulated measures of specific interfacial area, saturation, and capillary pressure are directly compared to up-scaled pore-scale simulation results. This research develops a numerical framework capable of capturing multiscale flow equations establishing a connection between pore-scale and continuum-scale simulations.

Data-driven uncertainty quantification in computational human head models