28#ifndef EWOMS_BLACK_OIL_SOLVENT_MODULE_HH
29#define EWOMS_BLACK_OIL_SOLVENT_MODULE_HH
31#include <dune/common/fvector.hh>
33#include <opm/common/Exceptions.hpp>
35#include <opm/material/common/MathToolbox.hpp>
36#include <opm/material/common/Valgrind.hpp>
37#include <opm/material/fluidsystems/blackoilpvt/SolventPvt.hpp>
66template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
81 using Toolbox = MathToolbox<Evaluation>;
90 static constexpr unsigned solventSaturationIdx = Indices::solventSaturationIdx;
91 static constexpr unsigned contiSolventEqIdx = Indices::contiSolventEqIdx;
92 static constexpr unsigned enableSolvent = enableSolventV;
93 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
94 static constexpr unsigned numPhases = FluidSystem::numPhases;
95 static constexpr bool blackoilConserveSurfaceVolume =
96 getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>();
97 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
102 { params_ = params; }
108 { params_.solventPvt_ = value; }
118 if constexpr (enableSolvent) {
127 Simulator& simulator)
129 if constexpr (enableSolvent) {
136 if constexpr (enableSolvent) {
137 return pvIdx == solventSaturationIdx;
148 return "saturation_solvent";
156 return static_cast<Scalar
>(1.0);
161 if constexpr (enableSolvent) {
162 return eqIdx == contiSolventEqIdx;
169 static std::string
eqName([[maybe_unused]]
unsigned eqIdx)
173 return "conti^solvent";
176 static Scalar
eqWeight([[maybe_unused]]
unsigned eqIdx)
181 return static_cast<Scalar
>(1.0);
184 template <
class LhsEval>
185 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
186 const IntensiveQuantities& intQuants)
188 if constexpr (enableSolvent) {
189 if constexpr (blackoilConserveSurfaceVolume) {
190 storage[contiSolventEqIdx] +=
191 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
192 Toolbox::template decay<LhsEval>(intQuants.solventSaturation()) *
193 Toolbox::template decay<LhsEval>(intQuants.solventInverseFormationVolumeFactor());
195 storage[contiSolventEqIdx] +=
196 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
197 Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(waterPhaseIdx)) *
198 Toolbox::template decay<LhsEval>(intQuants.fluidState().invB(waterPhaseIdx)) *
199 Toolbox::template decay<LhsEval>(intQuants.rsSolw());
203 storage[contiSolventEqIdx] +=
204 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
205 Toolbox::template decay<LhsEval>(intQuants.solventSaturation()) *
206 Toolbox::template decay<LhsEval>(intQuants.solventDensity());
208 storage[contiSolventEqIdx] +=
209 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
210 Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(waterPhaseIdx)) *
211 Toolbox::template decay<LhsEval>(intQuants.fluidState().density(waterPhaseIdx)) *
212 Toolbox::template decay<LhsEval>(intQuants.rsSolw());
219 [[maybe_unused]]
const ElementContext& elemCtx,
220 [[maybe_unused]]
unsigned scvfIdx,
221 [[maybe_unused]]
unsigned timeIdx)
224 if constexpr (enableSolvent) {
225 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
227 const unsigned upIdx = extQuants.solventUpstreamIndex();
228 const unsigned inIdx = extQuants.interiorIndex();
229 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
231 if constexpr (blackoilConserveSurfaceVolume) {
232 if (upIdx == inIdx) {
233 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() *
234 up.solventInverseFormationVolumeFactor();
237 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() *
238 decay<Scalar>(up.solventInverseFormationVolumeFactor());
243 if (upIdx == inIdx) {
244 flux[contiSolventEqIdx] +=
245 extQuants.volumeFlux(waterPhaseIdx) *
246 up.fluidState().invB(waterPhaseIdx) *
250 flux[contiSolventEqIdx] +=
251 extQuants.volumeFlux(waterPhaseIdx) *
252 decay<Scalar>(up.fluidState().invB(waterPhaseIdx)) *
253 decay<Scalar>(up.rsSolw());
258 if (upIdx == inIdx) {
259 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() * up.solventDensity();
262 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() * decay<Scalar>(up.solventDensity());
266 if (upIdx == inIdx) {
267 flux[contiSolventEqIdx] +=
268 extQuants.volumeFlux(waterPhaseIdx) *
269 up.fluidState().density(waterPhaseIdx) *
273 flux[contiSolventEqIdx] +=
274 extQuants.volumeFlux(waterPhaseIdx) *
275 decay<Scalar>(up.fluidState().density(waterPhaseIdx)) *
276 decay<Scalar>(up.rsSolw());
287 Scalar solventSaturation,
290 if constexpr (!enableSolvent) {
291 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Disabled);
296 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Ss);
297 priVars[solventSaturationIdx] = solventSaturation;
299 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Rsolw);
300 priVars[solventSaturationIdx] = solventRsw;
308 const PrimaryVariables& oldPv,
309 const EqVector& delta)
311 if constexpr (enableSolvent) {
313 newPv[solventSaturationIdx] = oldPv[solventSaturationIdx] - delta[solventSaturationIdx];
326 return static_cast<Scalar
>(0.0);
335 return std::abs(Toolbox::scalarValue(resid[contiSolventEqIdx]));
338 template <
class DofEntity>
339 static void serializeEntity(
const Model& model, std::ostream& outstream,
const DofEntity& dof)
341 if constexpr (enableSolvent) {
342 const unsigned dofIdx = model.dofMapper().index(dof);
344 const PrimaryVariables& priVars = model.solution(0)[dofIdx];
345 outstream << priVars[solventSaturationIdx];
349 template <
class DofEntity>
352 if constexpr (enableSolvent) {
353 const unsigned dofIdx = model.dofMapper().index(dof);
355 PrimaryVariables& priVars0 = model.solution(0)[dofIdx];
356 PrimaryVariables& priVars1 = model.solution(1)[dofIdx];
358 instream >> priVars0[solventSaturationIdx];
361 priVars1 = priVars0[solventSaturationIdx];
366 {
return params_.solventPvt_; }
369 {
return params_.co2GasPvt_; }
372 {
return params_.h2GasPvt_; }
375 {
return params_.brineCo2Pvt_; }
378 {
return params_.brineH2Pvt_; }
380 template <
class ElemContext>
381 static const TabulatedFunction&
ssfnKrg(
const ElemContext& elemCtx,
385 const unsigned satnumRegionIdx =
386 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
387 return params_.ssfnKrg_[satnumRegionIdx];
390 template <
class ElemContext>
391 static const TabulatedFunction&
ssfnKrs(
const ElemContext& elemCtx,
395 const unsigned satnumRegionIdx =
396 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
397 return params_.ssfnKrs_[satnumRegionIdx];
400 template <
class ElemContext>
401 static const TabulatedFunction&
sof2Krn(
const ElemContext& elemCtx,
405 const unsigned satnumRegionIdx =
406 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
407 return params_.sof2Krn_[satnumRegionIdx];
410 template <
class ElemContext>
411 static const TabulatedFunction&
misc(
const ElemContext& elemCtx,
415 const unsigned miscnumRegionIdx =
416 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
417 return params_.misc_[miscnumRegionIdx];
420 template <
class ElemContext>
421 static const TabulatedFunction&
pmisc(
const ElemContext& elemCtx,
425 const unsigned miscnumRegionIdx =
426 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
427 return params_.pmisc_[miscnumRegionIdx];
430 template <
class ElemContext>
431 static const TabulatedFunction&
msfnKrsg(
const ElemContext& elemCtx,
435 const unsigned satnumRegionIdx =
436 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
437 return params_.msfnKrsg_[satnumRegionIdx];
440 template <
class ElemContext>
441 static const TabulatedFunction&
msfnKro(
const ElemContext& elemCtx,
445 const unsigned satnumRegionIdx =
446 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
447 return params_.msfnKro_[satnumRegionIdx];
450 template <
class ElemContext>
451 static const TabulatedFunction&
sorwmis(
const ElemContext& elemCtx,
455 const unsigned miscnumRegionIdx =
456 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
457 return params_.sorwmis_[miscnumRegionIdx];
460 template <
class ElemContext>
461 static const TabulatedFunction&
sgcwmis(
const ElemContext& elemCtx,
465 const unsigned miscnumRegionIdx =
466 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
467 return params_.sgcwmis_[miscnumRegionIdx];
470 template <
class ElemContext>
471 static const TabulatedFunction&
tlPMixTable(
const ElemContext& elemCtx,
475 const unsigned miscnumRegionIdx =
476 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
477 return params_.tlPMixTable_[miscnumRegionIdx];
480 template <
class ElemContext>
485 const unsigned miscnumRegionIdx =
486 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
487 return params_.tlMixParamViscosity_[miscnumRegionIdx];
491 template <
class ElemContext>
496 const unsigned miscnumRegionIdx =
497 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
498 return params_.tlMixParamDensity_[miscnumRegionIdx];
502 {
return params_.isMiscible_; }
504 template <
class Value>
506 const Value& pressure,
const Value& saltConcentration)
514 return brineCo2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature,
515 pressure, saltConcentration);
518 return brineH2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature,
519 pressure, saltConcentration);
524 {
return params_.rsSolw_active_; }
527 {
return params_.co2sol_; }
530 {
return params_.h2sol_; }
536template <
class TypeTag,
bool enableSolventV>
537BlackOilSolventParams<typename BlackOilSolventModule<TypeTag, enableSolventV>::Scalar>
538BlackOilSolventModule<TypeTag, enableSolventV>::params_;
540template <
class TypeTag,
bool enableSolventV>
550template <
class TypeTag>
566 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
567 static constexpr int solventSaturationIdx = Indices::solventSaturationIdx;
568 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
569 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
570 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
571 static constexpr double cutOff = 1e-12;
584 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
585 this->solventPreSatFuncUpdate_(priVars, timeIdx, elemCtx.linearizationType());
589 const unsigned timeIdx,
592 auto& fs = asImp_().fluidState_;
593 solventSaturation_ = 0.0;
594 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
595 solventSaturation_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, linearizationType);
598 hydrocarbonSaturation_ = fs.saturation(gasPhaseIdx);
601 if (solventSaturation().value() < cutOff) {
607 fs.setSaturation(gasPhaseIdx, hydrocarbonSaturation_ + solventSaturation_);
621 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
622 this->solventPostSatFuncUpdate_(elemCtx.problem(), priVars, elemCtx.globalSpaceIndex(dofIdx, timeIdx), timeIdx, elemCtx.linearizationType());
627 template <
class Problem>
628 struct ProblemAndCellIndexOnlyContext
630 const Problem& problem_;
632 const Problem& problem()
const {
return problem_; }
633 unsigned int globalSpaceIndex([[maybe_unused]]
const unsigned int spaceIdx,
634 [[maybe_unused]]
const unsigned int timeIdx)
const
642 const PrimaryVariables& priVars,
643 const unsigned globalSpaceIdx,
644 const unsigned timeIdx,
649 auto& fs = asImp_().fluidState_;
650 fs.setSaturation(gasPhaseIdx, hydrocarbonSaturation_);
654 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
655 rsSolw_ = SolventModule::solubilityLimit(asImp_().pvtRegionIndex(),
656 fs.temperature(waterPhaseIdx),
657 fs.pressure(waterPhaseIdx),
658 fs.saltConcentration());
659 }
else if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Rsolw) {
660 rsSolw_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, linearizationType);
663 solventMobility_ = 0.0;
666 if (solventSaturation().value() < cutOff) {
670 ProblemAndCellIndexOnlyContext<Problem> elemCtx{problem, globalSpaceIdx};
674 if (SolventModule::isMiscible()) {
675 const Evaluation& p =
676 FluidSystem::phaseIsActive(oilPhaseIdx)
677 ? fs.pressure(oilPhaseIdx)
678 : fs.pressure(gasPhaseIdx);
679 const Evaluation pmisc =
680 SolventModule::pmisc(elemCtx, dofIdx, timeIdx).eval(p,
true);
681 const Evaluation& pgImisc = fs.pressure(gasPhaseIdx);
684 Evaluation pgMisc = 0.0;
685 std::array<Evaluation, numPhases> pC;
686 const auto& materialParams = problem.materialLawParams(elemCtx, dofIdx, timeIdx);
687 MaterialLaw::capillaryPressures(pC, materialParams, fs);
690 if (priVars.primaryVarsMeaningPressure() == PrimaryVariables::PressureMeaning::Pg) {
691 pgMisc = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx,
695 const Evaluation& po =
696 priVars.makeEvaluation(Indices::pressureSwitchIdx,
697 timeIdx, linearizationType);
698 pgMisc = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
701 fs.setPressure(gasPhaseIdx, pmisc * pgMisc + (1.0 - pmisc) * pgImisc);
704 const Evaluation gasSolventSat = hydrocarbonSaturation_ + solventSaturation_;
706 if (gasSolventSat.value() < cutOff) {
710 const Evaluation Fhydgas = hydrocarbonSaturation_ / gasSolventSat;
711 const Evaluation Fsolgas = solventSaturation_ / gasSolventSat;
714 if (SolventModule::isMiscible() && FluidSystem::phaseIsActive(oilPhaseIdx)) {
715 const auto& misc = SolventModule::misc(elemCtx, dofIdx, timeIdx);
716 const auto& pmisc = SolventModule::pmisc(elemCtx, dofIdx, timeIdx);
717 const Evaluation& p =
718 FluidSystem::phaseIsActive(oilPhaseIdx)
719 ? fs.pressure(oilPhaseIdx)
720 : fs.pressure(gasPhaseIdx);
721 const Evaluation miscibility =
722 misc.eval(Fsolgas,
true) *
726 const unsigned cellIdx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
727 const auto& materialLawManager = elemCtx.problem().materialLawManager();
728 const auto& scaledDrainageInfo =
729 materialLawManager->oilWaterScaledEpsInfoDrainage(cellIdx);
731 const Scalar sogcr = scaledDrainageInfo.Sogcr;
732 Evaluation sor = sogcr;
733 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
734 const Evaluation& sw = fs.saturation(waterPhaseIdx);
735 const auto& sorwmis = SolventModule::sorwmis(elemCtx, dofIdx, timeIdx);
737 sorwmis.eval(sw,
true) + (1.0 - miscibility) * sogcr;
739 const Scalar sgcr = scaledDrainageInfo.Sgcr;
740 Evaluation sgc = sgcr;
741 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
742 const Evaluation& sw = fs.saturation(waterPhaseIdx);
743 const auto& sgcwmis = SolventModule::sgcwmis(elemCtx, dofIdx, timeIdx);
744 sgc = miscibility * sgcwmis.eval(sw,
true) + (1.0 - miscibility) * sgcr;
747 Evaluation oilGasSolventSat = gasSolventSat;
748 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
749 oilGasSolventSat += fs.saturation(oilPhaseIdx);
751 const Evaluation zero = 0.0;
752 const Evaluation oilGasSolventEffSat = std::max(oilGasSolventSat - sor - sgc, zero);
754 Evaluation F_totalGas = 0.0;
755 if (oilGasSolventEffSat.value() > cutOff) {
756 const Evaluation gasSolventEffSat = std::max(gasSolventSat - sgc, zero);
757 F_totalGas = gasSolventEffSat / oilGasSolventEffSat;
759 const auto& msfnKro = SolventModule::msfnKro(elemCtx, dofIdx, timeIdx);
760 const auto& msfnKrsg = SolventModule::msfnKrsg(elemCtx, dofIdx, timeIdx);
761 const auto& sof2Krn = SolventModule::sof2Krn(elemCtx, dofIdx, timeIdx);
763 const Evaluation mkrgt = msfnKrsg.eval(F_totalGas,
true) *
764 sof2Krn.eval(oilGasSolventSat,
true);
765 const Evaluation mkro = msfnKro.eval(F_totalGas,
true) *
766 sof2Krn.eval(oilGasSolventSat,
true);
768 Evaluation& kro = asImp_().mobility_[oilPhaseIdx];
769 Evaluation& krg = asImp_().mobility_[gasPhaseIdx];
772 krg *= 1.0 - miscibility;
773 krg += miscibility * mkrgt;
774 kro *= 1.0 - miscibility;
775 kro += miscibility * mkro;
779 const auto& ssfnKrg = SolventModule::ssfnKrg(elemCtx, dofIdx, timeIdx);
780 const auto& ssfnKrs = SolventModule::ssfnKrs(elemCtx, dofIdx, timeIdx);
782 Evaluation& krg = asImp_().mobility_[gasPhaseIdx];
783 solventMobility_ = krg * ssfnKrs.eval(Fsolgas,
true);
784 krg *= ssfnKrg.eval(Fhydgas,
true);
797 const auto& iq = asImp_();
798 const unsigned pvtRegionIdx = iq.pvtRegionIndex();
799 const Evaluation& T = iq.fluidState().temperature(gasPhaseIdx);
800 const Evaluation& p = iq.fluidState().pressure(gasPhaseIdx);
802 const Evaluation rv = 0.0;
803 const Evaluation rvw = 0.0;
804 if (SolventModule::isCO2Sol() || SolventModule::isH2Sol() ){
805 if (SolventModule::isCO2Sol()) {
806 const auto& co2gasPvt = SolventModule::co2GasPvt();
807 solventInvFormationVolumeFactor_ =
808 co2gasPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rv, rvw);
809 solventRefDensity_ = co2gasPvt.gasReferenceDensity(pvtRegionIdx);
810 solventViscosity_ = co2gasPvt.viscosity(pvtRegionIdx, T, p, rv, rvw);
812 const auto& brineCo2Pvt = SolventModule::brineCo2Pvt();
813 auto& fs = asImp_().fluidState_;
814 const auto& bw = brineCo2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rsSolw());
816 const auto denw = bw * brineCo2Pvt.waterReferenceDensity(pvtRegionIdx) +
817 rsSolw() * bw * brineCo2Pvt.gasReferenceDensity(pvtRegionIdx);
818 fs.setDensity(waterPhaseIdx, denw);
819 fs.setInvB(waterPhaseIdx, bw);
820 Evaluation& mobw = asImp_().mobility_[waterPhaseIdx];
821 const auto& muWat = fs.viscosity(waterPhaseIdx);
822 const auto& muWatEff = brineCo2Pvt.viscosity(pvtRegionIdx, T, p, rsSolw());
823 mobw *= muWat / muWatEff;
825 const auto& h2gasPvt = SolventModule::h2GasPvt();
826 solventInvFormationVolumeFactor_ =
827 h2gasPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rv, rvw);
828 solventRefDensity_ = h2gasPvt.gasReferenceDensity(pvtRegionIdx);
829 solventViscosity_ = h2gasPvt.viscosity(pvtRegionIdx, T, p, rv, rvw);
831 const auto& brineH2Pvt = SolventModule::brineH2Pvt();
832 auto& fs = asImp_().fluidState_;
833 const auto& bw = brineH2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rsSolw());
835 const auto denw = bw * brineH2Pvt.waterReferenceDensity(pvtRegionIdx) +
836 rsSolw() * bw * brineH2Pvt.gasReferenceDensity(pvtRegionIdx);
837 fs.setDensity(waterPhaseIdx, denw);
838 fs.setInvB(waterPhaseIdx, bw);
839 Evaluation& mobw = asImp_().mobility_[waterPhaseIdx];
840 const auto& muWat = fs.viscosity(waterPhaseIdx);
841 const auto& muWatEff = brineH2Pvt.viscosity(pvtRegionIdx, T, p, rsSolw());
842 mobw *= muWat / muWatEff;
845 const auto& solventPvt = SolventModule::solventPvt();
846 solventInvFormationVolumeFactor_ = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
847 solventRefDensity_ = solventPvt.referenceDensity(pvtRegionIdx);
848 solventViscosity_ = solventPvt.viscosity(pvtRegionIdx, T, p);
851 solventDensity_ = solventInvFormationVolumeFactor_*solventRefDensity_;
852 effectiveProperties(elemCtx, scvIdx, timeIdx);
854 solventMobility_ /= solventViscosity_;
858 {
return solventSaturation_; }
864 {
return solventDensity_; }
867 {
return solventViscosity_; }
870 {
return solventMobility_; }
873 {
return solventInvFormationVolumeFactor_; }
877 {
return solventRefDensity_; }
882 void effectiveProperties(
const ElementContext& elemCtx,
886 if (!SolventModule::isMiscible()) {
892 if (solventSaturation() < cutOff) {
898 auto& fs = asImp_().fluidState_;
901 Evaluation oilEffSat = 0.0;
902 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
903 oilEffSat = fs.saturation(oilPhaseIdx);
905 Evaluation gasEffSat = fs.saturation(gasPhaseIdx);
906 Evaluation solventEffSat = solventSaturation();
907 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
908 const auto& sorwmis = SolventModule::sorwmis(elemCtx, scvIdx, timeIdx);
909 const auto& sgcwmis = SolventModule::sgcwmis(elemCtx, scvIdx, timeIdx);
910 const Evaluation zero = 0.0;
911 const Evaluation& sw = fs.saturation(waterPhaseIdx);
912 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
913 oilEffSat = std::max(oilEffSat - sorwmis.eval(sw,
true), zero);
915 gasEffSat = std::max(gasEffSat - sgcwmis.eval(sw,
true), zero);
916 solventEffSat = std::max(solventEffSat - sgcwmis.eval(sw,
true), zero);
918 const Evaluation oilGasSolventEffSat = oilEffSat + gasEffSat + solventEffSat;
919 const Evaluation oilSolventEffSat = oilEffSat + solventEffSat;
920 const Evaluation solventGasEffSat = solventEffSat + gasEffSat;
927 const Evaluation& p =
928 FluidSystem::phaseIsActive(oilPhaseIdx)
929 ? fs.pressure(oilPhaseIdx)
930 : fs.pressure(gasPhaseIdx);
932 const auto& pmiscTable = SolventModule::pmisc(elemCtx, scvIdx, timeIdx);
933 const Evaluation pmisc = pmiscTable.eval(p,
true);
934 const auto& tlPMixTable = SolventModule::tlPMixTable(elemCtx, scvIdx, timeIdx);
935 const Evaluation tlMixParamMu = SolventModule::tlMixParamViscosity(elemCtx, scvIdx, timeIdx) *
936 tlPMixTable.eval(p,
true);
938 const Evaluation& muGas = fs.viscosity(gasPhaseIdx);
939 const Evaluation& muSolvent = solventViscosity_;
941 assert(muGas.value() > 0);
942 assert(muSolvent.value() > 0);
943 const Evaluation muGasPow = pow(muGas, 0.25);
944 const Evaluation muSolventPow = pow(muSolvent, 0.25);
946 Evaluation muMixSolventGas = muGas;
947 if (solventGasEffSat > cutOff) {
948 muMixSolventGas *= muSolvent / pow(((gasEffSat / solventGasEffSat) * muSolventPow) +
949 ((solventEffSat / solventGasEffSat) * muGasPow), 4.0);
952 Evaluation muOil = 1.0;
953 Evaluation muOilPow = 1.0;
954 Evaluation muMixOilSolvent = 1.0;
955 Evaluation muOilEff = 1.0;
956 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
957 muOil = fs.viscosity(oilPhaseIdx);
958 assert(muOil.value() > 0);
959 muOilPow = pow(muOil, 0.25);
960 muMixOilSolvent = muOil;
961 if (oilSolventEffSat > cutOff) {
962 muMixOilSolvent *= muSolvent / pow(((oilEffSat / oilSolventEffSat) * muSolventPow) +
963 ((solventEffSat / oilSolventEffSat) * muOilPow), 4.0);
966 muOilEff = pow(muOil,1.0 - tlMixParamMu) * pow(muMixOilSolvent, tlMixParamMu);
968 Evaluation muMixSolventGasOil = muOil;
969 if (oilGasSolventEffSat > cutOff) {
970 muMixSolventGasOil *= muSolvent * muGas / pow(((oilEffSat / oilGasSolventEffSat) *
971 muSolventPow * muGasPow) +
972 ((solventEffSat / oilGasSolventEffSat) *
973 muOilPow * muGasPow) +
974 ((gasEffSat / oilGasSolventEffSat) *
975 muSolventPow * muOilPow), 4.0);
978 Evaluation muGasEff = pow(muGas,1.0 - tlMixParamMu) * pow(muMixSolventGas, tlMixParamMu);
979 const Evaluation muSolventEff = pow(muSolvent,1.0 - tlMixParamMu) * pow(muMixSolventGasOil, tlMixParamMu);
982 const Evaluation& rhoGas = fs.density(gasPhaseIdx);
983 Evaluation rhoOil = 0.0;
984 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
985 rhoOil = fs.density(oilPhaseIdx);
988 const Evaluation& rhoSolvent = solventDensity_;
994 const Evaluation tlMixParamRho = SolventModule::tlMixParamDensity(elemCtx, scvIdx, timeIdx) *
995 tlPMixTable.eval(p,
true);
999 const Evaluation muOilEffPow = pow(pow(muOil, 1.0 - tlMixParamRho) *
1000 pow(muMixOilSolvent, tlMixParamRho), 0.25);
1001 const Evaluation muGasEffPow = pow(pow(muGas, 1.0 - tlMixParamRho) *
1002 pow(muMixSolventGas, tlMixParamRho), 0.25);
1003 const Evaluation muSolventEffPow = pow(pow(muSolvent, 1.0 - tlMixParamRho) *
1004 pow(muMixSolventGasOil, tlMixParamRho), 0.25);
1006 const Evaluation oilGasEffSaturation = oilEffSat + gasEffSat;
1007 Evaluation sof = 0.0;
1008 Evaluation sgf = 0.0;
1009 if (oilGasEffSaturation.value() > cutOff) {
1010 sof = oilEffSat / oilGasEffSaturation;
1011 sgf = gasEffSat / oilGasEffSaturation;
1014 const Evaluation muSolventOilGasPow = muSolventPow * ((sgf * muOilPow) + (sof * muGasPow));
1016 Evaluation rhoMixSolventGasOil = 0.0;
1017 if (oilGasSolventEffSat.value() > cutOff) {
1018 rhoMixSolventGasOil = (rhoOil * oilEffSat / oilGasSolventEffSat) +
1019 (rhoGas * gasEffSat / oilGasSolventEffSat) +
1020 (rhoSolvent * solventEffSat / oilGasSolventEffSat);
1023 Evaluation rhoGasEff = 0.0;
1024 if (std::abs(muSolventPow.value() - muGasPow.value()) < cutOff) {
1025 rhoGasEff = ((1.0 - tlMixParamRho) * rhoGas) +
1026 (tlMixParamRho * rhoMixSolventGasOil);
1029 const Evaluation solventGasEffFraction = (muGasPow * (muSolventPow - muGasEffPow)) /
1030 (muGasEffPow * (muSolventPow - muGasPow));
1031 rhoGasEff = (rhoGas * solventGasEffFraction) + (rhoSolvent * (1.0 - solventGasEffFraction));
1034 Evaluation rhoSolventEff = 0.0;
1035 if (std::abs((muSolventOilGasPow.value() - (muOilPow.value() * muGasPow.value()))) < cutOff) {
1036 rhoSolventEff = ((1.0 - tlMixParamRho) * rhoSolvent) + (tlMixParamRho * rhoMixSolventGasOil);
1039 const Evaluation sfraction_se = (muSolventOilGasPow -
1040 muOilPow * muGasPow * muSolventPow / muSolventEffPow) /
1041 (muSolventOilGasPow - (muOilPow * muGasPow));
1042 rhoSolventEff = (rhoSolvent * sfraction_se) +
1043 (rhoGas * sgf * (1.0 - sfraction_se)) +
1044 (rhoOil * sof * (1.0 - sfraction_se));
1048 const unsigned pvtRegionIdx = asImp_().pvtRegionIndex();
1049 Evaluation bGasEff = rhoGasEff;
1050 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1051 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) +
1052 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) * fs.Rv());
1054 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1056 const Evaluation bSolventEff = rhoSolventEff / solventRefDensity();
1059 const Evaluation bg = fs.invB(gasPhaseIdx);
1060 const Evaluation bs = solventInverseFormationVolumeFactor();
1061 const Evaluation bg_eff = pmisc * bGasEff + (1.0 - pmisc) * bg;
1062 const Evaluation bs_eff = pmisc * bSolventEff + (1.0 - pmisc) * bs;
1065 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1066 fs.setDensity(gasPhaseIdx,
1068 (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) +
1069 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx)*fs.Rv()));
1071 fs.setDensity(gasPhaseIdx,
1072 bg_eff * FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1074 solventDensity_ = bs_eff * solventRefDensity();
1077 Evaluation& mobg = asImp_().mobility_[gasPhaseIdx];
1078 muGasEff = bg_eff / (pmisc * bGasEff / muGasEff + (1.0 - pmisc) * bg / muGas);
1079 mobg *= muGas / muGasEff;
1082 solventViscosity_ = bs_eff / (pmisc * bSolventEff / muSolventEff +
1083 (1.0 - pmisc) * bs / muSolvent);
1085 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1086 Evaluation rhoOilEff = 0.0;
1087 if (std::abs(muOilPow.value() - muSolventPow.value()) < cutOff) {
1088 rhoOilEff = ((1.0 - tlMixParamRho) * rhoOil) + (tlMixParamRho * rhoMixSolventGasOil);
1091 const Evaluation solventOilEffFraction = (muOilPow * (muOilEffPow - muSolventPow)) /
1092 (muOilEffPow * (muOilPow - muSolventPow));
1093 rhoOilEff = (rhoOil * solventOilEffFraction) + (rhoSolvent * (1.0 - solventOilEffFraction));
1095 const Evaluation bOilEff =
1096 rhoOilEff / (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) +
1097 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs());
1098 const Evaluation bo = fs.invB(oilPhaseIdx);
1099 const Evaluation bo_eff = pmisc * bOilEff + (1.0 - pmisc) * bo;
1100 fs.setDensity(oilPhaseIdx,
1101 bo_eff * (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) +
1102 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs()));
1105 Evaluation& mobo = asImp_().mobility_[oilPhaseIdx];
1106 muOilEff = bo_eff / (pmisc * bOilEff / muOilEff + (1.0 - pmisc) * bo / muOil);
1107 mobo *= muOil / muOilEff;
1113 {
return *
static_cast<Implementation*
>(
this); }
1126template <
class TypeTag>
1150 {
throw std::runtime_error(
"solventSaturation() called but solvents are disabled"); }
1153 {
throw std::runtime_error(
"rsSolw() called but solvents are disabled"); }
1156 {
throw std::runtime_error(
"solventDensity() called but solvents are disabled"); }
1159 {
throw std::runtime_error(
"solventViscosity() called but solvents are disabled"); }
1162 {
throw std::runtime_error(
"solventMobility() called but solvents are disabled"); }
1165 {
throw std::runtime_error(
"solventInverseFormationVolumeFactor() called but solvents are disabled"); }
1168 {
throw std::runtime_error(
"solventRefDensity() called but solvents are disabled"); }
1178template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
1191 using Toolbox = MathToolbox<Evaluation>;
1193 static constexpr unsigned gasPhaseIdx = FluidSystem::gasPhaseIdx;
1194 static constexpr int dimWorld = GridView::dimensionworld;
1196 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
1197 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
1204 template <
class Dummy =
bool>
1210 const auto& gradCalc = elemCtx.gradientCalculator();
1213 const auto& scvf = elemCtx.stencil(timeIdx).interiorFace(scvfIdx);
1214 const auto& faceNormal = scvf.normal();
1216 const unsigned i = scvf.interiorIndex();
1217 const unsigned j = scvf.exteriorIndex();
1220 DimEvalVector solventPGrad;
1222 gradCalc.calculateGradient(solventPGrad,
1226 Valgrind::CheckDefined(solventPGrad);
1229 if (Parameters::Get<Parameters::EnableGravity>()) {
1232 const auto& gIn = elemCtx.problem().gravity(elemCtx, i, timeIdx);
1233 const auto& gEx = elemCtx.problem().gravity(elemCtx, j, timeIdx);
1235 const auto& intQuantsIn = elemCtx.intensiveQuantities(i, timeIdx);
1236 const auto& intQuantsEx = elemCtx.intensiveQuantities(j, timeIdx);
1238 const auto& posIn = elemCtx.pos(i, timeIdx);
1239 const auto& posEx = elemCtx.pos(j, timeIdx);
1240 const auto& posFace = scvf.integrationPos();
1243 DimVector distVecIn(posIn);
1244 DimVector distVecEx(posEx);
1245 DimVector distVecTotal(posEx);
1247 distVecIn -= posFace;
1248 distVecEx -= posFace;
1249 distVecTotal -= posIn;
1250 const Scalar absDistTotalSquared = distVecTotal.two_norm2();
1253 auto rhoIn = intQuantsIn.solventDensity();
1254 auto pStatIn = - rhoIn * (gIn * distVecIn);
1258 Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1259 Scalar pStatEx = - rhoEx * (gEx * distVecEx);
1265 DimEvalVector f(distVecTotal);
1266 f *= (pStatEx - pStatIn) / absDistTotalSquared;
1269 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
1270 solventPGrad[dimIdx] += f[dimIdx];
1272 if (!isfinite(solventPGrad[dimIdx])) {
1273 throw NumericalProblem(
"Non-finite potential gradient for solvent 'phase'");
1279 Evaluation solventPGradNormal = 0.0;
1280 for (
unsigned dimIdx = 0; dimIdx < faceNormal.size(); ++dimIdx) {
1281 solventPGradNormal += solventPGrad[dimIdx]*faceNormal[dimIdx];
1284 if (solventPGradNormal > 0) {
1285 solventUpstreamDofIdx_ = j;
1286 solventDownstreamDofIdx_ = i;
1289 solventUpstreamDofIdx_ = i;
1290 solventDownstreamDofIdx_ = j;
1293 const auto& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1299 if (solventUpstreamDofIdx_ == i) {
1300 solventVolumeFlux_ = solventPGradNormal * up.solventMobility();
1303 solventVolumeFlux_ = solventPGradNormal * scalarValue(up.solventMobility());
1311 template <
class Dummy =
bool>
1317 const ExtensiveQuantities& extQuants = asImp_();
1319 const unsigned interiorDofIdx = extQuants.interiorIndex();
1320 const unsigned exteriorDofIdx = extQuants.exteriorIndex();
1321 assert(interiorDofIdx != exteriorDofIdx);
1323 const auto& intQuantsIn = elemCtx.intensiveQuantities(interiorDofIdx, timeIdx);
1324 const auto& intQuantsEx = elemCtx.intensiveQuantities(exteriorDofIdx, timeIdx);
1326 const unsigned I = elemCtx.globalSpaceIndex(interiorDofIdx, timeIdx);
1327 const unsigned J = elemCtx.globalSpaceIndex(exteriorDofIdx, timeIdx);
1329 const Scalar thpres = elemCtx.problem().thresholdPressure(I, J);
1330 const Scalar trans = elemCtx.problem().transmissibility(elemCtx, interiorDofIdx, exteriorDofIdx);
1331 const Scalar g = elemCtx.problem().gravity()[dimWorld - 1];
1333 const Scalar zIn = elemCtx.problem().dofCenterDepth(elemCtx, interiorDofIdx, timeIdx);
1334 const Scalar zEx = elemCtx.problem().dofCenterDepth(elemCtx, exteriorDofIdx, timeIdx);
1335 const Scalar distZ = zIn - zEx;
1337 const Evaluation& rhoIn = intQuantsIn.solventDensity();
1338 const Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1339 const Evaluation& rhoAvg = rhoIn * 0.5 + rhoEx * 0.5;
1341 const Evaluation& pressureInterior = intQuantsIn.fluidState().pressure(gasPhaseIdx);
1342 Evaluation pressureExterior = Toolbox::value(intQuantsEx.fluidState().pressure(gasPhaseIdx));
1343 pressureExterior += distZ * g * rhoAvg;
1345 Evaluation pressureDiffSolvent = pressureExterior - pressureInterior;
1346 if (std::abs(scalarValue(pressureDiffSolvent)) > thpres) {
1347 if (pressureDiffSolvent < 0.0) {
1348 pressureDiffSolvent += thpres;
1351 pressureDiffSolvent -= thpres;
1355 pressureDiffSolvent = 0.0;
1358 if (pressureDiffSolvent > 0.0) {
1359 solventUpstreamDofIdx_ = exteriorDofIdx;
1360 solventDownstreamDofIdx_ = interiorDofIdx;
1362 else if (pressureDiffSolvent < 0.0) {
1363 solventUpstreamDofIdx_ = interiorDofIdx;
1364 solventDownstreamDofIdx_ = exteriorDofIdx;
1370 solventUpstreamDofIdx_ = std::min(interiorDofIdx, exteriorDofIdx);
1371 solventDownstreamDofIdx_ = std::max(interiorDofIdx, exteriorDofIdx);
1372 solventVolumeFlux_ = 0.0;
1376 const Scalar faceArea = elemCtx.stencil(timeIdx).interiorFace(scvfIdx).area();
1377 const IntensiveQuantities& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1378 if (solventUpstreamDofIdx_ == interiorDofIdx) {
1379 solventVolumeFlux_ = up.solventMobility() * (-trans / faceArea) * pressureDiffSolvent;
1382 solventVolumeFlux_ =
1383 scalarValue(up.solventMobility()) * (-trans / faceArea) * pressureDiffSolvent;
1388 {
return solventUpstreamDofIdx_; }
1391 {
return solventDownstreamDofIdx_; }
1394 {
return solventVolumeFlux_; }
1400 Implementation& asImp_()
1401 {
return *
static_cast<Implementation*
>(
this); }
1403 Evaluation solventVolumeFlux_;
1404 unsigned solventUpstreamDofIdx_;
1405 unsigned solventDownstreamDofIdx_;
1408template <
class TypeTag>
1426 {
throw std::runtime_error(
"solventUpstreamIndex() called but solvents are disabled"); }
1429 {
throw std::runtime_error(
"solventDownstreamIndex() called but solvents are disabled"); }
1432 {
throw std::runtime_error(
"solventVolumeFlux() called but solvents are disabled"); }
1435 {
throw std::runtime_error(
"setSolventVolumeFlux() called but solvents are disabled"); }
Declares the properties required by the black oil model.
Contains the parameters required to extend the black-oil model by solvents.
unsigned solventUpstreamIndex() const
Definition: blackoilsolventmodules.hh:1425
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1428
void setSolventVolumeFlux(const Evaluation &)
Definition: blackoilsolventmodules.hh:1434
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1431
void updateVolumeFluxPerm(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1415
void updateVolumeFluxTrans(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1420
Provides the solvent specific extensive quantities to the generic black-oil module's extensive quanti...
Definition: blackoilsolventmodules.hh:1180
unsigned solventUpstreamIndex() const
Definition: blackoilsolventmodules.hh:1387
void updateVolumeFluxPerm(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Method which calculates the volume flux of the polymer "phase" using the pressure potential gradient ...
Definition: blackoilsolventmodules.hh:1206
void updateVolumeFluxTrans(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Method which calculates the volume flux of the polymer "phase" using the gas pressure potential diffe...
Definition: blackoilsolventmodules.hh:1313
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1390
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1393
void setSolventVolumeFlux(const Evaluation &solventVolumeFlux)
Definition: blackoilsolventmodules.hh:1396
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:1155
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:1164
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:1152
void solventPreSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1134
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:1158
Scalar solventRefDensity() const
Definition: blackoilsolventmodules.hh:1167
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:1161
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:1149
void solventPvtUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1144
void solventPostSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1139
void solventPreSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called before the saturation functions are doing their magic.
Definition: blackoilsolventmodules.hh:580
Implementation & asImp_()
Definition: blackoilsolventmodules.hh:1112
Evaluation rsSolw_
Definition: blackoilsolventmodules.hh:1117
Evaluation solventViscosity_
Definition: blackoilsolventmodules.hh:1119
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:866
void solventPvtUpdate_(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Update the intensive PVT properties needed to handle solvents from the primary variables.
Definition: blackoilsolventmodules.hh:793
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:860
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:872
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:863
void solventPostSatFuncUpdate_(const Problem &problem, const PrimaryVariables &priVars, const unsigned globalSpaceIdx, const unsigned timeIdx, const LinearizationType linearizationType)
Definition: blackoilsolventmodules.hh:641
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:869
Evaluation hydrocarbonSaturation_
Definition: blackoilsolventmodules.hh:1115
void solventPreSatFuncUpdate_(const PrimaryVariables &priVars, const unsigned timeIdx, const LinearizationType linearizationType)
Definition: blackoilsolventmodules.hh:588
void solventPostSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called after the saturation functions have been doing their magic.
Definition: blackoilsolventmodules.hh:617
Evaluation solventSaturation_
Definition: blackoilsolventmodules.hh:1116
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:857
Evaluation solventDensity_
Definition: blackoilsolventmodules.hh:1118
Evaluation solventInvFormationVolumeFactor_
Definition: blackoilsolventmodules.hh:1121
Scalar solventRefDensity_
Definition: blackoilsolventmodules.hh:1123
Scalar solventRefDensity() const
Definition: blackoilsolventmodules.hh:876
Evaluation solventMobility_
Definition: blackoilsolventmodules.hh:1120
Provides the volumetric quantities required for the equations needed by the solvents extension of the...
Definition: blackoilsolventmodules.hh:541
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilsolventmodules.hh:68
static void setSolventPvt(const SolventPvt &value)
Specify the solvent PVT of a all PVT regions.
Definition: blackoilsolventmodules.hh:107
static const TabulatedFunction & misc(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:411
static const TabulatedFunction & sof2Krn(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:401
static Scalar computeResidualError(const EqVector &resid)
Return how much a residual is considered an error.
Definition: blackoilsolventmodules.hh:332
static Scalar tlMixParamDensity(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:492
static const TabulatedFunction & ssfnKrg(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:381
static const BrineH2Pvt & brineH2Pvt()
Definition: blackoilsolventmodules.hh:377
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:151
static void setParams(BlackOilSolventParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilsolventmodules.hh:101
static std::string eqName(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:169
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:176
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilsolventmodules.hh:339
static bool isSolubleInWater()
Definition: blackoilsolventmodules.hh:523
static bool isCO2Sol()
Definition: blackoilsolventmodules.hh:526
static const TabulatedFunction & sgcwmis(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:461
static Scalar tlMixParamViscosity(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:481
static void registerOutputModules(Model &model, Simulator &simulator)
Register all solvent specific VTK and ECL output modules.
Definition: blackoilsolventmodules.hh:126
static const TabulatedFunction & tlPMixTable(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:471
static const TabulatedFunction & sorwmis(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:451
static Value solubilityLimit(unsigned pvtIdx, const Value &temperature, const Value &pressure, const Value &saltConcentration)
Definition: blackoilsolventmodules.hh:505
static const TabulatedFunction & pmisc(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:421
static void setIsMiscible(const bool isMiscible)
Definition: blackoilsolventmodules.hh:110
static void registerParameters()
Register all run-time parameters for the black-oil solvent module.
Definition: blackoilsolventmodules.hh:116
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilsolventmodules.hh:320
static const SolventPvt & solventPvt()
Definition: blackoilsolventmodules.hh:365
static bool isH2Sol()
Definition: blackoilsolventmodules.hh:529
static bool isMiscible()
Definition: blackoilsolventmodules.hh:501
static const H2GasPvt & h2GasPvt()
Definition: blackoilsolventmodules.hh:371
static void updatePrimaryVars(PrimaryVariables &newPv, const PrimaryVariables &oldPv, const EqVector &delta)
Do a Newton-Raphson update the primary variables of the solvents.
Definition: blackoilsolventmodules.hh:307
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar solventSaturation, Scalar solventRsw)
Assign the solvent specific primary variables to a PrimaryVariables object.
Definition: blackoilsolventmodules.hh:286
static const TabulatedFunction & msfnKrsg(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:431
static const TabulatedFunction & ssfnKrs(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:391
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:144
static bool eqApplies(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:159
static const TabulatedFunction & msfnKro(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:441
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilsolventmodules.hh:185
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilsolventmodules.hh:350
static const Co2GasPvt & co2GasPvt()
Definition: blackoilsolventmodules.hh:368
static const BrineCo2Pvt & brineCo2Pvt()
Definition: blackoilsolventmodules.hh:374
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:218
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:134
Callback class for a phase pressure.
Definition: quantitycallbacks.hh:85
void setPhaseIndex(unsigned phaseIdx)
Set the index of the fluid phase for which the pressure should be returned.
Definition: quantitycallbacks.hh:109
VTK output module for the black oil model's solvent related quantities.
Definition: vtkblackoilsolventmodule.hpp:54
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilsolventmodule.hpp:84
Defines the common parameters for the porous medium multi-phase models.
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
This method contains all callback classes for quantities that are required by some extensive quantiti...
Struct holding the parameters for the BlackOilSolventModule class.
Definition: blackoilsolventparams.hpp:49
::Opm::SolventPvt< Scalar > SolventPvt
Definition: blackoilsolventparams.hpp:54
::Opm::H2GasPvt< Scalar > H2GasPvt
Definition: blackoilsolventparams.hpp:53
::Opm::BrineCo2Pvt< Scalar > BrineCo2Pvt
Definition: blackoilsolventparams.hpp:50
::Opm::Co2GasPvt< Scalar > Co2GasPvt
Definition: blackoilsolventparams.hpp:52
Tabulated1DFunction< Scalar > TabulatedFunction
Definition: blackoilsolventparams.hpp:55
::Opm::BrineH2Pvt< Scalar > BrineH2Pvt
Definition: blackoilsolventparams.hpp:51
Definition: linearizationtype.hh:34