Opm::TransportSolverTwophaseCompressiblePolymer Class Reference

Implements a reordering transport solver for incompressible two-phase flow with polymer in the water phase. More...

#include <TransportSolverTwophaseCompressiblePolymer.hpp>

Inheritance diagram for Opm::TransportSolverTwophaseCompressiblePolymer:

Classes

class  ResCOnCurve
 
struct  ResidualC
 
class  ResidualCGrav
 
class  ResidualEquation
 
struct  ResidualS
 
class  ResidualSGrav
 
class  ResSOnCurve
 

Public Types

enum  SingleCellMethod { Bracketing, Newton, NewtonC, Gradient }
 
enum  GradientMethod { Analytic, FinDif }
 

Public Member Functions

 TransportSolverTwophaseCompressiblePolymer (const UnstructuredGrid &grid, const BlackoilPropertiesInterface &props, const PolymerProperties &polyprops, const SingleCellMethod method, const double tol, const int maxit)
 Construct solver. More...
 
void setPreferredMethod (SingleCellMethod method)
 Set the preferred method, Bracketing or Newton.
 
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)
 Solve for saturation, concentration and cmax at next timestep. More...
 
void initGravity (const double *grav)
 Initialise quantities needed by gravity solver. More...
 
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)
 Solve for gravity segregation. More...
 

Friends

class TransportSolverTwophaseCompressiblePolymer::ResidualEquation
 
class TransportSolverTwophaseCompressiblePolymer::ResSOnCurve
 
class TransportSolverTwophaseCompressiblePolymer::ResCOnCurve
 

Detailed Description

Implements a reordering transport solver for incompressible two-phase flow with polymer in the water phase.

Include permeability reduction effect.

Constructor & Destructor Documentation

◆ TransportSolverTwophaseCompressiblePolymer()

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]gridA 2d or 3d grid.
[in]propsRock and fluid properties.
[in]polypropsPolymer properties.
[in]rock_compRock compressibility properties
[in]methodBracketing: 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]tolTolerance used in the solver.
[in]maxitMaximum number of non-linear iterations used.

Member Function Documentation

◆ initGravity()

void Opm::TransportSolverTwophaseCompressiblePolymer::initGravity ( const double *  grav)

Initialise quantities needed by gravity solver.

Parameters
[in]gravGravity vector

◆ solve()

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]darcyfluxArray of signed face fluxes.
[in]initial_pressureArray with pressure at start of timestep.
[in]pressureArray with pressure.
[in]temperatureArray with temperature.
[in]porevolume0Array with pore volume at start of timestep.
[in]porevolumeArray with pore volume.
[in]sourceTransport 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_cArray of inflow polymer concentrations per cell.
[in]dtTime step.
[in,out]saturationPhase saturations.
[in,out]surfacevolSurface volumes.
[in,out]concentrationPolymer concentration.
[in,out]cmaxHighest concentration that has occured in a given cell.

◆ solveGravity()

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]columnsVector of cell-columns.
[in]dtTime step.
[in,out]saturationPhase saturations.
[in,out]surfacevolSurface volumes.
[in,out]concentrationPolymer concentration.
[in,out]cmaxHighest concentration that has occured in a given cell.

The documentation for this class was generated from the following files: