Numerical methods for coupled processes in fractured porous media

Numerical simulations have become essential in the planning and execution of operations in the subsurface, whether this is geothermal energy production or storage, carbon sequestration, petroleum production, or wastewater disposal. As the computational power increases, more complex models become feasible, not only in the form of more complicated physics, but also in the details of geometric constraints such as fractures, faults and wells. These features are often of interest as they can have a profound effect on different physical processes in the porous medium. This thesis focuses on modeling and simulations of fluid flow, transport and deformation of fractured porous media. The physical processes are formulated in a mixed-dimensional discrete fracture matrix model, where the rock matrix, fractures, and fracture intersections form a hierarchy of subdomains of different dimensions that are coupled through interface laws. A new discretization scheme for solving the deformation of a poroelastic rock coupled to a Coulomb friction law governing fracture deformation is presented. The novelty of this scheme comes from combining an existing finite-volume discretization for poroelasticity with a hybrid formulation that adds Lagrange multipliers on the fracture surface. This allows us to formulate the inequalities as complementary functions and solve the corresponding non-linear system using a semi-smooth Newton method. The mixed-dimensional framework is used to investigate non-linear coupled flow and transport. Here, we study how highly permeable fractures affect the viscous fingering in a porous medium and show that there is a complex interplay between the unstable viscous fingers and the fractures. The computer code of the above contributions of the thesis work has been implemented in the open-source framework PorePy. The introduction of fractures is a challenge to the discretization and the implementation of the governing equations, and the aim of this framework is to enable researchers to overcome many of the technical difficulties inherent to fractures, allowing them to easily develop models for fractured porous media. One of the large challenges for the mixed-dimensional discrete fracture matrix models is to create meshes that conform to the fractures, and we present a novel algorithm for constructing conforming Voronoi meshes. The proposed algorithm creates a mesh hierarchy, where the faces of the rock matrix mesh conform to the cells of the fractures, and the faces of the fracture mesh conform to the cells of the fracture intersections. The flexibility of the mixed-dimensional framework is exemplified by the wide range of applications and models studied within this thesis. While these physical processes might be fairly well known in a porous medium without fractures, the results of this thesis improves our understanding as well as the models and solution strategies for fractured porous media

Numerical methods for all-speed flows in fluid-dynamics and non-linear elasticity

In this thesis we are concerned with the numerical simulation of compressible materials flows, including gases, liquids and elastic solids. These materials are described by a monolithic Eulerian model of conservation laws, closed by an hyperelastic state law that includes the different behaviours of the considered materials. A novel implicit relaxation scheme to solve compressible flows at all speeds is proposed, with Mach numbers ranging from very small to the order of unity. The scheme is general and has the same formulation for all the considered materials, since a direct dependence on the state law is avoided via the relaxation. It is based on a fully implicit time discretization, easily implemented thanks to the linearity of the transport operator in the relaxation system. The spatial discretization is obtained by a combination of upwind and centered schemes in order to recover the correct numerical viscosity in different Mach regimes. The scheme is validated with one and two dimensional simulations of fluid flows and of deformations of compressible solids. We exploit the domain discretization through Cartesian grids, allowing for massively parallel computations (HPC) that drastically reduce the computational times on 2D test cases. Moreover, the scheme is adapted to the resolution on adaptive grids based on quadtrees, implementing adaptive mesh refinement techinques. The last part of the thesis is devoted to the numerical simulation of heterogeneous multi-material flows. A novel sharp interface method is proposed, with the derivation of implicit equilibrium conditions. The aim of the implicit framework is the solution of weakly compressible and low Mach flows, thus the proposed multi-material conditions are coupled with the implicit relaxation scheme that is solved in the bulk of the flow. Dans cette thèse on s’intéresse à la simulation numérique d’écoulements des matériaux compressibles, voir fluides et solides élastiques. Les matériaux considérés sont décrits avec un modèle monolithique eulérian, fermé avec une loi d’état hyperélastique qui considère les différents comportéments des matériaux. On propose un nouveau schéma de relaxation qui résout les écoulements compressibles dans des différents régimes, avec des nombres de Mach très petits jusqu’à l’ordre 1. Le schéma a une formulation générale qui est la même pour tous le matériaux considérés, parce que il ne dépend pas directement de la loi d’état. Il se base sur une discrétization complétement implicite, facile à implémenter grâce à la linearité de l’opérateur de transport du système de relaxation. La discrétization en éspace est donnée par la combinaison de flux upwind et centrés, pour retrouver la correcte viscosité numérique dans les différents régimes. L’utilisation de mailles cartésiennes pour les cas 2D s’adapte bien à une parallélisation massive, qui permet de réduire drastiquement le temps de calcul. De plus, le schéma a été adapté pour la résolution sur des mailles quadtree, pour implémenter l’adaptivité de la maille avec des critères entropiques. La dernière partie de la thèse concerne la simulation numérique d’écoulements multi-matériaux. On a proposé une nouvelle méthode d’interface “sharp”, en dérivant les conditions d’équilibre en implicite. L’objectif est la résolution d’interfaces physiques dans des régimes faiblement compressibles et avec un nombre de Mach faible, donc les conditions multi-matériaux sont couplées au schéma implicite de relaxation

