#include <IncompTpfa.hpp>

Public Member Functions

 IncompTpfa (const UnstructuredGrid &grid, const IncompPropertiesInterface &props, LinearSolverInterface &linsolver, const double *gravity, const Wells *wells, const std::vector< double > &src, const FlowBoundaryConditions *bcs)
 
 IncompTpfa (const UnstructuredGrid &grid, const IncompPropertiesInterface &props, const RockCompressibility *rock_comp_props, LinearSolverInterface &linsolver, const double residual_tol, const double change_tol, const int maxiter, const double *gravity, const Wells *wells, const std::vector< double > &src, const FlowBoundaryConditions *bcs)
 
virtual ~IncompTpfa ()
 Destructor. More...
 
void solve (const double dt, TwophaseState &state, WellState &well_state)
 
const std::vector< double > & getHalfTrans () const
 Expose read-only reference to internal half-transmissibility. More...
 

Protected Member Functions

void solveIncomp (const double dt, TwophaseState &state, WellState &well_state)
 
void solveRockComp (const double dt, TwophaseState &state, WellState &well_state)
 

Protected Attributes

const UnstructuredGridgrid_
 
const IncompPropertiesInterfaceprops_
 
const RockCompressibilityrock_comp_props_
 
const LinearSolverInterfacelinsolver_
 
const double residual_tol_
 
const double change_tol_
 
const int maxiter_
 
const double * gravity_
 
const Wellswells_
 
const std::vector< double > & src_
 
const FlowBoundaryConditionsbcs_
 
std::vector< double > htrans_
 
std::vector< double > gpress_
 
std::vector< int > allcells_
 
std::vector< double > trans_
 
std::vector< double > wdp_
 
std::vector< double > totmob_
 
std::vector< double > omega_
 
std::vector< double > gpress_omegaweighted_
 
std::vector< double > initial_porevol_
 
struct ifs_tpfa_forces forces_
 
std::vector< double > porevol_
 
std::vector< double > rock_comp_
 
std::vector< double > pressures_
 
struct ifs_tpfa_datah_
 

Detailed Description

Encapsulating a tpfa pressure solver for the incompressible-fluid case. Supports gravity, wells controlled by bhp or reservoir rates, boundary conditions and simple sources as driving forces. Rock compressibility can be included, and necessary nonlinear iterations are handled. Below we use the shortcuts D for the number of dimensions, N for the number of cells and F for the number of faces.

Constructor & Destructor Documentation

Opm::IncompTpfa::IncompTpfa ( const UnstructuredGrid grid,
const IncompPropertiesInterface props,
LinearSolverInterface linsolver,
const double *  gravity,
const Wells wells,
const std::vector< double > &  src,
const FlowBoundaryConditions bcs 
)

Construct solver for incompressible case.

Parameters
[in]gridA 2d or 3d grid.
[in]propsRock and fluid properties.
[in]linsolverLinear solver to use.
[in]gravityGravity vector. If non-null, the array should have D elements.
[in]wellsThe wells argument. Will be used in solution, is ignored if NULL. Note: this class observes the well object, and makes the assumption that the well topology and completions does not change during the run. However, controls (only) are allowed to change.
[in]srcSource terms. May be empty().
[in]bcsBoundary conditions, treat as all noflow if null.
Opm::IncompTpfa::IncompTpfa ( const UnstructuredGrid grid,
const IncompPropertiesInterface props,
const RockCompressibility rock_comp_props,
LinearSolverInterface linsolver,
const double  residual_tol,
const double  change_tol,
const int  maxiter,
const double *  gravity,
const Wells wells,
const std::vector< double > &  src,
const FlowBoundaryConditions bcs 
)

Construct solver, possibly with rock compressibility.

Parameters
[in]gridA 2d or 3d grid.
[in]propsRock and fluid properties.
[in]rock_comp_propsRock compressibility properties. May be null.
[in]linsolverLinear solver to use.
[in]residual_tolSolution accepted if inf-norm of residual is smaller.
[in]change_tolSolution accepted if inf-norm of change in pressure is smaller.
[in]maxiterMaximum acceptable number of iterations.
[in]gravityGravity vector. If non-null, the array should have D elements.
[in]wellsThe wells argument. Will be used in solution, is ignored if NULL. Note: this class observes the well object, and makes the assumption that the well topology and completions does not change during the run. However, controls (only) are allowed to change.
[in]srcSource terms. May be empty().
[in]bcsBoundary conditions, treat as all noflow if null.
virtual Opm::IncompTpfa::~IncompTpfa ( )
virtual

Destructor.

Member Function Documentation

const std::vector<double>& Opm::IncompTpfa::getHalfTrans ( ) const
inline

Expose read-only reference to internal half-transmissibility.

References htrans_.

void Opm::IncompTpfa::solve ( const double  dt,
TwophaseState state,
WellState well_state 
)

Solve the pressure equation. If there is no pressure dependency introduced by rock compressibility effects, the equation is linear, and it is solved directly. Otherwise, the nonlinear equations ares solved by a Newton-Raphson scheme. May throw an exception if the number of iterations exceed maxiter (set in constructor).

void Opm::IncompTpfa::solveIncomp ( const double  dt,
TwophaseState state,
WellState well_state 
)
protected
void Opm::IncompTpfa::solveRockComp ( const double  dt,
TwophaseState state,
WellState well_state 
)
protected

Member Data Documentation

std::vector<int> Opm::IncompTpfa::allcells_
protected
const FlowBoundaryConditions* Opm::IncompTpfa::bcs_
protected
const double Opm::IncompTpfa::change_tol_
protected
struct ifs_tpfa_forces Opm::IncompTpfa::forces_
protected
std::vector<double> Opm::IncompTpfa::gpress_
protected
std::vector<double> Opm::IncompTpfa::gpress_omegaweighted_
protected
const double* Opm::IncompTpfa::gravity_
protected
const UnstructuredGrid& Opm::IncompTpfa::grid_
protected
struct ifs_tpfa_data* Opm::IncompTpfa::h_
protected
std::vector<double> Opm::IncompTpfa::htrans_
protected

Referenced by getHalfTrans().

std::vector<double> Opm::IncompTpfa::initial_porevol_
protected
const LinearSolverInterface& Opm::IncompTpfa::linsolver_
protected
const int Opm::IncompTpfa::maxiter_
protected
std::vector<double> Opm::IncompTpfa::omega_
protected
std::vector<double> Opm::IncompTpfa::porevol_
protected
std::vector<double> Opm::IncompTpfa::pressures_
protected
const IncompPropertiesInterface& Opm::IncompTpfa::props_
protected
const double Opm::IncompTpfa::residual_tol_
protected
std::vector<double> Opm::IncompTpfa::rock_comp_
protected
const RockCompressibility* Opm::IncompTpfa::rock_comp_props_
protected
const std::vector<double>& Opm::IncompTpfa::src_
protected
std::vector<double> Opm::IncompTpfa::totmob_
protected
std::vector<double> Opm::IncompTpfa::trans_
protected
std::vector<double> Opm::IncompTpfa::wdp_
protected
const Wells* Opm::IncompTpfa::wells_
protected

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