31#ifndef OPM_FLOW_PROBLEM_BLACK_HPP
32#define OPM_FLOW_PROBLEM_BLACK_HPP
34#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
35#include <opm/material/fluidsystems/blackoilpvt/DryGasPvt.hpp>
36#include <opm/material/fluidsystems/blackoilpvt/WetGasPvt.hpp>
37#include <opm/material/fluidsystems/blackoilpvt/LiveOilPvt.hpp>
38#include <opm/material/fluidsystems/blackoilpvt/DeadOilPvt.hpp>
39#include <opm/material/fluidsystems/blackoilpvt/ConstantCompressibilityOilPvt.hpp>
40#include <opm/material/fluidsystems/blackoilpvt/ConstantCompressibilityWaterPvt.hpp>
44#include <opm/output/eclipse/EclipseIO.hpp>
72template <
class TypeTag>
157 DamarisWriterType::registerParameters();
170 simulator.vanguard().schedule(),
171 simulator.vanguard().actionState(),
172 simulator.vanguard().summaryState(),
174 simulator.vanguard().grid().comm())
179 const auto& vanguard = simulator.vanguard();
182 brineParams.template initFromState<enableBrine,
183 enableSaltPrecipitation>(vanguard.eclState());
186 DiffusionModule::initFromState(vanguard.eclState());
187 DispersionModule::initFromState(vanguard.eclState());
190 extboParams.template initFromState<enableExtbo>(vanguard.eclState());
194 foamParams.template initFromState<enableFoam>(vanguard.eclState());
198 micpParams.template initFromState<enableMICP>(vanguard.eclState());
202 polymerParams.template initFromState<enablePolymer, enablePolymerMolarWeight>(vanguard.eclState());
206 solventParams.template initFromState<enableSolvent>(vanguard.eclState(), vanguard.schedule());
210 eclWriter_ = std::make_unique<EclWriterType>(simulator);
214 damarisWriter_ = std::make_unique<DamarisWriterType>(simulator);
215 enableDamarisOutput_ = Parameters::Get<Parameters::EnableDamarisOutput>();
226 auto& simulator = this->simulator();
227 int episodeIdx = simulator.episodeIndex();
228 const auto& schedule = simulator.vanguard().schedule();
232 actionHandler_.evalUDQAssignments(episodeIdx, simulator.vanguard().udqState());
234 if (episodeIdx >= 0) {
235 const auto& oilVap = schedule[episodeIdx].oilvap();
236 if (oilVap.getType() == OilVaporizationProperties::OilVaporization::VAPPARS) {
237 FluidSystem::setVapPars(oilVap.vap1(), oilVap.vap2());
239 FluidSystem::setVapPars(0.0, 0.0);
243 ConvectiveMixingModule::beginEpisode(simulator.vanguard().eclState(), simulator.vanguard().schedule(), episodeIdx,
moduleParams_.convectiveMixingModuleParam);
254 FlowProblemType::finishInit();
255 auto& simulator = this->simulator();
257 auto finishTransmissibilities = [updated =
false,
this]()
mutable
259 if (updated) {
return; }
260 this->
transmissibilities_.finishInit([&vg = this->simulator().vanguard()](
const unsigned int it) {
261 return vg.gridIdxToEquilGridIdx(it);
271 if (simulator.vanguard().grid().comm().size() > 1) {
272 if (simulator.vanguard().grid().comm().rank() == 0)
273 eclWriter_->setTransmissibilities(&simulator.vanguard().globalTransmissibility());
275 finishTransmissibilities();
276 eclWriter_->setTransmissibilities(&simulator.problem().eclTransmissibilities());
279 std::function<
unsigned int(
unsigned int)> equilGridToGrid = [&simulator](
unsigned int i) {
280 return simulator.vanguard().gridEquilIdxToGridIdx(i);
282 eclWriter_->extractOutputTransAndNNC(equilGridToGrid);
284 simulator.vanguard().releaseGlobalTransmissibilities();
286 const auto& eclState = simulator.vanguard().eclState();
287 const auto& schedule = simulator.vanguard().schedule();
290 simulator.setStartTime(schedule.getStartTime());
291 simulator.setEndTime(schedule.simTime(schedule.size() - 1));
297 simulator.setEpisodeIndex(-1);
298 simulator.setEpisodeLength(0.0);
303 this->gravity_ = 0.0;
304 if (Parameters::Get<Parameters::EnableGravity>())
305 this->gravity_[dim - 1] = 9.80665;
306 if (!eclState.getInitConfig().hasGravity())
307 this->gravity_[dim - 1] = 0.0;
312 const auto& tuning = schedule[0].tuning();
321 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
322 FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
327 [&simulator](
const unsigned idx)
329 std::array<int,dim> coords;
330 simulator.vanguard().cartesianCoordinate(idx, coords);
331 for (auto& c : coords) {
344 finishTransmissibilities();
346 const auto& initconfig = eclState.getInitConfig();
348 if (initconfig.restartRequested())
357 if constexpr (getPropValue<TypeTag, Properties::EnablePolymer>()) {
358 const auto& vanguard = this->simulator().vanguard();
359 const auto& gridView = vanguard.gridView();
360 int numElements = gridView.size(0);
361 this->
polymer_.maxAdsorption.resize(numElements, 0.0);
370 this->
drift_.resize(this->model().numGridDof());
375 simulator.setTimeStepSize(0.0);
382 if (!initconfig.restartRequested()) {
383 simulator.startNextEpisode(schedule.seconds(1));
384 simulator.setEpisodeIndex(0);
385 simulator.setTimeStepIndex(0);
392 eclState.runspec().tabdims().getNumPVTTables());
400 FlowProblemType::endTimeStep();
402 const bool isSubStep = !this->simulator().episodeWillBeOver();
405 const auto& grid = this->simulator().vanguard().gridView().grid();
407 using GridType = std::remove_cv_t<std::remove_reference_t<
decltype(grid)>>;
408 constexpr bool isCpGrid = std::is_same_v<GridType, Dune::CpGrid>;
409 if (!isCpGrid || (grid.maxLevel() == 0)) {
410 this->eclWriter_->evalSummaryState(isSubStep);
414 OPM_TIMEBLOCK(applyActions);
416 const int episodeIdx = this->episodeIndex();
417 auto& simulator = this->simulator();
421 .applyActions(episodeIdx, simulator.time() + simulator.timeStepSize(),
422 [
this](
const bool global)
424 using TransUpdateQuantities = typename Vanguard::TransmissibilityType::TransUpdateQuantities;
425 this->transmissibilities_
426 .update(global, TransUpdateQuantities::All, [&vg = this->simulator().vanguard()]
427 (const unsigned int i)
429 return vg.gridIdxToEquilGridIdx(i);
435 if constexpr (enableMICP) {
436 auto& model = this->model();
437 const auto& residual = model.linearizer().residual();
439 for (
unsigned globalDofIdx = 0; globalDofIdx < residual.size(); ++globalDofIdx) {
440 auto& phi = this->referencePorosity_[1][globalDofIdx];
441 MICPModule::checkCloggingMICP(model, phi, globalDofIdx);
451 OPM_TIMEBLOCK(endEpisode);
452 const int episodeIdx = this->episodeIndex();
464 .evalUDQAssignments(episodeIdx, this->simulator().vanguard().udqState());
466 FlowProblemType::endEpisode();
470 if (enableEclOutput_){
471 eclWriter_->writeReports(timer);
482 FlowProblemType::writeOutput(verbose);
484 bool isSubStep = !this->simulator().episodeWillBeOver();
486 data::Solution localCellData = {};
490 if (enableDamarisOutput_) {
491 damarisWriter_->writeOutput(localCellData, isSubStep) ;
494 if (enableEclOutput_){
495 eclWriter_->writeOutput(std::move(localCellData), isSubStep);
501 OPM_TIMEBLOCK(finalizeOutput);
513 FlowProblemType::initialSolutionApplied();
518 this->thresholdPressures_.finishInit();
520 if (this->simulator().episodeIndex() == 0) {
521 eclWriter_->writeInitialFIPReport();
526 unsigned globalDofIdx,
527 unsigned timeIdx)
const override
529 this->aquiferModel_.addToSource(rate, globalDofIdx, timeIdx);
532 const auto& source = this->simulator().vanguard().schedule()[this->episodeIndex()].source();
533 std::array<int,3> ijk;
534 this->simulator().vanguard().cartesianCoordinate(globalDofIdx, ijk);
536 if (source.hasSource(ijk)) {
537 const int pvtRegionIdx = this->pvtRegionIndex(globalDofIdx);
538 static std::array<SourceComponent, 3> sc_map = {SourceComponent::WATER, SourceComponent::OIL, SourceComponent::GAS};
539 static std::array<int, 3> phidx_map = {FluidSystem::waterPhaseIdx, FluidSystem::oilPhaseIdx, FluidSystem::gasPhaseIdx};
540 static std::array<int, 3> cidx_map = {waterCompIdx, oilCompIdx, gasCompIdx};
542 for (
unsigned i = 0; i < phidx_map.size(); ++i) {
543 const auto phaseIdx = phidx_map[i];
544 const auto sourceComp = sc_map[i];
545 const auto compIdx = cidx_map[i];
546 if (!FluidSystem::phaseIsActive(phaseIdx)) {
549 Scalar mass_rate = source.rate({ijk, sourceComp}) / this->model().dofTotalVolume(globalDofIdx);
550 if constexpr (getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>()) {
551 mass_rate /= FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
553 rate[Indices::canonicalToActiveComponentIndex(compIdx)] += mass_rate;
556 if constexpr (enableSolvent) {
557 Scalar mass_rate = source.rate({ijk, SourceComponent::SOLVENT}) / this->model().dofTotalVolume(globalDofIdx);
558 if constexpr (getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>()) {
559 const auto& solventPvt = SolventModule::solventPvt();
560 mass_rate /= solventPvt.referenceDensity(pvtRegionIdx);
562 rate[Indices::contiSolventEqIdx] += mass_rate;
564 if constexpr (enablePolymer) {
565 rate[Indices::polymerConcentrationIdx] += source.rate({ijk, SourceComponent::POLYMER}) / this->model().dofTotalVolume(globalDofIdx);
567 if constexpr (enableEnergy) {
568 for (
unsigned i = 0; i < phidx_map.size(); ++i) {
569 const auto phaseIdx = phidx_map[i];
570 if (!FluidSystem::phaseIsActive(phaseIdx)) {
573 const auto sourceComp = sc_map[i];
574 if (source.hasHrate({ijk, sourceComp})) {
575 rate[Indices::contiEnergyEqIdx] += source.hrate({ijk, sourceComp}) / this->model().dofTotalVolume(globalDofIdx);
577 const auto& intQuants = this->simulator().model().intensiveQuantities(globalDofIdx, 0);
578 auto fs = intQuants.fluidState();
580 if (source.hasTemperature({ijk, sourceComp})) {
581 Scalar temperature = source.temperature({ijk, sourceComp});
582 fs.setTemperature(temperature);
584 const auto& h = FluidSystem::enthalpy(fs, phaseIdx, pvtRegionIdx);
585 Scalar mass_rate = source.rate({ijk, sourceComp})/ this->model().dofTotalVolume(globalDofIdx);
586 Scalar energy_rate = getValue(h)*mass_rate;
587 rate[Indices::contiEnergyEqIdx] += energy_rate;
595 if (this->enableDriftCompensation_) {
596 const auto& simulator = this->simulator();
597 const auto& model = this->model();
602 Scalar maxCompensation = model.newtonMethod().tolerance()/10;
603 Scalar poro = this->porosity(globalDofIdx, timeIdx);
604 Scalar dt = simulator.timeStepSize();
605 EqVector dofDriftRate = this->drift_[globalDofIdx];
606 dofDriftRate /= dt*model.dofTotalVolume(globalDofIdx);
609 for (
unsigned eqIdx = 0; eqIdx < numEq; ++ eqIdx) {
610 Scalar cnv = std::abs(dofDriftRate[eqIdx])*dt*model.eqWeight(globalDofIdx, eqIdx)/poro;
611 if (cnv > maxCompensation) {
612 dofDriftRate[eqIdx] *= maxCompensation/cnv;
616 for (
unsigned eqIdx = 0; eqIdx < numEq; ++ eqIdx)
617 rate[eqIdx] -= dofDriftRate[eqIdx];
626 template <
class LhsEval>
629 OPM_TIMEBLOCK_LOCAL(permFactTransMultiplier);
630 if (!enableSaltPrecipitation)
633 const auto& fs = intQuants.fluidState();
634 unsigned tableIdx = fs.pvtRegionIndex();
635 LhsEval porosityFactor = decay<LhsEval>(1. - fs.saltSaturation());
636 porosityFactor = min(porosityFactor, 1.0);
637 const auto& permfactTable = BrineModule::permfactTable(tableIdx);
638 return permfactTable.eval(porosityFactor,
true);
643 {
return initialFluidStates_[globalDofIdx]; }
646 {
return initialFluidStates_; }
649 {
return initialFluidStates_; }
652 {
return eclWriter_->eclIO(); }
655 {
return eclWriter_->setSubStepReport(report); }
658 {
return eclWriter_->setSimulationReport(report); }
662 OPM_TIMEBLOCK_LOCAL(boundaryFluidState);
663 const auto& bcprop = this->simulator().vanguard().schedule()[this->episodeIndex()].bcprop;
664 if (bcprop.size() > 0) {
665 FaceDir::DirEnum dir = FaceDir::FromIntersectionIndex(directionId);
669 if (this->bcindex_(dir)[globalDofIdx] == 0)
670 return initialFluidStates_[globalDofIdx];
672 const auto& bc = bcprop[this->bcindex_(dir)[globalDofIdx]];
673 if (bc.bctype == BCType::DIRICHLET )
675 InitialFluidState fluidState;
676 const int pvtRegionIdx = this->pvtRegionIndex(globalDofIdx);
677 fluidState.setPvtRegionIndex(pvtRegionIdx);
679 switch (bc.component) {
680 case BCComponent::OIL:
681 if (!FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx))
682 throw std::logic_error(
"oil is not active and you're trying to add oil BC");
684 fluidState.setSaturation(FluidSystem::oilPhaseIdx, 1.0);
686 case BCComponent::GAS:
687 if (!FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
688 throw std::logic_error(
"gas is not active and you're trying to add gas BC");
690 fluidState.setSaturation(FluidSystem::gasPhaseIdx, 1.0);
692 case BCComponent::WATER:
693 if (!FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx))
694 throw std::logic_error(
"water is not active and you're trying to add water BC");
696 fluidState.setSaturation(FluidSystem::waterPhaseIdx, 1.0);
698 case BCComponent::SOLVENT:
699 case BCComponent::POLYMER:
701 throw std::logic_error(
"you need to specify a valid component (OIL, WATER or GAS) when DIRICHLET type is set in BC");
703 fluidState.setTotalSaturation(1.0);
704 double pressure = initialFluidStates_[globalDofIdx].pressure(this->refPressurePhaseIdx_());
705 const auto pressure_input = bc.pressure;
706 if (pressure_input) {
707 pressure = *pressure_input;
710 std::array<Scalar, numPhases> pc = {0};
711 const auto& matParams = this->materialLawParams(globalDofIdx);
712 MaterialLaw::capillaryPressures(pc, matParams, fluidState);
713 Valgrind::CheckDefined(pressure);
714 Valgrind::CheckDefined(pc);
715 for (
unsigned activePhaseIdx = 0; activePhaseIdx < FluidSystem::numActivePhases(); ++activePhaseIdx) {
716 const auto phaseIdx = FluidSystem::activeToCanonicalPhaseIdx(activePhaseIdx);
718 fluidState.setPc(phaseIdx, pc[phaseIdx]);
719 if (Indices::oilEnabled)
720 fluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[oilPhaseIdx]));
721 else if (Indices::gasEnabled)
722 fluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[gasPhaseIdx]));
723 else if (Indices::waterEnabled)
725 fluidState.setPressure(phaseIdx, pressure);
728 double temperature = initialFluidStates_[globalDofIdx].temperature(0);
729 const auto temperature_input = bc.temperature;
730 if(temperature_input)
731 temperature = *temperature_input;
732 fluidState.setTemperature(temperature);
734 if (FluidSystem::enableDissolvedGas()) {
735 fluidState.setRs(0.0);
736 fluidState.setRv(0.0);
738 if (FluidSystem::enableDissolvedGasInWater()) {
739 fluidState.setRsw(0.0);
741 if (FluidSystem::enableVaporizedWater())
742 fluidState.setRvw(0.0);
744 for (
unsigned activePhaseIdx = 0; activePhaseIdx < FluidSystem::numActivePhases(); ++activePhaseIdx) {
745 const auto phaseIdx = FluidSystem::activeToCanonicalPhaseIdx(activePhaseIdx);
747 const auto& b = FluidSystem::inverseFormationVolumeFactor(fluidState, phaseIdx, pvtRegionIdx);
748 fluidState.setInvB(phaseIdx, b);
750 const auto& rho = FluidSystem::density(fluidState, phaseIdx, pvtRegionIdx);
751 fluidState.setDensity(phaseIdx, rho);
753 const auto& h = FluidSystem::enthalpy(fluidState, phaseIdx, pvtRegionIdx);
754 fluidState.setEnthalpy(phaseIdx, h);
757 fluidState.checkDefined();
761 return initialFluidStates_[globalDofIdx];
772 eclWriter_->mutableOutputModule().setCnvData(data);
781 return this->mixControls_.maxGasDissolutionFactor(timeIdx, globalDofIdx,
782 this->episodeIndex(),
783 this->pvtRegionIndex(globalDofIdx));
792 return this->mixControls_.maxOilVaporizationFactor(timeIdx, globalDofIdx,
793 this->episodeIndex(),
794 this->pvtRegionIndex(globalDofIdx));
807 int episodeIdx = this->episodeIndex();
808 return !this->mixControls_.drsdtActive(episodeIdx) &&
809 !this->mixControls_.drvdtActive(episodeIdx) &&
810 this->rockCompPoroMultWc_.empty() &&
811 this->rockCompPoroMult_.empty();
820 template <
class Context>
823 unsigned globalDofIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
825 values.setPvtRegionIndex(pvtRegionIndex(context, spaceIdx, timeIdx));
826 values.assignNaive(initialFluidStates_[globalDofIdx]);
828 SolventModule::assignPrimaryVars(values,
829 enableSolvent ? this->solventSaturation_[globalDofIdx] : 0.0,
830 enableSolvent ? this->solventRsw_[globalDofIdx] : 0.0);
832 if constexpr (enablePolymer)
833 values[Indices::polymerConcentrationIdx] = this->polymer_.concentration[globalDofIdx];
835 if constexpr (enablePolymerMolarWeight)
836 values[Indices::polymerMoleWeightIdx]= this->polymer_.moleWeight[globalDofIdx];
838 if constexpr (enableBrine) {
839 if (enableSaltPrecipitation && values.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) {
840 values[Indices::saltConcentrationIdx] = initialFluidStates_[globalDofIdx].saltSaturation();
843 values[Indices::saltConcentrationIdx] = initialFluidStates_[globalDofIdx].saltConcentration();
847 if constexpr (enableMICP){
848 values[Indices::microbialConcentrationIdx] = this->micp_.microbialConcentration[globalDofIdx];
849 values[Indices::oxygenConcentrationIdx]= this->micp_.oxygenConcentration[globalDofIdx];
850 values[Indices::ureaConcentrationIdx]= this->micp_.ureaConcentration[globalDofIdx];
851 values[Indices::calciteConcentrationIdx]= this->micp_.calciteConcentration[globalDofIdx];
852 values[Indices::biofilmConcentrationIdx]= this->micp_.biofilmConcentration[globalDofIdx];
855 values.checkDefined();
861 return this->mixControls_.drsdtcon(elemIdx, episodeIdx,
862 this->pvtRegionIndex(elemIdx));
870 template <
class Context>
872 const Context& context,
874 unsigned timeIdx)
const
876 OPM_TIMEBLOCK_LOCAL(eclProblemBoundary);
877 if (!context.intersection(spaceIdx).boundary())
880 if constexpr (!enableEnergy || !enableThermalFluxBoundaries)
888 unsigned interiorDofIdx = context.interiorScvIndex(spaceIdx, timeIdx);
889 unsigned globalDofIdx = context.globalSpaceIndex(interiorDofIdx, timeIdx);
890 values.setThermalFlow(context, spaceIdx, timeIdx, this->initialFluidStates_[globalDofIdx] );
893 if (this->nonTrivialBoundaryConditions()) {
894 unsigned indexInInside = context.intersection(spaceIdx).indexInInside();
895 unsigned interiorDofIdx = context.interiorScvIndex(spaceIdx, timeIdx);
896 unsigned globalDofIdx = context.globalSpaceIndex(interiorDofIdx, timeIdx);
897 unsigned pvtRegionIdx = pvtRegionIndex(context, spaceIdx, timeIdx);
898 const auto [type, massrate] = this->boundaryCondition(globalDofIdx, indexInInside);
899 if (type == BCType::THERMAL)
900 values.setThermalFlow(context, spaceIdx, timeIdx, this->boundaryFluidState(globalDofIdx, indexInInside));
901 else if (type == BCType::FREE || type == BCType::DIRICHLET)
902 values.setFreeFlow(context, spaceIdx, timeIdx, this->boundaryFluidState(globalDofIdx, indexInInside));
903 else if (type == BCType::RATE)
904 values.setMassRate(massrate, pvtRegionIdx);
913 {
return thresholdPressures_.thresholdPressure(elem1Idx, elem2Idx); }
916 {
return thresholdPressures_; }
919 {
return thresholdPressures_; }
923 return moduleParams_;
926 template<
class Serializer>
930 serializer(mixControls_);
931 serializer(*eclWriter_);
937 this->updateExplicitQuantities_(first_step_after_restart);
939 if constexpr (getPropValue<TypeTag, Properties::EnablePolymer>())
940 updateMaxPolymerAdsorption_();
942 mixControls_.updateExplicitQuantities(episodeIdx, timeStepSize);
948 this->updateProperty_(
"FlowProblemBlackoil::updateMaxPolymerAdsorption_() failed:",
951 this->updateMaxPolymerAdsorption_(compressedDofIdx,iq);
957 const Scalar pa = scalarValue(iq.polymerAdsorption());
958 auto& mpa = this->polymer_.maxAdsorption;
959 if (mpa[compressedDofIdx] < pa) {
960 mpa[compressedDofIdx] = pa;
969 std::vector<Scalar> sumInvB(numPhases, 0.0);
970 const auto& gridView = this->gridView();
972 for(
const auto& elem: elements(gridView, Dune::Partitions::interior)) {
973 elemCtx.updatePrimaryStencil(elem);
974 int elemIdx = elemCtx.globalSpaceIndex(0, 0);
975 const auto& dofFluidState = this->initialFluidStates_[elemIdx];
976 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
977 if (!FluidSystem::phaseIsActive(phaseIdx))
980 sumInvB[phaseIdx] += dofFluidState.invB(phaseIdx);
984 std::size_t numDof = this->model().numGridDof();
985 const auto& comm = this->simulator().vanguard().grid().comm();
986 comm.sum(sumInvB.data(),sumInvB.size());
987 Scalar numTotalDof = comm.sum(numDof);
989 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
990 if (!FluidSystem::phaseIsActive(phaseIdx))
993 Scalar avgB = numTotalDof / sumInvB[phaseIdx];
994 unsigned solventCompIdx = FluidSystem::solventComponentIndex(phaseIdx);
995 unsigned activeSolventCompIdx = Indices::canonicalToActiveComponentIndex(solventCompIdx);
996 this->model().setEqWeight(activeSolventCompIdx, avgB);
1003 OPM_TIMEBLOCK(updateCompositionChangeLimits);
1006 int episodeIdx = this->episodeIndex();
1007 std::array<bool,3> active{this->mixControls_.drsdtConvective(episodeIdx),
1008 this->mixControls_.drsdtActive(episodeIdx),
1009 this->mixControls_.drvdtActive(episodeIdx)};
1010 if (!active[0] && !active[1] && !active[2]) {
1014 this->updateProperty_(
"FlowProblemBlackoil::updateCompositionChangeLimits_()) failed:",
1015 [
this,episodeIdx,active](
unsigned compressedDofIdx,
1018 const DimMatrix& perm = this->intrinsicPermeability(compressedDofIdx);
1019 const Scalar distZ = active[0] ? this->simulator().vanguard().cellThickness(compressedDofIdx) : 0.0;
1020 const int pvtRegionIdx = this->pvtRegionIndex(compressedDofIdx);
1021 this->mixControls_.update(compressedDofIdx,
1024 this->gravity_[
dim - 1],
1038 if(this->simulator().vanguard().grid().maxLevel() > 0) {
1039 throw std::invalid_argument(
"Refined grids are not yet supported for restart ");
1043 auto& simulator = this->simulator();
1044 const auto& schedule = simulator.vanguard().schedule();
1045 const auto& eclState = simulator.vanguard().eclState();
1046 const auto& initconfig = eclState.getInitConfig();
1047 const int restart_step = initconfig.getRestartStep();
1049 simulator.setTime(schedule.seconds(restart_step));
1051 simulator.startNextEpisode(simulator.startTime() + simulator.time(),
1052 schedule.stepLength(restart_step));
1053 simulator.setEpisodeIndex(restart_step);
1055 this->eclWriter_->beginRestart();
1057 Scalar dt = std::min(this->eclWriter_->restartTimeStepSize(), simulator.episodeLength());
1058 simulator.setTimeStepSize(dt);
1060 std::size_t numElems = this->model().numGridDof();
1061 this->initialFluidStates_.resize(numElems);
1062 if constexpr (enableSolvent) {
1063 this->solventSaturation_.resize(numElems, 0.0);
1064 this->solventRsw_.resize(numElems, 0.0);
1067 if constexpr (enablePolymer)
1068 this->polymer_.concentration.resize(numElems, 0.0);
1070 if constexpr (enablePolymerMolarWeight) {
1071 const std::string msg {
"Support of the RESTART for polymer molecular weight "
1072 "is not implemented yet. The polymer weight value will be "
1073 "zero when RESTART begins"};
1074 OpmLog::warning(
"NO_POLYMW_RESTART", msg);
1075 this->polymer_.moleWeight.resize(numElems, 0.0);
1078 if constexpr (enableMICP) {
1079 this->micp_.resize(numElems);
1083 this->mixControls_.init(numElems, restart_step, eclState.runspec().tabdims().getNumPVTTables());
1085 for (std::size_t elemIdx = 0; elemIdx < numElems; ++elemIdx) {
1086 auto& elemFluidState = this->initialFluidStates_[elemIdx];
1087 elemFluidState.setPvtRegionIndex(pvtRegionIndex(elemIdx));
1088 this->eclWriter_->outputModule().initHysteresisParams(simulator, elemIdx);
1089 this->eclWriter_->outputModule().assignToFluidState(elemFluidState, elemIdx);
1098 auto ssol = enableSolvent
1099 ? this->eclWriter_->outputModule().getSolventSaturation(elemIdx)
1102 this->processRestartSaturations_(elemFluidState, ssol);
1104 if constexpr (enableSolvent) {
1105 this->solventSaturation_[elemIdx] = ssol;
1106 this->solventRsw_[elemIdx] = this->eclWriter_->outputModule().getSolventRsw(elemIdx);
1111 bool isThermal = eclState.getSimulationConfig().isThermal();
1112 bool needTemperature = (eclState.runspec().co2Storage() || eclState.runspec().h2Storage());
1113 if (!isThermal && needTemperature) {
1114 const auto& fp = simulator.vanguard().eclState().fieldProps();
1115 elemFluidState.setTemperature(fp.get_double(
"TEMPI")[elemIdx]);
1118 this->mixControls_.updateLastValues(elemIdx, elemFluidState.Rs(), elemFluidState.Rv());
1120 if constexpr (enablePolymer)
1121 this->polymer_.concentration[elemIdx] = this->eclWriter_->outputModule().getPolymerConcentration(elemIdx);
1122 if constexpr (enableMICP){
1123 this->micp_.microbialConcentration[elemIdx] = this->eclWriter_->outputModule().getMicrobialConcentration(elemIdx);
1124 this->micp_.oxygenConcentration[elemIdx] = this->eclWriter_->outputModule().getOxygenConcentration(elemIdx);
1125 this->micp_.ureaConcentration[elemIdx] = this->eclWriter_->outputModule().getUreaConcentration(elemIdx);
1126 this->micp_.biofilmConcentration[elemIdx] = this->eclWriter_->outputModule().getBiofilmConcentration(elemIdx);
1127 this->micp_.calciteConcentration[elemIdx] = this->eclWriter_->outputModule().getCalciteConcentration(elemIdx);
1132 const int episodeIdx = this->episodeIndex();
1133 this->mixControls_.updateMaxValues(episodeIdx, simulator.timeStepSize());
1138 auto& sol = this->model().solution(0);
1139 const auto& gridView = this->gridView();
1141 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
1142 elemCtx.updatePrimaryStencil(elem);
1143 int elemIdx = elemCtx.globalSpaceIndex(0, 0);
1144 this->initial(sol[elemIdx], elemCtx, 0, 0);
1152 this->model().syncOverlap();
1154 this->eclWriter_->endRestart();
1159 const auto& simulator = this->simulator();
1164 std::size_t numElems = this->model().numGridDof();
1165 this->initialFluidStates_.resize(numElems);
1166 for (std::size_t elemIdx = 0; elemIdx < numElems; ++elemIdx) {
1167 auto& elemFluidState = this->initialFluidStates_[elemIdx];
1174 const auto& simulator = this->simulator();
1175 const auto& vanguard = simulator.vanguard();
1176 const auto& eclState = vanguard.eclState();
1177 const auto& fp = eclState.fieldProps();
1178 bool has_swat = fp.has_double(
"SWAT");
1179 bool has_sgas = fp.has_double(
"SGAS");
1180 bool has_rs = fp.has_double(
"RS");
1181 bool has_rv = fp.has_double(
"RV");
1182 bool has_rvw = fp.has_double(
"RVW");
1183 bool has_pressure = fp.has_double(
"PRESSURE");
1184 bool has_salt = fp.has_double(
"SALT");
1185 bool has_saltp = fp.has_double(
"SALTP");
1188 if (Indices::numPhases > 1) {
1189 if (FluidSystem::phaseIsActive(waterPhaseIdx) && !has_swat)
1190 throw std::runtime_error(
"The ECL input file requires the presence of the SWAT keyword if "
1191 "the water phase is active");
1192 if (FluidSystem::phaseIsActive(gasPhaseIdx) && !has_sgas && FluidSystem::phaseIsActive(oilPhaseIdx))
1193 throw std::runtime_error(
"The ECL input file requires the presence of the SGAS keyword if "
1194 "the gas phase is active");
1197 throw std::runtime_error(
"The ECL input file requires the presence of the PRESSURE "
1198 "keyword if the model is initialized explicitly");
1199 if (FluidSystem::enableDissolvedGas() && !has_rs)
1200 throw std::runtime_error(
"The ECL input file requires the RS keyword to be present if"
1201 " dissolved gas is enabled");
1202 if (FluidSystem::enableVaporizedOil() && !has_rv)
1203 throw std::runtime_error(
"The ECL input file requires the RV keyword to be present if"
1204 " vaporized oil is enabled");
1205 if (FluidSystem::enableVaporizedWater() && !has_rvw)
1206 throw std::runtime_error(
"The ECL input file requires the RVW keyword to be present if"
1207 " vaporized water is enabled");
1208 if (enableBrine && !has_salt)
1209 throw std::runtime_error(
"The ECL input file requires the SALT keyword to be present if"
1210 " brine is enabled and the model is initialized explicitly");
1211 if (enableSaltPrecipitation && !has_saltp)
1212 throw std::runtime_error(
"The ECL input file requires the SALTP keyword to be present if"
1213 " salt precipitation is enabled and the model is initialized explicitly");
1215 std::size_t numDof = this->model().numGridDof();
1217 initialFluidStates_.resize(numDof);
1219 std::vector<double> waterSaturationData;
1220 std::vector<double> gasSaturationData;
1221 std::vector<double> pressureData;
1222 std::vector<double> rsData;
1223 std::vector<double> rvData;
1224 std::vector<double> rvwData;
1225 std::vector<double> tempiData;
1226 std::vector<double> saltData;
1227 std::vector<double> saltpData;
1229 if (FluidSystem::phaseIsActive(waterPhaseIdx) && Indices::numPhases > 1)
1230 waterSaturationData = fp.get_double(
"SWAT");
1232 waterSaturationData.resize(numDof);
1234 if (FluidSystem::phaseIsActive(gasPhaseIdx) && FluidSystem::phaseIsActive(oilPhaseIdx))
1235 gasSaturationData = fp.get_double(
"SGAS");
1237 gasSaturationData.resize(numDof);
1239 pressureData = fp.get_double(
"PRESSURE");
1240 if (FluidSystem::enableDissolvedGas())
1241 rsData = fp.get_double(
"RS");
1243 if (FluidSystem::enableVaporizedOil())
1244 rvData = fp.get_double(
"RV");
1246 if (FluidSystem::enableVaporizedWater())
1247 rvwData = fp.get_double(
"RVW");
1250 tempiData = fp.get_double(
"TEMPI");
1253 if constexpr (enableBrine)
1254 saltData = fp.get_double(
"SALT");
1257 if constexpr (enableSaltPrecipitation)
1258 saltpData = fp.get_double(
"SALTP");
1261 for (std::size_t dofIdx = 0; dofIdx < numDof; ++dofIdx) {
1262 auto& dofFluidState = initialFluidStates_[dofIdx];
1264 dofFluidState.setPvtRegionIndex(pvtRegionIndex(dofIdx));
1269 Scalar temperatureLoc = tempiData[dofIdx];
1270 if (!std::isfinite(temperatureLoc) || temperatureLoc <= 0)
1271 temperatureLoc = FluidSystem::surfaceTemperature;
1272 dofFluidState.setTemperature(temperatureLoc);
1277 if constexpr (enableBrine)
1278 dofFluidState.setSaltConcentration(saltData[dofIdx]);
1283 if constexpr (enableSaltPrecipitation)
1284 dofFluidState.setSaltSaturation(saltpData[dofIdx]);
1289 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx))
1290 dofFluidState.setSaturation(FluidSystem::waterPhaseIdx,
1291 waterSaturationData[dofIdx]);
1293 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)){
1294 if (!FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)){
1295 dofFluidState.setSaturation(FluidSystem::gasPhaseIdx,
1297 - waterSaturationData[dofIdx]);
1300 dofFluidState.setSaturation(FluidSystem::gasPhaseIdx,
1301 gasSaturationData[dofIdx]);
1303 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx))
1304 dofFluidState.setSaturation(FluidSystem::oilPhaseIdx,
1306 - waterSaturationData[dofIdx]
1307 - gasSaturationData[dofIdx]);
1312 Scalar pressure = pressureData[dofIdx];
1316 std::array<Scalar, numPhases> pc = {0};
1317 const auto& matParams = this->materialLawParams(dofIdx);
1318 MaterialLaw::capillaryPressures(pc, matParams, dofFluidState);
1319 Valgrind::CheckDefined(pressure);
1320 Valgrind::CheckDefined(pc);
1321 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1322 if (!FluidSystem::phaseIsActive(phaseIdx))
1325 if (Indices::oilEnabled)
1326 dofFluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[oilPhaseIdx]));
1327 else if (Indices::gasEnabled)
1328 dofFluidState.setPressure(phaseIdx, pressure + (pc[phaseIdx] - pc[gasPhaseIdx]));
1329 else if (Indices::waterEnabled)
1331 dofFluidState.setPressure(phaseIdx, pressure);
1334 if (FluidSystem::enableDissolvedGas())
1335 dofFluidState.setRs(rsData[dofIdx]);
1336 else if (Indices::gasEnabled && Indices::oilEnabled)
1337 dofFluidState.setRs(0.0);
1339 if (FluidSystem::enableVaporizedOil())
1340 dofFluidState.setRv(rvData[dofIdx]);
1341 else if (Indices::gasEnabled && Indices::oilEnabled)
1342 dofFluidState.setRv(0.0);
1344 if (FluidSystem::enableVaporizedWater())
1345 dofFluidState.setRvw(rvwData[dofIdx]);
1350 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1351 if (!FluidSystem::phaseIsActive(phaseIdx))
1354 const auto& b = FluidSystem::inverseFormationVolumeFactor(dofFluidState, phaseIdx, pvtRegionIndex(dofIdx));
1355 dofFluidState.setInvB(phaseIdx, b);
1357 const auto& rho = FluidSystem::density(dofFluidState, phaseIdx, pvtRegionIndex(dofIdx));
1358 dofFluidState.setDensity(phaseIdx, rho);
1369 const Scalar smallSaturationTolerance = 1.e-6;
1370 Scalar sumSaturation = 0.0;
1371 for (std::size_t phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1372 if (FluidSystem::phaseIsActive(phaseIdx)) {
1373 if (elemFluidState.saturation(phaseIdx) < smallSaturationTolerance)
1374 elemFluidState.setSaturation(phaseIdx, 0.0);
1376 sumSaturation += elemFluidState.saturation(phaseIdx);
1380 if constexpr (enableSolvent) {
1381 if (solventSaturation < smallSaturationTolerance)
1382 solventSaturation = 0.0;
1384 sumSaturation += solventSaturation;
1387 assert(sumSaturation > 0.0);
1389 for (std::size_t phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1390 if (FluidSystem::phaseIsActive(phaseIdx)) {
1391 const Scalar saturation = elemFluidState.saturation(phaseIdx) / sumSaturation;
1392 elemFluidState.setSaturation(phaseIdx, saturation);
1395 if constexpr (enableSolvent) {
1396 solventSaturation = solventSaturation / sumSaturation;
1402 FlowProblemType::readInitialCondition_();
1404 if constexpr (enableSolvent || enablePolymer || enablePolymerMolarWeight || enableMICP)
1405 this->readBlackoilExtentionsInitialConditions_(this->model().numGridDof(),
1408 enablePolymerMolarWeight,
1415 if constexpr (!enableSolvent)
1416 throw std::logic_error(
"solvent is disabled and you're trying to add solvent to BC");
1418 rate[Indices::solventSaturationIdx] = bc.rate;
1423 if constexpr (!enablePolymer)
1424 throw std::logic_error(
"polymer is disabled and you're trying to add polymer to BC");
1426 rate[Indices::polymerConcentrationIdx] = bc.rate;
1431 OPM_TIMEBLOCK(updateExplicitQuantities);
1432 const bool invalidateFromMaxWaterSat = this->updateMaxWaterSaturation_();
1433 const bool invalidateFromMinPressure = this->updateMinPressure_();
1436 const bool invalidateFromHyst = this->updateHysteresis_();
1437 const bool invalidateFromMaxOilSat = this->updateMaxOilSaturation_();
1440 const bool invalidateDRDT = !first_step_after_restart && this->updateCompositionChangeLimits_();
1443 const bool invalidateIntensiveQuantities
1444 = invalidateFromMaxWaterSat || invalidateFromMinPressure || invalidateFromHyst || invalidateFromMaxOilSat || invalidateDRDT;
1445 if (invalidateIntensiveQuantities) {
1446 OPM_TIMEBLOCK(beginTimeStepInvalidateIntensiveQuantities);
1447 this->model().invalidateAndUpdateIntensiveQuantities(0);
1450 this->updateRockCompTransMultVal_();
1460 bool enableDamarisOutput_ = false ;
1461 std::unique_ptr<DamarisWriterType> damarisWriter_;
Class handling Action support in simulator.
Definition: ActionHandler.hpp:51
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:50
static void setParams(BlackOilBrineParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilbrinemodules.hh:81
Definition: blackoilconvectivemixingmodule.hh:58
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: blackoildiffusionmodule.hh:48
Provides the auxiliary methods required for consideration of the dispersion equation.
Definition: blackoildispersionmodule.hh:48
Contains the high level supplements required to extend the black oil model.
Definition: blackoilextbomodules.hh:59
static void setParams(BlackOilExtboParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilextbomodules.hh:90
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition: blackoilfoammodules.hh:60
static void setParams(BlackOilFoamParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilfoammodules.hh:92
Contains the high level supplements required to extend the black oil model by MICP.
Definition: blackoilmicpmodules.hh:49
static void setParams(BlackOilMICPParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilmicpmodules.hh:83
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:54
static void setParams(BlackOilPolymerParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilpolymermodules.hh:88
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilsolventmodules.hh:58
static void setParams(BlackOilSolventParams< Scalar > &¶ms)
Set parameters.
Definition: blackoilsolventmodules.hh:90
Collects necessary output values and pass them to Damaris server processes.
Definition: DamarisWriter.hpp:90
Collects necessary output values and pass it to opm-common's ECL output.
Definition: EclWriter.hpp:101
static void registerParameters()
Definition: EclWriter.hpp:124
Computes the initial condition based on the EQUIL keyword from ECL.
Definition: EquilInitializer.hpp:58
const ScalarFluidState & initialFluidState(unsigned elemIdx) const
Return the initial thermodynamic state which should be used as the initial condition.
Definition: EquilInitializer.hpp:195
BlackOilFluidState< Scalar, FluidSystem, enableTemperature, enableEnergy, enableDissolution, enableVapwat, enableBrine, enableSaltPrecipitation, has_disgas_in_water, Indices::numPhases > ScalarFluidState
Definition: EquilInitializer.hpp:98
PolymerSolutionContainer< Scalar > polymer_
Definition: FlowGenericProblem.hpp:357
Scalar initialTimeStepSize_
Definition: FlowGenericProblem.hpp:368
bool enableTuning_
Definition: FlowGenericProblem.hpp:367
void readRockParameters_(const std::vector< Scalar > &cellCenterDepths, std::function< std::array< int, 3 >(const unsigned)> ijkIndex)
Definition: FlowGenericProblem_impl.hpp:138
std::vector< Scalar > maxOilSaturation_
Definition: FlowGenericProblem.hpp:358
void initFluidSystem_()
Definition: FlowGenericProblem_impl.hpp:482
Scalar maxTimeStepAfterWellEvent_
Definition: FlowGenericProblem.hpp:369
This problem simulates an input file given in the data format used by the commercial ECLiPSE simulato...
Definition: FlowProblemBlackoil.hpp:74
void updateExplicitQuantities_(int episodeIdx, int timeStepSize, const bool first_step_after_restart) override
Definition: FlowProblemBlackoil.hpp:935
bool updateMaxPolymerAdsorption_(unsigned compressedDofIdx, const IntensiveQuantities &iq)
Definition: FlowProblemBlackoil.hpp:955
void handlePolymerBC(const BCProp::BCFace &bc, RateVector &rate) const override
Definition: FlowProblemBlackoil.hpp:1421
void readInitialCondition_() override
Definition: FlowProblemBlackoil.hpp:1400
void readEquilInitialCondition_() override
Definition: FlowProblemBlackoil.hpp:1157
void handleSolventBC(const BCProp::BCFace &bc, RateVector &rate) const override
Definition: FlowProblemBlackoil.hpp:1413
Scalar maxGasDissolutionFactor(unsigned timeIdx, unsigned globalDofIdx) const
Returns the maximum value of the gas dissolution factor at the current time for a given degree of fre...
Definition: FlowProblemBlackoil.hpp:779
const std::vector< InitialFluidState > & initialFluidStates() const
Definition: FlowProblemBlackoil.hpp:648
void processRestartSaturations_(InitialFluidState &elemFluidState, Scalar &solventSaturation)
Definition: FlowProblemBlackoil.hpp:1365
std::vector< InitialFluidState > & initialFluidStates()
Definition: FlowProblemBlackoil.hpp:645
FlowProblemBlackoil(Simulator &simulator)
Definition: FlowProblemBlackoil.hpp:165
bool enableEclOutput_
Definition: FlowProblemBlackoil.hpp:1457
Scalar drsdtcon(unsigned elemIdx, int episodeIdx) const
Definition: FlowProblemBlackoil.hpp:859
Scalar maxOilVaporizationFactor(unsigned timeIdx, unsigned globalDofIdx) const
Returns the maximum value of the oil vaporization factor at the current time for a given degree of fr...
Definition: FlowProblemBlackoil.hpp:790
std::vector< InitialFluidState > initialFluidStates_
Definition: FlowProblemBlackoil.hpp:1455
void endTimeStep() override
Called by the simulator after each time integration.
Definition: FlowProblemBlackoil.hpp:398
void updateMaxPolymerAdsorption_()
Definition: FlowProblemBlackoil.hpp:945
const InitialFluidState & initialFluidState(unsigned globalDofIdx) const
Definition: FlowProblemBlackoil.hpp:642
void endEpisode() override
Called by the simulator after the end of an episode.
Definition: FlowProblemBlackoil.hpp:449
void setSubStepReport(const SimulatorReportSingle &report)
Definition: FlowProblemBlackoil.hpp:654
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Evaluate the initial value for a control volume.
Definition: FlowProblemBlackoil.hpp:821
void finishInit()
Called by the Opm::Simulator in order to initialize the problem.
Definition: FlowProblemBlackoil.hpp:249
const FlowThresholdPressure< TypeTag > & thresholdPressure() const
Definition: FlowProblemBlackoil.hpp:915
void finalizeOutput()
Definition: FlowProblemBlackoil.hpp:499
LhsEval permFactTransMultiplier(const IntensiveQuantities &intQuants) const
Calculate the transmissibility multiplier due to porosity reduction.
Definition: FlowProblemBlackoil.hpp:627
void boundary(BoundaryRateVector &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Evaluate the boundary conditions for a boundary segment.
Definition: FlowProblemBlackoil.hpp:871
InitialFluidState boundaryFluidState(unsigned globalDofIdx, const int directionId) const
Definition: FlowProblemBlackoil.hpp:660
std::unique_ptr< EclWriterType > eclWriter_
Definition: FlowProblemBlackoil.hpp:1458
void initialSolutionApplied() override
Callback used by the model to indicate that the initial solution has been determined for all degrees ...
Definition: FlowProblemBlackoil.hpp:511
const ModuleParams & moduleParams() const
Definition: FlowProblemBlackoil.hpp:921
void readEclRestartSolution_()
Definition: FlowProblemBlackoil.hpp:1035
FlowThresholdPressure< TypeTag > & thresholdPressure()
Definition: FlowProblemBlackoil.hpp:918
FlowThresholdPressure< TypeTag > thresholdPressures_
Definition: FlowProblemBlackoil.hpp:1453
void readExplicitInitialCondition_() override
Definition: FlowProblemBlackoil.hpp:1172
void beginEpisode() override
Called by the simulator before an episode begins.
Definition: FlowProblemBlackoil.hpp:222
bool recycleFirstIterationStorage() const
Return if the storage term of the first iteration is identical to the storage term for the solution o...
Definition: FlowProblemBlackoil.hpp:805
void serializeOp(Serializer &serializer)
Definition: FlowProblemBlackoil.hpp:927
MixingRateControls< FluidSystem > mixControls_
Definition: FlowProblemBlackoil.hpp:1463
void writeReports(const SimulatorTimer &timer)
Definition: FlowProblemBlackoil.hpp:469
ModuleParams moduleParams_
Definition: FlowProblemBlackoil.hpp:1467
void setSimulationReport(const SimulatorReport &report)
Definition: FlowProblemBlackoil.hpp:657
void addToSourceDense(RateVector &rate, unsigned globalDofIdx, unsigned timeIdx) const override
Definition: FlowProblemBlackoil.hpp:525
void computeAndSetEqWeights_()
Definition: FlowProblemBlackoil.hpp:967
ActionHandler< Scalar > actionHandler_
Definition: FlowProblemBlackoil.hpp:1465
void writeOutput(bool verbose=true)
Write the requested quantities of the current solution into the output files.
Definition: FlowProblemBlackoil.hpp:480
void updateExplicitQuantities_(const bool first_step_after_restart)
Definition: FlowProblemBlackoil.hpp:1429
static void registerParameters()
Registers all available parameters for the problem and the model.
Definition: FlowProblemBlackoil.hpp:151
const std::unique_ptr< EclWriterType > & eclWriter() const
Definition: FlowProblemBlackoil.hpp:765
const EclipseIO & eclIO() const
Definition: FlowProblemBlackoil.hpp:651
Scalar thresholdPressure(unsigned elem1Idx, unsigned elem2Idx) const
Definition: FlowProblemBlackoil.hpp:912
void setConvData(const std::vector< std::vector< int > > &data)
Definition: FlowProblemBlackoil.hpp:770
bool updateCompositionChangeLimits_()
Definition: FlowProblemBlackoil.hpp:1001
This problem simulates an input file given in the data format used by the commercial ECLiPSE simulato...
Definition: FlowProblem.hpp:91
@ numComponents
Definition: FlowProblem.hpp:113
@ enableSolvent
Definition: FlowProblem.hpp:127
@ enablePolymer
Definition: FlowProblem.hpp:124
unsigned pvtRegionIndex(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Returns the index of the relevant region for thermodynmic properties.
Definition: FlowProblem.hpp:813
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: FlowProblem.hpp:668
GetPropType< TypeTag, Properties::Vanguard > Vanguard
Definition: FlowProblem.hpp:104
@ enableMICP
Definition: FlowProblem.hpp:123
@ waterCompIdx
Definition: FlowProblem.hpp:139
@ enableDiffusion
Definition: FlowProblem.hpp:117
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: FlowProblem.hpp:98
GetPropType< TypeTag, Properties::EqVector > EqVector
Definition: FlowProblem.hpp:103
@ numEq
Definition: FlowProblem.hpp:111
@ enableEnergy
Definition: FlowProblem.hpp:119
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: FlowProblem.hpp:145
GlobalEqVector drift_
Definition: FlowProblem.hpp:1650
bool enableDriftCompensation_
Definition: FlowProblem.hpp:1649
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: FlowProblem.hpp:142
@ enableThermalFluxBoundaries
Definition: FlowProblem.hpp:129
@ dim
Definition: FlowProblem.hpp:107
Dune::FieldMatrix< Scalar, dimWorld, dimWorld > DimMatrix
Definition: FlowProblem.hpp:160
int episodeIndex() const
Definition: FlowProblem.hpp:297
GetPropType< TypeTag, Properties::Indices > Indices
Definition: FlowProblem.hpp:154
@ enableTemperature
Definition: FlowProblem.hpp:128
@ gasPhaseIdx
Definition: FlowProblem.hpp:131
@ oilCompIdx
Definition: FlowProblem.hpp:138
@ enableExperiments
Definition: FlowProblem.hpp:120
GetPropType< TypeTag, Properties::GlobalEqVector > GlobalEqVector
Definition: FlowProblem.hpp:102
@ dimWorld
Definition: FlowProblem.hpp:108
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: FlowProblem.hpp:143
@ enableSaltPrecipitation
Definition: FlowProblem.hpp:126
@ enablePolymerMolarWeight
Definition: FlowProblem.hpp:125
@ enableConvectiveMixing
Definition: FlowProblem.hpp:115
TracerModel tracerModel_
Definition: FlowProblem.hpp:1659
WellModel wellModel_
Definition: FlowProblem.hpp:1652
virtual void beginEpisode()
Called by the simulator before an episode begins.
Definition: FlowProblem.hpp:305
@ numPhases
Definition: FlowProblem.hpp:112
void readThermalParameters_()
Definition: FlowProblem.hpp:1349
@ enableFoam
Definition: FlowProblem.hpp:122
@ enableBrine
Definition: FlowProblem.hpp:116
@ waterPhaseIdx
Definition: FlowProblem.hpp:133
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: FlowProblem.hpp:155
bool enableVtkOutput_
Definition: FlowProblem.hpp:1655
GetPropType< TypeTag, Properties::GridView > GridView
Definition: FlowProblem.hpp:99
static void registerParameters()
Registers all available parameters for the problem and the model.
Definition: FlowProblem.hpp:178
void updatePffDofData_()
Definition: FlowProblem.hpp:1461
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition: FlowProblem.hpp:141
@ enableExtbo
Definition: FlowProblem.hpp:121
@ gasCompIdx
Definition: FlowProblem.hpp:137
void readBoundaryConditions_()
Definition: FlowProblem.hpp:1492
Vanguard::TransmissibilityType transmissibilities_
Definition: FlowProblem.hpp:1644
void writeOutput(bool verbose=true)
Write the requested quantities of the current solution into the output files.
Definition: FlowProblem.hpp:489
@ oilPhaseIdx
Definition: FlowProblem.hpp:132
@ enableDispersion
Definition: FlowProblem.hpp:118
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: FlowProblem.hpp:101
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: FlowProblem.hpp:151
void readMaterialParameters_()
Definition: FlowProblem.hpp:1315
This class calculates the threshold pressure for grid faces according to the Eclipse Reference Manual...
Definition: FlowThresholdPressure.hpp:59
Class handling mixing rate controls for a FlowProblemBlackoil.
Definition: MixingRateControls.hpp:46
Definition: SimulatorTimer.hpp:39
VTK output module for the tracer model's parameters.
Definition: VtkTracerModule.hpp:57
static void registerParameters()
Register all run-time parameters for the tracer VTK output module.
Definition: VtkTracerModule.hpp:83
@ NONE
Definition: DeferredLogger.hpp:46
static const int dim
Definition: structuredgridvanguard.hh:67
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:235
Struct holding the parameters for the BlackoilBrineModule class.
Definition: blackoilbrineparams.hpp:44
Struct holding the parameters for the BlackoilExtboModule class.
Definition: blackoilextboparams.hpp:49
Struct holding the parameters for the BlackoilFoamModule class.
Definition: blackoilfoamparams.hpp:46
Definition: blackoillocalresidualtpfa.hh:133
Struct holding the parameters for the BlackOilMICPModule class.
Definition: blackoilmicpparams.hpp:42
Struct holding the parameters for the BlackOilPolymerModule class.
Definition: blackoilpolymerparams.hpp:45
Struct holding the parameters for the BlackOilSolventModule class.
Definition: blackoilsolventparams.hpp:49
Definition: SimulatorReport.hpp:100
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34