Evolution of gas kick and overflow in wellbore and formation pressure inversion method under the condition of failure in well shut-in during a blowout

Guo-Shuai Ju，Tie Yan，Xiao-Feng Sun，Jing-Yu Qu，Qiao-Bo Hu

College of Petroleum Engineering，Northeast Petroleum University，Daqing，163318，Heilongjiang，People's Republic of China

Keywords:Gas kick Formation pressure Multiphase flow Computational model Long short-term memory

ABSTRACT With ongoing development of oil exploration and techniques，there is a significant need for improved well control strategies and formation pressure prediction methods.In this paper，a gas-liquid transient drift flow model was established according to the gas-liquid two-phase flow characteristics during the gas kick.A Roe scheme was used for numerical calculation based on the finite volume method.The changes of bottom-hole pressure，casing pressure，the development law of cross-sectional gas holdup，and gas velocity，along with the vertical well depth，were analyzed through simulation examples.The time-series characteristics of mud pit gain were obtained by adjusting the formation parameter.The complex nonlinear mapping relationship between the formation parameters and the mud pit gain was established.The long short-term memory network (LSTM) of deep learning was used to obtain a formation pressure inversion when the blowout is out of control and the well cannot be shut-in.Experimental data from a well were used to verify the gas-liquid two-phase transient drift flow model based on the finite volume method，demonstrating that this method is reliable，with greatly improved prediction accuracy.This approach provides theoretical support for the early monitoring of gas kick during drilling，and for well-killing design and construction after uncontrolled blowout.

1.Introduction

With the expansion of oil and gas exploration and increased construction of complex wells in deep，deepwater/ultra-deepwater，high-temperature，and high-pressure areas，there are increased requirements for safety and risk management for well drilling and well control(Bhandari et al.，2015;Sun et al.，2018;Xu et al.，2020).Therefore，it is important to study the development law of gasliquid two-phase transient flow characteristics in the wellbore during the gas kick and to predict the formation pressure when blowout is out of control and the well cannot be shut-in.This information can help ensure the safety of drilling，improve drilling efficiency，reduce drilling cost，and allow the restoration of control with improved well-killing design and construction in response to blowout.

The coupling process of gas-liquid two-phase flow plays an essential role in the process of gas kick.The gas-liquid two-phase flow in the wellbore has typically been studied using the drift flow model，which is composed of continuity equations of each phase and a mixing momentum equation.This model directly generates a set of equations to study the mathematical structure of the twophase flow，and has advantages of simplicity and transparency(Kulia et al.，2015，2016;Zeidan and Sekhar，2018;Shen，2020).Ishii and Hibiki(2011)studied the motion and constitutive equations of drift velocity under different two-phase flow patterns.Considering macroscopic factors such as the geometric parameters of the interface，the volume force field，the shear stress，and the momentum transfer at the interface，they established a constitutive equation describing the relative motion between phases in the drift flow model.Evaluating the drift flow model's practical significance for analysis of two-phase flow，Hibiki and Ishii (2003) determined the two-phase flow pattern distribution parameters and drift constitutive velocity equation in microgravity.Fjelde et al.(2016)established a transient calculation model for the numerical simulation of two-phase flow，effectively reducing numerical dissipation and discretization errors based on the drift flow model and the AUSMV scheme.Wang et al.(2016) proposed a new algorithm to solve the drift flow model by referring to the SIMPLE algorithm and the calculation results are in good agreement with those calculated by the Roe method.Ma et al.(2016)established a multiphase flow model combined with a well control hydraulic model of pressurecontrolled drilling based on the drift flow model.This model retained transient multiphase flow characteristics of fluid and gas in the wellbore and adopted an appropriate closure relationship.Wei et al.(2018) used the MUSCL solver framework to establish a gas-liquid two-phase transient flow model based on the drift flow model，described the gas-liquid phase relationship with Shi's relationship，and solved it through the second-order AUSMV numerical scheme.

For uncontrolled blowout caused by a severe gas kick，the conventional method of determining formation pressure by reading the casing pressure and standpipe pressure after shut-in the well will not work.Formation pressure is an essential key parameter in a well-killing design after a blowout is out of control(Anderson et al.，2011;Li et al.，2016;Meister et al.，2003).The interactive interpretation method can predict formation pressure in combination with the two-phase flow calculation model，continuously adjusting bottom hole parameters during forward modeling to obtain prediction results with high consistency with measured wellhead parameters.However，this method is a trial and error method based on manual optimization，requiring a large amount of calculation and with limited processing accuracy(Wang et al.，2021).In recent years，deep learning has attracted extensive attention for modeling.As a kind of recurrent neural network with memory ability，long short-term memory (LSTM) can intelligently process the timeseries evolution characteristics of data (Kratzert et al.，2018;Yildirim et al.，2019;Zhao et al.，2017).In this study，time-series evolution characteristics of wellhead mud pit gain were determined by adjusting the bottom hole parameters and using a transient calculation model.The corresponding training set was then established and input into the LSTM network.Measured time series evolution data of mud pit gain were then used with the trained LSTM network to realize the inversion and prediction of formation pressure.

The finite difference method has been widely used to solve the drift flow model.This method has low conservation，relies on high computational cost to satisfy the conservation relationship，and is accessible to divergence in the iterative calculation.There is no attempt to retrieve formation pressure in case of blowout out of control and unable to shut in the well.The earlier a gas kick is detected，the higher the possibility of successful well control.In this study，a gas-liquid two-phase transient drift flow model based on the finite volume method was established.The Roe scheme was used for numerical calculation to analyze gas kick and overflow development law in the wellbore.LSTM method based on deep learning was used for the inversion of formation pressure when blowout is out of control.This work can provide theoretical support for the early monitoring of gas kick during drilling and for wellkilling design and implementation after a blowout.

2.Mathematical model of gas-liquid two-phase transient flow

The gas-liquid two-phase transient flow in the wellbore after gas kick can be calculated as follows:(1) According to the continuity equation of gas phase and liquid phase and the mixture momentum equation of gas-liquid phase，a gas-liquid two-phase transient drift flow model can be established in matrix form.(2)Jacobian matrix transformation is then carried out for the equations.According to the chain rule，the sub-elements of the Jacobian matrix A can be derived by constructing intermediate functions.(3)The approximate linearization matrixof matrix A is obtained，and the complex nonlinear equations can be transformed into simpler linear equations.(4) The Roe scheme is then used to find the relationship between matrix A and.(5) The Roe flux difference scheme is used to solve the conserved variables，and these variables can then be updated on the time layer.

2.1.Governing equation of gas-liquid two-phase transient flow

The governing equations of gas-liquid two-phase transient flow in the wellbore include the continuity equations and the momentum equation.Continuity equations of gas phase and liquid phase are:

where αgand αlare the volume fraction of gas phase and liquid phase，respectively;ρgand ρlare the density of gas phase and liquid phase，respectively;vgand vlare the velocity of gas phase and liquid phase，respectively;p is the pressure;ρmis the mixing density of gas phase and liquid phase;Ffricis the friction term;g is the acceleration of gravity;θ is the inclination angle.

Because the interaction between the gas and liquid phase is ignored，it is necessary to introduce the gas-liquid drift relation Ffricto close the drift flow model.Shi's gas-liquid drift relationship has been widely used in the drilling industry and verified by many actual data (Shi et al.，2005).

The drift flow model equations of gas-liquid two-phase flow in vector form can then be obtained by simultaneous Eq.(1)-(3):

These governing equations are a class of typical hyperbolic nonlinear partial differential equations.The use of these equations may generate problems such as weak discontinuity，discontinuous derivatives of the solution，a discontinuous solution，and the strong discontinuity of the function itself.Additionally，a numerical method is required to solve the problem.In this paper，we use the finite volume method to solve the nonlinear governing equations based on the Roe scheme (Roe，1981)，including Jacobian matrix transformation and solving the approximate linearization matrix.The step-by-step calculations are detailed below.

2.2.Transformation of the Jacobian matrix

Without consideration of the source term S(W)，Eq.(4) can be written as a homogeneous equation by applying the chain rule:

The transformed Jacobian matrix A(W) is:

Due to the strong nonlinearity of the drift flow control equation's coefficient matrix A(W)，its eigenvalues and flow propagation direction cannot be determined accurately.Therefore，it is difficult to construct the upwind difference of the nonlinear vector equation.To address this，an intermediate function can be constructed and through chain rule，every submatrix of the Jacobian matrix can be obtained.Because the derivation process is lengthy，the derivation results of each sub-matrix are given directly here.

Among them，ε1，ε2，and ε3are parameters without physical meaning，and are used only to simplify the above formula's expression.

where C0is the gas distribution coefficient.

Using these equations，the expression of the Jacobian matrix A(W) is obtained.

2.3.Calculating the approximate linearization matrix of Jacobian matrix

According to the gas-liquid two-phase transient drift flow model，the splitters of the linearized matrix can be calculated from the continuity and momentum equations.Similarly，the split form of flux term F(W) can be obtained.

where L and R are the left and right sides of the flux.

Next，the Roe average values of the continuity equation，momentum equation，pressure term，void fraction，liquid fraction，and interphase slip relation can be calculated.Since there are many subelements in matrix splitting design，resulting in a cumbersome calculation process，the complex expression of splitting sub-terms is not presented here.Using the above steps，the approximate linearization matrix of the original Jacobian matrix can be obtained.

The relationship between gas density，pressure，and temperature is used as the auxiliary equation of the model.The auxiliary equation has three parts:the density relationship of compressible gas under different temperature and pressure conditions(ideal gas equation of state)，the relationship between pressure and viscosity of gas-phase(Dean-Stiel viscosity model(Tan et al.，2017))，and the general flow pattern discrimination relationship of gas-liquid twophase (including dispersed bubble flow，bubble flow，slug flow，stirred flow，annular fog flow，etc.).

2.4.Boundary conditions and the Roe algorithm to solve equations

The grid division of the wellbore is shown in Fig.1，with the bottom and wellhead set as the boundary conditions of the flow region.To simulate gas kick，the parameters of the bottom hole and wellhead will change with the change of gas-liquid two-phase coupling relationship in the wellbore，correlating to a change of boundary conditions.During drilling，when gas invasion occurs in the wellbore or if the wellhead back pressure is adjusted，the boundary conditions change dynamically，and the changes of these boundary conditions alter the flow parameters in the whole wellbore grid space.The boundary conditions can be obtained from the characteristic lines and compatibility equations.

The annular flow field inlet is set as the pressure-inlet boundary condition.The drilling fluid flows along the bottom hole inlet direction，and the inlet velocity is evenly distributed.The annulus outlet is set as the boundary condition of the pressure-outlet，with unidirectional flow back along the annulus.

According to the Roe flux difference scheme，the transient drift flow equations can be solved as:

Through the Jacobian matrix and its approximate linearization matrix obtained by the previous calculation，the conserved variables can then be calculated by the Roe iterative scheme according to Eq.(25).

The calculation scheme is as follows (see Fig.2):

(1) Input basic parameters，including well trajectory，wellbore structure，formation parameters，and physical parameters of drilling fluid and natural gas.

(2) Discretize the wellbore flow region and divide into onedimensional grids to determine the total simulation time，time step，and convergence conditions.

(3) Starting with each node parameter of the well before gas kick at initial values and using the bottom of the well at any time as the starting point (the first node as the boundary value)，the Roe scheme based on the finite volume method can be applied to get the 2nd，3rd，and nth recursively.

Fig.2.Flow chart to solve the gas-liquid two-phase transient drift flow model.

(4) The parameter values of each space node in the time layer are calculated step by step to obtain the velocity field and pressure field distribution of the whole wellbore space.All parameters of the time layer node are then used as the known data of the next layer，and the cycle continues until the preset total simulation time cycle is reached.

3.Formation pressure inversion method based on LSTM

Long short-term memory network (LSTM) is an improved algorithm based on recurrent neural networks (RNN).LSTM has the advantages of RNN for the processing of time-series data and overcomes the limitations of long-term dependence and easy gradient disappearance.LSTM retains the chain form of RNN and is composed of a series of recursively connected subnetworks of memory blocks.The LSTM fundamental structure consists of three gates in the interaction layer，i.e.，input gate，output gate，and forget gate (Greff et al.，2016).The cooperation of the three gates can effectively control the information coming into the model.The network and structure of LSTM are shown in Fig.3.

