Wellbore annulus water hammer pressure prediction based on transient multi-phase ﬂ ow characteristics

. Water hammer pressure has been known to cause formation fracture and well-control problems. Accurate prediction of water hammer pressure is crucially important to determine the selection of shut-in methods. In this study, the mathematic model of wellbore annulus transient water hammer has been established with the consideration of transient multi-phase ﬂ ow characteristics, and it has been solved by the Method Of Characteristic (MOC). Finally, this paper focused on the effects of gas cutting, shut-in time and friction on water hammer pressure, and gas kick time were also regarded to study on the in ﬂ uence of water hammer pressure. The results show that both the gas cutting and gas kick time have few in ﬂ uences on the maximum water hammer pressure, but intensi ﬁ ed the attenuation of water hammer pressure. Additionally, the peak value of water hammer pressure declines with the increase of the shut-in time, and the effect of friction loss on water hammer pressure became signi ﬁ cant with the increase of well depth. More importantly, both the additional water hammer pressure and Shut-In Casing Pressure (SICP) generated by the closure of BlowOut Preventer (BOP) are likely to cause formation at the shallow casing shoe damage.


Introduction
Gas kick takes place from time to time during the drilling process when shut-in is required. However, rapid shut-in can lead to a sudden change in flow velocity that will generate water hammer at the wellhead, posing serious threats to equipment Ouyang et al., 2009). Moreover, in gas and oil reservoir with narrow safety density window, where the reservoir pore pressure is close to fracturing pressure, the water hammer effects could result in a potential damage for formation, such as underground blowout and loss circulation (Jiang et al., 2014;Tang et al., 2014).
In the oil and gas industry, the study on water hammer mainly focuses on shut-in of water injection wells, raising the issues of sand production (Santarelli et al., 2000;Vaziri et al., 2008) and borehole stability (Han et al., 2002;Wang et al., 2008). Subsequently, measures such as alteration of valve installation position (Tang and Ouyang, 2010), adjusting the operating parameters (Choi and Huang, 2011) and controlling valve closure time (Zazovsky et al., 2014) have been taken to worsen the water hammer effects.
However, few studies have carried out to investigate drilling problems associated with water hammer. Jardine et al. (1993) firstly studied the hard or soft shut-in question. Nevertheless, this work neglected the impact of gas on water hammer wave speed. A similar work by Li and Zheng (1995), however, attracted greater attention. Li and Zheng (1995) were successful in calculating pressure for hard shutin with the consideration of gas void fraction effects. The influence of the unsteady friction, however, has not been taken into account in the model, so the water hammer wave attenuation cannot be exactly described. He et al. (2008) applied the method of Finite Element Method (FEM) to address the water hammer Partial Differential Equations (PDEs) and attempted to simulate water hammer by Automatic Dynamic Incremental Nonlinear Analysis (ADINA) finite element software. Nonetheless, this modeling method was formulated based on the assumption that the gas void fraction uniform distribution. Han et al. (2012 used the commercial software to represent the water hammer effects introduced during well startup and shut-in. Wang et al. (2016) investigated the water hammer effect caused by the sudden intrusion of formation fluid into the drilling process, while it does not refer to shut-in and only limited parameters such as gas influx rates was analyzed.
After reviewing the current research as mentioned above, it can be concluded that some unreasonable assumptions for calculating water hammer, such as flow parameters uniform distribution along the wellbore and water hammer speed without taking into account gas void fraction. In essence, the gas void fraction has a significant influence on water hammer speed (Esmaeilzadeh et al., 2009;Lee and Pejovic, 1996;Zhou et al., 2004). Additionally, when the formation gas entrance into the well, gas-liquid twophase flow will appear in the annulus, meaning that there was a single phase flow in the upper part of the annulus and gas-liquid two-phase flow in the lower part. And so, until the influx gas reached the ground, there was gas liquid two-phase flow in the entire annulus. Meanwhile, the influx gas will expand as it rises up the annulus, leading to the flow parameters (e.g., gas void fraction and mixture velocity) non-uniform distribution along the wellbore (Avelar et al., 2009;Sun et al., 2017;Yin et al., 2017).
In this paper, the mathematical model of water hammer fluctuations was established according to the transient multi-phase flow characteristics. The unsteady friction modeling was included. Then, the classical Method Of Characteristic (MOC) was also used to solve the governing equations for gas-liquid two-phase transient flow in the wellbore. Besides, both boundary and initial conditions of water hammer for well shut-in were determined. Finally, the effects of related parameters, such as gas cutting, shut-in time, well depth and gas kick time, on water hammer pressure were investigated.

Physical model
The physical model of the wellbore encountering gas kick during the drilling process is shown in Figure 1. Drilling fluid is pumped from the mud pits, down the drilling string, circulated down to the bit, through the drill bit nozzles, and back to the mud pits via the annulus. The annular BlowOut Preventer (BOP) is located at the top of the annulus and permits passage of drilling string. It is intended to seal the annulus space between the drilling string and the wellbore in a gas kick situation and to avoid an uncontrolled flow of gases or liquids from the well during drilling. When the gas kick is detected at the surface, the mud pump is turned off and BOP is closed. As the BOP is closed, the drilling fluid circulation path has to change to the choke line located below the annular BOP as depicted in Figure 1. It can allow fluids to flow across well control choke and bring it back into the mud pits during the well control operations.
During the normal drilling process, the bottom hole pressure is equal to or slightly higher than the formation pressure. In this case, the bottom hole pressure is sufficient to prevent any formation fluids invading into the wellbore. Thus, there will be only single fluid phase existed in the wellbore. Abnormal formation pressure, however, can be encountered during the drilling of a well in which problems can be unexpectedly occurred such as a gas kick. As a result, the natural gas will flow from the formation into the wellbore, which would form a gas-liquid two-phase flow in the annulus as shown in Figure 1. Besides, the leading edge of gas phase move forward gradually upward as the free gas is travelled up the wellbore. Once a gas kick is detected, the well should be shut-in timely through the BOP to stop the influx of formation fluids into the wellbore. The mathematic models were set up based on the following assumptions: 1. the fluid flow model in the wellbore is one-dimensional (1D) transient gas-liquid two-phase flow; 2. casing and drilling string are assumed to be linear elastic and the effects of cementing and formation are not taken into account; 3. the annulus fluid temperature profile is assumed equal to the formation temperature, and no heat transfer is accounted for; 4. drilling fluid and gas are regarded as compressible and the formation pressure is kept constant; 5. the influence of cuttings on water hammer wave speed is not considered; 6. the time for turning off the mud pump is not taken into account and the well control choke is closed before closing the annular BOP.
3 Mathematical model

Transient multi-phase model
In essence, transient multi-phase flow parameters are governed by the conservations of mass, momentum, and energy. In order to simplify the problem, the drilling fluid temperature profile was supposed to be linear, namely, it was equal to the formation temperature, and no heat transfer was accounted for. Consequently, for the 1D, the following mass and momentum conservation equations are applied: Conservation of gas-phase mass (gas-producing zone): where q g is the gas cutting speed, kg/(m s); q g is the gas density, kg/m 3 ; H g is the gas void fraction; u g is the gas velocity, m/s; A is the cross-sectional area of annulus, m 2 ; t is the time, s; z is the axial position, m. Conservation of gas-phase mass (non-gas producing zone): Conservation of liquid-phase mass: where q l is the drilling fluid density, kg/m 3 ; u l is the drilling fluid velocity, m/s; H l is the liquid holdup. Conservation of total momentum: where F r is the friction pressure drop, Pa; P is the pressure, Pa; g is the acceleration of gravity, m/s 2 .

Gas influx model
During drilling into a gas reservoir, formation gas begins to invade into the wellbore when the bottom hole pressure is lower than the formation pressure. The gas influx rate from the reservoir can be calculated from the binomial theorem equation (Huang and Ayoub, 2008): where P e is the formation pressure, MPa; P wf is the bottom hole pressure, MPa; K is the reservoir permeability, mD; h is the net-pay thickness, m; T is the average temperature,°C; l is the gas viscosity, mPa s; Z is the average compressibility factor; r e is the supply radius of gas reservoir, m; r w is the open hole radius, m; S k is the skin factor; q sc is the gas influx rate at standard condition, m 3 /s; r g is the relative density of the gas; b is the turbulent coefficient, m À1 .

Water hammer model
For one-dimensional unsteady fluid flow in the wellbore, the equation of motion and the continuity are the governing equations. The equation of motion and the continuity can be written as: where q is the mixture density, kg/m 3 ; u is the mixture velocity, m/s; k is the unsteady friction coefficient; D i is the inner diameter of annulus, mm; D o is the outer diameter of annulus, mm; s is the spatial coordinates, m; oz os ¼ À sin h, h is the angle between axis direction and horizontal direction. For a vertical well, sinh = 0. The hydroelasticity is defined as: where K is the liquid-phase volumetric elasticity modulus, MPa.
The pipe-wall elasticity is written as: where d 1 is the casing thickness, mm; E p is the casing elasticity modulus, MPa. Unsteady friction coefficient is calculated according to Vardy and Brown (2003) where Re is the Reynolds number. The water hammer wave speed considering gas content can be calculated as follows (Zhou et al., 2004): where d 2 is the drilling string thickness, mm; E l is the liquid-phase elasticity modulus, MPa; E i is the gasphase elasticity modulus, MPa; E d is the drilling string elasticity modulus, MPa; a m is the water hammer wave speed, m/s.

Solving method
In this paper, the classical MOC for water hammer PDEs was adopted. The MOC is used to transform the momentum and continuity PDEs into four ordinary differential equations. In order to apply the MOC, the equations (6) and (7) should be re-written as equations (12) and (13): By using an undetermined coefficient x, the linear combination of equations (12) and (13) in the form of L = L 1 + xL 2 can be expressed as: Typically, both the pressure P and velocity u are functions of distance s and time t, thus the total derivative describing both pressure P and velocity u are represented as follows: Compared equation (14) with equation (15), the undefined coefficient x can be determined as: Then, the substitution of equation (16) into equation (14) can lead to two sets of ordinary differential equations which are characterized by positive (C + : u + a m ) and negative (C À : u À a m ) equations (as shown in Fig. 2). In essence, the fluid velocity u is far less than the water hammer wave speed a m , thus the fluid velocity u both in the positive and negative characteristics equations can be ignored: According to Figure 2, the characteristic equations (17) and (18) depicted in are integrated along the positive (C + ) and negative (C À ) characteristic lines, respectively. The following characteristic equations are obtained:  Figure 3 shows the temporal and spatial mesh for the water hammer model. The abscissa is the space, and the ordinate is the time. Discretization of the partial differential equations (21) and (22), using a finite difference method, results in the following equations:

Boundary conditions
(1) Bottom hole boundary It is generally considered that the pressure boundary condition at i = 1 is consistent with the bottom hole pressure P wf : Substituting equation (23) into equation (6), the velocity boundary condition at i = 1 can be obtained as the following equation: (2) Wellhead boundary The variation of velocity at the wellhead is relevant to the closing law of BOP. Assuming that the outlet flow regularity conforms to the orifice discharge regularity, the velocity boundary condition at i = N is: Substituting equation (25) into equation (7), the pressure boundary condition at i = N can be calculated as follows: where s is the relative valve opening of BOP, dimensionless: In fact, there is no empirical formula available that can be used for calculating relative valve opening of BOP, so the relative valve opening of BOP was assumed to be in the same as in the characteristics of the valve opening in the process of change. For uniform valve closure arrangement, it can be expressed as (Karney and Ruus, 1985): s where T s is the BOP closing time, s.

Initial conditions
The annulus fluid velocity and pressure before shut-in can be obtained through simulation of gas kick process. The velocity and pressure of various nodes in wellbore annulus at the initial moment are respectively: where u 0 (i) and P 0 (i) are respectively the velocity and pressure of various nodes of wellbore annulus during gas kick.

Results and discussion
In order to study the variation of water hammer pressure resulting from the shutting-in of a kicking well, the distribution of flow parameters along well depth before shutting-in the well is needed to be determined by adopting multi-phase flow theory at first. Taking a vertical well in Tarim Basin as an example, the well was shut in when a gas kick was encountered during drilling. È311.1 mm borehole was drilled to 6300 m, and È244.5 mm casing was set at the depth of 6299.53 m; the gas kick occurred when È215.9 mm borehole was drilled to 6,436 m. È127 mm drill pipe was used. The density, plastic viscosity and yield point of drilling fluid were respectively, 1180 kg/m 3 , 24 mPa s and 8 Pa. The flow rate was 30 L/s. Other basic parameters are shown in Table 1.
Assuming that the bottom hole pressure is 0.5 MPa less than formation pressure, then formation fluids such as natural gas can flow into the well. The variations of gas void fraction and mixture velocity with gas kick time and well depth were shown in Figures 4 and 5 before shutting-in the well, respectively. It can be seen that the leading edge of gas-liquid two-phase flow continuously pushes forward to the wellhead as the gas kick time increases, and the influx gas migrated to the wellhead when the gas kick time was 45 min, which indicated that there was only liquid phase flow in the upper part of the wellbore annulus and gasliquid two-phase flow in the lower middle part. In the uncontaminated region, gas void fraction is 0 and liquid hold-up is 1. During the initial gas kick stage, the gas void fraction and mixture velocity have no significant changes; when increase sharply, that meant it is very close to the wellhead. Figure 6 illustrated that the variation of the well bottom hole pressure with time before shutting-in the well. It is observed that the bottom hole pressure is gradually decreases linearly in the early stage, but then decreases rapidly with the migration and expansion of gas, which supported by the gas void fraction variations in Figure 4. This happens due to reducing hydrostatic pressure in the annulus, noticing that gas density is considerably lower than the drilling fluid. Remarkably, during the drilling process, the annulus at the surface is open to the surface and its pressure is always equal to the atmospheric pressure.

Effects of gas cutting on water hammer pressure
When the difference between the formation pressure and the bottom hole pressure is 0.5 MPa and gas kick time was 14 min with 10 s of the shut-in time, respectively. The variation of water hammer pressure versus time with or without gas kick were shown in Figure 7. For this shut-in time, calculation results show that the maximum water hammer pressure with or without the gas kick is nearly the same. As can be seen from Figure 5, the mixture velocity near the surface whether to consider the impact of gas kick is approximately equal. It means that there was little change in mixture velocity within the same shut-in time, which can explain the phenomenon. While the BOP was totally closed, the attenuation trends and fluctuation amplitudes of water hammer pressure are quite different. The water hammer pressure is attenuated very fast with the increased of time under gas kick condition, and approached 0 MPa after 110 s. Nevertheless, the water hammer pressure decayed relatively slowly without consideration of gas kick, and approached 0 MPa after 150 s. This was mainly because the former taken the influence of both the free gas and the friction on water hammer wave attenuation into account, but the latter just the friction. Furthermore, the effect of free gas in the annulus on wave speed attenuation is greater on the order of magnitude than that of the friction.

Effects of shut-in time on water hammer pressure
When the difference between the formation pressure and the bottom hole pressure is 0.5 MPa and gas kick time was 14 min, respectively. The variation of water hammer pressure with time under different shut-in times was shown in Figure 8. It can be seen that the attenuation trends and fluctuation amplitudes of pressure are roughly similar. The changes of shut-in time have significant influence on water hammer pressure. The more the shut-in time is, the lower the maximum pressure generated by water hammer is. For instance, the maximum water hammer pressures are 1.45 MPa and 0.17 MPa at shut-in time 5 s and shut-in time 30 s, respectively. Obviously, when the shut-in time have changed from 5 s to 30 s, the maximum water hammer pressure is about 8.5 times higher than the time of 30 s. Additionally, it can be found that the maximum water hammer pressures are 0.21 MPa and 0.17 MPa at shut-in time 25 s and shut-in time 30 s, respectively. This implies that the peak value of water hammer pressure has no significant changes with further increase of the shut-in time when   it reaches a certain value, which in turn will cause additional influx from the formation.

Effects of well depth on water hammer pressure
When the difference between the formation pressure and the bottom hole pressure is 0.5 MPa and gas kick time was 14 min with 10 s of the shut-in time, respectively.
The variation of water hammer pressure with time at different well depths was shown in Figure 9. It is clear that the well depth is a negative correlation with water-hammer pressure. The wave amplitude decreases along the well depth from a maximum value at the surface to a minimum value at the bottom hole. The deeper the well depth is, the lower the water-hammer pressure is. The maximum water hammer pressure for the depth of 0 m and 6400 m are   approximately 0.68 MPa and 0.012 MPa, respectively. This happens because pressure wave travel and is attenuated due to energy loss through friction. Besides, the lower middle part of the wellbore is gas-liquid two-phase flow (shown in Fig. 4), where the existence of gas aggravates the attenuation of water hammer waves, further diminishing the water hammer pressure. However, it should be paid attention that, in most cases, the casing shoe is expected to be the weakest point in the open hole. The additional water hammer pressure generated by closing the BOP would fracture the shallowest exposed formation below the casing shoe, leading to an underground blowout.

Effects of gas kick time on water hammer pressure
When the difference between the formation pressure and the bottom hole pressure is 0.5 MPa and the shut-in time was 10 s, respectively. The variations of water hammer pressure with time under different gas kick times were shown in Figure 10. Obviously, the maximum water hammer pressure changed slightly but attenuated considerably with the increase of the gas kick time. Two reasons may account for the above for the phenomena. On the one hand, before the influx gas reaching the surface, there was a little variation in the mixture velocity (shown in Fig. 5). On the other hand, the formation gas entering the wellbore continues to migrate from bottom hole to surface. As a result, the leading edge of gas phase moves forward gradually upward with the increase of gas kick time. What is more, the gas void fraction increased due to volume expansion of invading gas in migration upwards, leading to the acceleration of the water hammer wave speed attenuation. Besides, during the shut-in stage, Shut-In Casing Pressure (SICP) refers to the difference between the formation pressure and hydrostatic pressure in the annulus. Despite that there is benefit to the attenuation of water hammer wave speed if the gas kick last for a long time, the bottom hole pressure has decreased by 0.67 MPa when the gas kick time increases from 14 min to 28 min (shown in Fig. 6). It means that the additional 0.67 MPa SICP would be exerted at the top of a wellbore. More seriously, if the SICP exceeds the maximum allowed SICP, it would be break the equipment or the formation. For the deep formation, however, the additional water hammer pressure applied to the formation caused by shutting in a well can be neglected.

Conclusion
In this work, considering the transient multi-phase flow characteristics and unsteady friction, MOC has been developed to model transient flow caused by shutting in the gas kick. The model was used to analyze effects of gas cutting, shut-in time, friction and gas kick time on water hammer pressure. The following are the results obtained: 1. It was found that the gas kick has no significant impact on maximum water hammer pressure, but the existing gas in the annulus remarkably reduces the water hammer wave speed. 2. The maximum water hammer pressure decreased as the shut-in time increased, but the shut-in time reaches a certain level, there is little effect of maximum water hammer pressure in case of further increasing shut-in time. 3. The compressibility of gas and the wall friction are beneficial to aggravating water hammer wave attenuation, demonstrating that additional water hammer pressure applied to the deep formation caused by shutting in a well can be neglected. Indeed, much more attention should be paid to the shallow casing shoe formation damage by the closure of BOP. 4. The SICP increase due to the influx gas is more serious; hence, it is advisable to shut in a gas kicking well timely upon occurrence of gas kick, if possible, to reduce the gas kick effects.