Synonym: = Volatile/Black Oil Reservoir Flow @model = Muskat - Leverett equation

Content

Expand: Content

Content
- Definition
Foundation
Mathematical Model – Integral Format
Mathematical Model – Divergence Format
Initial Conditions
Well Model
The list of dynamic flow properties and model parameters
Computational Model

Definition

Expand: Definition

Mathematical model of Modified Black Oil fluid ( = Volatile Oil ) reservoir flow predicts the temperature, pressure and flow speed distribution in reservoir with account for:

available historical data on surface flowrates and/or bottom hole pressure
available 3D geological model
PVT and SCAL model
specific wellbore designs
gravitational forces
heat propagation
adiabatic and Jole-Thomson heat effects

The Black Oil flow is specific type of the Volatile Oil flow with $\begin{array}{l}R_v=0\end{array}$ .

Foundation

Expand – Foundation

Consider Modified Black Oil fluid @model consisting of three chemical components $\begin{array}{l}C = \{ W, \, O, \, G \}\end{array}$ which may exists in three phases $\begin{array}{l}f = \{ w, \, o, \, g \}\end{array}$

The driving equations of Modified Black Oil Reservoir Flow Grid Computation @model are:

(1)	$\begin{array}{l}\displaystyle \frac{dm^*_C}{dt} \, = \, \, \mbox{const}\end{array}$

mass conservation of $\begin{array}{l}m^*_W, \, m^*_O, \, m^*_G\end{array}$ – stock tank mass of the chemical components $\begin{array}{l}C = \{ W, \, O, \, G \}\end{array}$

(2)	$\begin{array}{l}\displaystyle \bar u_f = - k_p \cdot \frac{k_{rf}}{\mu_f} \cdot \big( \hat k * \bar \nabla} (p + p_{cf}) - \rho_f \cdot \hat k * \bar g \big)\, \quad f = \{ w, \, o, \, g \}, \quad p_{cf} = p_{cf}(s_f)\end{array}$

(3)	$\begin{array}{l}\displaystyle k_p = \exp(c_k \cdot (p-p_0))\end{array}$

Darcy transport equation dor the velocity $\begin{array}{l}\bar u_f\end{array}$ of three phases: $\begin{array}{l}f = \{ w, \, o, \, g \}\end{array}$

(4)	$\begin{array}{l}\displaystyle \phi = \phi_0 \cdot \exp(c_\phi \cdot (p-p_0))\end{array}$

Pore compressibility equation on reservoir porosity $\begin{array}{l}\phi(p)\end{array}$

(5)	$\begin{array}{l}\displaystyle \rho_f = \rho_f(T,p) \, , \quad \mu_f = \mu_f(T,p) \, , \quad R_s = R_s(T,p) \, , \quad \quad R_v = R_v(T,p)\end{array}$

Equation of state

(6)

$\begin{array}{l}\displaystyle (\rho \,c_{pt})_p \frac{\partial T}{\partial t} - \ \phi \sum_{f = \{w,o,g \}} \rho_f \ c_{pf} \ \eta_{sf} \ \frac{\partial p_f}{\partial t} + \bigg( \sum_{f = \{w,o,g \}} \rho_f \ c_{pf} \ \epsilon_f \ \mathbf{u}_f \bigg) \nabla p + \bigg( \sum_{f = \{w,o,g \}} \rho_f \ c_{pf} \ \mathbf{u}_f \bigg) \ \nabla T - \nabla (\lambda_t \nabla T) = \rho^{\downarrow} ({\bf r}) \, c^{\downarrow} ({\bf r}) \, T^{\downarrow} ({\bf r}) \, q^{\downarrow}({\bf r})\, \delta({\bf r})\end{array}$

Heat balance and transfer equation

Mathematical Model – Integral Format

Expand – Integral Format

(7)	$\begin{array}{l}\displaystyle \frac{\partial m_W}{\partial t} + \int\rho_{Ww} \cdot\, \bar u_w \, d \bar \Sigma = \frac{d m^*_W}{dt}\end{array}$

(8)	$\begin{array}{l}\displaystyle \frac{\partial m_O}{\partial t} + \int \left( \rho_{Oo} \cdot \, \bar u_o + \rho_{Og} \cdot \, \bar u_g \right) d \bar \Sigma = \frac{d m^*_O}{dt}\end{array}$

(9)	$\begin{array}{l}\displaystyle \frac{\partial m_G}{\partial t} + \int \left( \rho_{Go} \cdot \, \bar u_o + \rho_{Gg} \cdot \, \bar u_g \right) d \bar \Sigma = \frac{d m^*_G}{dt}\end{array}$

