Skip to main content

Dynamic Kill Simulation

This page contains a conceptual description of the dynamic kill engine implemented in Oliasoft WellDesignβ„’#


Dynamic kill is a process to kill a blowout by pumping kill fluid through a relief well into the blowout well to suppress the flow. The model herein describes an advanced finite volume method to accurately calculate the process of dynamically killing a wellbore, with associated effects, such as necessary pump pressure.

Printable Version#

Oliasoft Technical Documentation - Dynamic Kill


The first input is the complete geometry of the flowing well and relief well. The geometry includes well trajectory, drill string, casing details, open hole dimensions, and depth of the intersection between the relief well and flowing well, technically speaking, inclinations, depths, inner and outer diameters of each section in wellbores.

In the flowing well, the fluid properties of crude oil and/or natural gas are required both at the standard conditions and in the reservoir. Users can choose the blowout rates for each scenario.

The kill mud has the potential to dynamically kill the blowout successfully by restoring the bottom hole pressure above the reservoir pressure, and thus, formation properties are also considered as inputs. The last required inputs are fluid properties and flow rate of kill fluid

Drift Flux Model#

The drift-flux model is widely used to predict two and three-phase flow in pipes and wellbores. The drift-flux model can handle co-current and counter current flow, and it allows the heavy and light phases to move in opposite directions when the mixture velocity is small. The drift-flux model relies on the slip relation, which requires empirical parameters.

Governing Equations#

The conservation equation of drift-flux model contains mass conservation equations for each phase, and a mixture momentum equation for gas gg and liquid ll ,

βˆ‚mlβˆ‚t+βˆ‚mlulβˆ‚x=Slβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(1)\frac{\partial m_l}{\partial t} + \frac{\partial m_l u_l}{\partial x} = S_{l} \;\;\;\;\;\;\;\;\;(1)

βˆ‚mgβˆ‚t+βˆ‚mgugβˆ‚x=Sgβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(2)\frac{\partial m_g}{\partial t} + \frac{\partial m_g u_g}{\partial x} = S_{g} \;\;\;\;\;\;\;\;\;(2)

βˆ‚moilβˆ‚t+βˆ‚moiluoilβˆ‚x=Soilβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(3)\frac{\partial m_{oil}}{\partial t} + \frac{\partial m_{oil} u_{oil}}{\partial x} = S_{oil} \;\;\;\;(3)

βˆ‚βˆ‚t(mgug+mlul)+βˆ‚βˆ‚x(mgug2+mlul2)+βˆ‚Pβˆ‚x=Gl+Gg+Fl+Fgβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(4)\frac{\partial}{\partial t}( m_g u_g + m_l u_l) + \frac{\partial}{\partial x}(m_g u_g^2 + m_l u_l^2) + \frac{\partial P}{\partial x} = G_l + G_g + F_l + F_g \;\;\;\;\;\;\;\;(4)

Equation of State#

In the equation of state (EOS), the density is a function of temperature and pressure in the thermal flow. In the isothermal flow, the EOS is expressed as a linear function of pressure for each phase,

ρ=ρ0+Pβˆ’P0c2β€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(5)\rho = \rho_0 + \frac{P - P_0}{c^2} \;\;\;\;\;\;\;\;(5)

where cc is the speed of sound.

Pressure Equation#

In multiphase flow, the volume conservation should be required,

Ξ±g+Ξ±l=1β€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(6)\alpha_g + \alpha_l = 1 \;\;\;\;\;\;\;\;(6)

which ensures the accuracy of the numerical solutions.

The pressure equation can be derived from the volume balance combined with EOS, or we can introduce the mass conservation equations additionally to obtain a more conservative balance equation for pressure.

Slip Relation#

The slip between the fluid phases is described by two basic parameters, the profile parameter C0C_0 and drift velocity UdU_d , where UdU_d is the difference of velocity between phases. And those parameters vary with different flow regimes.

ug=C0um+Udβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(7)u_g = C_0 u_m + U_d \;\;\;\;\;\;\;\;(7)

The average mixture velocity is the sum of superficial velocities usgu_{sg} and uslu_{sl} ,

um=usg+usl=Ξ±gug+Ξ±lulβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(8)u_m = u_{sg} + u_{sl} = \alpha_g u_g + \alpha_l u_l \;\;\;\;\;\;\;\;(8)

The Shi slip relation is employed in the present work, and a more continuous correlation is considered in the calculation to strengthen the robustness of the present model.

Numerical Scheme#

A flux splitting method, called advection upstream splitting method (AUSMV) is implemented to solve the compressible Euler equations, and it has been well-developed in recent years.

Flux Splitting Method#

The (AUSMV) scheme is developed based on the upwind concept as a flux function for solving a system of conservation equations, and it has been verified for its capability to solve a wide range of problems, including multiphase flow. The flux using (AUSMV) method is written as FAUSMV\mathcal{F}^{AUSMV}.

Numerical Discretazion#

Rewrite the conservation equations 1 to 4 in the Euler form,

βˆ‚wβˆ‚t+βˆ‚Fβˆ‚x=0β€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(9)\frac{\partial w}{\partial t} + \frac{\partial \mathcal{F}}{\partial x} = 0 \;\;\;\;\;\;\;\;(9)

The wellbore of the flowing well is divided into NN grid cells, where index jj represents the cell center and j+1/2j+1/2 locates at the cell interface. The conservative variable ww at time step n+1n+1 is,

wjn+1=wjnβˆ’Ξ”tΞ”z(Fj+1/2AUSMVβˆ’Fjβˆ’1/2AUSMV)β€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(10)w^{n+1}_j = w^{n}_j - \frac{\Delta t}{\Delta z}(\mathcal{F}^{AUSMV}_{j+1/2} - \mathcal{F}^{AUSMV}_{j-1/2}) \;\;\;\;\;\;\;\;(10)

Boundary Conditions#

The inlet and outlet boundaries are presented in this section. The pressure fluxes at inlet use extrapolation method, and the flow rates of gas and liquid at the inlet are given as initial values. For most of the scenarios, an open end at the outlet is considered, and thus, a prescribed pressure should be used. Otherwise, the pressure at the outlet can employ the extrapolation method if it is a closed-end situation.

Pressure Drop in the Relief Well#

Kill mud is typically pumped into the relief well annulus through kill/choke lines. The friction loss in the relief well plays an important role to calculate pump pressures during the dynamic kill operation. A practical hydraulic model, called the Herschel-Bulkley model is considered in for this,

Ο„=Ο„0+KΞ³Λ™nβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(11)\tau = \tau_0 + K \dot{\gamma}^n \;\;\;\;\;\;\;\;(11)

which is a relation between shear rate Ξ³Λ™\dot{\gamma} and shear stress Ο„\tau. The Herschel-Bulkley model is a rheological mathematical model for the flow of non-Newtonian fluids which can be developed for general use. The wall shear stress Ο„w\tau_w for a given flow rate qgivenq_{given} can be obtained using the above model. And Reynolds number ReRe for Herschel-Bulkley fluids is,

Re=12ρV2Ο„wRe = \frac{12\rho V^2}{\tau_w}

The transition interval Re1<Re<Re2Re_1 < Re < Re_2 is estimated using proposed correlations. If the flow type is laminar flow, Re<Re1Re < Re_1 ,

flam=24Ref_{lam} = \frac{24}{Re}

If it is turbulent flow, Re>Re2Re > Re_2,

1ftur=4N0.75log10[Reftur(1βˆ’N/2)]βˆ’0.4N1.2\frac{1}{\sqrt{f_{tur}}} = \frac{4}{N^{0.75}}log_{10}\bigg[Re f_{tur}^{(1- N/2)}\bigg] - \frac{0.4}{N^{1.2}}

For transitional flow,

flam=flam+(Reβˆ’Re1)(fturβˆ’flam)Re2βˆ’Re1f_{lam} = f_{lam} + \frac{(Re - Re_1)(f_{tur} - f{lam})}{Re_2 - Re_1}

The pressure drop in the relief well is given as,

dPfdl=2fρV2d2βˆ’d1β€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Šβ€…β€Š(12)\frac{dP_f}{dl}= \frac{2f\rho V^2}{d_2-d_1} \;\;\;\;\;\;\;\;(12)

A model from Founargiotakis et al. (2008) is used as an accurate methodology to predict the friction factors ff for laminar, transitional, and turbulent flow as well as the transition points from laminar to turbulent flow for the flow of Herschel-Bulkley fluids in a concentric annulus. Additionally, we introduce the effect of eccentricity and diameter ratio to meet actual conditions according to Erge et al. (2014)'s work.


The dynamic kill simulation is combining a multiphase flow model and a hydraulic model to estimate the kill duration, mud consumption, pump pressure, pump power, and other fluid quantities if it is necessary. kk represents gg, ll , or oiloil.

  1. The mass mkm_k is solved for each phase.
  2. Solve the momentum equation for mass flux mgug+mlulm_gu_g + m_lu_l.
  3. Use pressure equation to calculate pp.
  4. The density ρk\rho_k is yielded using EOS, and the volume fraction αk\alpha_k with known density.
  5. The velocity uku_k is obtained using the slip relation.
  6. If the kill mud is injected, the pump pressure is calculated by considering pressure drop in the blowout well and relief well.
  7. Advance to the next time step.


[1] H. Shi, J.A. Holmes, L.J. Durlofsky, K. Aziz, L.R. Diaz, B. Alkaya, and G. Oddie. Drift-flux modeling of multiphase flow in wellbores. In SPE Annual Technical Conference and Exhibition, Denver, Colorado, U.S.A., 2003. Society of Petroleum Engineers.

[2] A. R. Hasan and C. S. Kabir. A flow pattern independent drift flux model based void fraction correlation for a wide range of gasliquid two phase flow. International Journal of Multiphase Flow, 59:186–205, 2014.

[3] S. Evje and K. K. Fjelde. On a rough ausm scheme for a one-dimensional two-phase model. Computer & Fluids, 32:1497–1530, 2003.

[4] J. E. Udegbunam, K. K. Fjelde, S. Evje, and G. Nygaard. On the advection-upstream-splitting-method hybrid scheme: A simple transient-flow model for managed-pressure-drilling and underbalanced drilling applications. In SPE/IADC Managed Pressure Drilling & Underbalanced Operations Conference & Exhibition, Madrid, Spain, 2015. Society of Petroleum Engineers.

[5] J. E. Udegbunam, D. Sui, F. Moeinikia, A. C. V. M. Lage, K. K. Fjelde, Ø. Arild, H. A. Rabenjafimanantsoa, and F. H. Nygaard. A transient flow model for predicting pressure buildup in closed annuli. In Proceedings of the ASME 2017 36th International Conference on Ocean, Offshore and Arctic Engineering, Trondheim, Norway, 2017.

[6] K. Founargiotakis, V. C. Kelessidis, and R. Maglione. Laminar, transitional and turbulent flow of herschel-bulkley fluids in concentric annulus. The Canadian Journal of Chemical Engineering, 86:676–683, 2008

[7] O Erge, E. M. Ozbayoglu, S. Z. Miska, M. Yu, N. Takach, A. Saasen, R. May, and B. Hughes. The effects of drillstring eccentricity, rotation and buckling configurations on annular frictional pressure losses while circulating yield power law fluids. In SPE Drilling Conference and Exhibition, Fort Worth, Texas, USA, 2014