Classes | Functions
hybsys.h File Reference

Go to the source code of this file.

Classes

struct  hybsys
 
struct  hybsys_well
 

Functions

struct hybsys * hybsys_allocate_symm (int max_nconn, int nc, int nconn_tot)
 
struct hybsys * hybsys_allocate_unsymm (int max_nconn, int nc, int nconn_tot)
 
struct hybsys_well * hybsys_well_allocate_symm (int max_nconn, int nc, int *cwpos)
 
struct hybsys_well * hybsys_well_allocate_unsymm (int max_nconn, int nc, int *cwpos)
 
void hybsys_free (struct hybsys *sys)
 
void hybsys_well_free (struct hybsys_well *wsys)
 
void hybsys_init (int max_nconn, struct hybsys *sys)
 
void hybsys_schur_comp_symm (int nc, const int *pconn, const double *Binv, struct hybsys *sys)
 
void hybsys_schur_comp_unsymm (int nc, const int *pconn, const double *Binv, const double *BIV, const double *P, struct hybsys *sys)
 
void hybsys_schur_comp_gen (int nc, const int *pconn, const double *Binv, const double *C2, const double *P, struct hybsys *sys)
 
void hybsys_well_schur_comp_symm (int nc, const int *cwpos, double *WI, struct hybsys *sys, struct hybsys_well *wsys)
 
void hybsys_cellcontrib_symm (int c, int nconn, int p1, int p2, const double *gpress, const double *src, const double *Binv, struct hybsys *sys)
 
void hybsys_cellcontrib_unsymm (int c, int nconn, int p1, int p2, const double *gpress, const double *src, const double *Binv, struct hybsys *sys)
 
void hybsys_well_cellcontrib_symm (int c, int ngconn, int p1, const int *cwpos, const double *WI, const double *wdp, struct hybsys *sys, struct hybsys_well *wsys)
 
void hybsys_compute_press_flux (int nc, const int *pconn, const int *conn, const double *gpress, const double *Binv, const struct hybsys *sys, const double *pi, double *press, double *flux, double *work)
 
void hybsys_compute_press_flux_well (int nc, const int *pgconn, int nf, int nw, const int *pwconn, const int *wconn, const double *Binv, const double *WI, const double *wdp, const struct hybsys *sys, const struct hybsys_well *wsys, const double *pi, double *cpress, double *cflux, double *wpress, double *wflux, double *work)
 

Detailed Description

Routines and data structures to manage local contributions to a global system of simultaneous linear equations arising from a Schur complement reduction of an original block system.

Specifically, these data structures and related routines compute and store the elemental (cell-based) contributions of the Schur complement reduction of the block system of simultaneous linear equations

\[ \begin{pmatrix} B & C_1 & D \\ C_2^\mathsf{T} & P & 0 \\ D^\mathsf{T} & 0 & 0 \end{pmatrix} \begin{pmatrix} v \\ -p \\ \pi \end{pmatrix} = \begin{pmatrix} G \\ g \\ h \end{pmatrix} \]

in which $G$ accounts for effects of gravity. The traditional Schurcomplement reduction (block Gaussian elimination) then produces the equivalent system of simultaneous linear equations

\[ \begin{pmatrix} B & C_1 & D \\ 0 & -L & -F_2 \\ 0 & 0 & A \end{pmatrix} \begin{pmatrix} v \\ -p \\ \pi \end{pmatrix} = \begin{pmatrix} G \\ g - C_2^\mathsf{T}B^{-1}G \\ b \end{pmatrix}. \]

Here, the matrix $A$ and the right hand side vector $b$ are given by

\[ \begin{aligned} A &= D^\mathsf{T}B^{-1}D - F_1^\mathsf{T}L^{-1}F_2 \\ b &= D^\mathsf{T}B^{-1}G + F_1^\mathsf{T}L^{-1}(g - C_2^\mathsf{T}B^{-1}G) - h, \end{aligned} \]

and the component matrices $F_1$, $F_2$, and $L$ are given by

\[ F_1 = C_1^\mathsf{T}B^{-1}D, \quad F_2 = C_2^\mathsf{T}B^{-1}D, \quad L = C_2^\mathsf{T}B^{-1}C_1 - P. \]

