Continuous Galerkin Finite Element Method Research Articles

Numerical simulation of CO2-enhanced oil recovery in fractured shale reservoirs using discontinuous and continuous Galerkin finite element methods

Introduction: This study explores the potential of enhancing shale oil recovery and reducing CO2 emissions through CO2 injection in fractured shale reservoirs. The importance of this approach lies in its dual benefit: improving oil extraction efficiency and addressing environmental concerns associated with CO2 emissions.Method: We employed a discrete fracture-matrix model to simulate CO2 flooding in fractured shale reservoirs, utilizing both discontinuous Galerkin (DG) and continuous Galerkin (CG) finite element methods. The DG-CG FEM’s accuracy was validated against the McWhorter problem, ensuring the reliability of the simulation results. Our model also considered various factors, including reservoir heterogeneity, fracture permeability, CO2 injection volume, and gas injection patterns, to analyze their impact on shale oil recovery.Result: Our simulations revealed that fractured reservoirs significantly enhance shale oil production efficiency compared to homogeneous reservoirs, with an approximate 48.9% increase in production. A notable increase in shale oil production, by 15.8%, was observed when fracture permeability was increased by two orders of magnitude. Additionally, a fourfold increase in CO2 injection rate resulted in a 31.5% rise in shale oil production. Implementing a step-by-step reduction in injection volume while maintaining the total CO2 injection constant proved to be more effective than constant-rate injections.Discussion: The study demonstrates the effectiveness of CO2 flooding in fractured shale reservoirs for enhancing shale oil recovery.

FEMFFUSION and its verification using the C5G7 benchmark

FEMFFUSION is an open source code that solves the multigroup neutron transport equation using the diffusion and the SPN approximations. This code uses the continuous Galerkin finite element method to be able to deal with any type of geometry, structured or unstructured, and problem dimension. To verify the program, a research study of the two-dimensional C5G7 benchmark has been performed, including the static, time-dependent and neutron noise computations. Numerical results include the effective multiplication factor, keff, the spatial neutron fluxes, and neutron power distributions for the static computations. The total neutron power evolution of a given transient and the neutron flux distribution evolution for the different energy groups. In addition, neutron noise frequency domain calculation were completed for perturbations of some cross-section in the reactor core. All computation have been performed using different SPN approximations. FEMFFUSION results are verified against the official OECD/NEA benchmark results.

Fully coupled hydro-mechanical modeling of two-phase flow in deformable fractured porous media with discontinuous and continuous Galerkin method

Discrete fracture and matrix modeling of coupled flow and geomechanics is instrumental for understanding flow in fractured media for various geoengineering applications such as enhanced geothermal energy systems and groundwater remediation. We employ the governing equations for two-phase flow in deformable fractured porous media with a stress-dependent porosity and permeability model of matrix, and variable fracture aperture and the corresponding permeability. A finite element framework is presented, in which fractures are regarded as low-dimensional objects, to discretize the coupled two-phase flow and geomechanics in fractured porous media. A hybrid method, combining discontinuous Galerkin (DG) and continuous Galerkin (CG) finite element methods (FEM), is utilized to solve for the two-phase flow, while the solid deformation is approximated using a discontinuous Galerkin FEM approach. Several benchmark cases are utilized to examine the accuracy of the proposed hybrid DG-CG FEM method. Further validation is performed using more complex, realistic fracture configurations. As the mechanical characteristics of the fracture and the surrounding matrix differ, the simulation results demonstrate that displacement and stress are discontinuous on both sides of the fracture. While the case with low permeability fractures exhibits pressure jump across the fractures, the pressure changes are reasonably smooth for the conduit fractures.

Solution of Neutron Diffusion Problems by Discontinuous Galerkin Finite Element Method With Consideration of Discontinuity Factors