Two-phase flow properties upscaling in heterogeneous porous media

The groundwater specialists and the reservoir engineers share the same interest in simulating multiphase flow in soil with heterogeneous intrinsic properties. They also both face the challenge of going from a well-modeled micrometer scale to the reservoir scale with a controlled loss of information. This upscaling process is indeed worthy to make simulation over an entire reservoir manageable and stochastically repeatable. Two upscaling steps can be defined: one from the micrometer scale to the Darcy scale, and another from the Darcy scale to the reservoir scale. In this thesis, a new second upscaling multiscale algorithm Finite Volume Mixed Hybrid Multiscale Methods (Fv-MHMM) is investigated. Extension to a two-phase flow system is done by weakly and sequentially coupling saturation and pressure via IMPES-like method

Modelling multiple-material flows on adaptive unstructured meshes.

The ability to distinguish between regions with different material properties is essential when numerically modelling many physical systems. Using a dual control volume mesh that avoids the problem of corner coupling, the HyperC face value scheme is extended to multiple dimensions and applied as a device for material advection on unstructured simplex meshes. The new scheme performs well at maintaining sharp interfaces between materials and is shown to produce small advection errors, comparable to those of standard material advection methods on structured meshes. To further minimise numerical diffusion of material interfaces a total variation bounded flux limiter, UltraC, is defined using a normalised variable diagram. Combining the material tracking scheme with dynamically adapting meshes, the use of a minimally dissipative bounded projection algorithm for interpolating fields from the old mesh to the new, optimised mesh is demonstrated that conserves the mass of each material. More generally, material conservation during the advection process is ensured through the coupling of the material tracking scheme to the momentum and mass equations. A new element pair for the discretisation of velocity and pressure is proposed that maintains the stability of the system while conserving the mass of materials. When modelling multiple materials the use of independent advection algorithms for each material can lead to the problem of non-physical material overlap. A novel coupled flux limiter is developed to overcome this problem. This automatically generalises to arbitrary numbers of materials. Using the fully coupled (and rigorously verified) multi-material model, several geophysically relevant simulations are presented examining the generation of waves by landslides. The model is validated by demonstrating close agreement between model predictions and experimental results of wave generation, propagation and run-up. The simulations also showcase the powerful capabilities of an unstructured, adaptive multi-material model in real world scenarios

Complex flow and transport phenomena in porous media

This thesis analyzes partial differential equations related to the coupled surface and subsurface flows and develops efficient high order discontinuous Galerkin (DG) methods to solve them numerically. Specifically, the coupling of the Navier-Stokes and the Darcy's equations, which is encountered in the environmental problem of groundwater contamination through lakes and rivers, is considered. Predicting accurately the transport of contaminants by this coupled flow is of great importance for the remediation strategies. The first part of this thesis analyzes a weak formulation of the time-dependent Navier-Stokes equation coupled with the Darcy's equation through the Beavers-Joseph-Saffman condition. The analysis changes depending on whether the inertial forces are included in the interface conditions or not. The inclusion of the inertial forces (Model I) remedies the difficulty in the analysis caused by the nonlinear convection term; however, it does not reflect the physical interactions on the interface correctly. Hence, I also analyze the weak problem by omitting these forces (Model II) which complicates the analysis and necessitates an extra small data condition. For Model I, a fully discrete scheme based on the DG method and the Crank-Nicolson method is introduced. The convergence of the scheme is proven with optimal error estimates. The second part couples the surface flow and a convection-diffusion type transport with miscible displacement in the subsurface. Initially, I consider the coupled stationary Stokes and Darcy's equations for the flow and establish the existence of a weak solution. Next, imposing additional assumptions on the data, I extend the result to the nonlinear case where the surface flow is given by the Navier-Stokes equation. The analysis also applies to the particular case where the flow is loosely coupled to the transport, that is, the velocity field obtained from the flow is an input for the transport equation. The flow is discretized by combinations of the continuous finite element method and the DG method whereas the discretization of the transport is done by a combined DG and backward Euler methods. The scheme yields optimal error estimates and its robustness for fractured porous media is shown by a numerical example

