24#ifndef EWOMS_BLACK_OIL_PRIMARY_VARIABLES_HH
25#define EWOMS_BLACK_OIL_PRIMARY_VARIABLES_HH
38#include <dune/common/fvector.hh>
40#include <opm/material/constraintsolvers/NcpFlash.hpp>
41#include <opm/material/fluidstates/CompositionalFluidState.hpp>
42#include <opm/material/fluidstates/SimpleModularFluidState.hpp>
43#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
44#include <opm/material/common/Valgrind.hpp>
47 template<
class TypeTag,
class MyTypeTag>
56template <
class TypeTag,
bool enableSolvent>
59template <
class TypeTag,
bool enableExtbo>
62template <
class TypeTag,
bool enablePolymer>
65template <
class TypeTag,
bool enableBrine>
73template <
class TypeTag>
88 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
91 enum { waterSwitchIdx = Indices::waterSwitchIdx };
92 enum { pressureSwitchIdx = Indices::pressureSwitchIdx };
93 enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
94 enum { saltConcentrationIdx = Indices::saltConcentrationIdx };
95 enum { solventSaturationIdx = Indices::solventSaturationIdx };
97 static constexpr bool compositionSwitchEnabled = Indices::compositionSwitchIdx >= 0;
98 static constexpr bool waterEnabled = Indices::waterEnabled;
99 static constexpr bool gasEnabled = Indices::gasEnabled;
100 static constexpr bool oilEnabled = Indices::oilEnabled;
103 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
104 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
105 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
106 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
109 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
110 enum { enableSolvent = getPropValue<TypeTag, Properties::EnableSolvent>() };
111 enum { enableExtbo = getPropValue<TypeTag, Properties::EnableExtbo>() };
112 enum { enablePolymer = getPropValue<TypeTag, Properties::EnablePolymer>() };
113 enum { enableFoam = getPropValue<TypeTag, Properties::EnableFoam>() };
114 enum { enableBrine = getPropValue<TypeTag, Properties::EnableBrine>() };
115 enum { enableSaltPrecipitation = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>() };
116 enum { enableVapwat = getPropValue<TypeTag, Properties::EnableVapwat>() };
117 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
118 enum { enableTemperature = getPropValue<TypeTag, Properties::EnableTemperature>() };
119 enum { enableMICP = getPropValue<TypeTag, Properties::EnableMICP>() };
120 enum { gasCompIdx = FluidSystem::gasCompIdx };
121 enum { waterCompIdx = FluidSystem::waterCompIdx };
122 enum { oilCompIdx = FluidSystem::oilCompIdx };
124 using Toolbox = MathToolbox<Evaluation>;
125 using ComponentVector = Dune::FieldVector<Scalar, numComponents>;
134 static_assert(numPhases == 3,
"The black-oil model assumes three phases!");
135 static_assert(numComponents == 3,
"The black-oil model assumes three components!");
172 Valgrind::SetUndefined(*
this);
182 Valgrind::SetUndefined(primaryVarsMeaningWater_);
183 Valgrind::SetUndefined(primaryVarsMeaningGas_);
184 Valgrind::SetUndefined(primaryVarsMeaningPressure_);
185 Valgrind::SetUndefined(primaryVarsMeaningBrine_);
186 Valgrind::SetUndefined(primaryVarsMeaningSolvent_);
199 result.pvtRegionIdx_ = 1;
205 for (
size_t i = 0; i < result.size(); ++i) {
215 pressureScale_ = Parameters::get<TypeTag, Properties::PressureScale>();
220 Parameters::registerParam<TypeTag, Properties::PressureScale>
221 (
"Scaling of pressure primary variable");
226 pressureScale_ = val;
233 if (varIdx == pressureSwitchIdx) {
234 scale = this->pressureScale_;
236 if (std::is_same<Evaluation, Scalar>::value)
237 return (*
this)[varIdx] * scale;
240 if (timeIdx == linearizationType.time)
241 return Toolbox::createVariable((*
this)[varIdx], varIdx) * scale;
243 return Toolbox::createConstant((*
this)[varIdx]) * scale;
255 { pvtRegionIdx_ =
static_cast<unsigned short>(value); }
261 {
return pvtRegionIdx_; }
268 {
return primaryVarsMeaningWater_; }
275 { primaryVarsMeaningWater_ = newMeaning; }
282 {
return primaryVarsMeaningPressure_; }
289 { primaryVarsMeaningPressure_ = newMeaning; }
296 {
return primaryVarsMeaningGas_; }
303 { primaryVarsMeaningGas_ = newMeaning; }
306 {
return primaryVarsMeaningBrine_; }
314 { primaryVarsMeaningBrine_ = newMeaning; }
318 {
return primaryVarsMeaningSolvent_; }
326 { primaryVarsMeaningSolvent_ = newMeaning; }
331 template <
class Flu
idState>
333 const MaterialLawParams& matParams,
334 bool isInEquilibrium =
false)
336 using ConstEvaluation =
typename std::remove_reference<typename FluidState::Scalar>::type;
337 using FsEvaluation =
typename std::remove_const<ConstEvaluation>::type;
338 using FsToolbox = MathToolbox<FsEvaluation>;
342 for (
unsigned phaseIdx = 1; phaseIdx < numPhases; ++phaseIdx) {
343 Valgrind::CheckDefined(fluidState.temperature(0));
344 Valgrind::CheckDefined(fluidState.temperature(phaseIdx));
346 assert(fluidState.temperature(0) == fluidState.temperature(phaseIdx));
352 if (isInEquilibrium) {
359 typename FluidSystem::template ParameterCache<Scalar> paramCache;
360 paramCache.setRegionIndex(pvtRegionIdx_);
361 paramCache.setMaxOilSat(FsToolbox::value(fluidState.saturation(oilPhaseIdx)));
364 using NcpFlash = NcpFlash<Scalar, FluidSystem>;
365 using FlashFluidState = CompositionalFluidState<Scalar, FluidSystem>;
366 FlashFluidState fsFlash;
367 fsFlash.setTemperature(FsToolbox::value(fluidState.temperature(0)));
368 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
369 fsFlash.setPressure(phaseIdx, FsToolbox::value(fluidState.pressure(phaseIdx)));
370 fsFlash.setSaturation(phaseIdx, FsToolbox::value(fluidState.saturation(phaseIdx)));
371 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
372 fsFlash.setMoleFraction(phaseIdx, compIdx, FsToolbox::value(fluidState.moleFraction(phaseIdx, compIdx)));
375 paramCache.updateAll(fsFlash);
376 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
377 if (!FluidSystem::phaseIsActive(phaseIdx))
380 Scalar rho = FluidSystem::template density<FlashFluidState, Scalar>(fsFlash, paramCache, phaseIdx);
381 fsFlash.setDensity(phaseIdx, rho);
385 ComponentVector globalMolarities(0.0);
386 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
387 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
388 if (!FluidSystem::phaseIsActive(phaseIdx))
391 globalMolarities[compIdx] +=
392 fsFlash.saturation(phaseIdx) * fsFlash.molarity(phaseIdx, compIdx);
400 NcpFlash::template solve<MaterialLaw>(fsFlash, matParams, paramCache, globalMolarities);
409 template <
class Flu
idState>
412 using ConstEvaluation =
typename std::remove_reference<typename FluidState::Scalar>::type;
413 using FsEvaluation =
typename std::remove_const<ConstEvaluation>::type;
414 using FsToolbox = MathToolbox<FsEvaluation>;
416 bool gasPresent = FluidSystem::phaseIsActive(gasPhaseIdx)?(fluidState.saturation(gasPhaseIdx) > 0.0):
false;
417 bool oilPresent = FluidSystem::phaseIsActive(oilPhaseIdx)?(fluidState.saturation(oilPhaseIdx) > 0.0):
false;
418 bool waterPresent = FluidSystem::phaseIsActive(waterPhaseIdx)?(fluidState.saturation(waterPhaseIdx) > 0.0):
false;
419 const auto& saltSaturation = BlackOil::getSaltSaturation_<FluidSystem, FluidState, Scalar>(fluidState, pvtRegionIdx_);
420 bool precipitatedSaltPresent = enableSaltPrecipitation?(saltSaturation > 0.0):
false;
421 bool oneActivePhases = FluidSystem::numActivePhases() == 1;
428 if (gasPresent && FluidSystem::enableVaporizedOil() && !oilPresent){
430 }
else if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
432 }
else if ( waterPresent && FluidSystem::enableDissolvedGasInWater() && !gasPresent){
434 }
else if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
437 assert(FluidSystem::phaseIsActive(waterPhaseIdx));
445 if ( waterPresent && gasPresent ){
447 }
else if (gasPresent && FluidSystem::enableVaporizedWater()) {
449 }
else if (waterPresent && FluidSystem::enableDissolvedGasInWater()) {
451 }
else if (FluidSystem::phaseIsActive(waterPhaseIdx) && !oneActivePhases) {
462 if ( gasPresent && oilPresent ) {
464 }
else if (oilPresent && FluidSystem::enableDissolvedGas()) {
466 }
else if (gasPresent && FluidSystem::enableVaporizedOil()){
468 }
else if (FluidSystem::phaseIsActive(gasPhaseIdx) && FluidSystem::phaseIsActive(oilPhaseIdx)) {
475 if constexpr (enableSaltPrecipitation){
476 if (precipitatedSaltPresent)
487 this->setScaledPressure_(FsToolbox::value(fluidState.pressure(oilPhaseIdx)));
490 this->setScaledPressure_(FsToolbox::value(fluidState.pressure(gasPhaseIdx)));
493 this->setScaledPressure_(FsToolbox::value(fluidState.pressure(waterPhaseIdx)));
496 throw std::logic_error(
"No valid primary variable selected for pressure");
501 (*this)[waterSwitchIdx] = FsToolbox::value(fluidState.saturation(waterPhaseIdx));
506 const auto& rvw = BlackOil::getRvw_<FluidSystem, FluidState, Scalar>(fluidState, pvtRegionIdx_);
507 (*this)[waterSwitchIdx] = rvw;
512 const auto& Rsw = BlackOil::getRsw_<FluidSystem, FluidState, Scalar>(fluidState, pvtRegionIdx_);
513 (*this)[waterSwitchIdx] = Rsw;
521 throw std::logic_error(
"No valid primary variable selected for water");
526 (*this)[compositionSwitchIdx] = FsToolbox::value(fluidState.saturation(gasPhaseIdx));
531 const auto& rs = BlackOil::getRs_<FluidSystem, FluidState, Scalar>(fluidState, pvtRegionIdx_);
532 (*this)[compositionSwitchIdx] = rs;
537 const auto& rv = BlackOil::getRv_<FluidSystem, FluidState, Scalar>(fluidState, pvtRegionIdx_);
538 (*this)[compositionSwitchIdx] = rv;
546 throw std::logic_error(
"No valid primary variable selected for composision");
562 unsigned globalDofIdx,
563 [[maybe_unused]] Scalar swMaximum,
564 Scalar thresholdWaterFilledCell, Scalar eps = 0.0)
581 Scalar saltConcentration = 0.0;
582 const Scalar& T = asImp_().temperature_(problem, globalDofIdx);
584 sw = (*this)[waterSwitchIdx];
586 sg = (*this)[compositionSwitchIdx];
593 if constexpr (enableSaltPrecipitation) {
596 saltConcentration = saltSolubility;
597 Scalar saltSat = (*this)[saltConcentrationIdx];
600 (*this)[saltConcentrationIdx] = saltSolubility;
604 saltConcentration = (*this)[saltConcentrationIdx];
605 if (saltConcentration > saltSolubility + eps){
607 (*this)[saltConcentrationIdx] = 0.0;
615 if constexpr (enableSolvent) {
617 Scalar p = (*this)[pressureSwitchIdx];
620 Scalar solSat = (*this)[solventSaturationIdx];
623 (*this)[solventSaturationIdx] = solLimit;
628 Scalar rsolw = (*this)[solventSaturationIdx];
629 if (rsolw > solLimit + eps){
631 (*this)[solventSaturationIdx] = 0.0;
638 bool changed =
false;
646 if (sw >= thresholdWaterFilledCell && !FluidSystem::enableDissolvedGasInWater()) {
649 if constexpr (waterEnabled) {
650 (*this)[Indices::waterSwitchIdx] = std::min(swMaximum, sw);
654 if constexpr (compositionSwitchEnabled)
655 (*
this)[Indices::compositionSwitchIdx] = 0.0;
659 if constexpr (compositionSwitchEnabled)
668 unsigned satnumRegionIdx = problem.satnumRegionIndex(globalDofIdx);
669 Scalar Sp = saltConcentration_();
670 Scalar porosityFactor = min(1.0 - Sp, 1.0);
672 pcFactor_ = pcfactTable.eval(porosityFactor,
true);
682 if(sw < -eps && sg > eps && FluidSystem::enableVaporizedWater()) {
683 Scalar p = this->pressure_();
685 std::array<Scalar, numPhases> pC = { 0.0 };
686 const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
687 Scalar so = 1.0 - sg - solventSaturation_();
688 computeCapillaryPressures_(pC, so, sg + solventSaturation_(), 0.0, matParams);
689 p += pcFactor_ * (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
691 Scalar rvwSat = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(pvtRegionIdx_,
696 (*this)[Indices::waterSwitchIdx] = rvwSat;
702 if(sg < -eps && sw > eps && FluidSystem::enableDissolvedGasInWater()) {
703 const Scalar pg = this->pressure_();
705 std::array<Scalar, numPhases> pC = { 0.0 };
706 const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
707 Scalar so = 1.0 - sw - solventSaturation_();
708 computeCapillaryPressures_(pC, so, 0.0, sw, matParams);
709 Scalar pw = pg + pcFactor_ * (pC[waterPhaseIdx] - pC[gasPhaseIdx]);
710 Scalar rswSat = FluidSystem::waterPvt().saturatedGasDissolutionFactor(pvtRegionIdx_,
715 (*this)[Indices::waterSwitchIdx] = rswSat;
717 this->setScaledPressure_(pw);
725 const Scalar& rvw = (*this)[waterSwitchIdx];
726 Scalar p = this->pressure_();
728 std::array<Scalar, numPhases> pC = { 0.0 };
729 const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
730 Scalar so = 1.0 - sg - solventSaturation_();
731 computeCapillaryPressures_(pC, so, sg + solventSaturation_(), 0.0, matParams);
732 p += pcFactor_ * (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
734 Scalar rvwSat = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(pvtRegionIdx_,
739 if (rvw > rvwSat*(1.0 + eps)) {
741 (*this)[Indices::waterSwitchIdx] = 0.0;
751 const Scalar& pw = this->pressure_();
753 Scalar rswSat = FluidSystem::waterPvt().saturatedGasDissolutionFactor(pvtRegionIdx_,
758 Scalar rsw = (*this)[Indices::waterSwitchIdx];
762 (*this)[Indices::waterSwitchIdx] = 1.0;
764 std::array<Scalar, numPhases> pC = { 0.0 };
765 const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
766 computeCapillaryPressures_(pC, 0.0, 0.0, 1.0, matParams);
767 Scalar pg = pw + pcFactor_ * (pC[gasPhaseIdx] - pC[waterPhaseIdx]);
768 this->setScaledPressure_(pg);
778 throw std::logic_error(
"No valid primary variable selected for water");
792 Scalar s = 1.0 - sw - solventSaturation_();
793 if (sg < -eps && s > 0.0 && FluidSystem::enableDissolvedGas()) {
794 const Scalar po = this->pressure_();
796 Scalar soMax = std::max(s, problem.maxOilSaturation(globalDofIdx));
797 Scalar rsMax = problem.maxGasDissolutionFactor(0, globalDofIdx);
801 : FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_,
806 (*this)[Indices::compositionSwitchIdx] = std::min(rsMax, rsSat);
809 Scalar so = 1.0 - sw - solventSaturation_() - sg;
810 if (so < -eps && sg > 0.0 && FluidSystem::enableVaporizedOil()) {
815 const Scalar po = this->pressure_();
816 std::array<Scalar, numPhases> pC = { 0.0 };
817 const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
818 computeCapillaryPressures_(pC, 0.0, sg + solventSaturation_(), sw, matParams);
819 Scalar pg = po + pcFactor_ * (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
824 this->setScaledPressure_(pg);
825 Scalar soMax = problem.maxOilSaturation(globalDofIdx);
826 Scalar rvMax = problem.maxOilVaporizationFactor(0, globalDofIdx);
830 : FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_,
836 (*this)[Indices::compositionSwitchIdx] = std::min(rvMax, rvSat);
846 const Scalar po = this->pressure_();
847 Scalar so = 1.0 - sw - solventSaturation_();
848 Scalar soMax = std::max(so, problem.maxOilSaturation(globalDofIdx));
849 Scalar rsMax = problem.maxGasDissolutionFactor(0, globalDofIdx);
853 : FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_,
859 Scalar rs = (*this)[Indices::compositionSwitchIdx];
860 if (rs > std::min(rsMax, rsSat*(1.0 + eps))) {
863 (*this)[Indices::compositionSwitchIdx] = 0.0;
874 const Scalar pg = this->pressure_();
875 Scalar soMax = problem.maxOilSaturation(globalDofIdx);
876 Scalar rvMax = problem.maxOilVaporizationFactor(0, globalDofIdx);
880 : FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_,
886 Scalar rv = (*this)[Indices::compositionSwitchIdx];
887 if (rv > std::min(rvMax, rvSat*(1.0 + eps))) {
891 Scalar sg2 = 1.0 - sw - solventSaturation_();
892 std::array<Scalar, numPhases> pC = { 0.0 };
893 const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
894 computeCapillaryPressures_(pC,
896 sg2 + solventSaturation_(),
899 Scalar po = pg + pcFactor_ * (pC[oilPhaseIdx] - pC[gasPhaseIdx]);
903 this->setScaledPressure_(po);
904 (*this)[Indices::compositionSwitchIdx] = sg2;
914 throw std::logic_error(
"No valid primary variable selected for water");
926 sw = (*this)[Indices::waterSwitchIdx];
929 sg = (*this)[Indices::compositionSwitchIdx];
933 ssol =(*this) [Indices::solventSaturationIdx];
935 Scalar so = 1.0 - sw - sg - ssol;
936 sw = std::min(std::max(sw,0.0),1.0);
937 so = std::min(std::max(so,0.0),1.0);
938 sg = std::min(std::max(sg,0.0),1.0);
939 ssol = std::min(std::max(ssol,0.0),1.0);
940 Scalar st = sw + so + sg + ssol;
946 (*this)[Indices::waterSwitchIdx] = sw;
948 (*this)[Indices::compositionSwitchIdx] = sg;
950 (*this) [Indices::solventSaturationIdx] = ssol;
958 for (
unsigned i = 0; i < numEq; ++i)
975 for (
unsigned i = 0; i < this->size(); ++i)
976 Valgrind::CheckDefined((*
this)[i]);
979 Valgrind::CheckDefined(primaryVarsMeaningWater_);
980 Valgrind::CheckDefined(primaryVarsMeaningGas_);
981 Valgrind::CheckDefined(primaryVarsMeaningPressure_);
982 Valgrind::CheckDefined(primaryVarsMeaningBrine_);
983 Valgrind::CheckDefined(primaryVarsMeaningSolvent_);
985 Valgrind::CheckDefined(pvtRegionIdx_);
989 template<
class Serializer>
992 using FV = Dune::FieldVector<double,getPropValue<TypeTag, Properties::NumEq>()>;
993 serializer(
static_cast<FV&
>(*
this));
994 serializer(primaryVarsMeaningWater_);
995 serializer(primaryVarsMeaningPressure_);
996 serializer(primaryVarsMeaningGas_);
997 serializer(primaryVarsMeaningBrine_);
998 serializer(primaryVarsMeaningSolvent_);
999 serializer(pvtRegionIdx_);
1005 this->primaryVarsMeaningWater_ == rhs.primaryVarsMeaningWater_ &&
1006 this->primaryVarsMeaningPressure_ == rhs.primaryVarsMeaningPressure_ &&
1007 this->primaryVarsMeaningGas_ == rhs.primaryVarsMeaningGas_ &&
1008 this->primaryVarsMeaningBrine_ == rhs.primaryVarsMeaningBrine_ &&
1009 this->primaryVarsMeaningSolvent_ == rhs.primaryVarsMeaningSolvent_ &&
1010 this->pvtRegionIdx_ == rhs.pvtRegionIdx_;
1014 Implementation& asImp_()
1015 {
return *
static_cast<Implementation*
>(
this); }
1017 const Implementation& asImp_()
const
1018 {
return *
static_cast<const Implementation*
>(
this); }
1020 Scalar solventSaturation_()
const
1022 if constexpr (enableSolvent) {
1024 return (*
this)[Indices::solventSaturationIdx];
1029 Scalar zFraction_()
const
1031 if constexpr (enableExtbo)
1032 return (*
this)[Indices::zFractionIdx];
1037 Scalar polymerConcentration_()
const
1039 if constexpr (enablePolymer)
1040 return (*
this)[Indices::polymerConcentrationIdx];
1045 Scalar foamConcentration_()
const
1047 if constexpr (enableFoam)
1048 return (*
this)[Indices::foamConcentrationIdx];
1053 Scalar saltConcentration_()
const
1055 if constexpr (enableBrine)
1056 return (*
this)[Indices::saltConcentrationIdx];
1061 Scalar temperature_(
const Problem& problem, [[maybe_unused]]
unsigned globalDofIdx)
const
1063 if constexpr (enableEnergy)
1064 return (*
this)[Indices::temperatureIdx];
1065 else if constexpr( enableTemperature)
1066 return problem.temperature(globalDofIdx, 0);
1069 return FluidSystem::reservoirTemperature();
1072 Scalar microbialConcentration_()
const
1074 if constexpr (enableMICP)
1075 return (*
this)[Indices::microbialConcentrationIdx];
1080 Scalar oxygenConcentration_()
const
1082 if constexpr (enableMICP)
1083 return (*
this)[Indices::oxygenConcentrationIdx];
1088 Scalar ureaConcentration_()
const
1090 if constexpr (enableMICP)
1091 return (*
this)[Indices::ureaConcentrationIdx];
1096 Scalar biofilmConcentration_()
const
1098 if constexpr (enableMICP)
1099 return (*
this)[Indices::biofilmConcentrationIdx];
1104 Scalar calciteConcentration_()
const
1106 if constexpr (enableMICP)
1107 return (*
this)[Indices::calciteConcentrationIdx];
1112 template <
class Container>
1113 void computeCapillaryPressures_(Container& result,
1117 const MaterialLawParams& matParams)
const
1119 using SatOnlyFluidState = SimpleModularFluidState<Scalar,
1131 SatOnlyFluidState fluidState;
1132 fluidState.setSaturation(waterPhaseIdx, sw);
1133 fluidState.setSaturation(oilPhaseIdx, so);
1134 fluidState.setSaturation(gasPhaseIdx, sg);
1136 MaterialLaw::capillaryPressures(result, matParams, fluidState);
1139 Scalar pressure_()
const
1141 return (*
this)[Indices::pressureSwitchIdx] * this->pressureScale_;
1144 void setScaledPressure_(Scalar pressure)
1146 (*this)[Indices::pressureSwitchIdx] = pressure / (this->pressureScale_);
1154 unsigned short pvtRegionIdx_;
1156 static inline Scalar pressureScale_ = 1.0;
Contains the classes required to extend the black-oil model by brine.
Contains the classes required to extend the black-oil model by energy.
Contains the classes required to extend the black-oil model by solvent component. For details,...
Contains the classes required to extend the black-oil model to include the effects of foam.
Contains the classes required to extend the black-oil model by MICP.
Contains the classes required to extend the black-oil model by polymer.
Declares the properties required by the black oil model.
Contains the classes required to extend the black-oil model by solvents.
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:58
static Scalar saltSol(unsigned regionIdx)
Definition: blackoilbrinemodules.hh:394
static const TabulatedFunction & pcfactTable(unsigned satnumRegionIdx)
Definition: blackoilbrinemodules.hh:342
static bool hasPcfactTables()
Definition: blackoilbrinemodules.hh:386
Contains the high level supplements required to extend the black oil model by energy.
Definition: blackoilenergymodules.hh:52
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar)
Assign the energy specific primary variables to a PrimaryVariables object.
Definition: blackoilenergymodules.hh:247
Contains the high level supplements required to extend the black oil model.
Definition: blackoilextbomodules.hh:69
static Value rs(unsigned pvtRegionIdx, const Value &pressure, const Value &z)
Definition: blackoilextbomodules.hh:518
static Value rv(unsigned pvtRegionIdx, const Value &pressure, const Value &z)
Definition: blackoilextbomodules.hh:523
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition: blackoilfoammodules.hh:57
Contains the high level supplements required to extend the black oil model by MICP.
Definition: blackoilmicpmodules.hh:56
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:61
Represents the primary variables used by the black-oil model.
Definition: blackoilprimaryvariables.hh:75
bool operator==(const BlackOilPrimaryVariables &rhs) const
Definition: blackoilprimaryvariables.hh:1002
unsigned pvtRegionIndex() const
Return the index of the region which should be used for PVT properties.
Definition: blackoilprimaryvariables.hh:260
WaterMeaning
Definition: blackoilprimaryvariables.hh:138
BlackOilPrimaryVariables & operator=(const BlackOilPrimaryVariables &other)=default
static void registerParameters()
Definition: blackoilprimaryvariables.hh:218
GasMeaning primaryVarsMeaningGas() const
Return the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:295
PressureMeaning primaryVarsMeaningPressure() const
Return the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:281
void serializeOp(Serializer &serializer)
Definition: blackoilprimaryvariables.hh:990
static void init()
Definition: blackoilprimaryvariables.hh:212
void setPrimaryVarsMeaningSolvent(SolventMeaning newMeaning)
Set the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:325
BlackOilPrimaryVariables()
Definition: blackoilprimaryvariables.hh:169
static BlackOilPrimaryVariables serializationTestObject()
Definition: blackoilprimaryvariables.hh:196
SolventMeaning primaryVarsMeaningSolvent() const
Definition: blackoilprimaryvariables.hh:317
bool chopAndNormalizeSaturations()
Definition: blackoilprimaryvariables.hh:919
bool adaptPrimaryVariables(const Problem &problem, unsigned globalDofIdx, Scalar swMaximum, Scalar thresholdWaterFilledCell, Scalar eps=0.0)
Adapt the interpretation of the switching variables to be physically meaningful.
Definition: blackoilprimaryvariables.hh:561
void setPrimaryVarsMeaningPressure(PressureMeaning newMeaning)
Set the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:288
WaterMeaning primaryVarsMeaningWater() const
Return the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:267
PressureMeaning
Definition: blackoilprimaryvariables.hh:145
void setPrimaryVarsMeaningWater(WaterMeaning newMeaning)
Set the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:274
void assignNaive(const FluidState &fluidState)
Directly retrieve the primary variables from an arbitrary fluid state.
Definition: blackoilprimaryvariables.hh:410
Evaluation makeEvaluation(unsigned varIdx, unsigned timeIdx, LinearizationType linearizationType=LinearizationType()) const
Definition: blackoilprimaryvariables.hh:230
void assignMassConservative(const FluidState &fluidState, const MaterialLawParams &matParams, bool isInEquilibrium=false)
Set the primary variables from an arbitrary fluid state in a mass conservative way.
Definition: blackoilprimaryvariables.hh:332
BlackOilPrimaryVariables(Scalar value)
Constructor with assignment from scalar.
Definition: blackoilprimaryvariables.hh:179
void setPrimaryVarsMeaningBrine(BrineMeaning newMeaning)
Set the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:313
BlackOilPrimaryVariables(const BlackOilPrimaryVariables &value)=default
Copy constructor.
SolventMeaning
Definition: blackoilprimaryvariables.hh:163
BrineMeaning primaryVarsMeaningBrine() const
Definition: blackoilprimaryvariables.hh:305
BlackOilPrimaryVariables & operator=(Scalar value)
Definition: blackoilprimaryvariables.hh:956
GasMeaning
Definition: blackoilprimaryvariables.hh:150
BrineMeaning
Definition: blackoilprimaryvariables.hh:157
void setPressureScale(Scalar val)
Definition: blackoilprimaryvariables.hh:224
void setPrimaryVarsMeaningGas(GasMeaning newMeaning)
Set the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:302
void setPvtRegionIndex(unsigned value)
Set the index of the region which should be used for PVT properties.
Definition: blackoilprimaryvariables.hh:254
void checkDefined() const
Instruct valgrind to check the definedness of all attributes of this class.
Definition: blackoilprimaryvariables.hh:971
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilsolventmodules.hh:68
static bool isSolubleInWater()
Definition: blackoilsolventmodules.hh:747
static const Value solubilityLimit(unsigned pvtIdx, const Value &temperature, const Value &pressure, const Value &saltConcentration)
Definition: blackoilsolventmodules.hh:735
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:52
Definition: blackoilmodel.hh:72
Definition: blackoilboundaryratevector.hh:37
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:242
Definition: linearizationtype.hh:35
Definition: blackoilprimaryvariables.hh:48
static constexpr type value
Definition: blackoilprimaryvariables.hh:50
GetPropType< TypeTag, Scalar > type
Definition: blackoilprimaryvariables.hh:49