00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #ifndef OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
00022 #define OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
00023
00024 #include <opm/autodiff/BlackoilModelBase.hpp>
00025 #include <opm/autodiff/BlackoilModelParameters.hpp>
00026 #include <opm/polymer/PolymerProperties.hpp>
00027 #include <opm/polymer/fullyimplicit/PolymerPropsAd.hpp>
00028 #include <opm/polymer/PolymerBlackoilState.hpp>
00029 #include <opm/polymer/fullyimplicit/WellStateFullyImplicitBlackoilPolymer.hpp>
00030 #include <opm/autodiff/StandardWells.hpp>
00031 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
00032
00033 namespace Opm {
00034
00044 template<class Grid>
00045 class BlackoilPolymerModel : public BlackoilModelBase<Grid, StandardWells, BlackoilPolymerModel<Grid> >
00046 {
00047 public:
00048
00049
00050
00051 typedef BlackoilModelBase<Grid, StandardWells, BlackoilPolymerModel<Grid> > Base;
00052 typedef typename Base::ReservoirState ReservoirState;
00053 typedef typename Base::WellState WellState;
00054
00055
00056 friend class BlackoilModelBase<Grid, StandardWells, BlackoilPolymerModel<Grid> >;
00057
00077 BlackoilPolymerModel(const typename Base::ModelParameters& param,
00078 const Grid& grid,
00079 const BlackoilPropsAdFromDeck& fluid,
00080 const DerivedGeology& geo,
00081 const RockCompressibility* rock_comp_props,
00082 const PolymerPropsAd& polymer_props_ad,
00083 const StandardWells& well_model,
00084 const NewtonIterationBlackoilInterface& linsolver,
00085 std::shared_ptr< const EclipseState > eclipseState,
00086 const bool has_disgas,
00087 const bool has_vapoil,
00088 const bool has_polymer,
00089 const bool has_plyshlog,
00090 const bool has_shrate,
00091 const std::vector<double>& wells_rep_radius,
00092 const std::vector<double>& wells_perf_length,
00093 const std::vector<double>& wells_bore_diameter,
00094 const bool terminal_output);
00095
00100 void prepareStep(const SimulatorTimerInterface& timer,
00101 const ReservoirState& reservoir_state,
00102 const WellState& well_state);
00103
00108 void afterStep(const SimulatorTimerInterface& timer,
00109 ReservoirState& reservoir_state,
00110 WellState& well_state);
00111
00116 void updateState(const V& dx,
00117 ReservoirState& reservoir_state,
00118 WellState& well_state);
00119
00124 SimulatorReport
00125 assemble(const ReservoirState& reservoir_state,
00126 WellState& well_state,
00127 const bool initial_assembly);
00128
00129 using Base::wellModel;
00130
00131 protected:
00132
00133
00134
00135 typedef typename Base::SolutionState SolutionState;
00136 typedef typename Base::DataBlock DataBlock;
00137 enum { Concentration = CanonicalVariablePositions::Next };
00138
00139
00140
00141 const PolymerPropsAd& polymer_props_ad_;
00142 const bool has_polymer_;
00143 const bool has_plyshlog_;
00144 const bool has_shrate_;
00145 const int poly_pos_;
00146 V cmax_;
00147
00148
00149
00150 std::vector<double> wells_rep_radius_;
00151 std::vector<double> wells_perf_length_;
00152
00153 std::vector<double> wells_bore_diameter_;
00154
00155
00156 std::vector<double> shear_mult_faces_;
00157
00158 std::vector<double> shear_mult_wells_;
00159
00160
00161 using Base::grid_;
00162 using Base::fluid_;
00163 using Base::geo_;
00164 using Base::rock_comp_props_;
00165 using Base::linsolver_;
00166 using Base::active_;
00167 using Base::canph_;
00168 using Base::cells_;
00169 using Base::ops_;
00170 using Base::has_disgas_;
00171 using Base::has_vapoil_;
00172 using Base::param_;
00173 using Base::use_threshold_pressure_;
00174 using Base::threshold_pressures_by_connection_;
00175 using Base::sd_;
00176 using Base::phaseCondition_;
00177 using Base::residual_;
00178 using Base::terminal_output_;
00179 using Base::pvdt_;
00180 using Base::vfp_properties_;
00181
00182
00183
00184
00185 using Base::wells;
00186 using Base::wellsActive;
00187 using Base::variableState;
00188 using Base::computePressures;
00189 using Base::computeGasPressure;
00190 using Base::applyThresholdPressures;
00191 using Base::fluidViscosity;
00192 using Base::fluidReciprocFVF;
00193 using Base::fluidDensity;
00194 using Base::fluidRsSat;
00195 using Base::fluidRvSat;
00196 using Base::poroMult;
00197 using Base::transMult;
00198 using Base::updatePrimalVariableFromState;
00199 using Base::updatePhaseCondFromPrimalVariable;
00200 using Base::dpMaxRel;
00201 using Base::dsMax;
00202 using Base::drMaxRel;
00203 using Base::maxResidualAllowed;
00204
00205
00206
00207
00208 using Base::computeRelPerm;
00209
00210
00211 void
00212 makeConstantState(SolutionState& state) const;
00213
00214 std::vector<V>
00215 variableStateInitials(const ReservoirState& x,
00216 const WellState& xw) const;
00217
00218 std::vector<int>
00219 variableStateIndices() const;
00220
00221 SolutionState
00222 variableStateExtractVars(const ReservoirState& x,
00223 const std::vector<int>& indices,
00224 std::vector<ADB>& vars) const;
00225
00226 void
00227 computeAccum(const SolutionState& state,
00228 const int aix );
00229
00230 void
00231 computeInjectionMobility(const SolutionState& state,
00232 std::vector<ADB>& mob_perfcells);
00233
00234 void
00235 assembleMassBalanceEq(const SolutionState& state);
00236
00237 void
00238 addWellContributionToMassBalanceEq(const std::vector<ADB>& cq_s,
00239 const SolutionState& state,
00240 WellState& xw);
00241
00242 void updateEquationsScaling();
00243
00244 void
00245 computeMassFlux(const int actph ,
00246 const V& transi,
00247 const ADB& kr ,
00248 const ADB& mu ,
00249 const ADB& rho ,
00250 const ADB& p ,
00251 const SolutionState& state );
00252
00253 void
00254 computeCmax(ReservoirState& state);
00255
00256 ADB
00257 computeMc(const SolutionState& state) const;
00258
00259 const std::vector<PhasePresence>
00260 phaseCondition() const {return this->phaseCondition_;}
00261
00264 void computeWaterShearVelocityFaces(const V& transi,
00265 const std::vector<ADB>& phasePressure, const SolutionState& state,
00266 std::vector<double>& water_vel, std::vector<double>& visc_mult);
00267
00270 void computeWaterShearVelocityWells(const SolutionState& state, WellState& xw, const ADB& cq_sw,
00271 std::vector<double>& water_vel_wells, std::vector<double>& visc_mult_wells);
00272
00273 };
00274
00275
00276
00279 struct BlackoilPolymerSolutionState : public DefaultBlackoilSolutionState
00280 {
00281 explicit BlackoilPolymerSolutionState(const int np)
00282 : DefaultBlackoilSolutionState(np),
00283 concentration( ADB::null())
00284 {
00285 }
00286 ADB concentration;
00287 };
00288
00289
00290
00292 template <class Grid>
00293 struct ModelTraits< BlackoilPolymerModel<Grid> >
00294 {
00295 typedef PolymerBlackoilState ReservoirState;
00296 typedef WellStateFullyImplicitBlackoilPolymer WellState;
00297 typedef BlackoilModelParameters ModelParameters;
00298 typedef BlackoilPolymerSolutionState SolutionState;
00299 };
00300
00301 }
00302
00303 #include "BlackoilPolymerModel_impl.hpp"
00304
00305
00306 #endif // OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED