Multiphase flow behavior in deep water drilling: The influence of gas hydrate

Hydrate phase transition may bring hidden dangers for well control safety, even accidents of blowout happen in deep water drilling. Accurate calculation of all the drilling parameters is of great importance for well control safety. In this study, the mathematic models of transient heat transfer and multiphase flow have been established with the consideration of hydrate phase transition. Furthermore, the position where hydrate formation is most likely to occur was predicted with the proposed models, and the factors affecting hydrate formation region were also analyzed. Finally, the analysis of the effect of gas hydrate phase transition on some important parameters in well control, such as the pit gain, gas void fraction, and bottom hole pressure, is presented. The results show that hydrates form near the mud line during the drilling and shut‐in stage, but the hydrate formation region will be narrowed as well killing start. Additionally, the hydrate formation region decreased with the increase of flow rate and increased with the increasing shut‐in time, hydrate inhibitors, and riser insulation are presented to prevent hydrate blockage. More importantly, the gas void fraction and pit gain decrease due to the formation of hydrates, leading to the gas kick early detection not timely and gas kick more “hidden.”


| 1387
FU et al. decomposition of gas hydrates was also suggested to occur. Free gas will be formed in the wellbore annulus again, and it would bring hidden dangers for well control safety, even accidents of blowout happen [5][6][7][8][9] ; (b) hydrate blockages are the number one flow assurance and security problem. Hydrate blockages can form under various operation stages, such as circulation, shut-in, and well killing. 10,11 If hydrate blockages formed, wellbore and blowout preventer blocked, drilling fluid property changed and interruption of the normal operation, and it will bring serious disasters to the well control; (c) the concentration of hydrate and the shear rate of hydrate slurry have important effects on the morphology and rheological properties of hydrate slurry. [12][13][14] Besides, the rheological properties of drilling fluid in the riser can be also changed due to the low temperature environment of sea water, resulting an increase in both density and viscosity of drilling fluid. These will bring an extremely difficulty to the well control.
A number of mathematical models have been developed to simulate the well control process. The first well control model was developed by Leblanc and Lewis 15 ; however, the influence produced by gas-slip velocity and friction losses to annulus pressure is not taken into account. Subsequently, Records et al 16 improved the LeBlanc and Lewis model with the consideration of the annular friction losses, but not achieved a good agreement with the field data. Then, empirical correlations by experimental data were obtained to calculate the liquid holdup because of the importance for hydrostatic head, [17][18][19][20] which is only suitable for Newtonian fluid under steady flow condition. Due to its importance for the safety of well control, gas rise velocity in vertical and inclined well was experimentally investigated. [21][22][23] Meanwhile, significant progress has been made to acknowledge the mechanism of gas-fluid two-phase flow in wells. On the base of the one-dimensional mass, momentum, and energy conversation equations, several two-phase flow mechanistic models have been developed for simulating gas kick and were mainly classified into three classes: the homogenous model, the separated flow model, and the drift-flux model. The latter has been widely used in transient multiphase flow modeling for well control simulation due to its high accuracy. [24][25][26][27][28][29] Reviewing the previous studies, the available well control models can be used satisfactorily in well control simulation, but some actual situations cannot be properly evaluated yet due to model limitations. On the one hand, some well control models have ignored the effect of hydrate phase transition on hydrostatic head in deep water drilling. The gas kick will generally result in pit gain and bottom hole pressure less than those for conventional drilling; it could be more difficult to detect gas kick. Failure to detect a gas influx on entry to a well will make control of the kick when it reaches the surface more difficult. Besides, their models do not consider the effects of the drilling mud additives on hydrate formation. In essence, the recent studies showed that the drilling mud additives has a significant influence on hydrate formation in water-based mud. 30,31 However, this work is not the focus of the paper, which should be recommended in the future work. The main focus of this research is to develop the coupled multiphase flow modeling of a wellbore unsteady flow and a formation seepage with the consideration of the hydrate phase transition. On the other hand, many well control models ignore any temperature variations during the gas kick, but input a steady temperature variation with depth, namely the geothermal gradient. Additionally, effects of decomposition heat and heat of formation on wellbore temperature were seldom considered.
In this paper, a transient heat transfer model has been developed to simulate correctly and predict transient wellbore heat transmission with the consideration of heat of phase transformation under gas kick. Then, to predict pressure behavior inside the annulus as the influx gas moves from the formation into the wellbore and up to the wellhead, the coupled multiphase flow modeling of a wellbore unsteady flow and a formation seepage were established with the consideration of the hydrate phase transition. Finally, hydrate formation region was investigated under three operation stages, and several factors affecting the hydrate formation were analyzed. Moreover, the effects of hydrate phase transition on gas void fraction, pit gain, and bottom hole pressure are studied. And also, the influences of the parameter including formation permeability, initial differential pressure, and net pay thickness on gas kick development were analyzed.

