#include <TransportSolverTwophaseCompressiblePolymer.hpp>
|
| TransportSolverTwophaseCompressiblePolymer (const UnstructuredGrid &grid, const BlackoilPropertiesInterface &props, const PolymerProperties &polyprops, const SingleCellMethod method, const double tol, const int maxit) |
|
void | setPreferredMethod (SingleCellMethod method) |
| Set the preferred method, Bracketing or Newton. More...
|
|
void | solve (const double *darcyflux, const std::vector< double > &initial_pressure, const std::vector< double > &pressure, const std::vector< double > &temperature, const double *porevolume0, const double *porevolume, const double *source, const double *polymer_inflow_c, const double dt, std::vector< double > &saturation, std::vector< double > &surfacevol, std::vector< double > &concentration, std::vector< double > &cmax) |
|
void | initGravity (const double *grav) |
|
void | solveGravity (const std::vector< std::vector< int > > &columns, const double dt, std::vector< double > &saturation, std::vector< double > &surfacevol, std::vector< double > &concentration, std::vector< double > &cmax) |
|
Implements a reordering transport solver for incompressible two-phase flow with polymer in the water phase. Include permeability reduction effect.
Enumerator |
---|
Analytic |
|
FinDif |
|
Enumerator |
---|
Bracketing |
|
Newton |
|
NewtonC |
|
Gradient |
|
Opm::TransportSolverTwophaseCompressiblePolymer::TransportSolverTwophaseCompressiblePolymer |
( |
const UnstructuredGrid & |
grid, |
|
|
const BlackoilPropertiesInterface & |
props, |
|
|
const PolymerProperties & |
polyprops, |
|
|
const SingleCellMethod |
method, |
|
|
const double |
tol, |
|
|
const int |
maxit |
|
) |
| |
Construct solver. - Parameters
-
[in] | grid | A 2d or 3d grid. |
[in] | props | Rock and fluid properties. |
[in] | polyprops | Polymer properties. |
[in] | rock_comp | Rock compressibility properties |
[in] | method | Bracketing: solve for c in outer loop, s in inner loop, each solve being bracketed for robustness. Newton: solve simultaneously for c and s with Newton's method. (using gradient variant and bracketing as fallbacks). |
[in] | tol | Tolerance used in the solver. |
[in] | maxit | Maximum number of non-linear iterations used. |
void Opm::TransportSolverTwophaseCompressiblePolymer::initGravity |
( |
const double * |
grav | ) |
|
Initialise quantities needed by gravity solver. - Parameters
-
void Opm::TransportSolverTwophaseCompressiblePolymer::setPreferredMethod |
( |
SingleCellMethod |
method | ) |
|
Set the preferred method, Bracketing or Newton.
void Opm::TransportSolverTwophaseCompressiblePolymer::solve |
( |
const double * |
darcyflux, |
|
|
const std::vector< double > & |
initial_pressure, |
|
|
const std::vector< double > & |
pressure, |
|
|
const std::vector< double > & |
temperature, |
|
|
const double * |
porevolume0, |
|
|
const double * |
porevolume, |
|
|
const double * |
source, |
|
|
const double * |
polymer_inflow_c, |
|
|
const double |
dt, |
|
|
std::vector< double > & |
saturation, |
|
|
std::vector< double > & |
surfacevol, |
|
|
std::vector< double > & |
concentration, |
|
|
std::vector< double > & |
cmax |
|
) |
| |
Solve for saturation, concentration and cmax at next timestep. Using implicit Euler scheme, reordered. - Parameters
-
[in] | darcyflux | Array of signed face fluxes. |
[in] | initial_pressure | Array with pressure at start of timestep. |
[in] | pressure | Array with pressure. |
[in] | temperature | Array with temperature. |
[in] | porevolume0 | Array with pore volume at start of timestep. |
[in] | porevolume | Array with pore volume. |
[in] | source | Transport source term, to be interpreted by sign: (+) Inflow, value is first phase flow (water) per second, in surface volumes (unlike the incompressible version). (-) Outflow, value is total flow of all phases per second, in reservoir volumes. |
[in] | polymer_inflow_c | Array of inflow polymer concentrations per cell. |
[in] | dt | Time step. |
[in,out] | saturation | Phase saturations. |
[in,out] | surfacevol | Surface volumes. |
[in,out] | concentration | Polymer concentration. |
[in,out] | cmax | Highest concentration that has occured in a given cell. |
void Opm::TransportSolverTwophaseCompressiblePolymer::solveGravity |
( |
const std::vector< std::vector< int > > & |
columns, |
|
|
const double |
dt, |
|
|
std::vector< double > & |
saturation, |
|
|
std::vector< double > & |
surfacevol, |
|
|
std::vector< double > & |
concentration, |
|
|
std::vector< double > & |
cmax |
|
) |
| |
Solve for gravity segregation. This uses a column-wise nonlinear Gauss-Seidel approach. It assumes that the input columns contain cells in a single vertical stack, that do not interact with other columns (for gravity segregation. - Parameters
-
[in] | columns | Vector of cell-columns. |
[in] | dt | Time step. |
[in,out] | saturation | Phase saturations. |
[in,out] | surfacevol | Surface volumes. |
[in,out] | concentration | Polymer concentration. |
[in,out] | cmax | Highest concentration that has occured in a given cell. |
friend class TransportSolverTwophaseCompressiblePolymer::ResCOnCurve |
|
friend |
friend class TransportSolverTwophaseCompressiblePolymer::ResidualEquation |
|
friend |
friend class TransportSolverTwophaseCompressiblePolymer::ResSOnCurve |
|
friend |
The documentation for this class was generated from the following file:
|