Fig.1.The grid division of the wellbore.

The time series of wellhead mud pit gain of the input layer is x=（x1，x2，…，xt），the state of the hidden layer memory unit h=（h1，h2，…，ht）and output layer sequence y=（y1，y2，…，yt）can be calculated as follows:

The memory unit memorizes the historical information of the sequence data together with the hidden state.The information in the memory unit is controlled by three gating units.The forget gate is:

The input gate adds new information to the memory unit according to ht-1and xt:

The output gate determines htaccording to ht-1，xt，and ct:

where xtis the input at t time;ytis the output at t time;htis the hidden state at t time;Wxh，Whh，and Whyare the weights of input，hidden，and output，respectively;bhand byare the offset of hidden state and output，respectively;f(·) and g(·) are the activation functions of hidden layer and output layer，respectively;ht-1is the output of the last single state;Wxf，Whf，and Wcfare the corresponding weights of forgetting gate，respectively;bfis the offset of forget gate;σ is the sigmoid function;Wxi，Whi，and Wciare the corresponding weights of input gate;biis the offset of input gate;tanh is the hyperbolic tangent activation function;Wxcand Whcare the corresponding weights of memory units;bcis the bias term of memory unit weights;Wxo，Who，and Wcoare the corresponding weights of output gate;bois the bias of output gate;ctis the state value of cell structure at t time.

Fig.3.Structure of LSTM memory unit.

Fig.4.Photos of full-scale experimental well.

By changing the bottom hole parameters，including formation pressure and gas-liquid physical parameters，the model can be trained in real-time to meet the parameters of specific wells.The established gas-liquid two-phase transient calculation model can be used to obtain the time series data set of mud pit gain change to explore gas invasion and overflow development law under different formation parameters.These data are input into the LSTM neural network to learn the sequential operation law of mud pit gain fluctuation caused by a gas kick.LSTM training uses the mean relative error (MRE) and mean square error (MSE) as evaluation indexes.With a large amount of data in the time-series data set，the unique gating mechanism of LSTM and the use of a memory unit enable learning all historical information up to the current time and then allow modeling of this long-distance dependency when dealing with such a long input sequence.When training LSTM，the forward propagation is calculated using Eqs.(28)-(32)，and the coupling result of the current memory unit is calculated using Eq.(33)，combining the forget and input units.According to the error equation for gradient backpropagation training，the Adam algorithm can then be used as an adaptive momentum estimation algorithm for optimization training.

4.Case calculation

We tested this approach using a full-scale experimental well in the well control training center of the Daqing drilling engineering company (Fig.4).The gas-liquid two-phase transient drift flow model described above was used to predict the gas-liquid twophase development law in the wellbore after gas kick and thencalculate the gas holdup in the annulus at different times and the distribution of gas-phase velocity with well depth.

The test well has complete wellhead equipment(including bell mouth and BOP group)，experimental annulus and gas storage annulus，a 35 MPa full set of well control equipment，an energy storage air compressor，a simulated gas injection device，a measurement and control system，and a monitoring system.The specific parameters are as follows:vertical depth of 2000 m，continuous gas injection simulation realized at 2000 m borehole depth，hole diameter of 215.9 mm，and drill pipe outer diameter of 127 mm.The BHA is the experimental wellbore that consists of a bit，drill collar，drill pipe，variable thread joint，back pressure valve，annular pressure measuring nipple，and other components.Air is injected into the gas storage well through the air compressor until the pressure rises to the specified value.Drilling fluid enters the annulus through the surface manifold，drilling tools，and other equipment.The gas is injected into the annulus through a parasitic pipeline from the gas storage well.The parameters of the experimental well were used to generate the basic parameters of the simulation model，with consistent well depth，borehole diameter，and drill pipe diameter.The simulation model included 1000 grids，and the gas kick duration was set to 100 s(100-200 s after the start of the experiment).The total simulation time is 1000 s and the time step is 0.005 s.Other experimental parameters are listed in Table 1.

Table 1Well experimental parameters.

Table 2Evaluation indexes of LSTM model prediction results.