| PHYSICAL MODEL
The physic model of fluid flow under gas kick during deep water drilling is shown in Figure 1. With the mud line as F I G U R E 1 Physical model of the wellbore in deep water demarcation lines, the upper part is surrounded by sea water, and the lower part is surrounded by formation. The drilling fluid was pumped from the mud pit into the drilling string where it sprays out of nozzles on the drill bit. Then, the drilling fluid moves upward from the bottom hole to the ground through the annulus space between the drill string and the sides of the hole being drilled. The mathematic models were set up based on the following assumptions.
1. One-dimensional fluid flow and heat transfer in a wellbore; 2. Bingham plastic model is used to describe the rheological behavior of non-Newtonian fluids, and the drift-flux model is used in this study; 3. Free gas, gas hydrate, and drilling fluid at the same location have identical temperature and pressure, and all flow parameters are the same at the same cross section; 4. Gas is compressible, and drilling fluid is incompressible; 5. The formation temperature increases with the increase of depth at a constant geothermal temperature; 6. Hydrate phase transition process was considered, and no mass transfer occurs between the drilling fluid and the natural gas; 7. Formation and sea water properties were assumed to be constant; 8. The convection heat transfer in the radial direction was taken into account, but the axial heat conduction inside the fluid was ignored.

| Wellbore temperature model
In deep water drilling, thermal changes are of great importance, as their results may be more significant on gas expansion and hydrate phase transition. Simultaneously, in the process of phase transition, materials can absorb or release heat energy and accordingly can affect the annulus fluid temperature. Therefore, contribution of hydrate transition to temperature should be considered in the heat conservation equations. According to Figure  1, the annulus section above the mud line is a riser surrounded by sea water and that below the mud line is surrounded by the formation. The fluid in the drilling string travels downward and exchanges heat with both the fluid in the annulus and sea water. Meanwhile, the annulus fluid travels upward and exchanges heat with both the drilling string fluid and the surroundings. Take a small control cell of the wellbore annulus along the axial direction as example, similarly for drilling string (Figure 2A,B). In Figure 2, Q a (z) is the heat flowing into microunit the axial direction, Q a (z +Δz) is the heat flowing out of the microunit the axial direction, Q ap is heat transferring from drilling string fluid to annulus fluid in radium direction, Q sa is heat transferring from sea water to annulus fluid in radium direction, Q sa is heat transferring from formation to annulus fluid in radium direction, Q fa is the friction heat in the annulus, ±Q hyd is the heat of phase transformation, and note that "+" and "−" represents hydrates decomposition and hydrate formation, respectively.