In the case of incompressible flow, the matrix $C_2$ is the same as $C_1$ and $P=0$ whence the coefficient matrix $A$ of the Schur complement system $A\pi=b$ is symmetric.

A great deal of simplification arises from the simple characterisation of the $C_1$ and $D$ matrices. Specifically,

\[ (C_1)_{ij} = \begin{cases} 1, &\quad i\in\{\mathit{pconn}_j, \dots, \mathit{pconn}_{j+1}-1\}, \\ 0, &\quad \text{otherwise}, \end{cases} \]

and

\[ (D)_{ij} = \begin{cases} 1, &\quad \mathit{conn}_i = j, \\ 0, &\quad \text{otherwise}. \end{cases} \]

When viewed in the context of a single cell, then the $D$ matrix is, effectively, the identity with the $\mathit{conn}$ array simply affecting a finite-element style redistribution (assembly) of the local contributions. This module leverages that property extensively.

Function Documentation

struct hybsys* hybsys_allocate_symm ( int  max_nconn,
int  nc,
int  nconn_tot 
)

Allocate a hybrid system management structure suitable for discretising a symmetric (i.e., incompressible) flow problem on a grid model of given size.

Parameters
[in]max_nconnMaximum number of single cell faces.
[in]ncTotal number of grid cells.
[in]nconn_totAggregate number of cell faces for all cells.
Returns
Fully formed hybrid system management structure if successful or NULL in case of allocation failure.
struct hybsys* hybsys_allocate_unsymm ( int  max_nconn,
int  nc,
int  nconn_tot 
)

Allocate a hybrid system management structure suitable for discretising an unsymmetric (i.e., compressible) flow problem on a grid model of given size.

Parameters
[in]max_nconnMaximum number of single cell faces.
[in]ncTotal number of grid cells.
[in]nconn_totAggregate number of cell faces for all cells.
Returns
Fully formed hybrid system management structure if successful or NULL in case of allocation failure.
void hybsys_cellcontrib_symm ( int  c,
int  nconn,
int  p1,
int  p2,
const double *  gpress,
const double *  src,
const double *  Binv,
struct hybsys *  sys 
)

Compute final (symmetric) Schur complement contributions to global system of simultaneous linear equations.

This function forms the coefficient matrix

\[ S_c = D^\mathsf{T}B_c^{-1}D - F_c^\mathsf{T}L_c^{-1}F_c \]

and similar right-hand side $r_c$ elemental contributions. These values must be subsequently assembled into the global system using function hybsys_global_assemble_cell() after imposing any applicable boundary conditions.

This function overwrites the fields S and r of the hybrid system structure.

Parameters
[in]cCell for which to compute local contributions.
[in]nconnNumber of connections (faces) of cell c.
[in]p1Start address (into gpress) of the gravity contributions of cell c.
[in]p2Start address (into Binv) of the inverse inner product of cell c.
[in]gpressGravity contributions of all cells. Must include effects of multiple phases if applicable.
[in]srcExplicit source terms for all cells.
[in]BinvInverse inner products for all cells. Must include effects of multiple phases if applicable.
[in,out]sysHybrid system management structure allocated using hybsys_allocate_symm() and initialised using hybsys_init() and/or filled using function hybsys_schur_comp_symm() and hybsys_well_schur_comp_symm() if applicable.
void hybsys_cellcontrib_unsymm ( int  c,
int  nconn,
int  p1,
int  p2,
const double *  gpress,
const double *  src,
const double *  Binv,
struct hybsys *  sys 
)

Compute final (non-symmetric) Schur complement contributions to global system of simultaneous linear equations.

This function forms the coefficient matrix

\[ S_c = D^\mathsf{T}B_c^{-1}D - (F_1)_c^\mathsf{T}L_c^{-1}(F_2)_c \]

and similar right-hand side $r_c$ elemental contributions. These values must be subsequently assembled into the global system using function hybsys_global_assemble_cell() after imposing any applicable boundary conditions.

This function overwrites the fields S and r of the hybrid system structure.

