Compositional Multiphase Well Solver

Introduction

Here, we present a description of the well solvers. These solvers are designed to be coupled with the flow solvers. Their specific task is to implement the multi-segment well discretization using the fluid model used in the corresponding flow solver – i.e., either single-phase flow or compositional multiphase flow. In particular, the perforation terms computed by the well solvers are a source/sink term for the discrete reservoir equations assembled by the flow solvers.

In the present description, we focus on the compositional multiphase well solver. The structure of the single-phase well solver is analogous and will not be described here for brevity.

Theory

Here, we give an overview of the well formulation implemented in GEOSX. We review the set of Discrete Equations, and then we describe the Primary variables used the well solvers.

Discrete Equations

We assume that the well is discretized into segments denoted by the index i.

Mass Conservation

In well segment i, mass conservation for component c reads:

z^{upw}_{c,(i-1,i)} q_{(i-1,i)} - z^{upw}_{c,(i,i+1)} q_{(i,i+1)} + q^{perf}_{c,i} = 0,

where we have neglected the accumulation term. In the previous equation, z^{upw}_{c,(i,j)} is the upwinded mass fraction of component c at the interface between segments i and j, q_{(i,j)} is the total mass flux between segments i and j, and q^{perf}_{c,i} is the source/sink term generated by the perforations – i.e., the connections between the well and the reservoir. The upwinded mass fractions are computed as:

z^{upw}_{c,(i,j)} = \left\{ \begin{array}{cl} z_{c,i} & \text{if} \, \, q_{(i,j)} > 0 \\[10pt] z_{c,j} & \text{otherwise.} \end{array}\right.

The perforation terms are obtained with:

q^{perf}_{c,i} = \left\{ \begin{array}{cl} WI z_{c,i} \rho_{m,i} \lambda^{res}_T  \Delta \Phi & \text{if the well is upstream (i.e.,} \, \, \Delta \Phi > 0)  \\[10pt] WI x^{res}_{c,\ell} \rho^{res}_{\ell} \lambda^{res}_{\ell} \Delta \Phi & \text{otherwise,} \end{array}\right.

where \Delta \Phi = p_i - p^{res} + \rho_{m,i} g \Delta d_{i,perf} is the potential difference between the segment center and the reservoir center. In the expression of the potential difference, the mixture density is computed as \rho_{m,i} = \sum_{\ell} S_{\ell,i} \rho_{\ell,i}. The well index, WI, is currently an input of the simulator. The superscript res means that the variable is evaluated at the center of the reservoir element.

Volume Constraint Equation

As in the Compositional Multiphase Flow Solver, the system is closed with a volume constraint equation.

Pressure Relations

In the current implementation of the well solver, we assume a hydrostatic equilibrium:

p_{i+1} - p_i = \rho_{m,(i,i+1)} g \Delta d_{i,i+1},

where \rho_{m,(i,i+1)} is the arithmetic average of the mixture densities evaluated in segments i and i+1. Pressure drop components due to friction and acceleration are not implemented at the moment.

Pressure and Rate Controls

The well solver supports two types of control, namely, pressure control and rate control.

If pressure control is chosen, we add the following constraint for the pressure of the top segment of the well:

p_0 - p^{target} = 0

In this case, we check that at each iteration of the Newton solver, the rate at the top of the first segment is smaller than the maximum rate specified by the user. If this is not the case, we switch to rate control.

If rate control is used, we add the following constraint for the rate at the top of the first segment, denoted by q_{(-1,0)}:

q_{(-1,0)} - q^{target} = 0

If the pressure at the top segment becomes larger than the maximum pressure specified by the user, then we switch to pressure control.

To summarize, the compositional multiphase flow solver assembles a set of n_c+2 equations, i.e., n_c mass conservation equations and 1 volume constraint equation in each segment, plus 1 pressure relation at the interface between a segment and the next segment in the direction of the well head. For the top segment, the pressure relation is replaced with the control equation.

Number of equations Equation type
n_c Mass conservation equations
1 Pressure relation or control equation
1 Volume constraint

Primary variables

The well variable formulation is the same as that of the Compositional Multiphase Flow Solver. In a well segment, in addition to the n_c+1 primary variables of the Compositional Multiphase Flow Solver, namely, one pressure, p, and n_c component densities, \rho_c, we also treat the total mass flux at the interface with the next segment, denoted by q, as a primary variable.

Number of primary variables Variable type
1 Pressure
1 Total mass flux at the interface with next segment
n_c Component densities

Parameters

The following attributes are supported:

Name Type Default Description
allowLocalCompDensityChopping integer 1 Flag indicating whether local (cell-wise) chopping of negative compositions is allowed
cflFactor real64 0.5 Factor to apply to the CFL condition when calculating the maximum allowable time step. Values should be in the interval (0,1]
fluidNames string_array required Name of fluid constitutive object to use for this solver.
initialDt real64 1e+99 Initial time-step value required by the solver to the event manager.
logLevel integer 0 Log level
maxCompFractionChange real64 1 Maximum (absolute) change in a component fraction between two Newton iterations
name string required A name is required for any non-unique nodes
relPermNames string_array required Names of relative permeability constitutive models to use
targetRegions string_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager.
useMass integer 0 Use mass formulation instead of molar
wellTemperature real64 required Temperature
LinearSolverParameters node unique Element: LinearSolverParameters
NonlinearSolverParameters node unique Element: NonlinearSolverParameters
WellControls node   Element: WellControls

Example

    <CompositionalMultiphaseWell
      name="compositionalMultiphaseWell"
      logLevel="1"
      targetRegions="{ wellRegion1, wellRegion2, wellRegion3 }"
      fluidNames="{ fluid1 }"
      relPermNames="{ relperm }"
      wellTemperature="297.15"
      useMass="0">
      <WellControls
        name="wellControls1"
        type="producer"
        control="BHP"
        targetBHP="4e6"
        targetRate="1"/>
      <WellControls
        name="wellControls2"
        type="producer"
        control="BHP"
        targetBHP="4e6"
        targetRate="1"/>
      <WellControls
        name="wellControls3"
        type="injector"
        control="liquidRate"
        targetBHP="1e8"
        targetRate="5e-3"
        injectionStream="{ 0.1, 0.1, 0.1, 0.7 }"/>
    </CompositionalMultiphaseWell>