Two-point Flux Approximation Research Articles

Comparison of finite-volume schemes for diffusion problems

We present an abstract discretization framework and demonstrate that various cell-centered and hybrid finite-volume schemes fit into it. The different schemes considered in this work are then analyzed numerically for an elliptic model problem with respect to the properties consistency, coercivity, extremum principles, and sparsity. The test cases presented comprise of two- and three-dimensional setups, mildly and highly anisotropic tensors and grids of different complexities. The results show that all schemes show a similar convergence behavior, except for the two-point flux approximation scheme, and seem to be coercive. Furthermore, they confirm that linear schemes, in contrast to nonlinear schemes, are in general neither positivity-preserving nor satisfy discrete minimum or maximum principles.

A discrete modeling framework for reservoirs with complex fractured media: Theory, validation and case studies

An unstructured discrete fracture model (UDFM) is developed in this study to simulate the two-dimensional fluid flow in the complex fractured tight oil reservoirs. The model applies a control-volume finite-difference (CVFD) discretization scheme with the two-point flux-approximation (TPFA) to evaluate the transmissibility list. Discretization process of the mathematical model is greatly simplified by three differential operators that are defined to represent the expressions of divergence, gradient and average as their original forms, and thus, enable rapid prototyping of the numerical model. The discretized equation system is solved by a modified Newton's method, in which the Jacobian matrix is computed by the automatic differentiation (AD) technique. Our model is validated with the simulation results obtained by the commercial reservoir simulator, Eclipse. To test the model's potential in field-scale application, we implement it in five large-scale complex fractured systems that are radial domains with a uniform diameter of 1000 m but different fracture densities: 382/km2, 509/km2, 636/km2, 764/km2 and 1018/km2. Simulation results not only indicate the robustness of our model in simulating complex fractured systems but also highlight the significant effects of fracture network configuration and fracture density on the development features of fractured tight oil reservoirs.

Modeling and Simulation of Microbial Enhanced Oil Recovery Including Interfacial Area

The focus of this paper is the derivation of a nonstandard model for microbial enhanced oil recovery (MEOR) that includes the interfacial area (IFA) between oil and water. We consider the continuity equations for water and oil, a balance equation for the oil–water interfacial area, and advective–dispersive transport equations for bacteria, nutrients, and biosurfactants. Surfactants lower the interfacial tension (IFT), which improves oil recovery. Therefore, the parametrizations of the IFT reduction and residual oil saturation are included as a function of the surfactant concentration in the model. We consider for the first time in context of MEOR, the role of IFA in enhanced oil recovery. The motivation to include the IFA is to model the hysteresis in the capillary pressure–saturation relationship in a physically based manner, to include the effects of observed bacteria migration toward the oil–water interface and the production of biosurfactants at the oil–water interface. A comprehensive 2D implementation based on two-point flux approximation and backward Euler is proposed. An efficient and robust linearization scheme is used to solve the nonlinear systems at each time step. Illustrative numerical simulations are presented. The differences in the oil recovery profiles obtained with and without IFA are discussed. The presented model can also be used to design new experiments toward a better understanding and eventually optimization of MEOR.

A finite volume scheme with improved well modeling in subsurface flow simulation

We present the latest enhancement of the nonlinear monotone finite volume method for the near-well regions. The original nonlinear method is applicable for diffusion, advection-diffusion, and multiphase flow model equations with full anisotropic discontinuous permeability tensors on conformal polyhedral meshes. The approximation of the diffusive flux uses the nonlinear two-point stencil which reduces to the conventional two-point flux approximation (TPFA) on cubic meshes but has much better accuracy for the general case of non-orthogonal grids and anisotropic media. The latest modification of the nonlinear method takes into account the nonlinear (e.g., logarithmic) singularity of the pressure in the near-well region and introduces a correction to improve accuracy of the pressure and the flux calculation. In this paper, we consider a linear version of the nonlinear method waiving its monotonicity for sake of better accuracy. The new method is generalized for anisotropic media, polyhedral grids and nontrivial cases such as slanted, partially perforated wells or wells shifted from the cell center. Numerical experiments show noticeable reduction of numerical errors compared to the original monotone nonlinear FV scheme with the conventional Peaceman well model or with the given analytical well rate.

