Risk Assessment Method for Bullheading Killing Based on the Uncertainty of Formation Parameters

To analyze the risks caused by the uncertainty of formation parameters to bullheading killing, a method for quantitatively evaluating the bullheading killing risks is established. Firstly, considering the influence of gas invasion volume, formation fracture, and killing parameters, a bullheading killing model is established based on a gas-liquid two-phase flow. Then, the uncertainties of formation parameters (formation pressure and permeability) are quantified. Based on the shut-in wellhead information, the range of formation pressure is predicted with the gas column model and multiphase flow model. Considering the influence of formation fracture on the permeability, Monte Carlo random sampling is applied to predict the range of formation permeability. Based on industry standards, a safety pressure value is set up, and the wellbore pressure corresponding to all value combinations of formation parameters under the given killing parameters is obtained by killing model. Moreover, according to the probability and degree that wellbore pressure exceeds the safety value, the risks are rated to quantify the risk of bullheading killing. Under this circumstance, the feasibility and accuracy of this method are validated by practical cases, and it is found by simulation that flow rate can affect the risk of wellhead damage to the greatest extent, and there exists a critical rate. When the flow rate is greater than the critical rate, the increase in flow rate will greatly improve the risk probability. In such case, improving the density of kill fluid can reduce the risk of wellhead damage in a limited way, but it will greatly increase the risks of formation fracture and casing damage. Therefore, for bullheading killing, it is not advisable to employ high-density kill fluid. By this method, the bullheading killing risks can be fully assessed before actual construction, thus providing reference for determining reasonable construction parameters of bullheading killing.


Introduction
During the bullheading killing, the wellbore and formation will bear great pressure. If the pressure exceeds the safety value, more serious accidents such as wellbore damage and formation fracture will happen [1,2]. For example, in Dina Well 2-9 [3], the formation was fractured in the process of bullheading killing, which caused serious underground blowout. In Da Well 61-21 [4], when the kick was treated by bullheading, the flow rate was too large, resulting in formation fracture and failure of killing. In addition, it is difficult to acquire accurate formation parameters before site operation, and the uncertainty of formation parameters results in extra risks to the bullheading. Therefore, the high construction risks have imposed restrictions on the use of bullheading on site.
Scholars have conducted plenty of researches on drilling risks. Based on analytical hierarchy process (AHP) and fuzzy comprehensive evaluation model, Li [5] established a risk assessment model for wellbore integrity of natural gas hydrate according to the established risk acceptance criteria and risk assessment index system, determined the overall risk level of wellbore, and assessed the risks of wellbore integrity of natural gas hydrate. Based on FMEA (Failure Mode Effects Analysis) and with the help of fuzzy mathematics and analytic hierarchy process, Yan [6] established a fuzzy evaluation matrix aiming at the composition and equipment of well control system of jack-up drilling platform by fuzzy comprehensive evaluation method, confirmed the weighting coefficient of each factor according to analytic hierarchy process, and made a comprehensive evaluation of the risks of the whole well control system. Hu [7] analyzed the advantages of grey relational evaluation method in drilling well control risk assessment. After corresponding improvement by reference to the principle of TOPSIS (Technique for Order Preference by Similarity to an Ideal Solution), the method is more suitable for quantitative evaluation of single well control risks and drilling well control risks based on formation factors. Chang et al. [8] established the risk index evaluation system for drilling riser accidents by analysis on drilling riser accidents and identification of risk factors and determined the weight of each risk factor by AHP-EM method. In addition, he also defined a 9-tuple set for the modeling of drilling riser risks according to FPN theory and evaluated the risks of deepwater drilling riser based on fuzzy Petri net model. Liang et al. [9] proposed a fuzzy multilevel algorithm based on particle swarm optimization (PSO) to optimize the support vector regression (SVR), established a fuzzy multilevel drilling leakage risk assessment system, and studied the risk assessment results by PSO-SVR algorithm, thus realizing real-time risk assessment. By analyzing historical data and identifying unsafe factors with the help of expert experience, Wang et al. [10] divided 32 factors affecting the safety of drilling operation site into "unsafe behavior of human" and "unsafe state of objects" and conducted quantitative risk assessment on drilling operation site based on Bayesian network. Wu et al. [11] put forward a conditional probability method for determining dangerous events and their consequences based on deep belief network to evaluate offshore drilling risks. This method, based on the dynamic Bayesian network (DBN) theory, represents the uncertainty of the model and parameters by adding extra nodes. Tamim et al. [12] developed a variety of decision support algorithms by virtue of Bayesian network and real-time indexes to predict the drilling kick risk.
It is found through survey that analytic hierarchy process, fuzzy comprehensive evaluation method, and neural network method are generally employed for quantitative evaluation of drilling and well control risks, but there is still no method for quantitative evaluation of bullheading killing risks. Analytic hierarchy process is complex in structure and mainly relies on subjective judgment of experts. Therefore, it is more suitable for the calculation of index weight. For fuzzy comprehensive evaluation method, the mathematical model is complex, and the weight and membership function are generally determined subjectively, thus needing a large number of field statistical data. The neural network method is of low accuracy and needs a large number of training samples, with a large amount of calculation. Under this circumstance, a method for quantitatively evaluating the bullheading killing risks is studied in this paper for the first time. After the uncertainty of formation information is represented by value intervals, a killing model is established to calculate the wellbore pressure corresponding to different value combinations of formation parameters under the given killing parameters, and then, the killing risks are quantitatively evaluated according to the probability and degree that wellbore pressure exceeds the safety threshold.
In Section 2 of this paper, a bullheading killing model is established and a method for quantitatively evaluating the bullheading killing risks is elaborated. In Section 3, the bullheading killing model is validated using field data. And the risk assessment method proposed in this paper is validated by practical cases. Factors affecting the killing risks are also analyzed. At last, a conclusion is drawn.

