near-well modeling to reservoir heterogeneity and formation damage effects

Journal of Petroleum Science and Engineering 176 (2019) 640–652

Contents lists available at ScienceDirect

Journal of Petroleum Science and Engineering journal homepage: www.elsevier.com/locate/petrol

Applications of fully coupled well/near-well modeling to reservoir heterogeneity and formation damage effects

T

Jie Caoa,b, Nan Zhanga,b,∗, Thormod E. Johansenc a

School of Petroleum Engineering, Xi'an Shiyou University, Xi'an, 710065, China Xi'an Shiyou University Engineering Research Center of Development and Management for Low to Extra-Low Permeability Oil & Gas Reservoirs in West China, Ministry of Education, China c Upstream Petroleum Research and Consulting Inc., St. John's, Canada b

A R T I C LE I N FO

A B S T R A C T

Keywords: Well/near-well modeling Iterative coupling Heterogeneity Crossflow Formation damage

Applications of a new well/near-well fluid flow modeling scheme to heterogeneous reservoirs and formation damage are presented in this paper. This well/near-well model is based on an analytically solution of an axialradial productivity model for steady state flow in wells and near-well reservoirs, in which wellbore hydraulics and heterogeneity in the near-well region are incorporated. The numerical model includes a reservoir region describing both radial and axial reservoir flow. It can be used as a stand alone simulator, or it can be coupled with a reservoir simulator, using a new iterative coupling scheme described in this paper. When coupling with the reservoir simulator, the reservoir region in the well/near-well model is shared by and computed in both models, which ensures the stability and convergence of the iterative coupling scheme. A new coupling approach is presented without needing to determine any partial derivatives. The capabilities of the modeling methodology are demonstrated through applications in three challenging cases: highly heterogeneous reservoir, well/reservoir crossflow, and formation damaged reservoirs. In the highly heterogeneous near-well reservoir case, results show that the model is stable for a wide range of permeabilities with large contrasts. The detailed permeability information from well logging can be fully taken into account in this model. The model determines pressure distribution and volumetric flow rates including flow directions, therefore enabling the model to estimate crossflow as negative flow rates, for example when the well is producing from several non-communicating regions. Furthermore, the well/near-well model is iteratively coupled with a reservoir simulator such that the well completion effects are taken into consideration in large scale simulations. The transient case of this coupling to evaluate formation damage along a horizontal well is then presented. The model is not limited to single phase flow. A two phase flow example coupling the well/near-well model and a streamline reservoir model is finally presented.

1. Introduction Productivity estimation of horizontal/advanced wells is essential in reservoir developments. This has been investigated through analytical methods, semi-analytical methods and numerical methods. Most analytical models are easy to use and it is straightforward to calculate a productivity index (PI = q/Δp). They represent all details of the reservoir flow rigorously subject to the assumptions and boundary conditions for which they were derived. For example, Joshi (1988) presented a steady state inflow equation for a horizontal well draining symmetrically from an elliptical cylinder geometry, assuming constant outer boundary pressure and well bottom hole pressure. Babu and Odeh (1989) proposed a simplified well productivity equation for a

∗

horizontal well located arbitrarily in a closed rectangular shape at semisteady state flow conditions. Economides et al. (1996) presented a productivity index approximation for arbitrarily oriented wellbore using dimensionless pressure of a point source in a no flow boundary box shaped reservoir. Dikken (1990) presented an analytical method coupling single-phase turbulent well flow with stabilized reservoir flow. Pressure losses, especially those caused by friction in turbulent flow, along the horizontal well is considered a significant factor in productivity calculations. Schulkes and King (2002) presented the coupled hydraulic model for flow in an annulus and a liner. The model was solved analytically and numerically to investigate the simultaneous effects of the reservoir layer and the annulus boundary layer. The uniform pressure condition assumes the pressure is constant

Corresponding author. School of Petroleum Engineering, Xi'an Shiyou University, Xi'an, 710065, China. E-mail address: [email protected] (N. Zhang).

https://doi.org/10.1016/j.petrol.2019.01.091 Received 16 December 2017; Received in revised form 23 December 2018; Accepted 27 January 2019 Available online 01 February 2019 0920-4105/ © 2019 Elsevier B.V. All rights reserved.

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 1. Pressure and flow schematics in the near-well region.

pressure equations were fully coupled and solved implicitly. Holmes et al. (2010) also presented a comprehensive and flexible multi-segment well model which handles global phase-component partitioning for both isothermal and thermal modes. Wan et al. (2008) presented a coupled well and reservoir model for well completion optimization. The well model, representing completion details with a network of nodes and connections, determines flow rates and pressures in tubing and annulus given well block information, and passes pressure differences to a reservoir simulator as constant pressure constrains. This two-way coupling approach was called Tier-1 stage and extended by Grubert et al. (2009) to a Three-Tier approach. The Tier-2 stage is a dynamic one-way coupling, applying the well model results for a general depletion development. The Tier-3 stage is a static study of various completion methods run and optimized based on a steady state simulation, which is the simplest and fastest screening. Another one way coupling by feeding reservoir transient performance automatically to the wellbore model was applied in Livescu et al. (2010) for well completion screening and optimization. The Tier-1 and Tier-3 stages are also integrated by Jain et al. (2013) in the workflow of completion design for heavy oil production without explicit coupling of the well and reservoir simulators. For the iterative coupling approaches, the hydraulic well model and the reservoir model are solved independently at each time step and variables are exchanged until convergence, as discussed in Johansen and Khoriakov (2007) and Khoriakov et al. (2012). This approach offers more flexibility than the fully implicit solution. The iterative coupling approach was first presented by Brekke et al. (1993) through the connection points at the sand surface between the reservoir model and the well model. In these connection points, the transfer variables including ∂q phase flow rate (qα ), phase pressure ( pα ) and partial derivatives ( α ) are

