Interior boundary-aligned unstructured grid generation and cell-centered versus vertex-centered CVD-MPFA performance

Grid generation for reservoir simulation must honor classical key constraints and ensure boundary alignment such that control-volume boundaries are aligned with geological features including layers, shale barriers, fractures, faults, pinch-outs, and multilateral wells. Novel unstructured grid generation methods are proposed that automate control-volume and/or control point boundary alignment and yield perpendicular-bisector (PEBI) meshes both with respect to primal and dual (essentially PEBI) cells. In order to honor geological features in the primal configuration, we introduce the idea of protection circles that contain segments of key geological boundaries, while in order to generate a dual-cell feature aligned grid, we construct halos around key geological features. The grids generated are employed to study comparative performance of cell-centred versus cell-vertex flux-continuous control-volume distributed multi-point flux approximation (CVD-MPFA) finite-volume formulations using equivalent degrees of freedom and thus ensure application of the most efficient methods. The CVD-MPFA formulation (c.f. Edwards et al.) in cell-centred and cell-vertex modes is somewhat analogous and requires switching control-volume from primal to dual or vice versa, together with appropriate data structures and boundary conditions, however dual-cells are generated after primal grid generation. The relative benefits of both types of approximation, i.e., cell-centred versus vertex-centred, are contrasted in terms of flow resolution and degrees of freedom required

Quasi-positive families of flux continuous finite volumes schemes in two and three dimensions.

In this thesis, new families of full pressure support flux-continuous, locally conservative, finite-volume schemes are presented for solving the general geometry-permeability tensor pressure equation on structured and unstructured grids in two and three dimensions. The families of flux-continuous schemes have also been referred to in the literature as Multi-point Flux Approximation or MPFA schemes. The schemes are applicable to the general tensor pressure equation with discontinuous coefficients and remove the 0(1) errors introduced by standard reservoir simulation (two-point flux) schemes when applied to full, anisotropic and asymmetric permeability tensor flow approximation. Such tensors may arise when fine scale permeability distributions are upscaled to obtain gridblock-scale permeability distributions. In contrast to the previous MPFA schemes which assume point-wise pressure and flux continuity locally, the new families of schemes presented in the work recover full pressure continuity across the interface between neighboring subcells. The M-matrix conditions [1, 2] define the upper limits for ensuring a local maximum principle is obtained for full-tensor fields. A key condition is that the modulus of the off-diagonal tensor coefficients are bounded by the minimum of the diagonal coefficients. For higher anisotropic ratios, when the resulting discrete matrices violate these bounds these schemes can violate the maximum principle (as with more standard methods) and the numerical pressure solutions can consequently exhibit spurious oscillations. The new family of schemes yield improved performance for challenging problems where earlier flux-continuous schemes exhibit strong spurious oscillations. The M- matrix analysis leads to an optimal quadrature range for these methods. The degree of freedom within the family of full pressure continuity schemes presented is shown to maximise the quadrature range of the flux-continuous schemes. For strongly anisotropic full-tensor cases where M-matrix conditions are violated, it is shown that the earlier families of schemes cannot avoid decoupling of the solution which leads to severe spurious oscillations in the discrete solution. The full quadrature range of the new schemes permits use of quadrature points that were previously out of range for the earlier methods, and that the resulting schemes minimize spurious oscillations in discrete pressure solutions. The new formulation leads to a more robust quasi-positive family of flux-continuous schemes applicable to general discontinuous full-tensor fields. This work also extends the single parameter family of FPS schemes to double families of schemes with general flexibility in quadrature that allow different quadrature points to be used on different control-volume subfaces. The new schemes minimize spurious oscillations in discrete pressure solutions. The new formulation leads to more robust quasi-positive families of flux-continuous schemes applicable to general discontinuous full-tensor fields. The full pressure support flux continuous schemes also extend to 3D on structured and unstructured grids. Surface auxiliary control volume and volume auxiliary control volume are introduced to handle extra degrees of freedom which are required for full pressure continuity over neighboring subcell surface. The new schemes are shown to be beneficial in high anisotropic test cases while remaining comparable with previous tetrahedral pressure support (TPS) schemes in terms of convergence rate. Multi-family schemes in 3D are also presented in this work. This is the extension of 2D double family to 3D. Compared to single family FPS schemes, multi-family schemes are shown to be able to maximize the quadrature and have incomparable flexibility over previous schemes, leading to improved solutions