Material and Methods
2.1. Establishment of Bullheading Killing Model. The bullheading killing technology is a well control method [13,14]. After the wellhead is closed, the kill fluid is pumped through annulus or drill pipe to push the kick back to formation by high pressure at the wellhead. Annulus reverse bullheading is the main operation mode, and the change of wellbore pressure during forward bullheading is similar to that during reverse bullheading. Therefore, reverse bullheading is taken as an example in this paper for analysis.
Assumptions in the killing model established in this paper are made as follows: the well is vertical. After kick shut-in, there is a segment of pure gas column gathering at the bottom hole. Due to the difference in gas-liquid density, part of the gas slips upward in drilling fluid (water base) in the form of bubbles to form a segment of gas-liquid twophase flow, and there is a segment of uncontaminated drilling fluid at upper end of the gas-liquid two-phase flow ( Figure 1); the deformation of pipe column is not considered; the fluid is incompressible; the geothermal gradient is a constant or the geothermal gradient function is known; the formation in open hole section is homogeneous, and the fluid in the formation flows radially; and the gas holdup in gas-liquid two-phase flow is linearly distributed.
After shut-in, kill fluid will be pumped into the wellbore through kill manifold and flow downward to push the fluid in wellbore into the formation in turn. After the gas-liquid two-phase flow is completely pushed back into the formation, the killing process is finished (whether the drilling fluid needs to be pushed back into the formation should be determined as the case may be). Only when the bottom hole pressure is greater than formation pressure, the kick will flow to the formation. If the relationship between the recovery of formation pressure near the bottom hole and the build-up speed of hydrostatic fluid column pressure is taken into account, the bullheading killing can be divided into two processes, i.e., dynamic sealing stage and pushback stage.
2.1.1. Dynamic Sealing Stage. Dynamic sealing stage is as follows: if the formation pressure recovery speed at the bottom hole is higher than the build-up speed of hydrostatic fluid column pressure, the fluid in wellbore cannot be pushed back to the formation, the gas in wellbore will be compressed, and the wellhead and bottom hole pressure will be boosted until the bottom hole pressure is greater than formation pressure. However, if the build-up speed of hydrostatic fluid column 2 Lithosphere pressure is higher than the formation pressure recovery speed or there is no formation pressure recovery process, there will be no dynamic sealing stage. In dynamic sealing process, the kill fluid will push the fluid in wellbore to flow downwards, which must overcome the flow friction while flowing down. If the gas invading into the wellbore from the formation at this stage is ignored, the wellhead casing pressure can be acquired by the following formula: In the formula, P a is the wellhead pressure, MPa; P pt is the real-time formation pressure, MPa; P f a is the annulus flow friction, MPa; and P ha is the annulus hydrostatic fluid column pressure, MPa.
After the killing is started, bottom hole pressure will increase with the increase of fluid column pressure and can be obtained according to the following formula: In the formula, P w is the bottom hole pressure (MPa) and P as is the initial wellhead pressure (MPa). According to the formation pressure recovery characteristics, the real-time formation pressure P pt can be worked out through the following formula: In the formula, Pp is the formation pressure, MPa; P ws is the initial bottom hole pressure, MPa; a is the formation pressure recovery coefficient; t is the time, s; and a can be decided from nearby wells.
The maximum internal pressure that the casing can bear can be acquired according to internal pressure at casing shoe.
In the formula, P shoe is the pressure at casing shoe, MPa; P f shoe is the friction at casing shoe, MPa; and P hshoe is the hydrostatic fluid column pressure at casing shoe, MPa.
Annulus friction can be expressed by the Fanning-Darcy Formula [15]: In Formula (5), P f i is the fluid flow friction, MPa; f i is the Fanning friction factor of a certain type of fluid,  3 Lithosphere dimensionless; ρ i is the density of a certain type of fluid, kg/ m 3 ; v-flow velocity, m/s; D w is the inner casing diameter, mm; D z is the outer diameter of drill string, mm; and L i is the length of a certain type of fluid in wellbore, m.
The Fanning friction factor of each type of fluid can be worked out by Formula (6): In Formula (6), R ei is the Reynolds number of fluid, dimensionless; f g is the gas phase Fanning-Darcy friction factor; f l is the liquid phase Fanning-Darcy friction factor; f m is the Fanning-Darcy friction factor of gas liquid twophase flow, dimensionless; Re g is the gas phase Reynolds number, Re l is the liquid phase Reynolds number, Re m is the Reynolds number of gas liquid two-phase flow, dimensionless; and ε/D is the relative roughness factor, dimensionless.