Numerical modelling and analysis of reactive flow and wormhole formation in fractured carbonate rocks

Combining the two-scale continuum model and the discrete fracture network model, a continuum-based model is developed that calculates the reactive flow of acid in carbonate rock with a complex fracture network. The locations of fractures in the model are explicitly defined and the method can capture complex geometric relationships. The governing equations are discretized by the finite-volume method, where the diffusion term and convection term are discretized using the two-point flux approximation (TPFA) scheme and the upwind scheme, respectively. The physical domain is discretized by Delaunay triangulation. To keep the robustness and efficiency of the TPFA scheme, the optimization algorithm is used to move the centroid node of the control volume to improve the orthogonality of the grids. Numerical simulations of reactive flow in 2D fractured porous media, in cases with simple and complex fracture arrays, under linear and radial flow conditions, are presented. In particular, a sensitivity analysis of the dissolution process with respect to the presence of fractures, fracture aperture, fracture distribution, and acid injection rate, is conducted.

Improved $L^2$ estimate for gradient schemes and super-convergence of the TPFA finite volume scheme

The gradient discretisation method is a generic framework that is applicable to a number of schemes for diffusion equations, and provides in particular generic error estimates in $L^2$ and $H^1$-like norms. In this paper, we establish an improved $L^2$ error estimate for gradient schemes. This estimate is applied to a family of gradient schemes, namely, the Hybrid Mimetic Mixed (HMM) schemes, and yields an $\mathcal O(h^2)$ super-convergence rate in $L^2$ norm, provided local compensations occur between the cell points used to define the scheme and the centers of mass of the cells. To establish this result, a modified HMM method is designed by just changing the quadrature of the source term; this modified HMM enjoys a super-convergence result even on meshes without local compensations. Finally, the link between HMM and Two-Point Flux Approximation (TPFA) finite volume schemes is exploited to partially answer a long-standing conjecture on the super-convergence of TPFA schemes.

A two-step finite volume method for the simulation of multiphase fluid flow in heterogeneous and anisotropic reservoirs

Efficient and accurate simulation of multiphase flow and transport in porous media plays a key role in optimizing recovery of hydrocarbon resources. Two-Point Flux Approximation (TPFA) is widely used in the industry for discretizing the pressure equation because of its robustness and efficiency despite the fact that it introduces O(1) flux errors on non-K-orthogonal grids. Multi-Point Flux Approximation (MPFA) provides a consistent discretization of the flow equations at the cost of larger flux stencils which leads to a denser, asymmetric system of linear equations requiring more computer memory and work for linear solvers. This drawback of MPFA becomes more pronounced for generally unstructured triangular grids and is even worse in 3D applications. In this work, we develop a two-step finite volume method (TSFVM) that is consistent for the anisotropic pressure equation on triangular and tetrahedral grids with improved efficiency as compared to MPFA methods. During each time step, the pressure equation is solved implicitly and saturation is then updated explicitly (IMPES). To discretize the pressure equation with full permeability tensors, Galerkin finite element method (FEM) is employed in the first step to compute pressure at grid vertices. Pressures at grid vertices are then utilized in the second step to derive a two-point flux stencil with a constant for each grid face. Pressure at cell centers and flux across grid faces can then be solved using a TPFA approach. Finally, the flux field at the new time level is used to update saturations explicitly. Although during each time step our TSFVM needs to solve two systems of linear equations, results of the numerical experiments clearly demonstrate that it is still less expensive than the classical MPFA-O method while providing high-quality numerical solutions matching that of MPFA-O. The gained efficiency is especially significant for 3D tetrahedral grid. As the scale of the problem increases, MPFA becomes less appealing because of its large flux stencils that result in much denser matrices while our TSFVM still remains computationally tractable thanks to the flexibility of FEM in the first step and the two-point flux stencil in the second step. On the other hand, the first step of the TSFVM is reminiscent of the Control Volume Finite Element Method (CVFEM) but with important differences. Unlike CVFEM, our TSFVM does not need to construct a dual control volume mesh around grid vertices on the basis of the primary mesh and the saturation unknowns are not associated with node-centered control volumes but are piece-wise constant on grid cells in which reservoir properties are distributed. This ensures that a saturation value does not straddle across material discontinuities which CVFEM suffer from.