Abstract Neutronics calculations are the basis of reactor analysis and design. Finite element methods (FEM) have gained increasing attention in solving neutron transport problems for the rigorous mathematical formulation and the flexibility of handling complex geometric domains and boundary conditions. In order to reduce the computational errors caused by the homogenization of cross sections, this paper adopts the discontinuous Galerkin finite element method (DG-FEM) to solve the generalized eigenvalue problem formulated by the neutron diffusion theory and compensates the homogenization error by incorporating discontinuity factors. The results show that the discontinuous Galerkin finite element method can introduce the discontinuity factors with clear mathematical and physical meanings. The computational results of the discontinuous Galerkin finite element method are slightly better than those of the continuous Galerkin finite element method. However, the computation cost of the former is higher than that of the latter. Although good parallel efficiency can be achieved, the discontinuous Galerkin finite element method is not preferable for large-scale problems unless the effect of the discontinuity factors is significant.

A projection-based, semi-implicit time-stepping approach for the Cahn-Hilliard Navier-Stokes equations on adaptive octree meshes

The Cahn-Hilliard Navier-Stokes (CHNS) system provides a computationally tractable model that can be used to effectively capture interfacial dynamics in two-phase fluid flows. In this work we present a semi-implicit, projection-based finite element framework for solving the CHNS system. We use a projection-based semi-implicit time discretization for the Navier-Stokes equation and a fully-implicit time discretization for the Cahn-Hilliard equation. We use a conforming continuous Galerkin (cG) finite element method in space equipped with a residual-based variational multiscale (RBVMS) formulation. Pressure is decoupled using a projection step, which results in two linear positive semi-definite systems for velocity and pressure, instead of the saddle point system of a pressure-stabilized method. All the linear systems are solved using an efficient and scalable algebraic multigrid (AMG) method. We deploy this approach on a massively parallel numerical implementation using parallel octree-based adaptive meshes. The overall approach allows the use of relatively large time steps with much faster time-to-solve than similar fully-implicit methods. We present comprehensive numerical experiments showing detailed comparisons with results from the literature for canonical cases, including the single bubble rise and Rayleigh-Taylor instability.

A Priori and a Posteriori Error Analysis for Generic Linear Elliptic Problems

In this paper, a priori error analysis has been examined for the continuous Galerkin finite element method which is used for solving a generic scalar and a generic system of linear elliptic equations. We derived optimal order a priori error bounds in (energy) norm utilising standard a priori error analysis techniques and tools. Also, a posteriori error analysis is investigated for a generic scalar linear elliptic equation and for a generic system of linear elliptic equations. We derived optimal residual-based a posteriori error estimates energy technique in norm

Meet the Editors

Meet the Editors

A fully-coupled framework for solving Cahn-Hilliard Navier-Stokes equations: Second-order, energy-stable numerical methods on adaptive octree based meshes

We present a fully-coupled, implicit-in-time framework for solving a thermodynamically-consistent Cahn-Hilliard Navier-Stokes system that models two-phase flows. In this work, we extend the block iterative method presented in Khanwale et al. [Simulating two-phase flows with thermodynamically consistent energy stable Cahn-Hilliard Navier-Stokes equations on parallel adaptive octree based meshes, J. Comput. Phys. (2020)], to a fully-coupled, provably second-order accurate scheme in time, while maintaining energy-stability. The new method requires fewer matrix assemblies in each Newton iteration resulting in faster solution time. The method is based on a fully-implicit Crank-Nicolson scheme in time and a pressure stabilization for an equal order Galerkin formulation. That is, we use a conforming continuous Galerkin (cG) finite element method in space equipped with a residual-based variational multiscale (RBVMS) procedure to stabilize the pressure. We deploy this approach on a massively parallel numerical implementation using parallel octree-based adaptive meshes. We present comprehensive numerical experiments showing detailed comparisons with results from the literature for canonical cases, including the single bubble rise, Rayleigh-Taylor instability, and lid-driven cavity flow problems. We analyze in detail the scaling of our numerical implementation.

Dislocation transport using a time-explicit Runge–Kutta discontinuous Galerkin finite element approach

