Improving the robustness of the control volume finite element method with application to multiphase porous media flow

Control volume finite element methods (CVFEMs) have been proposed to simulate flow in heterogeneous porous media because they are better able to capture complex geometries using unstructured meshes. However, producing good quality meshes in such models is nontrivial and may sometimes be impossible, especially when all or parts of the domains have very large aspect ratio. A novel CVFEM is proposed here that uses a control volume representation for pressure and yields significant improvements in the quality of the pressure matrix. The method is initially evaluated and then applied to a series of test cases using unstructured (triangular/tetrahedral) meshes, and numerical results are in good agreement with semianalytically obtained solutions. The convergence of the pressure matrix is then studied using complex, heterogeneous example problems. The results demonstrate that the new formulation yields a pressure matrix than can be solved efficiently even on highly distorted, tetrahedral meshes in models of heterogeneous porous media with large permeability contrasts. The new approach allows effective application of CVFEM in such models.     

Monotone nonlinear finite‐volume method for nonisothermal two‐phase two‐component flow in porous media

This article presents a new nonlinear finite‐volume scheme for the nonisothermal two‐phase two‐component flow equations in porous media. The face fluxes are approximated by a nonlinear two‐point flux approximation, where transmissibilities nonlinearly depend on primary variables. Thereby, we mainly follow the ideas proposed by Le Potier combined with a harmonic averaging point interpolation strategy for the approximation of arbitrary heterogeneous permeability fields on polygonal grids. The behavior of this interpolation strategy is analyzed, and its limitation for highly anisotropic permeability tensors is demonstrated. Moreover, the condition numbers of occurring matrices are compared with linear finite‐volume schemes. Additionally, the convergence behavior of iterative solvers is investigated. Finally, it is shown that the nonlinear scheme is more efficient than its linear counterpart. Copyright © 2016 John Wiley & Sons, Ltd.     

Multiphase flow in heterogeneous porous media: A classical finite element method versus an implicit pressure–explicit saturation‐based mixed finite element–finite volume approach

Various discretization methods exist for the numerical simulation of multiphase flow in porous media. In this paper, two methods are introduced and analyzed—a full‐upwind Galerkin method which belongs to the classical finite element methods, and a mixed‐hybrid finite element method based on an implicit pressure–explicit saturation (IMPES) approach. Both methods are derived from the governing equations of two‐phase flow. Their discretization concepts are compared in detail. Their efficiency is discussed using several examples. Copyright © 1999 John Wiley & Sons, Ltd.     

A control volume finite element method for three-dimensional three-phase flows

A novel control volume finite element method with adaptive anisotropic unstructured meshes is presented for three-dimensional three-phase flows with interfacial tension. The numerical framework consists of a mixed control volume and finite element formulation with a new P₁DG-P₂ elements (linear discontinuous velocity between elements and quadratic continuous pressure between elements). A “volume of fluid” type method is used for the interface capturing, which is based on compressive control volume advection and second-order finite element methods. A force-balanced continuum surface force model is employed for the interfacial tension on unstructured meshes. The interfacial tension coefficient decomposition method is also used to deal with interfacial tension pairings between different phases. Numerical examples of benchmark tests and the dynamics of three-dimensional three-phase rising bubble, and droplet impact are presented. The results are compared with the analytical solutions and previously published experimental data, demonstrating the capability of the present method.     

An improved unstructured WENO method for compressible multi‐fluids flow

An improved high‐order accurate WENO finite volume method based on unstructured grids for compressible multi‐fluids flow is proposed in this paper. The third‐order accuracy WENO finite volume method based on triangle cell is used to discretize the governing equations. To have higher order of accuracy, the P1 polynomial is reconstructed firstly. After that, the P2 polynomial is reconstructed from the combination of the P1. The reconstructed coefficients are calculated by analytical form of inverse matrix rather than the numerical inversion. This greatly improved the efficiency and the robustness. Four examples are presented to examine this algorithm. Numerical results show that there is no spurious oscillation of velocity and pressure across the interface and high‐order accurate result can be achieved. Copyright © 2016 John Wiley & Sons, Ltd.     

A higher resolution edge‐based finite volume method for the simulation of the oil–water displacement in heterogeneous and anisotropic porous media using a modified IMPES method

