28#ifndef EWOMS_BLACK_OIL_POLYMER_MODULE_HH
29#define EWOMS_BLACK_OIL_POLYMER_MODULE_HH
31#include <dune/common/fvector.hh>
33#include <opm/common/OpmLog/OpmLog.hpp>
35#include <opm/material/common/MathToolbox.hpp>
62template <class TypeTag, bool enablePolymerV = getPropValue<TypeTag, Properties::EnablePolymer>()>
78 using Toolbox = MathToolbox<Evaluation>;
83 static constexpr unsigned polymerConcentrationIdx = Indices::polymerConcentrationIdx;
84 static constexpr unsigned polymerMoleWeightIdx = Indices::polymerMoleWeightIdx;
85 static constexpr unsigned contiPolymerEqIdx = Indices::contiPolymerEqIdx;
86 static constexpr unsigned contiPolymerMolarWeightEqIdx = Indices::contiPolymerMWEqIdx;
87 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
89 static constexpr unsigned enablePolymer = enablePolymerV;
90 static constexpr bool enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>();
92 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
93 static constexpr unsigned numPhases = FluidSystem::numPhases;
107 const auto iterTable = params_.plymwinjTables_.find(tableNumber);
108 if (iterTable != params_.plymwinjTables_.end()) {
109 return iterTable->second;
112 throw std::runtime_error(
" the PLYMWINJ table " +
std::to_string(tableNumber) +
" does not exist\n");
121 const auto iterTable = params_.skprwatTables_.find(tableNumber);
122 if (iterTable != params_.skprwatTables_.end()) {
123 return iterTable->second;
126 throw std::runtime_error(
" the SKPRWAT table " +
std::to_string(tableNumber) +
" does not exist\n");
136 const auto iterTable = params_.skprpolyTables_.find(tableNumber);
137 if (iterTable != params_.skprpolyTables_.end()) {
138 return iterTable->second;
141 throw std::runtime_error(
" the SKPRPOLY table " +
std::to_string(tableNumber) +
" does not exist\n");
150 if constexpr (enablePolymer) {
159 Simulator& simulator)
161 if constexpr (enablePolymer) {
168 if constexpr (enablePolymer) {
169 if constexpr (enablePolymerMolarWeight) {
170 return pvIdx == polymerConcentrationIdx || pvIdx == polymerMoleWeightIdx;
173 return pvIdx == polymerConcentrationIdx;
185 if (pvIdx == polymerConcentrationIdx) {
186 return "polymer_waterconcentration";
189 return "polymer_molecularweight";
198 return static_cast<Scalar
>(1.0);
203 if constexpr (enablePolymer) {
204 if constexpr (enablePolymerMolarWeight) {
205 return eqIdx == contiPolymerEqIdx || eqIdx == contiPolymerMolarWeightEqIdx;
208 return eqIdx == contiPolymerEqIdx;
216 static std::string
eqName(
unsigned eqIdx)
220 if (eqIdx == contiPolymerEqIdx) {
221 return "conti^polymer";
224 return "conti^polymer_molecularweight";
228 static Scalar
eqWeight([[maybe_unused]]
unsigned eqIdx)
233 return static_cast<Scalar
>(1.0);
237 template <
class LhsEval>
238 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
239 const IntensiveQuantities& intQuants)
241 if constexpr (enablePolymer) {
242 const auto& fs = intQuants.fluidState();
245 const LhsEval surfaceVolumeWater =
246 max(Toolbox::template decay<LhsEval>(fs.saturation(waterPhaseIdx)) *
247 Toolbox::template decay<LhsEval>(fs.invB(waterPhaseIdx)) *
248 Toolbox::template decay<LhsEval>(intQuants.porosity()),
252 const LhsEval massPolymer =
254 Toolbox::template decay<LhsEval>(intQuants.polymerConcentration()) *
255 (1.0 - Toolbox::template decay<LhsEval>(intQuants.polymerDeadPoreVolume()));
258 const LhsEval adsorptionPolymer =
259 Toolbox::template decay<LhsEval>(1.0 - intQuants.porosity()) *
260 Toolbox::template decay<LhsEval>(intQuants.polymerRockDensity()) *
261 Toolbox::template decay<LhsEval>(intQuants.polymerAdsorption());
263 LhsEval accumulationPolymer = massPolymer + adsorptionPolymer;
265 storage[contiPolymerEqIdx] += accumulationPolymer;
268 if constexpr (enablePolymerMolarWeight) {
269 accumulationPolymer = max(accumulationPolymer, 1e-10);
271 storage[contiPolymerMolarWeightEqIdx] +=
272 accumulationPolymer * Toolbox::template decay<LhsEval>(intQuants.polymerMoleWeight());
278 [[maybe_unused]]
const ElementContext& elemCtx,
279 [[maybe_unused]]
unsigned scvfIdx,
280 [[maybe_unused]]
unsigned timeIdx)
282 if constexpr (enablePolymer) {
283 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
285 const unsigned upIdx = extQuants.upstreamIndex(FluidSystem::waterPhaseIdx);
286 const unsigned inIdx = extQuants.interiorIndex();
287 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
288 const unsigned contiWaterEqIdx =
289 Indices::conti0EqIdx + Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
291 if (upIdx == inIdx) {
292 flux[contiPolymerEqIdx] =
293 extQuants.volumeFlux(waterPhaseIdx) *
294 up.fluidState().invB(waterPhaseIdx) *
295 up.polymerViscosityCorrection() /
296 extQuants.polymerShearFactor() *
297 up.polymerConcentration();
300 flux[contiWaterEqIdx] /= extQuants.waterShearFactor();
303 flux[contiPolymerEqIdx] =
304 extQuants.volumeFlux(waterPhaseIdx) *
305 decay<Scalar>(up.fluidState().invB(waterPhaseIdx)) *
306 decay<Scalar>(up.polymerViscosityCorrection()) /
307 decay<Scalar>(extQuants.polymerShearFactor()) *
308 decay<Scalar>(up.polymerConcentration());
311 flux[contiWaterEqIdx] /= decay<Scalar>(extQuants.waterShearFactor());
315 if constexpr (enablePolymerMolarWeight) {
316 if (upIdx == inIdx) {
317 flux[contiPolymerMolarWeightEqIdx] =
318 flux[contiPolymerEqIdx] * up.polymerMoleWeight();
321 flux[contiPolymerMolarWeightEqIdx] =
322 flux[contiPolymerEqIdx] * decay<Scalar>(up.polymerMoleWeight());
337 return static_cast<Scalar
>(0.0);
340 template <
class DofEntity>
341 static void serializeEntity(
const Model& model, std::ostream& outstream,
const DofEntity& dof)
343 if constexpr (enablePolymer) {
344 const unsigned dofIdx = model.dofMapper().index(dof);
345 const PrimaryVariables& priVars = model.solution(0)[dofIdx];
346 outstream << priVars[polymerConcentrationIdx];
347 outstream << priVars[polymerMoleWeightIdx];
351 template <
class DofEntity>
354 if constexpr (enablePolymer) {
355 const unsigned dofIdx = model.dofMapper().index(dof);
356 PrimaryVariables& priVars0 = model.solution(0)[dofIdx];
357 PrimaryVariables& priVars1 = model.solution(1)[dofIdx];
359 instream >> priVars0[polymerConcentrationIdx];
360 instream >> priVars0[polymerMoleWeightIdx];
363 priVars1[polymerConcentrationIdx] = priVars0[polymerConcentrationIdx];
364 priVars1[polymerMoleWeightIdx] = priVars0[polymerMoleWeightIdx];
372 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
373 return params_.plyrockDeadPoreVolume_[satnumRegionIdx];
380 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
381 return params_.plyrockResidualResistanceFactor_[satnumRegionIdx];
388 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
389 return params_.plyrockRockDensityFactor_[satnumRegionIdx];
396 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
397 return params_.plyrockAdsorbtionIndex_[satnumRegionIdx];
404 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
405 return params_.plyrockMaxAdsorbtion_[satnumRegionIdx];
408 static const TabulatedFunction&
413 const unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
414 return params_.plyadsAdsorbedPolymer_[satnumRegionIdx];
417 static const TabulatedFunction&
422 unsigned pvtnumRegionIdx =
423 elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
424 return params_.plyviscViscosityMultiplierTable_[pvtnumRegionIdx];
427 static const TabulatedFunction&
429 {
return params_.plyviscViscosityMultiplierTable_[pvtnumRegionIdx]; }
435 const unsigned polymerMixRegionIdx =
436 elemCtx.problem().plmixnumRegionIndex(elemCtx, scvIdx, timeIdx);
437 return params_.plymaxMaxConcentration_[polymerMixRegionIdx];
444 const unsigned polymerMixRegionIdx =
445 elemCtx.problem().plmixnumRegionIndex(elemCtx, scvIdx, timeIdx);
446 return params_.plymixparToddLongstaff_[polymerMixRegionIdx];
451 const unsigned scvIdx,
452 const unsigned timeIdx)
454 const unsigned polymerMixRegionIdx =
455 elemCtx.problem().plmixnumRegionIndex(elemCtx, scvIdx, timeIdx);
456 return params_.plyvmhCoefficients_[polymerMixRegionIdx];
460 {
return params_.hasPlyshlog_; }
463 {
return params_.hasShrate_; }
465 static Scalar
shrate(
unsigned pvtnumRegionIdx)
466 {
return params_.shrate_[pvtnumRegionIdx]; }
474 template <
class Evaluation>
476 unsigned pvtnumRegionIdx,
477 const Evaluation& v0)
479 using ToolboxLocal = MathToolbox<Evaluation>;
481 const auto& viscosityMultiplierTable = params_.plyviscViscosityMultiplierTable_[pvtnumRegionIdx];
482 const Scalar viscosityMultiplier =
483 viscosityMultiplierTable.eval(scalarValue(polymerConcentration),
true);
485 const Scalar eps = 1e-14;
487 if (std::abs((viscosityMultiplier - 1.0)) < eps) {
488 return ToolboxLocal::createConstant(v0, 1.0);
491 const std::vector<Scalar>& shearEffectRefLogVelocity =
492 params_.plyshlogShearEffectRefLogVelocity_[pvtnumRegionIdx];
493 const auto v0AbsLog = log(abs(v0));
495 if (v0AbsLog < shearEffectRefLogVelocity[0]) {
496 return ToolboxLocal::createConstant(v0, 1.0);
503 const std::vector<Scalar>& shearEffectRefMultiplier =
504 params_.plyshlogShearEffectRefMultiplier_[pvtnumRegionIdx];
505 const std::size_t numTableEntries = shearEffectRefLogVelocity.size();
506 assert(shearEffectRefMultiplier.size() == numTableEntries);
508 std::vector<Scalar> shearEffectMultiplier(numTableEntries, 1.0);
509 for (std::size_t i = 0; i < numTableEntries; ++i) {
510 shearEffectMultiplier[i] = (1.0 + (viscosityMultiplier - 1.0) *
511 shearEffectRefMultiplier[i]) / viscosityMultiplier;
512 shearEffectMultiplier[i] = log(shearEffectMultiplier[i]);
516 const TabulatedFunction logShearEffectMultiplier =
517 TabulatedFunction(numTableEntries, shearEffectRefLogVelocity,
518 shearEffectMultiplier,
false);
525 auto F = [&logShearEffectMultiplier, &v0AbsLog](
const Evaluation& u) {
526 return u + logShearEffectMultiplier.eval(u,
true) - v0AbsLog;
529 auto dF = [&logShearEffectMultiplier](
const Evaluation& u) {
530 return 1 + logShearEffectMultiplier.evalDerivative(u,
true);
536 bool converged =
false;
538 for (
int i = 0; i < 20; ++i) {
540 const auto df = dF(u);
542 if (std::abs(scalarValue(f)) < 1e-12) {
548 throw std::runtime_error(
"Not able to compute shear velocity. \n");
552 return exp(logShearEffectMultiplier.eval(u,
true));
562template <
class TypeTag,
bool enablePolymerV>
563BlackOilPolymerParams<typename BlackOilPolymerModule<TypeTag, enablePolymerV>::Scalar>
564BlackOilPolymerModule<TypeTag, enablePolymerV>::params_;
566template <
class TypeTag,
bool enablePolymerV>
576template <
class TypeTag>
591 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
592 static constexpr int polymerConcentrationIdx = Indices::polymerConcentrationIdx;
593 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
594 static constexpr bool enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>();
595 static constexpr int polymerMoleWeightIdx = Indices::polymerMoleWeightIdx;
607 const auto linearizationType = elemCtx.linearizationType();
608 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
609 polymerConcentration_ = priVars.makeEvaluation(polymerConcentrationIdx, timeIdx, linearizationType);
610 if constexpr (enablePolymerMolarWeight) {
611 polymerMoleWeight_ = priVars.makeEvaluation(polymerMoleWeightIdx, timeIdx, linearizationType);
615 const Scalar& maxAdsorbtion = PolymerModule::plyrockMaxAdsorbtion(elemCtx, dofIdx, timeIdx);
616 const auto& plyadsAdsorbedPolymer = PolymerModule::plyadsAdsorbedPolymer(elemCtx, dofIdx, timeIdx);
617 polymerAdsorption_ = plyadsAdsorbedPolymer.eval(polymerConcentration_,
true);
618 if (
static_cast<int>(PolymerModule::plyrockAdsorbtionIndex(elemCtx, dofIdx, timeIdx)) ==
621 const auto maxPolymerAdsorption =
622 elemCtx.problem().maxPolymerAdsorption(elemCtx, dofIdx, timeIdx);
623 polymerAdsorption_ = std::max(Evaluation(maxPolymerAdsorption), polymerAdsorption_);
627 const Scalar& residualResistanceFactor =
628 PolymerModule::plyrockResidualResistanceFactor(elemCtx, dofIdx, timeIdx);
629 const Evaluation resistanceFactor = 1.0 + (residualResistanceFactor - 1.0) *
630 polymerAdsorption_ / maxAdsorbtion;
633 if constexpr (!enablePolymerMolarWeight) {
634 const Scalar cmax = PolymerModule::plymaxMaxConcentration(elemCtx, dofIdx, timeIdx);
635 const auto& fs = asImp_().fluidState_;
636 const Evaluation& muWater = fs.viscosity(waterPhaseIdx);
637 const auto& viscosityMultiplier =
638 PolymerModule::plyviscViscosityMultiplierTable(elemCtx, dofIdx, timeIdx);
639 const Evaluation viscosityMixture =
640 viscosityMultiplier.eval(polymerConcentration_,
true) * muWater;
643 const Scalar plymixparToddLongstaff = PolymerModule::plymixparToddLongstaff(elemCtx, dofIdx, timeIdx);
644 const Evaluation viscosityPolymer = viscosityMultiplier.eval(cmax,
true) * muWater;
645 const Evaluation viscosityPolymerEffective =
646 pow(viscosityMixture, plymixparToddLongstaff) * pow(viscosityPolymer, 1.0 - plymixparToddLongstaff);
647 const Evaluation viscosityWaterEffective =
648 pow(viscosityMixture, plymixparToddLongstaff) * pow(muWater, 1.0 - plymixparToddLongstaff);
650 const Evaluation cbar = polymerConcentration_ / cmax;
652 waterViscosityCorrection_ = muWater * ((1.0 - cbar) / viscosityWaterEffective + cbar / viscosityPolymerEffective);
654 polymerViscosityCorrection_ = (muWater / waterViscosityCorrection_) / viscosityPolymerEffective;
657 const auto& plyvmhCoefficients = PolymerModule::plyvmhCoefficients(elemCtx, dofIdx, timeIdx);
658 const Scalar k_mh = plyvmhCoefficients.k_mh;
659 const Scalar a_mh = plyvmhCoefficients.a_mh;
660 const Scalar gamma = plyvmhCoefficients.gamma;
661 const Scalar kappa = plyvmhCoefficients.kappa;
665 const Evaluation intrinsicViscosity = k_mh * pow(polymerMoleWeight_ * 1000., a_mh);
666 const Evaluation x = polymerConcentration_ * intrinsicViscosity;
667 waterViscosityCorrection_ = 1.0 / (1.0 + gamma * (x + kappa * x * x));
668 polymerViscosityCorrection_ = 1.0;
672 asImp_().mobility_[waterPhaseIdx] *= waterViscosityCorrection_ / resistanceFactor;
675 polymerDeadPoreVolume_ = PolymerModule::plyrockDeadPoreVolume(elemCtx, dofIdx, timeIdx);
676 polymerRockDensity_ = PolymerModule::plyrockRockDensityFactor(elemCtx, dofIdx, timeIdx);
680 {
return polymerConcentration_; }
684 if constexpr (enablePolymerMolarWeight) {
685 return polymerMoleWeight_;
688 throw std::logic_error(
"polymerMoleWeight() is called but polymer milecular weight is disabled");
693 {
return polymerDeadPoreVolume_; }
696 {
return polymerAdsorption_; }
699 {
return polymerRockDensity_; }
703 {
return polymerViscosityCorrection_; }
707 {
return waterViscosityCorrection_; }
711 {
return *
static_cast<Implementation*
>(
this); }
723template <
class TypeTag>
737 {
throw std::logic_error(
"polymerMoleWeight() called but polymer molecular weight is disabled"); }
740 {
throw std::runtime_error(
"polymerConcentration() called but polymers are disabled"); }
743 {
throw std::runtime_error(
"polymerDeadPoreVolume() called but polymers are disabled"); }
746 {
throw std::runtime_error(
"polymerAdsorption() called but polymers are disabled"); }
749 {
throw std::runtime_error(
"polymerRockDensity() called but polymers are disabled"); }
752 {
throw std::runtime_error(
"polymerViscosityCorrection() called but polymers are disabled"); }
755 {
throw std::runtime_error(
"waterViscosityCorrection() called but polymers are disabled"); }
758template <
class TypeTag,
bool enablePolymerV>
768template <
class TypeTag>
781 static constexpr unsigned gasPhaseIdx = FluidSystem::gasPhaseIdx;
782 static constexpr int dimWorld = GridView::dimensionworld;
783 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
785 using Toolbox = MathToolbox<Evaluation>;
787 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
788 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
797 template <
class Dummy =
bool>
803 throw std::runtime_error(
"The extension of the blackoil model for polymers is not yet "
804 "implemented for problems specified using permeabilities.");
813 template <
class Dummy =
bool>
819 waterShearFactor_ = 1.0;
820 polymerShearFactor_ = 1.0;
822 if (!PolymerModule::hasPlyshlog()) {
826 const ExtensiveQuantities& extQuants = asImp_();
827 const unsigned upIdx = extQuants.upstreamIndex(waterPhaseIdx);
828 const unsigned interiorDofIdx = extQuants.interiorIndex();
829 const unsigned exteriorDofIdx = extQuants.exteriorIndex();
830 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
831 const auto& intQuantsIn = elemCtx.intensiveQuantities(interiorDofIdx, timeIdx);
832 const auto& intQuantsEx = elemCtx.intensiveQuantities(exteriorDofIdx, timeIdx);
835 const Evaluation poroAvg = intQuantsIn.porosity() * 0.5 + intQuantsEx.porosity() * 0.5;
836 const unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvfIdx, timeIdx);
837 const Evaluation& Sw = up.fluidState().saturation(waterPhaseIdx);
838 const unsigned cellIdx = elemCtx.globalSpaceIndex(scvfIdx, timeIdx);
839 const auto& materialLawManager = elemCtx.problem().materialLawManager();
840 const auto& scaledDrainageInfo =
841 materialLawManager->oilWaterScaledEpsInfoDrainage(cellIdx);
842 const Scalar& Swcr = scaledDrainageInfo.Swcr;
845 const Evaluation denom = max(poroAvg * (Sw - Swcr), 1e-12);
846 Evaluation waterVolumeVelocity = extQuants.volumeFlux(waterPhaseIdx) / denom;
849 if (PolymerModule::hasShrate()) {
850 const Evaluation& relWater = up.relativePermeability(waterPhaseIdx);
851 const Scalar trans = elemCtx.problem().transmissibility(elemCtx, interiorDofIdx, exteriorDofIdx);
853 const Scalar faceArea = elemCtx.stencil(timeIdx).interiorFace(scvfIdx).area();
854 const auto dist = elemCtx.pos(interiorDofIdx, timeIdx) - elemCtx.pos(exteriorDofIdx, timeIdx);
856 const Scalar absPerm = trans / faceArea * dist.two_norm();
857 waterVolumeVelocity *=
858 PolymerModule::shrate(pvtnumRegionIdx) * sqrt(poroAvg * Sw / (relWater * absPerm));
859 assert(isfinite(waterVolumeVelocity));
865 PolymerModule::computeShearFactor(up.polymerConcentration(),
867 waterVolumeVelocity);
868 polymerShearFactor_ =
869 PolymerModule::computeShearFactor(up.polymerConcentration(),
871 waterVolumeVelocity * up.polymerViscosityCorrection());
875 {
return polymerShearFactor_; }
878 {
return waterShearFactor_; }
881 Implementation& asImp_()
882 {
return *
static_cast<Implementation*
>(
this); }
884 Evaluation polymerShearFactor_;
885 Evaluation waterShearFactor_;
888template <
class TypeTag>
906 {
throw std::runtime_error(
"polymerShearFactor() called but polymers are disabled"); }
909 {
throw std::runtime_error(
"waterShearFactor() called but polymers are disabled"); }
Contains the parameters required to extend the black-oil model by polymer.
Declares the properties required by the black oil model.
const Evaluation & polymerShearFactor() const
Definition: blackoilpolymermodules.hh:905
void updateShearMultipliersPerm(const ElementContext &, unsigned, unsigned)
Definition: blackoilpolymermodules.hh:900
void updateShearMultipliers(const ElementContext &, unsigned, unsigned)
Definition: blackoilpolymermodules.hh:895
const Evaluation & waterShearFactor() const
Definition: blackoilpolymermodules.hh:908
const Evaluation & polymerShearFactor() const
Definition: blackoilpolymermodules.hh:874
void updateShearMultipliersPerm(const ElementContext &, unsigned, unsigned)
Method which calculates the shear factor based on flow velocity.
Definition: blackoilpolymermodules.hh:799
const Evaluation & waterShearFactor() const
Definition: blackoilpolymermodules.hh:877
void updateShearMultipliers(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Method which calculates the shear factor based on flow velocity.
Definition: blackoilpolymermodules.hh:815
Provides the polymer specific extensive quantities to the generic black-oil module's extensive quanti...
Definition: blackoilpolymermodules.hh:759
const Evaluation & polymerConcentration() const
Definition: blackoilpolymermodules.hh:739
const Evaluation & waterViscosityCorrection() const
Definition: blackoilpolymermodules.hh:754
const Evaluation & polymerRockDensity() const
Definition: blackoilpolymermodules.hh:748
const Evaluation & polymerDeadPoreVolume() const
Definition: blackoilpolymermodules.hh:742
const Evaluation & polymerMoleWeight() const
Definition: blackoilpolymermodules.hh:736
void polymerPropertiesUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilpolymermodules.hh:731
const Evaluation & polymerAdsorption() const
Definition: blackoilpolymermodules.hh:745
const Evaluation & polymerViscosityCorrection() const
Definition: blackoilpolymermodules.hh:751
Evaluation polymerViscosityCorrection_
Definition: blackoilpolymermodules.hh:719
Scalar polymerDeadPoreVolume_
Definition: blackoilpolymermodules.hh:716
void polymerPropertiesUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the intensive properties needed to handle polymers from the primary variables.
Definition: blackoilpolymermodules.hh:603
Implementation & asImp_()
Definition: blackoilpolymermodules.hh:710
Scalar polymerRockDensity_
Definition: blackoilpolymermodules.hh:717
const Evaluation & polymerMoleWeight() const
Definition: blackoilpolymermodules.hh:682
Evaluation polymerConcentration_
Definition: blackoilpolymermodules.hh:713
Scalar polymerDeadPoreVolume() const
Definition: blackoilpolymermodules.hh:692
Evaluation polymerMoleWeight_
Definition: blackoilpolymermodules.hh:715
const Evaluation & waterViscosityCorrection() const
Definition: blackoilpolymermodules.hh:706
const Evaluation & polymerAdsorption() const
Definition: blackoilpolymermodules.hh:695
const Evaluation & polymerConcentration() const
Definition: blackoilpolymermodules.hh:679
const Evaluation & polymerViscosityCorrection() const
Definition: blackoilpolymermodules.hh:702
Evaluation waterViscosityCorrection_
Definition: blackoilpolymermodules.hh:720
Scalar polymerRockDensity() const
Definition: blackoilpolymermodules.hh:698
Evaluation polymerAdsorption_
Definition: blackoilpolymermodules.hh:718
Provides the volumetric quantities required for the equations needed by the polymers extension of the...
Definition: blackoilpolymermodules.hh:567
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
static const TabulatedFunction & plyviscViscosityMultiplierTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:418
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilpolymermodules.hh:181
static const TabulatedFunction & plyadsAdsorbedPolymer(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:409
static TabulatedTwoDFunction & getPlymwinjTable(const int tableNumber)
get the PLYMWINJ table
Definition: blackoilpolymermodules.hh:105
static Scalar plyrockAdsorbtionIndex(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:392
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilpolymermodules.hh:341
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilpolymermodules.hh:228
static void registerParameters()
Register all run-time parameters for the black-oil polymer module.
Definition: blackoilpolymermodules.hh:148
static bool eqApplies(unsigned eqIdx)
Definition: blackoilpolymermodules.hh:201
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilpolymermodules.hh:193
static BlackOilPolymerParams< Scalar >::SkprpolyTable & getSkprpolyTable(const int tableNumber)
get the SKPRPOLY table
Definition: blackoilpolymermodules.hh:134
static bool hasShrate()
Definition: blackoilpolymermodules.hh:462
static Scalar plyrockRockDensityFactor(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:384
static Scalar plyrockDeadPoreVolume(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:368
static const TabulatedFunction & plyviscViscosityMultiplierTable(unsigned pvtnumRegionIdx)
Definition: blackoilpolymermodules.hh:428
static void registerOutputModules(Model &model, Simulator &simulator)
Register all polymer specific VTK and ECL output modules.
Definition: blackoilpolymermodules.hh:158
static Scalar plymaxMaxConcentration(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:431
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilpolymermodules.hh:331
Scalar molarMass() const
Definition: blackoilpolymermodules.hh:555
static const BlackOilPolymerParams< Scalar >::PlyvmhCoefficients & plyvmhCoefficients(const ElementContext &elemCtx, const unsigned scvIdx, const unsigned timeIdx)
Definition: blackoilpolymermodules.hh:450
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilpolymermodules.hh:166
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilpolymermodules.hh:238
static bool hasPlyshlog()
Definition: blackoilpolymermodules.hh:459
static Scalar shrate(unsigned pvtnumRegionIdx)
Definition: blackoilpolymermodules.hh:465
static std::string eqName(unsigned eqIdx)
Definition: blackoilpolymermodules.hh:216
static void setParams(BlackOilPolymerParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilpolymermodules.hh:97
static TabulatedTwoDFunction & getSkprwatTable(const int tableNumber)
get the SKPRWAT table
Definition: blackoilpolymermodules.hh:119
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilpolymermodules.hh:352
static Scalar plymixparToddLongstaff(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:440
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:277
static Scalar plyrockMaxAdsorbtion(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:400
static Scalar plyrockResidualResistanceFactor(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:376
static Evaluation computeShearFactor(const Evaluation &polymerConcentration, unsigned pvtnumRegionIdx, const Evaluation &v0)
Computes the shear factor.
Definition: blackoilpolymermodules.hh:475
VTK output module for the black oil model's polymer related quantities.
Definition: vtkblackoilpolymermodule.hpp:61
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilpolymermodule.hpp:91
Definition: blackoilboundaryratevector.hh:39
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
The Opm property system, traits with inheritance.
Definition: blackoilpolymerparams.hpp:85
Definition: blackoilpolymerparams.hpp:92
Struct holding the parameters for the BlackOilPolymerModule class.
Definition: blackoilpolymerparams.hpp:45
IntervalTabulated2DFunction< Scalar > TabulatedTwoDFunction
Definition: blackoilpolymerparams.hpp:47
Tabulated1DFunction< Scalar > TabulatedFunction
Definition: blackoilpolymerparams.hpp:46