95 using FluidState =
typename IntensiveQuantities::FluidState;
97 using Element =
typename GridView::template Codim<0>::Entity;
98 using ElementIterator =
typename GridView::template Codim<0>::Iterator;
101 using Dir = FaceDir::DirEnum;
108 static constexpr int conti0EqIdx = Indices::conti0EqIdx;
109 static constexpr int numPhases = FluidSystem::numPhases;
110 static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
111 static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
112 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
113 static constexpr int gasCompIdx = FluidSystem::gasCompIdx;
114 static constexpr int oilCompIdx = FluidSystem::oilCompIdx;
115 static constexpr int waterCompIdx = FluidSystem::waterCompIdx;
119 template<
int idx,
class VectorType>
120 static Scalar value_or_zero(
const VectorType&
v)
122 if constexpr (
idx == -1) {
125 return v.empty() ? 0.0 :
v[
idx];
133 :
BaseType(simulator.vanguard().eclState(),
134 simulator.vanguard().schedule(),
136 simulator.vanguard().summaryState(),
138 [
this](
const int idx)
139 {
return simulator_.problem().eclWriter().collectOnIORank().localIdxToGlobalIdx(
idx); },
140 simulator.vanguard().grid().comm(),
151 , simulator_(simulator)
152 , collectOnIORank_(collectOnIORank)
158 auto isCartIdxOnThisRank = [&collectOnIORank](
const int idx) {
159 return collectOnIORank.isCartIdxOnThisRank(
idx);
162 this->setupBlockData(isCartIdxOnThisRank);
164 if (! Parameters::Get<Parameters::OwnerCellsFirst>()) {
165 const std::string
msg =
"The output code does not support --owner-cells-first=false.";
166 if (collectOnIORank.isIORank()) {
173 auto rset = this->eclState_.fieldProps().fip_regions();
174 rset.push_back(
"PVTNUM");
179 this->regionAvgDensity_
180 .emplace(this->simulator_.gridView().comm(),
181 FluidSystem::numPhases,
rset,
182 [
fp = std::cref(
this->eclState_.fieldProps())]
183 (
const std::string&
rsetName) ->
decltype(
auto)
184 { return fp.get().get_int(rsetName); });
194 const unsigned reportStepNum,
197 const bool isRestart)
203 const auto& problem = this->simulator_.problem();
210 &problem.materialLawManager()->hysteresisConfig(),
211 problem.eclWriter().getOutputNnc().size());
216 const int reportStepNum)
218 this->setupElementExtractors_();
219 this->setupBlockExtractors_(
isSubStep, reportStepNum);
225 this->extractors_.clear();
226 this->blockExtractors_.clear();
227 this->extraBlockExtractors_.clear();
241 if (this->extractors_.empty()) {
245 const auto&
matLawManager = simulator_.problem().materialLawManager();
248 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
249 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
250 const auto& fs = intQuants.fluidState();
253 elemCtx.globalSpaceIndex(dofIdx, 0),
254 elemCtx.primaryVars(dofIdx, 0).pvtRegionIndex(),
262 if (FluidSystem::phaseIsActive(oilPhaseIdx) && FluidSystem::phaseIsActive(waterPhaseIdx)) {
268 if (FluidSystem::phaseIsActive(oilPhaseIdx) && FluidSystem::phaseIsActive(gasPhaseIdx)) {
280 void processElementBlockData(
const ElementContext& elemCtx)
287 if (this->blockExtractors_.empty() &&
this->extraBlockExtractors_.empty()) {
291 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
293 const auto globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
294 const auto cartesianIdx = elemCtx.simulator().vanguard().cartesianIndex(globalDofIdx);
298 if (
be_it == this->blockExtractors_.end() &&
304 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
305 const auto& fs = intQuants.fluidState();
307 const typename BlockExtractor::Context
ectx{
315 if (
be_it != this->blockExtractors_.end()) {
318 if (
bee_it != this->extraBlockExtractors_.end()) {
325 const std::size_t reportStepNum,
329 const Parallel::Communication& comm)
332 if (comm.rank() != 0) {
337 std::unique_ptr<FIPConfig>
fipSched;
338 if (reportStepNum > 0) {
339 const auto&
rpt = this->schedule_[reportStepNum-1].rpt_config.get();
342 const FIPConfig&
fipc = reportStepNum == 0 ? this->eclState_.getEclipseConfig().fip()
345 if (!
substep && !this->forceDisableFipOutput_ &&
fipc.output(FIPConfig::OutputField::FIELD)) {
347 this->logOutput_.timeStamp(
"BALANCE", elapsed, reportStepNum,
currentDate);
352 if (
fipc.output(FIPConfig::OutputField::FIPNUM)) {
355 if (
fipc.output(FIPConfig::OutputField::RESV))
356 this->logOutput_.fipResv(
inplace,
"FIPNUM");
359 if (
fipc.output(FIPConfig::OutputField::FIP)) {
360 for (
const auto&
reg :
this->regions_) {
361 if (
reg.first !=
"FIPNUM") {
362 std::ostringstream
ss;
363 ss <<
"BAL" <<
reg.first.substr(3);
364 this->logOutput_.timeStamp(
ss.str(), elapsed, reportStepNum,
currentDate);
367 if (
fipc.output(FIPConfig::OutputField::RESV))
403 template <
class ActiveIndex,
class CartesianIndex>
409 const auto identifyCell = [&activeIndex, &cartesianIndex](
const Element&
elem)
417 elem.partitionType() == Dune::InteriorEntity
422 const auto& stencil = elemCtx.stencil(
timeIdx);
423 const auto numInteriorFaces = elemCtx.numInteriorFaces(
timeIdx);
430 const auto rates = this->
445 this->interRegionFlows_.
clear();
461 return this->interRegionFlows_;
464 template <
class Flu
idState>
465 void assignToFluidState(FluidState& fs,
unsigned elemIdx)
const
468 if (this->saturation_[
phaseIdx].empty())
474 if (!this->fluidPressure_.empty()) {
477 std::array<Scalar, numPhases>
pc = {0};
478 const MaterialLawParams&
matParams = simulator_.problem().materialLawParams(elemIdx);
479 MaterialLaw::capillaryPressures(
pc,
matParams, fs);
480 Valgrind::CheckDefined(this->fluidPressure_[elemIdx]);
481 Valgrind::CheckDefined(
pc);
482 const auto& pressure = this->fluidPressure_[elemIdx];
484 if (!FluidSystem::phaseIsActive(
phaseIdx))
487 if (Indices::oilEnabled)
489 else if (Indices::gasEnabled)
491 else if (Indices::waterEnabled)
497 if (!this->temperature_.empty())
498 fs.setTemperature(this->temperature_[elemIdx]);
499 if (!this->rs_.empty())
500 fs.setRs(this->rs_[elemIdx]);
501 if (!this->rsw_.empty())
502 fs.setRsw(this->rsw_[elemIdx]);
503 if (!this->rv_.empty())
504 fs.setRv(this->rv_[elemIdx]);
505 if (!this->rvw_.empty())
506 fs.setRvw(this->rvw_[elemIdx]);
509 void initHysteresisParams(Simulator& simulator,
unsigned elemIdx)
const
511 if (!this->soMax_.empty())
512 simulator.problem().setMaxOilSaturation(elemIdx, this->soMax_[elemIdx]);
514 if (simulator.problem().materialLawManager()->enableHysteresis()) {
515 auto matLawManager = simulator.problem().materialLawManager();
517 if (FluidSystem::phaseIsActive(oilPhaseIdx)
518 && FluidSystem::phaseIsActive(waterPhaseIdx)) {
524 if (!this->soMax_.empty()) {
525 somax = this->soMax_[elemIdx];
529 if (!this->swMax_.empty()) {
530 swmax = this->swMax_[elemIdx];
534 if (!this->swmin_.empty()) {
535 swmin = this->swmin_[elemIdx];
539 somax, swmax, swmin, elemIdx);
541 if (FluidSystem::phaseIsActive(oilPhaseIdx)
542 && FluidSystem::phaseIsActive(gasPhaseIdx)) {
548 if (!this->sgmax_.empty()) {
549 sgmax = this->sgmax_[elemIdx];
553 if (!this->shmax_.empty()) {
554 shmax = this->shmax_[elemIdx];
558 if (!this->somin_.empty()) {
559 somin = this->somin_[elemIdx];
563 sgmax, shmax, somin, elemIdx);
568 if (simulator_.vanguard().eclState().fieldProps().has_double(
"SWATINIT")) {
569 simulator.problem().materialLawManager()
570 ->applyRestartSwatInit(elemIdx, this->ppcw_[elemIdx]);
574 void updateFluidInPlace(
const ElementContext& elemCtx)
576 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
577 updateFluidInPlace_(elemCtx, dofIdx);
581 void updateFluidInPlace(
const unsigned globalDofIdx,
582 const IntensiveQuantities& intQuants,
585 this->updateFluidInPlace_(globalDofIdx, intQuants,
totVolume);
589 template <
typename T>
590 using RemoveCVR = std::remove_cv_t<std::remove_reference_t<T>>;
592 template <
typename,
class =
void>
593 struct HasGeoMech :
public std::false_type {};
595 template <
typename Problem>
597 Problem, std::
void_t<decltype(std::declval<Problem>().geoMechModel())>
598 > :
public std::true_type {};
600 bool isDefunctParallelWell(std::string
wname)
const override
602 if (simulator_.gridView().comm().size() == 1)
604 const auto& parallelWells = simulator_.vanguard().parallelWells();
605 std::pair<std::string, bool> value {
wname,
true};
606 auto candidate = std::lower_bound(parallelWells.begin(), parallelWells.end(), value);
610 void updateFluidInPlace_(
const ElementContext& elemCtx,
const unsigned dofIdx)
612 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
613 const unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
614 const auto totVolume = elemCtx.simulator().model().dofTotalVolume(globalDofIdx);
616 this->updateFluidInPlace_(globalDofIdx, intQuants,
totVolume);
619 void updateFluidInPlace_(
const unsigned globalDofIdx,
620 const IntensiveQuantities& intQuants,
625 this->updateTotalVolumesAndPressures_(globalDofIdx, intQuants,
totVolume);
627 if (this->computeFip_) {
628 this->updatePhaseInplaceVolumes_(globalDofIdx, intQuants,
totVolume);
632 void createLocalRegion_(std::vector<int>& region)
634 std::size_t elemIdx = 0;
635 for (
const auto&
elem :
elements(simulator_.gridView())) {
636 if (
elem.partitionType() != Dune::InteriorEntity) {
644 template <
typename Flu
idState>
645 void aggregateAverageDensityContributions_(
const FluidState& fs,
646 const unsigned int globalDofIdx,
649 auto pvCellValue = RegionPhasePoreVolAverage::CellValue{};
652 for (
auto phaseIdx = 0*FluidSystem::numPhases;
655 if (! FluidSystem::phaseIsActive(
phaseIdx)) {
662 this->regionAvgDensity_
663 ->addCell(globalDofIdx,
664 RegionPhasePoreVolAverage::Phase {
phaseIdx },
684 data::InterRegFlowMap::FlowRates
685 getComponentSurfaceRates(
const ElementContext& elemCtx,
686 const Scalar faceArea,
688 const std::size_t
timeIdx)
const
690 using Component = data::InterRegFlowMap::Component;
692 auto rates = data::InterRegFlowMap::FlowRates {};
698 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
699 const auto&
up = elemCtx
702 const auto pvtReg =
up.pvtRegionIndex();
705 (
up.fluidState(), oilPhaseIdx,
pvtReg));
709 rates[Component::Oil] +=
qO;
711 if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
713 BlackOil::getRs_<FluidSystem, FluidState, Scalar>
716 rates[Component::Gas] +=
qO * Rs;
717 rates[Component::Disgas] +=
qO * Rs;
721 if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
722 const auto&
up = elemCtx
725 const auto pvtReg =
up.pvtRegionIndex();
728 (
up.fluidState(), gasPhaseIdx,
pvtReg));
732 rates[Component::Gas] +=
qG;
734 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
736 BlackOil::getRv_<FluidSystem, FluidState, Scalar>
739 rates[Component::Oil] +=
qG * Rv;
740 rates[Component::Vapoil] +=
qG * Rv;
744 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
745 const auto&
up = elemCtx
748 const auto pvtReg =
up.pvtRegionIndex();
751 (
up.fluidState(), waterPhaseIdx,
pvtReg));
753 rates[Component::Water] +=
760 template <
typename Flu
idState>
761 Scalar hydroCarbonFraction(
const FluidState& fs)
const
763 if (this->eclState_.runspec().co2Storage()) {
771 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
775 if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
782 void updateTotalVolumesAndPressures_(
const unsigned globalDofIdx,
783 const IntensiveQuantities& intQuants,
786 const auto& fs = intQuants.fluidState();
788 const double pv =
totVolume * intQuants.porosity().value();
789 const auto hydrocarbon = this->hydroCarbonFraction(fs);
791 if (! this->hydrocarbonPoreVolume_.empty()) {
792 this->fipC_.assignPoreVolume(globalDofIdx,
793 totVolume * intQuants.referencePorosity());
795 this->dynamicPoreVolume_[globalDofIdx] = pv;
796 this->hydrocarbonPoreVolume_[globalDofIdx] = pv *
hydrocarbon;
799 if (!this->pressureTimesHydrocarbonVolume_.empty() &&
800 !
this->pressureTimesPoreVolume_.empty())
802 assert(this->hydrocarbonPoreVolume_.size() ==
this->pressureTimesHydrocarbonVolume_.size());
803 assert(this->fipC_.get(Inplace::Phase::PoreVolume).size() ==
this->pressureTimesPoreVolume_.size());
805 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
806 this->pressureTimesPoreVolume_[globalDofIdx] =
807 getValue(fs.pressure(oilPhaseIdx)) * pv;
809 this->pressureTimesHydrocarbonVolume_[globalDofIdx] =
810 this->pressureTimesPoreVolume_[globalDofIdx] *
hydrocarbon;
812 else if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
813 this->pressureTimesPoreVolume_[globalDofIdx] =
814 getValue(fs.pressure(gasPhaseIdx)) * pv;
816 this->pressureTimesHydrocarbonVolume_[globalDofIdx] =
817 this->pressureTimesPoreVolume_[globalDofIdx] *
hydrocarbon;
819 else if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
820 this->pressureTimesPoreVolume_[globalDofIdx] =
821 getValue(fs.pressure(waterPhaseIdx)) * pv;
826 void updatePhaseInplaceVolumes_(
const unsigned globalDofIdx,
827 const IntensiveQuantities& intQuants,
830 std::array<Scalar, FluidSystem::numPhases> fip {};
831 std::array<Scalar, FluidSystem::numPhases>
fipr{};
833 const auto& fs = intQuants.fluidState();
834 const auto pv =
totVolume * intQuants.porosity().value();
837 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
848 this->fipC_.assignVolumesSurface(globalDofIdx, fip);
849 this->fipC_.assignVolumesReservoir(globalDofIdx,
850 fs.saltConcentration().value(),
853 if (FluidSystem::phaseIsActive(oilPhaseIdx) &&
854 FluidSystem::phaseIsActive(gasPhaseIdx))
856 this->updateOilGasDistribution(globalDofIdx, fs, fip);
859 if (FluidSystem::phaseIsActive(waterPhaseIdx) &&
860 FluidSystem::phaseIsActive(gasPhaseIdx))
862 this->updateGasWaterDistribution(globalDofIdx, fs, fip);
865 if (FluidSystem::phaseIsActive(gasPhaseIdx) &&
866 this->fipC_.hasCo2InGas())
868 this->updateCO2InGas(globalDofIdx, pv, intQuants);
871 if (this->fipC_.hasCo2InWater() &&
872 (FluidSystem::phaseIsActive(waterPhaseIdx) ||
873 FluidSystem::phaseIsActive(oilPhaseIdx)))
875 this->updateCO2InWater(globalDofIdx, pv, fs);
878 if constexpr(enableMICP) {
880 if (this->fipC_.hasMicrobialMass()) {
881 this->updateMicrobialMass(globalDofIdx, intQuants,
surfVolWat);
883 if (this->fipC_.hasOxygenMass()) {
884 this->updateOxygenMass(globalDofIdx, intQuants,
surfVolWat);
886 if (this->fipC_.hasUreaMass()) {
887 this->updateUreaMass(globalDofIdx, intQuants,
surfVolWat);
889 if (this->fipC_.hasBiofilmMass()) {
890 this->updateBiofilmMass(globalDofIdx, intQuants,
totVolume);
892 if (this->fipC_.hasCalciteMass()) {
893 this->updateCalciteMass(globalDofIdx, intQuants,
totVolume);
898 template <
typename Flu
idState,
typename FIPArray>
899 void updateOilGasDistribution(
const unsigned globalDofIdx,
900 const FluidState& fs,
910 template <
typename Flu
idState,
typename FIPArray>
911 void updateGasWaterDistribution(
const unsigned globalDofIdx,
912 const FluidState& fs,
922 template <
typename IntensiveQuantities>
923 void updateCO2InGas(
const unsigned globalDofIdx,
925 const IntensiveQuantities& intQuants)
928 ->oilWaterScaledEpsInfoDrainage(globalDofIdx);
930 const auto& fs = intQuants.fluidState();
932 if (this->simulator_.problem().materialLawManager()->enableHysteresis()) {
933 const auto&
matParams = simulator_.problem().materialLawParams(globalDofIdx);
934 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
938 if (this->fipC_.has(Inplace::Phase::CO2MassInGasPhaseMaximumTrapped) ||
939 this->fipC_.has(Inplace::Phase::CO2MassInGasPhaseMaximumUnTrapped))
941 if (this->simulator_.problem().materialLawManager()->enableHysteresis()) {
942 const auto&
matParams = simulator_.problem().materialLawParams(globalDofIdx);
948 const Scalar sg =
getValue(fs.saturation(gasPhaseIdx));
950 if (this->fipC_.has(Inplace::Phase::CO2MassInGasPhaseEffectiveTrapped) ||
951 this->fipC_.has(Inplace::Phase::CO2MassInGasPhaseEffectiveUnTrapped))
953 if (this->simulator_.problem().materialLawManager()->enableHysteresis()) {
954 const auto&
matParams = simulator_.problem().materialLawParams(globalDofIdx);
955 const double krg =
getValue(intQuants.relativePermeability(gasPhaseIdx));
960 const typename FIPContainer<FluidSystem>::Co2InGasInput
v{
965 FluidSystem::phaseIsActive(waterPhaseIdx)
966 ? FluidSystem::convertRvwToXgW(
getValue(fs.Rvw()), fs.pvtRegionIndex())
968 FluidSystem::molarMass(gasCompIdx, fs.pvtRegionIndex()),
973 this->fipC_.assignCo2InGas(globalDofIdx,
v);
976 template <
typename Flu
idState>
977 void updateCO2InWater(
const unsigned globalDofIdx,
979 const FluidState& fs)
981 const auto co2InWater = FluidSystem::phaseIsActive(oilPhaseIdx)
982 ? this->co2InWaterFromOil(fs, pv)
983 :
this->co2InWaterFromWater(fs, pv);
985 const Scalar mM = FluidSystem::molarMass(gasCompIdx, fs.pvtRegionIndex());
987 this->fipC_.assignCo2InWater(globalDofIdx,
co2InWater, mM);
990 template <
typename Flu
idState>
991 Scalar co2InWaterFromWater(
const FluidState& fs,
const double pv)
const
993 const double rhow =
getValue(fs.density(waterPhaseIdx));
994 const double sw =
getValue(fs.saturation(waterPhaseIdx));
995 const double xwG = FluidSystem::convertRswToXwG(
getValue(fs.Rsw()), fs.pvtRegionIndex());
997 const Scalar mM = FluidSystem::molarMass(gasCompIdx, fs.pvtRegionIndex());
1002 template <
typename Flu
idState>
1003 Scalar co2InWaterFromOil(
const FluidState& fs,
const double pv)
const
1006 const double so =
getValue(fs.saturation(oilPhaseIdx));
1007 const double xoG = FluidSystem::convertRsToXoG(
getValue(fs.Rs()), fs.pvtRegionIndex());
1009 const Scalar mM = FluidSystem::molarMass(gasCompIdx, fs.pvtRegionIndex());
1014 template <
typename IntensiveQuantities>
1015 void updateMicrobialMass(
const unsigned globalDofIdx,
1016 const IntensiveQuantities& intQuants,
1019 const Scalar
mass =
surfVolWat * intQuants.microbialConcentration().value();
1021 this->fipC_.assignMicrobialMass(globalDofIdx,
mass);
1024 template <
typename IntensiveQuantities>
1025 void updateOxygenMass(
const unsigned globalDofIdx,
1026 const IntensiveQuantities& intQuants,
1029 const Scalar
mass =
surfVolWat * intQuants.oxygenConcentration().value();
1031 this->fipC_.assignOxygenMass(globalDofIdx,
mass);
1034 template <
typename IntensiveQuantities>
1035 void updateUreaMass(
const unsigned globalDofIdx,
1036 const IntensiveQuantities& intQuants,
1039 const Scalar
mass =
surfVolWat * intQuants.ureaConcentration().value();
1041 this->fipC_.assignUreaMass(globalDofIdx,
mass);
1044 template <
typename IntensiveQuantities>
1045 void updateBiofilmMass(
const unsigned globalDofIdx,
1046 const IntensiveQuantities& intQuants,
1049 const Scalar
mass =
totVolume * intQuants.biofilmMass().value();
1051 this->fipC_.assignBiofilmMass(globalDofIdx,
mass);
1054 template <
typename IntensiveQuantities>
1055 void updateCalciteMass(
const unsigned globalDofIdx,
1056 const IntensiveQuantities& intQuants,
1059 const Scalar
mass =
totVolume * intQuants.calciteMass().value();
1061 this->fipC_.assignCalciteMass(globalDofIdx,
mass);
1065 void setupElementExtractors_()
1067 using Entry =
typename Extractor::Entry;
1068 using Context =
typename Extractor::Context;
1069 using ScalarEntry =
typename Extractor::ScalarEntry;
1070 using PhaseEntry =
typename Extractor::PhaseEntry;
1072 const bool hasResidual = simulator_.model().linearizer().residual().size() > 0;
1073 const auto&
hysteresisConfig = simulator_.problem().materialLawManager()->hysteresisConfig();
1076 Entry{PhaseEntry{&this->saturation_,
1077 [](
const unsigned phase,
const Context&
ectx)
1081 Entry{PhaseEntry{&this->invB_,
1082 [](
const unsigned phase,
const Context&
ectx)
1086 Entry{PhaseEntry{&this->density_,
1087 [](
const unsigned phase,
const Context&
ectx)
1091 Entry{PhaseEntry{&this->relativePermeability_,
1092 [](
const unsigned phase,
const Context&
ectx)
1093 {
return getValue(
ectx.intQuants.relativePermeability(phase)); }
1096 Entry{PhaseEntry{&this->viscosity_,
1099 if (this->extboC_.allocated() &&
phaseIdx == oilPhaseIdx) {
1102 else if (this->extboC_.allocated() &&
phaseIdx == gasPhaseIdx) {
1111 Entry{PhaseEntry{&this->residual_,
1112 [&
modelResid = this->simulator_.model().linearizer().residual()]
1115 const unsigned sIdx = FluidSystem::solventComponentIndex(
phaseIdx);
1122 Entry{ScalarEntry{&this->rockCompPorvMultiplier_,
1123 [&problem = this->simulator_.problem()](
const Context&
ectx)
1125 return problem.template
1131 Entry{ScalarEntry{&this->rockCompTransMultiplier_,
1132 [&problem = this->simulator_.problem()](
const Context&
ectx)
1139 Entry{ScalarEntry{&this->minimumOilPressure_,
1140 [&problem = this->simulator_.problem()](
const Context&
ectx)
1142 return std::min(
getValue(
ectx.fs.pressure(oilPhaseIdx)),
1143 problem.minOilPressure(
ectx.globalDofIdx));
1147 Entry{ScalarEntry{&this->bubblePointPressure_,
1149 &vanguard = this->simulator_.vanguard()](
const Context&
ectx)
1153 FluidSystem::bubblePointPressure(
ectx.fs,
1154 ectx.intQuants.pvtRegionIndex())
1164 Entry{ScalarEntry{&this->dewPointPressure_,
1166 &vanguard = this->simulator_.vanguard()](
const Context&
ectx)
1170 FluidSystem::dewPointPressure(
ectx.fs,
1171 ectx.intQuants.pvtRegionIndex())
1181 Entry{ScalarEntry{&this->overburdenPressure_,
1182 [&problem = simulator_.problem()](
const Context&
ectx)
1183 {
return problem.overburdenPressure(
ectx.globalDofIdx); }
1186 Entry{ScalarEntry{&this->temperature_,
1187 [](
const Context&
ectx)
1191 Entry{ScalarEntry{&this->sSol_,
1192 [](
const Context&
ectx)
1193 {
return getValue(
ectx.intQuants.solventSaturation()); }
1196 Entry{ScalarEntry{&this->rswSol_,
1197 [](
const Context&
ectx)
1201 Entry{ScalarEntry{&this->cPolymer_,
1202 [](
const Context&
ectx)
1203 {
return getValue(
ectx.intQuants.polymerConcentration()); }
1206 Entry{ScalarEntry{&this->cFoam_,
1207 [](
const Context&
ectx)
1208 {
return getValue(
ectx.intQuants.foamConcentration()); }
1211 Entry{ScalarEntry{&this->cSalt_,
1212 [](
const Context&
ectx)
1216 Entry{ScalarEntry{&this->pSalt_,
1217 [](
const Context&
ectx)
1221 Entry{ScalarEntry{&this->permFact_,
1222 [](
const Context&
ectx)
1226 Entry{ScalarEntry{&this->rPorV_,
1227 [&model = this->simulator_.model()](
const Context&
ectx)
1229 const auto totVolume = model.dofTotalVolume(
ectx.globalDofIdx);
1234 Entry{ScalarEntry{&this->rs_,
1235 [](
const Context&
ectx)
1239 Entry{ScalarEntry{&this->rv_,
1240 [](
const Context&
ectx)
1244 Entry{ScalarEntry{&this->rsw_,
1245 [](
const Context&
ectx)
1249 Entry{ScalarEntry{&this->rvw_,
1250 [](
const Context&
ectx)
1254 Entry{ScalarEntry{&this->ppcw_,
1255 [&
matLawManager = *this->simulator_.problem().materialLawManager()]
1256 (
const Context&
ectx)
1263 Entry{ScalarEntry{&this->drsdtcon_,
1264 [&problem = this->simulator_.problem()](
const Context&
ectx)
1266 return problem.drsdtcon(
ectx.globalDofIdx,
1271 Entry{ScalarEntry{&this->pcgw_,
1272 [](
const Context&
ectx)
1279 Entry{ScalarEntry{&this->pcow_,
1280 [](
const Context&
ectx)
1287 Entry{ScalarEntry{&this->pcog_,
1288 [](
const Context&
ectx)
1295 Entry{ScalarEntry{&this->fluidPressure_,
1296 [](
const Context&
ectx)
1298 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1302 else if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
1313 Entry{ScalarEntry{&this->gasDissolutionFactor_,
1314 [&problem = this->simulator_.problem()](
const Context&
ectx)
1316 const Scalar
SoMax = problem.maxOilSaturation(
ectx.globalDofIdx);
1317 return FluidSystem::template
1325 Entry{ScalarEntry{&this->oilVaporizationFactor_,
1326 [&problem = this->simulator_.problem()](
const Context&
ectx)
1328 const Scalar
SoMax = problem.maxOilSaturation(
ectx.globalDofIdx);
1329 return FluidSystem::template
1337 Entry{ScalarEntry{&this->gasDissolutionFactorInWater_,
1338 [&problem = this->simulator_.problem()](
const Context&
ectx)
1340 const Scalar
SwMax = problem.maxWaterSaturation(
ectx.globalDofIdx);
1341 return FluidSystem::template
1349 Entry{ScalarEntry{&this->waterVaporizationFactor_,
1350 [](
const Context&
ectx)
1352 return FluidSystem::template
1359 Entry{ScalarEntry{&this->gasFormationVolumeFactor_,
1360 [](
const Context&
ectx)
1362 return 1.0 / FluidSystem::template
1369 Entry{ScalarEntry{&this->saturatedOilFormationVolumeFactor_,
1370 [](
const Context&
ectx)
1372 return 1.0 / FluidSystem::template
1379 Entry{ScalarEntry{&this->oilSaturationPressure_,
1380 [](
const Context&
ectx)
1382 return FluidSystem::template
1389 Entry{ScalarEntry{&this->soMax_,
1390 [&problem = this->simulator_.problem()](
const Context&
ectx)
1392 return std::max(
getValue(
ectx.fs.saturation(oilPhaseIdx)),
1393 problem.maxOilSaturation(
ectx.globalDofIdx));
1398 Entry{ScalarEntry{&this->swMax_,
1399 [&problem = this->simulator_.problem()](
const Context&
ectx)
1401 return std::max(
getValue(
ectx.fs.saturation(waterPhaseIdx)),
1402 problem.maxWaterSaturation(
ectx.globalDofIdx));
1407 Entry{ScalarEntry{&this->soMax_,
1408 [](
const Context&
ectx)
1409 {
return ectx.hParams.somax; }
1413 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1414 FluidSystem::phaseIsActive(waterPhaseIdx)
1416 Entry{ScalarEntry{&this->swMax_,
1417 [](
const Context&
ectx)
1418 {
return ectx.hParams.swmax; }
1422 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1423 FluidSystem::phaseIsActive(waterPhaseIdx)
1425 Entry{ScalarEntry{&this->swmin_,
1426 [](
const Context&
ectx)
1427 {
return ectx.hParams.swmin; }
1431 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1432 FluidSystem::phaseIsActive(waterPhaseIdx)
1434 Entry{ScalarEntry{&this->sgmax_,
1435 [](
const Context&
ectx)
1436 {
return ectx.hParams.sgmax; }
1440 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1441 FluidSystem::phaseIsActive(gasPhaseIdx)
1443 Entry{ScalarEntry{&this->shmax_,
1444 [](
const Context&
ectx)
1445 {
return ectx.hParams.shmax; }
1449 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1450 FluidSystem::phaseIsActive(gasPhaseIdx)
1452 Entry{ScalarEntry{&this->somin_,
1453 [](
const Context&
ectx)
1454 {
return ectx.hParams.somin; }
1458 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1459 FluidSystem::phaseIsActive(gasPhaseIdx)
1461 Entry{[&model = this->simulator_.model(),
this](
const Context&
ectx)
1465 const auto porv =
ectx.intQuants.referencePorosity()
1466 * model.dofTotalVolume(
ectx.globalDofIdx);
1468 this->aggregateAverageDensityContributions_(
ectx.fs,
ectx.globalDofIdx,
1469 static_cast<double>(porv));
1470 }, this->regionAvgDensity_.has_value()
1472 Entry{[&
extboC = this->extboC_](
const Context&
ectx)
1475 ectx.intQuants.xVolume().value(),
1476 ectx.intQuants.yVolume().value());
1478 ectx.intQuants.zFraction().value());
1487 (1.0 -
ectx.intQuants.yVolume().value()) +
1491 (1.0 -
ectx.intQuants.xVolume().value());
1494 ectx.intQuants.yVolume().value() +
1498 ectx.intQuants.xVolume().value();
1499 const Scalar
rhoO = FluidSystem::referenceDensity(gasPhaseIdx,
ectx.pvtRegionIdx);
1500 const Scalar
rhoG = FluidSystem::referenceDensity(gasPhaseIdx,
ectx.pvtRegionIdx);
1501 const Scalar
rhoCO2 =
ectx.intQuants.zRefDensity();
1503 extboC.assignMassFractions(
ectx.globalDofIdx,
1507 }, this->extboC_.allocated()
1509 Entry{[&
micpC = this->micpC_](
const Context&
ectx)
1512 ectx.intQuants.microbialConcentration().value(),
1513 ectx.intQuants.oxygenConcentration().value(),
1514 ectx.intQuants.ureaConcentration().value(),
1515 ectx.intQuants.biofilmConcentration().value(),
1516 ectx.intQuants.calciteConcentration().value());
1517 }, this->micpC_.allocated()
1519 Entry{[&
rftC = this->rftC_,
1520 &vanguard = this->simulator_.vanguard()](
const Context&
ectx)
1524 [&fs =
ectx.fs]() {
return getValue(fs.pressure(oilPhaseIdx)); },
1525 [&fs =
ectx.fs]() {
return getValue(fs.saturation(waterPhaseIdx)); },
1526 [&fs =
ectx.fs]() {
return getValue(fs.saturation(gasPhaseIdx)); });
1529 Entry{[&
tC = this->tracerC_,
1530 &
tM = this->simulator_.problem().tracerModel()](
const Context&
ectx)
1532 tC.assignFreeConcentrations(
ectx.globalDofIdx,
1535 tC.assignSolConcentrations(
ectx.globalDofIdx,
1540 Entry{[&
flowsInf = this->simulator_.problem().model().linearizer().getFlowsInfo(),
1541 &
flowsC = this->flowsC_](
const Context&
ectx)
1543 constexpr auto gas_idx = Indices::gasEnabled ?
1544 conti0EqIdx + Indices::canonicalToActiveComponentIndex(gasCompIdx) : -1;
1545 constexpr auto oil_idx = Indices::oilEnabled ?
1546 conti0EqIdx + Indices::canonicalToActiveComponentIndex(oilCompIdx) : -1;
1547 constexpr auto water_idx = Indices::waterEnabled ?
1548 conti0EqIdx + Indices::canonicalToActiveComponentIndex(waterCompIdx) : -1;
1558 }, !this->simulator_.problem().model().linearizer().getFlowsInfo().empty()
1560 Entry{[&
floresInf = this->simulator_.problem().model().linearizer().getFloresInfo(),
1561 &
flowsC = this->flowsC_](
const Context&
ectx)
1563 constexpr auto gas_idx = Indices::gasEnabled ?
1564 conti0EqIdx + Indices::canonicalToActiveComponentIndex(gasCompIdx) : -1;
1565 constexpr auto oil_idx = Indices::oilEnabled ?
1566 conti0EqIdx + Indices::canonicalToActiveComponentIndex(oilCompIdx) : -1;
1567 constexpr auto water_idx = Indices::waterEnabled ?
1568 conti0EqIdx + Indices::canonicalToActiveComponentIndex(waterCompIdx) : -1;
1578 }, !this->simulator_.problem().model().linearizer().getFloresInfo().empty()
1585 Entry{ScalarEntry{&this->rv_,
1586 [&problem = this->simulator_.problem()](
const Context&
ectx)
1587 {
return problem.initialFluidState(
ectx.globalDofIdx).Rv(); }
1589 simulator_.episodeIndex() < 0 &&
1590 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1591 FluidSystem::phaseIsActive(gasPhaseIdx)
1593 Entry{ScalarEntry{&this->rs_,
1594 [&problem = this->simulator_.problem()](
const Context&
ectx)
1595 {
return problem.initialFluidState(
ectx.globalDofIdx).Rs(); }
1597 simulator_.episodeIndex() < 0 &&
1598 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1599 FluidSystem::phaseIsActive(gasPhaseIdx)
1601 Entry{ScalarEntry{&this->rsw_,
1602 [&problem = this->simulator_.problem()](
const Context&
ectx)
1603 {
return problem.initialFluidState(
ectx.globalDofIdx).Rsw(); }
1605 simulator_.episodeIndex() < 0 &&
1606 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1607 FluidSystem::phaseIsActive(gasPhaseIdx)
1609 Entry{ScalarEntry{&this->rvw_,
1610 [&problem = this->simulator_.problem()](
const Context&
ectx)
1611 {
return problem.initialFluidState(
ectx.globalDofIdx).Rvw(); }
1613 simulator_.episodeIndex() < 0 &&
1614 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1615 FluidSystem::phaseIsActive(gasPhaseIdx)
1618 Entry{PhaseEntry{&this->density_,
1619 [&problem = this->simulator_.problem()](
const unsigned phase,
1620 const Context&
ectx)
1622 const auto&
fsInitial = problem.initialFluidState(
ectx.globalDofIdx);
1625 ectx.intQuants.pvtRegionIndex());
1628 simulator_.episodeIndex() < 0 &&
1629 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1630 FluidSystem::phaseIsActive(gasPhaseIdx)
1632 Entry{PhaseEntry{&this->invB_,
1633 [&problem = this->simulator_.problem()](
const unsigned phase,
1634 const Context&
ectx)
1636 const auto&
fsInitial = problem.initialFluidState(
ectx.globalDofIdx);
1637 return FluidSystem::inverseFormationVolumeFactor(
fsInitial,
1639 ectx.intQuants.pvtRegionIndex());
1642 simulator_.episodeIndex() < 0 &&
1643 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1644 FluidSystem::phaseIsActive(gasPhaseIdx)
1646 Entry{PhaseEntry{&this->viscosity_,
1647 [&problem = this->simulator_.problem()](
const unsigned phase,
1648 const Context&
ectx)
1650 const auto&
fsInitial = problem.initialFluidState(
ectx.globalDofIdx);
1651 return FluidSystem::viscosity(
fsInitial,
1653 ectx.intQuants.pvtRegionIndex());
1656 simulator_.episodeIndex() < 0 &&
1657 FluidSystem::phaseIsActive(oilPhaseIdx) &&
1658 FluidSystem::phaseIsActive(gasPhaseIdx)
1666 if (this->mech_.allocated()) {
1667 this->extractors_.emplace(
1668 [&
mech = this->mech_,
1669 &model = simulator_.problem().geoMechModel()](
const Context&
ectx)
1671 mech.assignDelStress(
ectx.globalDofIdx,
1672 model.delstress(
ectx.globalDofIdx));
1674 mech.assignDisplacement(
ectx.globalDofIdx,
1675 model.disp(
ectx.globalDofIdx,
true));
1678 mech.assignFracStress(
ectx.globalDofIdx,
1679 model.fractureStress(
ectx.globalDofIdx));
1681 mech.assignLinStress(
ectx.globalDofIdx,
1682 model.linstress(
ectx.globalDofIdx));
1684 mech.assignPotentialForces(
ectx.globalDofIdx,
1685 model.mechPotentialForce(
ectx.globalDofIdx),
1686 model.mechPotentialPressForce(
ectx.globalDofIdx),
1687 model.mechPotentialTempForce(
ectx.globalDofIdx));
1689 mech.assignStrain(
ectx.globalDofIdx,
1690 model.strain(
ectx.globalDofIdx,
true));
1693 mech.assignStress(
ectx.globalDofIdx,
1694 model.stress(
ectx.globalDofIdx,
true));
1702 void setupBlockExtractors_(
const bool isSubStep,
1703 const int reportStepNum)
1706 using Context =
typename BlockExtractor::Context;
1707 using PhaseEntry =
typename BlockExtractor::PhaseEntry;
1708 using ScalarEntry =
typename BlockExtractor::ScalarEntry;
1710 using namespace std::string_view_literals;
1713 Entry{ScalarEntry{std::vector{
"BPR"sv,
"BPRESSUR"sv},
1714 [](
const Context&
ectx)
1716 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1719 else if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
1731 Entry{PhaseEntry{std::array{
1732 std::array{
"BWSAT"sv,
"BOSAT"sv,
"BGSAT"sv},
1733 std::array{
"BSWAT"sv,
"BSOIL"sv,
"BSGAS"sv}
1741 Entry{ScalarEntry{
"BNSAT",
1742 [](
const Context&
ectx)
1744 return ectx.intQuants.solventSaturation().value();
1748 Entry{ScalarEntry{std::vector{
"BTCNFHEA"sv,
"BTEMP"sv},
1749 [](
const Context&
ectx)
1751 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1754 else if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
1763 Entry{PhaseEntry{std::array{
1764 std::array{
"BWKR"sv,
"BOKR"sv,
"BGKR"sv},
1765 std::array{
"BKRW"sv,
"BKRO"sv,
"BKRG"sv}
1773 Entry{ScalarEntry{
"BKROG",
1774 [&problem = this->simulator_.problem()](
const Context&
ectx)
1777 problem.materialLawParams(
ectx.elemCtx,
1780 return getValue(MaterialLaw::template
1786 Entry{ScalarEntry{
"BKROW",
1787 [&problem = this->simulator_.problem()](
const Context&
ectx)
1792 return getValue(MaterialLaw::template
1798 Entry{ScalarEntry{
"BWPC",
1799 [](
const Context&
ectx)
1806 Entry{ScalarEntry{
"BGPC",
1807 [](
const Context&
ectx)
1814 Entry{ScalarEntry{
"BWPR",
1815 [](
const Context&
ectx)
1821 Entry{ScalarEntry{
"BGPR",
1822 [](
const Context&
ectx)
1828 Entry{PhaseEntry{std::array{
1829 std::array{
"BVWAT"sv,
"BVOIL"sv,
"BVGAS"sv},
1830 std::array{
"BWVIS"sv,
"BOVIS"sv,
"BGVIS"sv}
1838 Entry{PhaseEntry{std::array{
1839 std::array{
"BWDEN"sv,
"BODEN"sv,
"BGDEN"sv},
1840 std::array{
"BDENW"sv,
"BDENO"sv,
"BDENG"sv}
1848 Entry{ScalarEntry{
"BFLOWI",
1849 [&
flowsC = this->flowsC_](
const Context&
ectx)
1851 return flowsC.getFlow(
ectx.globalDofIdx, Dir::XPlus, waterCompIdx);
1855 Entry{ScalarEntry{
"BFLOWJ",
1856 [&
flowsC = this->flowsC_](
const Context&
ectx)
1858 return flowsC.getFlow(
ectx.globalDofIdx, Dir::YPlus, waterCompIdx);
1862 Entry{ScalarEntry{
"BFLOWK",
1863 [&
flowsC = this->flowsC_](
const Context&
ectx)
1865 return flowsC.getFlow(
ectx.globalDofIdx, Dir::ZPlus, waterCompIdx);
1869 Entry{ScalarEntry{
"BRPV",
1870 [&model = this->simulator_.model()](
const Context&
ectx)
1873 model.dofTotalVolume(
ectx.globalDofIdx);
1877 Entry{PhaseEntry{std::array{
"BWPV"sv,
"BOPV"sv,
"BGPV"sv},
1878 [&model = this->simulator_.model()](
const unsigned phaseIdx,
1879 const Context&
ectx)
1883 model.dofTotalVolume(
ectx.globalDofIdx);
1887 Entry{ScalarEntry{
"BRS",
1888 [](
const Context&
ectx)
1894 Entry{ScalarEntry{
"BRV",
1895 [](
const Context&
ectx)
1901 Entry{ScalarEntry{
"BOIP",
1902 [&model = this->simulator_.model()](
const Context&
ectx)
1909 model.dofTotalVolume(
ectx.globalDofIdx) *
1914 Entry{ScalarEntry{
"BGIP",
1915 [&model = this->simulator_.model()](
const Context&
ectx)
1920 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1932 model.dofTotalVolume(
ectx.globalDofIdx) *
1937 Entry{ScalarEntry{
"BWIP",
1938 [&model = this->simulator_.model()](
const Context&
ectx)
1942 model.dofTotalVolume(
ectx.globalDofIdx) *
1947 Entry{ScalarEntry{
"BOIPL",
1948 [&model = this->simulator_.model()](
const Context&
ectx)
1952 model.dofTotalVolume(
ectx.globalDofIdx) *
1957 Entry{ScalarEntry{
"BGIPL",
1958 [&model = this->simulator_.model()](
const Context&
ectx)
1961 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
1972 model.dofTotalVolume(
ectx.globalDofIdx) *
1977 Entry{ScalarEntry{
"BGIPG",
1978 [&model = this->simulator_.model()](
const Context&
ectx)
1982 model.dofTotalVolume(
ectx.globalDofIdx) *
1987 Entry{ScalarEntry{
"BOIPG",
1988 [&model = this->simulator_.model()](
const Context&
ectx)
1993 model.dofTotalVolume(
ectx.globalDofIdx) *
1998 Entry{PhaseEntry{std::array{
"BPPW"sv,
"BPPO"sv,
"BPPG"sv},
1999 [&
simConfig = this->eclState_.getSimulationConfig(),
2000 &
grav = this->simulator_.problem().gravity(),
2002 &problem = this->simulator_.problem(),
2005 auto phase = RegionPhasePoreVolAverage::Phase{};
2018 const auto region = RegionPhasePoreVolAverage::Region {
2019 ectx.elemCtx.primaryVars(
ectx.dofIdx, 0).pvtRegionIndex() + 1
2024 const auto press =
getValue(
ectx.fs.pressure(phase.ix));
2025 const auto dz = problem.dofCenterDepth(
ectx.globalDofIdx) - datum;
2026 return press - density*
dz*
grav[GridView::dimensionworld - 1];
2030 Entry{ScalarEntry{
"BMMIP",
2031 [&model = this->simulator_.model()](
const Context&
ectx)
2033 return getValue(
ectx.intQuants.microbialConcentration()) *
2035 model.dofTotalVolume(
ectx.globalDofIdx);
2039 Entry{ScalarEntry{
"BMOIP",
2040 [&model = this->simulator_.model()](
const Context&
ectx)
2042 return getValue(
ectx.intQuants.oxygenConcentration()) *
2044 model.dofTotalVolume(
ectx.globalDofIdx);
2048 Entry{ScalarEntry{
"BMUIP",
2049 [&model = this->simulator_.model()](
const Context&
ectx)
2053 model.dofTotalVolume(
ectx.globalDofIdx) * 1;
2057 Entry{ScalarEntry{
"BMBIP",
2058 [&model = this->simulator_.model()](
const Context&
ectx)
2060 return model.dofTotalVolume(
ectx.globalDofIdx) *
2065 Entry{ScalarEntry{
"BMCIP",
2066 [&model = this->simulator_.model()](
const Context&
ectx)
2068 return model.dofTotalVolume(
ectx.globalDofIdx) *
2073 Entry{ScalarEntry{
"BGMIP",
2074 [&model = this->simulator_.model()](
const Context&
ectx)
2079 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
2089 const Scalar
rhoG = FluidSystem::referenceDensity(gasPhaseIdx,
2090 ectx.intQuants.pvtRegionIndex());
2092 model.dofTotalVolume(
ectx.globalDofIdx) *
2098 Entry{ScalarEntry{
"BGMGP",
2099 [&model = this->simulator_.model()](
const Context&
ectx)
2101 const Scalar
rhoG = FluidSystem::referenceDensity(gasPhaseIdx,
2102 ectx.intQuants.pvtRegionIndex());
2105 model.dofTotalVolume(
ectx.globalDofIdx) *
2111 Entry{ScalarEntry{
"BGMDS",
2112 [&model = this->simulator_.model()](
const Context&
ectx)
2115 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
2125 const Scalar
rhoG = FluidSystem::referenceDensity(gasPhaseIdx,
2126 ectx.intQuants.pvtRegionIndex());
2128 model.dofTotalVolume(
ectx.globalDofIdx) *
2134 Entry{ScalarEntry{
"BGMST",
2135 [&model = this->simulator_.model(),
2136 &problem = this->simulator_.problem()](
const Context&
ectx)
2139 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2140 const Scalar sg =
getValue(
ectx.fs.saturation(gasPhaseIdx));
2142 if (problem.materialLawManager()->enableHysteresis()) {
2143 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2144 const Scalar
krg =
getValue(
ectx.intQuants.relativePermeability(gasPhaseIdx));
2145 strandedGas = MaterialLaw::strandedGasSaturation(
matParams, sg,
krg);
2147 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2148 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2150 return (1.0 - xgW) *
2151 model.dofTotalVolume(
ectx.globalDofIdx) *
2154 std::min(strandedGas, sg);
2158 Entry{ScalarEntry{
"BGMUS",
2159 [&model = this->simulator_.model(),
2160 &problem = this->simulator_.problem()](
const Context&
ectx)
2163 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2164 const Scalar sg =
getValue(
ectx.fs.saturation(gasPhaseIdx));
2166 if (problem.materialLawManager()->enableHysteresis()) {
2167 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2168 const Scalar
krg =
getValue(
ectx.intQuants.relativePermeability(gasPhaseIdx));
2169 strandedGas = MaterialLaw::strandedGasSaturation(
matParams, sg,
krg);
2171 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2172 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2174 return (1.0 - xgW) *
2175 model.dofTotalVolume(
ectx.globalDofIdx) *
2178 std::max(Scalar{0.0}, sg - strandedGas);
2182 Entry{ScalarEntry{
"BGMTR",
2183 [&model = this->simulator_.model(),
2184 &problem = this->simulator_.problem()](
const Context&
ectx)
2187 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2189 if (problem.materialLawManager()->enableHysteresis()) {
2190 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2191 trappedGas = MaterialLaw::trappedGasSaturation(
matParams,
true);
2193 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2194 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2196 return (1.0 - xgW) *
2197 model.dofTotalVolume(
ectx.globalDofIdx) *
2200 std::min(trappedGas,
getValue(
ectx.fs.saturation(gasPhaseIdx)));
2204 Entry{ScalarEntry{
"BGMMO",
2205 [&model = this->simulator_.model(),
2206 &problem = this->simulator_.problem()](
const Context&
ectx)
2209 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2211 if (problem.materialLawManager()->enableHysteresis()) {
2212 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2213 trappedGas = MaterialLaw::trappedGasSaturation(
matParams,
true);
2215 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2216 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2218 return (1.0 - xgW) *
2219 model.dofTotalVolume(
ectx.globalDofIdx) *
2222 std::max(Scalar{0.0},
getValue(
ectx.fs.saturation(gasPhaseIdx)) - trappedGas);
2226 Entry{ScalarEntry{
"BGKTR",
2227 [&model = this->simulator_.model(),
2228 &problem = this->simulator_.problem()](
const Context&
ectx)
2231 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2232 const Scalar sg =
getValue(
ectx.fs.saturation(gasPhaseIdx));
2234 if (problem.materialLawManager()->enableHysteresis()) {
2235 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2236 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
2242 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2243 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2245 return (1.0 - xgW) *
2246 model.dofTotalVolume(
ectx.globalDofIdx) *
2254 Entry{ScalarEntry{
"BGKMO",
2255 [&model = this->simulator_.model(),
2256 &problem = this->simulator_.problem()](
const Context&
ectx)
2259 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2260 const Scalar sg =
getValue(
ectx.fs.saturation(gasPhaseIdx));
2262 if (problem.materialLawManager()->enableHysteresis()) {
2263 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2264 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
2270 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2271 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2273 return (1.0 - xgW) *
2274 model.dofTotalVolume(
ectx.globalDofIdx) *
2282 Entry{ScalarEntry{
"BGCDI",
2283 [&model = this->simulator_.model(),
2284 &problem = this->simulator_.problem()](
const Context&
ectx)
2287 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2289 if (problem.materialLawManager()->enableHysteresis()) {
2290 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2291 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
2293 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2294 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2296 return (1.0 - xgW) *
2297 model.dofTotalVolume(
ectx.globalDofIdx) *
2300 std::min(sgcr,
getValue(
ectx.fs.saturation(gasPhaseIdx))) /
2301 FluidSystem::molarMass(gasCompIdx,
ectx.intQuants.pvtRegionIndex());
2305 Entry{ScalarEntry{
"BGCDM",
2306 [&model = this->simulator_.model(),
2307 &problem = this->simulator_.problem()](
const Context&
ectx)
2310 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2312 if (problem.materialLawManager()->enableHysteresis()) {
2313 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2314 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
2316 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2317 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2319 return (1.0 - xgW) *
2320 model.dofTotalVolume(
ectx.globalDofIdx) *
2323 std::max(Scalar{0.0},
getValue(
ectx.fs.saturation(gasPhaseIdx)) - sgcr) /
2324 FluidSystem::molarMass(gasCompIdx,
ectx.intQuants.pvtRegionIndex());
2328 Entry{ScalarEntry{
"BGKDI",
2329 [&model = this->simulator_.model(),
2330 &problem = this->simulator_.problem()](
const Context&
ectx)
2333 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2334 const Scalar sg =
getValue(
ectx.fs.saturation(gasPhaseIdx));
2336 if (problem.materialLawManager()->enableHysteresis()) {
2337 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2338 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
2344 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2345 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2347 return (1.0 - xgW) *
2348 model.dofTotalVolume(
ectx.globalDofIdx) *
2352 FluidSystem::molarMass(gasCompIdx,
ectx.intQuants.pvtRegionIndex());
2357 Entry{ScalarEntry{
"BGKDM",
2358 [&model = this->simulator_.model(),
2359 &problem = this->simulator_.problem()](
const Context&
ectx)
2362 ->oilWaterScaledEpsInfoDrainage(
ectx.dofIdx);
2363 const Scalar sg =
getValue(
ectx.fs.saturation(gasPhaseIdx));
2365 if (problem.materialLawManager()->enableHysteresis()) {
2366 const auto&
matParams = problem.materialLawParams(
ectx.dofIdx);
2367 sgcr = MaterialLaw::trappedGasSaturation(
matParams,
false);
2373 const Scalar xgW = FluidSystem::phaseIsActive(waterPhaseIdx) ?
2374 FluidSystem::convertRvwToXgW(
getValue(
ectx.fs.Rvw()),
ectx.intQuants.pvtRegionIndex())
2376 return (1.0 - xgW) *
2377 model.dofTotalVolume(
ectx.globalDofIdx) *
2381 FluidSystem::molarMass(gasCompIdx,
ectx.intQuants.pvtRegionIndex());
2386 Entry{ScalarEntry{
"BWCD",
2387 [&model = this->simulator_.model()](
const Context&
ectx)
2390 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
2400 const Scalar
rhoG = FluidSystem::referenceDensity(gasPhaseIdx,
2401 ectx.intQuants.pvtRegionIndex());
2403 model.dofTotalVolume(
ectx.globalDofIdx) *
2406 FluidSystem::molarMass(gasCompIdx,
ectx.intQuants.pvtRegionIndex());
2410 Entry{ScalarEntry{
"BWIPG",
2411 [&model = this->simulator_.model()](
const Context&
ectx)
2414 if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
2420 model.dofTotalVolume(
ectx.globalDofIdx) *
2425 Entry{ScalarEntry{
"BWIPL",
2426 [&model = this->simulator_.model()](
const Context&
ectx)
2430 model.dofTotalVolume(
ectx.globalDofIdx) *
2439 this->extraBlockData_.clear();
2442 const auto&
rpt = this->schedule_[reportStepNum - 1].rpt_config.get();
2443 if (
rpt.contains(
"WELLS") &&
rpt.at(
"WELLS") > 1) {
2444 this->setupExtraBlockData(reportStepNum,
2445 [&
c = this->collectOnIORank_](
const int idx)
2446 {
return c.isCartIdxOnThisRank(
idx); });
2457 const Simulator& simulator_;
2458 const CollectDataOnIORankType& collectOnIORank_;
2459 std::vector<typename Extractor::Entry> extractors_;