cfs_tpfa_residual.h File Reference
Include dependency graph for cfs_tpfa_residual.h:

Go to the source code of this file.

Classes

struct  cfs_tpfa_res_wells
 
struct  cfs_tpfa_res_forces
 
struct  cfs_tpfa_res_data
 

Functions

struct cfs_tpfa_res_datacfs_tpfa_res_construct (struct UnstructuredGrid *G, struct cfs_tpfa_res_wells *wells, int nphases)
 
void cfs_tpfa_res_destroy (struct cfs_tpfa_res_data *h)
 
int cfs_tpfa_res_assemble (struct UnstructuredGrid *G, double dt, struct cfs_tpfa_res_forces *forces, const double *zc, struct compr_quantities_gen *cq, const double *trans, const double *gravcap_f, const double *cpress, const double *wpress, const double *porevol, struct cfs_tpfa_res_data *h)
 
int cfs_tpfa_res_comprock_assemble (struct UnstructuredGrid *G, double dt, struct cfs_tpfa_res_forces *forces, const double *zc, struct compr_quantities_gen *cq, const double *trans, const double *gravcap_f, const double *cpress, const double *wpress, const double *porevol, const double *porevol0, const double *rock_comp, struct cfs_tpfa_res_data *h)
 
void cfs_tpfa_res_flux (struct UnstructuredGrid *G, struct cfs_tpfa_res_forces *forces, int np, const double *trans, const double *pmobc, const double *pmobf, const double *gravcap_f, const double *cpress, const double *wpress, double *fflux, double *wflux)
 
void cfs_tpfa_res_fpress (struct UnstructuredGrid *G, int np, const double *htrans, const double *pmobf, const double *gravcap_f, struct cfs_tpfa_res_data *h, const double *cpress, const double *fflux, double *fpress)
 

Detailed Description

Public interface to assembler for (compressible) discrete pressure system based on two-point flux approximation method. The assembler implements a residual formulation that for a single cell $i$ reads

\[ \mathit{pv}_i\cdot (1 - \sum_\alpha A_i^{-1}(p^{n+1},z^n)z^n) + \Delta t\sum_\alpha A_i^{-1}(p^{n+1},z^n) \Big(\sum_j A_{ij}v_{ij}^{n+1} - q_i^{n+1}\Big) = 0 \]

in which $\mathit{pv}_i$ is the (constant or pressure-dependent) pore-volume of cell $i$. Moreover, $\Delta t$ is the time step size, $n$ denotes the time level, $A$ is the pressure and mass-dependent fluid matrix that converts phase volumes at reservoir conditions into component volumes at surface conditions and $v_{ij}$ is the vector of outward (with respect to cell $i$) phase fluxes across the $ij$ cell-interface.

This module's usage model is intended to be

  1. Construct assembler
  2. for (each time step)
    1. while (pressure not converged)
      1. Assemble appropriate (Jacobian) system of linear equations
      2. Solve Jacobian system to derive pressure increments
      3. Include increments into current state
      4. Check convergence
    2. Derive fluxes and, optionally, interface pressures
    3. Solve transport by some means
  3. Destroy assembler

Function Documentation

int cfs_tpfa_res_assemble ( struct UnstructuredGrid G,
double  dt,
struct cfs_tpfa_res_forces forces,
const double *  zc,
struct compr_quantities_gen cq,
const double *  trans,
const double *  gravcap_f,
const double *  cpress,
const double *  wpress,
const double *  porevol,
struct cfs_tpfa_res_data h 
)

Assemble system of linear equations by linearising the residual around the current pressure point. Assume incompressible rock (i.e., that the pore-volume is independent of pressure).

The fully assembled system is presented in h->J and h->F and must be solved separately using external software.

Parameters
[in]GGrid.
[in]dtTime step size $\Delta t$.
[in]forcesDriving forces.
[in]zcComponent volumes, per pore-volume, at surface conditions for all components in all cells stored consecutively per cell. Array of size G->number_of_cells * cq->nphases.
[in]cqCompressible quantities describing the current fluid state. Fields Ac, dAc, Af, and phasemobf must be valid.
[in]transBackground transmissibilities as defined by function tpfa_trans_compute().
[in]gravcap_fDiscrete gravity and capillary forces.
[in]cpressCell pressures. One scalar value per grid cell.
[in]wpressWell (bottom-hole) pressures. One scalar value per well. NULL in case of no wells.
[in]porevolPore-volumes. One (positive) scalar value for each grid cell.
[in,out]hOn input-a valid (non-NULL) assembler obtained from a previous call to constructor function cfs_tpfa_res_construct(). On output-valid assembler that includes the new system of linear equations in its J and F fields.
Returns
1 if the assembled matrix was adjusted to remove a singularity. This happens if all fluids are incompressible and there are no pressure conditions on wells or boundaries. Otherwise return 0.
int cfs_tpfa_res_comprock_assemble ( struct UnstructuredGrid G,
double  dt,
struct cfs_tpfa_res_forces forces,
const double *  zc,
struct compr_quantities_gen cq,
const double *  trans,
const double *  gravcap_f,
const double *  cpress,
const double *  wpress,
const double *  porevol,
const double *  porevol0,
const double *  rock_comp,
struct cfs_tpfa_res_data h 
)

Assemble system of linear equations by linearising the residual around the current pressure point. Assume compressible rock (i.e., that the pore-volume depends on pressure).