Pushback Stage.
Pushback stage is as follows: when the bottom hole pressure is greater than formation pressure, the fluid in wellbore will be driven by kill fluid to flow to the formation in turn. This process is called pushback stage. After the fluid in wellbore is pushed back to the formation in turn, and gas-liquid two-phase flow is completely pushed into the formation, the killing process is basically finished, and whether the drilling fluid needs to be pushed back to the formation depends on the actual situation on site.
In this process, the fluid flow needs to overcome the filtrational resistance, and the change of wellhead casing pressure can be obtained according to Formula (7): where P si is the filtrational resistance of fluid, MPa. Bottom hole pressure: Internal pressure of the casing can still be worked out by Formula (4).
Based on the radial flow productivity formula [16], the d formula of oil gas entering the formation is as follows: In the formula, K φ is the solubility, μm 2 ; μ i is the fluid viscosity, mPa·s; h r is the reservoir thickness, m; r e is the radius of reservoir outer edge, m; r w is the borehole radius, m; and r e and outer boundary pressure can be assumed or calculated from the information of reservoir engineering.
The bubble slip velocity in killing process can be calculated by the Harmathy [17] Formula: In the formula, u ∞ is the bubble slip velocity, m/s; g is the gravitational acceleration, m/s 2 ; ρ l is the liquid phase density; and ρ g is the gas phase density, g/cm 3 .

Risk Assessment
Method considering the Uncertainty of Formation Parameters. If accurate formation parameters can be grasped, after killing parameters are worked out, it is only required to simulate the killing process by virtue of the physical and mathematical model established in the previous section to judge whether any accident will happen in killing process. However, it is difficult to accurately master the formation parameters in practical construction. Therefore, it is necessary to take into account the risks caused by the uncertainty of formation parameters.
First of all, ranges of formation parameters are calculated, and safety values of wellhead pressure, casing pressure, and bottom hole pressure are worked out according to industry standards. Then, different combinations of formation parameters are selected in value intervals to calculate the wellbore pressure in killing process. Finally, killing risks are evaluated according to the probability and degree that wellbore pressure exceeds the safety value.

Ranges of Uncertain Formation Parameters.
There are various formation factors that affect bullheading killing, and it will take a long time to carry out a full-scale simulation. Therefore, in this paper, two key parameters (formation pressure and formation permeability) that can impose the greatest influence on killing risks are selected to analyze the bullheading killing risks.
(1) Formation Pressure. Existing researches mostly assume that after kick shut-in, the gas accumulates at the bottom hole in the form of gas column, and formation pressure is predicted according to gas column model. In this paper, according to the literature [18][19][20][21][22][23][24][25][26], gas column model (assuming that the invading gas only exists at the bottom hole in the form of gas column) and multiphase flow model (assuming that the gas in wellbore is completely dispersed in the drilling fluid) are employed to predict the formation pressure based on 10 sets of data in the literature. The prediction results are presented in Figure 2.
Compared with the actual formation pressure, it is found that the gas column model is effective in predicting small kick and low wellhead pressure but ineffective for great kick 4 Lithosphere or excessively high wellhead pressure, with an error of -22.7%~18% between the predicted formation pressure and actual formation pressure. When multiphase flow model is adopted to predict formation pressure, the error is -20.1%~19.87%. Therefore, in this paper, both gas column model and multiphase flow model are employed to predict formation pressure at first, and value intervals are generated according to the corresponding error ranges. Then, a new value interval is generated based on the union of two intervals. Finally, the value range of the new interval is extended by 30% as the final value interval according to site operation experience, that is, P p ∈ min 0:7P pg min , 1:3P pm max À Á , P p min = min P pg min , P pm min À Á : (2) Formation permeability. For fractured reservoirs, the formation permeability is determined by formation matrix permeability and formation fracture [27]. According to the parallel plate theory, the natural fracture development area is treated by equivalent filtrational resistance method [28], and the influence of formation fracture width, length, and number on killing risks is expressed by permeability, which can be calculated by Formula (12).
In Formulas (12) and (13), K is the reservoir matrix permeability, μm 2 ; n l is the number of natural fractures, dimen-sionless; w f is the natural fracture width, mm; k f is the natural fracture permeability, mD; and l is the natural fracture length, m. Since a number of parameters are involved in the formula, relevant parameters can be randomly selected by Monte Carlo method [29][30][31][32] for simulation and calculation to acquire the value interval of formation permeability.

Risk
Rating. Risk level is determined by the probability and severity of the risk [33][34][35][36], and risk level evaluation value can be calculated by Formula (14).
In the formula, P i is the risk probability and j is the risk severity coefficient.
The risks are divided into 0-V levels, with an upward trend level by level. Risks at level II and below are acceptable, risks at level III should be paid close attention, and risks at level IV and above are unacceptable. Please refer to Table 1 for the specific scoring criteria.
For working conditions with high risk level, measures should be taken to reduce the risks before bullheading killing or killing by other means.
(1) Calculation of risk probability. Formation parameters are combined in value intervals to simulate the kick process, and kick volume is taken as the end condition and actual shut-in wellhead pressure is compared to judge whether the calculated wellhead pressure in simulation is reasonable (in this paper, it is reasonable as long as the error with the actual wellhead pressure is smaller than 30%). After noise suppression for formation parameter combinations, reasonable parameter combinations are brought into the killing model to work out the wellbore pressure. Then, the calculation 5 Lithosphere result is compared with the safety value. If the wellbore pressure is higher than the safety value, it should be recorded, and the range beyond the safety value should be calculated. This process should be repeated until all combinations of formation parameters in value intervals are taken. Finally, the risk probability P i can be calculated by Formula (15).
In the formula, N 1 refers to the number of times that wellbore pressure exceeds the safety value, and N E represents the number of effective simulations.
(2) Calculation of risk severity. Firstly, the coefficient α i is defined according to the degree that wellbore pressure exceeds the safety value, which represents the influence of the degree that wellbore pressure exceeds the safety value on killing risks under a single combination of formation parameters. The degree E of exceeding the risk assessment value can be worked out by Formula (16) [37,38].
According to the degree that wellbore pressure exceeds the safety value, after consulting the experts' opinions, different values are given to α i . Please refer to Table 2 for the specific values.
Risk severity coefficient j is determined by cumulative αi, as shown in Table 3.
(3) Calculation of safety value. According to industry requirements [39]: The maximum wellhead pressure shall not exceed 0.8 times of the grade of wellhead blowout preventer, In the formula, E AS is the wellbore risk assessment value and P F is the grade of blowout preventer (MPa).
In the killing process, the maximum casing pressure shall not be higher than the maximum internal pressure strength, In the formula, E TS is the wellbore risk assessment value and P T max is the maximum internal pressure strength of casing (MPa).
In the killing process, bottom hole pressure shall be lower than the formation fracture pressure, In the formula, E WS is the risk assessment value of the formation at bottom hole and P ps is the formation fracture pressure (MPa).