| Governing equations for heat transfer inside the annulus
According to the Figure 2A,B, the fluid temperature in the annulus was determined by the drilling string fluid and the heat generated by the friction, and the heat of phase transformation is also included. Notably, both the seawater and the formation rock get involved in the heat transmission with the annulus drilling fluid. Based on the first law of thermodynamics, the equation of energy conservation for the drilling fluid inside the annulus above and below the mud line are, respectively, described as 6,10 : (1) where d ri is the internal diameter of riser, m; d po is the outside diameter of drilling string, m; g is the gas density, kg/m 3 ; H g is the gas void fraction; l is the drilling fluid density, kg/m 3 ; H l is the liquid holdup; c l is the specific heat of drilling fluid, J/(kg °C); c g is the specific heat of gas, J/(kg °C); T p and T a are the temperature of drilling fluid in the drilling string and annulus, respectively; °C, T s , and T f are the temperature of sea water and formation, respectively; °C and Q fa are the heat quantity caused by friction during drilling fluid flows in the annulus, W/m; q l and q g are the mass flow of drilling fluid and gas, respectively, kg/s; r h is the formation and decomposition velocity of hydrates, kg/(m s); M H is the average molecular weight of hydrates, kg/mol; Δh is the decomposition heat of hydrates, J/mol; t is the time, s; z is the axial position, m; d w is the wellbore diameter, m; U s is the overall heat transfer coefficient between annulus and sea water, W/ (m 2 K); U a is the overall heat transfer coefficient between annulus and drilling string, W/(m 2 K); and U f is the overall heat transfer coefficient between annulus and formation, W/(m 2 K). The overall heat transfer coefficient U s is defined as (the details of the derivation are given in Appendix A): where h ri is the convective heat transfer coefficient of the inner wall for riser, W/(m 2°C ); h ro is the convective heat transfer coefficient of the outer wall for riser, W/(m 2°C ); d ro is the outside diameter of riser, m; and r is the thermal conductivity of the riser, W/(m °C).
The overall heat transfer coefficient U a is defined as: where h po is the convective heat transfer coefficient of the outer wall for drilling string, W/(m 2°C ); h pi is the convective heat transfer coefficient of the inner wall for drilling string, W/(m 2°C ); and dp is the thermal conductivity of the drilling string, W/(m °C). The overall heat transfer coefficient U f is defined as: where h f is the convective heat transfer coefficient of the borehole wall, W/(m 2°C ); n is the number of concentric casing and cement sheath around the annular fluid at a certain position in the wellbore; d jo is the outside diameter of the concentric medium j, m; d ji is the internal diameter of the concentric medium j, m; and j is the thermal conductivity of the concentric medium j, W/(m °C).

| Governing equations for heat transfer inside the drilling string
Similarly, the fluid temperature in the drilling string was determined by the annulus fluid and the heat generated by the friction. According to the first law of thermodynamics, the equation of energy conservation for the drilling fluid in the drilling string can be expressed as 32 : where d pi is the internal diameter of drilling string, m; Q fp is the heat quantity caused by friction during drilling fluid flows inside the drilling string, W/m.

| Formation heat transfer model
By the use of a cylindrical coordinate system, the thermal conduction differential equation in the formation can be described as 33 : where c f is the specific heat capacity of formation, m; f is the heat conductivity of formation, W/m; f is the formation density, kg/m 3 ; and r is the formation radius, m.

| Multiphase flow model
To predict pressure behavior at any time and position in the annulus when hydrate phase was taken into account, conservation of drilling fluid mass, mass of free gas, and mass of hydrate is developed separately. Concerning of momentum, only momentum conservation of the total mixture is considered. 6,10 The details of the derivation of the control equations are illustrated in Appendix B. Conservation of gas-phase mass (gas-producing zone): where q g is the gas cutting speed, kg/(m s); g is the gas density, kg/m 3 ; H g is the gas void fraction; u g is the gas velocity, m/s; is the cross-sectional area of annulus, m 2 ; and x g is the mass fraction of free gas. Conservation of gas-phase mass (non-gas-producing zone): Conservation of liquid-phase mass: where l is the drilling fluid density, kg/m 3 ; u l is the drilling fluid velocity, m/s; and H l is the liquid holdup. Conservation of hydrate-phase mass: where H is the hydrate density, kg/m 3 ; u H is the hydrate velocity, m/s; and H H is the volume fraction of hydrate.