A time-explicit Runge–Kutta discontinuous Galerkin (RKDG) finite element scheme is proposed to solve the dislocation transport initial boundary value problem in 3D. The dislocation density transport equation, which lies at the core of this problem, is a first-order unsteady-state advection–reaction-type hyperbolic partial differential equation; the DG approach is well suited to solve such equations that lack any diffusion terms. The development of the RKDG scheme follows the method of lines approach. First, a space semi-discretization is performed using the DG approach with upwinding to obtain a system of ordinary differential equations in time. Then, time discretization is performed using explicit RK schemes to solve this system. The 3D numerical implementation of the RKDG scheme is performed for the first-order (forward Euler), second-order and third-order RK methods using the strong stability preserving approach. These implementations provide (quasi-)optimal convergence rates for smooth solutions. A slope limiter is used to prevent spurious Gibbs oscillations arising from high-order space approximations (polynomial degree ⩾ 1) of rough solutions. A parametric study is performed to understand the influence of key parameters of the RKDG scheme on the stability of the solution predicted during a screw dislocation transport simulation. Then, annihilation of two oppositely signed screw dislocations and the expansion of a polygonal dislocation loop are simulated. The RKDG scheme is able to resolve the shock generated during dislocation annihilation without any spurious oscillations and predict the prismatic loop expansion with very low numerical diffusion. These results indicate that the proposed scheme is more robust and accurate in comparison to existing approaches based on the continuous Galerkin finite element method or the fast Fourier transform method.

A method of FE modeling multiphase compressible flow in hydrocarbon reservoirs

We propose a model and numerical method for multiphase multicomponent flow simulation in porous media for solving problems of hydrocarbon reservoir engineering. The method takes into account the complex structure of reservoirs, a large number of operating wells, various methods of specifying the dependences of the properties of phases on pressure and mass fractions of their components, the transition of components from one phase to another, as well as compressibility and gravity. The method provides a sufficiently high computational efficiency, local conservation of mass for each component and protects the requirement of constraint of phases saturations from violation due to the following aspects. To calculate the phase flows we use a combination of continuous Galerkin finite element method (CG FEM) with a new balancing method. This combination provides high accuracy of flows calculated from the pressure field on rather coarse meshes. We also propose to use the grouping of finite elements, which allows us to remove strict restrictions on the time step when calculating the phase flows between cells. Verification of the proposed method is performed by comparing with the results of SPE tests. The effectiveness of the method is demonstrated both on simple model problems and on problems of modeling a complex reservoir. For a real field of high-viscosity oil in Tatarstan (Russia), a good agreement between the computed and observed data is obtained for all high-rate wells over a fairly long period of field development.

Investigation of micro-shock waves in a planar magnetogasdynamic flow using the discontinuous Galerkin finite element method

The present work focuses on the numerical investigation of micro-shock wave propagation in a two-dimensional magnetogasdynamic flow in the framework of the Dis- continuous Galerkin-Finite Element Method (DG-FEM). The Lorentz force has been im- plemented in the compressible, viscous Navier–Stokes equations as a source term using first-order spatial and fourth-order temporal Runge–Kutta discretization schemes. To in- vestigate the effect of the electrical conductivity on the micro-shock wave propagation, a two-dimensional micro-shock channel problem with hydraulic diameter of 2.5 mm, length of 82 mm, and no-slip boundary conditions at the left and at the right wall is considered as a benchmark problem. In this case, acoustic waves are generated behind after the rupture of the membrane that separates two states of the same gas originally at different pressure and density and both initially at rest. The magnetic field is taken into account as uniform and stationary throughout the microchannel, and the numerical simulations are performed in a short physical time, before the reflection of the waves on the lateral wall. A detailed parametric study of the temperature, density, pressure, and u-velocity is carried out by a variation of the electrical conductivity of the magnetogasdynamic flow, under the assumption of low magnetic Reynolds numbers. It has been found that the jumps of the acoustic waves become significantly intensified when the electrical conductivity of the gas is increased. It has also been observed that the presence of the Lorentz force causes an acceleration in the gasflow towards the outlet section of the microchannel at the low Knudsen number of 0.05. The outcome of this research work could be relevant to biomedical applications where the ability to control the flow in a microchannel has a significant impact on the development of small devices aimed to deliver pharmaceutical drugs in specific locations.

A Phase-Field-Based Approach for Modeling Flow and Geomechanics in Fractured Reservoirs