Results and Discussion
3.1. Validation of the Bullheading Killing Model. The bullhead killing model is validated with field data. When Well YM4 was drilled to 7279 m in the third spud-in, it encountered abnormal high formation pressure and kick occurred. Under this circumstance, the well was shut in, the drilling tools were in good condition and at the bottom hole. Based on a large number of neighboring well data, the formation pressure coefficient was determined to be 1.56, and there might be high sulfur content. Therefore, bullheading was selected to kill the well. Please refer to Table 4 for other drilling or formation parameters.
The kill fluid with a density of 1.8 g/cm 3 was employed on site for killing operation with the flow rate of 20 L/s, and the same parameters were adopted in this paper for killing simulation. The changes of wellhead pressure in the killing process are as shown in Figure 3. It can be seen from Figure 3 that the actual maximum wellhead pressure is 21.5 MPa, and the maximum wellhead pressure calculated in simulation is 20.2 MPa, which is slightly lower than the actual wellhead pressure with an error of 7% and an average error of 3.5%. The reason for such error is that the height of gas-liquid two-phase flow and gas holdup are different from actual parameters and the filtrational resistance is a theoretical value. On the whole, the simulation kill curve is in good   6 Lithosphere agreement with the true kill curve, which validates the accuracy of this model.

Analysis on the Bottom Hole Pressure and Casing Shoe
Pressure. Assuming that a certain well suffers from kick during drilling, with a kick volume of 2 m 3 and invading gas of methane, the well is shut in immediately upon the kick is found. Drilling tools are at the bottom hole and in good condition. Please refer to Table 5 for other drilling and formation data.
The changes of wellhead pressure, bottom hole pressure, and formation pressure at casing shoe in killing process are acquired through killing simulation, as shown in Figure 4. As shown in Figure 4(a), the wellhead bears great pressure at the initial stage of killing process. Besides, the wellhead pressure first decreases and then increases at the stage when multiphase flow is pushed back to the formation, and the second peak point appears. Moreover, it can be seen from Figure 4(b) that the bottom hole pressure and casing shoe pressure increase gradually along with the killing, and they will bear great pressure in later killing period.