(10)	$\begin{array}{l}\displaystyle \rho_{Ww} = \frac{m_W}{V_w} = \frac{m_W}{V_W} \cdot \frac{V_W}{V_w} = \frac{\mathring \rho_W}{B_w}\end{array}$

(11)	$\begin{array}{l}\displaystyle \rho_{Oo} = \frac{m_{Oo}}{V_o} = \frac{m_{Oo}}{V_O} \cdot \frac{V_O}{V_o} = \frac{\mathring \rho_O}{B_o}\end{array}$

(12)	$\begin{array}{l}\displaystyle \rho_{Og} = \frac{m_{Og}}{V_g} = \frac{m_{Og}}{V_{Og}} \cdot \frac{V_{Og}}{V_G} \cdot \frac{V_G}{V_g} = \frac{\mathring \rho_O \cdot R_v}{B_g}\end{array}$

(13)	$\begin{array}{l}\displaystyle \rho_{Go} = \frac{m_{Go}}{V_o} = \frac{m_{Go}}{V_{Go}} \cdot \frac{V_{Go}}{V_O} \cdot \frac{V_O}{V_o} = \frac{\mathring \rho_G \cdot R_s}{B_o}\end{array}$

(14)	$\begin{array}{l}\displaystyle \rho_{Gg} = \frac{m_{Gg}}{V_g} = \frac{m_{Gg}}{V_G} \cdot \frac{V_G}{V_g} = \frac{\mathring \rho_G}{B_g}\end{array}$

(15)	$\begin{array}{l}\displaystyle \frac{\partial m_W}{\partial t} + \mathring \rho_W \int \frac{1}{B_w} \, \bar u_w \, d \bar \Sigma = \frac{d m^*_W}{dt}\end{array}$

(16)	$\begin{array}{l}\displaystyle \frac{\partial m_O}{\partial t} + \mathring \rho_O \int \left( \frac{1}{B_o} \cdot \, \bar u_o + \frac{R_v}{B_g} \cdot \, \bar u_g \right) d \bar \Sigma = \frac{d m^*_O}{dt}\end{array}$

(17)	$\begin{array}{l}\displaystyle \frac{\partial m_G}{\partial t} + \mathring \rho_G \int \left( \frac{R_s}{B_o} \cdot \, \bar u_o + \frac{1}{B_g} \cdot \, \bar u_g \right) d \bar \Sigma = \frac{d m^*_G}{dt}\end{array}$

(18)	$\begin{array}{l}\displaystyle V_W = \mathring \rho_W^{-1} \cdot m_W\end{array}$

(19)	$\begin{array}{l}\displaystyle V_O = \mathring \rho_O^{-1} \cdot m_O\end{array}$

(20)	$\begin{array}{l}\displaystyle V_G = \mathring \rho_G^{-1} \cdot m_G\end{array}$

(21)	$\begin{array}{l}\displaystyle q_W = \mathring \rho_W^{-1} \cdot \frac{d m^*_W}{dt}\end{array}$

(22)	$\begin{array}{l}\displaystyle q_O = \mathring \rho_O^{-1} \cdot \frac{d m^*_O}{dt}\end{array}$

(23)	$\begin{array}{l}\displaystyle q_G = \mathring \rho_G^{-1} \cdot \frac{d m^*_G}{dt}\end{array}$

(24)	$\begin{array}{l}\displaystyle \frac{\partial V_W}{\partial t} + \int \frac{1}{B_w} \, \bar u_w \, d \bar \Sigma = q_W\end{array}$

(25)	$\begin{array}{l}\displaystyle \frac{\partial V_O}{\partial t} + \int \left( \frac{1}{B_o} \cdot \, \bar u_o + \frac{R_v}{B_g} \cdot \, \bar u_g \right) d \bar \Sigma = q_O\end{array}$

(26)	$\begin{array}{l}\displaystyle \frac{\partial V_G}{\partial t} + \int \left( \frac{R_s}{B_o} \cdot \, \bar u_o + \frac{1}{B_g} \cdot \, \bar u_g \right) d \bar \Sigma = q_G\end{array}$

(27)	$\begin{array}{l}\displaystyle V_w = B_w \cdot V_W\end{array}$

(28)	$\begin{array}{l}\displaystyle V_o =\frac{B_o}{1-R_s \, R_v} \cdot (V_O - R_v \, V_G)\end{array}$

(29)	$\begin{array}{l}\displaystyle V_g =\frac{B_g}{1-R_s \, R_v} \cdot (V_G - R_s \, V_O)\end{array}$

(30)	$\begin{array}{l}\displaystyle V_\phi = V_w + V_o + V_g\end{array}$

(31)	$\begin{array}{l}\displaystyle s_w = \frac{V_w}{V_\phi}\end{array}$

