# Study of Homogeneous Reservoir Pressure Inversion Model Based on Permeability Mechanics and Interpretation Software Design

Authors: Zhongshuai Chen; Hongjian Ni; Zhiqi Sun; Shiping Zhang; Qisong Wang
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

Well test analysis is required during the extraction of oil and gas wells. The information on formation parameters can be inverted by measuring the change in wellbore pressure at production start-up or after well shutdown. In order to calculate the characteristic parameters of the well, this paper creates a well test interpretation model for homogeneous reservoirs based on the theory of seepage mechanics, uses the Stehfest–Laplace inversion numerical inversion algorithm, and builds the Gringarten–Bourdet logarithmic curves model. The model can be used to evaluate the homogeneous reservoir. We use this model to design the pressure inversion interpretation software to implement a pressure inversion method based on permeability mechanics theory by using computer. The software can obtain the reservoir characteristic parameters such as permeability (K), skin coefficient (S), and wellbore storage coefficient (C). The homogeneous formation Gringarten–Bourdet curves data are available at https://github.com/JXLiaoHIT/Study-of-homogeneous-reservoir-pressure-inversion-model.


---

## Body

## 1. Introduction

Well test analysis is an important part of evaluating oil and gas reservoirs [1]. It leverages the downhole pressure measured by pressure sensors to invert the stratigraphic parameters information. It is possible to obtain, e.g., skin coefficient and wellbore storage coefficient. And further information on formation properties, reservoir reserves, and oil production capacity can be inference. Figure 1 illustrates the full process of well test analysis. To establish the relationship between measured pressure and stratigraphic parameters, we need to develop a variety of stratigraphic models, and the homogeneous reservoir model is a typical model and is the subject of this paper.Figure 1
A pipeline of well test analysis.As a typical type of oil reservoir, the homogeneous oil reservoir has also been explored by many scholars. Muskat proposed a method of using pressure recovery curves to estimate formation parameters in 1937 [2]. He pioneered the analysis of the well test, but the method can only do some qualitative analysis because the compressibility of the fluid is ignored. In 1949, Van Everdingen et al. used the Laplace integral transform to obtain an analytical solution for unsteady flow in test well analysis, which laid a theoretical basis for test well interpretation [3]. In the early 1970s, Ramey, Agarwal, and others developed various test well interpretation curves for homogeneous reservoirs and created the curves fitting analysis method, by which the flow stages of fluids can be distinguished. However, the fitting analysis was not widely used because it was quite cumbersome [1]. In 1970, Gringarten et al. developed pressure interpretation curves for various models based on Ramey curves, using a combination of skin and tubular reservoir coefficients, and constructed a complete ““modern test well interpretation method” by integrating curves fitting with a conventional test well interpretation methods. This method encompasses the selection of test well models, estimation of formation parameters, and verification of interpretation results and is followed by the development of ““test well interpretation software.” Modern well test interpretation methods are beginning to be commonly used around the world [4]. In 1974, Gringarten et al. found that a large number of possible causes of transient pressure dynamic solutions for well and reservoir analysis could be constructed from transient Green’s functions [5]. In 1979, Gringarten et al. proposed dimensionless theoretical curve clusters of homogeneous reservoir considering reservoir and epidermal effects, defined dimensionless of pressure and time, and proposed to label different curves using dimensionless time combination parameters in typical curves [6]. In 1983, Bourdet et al. developed pressure derivative type curves and created a fitting interpretation method. This method developed well trial interpretation further [7]. The domestic well test analysis can date back to the 1940s; Karamay Oilfield popularized and applied the unstable well test interpretation method to research reservoir pressure systems. Well test expert Xianzhang Tong studied the conventional well test interpretation method and put forward many well test interpretation measures combining with the actual Chinese oilfield, which laid the foundation for well test work in China. According to the pressure change characteristics in the actual oilfield development, to master the stratum seepage dynamic information, the article obtained the stratum parameter information based on the fitting method of the pressure and pressure derivative map, which provides a reference for the better evaluation of the homogeneous oil and gas reservoir.In recent years, the application of test well interpretation software has improved the level of test well interpretation by a significant step. Full-featured well test interpretation software can be used not only to interpret data for various tests (e.g., pressure drop and pressure recovery tests, pressure drop and pressure recovery tests for variable production scenarios, capacity tests, interwell interference tests, and so on) in oil, gas, and water wells but also to conduct various well test studies and generate new well test interpretation plots [4]. The use of computer fitting eliminates the tedious steps of manual fitting while improving accuracy metrics. For the prediction of production, permeability, and porosity, various artificial intelligence algorithms based on big data from logging curves are proposed. Paper [8] summarised a variety of soft computing models in permeability predictions. According to [9], the initial well-known models for permeability prediction are ANN, fuzzy logic (FL), and neuro-fuzzy. Each of these models has its strength and weaknesses. Wang Fuhao et al. used gated recurrent unit (GRU) networks to predict reservoir production, and this model obtained a 5.68 RMSE (root mean square error) and a 4.16 MAE (mean absolute error) for the prediction of reservoir production [10]. Although artificial intelligence algorithms are the future of well testing analysis, predictions can still only be made for a limited number of formation parameters. At the same time, data-based prediction algorithms require a large amount of historical data, making it difficult to achieve satisfactory results for the analysis of new wells. Therefore, well test analysis methods based on the classical theory of permeability mechanics still have an important role to play.The rest of the paper organized as follows. We first give an overview of the theory of the radial flow model for homogeneous reservoirs, followed by simulation analysis with curve plotted after solving the formation equation. Then, we present the pressure inversion software designed by the above methods and give the inversion results of the well test data. The conclusions are drawn at last.


