28#ifndef OPM_ECL_WRITER_HPP
29#define OPM_ECL_WRITER_HPP
31#include <dune/grid/common/partitionset.hh>
33#include <opm/common/TimingMacros.hpp>
34#include <opm/common/OpmLog/OpmLog.hpp>
35#include <opm/input/eclipse/Schedule/RPTConfig.hpp>
37#include <opm/input/eclipse/Units/UnitSystem.hpp>
38#include <opm/input/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp>
40#include <opm/output/eclipse/Inplace.hpp>
41#include <opm/output/eclipse/RestartValue.hpp>
55#include <boost/date_time/posix_time/posix_time.hpp>
108template <
class TypeTag,
class OutputModule>
110 GetPropType<TypeTag, Properties::EquilGrid>,
111 GetPropType<TypeTag, Properties::GridView>,
112 GetPropType<TypeTag, Properties::ElementMapper>,
113 GetPropType<TypeTag, Properties::Scalar>>
123 using Element =
typename GridView::template Codim<0>::Entity;
125 using ElementIterator =
typename GridView::template Codim<0>::Iterator;
128 typedef Dune::MultipleCodimMultipleGeomTypeMapper< GridView > VertexMapper;
130 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
131 enum { enableMech = getPropValue<TypeTag, Properties::EnableMech>() };
132 enum { enableTemperature = getPropValue<TypeTag, Properties::EnableTemperature>() };
133 enum { enableSolvent = getPropValue<TypeTag, Properties::EnableSolvent>() };
139 OutputModule::registerParameters();
141 Parameters::Register<Parameters::EnableAsyncEclOutput>
142 (
"Write the ECL-formated results in a non-blocking way "
143 "(i.e., using a separate thread).");
144 Parameters::Register<Parameters::EnableEsmry>
145 (
"Write ESMRY file for fast loading of summary data.");
152 :
BaseType(simulator.vanguard().schedule(),
153 simulator.vanguard().eclState(),
154 simulator.vanguard().summaryConfig(),
155 simulator.vanguard().grid(),
156 ((simulator.vanguard().grid().comm().rank() == 0)
157 ? &simulator.vanguard().equilGrid()
159 simulator.vanguard().gridView(),
160 simulator.vanguard().cartesianIndexMapper(),
161 ((simulator.vanguard().grid().comm().rank() == 0)
162 ? &simulator.vanguard().equilCartesianIndexMapper()
164 Parameters::
Get<Parameters::EnableAsyncEclOutput>(),
165 Parameters::
Get<Parameters::EnableEsmry>())
166 , simulator_(simulator)
169 if (this->simulator_.vanguard().grid().comm().size() > 1) {
170 auto smryCfg = (this->simulator_.vanguard().grid().comm().rank() == 0)
171 ? this->
eclIO_->finalSummaryConfig()
174 eclBroadcast(this->simulator_.vanguard().grid().comm(), smryCfg);
176 this->outputModule_ = std::make_unique<OutputModule>
182 this->outputModule_ = std::make_unique<OutputModule>
183 (simulator, this->
eclIO_->finalSummaryConfig(), this->collectOnIORank_);
186 this->rank_ = this->simulator_.vanguard().grid().comm().rank();
188 this->simulator_.vanguard().eclState().computeFipRegionStatistics();
196 return simulator_.vanguard().equilGrid();
205 const int reportStepNum = simulator_.episodeIndex() + 1;
225 if (reportStepNum == 0)
228 const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
229 const Scalar totalCpuTime =
230 simulator_.executionTimer().realTimeElapsed() +
231 simulator_.setupTimer().realTimeElapsed() +
232 simulator_.vanguard().setupTime();
234 const auto localWellData = simulator_.problem().wellModel().wellData();
235 const auto localWBP = simulator_.problem().wellModel().wellBlockAveragePressures();
236 const auto localGroupAndNetworkData = simulator_.problem().wellModel()
237 .groupAndNetworkData(reportStepNum);
239 const auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
240 const auto localWellTestState = simulator_.problem().wellModel().wellTestState();
241 this->prepareLocalCellData(isSubStep, reportStepNum);
243 if (this->outputModule_->needInterfaceFluxes(isSubStep)) {
244 this->captureLocalFluxData();
250 std::map<std::pair<std::string,int>,
double> dummy;
252 outputModule_->getBlockData(),
256 localGroupAndNetworkData,
259 this->outputModule_->getInterRegFlows(),
266 if (! iregFlows.readIsConsistent()) {
267 throw std::runtime_error {
268 "Inconsistent inter-region flow "
269 "region set names in parallel"
277 this->simulator_.vanguard().grid().comm());
281 std::map<std::string, double> miscSummaryData;
282 std::map<std::string, std::vector<double>> regionData;
286 OPM_TIMEBLOCK(outputFipLogAndFipresvLog);
288 inplace = outputModule_->calc_inplace(miscSummaryData, regionData, simulator_.gridView().comm());
296 if (totalCpuTime != 0.0) {
297 miscSummaryData[
"TCPU"] = totalCpuTime;
326 : this->outputModule_->getBlockData();
330 : this->outputModule_->getInterRegFlows();
336 localGroupAndNetworkData,
342 this->outputModule_->initialInplace(),
344 this->summaryState(),
352 const auto& gridView = simulator_.vanguard().gridView();
356 this->outputModule_->
357 allocBuffers(num_interior, 0,
false,
false,
false);
360#pragma omp parallel for
362 for (
int dofIdx = 0; dofIdx < num_interior; ++dofIdx) {
363 const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, 0);
364 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
366 this->outputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
371 outputModule_->calc_initial_inplace(simulator_.gridView().comm());
374 const auto& fip = simulator_.vanguard().eclState().getEclipseConfig().fip();
375 if (fip.output(FIPConfig::OutputField::FIELD) ||
376 fip.output(FIPConfig::OutputField::RESV))
378 OPM_TIMEBLOCK(outputFipLogAndFipresvLog);
379 boost::posix_time::ptime start_time =
380 boost::posix_time::from_time_t(simulator_.vanguard().schedule().getStartTime());
383 inplace_ = outputModule_->initialInplace().value();
384 outputModule_->outputFipAndResvLog(inplace_, 0, 0.0, start_time,
385 false, simulator_.gridView().comm());
389 outputModule_->outputFipAndResvLogToCSV(0,
false, simulator_.gridView().comm());
406 const auto firstStep = this->initialStep();
410 const auto& rpt = this->
schedule_[simStep].rpt_config();
412 if (rpt.contains(
"WELSPECS") && (rpt.at(
"WELSPECS") > 0)) {
415 this->writeWellspecReport(timer);
423 if (rpt.contains(
"WELLS") && rpt.at(
"WELLS") > 0) {
424 this->writeWellflowReport(timer, simStep, rpt.at(
"WELLS"));
427 this->outputModule_->outputFipAndResvLog(this->inplace_,
432 simulator_.gridView().comm());
441 const int reportStepNum = simulator_.episodeIndex() + 1;
442 this->prepareLocalCellData(isSubStep, reportStepNum);
443 this->outputModule_->outputErrorLog(simulator_.gridView().comm());
446 auto localWellData = simulator_.problem().wellModel().wellData();
447 auto localGroupAndNetworkData = simulator_.problem().wellModel()
448 .groupAndNetworkData(reportStepNum);
450 auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
451 auto localWellTestState = simulator_.problem().wellModel().wellTestState();
453 const bool isFlowsn = this->outputModule_->getFlows().hasFlowsn();
454 auto flowsn = this->outputModule_->getFlows().getFlowsn();
456 const bool isFloresn = this->outputModule_->getFlows().hasFloresn();
457 auto floresn = this->outputModule_->getFlows().getFloresn();
459 if (! isSubStep || Parameters::Get<Parameters::EnableWriteAllSolutions>()) {
461 if (localCellData.empty()) {
462 this->outputModule_->assignToSolution(localCellData);
466 this->outputModule_->addRftDataToWells(localWellData,
468 simulator_.gridView().comm());
472 this->collectOnIORank_.doesNeedReordering())
480 this->outputModule_->getBlockData(),
481 this->outputModule_->getExtraBlockData(),
484 localGroupAndNetworkData,
494 this->outputModule_->assignGlobalFieldsToSolution(localCellData);
498 const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
499 const Scalar nextStepSize = simulator_.problem().nextTimeStepSize();
500 std::optional<int> timeStepIdx;
501 if (Parameters::Get<Parameters::EnableWriteAllSolutions>()) {
502 timeStepIdx = simulator_.timeStepIndex();
505 std::move(localCellData),
506 std::move(localWellData),
507 std::move(localGroupAndNetworkData),
508 std::move(localAquiferData),
509 std::move(localWellTestState),
512 this->summaryState(),
513 this->simulator_.problem().thresholdPressure().getRestartVector(),
514 curTime, nextStepSize,
515 Parameters::Get<Parameters::EclOutputDoublePrecision>(),
516 isFlowsn, std::move(flowsn),
517 isFloresn, std::move(floresn));
523 const auto enablePCHysteresis = simulator_.problem().materialLawManager()->enablePCHysteresis();
524 const auto enableNonWettingHysteresis = simulator_.problem().materialLawManager()->enableNonWettingHysteresis();
525 const auto enableWettingHysteresis = simulator_.problem().materialLawManager()->enableWettingHysteresis();
526 const auto oilActive = FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx);
527 const auto gasActive = FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx);
528 const auto waterActive = FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx);
529 const auto enableSwatinit = simulator_.vanguard().eclState().fieldProps().has_double(
"SWATINIT");
531 std::vector<RestartKey> solutionKeys {
532 {
"PRESSURE", UnitSystem::measure::pressure},
533 {
"SWAT", UnitSystem::measure::identity, waterActive},
534 {
"SGAS", UnitSystem::measure::identity, gasActive},
535 {
"TEMP", UnitSystem::measure::temperature, enableEnergy},
536 {
"SSOLVENT", UnitSystem::measure::identity, enableSolvent},
538 {
"RS", UnitSystem::measure::gas_oil_ratio, FluidSystem::enableDissolvedGas()},
539 {
"RV", UnitSystem::measure::oil_gas_ratio, FluidSystem::enableVaporizedOil()},
540 {
"RVW", UnitSystem::measure::oil_gas_ratio, FluidSystem::enableVaporizedWater()},
541 {
"RSW", UnitSystem::measure::gas_oil_ratio, FluidSystem::enableDissolvedGasInWater()},
543 {
"SGMAX", UnitSystem::measure::identity, enableNonWettingHysteresis && oilActive && gasActive},
544 {
"SHMAX", UnitSystem::measure::identity, enableWettingHysteresis && oilActive && gasActive},
546 {
"SOMAX", UnitSystem::measure::identity,
547 (enableNonWettingHysteresis && oilActive && waterActive)
548 || simulator_.problem().vapparsActive(simulator_.episodeIndex())},
550 {
"SOMIN", UnitSystem::measure::identity, enablePCHysteresis && oilActive && gasActive},
551 {
"SWHY1", UnitSystem::measure::identity, enablePCHysteresis && oilActive && waterActive},
552 {
"SWMAX", UnitSystem::measure::identity, enableWettingHysteresis && oilActive && waterActive},
554 {
"PPCW", UnitSystem::measure::pressure, enableSwatinit},
558 const auto& tracers = simulator_.vanguard().eclState().tracer();
560 for (
const auto& tracer : tracers) {
561 const auto enableSolTracer =
562 ((tracer.phase == Phase::GAS) && FluidSystem::enableDissolvedGas()) ||
563 ((tracer.phase == Phase::OIL) && FluidSystem::enableVaporizedOil());
565 solutionKeys.emplace_back(tracer.fname(), UnitSystem::measure::identity,
true);
566 solutionKeys.emplace_back(tracer.sname(), UnitSystem::measure::identity, enableSolTracer);
570 const auto& inputThpres = eclState().getSimulationConfig().getThresholdPressure();
571 const std::vector<RestartKey> extraKeys {
572 {
"OPMEXTRA", UnitSystem::measure::identity,
false},
573 {
"THRESHPR", UnitSystem::measure::pressure, inputThpres.active()},
576 const auto& gridView = this->simulator_.vanguard().gridView();
577 const auto numElements = gridView.size(0);
581 this->outputModule_->allocBuffers(numElements,
587 const auto restartSolution =
589 solutionKeys, gridView.comm(), 0);
591 if (!restartSolution.empty()) {
592 for (
auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
594 this->outputModule_->setRestart(restartSolution, elemIdx, globalIdx);
597 this->simulator_.problem().readSolutionFromOutputModule(0,
true);
598 ElementContext elemCtx(this->simulator_);
599 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
600 elemCtx.updatePrimaryStencil(elem);
601 elemCtx.updatePrimaryIntensiveQuantities(0);
603 this->outputModule_->updateFluidInPlace(elemCtx);
606 this->outputModule_->calc_initial_inplace(this->simulator_.gridView().comm());
614 const auto restartStepIdx = this->simulator_.vanguard()
615 .eclState().getInitConfig().getRestartStep();
617 this->outputModule_->allocBuffers(numElements,
625 const auto restartValues =
628 this->summaryState(),
629 solutionKeys, extraKeys, gridView.comm());
631 for (
auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
633 this->outputModule_->setRestart(restartValues.solution, elemIdx, globalIdx);
636 auto& tracer_model = simulator_.problem().tracerModel();
637 for (
int tracer_index = 0; tracer_index < tracer_model.numTracers(); ++tracer_index) {
640 const auto& free_tracer_name = tracer_model.fname(tracer_index);
641 const auto& free_tracer_solution = restartValues.solution
642 .template data<double>(free_tracer_name);
644 for (
auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
646 tracer_model.setFreeTracerConcentration
647 (tracer_index, elemIdx, free_tracer_solution[globalIdx]);
652 if ((tracer_model.phase(tracer_index) == Phase::GAS && FluidSystem::enableDissolvedGas()) ||
653 (tracer_model.phase(tracer_index) == Phase::OIL && FluidSystem::enableVaporizedOil()))
655 tracer_model.setEnableSolTracers(tracer_index,
true);
657 const auto& sol_tracer_name = tracer_model.sname(tracer_index);
658 const auto& sol_tracer_solution = restartValues.solution
659 .template data<double>(sol_tracer_name);
661 for (
auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
663 tracer_model.setSolTracerConcentration
664 (tracer_index, elemIdx, sol_tracer_solution[globalIdx]);
668 tracer_model.setEnableSolTracers(tracer_index,
false);
670 for (
auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
671 tracer_model.setSolTracerConcentration(tracer_index, elemIdx, 0.0);
676 if (inputThpres.active()) {
677 const_cast<Simulator&
>(this->simulator_)
678 .problem().thresholdPressure()
679 .setFromRestart(restartValues.getExtra(
"THRESHPR"));
682 restartTimeStepSize_ = restartValues.getExtra(
"OPMEXTRA")[0];
683 if (restartTimeStepSize_ <= 0) {
684 restartTimeStepSize_ = std::numeric_limits<double>::max();
688 this->simulator_.problem().wellModel()
689 .initFromRestartFile(restartValues);
691 if (!restartValues.aquifer.empty()) {
692 this->simulator_.problem().mutableAquiferModel()
693 .initFromRestart(restartValues.aquifer);
703 this->outputModule_->calc_initial_inplace(this->simulator_.gridView().comm());
706 if (this->outputModule_->initialInplace().has_value()) {
707 this->inplace_ = this->outputModule_->initialInplace().value();
713 {
return *outputModule_; }
716 {
return *outputModule_; }
719 {
return restartTimeStepSize_; }
721 template <
class Serializer>
724 serializer(*outputModule_);
728 static bool enableEclOutput_()
730 static bool enable = Parameters::Get<Parameters::EnableEclOutput>();
734 const EclipseState& eclState()
const
735 {
return simulator_.vanguard().eclState(); }
737 SummaryState& summaryState()
738 {
return simulator_.vanguard().summaryState(); }
740 Action::State& actionState()
741 {
return simulator_.vanguard().actionState(); }
744 {
return simulator_.vanguard().udqState(); }
746 const Schedule& schedule()
const
747 {
return simulator_.vanguard().schedule(); }
749 void prepareLocalCellData(
const bool isSubStep,
750 const int reportStepNum)
752 OPM_TIMEBLOCK(prepareLocalCellData);
754 if (this->outputModule_->localDataValid()) {
758 const auto& gridView = simulator_.vanguard().gridView();
763 this->outputModule_->
764 allocBuffers(num_interior, reportStepNum,
765 isSubStep && !Parameters::Get<Parameters::EnableWriteAllSolutions>(),
768 ElementContext elemCtx(simulator_);
773 OPM_TIMEBLOCK(prepareCellBasedData);
775 this->outputModule_->prepareDensityAccumulation();
776 this->outputModule_->setupExtractors(isSubStep, reportStepNum);
777 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
778 elemCtx.updatePrimaryStencil(elem);
779 elemCtx.updatePrimaryIntensiveQuantities(0);
781 this->outputModule_->processElement(elemCtx);
782 this->outputModule_->processElementBlockData(elemCtx);
784 this->outputModule_->clearExtractors();
786 this->outputModule_->accumulateDensityParallel();
790 OPM_TIMEBLOCK(prepareFluidInPlace);
793#pragma omp parallel for
795 for (
int dofIdx = 0; dofIdx < num_interior; ++dofIdx) {
796 const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, 0);
797 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
799 this->outputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
803 this->outputModule_->validateLocalData();
806 this->simulator_.vanguard().grid().comm());
809 void captureLocalFluxData()
811 OPM_TIMEBLOCK(captureLocalData);
813 const auto& gridView = this->simulator_.vanguard().gridView();
814 const auto timeIdx = 0u;
816 auto elemCtx = ElementContext { this->simulator_ };
818 const auto elemMapper = ElementMapper { gridView, Dune::mcmgElementLayout() };
819 const auto activeIndex = [&elemMapper](
const Element& e)
821 return elemMapper.index(e);
824 const auto cartesianIndex = [
this](
const int elemIndex)
826 return this->
cartMapper_.cartesianIndex(elemIndex);
829 this->outputModule_->initializeFluxData();
833 for (
const auto& elem : elements(gridView, Dune::Partitions::interiorBorder)) {
834 elemCtx.updateStencil(elem);
835 elemCtx.updateIntensiveQuantities(timeIdx);
836 elemCtx.updateExtensiveQuantities(timeIdx);
838 this->outputModule_->processFluxes(elemCtx, activeIndex, cartesianIndex);
842 this->simulator_.vanguard().grid().comm())
844 this->outputModule_->finalizeFluxData();
847 void writeWellspecReport(const SimulatorTimer& timer)
const
849 const auto changedWells = this->
schedule_
850 .changed_wells(timer.reportStepNum(), this->initialStep());
852 if (changedWells.empty()) {
856 this->outputModule_->outputWellspecReport(changedWells,
857 timer.reportStepNum(),
858 timer.simulationTimeElapsed(),
859 timer.currentDateTime());
862 void writeWellflowReport(
const SimulatorTimer& timer,
864 const int wellsRequest)
const
866 this->outputModule_->outputTimeStamp(
"WELLS",
867 timer.simulationTimeElapsed(),
868 timer.reportStepNum(),
869 timer.currentDateTime());
871 const auto wantConnData = wellsRequest > 1;
873 this->outputModule_->outputProdLog(simStep, wantConnData);
874 this->outputModule_->outputInjLog(simStep, wantConnData);
875 this->outputModule_->outputCumLog(simStep, wantConnData);
876 this->outputModule_->outputMSWLog(simStep);
879 int initialStep()
const
881 const auto& initConfig = this->eclState().cfg().init();
883 return initConfig.restartRequested()
884 ? initConfig.getRestartStep()
888 Simulator& simulator_;
889 std::unique_ptr<OutputModule> outputModule_;
890 Scalar restartTimeStepSize_;
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:192
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:158
Declares the properties required by the black oil model.
int localIdxToGlobalIdx(unsigned localIdx) const
Definition: CollectDataOnIORank_impl.hpp:1102
InterRegFlowMap & globalInterRegFlows()
Definition: CollectDataOnIORank.hpp:114
bool isParallel() const
Definition: CollectDataOnIORank.hpp:129
bool isIORank() const
Definition: CollectDataOnIORank.hpp:126
const std::map< std::pair< std::string, int >, double > & globalBlockData() const
Definition: CollectDataOnIORank.hpp:90
void collect(const data::Solution &localCellData, const std::map< std::pair< std::string, int >, double > &localBlockData, std::map< std::pair< std::string, int >, double > &localExtraBlockData, const data::Wells &localWellData, const data::WellBlockAveragePressures &localWBPData, const data::GroupAndNetworkValues &localGroupAndNetworkData, const data::Aquifers &localAquiferData, const WellTestState &localWellTestState, const InterRegFlowMap &interRegFlows, const std::array< FlowsData< double >, 3 > &localFlowsn, const std::array< FlowsData< double >, 3 > &localFloresn)
Definition: CollectDataOnIORank_impl.hpp:969
const data::Solution & globalCellData() const
Definition: CollectDataOnIORank.hpp:93
Definition: EclGenericWriter.hpp:65
void doWriteOutput(const int reportStepNum, const std::optional< int > timeStepNum, const bool isSubStep, data::Solution &&localCellData, data::Wells &&localWellData, data::GroupAndNetworkValues &&localGroupAndNetworkData, data::Aquifers &&localAquiferData, WellTestState &&localWTestState, const Action::State &actionState, const UDQState &udqState, const SummaryState &summaryState, const std::vector< GetPropType< TypeTag, Properties::Scalar > > &thresholdPressure, GetPropType< TypeTag, Properties::Scalar > curTime, GetPropType< TypeTag, Properties::Scalar > nextStepSize, bool doublePrecision, bool isFlowsn, std::array< FlowsData< double >, 3 > &&flowsn, bool isFloresn, std::array< FlowsData< double >, 3 > &&floresn)
Definition: EclGenericWriter_impl.hpp:533
CollectDataOnIORankType collectOnIORank_
Definition: EclGenericWriter.hpp:152
void evalSummary(int reportStepNum, GetPropType< TypeTag, Properties::Scalar > curTime, const data::Wells &localWellData, const data::WellBlockAveragePressures &localWBPData, const data::GroupAndNetworkValues &localGroupAndNetworkData, const std::map< int, data::AquiferData > &localAquiferData, const std::map< std::pair< std::string, int >, double > &blockData, const std::map< std::string, double > &miscSummaryData, const std::map< std::string, std::vector< double > > ®ionData, const Inplace &inplace, const std::optional< Inplace > &initialInPlace, const InterRegFlowMap &interRegFlows, SummaryState &summaryState, UDQState &udqState)
Definition: EclGenericWriter_impl.hpp:626
const Schedule & schedule_
Definition: EclGenericWriter.hpp:155
SimulatorReport simulation_report_
Definition: EclGenericWriter.hpp:165
SimulatorReportSingle sub_step_report_
Definition: EclGenericWriter.hpp:164
const Dune::CartesianIndexMapper< GetPropType< TypeTag, Properties::Grid > > & cartMapper_
Definition: EclGenericWriter.hpp:161
std::unique_ptr< EclipseIO > eclIO_
Definition: EclGenericWriter.hpp:157
Collects necessary output values and pass it to opm-common's ECL output.
Definition: EclWriter.hpp:114
OutputModule & mutableOutputModule() const
Definition: EclWriter.hpp:715
const OutputModule & outputModule() const
Definition: EclWriter.hpp:712
void evalSummaryState(bool isSubStep)
collect and pass data and pass it to eclIO writer
Definition: EclWriter.hpp:202
static void registerParameters()
Definition: EclWriter.hpp:137
void serializeOp(Serializer &serializer)
Definition: EclWriter.hpp:722
void writeInitialFIPReport()
Writes the initial FIP report as configured in RPTSOL.
Definition: EclWriter.hpp:350
void beginRestart()
Definition: EclWriter.hpp:521
EclWriter(Simulator &simulator)
Definition: EclWriter.hpp:151
void writeOutput(data::Solution &&localCellData, bool isSubStep)
Definition: EclWriter.hpp:437
void writeReports(const SimulatorTimer &timer)
Definition: EclWriter.hpp:392
void endRestart()
Definition: EclWriter.hpp:698
Scalar restartTimeStepSize() const
Definition: EclWriter.hpp:718
~EclWriter()
Definition: EclWriter.hpp:191
const EquilGrid & globalGrid() const
Definition: EclWriter.hpp:194
virtual int reportStepNum() const
Current report step number. This might differ from currentStepNum in case of sub stepping.
Definition: SimulatorTimerInterface.hpp:109
Definition: SimulatorTimer.hpp:39
virtual boost::posix_time::ptime currentDateTime() const
Return the current time as a posix time object.
double simulationTimeElapsed() const override
Defines the common properties required by the porous medium multi-phase models.
Definition: ActionHandler.hpp:34
Definition: blackoilnewtonmethodparams.hpp:31
auto Get(bool errorIfNotRegistered=true)
Retrieve a runtime parameter.
Definition: parametersystem.hpp:185
std::size_t countLocalInteriorCellsGridView(const GridView &gridView)
Get the number of local interior cells in a grid view.
Definition: countGlobalCells.hpp:47
Definition: blackoilboundaryratevector.hh:39
data::Solution loadParallelRestartSolution(const EclipseIO *eclIO, const std::vector< RestartKey > &solutionKeys, Parallel::Communication comm, const int step)
void eclBroadcast(Parallel::Communication, T &)
RestartValue loadParallelRestart(const EclipseIO *eclIO, Action::State &actionState, SummaryState &summaryState, const std::vector< RestartKey > &solutionKeys, const std::vector< RestartKey > &extraKeys, Parallel::Communication comm)
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
Definition: EclWriter.hpp:72
static constexpr bool value
Definition: EclWriter.hpp:72
Definition: EclWriter.hpp:69
static constexpr bool value
Definition: EclWriter.hpp:69
Definition: EclWriter.hpp:79
static constexpr bool value
Definition: EclWriter.hpp:79
Definition: EclWriter.hpp:76
static constexpr bool value
Definition: EclWriter.hpp:76
SimulatorReportSingle success
Definition: SimulatorReport.hpp:123
unsigned int min_linear_iterations
Definition: SimulatorReport.hpp:52
unsigned int total_newton_iterations
Definition: SimulatorReport.hpp:50
unsigned int max_linear_iterations
Definition: SimulatorReport.hpp:53
unsigned int total_linear_iterations
Definition: SimulatorReport.hpp:51