In this article, we present a higher‐order finite volume method with a ‘Modified Implicit Pressure Explicit Saturation’ (MIMPES) formulation to model the 2D incompressible and immiscible two‐phase flow of oil and water in heterogeneous and anisotropic porous media. We used a median‐dual vertex‐centered finite volume method with an edge‐based data structure to discretize both, the elliptic pressure and the hyperbolic saturation equations. In the classical IMPES approach, first, the pressure equation is solved implicitly from an initial saturation distribution; then, the velocity field is computed explicitly from the pressure field, and finally, the saturation equation is solved explicitly. This saturation field is then used to re‐compute the pressure field, and the process follows until the end of the simulation is reached. Because of the explicit solution of the saturation equation, severe time restrictions are imposed on the simulation. In order to circumvent this problem, an edge‐based implementation of the MIMPES method of Hurtado and co‐workers was developed. In the MIMPES approach, the pressure equation is solved, and the velocity field is computed less frequently than the saturation field, using the fact that, usually, the velocity field varies slowly throughout the simulation. The solution of the pressure equation is performed using a modification of Crumpton's two‐step approach, which was designed to handle material discontinuity properly. The saturation equation is solved explicitly using an edge‐based implementation of a modified second‐order monotonic upstream scheme for conservation laws type method. Some examples are presented in order to validate the proposed formulation. Our results match quite well with others found in literature. Copyright © 2016 John Wiley & Sons, Ltd.     

流体饱和两相多孔介质拟静态问题的混合有限元方法

针对基于混合物理论的两相多孔介质模型,采用Galerkin加权残值有限元法,导出求解所静态问题的基于us-uF-P变量的混合有限元方程,由于系统方程的系数矩阵非定,进而针对该方程组提出了一种失代求解方法,并由分片试验得出节点压力插值函数的阶须低于固体相节点的位移插值函数的阶的结论,算例结果表明,采用基于u2-uF-p变量的混合法计算所得的固体相和流体相速度以及固体相的有效应力与罚方法一致,而压力值的粗度高于罚方法。     

A CIP/multi‐moment finite volume method for shallow water equations with source terms

A novel finite volume method has been presented to solve the shallow water equations. In addition to the volume‐integrated average (VIA) for each mesh cell, the surface‐integrated average (SIA) is also treated as the model variable and is independently predicted. The numerical reconstruction is conducted based on both the VIA and the SIA. Different approaches are used to update VIA and SIA separately. The SIA is updated by a semi‐Lagrangian scheme in terms of the Riemann invariants of the shallow water equations, while the VIA is computed by a flux‐based finite volume formulation and is thus exactly conserved. Numerical oscillation can be effectively avoided through the use of a non‐oscillatory interpolation function. The numerical formulations for both SIA and VIA moments maintain exactly the balance between the fluxes and the source terms. 1D and 2D numerical formulations are validated with numerical experiments. Copyright © 2007 John Wiley & Sons, Ltd.     

A discontinuous Galerkin scheme for the numerical solution of flow problems with discontinuities

In this paper, a new discontinuous Galerkin finite element method for the numerical solution of flow problems with discontinuities is presented. The method is based on the limitation in every cell of the difference between the extrema values and the mean value of the numerical solution. The algorithm and technical details for the implementation of the method are presented in one‐and two‐dimensional problems. Numerical experiments for classical test problems are solved on unstructured triangulations to demonstrate the performance of the proposed method. Copyright © 2011 John Wiley & Sons, Ltd.     

A semi‐Lagrangian multi‐moment finite volume method with fourth‐order WENO projection

Numerical oscillation has been an open problem for high‐order numerical methods with increased local degrees of freedom (DOFs). Current strategies mainly follow the limiting projections derived originally for conventional finite volume methods and thus are not able to make full use of the sub‐cell information available in the local high‐order reconstructions. This paper presents a novel algorithm that introduces a nodal value‐based weighted essentially non‐oscillatory limiter for constrained interpolation profile/multi‐moment finite volume method (CIP/MM FVM) (Ii and Xiao, J. Comput. Phys., 222 (2007), 849–871) as an effort to pursue a better suited formulation to implement the limiting projection in schemes with local DOFs. The new scheme, CIP‐CSL‐WENO4 scheme, extends the CIP/MM FVM method by limiting the slope constraint in the interpolation function using the weighted essentially non‐oscillatory (WENO) reconstruction that makes use of the sub‐cell information available from the local DOFs and is built from the point values at the solution points within three neighboring cells, thus resulting a more compact WENO stencil. The proposed WENO limiter matches well the original CIP/MM FVM, which leads to a new scheme of high accuracy, algorithmic simplicity, and computational efficiency. We present the numerical results of benchmark tests for both scalar and Euler conservation laws to manifest the fourth‐order accuracy and oscillation‐suppressing property of the proposed scheme. Copyright © 2016 John Wiley & Sons, Ltd.     

A mass‐conserving Level‐Set method for modelling of multi‐phase flows