## 2. Theory of Radial Flow Model for Homogeneous Reservoirs

### 2.1. Theory of Homogeneous Formation Pressure Equation

On the basis of disregarding composite formations, reservoir models are usually classified into the following two categories:(1)
Homogeneous Oil Reservoir. The fluid involved in the reservoir is single, and the fluid flow coefficient () and the energy storage coefficient (ΦCth) are the same everywhere(2)
Heterogeneous Oil Reservoir. The flow coefficient (Kh/μ ) or energy storage coefficient (ΦCth) in the reservoir is quite differentThis paper mainly studies the homogeneous oil reservoir model and the pressure inversion methods.A typical model in homogeneous reservoirs is described as follows:(1)
Basic model is homogeneous reservoir(2)
Internal boundary conditions are straight wells with wellbore storage effect and skin effect(3)
Outer boundary conditions are as follows: infinity stratum-reservoirs are homogeneous and equal in thickness(4)
Fluids and rocks in the formation are slightly compressible(5)
Fluids in porous medium satisfy Darcy flow(6)
The fluid in the well is pumped at a constant flow rateq, and the formation pressure is the original formation pressure pi before the well is openedAfter considering the internal and external boundary conditions, the fixed solution problem in an infinite stratum can be represented by the following system of partial differential equations:The seepage control equation is as follows:(1)∂2pD∂rD2+2rD∂pD∂rD=∂pD∂tD.Initial conditions are as follows:(2)pDr,0=0.Internal boundary conditions are as follows:(3)CD∂psD∂tD−rD∂pD∂rD∣rD=1=1,psD=pD−S∂pD∂rD∣rD=1.Outer boundary conditions are as follows:(4)limrD⟶∞pDrD,tD=0.Form simultaneous equations from (1) to (4), and use the Laplace transform on them(5)d2p¯DdrD2+1rDdp¯DdrD=sp¯Dp¯wDs=p¯DrD,s−Sdp¯DdrDrD=1CDsp¯wD−dp¯DdrDrD=1=1slimrD⟶∞p¯DrD,s=0.s is a Laplace variable. The general solution of (5) is as follows:(6)p¯DrD,s=AI0rDs+BK0rDs.I0z represents the Bessel function with first class correction of order 0. K0z represents second-class of modified Bessel function. I0z⟶∞⟶∞. Using the properties of the function, the coefficient A=0 is determined by the external boundary condition, and the coefficient B is determined by the internal boundary condition as follows:(7)B=1s1CDsK0s+1+CDSssK1s.Therefore, the image function of the pressure distribution is solved in the Rasch space as follows:(8)p¯DrD,s=1sK0rDsCDsK0s+1+CDSssK1s.For well test analysis, the wellbore pressure must be obtained. After considering the epidermal factor, the wellbore pressure is as follows:(9)p¯wD=p¯DrD,s−Sdp¯DdrDrD=1.Substituting equation (8) into (9) while taking rD=1, the expression of the image function of the bottom pressure of the well can be found as follows:(10)p¯wDs=1sK0s+SsK1ssK1s+CDsK0s+sSK1s.The above equation establishes the expression of the wellbore pressure image function for the ““radial flow” model in homogeneous formations. In plotting pressure-time double logarithmic curves, we usually use a combination of parametersCDe2S rather than giving separate magnitudes for CD and S. Paper [11] gave the analytical solution of the expression as follows:(11)pwDtD=4π∫0∞1−e−α2tDdαα3α2CDS−1J1α+αCDJ0α2+α2CDS−1N1α+αCDN0α2.The analytical solution of the partial differential equations cannot be found in the case of complex boundary conditions. Therefore, the Laplace numerical inversion algorithm is chosen to find the numerical solution of the image function. The following equation defines the Stehfest numerical inversion algorithm [12]:(12)ft=ln2t∑i=1NVif¯s,where N is even, s is iln2/t, and the expression of Vi is as follows:(13)Vi=−1N/2+i∑k=i+1/2mini,N/2kN/22k!N/2−k!k!k−1!i−k!2k−i!.Note that the wrong equation of the Stehfest numerical inversion algorithm given in paper [12] and paper [13] gives the errata statement. According to (13), given N and i values, find the value and then reverse and perform the original function by the image function according to (12). In the equation, it is stipulated that N is even, and the value of N has a great impact on the calculation accuracy. For different types of functions, the N value should be adjusted in practice. By comparison with the standard type curves, N=8 is selected for numerical inversion, in which case the error is minimal.


### 2.2. Pressure Inversion Theory

The key step of pressure inversion is to fit the well test pressure data to the Gringarten and Bourdet type curves. When we find the best G-B curve consistent with the pressure data trend, select any point in the double logarithmic curve graph and record the horizontal and vertical coordinates of the two curves. As shown in Figure2, mark the coordinates of the fitted point as pD/ΔpM and tD/CD/tM.Figure 2
The fitting step in pressure inversion. This is a screenshot of the well test pressure inversion software we have made. The screenshot shows the fine fitting process of the raw data and the G-B model we built, and software records the coordinate values of the following two points: (tDCD, pD) and (Δt, Δp).The inversion calculation equation is as follows:(14)K=141.2qμBhpDΔpMKhμ=141.2qBpDΔpMKμ=141.2qBhpDΔpMKh=141.2qμBpDΔpM.In equation (14), K,Kh/μ,K/μ,andKh are called permeability, flow coefficient, mobility, and formation coefficient, respectively.Wellbore storage coefficient is calculated from time fitted value as follows:(15)C=0.000295Khμ1tD/CD/tM.Wellbore storage coefficient (C) is used to calculate dimensionless wellbore storage coefficient (CD) as follows:(16)CD=0.8936CϕCthrw2.Finally, the skin coefficientS is calculated by the curve fitting value CDe2SM:(17)S=12lnCDe2SMCD.


## 2.1. Theory of Homogeneous Formation Pressure Equation

On the basis of disregarding composite formations, reservoir models are usually classified into the following two categories:(1)
Homogeneous Oil Reservoir. The fluid involved in the reservoir is single, and the fluid flow coefficient () and the energy storage coefficient (ΦCth) are the same everywhere(2)
Heterogeneous Oil Reservoir. The flow coefficient (Kh/μ ) or energy storage coefficient (ΦCth) in the reservoir is quite differentThis paper mainly studies the homogeneous oil reservoir model and the pressure inversion methods.A typical model in homogeneous reservoirs is described as follows:(1)
Basic model is homogeneous reservoir(2)
Internal boundary conditions are straight wells with wellbore storage effect and skin effect(3)
Outer boundary conditions are as follows: infinity stratum-reservoirs are homogeneous and equal in thickness(4)
Fluids and rocks in the formation are slightly compressible(5)
Fluids in porous medium satisfy Darcy flow(6)
The fluid in the well is pumped at a constant flow rateq, and the formation pressure is the original formation pressure pi before the well is openedAfter considering the internal and external boundary conditions, the fixed solution problem in an infinite stratum can be represented by the following system of partial differential equations:The seepage control equation is as follows:(1)∂2pD∂rD2+2rD∂pD∂rD=∂pD∂tD.Initial conditions are as follows:(2)pDr,0=0.Internal boundary conditions are as follows:(3)CD∂psD∂tD−rD∂pD∂rD∣rD=1=1,psD=pD−S∂pD∂rD∣rD=1.Outer boundary conditions are as follows:(4)limrD⟶∞pDrD,tD=0.Form simultaneous equations from (1) to (4), and use the Laplace transform on them(5)d2p¯DdrD2+1rDdp¯DdrD=sp¯Dp¯wDs=p¯DrD,s−Sdp¯DdrDrD=1CDsp¯wD−dp¯DdrDrD=1=1slimrD⟶∞p¯DrD,s=0.s is a Laplace variable. The general solution of (5) is as follows:(6)p¯DrD,s=AI0rDs+BK0rDs.I0z represents the Bessel function with first class correction of order 0. K0z represents second-class of modified Bessel function. I0z⟶∞⟶∞. Using the properties of the function, the coefficient A=0 is determined by the external boundary condition, and the coefficient B is determined by the internal boundary condition as follows:(7)B=1s1CDsK0s+1+CDSssK1s.Therefore, the image function of the pressure distribution is solved in the Rasch space as follows:(8)p¯DrD,s=1sK0rDsCDsK0s+1+CDSssK1s.For well test analysis, the wellbore pressure must be obtained. After considering the epidermal factor, the wellbore pressure is as follows:(9)p¯wD=p¯DrD,s−Sdp¯DdrDrD=1.Substituting equation (8) into (9) while taking rD=1, the expression of the image function of the bottom pressure of the well can be found as follows:(10)p¯wDs=1sK0s+SsK1ssK1s+CDsK0s+sSK1s.The above equation establishes the expression of the wellbore pressure image function for the ““radial flow” model in homogeneous formations. In plotting pressure-time double logarithmic curves, we usually use a combination of parametersCDe2S rather than giving separate magnitudes for CD and S. Paper [11] gave the analytical solution of the expression as follows:(11)pwDtD=4π∫0∞1−e−α2tDdαα3α2CDS−1J1α+αCDJ0α2+α2CDS−1N1α+αCDN0α2.The analytical solution of the partial differential equations cannot be found in the case of complex boundary conditions. Therefore, the Laplace numerical inversion algorithm is chosen to find the numerical solution of the image function. The following equation defines the Stehfest numerical inversion algorithm [12]:(12)ft=ln2t∑i=1NVif¯s,where N is even, s is iln2/t, and the expression of Vi is as follows:(13)Vi=−1N/2+i∑k=i+1/2mini,N/2kN/22k!N/2−k!k!k−1!i−k!2k−i!.Note that the wrong equation of the Stehfest numerical inversion algorithm given in paper [12] and paper [13] gives the errata statement. According to (13), given N and i values, find the value and then reverse and perform the original function by the image function according to (12). In the equation, it is stipulated that N is even, and the value of N has a great impact on the calculation accuracy. For different types of functions, the N value should be adjusted in practice. By comparison with the standard type curves, N=8 is selected for numerical inversion, in which case the error is minimal.


## 2.2. Pressure Inversion Theory

The key step of pressure inversion is to fit the well test pressure data to the Gringarten and Bourdet type curves. When we find the best G-B curve consistent with the pressure data trend, select any point in the double logarithmic curve graph and record the horizontal and vertical coordinates of the two curves. As shown in Figure2, mark the coordinates of the fitted point as pD/ΔpM and tD/CD/tM.Figure 2
The fitting step in pressure inversion. This is a screenshot of the well test pressure inversion software we have made. The screenshot shows the fine fitting process of the raw data and the G-B model we built, and software records the coordinate values of the following two points: (tDCD, pD) and (Δt, Δp).The inversion calculation equation is as follows:(14)K=141.2qμBhpDΔpMKhμ=141.2qBpDΔpMKμ=141.2qBhpDΔpMKh=141.2qμBpDΔpM.In equation (14), K,Kh/μ,K/μ,andKh are called permeability, flow coefficient, mobility, and formation coefficient, respectively.Wellbore storage coefficient is calculated from time fitted value as follows:(15)C=0.000295Khμ1tD/CD/tM.Wellbore storage coefficient (C) is used to calculate dimensionless wellbore storage coefficient (CD) as follows:(16)CD=0.8936CϕCthrw2.Finally, the skin coefficientS is calculated by the curve fitting value CDe2SM:(17)S=12lnCDe2SMCD.