Multidimensional upwind schemes and higher resolution methods for three-component two-phase systems including gravity driven flow in porous media on unstructured grids

Standard reservoir simulation schemes employ single-point upstream weighting for approximation of the convective fluxes when multiple phases or components are present. These schemes introduce both coordinate-line numerical diffusion and crosswind diffusion into the solution that is grid and geometry dependent.Families of locally conservative multidimensional upwind schemes are presented for essentially hyperbolic three-component two-phase flow systems of conservation laws in porous media including counter current gravity flow on unstructured grids. The multidimensional methods employ cell-based tracing, which involves tracing characteristic wave directions over each control-volume subquadrant. The multidimensional methods reduce crosswind diffusion inherent in standard methods for convective flow approximation in porous media. The schemes are coupled with continuous Darcy-flux approximations resulting from the elliptic pressure equation on unstructured grids.Characteristic upwind approximations are proposed and compared with the classical upstream weighting schemes for cases including gravity segregated flow. When dealing with systems of hyperbolic equations, upwind characteristic wave decomposition is used for wave tracing. The multidimensional upwind cell-based tracing formulations are designed for unstructured grids (and include structured grids by default) and are stable subject to conditions on the tracing direction and CFL number and satisfy a local maximum principle that ensures solutions are free of spurious oscillations.Benefits of the resulting schemes are demonstrated for two-phase flow and a three-component two-phase flow system including gravity segregated flow. The multidimensional cell based schemes are shown to reduce crosswind diffusion induced by standard upwind methods, and prove to be particularly effective when flow is strongly non-aligned with the grid, leading to improved resolution of numerical saturation and concentration fronts. Extension of higher order schemes to a three-component two-phase flow systems of conservation laws on unstructured grids is also presented, which provides a significant improvement in flow resolution for the system cases. Comparison is drawn between the methods

Simulation of multi-component multi-phase fluid flow in two-dimensional anisotropic heterogeneous porous media using high-order control volume distributed methods

Peer reviewedPostprin

Multi-dimensional higher resolution methods for flow in porous media.

Currently standard first order single-point upstream weighting methods are employed in reservoir simulation for integrating the essentially hyperbolic system components. These methods introduce both coordinate-line numerical diffusion (even in 1-D) and cross-wind diffusion into the solution that is grid and geometry dependent. These effects are particularly important when steep fronts and shocks are present and for cases where flow is across grid coordinate lines. In this thesis, families of novel edge-based and cell-based truly multidimensional upwind formulations that upwind in the direction of the wave paths in order to minimise crosswind diffusion are presented for hyperbolic conservation laws on structured and unstructured triangular and quadrilateral grids in two dimensions. Higher resolution as well as higher order multidimensional formulations are also developed for general structured and unstructured grids. The schemes are coupled with existing consistent and efficient continuous CVD (MPFA) Darcy flux approximations. They are formulated using an IMPES (Implicit in Pressure Explicit in Saturation) strategy for solving the coupled elliptic (pressure) and hyperbolic (saturation) system of equations governing the multi-phase multi-component flow in porous media. The new methods are compared with single point upstream weighting for two-phase and three-component two-phase flow problems. The tests arc conducted on both structured and unstructured grids and involve full-tensor coefficient velocity fields in homogeneous and heterogeneous domains. The comparisons demonstrate the benefits of multidimensional and higher order multidimensional schemes in terms of improved front resolution together with significant reduction in cross-wind diffusion

Control-volume distributed multi-point flux approximation coupled with a lower-dimensional fracture model