Conservation of total momentum:
where F r is the friction pressure loss, MPa; P is the pressure, MPa; and g is the acceleration of gravity, m/s 2 .
The drift-flux model (DFM) for two-phase flow is used in this study, and the relationship between the velocity of the two phases is expressed by 34 : where u m is the mixture average velocity, m/s. The parameters C o and u gr are the distribution coefficient and drift velocity of gas relative to liquid, which are determined according to the flow pattern. [35][36][37][38] The Bingham plastic model is adopted for describing flow characteristics of drilling fluid, and the frictional pressure loss (F r term in Equation (12)) for both single-phase and twophase flow is giver by: where m is the average density of the mixture, kg/m 3 ; D out is the borehole diameter, m; D in is the outer diameter of the drilling string, m; f r is the friction factor. For single-phase flow, the friction factor f r will have different expression depending on whether the flow is laminar flow or turbulent flow. Thus, the single-phase friction factors are determined by f r = 24 38 where Re g is the generalized Reynolds number; Δ is the absolute roughness, m; D eff is the equivalent diameter, m; and n ′ is the generalized flow behavior index. However, the frictional pressure loss calculation of the two-phase flow pattern may refer to the models of Hasan and Kabir. 35,37 The Minton-Bern model is used to calculate the apparent viscosity 39 : where o is the viscosity at P o MPa and T o °C, mPa s； and are empirical constants.
The gas influx rate from the reservoir can be calculated from the binomial theorem Equation 40 : where P e is the formation pressure, MPa; P wf is the bottom hole pressure, MPa; K is the reservoir permeability, mD; h e is the net pay thickness, m; T is the average temperature, °C; 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; and is the turbulent coefficient, m −1 .

| Hydrate phase equilibrium and decomposition rate equations
Reliable hydrate phase behavior predictions are critical to control the formation of natural gas hydrates. In this paper, the classic vdW-P model proposed by van der Waals and Platteeuw is adopted to predict the hydrate phase equilibrium condition. 41 where Δ 0 is the chemical potential difference of empty hydrate lattice and pure water under standard condition, J/mol; R is the gas constant, J/(mol K); T o is the temperature under standard condition, K; T H is the phase temperature of hydrate formation, K; ΔC K is the specific heat capacity difference between empty hydrate lattice and pure water, J/(mol K); P H is the phase pressure of hydrate formation, Pa; P o is the pressure under standard condition, Pa; ΔH o is the enthalpy difference between pure water and empty hydrate lattice, J/kg; ΔV is the specific volume difference between pure water and empty hydrate lattice, m 3 /kg; f w is the water fugacity in water-rich phase, Pa; f wr is the water fugacity in reference state T H and P H , Pa; l is the number of hydrate types; M i is the cavity number of i-type hydrate in unit water molecule; L is the component number of hydrate; and ij is the occupancy of j-type guest molecule in i-type cavity.
If inhibitors are added into the drilling fluid, the relationship between f w and f wr is expressed as: where x w is the water molar fraction in the water-rich phase and y w is the water activity coefficient in the water-rich phase.
Similarly, Kim-Bishnoi model is adopted to calculate hydrate dissociation rate. 42 where n H is the total moles of CH 4 contained in hydrate particles; A s is the surface area of one hydrate particle, m 2 ; K d is the hydrate decomposition rate constant, mol/(m 2 Pa s); f eq is the three-phase equilibrium fugacity, Pa; and f is the gas fugacity, Pa.
The classic V-B model proposed by van der Vysniauskas and Bishnoi is adopted to calculate hydrate formation rate. 43,44 where r is the rate of methane consumption during hydrate formation, m 3 /s; A c is a lumped pre-exponential constant; a s is the total surface area of the gas-water interface, m 2 ; ΔE a is the energy of activation, KJ/(g mol); T is the system temperature, °C; a is the arbitrary constant; b is the arbitrary constant; ΔT is the super-cooling, °C; and is the overall order of reaction with respect to pressure.

| boundary conditions
The inlet temperature of the drilling fluid in the drilling string can be measured, and the boundary condition can be described as: where T in is the inlet fluid temperature in the drilling string, °C.
The fluid temperature in the drilling string is the same as that of the annulus at the bottom hole, and it can be defined as Equation (23).
where H is the well depth, m.
The formation temperature profile is linear based on the geothermal gradient, and it is also assumed that at the outer boundary of the formation temperature does not change. The initial conditions are given in Equations (24) and (25), respectively.
where h sea is the depth of the sea water, m and G f is the geothermal gradient, °C/m.
At the interface between the wellbore and formation, the heat of flux out of the formation is equal to the heat of flux into the annulus, which can be defined as Equation (26).
During the drilling process, the annulus at the surface is open and its pressure is always equal to the atmospheric pressure, which is expressed as Equation (27).
At the depth of gas kick, the boundary conditions are defined as: where u sg is the superficial velocity of gas phase, m/s; u sl is the superficial velocity of liquid phase, m/s; h is the depth of some point in wellbore, m; and Q mud is the flow rate, m 3 /s.

