Modelling of Formation Pore Pressure Inversion during Tight Reservoir Drilling

Identifying and controlling a kicking well hinge on quickly obtaining reliable and accurate formation pore pressure. In this study, we derive an analytical model for estimating formation pore pressure when a gas kick occurs during tight reservoir drilling. The model considers the variations in gas volume and pressures in the annulus affected by mutual coupling between the wellbore and formation, as well as bubble migration and expansion in the annulus. Additionally, a numerical computation method that reduces the effect of measurement noise using the Hooke-Jeeves algorithm is proposed. The method is capable of estimating pore pressure during the early stage of a kick in real time, is robust to the inherit noise of the measurements, and can be applied in scenarios when a well shut-in process cannot be performed. The simulation results demonstrate that both kick simulation and formation pore pressure inversion can be conducted via the proposed methodology. The errors of the pore pressure estimating results are less than 2.03% compared to the field data of seven wells. The method is tested and validated to be robust to noise and maintain good convergence performance, thereby providing drilling engineers with a simple and quick way to estimate pore pressure during a kick.


Introduction
A well kick occurs when the formation pressure is higher than the wellbore pressure, causing fluid to flow from the formation into the wellbore. This can be especially dangerous when the fluid is gas, as blowouts can cause damage to equipment and personnel, even death. To reduce this risk, quick and reliable formation pore pressure measurements are needed to select an appropriate well killing method and design optimal hydraulic parameters (e.g., displacement and density of weighted kill mud) to handle the gas influx rapidly.
Currently, formation pore pressure models are commonly based on seismic and logging data, and the uncertainty within these models is influenced primarily by widely used model coefficients and measurement accuracy [1][2][3][4]. Inaccurate prediction and underestimation of formation pore pressure typically lead to pressure underbalance in an open hole section of the well, causing fluid to enter the wellbore. Once a gas kick occurs, appropriate operations, which are designed based on specific pore pressures, should be per-formed as quickly as possible to stop the influx and circulate out the formation fluid; therefore, quick estimation of the formation pore pressure is needed to effectively respond to a change in the drilling environment.
As an inverse problem, the pore pressure can be evaluated from the evolution process of gas influx, similar to pressure transient analysis during well testing [5][6][7][8]. To combat this problem, Samuel [9] and Miska et al. [10] proposed a method to estimate formation pore pressure and permeability based on buildup data collected upon shutting in a kicking well. In underbalanced pressure drilling (UBD) and managed pressure drilling (MPD), the data of pressure buildup and fluid flow metering on the surface can facilitate real-time reservoir characterization [11]. Using annular volumetric outflow data and bottomhole transient pressure, Kneissl [12] proposed a model for reservoir characterization during a gas kick in UBD. Likewise, Vefring et al. [13] proposed a novel method during UBD, which utilizes a dynamic well fluid-flow model coupled with a transient reservoir model. It adopts the Levenberg-Marquardt method and an ensemble Kalman filter to estimate reservoir properties. For managed pressure drilling, Gravdal et al. [14] proposed a method that performs a polynomial curve fitting algorithm. The model relies on the characteristics of a pressure-buildup curve during a kick and was verified as being strongly antinoise. Employing the method of rate-integral productivity index (RIPI) analysis, Suryanarayana et al. [15] and Shayegi et al. [16] evaluated the reservoir pressure and productivity during UBD. Besides, a genetic algorithm in conjunction with a transient two-phase reservoir simulator was developed by Biswas et al. [17] to analyse the variation of formation permeability along the wellbore. Generally, these models commonly rely on a complicated wellbore multiphase flow model, which is difficult to be incorporated into the numerical parameter inversion algorithm. Furthermore, the buildup data of bottomhole transient pressures and surface shut-in pressure measurements (e.g., drill pipe pressure and casing pressure) are necessary in the traditional models. However, they can be restricted or unavailable to some extent in several conditions: (i) there exist no bottomhole pressure measurements. (ii) The surface pressures during the shut-in process of the kicking well are difficult to collect because an undesirable well shut-in operation can lead to fracturing of the formation and well blowout [14] when a kick size is larger than the predetermined kick tolerance or the pressure-bearing capacities of the well control device and formations are low. (iii) There e is a risk of gas entering the drill string during shut-in when the float valves cannot work. In this scenario, the pore pressure estimated via shut-in drill pipe pressure can be larger than the actual value. (iv) Under some conditions, such as low formation permeability or large gas influx size, the time for stabilizing the standpipe pressure and casing pressure can be more than 2-3 hours, which may lead to a delay in controlling the kicking well quickly and an increase in nonproductive time.
Here, we derive an analytical model for estimating pore pressure during a gas kick in tight reservoir. The model uses the flow rate measurements taken as gas enters, rises, and expands in the annulus, rather than relying on the pressure-buildup curve. The analytical model, combined with a fast search algorithm, can ensure rapid convergence and robustness to noise. Additionally, the model is independent of the pressure measurements and capable of estimating pore pressure as soon as a kick occurs. Applied to the drilled gas-kick scenarios, the simulated results indicate that the model is in good quantitative agreement compared to commercial software and field data.