(32)	$\begin{array}{l}\displaystyle s_o = \frac{V_o}{V_\phi}\end{array}$

(33)	$\begin{array}{l}\displaystyle s_g = \frac{V_g}{V_\phi}\end{array}$

(34)	$\begin{array}{l}\displaystyle s_w + s_o + s_g = 1\end{array}$

(35)	$\begin{array}{l}\displaystyle V_\phi = V_w + V_o + V_g = B_w \cdot V_W + \frac{B_o}{1-R_s \, R_v} \cdot \left( V_O - R_v \, V_G \right) + \frac{B_g}{1-R_s \, R_v} \cdot \left( V_G - R_s \, V_O \right)\end{array}$

(36)	$\begin{array}{l}\displaystyle V_\phi = V_w + V_o + V_g = B_w \cdot V_W + \frac{B_o - R_s \cdot B_g}{1-R_s \, R_v} \cdot V_O + \frac{B_g - R_v \cdot B_o}{1-R_s \, R_v} \cdot V_G\end{array}$

(37)	$\begin{array}{l}\displaystyle \phi = \frac{V_\phi}{V} = \frac{V_w}{V} + \frac{V_o}{V} + \frac{V_g}{V}\end{array}$

(38)	$\begin{array}{l}\displaystyle \phi = B_w \cdot \frac{V_W}{V} + \frac{B_o - R_s \cdot B_g}{1-R_s \, R_v} \cdot \frac{V_O}{V} + \frac{B_g - R_v \cdot B_o}{1-R_s \, R_v} \cdot \frac{V_G}{V}\end{array}$

(39)	$\begin{array}{l}\displaystyle \phi = \phi_0 \cdot \exp \big( c_\phi \cdot (p-p_0) \big)\end{array}$

(40)	$\begin{array}{l}\displaystyle B_w \cdot \frac{V_W}{V} + \frac{B_o - R_s \cdot B_g}{1-R_s \, R_v} \cdot \frac{V_O}{V} + \frac{B_g - R_v \cdot B_o}{1-R_s \, R_v} \cdot \frac{V_G}{V} = \phi_0 \cdot \exp \big( c_\phi \cdot (p-p_0) \big)\end{array}$

(41)	$\begin{array}{l}\displaystyle B_w \cdot \frac{V_W}{V_{\phi0}} + \frac{B_o - R_s \cdot B_g}{1-R_s \, R_v} \cdot \frac{V_O}{V_{\phi0}} + \frac{B_g - R_v \cdot B_o}{1-R_s \, R_v} \cdot \frac{V_G}{V_{\phi0}} = \exp \big( c_\phi \cdot (p-p_0) \big), \quad V_{\phi0} = V \cdot \phi_0\end{array}$

(42)	$\begin{array}{l}\displaystyle \tilde B_w = \tilde B_w(T,p) = B_w\end{array}$

(43)	$\begin{array}{l}\displaystyle \tilde B_o = \tilde B_o(T,p) = \frac{B_o - R_s \cdot B_g}{1-R_s \, R_v}\end{array}$

(44)	$\begin{array}{l}\displaystyle \tilde B_g = \tilde B_g(T,p) = \frac{B_g - R_v \cdot B_o}{1-R_s \, R_v}\end{array}$

(45)	$\begin{array}{l}\displaystyle \tilde B_w \cdot \frac{V_W}{V_{\phi0}} + \tilde B_o \cdot \frac{V_O}{V_{\phi0}} + \tilde B_g \cdot \frac{V_G}{V_{\phi0}} = \exp \big( c_\phi \cdot (p-p_0) \big), \quad V_{\phi0} = V \cdot \phi_0\end{array}$

(46)	$\begin{array}{l}\displaystyle \bar u_w = - k_p \cdot \frac{k_{rw}}{\mu_w} \cdot \hat k * \big( \bar \nabla} p_w - \rho_w \cdot \bar g \big)\end{array}$

(47)	$\begin{array}{l}\displaystyle \bar u_o = - k_p \cdot \frac{k_{ro}}{\mu_o} \cdot \hat k * \big( \bar \nabla p_o - \rho_o \cdot \bar g \big)\end{array}$

(48)	$\begin{array}{l}\displaystyle \bar u_g = - k_p \cdot \frac{k_{rg}}{\mu_g} \cdot \hat k * \big( \bar \nabla p_g - \rho_g \cdot \bar g \big)\end{array}$

(49)	$\begin{array}{l}\displaystyle k_p = \exp \big[ c_k \, (\,p - p_{\mathrm{ref}}\,) \big]\end{array}$

Given arbitrary coordinate system: $\begin{array}{l}\bold e = \{ \, \bar e_1, \, \bar e_2, \, \bar e_3 \, \}\end{array}$ :