inside the wellbore and the well has an infinite conductivity. Under the uniform pressure condition, shown in Fig. 1 a), the flow rate increases linearly from the toe to the heel. However, the well inflow including wellbore hydraulics is not normal to the well trajectory, as shown in Fig. 1 b). The color map in Fig. 1 represents the pressure while the arrow length represents the flux. The axial component of the inflow, i.e. the axial reservoir flow along the well trajectory used in this paper, was solved using the coupled axial-radial productivity model for steady state in Johansen et al. (2015), and semi-steady state flow in Johansen et al. (2016). Semi-analytical methods for horizontal and multi-lateral wells aim to maintain higher accuracy through combining analytical well inflow and wellbore hydraulics. Penmatcha and Aziz (1999) presented a semianalytical model for transient reservoir/wellbore coupled models for infinite and finite conductivity wellbores using superposition in space and time for each well segment. Kurtoglu et al. (2008) proposed a semianalytical representation method to determine the well and near-well flow convergence problems in reservoir simulation. It used the analytical solution from potential theory to determine the pressure in the well and near-well domain treating the well as a source term. Similarly, Tabatabaei and Ghalambor (2011) coupled the reservoir inflow to the wellbore with a box-shaped drainage domain, and the flow within the wellbore used a semi-analytical model. Besides Cartesian and hybrid grids, Karimi-Fard and Durlofsky (2012) investigated using unstructured grids to improve the prediction of fluid saturation and pressure distribution. Cao et al. (2015) applied the analytical coupled axial-radial near-well productivity models in a numerical simulation scheme, which was proved to be of higher order than a standard finite difference method. In conventional finite difference methods, the well is represented by point- or line-sources, and local grid refinements are usually needed in the near-well region, where pressure and saturation change sharply. The hybrid grid method was proposed by Pedrosa and Aziz (1986) as a better way to represent wells in reservoir simulation and to apply local grid refinement. Palagi and Aziz (1994) investigated the Voronoi grids, a flexible gridding technique also known as perpendicular bisection grid for complicated boundary conditions and near-well refinement. Krogstad and Durlofsky (2007) proposed a multiscale mixed finite element method to couple the wellbore and near-well flow, which included fine grid blocks resolving the well trajectory and coarse grid blocks representing the reservoir. To consider both remote reservoir effects and well completion and hydraulic effects, reservoir simulation employs either fully coupled and full implicit approach or an iterative approach. Holmes et al. (1998) presented the multisegment well model for advanced horizontal wells, representing the wellbore by a number of segments. The flow rate and

∂pα

first generated from the reservoir model by determining the inflow performance relationship curves, and subsequently transferred to the well model. The process continues to update the transfer variables using ∂q the partial derivatives α until the convergence criteria are satisfied. ∂pα

One challenge in this is the determination of the partial derivatives. In Brekke et al. (1993), the secant method was chosen, which needs two initial iteration values. Advanced techniques for faster and more accurate well and nearwell modelling are in great demand. The stable and efficient dynamic coupling schemes between the well/near-well flow and the reservoir flow are essential for evaluation of advanced wells. For example, as intelligent well completions are being used increasingly, the dynamic simulation of advanced wells in the reservoir is needed for the short and long term production optimization (Ranjith et al., 2017). Furthermore, The near-well heterogeneity (Fokker et al., 2012; Tamayo-Mas et al., 2016; Sun et al., 2017), crossflow (Holmes, 2001) and formation 641

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

surrounded by a circular cylindrical reservoir domain, which for a given segment is assumed to be homogeneous and isotropic. Each segment has a radial permeability kr in the plane perpendicular to the well trajectory and an axial permeability k x along the well trajectory. We assume the model is two dimensional with flow in the x − and r − directions only, ignoring flow in the angular direction. The radial diffusivity equation was solved analytically for steady state flow (Johansen et al., 2015) and semi-steady state flow (Johansen et al., 2016). It is emphasized that these are fairly complex analytical solutions of partial differential equations with specific boundary conditions. For a full understanding of these models, the detailed derivations are required. For steady state flow, the analytical solution to the radial and axial flow rates are given by

damage effects (Renard and Dupuy, 1991; Parn-anurak and Engler, 2005) are by experience challenging cases and significant for well productivity estimation and reservoir management. In this research, the numerical well/near-well model represents near-well flow considering formation damage effects and well hydraulics and completion effects. The main difference between the new numerical model and most previous models is that the numerical scheme used here is based on an analytical solution of coupled axialradial productivity model, in which axial near-well reservoir flow along well trajectory and radial inflow is solved simultaneously. The nearwell reservoir region is shared by and computed in both the well model and the reservoir model. We also present a new coupling scheme using an iterative approach without needing to determine any partial derivatives. This will be described in details in Section 3. The numerical well/near-well model is then applied to study near-well heterogeneity effects, cases where crossflow occurs, and formation damage effects. As a stand alone simulator, it is demonstrated that large permeability contrasts can be handled without smoothing of contrasts. The model determines pressure distribution and volumetric flow rates including flow directions; hence crossflow between the well and the near-well reservoir will be represented at any location. Experience has shown that large permeability contrasts and crossflow represent numerical challenges for lower order methods, while the method in this paper is higher order. This is believed to be the reason why these difficulties and not encountered in this paper. The new well/near-well model can also be coupled with standard finite difference reservoir simulators such that both the well completion effects and remote reservoir effects are taken into account. During the drilling process, the mud invasion varies along the well trajectory, resulting in a non-uniform skin distribution. The skin distribution is usually assumed to be a truncated elliptical cone (Frick and Economides, 1993). Calculated examples considering unevenly distributed formation damage effects using this coupling also demonstrate the usefulness of the new well/near-well model. In a purely Cartesian model, the effect of formation damage can only be represented by a skin factor, or by a lower permeability in the inner Cartesian grid blocks. In this paper, such damage is accurately represented in a cylindrical geometry. Furthermore, the well/near-well model is capable for multiphase flow. An example coupling the well/near-well model and a streamline reservoir simulator is presented at the end, demonstrating a fast concept for modeling two phase flow in both the remote reservoir and the well/near-well domain.

Qr = Tr (Δr pL + Δr p0 );

(1)

Qx = Tx1 (Δr p0 − Δr pL ) + Tx 2 Δx pw ,

(2)

where

Tr =

πkr L ; μ [ln(re / rw ) + ST ]

(3)

Tx1 =

πk x ⎡r 2ln ⎛ re ⎞ + 1 r 2 − 1 r 2 ⎤; e w e μLln(re / rw ) ⎢ 2 2 ⎥ ⎝ rw ⎠ ⎦ ⎣

(4)

Tx 2 =

πk x (re2 − rw2 ) , μL

(5)

⎜

⎟

where Δr pL , Δr p0 , Δr p0 and Δr PL are pressure differences in radial (Δr ) and axial (Δx ) directions, respectively; L is the segment length; re is the external radius; rw is the wellbore radius; μ is the fluid viscosity; and ST is the total skin including mechanical skin like perforation skin and formation damage skin. For multi-phase flow, the transmissibility terms (Tr , Tr and Tr ) should be multiplied by the relative permeability of each phase, which is a function of the phase saturations. For each segment in the near-well region and in the well, conservation equations for momentum and mass are used to derive the above equations. The numerical model in this paper is using a discretization of a circular cylindrical well segment. The radial circulars from the innermost to the outermost are the wellbore, the annulus, the near-well reservoir and the external reservoir. In Fig. 3, this is depicted for the case of two reservoir tori and one well annulus. The reservoir part of each segment torus is assumed to be homogeneous and isotropic, consisting of four pressure nodes located at the segment “corners” see Fig. 3. Two of them are located at the outer torus surface and two on the sand surface. Recall that we assume no angular flow, i.e. the pressure is constant on a concentric circle and can thus be represented by a single node. In a finite difference simulator, these two outermost nodes correspond to the pressure nodes of adjacent grid blocks along the well trajectory. Based on the solution of the analytical model, the pressure distribution is then piece-wise linear in the axial direction and piecewise logarithmic in the radial direction. This is in contrast to the piecewise constant distribution used in a standard finite difference method. Therefore, the model presented in this paper is of higher order accuracy, which is proved in Cao et al. (2015). The flow direction in the well/near-well model is in general not known. Hence, it must be solved for in the numerical modeling. The positive direction in this system is defined as (shown in Figs. 2 and 3): flow from the heel to the toe in axial direction; from reservoir to the wellbore in radial direction. The flow rate may be negative. For example, if the well is completed through several non-communicating parts of the reservoir (Holmes, 2001), crossflow from the well to the reservoir may occur. It is also possible that the axial flow is not from the heel to the toe, for example with high near-well reservoir pressure near the heel or with inflow control devices in the annulus. Two such cases are considered in the following sections for various external reservoir pressure distributions.

2. The well/near-well numerical model The coupled axial-radial near-well flow model was presented in previous research and is briefly summarized below. A given well trajectory is approximated by a curve composed of piece-wise linear segments. Consider one such segment as shown in Fig. 2. The well is

Fig. 2. The well segment model with coupled axial-radial flow. 642

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 3. The well/near-well model cross section view.

well/near-well model. The well/near-well model (below denoted by Well Model) can consist of any number of concentric (circular) cylinders representing the near well reservoir and well completion details (by using completion skins) with possibly multiple annuli. Here, we assume only one annulus such as for a liner in an open hole or a bare foot open hole completion, and one near-well reservoir where formation damage may exist, as shown in Fig. 3. Also in Fig. 4, we notice that in the reservoir simulation model (briefly denoted Reservoir Model), the well block pressure nodes are shifted relative to the well segment external pressure nodes in the Well Model. The two simulation grids are staggered; and the external Well Model pressures are taken from the simulation grid block pressures. In cases where the Reservoir Model uses a block centered grid system, the pressure at the edges are determined based on the flux continuity at those edges. For example as shown in Fig. 4,

The horizontal well completion is along the axis of a cylindrical reservoir. The heel and the toe ends of the near-well reservoir are treated as (but not limited to) no-flow boundaries. The volumetric flow rate or wellbore pressure at the heel segment external node together with the external reservoir pressures are imposed as boundary conditions. Therefore, any combination thereoff such as target rate with limiting bottom hole pressure can also be used. The momentum balance equation for flow between adjacent reservoir segments is described by Eqs. (2)–(5). For the frictional flow between adjacent well segments in the tubing and in the annulus, the momentum balance equation is expressed as:

∂pf ∂x

=−

fρv 2 . 2D

(6)

where f is the Darcy friction factor (Colebrook, 1939; Moody, 1947), ρ is the fluid density, v is the volumetric flux and D is the pipe diameter. For 0.3164 turbulent flow in smooth pipes, f = 4 , where Re is the Reynolds Re number. Considering the friction pressure loss ΔPf , the frictional flow in the pipe (Qf ) can be written as

Qf =

4 sign (ΔP ) Tf ΔP f7 ,

Pi − 1/2 =

D19/7 , ρ3/7 μ1/7 Δx 4/7

(9)

where i denotes the grid block number in the axial direction. In coupled simulations, the communication between the two models takes place at two levels. The well block pressures are used in the Well Model as the external pressures. In addition, the annulus pressures from the Well Model are used as flowing wellbore (connection) pressures in the Reservoir Model. When running the reservoir simulator as standalone, a connection factor for each well grid block is used to calculate the well inflow, such as the Peaceman model (1983). As we shall see, this is needed only on the first time step in the coupled model. In this way, the coupled model is forced to incorporate detailed information which would not be possible if the two models were used as stand-alone models; for example the grid block pressure is influenced by large scale