Parameters
[in]cCell for which to compute local contributions.
[in]nconnNumber of connections (faces) of cell c.
[in]p1Start address (into gpress) of the gravity contributions of cell c.
[in]p2Start address (into Binv) of the inverse inner product of cell c.
[in]gpressGravity contributions of all cells. Must include effects of multiple phases if applicable.
[in]srcExplicit source terms for all cells.
[in]BinvInverse inner products for all cells. Must include effects of multiple phases if applicable.
[in,out]sysHybrid system management structure allocated using hybsys_allocate_symm() and initialised using hybsys_init() and/or filled using functions hybsys_schur_comp_unsymm() or hybsys_schur_comp_gen().
void hybsys_compute_press_flux ( int  nc,
const int *  pconn,
const int *  conn,
const double *  gpress,
const double *  Binv,
const struct hybsys *  sys,
const double *  pi,
double *  press,
double *  flux,
double *  work 
)

Recover cell pressures and outward fluxes (with respect to cells–i.e., the ``half-face fluxes'') through back substitution after solving a symmetric (i.e., incompressible) Schur complement system of simultaneous linear equations.

Specifically, given the solution $\pi$ to the global system of simultaneous linear equations, $A\pi=b$, that arises as a result of the Schur complement analysis, this function recovers the cell pressures $p$ and outward fluxes $v$ defined by

\[ \begin{aligned} Lp &= g - C_2^\mathsf{T}B^{-1}G + F_2\pi \\ Bv &= G + C_1p - D\pi \end{aligned}. \]

Parameters
[in]ncTotal number of grid cells.
[in]pconnCell-to-face start pointers.
[in]connCell-to-face mapping.
[in]gpressGravity contributions of all cells. Must coincide with equally named parameter in calls to cell contribution functions such as hybsys_cellcontrib_symm().
[in]BinvInverse inner products for all cells. Must coincide with equally named parameter in calls to contribution functions such as hybsys_cellcontrib_symm().
[in]sysHybrid system management structure coinciding with equally named parameter in contribution functions such as hybsys_cellcontrib_symm() or hybsys_cellcontrib_unsymm().
[in]piSolution (interface/contact pressure) obtained from solving the global system $A\pi = b$.
[out]pressCell pressures, $p$. Array of size nc.
[out]fluxOutward interface fluxes, $v$. Array of size pconn[nc].
[in,out]workScratch array for temporary results. Array of size at least $\max_c \{ \mathit{pconn}_{c + 1} - \mathit{pconn}_c \} $.
void hybsys_compute_press_flux_well ( int  nc,
const int *  pgconn,
int  nf,
int  nw,
const int *  pwconn,
const int *  wconn,
const double *  Binv,
const double *  WI,
const double *  wdp,
const struct hybsys *  sys,
const struct hybsys_well *  wsys,
const double *  pi,
double *  cpress,
double *  cflux,
double *  wpress,
double *  wflux,
double *  work 
)

Recover well pressures (i.e., bottom-hole pressure values) and well connection (perforation) fluxes.

Specifically, this function performs the same role (i.e., back-substitution) for wells as function hybsys_compute_press_flux() does for grid cells and grid contacts (interfaces).

Parameters
[in]ncTotal number of grid cells.
[in]pgconnCell-to-face start pointers.
[in]nfTotal number of grid faces.
[in]nwTotal number of wells.
[in]pwconnCell-to-well start pointers. If nw > 0, then this parameter must coincide with the cwpos array used in call to hybsys_well_schur_comp_symm().
[in]wconnCell-to-well mapping.
[in]BinvInverse inner products for all cells. Must coincide with equally named parameter in calls to contribution functions such as hybsys_well_cellcontrib_symm().
[in]WIPeaceman well connection indices. Array of size pwconn[nc]. Must coincide with equally named parameter of contribution function hybsys_well_cellcontrib_symm().
[in]wdpWell connection gravity pressure adjustments.
[in]sysHybrid system management structure coinciding with equally named parameter in contribution functions such as hybsys_cellcontrib_symm() and hybsys_well_cellcontrib_symm().
[in]wsysHybrid well-system management structure. Must coincide with equally named paramter of contribution function hybsys_well_cellcontrib_symm().
[in]piSolution (interface/contact pressure and well BHPs) obtained from solving the global system $A\pi = b$.
[in]cpressCell pressures, $p$, obtained from a previous call to function hybsys_compute_press_flux().
[in]cfluxOutward fluxes, $v$, obtained from a previous call to function hybsys_compute_press_flux().
[out]wpressWell (i.e., bottom-hole) pressures. Array of size nw.
[out]wfluxWell connection (perforation) fluxes. Array of size pwconn[nw].
[in,out]workScratch array for storing intermediate results. Array of size at least $\max_w \{ \mathit{pwconn}_{w + 1} - \mathit{pwconn}_w\}$.
void hybsys_free ( struct hybsys *  sys)

