Polyhedral Grids Research Articles

Numerical modeling of non-Newtonian fluid flow in fractures and porous media

Non-Newtonian fluids having Bingham or power-law rheology are common in many applications within drilling and reservoir engineering. Examples of such fluids are drilling muds, foams, heavy oil, hydraulic-fracturing and other stimulation fluids, and cement slurries. Despite the importance of non-Newtonian rheology, it is rarely used in reservoir simulators and fracture flow simulations. We study two types of non-Newtonian rheology: the truncated power-law (Ostwald-de Waele) fluid and the Bingham fluid. For either of the two types of non-Newtonian rheology, we construct relationships between the superficial fluid velocity and the pressure gradient in fractures and porous media. The Bingham fluid is regularized by means of Papanastasiou-type regularization for porous media and by means of a simple hyperbolic function for fracture flow. Approximation by Taylor expansion is used to evaluate the fluid velocity for small pressure gradients to reduce rounding errors. We report simulations of flow in rough-walled fractures for different rheologies and study the effect of fluid parameters on the flow channelization in rough-walled fractures. This effect is known from previous studies. We demonstrate how rheologies on different domains can be included in a fully-unstructured reservoir simulation that incorporates discrete fracture modeling (DFM). The above formulation was implemented in the open-source MATLAB Reservoir Simulation Toolbox (MRST), which uses fully implicit discretization on general polyhedral grids, including industry standard grids with DFM. This robust implementation is an important step towards hydro-mechanically coupled simulation of hydraulic fracturing with realistic non-Newtonian fluid rheology since most hydraulic fracturing models implemented so far make use of oversimplified rheological models (e.g., Newtonian or pure power-law).

Virtual element method for geomechanical simulations of reservoir models

In this paper, we study the use of virtual element method (VEM) for geomechanics. Our emphasis is on applications to reservoir simulations. The physical processes behind the formation of the reservoirs, such as sedimentation, erosion, and faulting, lead to complex geometrical structures. A minimal representation, with respect to the physical parameters of the system, then naturally leads to general polyhedral grids. Numerical methods which can directly handle this representation will be highly favorable, in particular in the setting of advanced work-flows. The virtual element method is a promising candidate to solve the linear elasticity equations on such models. In this paper, we investigate some of the limits of the VEM method when used on reservoir models. First, we demonstrate that care must be taken to make the method robust for highly elongated cells, which is common in these applications, and show the importance of calculating forces in terms of traction on the boundary of the elements for elongated distorted cells. Second, we study the effect of triangulations on the surfaces of curved faces, which also naturally occur in subsurface models. We also demonstrate how a more stable stabilization term for reservoir application can be derived.

Stochastic finite integration technique for magnetostatic problems

AbstractBoth nonintrusive and intrusive stochastic approaches on the basis of polynomial chaos expansion are presented for the finite integration technique over generic polyhedral grids for 3D magnetostatic linear problems. Such algorithms outperform Monte Carlo methods (when the number of random parameters is small), both in accuracy and efficiency. A novel algorithm for the intrusive approach is also provided, by which the intrusive approach becomes less computationally expensive than the nonintrusive approach. Validation is performed by solving a magnetic circuit where the reluctivity is uncertain.

A Multiscale Restriction-Smoothed Basis Method for Compressible Black-Oil Models

Summary Simulation problems encountered in reservoir management are often computationally expensive because of the complex fluid physics for multiphase flow and the large number of grid cells required to honor geological heterogeneity. Multiscale methods have been proposed as a computationally inexpensive alternative to traditional fine-scale solvers for computing conservative approximations of the pressure and velocity fields on high-resolution geocellular models. Although a wide variety of such multiscale methods have been discussed in the literature, these methods have not yet seen widespread use in industry. One reason may be that no method has been presented so far that handles the combination of realistic flow physics and industry-standard grid formats in their full complexity. Herein, we present a multiscale method that handles both the most widespread type of flow physics (black-oil-type models) and standard grid formats such as corner-point, stair-stepped, and perpendicular bisector (PEBI), as well as general unstructured, polyhedral grids. Our approach is derived from a finite-volume formulation in which the basis functions are constructed by use of restricted smoothing to effectively capture the local features of the permeability. The method can also be formulated easily for other types of flow models, provided that one has a reliable (iterative) solution strategy that computes flow and transport in separate steps. The proposed method is implemented as open-source software and validated on a number of two- and three-phase test cases with significant compressibility and gas dissolution. The test cases include both synthetic models and models of real fields with complex wells, faults, and inactive and degenerate cells. Through a prescribed tolerance, the solver can be set to either converge to a sequential solution or the fully implicit solution, in both cases with a significant speedup compared with a fine-scale multigrid solver. Altogether, this ensures that one can easily and systematically trade accuracy for efficiency, or vice versa.

