21 #ifndef OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
22 #define OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
24 #include <opm/autodiff/BlackoilModelBase.hpp>
25 #include <opm/autodiff/BlackoilModelParameters.hpp>
26 #include <opm/polymer/PolymerProperties.hpp>
27 #include <opm/polymer/fullyimplicit/PolymerPropsAd.hpp>
28 #include <opm/polymer/PolymerBlackoilState.hpp>
29 #include <opm/polymer/fullyimplicit/WellStateFullyImplicitBlackoilPolymer.hpp>
30 #include <opm/autodiff/StandardWells.hpp>
31 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
52 typedef typename Base::ReservoirState ReservoirState;
53 typedef typename Base::WellState WellState;
81 const RockCompressibility* rock_comp_props,
83 const StandardWells& well_model,
85 std::shared_ptr< const EclipseState > eclipseState,
86 const bool has_disgas,
87 const bool has_vapoil,
88 const bool has_polymer,
89 const bool has_plyshlog,
90 const bool has_shrate,
91 const std::vector<double>& wells_rep_radius,
92 const std::vector<double>& wells_perf_length,
93 const std::vector<double>& wells_bore_diameter,
94 const bool terminal_output);
101 const ReservoirState& reservoir_state,
102 const WellState& well_state);
109 ReservoirState& reservoir_state,
110 WellState& well_state);
117 ReservoirState& reservoir_state,
118 WellState& well_state);
125 assemble(
const ReservoirState& reservoir_state,
126 WellState& well_state,
127 const bool initial_assembly);
135 typedef typename Base::SolutionState SolutionState;
136 typedef typename Base::DataBlock DataBlock;
137 enum { Concentration = CanonicalVariablePositions::Next };
142 const bool has_polymer_;
143 const bool has_plyshlog_;
144 const bool has_shrate_;
150 std::vector<double> wells_rep_radius_;
151 std::vector<double> wells_perf_length_;
153 std::vector<double> wells_bore_diameter_;
156 std::vector<double> shear_mult_faces_;
158 std::vector<double> shear_mult_wells_;
164 using Base::rock_comp_props_;
165 using Base::linsolver_;
170 using Base::has_disgas_;
171 using Base::has_vapoil_;
173 using Base::use_threshold_pressure_;
174 using Base::threshold_pressures_by_connection_;
176 using Base::phaseCondition_;
177 using Base::residual_;
180 using Base::vfp_properties_;
187 using Base::variableState;
188 using Base::computePressures;
189 using Base::computeGasPressure;
190 using Base::applyThresholdPressures;
191 using Base::fluidViscosity;
192 using Base::fluidReciprocFVF;
193 using Base::fluidDensity;
194 using Base::fluidRsSat;
195 using Base::fluidRvSat;
196 using Base::poroMult;
197 using Base::transMult;
200 using Base::dpMaxRel;
202 using Base::drMaxRel;
203 using Base::maxResidualAllowed;
208 using Base::computeRelPerm;
212 makeConstantState(SolutionState& state)
const;
215 variableStateInitials(
const ReservoirState& x,
216 const WellState& xw)
const;
219 variableStateIndices()
const;
222 variableStateExtractVars(
const ReservoirState& x,
223 const std::vector<int>& indices,
224 std::vector<ADB>& vars)
const;
227 computeAccum(
const SolutionState& state,
231 computeInjectionMobility(
const SolutionState& state,
232 std::vector<ADB>& mob_perfcells);
235 assembleMassBalanceEq(
const SolutionState& state);
238 addWellContributionToMassBalanceEq(
const std::vector<ADB>& cq_s,
239 const SolutionState& state,
242 void updateEquationsScaling();
245 computeMassFlux(
const int actph ,
251 const SolutionState& state );
254 computeCmax(ReservoirState& state);
257 computeMc(
const SolutionState& state)
const;
259 const std::vector<PhasePresence>
260 phaseCondition()
const {
return this->phaseCondition_;}
265 const std::vector<ADB>& phasePressure,
const SolutionState& state,
266 std::vector<double>& water_vel, std::vector<double>& visc_mult);
271 std::vector<double>& water_vel_wells, std::vector<double>& visc_mult_wells);
292 template <
class Gr
id>
303 #include "BlackoilPolymerModel_impl.hpp"
306 #endif // OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
static AutoDiffBlock null()
Construct an empty AutoDiffBlock.
Definition: AutoDiffBlock.hpp:104
Need to include concentration in our state variables, otherwise all is as the default blackoil model...
Definition: BlackoilPolymerModel.hpp:279
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.
Definition: BlackoilPolymerModel_impl.hpp:599
Definition: PolymerPropsAd.hpp:32
SimulatorReport assemble(const ReservoirState &reservoir_state, WellState &well_state, const bool initial_assembly)
Assemble the residual and Jacobian of the nonlinear system.
Definition: BlackoilPolymerModel_impl.hpp:496
Definition: WellStateFullyImplicitBlackoilPolymer.hpp:28
A model implementation for three-phase black oil.
Definition: BlackoilModelBase.hpp:76
const Wells & wells() const
return the Well struct in the WellModel
Definition: BlackoilModelBase.hpp:383
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
void updateState(const V &dx, ReservoirState &reservoir_state, WellState &well_state)
Apply an update to the primary variables, chopped if appropriate.
Definition: BlackoilPolymerModel_impl.hpp:409
void prepareStep(const SimulatorTimerInterface &timer, const ReservoirState &reservoir_state, const WellState &well_state)
Called once before each time step.
Definition: BlackoilPolymerModel_impl.hpp:124
void afterStep(const SimulatorTimerInterface &timer, ReservoirState &reservoir_state, WellState &well_state)
Called once after each time step.
Definition: BlackoilPolymerModel_impl.hpp:141
void updatePrimalVariableFromState(const ReservoirState &state)
update the primal variable for Sg, Rv or Rs.
Simulator state for a compressible two-phase simulator with polymer.
Definition: PolymerBlackoilState.hpp:33
Traits to encapsulate the types used by classes using or extending this model.
Definition: BlackoilModelBase.hpp:60
StandardWells & wellModel()
return the WellModel object
Definition: BlackoilModelBase.hpp:268
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.
Definition: BlackoilPolymerModel_impl.hpp:77
A model implementation for three-phase black oil with polymer.
Definition: BlackoilPolymerModel.hpp:45
bool terminal_output_
Whether we print something to std::cout.
Definition: BlackoilModelBase.hpp:353
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
Struct for containing iteration variables.
Definition: DefaultBlackoilSolutionState.hpp:29
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
void updatePhaseCondFromPrimalVariable(const ReservoirState &state)
Update the phaseCondition_ member based on the primalVariable_ member.
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilModelBase.hpp:386
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...
Definition: BlackoilPolymerModel_impl.hpp:715