The transient simulation analysis of gas-liquid two-phase flow was carried out using the gas-liquid two-phase transient drift flow model described in this work.And the specific results are detailed below.Fig.5 shows the variation of bottom hole pressure with time within 110 s before the experiment (10 s after gas kick) and the distribution of cross-sectional gas holdup and gas velocity with well depth at t=110 s.When the bottom hole gas first intrudes into the wellbore，the velocity of the liquid column in the wellbore increases because the bottom hole gas increases the annular circulation friction resistance.At the same time，because the invaded volume occupies part of the liquid column volume，the hydrostatic fluid column pressure decreases.However，the gas does not expand significantly due to the small air intake.The hydrostatic fluid column pressure is less，and the increase of annular circulating friction is more significant than the loss of hydrostatic fluid column pressure，so the bottom hole pressure increases.There is some free gas at the bottom of the well，and the cross-section gas holdup is no longer zero.

Fig.6 shows the change of bottom hole pressure within the first 210 s of the experiment (10 s after the end of gas kick) and the distribution of cross-sectional gas holdup and gas velocity with well depth at t=210 s.At t=210 s，all the gas has been completely injected into the wellbore.With increased air intake，the bottom hole pressure decreases.The larger the bottom hole pressure，the smaller the volume occupied by gas and the smaller the pressure loss of the hydrostatic column，resulting in a minor reduction of the bottom hole pressure.With the constant gas intrusion rate，the pressure fluctuates widely at critical points during the initial and completion stages of the intake，but quickly returns to normal.

Fig.7 shows the variation of bottom hole pressure within the first 680 s of the experiment and the distribution of gas holdup and gas velocity with well depth at t=680 s.When the gas column migrates to the middle of the wellbore，as the gas column moves upward along the wellbore，the wellbore pressure becomes smaller，and the gas expands continuously.The maximum cross-sectional gas holdup increases from 0.2 at the bottom of the well to 0.3.The gas velocity increases with the decrease of the well depth，but the velocity change is not apparent.The gas expansion increases the velocity of the liquid column which increases the circulating friction pressure drop，but the gas expansion also reduces the hydrostatic pressure.The hydrostatic fluid column pressure loss is dominant，so the bottom hole pressure continues to decrease，but the decrease rate is relatively slow.

Fig.5.Development law of gas kick during the first 110 s (10 s after gas kick).

Fig.6.Development law of gas kick during the first 210 s (10 s after the end of gas kick).

Fig.8 shows the variation of bottom hole pressure within the first 980 s before the experiment and the distribution of gas holdup and gas velocity with well depth at t=980 s.When the gas column moves close to the wellhead，the gas in the annulus expands rapidly.The closer the wellhead is，the higher the gas holdup.The wellhead is almost full of gas，and the upper liquid column is nearly emptied，resulting in a sharp decrease in hydrostatic fluid column pressure.Although the gas velocity increases sharply，the bottom hole pressure decreases sharply because the gas density is far less than the drilling fluid density.The increase in annular friction loss is negligible compared with the hydrostatic fluid column pressure.The variation of bottom hole pressure and casing pressure as measured are shown in Fig.8，where the dotted lines show the measured pressure fluctuation curves.We found that the gas-liquid two-phase transient calculation model can better predict the development law of bottom hole pressure and casing pressure compared with the model，with average relative errors of prediction of 11.58% and 8.41%，respectively.The analysis error may be caused by the difference between the discrimination formula of the gas-liquid two-phase flow pattern and the actual situation as well as by inconsistency between the drill collar and drill pipe.

The gas-liquid two-phase transient calculation model can effectively simulate the fluid changes in the wellbore，allowing prediction of the gas-liquid two-phase transient flow development and the gas holdup and gas velocity distribution law at different well depths.However，it is difficult to predict the formation pressure when blowout is out of control due to a severe gas kick.The use of the LSTM method with deep learning can determine the nonlinear mapping relationship between the mud pit gain and formation pressure.Studying the time series characteristics and evolution mechanism of mud pit gain fluctuation caused by the gas kick and using measured early wellhead mud pit gain fluctuation data allows determination of inversion and prediction of formation pressure.

Fig.7.Development law of gas kick overflow during the first 680 s.