Mathematical Model
2.1. Model Formulation. Gas flows into the wellbore, accompanied by anterior bubbles rising and expanding gradually in the annulus, during a gas influx. Thus, the bottomhole pressure drops and the underbalance pressure increases; meanwhile, the thickness of the exposed reservoir increases dynamically, which leads to enlargement of the mass transfer between the wellbore and reservoir. To express the complex process, some necessary assumptions are adopted.
(1) The classical gas distribution assumption is adopted, in which gas enters the wellbore as individual bubbles and the bubbles rise in the annulus [18][19][20]. Actually, it can be reasonable because the gas influx rate in tight reservoir is generally low. To simplify, we assume that the bubbles are spherical and the time interval of adjacent bubbles is Δt. This assumption can be suitable for the early stage of a well kick [21,22] (2) The dissolution of gas in the drilling fluid is neglected. Furthermore, the density and rheological property of the drilling fluid are assumed to be constant [23] (3) According to the PVT equation of the gas state, PV ∝ ZT. In the early stage of a kick, change in the product of the compressibility factor and fluid temperature near the bottom hole is negligible and the value of PV for a single bubble is assumed to be constant [18] (4) Consider the change in thickness of the exposed reservoir when a kick occurs during drilling Based on the assumptions above, the proposed model is suitable for the early stage of a kick, namely, a small kick whose size is generally defined as less than 1 m 3 (pit gain). According to the well control handbook, these wells are commonly suggested to be shut-in when a detected pit gain is larger than 1 m 3 .
2.1.1. Model of Gas Volume in the Wellbore. The process of gas bubble formation, rising, and expansion in the annulus during a kick is shown in Figure 1. Based on that, we analyse the relationship between gas volume and annulus pressures at different times. Initially, the bottomhole pressure is equal to the sum of the hydrostatic pressure and friction pressure without gas entering into the wellbore [19], i.e., where p wf ,0 is the bottomhole pressure at the initial time, Pa; ρ w is the density of the drilling fluid, kg/m 3 ; g is the 2 Geofluids gravitational acceleration, m/s 2 ; L is the initial well depth, m; f is the frictional resistance coefficient [24]; v is the fluid velocity in the annulus, m/s; and D is the equivalent diameter of the annulus, m. When gas flows into the wellbore, every bubble forms within Δt time. At time t, there exist N bubbles in the annulus, N = t/Δt. The pit gain, bottomhole pressure, and position of the top bubble are as follows [19,21]: where V pit,t is the pit gain at time t which can be measured on the surface, m 3 ; V g,t is the total gas volume in the annulus measured at the surface, m 3 ; A is the crosssectional area of the annulus, m 3 ; R P is the rate of penetration (ROP), m/s; t is time, s; p wf ,t is the bottomhole pressure, Pa; ρ g is the gas density, kg/m 3 ; l g,t is the position of the top bubble at time t, m; λ 0 is the gas distribution parameter [25]; v ∞ is the rising velocity of the bubble, m/s; σ is the surface tension, N/m; v g is the gas velocity in the annulus, m/s; and N is the number of bubbles at time t.
From time t to time t + Δt, the former bubbles (1 − N) rise and expand in the annulus; meanwhile, gas enters the wellbore from the formation and forms a new bubble (N + 1) at the bottom of the wellbore. The pressure on the top bubble in the annulus is then The pressure on the bubble when it formed initially is expressed as The product of the pressure and volume of a single bubble is assumed to be constant in the early stage of a kick, where p 1,0 and p 1,t+Δt are the pressures on the first bubble at the initial time and time t + Δt, respectively, Pa; V 1,0 and V 1,t+Δt are the volumes of the first bubble at the initial time and time t + Δt, respectively, m 3 .
From time t to time t + Δt, the volume change of the first bubble is By analogy, the volume change of the i th bubble is where V i,t and V i,t+Δt are the volumes of the i th bubble at time t and time t + Δt, respectively, m 3 . For the whole annulus, the total change of gas volume is equal to the sum of the expansion volume of the former (1 − N) th bubbles and volume of the (N + 1) th bubble, i.e., The rate of gas influx at the bottom hole is [26] where C is a variable related to time, i.e., where q is the rate of the gas influx, m 3 /s and △p is the underbalance pressure at time t, Δp = p 0 + p, Pa. According to equation (11), the initial volume of the i th bubble is where Equation (13) is introduced into equation (10) Rewriting equation (15) into an integral form, Equation (16) is the derivative of gas volume, which is a function of underbalance pressure, kicking time, ROP, and reservoir properties.
2.1.2. Model of Pore Pressure Inversion. The formation pressure can be expressed as [10] where p pore and p 0 are the formation pressure and initial underbalance pressure, respectively, Pa. p is the variation in underbalance pressure, Pa. At the initial time, there exists only single phase liquid flow in the annulus, p = 0, and the bottomhole pressure (p wf ,0 ) can be calculated easily. Therefore, the estimation of the formation pore pressure (p pore ) mainly lies in the inversion of the initial underbalance pressure (p 0 ). Through derivation of equation (3) and (17), The variation in total gas volume consists of two parts. The first part is volume increment caused by gas rise and expansion in the annulus; the second is the latest volume of gas entering from the formation into the wellbore. Combining equation (16) with equation (18), the equation for initial underbalance pressure inversion during a well kick can be obtained If the gas kick occurs when drilling is stopped, i.e., R P = 0, the equation can be written as where h loss is the thickness of the exposed reservoir, m. If the formation pore pressure and permeability are known, kick prediction can be conducted using the proposed model. Furthermore, either formation pore pressure or permeability can be estimated based on the flow rate data, although they cannot be calculated simultaneously in most cases [13]. For pore pressure inversion, we assume formation permeability, which is usually obtained using data from the core, well log, initial production, and drill cuttings [27][28][29] to be known in the model.
Utilizing the model in an analytical form, one fundamental feature of our method is that it is simple and easy to perform; furthermore, it can be used to estimate the formation pore pressure as soon as a gas kick is detected, not relying on in situ pressure measurements.