(7)

where Tf is the frictional flow transmissibility

Tf = cf

ki − 1 Δx i Pi − 1 + ki Δx i Pi , ki − 1 Δx i + ki Δx i − 1

(8)

and cf is a constant accounting for units conversion; for example cf = 1.4006 × 108 in SI units. The frictional flow in the annulus is usually treated in the same manner as flow in a pipe, except the hydraulic diameter Da − Dp being used for D, where Da and Dp is the inner annulus diameter and outer liner diameter, respectively. The complete numerical model is given by Eqs. (1), (2) and (7) for single phase steady state flow, together with the associated parameters described above. For multiphase flow, the momentum balance equations are applied to each phase in the same manner, as described in details by Johansen and Khoriakov (2007). Here, the unknowns consist of 4N flow rates and 2N + 2 pressures. The equations are 4N momentum balances, 2N material balances and 2 boundary conditions. Therefore, the problem is closed. Since most of the equations are nonlinear, an iterative method must be used. For this, the Newton-Raphson method was chosen. 3. Coupling of well/near-well model and reservoir model In this section, the iterative coupling scheme between well/nearwell model and the reservoir model is presented in detail. In general, flow transients are represented through the coupling of the two models on each time step, notwithstanding the steady state assumption in the

Fig. 4. Coupling scheme with Reservoir Model and Well Model. 643

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

worthwhile to note that in the coupled simulation, the only time when conventional connection factors are used is for n = 0 , i.e. the very first time step. In this time step, the Reservoir Model is using the heel bottom hole pressure as flowing well connection pressure in all connection, i.e. the detailed completion data resulting from local pressure losses (for example friction) in the well and near-well region are ignored. However, through the iterations these details are represented by updated connection factors ω (Eq. (11)) and well connection pressures. In the coupling process, transients are present through the time stepping of the Reservoir Model. To allow interpretation of transient well tests in long, horizontal well completions, it requires wellbore transients for wellbore storage calculations as well as reservoir transients. Hence, the well model of this paper must be extended to transient flow of compressible fluids in wells, as presented by Khoriakov et al. (2012). This constitutes as further development of this model. The above coupling was described for single phase flow. For multiphase flow, the saturations at well completion grid blocks are calculated by the Reservoir Model and can therefore be included in the iterations by using reservoir flow rates for each phase as kr (S ) ∗ qB , where kr is phase relative permeability, and well flow rates as α ∗ q where α is phase hold up. This was described in detail in Johansen and Khoriakov (2007). In the following Sections, three main applications with calculated examples are described. The first part in Section 4 presents the application of the model in heterogeneous reservoirs. Various heterogeneity distributions in the near-well region are used in the model and results show that the model can handle large permeabilities contrast without any smoothing treatment. Next, in Section 5, various external reservoir pressure profiles are specified. In a poorly communicating reservoir, the well inflow rates are negative in some parts, which means flow direction is from well to external reservoir, i.e. crossflow occurs. In this case, the crossflow is captured by the new model. In Section 6, the coupling of reservoir models with the well/near-well model are demonstrated and compared with coupling reservoir models with a friction well model. An application is also presented for a horizontal well with uneven distribution of formation damage.

(remote) reservoir effects (faults, other wells ect.) and the well annulus pressures is influenced by completion details (inflow devices, perforations, any well flow restriction). The coupling between the two models is iterative on each reservoir simulation time step. It is assumed that only slow transients are important. In and near the well the flow is fast, and a steady state well model is therefore used. In Fig. 4, the iterative coupling is illustrated for one well segment. Note that in Fig. 4, the liner is not shown since it is not being used in the coupling. However, it is being used in well hydraulics calculation, which influences the annulus pressure. This iteration is performed on all segments simultaneously, and on each time step. The grid block pressures PB from the Reservoir Model are used as external pressures in the Well Model, and the annulus pressures Pw from the Well Model are used as connection wellbore pressure in the Reservoir Model. We describe an iteration cycle k → k + 1 which starts (k = 1) at the end of a reservoir simulation time step nΔt , where Δt is the time step and n is the number of time steps. First, the Reservoir Model is started at time nΔt using a constant connection wellbore pressure over the time step. This produces estimates of reservoir grid block pressure PB(k ) and well connection inflow rate qB(k ) at time (n + 1)Δt . Second, the Well Model is using PB(k ) as external node pressures and is run to produce annulus pressure Pw(k ) and flow rates qw(k ) from external nodes to adjacent annulus nodes. The connection factor ω(k ) for grid blocks are then calculated using

ω(k ) =

qw(k ) PB(k )

− Pw(k )

.

(10)

If for a given tolerance ε and reference flow rate qo , N

1/2

⎡ (k ) (k ) 2⎤ ⎢∑ (qB, i − qw, i ) ⎥ 1 i = ⎦ ⎣

< εqo ,

(11)

the iterations are stopped and values at time (n + 1)Δt are: Pw(k ) is flowing well connection bottom hole pressure, PBk is grid block pressure and ω(k ) is well connection factor; k here being the last iteration count. In Eq. (11) the sum is over all the segments. If Eq. (11) is not satisfied, k : =k + 1; and the procedure is repeated using values from the previous iteration. The process is summarized in the flow chart as shown in Fig. 5. In the above iteration process, both the Reservoir Model and the Well Model are run with the same overall well operating conditions, such as heel bottom hole pressure or total flow rate. It is also