| Initial conditions
Before the normal drilling, the wellbore temperature is equal to the surrounding seawater temperature above the (19) ln mud line; the wellbore temperature is equal to the surrounding formation temperature below the mud line. The initial conditions are described in Equations (29) and (30), respectively.
The end temperature profiles during the normal circulation stage are defined as the initial condition of temperature profiles during the gas kick stage. The initial conditions were defined as: where t end is the time of normal circulation stage, s.
Initially, there will be only single fluid phase existed in the wellbore. Thus, the primary state variables can be described as:

| RESULTS AND DISCUSSIONS
To analyze the influence of the hydrate phase transition on the annular flow parameters when gas kick occurs during deep water drilling, thus the L1 well located in the South China Sea is taken as an example, and here assume that the measured depth of the well was 3715 m at the time of the gas kick. The wellbore data and casing program were used as the input parameters of our model, which are shown in Tables 1 and 2, respectively. The thermophysical parameters of the heat transfer mediums are summarized in Table 3. The sea water temperature profile used correlations in literature were developed by Xu et al. 45

| Hydrate formation prediction
Before the effect of hydrae phase transition on the behavior of multiphase flow was investigated, it is necessary to determine hydrate formation region in the annulus at first. The hydrate formation region in the annulus can be obtained by solving thermodynamic balance equation, temperature governing equations and multiphase flow governing equations; that is, the intersecting region between annular temperature curve and phase state curve is the hydrate formation region.

| Different drilling conditions
The hydrate formation regions in the annulus under different drilling conditions are shown in Figure 3. During circulation stages, no formation gas entrance into the well under non-gas kick conditions, resulting in no formation of hydrate. Conversely, under gas kick conditions, the invading gas will reduce the specific heat and thermal conductivity of the drilling fluid and release heat after gas hydrate formation, leading to the annulus temperature changes in comparison with that of non-gas kick. Hence, there is closure region between the annulus temperature curve and phase state curve.
During the shutting-in stages, due to the annulus fluid was in static state, heat transfer between surroundings and the fluid in the annulus was mainly through heat conduction, which eventually make the annulus fluid temperature to get close to the ambient temperature. Therefore, the closure region between phase state curve and annular temperature curve is further expansion. During the well killing stages, the fluid in the annulus is again in flowing state. With the migration of the fluid in the annulus, the annulus fluid temperature in the open-hole section was lower than the annulus fluid temperature at the shutting-in stages. However, the heat was continuously carried from borehole to ground when the fluid moves upward from the bottom hole to the ground, leading to the elevation of annulus fluid temperature. Consequently, compared with the results of the shutting-in stages, the hydrate formation region in the annulus decreased.

| Sensitivity analysis
In order to investigate the hydrate formation region under different parameters, and the influence degree of each parameter on the annulus fluid temperature, we performed a sensitivity analysis in this section based on the developed model.

Flow rate
The hydrate formation regions at different flow rates during drilling process are shown in Figure 4. Obviously, the annular fluid temperature above mud line increased with the increasing of flow rate gradually. In contrast, the annulus fluid temperature in the open-hole section gradually decreased with the flow rate increasing. This is attributed to the fluid in the annulus being not sufficiently heat exchange with seawater and formation due to the increase of flow rate. Hence, the annulus fluid temperature has a tendency to get close to the ambient temperature in case of low flow rate. Noticeably, for Φ311.1 mm borehole, the normal flow rate is 55-60 L/s, while the well killing rate is generally 1/2-1/3 of the normal flow rate, that is, 20-30 L/s. As seen from Figure  4, it can be concluded that the higher the flow rate is, the smaller the hydrate formation region will become. This means that an increasing flow rate is conducive to inhibiting hydrate formation.