Summary Modeling the geomechanical deformations of fracture networks has become an integral part of designing enhanced geothermal systems and recovery mechanisms for unconventional reservoirs. Stress changes in the reservoir can cause variations in the apertures of fractures resulting in large changes in their transmissivities. At the same time, sustained high-injection pressures can induce shear slipping along existing fractures and faults and trigger seismic activity. In this work, we extend the phase-field method to solve for flow and geomechanical deformations in naturally fractured reservoirs. The framework can predict the opening/closure of fractures as well as their shear slipping because of induced stresses and poromechanical effects. The flow through fractures is modeled on spatially nonconforming grids using the enhanced velocity mixed finite element method. The geomechanics equations are discretized using the continuous Galerkin (CG) finite element method. The flow and mechanics equations are decoupled using the fixed stress iterative scheme. The implementation is validated against the analytical solutions of Mandel’s problem and Sneddon’s benchmark test. Two synthetic examples are presented to illustrate the impact of poroelastic deformations and the accompanying dynamic behavior of fractures on the safety and productivity of subsurface projects.

Energy stable and accurate coupling of finite element methods and finite difference methods

We introduce a hybrid method to couple continuous Galerkin finite element methods and high-order finite difference methods in a nonconforming multiblock fashion. The aim is to optimize computational efficiency when complex geometries are present. The proposed coupling technique requires minimal changes in the existing schemes while maintaining strict stability, accuracy, and energy conservation. Results are demonstrated on linear and nonlinear scalar conservation laws in two spatial dimensions.

Theory and Implementation of Coupled Port-Hamiltonian Continuum and Lumped Parameter Models

A continuous Galerkin finite element method that allows mixed boundary conditions without the need for Lagrange multipliers or user-defined parameters is developed. A mixed coupling of Lagrange and Raviart-Thomas basis functions are used. The method is proven to have a Hamiltonian-conserving spatial discretisation and a symplectic time discretisation. The energy residual is therefore guaranteed to be bounded for general problems and exactly conserved for linear problems. The linear 2D wave equation is discretised and modelled by making use of a port-Hamiltonian framework. This model is verified against an analytic solution and shown to have standard order of convergence for the temporal and spatial discretisation. The error growth over time is shown to grow linearly for this symplectic method, which agrees with theoretical results. A modal analysis is performed which verifies that the eigenvalues of the model accurately converge to the exact eigenvalues, as the mesh is refined. The port-Hamiltonian framework allows boundary coupling with bond-graph or, more generally, lumped parameter models, therefore unifying the two fields of lumped parameter modelling and continuum modelling of Hamiltonian systems. The wave domain discretisation is shown to be equivalent to a coupling of canonical port-Hamiltonian forms. This feature allows the model to have mixed boundary conditions as well as to have mixed causality interconnections with other port-Hamiltonian models. A model of the 2D wave equation is coupled, in a monolithic manner, with a lumped parameter model of an electromechanical linear actuator. The combined model is also verified to conserve energy exactly.

A weak Galerkin finite element method for time fractional reaction-diffusion-convection problems with variable coefficients

In this paper, a weak Galerkin finite element method for solving the time fractional reaction-convection diffusion problem is proposed. We use the well known L1 discretization in time and a weak Galerkin finite element method on uniform mesh in space. Both continuous and discrete time weak Galerkin finite element method are considered and analyzed. The stability of the discrete time scheme is proved. The error estimates for both schemes are given. Finally, we give some numerical experiments to show the efficiency of the proposed method.

A CG–DG method for Maxwell’s equations in Cole–Cole dispersive media

The numerical approximation of Maxwell’s equations in Cole–Cole dispersive media is studied for two and three-dimensional cases. We construct a combined scheme using the standard continuous Galerkin (CG) finite element method in time and discontinuous Galerkin (DG) method in space to discretize the Cole–Cole dispersive model. The L2-stability and error estimate of the proposed scheme are analyzed. Furthermore, we use sum-of-exponential approximation for the convolution kernel to speed up the evaluation of the Caputo derivative. Numerical examples are provided to demonstrate the performance of the proposed numerical methods.

Achieving local mass conservation when using continuous Galerkin finite element methods to solve solute transport equations with spatially variable coefficients in a transient state