Dispose of memory resources previously obtained through one of the allocation functions, hybsys_allocate_symm() or hybsys_allocate_unsymm().

Following a call to hybsys_free(), the input pointer is no longer valid. hybsys_free(NULL) does nothing.

Parameters
[in,out]sysPreviously allocated hybrid system management structure (or NULL).
void hybsys_init ( int  max_nconn,
struct hybsys *  sys 
)

Perform post-construction dynamic initialisation of system structure obtained from function hybsys_allocate_symm() or hybsys_allocate_unsymm().

Parameters
[in]max_nconnMaximum number of single cell faces. Must coincide with the equally named parameter of functions hybsys_allocate_symm() or hybsys_allocate_unsymm().
[in,out]sysPreviously allocated hybrid system management structure.
void hybsys_schur_comp_gen ( int  nc,
const int *  pconn,
const double *  Binv,
const double *  C2,
const double *  P,
struct hybsys *  sys 
)

Compute elemental (per-cell) contributions to unsymmetric Schur system of simultaneous linear equations.

This function assumes that the coefficient matrix of the hybrid system of linear equations is that of the introduction with no additional provisions. In other words, this function assumes that the coefficient matrix is of the form

\[ \begin{pmatrix} B & C_1 & D \\ C_2^\mathsf{T} & P & 0 \\ D^\mathsf{T} & 0 & 0 \end{pmatrix}. \]

This function fills the F1, F2 and L fields of the management structure.

Parameters
[in]ncTotal number of grid cells.
[in]pconnCell-to-face start pointers.
[in]BinvInverse inner product results, usually computed using mim_ip_simple_all() and mim_ip_mobility_update().
[in]C2Explicit representation of the $C_2$ matrix as a linear array. Assumed to only contain the (structurally) non-zero matrix elements (that correspond to the non-zero structure of $C_1$).
[in]PPer cell compressible accumulation term. One scalar per cell.
[in,out]sysHybrid system management structure allocated using hybsys_allocate_symm() and initialised using hybsys_init().
void hybsys_schur_comp_symm ( int  nc,
const int *  pconn,
const double *  Binv,
struct hybsys *  sys 
)

Compute elemental (per-cell) contributions to symmetric Schur system of simultaneous linear equations.

This function assumes that the coefficient matrix of the hybrid system of linear equations is that of the introduction with the additional provision that $C_1=C_2=C$ and that $P=0$. In other words, this function assumes that the coefficient matrix is of the form

\[ \begin{pmatrix} B & C & D \\ C^\mathsf{T} & 0 & 0 \\ D^\mathsf{T} & 0 & 0 \end{pmatrix}. \]

This function fills the F1 and L fields of the management structure.

Parameters
[in]ncTotal number of grid cells.
[in]pconnCell-to-face start pointers.
[in]BinvInverse inner product results, usually computed using mim_ip_simple_all() and mim_ip_mobility_update().
[in,out]sysHybrid system management structure allocated using hybsys_allocate_symm() and initialised using hybsys_init().
void hybsys_schur_comp_unsymm ( int  nc,
const int *  pconn,
const double *  Binv,
const double *  BIV,
const double *  P,
struct hybsys *  sys 
)

Compute elemental (per-cell) contributions to unsymmetric Schur system of simultaneous linear equations.

This function assumes that the coefficient matrix of the hybrid system of linear equations is that of the introduction with the additional provision that $C_2=C_1-V$. In other words, this function assumes that the coefficient matrix is of the form

\[ \begin{pmatrix} B & C & D \\ (C-V)^\mathsf{T} & P & 0 \\ D^\mathsf{T} & 0 & 0 \end{pmatrix}. \]

This matrix arises in the `` $v^2$'' phase compressibility formulation of the compressible black-oil model. This function fills the F1, F2 and L fields of the management structure.

