87 using typename FlowProblemType::Scalar;
88 using typename FlowProblemType::Simulator;
89 using typename FlowProblemType::GridView;
90 using typename FlowProblemType::FluidSystem;
91 using typename FlowProblemType::Vanguard;
92 using typename FlowProblemType::GlobalEqVector;
93 using typename FlowProblemType::EqVector;
94 using FlowProblemType::dim;
95 using FlowProblemType::dimWorld;
96 using FlowProblemType::numEq;
97 using FlowProblemType::numPhases;
98 using FlowProblemType::numComponents;
101 using FlowProblemType::enableConvectiveMixing;
102 using FlowProblemType::enableBrine;
103 using FlowProblemType::enableDiffusion;
104 using FlowProblemType::enableDispersion;
105 using FlowProblemType::enableEnergy;
106 using FlowProblemType::enableExperiments;
107 using FlowProblemType::enableExtbo;
108 using FlowProblemType::enableFoam;
109 using FlowProblemType::enableMICP;
110 using FlowProblemType::enablePolymer;
111 using FlowProblemType::enablePolymerMolarWeight;
112 using FlowProblemType::enableSaltPrecipitation;
113 using FlowProblemType::enableSolvent;
114 using FlowProblemType::enableTemperature;
115 using FlowProblemType::enableThermalFluxBoundaries;
117 using FlowProblemType::gasPhaseIdx;
118 using FlowProblemType::oilPhaseIdx;
119 using FlowProblemType::waterPhaseIdx;
121 using FlowProblemType::waterCompIdx;
122 using FlowProblemType::oilCompIdx;
123 using FlowProblemType::gasCompIdx;
126 using typename FlowProblemType::RateVector;
127 using typename FlowProblemType::PrimaryVariables;
128 using typename FlowProblemType::Indices;
129 using typename FlowProblemType::IntensiveQuantities;
130 using typename FlowProblemType::ElementContext;
132 using typename FlowProblemType::MaterialLaw;
133 using typename FlowProblemType::DimMatrix;
146 using InitialFluidState =
typename EquilInitializer<TypeTag>::ScalarFluidState;
164 EclWriterType::registerParameters();
166 DamarisWriterType::registerParameters();
176 , thresholdPressures_(simulator)
177 , mixControls_(simulator.vanguard().schedule())
178 , actionHandler_(simulator.vanguard().eclState(),
179 simulator.vanguard().schedule(),
180 simulator.vanguard().actionState(),
181 simulator.vanguard().summaryState(),
183 simulator.vanguard().grid().comm())
188 const auto& vanguard = simulator.vanguard();
192 enableSaltPrecipitation>(vanguard.eclState());
195 DiffusionModule::initFromState(vanguard.eclState());
196 DispersionModule::initFromState(vanguard.eclState());
219 eclWriter_ = std::make_unique<EclWriterType>(simulator);
220 enableEclOutput_ = Parameters::Get<Parameters::EnableEclOutput>();
225 enableDamarisOutput_ = Parameters::Get<Parameters::EnableDamarisOutput>();
236 auto& simulator = this->simulator();
238 const int episodeIdx = simulator.episodeIndex();
239 const auto& schedule = simulator.vanguard().schedule();
244 .evalUDQAssignments(
episodeIdx, simulator.vanguard().udqState());
248 if (
oilVap.getType() == OilVaporizationProperties::OilVaporization::VAPPARS) {
252 FluidSystem::setVapPars(0.0, 0.0);
256 ConvectiveMixingModule::beginEpisode(simulator.vanguard().eclState(), schedule,
episodeIdx,
257 this->moduleParams_.convectiveMixingModuleParam);
268 FlowProblemType::finishInit();
270 auto& simulator = this->simulator();
276 this->transmissibilities_.finishInit([&
vg = this->simulator().vanguard()](
const unsigned int it) {
277 return vg.gridIdxToEquilGridIdx(it);
287 if (enableEclOutput_) {
288 if (simulator.vanguard().grid().comm().size() > 1) {
289 if (simulator.vanguard().grid().comm().rank() == 0)
290 eclWriter_->setTransmissibilities(&simulator.vanguard().globalTransmissibility());
293 eclWriter_->setTransmissibilities(&simulator.problem().eclTransmissibilities());
296 std::function<
unsigned int(
unsigned int)>
equilGridToGrid = [&simulator](
unsigned int i) {
297 return simulator.vanguard().gridEquilIdxToGridIdx(i);
302 simulator.vanguard().releaseGlobalTransmissibilities();
304 const auto& eclState = simulator.vanguard().eclState();
305 const auto& schedule = simulator.vanguard().schedule();
308 simulator.setStartTime(schedule.getStartTime());
309 simulator.setEndTime(schedule.simTime(schedule.size() - 1));
315 simulator.setEpisodeIndex(-1);
316 simulator.setEpisodeLength(0.0);
321 this->gravity_ = 0.0;
322 if (Parameters::Get<Parameters::EnableGravity>() &&
323 eclState.getInitConfig().hasGravity())
326 this->gravity_[dim - 1] = unit::gravity;
329 if (this->enableTuning_) {
332 const auto&
tuning = schedule[0].tuning();
333 this->initialTimeStepSize_ =
tuning.TSINIT.has_value() ?
tuning.TSINIT.value() : -1.0;
334 this->maxTimeStepAfterWellEvent_ =
tuning.TMAXWC;
337 this->initFluidSystem_();
339 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
340 FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
341 this->maxOilSaturation_.resize(this->model().numGridDof(), 0.0);
344 this->readRockParameters_(simulator.vanguard().cellCenterDepths(),
345 [&simulator](
const unsigned idx)
347 std::array<int,dim> coords;
348 simulator.vanguard().cartesianCoordinate(idx, coords);
349 for (auto& c : coords) {
355 this->readMaterialParameters_();
356 this->readThermalParameters_();
359 if (enableEclOutput_) {
360 this->eclWriter_->writeInit();
365 const auto&
initconfig = eclState.getInitConfig();
366 this->tracerModel_.init(
initconfig.restartRequested());
368 this->readEclRestartSolution_();
371 this->readInitialCondition_();
374 this->tracerModel_.prepareTracerBatches();
376 this->updatePffDofData_();
379 const auto& vanguard = this->simulator().vanguard();
380 const auto& gridView = vanguard.gridView();
382 this->polymer_.maxAdsorption.resize(
numElements, 0.0);
385 this->readBoundaryConditions_();
388 this->computeAndSetEqWeights_();
390 if (this->enableDriftCompensation_) {
391 this->drift_.resize(this->model().numGridDof());
398 if (!
initconfig.restartRequested() && !eclState.getIOConfig().initOnly()) {
399 simulator.startNextEpisode(schedule.seconds(1));
400 simulator.setEpisodeIndex(0);
401 simulator.setTimeStepIndex(0);
404 if (Parameters::Get<Parameters::CheckSatfuncConsistency>() &&
405 ! this->satfuncConsistencyRequirementsMet())
414 this->simulator().vanguard().grid().comm().barrier();
416 throw std::domain_error {
417 "Saturation function end-points do not "
418 "meet requisite consistency conditions"
425 this->mixControls_.init(this->model().numGridDof(),
426 this->episodeIndex(),
427 eclState.runspec().tabdims().getNumPVTTables());
429 if (this->enableVtkOutput_() && eclState.getIOConfig().initOnly()) {
430 simulator.setTimeStepSize(0.0);
431 simulator.model().applyInitialSolution();
441 FlowProblemType::endTimeStep();
442 this->endStepApplyAction();
445 void endStepApplyAction()
449 this->eclWriter().mutableOutputModule().invalidateLocalData();
452 const auto& grid = this->simulator().vanguard().gridView().grid();
454 using GridType = std::remove_cv_t<std::remove_reference_t<
decltype(grid)>>;
455 constexpr bool isCpGrid = std::is_same_v<GridType, Dune::CpGrid>;
456 if (!
isCpGrid || (grid.maxLevel() == 0)) {
457 const bool isSubStep = !this->simulator().episodeWillBeOver();
459 this->eclWriter_->evalSummaryState(
isSubStep);
466 auto& simulator = this->simulator();
470 this->simulator().vanguard().schedule().clearEvents(
episodeIdx);
474 .applyActions(
episodeIdx, simulator.time() + simulator.timeStepSize(),
475 [
this](
const bool global)
477 using TransUpdateQuantities = typename
478 Vanguard::TransmissibilityType::TransUpdateQuantities;
480 this->transmissibilities_
481 .update(global, TransUpdateQuantities::All,
482 [&vg = this->simulator().vanguard()]
483 (const unsigned int i)
485 return vg.gridIdxToEquilGridIdx(i);
496 OPM_TIMEBLOCK(endEpisode);
509 .evalUDQAssignments(this->episodeIndex(), this->simulator().vanguard().udqState());
511 FlowProblemType::endEpisode();
516 if (this->enableEclOutput_) {
517 this->eclWriter_->writeReports(timer);
528 FlowProblemType::writeOutput(verbose);
530 bool isSubStep = !this->simulator().episodeWillBeOver();
532 data::Solution localCellData = {};
536 if (enableDamarisOutput_) {
537 damarisWriter_->writeOutput(localCellData, isSubStep) ;
540 if (enableEclOutput_){
541 eclWriter_->writeOutput(std::move(localCellData), isSubStep);
545 void finalizeOutput()
547 OPM_TIMEBLOCK(finalizeOutput);
559 FlowProblemType::initialSolutionApplied();
564 this->thresholdPressures_.finishInit();
567 const auto& grid = this->simulator().vanguard().gridView().grid();
569 using GridType = std::remove_cv_t<std::remove_reference_t<
decltype(grid)>>;
570 constexpr bool isCpGrid = std::is_same_v<GridType, Dune::CpGrid>;
572 if (!isCpGrid || (grid.maxLevel() == 0)) {
573 if (this->simulator().episodeIndex() == 0) {
574 eclWriter_->writeInitialFIPReport();
579 void addToSourceDense(RateVector& rate,
580 unsigned globalDofIdx,
581 unsigned timeIdx)
const override
583 this->aquiferModel_.addToSource(rate, globalDofIdx, timeIdx);
586 const auto& source = this->simulator().vanguard().schedule()[this->episodeIndex()].source();
587 std::array<int,3> ijk;
588 this->simulator().vanguard().cartesianCoordinate(globalDofIdx, ijk);
590 if (source.hasSource(ijk)) {
591 const int pvtRegionIdx = this->pvtRegionIndex(globalDofIdx);
592 static std::array<SourceComponent, 3> sc_map = {SourceComponent::WATER, SourceComponent::OIL, SourceComponent::GAS};
593 static std::array<int, 3> phidx_map = {FluidSystem::waterPhaseIdx, FluidSystem::oilPhaseIdx, FluidSystem::gasPhaseIdx};
594 static std::array<int, 3> cidx_map = {waterCompIdx, oilCompIdx, gasCompIdx};
596 for (
unsigned i = 0; i < phidx_map.size(); ++i) {
597 const auto phaseIdx = phidx_map[i];
598 const auto sourceComp = sc_map[i];
599 const auto compIdx = cidx_map[i];
600 if (!FluidSystem::phaseIsActive(phaseIdx)) {
603 Scalar mass_rate = source.rate(ijk, sourceComp) / this->model().dofTotalVolume(globalDofIdx);
604 if constexpr (getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>()) {
605 mass_rate /= FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
607 rate[Indices::canonicalToActiveComponentIndex(compIdx)] += mass_rate;
610 if constexpr (enableSolvent) {
611 Scalar mass_rate = source.rate(ijk, SourceComponent::SOLVENT) / this->model().dofTotalVolume(globalDofIdx);
612 if constexpr (getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>()) {
613 const auto& solventPvt = SolventModule::solventPvt();
614 mass_rate /= solventPvt.referenceDensity(pvtRegionIdx);
616 rate[Indices::contiSolventEqIdx] += mass_rate;
618 if constexpr (enablePolymer) {
619 rate[Indices::polymerConcentrationIdx] += source.rate(ijk, SourceComponent::POLYMER) / this->model().dofTotalVolume(globalDofIdx);
621 if constexpr (enableMICP) {
622 rate[Indices::microbialConcentrationIdx] += source.rate(ijk, SourceComponent::MICR) / this->model().dofTotalVolume(globalDofIdx);
623 rate[Indices::oxygenConcentrationIdx] += source.rate(ijk, SourceComponent::OXYG) / this->model().dofTotalVolume(globalDofIdx);
624 rate[Indices::ureaConcentrationIdx] += source.rate(ijk, SourceComponent::UREA) / (this->model().dofTotalVolume(globalDofIdx));
626 if constexpr (enableEnergy) {
627 for (
unsigned i = 0; i < phidx_map.size(); ++i) {
628 const auto phaseIdx = phidx_map[i];
629 if (!FluidSystem::phaseIsActive(phaseIdx)) {
632 const auto sourceComp = sc_map[i];
633 const auto source_hrate = source.hrate(ijk, sourceComp);
635 rate[Indices::contiEnergyEqIdx] += source_hrate.value() / this->model().dofTotalVolume(globalDofIdx);
637 const auto& intQuants = this->simulator().model().intensiveQuantities(globalDofIdx, 0);
638 auto fs = intQuants.fluidState();
640 const auto source_temp = source.temperature(ijk, sourceComp);
642 Scalar temperature = source_temp.value();
643 fs.setTemperature(temperature);
645 const auto& h = FluidSystem::enthalpy(fs, phaseIdx, pvtRegionIdx);
646 Scalar mass_rate = source.rate(ijk, sourceComp)/ this->model().dofTotalVolume(globalDofIdx);
647 Scalar energy_rate = getValue(h)*mass_rate;
648 rate[Indices::contiEnergyEqIdx] += energy_rate;
656 if (this->enableDriftCompensation_) {
657 const auto& simulator = this->simulator();
658 const auto& model = this->model();
663 Scalar maxCompensation = model.newtonMethod().tolerance()/10;
664 Scalar poro = this->porosity(globalDofIdx, timeIdx);
665 Scalar dt = simulator.timeStepSize();
666 EqVector dofDriftRate = this->drift_[globalDofIdx];
667 dofDriftRate /= dt*model.dofTotalVolume(globalDofIdx);
670 for (
unsigned eqIdx = 0; eqIdx < numEq; ++ eqIdx) {
671 Scalar cnv = std::abs(dofDriftRate[eqIdx])*dt*model.eqWeight(globalDofIdx, eqIdx)/poro;
672 if (cnv > maxCompensation) {
673 dofDriftRate[eqIdx] *= maxCompensation/cnv;
677 for (
unsigned eqIdx = 0; eqIdx < numEq; ++ eqIdx)
678 rate[eqIdx] -= dofDriftRate[eqIdx];
687 template <
class LhsEval>
690 OPM_TIMEBLOCK_LOCAL(permFactTransMultiplier);
691 if constexpr (enableSaltPrecipitation) {
692 const auto& fs = intQuants.fluidState();
693 unsigned tableIdx = this->simulator().problem().satnumRegionIndex(elementIdx);
694 LhsEval porosityFactor = decay<LhsEval>(1. - fs.saltSaturation());
695 porosityFactor = min(porosityFactor, 1.0);
696 const auto& permfactTable = BrineModule::permfactTable(tableIdx);
697 return permfactTable.eval(porosityFactor,
true);
699 else if constexpr (enableMICP) {
700 return intQuants.permFactor().value();
708 const InitialFluidState& initialFluidState(
unsigned globalDofIdx)
const
709 {
return initialFluidStates_[globalDofIdx]; }
711 std::vector<InitialFluidState>& initialFluidStates()
712 {
return initialFluidStates_; }
714 const std::vector<InitialFluidState>& initialFluidStates()
const
715 {
return initialFluidStates_; }
717 const EclipseIO& eclIO()
const
718 {
return eclWriter_->eclIO(); }
720 void setSubStepReport(
const SimulatorReportSingle& report)
721 {
return eclWriter_->setSubStepReport(report); }
723 void setSimulationReport(
const SimulatorReport& report)
724 {
return eclWriter_->setSimulationReport(report); }
726 InitialFluidState boundaryFluidState(
unsigned globalDofIdx,
const int directionId)
const
728 OPM_TIMEBLOCK_LOCAL(boundaryFluidState);
729 const auto& bcprop = this->simulator().vanguard().schedule()[this->episodeIndex()].bcprop;
730 if (bcprop.size() > 0) {
731 FaceDir::DirEnum dir = FaceDir::FromIntersectionIndex(directionId);
735 if (this->bcindex_(dir)[globalDofIdx] == 0)
736 return initialFluidStates_[globalDofIdx];
738 const auto& bc = bcprop[this->bcindex_(dir)[globalDofIdx]];
739 if (bc.bctype == BCType::DIRICHLET )
741 InitialFluidState fluidState;
742 const int pvtRegionIdx = this->pvtRegionIndex(globalDofIdx);
743 fluidState.setPvtRegionIndex(pvtRegionIdx);
745 switch (bc.component) {
746 case BCComponent::OIL:
747 if (!FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx))
748 throw std::logic_error(
"oil is not active and you're trying to add oil BC");
750 fluidState.setSaturation(FluidSystem::oilPhaseIdx, 1.0);
752 case BCComponent::GAS:
753 if (!FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
754 throw std::logic_error(
"gas is not active and you're trying to add gas BC");
756 fluidState.setSaturation(FluidSystem::gasPhaseIdx, 1.0);
758 case BCComponent::WATER:
759 if (!FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx))
760 throw std::logic_error(
"water is not active and you're trying to add water BC");
762 fluidState.setSaturation(FluidSystem::waterPhaseIdx, 1.0);
764 case BCComponent::SOLVENT:
765 case BCComponent::POLYMER:
766 case BCComponent::MICR:
767 case BCComponent::OXYG:
768 case BCComponent::UREA:
769 case BCComponent::NONE:
770 throw std::logic_error(
"you need to specify a valid component (OIL, WATER or GAS) when DIRICHLET type is set in BC");
772 fluidState.setTotalSaturation(1.0);
773 double pressure = initialFluidStates_[globalDofIdx].pressure(this->refPressurePhaseIdx_());
774 const auto pressure_input = bc.pressure;
775 if (pressure_input) {
776 pressure = *pressure_input;
779 std::array<Scalar, numPhases> pc = {0};
780 const auto& matParams = this->materialLawParams(globalDofIdx);
781 MaterialLaw::capillaryPressures(pc, matParams, fluidState);
782 Valgrind::CheckDefined(pressure);
783 Valgrind::CheckDefined(pc);
784 for (
unsigned activePhaseIdx = 0; activePhaseIdx < FluidSystem::numActivePhases(); ++activePhaseIdx) {
785 const auto phaseIdx = FluidSystem::activeToCanonicalPhaseIdx(activePhaseIdx);
787 fluidState.setPc(phaseIdx, pc[phaseIdx]);
788 if (Indices::oilEnabled)
789 fluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[oilPhaseIdx]));
790 else if (Indices::gasEnabled)
791 fluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[gasPhaseIdx]));
792 else if (Indices::waterEnabled)
794 fluidState.setPressure(phaseIdx, pressure);
797 double temperature = initialFluidStates_[globalDofIdx].temperature(0);
798 const auto temperature_input = bc.temperature;
799 if(temperature_input)
800 temperature = *temperature_input;
801 fluidState.setTemperature(temperature);
803 if (FluidSystem::enableDissolvedGas()) {
804 fluidState.setRs(0.0);
805 fluidState.setRv(0.0);
807 if (FluidSystem::enableDissolvedGasInWater()) {
808 fluidState.setRsw(0.0);
810 if (FluidSystem::enableVaporizedWater())
811 fluidState.setRvw(0.0);
813 for (
unsigned activePhaseIdx = 0; activePhaseIdx < FluidSystem::numActivePhases(); ++activePhaseIdx) {
814 const auto phaseIdx = FluidSystem::activeToCanonicalPhaseIdx(activePhaseIdx);
816 const auto& b = FluidSystem::inverseFormationVolumeFactor(fluidState, phaseIdx, pvtRegionIdx);
817 fluidState.setInvB(phaseIdx, b);
819 const auto& rho = FluidSystem::density(fluidState, phaseIdx, pvtRegionIdx);
820 fluidState.setDensity(phaseIdx, rho);
822 const auto& h = FluidSystem::enthalpy(fluidState, phaseIdx, pvtRegionIdx);
823 fluidState.setEnthalpy(phaseIdx, h);
826 fluidState.checkDefined();
830 return initialFluidStates_[globalDofIdx];
834 const EclWriterType& eclWriter()
const
835 {
return *eclWriter_; }
837 EclWriterType& eclWriter()
838 {
return *eclWriter_; }
846 return this->mixControls_.maxGasDissolutionFactor(timeIdx, globalDofIdx,
847 this->episodeIndex(),
848 this->pvtRegionIndex(globalDofIdx));
857 return this->mixControls_.maxOilVaporizationFactor(timeIdx, globalDofIdx,
858 this->episodeIndex(),
859 this->pvtRegionIndex(globalDofIdx));
872 int episodeIdx = this->episodeIndex();
873 return !this->mixControls_.drsdtActive(episodeIdx) &&
874 !this->mixControls_.drvdtActive(episodeIdx) &&
875 this->rockCompPoroMultWc_.empty() &&
876 this->rockCompPoroMult_.empty();
885 template <
class Context>
886 void initial(PrimaryVariables& values,
const Context& context,
unsigned spaceIdx,
unsigned timeIdx)
const
888 unsigned globalDofIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
890 values.setPvtRegionIndex(pvtRegionIndex(context, spaceIdx, timeIdx));
891 values.assignNaive(initialFluidStates_[globalDofIdx]);
893 SolventModule::assignPrimaryVars(values,
894 enableSolvent ? this->solventSaturation_[globalDofIdx] : 0.0,
895 enableSolvent ? this->solventRsw_[globalDofIdx] : 0.0);
897 if constexpr (enablePolymer)
898 values[Indices::polymerConcentrationIdx] = this->polymer_.concentration[globalDofIdx];
900 if constexpr (enablePolymerMolarWeight)
901 values[Indices::polymerMoleWeightIdx]= this->polymer_.moleWeight[globalDofIdx];
903 if constexpr (enableBrine) {
904 if (enableSaltPrecipitation && values.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) {
905 values[Indices::saltConcentrationIdx] = initialFluidStates_[globalDofIdx].saltSaturation();
908 values[Indices::saltConcentrationIdx] = initialFluidStates_[globalDofIdx].saltConcentration();
912 if constexpr (enableMICP){
913 values[Indices::microbialConcentrationIdx] = this->micp_.microbialConcentration[globalDofIdx];
914 values[Indices::oxygenConcentrationIdx]= this->micp_.oxygenConcentration[globalDofIdx];
915 values[Indices::ureaConcentrationIdx]= this->micp_.ureaConcentration[globalDofIdx];
916 values[Indices::calciteConcentrationIdx]= this->micp_.calciteConcentration[globalDofIdx];
917 values[Indices::biofilmConcentrationIdx]= this->micp_.biofilmConcentration[globalDofIdx];
920 values.checkDefined();
924 Scalar drsdtcon(
unsigned elemIdx,
int episodeIdx)
const
926 return this->mixControls_.drsdtcon(elemIdx, episodeIdx,
927 this->pvtRegionIndex(elemIdx));
930 bool drsdtconIsActive(
unsigned elemIdx,
int episodeIdx)
const
932 return this->mixControls_.drsdtConvective(episodeIdx, this->pvtRegionIndex(elemIdx));
940 template <
class Context>
942 const Context& context,
944 unsigned timeIdx)
const
946 OPM_TIMEBLOCK_LOCAL(eclProblemBoundary);
947 if (!context.intersection(spaceIdx).boundary())
950 if constexpr (!enableEnergy || !enableThermalFluxBoundaries)
958 unsigned interiorDofIdx = context.interiorScvIndex(spaceIdx, timeIdx);
959 unsigned globalDofIdx = context.globalSpaceIndex(interiorDofIdx, timeIdx);
960 values.setThermalFlow(context, spaceIdx, timeIdx, this->initialFluidStates_[globalDofIdx] );
963 if (this->nonTrivialBoundaryConditions()) {
964 unsigned indexInInside = context.intersection(spaceIdx).indexInInside();
965 unsigned interiorDofIdx = context.interiorScvIndex(spaceIdx, timeIdx);
966 unsigned globalDofIdx = context.globalSpaceIndex(interiorDofIdx, timeIdx);
967 unsigned pvtRegionIdx = pvtRegionIndex(context, spaceIdx, timeIdx);
968 const auto [type, massrate] = this->boundaryCondition(globalDofIdx, indexInInside);
969 if (type == BCType::THERMAL)
970 values.setThermalFlow(context, spaceIdx, timeIdx, this->boundaryFluidState(globalDofIdx, indexInInside));
971 else if (type == BCType::FREE || type == BCType::DIRICHLET)
972 values.setFreeFlow(context, spaceIdx, timeIdx, this->boundaryFluidState(globalDofIdx, indexInInside));
973 else if (type == BCType::RATE)
974 values.setMassRate(massrate, pvtRegionIdx);
984 auto& simulator = this->simulator();
985 const auto& eclState = simulator.vanguard().eclState();
987 std::size_t numElems = this->model().numGridDof();
988 this->initialFluidStates_.resize(numElems);
989 if constexpr (enableSolvent) {
990 this->solventSaturation_.resize(numElems, 0.0);
991 this->solventRsw_.resize(numElems, 0.0);
994 if constexpr (enablePolymer)
995 this->polymer_.concentration.resize(numElems, 0.0);
997 if constexpr (enablePolymerMolarWeight) {
998 const std::string msg {
"Support of the RESTART for polymer molecular weight "
999 "is not implemented yet. The polymer weight value will be "
1000 "zero when RESTART begins"};
1001 OpmLog::warning(
"NO_POLYMW_RESTART", msg);
1002 this->polymer_.moleWeight.resize(numElems, 0.0);
1005 if constexpr (enableMICP) {
1006 this->micp_.resize(numElems);
1010 this->mixControls_.init(numElems, restart_step, eclState.runspec().tabdims().getNumPVTTables());
1012 if constexpr (enableMICP) {
1013 this->micp_ = this->eclWriter_->outputModule().getMICP().getSolution();
1016 for (std::size_t elemIdx = 0; elemIdx < numElems; ++elemIdx) {
1017 auto& elemFluidState = this->initialFluidStates_[elemIdx];
1018 elemFluidState.setPvtRegionIndex(pvtRegionIndex(elemIdx));
1019 this->eclWriter_->outputModule().initHysteresisParams(simulator, elemIdx);
1020 this->eclWriter_->outputModule().assignToFluidState(elemFluidState, elemIdx);
1029 auto ssol = enableSolvent
1030 ? this->eclWriter_->outputModule().getSolventSaturation(elemIdx)
1033 this->processRestartSaturations_(elemFluidState, ssol);
1035 if constexpr (enableSolvent) {
1036 this->solventSaturation_[elemIdx] = ssol;
1037 this->solventRsw_[elemIdx] = this->eclWriter_->outputModule().getSolventRsw(elemIdx);
1042 bool isThermal = eclState.getSimulationConfig().isThermal();
1043 bool needTemperature = (eclState.runspec().co2Storage() || eclState.runspec().h2Storage());
1044 if (!isThermal && needTemperature) {
1045 const auto& fp = simulator.vanguard().eclState().fieldProps();
1046 elemFluidState.setTemperature(fp.get_double(
"TEMPI")[elemIdx]);
1049 this->mixControls_.updateLastValues(elemIdx, elemFluidState.Rs(), elemFluidState.Rv());
1051 if constexpr (enablePolymer)
1052 this->polymer_.concentration[elemIdx] = this->eclWriter_->outputModule().getPolymerConcentration(elemIdx);
1056 const int episodeIdx = this->episodeIndex();
1057 this->mixControls_.updateMaxValues(episodeIdx, simulator.timeStepSize());
1062 auto& sol = this->model().solution(0);
1063 const auto& gridView = this->gridView();
1064 ElementContext elemCtx(simulator);
1065 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
1066 elemCtx.updatePrimaryStencil(elem);
1067 int elemIdx = elemCtx.globalSpaceIndex(0, 0);
1068 this->initial(sol[elemIdx], elemCtx, 0, 0);
1076 this->model().syncOverlap();
1079 this->updateReferencePorosity_();
1080 this->mixControls_.init(this->model().numGridDof(),
1081 this->episodeIndex(),
1082 eclState.runspec().tabdims().getNumPVTTables());
1090 {
return thresholdPressures_.thresholdPressure(elem1Idx, elem2Idx); }
1093 {
return thresholdPressures_; }
1095 FlowThresholdPressure<TypeTag>& thresholdPressure()
1096 {
return thresholdPressures_; }
1098 const ModuleParams& moduleParams()
const
1100 return moduleParams_;
1103 template<
class Serializer>
1104 void serializeOp(Serializer& serializer)
1106 serializer(
static_cast<FlowProblemType&
>(*
this));
1107 serializer(mixControls_);
1108 serializer(*eclWriter_);
1112 void updateExplicitQuantities_(
int episodeIdx,
int timeStepSize,
const bool first_step_after_restart)
override
1114 this->updateExplicitQuantities_(first_step_after_restart);
1116 if constexpr (getPropValue<TypeTag, Properties::EnablePolymer>())
1117 updateMaxPolymerAdsorption_();
1119 mixControls_.updateExplicitQuantities(episodeIdx, timeStepSize);
1122 void updateMaxPolymerAdsorption_()
1125 this->updateProperty_(
"FlowProblemBlackoil::updateMaxPolymerAdsorption_() failed:",
1126 [
this](
unsigned compressedDofIdx,
const IntensiveQuantities& iq)
1128 this->updateMaxPolymerAdsorption_(compressedDofIdx,iq);
1132 bool updateMaxPolymerAdsorption_(
unsigned compressedDofIdx,
const IntensiveQuantities& iq)
1134 const Scalar pa = scalarValue(iq.polymerAdsorption());
1135 auto& mpa = this->polymer_.maxAdsorption;
1136 if (mpa[compressedDofIdx] < pa) {
1137 mpa[compressedDofIdx] = pa;
1144 void computeAndSetEqWeights_()
1146 std::vector<Scalar> sumInvB(numPhases, 0.0);
1147 const auto& gridView = this->gridView();
1148 ElementContext elemCtx(this->simulator());
1149 for(
const auto& elem: elements(gridView, Dune::Partitions::interior)) {
1150 elemCtx.updatePrimaryStencil(elem);
1151 int elemIdx = elemCtx.globalSpaceIndex(0, 0);
1152 const auto& dofFluidState = this->initialFluidStates_[elemIdx];
1153 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1154 if (!FluidSystem::phaseIsActive(phaseIdx))
1157 sumInvB[phaseIdx] += dofFluidState.invB(phaseIdx);
1161 std::size_t numDof = this->model().numGridDof();
1162 const auto& comm = this->simulator().vanguard().grid().comm();
1163 comm.sum(sumInvB.data(),sumInvB.size());
1164 Scalar numTotalDof = comm.sum(numDof);
1166 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1167 if (!FluidSystem::phaseIsActive(phaseIdx))
1170 Scalar avgB = numTotalDof / sumInvB[phaseIdx];
1171 unsigned solventCompIdx = FluidSystem::solventComponentIndex(phaseIdx);
1172 unsigned activeSolventCompIdx = Indices::canonicalToActiveComponentIndex(solventCompIdx);
1173 this->model().setEqWeight(activeSolventCompIdx, avgB);
1178 bool updateCompositionChangeLimits_()
1180 OPM_TIMEBLOCK(updateCompositionChangeLimits);
1183 int episodeIdx = this->episodeIndex();
1184 std::array<bool,3> active{this->mixControls_.drsdtConvective(episodeIdx),
1185 this->mixControls_.drsdtActive(episodeIdx),
1186 this->mixControls_.drvdtActive(episodeIdx)};
1187 if (!active[0] && !active[1] && !active[2]) {
1191 this->updateProperty_(
"FlowProblemBlackoil::updateCompositionChangeLimits_()) failed:",
1192 [
this,episodeIdx,active](
unsigned compressedDofIdx,
1193 const IntensiveQuantities& iq)
1195 const DimMatrix& perm = this->intrinsicPermeability(compressedDofIdx);
1196 const Scalar distZ = active[0] ? this->simulator().vanguard().cellThickness(compressedDofIdx) : 0.0;
1197 const int pvtRegionIdx = this->pvtRegionIndex(compressedDofIdx);
1198 this->mixControls_.update(compressedDofIdx,
1201 this->gravity_[dim - 1],
1202 perm[dim - 1][dim - 1],
1211 void readEclRestartSolution_()
1214 if(this->simulator().vanguard().grid().maxLevel() > 0) {
1215 throw std::invalid_argument(
"Refined grids are not yet supported for restart ");
1219 auto& simulator = this->simulator();
1220 const auto& schedule = simulator.vanguard().schedule();
1221 const auto& eclState = simulator.vanguard().eclState();
1222 const auto& initconfig = eclState.getInitConfig();
1223 const int restart_step = initconfig.getRestartStep();
1225 simulator.setTime(schedule.seconds(restart_step));
1227 simulator.startNextEpisode(simulator.startTime() + simulator.time(),
1228 schedule.stepLength(restart_step));
1229 simulator.setEpisodeIndex(restart_step);
1231 this->eclWriter_->beginRestart();
1233 Scalar dt = std::min(this->eclWriter_->restartTimeStepSize(), simulator.episodeLength());
1234 simulator.setTimeStepSize(dt);
1236 this->readSolutionFromOutputModule(restart_step,
false);
1238 this->eclWriter_->endRestart();
1241 void readEquilInitialCondition_()
override
1243 const auto& simulator = this->simulator();
1246 EquilInitializer<TypeTag> equilInitializer(simulator, *(this->materialLawManager_));
1248 std::size_t numElems = this->model().numGridDof();
1249 this->initialFluidStates_.resize(numElems);
1250 for (std::size_t elemIdx = 0; elemIdx < numElems; ++elemIdx) {
1251 auto& elemFluidState = this->initialFluidStates_[elemIdx];
1252 elemFluidState.assign(equilInitializer.initialFluidState(elemIdx));
1256 void readExplicitInitialCondition_()
override
1258 const auto& simulator = this->simulator();
1259 const auto& vanguard = simulator.vanguard();
1260 const auto& eclState = vanguard.eclState();
1261 const auto& fp = eclState.fieldProps();
1262 bool has_swat = fp.has_double(
"SWAT");
1263 bool has_sgas = fp.has_double(
"SGAS");
1264 bool has_rs = fp.has_double(
"RS");
1265 bool has_rsw = fp.has_double(
"RSW");
1266 bool has_rv = fp.has_double(
"RV");
1267 bool has_rvw = fp.has_double(
"RVW");
1268 bool has_pressure = fp.has_double(
"PRESSURE");
1269 bool has_salt = fp.has_double(
"SALT");
1270 bool has_saltp = fp.has_double(
"SALTP");
1273 if (Indices::numPhases > 1) {
1274 if (FluidSystem::phaseIsActive(waterPhaseIdx) && !has_swat)
1275 throw std::runtime_error(
"The ECL input file requires the presence of the SWAT keyword if "
1276 "the water phase is active");
1277 if (FluidSystem::phaseIsActive(gasPhaseIdx) && !has_sgas && FluidSystem::phaseIsActive(oilPhaseIdx))
1278 throw std::runtime_error(
"The ECL input file requires the presence of the SGAS keyword if "
1279 "the gas phase is active");
1282 throw std::runtime_error(
"The ECL input file requires the presence of the PRESSURE "
1283 "keyword if the model is initialized explicitly");
1284 if (FluidSystem::enableDissolvedGas() && !has_rs)
1285 throw std::runtime_error(
"The ECL input file requires the RS keyword to be present if"
1286 " dissolved gas is enabled and the model is initialized explicitly");
1287 if (FluidSystem::enableDissolvedGasInWater() && !has_rsw)
1288 OpmLog::warning(
"The model is initialized explicitly and the RSW keyword is not present in the"
1289 " ECL input file. The RSW values are set equal to 0");
1290 if (FluidSystem::enableVaporizedOil() && !has_rv)
1291 throw std::runtime_error(
"The ECL input file requires the RV keyword to be present if"
1292 " vaporized oil is enabled and the model is initialized explicitly");
1293 if (FluidSystem::enableVaporizedWater() && !has_rvw)
1294 throw std::runtime_error(
"The ECL input file requires the RVW keyword to be present if"
1295 " vaporized water is enabled and the model is initialized explicitly");
1296 if (enableBrine && !has_salt)
1297 throw std::runtime_error(
"The ECL input file requires the SALT keyword to be present if"
1298 " brine is enabled and the model is initialized explicitly");
1299 if (enableSaltPrecipitation && !has_saltp)
1300 throw std::runtime_error(
"The ECL input file requires the SALTP keyword to be present if"
1301 " salt precipitation is enabled and the model is initialized explicitly");
1303 std::size_t numDof = this->model().numGridDof();
1305 initialFluidStates_.resize(numDof);
1307 std::vector<double> waterSaturationData;
1308 std::vector<double> gasSaturationData;
1309 std::vector<double> pressureData;
1310 std::vector<double> rsData;
1311 std::vector<double> rswData;
1312 std::vector<double> rvData;
1313 std::vector<double> rvwData;
1314 std::vector<double> tempiData;
1315 std::vector<double> saltData;
1316 std::vector<double> saltpData;
1318 if (FluidSystem::phaseIsActive(waterPhaseIdx) && Indices::numPhases > 1)
1319 waterSaturationData = fp.get_double(
"SWAT");
1321 waterSaturationData.resize(numDof);
1323 if (FluidSystem::phaseIsActive(gasPhaseIdx) && FluidSystem::phaseIsActive(oilPhaseIdx))
1324 gasSaturationData = fp.get_double(
"SGAS");
1326 gasSaturationData.resize(numDof);
1328 pressureData = fp.get_double(
"PRESSURE");
1329 if (FluidSystem::enableDissolvedGas())
1330 rsData = fp.get_double(
"RS");
1332 if (FluidSystem::enableDissolvedGasInWater() && has_rsw)
1333 rswData = fp.get_double(
"RSW");
1335 if (FluidSystem::enableVaporizedOil())
1336 rvData = fp.get_double(
"RV");
1338 if (FluidSystem::enableVaporizedWater())
1339 rvwData = fp.get_double(
"RVW");
1342 tempiData = fp.get_double(
"TEMPI");
1345 if constexpr (enableBrine)
1346 saltData = fp.get_double(
"SALT");
1349 if constexpr (enableSaltPrecipitation)
1350 saltpData = fp.get_double(
"SALTP");
1353 for (std::size_t dofIdx = 0; dofIdx < numDof; ++dofIdx) {
1354 auto& dofFluidState = initialFluidStates_[dofIdx];
1356 dofFluidState.setPvtRegionIndex(pvtRegionIndex(dofIdx));
1361 Scalar temperatureLoc = tempiData[dofIdx];
1362 if (!std::isfinite(temperatureLoc) || temperatureLoc <= 0)
1363 temperatureLoc = FluidSystem::surfaceTemperature;
1364 dofFluidState.setTemperature(temperatureLoc);
1369 if constexpr (enableBrine)
1370 dofFluidState.setSaltConcentration(saltData[dofIdx]);
1375 if constexpr (enableSaltPrecipitation)
1376 dofFluidState.setSaltSaturation(saltpData[dofIdx]);
1381 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx))
1382 dofFluidState.setSaturation(FluidSystem::waterPhaseIdx,
1383 waterSaturationData[dofIdx]);
1385 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)){
1386 if (!FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)){
1387 dofFluidState.setSaturation(FluidSystem::gasPhaseIdx,
1389 - waterSaturationData[dofIdx]);
1392 dofFluidState.setSaturation(FluidSystem::gasPhaseIdx,
1393 gasSaturationData[dofIdx]);
1395 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
1396 const Scalar soil = 1.0 - waterSaturationData[dofIdx] - gasSaturationData[dofIdx];
1397 if (soil < smallSaturationTolerance_) {
1398 dofFluidState.setSaturation(FluidSystem::oilPhaseIdx, 0.0);
1401 dofFluidState.setSaturation(FluidSystem::oilPhaseIdx, soil);
1408 Scalar pressure = pressureData[dofIdx];
1412 std::array<Scalar, numPhases> pc = {0};
1413 const auto& matParams = this->materialLawParams(dofIdx);
1414 MaterialLaw::capillaryPressures(pc, matParams, dofFluidState);
1415 Valgrind::CheckDefined(pressure);
1416 Valgrind::CheckDefined(pc);
1417 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1418 if (!FluidSystem::phaseIsActive(phaseIdx))
1421 if (Indices::oilEnabled)
1422 dofFluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[oilPhaseIdx]));
1423 else if (Indices::gasEnabled)
1424 dofFluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[gasPhaseIdx]));
1425 else if (Indices::waterEnabled)
1427 dofFluidState.setPressure(phaseIdx, pressure);
1430 if (FluidSystem::enableDissolvedGas())
1431 dofFluidState.setRs(rsData[dofIdx]);
1432 else if (Indices::gasEnabled && Indices::oilEnabled)
1433 dofFluidState.setRs(0.0);
1435 if (FluidSystem::enableDissolvedGasInWater() && has_rsw)
1436 dofFluidState.setRsw(rswData[dofIdx]);
1438 if (FluidSystem::enableVaporizedOil())
1439 dofFluidState.setRv(rvData[dofIdx]);
1440 else if (Indices::gasEnabled && Indices::oilEnabled)
1441 dofFluidState.setRv(0.0);
1443 if (FluidSystem::enableVaporizedWater())
1444 dofFluidState.setRvw(rvwData[dofIdx]);
1449 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1450 if (!FluidSystem::phaseIsActive(phaseIdx))
1453 const auto& b = FluidSystem::inverseFormationVolumeFactor(dofFluidState, phaseIdx, pvtRegionIndex(dofIdx));
1454 dofFluidState.setInvB(phaseIdx, b);
1456 const auto& rho = FluidSystem::density(dofFluidState, phaseIdx, pvtRegionIndex(dofIdx));
1457 dofFluidState.setDensity(phaseIdx, rho);
1464 void processRestartSaturations_(InitialFluidState& elemFluidState, Scalar& solventSaturation)
1468 Scalar sumSaturation = 0.0;
1469 for (std::size_t phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1470 if (FluidSystem::phaseIsActive(phaseIdx)) {
1471 if (elemFluidState.saturation(phaseIdx) < smallSaturationTolerance_)
1472 elemFluidState.setSaturation(phaseIdx, 0.0);
1474 sumSaturation += elemFluidState.saturation(phaseIdx);
1478 if constexpr (enableSolvent) {
1479 if (solventSaturation < smallSaturationTolerance_)
1480 solventSaturation = 0.0;
1482 sumSaturation += solventSaturation;
1485 assert(sumSaturation > 0.0);
1487 for (std::size_t phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1488 if (FluidSystem::phaseIsActive(phaseIdx)) {
1489 const Scalar saturation = elemFluidState.saturation(phaseIdx) / sumSaturation;
1490 elemFluidState.setSaturation(phaseIdx, saturation);
1493 if constexpr (enableSolvent) {
1494 solventSaturation = solventSaturation / sumSaturation;
1498 void readInitialCondition_()
override
1500 FlowProblemType::readInitialCondition_();
1502 if constexpr (enableSolvent || enablePolymer || enablePolymerMolarWeight || enableMICP)
1503 this->readBlackoilExtentionsInitialConditions_(this->model().numGridDof(),
1506 enablePolymerMolarWeight,
1511 void handleSolventBC(
const BCProp::BCFace& bc, RateVector& rate)
const override
1513 if constexpr (!enableSolvent)
1514 throw std::logic_error(
"solvent is disabled and you're trying to add solvent to BC");
1516 rate[Indices::solventSaturationIdx] = bc.rate;
1519 void handlePolymerBC(
const BCProp::BCFace& bc, RateVector& rate)
const override
1521 if constexpr (!enablePolymer)
1522 throw std::logic_error(
"polymer is disabled and you're trying to add polymer to BC");
1524 rate[Indices::polymerConcentrationIdx] = bc.rate;
1527 void handleMicrBC(
const BCProp::BCFace& bc, RateVector& rate)
const override
1529 if constexpr (!enableMICP)
1530 throw std::logic_error(
"MICP is disabled and you're trying to add microbes to BC");
1532 rate[Indices::microbialConcentrationIdx] = bc.rate;
1535 void handleOxygBC(
const BCProp::BCFace& bc, RateVector& rate)
const override
1537 if constexpr (!enableMICP)
1538 throw std::logic_error(
"MICP is disabled and you're trying to add oxygen to BC");
1540 rate[Indices::oxygenConcentrationIdx] = bc.rate;
1543 void handleUreaBC(
const BCProp::BCFace& bc, RateVector& rate)
const override
1545 if constexpr (!enableMICP)
1546 throw std::logic_error(
"MICP is disabled and you're trying to add urea to BC");
1548 rate[Indices::ureaConcentrationIdx] = bc.rate;
1550 rate[Indices::ureaConcentrationIdx] *= getPropValue<TypeTag, Properties::BlackOilUreaScalingFactor>();
1553 void updateExplicitQuantities_(
const bool first_step_after_restart)
1555 OPM_TIMEBLOCK(updateExplicitQuantities);
1556 const bool invalidateFromMaxWaterSat = this->updateMaxWaterSaturation_();
1557 const bool invalidateFromMinPressure = this->updateMinPressure_();
1560 const bool invalidateFromHyst = this->updateHysteresis_();
1561 const bool invalidateFromMaxOilSat = this->updateMaxOilSaturation_();
1564 const bool invalidateDRDT = !first_step_after_restart && this->updateCompositionChangeLimits_();
1567 const bool invalidateIntensiveQuantities
1568 = invalidateFromMaxWaterSat || invalidateFromMinPressure || invalidateFromHyst || invalidateFromMaxOilSat || invalidateDRDT;
1569 if (invalidateIntensiveQuantities) {
1570 OPM_TIMEBLOCK(beginTimeStepInvalidateIntensiveQuantities);
1571 this->model().invalidateAndUpdateIntensiveQuantities(0);
1574 this->updateRockCompTransMultVal_();
1577 bool satfuncConsistencyRequirementsMet()
const
1579 if (
const auto nph = FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)
1580 + FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)
1581 + FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx);
1590 const auto numSamplePoints =
static_cast<std::size_t
>
1591 (Parameters::Get<Parameters::NumSatfuncConsistencySamplePoints>());
1593 auto sfuncConsistencyChecks =
1594 Satfunc::PhaseChecks::SatfuncConsistencyCheckManager<Scalar> {
1595 numSamplePoints, this->simulator().vanguard().eclState(),
1596 [&cmap = this->simulator().vanguard().cartesianIndexMapper()](
const int elemIdx)
1597 {
return cmap.cartesianIndex(elemIdx); }
1600 const auto ioRank = 0;
1601 const auto isIoRank = this->simulator().vanguard()
1602 .grid().comm().rank() == ioRank;
1604 sfuncConsistencyChecks.collectFailuresTo(ioRank)
1605 .run(this->simulator().vanguard().grid().leafGridView(),
1606 [&vg = this->simulator().vanguard(),
1607 &emap = this->simulator().model().elementMapper()]
1609 {
return vg.gridIdxToEquilGridIdx(emap.index(elem)); });
1611 using ViolationLevel =
typename Satfunc::PhaseChecks::
1612 SatfuncConsistencyCheckManager<Scalar>::ViolationLevel;
1614 auto reportFailures = [&sfuncConsistencyChecks]
1615 (
const ViolationLevel level)
1617 sfuncConsistencyChecks.reportFailures
1618 (level, [](std::string_view record)
1619 { OpmLog::info(std::string { record }); });
1622 if (sfuncConsistencyChecks.anyFailedStandardChecks()) {
1624 OpmLog::warning(
"Saturation Function "
1625 "End-point Consistency Problems");
1627 reportFailures(ViolationLevel::Standard);
1631 if (sfuncConsistencyChecks.anyFailedCriticalChecks()) {
1633 OpmLog::error(
"Saturation Function "
1634 "End-point Consistency Failures");
1636 reportFailures(ViolationLevel::Critical);
1649 FlowThresholdPressure<TypeTag> thresholdPressures_;
1651 std::vector<InitialFluidState> initialFluidStates_;
1653 bool enableEclOutput_;
1654 std::unique_ptr<EclWriterType> eclWriter_;
1656 const Scalar smallSaturationTolerance_ = 1.e-6;
1658 bool enableDamarisOutput_ = false ;
1659 std::unique_ptr<DamarisWriterType> damarisWriter_;
1661 MixingRateControls<FluidSystem> mixControls_;
1663 ActionHandler<Scalar> actionHandler_;
1665 ModuleParams moduleParams_;