Computational framework for Launch, Ascent, and Vehicle Aerodynamics (LAVA)

The Launch Ascent and Vehicle Aerodynamics (LAVA) framework, developed at NASA Ames Research Center, is introduced. This technology originated from addressing some of the key challenges that were present during the re-design of the launch infrastructure at Kennedy Space Center. The solver combines Computational Fluid Dynamics (CFD) capabilities with auxiliary modules, such as Conjugate Heat Transfer (CHT) and Computational Aero-Acoustics (CAA). LAVA is designed to be grid agnostic, i.e., it can handle block-structured Cartesian, generalized curvilinear overset and unstructured polyhedral grids either as stand-alone mode or by coupling different grid types through an overset interface. A description of the spatial discretizations utilized for each grid type, along with the available explicit and implicit time-stepping schemes, is provided. The overset grid coupling procedure for Cartesian and unstructured mesh types, as well as the CHT and CAA capabilities is discussed in some detail. Several NASA mission related applications are also presented to demonstrate the capabilities for large-scale applications, such as pressure, thermal and acoustic analyses of the geometrically complex launch environment, steady and unsteady ascent aerodynamics, plume-induced flow separation analyses of heavy lift launch vehicles and aeroacoustic applications. In addition, two validation cases related to NASA aeronautics applications are presented: the 1st AIAA Sonic Boom Prediction Workshop test cases and a computational study of slat noise.

Stochastic Finite Integration Technique for Eddy-Current Problems

A stochastic formulation of the finite integration technique over generic polyhedral grids for 3-D eddy-current problems is derived using polynomial chaos expansion. The formulation is then applied to a test case representing an ac resistance welding system. Validation of the formulation is carried out by comparison with Monte Carlo simulations.

Quasi-DNS capabilities of OpenFOAM for different mesh types

Experimental limitations for certain nuclear reactor safety applications have pushed forward the demand for high fidelity DNS reference solutions for complex geometric configurations such as a T-junction or a spherical pebble bed. The application of traditional high-order DNS codes is limited to simple flow domains such as a periodic box or channel. As a possibility to create reference DNS solutions for more complex geometries, we have assessed the (quasi-)DNS capabilities of the OpenFOAM finite volume CFD solver for both structured hexahedral meshes and arbitrary polyhedral meshes. The feasibility of (quasi-)DNS analyses on polyhedral grids is of main interest, since this may offer the possibility to significantly expand the availability of (quasi-)DNS-quality data on arbitrarily complex geometries. In order to have a basis for the considered assessment, the mutual differences between generally recognized reference DNS data bases for turbulent channel and pipe flows are determined first. Subsequently, the differences between these reference DNS solutions and the present OpenFOAM (quasi-)DNS solutions are quantified for the considered mesh types. We use an existing finite volume CFD method and well known turbulent channel and pipe flow DNS reference cases for the assessments in this paper. New in this paper are the application of this CFD method to (quasi-)DNS analyses using arbitrary polyhedral meshes, and the quantification of respectively the mutual differences between generally recognized reference DNS data bases and the differences between the obtained OpenFOAM (quasi-)DNS data and these reference DNS data bases.Based on the presented assessment, it is observed that the differences between the OpenFOAM solutions and the considered reference DNS solutions are practically the same as the mutual differences between these reference DNS solutions when structured hexahedral meshes are used. Furthermore, it is observed that the differences as obtained by OpenFOAM on extruded polyhedral meshes are practically the same as those obtained for the structured hexahedral meshes. In contrast, the full polyhedral mesh shows somewhat larger differences near the peaks in the rms velocity profiles, whereas the differences in the bulk flow are again practically the same as those for the hexahedral grids.

