Efficient Field-Scale Simulation of Black Oil in a Naturally Fractured Reservoir Through Discrete Fracture Networks and Homogenized Media

Summary This paper describes a hybrid finite volume method, designed to simulate multiphase flow in a field-scale naturally fractured reservoir. Lee et al. (WRR 37:443-455, 2001) developed a hierarchical approach in which the permeability contribution from short fractures is derived in an analytical expression that from medium fractures is numerically solved using a boundary element method. The long fractures are modeled explicitly as major fluid conduits. Reservoirs with well-developed natural fractures include many complex fracture networks that cannot be easily modeled by simple long fracture formulation and/or homogenized single continuity model. We thus propose a numerically efficient hybrid method in which small and medium fractures are modeled by effective permeability, and large fractures are modeled by discrete fracture networks. A simple, systematic way is devised to calculate transport parameters between fracture networks and discretized, homogenized media. An efficient numerical algorithm is also devised to solve the dual system of fracture network and finite volume grid. Black oil formulation is implemented in the simulator to demonstrate practical applications of this hybrid finite volume method. Introduction Many reservoirs are highly fractured due to the complex tectonic movement and sedimentation process the formation has experienced. The permeability of a fracture is usually much larger than that of the rock matrix; as a result, the fluid will flow mostly through the fracture network, if the fractures are connected. This implies that the fracture connectivities and their distribution will determine fluid transport in a naturally fractured reservoir (Long and Witherspoon 1985). Because of statistically complex distribution of geological heterogeneity and multiple length and time scales in natural porous media, three approaches (Smith and Schwartz 1993) are commonly used in describing fluid flow and solute transport in naturally fractured formations:discrete fracture models;continuum models using effective properties for discrete grids; andhybrid models that combine discrete, large features and equivalent continuum. Currently, most reservoir simulators use dual continuum formulations (i.e., dual porosity/permeability) for naturally fractured reservoirs in which matrix blocks are divided by very regular fracture patterns (Kazemi et al. 1976, Van Golf-Racht 1982). Part of the primary input into these simulation models is the permeability of the fracture system assigned to the individual grid-blocks. This value can only be reasonably calculated if the fracture systems are regular and well connected. Field characterization studies have shown, however, that fracture systems are very irregular, often disconnected, and occur in swarms (Laubach 1991, Lorenz and Hill 1991, Narr et al. 2003). Most naturally fractured reservoirs include fractures of multiple- length scales. The effective grid-block permeability calculated by the boundary element method becomes expensive as the number of fractures increases. The calculated effective properties for grid-blocks also underestimates the properties for long fractures whose length scale is much larger than the grid-block size. Lee et al. (2001) proposed a hierarchical method to model fluid flow in a reservoir with multiple-length scaled fractures. They assumed that short fractures are randomly distributed and contribute to increasing the effective matrix permeability. An asymptotic solution representing the permeability contribution from short fractures was derived. With the short fracture contribution to permeability, the effective matrix permeability can be expressed in a general tensor form. Thus, they also developed a boundary element method for Darcy's equation with tensor permeability. For medium-length fractures in a grid-block, a coupled system of Poisson equations with tensor permeability was solved numerically using a boundary element method. The grid-block effective permeabilities were used with a finite difference simulator to compute flow through the fracture system. The simulator was enhanced to use a control-volume finite difference formulation (Lee et al. 1998, 2002) for general tensor permeability input (i.e., 9-point stencil for 2-D and 27-point stencil for 3-D). In addition, long fractures were explicitly modeled by using the transport index between fracture and matrix in a gridblock. In this paper we adopt their transport index concept and extend the hierarchical method:to include networks of long fractures;to model fracture as a two-dimensional plane; andto allow fractures to intersect with well bore. This generalization allows us to model a more realistic and complex fracture network that can be found in naturally fractured reservoirs. To demonstrate this new method for practical reservoir applications, we furthermore implement a black oil formulation in the simulator. We explicitly model long fractures as a two-dimensional plane that can penetrate several layers. The method, presented in this paper, allows a general description of fracture orientation in space. For simplicity of computational implementation however, both the medium-length and long fractures considered in this paper are assumed to be perpendicular to bedding boundaries. In addition, we derive a source/sink term to model the flux between matrix and long fracture networks. This source/sink allows for coupling multiphase flow equations in long fractures and matrix. The paper is organized as follows. In Section 2 black oil formulation is briefly summarized and the transport equations for three phase flow are presented. The fracture characterization and hierarchical modeling approach, based on fracture length, are discussed in Section 3. In Section 4 we review homogenization of short and medium fractures, which is part of our hierarchical approach to modeling flow in porous media with multiple length-scale fractures. In Section 5 we discuss a discrete network model of long fractures. We also derive transfer indices between fracture and effective matrix blocks. In Section 6 we present numerical examples for tracer transport in a model with simple fracture network, interactions of fractures and wells, and black oil production in a reservoir with a complex fracture network system. Finally, the summary of our main results and conclusion follows.