Public interface to assembler for (compressible) discrete pressure system based on two-point flux approximation method. More...
#include <opm/core/grid.h>
#include <opm/core/wells.h>
#include <opm/core/pressure/tpfa/compr_source.h>
Go to the source code of this file.
Classes | |
struct | cfs_tpfa_res_wells |
Type encapsulating well topology and completion data (e.g., phase mobilities per connection (perforation)). More... | |
struct | cfs_tpfa_res_forces |
Type encapsulating all driving forces affecting the discrete pressure system. More... | |
struct | cfs_tpfa_res_data |
Result structure that presents the fully assembled system of linear equations, linearised around the current pressure point. More... | |
Functions | |
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. More... | |
void | cfs_tpfa_res_destroy (struct cfs_tpfa_res_data *h) |
Destroy assembler for system of linear equations. More... | |
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. More... | |
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. More... | |
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. More... | |
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. More... | |
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 reads
in which is the (constant or pressure-dependent) pore-volume of cell
. Moreover,
is the time step size,
denotes the time level,
is the pressure and mass-dependent fluid matrix that converts phase volumes at reservoir conditions into component volumes at surface conditions and
is the vector of outward (with respect to cell
) phase fluxes across the
cell-interface.
This module's usage model is intended to be
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.
[in] | G | Grid. |
[in] | dt | Time step size ![]() |
[in] | forces | Driving forces. |
[in] | zc | Component 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] | cq | Compressible quantities describing the current fluid state. Fields Ac , dAc , Af , and phasemobf must be valid. |
[in] | trans | Background transmissibilities as defined by function tpfa_trans_compute(). |
[in] | gravcap_f | Discrete gravity and capillary forces. |
[in] | cpress | Cell pressures. One scalar value per grid cell. |
[in] | wpress | Well (bottom-hole) pressures. One scalar value per well. NULL in case of no wells. |
[in] | porevol | Pore-volumes. One (positive) scalar value for each grid cell. |
[in,out] | h | On 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. |
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.
[in] | G | Grid. |
[in] | dt | Time step size ![]() |
[in] | forces | Driving forces. |
[in] | zc | Component 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] | cq | Compressible quantities describing the current fluid state. Fields Ac , dAc , Af , and phasemobf must be valid. |
[in] | trans | Background transmissibilities as defined by function tpfa_trans_compute(). |
[in] | gravcap_f | Discrete gravity and capillary forces. |
[in] | cpress | Cell pressures. One scalar value per grid cell. |
[in] | wpress | Well (bottom-hole) pressures. One scalar value per well. NULL in case of no wells. |
[in] | porevol | Pore-volumes. One (positive) scalar value for each grid cell. |
[in] | porevol0 | Pore-volumes at start of time step (i.e., at time level ![]() |
[in] | rock_comp | Rock compressibility. One non-negative scalar for each grid cell. |
[in,out] | h | On 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. |
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.
[in] | G | Grid |
[in] | wells | Well 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] | nphases | Number of active fluid phases in this simulation run. Needed to correctly size various internal work arrays. |
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)
).
[in,out] | h | On 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.
[in] | G | Grid |
[in] | forces | Driving 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] | np | Number of fluid phases (and components). |
[in] | trans | Background transmissibilities as defined by function tpfa_trans_compute(). Must correspond to equally named parameter of the assembly functions. |
[in] | pmobc | Phase mobilities stored consecutively per cell with phase index cycling the most rapidly. Array of size G->number_of_cells * np . |
[in] | pmobf | Phase mobilities stored consecutively per interface with phase index cycling the most rapidly. Array of size G->number_of_faces * np . |
[in] | gravcap_f | Discrete gravity and capillary forces. |
[in] | cpress | Cell pressure values. One (positive) scalar for each grid cell. |
[in] | wpress | Well (bottom-hole) pressure values. One (positive) scalar value for each well. NULL in case of no wells. |
[out] | fflux | Total Darcy interface fluxes. One scalar value for each interface (inter-cell connection). Array of size G->number_of_faces . |
[out] | wflux | Total 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.
[in] | G | Grid |
[in] | np | Number of fluid phases (and components). |
[in] | htrans | Background one-sided ("half") transmissibilities as defined by function tpfa_htrans_compute(). |
[in] | pmobf | Phase mobilities stored consecutively per interface with phase index cycling the most rapidly. Array of size G->number_of_faces * np . |
[in] | gravcap_f | Discrete gravity and capillary forces. |
[in] | h | System 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] | cpress | Cell pressure values. One (positive) scalar for each grid cell. |
[in] | fflux | Total 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] | fpress | Interface pressure values. One (positive) scalar for each interface. Array of size G->number_of_faces . |