Contrasting phase field method and pairwise force smoothed particle hydrodynamics method in simulating multiphase flow through fracture-vug medium

Highlights

• The PFM and the PF-SPH methods are utilized for simulating multiphase flow in artificial fracture-vug medium.

• Comparing the results of simulation and experiment, PF-SPH is computationally more efficient and accurate than PFM.

• It is a stable stratified flow for LCN≤ −4 (gas) or −2.3 (water) and an unstable jet flow for LCN≥ −3 (gas) or −1.3 (water).

Abstract

Strong heterogeneity and anisotropy result in complex flow and relatively low oil recovery in fractured-vuggy reservoirs. In this work, we assess the Cahn-Hilliard phase field method (PFM) and the pairwise force smoothed particle hydrodynamics (PF-SPH) method in modeling 2D multiphase flow through fracture-vug medium, in terms of different capillary numbers, mobility ratios, gravity, and wettability. Comparing the results of numerical simulations and physical experiments, the results of the PF-SPH method are closer to the physical experimental displacement process than those of the PFM, and the PF-SPH method is found to be computationally more efficient and accurate than the PFM. The experimental and numerical results reveal that when the injection velocity, represented via logarithm of the capillary number (LCN), is not greater than −4 (gas phase case) or −2.3 (water phase case) in the fracture-vug medium, the gravitative differentiation of the oil-gas or oil-water is significant, and the oil-gas or oil-water flow is a stable stratified flow. When the injection velocity LCN exceeds −3 (gas phase case) or −1.3 (water phase case) in the vertical section, the gravitative differentiation of the oil-gas or oil-water is not obvious, and the fluid flow pattern appears as an unstable jet flow.

Introduction

In general, fractured-vuggy reservoirs exhibit strong heterogeneity and anisotropy. Fractured-vuggy reservoirs consist of a complex system of media, such as matrix, fractures and vugs, which vary significantly in scale, from micrometers to meters (Yang et al., 2020). Fluid exchange occurs between different media. The connectivity and distribution of these media are difficult to conceptualize. Therefore, it is a challenge to find a representative elementary volume (REV) in fractured-vuggy reservoirs. In addition, fluid flow in the fracture-vug medium is complex, in that it constitutes free flow in unfilled vugs and large fractures, but Darcy flow in micro-fractures and matrix. Meanwhile, the fluid flow regimes may range from laminar flow to turbulent flow in fracture-vug medium (Yao and Huang, 2017).

There are two types of approaches to simulate complex fluid flow through fractured-vuggy media: the discontinuum approach (Fadlelmula et al., 2015) and the continuum approach (Camacho Velazquez et al., 2002; Liu et al., 2003; Sitharam et al., 2001; Wu et al., 2006, 2011). The following are some commonly used simulation methods to describe fluid flow in fracture-vug medium: (a) Darcy-Stokes methods, which are Stokes equations in the vugs, Darcy's law in porous rock, and the Beavers-Joseph-Saffman boundary condition on the interface between the two regions (Arbogast and Brunson, 2007; Arbogast and Lehr, 2006; Huang et al., 2010; Peng et al., 2009; Popov et al., 2009); (b) Stokes-Brinkman methods (Brinkman, 1949), which continuously vary from a Darcy dominated flow to a Stokes dominated flow, which avoid the explicit formulation of the boundary conditions at the fluid-porous interfaces (Gulbransen et al., 2009; He et al., 2015; Popov et al., 2007); and (c) the discrete fracture-vug network model (DFVN) (Huang et al., 2010; Yao and Huang, 2017), which constitutes a discontinuum approach based on a two-scale homogenization limit theory, vugs and fractures as discrete units are embedded into rock, together with a microscopic Darcy's law governing the porous flow in the fracture and matrix, Navier–Stokes equation governing the free flow in the vug, and the coupling of free flow and porous flow. Since the permeability of the matrix system is much smaller than the other two systems, several studies indicate that fluid flow in the matrix system is negligible (Zhang et al., 2011; He et al., 2017; Hui et al., 2015; Lyu et al., 2017).