## 3. Well Test Curve Analysis

### 3.1. Gringarten Type Curves

Gringarten type curves were proposed by Gringarten in 1979 [6]. It consists of a set of dimensionless time vs. pressure double logarithmic curves. The Stehfest numerical inversion is conducted to equation (10) to draw the Gringarten type curves. Equation (10) has two parameters S and CD, but we use one combined parameter CDe2S to replace these two parameters so that the pressure equation only has one variable parameter. The horizontal and vertical coordinates are tD/CD and pD.Different wellbore conditions are determined according to the range of values of the combination parameters.S>0 corresponds to contaminated or uncontaminated wells, and S<0 corresponds to acidification effective wells or fracking effective wells. In the numerical inversion of the theoretical model, the positive and negative skin coefficient needs to be solved using different equations:Contaminated well:103<=CDe2SUncontaminated well:5<=CDe2S<103Acidification effective well:0.5<=CDe2S<5Fracking effective well:0.5<CDe2SIn equation (10), Cd and S are two independent variables, not in the form of combined parameters CDe2S. As a result, one combined parameter has multiple sets and pairs of parameters of S and CD, which does not ensure inversion uniqueness. According to the conclusion of paper [1], when S>0, the deviations are taken to be lower than 0.1% whether we take S=0 or S≠0. Therefore, for all combinations of the parameters CDe2S>5, take S=0 and substitute into expression (10).WhenS>0, the effective wellbore concept defined by (14) regards the dredging action around the well wall as the larger wellbore radius.(18)rwe=rwe−S.rw represents the actual radius of the wellbore, which is equivalent to expanding the effective wellbore radius under the effect of the negative skin factor. Specifically, it is let S=0 in equation (10); the dimensionless time tD and dimensionless pipe storage constants CD are replaced with tDe2S and CDe2S. Equation (10) is reduced to the following:(19)p¯wDs=1sK0ssK1s+sCDe2SK0s.In this form, the value of the combination parameterCDe2S can be brought directly into equation (19). Due to the change in the reference length, at the numerical inversion, it is taken as follows:(20)s=iln2tD/CDCDe2S.The Gringarten pressure type curves established by the above method are shown in Figure3. The unconstrained pressure increased with the increase in combined parameters; the value we take is [1e-3, 5e-3, 1e-2, 5e-2, 1e-1, 5e-1, 1, 5, 10, 100, 1000, 1e4, 1e6, 1e8, 1e10, 1e15, 1e20, 1e30]. In the early stage, the angle between the curve and the horizontal axis is 45°when the logarithmic periods of the horizontal and vertical coordinates are equal. In the middle stage, the curve appears as a straight line nearly parallel to the horizontal axis, with different degrees of ““upward curvature” depending on the values of the combination parameters.Figure 3
The model of homogeneous formation radial flow Gringarten curve clusters.


### 3.2. Bourdet Type Curves

Bourdet calculated the dimensionless pressure derivative based on the Gringarten pressure curves [7]. A cluster of curves drawn with the same combined parameter value is called the Bourdet type curves, usually combined with the Gringarten curves to form the composite version for pressure inversion. The advantage of the pressure derivative curves is that it can distinguish areas that are too densely populated in the Gringarten curves and can be used to identify flow stages because of the clear distinction between the pressure derivative curves for different flow stages. The Gringarten–Bourdet curves have become the standard method for interpreting pressure inversions in well testing.In the interpretation of the pressure derivatives, the rate of pressure change relative to time during the test period is considered. To enhance the identification of the radial flow period, the pressure derivative is derived against the natural logarithm of time. Using the natural logarithm, the pressure derivative is expressed as the time derivative multiplied by the time takenΔt.(21)Δp′=dpdlnΔt=Δtdpdt.According to this method, a radial flow model of the homogeneous reservoir with infinite formation straight well with wellbore storage effect and skin effect is developed. As shown in Figure4, the solid and dashed lines correspond to the dimensionless pressure plot curves and pressure derivative plot curves for different values of the combined parameters, respectively. The combination parameters take values in the range of 0.1 to 1e60.Figure 4
Gringarten–Bourdet model of homogeneous oil reservoir.It can be seen that the pressure derivative curves are more significant than the pressure curve features. Figure5 shows a set of Gringarten–Bourdet plots. In the early stage, the Gringarten and Bourdet curves overlap. At this period, the formation pressure is influenced only by the wellbore storage effect, also known as the pure wellbore storage stage. And in the transition stage, the Bourdet curves are peaked with different combinations of parameter values, and the apex of the ““peak” and the magnitude of the change are significantly different. During this stage, the formation pressure is influenced by a combination of the wellbore storage effect and the skin effect. This stage can be used to determine which sample curve the measured wellbore pressure belongs to. In the last stage of the curve, the Bourdet curves show a ““0.5 line” parallel to the horizontal axis. This stage is known as the radial flow stage. Fitting points using derivatives are determined by early stage (straight lines with a slope of 1) and radial flow stage (0.5 lines) [14].Figure 5
Drawing edition of radial flow pressure and pressure derivations of homogeneous oil reservoir.