A numerical simulation model for multi-scale flow in tight oil reservoirs

Abstract A discrete fracture model for multi-scale flow in large-scale fractured tight oil reservoirs is proposed considering the compressibility of reservoir rock and fluid, and the non-linear flow in the tight matrix. Validation of the model is performed, followed by the field application of the model. The two-point flux-approximation scheme is adopted in the model to calculate conductivity, and small grids at the fracture intersections are eliminated by the “star-delta” transformation method to improve the computational stability. The fully implicit discretization scheme is performed on the temporal domain. Automatic differentiation technique which can improve model establishment efficiency and computational accuracy is applied in the model to solve the numerical model. The model is validated with the simulation results of Eclipse and the historical production data of a long fractured horizontal well in a tight oil reservoir in Xinjiang oilfield. Simulation results of a field-scale reservoir show that the model proposed can simulate reservoirs with large-scale complex fracture systems; well productivity is positively correlated with the scale of the stimulated reservoir volume, and the difference in planar fracture density and fracture connectivity are proved to be the key factors that lead to the heterogeneous distribution of remaining oil in tight oil reservoirs.

A vertex-centered and positivity-preserving scheme for anisotropic diffusion problems on arbitrary polygonal grids

We suggest a new positivity-preserving finite volume scheme for anisotropic diffusion problems on arbitrary polygonal grids. The scheme has vertex-centered, edge-midpoint and cell-centered unknowns. The vertex-centered unknowns are primary and have finite volume equations associated with them. The edge-midpoint and cell-centered unknowns are treated as auxiliary ones and are interpolated by the primary unknowns, which makes the final scheme a pure vertex-centered one. Unlike most existing positivity-preserving schemes, the construction of the scheme is based on a special nonlinear two-point flux approximation that has a fixed stencil and does not require the convex decomposition of the co-normal. In order to solve efficiently the nonlinear systems resulting from the nonlinear scheme, Picard method and its Anderson acceleration are discussed. Numerical experiments demonstrate the second-order accuracy and well positivity of the solution for heterogeneous and anisotropic problems on severely distorted grids. The high efficiency of the Anderson acceleration is also shown on reduction of the number of nonlinear iterations. Moreover, the proposed scheme does not have the so-called numerical heat-barrier issue suffered by most existing cell-centered and hybrid schemes. However, further improvements have to be made if the solution is very close to the machine precision and the mesh distortion is very severe.

A simplified enhanced MPFA formulation for the elliptic equation on general grids

Multi-point flux approximation (MPFA) has proven to be a powerful tool for discretizing the diffusion equation on general grids with heterogeneous anisotropic permeability tensors and hence, removing the O(1) error introduced by two-point flux approximation (TPFA) for non-K-orthogonal grids. However, it is well known that the classical MPFA-O suffers from monotonicity issues and strong unphysical oscillations can be present for highly anisotropic media. Enriched MPFA (EMPFA) and MPFA with full pressure support (FPS) have been proposed in the literature to reduce the strength of the oscillations. In this work, we present a simplified enhanced MPFA formulation (eMPFA) on general grids in 2D. Similar to the MPFA-O method, our formulation starts with an interaction region formed around the vertex of control volumes. Potential at the vertex is introduced as an auxiliary unknown in addition to potential at the centroids of the control volume faces. The original EMPFA and FPS solve the diffusion equation on a small volume around the vertex to close the local system of equations. To simplify the formulation, we propose an average technique to approximate potential at the vertex. Utilizing the potential value at the vertex, full potential continuity can be imposed on control volume faces to construct a more accurate flux across control volume faces. Extensive numerical experiments are conducted to test the eMPFA formulation. Specifically, we compare eMPFA with MPFA-O and two variants of EMPFA in detail. The results show that our average technique is quite robust even for highly distorted quadrilateral grids with large permeability anisotropy ratios and solutions of our formulation are in excellent agreement with EMPFA using bilinear pressure support. The new scheme also reproduces linear potential solutions and has comparable, and in some cases even better, convergence properties to those of the MPFA-O method. Finally, we generalize the formulation to include structured and unstructured triangular and polygonal grids and present the results.