4. Application in heterogeneous reservoirs Reservoir permeability is mostly heterogeneous. This heterogeneity results in changes of transmissibility between adjacent reservoir grid blocks, and between well blocks and wellbores. In the new well/near well model, calculation of flow rates and pressures are segment-based, which means that each segment of well and near-well reservoir are given a unique set of properties such as permeability. In this subsection, various distributions of permeability are used in the model to verify the convergence of the simulation. Basic parameters used here are given in Table 1. We start with a homogeneous case, where the permeability is 50 mD for comparison purposes. The permeability distribution and the solution of pressure distribution are given in Figs. 6 and 7. Next, a step wise permeability distribution is imposed with a permeability in the left half of the well of 5 mD and in the right half of 1000 mD, as shown in Fig. 8. Table 1 Parameters used in the coupling model.

Fig. 5. The iterative coupling scheme with Well and Reservoir models. 644

Parameter

Unit

Value

Wellbore length Wellbore radius External reservoir radius Oil viscosity Oil formation volume factor Oil Density

Lw rw re μo Bo ρo

[m] [m] [m] [cP ]

External pressure Wellbore pressure

Pe Pw

[kg / m3 ] [bar ] [bar ]

1000 0.1 20 1.18 1.3 800 320 300

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 9. Pressure distribution with step wise permeability.

Fig. 6. Permeability distribution in homogeneous case.

Fig. 10. Well inflow rate with step wise permeability.

Fig. 7. Pressure distribution in homogeneous case.

Fig. 11. Wellbore pipe flow rate with step wise permeability. Fig. 8. Step wise permeability distribution.

in Fig. 10. Next, a randomized permeability distribution within the range from 20 to 200 mD is used in the near-well model, as for example often is experienced in well logs in shaly reservoirs. In real reservoirs, the permeability often varies in a more smooth way, which is easier to handle numerically than the permeability shown in Fig. 12. In this case, the pressure is solved without smoothing of the permeability field and the results show a set of smooth curves in the wellbore and the near-

The pressure profiles (Fig. 9), well specific inflow rates and the flow inside the wellbore are shown in Fig. 10 and Fig. 11, respectively. The specific well inflow first increases from toe to the middle due to increasing pressure difference between the wellbore and near-well reservoir. The specific well inflow still increases from the middle to the heel although much less significantly since permeability is much lower in this part. The corresponding well inflow rate per length unit is shown

645

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 15. Wellbore pipe flow rate with randomized permeability. Fig. 12. Randomized permeability distribution.

5. Application to crossflow situations In the well/near-well model, the flow rate is defined with directions. The positive direction in well inflow means flow from near-well reservoir toward the well, and in axial near-well reservoir flow is positive from toe to heel. Depending on the external reservoir pressure, the flow rate may result in negative values. Two representative types of special external reservoir pressure distributions are used in this Section. Basic parameters used are the same as in the previous Section, except external reservoir pressure and the reservoir permeability, the latter being 100 mD . In the first case, the external reservoir pressure is decreasing from the heel to the toe. The solutions of wellbore pressure and near-well reservoir pressure in this case are plotted in Fig. 16. The well inflow in Fig. 17 is generally increasing from toe to the heel while the axial reservoir flow in Fig. 18 is partly negative. This means that in the section from the toe to around 300 m from the heel, the axial reservoir flow is toward the toe side, as shown in Fig. 19. In the second case, the external reservoir has a low pressure region due to lateral poor communication within the compartmentalized reservoir. The external reservoir pressure is 330 bar at the heel and the toe ends, separated by the low pressure region of 315 bar . The bottom hole pressure is 310 bar . The calculated pressure distribution is shown in Fig. 20, where in the low pressure region, the wellbore pressure is larger than the near-well reservoir pressure, which is larger than external reservoir pressure. In this case, crossflow from the wellbore to the near-well reservoir occurs, which represents loss of fluids. As shown in Fig. 21, the well inflow is negative in the low pressure region. In the

Fig. 13. Pressure distribution with randomized permeability.

Fig. 14. Well inflow rate with randomized permeability.

well reservoir (Fig. 13). The well inflow in Fig. 14 oscillates from the toe to the heel due to permeability heterogeneity, but generally reflects an increasing trend. The flow rate inside the wellbore generally increases from toe to heel however in this case in a less smooth fashion due to varying inflow (see Fig. 15).

Fig. 16. Pressure profiles with linear external reservoir pressure. 646

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 17. Well inflow rate with linear external reservoir pressure.

Fig. 20. Pressure profiles with a low external reservoir region.

Fig. 21. Well inflow rate with a low external reservoir region.

Fig. 18. Axial near-well reservoir rate with linear external reservoir pressure.

Fig. 19. The flow direction in the well/near-well model given linear external reservoir pressure.

Fig. 22. Axial near-well reservoir rate with a low external reservoir region.

6. Coupled simulation, formation damage

axial near-well reservoir flow in Fig. 22, there are two jumps in flow rates due to the sharp change of pressure. The flow directions are shown in Fig. 23 to represent the crossflow. Experience has shown that cross flow may be difficult to capture by low order methods, as opposed to the approach used in this paper.