## 3.1. Gringarten Type Curves

Gringarten type curves were proposed by Gringarten in 1979 [6]. It consists of a set of dimensionless time vs. pressure double logarithmic curves. The Stehfest numerical inversion is conducted to equation (10) to draw the Gringarten type curves. Equation (10) has two parameters S and CD, but we use one combined parameter CDe2S to replace these two parameters so that the pressure equation only has one variable parameter. The horizontal and vertical coordinates are tD/CD and pD.Different wellbore conditions are determined according to the range of values of the combination parameters.S>0 corresponds to contaminated or uncontaminated wells, and S<0 corresponds to acidification effective wells or fracking effective wells. In the numerical inversion of the theoretical model, the positive and negative skin coefficient needs to be solved using different equations:Contaminated well:103<=CDe2SUncontaminated well:5<=CDe2S<103Acidification effective well:0.5<=CDe2S<5Fracking effective well:0.5<CDe2SIn equation (10), Cd and S are two independent variables, not in the form of combined parameters CDe2S. As a result, one combined parameter has multiple sets and pairs of parameters of S and CD, which does not ensure inversion uniqueness. According to the conclusion of paper [1], when S>0, the deviations are taken to be lower than 0.1% whether we take S=0 or S≠0. Therefore, for all combinations of the parameters CDe2S>5, take S=0 and substitute into expression (10).WhenS>0, the effective wellbore concept defined by (14) regards the dredging action around the well wall as the larger wellbore radius.(18)rwe=rwe−S.rw represents the actual radius of the wellbore, which is equivalent to expanding the effective wellbore radius under the effect of the negative skin factor. Specifically, it is let S=0 in equation (10); the dimensionless time tD and dimensionless pipe storage constants CD are replaced with tDe2S and CDe2S. Equation (10) is reduced to the following:(19)p¯wDs=1sK0ssK1s+sCDe2SK0s.In this form, the value of the combination parameterCDe2S can be brought directly into equation (19). Due to the change in the reference length, at the numerical inversion, it is taken as follows:(20)s=iln2tD/CDCDe2S.The Gringarten pressure type curves established by the above method are shown in Figure3. The unconstrained pressure increased with the increase in combined parameters; the value we take is [1e-3, 5e-3, 1e-2, 5e-2, 1e-1, 5e-1, 1, 5, 10, 100, 1000, 1e4, 1e6, 1e8, 1e10, 1e15, 1e20, 1e30]. In the early stage, the angle between the curve and the horizontal axis is 45°when the logarithmic periods of the horizontal and vertical coordinates are equal. In the middle stage, the curve appears as a straight line nearly parallel to the horizontal axis, with different degrees of ““upward curvature” depending on the values of the combination parameters.Figure 3
The model of homogeneous formation radial flow Gringarten curve clusters.


## 3.2. Bourdet Type Curves

Bourdet calculated the dimensionless pressure derivative based on the Gringarten pressure curves [7]. A cluster of curves drawn with the same combined parameter value is called the Bourdet type curves, usually combined with the Gringarten curves to form the composite version for pressure inversion. The advantage of the pressure derivative curves is that it can distinguish areas that are too densely populated in the Gringarten curves and can be used to identify flow stages because of the clear distinction between the pressure derivative curves for different flow stages. The Gringarten–Bourdet curves have become the standard method for interpreting pressure inversions in well testing.In the interpretation of the pressure derivatives, the rate of pressure change relative to time during the test period is considered. To enhance the identification of the radial flow period, the pressure derivative is derived against the natural logarithm of time. Using the natural logarithm, the pressure derivative is expressed as the time derivative multiplied by the time takenΔt.(21)Δp′=dpdlnΔt=Δtdpdt.According to this method, a radial flow model of the homogeneous reservoir with infinite formation straight well with wellbore storage effect and skin effect is developed. As shown in Figure4, the solid and dashed lines correspond to the dimensionless pressure plot curves and pressure derivative plot curves for different values of the combined parameters, respectively. The combination parameters take values in the range of 0.1 to 1e60.Figure 4
Gringarten–Bourdet model of homogeneous oil reservoir.It can be seen that the pressure derivative curves are more significant than the pressure curve features. Figure5 shows a set of Gringarten–Bourdet plots. In the early stage, the Gringarten and Bourdet curves overlap. At this period, the formation pressure is influenced only by the wellbore storage effect, also known as the pure wellbore storage stage. And in the transition stage, the Bourdet curves are peaked with different combinations of parameter values, and the apex of the ““peak” and the magnitude of the change are significantly different. During this stage, the formation pressure is influenced by a combination of the wellbore storage effect and the skin effect. This stage can be used to determine which sample curve the measured wellbore pressure belongs to. In the last stage of the curve, the Bourdet curves show a ““0.5 line” parallel to the horizontal axis. This stage is known as the radial flow stage. Fitting points using derivatives are determined by early stage (straight lines with a slope of 1) and radial flow stage (0.5 lines) [14].Figure 5
Drawing edition of radial flow pressure and pressure derivations of homogeneous oil reservoir.


## 4. Development of the Pressure Inversion Interpretation Software

### 4.1. Pressure Inversion Software Algorithms and Execution Process

The pressure inversion interpretation software replaces the artificial fitting process and integrates the inversion algorithm. By entering the basic formation parameters and obtaining the measured data of the wellbore pressure drawdown (build-up), the interpretation version is selected and fitted and the formation parameters information is reversed. The operation interface is simple and easy to use. The user enters the initial parameters at the initial interface, distinguishes the pressure drawdown data or builds up through the production time before closing the well, and calls up different data processing algorithms.


