- Email: [email protected]
Method for calculating production indices of multilayer water drive reservoirs
Journal of Petroleum Science and Engineering 75 (2010) 66–70
Contents lists available at ScienceDirect
Journal of Petroleum Science and Engineering j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / p e t r o l
Method for calculating production indices of multilayer water drive reservoirs Chuanzhi Cui ⁎, Xiaoyan Zhao College of Petroleum Engineering, China University of Petroleum, Dongying, Shandong 257061, China
a r t i c l e
i n f o
Article history: Received 5 May 2009 Accepted 11 October 2010 Available online 21 October 2010 Keywords: water drive reservoir multilayer reservoir production indices mathematical model calculation method
a b s t r a c t Based on the dual-phase unstable flow theory of oil and water, the production indices calculating method of multilayer water drive reservoirs are established, which takes the changes of rock permeability and oil viscosity into consideration. The flowing resistance and production pressure drop of every sublayer is taken as the weighting coefficient in the allocation of well productivity. The solution of pressure equation is obtained by the implicit difference equation. And the water saturation is calculated by the streamline model, which enhanced the stability of solution. Using this model the production indices such as water cut and the degree of reserve recovery can be calculated for the development unit of multilayer reservoir at different development stage. The calculation method is applied into field example and the result shows the good agreement with the actual case. This research provides an important measure for the optimization of layers recombination modes. © 2010 Elsevier B.V. All rights reserved.
1. Introduction Continental sedimentary reservoirs are characterized by many longitudinal layer divisions and strong inter-layer heterogeneity. Due to the limitation of oil production technology, we generally consider several or dozens of adjacent layers as one development unit and share the same well pattern for production, which causes prominent conflict between layers and makes it difficult to further improve production efficiency. As production proceeds, it becomes difficult to develop the potentiality of remaining oil effectively by applying conventional well-pattern adjusting technology to the development unit during later period of production, so subdivision or recombination of layers to new development units becomes an in-depth production technology to further improve oil recovery. Different oil fields follow different standards for recombination and subdivision in accordance with geological conditions of their reservoirs and water flooding characteristics. Besides, recombination modes of layers is an optimization process during later period of production. Apart from original physical parameters of reservoir, the variation of reservoir parameters during production, the distribution rules of the remaining oil, production potentiality at each layer and the well-pattern conditions should be considered. Based on those factors, various recombination modes of layers were established, and the optimal mode was determined based on the predicted production indices. Anbarci et al. (1989) developed an analytical solution of a multilayered composite model to analyze front locations in different
⁎ Corresponding author. E-mail addresses: [email protected] (C. Cui), [email protected] (X. Zhao). 0920-4105/$ – see front matter © 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.petrol.2010.10.003
layers when secondary and tertiary recovery processes are performed. In their model a pseudo-steady state approximation for interlayer flow was used, which transforms a two-dimensional problem into a onedimensional problem. Gomes and Ambastha (1993) extended the model and demonstrated that it can be used for a general n-layer composite reservoir. Atsushi Sagawa et al. (2000) adopted the model developed by Gomes and Ambastha (1993) to calculate the pressure response of reservoirs including the high permeability lens. E.P. Lolon et al. (2008) validates the concept of developing explicit solutions in the real domain for modeling the performance of a multilayer reservoir system. This process provides solutions for multilayer reservoir systems that can be used to model well test or production data. Ahmed Aly and W.J. Lee (1995) present a description and verification of a new semianalytical simulator used to model the wellbore performance of multilayered reservoirs with unequal initial pressures. The semi-analytical simulator proved to be faster than a conventional three-dimensional, finite-difference commercial simulator. In the above literatures the single phase fluid flow is assumed and only the pressure performance is analyzed. In the calculation of production indices of oil–water two-phase flow in multilayer reservoir system, generally a reservoir is simplified, and oil well rows and water well rows are simplified as outlet sap and water injector line. And then, one-dimensional flow formula is applied for calculation. Some presumptions required in the earliest equivalent flowing resistance method do not fit with real conditions, and calculation methods based on single-layer one-dimensional dualphase and multi-layer one-dimensional dual-phase have been proposed in 1960s. In single-layer one-dimensional dual-phase calculation, multiple layers are combined into one (Jiang et al., 2001), so it cannot reflect the impact of inter-layer heterogeneity. Lang (1991) showed that
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
the Stiles method (1949) is frequently applied in multi-layer onedimensional dual-phase. However, the method is based on the presumption that the fluid flow in each sublayer is piston-type flow, which is inconsistent with the actual situation. At present the technology of reservoir numerical simulation is often used to calculate the production indices after layers recombination (Sun, 2005). However, as a reservoir has many layers, it is difficult to build geological and dynamic reservoir models, and reservoir simulation costs much time and labor. Besides, the variation of reservoir parameters was ignored in present reservoir numericalsimulation software. Even if the variation of reservoir parameters is considered in the calculation by means of restarting method, the operation is still complicated, and some error may exist in calculation results. In this paper, a method for calculating production indices of multilayer reservoirs under the unstable flow of oil–water phases are introduced. The compressibility of rock and fluid as well as the changes of permeability and crude-oil viscosity are taken into account.
67
3. Building and solving of mathematical model based on one-dimensional one-way flow 3.1. Mathematical model for one-dimensional one-way flow The continuity equations of oil–water phases under one-dimensional one-way flow are: ∂ kk ∂P ∂ ðϕρo So Þ ρo ro + qo = μ o ∂x ∂x ∂t
ð1Þ
∂ kk ∂P ∂ ðϕρw Sw Þ ρw rw + qw = μ w ∂x ∂x ∂t
ð2Þ
Auxiliary equation : So + Sw = 1
ð3Þ
Initial conditions : pðxÞj t = 0 = pi ; Sw ðxÞj t = 0 = Swc
ð4Þ
Outer boundary conditions: the outer boundary of one-dimensional profile model is regarded as closed, i.e., 2. Reservoir model building and calculating idea It is presumed that a production unit is composed of number n of sublayers, and each sublayer has its own parameters such as permeability, thickness, porosity and formation crude-oil viscosity. Thus, an injection–production profile model is built, as shown in Fig. 1. For the convenience of production indices calculation after combination of different sublayers, the fluid flow in multilayer profile model is regarded as the combination of one-dimensional one-way flows in n sublayers. The pressure, saturation distribution and production indices of single layer were calculated. The production indices of each layer in one same production unit at the same moment were added up. Thus, the production indices of one same production unit can be obtained. It is presumed that dual-phase flow of oil and water exists in layers, and fluid flow follows Darcy linear flow law. Both rock and fluid are compressible, and no vertical convection occurs among layers. On this basis, mathematical model was built to calculate production indices, with the influence of capillary force and gravity being ignored.
∂P ∂x
j
=0
ð5Þ
=0
ð6Þ
x=0
and ∂P ∂x
j
x=L
Inner boundary conditions can be classified into two types: one is that the well productivity is fixed, and the other is that the well bottom-hole pressure is fixed. In the above formulas, p, S,ρ, μ, q, k and kr represent pressure, saturation, density, viscosity, productivity, absolute permeability and relative permeability, respectively (subscript o and w represent oil phase and water phase, respectively). ϕ is the porosity. 3.2. Calculation of pressure and saturation distribution 3.2.1. Calculation of pressure distribution With the item on the right of oil and water mass conservation equation being expanded, the following expressions are obtained: ∂ kk ∂P ∂P ∂S ρo ro + qo = βo + ϕρo o μ o ∂x ∂t ∂x ∂t
ð7Þ
∂ kk ∂P ∂P ∂S ρw rw + qw = βw + ϕρw w μ w ∂x ∂t ∂x ∂t
ð8Þ
where βo = ϕρoSo(co + cϕ) and βw = ϕρwSw(cw + cϕ) If A = ρo/ρw and (7) + (8) × A, So and Sw are eliminated to obtain a pressure equation containing only variable P: ∂ ∂P ∂ ∂P ∂P λo +A λw + qo + Aqw = ðβo + Aβw Þ ∂x ∂x ∂x ∂x ∂t
ð9Þ
where λo = ρo kkμ oro and λw = ρw kkμwrw : For Eq. (9), block-centered grids are employed for differentiation to obtain implicit difference equation, with which and inner and outer boundary conditions, equation set are obtained. Chasing algorithm is used for solution, and pressure distribution can be obtained.
Fig. 1. Multi-layer one-dimension flow model.
3.2.2. Calculation of water saturation distribution The streamline method is often used in reservoir simulation (Cheng et al., 2004; Mallison, et al., 2004; Bhambri and Mohanty, 2008; Maschio and Schiozer, 2003), which enhanced the stability and convergence of saturation calculation. Here the streamline model is
68
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
applied for saturation calculation. It is assumed that a streamline exists between oil and water wells. The injected water starts from the water well and flows into the oil well. And the streamline has equivalent section area A and at a site l length away from the water injector the porosity is ϕ(l), the water saturation is Sw(l), the flow velocity is v(l), and the water cut is fw(Sw). The mass conservation equation on stream line (l,l + Δl) within (t,t + Δt) is that water accumulation equals to the difference of water influx with water outflux, i.e.: 0 h 0 i 0 h 0 0 i 0 = ϕ l AΔl Sw l ; t + Δt −Sw l ; t ΔtAv l fw l; t −fw l + Δl; t
ð10Þ where l ' ∈ [l, l + Δl] and t ' ∈ [t, t + Δt] are the mean value obtained according to mean value theorem. With limit value being taken from formula (10), Buckley–Levertt equation is obtained: ∂Sw V ðlÞ ∂fw + =0 ϕðlÞ ∂l ∂t .
l ϕðlÞ
If τ = ∫
l0
V ðl Þ
dl; then there exists :
ð11Þ
∂fw ðt−t0 Þ dSw
∂Sw ∂f + w =0 ∂t ∂τ
ð12Þ
ð13Þ
For τ at any position and any time t, Sw(τ, t)can be obtained through formula (13). 4. Introduction of parameter variation in water drive reservoirs Wu (2006) showed the variation rule of rock parameters washed by long-term injected water and its impact on development effect. Cui and Zhao (2004) developed the reservoir numerical simulation with the variation of rock parameters. Present research work (Wu, 2006; Deng and Wu, 1996; You et al., 2007) concerning the variation laws of reservoir parameters involves three types, i.e., reservoir parameters vary with pressure, water-passing multiple or water cut. Here, the variation of reservoir parameters following water cut is introduced for analysis. Water cut is obtained by means of calculating based on the water saturation at each grid in every layer, and the parameters of rock and fluid after variation are obtained through calculation based on the relation of rock & fluid parameters with water cut. The relationships of rock and fluid parameters with water cut in the second member of the Shahejie formation in the second district of Shengtuo Oilfield were adopted (Cui and Zhao, 2004). The variation multiple of permeability vs. water cut is expressed as: kc = 1:0733 + 0:0034fw
ð14Þ
The variation law of formation crude-oil viscosity with different water cut is expressed as: 0:0122fw
μo = μoi e
5.1. Allocation of well productivity at each sublayer In the method to treat well item in reservoir numerical simulation (Chen, 1989; Fanchi, 2006), the allocation of well productivity is, under fixed total well liquid production rate and water injection rate, treated based on the parameters of grid where a well is in, and the impact of the flow resistance in the whole formation cannot be considered. Here, the flow resistance and pressure drop of a whole sublayer between injection and production wells are taken into account to allocate the well productivity at each sublayer in one same production unit. Suppose there are number N of sublayers in one production unit, and each sublayer is divided into number m of grids. m
If PIk = ∑
j=1
∂fw = dS . It is presumed that its solution (i.e., As fw = fw ðSw Þ, dτ dt w characteristic line) is τ = τc ðt Þ. On the characteristic line, Sw ðτc ðt Þ; t Þ is unrelated to time t, so ∂ fw/d Sw is also unrelated to t on the characteristic line and is a function of initial condition Sw ðτc ðt0 Þ; t0 Þ. Then
τ = τc ðt0 Þ +
5. Treatment of well productivity and bottom-hole flowing pressure
ð15Þ
where fw is the water cut (%); kc is the variation multiple of permeability; μoi is formation crude-oil viscosity at the beginning of production (mPa s) and μoi = 18 mPa s in the second member of the Shahejie formation in the second district of Shengtuo Oilfield.
kkro kk + rw hkj Pwi −Pwf kj μo μ w kj
ð16Þ
If Q l is the total single-well liquid productivity in one production unit, then the liquid productivity at sublayer k is: Q lk =
PIk
. N
∑ PIj
⋅Q l
ð17Þ
j=1
where Pwi and Pwf represent the bottom-hole flowing pressure of water injection well and production well at sublayer No. k. The oil production rate and water production rate at sublayer k of an oil well is distributed based on the mobility ratio of the grid where the well is in, i.e.: Q ok =
λo λo + λw
Q lk
ð18Þ
k
Q wk = Q lk −Q ok
ð19Þ
Usually one set of well pattern corresponds to one flow unit. No matter how many sublayers that one flow unit possesses, the above method about the allocation of historical production is still applicable. In the condition of complication with commingling different flow units, each flow unit has independent well pattern. Then, a pair of injection–production wells can be simulated for each flow unit, and the allocation of well productivity is calculated based on the rock & fluid parameters at each sublayer in one same flow unit. 5.2. Calculation of bottom-hole flowing pressure at each sublayer Under given bottom-hole flowing pressure, it is possible to combine layers with larger vertical intervals since layers recombination technique breaks the prior principle of adjacency combination. Therefore, the influence of pressure difference in the shaft due to inter-layer depth difference on production should be considered in calculation. According to the pressure drop expression of vertical pipe flow, the formula for calculating inter-layer pressure drop can be obtained:
Δpvi = ρm gΔzi + λm
ρm vm i 2 Δzi 2D
ð20Þ
where Δpvi is the pressure drop between two layers; Δzi is the depth 64 difference between two layers; λm = Re ; ρm is the density of fluid m mixture, and ρm = ρo fo + ρwfw, where fo is the volume oil cut, fw is the
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
69
Table 1 Reservoir parameters of the 4th–6th sand groups. Parameters
Value
Bearing area (km2) Geological reserves (104 t) Average effective thickness (m) Reservoir depth (m) Average effective permeability (μm2) Porosity (%) Original formation pressure (MPa) Bubble pressure (MPa) Formation oil viscosity (mPa s) Oil formation volume factor Surface density of crude oil (g/cm3)
7.4 1204 10.8 1879–2060 1.6 28 20.7 11.4 18 1.12 0.901
volume water cut; νmi is the mean flow velocity of fluid mixture when passing through cross section, and νmi = Q li = 14 πD2 ; Q li is the volume flow rate of oil–water mixture; D is the shaft diameter at vertical section. 6. Field example With the above calculation method being programmed, calculation of production indices was performed with the 4th–6th sand groups in the second member of the Shahejie formation as the example. The second member of the Shahejie Formation in the second district of Shengtuo Oilfield is in the east of Shengtuo oilfiled. The 4th– 6th sand groups in the second member of the Shahejie Formation are the delta plain facies reservoir, and the reservoir parameters are listed in Table 1. The 4th–6th sand groups consists of 18 sublayers, of which the 4th sand group is divided into five sublayers, the 5th sand group is divided into six sublayers, and the 6th sand group is divided into seven sublayers. In 1975, as an independent development unit, the 4th–6th sand groups was started to develop. Currently, the degree of reserve recovery has reached 35.49%, and comprehensive water cut reached 96.5%. A profile model was built according to the physical properties of eighteen sublayers and the producer-injector well spacing is 300 m. The physical parameters of each sublayer are shown in Table 2. Two curves of water cut vs. degree of reserve recovery in commingled production of the 4th–6th sand groups are showed in Fig. 2. In one curve the variation of permeability and crude-oil viscosity is considered, which is calculated based on the model
Fig. 2. Influence of variation in permeability and crude-oil viscosity on calculation results.
mentioned, and the variation of permeability and crude-oil viscosity vs. water cut was calculated according to formulas (14) and (15). In the other curve the variation of permeability and crude-oil viscosity is ignored. When the variation of permeability and crude-oil viscosity are ignored, the calculated degree of reserve recovery amounted to 42.97% and was 7.48 more than the actual degree of reserve recovery (35.49%), which is a large error. With the variation of permeability and crude-oil viscosity being considered, current comprehensive water cut was calculated as 95.8%, and degree of reserve recovery as 35.36%, which accords with current production situation, indicating that the method can well reflect the influence of inter-layer heterogeneity on production and achieves favorable accuracy and reliability. And this result also demonstrates that the variation of the rock parameters and fluid parameters should be considered during the calculation. Current water cut and degree of reserve recovery of each sublayer calculated based on the model are listed in Table 2. Based on the calculation results, the sublayers were recombined to two flow units. Main productive layers 41, 42, 52, 55 and 65 with better physical properties as well as higher water cut and degree of reserve recovery were combined into one flow unit, while others formed another flow unit. Production indices after layers recombination are shown in Fig. 3. As contrast, the production indices of regarding eighteen sublayers as one flow unit is also predicted, which is also showed in Fig. 3. It can be seen that the water cut began to decline after layers recombination. It was predicted that the degree of reserve recovery will increase by 2.85% in ten years after layers recombination. This calculation result also demonstrates that whether at start or in the late of a field development the method is applicable to the calculating production indices of multilayer water drive reservoirs.
Table 2 Physical parameters and calculated case of sublayers. Sublayer name 1
4 42 43 44 45 51 52 53 54 55 56 61 62 63 64 65 66 67
Depth (m)
Net thickness (m)
Permeability (10−3um2)
Porosity (f)
Geological reserves (104 t)
Water cut (%)
Degree of reserve recovery (%)
1879.6 1899.5 1912.7 1919.1 1930.9 1949.6 1955.2 1961.7 1976.4 1989.2 2003.5 2007.4 2016.3 2025.4 2029.0 2034.1 2045.8 2056.1
3.440 3.088 1.368 1.968 2.640 0.864 1.080 2.064 2.376 3.360 1.104 1.440 1.728 0.888 1.344 1.608 1.488 2.544
3041.5 3796.0 1084.0 925.0 1427.9 1100.9 2790.4 1329.7 1445.6 3512.0 1735.1 1581.6 615.8 903.7 1276.8 2767.1 1903.7 1443.0
0.274 0.286 0.278 0.271 0.261 0.276 0.267 0.277 0.270 0.288 0.270 0.288 0.230 0.255 0.279 0.275 0.284 0.268
198.9 160.0 95.3 69.5 60.4 19.1 49.5 82.3 46.4 97.6 38.8 48.3 67.6 7.6 28.6 38.3 39.5 56.5
97.12 97.78 87.2 85.93 91.43 87.38 96.75 90.12 91.28 97.58 92.43 91.46 81.39 86.37 89.47 96.67 93.11 91.32
42.44 45.49 25.77 24.04 30.86 26.09 41.38 28.81 30.51 44.04 33.73 30.92 21.23 24.5 28.09 40.72 34.53 30.60
70
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
ρ μ ϕ
Density, g/cm3 Viscosity, mPa s Porosity, fractional
Subscripts k Layer order number l Liquid phase or distance m Mixture or layer number o Oil phase w Water phase Fig. 3. Predicted production curve after layers recombination to two units.
References 7. Conclusion The following conclusions are drawn from this study: (1) Mathematical model and method for calculating production indices were developed for multilayer water drive reservoirs. During the process of solving the model, implicit difference equation was adopted for pressure calculation, and streamline model was adopted to calculate water saturation. The flow resistance and pressure drop of one sublayer was considered comprehensively during single-well productivity allocation. And the variation in permeability and formation crude-oil viscosity was considered so that calculation accuracy and reliability can be improved. (2) The model was implemented based on one-dimension multilayer profile model, without the influence of plane heterogeneity being taken into account. Improvement can be made in future to fit the calculation results of the model better with real conditions. And the method is only applicable for oil/water systems, so as soon as gas is liberated (due to pressure reduction) the method is unlikely to be applicable. Generally the reservoir pressure is higher than the bubble point pressure in water drive reservoirs, and the dominated flow is oil–water two phase flow, so the method can be applicable for most water drive reservoirs. Nomenclature fw h k kc kr L p Pwi Pwf q S t μoi v
Water cut, fractional Net thickness, cm Absolute permeability, um2 Variation multiple of permeability, dimensionless Relative permeability, fractional Distance from water well to oil well, cm Pressure, 10−1 MPa Water well bottom-hole flowing pressure, 10−1 MPa Oil well bottom-hole flowing pressure, 10−1 MPa Production rate, cm3/s Saturation, fractional Time, s Initial formation crude-oil viscosity, mPa s Velocity, cm/s
Aly, Ahmed, Lee, W.J., 1995. Computational modeling of multi-layered reservoirs with unequal initial pressures: development of a new pre-production well test. Paper SPE 29586 prepared for presentation at the SPE Rocky Mountain Regional/LowPermeability Reservoirs Symposium held in Denver, CO, U. S.A., March 20–22. Anbarci, K., Grader, A.S., Ertekin, T., 1989. Determination of front locations in a multilayer composite reservoir. Paper SPE 19799 presented at the 64th Annual Technical Conference and Exhibition, San Antonio, TX, Oct. 8–11. Bhambri, P., Mohanty, K.K., 2008. Two- and three-hydrocarbon phase streamline-based compositional simulation of gas injections. J. Petrol. Sci. Eng. 62, 16–27. Chen, Yue-ming, 1989. Base of reservoir numerical simulation. Petroleum University Press, Dongying. (In Chinese). Cheng, H., Kharghoria, A., He, Z., Datta-Gupta, A., 2004. Fast history matching of finite difference models using streamline-derived sensitivities. Paper SPE 89447 presented at SPE/DOE Symposium on Improved Oil Recovery, Tulsa, OK., April 17–21. Cui, Chuan-zhi, Zhao, Xiao-yan, 2004. The reservoir numerical simulation study with the variety of reservoir parameters. J. Hydrodynamics 19, 912–915 (Supplement, In Chinese). Deng, Yu-zhen, Wu, Su-ying, 1996. Study on the variation rules of reservoir parameters during water flooding. PGRE 3 (4), 44–52 (In Chinese). Fanchi, John R., 2006. Principles of applied reservoir simulation, Third edition. Gulf Professional Publishing. Gomes, E., Ambastha, A.K., 1993. An analytical pressure-transient model for multilayered composite reservoirs with pseudo steady-state formation crossflow. Paper SPE 26049 presented at the Western Regional Meeting Anchorage, AK, May 26–28. Jiang, Han-qiao, Yao, Jun, Jiang, Rui-zhong, 2001. Theory and method of petroleum engineer. Petroleum University Press, Dongying. (In Chinese). Lang, Zhao-xin, 1991. Base of reservoir engineering. Petroleum Industry Press, Beijing. (In Chinese). Lolon, E.P., Archer, R.A., llk, D., Blasingame, T.A., 2008. New semi-analytical solutions for multilayer reservoirs. Paper SPE 114946 presented at the CIPC/SPE Gas Technology Symposium 2008 Joint Conference Held in Calgary, Alberta, Canada, 16–19, Jun. Mallison, B.T., Gerritsen, M.G., Matringe, S.F., 2004. Improved mappings for streamline based simulation. Paper SPE 89352 presented at SPE/DOE Symposium on Improved Oil Recovery, Tulsa, OK, April 17–21. Maschio, Ce´lio, Schiozer, Denis Jose, 2003. A new upscaling technique based on Dykstra–Parsons coefficient: evaluation with streamline reservoir simulation. J. Petrol. Sci. Eng. 40, 27–36. Sagawa, Atsushi, Corbett, Patrick W.M., Davies, David R., 2000. Pressure transient analysis of reservoirs with a high permeability lens intersected by the wellbore. J. Petrol. Sci. Eng. 27, 165–177. Stiles, W.E., 1949. Use of permeability distribution in water flood calculations. Petrol. Trans. AIME 9–13 (Jan.). Sun, Guo, 2005. Investigation on optimized well pattern rearrangement in Shengtuo oil field during the ultra high water cut Period. PGRE 12 (3), 48–50. Wu, Su-ying, 2006. Variation rule of oil layer parameters washed by long-term injected water and its impact on development effect. Petrol. Geol. Oil Field Dev. Daqing 25 (4), 35–37 (In Chinese). You, Qi-dong, Zhou, Fang-xi, Zhang, Jian-liang, et al., 2007. Law and mechanism of parameters change in salutiferous reservoir. J. China Univ. Petrol. 31 (2), 79–82 (In Chinese).
Contents lists available at ScienceDirect
Journal of Petroleum Science and Engineering j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / p e t r o l
Method for calculating production indices of multilayer water drive reservoirs Chuanzhi Cui ⁎, Xiaoyan Zhao College of Petroleum Engineering, China University of Petroleum, Dongying, Shandong 257061, China
a r t i c l e
i n f o
Article history: Received 5 May 2009 Accepted 11 October 2010 Available online 21 October 2010 Keywords: water drive reservoir multilayer reservoir production indices mathematical model calculation method
a b s t r a c t Based on the dual-phase unstable flow theory of oil and water, the production indices calculating method of multilayer water drive reservoirs are established, which takes the changes of rock permeability and oil viscosity into consideration. The flowing resistance and production pressure drop of every sublayer is taken as the weighting coefficient in the allocation of well productivity. The solution of pressure equation is obtained by the implicit difference equation. And the water saturation is calculated by the streamline model, which enhanced the stability of solution. Using this model the production indices such as water cut and the degree of reserve recovery can be calculated for the development unit of multilayer reservoir at different development stage. The calculation method is applied into field example and the result shows the good agreement with the actual case. This research provides an important measure for the optimization of layers recombination modes. © 2010 Elsevier B.V. All rights reserved.
1. Introduction Continental sedimentary reservoirs are characterized by many longitudinal layer divisions and strong inter-layer heterogeneity. Due to the limitation of oil production technology, we generally consider several or dozens of adjacent layers as one development unit and share the same well pattern for production, which causes prominent conflict between layers and makes it difficult to further improve production efficiency. As production proceeds, it becomes difficult to develop the potentiality of remaining oil effectively by applying conventional well-pattern adjusting technology to the development unit during later period of production, so subdivision or recombination of layers to new development units becomes an in-depth production technology to further improve oil recovery. Different oil fields follow different standards for recombination and subdivision in accordance with geological conditions of their reservoirs and water flooding characteristics. Besides, recombination modes of layers is an optimization process during later period of production. Apart from original physical parameters of reservoir, the variation of reservoir parameters during production, the distribution rules of the remaining oil, production potentiality at each layer and the well-pattern conditions should be considered. Based on those factors, various recombination modes of layers were established, and the optimal mode was determined based on the predicted production indices. Anbarci et al. (1989) developed an analytical solution of a multilayered composite model to analyze front locations in different
⁎ Corresponding author. E-mail addresses: [email protected] (C. Cui), [email protected] (X. Zhao). 0920-4105/$ – see front matter © 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.petrol.2010.10.003
layers when secondary and tertiary recovery processes are performed. In their model a pseudo-steady state approximation for interlayer flow was used, which transforms a two-dimensional problem into a onedimensional problem. Gomes and Ambastha (1993) extended the model and demonstrated that it can be used for a general n-layer composite reservoir. Atsushi Sagawa et al. (2000) adopted the model developed by Gomes and Ambastha (1993) to calculate the pressure response of reservoirs including the high permeability lens. E.P. Lolon et al. (2008) validates the concept of developing explicit solutions in the real domain for modeling the performance of a multilayer reservoir system. This process provides solutions for multilayer reservoir systems that can be used to model well test or production data. Ahmed Aly and W.J. Lee (1995) present a description and verification of a new semianalytical simulator used to model the wellbore performance of multilayered reservoirs with unequal initial pressures. The semi-analytical simulator proved to be faster than a conventional three-dimensional, finite-difference commercial simulator. In the above literatures the single phase fluid flow is assumed and only the pressure performance is analyzed. In the calculation of production indices of oil–water two-phase flow in multilayer reservoir system, generally a reservoir is simplified, and oil well rows and water well rows are simplified as outlet sap and water injector line. And then, one-dimensional flow formula is applied for calculation. Some presumptions required in the earliest equivalent flowing resistance method do not fit with real conditions, and calculation methods based on single-layer one-dimensional dualphase and multi-layer one-dimensional dual-phase have been proposed in 1960s. In single-layer one-dimensional dual-phase calculation, multiple layers are combined into one (Jiang et al., 2001), so it cannot reflect the impact of inter-layer heterogeneity. Lang (1991) showed that
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
the Stiles method (1949) is frequently applied in multi-layer onedimensional dual-phase. However, the method is based on the presumption that the fluid flow in each sublayer is piston-type flow, which is inconsistent with the actual situation. At present the technology of reservoir numerical simulation is often used to calculate the production indices after layers recombination (Sun, 2005). However, as a reservoir has many layers, it is difficult to build geological and dynamic reservoir models, and reservoir simulation costs much time and labor. Besides, the variation of reservoir parameters was ignored in present reservoir numericalsimulation software. Even if the variation of reservoir parameters is considered in the calculation by means of restarting method, the operation is still complicated, and some error may exist in calculation results. In this paper, a method for calculating production indices of multilayer reservoirs under the unstable flow of oil–water phases are introduced. The compressibility of rock and fluid as well as the changes of permeability and crude-oil viscosity are taken into account.
67
3. Building and solving of mathematical model based on one-dimensional one-way flow 3.1. Mathematical model for one-dimensional one-way flow The continuity equations of oil–water phases under one-dimensional one-way flow are: ∂ kk ∂P ∂ ðϕρo So Þ ρo ro + qo = μ o ∂x ∂x ∂t
ð1Þ
∂ kk ∂P ∂ ðϕρw Sw Þ ρw rw + qw = μ w ∂x ∂x ∂t
ð2Þ
Auxiliary equation : So + Sw = 1
ð3Þ
Initial conditions : pðxÞj t = 0 = pi ; Sw ðxÞj t = 0 = Swc
ð4Þ
Outer boundary conditions: the outer boundary of one-dimensional profile model is regarded as closed, i.e., 2. Reservoir model building and calculating idea It is presumed that a production unit is composed of number n of sublayers, and each sublayer has its own parameters such as permeability, thickness, porosity and formation crude-oil viscosity. Thus, an injection–production profile model is built, as shown in Fig. 1. For the convenience of production indices calculation after combination of different sublayers, the fluid flow in multilayer profile model is regarded as the combination of one-dimensional one-way flows in n sublayers. The pressure, saturation distribution and production indices of single layer were calculated. The production indices of each layer in one same production unit at the same moment were added up. Thus, the production indices of one same production unit can be obtained. It is presumed that dual-phase flow of oil and water exists in layers, and fluid flow follows Darcy linear flow law. Both rock and fluid are compressible, and no vertical convection occurs among layers. On this basis, mathematical model was built to calculate production indices, with the influence of capillary force and gravity being ignored.
∂P ∂x
j
=0
ð5Þ
=0
ð6Þ
x=0
and ∂P ∂x
j
x=L
Inner boundary conditions can be classified into two types: one is that the well productivity is fixed, and the other is that the well bottom-hole pressure is fixed. In the above formulas, p, S,ρ, μ, q, k and kr represent pressure, saturation, density, viscosity, productivity, absolute permeability and relative permeability, respectively (subscript o and w represent oil phase and water phase, respectively). ϕ is the porosity. 3.2. Calculation of pressure and saturation distribution 3.2.1. Calculation of pressure distribution With the item on the right of oil and water mass conservation equation being expanded, the following expressions are obtained: ∂ kk ∂P ∂P ∂S ρo ro + qo = βo + ϕρo o μ o ∂x ∂t ∂x ∂t
ð7Þ
∂ kk ∂P ∂P ∂S ρw rw + qw = βw + ϕρw w μ w ∂x ∂t ∂x ∂t
ð8Þ
where βo = ϕρoSo(co + cϕ) and βw = ϕρwSw(cw + cϕ) If A = ρo/ρw and (7) + (8) × A, So and Sw are eliminated to obtain a pressure equation containing only variable P: ∂ ∂P ∂ ∂P ∂P λo +A λw + qo + Aqw = ðβo + Aβw Þ ∂x ∂x ∂x ∂x ∂t
ð9Þ
where λo = ρo kkμ oro and λw = ρw kkμwrw : For Eq. (9), block-centered grids are employed for differentiation to obtain implicit difference equation, with which and inner and outer boundary conditions, equation set are obtained. Chasing algorithm is used for solution, and pressure distribution can be obtained.
Fig. 1. Multi-layer one-dimension flow model.
3.2.2. Calculation of water saturation distribution The streamline method is often used in reservoir simulation (Cheng et al., 2004; Mallison, et al., 2004; Bhambri and Mohanty, 2008; Maschio and Schiozer, 2003), which enhanced the stability and convergence of saturation calculation. Here the streamline model is
68
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
applied for saturation calculation. It is assumed that a streamline exists between oil and water wells. The injected water starts from the water well and flows into the oil well. And the streamline has equivalent section area A and at a site l length away from the water injector the porosity is ϕ(l), the water saturation is Sw(l), the flow velocity is v(l), and the water cut is fw(Sw). The mass conservation equation on stream line (l,l + Δl) within (t,t + Δt) is that water accumulation equals to the difference of water influx with water outflux, i.e.: 0 h 0 i 0 h 0 0 i 0 = ϕ l AΔl Sw l ; t + Δt −Sw l ; t ΔtAv l fw l; t −fw l + Δl; t
ð10Þ where l ' ∈ [l, l + Δl] and t ' ∈ [t, t + Δt] are the mean value obtained according to mean value theorem. With limit value being taken from formula (10), Buckley–Levertt equation is obtained: ∂Sw V ðlÞ ∂fw + =0 ϕðlÞ ∂l ∂t .
l ϕðlÞ
If τ = ∫
l0
V ðl Þ
dl; then there exists :
ð11Þ
∂fw ðt−t0 Þ dSw
∂Sw ∂f + w =0 ∂t ∂τ
ð12Þ
ð13Þ
For τ at any position and any time t, Sw(τ, t)can be obtained through formula (13). 4. Introduction of parameter variation in water drive reservoirs Wu (2006) showed the variation rule of rock parameters washed by long-term injected water and its impact on development effect. Cui and Zhao (2004) developed the reservoir numerical simulation with the variation of rock parameters. Present research work (Wu, 2006; Deng and Wu, 1996; You et al., 2007) concerning the variation laws of reservoir parameters involves three types, i.e., reservoir parameters vary with pressure, water-passing multiple or water cut. Here, the variation of reservoir parameters following water cut is introduced for analysis. Water cut is obtained by means of calculating based on the water saturation at each grid in every layer, and the parameters of rock and fluid after variation are obtained through calculation based on the relation of rock & fluid parameters with water cut. The relationships of rock and fluid parameters with water cut in the second member of the Shahejie formation in the second district of Shengtuo Oilfield were adopted (Cui and Zhao, 2004). The variation multiple of permeability vs. water cut is expressed as: kc = 1:0733 + 0:0034fw
ð14Þ
The variation law of formation crude-oil viscosity with different water cut is expressed as: 0:0122fw
μo = μoi e
5.1. Allocation of well productivity at each sublayer In the method to treat well item in reservoir numerical simulation (Chen, 1989; Fanchi, 2006), the allocation of well productivity is, under fixed total well liquid production rate and water injection rate, treated based on the parameters of grid where a well is in, and the impact of the flow resistance in the whole formation cannot be considered. Here, the flow resistance and pressure drop of a whole sublayer between injection and production wells are taken into account to allocate the well productivity at each sublayer in one same production unit. Suppose there are number N of sublayers in one production unit, and each sublayer is divided into number m of grids. m
If PIk = ∑
j=1
∂fw = dS . It is presumed that its solution (i.e., As fw = fw ðSw Þ, dτ dt w characteristic line) is τ = τc ðt Þ. On the characteristic line, Sw ðτc ðt Þ; t Þ is unrelated to time t, so ∂ fw/d Sw is also unrelated to t on the characteristic line and is a function of initial condition Sw ðτc ðt0 Þ; t0 Þ. Then
τ = τc ðt0 Þ +
5. Treatment of well productivity and bottom-hole flowing pressure
ð15Þ
where fw is the water cut (%); kc is the variation multiple of permeability; μoi is formation crude-oil viscosity at the beginning of production (mPa s) and μoi = 18 mPa s in the second member of the Shahejie formation in the second district of Shengtuo Oilfield.
kkro kk + rw hkj Pwi −Pwf kj μo μ w kj
ð16Þ
If Q l is the total single-well liquid productivity in one production unit, then the liquid productivity at sublayer k is: Q lk =
PIk
. N
∑ PIj
⋅Q l
ð17Þ
j=1
where Pwi and Pwf represent the bottom-hole flowing pressure of water injection well and production well at sublayer No. k. The oil production rate and water production rate at sublayer k of an oil well is distributed based on the mobility ratio of the grid where the well is in, i.e.: Q ok =
λo λo + λw
Q lk
ð18Þ
k
Q wk = Q lk −Q ok
ð19Þ
Usually one set of well pattern corresponds to one flow unit. No matter how many sublayers that one flow unit possesses, the above method about the allocation of historical production is still applicable. In the condition of complication with commingling different flow units, each flow unit has independent well pattern. Then, a pair of injection–production wells can be simulated for each flow unit, and the allocation of well productivity is calculated based on the rock & fluid parameters at each sublayer in one same flow unit. 5.2. Calculation of bottom-hole flowing pressure at each sublayer Under given bottom-hole flowing pressure, it is possible to combine layers with larger vertical intervals since layers recombination technique breaks the prior principle of adjacency combination. Therefore, the influence of pressure difference in the shaft due to inter-layer depth difference on production should be considered in calculation. According to the pressure drop expression of vertical pipe flow, the formula for calculating inter-layer pressure drop can be obtained:
Δpvi = ρm gΔzi + λm
ρm vm i 2 Δzi 2D
ð20Þ
where Δpvi is the pressure drop between two layers; Δzi is the depth 64 difference between two layers; λm = Re ; ρm is the density of fluid m mixture, and ρm = ρo fo + ρwfw, where fo is the volume oil cut, fw is the
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
69
Table 1 Reservoir parameters of the 4th–6th sand groups. Parameters
Value
Bearing area (km2) Geological reserves (104 t) Average effective thickness (m) Reservoir depth (m) Average effective permeability (μm2) Porosity (%) Original formation pressure (MPa) Bubble pressure (MPa) Formation oil viscosity (mPa s) Oil formation volume factor Surface density of crude oil (g/cm3)
7.4 1204 10.8 1879–2060 1.6 28 20.7 11.4 18 1.12 0.901
volume water cut; νmi is the mean flow velocity of fluid mixture when passing through cross section, and νmi = Q li = 14 πD2 ; Q li is the volume flow rate of oil–water mixture; D is the shaft diameter at vertical section. 6. Field example With the above calculation method being programmed, calculation of production indices was performed with the 4th–6th sand groups in the second member of the Shahejie formation as the example. The second member of the Shahejie Formation in the second district of Shengtuo Oilfield is in the east of Shengtuo oilfiled. The 4th– 6th sand groups in the second member of the Shahejie Formation are the delta plain facies reservoir, and the reservoir parameters are listed in Table 1. The 4th–6th sand groups consists of 18 sublayers, of which the 4th sand group is divided into five sublayers, the 5th sand group is divided into six sublayers, and the 6th sand group is divided into seven sublayers. In 1975, as an independent development unit, the 4th–6th sand groups was started to develop. Currently, the degree of reserve recovery has reached 35.49%, and comprehensive water cut reached 96.5%. A profile model was built according to the physical properties of eighteen sublayers and the producer-injector well spacing is 300 m. The physical parameters of each sublayer are shown in Table 2. Two curves of water cut vs. degree of reserve recovery in commingled production of the 4th–6th sand groups are showed in Fig. 2. In one curve the variation of permeability and crude-oil viscosity is considered, which is calculated based on the model
Fig. 2. Influence of variation in permeability and crude-oil viscosity on calculation results.
mentioned, and the variation of permeability and crude-oil viscosity vs. water cut was calculated according to formulas (14) and (15). In the other curve the variation of permeability and crude-oil viscosity is ignored. When the variation of permeability and crude-oil viscosity are ignored, the calculated degree of reserve recovery amounted to 42.97% and was 7.48 more than the actual degree of reserve recovery (35.49%), which is a large error. With the variation of permeability and crude-oil viscosity being considered, current comprehensive water cut was calculated as 95.8%, and degree of reserve recovery as 35.36%, which accords with current production situation, indicating that the method can well reflect the influence of inter-layer heterogeneity on production and achieves favorable accuracy and reliability. And this result also demonstrates that the variation of the rock parameters and fluid parameters should be considered during the calculation. Current water cut and degree of reserve recovery of each sublayer calculated based on the model are listed in Table 2. Based on the calculation results, the sublayers were recombined to two flow units. Main productive layers 41, 42, 52, 55 and 65 with better physical properties as well as higher water cut and degree of reserve recovery were combined into one flow unit, while others formed another flow unit. Production indices after layers recombination are shown in Fig. 3. As contrast, the production indices of regarding eighteen sublayers as one flow unit is also predicted, which is also showed in Fig. 3. It can be seen that the water cut began to decline after layers recombination. It was predicted that the degree of reserve recovery will increase by 2.85% in ten years after layers recombination. This calculation result also demonstrates that whether at start or in the late of a field development the method is applicable to the calculating production indices of multilayer water drive reservoirs.
Table 2 Physical parameters and calculated case of sublayers. Sublayer name 1
4 42 43 44 45 51 52 53 54 55 56 61 62 63 64 65 66 67
Depth (m)
Net thickness (m)
Permeability (10−3um2)
Porosity (f)
Geological reserves (104 t)
Water cut (%)
Degree of reserve recovery (%)
1879.6 1899.5 1912.7 1919.1 1930.9 1949.6 1955.2 1961.7 1976.4 1989.2 2003.5 2007.4 2016.3 2025.4 2029.0 2034.1 2045.8 2056.1
3.440 3.088 1.368 1.968 2.640 0.864 1.080 2.064 2.376 3.360 1.104 1.440 1.728 0.888 1.344 1.608 1.488 2.544
3041.5 3796.0 1084.0 925.0 1427.9 1100.9 2790.4 1329.7 1445.6 3512.0 1735.1 1581.6 615.8 903.7 1276.8 2767.1 1903.7 1443.0
0.274 0.286 0.278 0.271 0.261 0.276 0.267 0.277 0.270 0.288 0.270 0.288 0.230 0.255 0.279 0.275 0.284 0.268
198.9 160.0 95.3 69.5 60.4 19.1 49.5 82.3 46.4 97.6 38.8 48.3 67.6 7.6 28.6 38.3 39.5 56.5
97.12 97.78 87.2 85.93 91.43 87.38 96.75 90.12 91.28 97.58 92.43 91.46 81.39 86.37 89.47 96.67 93.11 91.32
42.44 45.49 25.77 24.04 30.86 26.09 41.38 28.81 30.51 44.04 33.73 30.92 21.23 24.5 28.09 40.72 34.53 30.60
70
C. Cui, X. Zhao / Journal of Petroleum Science and Engineering 75 (2010) 66–70
ρ μ ϕ
Density, g/cm3 Viscosity, mPa s Porosity, fractional
Subscripts k Layer order number l Liquid phase or distance m Mixture or layer number o Oil phase w Water phase Fig. 3. Predicted production curve after layers recombination to two units.
References 7. Conclusion The following conclusions are drawn from this study: (1) Mathematical model and method for calculating production indices were developed for multilayer water drive reservoirs. During the process of solving the model, implicit difference equation was adopted for pressure calculation, and streamline model was adopted to calculate water saturation. The flow resistance and pressure drop of one sublayer was considered comprehensively during single-well productivity allocation. And the variation in permeability and formation crude-oil viscosity was considered so that calculation accuracy and reliability can be improved. (2) The model was implemented based on one-dimension multilayer profile model, without the influence of plane heterogeneity being taken into account. Improvement can be made in future to fit the calculation results of the model better with real conditions. And the method is only applicable for oil/water systems, so as soon as gas is liberated (due to pressure reduction) the method is unlikely to be applicable. Generally the reservoir pressure is higher than the bubble point pressure in water drive reservoirs, and the dominated flow is oil–water two phase flow, so the method can be applicable for most water drive reservoirs. Nomenclature fw h k kc kr L p Pwi Pwf q S t μoi v
Water cut, fractional Net thickness, cm Absolute permeability, um2 Variation multiple of permeability, dimensionless Relative permeability, fractional Distance from water well to oil well, cm Pressure, 10−1 MPa Water well bottom-hole flowing pressure, 10−1 MPa Oil well bottom-hole flowing pressure, 10−1 MPa Production rate, cm3/s Saturation, fractional Time, s Initial formation crude-oil viscosity, mPa s Velocity, cm/s
Aly, Ahmed, Lee, W.J., 1995. Computational modeling of multi-layered reservoirs with unequal initial pressures: development of a new pre-production well test. Paper SPE 29586 prepared for presentation at the SPE Rocky Mountain Regional/LowPermeability Reservoirs Symposium held in Denver, CO, U. S.A., March 20–22. Anbarci, K., Grader, A.S., Ertekin, T., 1989. Determination of front locations in a multilayer composite reservoir. Paper SPE 19799 presented at the 64th Annual Technical Conference and Exhibition, San Antonio, TX, Oct. 8–11. Bhambri, P., Mohanty, K.K., 2008. Two- and three-hydrocarbon phase streamline-based compositional simulation of gas injections. J. Petrol. Sci. Eng. 62, 16–27. Chen, Yue-ming, 1989. Base of reservoir numerical simulation. Petroleum University Press, Dongying. (In Chinese). Cheng, H., Kharghoria, A., He, Z., Datta-Gupta, A., 2004. Fast history matching of finite difference models using streamline-derived sensitivities. Paper SPE 89447 presented at SPE/DOE Symposium on Improved Oil Recovery, Tulsa, OK., April 17–21. Cui, Chuan-zhi, Zhao, Xiao-yan, 2004. The reservoir numerical simulation study with the variety of reservoir parameters. J. Hydrodynamics 19, 912–915 (Supplement, In Chinese). Deng, Yu-zhen, Wu, Su-ying, 1996. Study on the variation rules of reservoir parameters during water flooding. PGRE 3 (4), 44–52 (In Chinese). Fanchi, John R., 2006. Principles of applied reservoir simulation, Third edition. Gulf Professional Publishing. Gomes, E., Ambastha, A.K., 1993. An analytical pressure-transient model for multilayered composite reservoirs with pseudo steady-state formation crossflow. Paper SPE 26049 presented at the Western Regional Meeting Anchorage, AK, May 26–28. Jiang, Han-qiao, Yao, Jun, Jiang, Rui-zhong, 2001. Theory and method of petroleum engineer. Petroleum University Press, Dongying. (In Chinese). Lang, Zhao-xin, 1991. Base of reservoir engineering. Petroleum Industry Press, Beijing. (In Chinese). Lolon, E.P., Archer, R.A., llk, D., Blasingame, T.A., 2008. New semi-analytical solutions for multilayer reservoirs. Paper SPE 114946 presented at the CIPC/SPE Gas Technology Symposium 2008 Joint Conference Held in Calgary, Alberta, Canada, 16–19, Jun. Mallison, B.T., Gerritsen, M.G., Matringe, S.F., 2004. Improved mappings for streamline based simulation. Paper SPE 89352 presented at SPE/DOE Symposium on Improved Oil Recovery, Tulsa, OK, April 17–21. Maschio, Ce´lio, Schiozer, Denis Jose, 2003. A new upscaling technique based on Dykstra–Parsons coefficient: evaluation with streamline reservoir simulation. J. Petrol. Sci. Eng. 40, 27–36. Sagawa, Atsushi, Corbett, Patrick W.M., Davies, David R., 2000. Pressure transient analysis of reservoirs with a high permeability lens intersected by the wellbore. J. Petrol. Sci. Eng. 27, 165–177. Stiles, W.E., 1949. Use of permeability distribution in water flood calculations. Petrol. Trans. AIME 9–13 (Jan.). Sun, Guo, 2005. Investigation on optimized well pattern rearrangement in Shengtuo oil field during the ultra high water cut Period. PGRE 12 (3), 48–50. Wu, Su-ying, 2006. Variation rule of oil layer parameters washed by long-term injected water and its impact on development effect. Petrol. Geol. Oil Field Dev. Daqing 25 (4), 35–37 (In Chinese). You, Qi-dong, Zhou, Fang-xi, Zhang, Jian-liang, et al., 2007. Law and mechanism of parameters change in salutiferous reservoir. J. China Univ. Petrol. 31 (2), 79–82 (In Chinese).