Shut-in time
The hydrate formation regions under different shutting time during shut-in stage are shown in Figure 5. Above the mud line, the annulus fluid temperature decreased as the shut-in time increased. This occurred because the annulus fluid temperature was higher than the seawater temperature during the circulation stage, as shown in Figure 4. Seawater absorbs heat from the annulus fluid by heat conduction, reducing the annulus fluid temperature. In the open-hole section, the annulus fluid temperature increased with the increase of shut-in time. That is because the annulus fluid temperature was lower than the formation temperature. The increase in annulus fluid temperature was the result of the heat quantity being transferred from the annulus to the formation by heat conduction. As show in Figure 5, it can be concluded that the longer the shut-in time is, the larger the hydrate formation region will become. It is easier to accumulate hydrates near the mud line, leading to the blowout preventer blocked. Consequently, in order to avoid hydrate formation during deep water drilling, shut-in time shall be shortened as much as possible.

Well killing time
The hydrate formation regions under different well killing time when the well kill rate is 20 L/s after shut-in for 4 hours are shown in Figure 6. The annulus fluid temperature gradually deviates from the ambient temperature with the increasing of well killing time. As illustrated in Figure 5, the annulus fluid temperature recovered to its initial ambient temperature after a long period of shut-in operation. After carrying out well killing operations, the annular fluid changing from static to flowing state, the heat quantity was carried to the surface because of the annulus fluid flowing up when the well killing time was increased, leading to an increase in the annulus fluid temperature above the mud line. As the same reason, the annulus fluid temperature in the open-hole section was decreased with the increase of time. Combined with the results of Figure 5, the hydrate formation region decreased during the well killing stages. It means that the hydrate decomposing, which led to potential safety hazards.

Hydrate inhibitors
The hydrate formation regions, using different sodium chloride concentration (C Nacl ) and ethanol concentration (C CH 3 CH 2 OH ) at 60 L/s flow rate, are presented in Figures 7 and 8. The phase state curve moves downward as both the sodium chloride concentration and ethanediol concentration increased, resulting in the decrease of hydrate formation region. Remarkably, when the sodium chloride concentration and ethanediol concentration are up to 5% and 3% respectively, the hydrate cannot form in the wellbore, which demonstrates that among these common inhibitors, the effect of ethanediol is the best. Therefore, a certain concentration of hydrate inhibitors can be added to avoid hydrate formation in the wellbore.

Riser insulating layer
The hydrate formation regions with or without riser insulating layer during shut-in stage after shut-in for 4 hours are shown in Figure 9. If the riser was not covered by the insulating layer, the minimum annular fluid temperature above the mud line is approximately 5°C. Conversely, the temperature will reach to 16°C, and the temperature difference can be 11°C. This is due to the thermal conductivity coefficient of riser insulating layer was 1902 times lower than that of riser (Table 3), which led to the heat quantity being transferred from the annulus to the seawater by heat conduction was little. Hence, it is observed that there was little variation in annulus fluid temperature when the rise is covered by insulating layer. Moreover, the closure region between the phase state curve and the annular fluid temperature curve significantly narrowed, indicating that the hydrate formation region decreased. Consequently, it will prevent hydrate from forming more effectively with the use of heat-insulating layer.

| Comparisons of different multiphase flow models
The work of Wang et al is generally considered the first to propose a comprehensive model for well control in deep water drilling. Wang et al 6 established the basic hydrodynamic models, including mass, momentum, and energy conservation equations, for annulus flow with gas hydrate phase transition during gas kick. The basic calculation data derived from the manuscript are described in detail in Section 4. Other important parameters in this case are as follows: the initial difference between the formation pressure and the bottom hole pressure is 0.5 MPa, the formation permeability is 30 mD, the net pay thickness is 2 m, and the constant influx gas velocity is 0.01 m 3 Figure 10 shows the comparisons of influx gas velocity when the front of the gas reaches the wellhead. It could be found that the influx gas velocity between the proposed model and the Wang et al model was completely different. The influx gas velocity was not a constant value but increased exponentially with the increased of the gas kick time in the proposed model. This phenomenon can be explained as follows: During the initial stage of gas kick, there is a little difference between the formation pressure and the bottom hole pressure, and thus, the influx gas velocity is small. Meanwhile, the gas will continue to expand when it was circulated from the bottom hole to the wellhead, which created a reduction in hydrostatic pressure because the gas density is considerably lower than the drilling fluid. As a result, the difference between the formation pressure and the bottom hole pressure increased, and then the influx gas velocity becomes higher and higher.  Figure 11 shows the comparisons of bottom hole pressure when the front of the gas reaches the wellhead. It can be observed that the variation regularities of bottom hole pressure are basically similar in cases with and without consideration of the wellbore unsteady flow not being coupled with formation seepage. During the initial stage of gas kick, there is no obvious difference between these two models. Note that the lower bottom hole pressure calculated by the Wang et al model was obtained because the given influx gas velocity is higher than the results of the proposed model. Subsequently, the trend of the bottom hole pressure calculated by two models showed in opposite direction. The reason is that the gas began to expand as it is circulated up the annulus, leading to the decrease of hydrostatic pressure. As a result, the difference between the formation pressure and the bottom hole pressure will be increased, which cause increase in influx gas velocity as showed in Figure 10. The higher the influx gas velocity is, the more the bottom hole pressure is reduced, and further accelerate the increase of influx gas velocity, and thus form a "vicious" cycle. As shown in Figure 11, the difference between the two results which are calculated by the two models, respectively, had been widening with the increase of gas kick time. During the final stage of gas kick, however, the difference between the two results decreased. This was mainly because the initial influx gas velocity given by Wang et al model is higher than that of our model (shown in Figure 10), meaning that more formation gas entrance into wellbore. Once it was circulated up near the wellhead, the gas volume expansion will be larger and then created a larger reduction in the bottom hole pressure.