Second-order accurate finite volume schemes with the discrete maximum principle for solving Richards’ equation on unstructured meshes

Richards’s equation describes steady-state or transient flow in a variably saturated medium. For a medium having multiple layers of soils that are not aligned with coordinate axes, a mesh fitted to these layers is no longer orthogonal and the classical two-point flux approximation finite volume scheme is no longer accurate. We propose new second-order accurate nonlinear finite volume (NFV) schemes for the head and pressure formulations of Richards’ equation. We prove that the discrete maximum principles hold for both formulations at steady-state which mimics similar properties of the continuum solution. The second-order accuracy is achieved using high-order upwind algorithms for the relative permeability. Numerical simulations of water infiltration into a dry soil show significant advantage of the second-order NFV schemes over the first-order NFV schemes even on coarse meshes. Since explicit calculation of the Jacobian matrix becomes prohibitively expensive for high-order schemes due to build-in reconstruction and slope limiting algorithms, we study numerically the preconditioning strategy introduced recently in Lipnikov et al. (2016) that uses a stable approximation of the continuum Jacobian. Numerical simulations show that the new preconditioner reduces computational cost up to 2–3 times in comparison with the conventional preconditioners.

A generalized Mimetic Finite Difference method and Two-Point Flux schemes over Voronoi diagrams

We develop a generalization of the mimetic finite difference (MFD) method for second order elliptic problems that extends the family of convergent schemes to include two-point flux approximation (TPFA) methods over general Voronoi meshes, which are known to satisfy the discrete maximum principle. The method satisfies a modified consistency condition, which utilizes element and face weighting functions. This results in shifting the points on the elements and faces where the pressure and the flux are most accurately approximated. The flux bilinear form is non-symmetric in general, although it reduces to a symmetric form in the case of TPFA. It can be defined as the L 2 -inner product of vectors in two H ( Ω ;div) discrete spaces, which are constructed via suitable lifting operators. A specific construction of such lifting operators is presented on rectangles. We note that a different choice is made for test and trial spaces, therefore the method can be viewed as a H ( Ω ;div)-conforming Petrov–Galerkin Mixed Finite Element method. We prove first-order convergence in pressure and flux, and superconvergence of the pressure under further restrictions. We present numerical results that support the theory.

Implementation of a Vertex-Centered Method Inside an Industrial Reservoir Simulator: Practical Issues and Comprehensive Comparison With Corner-Point Grids and Perpendicular-Bisector-Grid Models on a Field Case

Summary Corner-point grids (CPGs) and pillar-based unstructured grids do not provide an effective work flow for translating Earth models into simulation models. Such a work flow requires grids that allow an accurate representation of the near-well flow, preserve geological accuracy, and offer flexible resolution control. Hence, a 3D unstructured approach is required. Significant work has been performed for generating unstructured grids, and modeling hydraulic-fracture flow for gas-shale simulation has given a new impulse for unstructured gridding. Recent methods such as vertex-approximate gradient (VAG) or more-mature ones such as multipoint flux approximation (MPFA) provide a numerical scheme dependent on multipoint stencil more physical than two-point flux-approximation (TPFA) methods. This paper presents the implementation of VAG and MPFA schemes inside a next-generation reservoir simulator starting from a source code calculating multipoint flux nonneighbor connections (NNCs) for any polygonal-shaped control volume. The unstructured-scheme approach has been developed as an in-house extension to a next-generation multicompany collaborative reservoir simulator that is designed for handling unstructured grids. The main issues addressed are the introduction of vertices unknowns among the usual cell-center variables, the assignment of vertices properties in the reservoir-simulator model, and the link with the well model. Incidentally, the definition of an exchange format to describe the unstructured geometry (vertices, edges, faces, control volumes) on a large reservoir-simulation model is proposed. Three simulation examples are presented, and we compare results, accuracy, and performance of multipoint-scheme methods such as VAG and MPFA on unstructured grids. Results are compared with TPFA methods on refined structured CPGs and TPFA methods on unstructured Voronoi grids. The two first test cases are academic models, and the third one is a field model.

