The Well Model


The well is one of the most important components in a dynamic reservoir model, but also one of the most challenging to implement due to the different nature of the fluid flows occurring in the well and in the rock. While the fluid flow in the porous media is described by the Darcy’s flow regime, in the well pipes a variety of flow regimes (mist, churn, bubble, annual) exist, with a transition zone at the rock face.
Following the most common approach used in reservoir engineering, the section of the well connecting the reservoir and the surface facilities is modeled as one-dimensional, and its flow is considered to be much faster than the Darcy flow taking place in the reservoir, i.e. characterized by much smaller time scales. Therefore, the well response is modeled as instantaneous. The fluid flow in the tubing string to surface can often be described using tubing performance or ’lift curves’. In the well section open to the reservoir, the dynamic effects (e.g. inertia, friction, etc) can be neglected, and only the pressure variations due to hydrostatic effects are modeled.
Currently, in PFLOTRAN only the well section open to the reservoir is modeled, and both the coupling with tubing stream performance tables and the surface facilities are not treated.
In this brief description, “well connection” is the fluid flow connection between the well bore and its perforated grid block.
The figure below shows a schematic representation of a well perforating several grid blocks:
Neglecting the dynamic effects, the well pressure can be described by one additional variable, the well pressure at a specific reference depth, commonly known as the well BHP (bottom hole pressure). The pressure at various depths in the well is obtained by adding the hydrostatic corrections (\(H_j\)).

Inflow Performance Relationship

The equation governing the fluid exchange in each well connection is modelled by an inflow performance relation (IPR):

(1)\[q^w_{ α, j} =\Gamma_{w,j} \lambda_{\alpha,j} ( p_{\alpha,j} − p^{\mbox{ref}}_w −H_{w,j} )\]


\[\begin{split}\begin{split} &\alpha&=\mbox{fluid phase index} \\ &j&=\mbox{well connection index} \\ &q_{\alpha,j}^w &=\mbox{volumetric−rate for the } \alpha \mbox{ phase entering the well through the j connection} \\ &p_{α,j} &=\mbox{grid block pressure of the } \alpha \mbox{ phase for the j th connection} \\ &H_{w,j} &=\mbox{hydrostatic pressure correction for j th connection} \\ &\lambda_{m,j} &=\mbox{mobilty of the } \alpha \mbox{ phase for the j th connection} \\ &\Gamma_{w,j} &=\mbox{well connection factor (see below)} \\ \end{split}\end{split}\]

The well connection factor, which relates the pressure drop between the well and a perforated grid block, with the volumetric fluid rate, is modeled with an expression proposed by Peaceman [Pea]:

(2)\[\Gamma_{w,j} = \frac{2 \pi k h }{ | \ln(r_0/r_w) + S | },\]


\[\begin{split}\begin{split} &k &=( K_i K_j )^{(1/2)} \\ &i,j&=\mbox{directions defining the plane perpendicularly drilled by the well} \\ &K_i , K_j &= \mbox{perforated grid bock permeeabilites} \\ &r_w &=\mbox{well radius} \\ &h&= \mbox{perforated grid block tickness} \\ &S&=\mbox{skin factor} \\ &r_0 &=\mbox{grid block pressure equivalent radius (see below)} \end{split}\end{split}\]

and the grid block pressure equation radius \(r_0\) is given by:

(3)\[r_0 = \frac{ \left[ \Delta x_i \left( \frac{K_j}{K_i} \right)^{(1/2)} \Delta x_j \left( \frac{K_i}{K_j} \right)^{(1/2)} \right] }{ \left( \frac{K_j}{K_i} \right)^{(1/4)} + \left( \frac{K_i}{K_j} \right)^{(1/4)} }\]

Note the above is for a well completion in the z-direction (the default). If the direction of completion were the x-direction, then h in the deminator would become Δx for the cell, and the dimensions appearing in (3) would be \(\Delta y\) and \(\Delta y\).

The WELL_DATA well model