Numerical Simulation.
Formation pore pressure estimation is considered to be an inverse problem of kick simula-tion. For multiphase flow models, the solution of the inverse problem is complex because a set of partial differential equations and auxiliary equations [30][31][32] needs to be coupled and solved under specific initial and boundary conditions [33]. In this work, using an analytical model that expresses the relationship of the model parameters and output values directly, it is simple to estimate formation pore pressure.
Mathematically, the inversion of the formation pore pressure seeks the optimum value (p) to minimize the objective function, JðpÞ. Regarding the numerical calculation, we adopt the Hooke-Jeeves search algorithm, which was proposed by Hooke and Jeeves [34] and proved as convergent by Torczon [35]. Here, the detailed numerical simulation procedure is designed and presented in Appendix A.
If the value of JðpÞ is minimum when p = x i , the system deviation of the proposed model is [33] J sys = e Here, J sys is the system deviation generated by the model are the calculated data series and "ideal" measurements of the pit gain, respectively.
In the field, the measurements usually have random noises because of environmental factors (e.g., wind, current, and waves), negligence of operators, instrument error, etc. The value of JðpÞ, influenced by both system deviation and noise, is Here, J rel is the real deviation generated by the model and measurement noise; E ! noise is the data series of measurement errors at different times.
In equation (22), the second item on the right-hand side represents the influence of the noise on the solution. The value of J rel is also nearly minimal when p = x i , and the effect of noise on the search algorithm is negligible. The second item is small in most cases because the noise is commonly uncertain and random. noise vary considerably. Therefore, a search algorithm is used rather than a direct solution to create a method that is robust to noise. Moreover, the analytical form of the pore pressure inversion model and the fast search algorithm can ensure a quick solution process.

Case Study
In this section, analyses of kick simulation and formation pore pressure inversion are conducted using the proposed model. For kick simulation, we compare the results of this model with the "kick" module of the Drillbench software; moreover, the influence of the formation permeability and ROP on the kick development process is analysed. Subsequently, we estimate the formation pore pressure based on pit gain data and study the influence of noise and underbalance pressure on the model performance. Major parameters are shown in Table 1.