Cell-centered nonlinear finite-volume methods for the heterogeneous anisotropic diffusion problem

We present two new cell-centered nonlinear finite-volume methods for the heterogeneous, anisotropic diffusion problem. The schemes split the interfacial flux into harmonic and transversal components. Specifically, linear combinations of the transversal vector and the co-normal are used that lead to significant improvements in terms of the mesh-locking effects. The harmonic component of the flux is represented using a conventional monotone two-point flux approximation; the component along the parameterized direction is treated nonlinearly to satisfy either positivity of the solution as in [29], or the discrete maximum principle as in [9]. In order to make the method purely cell-centered, we derive a homogenization function that allows for seamless interpolation in the presence of heterogeneity following a strategy similar to [46]. The performance of the new schemes is compared with existing multi-point flux approximation methods [3,5]. The robustness of the scheme with respect to the mesh-locking problem is demonstrated using several challenging test cases.

A Globally Coupled Pressure Method for the Discretization of the Tensor-Pressure Equation on Non-K-orthogonal Grids

Summary Complex permeability tensors together with general nonorthogonal and unstructured grids pose great challenges to reservoir simulation. The widely used two-point flux approximation (TPFA) is inadequate for a rigorous discretization of the flow equations on such challenging grids. Multipoint flux approximation (MPFA) methods have been proposed to meet the challenges and are currently being deployed in next-generation simulators. In this work, we propose an alternative flux-continuous cell-centered finite-volume method called the globally coupled pressure (GCP) method to discretize the pressure equation on general grids with full permeability tensors. To accurately construct fluxes through control-volume interfaces, pressure at the centroid of all interfaces is introduced as auxiliary unknowns. Flux continuity across each interface gives one equation. Assembling all the flux-continuity equations together gives a system of linear equations that can be solved simultaneously for all the auxiliary unknowns. Flux across control-volume interfaces can then be approximated with the pressure values at control-volume centers only. The fundamental difference between the GCP method and MPFA methods is that, in the latter, auxiliary unknowns are locally coupled within an interaction region and then eliminated in a local stencil by imposing flux-continuity conditions, whereas in the former, all the auxiliary unknowns are globally coupled and can only be eliminated in a global stencil. Consequently, control volumes in our GCP method are directly associated with the edges of the original grid and not by means of a dual grid overlaid and allied with the centers of the grid. Two variants of the GCP method are presented here, and extensive numerical experiments are conducted to test the performance of the GCP methods. The results show that both variants of our GCP method are in good agreement with the classical MPFA-O method on non-K-orthogonal grids for less-challenging problems. Convergence studies reveal that the first variant of our GCP method has slower convergence rates than the MPFA-O method for some problems. However, the second variant of GCP has comparable, and in some cases, better convergence properties compared with the MPFA-O method. With numerical experiments, we further investigate monotonicity properties of our GCP method on highly anisotropic media. For Dirichlet boundary conditions, our GCP methods also suffer from nonphysical oscillations, with some degrees of improvement over the MFPA-O method. When no-flow boundary conditions are used, our GCP method is much more robust and does not produce spurious boundary extrema as MPFA methods do. Finally, we extend our GCP method to fully unstructured grids, and the results show that it is also more robust than the MPFA-O method on unstructured grids.

An exponential integrator for finite volume discretization of a reaction–advection–diffusion equation

We consider the numerical approximation of a general second order semi-linear parabolic partial differential equation. Equations of this type arise in many contexts, such as transport in porous media. Using the finite volume with two-point flux approximation on regular mesh combined with exponential time differencing of order one (ETD1) for temporal discretization, we derive the L2 estimate under the assumption that the non linear term is locally Lipschitz. Numerical simulations to sustain the theoretical results are provided.

Monotone multiscale finite volume method