A mass‐conserving Level‐Set method to model bubbly flows is presented. The method can handle high density‐ratio flows with complex interface topologies, such as flows with simultaneous occurrence of bubbles and droplets. Aspects taken into account are: a sharp front (density changes abruptly), arbitrarily shaped interfaces, surface tension, buoyancy and coalescence of droplets/bubbles. Attention is paid to mass‐conservation and integrity of the interface. The proposed computational method is a Level‐Set method, where a Volume‐of‐Fluid function is used to conserve mass when the interface is advected. The aim of the method is to combine the advantages of the Level‐Set and Volume‐of‐Fluid methods without the disadvantages. The flow is computed with a pressure correction method with the Marker‐and‐Cell scheme. Interface conditions are satisfied by means of the continuous surface force methodology and the jump in the density field is maintained similar to the ghost fluid method for incompressible flows. Copyright © 2005 John Wiley & Sons, Ltd.     

Improving Eulerian two‐phase flow finite element approximation with discontinuous gradient pressure shape functions

In this paper we present a problem we have encountered using a stabilized finite element method on fixed grids for flows with interfaces modelled with the level set approach. We propose a solution based on enriching the pressure shape functions on the elements cut by the interface. The enrichment is used to enable the pressure gradient to be discontinuous at the interface, thus improving the ability to simulate the behaviour of fluids with different density under a gravitational force. The additional shape function used is local to each element and the corresponding degree of freedom can therefore be condensed prior to assembly, making the implementation quite simple on any existing finite element code. Copyright © 2005 John Wiley & Sons, Ltd.     

Well‐balanced discontinuous Galerkin method and finite volume WENO scheme based on hydrostatic reconstruction for blood flow model in arteries

The blood flow model in arteries admits the steady state solutions, for which the flux gradient is nonzero, and is exactly balanced by the source term. In this paper, by means of hydrostatic reconstruction, we construct a high order discontinuous Galerkin method, which exactly preserves the dead‐man steady state, which is characterized by a discharge equal to zero (analogue to hydrostatic equilibrium). Moreover, the method maintains genuine high order of accuracy. Subsequently, we apply the key idea to finite volume weighted essentially non‐oscillatory schemes and obtain a well‐balanced finite volume weighted essentially non‐oscillatory scheme. Extensive numerical experiments are performed to verify the well‐balanced property, high order accuracy, as well as good resolution for smooth and discontinuous solutions.     

Thermo‐hydrodynamic‐mechanical multiphase flow model for CO2 sequestration in fracturing porous media

In this paper, we introduce a fully coupled thermo‐hydrodynamic‐mechanical computational model for multiphase flow in a deformable porous solid, exhibiting crack propagation due to fluid dynamics, with focus on CO₂ geosequestration. The geometry is described by a matrix domain, a fracture domain, and a matrix‐fracture domain. The fluid flow in the matrix domain is governed by Darcy's law and that in the crack is governed by the Navier–Stokes equations. At the matrix‐fracture domain, the fluid flow is governed by a leakage term derived from Darcy's law. Upon crack propagation, the conservation of mass and energy of the crack fluid is constrained by the isentropic process. We utilize the representative elementary volume‐averaging theory to formulate the mathematical model of the porous matrix, and the drift flux model to formulate the fluid dynamics in the fracture. The numerical solution is conducted using a mixed finite element discretization scheme. The standard Galerkin finite element method is utilized to discretize the diffusive dominant field equations, and the extended finite element method is utilized to discretize the crack propagation, and the fluid leakage at the boundaries between layers of different physical properties. A numerical example is given to demonstrate the computational capability of the model. It shows that the model, despite the relatively large number of degrees of freedom of different physical nature per node, is computationally efficient, and geometry and effectively mesh independent. Copyright © 2017 John Wiley & Sons, Ltd.     

Finite‐volume multi‐stage schemes for shallow‐water flow simulations

A finite‐volume multi‐stage (FMUSTA) scheme is proposed for simulating the free‐surface shallow‐water flows with the hydraulic shocks. On the basis of the multi‐stage (MUSTA) method, the original Riemann problem is transformed to an independent MUSTA mesh. The local Lax–Friedrichs scheme is then adopted for solving the solution of the Riemann problem at the cell interface on the MUSTA mesh. The resulting first‐order monotonic FMUSTA scheme, which does not require the use of the eigenstructure and the special treatment of entropy fixes, has the generality as well as simplicity. In order to achieve the high‐resolution property, the monotonic upstream schemes for conservation laws (MUSCL) method are used. For modeling shallow‐water flows with source terms, the surface gradient method (SGM) is adopted. The proposed schemes are verified using the simulations of six shallow‐water problems, including the 1D idealized dam breaking, the steady transcritical flow over a hump, the 2D oblique hydraulic jump, the circular dam breaking and two dam‐break experiments. The simulated results by the proposed schemes are in satisfactory agreement with the exact solutions and experimental data. It is demonstrated that the proposed FMUSTA schemes have superior overall numerical accuracy among the schemes tested such as the commonly adopted Roe and HLL schemes. Copyright © 2007 John Wiley & Sons, Ltd.     

