A model implementation for three-phase black oil. More...
#include <BlackoilModelEbos.hpp>
Classes | |
class | WellModelMatrixAdapter |
Adapter to turn a matrix into a linear operator. More... | |
Public Types | |
typedef BlackoilState | ReservoirState |
typedef WellStateFullyImplicitBlackoil | WellState |
typedef BlackoilModelParameters | ModelParameters |
typedef double | Scalar |
typedef Dune::FieldVector < Scalar, numEq > | VectorBlockType |
typedef Dune::FieldMatrix < Scalar, numEq, numEq > | MatrixBlockType |
typedef Dune::BCRSMatrix < MatrixBlockType > | Mat |
typedef Dune::BlockVector < VectorBlockType > | BVector |
typedef ISTLSolver < MatrixBlockType, VectorBlockType, Indices::pressureSwitchIdx > | ISTLSolverType |
typedef Opm::FIPData | FIPDataType |
typedef FluidSystems::BlackOil < double > | FluidSystem |
Public Member Functions | |
typedef | GET_PROP_TYPE (TypeTag, Simulator) Simulator |
typedef | GET_PROP_TYPE (TypeTag, Grid) Grid |
typedef | GET_PROP_TYPE (TypeTag, ElementContext) ElementContext |
typedef | GET_PROP_TYPE (TypeTag, SolutionVector) SolutionVector |
typedef | GET_PROP_TYPE (TypeTag, PrimaryVariables) PrimaryVariables |
typedef | GET_PROP_TYPE (TypeTag, FluidSystem) FluidSystem |
typedef | GET_PROP_TYPE (TypeTag, Indices) Indices |
typedef | GET_PROP_TYPE (TypeTag, MaterialLaw) MaterialLaw |
typedef | GET_PROP_TYPE (TypeTag, MaterialLawParams) MaterialLawParams |
BlackoilModelEbos (Simulator &ebosSimulator, const ModelParameters ¶m, BlackoilWellModel< TypeTag > &well_model, RateConverterType &rate_converter, const NewtonIterationBlackoilInterface &linsolver, const bool terminal_output) | |
Construct the model. | |
bool | isParallel () const |
const EclipseState & | eclState () const |
void | prepareStep (const SimulatorTimerInterface &timer, const ReservoirState &, const WellState &) |
Called once before each time step. | |
template<class NonlinearSolverType > | |
SimulatorReport | nonlinearIteration (const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &nonlinear_solver, ReservoirState &, WellState &well_state) |
Called once per nonlinear iteration. | |
void | printIf (int c, double x, double y, double eps, std::string type) |
void | afterStep (const SimulatorTimerInterface &timer, const ReservoirState &reservoir_state, WellState &well_state) |
Called once after each time step. | |
SimulatorReport | assemble (const SimulatorTimerInterface &timer, const int iterationIdx, WellState &well_state) |
Assemble the residual and Jacobian of the nonlinear system. | |
template<class Dummy > | |
double | relativeChange (const Dummy &, const Dummy &) const |
int | sizeNonLinear () const |
The size (number of unknowns) of the nonlinear system of equations. | |
int | linearIterationsLastSolve () const |
Number of linear iterations used in last call to solveJacobianSystem(). | |
void | solveJacobianSystem (BVector &x) const |
Solve the Jacobian system Jx = r where J is the Jacobian and r is the residual. | |
void | updateState (const BVector &dx) |
Apply an update to the primary variables, chopped if appropriate. | |
bool | terminalOutputEnabled () const |
Return true if output to cout is wanted. | |
template<class CollectiveCommunication > | |
double | convergenceReduction (const CollectiveCommunication &comm, const double pvSumLocal, std::vector< Scalar > &R_sum, std::vector< Scalar > &maxCoeff, std::vector< Scalar > &B_avg) |
bool | getConvergence (const SimulatorTimerInterface &timer, const int iteration, std::vector< double > &residual_norms) |
Compute convergence based on total mass balance (tol_mb) and maximum residual mass balance (tol_cnv). | |
int | numPhases () const |
The number of active fluid phases in the model. | |
int | numComponents () const |
template<class T > | |
std::vector< std::vector < double > > | computeFluidInPlace (const T &, const std::vector< int > &fipnum) const |
Wrapper required due to not following generic API. | |
std::vector< std::vector < double > > | computeFluidInPlace (const std::vector< int > &fipnum) const |
SimulationDataContainer | getSimulatorData (const SimulationDataContainer &localState) const |
const FIPDataType & | getFIPData () const |
const Simulator & | ebosSimulator () const |
const SimulatorReport & | failureReport () const |
return the statistics if the nonlinearIteration() method failed | |
BlackoilWellModel< TypeTag > & | wellModel () |
return the StandardWells object | |
const BlackoilWellModel < TypeTag > & | wellModel () const |
int | numWells () const |
bool | localWellsActive () const |
return true if wells are available on this process | |
int | flowPhaseToEbosCompIdx (const int phaseIdx) const |
int | flowPhaseToEbosPhaseIdx (const int phaseIdx) const |
void | beginReportStep () |
void | endReportStep () |
Public Attributes | |
bool | isBeginReportStep_ |
std::vector< bool > | wasSwitched_ |
Static Public Attributes | |
static const int | numEq = Indices::numEq |
static const int | contiSolventEqIdx = Indices::contiSolventEqIdx |
static const int | contiPolymerEqIdx = Indices::contiPolymerEqIdx |
static const int | solventSaturationIdx = Indices::solventSaturationIdx |
static const int | polymerConcentrationIdx = Indices::polymerConcentrationIdx |
Protected Member Functions | |
const ISTLSolverType & | istlSolver () const |
Protected Attributes | |
Simulator & | ebosSimulator_ |
const Grid & | grid_ |
const ISTLSolverType * | istlSolver_ |
const PhaseUsage | phaseUsage_ |
VFPProperties | vfp_properties_ |
const std::vector< bool > | active_ |
const std::vector< int > | cells_ |
const bool | has_disgas_ |
const bool | has_vapoil_ |
const bool | has_solvent_ |
const bool | has_polymer_ |
ModelParameters | param_ |
SimulatorReport | failureReport_ |
BlackoilWellModel< TypeTag > & | well_model_ |
bool | terminal_output_ |
Whether we print something to std::cout. | |
long int | global_nc_ |
The number of cells of the global grid. | |
RateConverterType & | rate_converter_ |
std::vector< std::vector < double > > | residual_norms_history_ |
double | current_relaxation_ |
BVector | dx_old_ |
FIPDataType | fip_ |
A model implementation for three-phase black oil.
The simulator is capable of handling three-phase problems where gas can be dissolved in oil and vice versa. It uses an industry-standard TPFA discretization with per-phase upwind weighting of mobilities.
Opm::BlackoilModelEbos< TypeTag >::BlackoilModelEbos | ( | Simulator & | ebosSimulator, | |
const ModelParameters & | param, | |||
BlackoilWellModel< TypeTag > & | well_model, | |||
RateConverterType & | rate_converter, | |||
const NewtonIterationBlackoilInterface & | linsolver, | |||
const bool | terminal_output | |||
) | [inline] |
Construct the model.
It will retain references to the arguments of this functions, and they are expected to remain in scope for the lifetime of the solver.
[in] | param | parameters |
[in] | grid | grid data structure |
[in] | wells | well structure |
[in] | vfp_properties | Vertical flow performance tables |
[in] | linsolver | linear solver |
[in] | eclState | eclipse state |
[in] | terminal_output | request output to cout/cerr |
void Opm::BlackoilModelEbos< TypeTag >::afterStep | ( | const SimulatorTimerInterface & | timer, | |
const ReservoirState & | reservoir_state, | |||
WellState & | well_state | |||
) | [inline] |
Called once after each time step.
In this class, this function does nothing.
[in] | timer | simulation timer |
[in,out] | reservoir_state | reservoir state variables |
[in,out] | well_state | well state variables |
SimulatorReport Opm::BlackoilModelEbos< TypeTag >::assemble | ( | const SimulatorTimerInterface & | timer, | |
const int | iterationIdx, | |||
WellState & | well_state | |||
) | [inline] |
Assemble the residual and Jacobian of the nonlinear system.
[in] | reservoir_state | reservoir state variables |
[in,out] | well_state | well state variables |
[in] | initial_assembly | pass true if this is the first call to assemble() in this timestep |
bool Opm::BlackoilModelEbos< TypeTag >::getConvergence | ( | const SimulatorTimerInterface & | timer, | |
const int | iteration, | |||
std::vector< double > & | residual_norms | |||
) | [inline] |
Compute convergence based on total mass balance (tol_mb) and maximum residual mass balance (tol_cnv).
[in] | timer | simulation timer |
[in] | dt | timestep length |
[in] | iteration | current iteration number |
SimulatorReport Opm::BlackoilModelEbos< TypeTag >::nonlinearIteration | ( | const int | iteration, | |
const SimulatorTimerInterface & | timer, | |||
NonlinearSolverType & | nonlinear_solver, | |||
ReservoirState & | , | |||
WellState & | well_state | |||
) | [inline] |
Called once per nonlinear iteration.
This model will perform a Newton-Raphson update, changing reservoir_state and well_state. It will also use the nonlinear_solver to do relaxation of updates if necessary.
[in] | iteration | should be 0 for the first call of a new timestep |
[in] | timer | simulation timer |
[in] | nonlinear_solver | nonlinear solver used (for oscillation/relaxation control) |
[in,out] | reservoir_state | reservoir state variables |
[in,out] | well_state | well state variables |
void Opm::BlackoilModelEbos< TypeTag >::prepareStep | ( | const SimulatorTimerInterface & | timer, | |
const ReservoirState & | , | |||
const WellState & | ||||
) | [inline] |
Called once before each time step.
[in] | timer | simulation timer |
[in,out] | reservoir_state | reservoir state variables |
[in,out] | well_state | well state variables |
void Opm::BlackoilModelEbos< TypeTag >::solveJacobianSystem | ( | BVector & | x | ) | const [inline] |
Solve the Jacobian system Jx = r where J is the Jacobian and r is the residual.
void Opm::BlackoilModelEbos< TypeTag >::updateState | ( | const BVector & | dx | ) | [inline] |
Apply an update to the primary variables, chopped if appropriate.
[in] | dx | updates to apply to primary variables |
[in,out] | reservoir_state | reservoir state variables |
[in,out] | well_state | well state variables |