hp-FEM for Two-component Flows with Applications in Optofluidics

This thesis is concerned with the application of hp-adaptive finite element methods to a mathematical model of immiscible two-component flows. With the aim of simulating the flow processes in microfluidic optical devices based on liquid-liquid interfaces, we couple the time-dependent incompressible Navier-Stokes equations with a level set method to describe the flow of the fluids and the evolution of the interface between them

A finite volume approach for the numerical analysis and solution of the Buckley-Leverett equation including capillary pressure

The study of petroleum recovery is significant for reservoir engineers. Mathematical models of the immiscible displacement process contain various assumptions and parameters, resulting in nonlinear governing equations which are tough to solve. The Buckley-Leverett equation is one such model, where controlling forces like gravity and capillary forces directly act on saturation profiles. These saturation profiles have important features during oil recovery. In this thesis, the Buckley-Leverett equation is solved through a finite volume scheme, and capillary forces are considered during this calculation. The detailed derivation and calculation are also illustrated here. First, the method of characteristics is used to calculate the shock speed and characteristics curve behaviour of the Buckley-Leverett equation without capillary forces. After that, the local Lax-Friedrichs finite-volume scheme is applied to the governing equation (assuming there are no capillary and gravity forces). This mathematical formulation is used for the next calculation, where the cell-centred finite volume scheme is applied to the Buckley- Leverett equation including capillary forces. All calculations are performed in MATLAB. The fidelity is also checked when the finite-volume scheme is computed in the case where an analytical solution is known. Without capillary pressure, all numerical solutions are calculated using explicit methods and smaller time steps are used for stability. Later, the fixed-point iteration method is followed to enable the stability of the local Lax-Friedrichs and Cell-centred finite volume schemes using an implicit formulation. Here, we capture the number of iterations per time-steps (including maximum and average iterations per time-step) to get the solution of water saturation for a new time-step and obtain the saturation profile. The cumulative oil production is calculated for this study and illustrates capillary effects. The influence of viscosity ratio and permeability in capillary effects is also tested in this study. Finally, we run a case study with valid field data and check every calculation to highlight that our proposed numerical schemes can capture capillary pressure effects by generating shock waves and providing single-valued saturation at each position. These saturation profiles help find the amount of water needed in an injection well to displace oil through a production well and obtains good recovery using the water flooding technique

Evolution equations in physical chemistry

textWe analyze a number of systems of evolution equations that arise in the study of physical chemistry. First we discuss the well-posedness of a system of mixing compressible barotropic multicomponent flows. We discuss the regularity of these variational solutions, their existence and uniqueness, and we analyze the emergence of a novel type of entropy that is derived for the system of equations. Next we present a numerical scheme, in the form of a discontinuous Galerkin (DG) finite element method, to model this compressible barotropic multifluid. We find that the DG method provides stable and accurate solutions to our system, and that further, these solutions are energy consistent; which is to say that they satisfy the classical entropy of the system in addition to an additional integral inequality. We discuss the initial-boundary problem and the existence of weak entropy at the boundaries. Next we extend these results to include more complicated transport properties (i.e. mass diffusion), where exotic acoustic and chemical inlets are explicitly shown. We continue by developing a mixed method discontinuous Galerkin finite element method to model quantum hydrodynamic fluids, which emerge in the study of chemical and molecular dynamics. These solutions are solved in the conservation form, or Eulerian frame, and show a notable scale invariance which makes them particularly attractive for high dimensional calculations. Finally we implement a wide class of chemical reactors using an adapted discontinuous Galerkin finite element scheme, where reaction terms are analytically integrated locally in time. We show that these solutions, both in stationary and in flow reactors, show remarkable stability, accuracy and consistency.Chemistry and Biochemistr