Parameters
[in]ncTotal number of grid cells.
[in]pconnCell-to-face start pointers.
[in]BinvInverse inner product results, usually computed using mim_ip_simple_all() and mim_ip_mobility_update().
[in]BIV$B^{-1}v$ in which $v$ is the flux field of a previous time step or non-linear iteration.
[in]PPer cell compressible accumulation term. One scalar per cell.
[in,out]sysHybrid system management structure allocated using hybsys_allocate_symm() and initialised using hybsys_init().
struct hybsys_well* hybsys_well_allocate_symm ( int  max_nconn,
int  nc,
int *  cwpos 
)

Allocate a hybrid system management structure suitable for discretising an incompressible (i.e., symmetric) well flow problem on a grid model of given size.

Parameters
[in]max_nconnMaximum number of single cell faces.
[in]ncTotal number of grid cells.
[in]cwposIndirection array that defines each cell's connecting wells. Values typically computed using function derive_cell_wells().
Returns
Fully formed hybrid system management structure if successful or NULL in case of allocation failure.
struct hybsys_well* hybsys_well_allocate_unsymm ( int  max_nconn,
int  nc,
int *  cwpos 
)

Allocate a hybrid system management structure suitable for discretising a compressible (i.e., unsymmetric) well flow problem on a grid model of given size.

Parameters
[in]max_nconnMaximum number of single cell faces.
[in]ncTotal number of grid cells.
[in]cwposIndirection array that defines each cell's connecting wells. Values typically computed using function derive_cell_wells().
Returns
Fully formed hybrid system management structure if successful or NULL in case of allocation failure.
void hybsys_well_cellcontrib_symm ( int  c,
int  ngconn,
int  p1,
const int *  cwpos,
const double *  WI,
const double *  wdp,
struct hybsys *  sys,
struct hybsys_well *  wsys 
)

Form elemental direct contributions to global system of simultaneous linear equations from cell<->well interactions.

Plays a role similar to function hybsys_cellcontrib_symm(), but for wells.

Parameters
[in]cCell for which to compute cell<->well Schur complement
[in]ngconnNumber of inter-cell connections (faces) of cell c.
[in]p1Start index (into sys->F1) of cell c.
[in]cwposIndirection array that defines each cell's connecting wells. Must coincide with equally named parameter of function hybsys_well_schur_comp_symm().
[in]WIPeaceman well connection indices. Array of size pwconn[nc]. Must coincide with equally named parameter of contribution function hybsys_well_schur_comp_symm().
[in]wdpWell connection gravity pressure adjustments. One scalar for each well connection in an array of size pwconn[nc].
[in,out]sysHybrid system management structure filled using functions hybsys_schur_comp_unsymm() or hybsys_schur_comp_gen().
[in,out]wsysHybrid well-system management structure filled using function hybsys_well_schur_comp_symm().
void hybsys_well_free ( struct hybsys_well *  wsys)

Dispose of memory resources previously obtained through one of the allocation functions, hybsys_well_allocate_symm() or hybsys_well_allocate_unsymm().

Following a call to hybsys_well_free(), the input pointer is no longer valid. hybsys_well_free(NULL) does nothing.

Parameters
[in,out]wsysPreviously allocated hybrid system management structure (or NULL).
void hybsys_well_schur_comp_symm ( int  nc,
const int *  cwpos,
double *  WI,
struct hybsys *  sys,
struct hybsys_well *  wsys 
)

Compute elemental contributions to global, symmetric system of simultaneous linear equations from cell<->well connections.

Specifically, for a well w intersecting a cell c, this function computes the elemental contributions

\[ (F_1)_{wc} = C_{wc}^\mathsf{T} B_{wc}^{-1} D_{wc} = \mathit{WI}_{wc} \]

and

\[ L_{wc} = C_{wc}^\mathsf{T} B_{wc}^{-1} C_{wc} = \mathit{WI}_{wc} \]

and incorporates the contributions into the global system quantities as appropriate.

This function modifies sys->L and wsys->F1.

Parameters
[in]ncTotal number of grid cells.
[in]cwposIndirection array that defines each cell's connecting wells. Values typically computed using function derive_cell_wells().
[in]WIPeaceman well connection indices. Array of size cwpos[nc]. Must incorporate effects of multiple phases (i.e., total mobility) if applicable.
[in,out]sysHybrid system management structure allocated using hybsys_allocate_symm() and initialised using hybsys_init() and/or filled using function hybsys_schur_comp_symm().
[in,out]wsysHybrid well-system management structure obtained from function hybsys_well_allocate_symm().