(50)	$\begin{array}{l}\displaystyle \hat k * \bar v = \bar e_1 \cdot (k_{11} \, v_1 + k_{12} \, v_2 + k_{13} \, v_3) + \bar e_2 \cdot (k_{21} \, v_1 + k_{22} \, v_2 + k_{23} \, v_3) + \bar e_3 \cdot (k_{31} \, v_1 + k_{32} \, v_2 + k_{33} \, v_3)\end{array}$

(51)

$\begin{array}{l}\displaystyle \hat k * \bar \nabla = \bar e_1 \, (k_{11} \, \partial_1 + k_{12} \, \partial_2 + k_{13} \, \partial_3) + \bar e_2 \,(k_{21} \, \partial_2 + k_{22} \, \partial_2 + k_{23} \, \partial_3) + \bar e_3 \, (k_{31} \, \partial_1 + k_{32} \, \partial_2 + k_{33} \, \partial_3)\end{array}$

(52)	$\begin{array}{l}\displaystyle \hat k * \bar g = \bar e_1 \, (k_{11} \, g_1 + k_{12} \, g_2 + k_{13} \, g_3) + \bar e_2 \,(k_{21} \, g_1 + k_{22} \, g_2 + k_{23} \, g_3) + \bar e_3 \, (k_{31} \, g_1 + k_{32} \, g_2 + k_{33} \, g_3)\end{array}$

(53)	$\begin{array}{l}\displaystyle \rho_w = \frac{ \mathring \rho_W}{B_w}\end{array}$

(54)	$\begin{array}{l}\displaystyle \rho_o = \frac{ \mathring \rho_O + \mathring \rho_G \cdot R_s}{B_o}\end{array}$

(55)	$\begin{array}{l}\displaystyle \rho_g = \frac{ \mathring \rho_G + \mathring \rho_O \cdot R_v}{B_g}\end{array}$

(56)	$\begin{array}{l}\displaystyle \bar u_w = - k_p \cdot \frac{k_{rw}}{\mu_w} \cdot \big( \hat k * \bar \nabla} p_w - \rho_w \cdot \hat k * \bar g \big)\end{array}$

(57)	$\begin{array}{l}\displaystyle \bar u_o = - k_p \cdot \frac{k_{ro}}{\mu_o} \cdot \big( \hat k * \bar \nabla p_o - \rho_o \cdot \hat k * \bar g \big)\end{array}$

(58)	$\begin{array}{l}\displaystyle \bar u_g = - k_p \cdot \frac{k_{rg}}{\mu_g} \cdot \big( \hat k * \bar \nabla p_g - \rho_g \cdot \hat k * \bar g \big)\end{array}$

(59)	$\begin{array}{l}\displaystyle \bar u_w = - k_p \cdot \frac{k_{rw}}{\mu_w} \cdot \big( \hat k * \bar \nabla} (p + p_{cw}) - \rho_w \cdot \hat k * \bar g \big)\end{array}$

(60)	$\begin{array}{l}\displaystyle \bar u_o = - k_p \cdot \frac{k_{ro}}{\mu_o} \cdot \big( \hat k * \bar \nabla (p + p_{co}) - \rho_o \cdot \hat k * \bar g \big)\end{array}$

(61)	$\begin{array}{l}\displaystyle \bar u_g = - k_p \cdot \frac{k_{rg}}{\mu_g} \cdot \big( \hat k * \bar \nabla (p + p_{cg}) - \rho_g \cdot \hat k * \bar g \big)\end{array}$

(62)	$\begin{array}{l}\displaystyle p_w = p + p_{cw}, \quad p_{cw} = \frac{1}{3} \cdot ( -2 \cdot p_{cow} + p_{cog} )\end{array}$

(63)	$\begin{array}{l}\displaystyle p_o = p + p_{co}, \quad p_{co} = \frac{1}{3} \cdot (p_{cow} + p_{cog})\end{array}$

(64)	$\begin{array}{l}\displaystyle p_g = p + p_{cg}, \quad p_{cg} = \frac{1}{3} \cdot (p_{cow} - 2 \cdot p_{cog} )\end{array}$

(65)	$\begin{array}{l}\displaystyle p_{cw} = p_{cw}(s)\end{array}$

(66)	$\begin{array}{l}\displaystyle p_{cow} = p_{co}(s)\end{array}$

(67)	$\begin{array}{l}\displaystyle p_{cog} = p_{cg}(s)\end{array}$

(68)	$\begin{array}{l}\displaystyle p = \frac{1}{3} \left( p_o + p_g + p_w \right)\end{array}$

(69)	$\begin{array}{l}\displaystyle p_f = p + p_{cf}(s)\end{array}$