The MultiScale Finite Volume (MSFV) method is known to produce non-monotone solutions. The causes of the non-monotone solutions are identified and connected to the local flux across the boundaries of primal coarse cells induced by the basis functions. We propose a monotone MSFV (m-MSFV) method based on a local stencil-fix that guarantees monotonicity of the coarse-scale operator, and thus, the resulting approximate fine-scale solution. Detection of non-physical transmissibility coefficients that lead to non-monotone solutions is achieved using local information only and is performed algebraically. For these ‘critical’ primal coarse-grid interfaces, a monotone local flux approximation, specifically, a Two-Point Flux Approximation (TPFA), is employed. Alternatively, a local linear boundary condition can be used for the dual basis functions to reduce the degree of non-monotonicity. The local nature of the two strategies allows for ensuring monotonicity in local sub-regions, where the non-physical transmissibility occurs. For practical applications, an adaptive approach based on normalized positive off-diagonal coarse-scale transmissibility coefficients is developed. Based on the histogram of these normalized coefficients, one can remove the large peaks by applying the proposed modifications only for a small fraction of the primal coarse grids. Though the m-MSFV approach can guarantee monotonicity of the solutions to any desired level, numerical results illustrate that employing the m-MSFV modifications only for a small fraction of the domain can significantly reduce the non-monotonicity of the conservative MSFV solutions.

The single-cell transport problem for two-phase flow with polymer

Polymer injection is a widespread strategy in enhanced oil recovery. Polymer increases the water viscosity and creates a more favorable mobility ratio between the injected water and the displaced oil. The computational cost of simulating polymer injection can be significantly reduced if one splits the governing system of two-phase equations into a pressure equation and a set of saturation/component equations and use a Gauss–Seidel algorithm with optimal cell ordering to solve the nonlinear systems arising from an implicit discretization of the saturation/component equations. This approach relies on a robust single-cell solver that computes the saturation and polymer concentration of a cell, given the total flux and the saturation and polymer concentration of the neighboring cells. In this paper, we consider a relatively comprehensive polymer model used in an industry-standard simulator, and show that, in the case of a discretization using a two-point flux approximation, the single-cell problem always admits a solution that is also unique.

A Second-Order Positivity-Preserving Finite Volume Scheme for Diffusion Equations on General Meshes

We propose a new positivity-preserving finite volume scheme for the anisotropic diffusion problems on general polygonal meshes based on a new nonlinear two-point flux approximation. The scheme uses both cell-centered and cell-vertex unknowns. The cell-vertex unknowns are treated as auxiliary ones and are computed by two second-order interpolation algorithms. Due to the new nonlinear two-point flux formulation, it is not required to replace the interpolation algorithm with positivity-preserving but usually low-order accurate ones whenever negative interpolation weights occur and it is also unnecessary to require the decomposition of the conormal vector to be a convex one. Moreover, the new nonlinear two-point flux approximation has a fixed stencil. These features make our scheme more flexible, easy for implementation, and different from other existing nonlinear schemes based on Le Potier's two-point flux approximation. The positivity-preserving property of our scheme is proved theoretically and numerical r...

Interpolation-based second-order monotone finite volume schemes for anisotropic diffusion equations on general grids

We propose two interpolation-based monotone schemes for the anisotropic diffusion problems on unstructured polygonal meshes through the linearity-preserving approach. The new schemes are characterized by their nonlinear two-point flux approximation, which is different from the existing ones and has no constraint on the associated interpolation algorithm for auxiliary unknowns. Thanks to the new nonlinear two-point flux formulation, it is no longer required that the interpolation algorithm should be a positivity-preserving one. The first scheme employs vertex unknowns as the auxiliary ones, and a second-order but not positivity-preserving interpolation algorithm is utilized. The second scheme uses the so-called harmonic averaging points located on cell edges to define the auxiliary unknowns, and a second-order positivity-preserving interpolation method is employed. Both schemes have nearly the same convergence rates as compared with their related second-order linear schemes. Numerical results demonstrate that the new schemes are monotone, and have the second-order accuracy for the solution and first-order for its gradient on severely distorted meshes.