A variety of multiphase flow simulation methods that can be implemented at pore-scales have been suggested, e.g., the level set method (LSM) (Amiri and Hamouda, 2013; Amiri et al., 2015; Osher and Sethian, 1988; Prodanovic et al., 2010; Sethian and Smereka, 2003), the phase field method (PFM) (Amiri and Hamouda, 2013; Badalassi et al., 2003; Jacqmin, 1999), the pairwise force smoothed particle hydrodynamics (PF-SPH) method (Bandara et al., 2013; Tartakovsky and Panchenko, 2016), etc. In the LSM (Ahmadi and Shadizadeh, 2013), re-initialization is introduced to ensure sharp interfaces. However, contact discontinuity at the interface is dissipated, the algorithm is not conserved, and the interface is distorted. The PFM is derived from the theory proposed by Der Waal (Anderson et al., 1998). The interface of the PFM is formed by two-phase interdiffusion. The PFM possesses certain advantages. Firstly, the evolution of the interface is semi-continuous and does not require initialization like the LSM (Zhou et al., 2010a). Secondly, the energy equation ensures stable calculation, and it can calculate the phenomenon of breaks or fusions of fluid flow. For isothermal flow, the interfacial tension is equal to the sum of the free energy densities through the interface (Yue et al., 2006).

Both the LSM and the PFM compute multiphase flow on a fixed Eulerian grid with an implicit front tracking scheme (Amiri and Hamouda, 2013). Due to the high computational cost, the LSM and the PFM are suitable for simple pore geometries (Bandara et al., 2013). In addition, based on Euler grids, due to mesh adaptability and connectivity issues, these techniques often encounter difficulties in dealing with complex merging and fragmentation of dynamic interfaces (Xu et al., 2015). In contrast, as a meshless Lagrangian method, smoothed particle hydrodynamics (SPH) offers significant advantages in this respect. In the SPH method, the continuum is decomposed into a collection of discrete particles associated with its own physical properties, which constantly update the density and velocity through explicit integration of the continuity and momentum equations (Xu et al., 2015). Therefore, the interface shape can be easily identified and tracked based on the natural distribution of the particles in each phase.

The pairwise force smoothed particle hydrodynamics (PF-SPH) method is used to simulate multiphase flow in porous media (Tartakovsky and Meakin, 2006; Tartakovsky and Panchenko, 2016). Kazemi et al. (2019) and Khayyer et al. (2018) solved the problem of fluid flow interaction at the interface between porous media flow and free flow, based on the weakly compressible SPH (WCSPH) and the incompressible SPH (ISPH) model, respectively. Zhou et al. (2010b) simulated the flow mechanism of water flooding in fractured-vuggy media based on the SPH method, and studied the effects of water injection velocity, wettability, fracture opening, and other factors on the recovery. As a fully Lagrangian method, the PF-SPH does not suffer from numerical dispersion, and interface dynamics can be modeled without complex interface tracking schemes. Moreover, due to its molecular dynamics isomorphism, complex dynamic wetting behaviors can be calculated (Bandara et al., 2013). The surface tension and wetting behavior of fluids are modeled by pairs of molecular forces in the PF-SPH method. The accuracy and consistency of the PF-SPH model has been demonstrated for simulating interface tension and (static or dynamic) contact angles in regions with simple geometries, such as flat surfaces and fractures with uniform apertures (Kordilla et al., 2013; Tartakovsky and Meakin, 2006).

In this paper, we contrast the PFM (Amiri and Hamouda, 2013; Badalassi et al., 2003; Jacqmin, 1999) and PS-SPH (Tartakovsky and Meakin, 2006; Tartakovsky and Panchenko, 2016) methods in modeling multiphase flow through the fracture-vug medium. Firstly, based on a series of displacement experiments in a micromodel that was created with three-dimensional (3D) printing technology representing a uniform two-dimensional (2D) fracture-vug medium by Yang et al. (2020), we utilize the PFM and PF-SPH methods to simulate displacement of the oil phase that originally occupies the micromodel, while the other fluid (water or gas) is injected into the media at a constant flow rate under different capillary numbers, mobility ratios, gravity, and wettability. Through comparing the results of the simulation methods with those of physical experiments, we study the characteristics of PFM and PF-SPH methods for multiphase flow in fracture-vug medium in displacement process, accuracy, computation time, etc. In the next section, the numerical simulation methods are introduced. Subsequently, the following section mainly focuses on the comparative results of physical experiments and numerical simulations, and discusses observations in detail. The paper concludes with a summary of the main findings from this study.