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>
64template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
79 using Toolbox = MathToolbox<Evaluation>;
88 static constexpr unsigned solventSaturationIdx = Indices::solventSaturationIdx;
89 static constexpr unsigned contiSolventEqIdx = Indices::contiSolventEqIdx;
90 static constexpr unsigned enableSolvent = enableSolventV;
91 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
92 static constexpr unsigned numPhases = FluidSystem::numPhases;
93 static constexpr bool blackoilConserveSurfaceVolume =
94 getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>();
95 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
100 { params_ = params; }
106 { params_.solventPvt_ = value; }
116 if constexpr (enableSolvent) {
125 Simulator& simulator)
127 if constexpr (enableSolvent) {
134 if constexpr (enableSolvent) {
135 return pvIdx == solventSaturationIdx;
146 return "saturation_solvent";
154 return static_cast<Scalar
>(1.0);
159 if constexpr (enableSolvent) {
160 return eqIdx == contiSolventEqIdx;
167 static std::string
eqName([[maybe_unused]]
unsigned eqIdx)
171 return "conti^solvent";
174 static Scalar
eqWeight([[maybe_unused]]
unsigned eqIdx)
179 return static_cast<Scalar
>(1.0);
182 template <
class LhsEval>
183 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
184 const IntensiveQuantities& intQuants)
186 if constexpr (enableSolvent) {
187 if constexpr (blackoilConserveSurfaceVolume) {
188 storage[contiSolventEqIdx] +=
189 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
190 Toolbox::template decay<LhsEval>(intQuants.solventSaturation()) *
191 Toolbox::template decay<LhsEval>(intQuants.solventInverseFormationVolumeFactor());
193 storage[contiSolventEqIdx] +=
194 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
195 Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(waterPhaseIdx)) *
196 Toolbox::template decay<LhsEval>(intQuants.fluidState().invB(waterPhaseIdx)) *
197 Toolbox::template decay<LhsEval>(intQuants.rsSolw());
201 storage[contiSolventEqIdx] +=
202 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
203 Toolbox::template decay<LhsEval>(intQuants.solventSaturation()) *
204 Toolbox::template decay<LhsEval>(intQuants.solventDensity());
206 storage[contiSolventEqIdx] +=
207 Toolbox::template decay<LhsEval>(intQuants.porosity()) *
208 Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(waterPhaseIdx)) *
209 Toolbox::template decay<LhsEval>(intQuants.fluidState().density(waterPhaseIdx)) *
210 Toolbox::template decay<LhsEval>(intQuants.rsSolw());
217 [[maybe_unused]]
const ElementContext& elemCtx,
218 [[maybe_unused]]
unsigned scvfIdx,
219 [[maybe_unused]]
unsigned timeIdx)
222 if constexpr (enableSolvent) {
223 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
225 const unsigned upIdx = extQuants.solventUpstreamIndex();
226 const unsigned inIdx = extQuants.interiorIndex();
227 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
229 if constexpr (blackoilConserveSurfaceVolume) {
230 if (upIdx == inIdx) {
231 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() *
232 up.solventInverseFormationVolumeFactor();
235 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() *
236 decay<Scalar>(up.solventInverseFormationVolumeFactor());
241 if (upIdx == inIdx) {
242 flux[contiSolventEqIdx] +=
243 extQuants.volumeFlux(waterPhaseIdx) *
244 up.fluidState().invB(waterPhaseIdx) *
248 flux[contiSolventEqIdx] +=
249 extQuants.volumeFlux(waterPhaseIdx) *
250 decay<Scalar>(up.fluidState().invB(waterPhaseIdx)) *
251 decay<Scalar>(up.rsSolw());
256 if (upIdx == inIdx) {
257 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() * up.solventDensity();
260 flux[contiSolventEqIdx] = extQuants.solventVolumeFlux() * decay<Scalar>(up.solventDensity());
264 if (upIdx == inIdx) {
265 flux[contiSolventEqIdx] +=
266 extQuants.volumeFlux(waterPhaseIdx) *
267 up.fluidState().density(waterPhaseIdx) *
271 flux[contiSolventEqIdx] +=
272 extQuants.volumeFlux(waterPhaseIdx) *
273 decay<Scalar>(up.fluidState().density(waterPhaseIdx)) *
274 decay<Scalar>(up.rsSolw());
285 Scalar solventSaturation,
288 if constexpr (!enableSolvent) {
289 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Disabled);
294 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Ss);
295 priVars[solventSaturationIdx] = solventSaturation;
297 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Rsolw);
298 priVars[solventSaturationIdx] = solventRsw;
306 const PrimaryVariables& oldPv,
307 const EqVector& delta)
309 if constexpr (enableSolvent) {
311 newPv[solventSaturationIdx] = oldPv[solventSaturationIdx] - delta[solventSaturationIdx];
324 return static_cast<Scalar
>(0.0);
333 return std::abs(Toolbox::scalarValue(resid[contiSolventEqIdx]));
336 template <
class DofEntity>
337 static void serializeEntity(
const Model& model, std::ostream& outstream,
const DofEntity& dof)
339 if constexpr (enableSolvent) {
340 const unsigned dofIdx = model.dofMapper().index(dof);
342 const PrimaryVariables& priVars = model.solution(0)[dofIdx];
343 outstream << priVars[solventSaturationIdx];
347 template <
class DofEntity>
350 if constexpr (enableSolvent) {
351 const unsigned dofIdx = model.dofMapper().index(dof);
353 PrimaryVariables& priVars0 = model.solution(0)[dofIdx];
354 PrimaryVariables& priVars1 = model.solution(1)[dofIdx];
356 instream >> priVars0[solventSaturationIdx];
359 priVars1 = priVars0[solventSaturationIdx];
364 {
return params_.solventPvt_; }
367 {
return params_.co2GasPvt_; }
370 {
return params_.h2GasPvt_; }
373 {
return params_.brineCo2Pvt_; }
376 {
return params_.brineH2Pvt_; }
378 static const TabulatedFunction&
ssfnKrg(
const ElementContext& elemCtx,
382 const unsigned satnumRegionIdx =
383 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
384 return params_.ssfnKrg_[satnumRegionIdx];
387 static const TabulatedFunction&
ssfnKrs(
const ElementContext& elemCtx,
391 const unsigned satnumRegionIdx =
392 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
393 return params_.ssfnKrs_[satnumRegionIdx];
396 static const TabulatedFunction&
sof2Krn(
const ElementContext& elemCtx,
400 const unsigned satnumRegionIdx =
401 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
402 return params_.sof2Krn_[satnumRegionIdx];
405 static const TabulatedFunction&
misc(
const ElementContext& elemCtx,
409 const unsigned miscnumRegionIdx =
410 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
411 return params_.misc_[miscnumRegionIdx];
414 static const TabulatedFunction&
pmisc(
const ElementContext& elemCtx,
418 const unsigned miscnumRegionIdx =
419 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
420 return params_.pmisc_[miscnumRegionIdx];
423 static const TabulatedFunction&
msfnKrsg(
const ElementContext& elemCtx,
427 const unsigned satnumRegionIdx =
428 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
429 return params_.msfnKrsg_[satnumRegionIdx];
432 static const TabulatedFunction&
msfnKro(
const ElementContext& elemCtx,
436 const unsigned satnumRegionIdx =
437 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
438 return params_.msfnKro_[satnumRegionIdx];
441 static const TabulatedFunction&
sorwmis(
const ElementContext& elemCtx,
445 const unsigned miscnumRegionIdx =
446 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
447 return params_.sorwmis_[miscnumRegionIdx];
450 static const TabulatedFunction&
sgcwmis(
const ElementContext& elemCtx,
454 const unsigned miscnumRegionIdx =
455 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
456 return params_.sgcwmis_[miscnumRegionIdx];
459 static const TabulatedFunction&
tlPMixTable(
const ElementContext& elemCtx,
463 const unsigned miscnumRegionIdx =
464 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
465 return params_.tlPMixTable_[miscnumRegionIdx];
472 const unsigned miscnumRegionIdx =
473 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
474 return params_.tlMixParamViscosity_[miscnumRegionIdx];
481 const unsigned miscnumRegionIdx =
482 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
483 return params_.tlMixParamDensity_[miscnumRegionIdx];
487 {
return params_.isMiscible_; }
489 template <
class Value>
491 const Value& pressure,
const Value& saltConcentration)
499 return brineCo2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature,
500 pressure, saltConcentration);
503 return brineH2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature,
504 pressure, saltConcentration);
509 {
return params_.rsSolw_active_; }
512 {
return params_.co2sol_; }
515 {
return params_.h2sol_; }
521template <
class TypeTag,
bool enableSolventV>
522BlackOilSolventParams<typename BlackOilSolventModule<TypeTag, enableSolventV>::Scalar>
523BlackOilSolventModule<TypeTag, enableSolventV>::params_;
525template <
class TypeTag,
bool enableSolventV>
535template <
class TypeTag>
550 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
551 static constexpr int solventSaturationIdx = Indices::solventSaturationIdx;
552 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
553 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
554 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
555 static constexpr double cutOff = 1e-12;
568 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
570 auto& fs = asImp_().fluidState_;
571 solventSaturation_ = 0.0;
572 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
573 solventSaturation_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx,
574 elemCtx.linearizationType());
577 hydrocarbonSaturation_ = fs.saturation(gasPhaseIdx);
580 if (solventSaturation().value() < cutOff) {
586 fs.setSaturation(gasPhaseIdx, hydrocarbonSaturation_ + solventSaturation_);
602 auto& fs = asImp_().fluidState_;
603 fs.setSaturation(gasPhaseIdx, hydrocarbonSaturation_);
607 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
608 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
609 rsSolw_ = SolventModule::solubilityLimit(asImp_().pvtRegionIndex(),
610 fs.temperature(waterPhaseIdx),
611 fs.pressure(waterPhaseIdx),
612 fs.saltConcentration());
613 }
else if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Rsolw) {
614 rsSolw_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx,
615 elemCtx.linearizationType());
618 solventMobility_ = 0.0;
621 if (solventSaturation().value() < cutOff) {
626 if (SolventModule::isMiscible()) {
627 const Evaluation& p =
628 FluidSystem::phaseIsActive(oilPhaseIdx)
629 ? fs.pressure(oilPhaseIdx)
630 : fs.pressure(gasPhaseIdx);
631 const Evaluation pmisc =
632 SolventModule::pmisc(elemCtx, dofIdx, timeIdx).eval(p,
true);
633 const Evaluation& pgImisc = fs.pressure(gasPhaseIdx);
636 const auto& problem = elemCtx.problem();
637 Evaluation pgMisc = 0.0;
638 std::array<Evaluation, numPhases> pC;
639 const auto& materialParams = problem.materialLawParams(elemCtx, dofIdx, timeIdx);
640 MaterialLaw::capillaryPressures(pC, materialParams, fs);
643 const auto linearizationType = elemCtx.linearizationType();
644 if (priVars.primaryVarsMeaningPressure() == PrimaryVariables::PressureMeaning::Pg) {
645 pgMisc = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx,
649 const Evaluation& po =
650 priVars.makeEvaluation(Indices::pressureSwitchIdx,
651 timeIdx, linearizationType);
652 pgMisc = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
655 fs.setPressure(gasPhaseIdx, pmisc * pgMisc + (1.0 - pmisc) * pgImisc);
658 const Evaluation gasSolventSat = hydrocarbonSaturation_ + solventSaturation_;
660 if (gasSolventSat.value() < cutOff) {
664 const Evaluation Fhydgas = hydrocarbonSaturation_ / gasSolventSat;
665 const Evaluation Fsolgas = solventSaturation_ / gasSolventSat;
668 if (SolventModule::isMiscible() && FluidSystem::phaseIsActive(oilPhaseIdx)) {
669 const auto& misc = SolventModule::misc(elemCtx, dofIdx, timeIdx);
670 const auto& pmisc = SolventModule::pmisc(elemCtx, dofIdx, timeIdx);
671 const Evaluation& p =
672 FluidSystem::phaseIsActive(oilPhaseIdx)
673 ? fs.pressure(oilPhaseIdx)
674 : fs.pressure(gasPhaseIdx);
675 const Evaluation miscibility =
676 misc.eval(Fsolgas,
true) *
680 const unsigned cellIdx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
681 const auto& materialLawManager = elemCtx.problem().materialLawManager();
682 const auto& scaledDrainageInfo =
683 materialLawManager->oilWaterScaledEpsInfoDrainage(cellIdx);
685 const Scalar sogcr = scaledDrainageInfo.Sogcr;
686 Evaluation sor = sogcr;
687 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
688 const Evaluation& sw = fs.saturation(waterPhaseIdx);
689 const auto& sorwmis = SolventModule::sorwmis(elemCtx, dofIdx, timeIdx);
691 sorwmis.eval(sw,
true) + (1.0 - miscibility) * sogcr;
693 const Scalar sgcr = scaledDrainageInfo.Sgcr;
694 Evaluation sgc = sgcr;
695 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
696 const Evaluation& sw = fs.saturation(waterPhaseIdx);
697 const auto& sgcwmis = SolventModule::sgcwmis(elemCtx, dofIdx, timeIdx);
698 sgc = miscibility * sgcwmis.eval(sw,
true) + (1.0 - miscibility) * sgcr;
701 Evaluation oilGasSolventSat = gasSolventSat;
702 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
703 oilGasSolventSat += fs.saturation(oilPhaseIdx);
705 const Evaluation zero = 0.0;
706 const Evaluation oilGasSolventEffSat = std::max(oilGasSolventSat - sor - sgc, zero);
708 Evaluation F_totalGas = 0.0;
709 if (oilGasSolventEffSat.value() > cutOff) {
710 const Evaluation gasSolventEffSat = std::max(gasSolventSat - sgc, zero);
711 F_totalGas = gasSolventEffSat / oilGasSolventEffSat;
713 const auto& msfnKro = SolventModule::msfnKro(elemCtx, dofIdx, timeIdx);
714 const auto& msfnKrsg = SolventModule::msfnKrsg(elemCtx, dofIdx, timeIdx);
715 const auto& sof2Krn = SolventModule::sof2Krn(elemCtx, dofIdx, timeIdx);
717 const Evaluation mkrgt = msfnKrsg.eval(F_totalGas,
true) *
718 sof2Krn.eval(oilGasSolventSat,
true);
719 const Evaluation mkro = msfnKro.eval(F_totalGas,
true) *
720 sof2Krn.eval(oilGasSolventSat,
true);
722 Evaluation& kro = asImp_().mobility_[oilPhaseIdx];
723 Evaluation& krg = asImp_().mobility_[gasPhaseIdx];
726 krg *= 1.0 - miscibility;
727 krg += miscibility * mkrgt;
728 kro *= 1.0 - miscibility;
729 kro += miscibility * mkro;
733 const auto& ssfnKrg = SolventModule::ssfnKrg(elemCtx, dofIdx, timeIdx);
734 const auto& ssfnKrs = SolventModule::ssfnKrs(elemCtx, dofIdx, timeIdx);
736 Evaluation& krg = asImp_().mobility_[gasPhaseIdx];
737 solventMobility_ = krg * ssfnKrs.eval(Fsolgas,
true);
738 krg *= ssfnKrg.eval(Fhydgas,
true);
751 const auto& iq = asImp_();
752 const unsigned pvtRegionIdx = iq.pvtRegionIndex();
753 const Evaluation& T = iq.fluidState().temperature(gasPhaseIdx);
754 const Evaluation& p = iq.fluidState().pressure(gasPhaseIdx);
756 const Evaluation rv = 0.0;
757 const Evaluation rvw = 0.0;
758 if (SolventModule::isCO2Sol() || SolventModule::isH2Sol() ){
759 if (SolventModule::isCO2Sol()) {
760 const auto& co2gasPvt = SolventModule::co2GasPvt();
761 solventInvFormationVolumeFactor_ =
762 co2gasPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rv, rvw);
763 solventRefDensity_ = co2gasPvt.gasReferenceDensity(pvtRegionIdx);
764 solventViscosity_ = co2gasPvt.viscosity(pvtRegionIdx, T, p, rv, rvw);
766 const auto& brineCo2Pvt = SolventModule::brineCo2Pvt();
767 auto& fs = asImp_().fluidState_;
768 const auto& bw = brineCo2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rsSolw());
770 const auto denw = bw * brineCo2Pvt.waterReferenceDensity(pvtRegionIdx) +
771 rsSolw() * bw * brineCo2Pvt.gasReferenceDensity(pvtRegionIdx);
772 fs.setDensity(waterPhaseIdx, denw);
773 fs.setInvB(waterPhaseIdx, bw);
774 Evaluation& mobw = asImp_().mobility_[waterPhaseIdx];
775 const auto& muWat = fs.viscosity(waterPhaseIdx);
776 const auto& muWatEff = brineCo2Pvt.viscosity(pvtRegionIdx, T, p, rsSolw());
777 mobw *= muWat / muWatEff;
779 const auto& h2gasPvt = SolventModule::h2GasPvt();
780 solventInvFormationVolumeFactor_ =
781 h2gasPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rv, rvw);
782 solventRefDensity_ = h2gasPvt.gasReferenceDensity(pvtRegionIdx);
783 solventViscosity_ = h2gasPvt.viscosity(pvtRegionIdx, T, p, rv, rvw);
785 const auto& brineH2Pvt = SolventModule::brineH2Pvt();
786 auto& fs = asImp_().fluidState_;
787 const auto& bw = brineH2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rsSolw());
789 const auto denw = bw * brineH2Pvt.waterReferenceDensity(pvtRegionIdx) +
790 rsSolw() * bw * brineH2Pvt.gasReferenceDensity(pvtRegionIdx);
791 fs.setDensity(waterPhaseIdx, denw);
792 fs.setInvB(waterPhaseIdx, bw);
793 Evaluation& mobw = asImp_().mobility_[waterPhaseIdx];
794 const auto& muWat = fs.viscosity(waterPhaseIdx);
795 const auto& muWatEff = brineH2Pvt.viscosity(pvtRegionIdx, T, p, rsSolw());
796 mobw *= muWat / muWatEff;
799 const auto& solventPvt = SolventModule::solventPvt();
800 solventInvFormationVolumeFactor_ = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
801 solventRefDensity_ = solventPvt.referenceDensity(pvtRegionIdx);
802 solventViscosity_ = solventPvt.viscosity(pvtRegionIdx, T, p);
805 solventDensity_ = solventInvFormationVolumeFactor_*solventRefDensity_;
806 effectiveProperties(elemCtx, scvIdx, timeIdx);
808 solventMobility_ /= solventViscosity_;
812 {
return solventSaturation_; }
818 {
return solventDensity_; }
821 {
return solventViscosity_; }
824 {
return solventMobility_; }
827 {
return solventInvFormationVolumeFactor_; }
831 {
return solventRefDensity_; }
836 void effectiveProperties(
const ElementContext& elemCtx,
840 if (!SolventModule::isMiscible()) {
846 if (solventSaturation() < cutOff) {
852 auto& fs = asImp_().fluidState_;
855 Evaluation oilEffSat = 0.0;
856 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
857 oilEffSat = fs.saturation(oilPhaseIdx);
859 Evaluation gasEffSat = fs.saturation(gasPhaseIdx);
860 Evaluation solventEffSat = solventSaturation();
861 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
862 const auto& sorwmis = SolventModule::sorwmis(elemCtx, scvIdx, timeIdx);
863 const auto& sgcwmis = SolventModule::sgcwmis(elemCtx, scvIdx, timeIdx);
864 const Evaluation zero = 0.0;
865 const Evaluation& sw = fs.saturation(waterPhaseIdx);
866 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
867 oilEffSat = std::max(oilEffSat - sorwmis.eval(sw,
true), zero);
869 gasEffSat = std::max(gasEffSat - sgcwmis.eval(sw,
true), zero);
870 solventEffSat = std::max(solventEffSat - sgcwmis.eval(sw,
true), zero);
872 const Evaluation oilGasSolventEffSat = oilEffSat + gasEffSat + solventEffSat;
873 const Evaluation oilSolventEffSat = oilEffSat + solventEffSat;
874 const Evaluation solventGasEffSat = solventEffSat + gasEffSat;
881 const Evaluation& p =
882 FluidSystem::phaseIsActive(oilPhaseIdx)
883 ? fs.pressure(oilPhaseIdx)
884 : fs.pressure(gasPhaseIdx);
886 const auto& pmiscTable = SolventModule::pmisc(elemCtx, scvIdx, timeIdx);
887 const Evaluation pmisc = pmiscTable.eval(p,
true);
888 const auto& tlPMixTable = SolventModule::tlPMixTable(elemCtx, scvIdx, timeIdx);
889 const Evaluation tlMixParamMu = SolventModule::tlMixParamViscosity(elemCtx, scvIdx, timeIdx) *
890 tlPMixTable.eval(p,
true);
892 const Evaluation& muGas = fs.viscosity(gasPhaseIdx);
893 const Evaluation& muSolvent = solventViscosity_;
895 assert(muGas.value() > 0);
896 assert(muSolvent.value() > 0);
897 const Evaluation muGasPow = pow(muGas, 0.25);
898 const Evaluation muSolventPow = pow(muSolvent, 0.25);
900 Evaluation muMixSolventGas = muGas;
901 if (solventGasEffSat > cutOff) {
902 muMixSolventGas *= muSolvent / pow(((gasEffSat / solventGasEffSat) * muSolventPow) +
903 ((solventEffSat / solventGasEffSat) * muGasPow), 4.0);
906 Evaluation muOil = 1.0;
907 Evaluation muOilPow = 1.0;
908 Evaluation muMixOilSolvent = 1.0;
909 Evaluation muOilEff = 1.0;
910 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
911 muOil = fs.viscosity(oilPhaseIdx);
912 assert(muOil.value() > 0);
913 muOilPow = pow(muOil, 0.25);
914 muMixOilSolvent = muOil;
915 if (oilSolventEffSat > cutOff) {
916 muMixOilSolvent *= muSolvent / pow(((oilEffSat / oilSolventEffSat) * muSolventPow) +
917 ((solventEffSat / oilSolventEffSat) * muOilPow), 4.0);
920 muOilEff = pow(muOil,1.0 - tlMixParamMu) * pow(muMixOilSolvent, tlMixParamMu);
922 Evaluation muMixSolventGasOil = muOil;
923 if (oilGasSolventEffSat > cutOff) {
924 muMixSolventGasOil *= muSolvent * muGas / pow(((oilEffSat / oilGasSolventEffSat) *
925 muSolventPow * muGasPow) +
926 ((solventEffSat / oilGasSolventEffSat) *
927 muOilPow * muGasPow) +
928 ((gasEffSat / oilGasSolventEffSat) *
929 muSolventPow * muOilPow), 4.0);
932 Evaluation muGasEff = pow(muGas,1.0 - tlMixParamMu) * pow(muMixSolventGas, tlMixParamMu);
933 const Evaluation muSolventEff = pow(muSolvent,1.0 - tlMixParamMu) * pow(muMixSolventGasOil, tlMixParamMu);
936 const Evaluation& rhoGas = fs.density(gasPhaseIdx);
937 Evaluation rhoOil = 0.0;
938 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
939 rhoOil = fs.density(oilPhaseIdx);
942 const Evaluation& rhoSolvent = solventDensity_;
948 const Evaluation tlMixParamRho = SolventModule::tlMixParamDensity(elemCtx, scvIdx, timeIdx) *
949 tlPMixTable.eval(p,
true);
953 const Evaluation muOilEffPow = pow(pow(muOil, 1.0 - tlMixParamRho) *
954 pow(muMixOilSolvent, tlMixParamRho), 0.25);
955 const Evaluation muGasEffPow = pow(pow(muGas, 1.0 - tlMixParamRho) *
956 pow(muMixSolventGas, tlMixParamRho), 0.25);
957 const Evaluation muSolventEffPow = pow(pow(muSolvent, 1.0 - tlMixParamRho) *
958 pow(muMixSolventGasOil, tlMixParamRho), 0.25);
960 const Evaluation oilGasEffSaturation = oilEffSat + gasEffSat;
961 Evaluation sof = 0.0;
962 Evaluation sgf = 0.0;
963 if (oilGasEffSaturation.value() > cutOff) {
964 sof = oilEffSat / oilGasEffSaturation;
965 sgf = gasEffSat / oilGasEffSaturation;
968 const Evaluation muSolventOilGasPow = muSolventPow * ((sgf * muOilPow) + (sof * muGasPow));
970 Evaluation rhoMixSolventGasOil = 0.0;
971 if (oilGasSolventEffSat.value() > cutOff) {
972 rhoMixSolventGasOil = (rhoOil * oilEffSat / oilGasSolventEffSat) +
973 (rhoGas * gasEffSat / oilGasSolventEffSat) +
974 (rhoSolvent * solventEffSat / oilGasSolventEffSat);
977 Evaluation rhoGasEff = 0.0;
978 if (std::abs(muSolventPow.value() - muGasPow.value()) < cutOff) {
979 rhoGasEff = ((1.0 - tlMixParamRho) * rhoGas) +
980 (tlMixParamRho * rhoMixSolventGasOil);
983 const Evaluation solventGasEffFraction = (muGasPow * (muSolventPow - muGasEffPow)) /
984 (muGasEffPow * (muSolventPow - muGasPow));
985 rhoGasEff = (rhoGas * solventGasEffFraction) + (rhoSolvent * (1.0 - solventGasEffFraction));
988 Evaluation rhoSolventEff = 0.0;
989 if (std::abs((muSolventOilGasPow.value() - (muOilPow.value() * muGasPow.value()))) < cutOff) {
990 rhoSolventEff = ((1.0 - tlMixParamRho) * rhoSolvent) + (tlMixParamRho * rhoMixSolventGasOil);
993 const Evaluation sfraction_se = (muSolventOilGasPow -
994 muOilPow * muGasPow * muSolventPow / muSolventEffPow) /
995 (muSolventOilGasPow - (muOilPow * muGasPow));
996 rhoSolventEff = (rhoSolvent * sfraction_se) +
997 (rhoGas * sgf * (1.0 - sfraction_se)) +
998 (rhoOil * sof * (1.0 - sfraction_se));
1002 const unsigned pvtRegionIdx = asImp_().pvtRegionIndex();
1003 Evaluation bGasEff = rhoGasEff;
1004 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1005 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) +
1006 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) * fs.Rv());
1008 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1010 const Evaluation bSolventEff = rhoSolventEff / solventRefDensity();
1013 const Evaluation bg = fs.invB(gasPhaseIdx);
1014 const Evaluation bs = solventInverseFormationVolumeFactor();
1015 const Evaluation bg_eff = pmisc * bGasEff + (1.0 - pmisc) * bg;
1016 const Evaluation bs_eff = pmisc * bSolventEff + (1.0 - pmisc) * bs;
1019 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1020 fs.setDensity(gasPhaseIdx,
1022 (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) +
1023 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx)*fs.Rv()));
1025 fs.setDensity(gasPhaseIdx,
1026 bg_eff * FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1028 solventDensity_ = bs_eff * solventRefDensity();
1031 Evaluation& mobg = asImp_().mobility_[gasPhaseIdx];
1032 muGasEff = bg_eff / (pmisc * bGasEff / muGasEff + (1.0 - pmisc) * bg / muGas);
1033 mobg *= muGas / muGasEff;
1036 solventViscosity_ = bs_eff / (pmisc * bSolventEff / muSolventEff +
1037 (1.0 - pmisc) * bs / muSolvent);
1039 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1040 Evaluation rhoOilEff = 0.0;
1041 if (std::abs(muOilPow.value() - muSolventPow.value()) < cutOff) {
1042 rhoOilEff = ((1.0 - tlMixParamRho) * rhoOil) + (tlMixParamRho * rhoMixSolventGasOil);
1045 const Evaluation solventOilEffFraction = (muOilPow * (muOilEffPow - muSolventPow)) /
1046 (muOilEffPow * (muOilPow - muSolventPow));
1047 rhoOilEff = (rhoOil * solventOilEffFraction) + (rhoSolvent * (1.0 - solventOilEffFraction));
1049 const Evaluation bOilEff =
1050 rhoOilEff / (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) +
1051 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs());
1052 const Evaluation bo = fs.invB(oilPhaseIdx);
1053 const Evaluation bo_eff = pmisc * bOilEff + (1.0 - pmisc) * bo;
1054 fs.setDensity(oilPhaseIdx,
1055 bo_eff * (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) +
1056 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs()));
1059 Evaluation& mobo = asImp_().mobility_[oilPhaseIdx];
1060 muOilEff = bo_eff / (pmisc * bOilEff / muOilEff + (1.0 - pmisc) * bo / muOil);
1061 mobo *= muOil / muOilEff;
1067 {
return *
static_cast<Implementation*
>(
this); }
1080template <
class TypeTag>
1104 {
throw std::runtime_error(
"solventSaturation() called but solvents are disabled"); }
1107 {
throw std::runtime_error(
"rsSolw() called but solvents are disabled"); }
1110 {
throw std::runtime_error(
"solventDensity() called but solvents are disabled"); }
1113 {
throw std::runtime_error(
"solventViscosity() called but solvents are disabled"); }
1116 {
throw std::runtime_error(
"solventMobility() called but solvents are disabled"); }
1119 {
throw std::runtime_error(
"solventInverseFormationVolumeFactor() called but solvents are disabled"); }
1122 {
throw std::runtime_error(
"solventRefDensity() called but solvents are disabled"); }
1132template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
1145 using Toolbox = MathToolbox<Evaluation>;
1147 static constexpr unsigned gasPhaseIdx = FluidSystem::gasPhaseIdx;
1148 static constexpr int dimWorld = GridView::dimensionworld;
1150 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
1151 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
1158 template <
class Dummy =
bool>
1164 const auto& gradCalc = elemCtx.gradientCalculator();
1167 const auto& scvf = elemCtx.stencil(timeIdx).interiorFace(scvfIdx);
1168 const auto& faceNormal = scvf.normal();
1170 const unsigned i = scvf.interiorIndex();
1171 const unsigned j = scvf.exteriorIndex();
1174 DimEvalVector solventPGrad;
1176 gradCalc.calculateGradient(solventPGrad,
1180 Valgrind::CheckDefined(solventPGrad);
1183 if (Parameters::Get<Parameters::EnableGravity>()) {
1186 const auto& gIn = elemCtx.problem().gravity(elemCtx, i, timeIdx);
1187 const auto& gEx = elemCtx.problem().gravity(elemCtx, j, timeIdx);
1189 const auto& intQuantsIn = elemCtx.intensiveQuantities(i, timeIdx);
1190 const auto& intQuantsEx = elemCtx.intensiveQuantities(j, timeIdx);
1192 const auto& posIn = elemCtx.pos(i, timeIdx);
1193 const auto& posEx = elemCtx.pos(j, timeIdx);
1194 const auto& posFace = scvf.integrationPos();
1197 DimVector distVecIn(posIn);
1198 DimVector distVecEx(posEx);
1199 DimVector distVecTotal(posEx);
1201 distVecIn -= posFace;
1202 distVecEx -= posFace;
1203 distVecTotal -= posIn;
1204 const Scalar absDistTotalSquared = distVecTotal.two_norm2();
1207 auto rhoIn = intQuantsIn.solventDensity();
1208 auto pStatIn = - rhoIn * (gIn * distVecIn);
1212 Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1213 Scalar pStatEx = - rhoEx * (gEx * distVecEx);
1219 DimEvalVector f(distVecTotal);
1220 f *= (pStatEx - pStatIn) / absDistTotalSquared;
1223 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
1224 solventPGrad[dimIdx] += f[dimIdx];
1226 if (!isfinite(solventPGrad[dimIdx])) {
1227 throw NumericalProblem(
"Non-finite potential gradient for solvent 'phase'");
1233 Evaluation solventPGradNormal = 0.0;
1234 for (
unsigned dimIdx = 0; dimIdx < faceNormal.size(); ++dimIdx) {
1235 solventPGradNormal += solventPGrad[dimIdx]*faceNormal[dimIdx];
1238 if (solventPGradNormal > 0) {
1239 solventUpstreamDofIdx_ = j;
1240 solventDownstreamDofIdx_ = i;
1243 solventUpstreamDofIdx_ = i;
1244 solventDownstreamDofIdx_ = j;
1247 const auto& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1253 if (solventUpstreamDofIdx_ == i) {
1254 solventVolumeFlux_ = solventPGradNormal * up.solventMobility();
1257 solventVolumeFlux_ = solventPGradNormal * scalarValue(up.solventMobility());
1265 template <
class Dummy =
bool>
1271 const ExtensiveQuantities& extQuants = asImp_();
1273 const unsigned interiorDofIdx = extQuants.interiorIndex();
1274 const unsigned exteriorDofIdx = extQuants.exteriorIndex();
1275 assert(interiorDofIdx != exteriorDofIdx);
1277 const auto& intQuantsIn = elemCtx.intensiveQuantities(interiorDofIdx, timeIdx);
1278 const auto& intQuantsEx = elemCtx.intensiveQuantities(exteriorDofIdx, timeIdx);
1280 const unsigned I = elemCtx.globalSpaceIndex(interiorDofIdx, timeIdx);
1281 const unsigned J = elemCtx.globalSpaceIndex(exteriorDofIdx, timeIdx);
1283 const Scalar thpres = elemCtx.problem().thresholdPressure(I, J);
1284 const Scalar trans = elemCtx.problem().transmissibility(elemCtx, interiorDofIdx, exteriorDofIdx);
1285 const Scalar g = elemCtx.problem().gravity()[dimWorld - 1];
1287 const Scalar zIn = elemCtx.problem().dofCenterDepth(elemCtx, interiorDofIdx, timeIdx);
1288 const Scalar zEx = elemCtx.problem().dofCenterDepth(elemCtx, exteriorDofIdx, timeIdx);
1289 const Scalar distZ = zIn - zEx;
1291 const Evaluation& rhoIn = intQuantsIn.solventDensity();
1292 const Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1293 const Evaluation& rhoAvg = rhoIn * 0.5 + rhoEx * 0.5;
1295 const Evaluation& pressureInterior = intQuantsIn.fluidState().pressure(gasPhaseIdx);
1296 Evaluation pressureExterior = Toolbox::value(intQuantsEx.fluidState().pressure(gasPhaseIdx));
1297 pressureExterior += distZ * g * rhoAvg;
1299 Evaluation pressureDiffSolvent = pressureExterior - pressureInterior;
1300 if (std::abs(scalarValue(pressureDiffSolvent)) > thpres) {
1301 if (pressureDiffSolvent < 0.0) {
1302 pressureDiffSolvent += thpres;
1305 pressureDiffSolvent -= thpres;
1309 pressureDiffSolvent = 0.0;
1312 if (pressureDiffSolvent > 0.0) {
1313 solventUpstreamDofIdx_ = exteriorDofIdx;
1314 solventDownstreamDofIdx_ = interiorDofIdx;
1316 else if (pressureDiffSolvent < 0.0) {
1317 solventUpstreamDofIdx_ = interiorDofIdx;
1318 solventDownstreamDofIdx_ = exteriorDofIdx;
1324 solventUpstreamDofIdx_ = std::min(interiorDofIdx, exteriorDofIdx);
1325 solventDownstreamDofIdx_ = std::max(interiorDofIdx, exteriorDofIdx);
1326 solventVolumeFlux_ = 0.0;
1330 const Scalar faceArea = elemCtx.stencil(timeIdx).interiorFace(scvfIdx).area();
1331 const IntensiveQuantities& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1332 if (solventUpstreamDofIdx_ == interiorDofIdx) {
1333 solventVolumeFlux_ = up.solventMobility() * (-trans / faceArea) * pressureDiffSolvent;
1336 solventVolumeFlux_ =
1337 scalarValue(up.solventMobility()) * (-trans / faceArea) * pressureDiffSolvent;
1342 {
return solventUpstreamDofIdx_; }
1345 {
return solventDownstreamDofIdx_; }
1348 {
return solventVolumeFlux_; }
1354 Implementation& asImp_()
1355 {
return *
static_cast<Implementation*
>(
this); }
1357 Evaluation solventVolumeFlux_;
1358 unsigned solventUpstreamDofIdx_;
1359 unsigned solventDownstreamDofIdx_;
1362template <
class TypeTag>
1380 {
throw std::runtime_error(
"solventUpstreamIndex() called but solvents are disabled"); }
1383 {
throw std::runtime_error(
"solventDownstreamIndex() called but solvents are disabled"); }
1386 {
throw std::runtime_error(
"solventVolumeFlux() called but solvents are disabled"); }
1389 {
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:1379
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1382
void setSolventVolumeFlux(const Evaluation &)
Definition: blackoilsolventmodules.hh:1388
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1385
void updateVolumeFluxPerm(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1369
void updateVolumeFluxTrans(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1374
Provides the solvent specific extensive quantities to the generic black-oil module's extensive quanti...
Definition: blackoilsolventmodules.hh:1134
unsigned solventUpstreamIndex() const
Definition: blackoilsolventmodules.hh:1341
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:1160
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:1267
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1344
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1347
void setSolventVolumeFlux(const Evaluation &solventVolumeFlux)
Definition: blackoilsolventmodules.hh:1350
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:1109
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:1118
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:1106
void solventPreSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1088
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:1112
Scalar solventRefDensity() const
Definition: blackoilsolventmodules.hh:1121
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:1115
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:1103
void solventPvtUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1098
void solventPostSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1093
void solventPreSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called before the saturation functions are doing their magic.
Definition: blackoilsolventmodules.hh:564
Implementation & asImp_()
Definition: blackoilsolventmodules.hh:1066
Evaluation rsSolw_
Definition: blackoilsolventmodules.hh:1071
Evaluation solventViscosity_
Definition: blackoilsolventmodules.hh:1073
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:820
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:747
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:814
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:826
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:817
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:823
Evaluation hydrocarbonSaturation_
Definition: blackoilsolventmodules.hh:1069
void solventPostSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called after the saturation functions have been doing their magic.
Definition: blackoilsolventmodules.hh:596
Evaluation solventSaturation_
Definition: blackoilsolventmodules.hh:1070
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:811
Evaluation solventDensity_
Definition: blackoilsolventmodules.hh:1072
Evaluation solventInvFormationVolumeFactor_
Definition: blackoilsolventmodules.hh:1075
Scalar solventRefDensity_
Definition: blackoilsolventmodules.hh:1077
Scalar solventRefDensity() const
Definition: blackoilsolventmodules.hh:830
Evaluation solventMobility_
Definition: blackoilsolventmodules.hh:1074
Provides the volumetric quantities required for the equations needed by the solvents extension of the...
Definition: blackoilsolventmodules.hh:526
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilsolventmodules.hh:66
static void setSolventPvt(const SolventPvt &value)
Specify the solvent PVT of a all PVT regions.
Definition: blackoilsolventmodules.hh:105
static const TabulatedFunction & sof2Krn(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:396
static const TabulatedFunction & ssfnKrs(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:387
static const TabulatedFunction & pmisc(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:414
static const TabulatedFunction & msfnKrsg(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:423
static const TabulatedFunction & sgcwmis(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:450
static Scalar computeResidualError(const EqVector &resid)
Return how much a residual is considered an error.
Definition: blackoilsolventmodules.hh:330
static const BrineH2Pvt & brineH2Pvt()
Definition: blackoilsolventmodules.hh:375
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:149
static void setParams(BlackOilSolventParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilsolventmodules.hh:99
static std::string eqName(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:167
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:174
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilsolventmodules.hh:337
static bool isSolubleInWater()
Definition: blackoilsolventmodules.hh:508
static bool isCO2Sol()
Definition: blackoilsolventmodules.hh:511
static const TabulatedFunction & msfnKro(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:432
static void registerOutputModules(Model &model, Simulator &simulator)
Register all solvent specific VTK and ECL output modules.
Definition: blackoilsolventmodules.hh:124
static Scalar tlMixParamDensity(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:477
static Value solubilityLimit(unsigned pvtIdx, const Value &temperature, const Value &pressure, const Value &saltConcentration)
Definition: blackoilsolventmodules.hh:490
static const TabulatedFunction & ssfnKrg(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:378
static void setIsMiscible(const bool isMiscible)
Definition: blackoilsolventmodules.hh:108
static void registerParameters()
Register all run-time parameters for the black-oil solvent module.
Definition: blackoilsolventmodules.hh:114
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilsolventmodules.hh:318
static const SolventPvt & solventPvt()
Definition: blackoilsolventmodules.hh:363
static bool isH2Sol()
Definition: blackoilsolventmodules.hh:514
static bool isMiscible()
Definition: blackoilsolventmodules.hh:486
static const H2GasPvt & h2GasPvt()
Definition: blackoilsolventmodules.hh:369
static Scalar tlMixParamViscosity(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:468
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:305
static const TabulatedFunction & tlPMixTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:459
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar solventSaturation, Scalar solventRsw)
Assign the solvent specific primary variables to a PrimaryVariables object.
Definition: blackoilsolventmodules.hh:284
static const TabulatedFunction & sorwmis(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:441
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:142
static bool eqApplies(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:157
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilsolventmodules.hh:183
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilsolventmodules.hh:348
static const Co2GasPvt & co2GasPvt()
Definition: blackoilsolventmodules.hh:366
static const BrineCo2Pvt & brineCo2Pvt()
Definition: blackoilsolventmodules.hh:372
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:216
static const TabulatedFunction & misc(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:405
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:132
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