Mathematical Model – Divergence Format

Expand – Divergence Format

The Volatile Oil flow dynamics is defined by the following set of 3D equations:

(70)	$\begin{array}{l}\displaystyle \partial_t \bigg [ \phi \ \bigg ( \frac{s_w}{B_w} \bigg ) \bigg ] + \nabla \bigg ( \frac{1}{B_w} \ \mathbf{u}_w \bigg ) = q_W (\mathbf{r})\end{array}$

(71)	$\begin{array}{l}\displaystyle \partial_t \bigg [ \phi \ \bigg ( \frac{s_o}{B_o} + \frac{R_v \ s_g} {B_g} \bigg ) \bigg ] + \nabla \bigg ( \frac{1}{B_o} \ \mathbf{u}_o + \frac{R_v}{B_g} \ \mathbf{u}_g \bigg ) = q_O(\mathbf{r})\end{array}$

(72)	$\begin{array}{l}\displaystyle \partial_t \bigg [ \phi \ \bigg ( \frac{s_g}{B_g} + \frac{R_s \ s_o} {B_o} \bigg ) \bigg ] + \nabla \bigg ( \frac{1}{B_g} \ \mathbf{u}_g + \frac{R_s}{B_o} \ \mathbf{u}_o \bigg ) = q_G (\mathbf{r})\end{array}$

(73)	$\begin{array}{l}\displaystyle \mathbf{u}_w = - k_a \ \frac{k_{rw}(s_w, s_g)}{\mu_w} \ ( \nabla P_w - \rho_w \ \mathbf{g} )\end{array}$

(74)	$\begin{array}{l}\displaystyle \mathbf{u}_o = - k_a \ \frac{k_{ro}(s_w, s_g)}{\mu_o} \ ( \nabla P_o - \rho_o \ \mathbf{g} )\end{array}$

(75)	$\begin{array}{l}\displaystyle \mathbf{u}_g = - k_a \ \frac{k_{rg}(s_w, s_g)}{\mu_g} \ ( \nabla P_g - \rho_g \ \mathbf{g} )\end{array}$

(76)	$\begin{array}{l}\displaystyle P_o - P_w = P_{cow}(s_w)\end{array}$

(77)	$\begin{array}{l}\displaystyle P_o - P_g = P_{cog}(s_g)\end{array}$

(78)	$\begin{array}{l}\displaystyle s_w + s_o + s_g = 1\end{array}$

(79)

$\begin{array}{l}\displaystyle (\rho \,c_{pt})_p \frac{\partial T}{\partial t} - \ \phi \sum_{a = \{w,o,g \}} \rho_\alpha \ c_{p \alpha} \ \eta_{s \alpha} \ \frac{\partial P_\alpha}{\partial t} + \bigg( \sum_{a = \{w,o,g \}} \rho_\alpha \ c_{p \alpha} \ \epsilon_\alpha \ \mathbf{u}_\alpha \bigg) \nabla P + \bigg( \sum_{a = \{w,o,g \}} \rho_\alpha \ c_{p \alpha} \ \mathbf{u}_\alpha \bigg) \ \nabla T - \nabla (\lambda_t \nabla T) = \rho_{\rm inj} \, c_{\rm inj} \, T_{\rm inj} \, q_{\rm inj}({\bf r})\, \delta({\bf r})\end{array}$

The disambiguation of the properties in the above equation is brought in The list of dynamic flow properties and model parameters.

The right sides of equations (70) – (72) suggest no sources of flow except the contacts between wells and reservoir which is specified by well models as boundary conditions (see below).

Initial Conditions

Expand: Initial Conditions

Initial temperature distribution is set as input:

$\begin{array}{l}\displaystyle T(0, \mathbf{r}) = T_0(\mathbf{r})\end{array}$

In case the simulation is performed over the undisturbed reservoir then initial temperature distribution is geothermal.

The initial condition on phase pressure, phase velocities and phase saturations is set by one of the following options: Equilibrium Start and Non-equilibrium Start.

Condition I – Equilibrium Start

Equilibrium Start means that flow was not happening before the start: $\begin{array}{l}\{ \mathbf{u}_w = 0, \ \mathbf{u}_o = 0, \ \mathbf{u}_g =0 \}\end{array}$ and correspondingly phase pressure $\begin{array}{l}\{ P_w, \ P_o, \ P_g \}\end{array}$ and phase saturations $\begin{array}{l}\{ s_w, \ s_o, \ s_g, \}\end{array}$ were in stationary (not varying in time) conditions:

(80)	$\begin{array}{l}\displaystyle \nabla \cdot \bigg ( \frac{1}{B_w} \ \mathbf{u}_w \bigg )_{t=0} = 0\end{array}$