#### 4.1.1. Data Preprocessing

The measured data will be treated differently depending on the pressure build-up and pressure drawdown because G-B type curves are based on the pressure drawdown method. So, in the event of a pressure drawdown, i.e., when measuring the pressure at the start-up production of the well, we can directly use the measured pressure double logarithmic curve to fit G-B type curves.In the pressure build-up test after well shutdown, the pressure build-up value at momentΔt of well shutdown is as follows:(22)Δpbuild_upΔt=141.2qμBKhpDtpD−pDtp+ΔtD+pDΔtD,where tp is the well shutdown time. In this case, the pressure derivative changes from equation (21) to the following equation:(23)pDbuild_up′=dpDuild_updln△tD/CD⋅tp+△tDtpD=dpDuild_uptDd△tD/CD⋅△tDCD⋅tp+△tDtpD.For pressure build-up analysis, the pressure derivative of the raw pressure data is calculated by multiplying equation (21) by the well shutdown time factor tp+△tD/tpD. The derivative curves for the pressure build-up test and the pressure drawdown test are shown in Figure 6.Figure 6
Derivative of pressure build-up test calculation and derivative of pressure drawdown test calculation. We use equation (25) to smooth the curves.Then, consider the calculation of the pressure derivative. Because the actual well testing pressure data are not equally spaced, the derivatives are discontinuous when the pressure derivatives are calculated by the differential method. We use the ““moving window method” to smooth the pressure derivative. To find the derivative of the pressurep at point ti, the following equation is used:(24)pdtti=pti−pt1ti−t1t2−ti+pt2−ptit2−titi−t1,where t1,p1andt2,p2 are the furthest data point that satisfies t−ti⩽d.Derivative for time is as follows:(25)p′ti=dpdlntti=pti−pt1/lnti−lnt1lnt2−lnti+pt2−pti/lnt2−lntilnti−lnt1lnt2−lnti.In Figure7, we compared the raw data pressure derivatives calculated using equation (25) with the directly calculated pressure derivatives. As can be seen from the shape of the curve, the curve profile is effectively retained by smoothing, which helps to fit.Figure 7
Calculation of derivatives from well log pressures. The red dots are the double logarithmic curve plotted after calculating the pressure differenceΔp and time difference Δt from the raw pressure data p−t; the black dots are its derivative Δp′.
(a)(b)


#### 4.1.2. Program Execution Workflow

Flow chart in Figure8 shows the process of software. To enhance the accuracy of the fit, we divided the fitting step into two steps. After raw data preprocessing, a crude fit is first performed; this stage mainly matches the Gringarten curves. Because of the close trend of adjacent Gringarten curves in the interval with high combination parameter CDe2S, it is difficult to select the curve with the highest degree of matching. The software records the parameters of the 3 Gringarten curves adjacent to the matching curve during the crude fitting stage, which is used to generate the curve during the fine fitting stage. Then, in the fine fitting stage, the software generates 3 G-B curves based on the parameters from the crude fitting stage, which point the Bourdet curve with the closest fitting trend. The processing of the fine fit is shown in Figure 2.Figure 8
Flow chart of homogeneous reservoir pressure inversion software.After the fitting process, the software records the fitted parameters (pD), (Δp), (tD/CD), and (Δt) and pressure inversion calculations using equations (14)–(17). At last, parameters such as stratum permeability (K), epidermis coefficient (S), pipe storage coefficient (C), and so on can be obtained.


### 4.2. Verification of Homogeneous Oil Reservoir Pressure Recovery Data

Two sets of homogeneous reservoir logging data were used to verify the accuracy of the radial flow model and the inversion software. The well logging data are difficult to obtain; some data we used are from paper [7]; paper data are English units; to avoid the error caused by unit conversion, the software inversion algorithm is also based on English units.Here, we have performed two sets of data calculations, which is the pressure build-up data measured after well shutdown. The initial parameters at the wellbore are shown in Tables1 and 2. In homogeneous strata, the parameters obtained by our inversion are the formation coefficients, skin coefficients, and pipe storage coefficients, which are shown in Tables 3 and 4. From equation (14), there are four sets of parameters related to permeability, but the parameters are calculated in the same way, so the results are given for comparison for only one of them.Table 1
Initial parameters of the well 1.
Parameter nameParameter valueq crude oil output (barrel per daily, bbl/d)174tp production time before well closing (hours, h)15.33B formation volume coefficient1.06Ct total compression factor4.2×10−6h formation thickness (foot, ft)107rw wellbore radius (foot, ft)0.29μ crude oil viscosity (Bali, cp)2.5Φ porosity0.25Table 2
Initial parameters of the well 2.
Parameter nameParameter valueq crude oil output (barrel per daily, bbl/d)1500tp production time before well closing (hours, h)18.04B formation volume coefficient1.3Ct total compression factor1×10−5h formation thickness (foot, ft)73rw wellbore radius (foot, ft)0.401μ crude oil viscosity (Bali, cp)0.5Φ porosity0.2Table 3
Well 1 inversion results.
Inversion resultsSoftware calculation resultsBaseline resultsRelative errorFormation coefficientKh (mD ft)1182.6716511.52%Skin coefficientS7.677.70.39%Wellbore storage coefficientC9.15×10−39.3×10−31.61%Table 4
Well 2 inversion results.
Inversion resultsSoftware calculation resultsBaseline resultsRelative errorFormation coefficientKh (mD ft)151516307.06%Skin coefficientS−3.573.51.96%Wellbore storage coefficientC7.58×10−27.3×10−23.84%The results show that the difference between the software inversion and the paper calculation is not large; however, the numerical inversion process still requires manual participation, and the curves may have deviations depending on the fitted coordinate points or different subjective human choices, thus leading to large errors in the inversion results, so the pressure inversion should be experimented several times until a set of inversion results with less data dispersion is found.


## 4.1. Pressure Inversion Software Algorithms and Execution Process

The pressure inversion interpretation software replaces the artificial fitting process and integrates the inversion algorithm. By entering the basic formation parameters and obtaining the measured data of the wellbore pressure drawdown (build-up), the interpretation version is selected and fitted and the formation parameters information is reversed. The operation interface is simple and easy to use. The user enters the initial parameters at the initial interface, distinguishes the pressure drawdown data or builds up through the production time before closing the well, and calls up different data processing algorithms.


### 4.1.1. Data Preprocessing

The measured data will be treated differently depending on the pressure build-up and pressure drawdown because G-B type curves are based on the pressure drawdown method. So, in the event of a pressure drawdown, i.e., when measuring the pressure at the start-up production of the well, we can directly use the measured pressure double logarithmic curve to fit G-B type curves.In the pressure build-up test after well shutdown, the pressure build-up value at momentΔt of well shutdown is as follows:(22)Δpbuild_upΔt=141.2qμBKhpDtpD−pDtp+ΔtD+pDΔtD,where tp is the well shutdown time. In this case, the pressure derivative changes from equation (21) to the following equation:(23)pDbuild_up′=dpDuild_updln△tD/CD⋅tp+△tDtpD=dpDuild_uptDd△tD/CD⋅△tDCD⋅tp+△tDtpD.For pressure build-up analysis, the pressure derivative of the raw pressure data is calculated by multiplying equation (21) by the well shutdown time factor tp+△tD/tpD. The derivative curves for the pressure build-up test and the pressure drawdown test are shown in Figure 6.Figure 6
Derivative of pressure build-up test calculation and derivative of pressure drawdown test calculation. We use equation (25) to smooth the curves.Then, consider the calculation of the pressure derivative. Because the actual well testing pressure data are not equally spaced, the derivatives are discontinuous when the pressure derivatives are calculated by the differential method. We use the ““moving window method” to smooth the pressure derivative. To find the derivative of the pressurep at point ti, the following equation is used:(24)pdtti=pti−pt1ti−t1t2−ti+pt2−ptit2−titi−t1,where t1,p1andt2,p2 are the furthest data point that satisfies t−ti⩽d.Derivative for time is as follows:(25)p′ti=dpdlntti=pti−pt1/lnti−lnt1lnt2−lnti+pt2−pti/lnt2−lntilnti−lnt1lnt2−lnti.In Figure7, we compared the raw data pressure derivatives calculated using equation (25) with the directly calculated pressure derivatives. As can be seen from the shape of the curve, the curve profile is effectively retained by smoothing, which helps to fit.Figure 7
Calculation of derivatives from well log pressures. The red dots are the double logarithmic curve plotted after calculating the pressure differenceΔp and time difference Δt from the raw pressure data p−t; the black dots are its derivative Δp′.
(a)(b)


### 4.1.2. Program Execution Workflow

Flow chart in Figure8 shows the process of software. To enhance the accuracy of the fit, we divided the fitting step into two steps. After raw data preprocessing, a crude fit is first performed; this stage mainly matches the Gringarten curves. Because of the close trend of adjacent Gringarten curves in the interval with high combination parameter CDe2S, it is difficult to select the curve with the highest degree of matching. The software records the parameters of the 3 Gringarten curves adjacent to the matching curve during the crude fitting stage, which is used to generate the curve during the fine fitting stage. Then, in the fine fitting stage, the software generates 3 G-B curves based on the parameters from the crude fitting stage, which point the Bourdet curve with the closest fitting trend. The processing of the fine fit is shown in Figure 2.Figure 8
Flow chart of homogeneous reservoir pressure inversion software.After the fitting process, the software records the fitted parameters (pD), (Δp), (tD/CD), and (Δt) and pressure inversion calculations using equations (14)–(17). At last, parameters such as stratum permeability (K), epidermis coefficient (S), pipe storage coefficient (C), and so on can be obtained.


## 4.1.1. Data Preprocessing

The measured data will be treated differently depending on the pressure build-up and pressure drawdown because G-B type curves are based on the pressure drawdown method. So, in the event of a pressure drawdown, i.e., when measuring the pressure at the start-up production of the well, we can directly use the measured pressure double logarithmic curve to fit G-B type curves.In the pressure build-up test after well shutdown, the pressure build-up value at momentΔt of well shutdown is as follows:(22)Δpbuild_upΔt=141.2qμBKhpDtpD−pDtp+ΔtD+pDΔtD,where tp is the well shutdown time. In this case, the pressure derivative changes from equation (21) to the following equation:(23)pDbuild_up′=dpDuild_updln△tD/CD⋅tp+△tDtpD=dpDuild_uptDd△tD/CD⋅△tDCD⋅tp+△tDtpD.For pressure build-up analysis, the pressure derivative of the raw pressure data is calculated by multiplying equation (21) by the well shutdown time factor tp+△tD/tpD. The derivative curves for the pressure build-up test and the pressure drawdown test are shown in Figure 6.Figure 6
Derivative of pressure build-up test calculation and derivative of pressure drawdown test calculation. We use equation (25) to smooth the curves.Then, consider the calculation of the pressure derivative. Because the actual well testing pressure data are not equally spaced, the derivatives are discontinuous when the pressure derivatives are calculated by the differential method. We use the ““moving window method” to smooth the pressure derivative. To find the derivative of the pressurep at point ti, the following equation is used:(24)pdtti=pti−pt1ti−t1t2−ti+pt2−ptit2−titi−t1,where t1,p1andt2,p2 are the furthest data point that satisfies t−ti⩽d.Derivative for time is as follows:(25)p′ti=dpdlntti=pti−pt1/lnti−lnt1lnt2−lnti+pt2−pti/lnt2−lntilnti−lnt1lnt2−lnti.In Figure7, we compared the raw data pressure derivatives calculated using equation (25) with the directly calculated pressure derivatives. As can be seen from the shape of the curve, the curve profile is effectively retained by smoothing, which helps to fit.Figure 7
Calculation of derivatives from well log pressures. The red dots are the double logarithmic curve plotted after calculating the pressure differenceΔp and time difference Δt from the raw pressure data p−t; the black dots are its derivative Δp′.
(a)(b)


## 4.1.2. Program Execution Workflow

Flow chart in Figure8 shows the process of software. To enhance the accuracy of the fit, we divided the fitting step into two steps. After raw data preprocessing, a crude fit is first performed; this stage mainly matches the Gringarten curves. Because of the close trend of adjacent Gringarten curves in the interval with high combination parameter CDe2S, it is difficult to select the curve with the highest degree of matching. The software records the parameters of the 3 Gringarten curves adjacent to the matching curve during the crude fitting stage, which is used to generate the curve during the fine fitting stage. Then, in the fine fitting stage, the software generates 3 G-B curves based on the parameters from the crude fitting stage, which point the Bourdet curve with the closest fitting trend. The processing of the fine fit is shown in Figure 2.Figure 8
Flow chart of homogeneous reservoir pressure inversion software.After the fitting process, the software records the fitted parameters (pD), (Δp), (tD/CD), and (Δt) and pressure inversion calculations using equations (14)–(17). At last, parameters such as stratum permeability (K), epidermis coefficient (S), pipe storage coefficient (C), and so on can be obtained.


## 4.2. Verification of Homogeneous Oil Reservoir Pressure Recovery Data

Two sets of homogeneous reservoir logging data were used to verify the accuracy of the radial flow model and the inversion software. The well logging data are difficult to obtain; some data we used are from paper [7]; paper data are English units; to avoid the error caused by unit conversion, the software inversion algorithm is also based on English units.Here, we have performed two sets of data calculations, which is the pressure build-up data measured after well shutdown. The initial parameters at the wellbore are shown in Tables1 and 2. In homogeneous strata, the parameters obtained by our inversion are the formation coefficients, skin coefficients, and pipe storage coefficients, which are shown in Tables 3 and 4. From equation (14), there are four sets of parameters related to permeability, but the parameters are calculated in the same way, so the results are given for comparison for only one of them.Table 1
Initial parameters of the well 1.
Parameter nameParameter valueq crude oil output (barrel per daily, bbl/d)174tp production time before well closing (hours, h)15.33B formation volume coefficient1.06Ct total compression factor4.2×10−6h formation thickness (foot, ft)107rw wellbore radius (foot, ft)0.29μ crude oil viscosity (Bali, cp)2.5Φ porosity0.25Table 2
Initial parameters of the well 2.
Parameter nameParameter valueq crude oil output (barrel per daily, bbl/d)1500tp production time before well closing (hours, h)18.04B formation volume coefficient1.3Ct total compression factor1×10−5h formation thickness (foot, ft)73rw wellbore radius (foot, ft)0.401μ crude oil viscosity (Bali, cp)0.5Φ porosity0.2Table 3
Well 1 inversion results.
Inversion resultsSoftware calculation resultsBaseline resultsRelative errorFormation coefficientKh (mD ft)1182.6716511.52%Skin coefficientS7.677.70.39%Wellbore storage coefficientC9.15×10−39.3×10−31.61%Table 4
Well 2 inversion results.
Inversion resultsSoftware calculation resultsBaseline resultsRelative errorFormation coefficientKh (mD ft)151516307.06%Skin coefficientS−3.573.51.96%Wellbore storage coefficientC7.58×10−27.3×10−23.84%The results show that the difference between the software inversion and the paper calculation is not large; however, the numerical inversion process still requires manual participation, and the curves may have deviations depending on the fitted coordinate points or different subjective human choices, thus leading to large errors in the inversion results, so the pressure inversion should be experimented several times until a set of inversion results with less data dispersion is found.


## 5. Conclusions

(1)
Homogeneous reservoir pressure derivative curves have distinct characteristics, and the characteristics of the transition section curves can be used to determine which curve the measured wellbore pressure belongs to. The radial flow stage pressure derivative plot curve shows a “0.5 line” parallel to the horizontal axis, which is used to classify the flow stage of the wellbore fluid.(2)
The pressure inversion interpretation software was developed through the study of the curves fitting method, which can invert the values of formation parameters such as permeability (k), skin coefficient (S), and storage coefficient (C). We have designed the pressure inversion fitting procedure so that it can be done quickly for pressure inversion of log data. The inversion results were calculated from the wellbore measured pressure build-up data with an accuracy of more than 90%.(3)
The numerical inversion process of software inversion still requires manual participation, and the inversion results will show different degrees of errors with different choices of fitting points.



---

*Source: 4494678-2021-12-15.xml*