A model implementation for three-phase black oil with one extra component. More...
#include <BlackoilSolventModel.hpp>
Public Types | |
typedef BlackoilModelBase < Grid, StandardWellsSolvent, BlackoilSolventModel< Grid > > | Base |
typedef Base::ReservoirState | ReservoirState |
typedef Base::WellState | WellState |
Public Member Functions | |
BlackoilSolventModel (const typename Base::ModelParameters ¶m, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const SolventPropsAdFromDeck &solvent_props, const StandardWellsSolvent &well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclState, const bool has_disgas, const bool has_vapoil, const bool terminal_output, const bool has_solvent, const bool is_miscible) | |
Construct the model. | |
void | updateState (const V &dx, ReservoirState &reservoir_state, WellState &well_state) |
Apply an update to the primary variables, chopped if appropriate. | |
std::vector< std::vector < double > > | computeFluidInPlace (const ReservoirState &x, const std::vector< int > &fipnum) |
Protected Types | |
enum | { Solvent = CanonicalVariablePositions::Next } |
typedef Base::SolutionState | SolutionState |
typedef Base::DataBlock | DataBlock |
Protected Member Functions | |
std::vector< ADB > | computeRelPerm (const SolutionState &state) const |
ADB | fluidViscosity (const int phase, const ADB &p, const ADB &temp, const ADB &rs, const ADB &rv, const std::vector< PhasePresence > &cond) const |
ADB | fluidReciprocFVF (const int phase, const ADB &p, const ADB &temp, const ADB &rs, const ADB &rv, const std::vector< PhasePresence > &cond) const |
ADB | fluidDensity (const int phase, const ADB &b, const ADB &rs, const ADB &rv) const |
void | makeConstantState (SolutionState &state) const |
std::vector< V > | variableStateInitials (const ReservoirState &x, const WellState &xw) const |
std::vector< int > | variableStateIndices () const |
SolutionState | variableStateExtractVars (const ReservoirState &x, const std::vector< int > &indices, std::vector< ADB > &vars) const |
void | computeAccum (const SolutionState &state, const int aix) |
void | assembleMassBalanceEq (const SolutionState &state) |
void | addWellContributionToMassBalanceEq (const std::vector< ADB > &cq_s, const SolutionState &state, WellState &xw) |
void | updateEquationsScaling () |
Update the scaling factors for mass balance equations. | |
void | computeMassFlux (const int actph, const V &transi, const ADB &kr, const ADB &mu, const ADB &rho, const ADB &p, const SolutionState &state) |
const std::vector< PhasePresence > | phaseCondition () const |
void | computeEffectiveProperties (const SolutionState &state) |
void | computeToddLongstaffMixing (std::vector< ADB > &viscosity, std::vector< ADB > &density, const std::vector< ADB > &saturations, const ADB po, const Opm::PhaseUsage pu) |
std::vector< ADB > | computePressures (const ADB &po, const ADB &sw, const ADB &so, const ADB &sg, const ADB &ss) const |
Protected Attributes | |
const bool | has_solvent_ |
const int | solvent_pos_ |
const SolventPropsAdFromDeck & | solvent_props_ |
const bool | is_miscible_ |
std::vector< ADB > | mu_eff_ |
std::vector< ADB > | b_eff_ |
Friends | |
class | BlackoilModelBase< Grid, StandardWellsSolvent, BlackoilSolventModel< Grid > > |
A model implementation for three-phase black oil with one extra component.
It uses automatic differentiation via the class AutoDiffBlock to simplify assembly of the jacobian matrix.
Opm::BlackoilSolventModel< Grid >::BlackoilSolventModel | ( | const typename Base::ModelParameters & | param, | |
const Grid & | grid, | |||
const BlackoilPropsAdFromDeck & | fluid, | |||
const DerivedGeology & | geo, | |||
const RockCompressibility * | rock_comp_props, | |||
const SolventPropsAdFromDeck & | solvent_props, | |||
const StandardWellsSolvent & | well_model, | |||
const NewtonIterationBlackoilInterface & | linsolver, | |||
std::shared_ptr< const EclipseState > | eclState, | |||
const bool | has_disgas, | |||
const bool | has_vapoil, | |||
const bool | terminal_output, | |||
const bool | has_solvent, | |||
const bool | is_miscible | |||
) | [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] | fluid | fluid properties |
[in] | geo | rock properties |
[in] | rock_comp_props | if non-null, rock compressibility properties |
[in] | solvent_props | solvent properties |
[in] | wells | well structure |
[in] | linsolver | linear solver |
[in] | has_disgas | turn on dissolved gas |
[in] | has_vapoil | turn on vaporized oil feature |
[in] | terminal_output | request output to cout/cerr |
[in] | has_solvent | turn on solvent feature |
[in] | is_miscible | turn on miscible feature |
void Opm::BlackoilSolventModel< Grid >::updateState | ( | const V & | dx, | |
ReservoirState & | reservoir_state, | |||
WellState & | well_state | |||
) | [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 |