Opm::TransportSolverTwophaseReorder Class Reference

Implements a reordering transport solver for incompressible two-phase flow. More...

#include <TransportSolverTwophaseReorder.hpp>

Inheritance diagram for Opm::TransportSolverTwophaseReorder:
Inheritance graph

Public Member Functions

 TransportSolverTwophaseReorder (const UnstructuredGrid &grid, const Opm::IncompPropertiesInterface &props, const double *gravity, const double tol, const int maxit)
 
virtual ~TransportSolverTwophaseReorder ()
 
virtual void solve (const double *porevolume, const double *source, const double dt, TwophaseState &state)
 
void solveGravity (const double *porevolume, const double dt, TwophaseState &state)
 
const std::vector< int > & getReorderIterations () const
 

Private Member Functions

void reorderAndTransport (const UnstructuredGrid &grid, const double *darcyflux)
 
const std::vector< int > & sequence () const
 
const std::vector< int > & components () const
 

Detailed Description

Implements a reordering transport solver for incompressible two-phase flow.

Constructor & Destructor Documentation

Opm::TransportSolverTwophaseReorder::TransportSolverTwophaseReorder ( const UnstructuredGrid grid,
const Opm::IncompPropertiesInterface props,
const double *  gravity,
const double  tol,
const int  maxit 
)

Construct solver.

Parameters
[in]gridA 2d or 3d grid.
[in]propsRock and fluid properties.
[in]gravityGravity vector (null for no gravity).
[in]tolTolerance used in the solver.
[in]maxitMaximum number of non-linear iterations used.
virtual Opm::TransportSolverTwophaseReorder::~TransportSolverTwophaseReorder ( )
virtual

Member Function Documentation

const std::vector<int>& Opm::TransportSolverTwophaseReorder::getReorderIterations ( ) const
virtual void Opm::TransportSolverTwophaseReorder::solve ( const double *  porevolume,
const double *  source,
const double  dt,
TwophaseState state 
)
virtual

Solve for saturation at next timestep. Note that this only performs advection by total velocity, and no gravity segregation.

Parameters
[in]porevolumeArray of pore volumes.
[in]sourceTransport source term. For interpretation see Opm::computeTransportSource().
[in]dtTime step.
[in,out]stateReservoir state. Calling solve() will read state.faceflux() and read and write state.saturation().

Implements Opm::TransportSolverTwophaseInterface.

void Opm::TransportSolverTwophaseReorder::solveGravity ( const double *  porevolume,
const double  dt,
TwophaseState state 
)

Solve for gravity segregation. This uses a column-wise nonlinear Gauss-Seidel approach. It assumes that the grid can be divided into vertical columns that do not interact with each other (for gravity segregation).

Parameters
[in]porevolumeArray of pore volumes.
[in]dtTime step.
[in,out]stateReservoir state. Calling solveGravity() will read state.faceflux() and read and write state.saturation().

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