Kick Simulation.
The gas kick incident may occur under drilling and cessation conditions. In simulation, the main difference lies in the exposed reservoir thickness, which is assumed as variable and constant for drilling and cessation conditions, respectively. Applying the proposed analytical    In Figure 2, the curves of the two methods show the same trend under the same kicking condition. For a kick while drilling, the time until the pit gain reaches 1 m 3 is 336.9 s, as calculated via the analytical method, and 342.2 s, as calculated via the numerical method (relative error: 1.55%); similarly, they are 163.3 s and 168.5 s, respectively, for a kick when drilling is stopped (relative error: 3.09%). The main reason for the slight difference is that gas and liquid compressibility is neglected in the analytical model. For kick prediction in the early stage, the calculation error is small, at least within acceptable bounds.
With different ROP and reservoir permeabilities, the sensitivity analysis of pit gains during kicks is as shown in Figure 3.
In Figure 3, the pit gains increase in a nonlinear manner as time passes and the slopes for the curves increase gradually. Furthermore, when the ROP increases, the pit gain concurrently increases significantly. This is because the thickness of the exposed reservoir increases, accompanied by a decrease in bottomhole pressure, as the time or ROP increases, which leads to a latent and abrupt kick.
Similarly, the pit gain rises when the formation permeability increases. Here, as the formation permeability increases from 100 md to 200 md, the pit gain rises from 2.65 m 3 to 7.10 m 3 (2.68 times larger). For drilled formations with high pore pressure and permeability, low ROP and accurate kick detection need to be maintained.

Formation Pressure Inversion.
In this section, the model is applied to a drilled gas-kick situation. The true pit gain data    6 Geofluids yielded while kicking follows the same scenario as that in Figure 2. However, noise commonly exists in a true data set, and the influx rates are varied due to different pore pressures. First, the effect of noise within the dataset is analysed to test the accuracy and noise control of the model. The background research on quantifying noise within a dataset is a broad topic that is handled elsewhere [36][37][38]. Here, we add random noise to the simulated dataset to imitate a more realistic data set based on the data of the red solid curve in Figure 2 with a time interval of 1 s (i.e., sampling at 1.0 Hz). Using the generated dataset, the pressure inversion is then conducted.
For example, if the noise level is 0.1 m 3 , we add random noises that range from -0.1 m 3 to 0.1 m 3 to each value of the base dataset and obtain 1000 sets of different data series after sampling 1000 times. One sample of the obtained data series is shown in Figure 4.
The distribution of the estimated initial underbalance pressure, with a noise level of 0.1 m 3 , is shown in Figure 5. The red bars refer to the probabilities of inversion results and the blue curve refers to the fitting form of the probability distribution.
The estimated values are close to a Gaussian distribution with a mean of 1.5029 MPa (real value: 1.5 MPa). The 95% confidence interval of its distribution is 1.445 MPa to 1.561 MPa. The error of the model is less than 3.67%, with a 95% confidence interval, when the added noise level is 0.1 m 3 .
Similar to the sampling method in Figure 5, noises at different levels are added to the base dataset to test the proposed method. The noise levels change from 0 m 3 to 0.1 m 3 , and we sample 1000 times at each noise level, creating 100 × 1000 sets of "real" data series to estimate the underbalance pressure. Figure 6 shows the distribution of inversion results of the initial underbalance pressure at different noise levels. The red curve refers to the means of the samples (1000 times at a noise level); the blue points refer to the values with a 95% confidence interval.
These results suggest that the means of the estimation results in the samples are almost 1.5 MPa, with slight fluctuations as the noise level increases. With a 95% confidence interval, the error of the estimated value increases as noise increases. Overall, the proposed method has good noise can-celling performance and its errors are acceptable, even within very noisy datasets.
Actually, the influx rates vary due to different underbalance pressures. Under the condition of different initial underbalance pressures (0.1 MPa-2.5 MPa), the base data series of the pit gain, changing from 0 m 3 to 1 m 3 , are simulated; furthermore, random noise (noise level: 0.1 m 3 ) is added to each data series, and 60 × 1000 sets of data series are obtained. Figure 7 shows the distribution of the inversion results with different initial underbalance pressures. The red curve refers to the means of the samples; the blue points refer to the 95% confidence interval.
In Figure 7, the means of the inversion results are nearly equal to the actual values. With a 95% confidence level, the absolute errors with different initial underbalance pressures are nearly the same, and the relative error of the inversion results decreases as the initial underbalance pressure increases. For example, the relative error of estimated results when p 0 equals 0.5 MPa is 11.7%, while it is 2.83% when p 0 equals 2.0 MPa.