可压缩流体流动多尺度混合有限元数值方法研究

可压缩流体是天然油藏中广泛存在的一种流体,研究其在多孔介质中的渗流规律对于油藏开发具有重要意义。本文采用多尺度混合有限元方法,对可压缩流体渗流问题进行了研究。考虑流体的可压缩性以及介质形变,推导得到了可压缩流体渗流问题的多尺度计算格式。数值计算结果表明,多尺度混合有限元适于求解非均质性和可压缩流问题,具有节省计算量、计算精度高等优势,对于实际大规模油藏模拟具有重要意义。     

Geometrical interpretation of the multi‐point flux approximation L‐method

In this paper, we first investigate the influence of different Dirichlet boundary discretizations on the convergence rate of the multi‐point flux approximation (MPFA) L‐method by the numerical comparisons between the MPFA O‐ and L‐method, and show how important it is for this new method to handle Dirichlet boundary conditions in a suitable way. A new Dirichlet boundary strategy is proposed, which in some sense can well recover the superconvergence rate of the normal velocity. In the second part of the work, the MPFA L‐method with homogeneous media is studied. A systematic concept and geometrical interpretations of the L‐method are given and illustrated, which yield more insight into the L‐method. Finally, we apply the MPFA L‐method for two‐phase flow in porous media on different quadrilateral grids and compare its numerical results for the pressure and saturation with the results of the two‐point flux approximation method. Copyright © 2008 John Wiley & Sons, Ltd.     

A finite element thermally coupled flow formulation for phase‐change problems

A finite element, thermally coupled incompressible flow formulation considering phase‐change effects is presented. This formulation accounts for natural convection, temperature‐dependent material properties and isothermal and non‐isothermal phase‐change models. In this context, the full Navier–Stokes equations are solved using a generalized streamline operator (GSO) technique. The highly non‐linear phase‐change effects are treated with a temperature‐based algorithm, which provides stability and convergence of the numerical solution. The Boussinesq approximation is used in order to consider the temperature‐dependent density variation. Furthermore, the numerical solution of the coupled problem is approached with a staggered incremental‐iterative solution scheme, such that the convergence criteria are written in terms of the residual vectors. Finally, this formulation is used for the solutions of solidification and melting problems validating some numerical results with other existing solutions obtained with different methodologies. Copyright © 2000 John Wiley & Sons, Ltd.     

A hybrid FVM–LBM method for single and multi‐fluid compressible flow problems

The lattice Boltzmann method (LBM) has established itself as an alternative approach to solve the fluid flow equations. In this work we combine LBM with the conventional finite volume method (FVM), and propose a non‐iterative hybrid method for the simulation of compressible flows. LBM is used to calculate the inter‐cell face fluxes and FVM is used to calculate the node parameters. The hybrid method is benchmarked for several one‐dimensional and two‐dimensional test cases. The results obtained by the hybrid method show a steeper and more accurate shock profile as compared with the results obtained by the widely used Godunov scheme or by a representative flux vector splitting scheme. Additional features of the proposed scheme are that it can be implemented on a non‐uniform grid, study of multi‐fluid problems is possible, and it is easily extendable to multi‐dimensions. These features have been demonstrated in this work. The proposed method is therefore robust and can possibly be applied to a variety of compressible flow situations. Copyright © 2009 John Wiley & Sons, Ltd.     

A higher‐order unsplit 2D direct Eulerian finite volume method for two‐material compressible flows based on the MOOD paradigms

A higher‐order unsplit multi‐dimensional discretization of the diffuse interface model for two‐material compressible flows proposed by R. Saurel, F. Petitpas and R. A. Berry in 2009 is developed. The proposed higher‐order method is based on the concepts of the Multidimensional Optimal Order Detection (MOOD) method introduced in three recent papers for single‐material flows. The first‐order unsplit multi‐dimensional Finite Volume discretization presented by SPB serves as foundation for the development of the higher‐order unlimited schemes. Specific detection criteria along with a novel decrementing algorithm for the MOOD method are designed in order to deal with the complexity of multi‐material flows. Numerically, we compare errors and computational times on several 1D problems (stringent shock tube and cavitation problems) computed on 2D meshes with the second‐ and fourth‐order MOOD methods using a classical MUSCL method as reference. Several simulations of a 2D shocked R22 bubble in the air are also presented on Cartesian and unstructured meshes with the second‐ and fourth‐order MOOD methods, and qualitative comparisons confirm the conclusions obtained with 1D problems. These numerical results demonstrate the robustness of the MOOD approach and the interest of using more than second‐order methods even for locally singular solutions of complex physics models. Copyright © 2014 John Wiley & Sons, Ltd.