Fig.8.Development law of gas kick overflow in 980 s (comparison between experimental (dotted lines) and predicted (solid lines) values).

Fig.9.Comparison of measured and LSTM predicted values of time series evolution of mud pit gain under four different formation pressures:(a) 20.2 MPa;(b) 20.5 MPa;(c)21.0 MPa;(d) 21.4 MPa.

The fluctuation time-series data set of mud pit gain is calculated by adjusting formation parameters according to the established gas-liquid two-phase transient drift flow model to meet the specific LSTM model training needs.Here，the sampling interval was 0.5 s，and 1045 training sets of time series evolution samples were obtained.These data were divided into the training set and test set at a ratio of 8:2，and then the data training set were further divided into the training set and verification set at a ratio of 8:2 for crosstraining.The data set was then normalized and preprocessed to enhance robustness and generalization.The LSTM algorithm establishes the inversion prediction model by combining the entire connection layer and nonlinear activation function for intelligent training.This was done as follows.Using double-layer LSTM and the entire connection layer in series，the neuron parameters of the hidden layer were 256，512，and 1000，the learning rate range was 0.005-0.00001，the batch size was 128，and the number of epochs was 4000.As the number of iterations increased，the error of training samples decreased.After 2750 iterations，the model converged，and the final model was obtained.The prediction accuracy of the model was 87.3%.The time series evolution data sets of four groups of mud pit gain data under different formation pressures were simulated using this model and also obtained through full-scale experiments，as shown in Fig.9.The inversion values of different formation pressures and the model's prediction indexes are listed in Table 2.

It can be seen from Fig.9 that the prediction results of the LSTM model at the local position deviate significantly from the measured values.For example，there is significant prediction error of gas invasion overflow in the initial stage.Still，the overall prediction results are similar，which can better reduce the fluctuation of mud pit gain caused by the gas kick.As shown in Table 2，the inversion accuracy of formation pressure is high and the mean relative error in the testing process is about 11.5%，indicating that the LSTM model can be used for effective inversion and prediction of formation pressure.

The inversion methods of formation pressure when blowout is out of control can be summarized as follows:(1)the gas-liquid twophase transient calculation model based on finite volume method is first used to calculate the time-series variation characteristics of mud pit gain under field conditions;(2) the training set is established and input into the algorithm to train the LSTM neural network;(3) the time series characteristics of mud pit gain measured in the early stage of gas invasion are substituted into the training set to invert the formation pressure.

5.Conclusions

In this paper，we established a drift flow model of gas-liquid two-phase flow in the wellbore.In this approach，the governing equations can be solved using the finite volume method and the Roe scheme.The corresponding calculation conditions are then used to obtain the transient development law of gas-liquid twophase flow in the wellbore after gas kick.

When gas intrudes into the wellbore，the annular circulation friction resistance increases，and the bottom hole pressure increases.When the gas column moves to the middle of the wellbore，the gas expands continuously.The circulation friction resistance increases，and the hydrostatic fluid column pressure decreases，so the bottom hole pressure decreases continuously.When the gas column moves close to the wellhead，the gas expands rapidly and the hydrostatic fluid column pressure decreases sharply，so the bottom hole pressure decreases sharply.

Based on the deep learning algorithm framework，the timeseries variation characteristics of mud pit gain were obtained by the adjusting bottom hole parameters according to the transient calculation model.The corresponding training set was then established and input into the algorithm to train the LSTM neural network，extract the time relationship between formation pressure and mud pit gain，describe the evolution characteristics of this relationship，and establish an accurate inversion prediction model.The measured time series evolution data of mud pit gain were entered into the trained LSTM network to realize formation pressure inversion when the blowout is out of control，and the well cannot be shut-in.

Experimental data were collected from the well control training center of a Daqing Drilling Engineering Company and used to confirm the reliability of prediction using the gas-liquid two-phase transient drift flow model.The formation pressure inversion and prediction results based on the LSTM algorithm exhibited high accuracy，providing theoretical support for the early monitoring of gas kick and overflow and prediction of formation pressure after uncontrolled blowout.

Acknowledgements

This work was financially supported by the National Natural Science Foundation of China (Grant No.51974090，51474073).