A FIT Formulation of Bianisotropic Materials Over Polyhedral Grids

Modeling bianisotropic constitutive equations, i.e., magnetoelectric coupling, in electromagnetics simulation is increasingly important, in particular in metamaterials applications. This paper introduces for the first time such constitutive relationships in the framework of finite integration technique and, furthermore, does so by allowing full generality in the discretization through arbitrary polyhedral grids. The resulting formulation is consistent, stable, and preserves the thermodynamic properties of the bianisotropic constitutive equations because of the energetic approach used to construct the interpolating functions.

A family of linearity-preserving schemes for anisotropic diffusion problems on arbitrary polyhedral grids

A family of cell-centered finite volume schemes are proposed for anisotropic diffusion problems on arbitrary polyhedral grids with planar facets. The derivation of the schemes is done under a general framework through a certain linearity-preserving approach. The key ingredient of our algorithm is to employ solely the so-called harmonic averaging points located at the cell interfaces to define the auxiliary unknowns, which not only makes the interpolation procedure for auxiliary unknowns simple and positivity-preserving, but also reduces the stencil of the schemes. The final schemes are cell-centered with a small stencil of 25-point on the structured hexahedral grids. Moreover, the schemes satisfy the local conservation condition, treat discontinuity exactly and allow for a simple stability analysis. A second-order accuracy in the L2 norm and a first-order accuracy in the H1 norm are observed numerically on general distorted meshes in case that the diffusion tensor is anisotropic and discontinuous.

Discontinuous Galerkin and mimetic finite difference methods for coupled Stokes–Darcy flows on polygonal and polyhedral grids

We study locally mass conservative approximations of coupled Darcy and Stokes flows on polygonal and polyhedral meshes. The discontinuous Galerkin (DG) finite element method is used in the Stokes region and the mimetic finite difference method is used in the Darcy region. DG finite element spaces are defined on polygonal and polyhedral grids by introducing lifting operators mapping mimetic degrees of freedom to functional spaces. Optimal convergence estimates for the numerical scheme are derived. Results from computational experiments supporting the theory are presented.

An accurate and efficient method for automatic deformation of unstructured polyhedral grids to simulate the flow induced by the motion of solid objects

In the present study, an efficient numerical method to deform an unstructured polyhedral grid for accurate simulations of a flow induced by a moving solid object is proposed. In order to overcome inefficiency of the existing methods based on the vertices of a computational cell, the present approach deforms the polyhedral grids adaptively using a network based on the centre of each computational cell. It is shown that the proposed method provides with an improved efficiency for an arbitrary Lagrangian-Eulerian (ALE) simulation of a flow with a moving solid object compared to the existing methods based on the vertices. It is also shown that a polyhedral grid results in a better numerical accuracy than a triangular grid. The present numerical methods were applied to flows around a fixed and oscillating circular cylinder with various Reynolds numbers and oscillating frequencies. The produced numerical results were verified against aerodynamic characteristics of the previous numerical studies.

Two-phase water flooding simulations on dynamic adaptive octree grids with two-point nonlinear fluxes

We present a method for numerical simulation of the two-phase water flooding problem on general polyhedral grids not aligned with permeability tensor K (K-nonorthogonal grids) and dynamic octree grids adapted to the front between the phases. The discretization is based on the cell-centered monotone finite volume (FV) method with the nonlinear two-point flux approximation (TPFA) applicable to general K-non-orthogonal polyhedral grids. We use fully implicit discretization in time to avoid the restriction on the time step caused by the minimal mesh size. In our numerical experiments we demonstrate the superiority of the nonlinear TPFA on a K-non-orthogonal grid over linear TPFA and considerable speed-up of the simulation on a dynamically adapted octree grid with minimal loss in accuracy compared to the simulation on a fine regular grid.

ParNCL and ParGAL: Data-parallel Tools for Postprocessing of Large-scale Earth Science Data