The well model is an implicit multi-mode which allows automatic switching between multiple well targets.
As mentioned in the previous section, a well is a collection of completions which provide connections to the reservoir. A connection is generally a section of a well trajectory which is completed (i.e. is connected to the reservoir by shooting through the well lining or by being unlined) and which lies within a reservoir simulation cell. Suppose a well has \(N_k\) completions, each with a completion connection factor of \(\Gamma_k =ccf_k\) . The use of the completion connection factor follows Eclipse© practice. Let the depth at which the well bottom hole pressure, \(P_w\), is defined be \(d_w\). Suppose the reservoir model contains \(N_c\) components, flowing in \(N_p\) phases.
The total well flow rate of component c is obtained by a sum over completions:
(4)\[f_{cw} = \sum_k{ f_{ck}}\]
where \(f_{ck}\) is the flow rate of component \(c\) through completion \(k\). A slight generalisation of equation (4) to allow for multiple components in a phase (for example dissolved gas in the black oil mode) yields:
(5)\[f_{ck} =ccf_k \sum_p{ M cp P_{pddk}},\]
where \(M_cp\) is the mobility of component \(c\), and \(P_{pddk}\) is the drawdown pressure for phase \(p\) from completion \(k\), given by:
\[P_{pddk} =P_{pa} -P_w +P_{hh}\]
\(P_{pa}\) is the phase pressure in the connected active cell \(a\), and \(P_{hh}\) is a pressure correction from \(d_w\) to the completion depth \(d_k\).
The term \(P_{hh}\) will generally include terms for gravity head, friction and fluid acceleration. However, at present we will consider the gravity correction only, and \(P_{h}\) is the ‘hydrostatic head’:
\[P_{hh} =g \int{\mathbf{e}_z \rho_w (z) d\mathbf{t}}\]
The integral is taken along the well trajectory from measured depth \(md_w\) to measured depth \(md_k\) . Measured depth is simply the distance along the well trajectory. If we assume a single value for the fluid density in the well we get:
\[P_{hh} =g.(z_k -z_w )\rho_w\]
Where \(rho_w\) is the mass density of the wellbore fluid, and \(g\) is the gravity constant.
Note that, in this work, production is regarded as a positive flow, and depth is measured upwards (following PFLOTRAN and not most simulation codes). So, if the well reference depth is above the completion \(k\), \(z_w >z_k\) , then \(P_{hh}\) will be negative and the pressure in the well at depth \(z_k\), \((P_w -P_{hh} )\), will be greater that the pressure at \(z_w\).
For a well operating on a surface flow rate target, \(T\), we require that the well residual condition is satisfied:
\[R_w =T(P_w )-T_{req} = 0\]
where \(T_{req}\) is the required target rate – normally a surface rate condition such as surface oil rate. This is a condition on \(P_w\).
In general, there will be a number of different target conditions on a well (bottom hole pressure, oil rate, liquid rate, gas rate etc.) and the well model will cycle through these until it finds the most restrictive one. When the correct target condition is found it is used to set the well operating point, and all the other conditions should be satisfied.
Once the \(P_w\) value needed to satisfy the target has been found, this can be substituted back into the well flow equations. This will yield the well flows as a function of \(P_w\) and the cell solution variables \(\{X_a \}\), for all the active cells with completions on the well. In addition, we require the derivatives \(df/dX\) of the flows with respect to the cell variables.
It is possible to convert derivatives with respect to \(P_w\) into further derivatives with respect to \(\{X_a \}\) using the implicit function theorem:
\[\frac{dP_w}{dX_a} = -\left(\frac{\partial R_w}{\partial P_w}\right)^{-1} \left(\frac{\partial R_w}{\partial X_a} \right)\]
This yields a fully coupled system – as all the flow rates depend on \(P_w\) , and \(P_w\) depends on all the variables \(\{X_a \}\) in the cells coupled to the well, all the flows finally depend on all the coupled variables. This treatment is used in the first PFLOTRAN well model. In modes which use numerical differentiation to find the Jacobian, the implicit function theorem is not required, as the SNES differencing in PFLOTRAN will find all the cell variable dependencies.

Crossflow and non-crossflow models

In a cross-flow model, we allow some of the completions of a producer to inject, and some of the completions of an injector to produce. This is a common occurrence in any situation in which the fluid in the wellbore has a different density to that in the reservoir, so that the hydrostatic gradient is different inside and outside the well. The extreme case of crossflow is, or course, a ‘stopped well’, in which the top of the well is closed, so that well is at zero surface rate, but the completions are left open. This is common when a well is shut, particularly for a short period. The alternative is to completely close the well by cementing the existing completions.
In the OGS_1.3 release of PFLOTRAN for reservoir engineering, a cross-flow model has been implemented.

Transient and non-transient models

In a transient model, we would include the effect of the finite volume of the well. The change in the amount of fluid in the wellbore is related to the flows over the time step. In a non-transient model, the well volume is regarded as negligible and the well model is assumed to be in steady state at the implicit (end of time-step) conditions.
For the present, we consider a non-transient model, as this is less likely to impose time-step limitations. In a transient well, the well is effectively being modelled as a low-pore-volume cell or set of cells, and the fast changes in the well saturations can limit the step-length that the non-linear solver can handle. The non-transient model has few disadvantages but will be unable to model wellbore storage effects in well test analysis. Wellbore volumes are, of course, generally negligible compared to reservoir volumes.

Single variable and multi-variable well models

The well bhp, \(P_w\), (bhp being the bottom hole pressure) is commonly used as a solution variable for the well. It is possible to solve the well model using just this variable. The big advantage is then that single-variable non-linear solvers exists which will converge unconditionally to the correct result for a continuous function for which a solution exists – an example is the bisection algorithm. For more than one variable a Newton method is generally used but is not as reliable.
In a non-crossflow well model a simple single variable well model is quite sufficient. The composition in the wellbore is known, either as the sum of the producting completion flows, or from the nature of the well injection fluid. In the OGS 1.3 cross-flow model, the wellbore composition is not initially known, and for any given \(P_w\), an inner consistency condition is converged-out to determine the wellbore density.