Basic Data and Risk Assessment. Well
LG34 is an exploration well [40]. Well type: vertical well; design well depth: 6,950 m; target stratum: the 1st member of Ordovician Yingshan Formation, as well as Yijianfang Formation and Lianglitage Formation; drilling finish stratum position: Ordovician Yingshan Formation. The wellbore structure of Well LG34 is as shown in Figure 5.
Kick happened when Well LG34 was drilled to the depth of 6694.00 m with Φ152.4 mm drill bit. At this time, the shut-in standpipe pressure was 4 MPa, casing pressure was 0, and the drilling assembly used was composed of one Φ152.4 mm PDC bit, one Φ330 mm NC 350 joint, 18 Φ120.7 mm drill collars, 291 Φ114.3 mm heavy weight drill pipes, and 365 Φ127 mm drill pipes. The wellhead casing pressure increased to 41 MPa after choke circulating killing failed for three times. The initial conditions of Well LG34 before bullheading killing are presented in Table 6.
The grade of wellhead blowout preventer for Well LG34 is 70 MPa, the formation fracture pressure coefficient is 1.3, and the maximum internal pressure of casing is 78 MPa. Please refer to Table 7 for value intervals of formation parameters.
When the wellhead pressure is 41 MPa, the kill fluid with a density of 1.18~1.22 g/cm 3 is employed on site for killing at the flow rate of 3.9~10 L/s. In this paper, the kill fluids with a density of 1.18 g/cm 3 and 1.22 g/cm 3 are used, respectively, to simulate the killing at the flow rate of 10 L/s for risk assessment.   It is found through testing that accuracy requirements can be met when the step size of formation pressure value is 0.1 MPa and the step size of formation permeability value is 1.0 mD. Based on the above parameters, the maximum wellbore pressure under a single set of killing parameters is as shown in Figure 6, in which Figure 6(a) presents the maximum wellhead pressure distribution, Figure 6(b) shows the maximum bottom hole pressure distribution, and Figure 6(c) indicates the distribution of maximum internal pressure of casing at the casing shoe. As shown in Figure 6(a), when the formation pressure is 74 MPa and the permeability is 50 mD, the minimum wellhead casing pressure calculated is 49.6 MPa. When the formation pressure is 83.5 MPa and the permeability is 5 mD, the maxi-mum wellhead casing pressure worked out is 63.2 MPa. In this way, the minimum wellhead casing pressure differs from the maximum wellhead casing pressure by 13.6 MPa. In addition, it can be seen from Figure 6(b) and Figure 6(c) that since the open hole section of Well LG34 is short, the internal pressure of casing at the casing shoe and the pressure of bottom hole formation are not quite different. When the formation pressure is 74 MPa and the permeability is 50 mD, the minimum bottom hole pressure and minimum internal pressure of casing calculated are 78.9 MPa and 77.9 MPa, respectively. However, when the formation pressure is 83.5 MPa and the permeability is 5 mD, the maximum bottom hole pressure and minimum internal pressure of casing are calculated to be 84.74 MPa and 83.75 MPa, respectively. When bullheading is adopted for killing with kill fluid of different densities, the risk assessment results are as shown in Table 8. According to Table 8, when kill fluids with a density of 1.18 g/cm 3 and 1.22 g/cm 3 are, respectively, used for bullheading killing at the flow rate of 10 L/s, the risks of both wellhead damage and formation fracture belong to level I, and the risk of wellbore damage belongs to level 0, which is benefited from the high grade of wellhead blowout preventer and the good pressure bearing capacity of casing and formation.

Parametric Analysis
3.4.1. Killing Flow Rate. When initial wellhead pressure is 41 MPa, the kill fluid with a density of 1.18 g/cm 3 is employed for killing at different flow rates, and the probability of various risks is as shown in Figure 7.
It can be seen from Figure 7 that with the increase of flow rate, the probability of various risks is also improved. The increase of flow rate imposes great impact on the risk of wellhead damage. When the flow rate increases from 10 L/s to 12 L/s, the risk probability of wellhead damage increases from 1% to 13%, with small increasing range. However, upon the flow rate exceeds 14 L/s, the risk probability of wellhead damage is improved from 14% to 54%, which means that for the increase of flow rate by every 1 L/ s, the risk probability is improved by 7%, with great increasing range and rate. Thus, it can be seen that there exists a critical rate for the risk of wellhead damage. Before the flow rate reaches the critical rate, the increase of flow rate will not impose great impact on the risk of wellhead damage. However, upon the flow rate exceeds the critical rate, once the flow rate increases, the risk probability will be improved to a great extent. However, the increase of flow rate has little effect on the risk of casing damage. When the flow rate is smaller than 12 L/s, there is no risk of casing damage. When the flow rate increases to 20 L/s, the risk probability is improved to 17%, with little change in both the increasing range and rate. Furthermore, the increase of flow rate also imposes general influence on the risk of formation fracture.  Shut-in standpipe pressure (MPa) Drilling fluid viscosity (mPa·s) 0 3 0 Drilling fluid density (g/cm 3 ) Kick volume (m 3 ) 1. 17 6 Drilling fluid density at the outlet (g/cm 3 ) Grade of blowout preventer 1.09 70 When the flow rate is smaller than 15 L/s, there is no risk of formation fracture. When the flow rate is increased to 20 L/s, the risk probability is improved to 7%, with little change in both increasing range and rate. Please refer to Table 9 for risk level assessment results of bullheading killing at different flow rates.
As shown in Table 9, the risk level of wellhead damage is increased from level 1 to level 2 after the flow rate exceeds 14 L/s and to level 3 when the flow rate reaches 20 L/s. Thus, it can be seen that after the flow rate exceeds 14 L/s, the risk probability is improved to a great extent, but the wellhead casing pressure is not much higher than the safety value, and the risk level coefficient j is not large. Therefore, even if the risk probability is as high as 54%, the final risk level still does not exceed level 3.