The iterative coupling scheme described in Section 3 is demonstrated through two examples. For convenience of comparison, two Well Models are used; one is the new well/near-well model and one is a simple well model that considers friction pressure loss in a pipe with radial inflow (using Peaceman's model), denoted as Friction Model. One major difference in the two models is that axial reservoir flow is ignored 647

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 23. The flow direction in the well/near-well model given a low external reservoir region. Fig. 24. The wellbore pressure in different iterations for the first time step. Table 2 Parameters used in the coupling model. Parameter

Unit

Value

Absolute permeability Oil relative permeability Connate water saturation Grid length in x − direction Grid length in y− direction Grid length in z − direction Time step Porosity Initial pressure Wellbore pressure at the heel Tolerance

K kro (Swc ) Swc dx dy dz dt ϕ Pinit Pheel ε

[mD ]

m m m Day [bar ] [bar ] %

20 0.7 0.1 50 50 50 0.5 0.2 300 250 1

in Friction Model, which is usually insignificant in low permeability reservoir and low productivity wells while significant for high productivity wells. In the first example, the Reservoir Model (a commercially available finite difference simulator) is coupled with two Well Models using the coupling scheme presented in Section 3. This is used to validate the methodology in a low permeability reservoir. In the second example, it demonstrates the applicability of this coupling by modeling a horizontal well with monotonic decreased skin from the heel to the toe. The basic information used in both the Reservoir Model and the Well Mode are shown in Table 2. In the Reservoir Model, the number of grid blocks is 20 × 10 × 5 in x , y and z directions, and the horizontal well is located in the middle layer of the reservoir with all 20 well blocks completed open hole. Since the Cartesian and Cylindrical grids are used in the Reservoir Model and the Well Model, respectively, the volume of the near-well segment is made the same with the volume of the well block. In the beginning, initial reservoir pressure is pinit = 300 bar and wellbore pressure is 250 bar . Using the coupling scheme in Section 3, the reservoir model is coupled with the two well models on each time step until convergence tolerance, ε < 1%, is satisfied. Then the results of pressure and saturation are recorded and used as the starting point for the next time step. In the first example, the iteration process for the first time step is shown in Figs. 24 and 25 for wellbore pressure and well inflow per segment, respectively. At the first time step, the number of iteration steps is usually larger (typically 50) than for the later time steps (typically 5) since the initial wellbore pressure is constant as frictional loss is ignored on the first time step. However, not many iterations steps are need to achieve a relative quick convergence at later times. The results in the first time step from both the Well/near-well model and the Friction model are compared in terms of wellbore pressure (Fig. 26), specific well inflow per segment (Fig. 27) and well block pressure (Fig. 28). The difference between the two well models are indistinguishable, with largest relative error being less than 1%. However, it

Fig. 25. The well inflow in different iterations for the first time step.

Fig. 26. Comparison of wellbore pressure between well models for the first time step.

also means that sharing information between the two models without iteration can hardly yield accurate results. Continue the coupling for 5 time steps, the total well productivity from both models are compared in Fig. 29. The two coupling results show a good agreement with each other. This verifies the effectiveness of the iterative coupling scheme and also the accuracy of the Well/near-well model. It also indicates that in low permeability reservoirs, the axial flow rate is small and can usually be ignored; i.e. the well inflow is the dominant flow in the near648

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 27. Comparison of well inflow between well models for the first time step. Fig. 30. Skin factors and axial permeabilities along the wellbore.

damage skin factor by Furui et al. (2003) is applied; i.e.

⎧ 1 ⎡ rdh (x ) k s (x ) = ⎛ − 1⎞ ln + ⎢ ⎨ + k ( x ) 1 I rw ani ⎝ d ⎠ ⎩ ⎣ ⎜

⎟

2

⎫ ⎛ rdh (x ) ⎞ + I 2 − 1 ⎤ ani ⎥⎬ r ⎝ w ⎠ ⎦⎭

⎜

⎟

(12) k

Here, given k = 3; maximum damaged zone radius rdh = 5 inch; d anisotropy ratio Iani = 1; permeability K = 200 mD , the skin factor distribution displayed in Fig. 30. The skin factors are used in the reservoir model and well models. With the same procedures as used in the first example, the results from coupling the Reservoir Model with two Well Models are compared for the first time step in terms of wellbore pressure (Fig. 31), well inflow per segment (Fig. 32) and well block pressure (Fig. 33). Continuing both couplings for 5 time steps, the total well productivity from both models are compared in Fig. 34. The total productivity shows clear difference with around 7 − 8% relative errors in each time step. The new well model results in a lower total productivity due to extra nearwell reservoir pressure drop in axial direction for formation damaged reservoirs. Although the above examples are for single phase flow, the concept is not limited to this. More specifically, multi-phase flow can be implemented in the well/near-well modelling using relative permeability in Eqs. (1)–(5) and Drift-Flux model (Shi et al., 2005) in Eq. (8). Any numerical reservoir simulator methods can be implemented in the iterative coupling scheme presented in this paper. As an example of

Fig. 28. Comparison of well block pressure between well models for the first time step.

Fig. 29. Comparison of total productivity between well models for various time steps.

well region. In the second example, unevenly distributed skin factors are introduced in the near-well region. The non-uniform skin distribution is usually caused by varying mud invasion along the wellbore, and is usually assumed to be cone shaped (Frick and Economides (1993) and Furui et al. (2003)). The physical process causing formation damage can be simulated by Lohne et al. (2010). Here, the final formation

Fig. 31. Comparison of wellbore pressure between well models with unevenly distributed skin factors. 649

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Table 3 Reservoir, wellbore and fluids properties. Parameter

Unit

Value

Reservoir length, Lres Reservoir width, Wres Reservoir thickness, Dres Well length, L w Oil viscosity, μo Water viscosity, μ w Wellbore radius, rw Oil density, ρ

m m m m cp cp m

1000 1000 1000 500 1.0 1.0 0.1 800

kg / m3

Water density, ρ Injector pressure, Pinj Producer bottom hole pressure, PBH Connate water saturation, Swc Residual oil saturation, Swc Maximum water relative permeability, aw Maximum oil relative permeability, aw Corey's model exponent for water, nw Corey's model exponent for oil, no

Fig. 32. Comparison of well inflow rates between well models with unevenly distributed skin factors.

kg / m3 bar

800

bar

100 0.15 0.15 1 1 2 2

350

flow, the Corey's relative permeability model is used, i.e. n

n

s w − Swc ⎞ w 1 − s w − Sor ⎞ o krw (s w ) = a w ⎛ ; kro (s w ) = ao ⎛ ; 1 − S − S 1 wc or ⎝ ⎠ ⎝ − Swc − Sor ⎠ ⎜

⎟

⎜

⎟

(13)

where, Swc and Sor are the connate water and residual oil saturation, respectively; a w and ao are maximum relative permeabilities for water and oil; n w and no are exponents for water and oil, as shown in Table 3. Using the coupling scheme described in Section 3, the pressure in the reservoir and well/near-well domain is first solved. Then, the streamlines are generated based on the pressure solution using a higher order velocity interpretation method (Zhang et al., 2018) than the classical Pollock's method (Pollock, 1988). Subsequently, the two phase flow problem are solved analytically along the streamlines, i.e. the fluid saturations are moved along the streamlines. The streamline simulation methods have been proved with high efficiency in heterogeneous reservoirs (Datta-Gupta and King, 2007) and widely applied as an alternative to the conventional finite difference reservoir simulators (Thiele et al., 2010). The results for a time step before water breakthrough are shown in Figs. 35 and 36, for heterogeneous reservoirs with and without a shale layer. The shale layer is impermeable, representing another type of reservoir heterogeneity. The horizontal well is producing at the top of the heterogeneous reservoir, with a vertical injector at the bottom corner. The streamlines clearly show the flow path from the injector to the horizontal producer. The color of the streamline represents the water saturation, where in the middle the sharp jump represents the front water saturation. For the case with the shale layer, the streamlines indicate the early water breakthrough near the toe. For the case with impermeable shale layer, the injected water bypasses the shale layer into the producer.

Fig. 33. Comparison of well block pressure between well models with unevenly distributed skin factors.

7. Conclusion The numerical well/near-well model is presented together with applications in heterogeneous reservoirs. It is applied both as a standalone simulator, and in a coupled approach with a reservoir simulator. A strong iterative coupling schematics avoiding the use of partial derivatives is also presented for the well/near-well model and a reservoir model to represent transient single phase and two phase flows. To account for reservoir heterogeneity and formation damage effects, the well/near-well model is applied for heterogeneous reservoirs, cross flow cases, and non-uniform formation damage effects along horizontal wellbore. The results show that the model can handle a large range of permeability heterogeneity, detect and represent crossflow between well and near-well reservoirs, simulate horizontal well

Fig. 34. Comparison of total well productivity with unevenly distributed skin factors in various time steps.

such, a coupling using a streamline model for the two phase flow problem (Cao et al., 2018) in heterogeneous reservoir is demonstrated. The basic information about the reservoir, the wellbore and the fluids is shown in Table 3. The reservoir permeability is heterogeneous, using randomly generated values between 5 and 500 mD . For the two phase 650

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Fig. 35. Coupled well/near-well model and streamline model for a heterogeneous reservoir.

Fig. 36. Coupled well/near-well model and streamline model for a heterogeneous reservoir with a shale layer.

formation damage accurately with unevenly distribution local skin factors. The results demonstrate the stability and accuracy of the iterative coupling scheme between the well/near-well reservoir model and a reservoir model. They also indicate that the new iterative coupling between well and reservoir model is necessary to achieve fast convergence for transient cases. The well/near-well model has not been applied for multi-stage fractured wells, in which the fractures dominate the flow with several orders of magnitude larger permeability than the formations. This, together with more field applications with complicated reservoirs and wells constitute parts of the future research work.

Cao, J., James, L.A., Johansen, T.E., 2015. A new coupled axial-radial productivity model with application to high order numerical well modeling. In: Presented at SPE Reservoir Characterisation and Simulation Conference and Exhibition, Abu Dhabi, UAE, September. Cao, J., Zhang, N., James, L.A., Johansen, T.E., 2018. Modeling semi-steady state nearwell flow performance for horizontal wells in anisotropic reservoirs. Comput. Geosci. 22 (3), 725–744. Colebrook, C.F., 1939. Turbulent flow in pipes, with particular reference to the transition region between the smooth and rough pipe laws. J. ICE 11 (4), 133–156. Datta-Gupta, A., King, M.J., 2007. Streamline Simulation: Theory and Practice. SPE Textbook Series, vol 11 Society of Petroleum Engineers, Richardson, TX. Dikken, B., 1990. Pressure drop in horizontal wells and its effect on production performance. J. Petrol. Technol. 42 (11), 1426–1433. Economides, M.J., Brand, C.W., Frick, T.P., 1996. Well configurations in anisotropic reservoirs. SPE Form. Eval. 11, 257–262. Fokker, P.A., Borello, E.S., Serazio, C., Verga, F., 2012. Estimating reservoir heterogeneities from pulse testing. J. Petrol. Sci. Eng. 86, 15–26. Frick, T.P., Economides, M.J., 1993. Horizontal well damage characterization and removal. SPE Prod. Facil. 8 (01), 15–22. Furui, K., Zhu, D., Hill, A.D., 2003. A rigorous formation damage skin factor and reservoir inflow model for a horizontal well (includes associated papers 88817 and 88818). SPE Prod. Facil. 18 (03), 151–157. Grubert, M.A., Wan, J., Ghai, S.S., Livescu, S., Brown, W.P., Long, T.A., 2009. Coupled completion and reservoir simulation technology for well performance optimization. In: SPE/EAGE Reservoir Characterization & Simulation Conference. Holmes, J.A., 2001. Modeling advanced wells in reservoir simulation. J. Petrol. Technol. 53 (11), 54–66. Holmes, J.A., Barkve, T., Lund, O., 1998. Application of a multisegment well model to simulate flow in advanced wells. In: Paper SPE 50646 Presented at the European

Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.petrol.2019.01.091. References Babu, D.K., Odeh, A., 1989. Productivity of a horizontal well. SPE Reservoir Eng. 4 (4), 417–421. Brekke, K., Johansen, T.E., Olufsen, R., 1993. A new modular approach to comprehensive simulation of horizontal wells. In: Paper SPE 26518 Presented at SPE Annual Technical Conference and Exhibition, 3-6 October, Houston, Texas.

651

Journal of Petroleum Science and Engineering 176 (2019) 640–652

J. Cao et al.

Parn-anurak, S., Engler, T.W., 2005. Modeling of fluid filtration and near-wellbore damage along a horizontal well. J. Petrol. Sci. Eng. 46 (3), 149–160. Peaceman, D.W., 1983. Interpretation of well-block pressures in numerical reservoir simulation with nonsquare grid blocks and anisotropic permeability. SPE J. 23 (03), 531–543. Pedrosa Jr., O.A., Aziz, K., 1986. Use of a hybrid grid in reservoir simulation. SPE Reservoir Eng. 1 (06), 611–621. Penmatcha, V., Aziz, K., 1999. Comprehensive reservoir/wellbore model for horizontal wells. SPE J. 4 (3), 224–234. Pollock, D.W., 1988. Semianalytical computation of path lines for finite difference models. Gr. Water 26 (6), 743–750. Ranjith, R., Suhag, A., Balaji, K., Putra, D., Dhannoon, D., Saracoglu, O., et al., 2017. Production optimization through utilization of smart wells in intelligent fields. In: SPE Western Regional Meeting. Renard, G., Dupuy, J.M., 1991. formation damage effects on horizontal-well flow efficiency (includes associated papers 23526 and 23833 and 23839). J. Petrol. Technol. 43 (07), 786–869. Schulkes, R.M., King, A.C., 2002. Boundary-layer behaviour in a hydraulic theory of horizontal oil wells. J. Eng. Math. 42 (1), 23–44. Shi, H., Holmes, J.A., Durlofsky, L.J., Aziz, K., Diaz, L., Alkaya, B., Oddie, G., 2005. Driftflux modeling of two-phase flow in wellbores. SPE J. 10 (01), 24–33. Sun, H., Vega, S., Tao, G., 2017. Analysis of heterogeneity and permeability anisotropy in carbonate rock samples using digital rock physics. J. Petrol. Sci. Eng. 156, 419–429. Tabatabaei, M., Ghalambor, A., 2011. A new method to predict performance of horizontal and multilateral wells. SPE Prod. Oper. 26 (01), 75–87. Tamayo-Mas, E., Mustapha, H., Dimitrakopoulos, R., 2016. Testing geological heterogeneity representations for enhanced oil recovery techniques. J. Petrol. Sci. Eng. 146, 222–240. Thiele, M.R., Batycky, R.P., Fenwick, D.H., 2010. Streamline simulation for modern reservoir-engineering workflows. J. Petrol. Technol. 62 (1), 64–70. Wan, J., Dale, B.A., Ellison, T.K., Benish, T.G., Grubert, M., 2008. Coupled well and reservoir simulation models to optimize completions design and operations for subsurface control. In: Europec/EAGE Conference and Exhibition, February. Society of Petroleum Engineers. Zhang, N., Cao, J., James, L.A., Johansen, T.E., 2018. Streamline tracing methods based on piecewise polynomial pressure approximations. Water Resour. Res. 54 (9), 7029–7045.

Petroleum Conference, 20-22 October, The Hague, Netherlands. Holmes, J.A., Byer, T.J., Edwards, D.A., Stone, T.W., Pallister, I., Shaw, G.J., Walsh, D., 2010. A unified wellbore model for reservoir simulation. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Jain, R., Syal, S., Long, T.A., Wattenbarger, R.C., Kosik, I.J., 2013. An integrated approach to design completions for horizontal wells for unconventional reservoirs. SPE J. 18 (06), 1–026. Johansen, T.E., Khoriakov, V., 2007. Iterative techniques in modeling of multi-phase flow in advanced wells and the near well region. J. Petrol. Sci. Eng. 58 (1), 49–67. Johansen, T.E., James, L.A., Cao, J., 2015. Analytical coupled axial and radial productivity model for steady-state flow in horizontal wells. Int. J. Petrol. Eng. 1 (4), 290–307. Johansen, T.E., Cao, J., James, L.A., 2016. Analytical coupled axial-radial semi-steady state productivity model for horizontal wells in anisotropic medium. Int. J. Petrol. Eng. 2 (3), 225–245 2016. Joshi, S.D., 1988. Augmentation of well productivity with slant and horizontal wells (includes associated papers 24547 and 25308). J. Petrol. Technol. 40 (06), 729–739. Karimi-Fard, M., Durlofsky, L., 2012. Accurate resolution of near-well effects in upscaled models using flow-based unstructured local grid refinement. SPE J. 17 (04), 1–084. Khoriakov, V., Johansen, A.C., Johansen, T.E., 2012. Transient flow modeling of advanced wells. J. Petrol. Sci. Eng. 86, 99–110. Krogstad, S., Durlofsky, L.J., 2007. Multiscale mixed-finite-element modeling of coupled wellbore/near-well flow. In: Paper SPE 106179, presented at the SPE Reservoir Simulation Symposium held in Houston, Texas, U.S.A., 26–28 February. Kurtoglu, B., Medeiros, F., Ozkan, E., Kazemi, H., 2008. Semianalytical representation of wells and near well flow convergence in numerical reservoir simulation. In: SPE Annual Technical Conference and Exhibition. Livescu, S., Brown, W.P., Jain, R., Grubert, M., Ghai, S.S., Lee, L.B.W., Long, T., 2010. Application of a coupled wellbore/reservoir simulator to well-performance optimization. In: SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers. Lohne, A., Han, L., Van Der Zwaag, C., van Velzen, H., Mathisen, A.M., Twyman, A., Hatzignatiou, D.G., 2010. Formation-damage and well-productivity simulation. SPE J. 15 (03), 751–769. Moody, L.F., 1947. An approximate formula for pipe friction factors. Trans. ASME 69, 1005–1006. Palagi, C.L., Aziz, K., 1994. Use of Voronoi grid in reservoir simulation. SPE Adv. Technol. 2 (02), 69–77.

652