Earth science high-performance applications often require extensive analysis of their output in order to complete the scien- tific goals or produce a visual image or animation. Often this analysis cannot be done in situ because it requires calculating time-series statistics from state sampled over the entire length of the run or analyzing the relationship between similar time series from previous simulations or observations. Many of the tools used for this postprocessing are not themselves high- performance applications, but the new Parallel Gridded Analysis Library (ParGAL) provides high-performance data-parallel versions of several common analysis algorithms for data from a structured or unstructured grid simulation. The library builds on several scalable systems, including the Mesh Oriented DataBase (MOAB), a library for representing mesh data that sup- ports structured, unstructured finite element, and polyhedral grids; the Parallel-NetCDF (PNetCDF) library; and Intrepid, an extensible library for computing operators (such as gradient, curl, and divergence) acting on discretized fields. We have used ParGAL to implement a parallel version of the NCAR Command Language (NCL) a scripting language widely used in the climate community for analysis and visualization. The data-parallel algorithms in ParGAL/ParNCL are both higher performing and more flexible than their serial counterparts.

An efficient algorithm for mapping imaging data to 3D unstructured grids in computational biomechanics

Geometries for organ scale and multiscale simulations of organ function are now routinely derived from imaging data. However, medical images may also contain spatially heterogeneous information other than geometry that are relevant to such simulations either as initial conditions or in the form of model parameters. In this manuscript, we present an algorithm for the efficient and robust mapping of such data to imaging-based unstructured polyhedral grids in parallel. We then illustrate the application of our mapping algorithm to three different mapping problems: (i) the mapping of MRI diffusion tensor data to an unstructured ventricular grid; (ii) the mapping of serial cyrosection histology data to an unstructured mouse brain grid; and (iii) the mapping of computed tomography-derived volumetric strain data to an unstructured multiscale lung grid. Execution times and parallel performance are reported for each case.

Accurate and robust CFD algorithms applied to 3D arbitrary polyhedral grids

A parallel, accurate, robust and grid-transparent CFD solver was developed to solve compressible flow on 3D arbitrary polyhedral grids. To improve spacial accurate, a constrained least-squares reconstruction method is developed. To accelerate convergence, a matrix free implicit method GMRES+LU-SGS, and a parallel method using OpenMP are presented on shared-memory parallel systems with help of a special grid reordering technique. Several typical test cases, including subsonic, transonic and hypersonic flows, prove that the resulting solver performs fast, accurate and robust.

Interactive Volume Visualization of General Polyhedral Grids

This paper presents a novel framework for visualizing volumetric data specified on complex polyhedral grids, without the need to perform any kind of a priori tetrahedralization. These grids are composed of polyhedra that often are non-convex and have an arbitrary number of faces, where the faces can be non-planar with an arbitrary number of vertices. The importance of such grids in state-of-the-art simulation packages is increasing rapidly. We propose a very compact, face-based data structure for representing such meshes for visualization, called two-sided face sequence lists (TSFSL), as well as an algorithm for direct GPU-based ray-casting using this representation. The TSFSL data structure is able to represent the entire mesh topology in a 1D TSFSL data array of face records, which facilitates the use of efficient 1D texture accesses for visualization. In order to scale to large data sizes, we employ a mesh decomposition into bricks that can be handled independently, where each brick is then composed of its own TSFSL array. This bricking enables memory savings and performance improvements for large meshes. We illustrate the feasibility of our approach with real-world application results, by visualizing highly complex polyhedral data from commercial state-of-the-art simulation packages.

Velocity interpolation and streamline tracing on irregular geometries

Streamline tracing on irregular grids requires reliable interpolation of velocity fields. We propose a new method for direct streamline tracing on polygon and polytope cells. While some numerical methods provide a basis function that can be used for interpolation, other methods provide only the fluxes at the faces of the elements. We introduce the concept of full- and raw-field methods. Full-field methods have built-in interpolation but are often not defined on general grids such as polygonal and polyhedral grids which we examine here. Also, reliability issues may arise on non-simplicial meshes in terms of not being able to reproduce constant velocity fields. We propose an interpolation in H(div) and H(curl) valid on general grids that is based on barycentric coordinates and that reproduces uniform flow. The interpolation can be used to compute the streamline directly on the complex cell geometry. The method generalizes to convex polytopes in 3D, with a restriction on the polytope topology near corners that is shown to be satisfied by several popular grid types. Numerical results confirm that the method is applicable to general grids and preserves uniform flow.

On the flexibility of agglomeration based physical space discontinuous Galerkin discretizations

In this work we show that the flexibility of the discontinuous Galerkin (dG) discretization can be fruitfully exploited to implement numerical solution strategies based on the use of elements with very general shapes. Thanks to the freedom in defining the mesh topology, we propose a new h-adaptive technique based on agglomeration coarsening of a fine mesh. The possibility to enhance the error distribution over the computational domain is investigated on a Poisson problem with the goal of obtaining a mesh independent discretization. The main building block of our dG method consists of defining discrete polynomial spaces directly on physical frame elements. For this purpose we orthonormalize with respect to the L 2-product a set of monomials relocated in a specific element frame and we introduce an easy way to reduce the cost related to numerical integration on agglomerated meshes. To complete the dG formulation for second order problems, two extensions of the BR2 scheme to arbitrary polyhedral grids, including an estimate of the stabilization parameter ensuring the coercivity property, are here proposed.

Open-source MATLAB implementation of consistent discretisations on complex grids

Accurate geological modelling of features such as faults, fractures or erosion requires grids that are flexible with respect to geometry. Such grids generally contain polyhedral cells and complex grid-cell connectivities. The grid representation for polyhedral grids in turn affects the efficient implementation of numerical methods for subsurface flow simulations. It is well known that conventional two-point flux-approximation methods are only consistent for K-orthogonal grids and will, therefore, not converge in the general case. In recent years, there has been significant research into consistent and convergent methods, including mixed, multipoint and mimetic discretisation methods. Likewise, the so-called multiscale methods based upon hierarchically coarsened grids have received a lot of attention. The paper does not propose novel mathematical methods but instead presents an open-source Matlab® toolkit that can be used as an efficient test platform for (new) discretisation and solution methods in reservoir simulation. The aim of the toolkit is to support reproducible research and simplify the development, verification and validation and testing and comparison of new discretisation and solution methods on general unstructured grids, including in particular corner point and 2.5D PEBI grids. The toolkit consists of a set of data structures and routines for creating, manipulating and visualising petrophysical data, fluid models and (unstructured) grids, including support for industry standard input formats, as well as routines for computing single and multiphase (incompressible) flow. We review key features of the toolkit and discuss a generic mimetic formulation that includes many known discretisation methods, including both the standard two-point method as well as consistent and convergent multipoint and mimetic methods. Apart from the core routines and data structures, the toolkit contains add-on modules that implement more advanced solvers and functionality. Herein, we show examples of multiscale methods and adjoint methods for use in optimisation of rates and placement of wells.

Topics in structure-preserving discretization

In the last few decades the concepts of structure-preserving discretization, geometric integration and compatible discretizations have emerged as subfields in the numerical approximation of ordinary and partial differential equations. The article discusses certain selected topics within these areas; discretization techniques both in space and time are considered. Lie group integrators are discussed with particular focus on the application to partial differential equations, followed by a discussion of how time integrators can be designed to preserve first integrals in the differential equation using discrete gradients and discrete variational derivatives.Lie group integrators depend crucially on fast and structure-preserving algorithms for computing matrix exponentials. Preservation of domain symmetries is of particular interest in the application of Lie group integrators to PDEs. The equivariance of linear operators and Fourier transforms on non-commutative groups is used to construct fast structure-preserving algorithms for computing exponentials. The theory of Weyl groups is employed in the construction of high-order spectral element discretizations, based on multivariate Chebyshev polynomials on triangles, simplexes and simplicial complexes.The theory of mixed finite elements is developed in terms of special inverse systems of complexes of differential forms, where the inclusion of cells corresponds to pullback of forms. The theory covers, for instance, composite piecewise polynomial finite elements of variable order over polyhedral grids. Under natural algebraic and metric conditions, interpolators and smoothers are constructed, which commute with the exterior derivative and whose product is uniformly stable in Lebesgue spaces. As a consequence we obtain not only eigenpair approximation for the Hodge–Laplacian in mixed form, but also variants of Sobolev injections and translation estimates adapted to variational discretizations.