A novel cell-centered control-volume distributed multi-point flux approximation (CVD-MPFA) finite-volume formulation is presented for discrete fracture-(rock)matrix flow simulations. The grid is aligned with the fractures and barriers which are then modeled by lower-dimensional interfaces located between rock matrix cells in the physical domain. The n D (n-dimension) pressure equation in the rock matrix is coupled with the (n−1)D pressure equation in the fractures, leading to the first reduced dimensional flux-continuous CVD-MPFA formulation. This formulation naturally handles fractures efficiently on unstructured grids. Matrix-fracture fluxes are expressed in terms of matrix and fracture pressures, resulting in a transfer function, which is added to the lower-dimensional flow equation. An additional transmission condition is used between matrix cells separated by low permeable fractures to couple the velocity and pressure jump across the fractures. Numerical tests serve to assess the convergence and accuracy of the lower-dimensional fracture model for lower anisotropic fractures having a range of apertures and permeability tensors. A tracer flow transport equation is solved for problems with single and intersecting fractures. A lower-dimensional mass balance for intersecting fracture cells circumvents the more restrictive CFL condition corresponding to standard equi-dimensional approximation with explicit time discretization. Lower-dimensional fracture model results are compared with hybrid-grid and equi-dimensional model results. Fractures and barriers are efficiently modeled by lower-dimensional interfaces which yield comparable results to those of the equi-dimensional model. Highly conductive fractures are modeled as lower-dimensional entities without the use of locally refined grids that are required by the equi-dimensional model, while pressure continuity across fractures is built into the model, without depending on the extra degrees of freedom which must be added locally by the hybrid-grid method. The lower-dimensional fracture model also yields improved results when compared to those of the hybrid-grid model for fractures with low-permeability in the normal direction to the fracture where pressure is discontinuous. In addition, transient pressure simulation involving geologically representative complex fracture networks is presented

Three-dimensional unstructured gridding for complex wells and geological features in subsurface reservoirs, with CVD-MPFA discretization performance

Grid generation for reservoir simulation, must honour classical key geological features and multilateral wells. The features to be honoured are classified into two groups; (1) involving layers, faults, pinchouts and fractures, and (2) involving well distributions. In the former, control-volume boundary aligned grids (BAGs) are required, while in the latter, control-point (defined as the centroid of the control-volume) well aligned grids (WAGs) are required. Depending on discretization method type and formulation, a choice of control-point and control-volume type is made, i.e. for a cell-centered method the primal grid cells act as control-volumes, otherwise for a vertex-centered method the dual-grid cells act as control-volumes. Novel three-dimensional unstructured grid generation methods are proposed that automate control-volume boundary alignment to geological features and control point alignment to complex wells, yielding essentially perpendicular bisector (PEBI) meshes either with respect to primal or dual-cells depending on grid type. Both grid types use tetrahedra, pyramids, prisms and hexahedra as grid elements. Primal-cell feature aligned grids are generated using special boundary surface protection techniques together with constrained cell-centered well trajectory alignment. Dual-cell feature aligned grids are generated from underlying primal-meshes, whereby features are protected such that dual-cell control-volume faces are aligned with interior feature boundaries, together with protected vertex-centered (control point) well trajectory alignment. The novel methods of grid generation presented enable practical application of both method types in 3-D for the first time. The primal and dual grids generated here demonstrate the gridding methods, and enable the first comparative performance study of cell-vertex versus cell-centered control-volume distributed multi-point flux approximation (CVD-MPFA) finite-volume formulations using equivalent mesh resolution on challenging problems in 3-D. Pressure fields computed by the cell-centered and vertex-centered CVD-MPFA schemes are compared and contrasted relative to the respective degrees of freedom employed, and demonstrate the relative benefits of each approximation type. Stability limits of the methods are also explored. For a given mesh the cell-vertex method uses approximately a fifth of the unknowns used by a cell-centered method and proves to be the most beneficial with respect to accuracy and efficiency. Numerical results show that vertex-centered CVD-MPFA methods outperform cell-centered CVD-MPFA method