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>()>
80 using Toolbox = MathToolbox<Evaluation>;
89 static constexpr unsigned solventSaturationIdx = Indices::solventSaturationIdx;
90 static constexpr unsigned contiSolventEqIdx = Indices::contiSolventEqIdx;
91 static constexpr unsigned enableSolvent = enableSolventV;
92 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
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 template <
class ElemContext>
379 static const TabulatedFunction&
ssfnKrg(
const ElemContext& elemCtx,
383 const unsigned satnumRegionIdx =
384 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
385 return params_.ssfnKrg_[satnumRegionIdx];
388 template <
class ElemContext>
389 static const TabulatedFunction&
ssfnKrs(
const ElemContext& elemCtx,
393 const unsigned satnumRegionIdx =
394 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
395 return params_.ssfnKrs_[satnumRegionIdx];
398 template <
class ElemContext>
399 static const TabulatedFunction&
sof2Krn(
const ElemContext& elemCtx,
403 const unsigned satnumRegionIdx =
404 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
405 return params_.sof2Krn_[satnumRegionIdx];
408 template <
class ElemContext>
409 static const TabulatedFunction&
misc(
const ElemContext& elemCtx,
413 const unsigned miscnumRegionIdx =
414 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
415 return params_.misc_[miscnumRegionIdx];
418 template <
class ElemContext>
419 static const TabulatedFunction&
pmisc(
const ElemContext& elemCtx,
423 const unsigned miscnumRegionIdx =
424 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
425 return params_.pmisc_[miscnumRegionIdx];
428 template <
class ElemContext>
429 static const TabulatedFunction&
msfnKrsg(
const ElemContext& elemCtx,
433 const unsigned satnumRegionIdx =
434 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
435 return params_.msfnKrsg_[satnumRegionIdx];
438 template <
class ElemContext>
439 static const TabulatedFunction&
msfnKro(
const ElemContext& elemCtx,
443 const unsigned satnumRegionIdx =
444 elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
445 return params_.msfnKro_[satnumRegionIdx];
448 template <
class ElemContext>
449 static const TabulatedFunction&
sorwmis(
const ElemContext& elemCtx,
453 const unsigned miscnumRegionIdx =
454 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
455 return params_.sorwmis_[miscnumRegionIdx];
458 template <
class ElemContext>
459 static const TabulatedFunction&
sgcwmis(
const ElemContext& elemCtx,
463 const unsigned miscnumRegionIdx =
464 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
465 return params_.sgcwmis_[miscnumRegionIdx];
468 template <
class ElemContext>
469 static const TabulatedFunction&
tlPMixTable(
const ElemContext& elemCtx,
473 const unsigned miscnumRegionIdx =
474 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
475 return params_.tlPMixTable_[miscnumRegionIdx];
478 template <
class ElemContext>
483 const unsigned miscnumRegionIdx =
484 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
485 return params_.tlMixParamViscosity_[miscnumRegionIdx];
489 template <
class ElemContext>
494 const unsigned miscnumRegionIdx =
495 elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
496 return params_.tlMixParamDensity_[miscnumRegionIdx];
500 {
return params_.isMiscible_; }
502 template <
class Value>
504 const Value& pressure,
const Value& saltConcentration)
512 return brineCo2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature,
513 pressure, saltConcentration);
516 return brineH2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature,
517 pressure, saltConcentration);
522 {
return params_.rsSolw_active_; }
525 {
return params_.co2sol_; }
528 {
return params_.h2sol_; }
534template <
class TypeTag,
bool enableSolventV>
535BlackOilSolventParams<typename BlackOilSolventModule<TypeTag, enableSolventV>::Scalar>
536BlackOilSolventModule<TypeTag, enableSolventV>::params_;
538template <
class TypeTag,
bool enableSolventV>
548template <
class TypeTag>
564 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
565 static constexpr int solventSaturationIdx = Indices::solventSaturationIdx;
566 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
567 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
568 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
569 static constexpr double cutOff = 1e-12;
582 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
583 this->solventPreSatFuncUpdate_(priVars, timeIdx, elemCtx.linearizationType());
587 const unsigned timeIdx,
590 auto& fs = asImp_().fluidState_;
591 solventSaturation_ = 0.0;
592 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
593 solventSaturation_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, linearizationType);
596 hydrocarbonSaturation_ = fs.saturation(gasPhaseIdx);
599 if (solventSaturation().value() < cutOff) {
605 fs.setSaturation(gasPhaseIdx, hydrocarbonSaturation_ + solventSaturation_);
619 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
620 this->solventPostSatFuncUpdate_(elemCtx.problem(), priVars, elemCtx.globalSpaceIndex(dofIdx, timeIdx), timeIdx, elemCtx.linearizationType());
625 template <
class Problem>
626 struct ProblemAndCellIndexOnlyContext
628 const Problem& problem_;
630 const Problem& problem()
const {
return problem_; }
631 unsigned int globalSpaceIndex([[maybe_unused]]
const unsigned int spaceIdx,
632 [[maybe_unused]]
const unsigned int timeIdx)
const
640 const PrimaryVariables& priVars,
641 const unsigned globalSpaceIdx,
642 const unsigned timeIdx,
647 auto& fs = asImp_().fluidState_;
648 fs.setSaturation(gasPhaseIdx, hydrocarbonSaturation_);
652 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
653 rsSolw_ = SolventModule::solubilityLimit(asImp_().pvtRegionIndex(),
654 fs.temperature(waterPhaseIdx),
655 fs.pressure(waterPhaseIdx),
656 fs.saltConcentration());
657 }
else if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Rsolw) {
658 rsSolw_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, linearizationType);
661 solventMobility_ = 0.0;
664 if (solventSaturation().value() < cutOff) {
668 ProblemAndCellIndexOnlyContext<Problem> elemCtx{problem, globalSpaceIdx};
672 if (SolventModule::isMiscible()) {
673 const Evaluation& p =
674 FluidSystem::phaseIsActive(oilPhaseIdx)
675 ? fs.pressure(oilPhaseIdx)
676 : fs.pressure(gasPhaseIdx);
677 const Evaluation pmisc =
678 SolventModule::pmisc(elemCtx, dofIdx, timeIdx).eval(p,
true);
679 const Evaluation& pgImisc = fs.pressure(gasPhaseIdx);
682 Evaluation pgMisc = 0.0;
683 std::array<Evaluation, numPhases> pC;
684 const auto& materialParams = problem.materialLawParams(elemCtx, dofIdx, timeIdx);
685 MaterialLaw::capillaryPressures(pC, materialParams, fs);
688 if (priVars.primaryVarsMeaningPressure() == PrimaryVariables::PressureMeaning::Pg) {
689 pgMisc = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx,
693 const Evaluation& po =
694 priVars.makeEvaluation(Indices::pressureSwitchIdx,
695 timeIdx, linearizationType);
696 pgMisc = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
699 fs.setPressure(gasPhaseIdx, pmisc * pgMisc + (1.0 - pmisc) * pgImisc);
702 const Evaluation gasSolventSat = hydrocarbonSaturation_ + solventSaturation_;
704 if (gasSolventSat.value() < cutOff) {
708 const Evaluation Fhydgas = hydrocarbonSaturation_ / gasSolventSat;
709 const Evaluation Fsolgas = solventSaturation_ / gasSolventSat;
712 if (SolventModule::isMiscible() && FluidSystem::phaseIsActive(oilPhaseIdx)) {
713 const auto& misc = SolventModule::misc(elemCtx, dofIdx, timeIdx);
714 const auto& pmisc = SolventModule::pmisc(elemCtx, dofIdx, timeIdx);
715 const Evaluation& p =
716 FluidSystem::phaseIsActive(oilPhaseIdx)
717 ? fs.pressure(oilPhaseIdx)
718 : fs.pressure(gasPhaseIdx);
719 const Evaluation miscibility =
720 misc.eval(Fsolgas,
true) *
724 const unsigned cellIdx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
725 const auto& materialLawManager = elemCtx.problem().materialLawManager();
726 const auto& scaledDrainageInfo =
727 materialLawManager->oilWaterScaledEpsInfoDrainage(cellIdx);
729 const Scalar sogcr = scaledDrainageInfo.Sogcr;
730 Evaluation sor = sogcr;
731 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
732 const Evaluation& sw = fs.saturation(waterPhaseIdx);
733 const auto& sorwmis = SolventModule::sorwmis(elemCtx, dofIdx, timeIdx);
735 sorwmis.eval(sw,
true) + (1.0 - miscibility) * sogcr;
737 const Scalar sgcr = scaledDrainageInfo.Sgcr;
738 Evaluation sgc = sgcr;
739 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
740 const Evaluation& sw = fs.saturation(waterPhaseIdx);
741 const auto& sgcwmis = SolventModule::sgcwmis(elemCtx, dofIdx, timeIdx);
742 sgc = miscibility * sgcwmis.eval(sw,
true) + (1.0 - miscibility) * sgcr;
745 Evaluation oilGasSolventSat = gasSolventSat;
746 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
747 oilGasSolventSat += fs.saturation(oilPhaseIdx);
749 const Evaluation zero = 0.0;
750 const Evaluation oilGasSolventEffSat = std::max(oilGasSolventSat - sor - sgc, zero);
752 Evaluation F_totalGas = 0.0;
753 if (oilGasSolventEffSat.value() > cutOff) {
754 const Evaluation gasSolventEffSat = std::max(gasSolventSat - sgc, zero);
755 F_totalGas = gasSolventEffSat / oilGasSolventEffSat;
757 const auto& msfnKro = SolventModule::msfnKro(elemCtx, dofIdx, timeIdx);
758 const auto& msfnKrsg = SolventModule::msfnKrsg(elemCtx, dofIdx, timeIdx);
759 const auto& sof2Krn = SolventModule::sof2Krn(elemCtx, dofIdx, timeIdx);
761 const Evaluation mkrgt = msfnKrsg.eval(F_totalGas,
true) *
762 sof2Krn.eval(oilGasSolventSat,
true);
763 const Evaluation mkro = msfnKro.eval(F_totalGas,
true) *
764 sof2Krn.eval(oilGasSolventSat,
true);
766 Evaluation& kro = asImp_().mobility_[oilPhaseIdx];
767 Evaluation& krg = asImp_().mobility_[gasPhaseIdx];
770 krg *= 1.0 - miscibility;
771 krg += miscibility * mkrgt;
772 kro *= 1.0 - miscibility;
773 kro += miscibility * mkro;
777 const auto& ssfnKrg = SolventModule::ssfnKrg(elemCtx, dofIdx, timeIdx);
778 const auto& ssfnKrs = SolventModule::ssfnKrs(elemCtx, dofIdx, timeIdx);
780 Evaluation& krg = asImp_().mobility_[gasPhaseIdx];
781 solventMobility_ = krg * ssfnKrs.eval(Fsolgas,
true);
782 krg *= ssfnKrg.eval(Fhydgas,
true);
795 const auto& iq = asImp_();
796 const unsigned pvtRegionIdx = iq.pvtRegionIndex();
797 const Evaluation& T = iq.fluidState().temperature(gasPhaseIdx);
798 const Evaluation& p = iq.fluidState().pressure(gasPhaseIdx);
800 const Evaluation rv = 0.0;
801 const Evaluation rvw = 0.0;
802 if (SolventModule::isCO2Sol() || SolventModule::isH2Sol() ){
803 if (SolventModule::isCO2Sol()) {
804 const auto& co2gasPvt = SolventModule::co2GasPvt();
805 solventInvFormationVolumeFactor_ =
806 co2gasPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rv, rvw);
807 solventRefDensity_ = co2gasPvt.gasReferenceDensity(pvtRegionIdx);
808 solventViscosity_ = co2gasPvt.viscosity(pvtRegionIdx, T, p, rv, rvw);
810 const auto& brineCo2Pvt = SolventModule::brineCo2Pvt();
811 auto& fs = asImp_().fluidState_;
812 const auto& bw = brineCo2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rsSolw());
814 const auto denw = bw * brineCo2Pvt.waterReferenceDensity(pvtRegionIdx) +
815 rsSolw() * bw * brineCo2Pvt.gasReferenceDensity(pvtRegionIdx);
816 fs.setDensity(waterPhaseIdx, denw);
817 fs.setInvB(waterPhaseIdx, bw);
818 Evaluation& mobw = asImp_().mobility_[waterPhaseIdx];
819 const auto& muWat = fs.viscosity(waterPhaseIdx);
820 const auto& muWatEff = brineCo2Pvt.viscosity(pvtRegionIdx, T, p, rsSolw());
821 mobw *= muWat / muWatEff;
823 const auto& h2gasPvt = SolventModule::h2GasPvt();
824 solventInvFormationVolumeFactor_ =
825 h2gasPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rv, rvw);
826 solventRefDensity_ = h2gasPvt.gasReferenceDensity(pvtRegionIdx);
827 solventViscosity_ = h2gasPvt.viscosity(pvtRegionIdx, T, p, rv, rvw);
829 const auto& brineH2Pvt = SolventModule::brineH2Pvt();
830 auto& fs = asImp_().fluidState_;
831 const auto& bw = brineH2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p, rsSolw());
833 const auto denw = bw * brineH2Pvt.waterReferenceDensity(pvtRegionIdx) +
834 rsSolw() * bw * brineH2Pvt.gasReferenceDensity(pvtRegionIdx);
835 fs.setDensity(waterPhaseIdx, denw);
836 fs.setInvB(waterPhaseIdx, bw);
837 Evaluation& mobw = asImp_().mobility_[waterPhaseIdx];
838 const auto& muWat = fs.viscosity(waterPhaseIdx);
839 const auto& muWatEff = brineH2Pvt.viscosity(pvtRegionIdx, T, p, rsSolw());
840 mobw *= muWat / muWatEff;
843 const auto& solventPvt = SolventModule::solventPvt();
844 solventInvFormationVolumeFactor_ = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
845 solventRefDensity_ = solventPvt.referenceDensity(pvtRegionIdx);
846 solventViscosity_ = solventPvt.viscosity(pvtRegionIdx, T, p);
849 solventDensity_ = solventInvFormationVolumeFactor_*solventRefDensity_;
850 effectiveProperties(elemCtx, scvIdx, timeIdx);
852 solventMobility_ /= solventViscosity_;
856 {
return solventSaturation_; }
862 {
return solventDensity_; }
865 {
return solventViscosity_; }
868 {
return solventMobility_; }
871 {
return solventInvFormationVolumeFactor_; }
875 {
return solventRefDensity_; }
880 void effectiveProperties(
const ElementContext& elemCtx,
884 if (!SolventModule::isMiscible()) {
890 if (solventSaturation() < cutOff) {
896 auto& fs = asImp_().fluidState_;
899 Evaluation oilEffSat = 0.0;
900 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
901 oilEffSat = fs.saturation(oilPhaseIdx);
903 Evaluation gasEffSat = fs.saturation(gasPhaseIdx);
904 Evaluation solventEffSat = solventSaturation();
905 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
906 const auto& sorwmis = SolventModule::sorwmis(elemCtx, scvIdx, timeIdx);
907 const auto& sgcwmis = SolventModule::sgcwmis(elemCtx, scvIdx, timeIdx);
908 const Evaluation zero = 0.0;
909 const Evaluation& sw = fs.saturation(waterPhaseIdx);
910 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
911 oilEffSat = std::max(oilEffSat - sorwmis.eval(sw,
true), zero);
913 gasEffSat = std::max(gasEffSat - sgcwmis.eval(sw,
true), zero);
914 solventEffSat = std::max(solventEffSat - sgcwmis.eval(sw,
true), zero);
916 const Evaluation oilGasSolventEffSat = oilEffSat + gasEffSat + solventEffSat;
917 const Evaluation oilSolventEffSat = oilEffSat + solventEffSat;
918 const Evaluation solventGasEffSat = solventEffSat + gasEffSat;
925 const Evaluation& p =
926 FluidSystem::phaseIsActive(oilPhaseIdx)
927 ? fs.pressure(oilPhaseIdx)
928 : fs.pressure(gasPhaseIdx);
930 const auto& pmiscTable = SolventModule::pmisc(elemCtx, scvIdx, timeIdx);
931 const Evaluation pmisc = pmiscTable.eval(p,
true);
932 const auto& tlPMixTable = SolventModule::tlPMixTable(elemCtx, scvIdx, timeIdx);
933 const Evaluation tlMixParamMu = SolventModule::tlMixParamViscosity(elemCtx, scvIdx, timeIdx) *
934 tlPMixTable.eval(p,
true);
936 const Evaluation& muGas = fs.viscosity(gasPhaseIdx);
937 const Evaluation& muSolvent = solventViscosity_;
939 assert(muGas.value() > 0);
940 assert(muSolvent.value() > 0);
941 const Evaluation muGasPow = pow(muGas, 0.25);
942 const Evaluation muSolventPow = pow(muSolvent, 0.25);
944 Evaluation muMixSolventGas = muGas;
945 if (solventGasEffSat > cutOff) {
946 muMixSolventGas *= muSolvent / pow(((gasEffSat / solventGasEffSat) * muSolventPow) +
947 ((solventEffSat / solventGasEffSat) * muGasPow), 4.0);
950 Evaluation muOil = 1.0;
951 Evaluation muOilPow = 1.0;
952 Evaluation muMixOilSolvent = 1.0;
953 Evaluation muOilEff = 1.0;
954 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
955 muOil = fs.viscosity(oilPhaseIdx);
956 assert(muOil.value() > 0);
957 muOilPow = pow(muOil, 0.25);
958 muMixOilSolvent = muOil;
959 if (oilSolventEffSat > cutOff) {
960 muMixOilSolvent *= muSolvent / pow(((oilEffSat / oilSolventEffSat) * muSolventPow) +
961 ((solventEffSat / oilSolventEffSat) * muOilPow), 4.0);
964 muOilEff = pow(muOil,1.0 - tlMixParamMu) * pow(muMixOilSolvent, tlMixParamMu);
966 Evaluation muMixSolventGasOil = muOil;
967 if (oilGasSolventEffSat > cutOff) {
968 muMixSolventGasOil *= muSolvent * muGas / pow(((oilEffSat / oilGasSolventEffSat) *
969 muSolventPow * muGasPow) +
970 ((solventEffSat / oilGasSolventEffSat) *
971 muOilPow * muGasPow) +
972 ((gasEffSat / oilGasSolventEffSat) *
973 muSolventPow * muOilPow), 4.0);
976 Evaluation muGasEff = pow(muGas,1.0 - tlMixParamMu) * pow(muMixSolventGas, tlMixParamMu);
977 const Evaluation muSolventEff = pow(muSolvent,1.0 - tlMixParamMu) * pow(muMixSolventGasOil, tlMixParamMu);
980 const Evaluation& rhoGas = fs.density(gasPhaseIdx);
981 Evaluation rhoOil = 0.0;
982 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
983 rhoOil = fs.density(oilPhaseIdx);
986 const Evaluation& rhoSolvent = solventDensity_;
992 const Evaluation tlMixParamRho = SolventModule::tlMixParamDensity(elemCtx, scvIdx, timeIdx) *
993 tlPMixTable.eval(p,
true);
997 const Evaluation muOilEffPow = pow(pow(muOil, 1.0 - tlMixParamRho) *
998 pow(muMixOilSolvent, tlMixParamRho), 0.25);
999 const Evaluation muGasEffPow = pow(pow(muGas, 1.0 - tlMixParamRho) *
1000 pow(muMixSolventGas, tlMixParamRho), 0.25);
1001 const Evaluation muSolventEffPow = pow(pow(muSolvent, 1.0 - tlMixParamRho) *
1002 pow(muMixSolventGasOil, tlMixParamRho), 0.25);
1004 const Evaluation oilGasEffSaturation = oilEffSat + gasEffSat;
1005 Evaluation sof = 0.0;
1006 Evaluation sgf = 0.0;
1007 if (oilGasEffSaturation.value() > cutOff) {
1008 sof = oilEffSat / oilGasEffSaturation;
1009 sgf = gasEffSat / oilGasEffSaturation;
1012 const Evaluation muSolventOilGasPow = muSolventPow * ((sgf * muOilPow) + (sof * muGasPow));
1014 Evaluation rhoMixSolventGasOil = 0.0;
1015 if (oilGasSolventEffSat.value() > cutOff) {
1016 rhoMixSolventGasOil = (rhoOil * oilEffSat / oilGasSolventEffSat) +
1017 (rhoGas * gasEffSat / oilGasSolventEffSat) +
1018 (rhoSolvent * solventEffSat / oilGasSolventEffSat);
1021 Evaluation rhoGasEff = 0.0;
1022 if (std::abs(muSolventPow.value() - muGasPow.value()) < cutOff) {
1023 rhoGasEff = ((1.0 - tlMixParamRho) * rhoGas) +
1024 (tlMixParamRho * rhoMixSolventGasOil);
1027 const Evaluation solventGasEffFraction = (muGasPow * (muSolventPow - muGasEffPow)) /
1028 (muGasEffPow * (muSolventPow - muGasPow));
1029 rhoGasEff = (rhoGas * solventGasEffFraction) + (rhoSolvent * (1.0 - solventGasEffFraction));
1032 Evaluation rhoSolventEff = 0.0;
1033 if (std::abs((muSolventOilGasPow.value() - (muOilPow.value() * muGasPow.value()))) < cutOff) {
1034 rhoSolventEff = ((1.0 - tlMixParamRho) * rhoSolvent) + (tlMixParamRho * rhoMixSolventGasOil);
1037 const Evaluation sfraction_se = (muSolventOilGasPow -
1038 muOilPow * muGasPow * muSolventPow / muSolventEffPow) /
1039 (muSolventOilGasPow - (muOilPow * muGasPow));
1040 rhoSolventEff = (rhoSolvent * sfraction_se) +
1041 (rhoGas * sgf * (1.0 - sfraction_se)) +
1042 (rhoOil * sof * (1.0 - sfraction_se));
1046 const unsigned pvtRegionIdx = asImp_().pvtRegionIndex();
1047 Evaluation bGasEff = rhoGasEff;
1048 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1049 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) +
1050 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) * fs.Rv());
1052 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1054 const Evaluation bSolventEff = rhoSolventEff / solventRefDensity();
1057 const Evaluation bg = fs.invB(gasPhaseIdx);
1058 const Evaluation bs = solventInverseFormationVolumeFactor();
1059 const Evaluation bg_eff = pmisc * bGasEff + (1.0 - pmisc) * bg;
1060 const Evaluation bs_eff = pmisc * bSolventEff + (1.0 - pmisc) * bs;
1063 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1064 fs.setDensity(gasPhaseIdx,
1066 (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) +
1067 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx)*fs.Rv()));
1069 fs.setDensity(gasPhaseIdx,
1070 bg_eff * FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1072 solventDensity_ = bs_eff * solventRefDensity();
1075 Evaluation& mobg = asImp_().mobility_[gasPhaseIdx];
1076 muGasEff = bg_eff / (pmisc * bGasEff / muGasEff + (1.0 - pmisc) * bg / muGas);
1077 mobg *= muGas / muGasEff;
1080 solventViscosity_ = bs_eff / (pmisc * bSolventEff / muSolventEff +
1081 (1.0 - pmisc) * bs / muSolvent);
1083 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1084 Evaluation rhoOilEff = 0.0;
1085 if (std::abs(muOilPow.value() - muSolventPow.value()) < cutOff) {
1086 rhoOilEff = ((1.0 - tlMixParamRho) * rhoOil) + (tlMixParamRho * rhoMixSolventGasOil);
1089 const Evaluation solventOilEffFraction = (muOilPow * (muOilEffPow - muSolventPow)) /
1090 (muOilEffPow * (muOilPow - muSolventPow));
1091 rhoOilEff = (rhoOil * solventOilEffFraction) + (rhoSolvent * (1.0 - solventOilEffFraction));
1093 const Evaluation bOilEff =
1094 rhoOilEff / (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) +
1095 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs());
1096 const Evaluation bo = fs.invB(oilPhaseIdx);
1097 const Evaluation bo_eff = pmisc * bOilEff + (1.0 - pmisc) * bo;
1098 fs.setDensity(oilPhaseIdx,
1099 bo_eff * (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) +
1100 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs()));
1103 Evaluation& mobo = asImp_().mobility_[oilPhaseIdx];
1104 muOilEff = bo_eff / (pmisc * bOilEff / muOilEff + (1.0 - pmisc) * bo / muOil);
1105 mobo *= muOil / muOilEff;
1111 {
return *
static_cast<Implementation*
>(
this); }
1124template <
class TypeTag>
1148 {
throw std::runtime_error(
"solventSaturation() called but solvents are disabled"); }
1151 {
throw std::runtime_error(
"rsSolw() called but solvents are disabled"); }
1154 {
throw std::runtime_error(
"solventDensity() called but solvents are disabled"); }
1157 {
throw std::runtime_error(
"solventViscosity() called but solvents are disabled"); }
1160 {
throw std::runtime_error(
"solventMobility() called but solvents are disabled"); }
1163 {
throw std::runtime_error(
"solventInverseFormationVolumeFactor() called but solvents are disabled"); }
1166 {
throw std::runtime_error(
"solventRefDensity() called but solvents are disabled"); }
1176template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
1189 using Toolbox = MathToolbox<Evaluation>;
1191 static constexpr unsigned gasPhaseIdx = FluidSystem::gasPhaseIdx;
1192 static constexpr int dimWorld = GridView::dimensionworld;
1194 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
1195 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
1202 template <
class Dummy =
bool>
1208 const auto& gradCalc = elemCtx.gradientCalculator();
1211 const auto& scvf = elemCtx.stencil(timeIdx).interiorFace(scvfIdx);
1212 const auto& faceNormal = scvf.normal();
1214 const unsigned i = scvf.interiorIndex();
1215 const unsigned j = scvf.exteriorIndex();
1218 DimEvalVector solventPGrad;
1220 gradCalc.calculateGradient(solventPGrad,
1224 Valgrind::CheckDefined(solventPGrad);
1227 if (Parameters::Get<Parameters::EnableGravity>()) {
1230 const auto& gIn = elemCtx.problem().gravity(elemCtx, i, timeIdx);
1231 const auto& gEx = elemCtx.problem().gravity(elemCtx, j, timeIdx);
1233 const auto& intQuantsIn = elemCtx.intensiveQuantities(i, timeIdx);
1234 const auto& intQuantsEx = elemCtx.intensiveQuantities(j, timeIdx);
1236 const auto& posIn = elemCtx.pos(i, timeIdx);
1237 const auto& posEx = elemCtx.pos(j, timeIdx);
1238 const auto& posFace = scvf.integrationPos();
1241 DimVector distVecIn(posIn);
1242 DimVector distVecEx(posEx);
1243 DimVector distVecTotal(posEx);
1245 distVecIn -= posFace;
1246 distVecEx -= posFace;
1247 distVecTotal -= posIn;
1248 const Scalar absDistTotalSquared = distVecTotal.two_norm2();
1251 auto rhoIn = intQuantsIn.solventDensity();
1252 auto pStatIn = - rhoIn * (gIn * distVecIn);
1256 Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1257 Scalar pStatEx = - rhoEx * (gEx * distVecEx);
1263 DimEvalVector f(distVecTotal);
1264 f *= (pStatEx - pStatIn) / absDistTotalSquared;
1267 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
1268 solventPGrad[dimIdx] += f[dimIdx];
1270 if (!isfinite(solventPGrad[dimIdx])) {
1271 throw NumericalProblem(
"Non-finite potential gradient for solvent 'phase'");
1277 Evaluation solventPGradNormal = 0.0;
1278 for (
unsigned dimIdx = 0; dimIdx < faceNormal.size(); ++dimIdx) {
1279 solventPGradNormal += solventPGrad[dimIdx]*faceNormal[dimIdx];
1282 if (solventPGradNormal > 0) {
1283 solventUpstreamDofIdx_ = j;
1284 solventDownstreamDofIdx_ = i;
1287 solventUpstreamDofIdx_ = i;
1288 solventDownstreamDofIdx_ = j;
1291 const auto& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1297 if (solventUpstreamDofIdx_ == i) {
1298 solventVolumeFlux_ = solventPGradNormal * up.solventMobility();
1301 solventVolumeFlux_ = solventPGradNormal * scalarValue(up.solventMobility());
1309 template <
class Dummy =
bool>
1315 const ExtensiveQuantities& extQuants = asImp_();
1317 const unsigned interiorDofIdx = extQuants.interiorIndex();
1318 const unsigned exteriorDofIdx = extQuants.exteriorIndex();
1319 assert(interiorDofIdx != exteriorDofIdx);
1321 const auto& intQuantsIn = elemCtx.intensiveQuantities(interiorDofIdx, timeIdx);
1322 const auto& intQuantsEx = elemCtx.intensiveQuantities(exteriorDofIdx, timeIdx);
1324 const unsigned I = elemCtx.globalSpaceIndex(interiorDofIdx, timeIdx);
1325 const unsigned J = elemCtx.globalSpaceIndex(exteriorDofIdx, timeIdx);
1327 const Scalar thpres = elemCtx.problem().thresholdPressure(I, J);
1328 const Scalar trans = elemCtx.problem().transmissibility(elemCtx, interiorDofIdx, exteriorDofIdx);
1329 const Scalar g = elemCtx.problem().gravity()[dimWorld - 1];
1331 const Scalar zIn = elemCtx.problem().dofCenterDepth(elemCtx, interiorDofIdx, timeIdx);
1332 const Scalar zEx = elemCtx.problem().dofCenterDepth(elemCtx, exteriorDofIdx, timeIdx);
1333 const Scalar distZ = zIn - zEx;
1335 const Evaluation& rhoIn = intQuantsIn.solventDensity();
1336 const Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1337 const Evaluation& rhoAvg = rhoIn * 0.5 + rhoEx * 0.5;
1339 const Evaluation& pressureInterior = intQuantsIn.fluidState().pressure(gasPhaseIdx);
1340 Evaluation pressureExterior = Toolbox::value(intQuantsEx.fluidState().pressure(gasPhaseIdx));
1341 pressureExterior += distZ * g * rhoAvg;
1343 Evaluation pressureDiffSolvent = pressureExterior - pressureInterior;
1344 if (std::abs(scalarValue(pressureDiffSolvent)) > thpres) {
1345 if (pressureDiffSolvent < 0.0) {
1346 pressureDiffSolvent += thpres;
1349 pressureDiffSolvent -= thpres;
1353 pressureDiffSolvent = 0.0;
1356 if (pressureDiffSolvent > 0.0) {
1357 solventUpstreamDofIdx_ = exteriorDofIdx;
1358 solventDownstreamDofIdx_ = interiorDofIdx;
1360 else if (pressureDiffSolvent < 0.0) {
1361 solventUpstreamDofIdx_ = interiorDofIdx;
1362 solventDownstreamDofIdx_ = exteriorDofIdx;
1368 solventUpstreamDofIdx_ = std::min(interiorDofIdx, exteriorDofIdx);
1369 solventDownstreamDofIdx_ = std::max(interiorDofIdx, exteriorDofIdx);
1370 solventVolumeFlux_ = 0.0;
1374 const Scalar faceArea = elemCtx.stencil(timeIdx).interiorFace(scvfIdx).area();
1375 const IntensiveQuantities& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1376 if (solventUpstreamDofIdx_ == interiorDofIdx) {
1377 solventVolumeFlux_ = up.solventMobility() * (-trans / faceArea) * pressureDiffSolvent;
1380 solventVolumeFlux_ =
1381 scalarValue(up.solventMobility()) * (-trans / faceArea) * pressureDiffSolvent;
1386 {
return solventUpstreamDofIdx_; }
1389 {
return solventDownstreamDofIdx_; }
1392 {
return solventVolumeFlux_; }
1398 Implementation& asImp_()
1399 {
return *
static_cast<Implementation*
>(
this); }
1401 Evaluation solventVolumeFlux_;
1402 unsigned solventUpstreamDofIdx_;
1403 unsigned solventDownstreamDofIdx_;
1406template <
class TypeTag>
1424 {
throw std::runtime_error(
"solventUpstreamIndex() called but solvents are disabled"); }
1427 {
throw std::runtime_error(
"solventDownstreamIndex() called but solvents are disabled"); }
1430 {
throw std::runtime_error(
"solventVolumeFlux() called but solvents are disabled"); }
1433 {
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:1423
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1426
void setSolventVolumeFlux(const Evaluation &)
Definition: blackoilsolventmodules.hh:1432
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1429
void updateVolumeFluxPerm(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1413
void updateVolumeFluxTrans(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1418
Provides the solvent specific extensive quantities to the generic black-oil module's extensive quanti...
Definition: blackoilsolventmodules.hh:1178
unsigned solventUpstreamIndex() const
Definition: blackoilsolventmodules.hh:1385
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:1204
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:1311
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1388
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1391
void setSolventVolumeFlux(const Evaluation &solventVolumeFlux)
Definition: blackoilsolventmodules.hh:1394
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:1153
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:1162
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:1150
void solventPreSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1132
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:1156
Scalar solventRefDensity() const
Definition: blackoilsolventmodules.hh:1165
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:1159
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:1147
void solventPvtUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1142
void solventPostSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1137
void solventPreSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called before the saturation functions are doing their magic.
Definition: blackoilsolventmodules.hh:578
Implementation & asImp_()
Definition: blackoilsolventmodules.hh:1110
Evaluation rsSolw_
Definition: blackoilsolventmodules.hh:1115
Evaluation solventViscosity_
Definition: blackoilsolventmodules.hh:1117
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:864
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:791
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:858
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:870
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:861
void solventPostSatFuncUpdate_(const Problem &problem, const PrimaryVariables &priVars, const unsigned globalSpaceIdx, const unsigned timeIdx, const LinearizationType linearizationType)
Definition: blackoilsolventmodules.hh:639
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:867
Evaluation hydrocarbonSaturation_
Definition: blackoilsolventmodules.hh:1113
void solventPreSatFuncUpdate_(const PrimaryVariables &priVars, const unsigned timeIdx, const LinearizationType linearizationType)
Definition: blackoilsolventmodules.hh:586
void solventPostSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called after the saturation functions have been doing their magic.
Definition: blackoilsolventmodules.hh:615
Evaluation solventSaturation_
Definition: blackoilsolventmodules.hh:1114
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:855
Evaluation solventDensity_
Definition: blackoilsolventmodules.hh:1116
Evaluation solventInvFormationVolumeFactor_
Definition: blackoilsolventmodules.hh:1119
Scalar solventRefDensity_
Definition: blackoilsolventmodules.hh:1121
Scalar solventRefDensity() const
Definition: blackoilsolventmodules.hh:874
Evaluation solventMobility_
Definition: blackoilsolventmodules.hh:1118
Provides the volumetric quantities required for the equations needed by the solvents extension of the...
Definition: blackoilsolventmodules.hh:539
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:105
static const TabulatedFunction & misc(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:409
static const TabulatedFunction & sof2Krn(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:399
static Scalar computeResidualError(const EqVector &resid)
Return how much a residual is considered an error.
Definition: blackoilsolventmodules.hh:330
static Scalar tlMixParamDensity(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:490
static const TabulatedFunction & ssfnKrg(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:379
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:521
static bool isCO2Sol()
Definition: blackoilsolventmodules.hh:524
static const TabulatedFunction & sgcwmis(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:459
static Scalar tlMixParamViscosity(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:479
static void registerOutputModules(Model &model, Simulator &simulator)
Register all solvent specific VTK and ECL output modules.
Definition: blackoilsolventmodules.hh:124
static const TabulatedFunction & tlPMixTable(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:469
static const TabulatedFunction & sorwmis(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:449
static Value solubilityLimit(unsigned pvtIdx, const Value &temperature, const Value &pressure, const Value &saltConcentration)
Definition: blackoilsolventmodules.hh:503
static const TabulatedFunction & pmisc(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:419
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:527
static bool isMiscible()
Definition: blackoilsolventmodules.hh:499
static const H2GasPvt & h2GasPvt()
Definition: blackoilsolventmodules.hh:369
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 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 & msfnKrsg(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:429
static const TabulatedFunction & ssfnKrs(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:389
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:142
static bool eqApplies(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:157
static const TabulatedFunction & msfnKro(const ElemContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:439
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 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: blackoilbioeffectsmodules.hh:43
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