Implements a reordering transport solver for incompressible two-phase flow with polymer in the water phase. More...
#include <TransportSolverTwophasePolymer.hpp>
Classes | |
struct | ResidualC |
class | ResidualCGrav |
class | ResidualEquation |
struct | ResidualS |
class | ResidualSGrav |
Public Types | |
enum | SingleCellMethod { Bracketing, Newton, Gradient, NewtonSimpleSC, NewtonSimpleC } |
enum | GradientMethod { Analytic, FinDif } |
Public Member Functions | |
TransportSolverTwophasePolymer (const UnstructuredGrid &grid, const IncompPropertiesInterface &props, const PolymerProperties &polyprops, const SingleCellMethod method, const double tol, const int maxit) | |
Construct solver. | |
void | setPreferredMethod (SingleCellMethod method) |
Set the preferred method, Bracketing or Newton. | |
void | solve (const double *darcyflux, const double *porevolume, const double *source, const double *polymer_inflow_c, const double dt, std::vector< double > &saturation, std::vector< double > &concentration, std::vector< double > &cmax) |
Solve for saturation, concentration and cmax at next timestep. | |
void | solveGravity (const std::vector< std::vector< int > > &columns, const double *porevolume, const double dt, std::vector< double > &saturation, std::vector< double > &concentration, std::vector< double > &cmax) |
Solve for gravity segregation. | |
virtual void | solveSingleCell (const int cell) |
virtual void | solveMultiCell (const int num_cells, const int *cells) |
void | solveSingleCellBracketing (int cell) |
void | solveSingleCellNewton (int cell) |
void | solveSingleCellGradient (int cell) |
void | solveSingleCellNewtonSimple (int cell, bool use_sc) |
void | initGravity (const double *grav) |
void | solveSingleCellGravity (const std::vector< int > &cells, const int pos, const double *gravflux) |
int | solveGravityColumn (const std::vector< int > &cells) |
void | scToc (const double *x, double *x_c) const |
Implements a reordering transport solver for incompressible two-phase flow with polymer in the water phase.
Include permeability reduction effect.
Opm::TransportSolverTwophasePolymer::TransportSolverTwophasePolymer | ( | const UnstructuredGrid & | grid, | |
const IncompPropertiesInterface & | props, | |||
const PolymerProperties & | polyprops, | |||
const SingleCellMethod | method, | |||
const double | tol, | |||
const int | maxit | |||
) |
Construct solver.
[in] | grid | A 2d or 3d grid. |
[in] | props | Rock and fluid properties. |
[in] | polyprops | Polymer 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::TransportSolverTwophasePolymer::solve | ( | const double * | darcyflux, | |
const double * | porevolume, | |||
const double * | source, | |||
const double * | polymer_inflow_c, | |||
const double | dt, | |||
std::vector< double > & | saturation, | |||
std::vector< double > & | concentration, | |||
std::vector< double > & | cmax | |||
) |
Solve for saturation, concentration and cmax at next timestep.
Using implicit Euler scheme, reordered.
[in] | darcyflux | Array of signed face fluxes. |
[in] | porevolume | Array of pore volumes. |
[in] | source | Transport source term, to be interpreted by sign: (+) Inflow, value is first phase flow (water) per second, in reservoir volumes. (-) 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] | concentration | Polymer concentration. |
[in,out] | cmax | Highest concentration that has occured in a given cell. |
void Opm::TransportSolverTwophasePolymer::solveGravity | ( | const std::vector< std::vector< int > > & | columns, | |
const double * | porevolume, | |||
const double | dt, | |||
std::vector< double > & | saturation, | |||
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.
[in] | columns | Vector of cell-columns. |
[in] | porevolume | Array of pore volumes. |
[in] | dt | Time step. |
[in,out] | saturation | Phase saturations. |
[in,out] | concentration | Polymer concentration. |
[in,out] | cmax | Highest concentration that has occured in a given cell. |