| Effects of hydrate phase transition on gas void fraction
The variations of gas void fraction with gas kick time and well depth are shown in Figure 12. The sharp reduction in the gas void fraction at 2900 m and 1454 m was attributed to the variation of cross-sectional area. When the well depth is less than 1454 m, the annulus was formed by the Φ508 mm riser and the Φ149.2 mm drilling string; when the well depth was 1454-2900 m, the annulus was formed by the Φ339.7 mm casing and the Φ149. . Therefore, the sudden increase of the cross-sectional area reduced the gas void fraction. As seen from Figure 12, in 30 minutes, the variations of gas void fraction are basically the same in cases with and without consideration of the hydrate phase transition. That was because the influx gas migration time to the hydrate formation region was 30 minutes. In this situation, no gas hydrate formed and the gas void fraction distribution was the same. Then, the gas void fraction decreased gradually because of hydrate formation compared with the case of no hydrate phase, which would make it difficult to the early detection of gas kick. Subsequently, the gas will leave the hydrate formation region as it migrated to the wellhead, so the gas void fraction increased due to the gas expansion. Besides, when the gas kick time was 60 minutes, the influx gas just migrated to the wellhead and the gas void fraction exhibited little difference between the two situations. Two reasons might contribute to the phenomenon. On one hand, during the drilling process, the wellhead is open and its pressure is always equal to the atmospheric pressure. On the other hand, the annulus fluid temperature gradually increased near the wellhead, as shown in Figure 3.

| Effects of hydrate phase transition on pit gain
The pit gain vs gas kick time is shown in Figure 13. The pit gain gradually increased with the increase of gas kick time, showing an exponential function of time. It should be noted that the differences begin between the two curves at the gas kick time of 30 minutes, primarily because the gas enters the hydrate formation region, as shown in Figure  12. Subsequently, the pit gain difference became larger over time on account of hydrate formation consumption of different amounts of free gas. Furthermore, pit gain was the monitoring index of gas kick, and the hydrate phase transition delayed the detection of the gas kick. Taking the pit gain 1 m 3 as an example, the time when pit gains reached 1 m 3 was 52 minutes and 55 minutes, respectively, with and without the consideration of the hydrate phase transition, as observed in Figure 13. Therefore, the detection of the gas kick was delayed by 3 minutes, which led to the early gas kick concealed and more difficult for detection during deep water drilling.