(81)	$\begin{array}{l}\displaystyle \nabla \cdot \bigg ( \frac{1}{B_o} \ \mathbf{u}_o + \frac{R_v}{B_g} \ \mathbf{u}_g \bigg )_{t=0} = 0\end{array}$

(82)	$\begin{array}{l}\displaystyle \nabla \cdot \bigg ( \frac{1}{B_g} \ \mathbf{u}_g + \frac{R_s}{B_o} \ \mathbf{u}_o \bigg )_{t=0} = 0\end{array}$

(83)	$\begin{array}{l}\displaystyle \mathbf{u}_w(0, \mathbf{r}) = - k_a \ \frac{k_{rw}(s_w, s_g)}{\mu_w} \ (\nabla P_w(0, \mathbf{r}) - \rho_w \ \mathbf{g} ) = 0\end{array}$

(84)	$\begin{array}{l}\displaystyle \mathbf{u}_o(0, \mathbf{r}) = - k_a \ \frac{k_{ro}(s_w, s_g)}{\mu_o} \ ( \nabla P_o(0, \mathbf{r}) - \rho_o \ \mathbf{g} ) = 0\end{array}$

(85)	$\begin{array}{l}\displaystyle \mathbf{u}_g(0, \mathbf{r}) = - k_a \ \frac{k_{rg}(s_w, s_g)}{\mu_g} \ ( \nabla P_g(0, \mathbf{r}) - \rho_g \ \mathbf{g} ) = 0\end{array}$

(86)	$\begin{array}{l}\displaystyle P_o(0, \mathbf{r}) - P_w(0, \mathbf{r}) = P_{cow}(s_w)\end{array}$

(87)	$\begin{array}{l}\displaystyle P_o(0, \mathbf{r}) - P_g(0, \mathbf{r}) = P_{cog}(s_g)\end{array}$

(88)	$\begin{array}{l}\displaystyle s_w + s_o + s_g = 1\end{array}$

Condition II – Non-equilibrium Start

Non-equilibrium Start means that flow happening before the start: $\begin{array}{l}\mathbf{u}_w^2 + \mathbf{u}_o^2 + \mathbf{u}_g^2 > 0\end{array}$ and correspondingly phase pressure $\begin{array}{l}\{ P_w, \ P_o, \ P_g \}\end{array}$ and phase saturations $\begin{array}{l}\{ s_w, \ s_o, \ s_g, \}\end{array}$ were in not in equilibrium:

(89)	$\begin{array}{l}\displaystyle s_w(0, \mathbf{r}) + s_o(0, \mathbf{r}) + s_g(0, \mathbf{r}) = 1\end{array}$

pressure distribution $\begin{array}{l}\{ P_w, \ P_o, \ P_g \}\end{array}$ could be arbitrary providing the capillary constraints:

(90)	$\begin{array}{l}\displaystyle P_o(0, \mathbf{r}) - P_w(0, \mathbf{r}) = P_{cow}(s_w)\end{array}$

(91)	$\begin{array}{l}\displaystyle P_o(0, \mathbf{r}) - P_g(0, \mathbf{r}) = P_{cog}(s_g).\end{array}$

The phase velocities are initialized as:

(92)	$\begin{array}{l}\displaystyle \mathbf{u}_w(0, \mathbf{r}) = - k_a \ \frac{k_{rw}(s_w, s_g)}{\mu_w} \ ( \nabla P_w(0, \mathbf{r}) - \rho_w \ \mathbf{g} )\end{array}$

(93)	$\begin{array}{l}\displaystyle \mathbf{u}_o(0, \mathbf{r}) = - k_a \ \frac{k_{ro}(s_w, s_g)}{\mu_o} \ ( \nabla P_o(0, \mathbf{r}) - \rho_o \ \mathbf{g} )\end{array}$

(94)	$\begin{array}{l}\displaystyle \mathbf{u}_g(0, \mathbf{r}) = - k_a \ \frac{k_{rg}(s_w, s_g)}{\mu_g} \ ( \nabla P_g(0, \mathbf{r}) - \rho_g \ \mathbf{g} )\end{array}$

In practice, the non-equilibrium conditions before the start is usually a result of previous flow simulations for the same reservoir, sometimes using a different grid-structure.

External Boundary Condition

Expand: External Boundary Condition

The external boundary condition for the temperature is usually set by one if the two options:

External Temperature Boundary Condition I – Fixed Temperature

$\begin{array}{l}\displaystyle T(t, \mathbf{r}) |_{\Gamma_e} = T_e( \mathbf{r})\end{array}$

External Temperature Boundary Condition II – Fixed Heat Exchange