Applications
In hydrocarbon drilling, kicks typically occur when highpressure formations are drilled using a drilling fluid with an improper density to compensate for the additional pressure. We collect data from seven drilled gas-kick wells that were mainly exploration and wildcat wells in the Sichuan and Xinjiang provinces, China. There were obvious signs when the gas influx began, e.g., increase in ROP, pit gain rising. However, due to the weakness of the well control device or formations, it was difficult to gain the surface pressures during the well shut-in processes, leading to difficulty in killing the wells quickly after the kicks occurred. For example, regarding well 1, after shut-in of 197 min, the standpipe pressure still increased and did not ultimately stabilize. Regarding well 3, the shut-in operation led to fracturing of the formation. Its casing pressure increased to 20 MPa and then dropped to 0 MPa; therefore, no reliable pore pressure was obtained and the well was not killed successfully.
Based on the recorded measurements (e.g., ROP and pit gains) during the kick progression, pore pressure inversions were performed, as shown in Table 2. The main parameters in the model were estimated as follows: (1) the thickness of the exposed reservoir was obtained based on bit footage during high ROP. For example, regarding well 1: at 20 : 32, ROP increased from 0.00027 m/s to 0.0008 m/s at 4481.88 m. At 20 : 45, drilling stopped; (2) the formation permeability was obtained through cores, well log, and information of adjacent wells. For example, the permeability of well 6 was obtained from information of an adjacent well. The other wells were cored in the kicking formations, with core recovery rates greater than 95%; (3) actual pore pressures were mainly obtained by well killing and well shut-in, along with productivity testing. For instance, because of uncertain pore pressure in well 1, the circulate-and-weight method was used to control the well; the increment was 0.03 g/cm 3 for each circle. At 9 : 00 17 August, the outflow mud density was 2.05 g/cm 3 , and the pump stopped with a slight kick. At 9 : 40 18 August, the outflow mud density was 2.12 g/cm 3 , and the pump stopped without a kick. Regarding well 4, 140 m 3 heavy mud (1.79 g/cm 3 ) was injected. After shut-in of 25 min, surface pressure stabilized: casing pressure was 8.5 MPa and standpipe pressure was 2.0 MPa. The detailed information of the kicking wells is given in Table 3 of Appendix B.
In Table 2, the errors of estimating pore pressures are less than 2.03%. The errors are caused mainly by the deviations of parameter estimations (e.g., true pore pressure, thickness of exposed reservoir, and permeability), as well as the model itself. Likewise, in other scenarios where the pore pressure is known, the formation permeability can be obtained using the proposed method.

Conclusions
In this study, we propose an analytical model for estimating the pore pressure on a kicking well. This model considers the variations in gas volume and pressures in the annulus influenced by mutual coupling between the wellbore and formation and bubble migration and expansion in the wellbore. Furthermore, a numerical calculation method using the Hooke-Jeeves algorithm is proposed. The methodology has the following merits: (i) It is fast and simple to use, especially for drilling into a tight reservoir when the gas influx rate is low, a wellbore bubble flow is formed and the required shut-in time is large (ii) It is applicable to scenarios where the pore pressure cannot be easily obtained via conventional methods because it capitalizes on the measurements of flow rate and does not rely on pressure measurements (iii) It estimates the formation pore pressure during the early stage of kicks, which can provide invaluable time and information for controlling a kicking well (iv) It is not sensitive to the inherit noise of measurements (v) It estimates the formation permeability if the pore pressure can be obtained by well shut-in Applied to inverse the formation pore pressure based on pit gain data, the analytical model is found to quickly converge and is both reliable and robust with different noise levels and underbalance pressures. Considering a noise level of 0.1 m 3 , we found that the distribution of the estimated underbalance pressures is approximately a Gaussian distribution and the relative errors are less than 3.67%, with a 95% confidence interval.
Additionally, the results of the proposed methodology produce very consistent trends and are in good quantitative agreement compared to the results of numerical software and field data.

B. The Detailed Information of the Kicking Wells
(1) The thicknesses of the exposed reservoir were obtained based on bit footage during high ROP (2) The permeabilities of formations were obtained based on adjacent wells, core or well log (i) For example, regarding well 1, the core recovery rate was 98.13%, with coring in the kicking formation. As for the other wells, except for well 6 (ii) The permeability of well 6 was obtained based on the information of its adjacent well (3) The real pore pressures were obtained via well killing, well shut-in and productivity testing (i) Well 1: the circulate-and-weight method was used and the increment was 0.03 g/cm 3 for each circle. At 9 : 00 17 August, outflow mud density was 2.05 g/cm 3 , and the pump stopped with a slight kick. At 9 : 40 18 August, outflow mud density was 2.12 g/cm 3 , and the pump stopped without a kick (ii) Well 2: after shut-in of 12.8 hours, pressure stabilized, with a casing pressure of 29 MPa and a standpipe pressure of 8 MPa Gravitational acceleration, m/s 2 J: Objective function for numerical simulation K: Formation permeability, m 2 L: Initial well depth, m l g,t : Position of the top bubble at time t, m M: Number of sampling points