Mass balance checks are often used to evaluate the abilities of numerical simulators and assess the accuracy of the solutions but are often insufficient. In addition, there is a need to calculate the difference in mass fluxes, through downstream subdomain analysis, to determine whether remediation technologies are suitable for upstream contaminated sites. The Galerkin finite element method (GFEM) has been commonly considered locally nonconservative, whereas the finite difference method (FDM) is a local mass conservation method. In this study, we proved that the GFEM has local mass balance equations on “mass conserved elements” (constructed by connecting lines that perpendicularly cross element faces at midpoints) in solute transport problems with spatially variable velocities and dispersion coefficients, identical to discretized balance equations in the FDM. Furthermore, we found that locally conservative GFEMs with two types of approaches involving spatially variable coefficients (i.e., a linear combination approach and an average approach) give rise to two types of postprocessing techniques for calculating local mass fluxes. These are equivalent to an area-weighted method and first-order Taylor series method, respectively. Finally, to illustrate the practical applicability of these GFEMs, we used the proposed postprocessing approaches to resolve practical problems involving computation of local mass fluxes, based on information directly obtainable from finite element solutions of concentration distribution. Our findings revealed that the local mass balance errors were negligibly small, regardless of numerical concentration accuracy.

Enriched Galerkin discretization for modeling poroelasticity and permeability alteration in heterogeneous porous media

In this paper, we utilize the enriched Galerkin (EG) finite element method for the flow equation in Biot's system, which provides a robust locally conservative flux in heterogeneous porous media. The computational algorithm to solve the coupled system with the permeability alteration is presented with the linearization and Picard's iterative scheme. The block structure is utilized for the linear system in numerical discretization, and the computer code is shared in the open-source platform. In the numerical experiments, we compare the proposed EG method with the classical continuous Galerkin (CG) and discontinuous Galerkin (DG) finite element methods in different scenarios, including the North sea reservoirs setup. While DG and EG methods provide similar approximations for the pressure solutions, the CG method produces spurious oscillations in fluid pressure and volumetric strain solutions near the material interfaces, especially for the soft materials. The difference of flux approximation between EG and DG methods is insignificant; still, the EG method demands approximately two and three times fewer degrees of freedom than the DG method for two- and three-dimensional geometries.

Mesh quality metrics for isogeometric Bernstein–Bézier discretizations

High-order finite element methods harbor the potential to deliver improved accuracy per degree of freedom versus low-order methods. Their success, however, hinges upon the use of a curvilinear mesh of not only sufficiently high accuracy but also sufficiently high quality. In this paper, theoretical results are presented quantifying the impact of mesh parameterization on the accuracy of a high-order finite element approximation, and a formal definition of shape regularity is introduced for curvilinear meshes based on these results. This formal definition of shape regularity in turn inspires a new set of quality metrics for curvilinear finite elements. Computable bounds are established for these quality metrics using the Bernstein–Bézier form, and a new curvilinear mesh optimization procedure is proposed based on these bounds. Numerical results confirming the importance of shape regularity in the context of high-order finite element methods are presented, and numerical results demonstrating the promise of the proposed curvilinear mesh optimization procedure are also provided. The theoretical results in this paper apply to any piecewise-polynomial or piecewise-rational finite element method posed on a mesh of polynomial or rational mapped simplices and hypercubes. As such, they apply not only to classical continuous Galerkin finite element methods but also to discontinuous Galerkin finite element methods and even isogeometric methods based on NURBS, T-splines, or hierarchical B-splines.

Modeling interactions of natural and two-phase fluid-filled fracture propagation in porous media

In this paper, a novel computational framework is introduced for simulation of multiphase flow, geomechanics, and fracture propagation in porous media based on Biot’s model for poroelasticity by focusing on interactions between hydraulic and natural fractures. Since realistic porous media contain many natural fractures, it is important not only to stimulate hydraulic fractures but also to study the interaction between natural and hydraulic fractures. Here, state-of-the-art numerical modeling of natural and hydraulic fractures using a diffusive adaptive finite element phase field approach is employed. The locally mass conservative enriched Galerkin finite element methods (EG) are utilized to model two-phase flow in propagating fractures with relative permeability and capillary pressure. Geomechanics approximated by a continuous Galerkin finite element method is coupled to multiphase flow by applying an iteratively coupled scheme. Numerical examples are presented that demonstrate the effectiveness of this framework for different propagation scenarios by varying the degrees of physics. In addition, the capabilities to perform high-fidelity simulations on complex fracture networks, with randomly joined diffusive natural fractures, are illustrated.