Kill Fluid Density.
When the initial wellhead pressure is 41 MPa, at the flow rate of 15 L/s, the kill fluids with a density of 1.18~1.26 g/cm 3 are employed, respectively, for killing at different flow rates, and relevant risks are assessed. The risk probability is as shown in Figure 8.
It can be seen from Figure 8 that the risk probability of casing damage will greatly increase with the increase of kill fluid density. When the density increases from 1.18 g/cm 3 to 1.2 g/cm 3 , the risk probability is improved from 6.3% to 15%. As the kill fluid density increases continuously,      10 Lithosphere the risk probability will increase at a higher rate and with a wider range. When the kill fluid density reaches 1.24 g/ cm 3 , the risk probability will be as high as 68%. In addition, the risk probability of formation fracture is also improved to a great extent, but the risk of wellhead damage basically remains unchanged, which is because the drilling fluid density has limited effect on reducing the maximum wellhead pressure. As Table 10 shows, when the kill fluid density reaches 1.22 g/cm 3 , the risk of casing damage has been up to level III, and when the density reaches 1.24 g/cm 3 , the risk level has increased to level V. Therefore, when designing the kill fluid density before killing, it is required to fully consider the influence of casing and formation and avoid using the kill fluid of a high density for killing.

Conclusions
Based on the gas-liquid two-phase flow theory, a wellbore pressure prediction model for bullheading killing is established in consideration with the influence of gas invasion volume, formation fracture and killing parameters. In addition, the uncertainty of formation parameters is represented by establishing value intervals, the formation pressure and   11 Lithosphere permeability are selected as key formation parameters, and their value intervals are predicted, respectively. Moreover, safety values are set based on industry standards, and wellbore pressure corresponding to all different value combinations of formation parameters under the given killing parameters is worked out by virtue of the killing model. Finally, risks are rated according to the probability and degree that wellbore pressure exceeds the safety value, and the bullheading killing risks are quantitatively evaluated. Based on the results, the following conclusions are drawn: (1) The bullheading killing model is established and validated by actual killing data of Well YM4. The average error between the wellhead pressure calculated by the model and the wellhead pressure measured in practical killing process is only 3.5%. The maximum wellhead pressure calculated by the model is 20.2 MPa, and the maximum wellhead pressure measured in practical killing is 21.5 MPa, with an error of smaller than 10%. Therefore, the killing model is accurate and reliable and can be used to predict wellbore pressure before site operation (2) In this paper, the risk theory is introduced into killing process for the first time, the method of quantitatively evaluating the bullheading killing risks is put forward, and the feasibility and accuracy of this method are validated according to the killing case of Well LG34. By this method, the risks of killing the Well LG34 by bullheading with the kill fluid with a density of 1.22 g/cm 3 at the flow rate of 10 L/s are evaluated. The evaluation results show that even if the wellhead pressure before killing is as high as 41 MPa, the risks during the bullheading killing are still extremely low, which breaks the inherent understanding that there are high risks in the bullheading killing under high wellhead pressure (3) In the process of bullheading killing, the risks of formation fracture, casing damage, and wellhead damage will all increase with the increase of flow rate.