A model implementation for three-phase black oil with polymer. More...
#include <BlackoilPolymerModel.hpp>
Public Types | |
typedef BlackoilModelBase < Grid, StandardWells, BlackoilPolymerModel< Grid > > | Base |
typedef Base::ReservoirState | ReservoirState |
typedef Base::WellState | WellState |
Public Member Functions | |
BlackoilPolymerModel (const typename Base::ModelParameters ¶m, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const PolymerPropsAd &polymer_props_ad, const StandardWells &well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclipseState, const bool has_disgas, const bool has_vapoil, const bool has_polymer, const bool has_plyshlog, const bool has_shrate, const std::vector< double > &wells_rep_radius, const std::vector< double > &wells_perf_length, const std::vector< double > &wells_bore_diameter, const bool terminal_output) | |
Construct the model. | |
void | prepareStep (const SimulatorTimerInterface &timer, const ReservoirState &reservoir_state, const WellState &well_state) |
Called once before each time step. | |
void | afterStep (const SimulatorTimerInterface &timer, ReservoirState &reservoir_state, WellState &well_state) |
Called once after each time step. | |
void | updateState (const V &dx, ReservoirState &reservoir_state, WellState &well_state) |
Apply an update to the primary variables, chopped if appropriate. | |
SimulatorReport | assemble (const ReservoirState &reservoir_state, WellState &well_state, const bool initial_assembly) |
Assemble the residual and Jacobian of the nonlinear system. | |
Protected Types | |
enum | { Concentration = CanonicalVariablePositions::Next } |
typedef Base::SolutionState | SolutionState |
typedef Base::DataBlock | DataBlock |
Protected Member Functions | |
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 | computeInjectionMobility (const SolutionState &state, std::vector< ADB > &mob_perfcells) |
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) |
void | computeCmax (ReservoirState &state) |
ADB | computeMc (const SolutionState &state) const |
const std::vector< PhasePresence > | phaseCondition () const |
void | computeWaterShearVelocityFaces (const V &transi, const std::vector< ADB > &phasePressure, const SolutionState &state, std::vector< double > &water_vel, std::vector< double > &visc_mult) |
Computing the water velocity without shear-thinning for the cell faces. | |
void | computeWaterShearVelocityWells (const SolutionState &state, WellState &xw, const ADB &cq_sw, std::vector< double > &water_vel_wells, std::vector< double > &visc_mult_wells) |
Computing the water velocity without shear-thinning for the well perforations based on the water flux rate. | |
Protected Attributes | |
const PolymerPropsAd & | polymer_props_ad_ |
const bool | has_polymer_ |
const bool | has_plyshlog_ |
const bool | has_shrate_ |
const int | poly_pos_ |
V | cmax_ |
std::vector< double > | wells_rep_radius_ |
std::vector< double > | wells_perf_length_ |
std::vector< double > | wells_bore_diameter_ |
std::vector< double > | shear_mult_faces_ |
std::vector< double > | shear_mult_wells_ |
Friends | |
class | BlackoilModelBase< Grid, StandardWells, BlackoilPolymerModel< Grid > > |
A model implementation for three-phase black oil with polymer.
The simulator is capable of handling three-phase problems where gas can be dissolved in oil and vice versa, with polymer in the water phase. It uses an industry-standard TPFA discretization with per-phase upwind weighting of mobilities.
It uses automatic differentiation via the class AutoDiffBlock to simplify assembly of the jacobian matrix.
Opm::BlackoilPolymerModel< Grid >::BlackoilPolymerModel | ( | const typename Base::ModelParameters & | param, | |
const Grid & | grid, | |||
const BlackoilPropsAdFromDeck & | fluid, | |||
const DerivedGeology & | geo, | |||
const RockCompressibility * | rock_comp_props, | |||
const PolymerPropsAd & | polymer_props_ad, | |||
const StandardWells & | well_model, | |||
const NewtonIterationBlackoilInterface & | linsolver, | |||
std::shared_ptr< const EclipseState > | eclipseState, | |||
const bool | has_disgas, | |||
const bool | has_vapoil, | |||
const bool | has_polymer, | |||
const bool | has_plyshlog, | |||
const bool | has_shrate, | |||
const std::vector< double > & | wells_rep_radius, | |||
const std::vector< double > & | wells_perf_length, | |||
const std::vector< double > & | wells_bore_diameter, | |||
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] | fluid | fluid properties |
[in] | geo | rock properties |
[in] | rock_comp_props | if non-null, rock compressibility 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] | has_polymer | turn on polymer feature |
[in] | has_plyshlog | true when PLYSHLOG keyword available |
[in] | has_shrate | true when PLYSHLOG keyword available |
[in] | wells_rep_radius | representative radius of well perforations during shear effects calculation |
[in] | wells_perf_length | perforation length for well perforations |
[in] | wells_bore_diameter | wellbore diameters for well performations |
[in] | terminal_output | request output to cout/cerr |
void Opm::BlackoilPolymerModel< Grid >::afterStep | ( | const SimulatorTimerInterface & | timer, | |
ReservoirState & | reservoir_state, | |||
WellState & | well_state | |||
) | [inline] |
Called once after each time step.
[in] | timer | simulation timer |
[in,out] | reservoir_state | reservoir state variables |
[in,out] | well_state | well state variables |
SimulatorReport Opm::BlackoilPolymerModel< Grid >::assemble | ( | const ReservoirState & | reservoir_state, | |
WellState & | well_state, | |||
const bool | initial_assembly | |||
) | [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 |
void Opm::BlackoilPolymerModel< Grid >::computeWaterShearVelocityFaces | ( | const V & | transi, | |
const std::vector< ADB > & | phasePressure, | |||
const SolutionState & | state, | |||
std::vector< double > & | water_vel, | |||
std::vector< double > & | visc_mult | |||
) | [inline, protected] |
Computing the water velocity without shear-thinning for the cell faces.
The water velocity will be used for shear-thinning calculation.
void Opm::BlackoilPolymerModel< Grid >::computeWaterShearVelocityWells | ( | const SolutionState & | state, | |
WellState & | xw, | |||
const ADB & | cq_sw, | |||
std::vector< double > & | water_vel_wells, | |||
std::vector< double > & | visc_mult_wells | |||
) | [inline, protected] |
Computing the water velocity without shear-thinning for the well perforations based on the water flux rate.
The water velocity will be used for shear-thinning calculation.
void Opm::BlackoilPolymerModel< Grid >::prepareStep | ( | const SimulatorTimerInterface & | timer, | |
const ReservoirState & | reservoir_state, | |||
const WellState & | well_state | |||
) | [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::BlackoilPolymerModel< 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 |