28#ifndef EWOMS_BLACK_OIL_SOLVENT_MODULE_HH
29#define EWOMS_BLACK_OIL_SOLVENT_MODULE_HH
33#include <opm/common/Exceptions.hpp>
39#include <opm/material/fluidsystems/blackoilpvt/SolventPvt.hpp>
42#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
43#include <opm/input/eclipse/EclipseState/Tables/SsfnTable.hpp>
44#include <opm/input/eclipse/EclipseState/Tables/Sof2Table.hpp>
45#include <opm/input/eclipse/EclipseState/Tables/MsfnTable.hpp>
46#include <opm/input/eclipse/EclipseState/Tables/PmiscTable.hpp>
47#include <opm/input/eclipse/EclipseState/Tables/MiscTable.hpp>
48#include <opm/input/eclipse/EclipseState/Tables/SorwmisTable.hpp>
49#include <opm/input/eclipse/EclipseState/Tables/SgcwmisTable.hpp>
50#include <opm/input/eclipse/EclipseState/Tables/TlpmixpaTable.hpp>
53#include <opm/material/common/Valgrind.hpp>
55#include <dune/common/fvector.hh>
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 = getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>();
96 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
104 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
108 if (enableSolvent && !eclState.runspec().phases().active(Phase::SOLVENT))
109 throw std::runtime_error(
"Non-trivial solvent treatment requested at compile "
110 "time, but the deck does not contain the SOLVENT keyword");
111 else if (!enableSolvent && eclState.runspec().phases().active(Phase::SOLVENT))
112 throw std::runtime_error(
"Solvent treatment disabled at compile time, but the deck "
113 "contains the SOLVENT keyword");
115 if (!eclState.runspec().phases().active(Phase::SOLVENT))
118 params_.co2sol_ = eclState.runspec().co2Sol();
119 params_.h2sol_ = eclState.runspec().h2Sol();
122 throw std::runtime_error(
"CO2SOL and H2SOL can not be used together");
127 params_.co2GasPvt_.initFromState(eclState, schedule);
128 params_.brineCo2Pvt_.initFromState(eclState, schedule);
130 params_.h2GasPvt_.initFromState(eclState, schedule);
131 params_.brineH2Pvt_.initFromState(eclState, schedule);
133 if (eclState.getSimulationConfig().hasDISGASW()) {
134 params_.rsSolw_active_ =
true;
137 params_.solventPvt_.initFromState(eclState, schedule);
139 const auto& tableManager = eclState.getTableManager();
141 const auto& ssfnTables = tableManager.getSsfnTables();
142 unsigned numSatRegions = tableManager.getTabdims().getNumSatTables();
143 params_.setNumSatRegions(numSatRegions);
144 for (
unsigned satRegionIdx = 0; satRegionIdx < numSatRegions; ++ satRegionIdx) {
145 const auto& ssfnTable = ssfnTables.template getTable<SsfnTable>(satRegionIdx);
146 params_.ssfnKrg_[satRegionIdx].setXYContainers(ssfnTable.getSolventFractionColumn(),
147 ssfnTable.getGasRelPermMultiplierColumn(),
149 params_.ssfnKrs_[satRegionIdx].setXYContainers(ssfnTable.getSolventFractionColumn(),
150 ssfnTable.getSolventRelPermMultiplierColumn(),
155 params_.isMiscible_ =
false;
156 if (!eclState.getTableManager().getMiscTables().empty()) {
157 params_.isMiscible_ =
true;
159 unsigned numMiscRegions = 1;
162 const auto& sof2Tables = tableManager.getSof2Tables();
163 if (!sof2Tables.empty()) {
165 params_.sof2Krn_.resize(numSatRegions);
166 for (
unsigned satRegionIdx = 0; satRegionIdx < numSatRegions; ++ satRegionIdx) {
167 const auto& sof2Table = sof2Tables.template getTable<Sof2Table>(satRegionIdx);
168 params_.sof2Krn_[satRegionIdx].setXYContainers(sof2Table.getSoColumn(),
169 sof2Table.getKroColumn(),
173 else if(eclState.runspec().phases().active(Phase::OIL))
174 throw std::runtime_error(
"SOF2 must be specified in MISCIBLE (SOLVENT and OIL) runs\n");
176 const auto& miscTables = tableManager.getMiscTables();
177 if (!miscTables.empty()) {
178 assert(numMiscRegions == miscTables.size());
181 params_.misc_.resize(numMiscRegions);
182 for (
unsigned miscRegionIdx = 0; miscRegionIdx < numMiscRegions; ++miscRegionIdx) {
183 const auto& miscTable = miscTables.template getTable<MiscTable>(miscRegionIdx);
186 const auto& solventFraction = miscTable.getSolventFractionColumn();
187 const auto&
misc = miscTable.getMiscibilityColumn();
188 params_.misc_[miscRegionIdx].setXYContainers(solventFraction,
misc);
192 throw std::runtime_error(
"MISC must be specified in MISCIBLE (SOLVENT) runs\n");
195 params_.pmisc_.resize(numMiscRegions);
196 const auto& pmiscTables = tableManager.getPmiscTables();
197 if (!pmiscTables.empty()) {
198 assert(numMiscRegions == pmiscTables.size());
200 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
201 const auto& pmiscTable = pmiscTables.template getTable<PmiscTable>(regionIdx);
204 const auto& po = pmiscTable.getOilPhasePressureColumn();
205 const auto&
pmisc = pmiscTable.getMiscibilityColumn();
207 params_.pmisc_[regionIdx].setXYContainers(po,
pmisc);
211 std::vector<double> x = {0.0,1.0e20};
212 std::vector<double> y = {1.0,1.0};
213 TabulatedFunction constant = TabulatedFunction(2, x, y);
214 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
215 params_.pmisc_[regionIdx] = constant;
220 params_.msfnKrsg_.resize(numSatRegions);
221 params_.msfnKro_.resize(numSatRegions);
222 const auto& msfnTables = tableManager.getMsfnTables();
223 if (!msfnTables.empty()) {
224 assert(numSatRegions == msfnTables.size());
226 for (
unsigned regionIdx = 0; regionIdx < numSatRegions; ++regionIdx) {
227 const MsfnTable& msfnTable = msfnTables.template getTable<MsfnTable>(regionIdx);
231 const auto& Ssg = msfnTable.getGasPhaseFractionColumn();
232 const auto& krsg = msfnTable.getGasSolventRelpermMultiplierColumn();
233 const auto& kro = msfnTable.getOilRelpermMultiplierColumn();
235 params_.msfnKrsg_[regionIdx].setXYContainers(Ssg, krsg);
236 params_.msfnKro_[regionIdx].setXYContainers(Ssg, kro);
240 std::vector<double> x = {0.0,1.0};
241 std::vector<double> y = {1.0,0.0};
242 TabulatedFunction unit = TabulatedFunction(2, x, x);
243 TabulatedFunction invUnit = TabulatedFunction(2, x, y);
245 for (
unsigned regionIdx = 0; regionIdx < numSatRegions; ++regionIdx) {
246 params_.setMsfn(regionIdx, unit, invUnit);
250 params_.sorwmis_.resize(numMiscRegions);
251 const auto& sorwmisTables = tableManager.getSorwmisTables();
252 if (!sorwmisTables.empty()) {
253 assert(numMiscRegions == sorwmisTables.size());
255 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
256 const auto& sorwmisTable = sorwmisTables.template getTable<SorwmisTable>(regionIdx);
259 const auto& sw = sorwmisTable.getWaterSaturationColumn();
260 const auto&
sorwmis = sorwmisTable.getMiscibleResidualOilColumn();
262 params_.sorwmis_[regionIdx].setXYContainers(sw,
sorwmis);
267 std::vector<double> x = {0.0,1.0};
268 std::vector<double> y = {0.0,0.0};
269 TabulatedFunction zero = TabulatedFunction(2, x, y);
270 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
271 params_.sorwmis_[regionIdx] = zero;
276 params_.sgcwmis_.resize(numMiscRegions);
277 const auto& sgcwmisTables = tableManager.getSgcwmisTables();
278 if (!sgcwmisTables.empty()) {
279 assert(numMiscRegions == sgcwmisTables.size());
281 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
282 const auto& sgcwmisTable = sgcwmisTables.template getTable<SgcwmisTable>(regionIdx);
285 const auto& sw = sgcwmisTable.getWaterSaturationColumn();
286 const auto&
sgcwmis = sgcwmisTable.getMiscibleResidualGasColumn();
288 params_.sgcwmis_[regionIdx].setXYContainers(sw,
sgcwmis);
293 std::vector<double> x = {0.0,1.0};
294 std::vector<double> y = {0.0,0.0};
295 TabulatedFunction zero = TabulatedFunction(2, x, y);
296 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx)
297 params_.sgcwmis_[regionIdx] = zero;
300 const auto& tlmixpar = eclState.getTableManager().getTLMixpar();
301 if (!tlmixpar.empty()) {
303 params_.tlMixParamViscosity_.resize(numMiscRegions);
304 params_.tlMixParamDensity_.resize(numMiscRegions);
306 assert(numMiscRegions == tlmixpar.size());
307 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
308 const auto& tlp = tlmixpar[regionIdx];
309 params_.tlMixParamViscosity_[regionIdx] = tlp.viscosity_parameter;
310 params_.tlMixParamDensity_[regionIdx] = tlp.density_parameter;
314 throw std::runtime_error(
"TLMIXPAR must be specified in MISCIBLE (SOLVENT) runs\n");
317 params_.tlPMixTable_.resize(numMiscRegions);
318 if (!eclState.getTableManager().getTlpmixpaTables().empty()) {
319 const auto& tlpmixparTables = tableManager.getTlpmixpaTables();
320 if (!tlpmixparTables.empty()) {
321 assert(numMiscRegions == tlpmixparTables.size());
322 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx) {
323 const auto& tlpmixparTable = tlpmixparTables.template getTable<TlpmixpaTable>(regionIdx);
326 const auto& po = tlpmixparTable.getOilPhasePressureColumn();
327 const auto& tlpmixpa = tlpmixparTable.getMiscibilityColumn();
329 params_.tlPMixTable_[regionIdx].setXYContainers(po, tlpmixpa);
334 if (params_.pmisc_.size() > 0)
335 params_.tlPMixTable_ = params_.pmisc_;
337 throw std::invalid_argument(
"If the pressure dependent TL values in "
338 "TLPMIXPA is defaulted (no entries), then "
339 "the PMISC tables must be specified.");
344 std::vector<double> x = {0.0,1.0e20};
345 std::vector<double> y = {1.0,1.0};
346 TabulatedFunction ones = TabulatedFunction(2, x, y);
347 for (
unsigned regionIdx = 0; regionIdx < numMiscRegions; ++regionIdx)
348 params_.tlPMixTable_[regionIdx] = ones;
358 { params_.solventPvt_ = value; }
369 if constexpr (enableSolvent)
377 Simulator& simulator)
379 if constexpr (enableSolvent)
385 if constexpr (enableSolvent)
386 return pvIdx == solventSaturationIdx;
395 return "saturation_solvent";
403 return static_cast<Scalar
>(1.0);
408 if constexpr (enableSolvent)
409 return eqIdx == contiSolventEqIdx;
414 static std::string
eqName([[maybe_unused]]
unsigned eqIdx)
418 return "conti^solvent";
421 static Scalar
eqWeight([[maybe_unused]]
unsigned eqIdx)
426 return static_cast<Scalar
>(1.0);
429 template <
class LhsEval>
430 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
431 const IntensiveQuantities& intQuants)
433 if constexpr (enableSolvent) {
434 if constexpr (blackoilConserveSurfaceVolume) {
435 storage[contiSolventEqIdx] +=
436 Toolbox::template decay<LhsEval>(intQuants.porosity())
437 * Toolbox::template decay<LhsEval>(intQuants.solventSaturation())
438 * Toolbox::template decay<LhsEval>(intQuants.solventInverseFormationVolumeFactor());
440 storage[contiSolventEqIdx] += Toolbox::template decay<LhsEval>(intQuants.porosity())
441 * Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(waterPhaseIdx))
442 * Toolbox::template decay<LhsEval>(intQuants.fluidState().invB(waterPhaseIdx))
443 * Toolbox::template decay<LhsEval>(intQuants.rsSolw());
447 storage[contiSolventEqIdx] +=
448 Toolbox::template decay<LhsEval>(intQuants.porosity())
449 * Toolbox::template decay<LhsEval>(intQuants.solventSaturation())
450 * Toolbox::template decay<LhsEval>(intQuants.solventDensity());
452 storage[contiSolventEqIdx] += Toolbox::template decay<LhsEval>(intQuants.porosity())
453 * Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(waterPhaseIdx))
454 * Toolbox::template decay<LhsEval>(intQuants.fluidState().density(waterPhaseIdx))
455 * Toolbox::template decay<LhsEval>(intQuants.rsSolw());
463 [[maybe_unused]]
const ElementContext& elemCtx,
464 [[maybe_unused]]
unsigned scvfIdx,
465 [[maybe_unused]]
unsigned timeIdx)
468 if constexpr (enableSolvent) {
469 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
471 unsigned upIdx = extQuants.solventUpstreamIndex();
472 unsigned inIdx = extQuants.interiorIndex();
473 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
475 if constexpr (blackoilConserveSurfaceVolume) {
477 flux[contiSolventEqIdx] =
478 extQuants.solventVolumeFlux()
479 *up.solventInverseFormationVolumeFactor();
481 flux[contiSolventEqIdx] =
482 extQuants.solventVolumeFlux()
483 *decay<Scalar>(up.solventInverseFormationVolumeFactor());
488 flux[contiSolventEqIdx] +=
489 extQuants.volumeFlux(waterPhaseIdx)
490 * up.fluidState().invB(waterPhaseIdx)
493 flux[contiSolventEqIdx] +=
494 extQuants.volumeFlux(waterPhaseIdx)
495 *decay<Scalar>(up.fluidState().invB(waterPhaseIdx))
496 *decay<Scalar>(up.rsSolw());
501 flux[contiSolventEqIdx] =
502 extQuants.solventVolumeFlux()
503 *up.solventDensity();
505 flux[contiSolventEqIdx] =
506 extQuants.solventVolumeFlux()
507 *decay<Scalar>(up.solventDensity());
512 flux[contiSolventEqIdx] +=
513 extQuants.volumeFlux(waterPhaseIdx)
514 * up.fluidState().density(waterPhaseIdx)
517 flux[contiSolventEqIdx] +=
518 extQuants.volumeFlux(waterPhaseIdx)
519 *decay<Scalar>(up.fluidState().density(waterPhaseIdx))
520 *decay<Scalar>(up.rsSolw());
530 Scalar solventSaturation,
533 if constexpr (!enableSolvent) {
534 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Disabled);
539 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Ss);
540 priVars[solventSaturationIdx] = solventSaturation;
542 priVars.setPrimaryVarsMeaningSolvent(PrimaryVariables::SolventMeaning::Rsolw);
543 priVars[solventSaturationIdx] = solventRsw;
551 const PrimaryVariables& oldPv,
552 const EqVector& delta)
554 if constexpr (enableSolvent)
556 newPv[solventSaturationIdx] = oldPv[solventSaturationIdx] - delta[solventSaturationIdx];
568 return static_cast<Scalar
>(0.0);
577 return std::abs(Toolbox::scalarValue(resid[contiSolventEqIdx]));
580 template <
class DofEntity>
581 static void serializeEntity(
const Model& model, std::ostream& outstream,
const DofEntity& dof)
583 if constexpr (enableSolvent) {
584 unsigned dofIdx = model.dofMapper().index(dof);
586 const PrimaryVariables& priVars = model.solution(0)[dofIdx];
587 outstream << priVars[solventSaturationIdx];
591 template <
class DofEntity>
594 if constexpr (enableSolvent) {
595 unsigned dofIdx = model.dofMapper().index(dof);
597 PrimaryVariables& priVars0 = model.solution(0)[dofIdx];
598 PrimaryVariables& priVars1 = model.solution(1)[dofIdx];
600 instream >> priVars0[solventSaturationIdx];
603 priVars1 = priVars0[solventSaturationIdx];
609 return params_.solventPvt_;
615 return params_.co2GasPvt_;
620 return params_.h2GasPvt_;
625 return params_.brineCo2Pvt_;
630 return params_.brineH2Pvt_;
633 static const TabulatedFunction&
ssfnKrg(
const ElementContext& elemCtx,
637 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
638 return params_.ssfnKrg_[satnumRegionIdx];
641 static const TabulatedFunction&
ssfnKrs(
const ElementContext& elemCtx,
645 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
646 return params_.ssfnKrs_[satnumRegionIdx];
649 static const TabulatedFunction&
sof2Krn(
const ElementContext& elemCtx,
653 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
654 return params_.sof2Krn_[satnumRegionIdx];
657 static const TabulatedFunction&
misc(
const ElementContext& elemCtx,
661 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
662 return params_.misc_[miscnumRegionIdx];
665 static const TabulatedFunction&
pmisc(
const ElementContext& elemCtx,
669 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
670 return params_.pmisc_[miscnumRegionIdx];
673 static const TabulatedFunction&
msfnKrsg(
const ElementContext& elemCtx,
677 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
678 return params_.msfnKrsg_[satnumRegionIdx];
681 static const TabulatedFunction&
msfnKro(
const ElementContext& elemCtx,
685 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
686 return params_.msfnKro_[satnumRegionIdx];
689 static const TabulatedFunction&
sorwmis(
const ElementContext& elemCtx,
693 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
694 return params_.sorwmis_[miscnumRegionIdx];
697 static const TabulatedFunction&
sgcwmis(
const ElementContext& elemCtx,
701 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
702 return params_.sgcwmis_[miscnumRegionIdx];
705 static const TabulatedFunction&
tlPMixTable(
const ElementContext& elemCtx,
709 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
710 return params_.tlPMixTable_[miscnumRegionIdx];
717 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
718 return params_.tlMixParamViscosity_[miscnumRegionIdx];
725 unsigned miscnumRegionIdx = elemCtx.problem().miscnumRegionIndex(elemCtx, scvIdx, timeIdx);
726 return params_.tlMixParamDensity_[miscnumRegionIdx];
731 return params_.isMiscible_;
734 template <
class Value>
735 static const Value
solubilityLimit(
unsigned pvtIdx,
const Value& temperature,
const Value& pressure,
const Value& saltConcentration)
742 return brineCo2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature, pressure, saltConcentration);
744 return brineH2Pvt().saturatedGasDissolutionFactor(pvtIdx, temperature, pressure, saltConcentration);
749 return params_.rsSolw_active_;
754 return params_.co2sol_;
759 return params_.h2sol_;
766template <
class TypeTag,
bool enableSolventV>
767BlackOilSolventParams<typename BlackOilSolventModule<TypeTag, enableSolventV>::Scalar>
768BlackOilSolventModule<TypeTag, enableSolventV>::params_;
777template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
792 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
793 static constexpr int solventSaturationIdx = Indices::solventSaturationIdx;
794 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
795 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
796 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
797 static constexpr double cutOff = 1e-12;
811 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
813 auto& fs =
asImp_().fluidState_;
815 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
816 solventSaturation_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, elemCtx.linearizationType());
843 auto& fs =
asImp_().fluidState_;
849 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
850 if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Ss) {
852 }
else if (priVars.primaryVarsMeaningSolvent() == PrimaryVariables::SolventMeaning::Rsolw) {
853 rsSolw_ = priVars.makeEvaluation(solventSaturationIdx, timeIdx, elemCtx.linearizationType());
864 const Evaluation& p = FluidSystem::phaseIsActive(oilPhaseIdx)? fs.pressure(oilPhaseIdx) : fs.pressure(gasPhaseIdx);
866 const Evaluation& pgImisc = fs.pressure(gasPhaseIdx);
869 const auto& problem = elemCtx.problem();
870 Evaluation pgMisc = 0.0;
871 std::array<Evaluation, numPhases> pC;
872 const auto& materialParams = problem.materialLawParams(elemCtx, dofIdx, timeIdx);
873 MaterialLaw::capillaryPressures(pC, materialParams, fs);
876 const auto linearizationType = elemCtx.linearizationType();
877 if (priVars.primaryVarsMeaningPressure() == PrimaryVariables::PressureMeaning::Pg)
878 pgMisc = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx, linearizationType);
880 const Evaluation& po = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx, linearizationType);
881 pgMisc = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
884 fs.setPressure(gasPhaseIdx, pmisc * pgMisc + (1.0 - pmisc) * pgImisc);
890 if (gasSolventSat.value() < cutOff)
900 const Evaluation& p = FluidSystem::phaseIsActive(oilPhaseIdx)? fs.pressure(oilPhaseIdx) : fs.pressure(gasPhaseIdx);
901 const Evaluation miscibility = misc.eval(Fsolgas,
true) * pmisc.eval(p,
true);
904 unsigned cellIdx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
905 const auto& materialLawManager = elemCtx.problem().materialLawManager();
906 const auto& scaledDrainageInfo =
907 materialLawManager->oilWaterScaledEpsInfoDrainage(cellIdx);
909 const Scalar& sogcr = scaledDrainageInfo.Sogcr;
910 Evaluation sor = sogcr;
911 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
912 const Evaluation& sw = fs.saturation(waterPhaseIdx);
914 sor = miscibility * sorwmis.eval(sw,
true) + (1.0 - miscibility) * sogcr;
916 const Scalar& sgcr = scaledDrainageInfo.Sgcr;
917 Evaluation sgc = sgcr;
918 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
919 const Evaluation& sw = fs.saturation(waterPhaseIdx);
921 sgc = miscibility * sgcwmis.eval(sw,
true) + (1.0 - miscibility) * sgcr;
924 Evaluation oilGasSolventSat = gasSolventSat;
925 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
926 oilGasSolventSat += fs.saturation(oilPhaseIdx);
928 const Evaluation zero = 0.0;
929 const Evaluation oilGasSolventEffSat = std::max(oilGasSolventSat - sor - sgc, zero);
931 Evaluation F_totalGas = 0.0;
932 if (oilGasSolventEffSat.value() > cutOff) {
933 const Evaluation gasSolventEffSat = std::max(gasSolventSat - sgc, zero);
934 F_totalGas = gasSolventEffSat / oilGasSolventEffSat;
940 const Evaluation mkrgt = msfnKrsg.eval(F_totalGas,
true) * sof2Krn.eval(oilGasSolventSat,
true);
941 const Evaluation mkro = msfnKro.eval(F_totalGas,
true) * sof2Krn.eval(oilGasSolventSat,
true);
943 Evaluation& kro =
asImp_().mobility_[oilPhaseIdx];
944 Evaluation& krg =
asImp_().mobility_[gasPhaseIdx];
947 krg *= (1.0 - miscibility);
948 krg += miscibility * mkrgt;
949 kro *= (1.0 - miscibility);
950 kro += miscibility * mkro;
957 Evaluation& krg =
asImp_().mobility_[gasPhaseIdx];
959 krg *= ssfnKrg.eval(Fhydgas,
true);
973 const auto& iq =
asImp_();
974 unsigned pvtRegionIdx = iq.pvtRegionIndex();
975 const Evaluation& T = iq.fluidState().temperature(gasPhaseIdx);
976 const Evaluation& p = iq.fluidState().pressure(gasPhaseIdx);
978 const Evaluation rv = 0.0;
979 const Evaluation rvw = 0.0;
988 auto& fs =
asImp_().fluidState_;
989 const auto& bw = brineCo2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p,
rsSolw());
991 const auto denw = bw*brineCo2Pvt.waterReferenceDensity(pvtRegionIdx)
992 +
rsSolw()*bw*brineCo2Pvt.gasReferenceDensity(pvtRegionIdx);
993 fs.setDensity(waterPhaseIdx, denw);
994 fs.setInvB(waterPhaseIdx, bw);
995 Evaluation& mobw =
asImp_().mobility_[waterPhaseIdx];
996 const auto& muWat = fs.viscosity(waterPhaseIdx);
997 const auto& muWatEff = brineCo2Pvt.viscosity(pvtRegionIdx, T, p,
rsSolw());
998 mobw *= muWat / muWatEff;
1006 auto& fs =
asImp_().fluidState_;
1007 const auto& bw = brineH2Pvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p,
rsSolw());
1009 const auto denw = bw*brineH2Pvt.waterReferenceDensity(pvtRegionIdx)
1010 +
rsSolw()*bw*brineH2Pvt.gasReferenceDensity(pvtRegionIdx);
1011 fs.setDensity(waterPhaseIdx, denw);
1012 fs.setInvB(waterPhaseIdx, bw);
1013 Evaluation& mobw =
asImp_().mobility_[waterPhaseIdx];
1014 const auto& muWat = fs.viscosity(waterPhaseIdx);
1015 const auto& muWatEff = brineH2Pvt.viscosity(pvtRegionIdx, T, p,
rsSolw());
1016 mobw *= muWat / muWatEff;
1026 effectiveProperties(elemCtx, scvIdx, timeIdx);
1058 void effectiveProperties(
const ElementContext& elemCtx,
1072 auto& fs =
asImp_().fluidState_;
1075 Evaluation oilEffSat = 0.0;
1076 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1077 oilEffSat = fs.saturation(oilPhaseIdx);
1079 Evaluation gasEffSat = fs.saturation(gasPhaseIdx);
1081 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
1084 const Evaluation zero = 0.0;
1085 const Evaluation& sw = fs.saturation(waterPhaseIdx);
1086 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1087 oilEffSat = std::max(oilEffSat - sorwmis.eval(sw,
true), zero);
1089 gasEffSat = std::max(gasEffSat - sgcwmis.eval(sw,
true), zero);
1090 solventEffSat = std::max(solventEffSat - sgcwmis.eval(sw,
true), zero);
1092 const Evaluation oilGasSolventEffSat = oilEffSat + gasEffSat + solventEffSat;
1093 const Evaluation oilSolventEffSat = oilEffSat + solventEffSat;
1094 const Evaluation solventGasEffSat = solventEffSat + gasEffSat;
1101 const Evaluation& p = FluidSystem::phaseIsActive(oilPhaseIdx)? fs.pressure(oilPhaseIdx) : fs.pressure(gasPhaseIdx);
1104 const Evaluation pmisc = pmiscTable.eval(p,
true);
1108 const Evaluation& muGas = fs.viscosity(gasPhaseIdx);
1111 assert(muGas.value() > 0);
1112 assert(muSolvent.value() > 0);
1113 const Evaluation muGasPow = pow(muGas, 0.25);
1114 const Evaluation muSolventPow = pow(muSolvent, 0.25);
1116 Evaluation muMixSolventGas = muGas;
1117 if (solventGasEffSat > cutOff)
1118 muMixSolventGas *= muSolvent / pow(((gasEffSat / solventGasEffSat) * muSolventPow) + ((solventEffSat / solventGasEffSat) * muGasPow) , 4.0);
1120 Evaluation muOil = 1.0;
1121 Evaluation muOilPow = 1.0;
1122 Evaluation muMixOilSolvent = 1.0;
1123 Evaluation muOilEff = 1.0;
1124 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1125 muOil = fs.viscosity(oilPhaseIdx);
1126 assert(muOil.value() > 0);
1127 muOilPow = pow(muOil, 0.25);
1128 muMixOilSolvent = muOil;
1129 if (oilSolventEffSat > cutOff)
1130 muMixOilSolvent *= muSolvent / pow(((oilEffSat / oilSolventEffSat) * muSolventPow) + ((solventEffSat / oilSolventEffSat) * muOilPow) , 4.0);
1132 muOilEff = pow(muOil,1.0 - tlMixParamMu) * pow(muMixOilSolvent, tlMixParamMu);
1134 Evaluation muMixSolventGasOil = muOil;
1135 if (oilGasSolventEffSat > cutOff)
1136 muMixSolventGasOil *= muSolvent * muGas / pow(((oilEffSat / oilGasSolventEffSat) * muSolventPow * muGasPow)
1137 + ((solventEffSat / oilGasSolventEffSat) * muOilPow * muGasPow) + ((gasEffSat / oilGasSolventEffSat) * muSolventPow * muOilPow), 4.0);
1139 Evaluation muGasEff = pow(muGas,1.0 - tlMixParamMu) * pow(muMixSolventGas, tlMixParamMu);
1140 Evaluation muSolventEff = pow(muSolvent,1.0 - tlMixParamMu) * pow(muMixSolventGasOil, tlMixParamMu);
1143 const Evaluation& rhoGas = fs.density(gasPhaseIdx);
1144 Evaluation rhoOil = 0.0;
1145 if (FluidSystem::phaseIsActive(oilPhaseIdx))
1146 rhoOil = fs.density(oilPhaseIdx);
1157 const Evaluation muOilEffPow = pow(pow(muOil, 1.0 - tlMixParamRho) * pow(muMixOilSolvent, tlMixParamRho), 0.25);
1158 const Evaluation muGasEffPow = pow(pow(muGas, 1.0 - tlMixParamRho) * pow(muMixSolventGas, tlMixParamRho), 0.25);
1159 const Evaluation muSolventEffPow = pow(pow(muSolvent, 1.0 - tlMixParamRho) * pow(muMixSolventGasOil, tlMixParamRho), 0.25);
1161 const Evaluation oilGasEffSaturation = oilEffSat + gasEffSat;
1162 Evaluation sof = 0.0;
1163 Evaluation sgf = 0.0;
1164 if (oilGasEffSaturation.value() > cutOff) {
1165 sof = oilEffSat / oilGasEffSaturation;
1166 sgf = gasEffSat / oilGasEffSaturation;
1169 const Evaluation muSolventOilGasPow = muSolventPow * ((sgf * muOilPow) + (sof * muGasPow));
1171 Evaluation rhoMixSolventGasOil = 0.0;
1172 if (oilGasSolventEffSat.value() > cutOff)
1173 rhoMixSolventGasOil = (rhoOil * oilEffSat / oilGasSolventEffSat) + (rhoGas * gasEffSat / oilGasSolventEffSat) + (rhoSolvent * solventEffSat / oilGasSolventEffSat);
1175 Evaluation rhoGasEff = 0.0;
1176 if (std::abs(muSolventPow.value() - muGasPow.value()) < cutOff)
1177 rhoGasEff = ((1.0 - tlMixParamRho) * rhoGas) + (tlMixParamRho * rhoMixSolventGasOil);
1179 const Evaluation solventGasEffFraction = (muGasPow * (muSolventPow - muGasEffPow)) / (muGasEffPow * (muSolventPow - muGasPow));
1180 rhoGasEff = (rhoGas * solventGasEffFraction) + (rhoSolvent * (1.0 - solventGasEffFraction));
1183 Evaluation rhoSolventEff = 0.0;
1184 if (std::abs((muSolventOilGasPow.value() - (muOilPow.value() * muGasPow.value()))) < cutOff)
1185 rhoSolventEff = ((1.0 - tlMixParamRho) * rhoSolvent) + (tlMixParamRho * rhoMixSolventGasOil);
1187 const Evaluation sfraction_se = (muSolventOilGasPow - (muOilPow * muGasPow * muSolventPow / muSolventEffPow)) / (muSolventOilGasPow - (muOilPow * muGasPow));
1188 rhoSolventEff = (rhoSolvent * sfraction_se) + (rhoGas * sgf * (1.0 - sfraction_se)) + (rhoOil * sof * (1.0 - sfraction_se));
1192 unsigned pvtRegionIdx =
asImp_().pvtRegionIndex();
1193 Evaluation bGasEff = rhoGasEff;
1194 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1195 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) + FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) * fs.Rv());
1197 bGasEff /= (FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1202 const Evaluation bg = fs.invB(gasPhaseIdx);
1204 const Evaluation bg_eff = pmisc * bGasEff + (1.0 - pmisc) * bg;
1205 const Evaluation bs_eff = pmisc * bSolventEff + (1.0 - pmisc) * bs;
1208 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1209 fs.setDensity(gasPhaseIdx,
1211 *(FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx)
1212 + FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx)*fs.Rv()));
1214 fs.setDensity(gasPhaseIdx,
1216 *FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx));
1221 Evaluation& mobg =
asImp_().mobility_[gasPhaseIdx];
1222 muGasEff = bg_eff / (pmisc * bGasEff / muGasEff + (1.0 - pmisc) * bg / muGas);
1223 mobg *= muGas / muGasEff;
1226 solventViscosity_ = bs_eff / (pmisc * bSolventEff / muSolventEff + (1.0 - pmisc) * bs / muSolvent);
1228 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1229 Evaluation rhoOilEff = 0.0;
1230 if (std::abs(muOilPow.value() - muSolventPow.value()) < cutOff) {
1231 rhoOilEff = ((1.0 - tlMixParamRho) * rhoOil) + (tlMixParamRho * rhoMixSolventGasOil);
1234 const Evaluation solventOilEffFraction = (muOilPow * (muOilEffPow - muSolventPow)) / (muOilEffPow * (muOilPow - muSolventPow));
1235 rhoOilEff = (rhoOil * solventOilEffFraction) + (rhoSolvent * (1.0 - solventOilEffFraction));
1237 const Evaluation bOilEff = rhoOilEff / (FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx) + FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx) * fs.Rs());
1238 const Evaluation bo = fs.invB(oilPhaseIdx);
1239 const Evaluation bo_eff = pmisc * bOilEff + (1.0 - pmisc) * bo;
1240 fs.setDensity(oilPhaseIdx,
1242 *(FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx)
1243 + FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx)*fs.Rs()));
1246 Evaluation& mobo =
asImp_().mobility_[oilPhaseIdx];
1247 muOilEff = bo_eff / (pmisc * bOilEff / muOilEff + (1.0 - pmisc) * bo / muOil);
1248 mobo *= muOil / muOilEff;
1254 {
return *
static_cast<Implementation*
>(
this); }
1267template <
class TypeTag>
1292 {
throw std::runtime_error(
"solventSaturation() called but solvents are disabled"); }
1295 {
throw std::runtime_error(
"rsSolw() called but solvents are disabled"); }
1298 {
throw std::runtime_error(
"solventDensity() called but solvents are disabled"); }
1301 {
throw std::runtime_error(
"solventViscosity() called but solvents are disabled"); }
1304 {
throw std::runtime_error(
"solventMobility() called but solvents are disabled"); }
1307 {
throw std::runtime_error(
"solventInverseFormationVolumeFactor() called but solvents are disabled"); }
1310 {
throw std::runtime_error(
"solventRefDensity() called but solvents are disabled"); }
1320template <class TypeTag, bool enableSolventV = getPropValue<TypeTag, Properties::EnableSolvent>()>
1333 using Toolbox = MathToolbox<Evaluation>;
1335 static constexpr unsigned gasPhaseIdx = FluidSystem::gasPhaseIdx;
1336 static constexpr int dimWorld = GridView::dimensionworld;
1338 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
1339 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
1346 template <
class Dummy =
bool>
1352 const auto& gradCalc = elemCtx.gradientCalculator();
1355 const auto& scvf = elemCtx.stencil(timeIdx).interiorFace(scvfIdx);
1356 const auto& faceNormal = scvf.normal();
1358 unsigned i = scvf.interiorIndex();
1359 unsigned j = scvf.exteriorIndex();
1362 DimEvalVector solventPGrad;
1364 gradCalc.calculateGradient(solventPGrad,
1368 Valgrind::CheckDefined(solventPGrad);
1371 if (Parameters::get<TypeTag, Properties::EnableGravity>()) {
1374 const auto& gIn = elemCtx.problem().gravity(elemCtx, i, timeIdx);
1375 const auto& gEx = elemCtx.problem().gravity(elemCtx, j, timeIdx);
1377 const auto& intQuantsIn = elemCtx.intensiveQuantities(i, timeIdx);
1378 const auto& intQuantsEx = elemCtx.intensiveQuantities(j, timeIdx);
1380 const auto& posIn = elemCtx.pos(i, timeIdx);
1381 const auto& posEx = elemCtx.pos(j, timeIdx);
1382 const auto& posFace = scvf.integrationPos();
1385 DimVector distVecIn(posIn);
1386 DimVector distVecEx(posEx);
1387 DimVector distVecTotal(posEx);
1389 distVecIn -= posFace;
1390 distVecEx -= posFace;
1391 distVecTotal -= posIn;
1392 Scalar absDistTotalSquared = distVecTotal.two_norm2();
1395 auto rhoIn = intQuantsIn.solventDensity();
1396 auto pStatIn = - rhoIn*(gIn*distVecIn);
1400 Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1401 Scalar pStatEx = - rhoEx*(gEx*distVecEx);
1407 DimEvalVector f(distVecTotal);
1408 f *= (pStatEx - pStatIn)/absDistTotalSquared;
1411 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
1412 solventPGrad[dimIdx] += f[dimIdx];
1414 if (!isfinite(solventPGrad[dimIdx]))
1415 throw NumericalProblem(
"Non-finite potential gradient for solvent 'phase'");
1420 Evaluation solventPGradNormal = 0.0;
1421 for (
unsigned dimIdx = 0; dimIdx < faceNormal.size(); ++dimIdx)
1422 solventPGradNormal += solventPGrad[dimIdx]*faceNormal[dimIdx];
1424 if (solventPGradNormal > 0) {
1425 solventUpstreamDofIdx_ = j;
1426 solventDownstreamDofIdx_ = i;
1429 solventUpstreamDofIdx_ = i;
1430 solventDownstreamDofIdx_ = j;
1433 const auto& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1439 if (solventUpstreamDofIdx_ == i)
1440 solventVolumeFlux_ = solventPGradNormal*up.solventMobility();
1442 solventVolumeFlux_ = solventPGradNormal*scalarValue(up.solventMobility());
1449 template <
class Dummy =
bool>
1455 const ExtensiveQuantities& extQuants = asImp_();
1457 unsigned interiorDofIdx = extQuants.interiorIndex();
1458 unsigned exteriorDofIdx = extQuants.exteriorIndex();
1459 assert(interiorDofIdx != exteriorDofIdx);
1461 const auto& intQuantsIn = elemCtx.intensiveQuantities(interiorDofIdx, timeIdx);
1462 const auto& intQuantsEx = elemCtx.intensiveQuantities(exteriorDofIdx, timeIdx);
1464 unsigned I = elemCtx.globalSpaceIndex(interiorDofIdx, timeIdx);
1465 unsigned J = elemCtx.globalSpaceIndex(exteriorDofIdx, timeIdx);
1467 Scalar thpres = elemCtx.problem().thresholdPressure(I, J);
1468 Scalar trans = elemCtx.problem().transmissibility(elemCtx, interiorDofIdx, exteriorDofIdx);
1469 Scalar g = elemCtx.problem().gravity()[dimWorld - 1];
1471 Scalar zIn = elemCtx.problem().dofCenterDepth(elemCtx, interiorDofIdx, timeIdx);
1472 Scalar zEx = elemCtx.problem().dofCenterDepth(elemCtx, exteriorDofIdx, timeIdx);
1473 Scalar distZ = zIn - zEx;
1475 const Evaluation& rhoIn = intQuantsIn.solventDensity();
1476 Scalar rhoEx = Toolbox::value(intQuantsEx.solventDensity());
1477 const Evaluation& rhoAvg = rhoIn*0.5 + rhoEx*0.5;
1479 const Evaluation& pressureInterior = intQuantsIn.fluidState().pressure(gasPhaseIdx);
1480 Evaluation pressureExterior = Toolbox::value(intQuantsEx.fluidState().pressure(gasPhaseIdx));
1481 pressureExterior += distZ*g*rhoAvg;
1483 Evaluation pressureDiffSolvent = pressureExterior - pressureInterior;
1484 if (std::abs(scalarValue(pressureDiffSolvent)) > thpres) {
1485 if (pressureDiffSolvent < 0.0)
1486 pressureDiffSolvent += thpres;
1488 pressureDiffSolvent -= thpres;
1491 pressureDiffSolvent = 0.0;
1493 if (pressureDiffSolvent > 0.0) {
1494 solventUpstreamDofIdx_ = exteriorDofIdx;
1495 solventDownstreamDofIdx_ = interiorDofIdx;
1497 else if (pressureDiffSolvent < 0.0) {
1498 solventUpstreamDofIdx_ = interiorDofIdx;
1499 solventDownstreamDofIdx_ = exteriorDofIdx;
1505 solventUpstreamDofIdx_ = std::min(interiorDofIdx, exteriorDofIdx);
1506 solventDownstreamDofIdx_ = std::max(interiorDofIdx, exteriorDofIdx);
1507 solventVolumeFlux_ = 0.0;
1511 Scalar faceArea = elemCtx.stencil(timeIdx).interiorFace(scvfIdx).area();
1512 const IntensiveQuantities& up = elemCtx.intensiveQuantities(solventUpstreamDofIdx_, timeIdx);
1513 if (solventUpstreamDofIdx_ == interiorDofIdx)
1514 solventVolumeFlux_ =
1515 up.solventMobility()
1517 *pressureDiffSolvent;
1519 solventVolumeFlux_ =
1520 scalarValue(up.solventMobility())
1522 *pressureDiffSolvent;
1526 {
return solventUpstreamDofIdx_; }
1529 {
return solventDownstreamDofIdx_; }
1532 {
return solventVolumeFlux_; }
1538 Implementation& asImp_()
1539 {
return *
static_cast<Implementation*
>(
this); }
1541 Evaluation solventVolumeFlux_;
1542 unsigned solventUpstreamDofIdx_;
1543 unsigned solventDownstreamDofIdx_;
1546template <
class TypeTag>
1564 {
throw std::runtime_error(
"solventUpstreamIndex() called but solvents are disabled"); }
1567 {
throw std::runtime_error(
"solventDownstreamIndex() called but solvents are disabled"); }
1570 {
throw std::runtime_error(
"solventVolumeFlux() called but solvents are disabled"); }
1573 {
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:1563
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1566
void setSolventVolumeFlux(const Evaluation &)
Definition: blackoilsolventmodules.hh:1572
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1569
void updateVolumeFluxPerm(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1553
void updateVolumeFluxTrans(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1558
Provides the solvent specific extensive quantities to the generic black-oil module's extensive quanti...
Definition: blackoilsolventmodules.hh:1322
unsigned solventUpstreamIndex() const
Definition: blackoilsolventmodules.hh:1525
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:1348
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:1451
unsigned solventDownstreamIndex() const
Definition: blackoilsolventmodules.hh:1528
const Evaluation & solventVolumeFlux() const
Definition: blackoilsolventmodules.hh:1531
void setSolventVolumeFlux(const Evaluation &solventVolumeFlux)
Definition: blackoilsolventmodules.hh:1534
const Scalar & solventRefDensity() const
Definition: blackoilsolventmodules.hh:1309
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:1297
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:1306
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:1294
void solventPreSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1276
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:1300
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:1303
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:1291
void solventPvtUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1286
void solventPostSatFuncUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilsolventmodules.hh:1281
Provides the volumetric quantities required for the equations needed by the solvents extension of the...
Definition: blackoilsolventmodules.hh:779
void solventPreSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called before the saturation functions are doing their magic.
Definition: blackoilsolventmodules.hh:807
const Evaluation & solventMobility() const
Definition: blackoilsolventmodules.hh:1045
Evaluation solventViscosity_
Definition: blackoilsolventmodules.hh:1260
Evaluation solventMobility_
Definition: blackoilsolventmodules.hh:1261
Scalar solventRefDensity_
Definition: blackoilsolventmodules.hh:1264
Evaluation solventDensity_
Definition: blackoilsolventmodules.hh:1259
Evaluation rsSolw_
Definition: blackoilsolventmodules.hh:1258
Implementation & asImp_()
Definition: blackoilsolventmodules.hh:1253
Evaluation solventSaturation_
Definition: blackoilsolventmodules.hh:1257
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:969
const Evaluation & solventDensity() const
Definition: blackoilsolventmodules.hh:1039
Evaluation solventInvFormationVolumeFactor_
Definition: blackoilsolventmodules.hh:1262
const Evaluation & rsSolw() const
Definition: blackoilsolventmodules.hh:1036
const Evaluation & solventViscosity() const
Definition: blackoilsolventmodules.hh:1042
void solventPostSatFuncUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Called after the saturation functions have been doing their magic.
Definition: blackoilsolventmodules.hh:837
const Evaluation & solventInverseFormationVolumeFactor() const
Definition: blackoilsolventmodules.hh:1048
const Evaluation & solventSaturation() const
Definition: blackoilsolventmodules.hh:1033
const Scalar & solventRefDensity() const
Definition: blackoilsolventmodules.hh:1052
Evaluation hydrocarbonSaturation_
Definition: blackoilsolventmodules.hh:1256
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:357
static const TabulatedFunction & sof2Krn(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:649
static const TabulatedFunction & ssfnKrs(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:641
static const TabulatedFunction & pmisc(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:665
static const TabulatedFunction & msfnKrsg(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:673
static const TabulatedFunction & sgcwmis(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:697
static Scalar computeResidualError(const EqVector &resid)
Return how much a residual is considered an error.
Definition: blackoilsolventmodules.hh:574
static const BrineH2Pvt & brineH2Pvt()
Definition: blackoilsolventmodules.hh:628
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:398
static std::string eqName(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:414
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:421
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilsolventmodules.hh:581
static bool isSolubleInWater()
Definition: blackoilsolventmodules.hh:747
static bool isCO2Sol()
Definition: blackoilsolventmodules.hh:752
static const TabulatedFunction & msfnKro(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:681
static void registerOutputModules(Model &model, Simulator &simulator)
Register all solvent specific VTK and ECL output modules.
Definition: blackoilsolventmodules.hh:376
static const TabulatedFunction & ssfnKrg(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:633
static void setIsMiscible(const bool isMiscible)
Definition: blackoilsolventmodules.hh:361
static void registerParameters()
Register all run-time parameters for the black-oil solvent module.
Definition: blackoilsolventmodules.hh:367
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilsolventmodules.hh:562
static const SolventPvt & solventPvt()
Definition: blackoilsolventmodules.hh:607
static bool isH2Sol()
Definition: blackoilsolventmodules.hh:757
static bool isMiscible()
Definition: blackoilsolventmodules.hh:729
static const H2GasPvt & h2GasPvt()
Definition: blackoilsolventmodules.hh:618
static const Scalar & tlMixParamDensity(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:721
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:550
static const TabulatedFunction & tlPMixTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:705
static const Scalar & tlMixParamViscosity(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:713
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar solventSaturation, Scalar solventRsw)
Assign the solvent specific primary variables to a PrimaryVariables object.
Definition: blackoilsolventmodules.hh:529
static const TabulatedFunction & sorwmis(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:689
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:391
static bool eqApplies(unsigned eqIdx)
Definition: blackoilsolventmodules.hh:406
static const Value solubilityLimit(unsigned pvtIdx, const Value &temperature, const Value &pressure, const Value &saltConcentration)
Definition: blackoilsolventmodules.hh:735
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilsolventmodules.hh:430
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilsolventmodules.hh:592
static const Co2GasPvt & co2GasPvt()
Definition: blackoilsolventmodules.hh:613
static const BrineCo2Pvt & brineCo2Pvt()
Definition: blackoilsolventmodules.hh:623
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:462
static const TabulatedFunction & misc(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:657
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilsolventmodules.hh:383
Callback class for a phase pressure.
Definition: quantitycallbacks.hh:84
void setPhaseIndex(unsigned phaseIdx)
Set the index of the fluid phase for which the pressure should be returned.
Definition: quantitycallbacks.hh:108
VTK output module for the black oil model's solvent related quantities.
Definition: vtkblackoilsolventmodule.hh:86
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:242
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.hh:43
::Opm::SolventPvt< Scalar > SolventPvt
Definition: blackoilsolventparams.hh:46
::Opm::H2GasPvt< Scalar > H2GasPvt
Definition: blackoilsolventparams.hh:52
::Opm::BrineCo2Pvt< Scalar > BrineCo2Pvt
Definition: blackoilsolventparams.hh:55
::Opm::Co2GasPvt< Scalar > Co2GasPvt
Definition: blackoilsolventparams.hh:49
Tabulated1DFunction< Scalar > TabulatedFunction
Definition: blackoilsolventparams.hh:44
::Opm::BrineH2Pvt< Scalar > BrineH2Pvt
Definition: blackoilsolventparams.hh:58