The fully assembled system is presented in h->J and h->F and must be solved separately using external software.

Parameters
[in]GGrid.
[in]dtTime step size $\Delta t$.
[in]forcesDriving forces.
[in]zcComponent volumes, per pore-volume, at surface conditions for all components in all cells stored consecutively per cell. Array of size G->number_of_cells * cq->nphases.
[in]cqCompressible quantities describing the current fluid state. Fields Ac, dAc, Af, and phasemobf must be valid.
[in]transBackground transmissibilities as defined by function tpfa_trans_compute().
[in]gravcap_fDiscrete gravity and capillary forces.
[in]cpressCell pressures. One scalar value per grid cell.
[in]wpressWell (bottom-hole) pressures. One scalar value per well. NULL in case of no wells.
[in]porevolPore-volumes. One (positive) scalar value for each grid cell.
[in]porevol0Pore-volumes at start of time step (i.e., at time level $n$). One (positive) scalar value for each grid cell.
[in]rock_compRock compressibility. One non-negative scalar for each grid cell.
[in,out]hOn input-a valid (non-NULL) assembler obtained from a previous call to constructor function cfs_tpfa_res_construct(). On output-valid assembler that includes the new system of linear equations in its J and F fields.
Returns
1 if the assembled matrix was adjusted to remove a singularity. This happens if all fluids are incompressible, the rock is incompressible, and there are no pressure conditions on wells or boundaries. Otherwise return 0.
struct cfs_tpfa_res_data* cfs_tpfa_res_construct ( struct UnstructuredGrid G,
struct cfs_tpfa_res_wells wells,
int  nphases 
)

Construct assembler for system of linear equations.

Parameters
[in]GGrid
[in]wellsWell description. NULL in case of no wells. For backwards compatibility, the constructor also interprets (wells != NULL) && (wells->W == NULL) as "no wells present", but new code should use wells == NULL to signify "no wells".
[in]nphasesNumber of active fluid phases in this simulation run. Needed to correctly size various internal work arrays.
Returns
Fully formed assembler structure suitable for forming systems of linear equations using, e.g., function cfs_tpfa_res_assemble(). NULL in case of allocation failure. Must be destroyed using function cfs_tpfa_res_destroy().
void cfs_tpfa_res_destroy ( struct cfs_tpfa_res_data h)

Destroy assembler for system of linear equations.

Disposes of all resources acquired in a previous call to construction function cfs_tpfa_res_construct(). Note that the statement cfs_tpfa_res_destroy(NULL) is supported and benign (i.e., behaves like free(NULL)).

Parameters
[in,out]hOn input - assembler obtained through a previous call to construction function cfs_tpfa_res_construct(). On output - invalid pointer.
void cfs_tpfa_res_flux ( struct UnstructuredGrid G,
struct cfs_tpfa_res_forces forces,
int  np,
const double *  trans,
const double *  pmobc,
const double *  pmobf,
const double *  gravcap_f,
const double *  cpress,
const double *  wpress,
double *  fflux,
double *  wflux 
)

Derive interface (total) Darcy fluxes from (converged) pressure solution.

Parameters
[in]GGrid
[in]forcesDriving forces. Must correspond to those used when forming the system of linear equations, e.g., in the call to function cfs_tpfa_res_assemble().
[in]npNumber of fluid phases (and components).
[in]transBackground transmissibilities as defined by function tpfa_trans_compute(). Must correspond to equally named parameter of the assembly functions.
[in]pmobcPhase mobilities stored consecutively per cell with phase index cycling the most rapidly. Array of size G->number_of_cells * np.
[in]pmobfPhase mobilities stored consecutively per interface with phase index cycling the most rapidly. Array of size G->number_of_faces * np.
[in]gravcap_fDiscrete gravity and capillary forces.
[in]cpressCell pressure values. One (positive) scalar for each grid cell.
[in]wpressWell (bottom-hole) pressure values. One (positive) scalar value for each well. NULL in case of no wells.
[out]ffluxTotal Darcy interface fluxes. One scalar value for each interface (inter-cell connection). Array of size G->number_of_faces.
[out]wfluxTotal Darcy well connection fluxes. One scalar value for each well connection (perforation). Array of size forces->wells->W->well_connpos [forces->wells->W->number_of_wells].
void cfs_tpfa_res_fpress ( struct UnstructuredGrid G,
int  np,
const double *  htrans,
const double *  pmobf,
const double *  gravcap_f,
struct cfs_tpfa_res_data h,
const double *  cpress,
const double *  fflux,
double *  fpress 
)

Derive interface pressures from (converged) pressure solution.

Parameters
[in]GGrid
[in]npNumber of fluid phases (and components).
[in]htransBackground one-sided ("half") transmissibilities as defined by function tpfa_htrans_compute().
[in]pmobfPhase mobilities stored consecutively per interface with phase index cycling the most rapidly. Array of size G->number_of_faces * np.
[in]gravcap_fDiscrete gravity and capillary forces.
[in]hSystem assembler. Must correspond to the assembler state used to form the final system of linear equations from which the converged pressure solution was derived.
[in]cpressCell pressure values. One (positive) scalar for each grid cell.
[in]ffluxTotal Darcy interface fluxes. One scalar value for each interface (inter-cell connection). Array of size G->number_of_faces. Typically computed using function cfs_tpfa_res_flux().
[out]fpressInterface pressure values. One (positive) scalar for each interface. Array of size G->number_of_faces.