| Effects of hydrate phase transition on bottom hole pressure
The bottom hole pressure vs gas kick time is shown in Figure  14. It is observed that the bottom hole pressure gradually decreases linearly at the beginning of gas kick, but then decreases rapidly over gas kick time. The bottom hole pressure was composed of hydrostatic pressure and friction pressure caused by drilling fluid movement through the annulus. The influx gas will expand as it flowing from well bottom to the ground, resulting in a reduction in hydrostatic pressure. Simultaneously, the mixture average velocity was accelerated by the influx gas as it rises up the annulus, which indirectly led to the increase of friction loss. However, it is generally acknowledged that the impact of gas expansion on bottom hole pressure is greater than the latter. So the bottom hole pressure will continue to decrease over gas kick time. Remarkably, from the results in Figures 12 and 13, the influx gas flowing from well bottom to hydrae formation region was 30 minutes, so the bottom hole pressure is similar. When the gas kick time is more than 30 min, it is easy to discern that the hydrate formation would consume a part of gas that enters the hydrate formation region near Besides, the gas consumed by phase transition is relatively small vs the riser with a large inner diameter. Therefore, the difference of the bottom hole pressure was only 0.1 MPa at the gas kick time of 60 minutes.

| Effects of the key parameters on well control
The gas void fraction distribution along the wellbore under different formation permeabilities was shown in Figure 15 when the initial difference between the formation pressure and the bottom hole pressure is 0.5 MPa and net pay thickness is 2 m, respectively. It was displayed that the gas void fraction distribution along the wellbore with different formation permeabilities got very large differences. The gas void fraction at a given depth increased with the increase of formation permeability. Because there'll probably be more formation gas entering into the wellbore within the same time with the increase of formation permeability. Additionally, the bottom hole pressure will be reduced more rapidly with large influx gas due to gas expansion, which, in turn, lead to the higher influx gas velocity. For instance, the gas void fractions at bottom hole are 3.1%, 4.5%, and 6.5% at formation permeability 10 mD, formation permeability 20 mD, and formation permeability 30 mD, respectively. Simultaneously, the pit gain will be risen during the phase of gas circulating form bottom hole to the wellhead as shown in Figure 16. In addition, as Figure  16 shows that the peak value of pit gain was obtained when the gas front reaches to the surface, mainly due to the effect of gas expansion. It was also found that the trends of pit gain were almost the same under different formation permeabilities, but increased obviously with formation permeability. It likely indicates that the effect of formation permeability on pit gain is significant. Figure 17 shows the variation of bottom hole pressure with different formation permeabilities. The bottom hole pressure shows an opposite trend compared to the pit gain. The bottom hole pressure increased more rapidly in high formation permeability than that in low formation permeability. The changes in gas void fraction distribution support the bottom hole pressure variations in Figure 17. The higher the gas void fraction distribution along the wellbore, the more is the bottom hole pressure decreased. Therefore, the gas kick development will be accelerating if high permeability formations are encountered, which make the bottom hole pressure drop rapidly and pose more challenges to well control. Also, the study of the effects of differential pressures and net pay thickness on well control was presented in Figures 18 and 19, respectively. Obviously, it can be seen that both the differential pressure and the net pay thickness have a similar effect as the formation permeability. The pit gain increased gradually with an increase both in differential pressure and net pay thickness, whereas the bottom hole pressure decreased.
In these cases, when drilling an oil/gas well into a formation, well control operations should be done to handle a gas influx as quickly as possible if a gas kick occurs. It also indicates that the early kick detection should be enhanced to prevent the risk of well blowout during deep water drilling.

| CONCLUSIONS
In this work, considering the transient multiphase flow characteristics of deep water drilling, the mathematical model was built for deep water gas kick with consideration of the hydrate phase transition. The model was used to analyze effects of flow rate, shut-in time, well killing time, hydrate inhibitors, and riser insulating layer on hydrae formation region. The influences of hydrate transition on gas void faction, pit gain, and bottom hole pressure have been investigated. The following are the results obtained: 1. The different drilling conditions will lead to the changes of annulus fluid temperature, which results in a variation of the hydrate formation region. During circulation stage, there is hydrate formation region under the gas kick condition. The hydrate formation region increased during the shut-in stage and decreased during the well killing stage. would be reduced compared with that of no hydrate transition, which made the gas kick became "hidden," resulting in the gas kick not being able to be detected in time. 4. The bottom hole pressure gradually decreased as the gas kick time increased, whereas there were few variations of the bottom hole pressure in cases with and without consideration of the hydrate phase transition. 5. It is found that, with an increase in formation permeability, initial differential pressure, and net pay thickness, the gas void fraction and pit gain increased, whereas the bottom hole pressure decreased. Thus, early kick detection should be further strengthened.