(95)	$\begin{array}{l}\displaystyle \big( \mathbf{n}, \nabla T(t, \mathbf{r} \big) \big \|_{\Gamma_e} = \zeta \cdot \big( T(t, \mathbf{r}) - T_e( \mathbf{r}) \big)\end{array}$

where $\begin{array}{l}\zeta\end{array}$ – heat exchange coefficient at model boundary.

The external boundary condition for phase pressure, phase velocities and phase saturations is set by one of the two popular options:

External Pressure Boundary Condition I – Non-permeable boundary $\begin{array}{l}\Gamma_e\end{array}$

(96)	$\begin{array}{l}\displaystyle \big( \mathbf{n}, \ (\nabla P_\alpha(t, \mathbf{r}) - \rho_\alpha \mathbf{r}) \big) \big\|_{\Gamma_e} = 0\end{array}$

where $\begin{array}{l}\mathbf{n}\end{array}$ – normal vector to the boundary $\begin{array}{l}\Gamma_e\end{array}$ and $\begin{array}{l}\alpha = \{ w, o, g \}\end{array}$ .

External Pressure Boundary Condition II – constant-pressure boundary $\begin{array}{l}\Gamma_e\end{array}$

(97)	$\begin{array}{l}\displaystyle P_\alpha(t, \mathbf{r}) \big\|_{\Gamma_e} = P_i = const\end{array}$

where $\begin{array}{l}\alpha = \{ w, o, g \}\end{array}$ .

Well Model

Expand: Well Model

Well flow model (don't get confused with Wellbore Flow Model) simulates the flow at the contact between well and reservoir thus relating the sandface flow rates and pressure distribution in reservoir around the well:

(98)	$\begin{array}{l}\displaystyle P_w(t, \mathbf{r}) \big\|_{\Gamma_{WRC}} = P_o(t, \mathbf{r}) \big\|_{\Gamma_{WRC}} = P_g(t, \mathbf{r}) \big\|_{\Gamma_{WRC}}= P_{wf}(t) \big\|_{\Gamma_{WRC}}\end{array}$

where $\begin{array}{l}\Gamma_{WRC}\end{array}$ is a well-reservoir.

Bottom-hole pressure $\begin{array}{l}P_{wf}(t) \big|_{\Gamma_{WFC}} = P_{wf}(t,h)\end{array}$ at the contact (at depth $\begin{array}{l}h\end{array}$ along-hole) is set by one of the three popular conditions (traditionally called "Controls"):

Well Condition I – Pressure Control
Well Condition II – Liquid Control
Well Condition III – Oil Control

The list of dynamic flow properties and model parameters

Expand: The list of dynamic flow properties and model parameters

$\begin{array}{l}(t,x,y,z)\end{array}$	time and space corrdinates , $\begin{array}{l}z\end{array}$ -axis is orientated towards the Earth centre, $\begin{array}{l}(x,y)\end{array}$ define transversal plane to the $\begin{array}{l}z\end{array}$ -axis
$\begin{array}{l}\mathbf{r} = (x, \ y, \ z)\end{array}$	position vector at which the flow equations are set
$\begin{array}{l}q_{mW} = \frac{d m_W}{dt}\end{array}$	speed of water-component mass change in wellbore draining points
$\begin{array}{l}q_{mO} = \frac{d m_O}{dt}\end{array}$	speed of oil-component mass change in wellbore draining points
$\begin{array}{l}q_{mG} = \frac{d m_G}{dt}\end{array}$	speed of gas-component mass change in wellbore draining points
$\begin{array}{l}q_W = \frac{1}{\rho_W^{\LARGE \circ}} \frac{d m_W}{dt} = \frac{d V_{Ww}^{\LARGE \circ}}{dt} = \frac{1}{B_w} q_w\end{array}$	volumetric water-component flow rate in wellbore draining points recalculated to standard surface conditions
$\begin{array}{l}q_O = \frac{1}{\rho_O^{\LARGE \circ}} \frac{d m_O}{dt} = \frac{d V_{Oo}^{\LARGE \circ}}{dt} + \frac{d V_{Og}^{\LARGE \circ}}{dt} = \frac{1}{B_o} q_o + \frac{R_v}{B_g} q_g\end{array}$	volumetric oil-component flow rate in wellbore draining points recalculated to standard surface conditions
$\begin{array}{l}q_G = \frac{1}{\rho_G^{\LARGE \circ}} \frac{d m_G}{dt} = \frac{d V_{Gg}^{\LARGE \circ}}{dt} + \frac{d V_{Go}^{\LARGE \circ}}{dt} = \frac{1}{B_g} q_g + \frac{R_s}{B_o} q_o\end{array}$	volumetric gas-component flow rate in wellbore draining points recalculated to standard surface conditions
$\begin{array}{l}q_w = \frac{d V_w}{dt}\end{array}$	volumetric water-phase flow rate in wellbore draining points
$\begin{array}{l}q_o = \frac{d V_o}{dt}\end{array}$	volumetric oil-phase flow rate in wellbore draining points
$\begin{array}{l}q_g = \frac{d V_g}{dt}\end{array}$	volumetric gas-phase flow rate in wellbore draining points
$\begin{array}{l}q^S_W =\frac{dV_{Ww}^S}{dt}\end{array}$	total well volumetric water-component flow rate
$\begin{array}{l}q^S_O = \frac{d (V_{Oo}^S + V_{Og}^S )}{dt}\end{array}$	total well volumetric oil-component flow rate
$\begin{array}{l}q^S_G = \frac{d (V_{Gg}^S + V_{Go}^S )}{dt}\end{array}$	total well volumetric gas-component flow rate
$\begin{array}{l}q^S_L = q^S_W + q^S_O\end{array}$	total well volumetric liquid-component flow rate
$\begin{array}{l}P_w = P_w (t, \vec r)\end{array}$	water-phase pressure pressure distribution and dynamics
$\begin{array}{l}P_o = P_o (t, \vec r)\end{array}$	oil-phase pressure pressure distribution and dynamics
$\begin{array}{l}P_g = P_g (t, \vec r)\end{array}$	gas-phase pressure pressure distribution and dynamics
$\begin{array}{l}\vec u_w = \vec u_w (t, \vec r)\end{array}$	water-phase flow speed distribution and dynamics
$\begin{array}{l}\vec u_o = \vec u_o (t, \vec r)\end{array}$	oil-phase flow speed distribution and dynamics
$\begin{array}{l}\vec u_g = \vec u_g (t, \vec r)\end{array}$	gas-phase flow speed distribution and dynamics
$\begin{array}{l}P_{cow} = P_{cow} (s_w)\end{array}$	capillary pressure at the oil-water phase contact as function of water saturation
$\begin{array}{l}P_{cog} = P_{cog} (s_ g)\end{array}$	capillary pressure at the oil-gas phase contact as function of gas saturation
$\begin{array}{l}k_{rw} = k_{rw}(s_w, \ s_g)\end{array}$	relative formation permeability to water flow as function of water and gas saturation
$\begin{array}{l}k_{ro} = k_{ro}(s_w, \ s_g)\end{array}$	relative formation permeability to oil flow as function of water and gas saturation
$\begin{array}{l}k_{rg} = k_{rg}(s_w, \ s_g)\end{array}$	relative formation permeability to gas flow as function of water and gas saturation
$\begin{array}{l}\phi = \phi(P)\end{array}$	porosity as function of formation pressure
$\begin{array}{l}k_a = k_a(P)\end{array}$	absolute formation permeability to air
$\begin{array}{l}\vec g = (0, \ 0, \ g)\end{array}$	gravitational acceleration vector
$\begin{array}{l}g = 9.81 \ \rm m/s^2\end{array}$	gravitational acceleration constant
$\begin{array}{l}\rho_\alpha(P,T)\end{array}$	mass density of $\begin{array}{l}\alpha\end{array}$ -phase fluid
$\begin{array}{l}\mu_\alpha(P,T)\end{array}$	viscosity of $\begin{array}{l}\alpha\end{array}$ -phase fluid
$\begin{array}{l}\lambda_t(P,T,s_w, s_o, s_g)\end{array}$	effective thermal conductivity of the rocks with account for multiphase fluid saturation
$\begin{array}{l}\lambda_r(P,T)\end{array}$	rock matrix thermal conductivity
$\begin{array}{l}\lambda_\alpha(P,T)\end{array}$	thermal conductivity of $\begin{array}{l}\alpha\end{array}$ -phase fluid
$\begin{array}{l}\rho_r(P,T)\end{array}$	rock matrix mass density
$\begin{array}{l}\eta_{s \alpha}(P,T)\end{array}$	differential adiabatic coefficient of $\begin{array}{l}\alpha\end{array}$ -phase fluid
$\begin{array}{l}c_{pr}(P,T)\end{array}$	specific isobaric heat capacity of the rock matrix
$\begin{array}{l}c_{p\alpha}(P,T)\end{array}$	specific isobaric heat capacity of $\begin{array}{l}\alpha\end{array}$ -phase fluid
$\begin{array}{l}\epsilon_\alpha (P, T)\end{array}$	differential Joule–Thomson coefficient of $\begin{array}{l}\alpha\end{array}$ -phase fluid

Computational Model

Modified Black Oil Reservoir Flow Grid Computation @model

Page tree

Muskat - Leverett equation