28#ifndef OPM_ECL_WRITER_HPP
29#define OPM_ECL_WRITER_HPP
31#include <dune/grid/common/partitionset.hh>
33#include <opm/common/OpmLog/OpmLog.hpp>
35#include <opm/input/eclipse/Units/UnitSystem.hpp>
37#include <opm/output/eclipse/RestartValue.hpp>
55template<
class TypeTag,
class MyTypeTag>
57 using type = UndefinedProperty;
59template<
class TypeTag,
class MyTypeTag>
61 using type = UndefinedProperty;
63template<
class TypeTag,
class MyTypeTag>
65 using type = UndefinedProperty;
67template<
class TypeTag,
class MyTypeTag>
69 using type = UndefinedProperty;
75namespace Action {
class State; }
94template <
class TypeTag>
96 GetPropType<TypeTag, Properties::EquilGrid>,
97 GetPropType<TypeTag, Properties::GridView>,
98 GetPropType<TypeTag, Properties::ElementMapper>,
99 GetPropType<TypeTag, Properties::Scalar>>
101 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
102 using Vanguard = GetPropType<TypeTag, Properties::Vanguard>;
103 using GridView = GetPropType<TypeTag, Properties::GridView>;
104 using Grid = GetPropType<TypeTag, Properties::Grid>;
105 using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
106 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
107 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
108 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
109 using Element =
typename GridView::template Codim<0>::Entity;
110 using ElementMapper = GetPropType<TypeTag, Properties::ElementMapper>;
111 using ElementIterator =
typename GridView::template Codim<0>::Iterator;
114 typedef Dune::MultipleCodimMultipleGeomTypeMapper< GridView > VertexMapper;
116 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
117 enum { enableMech = getPropValue<TypeTag, Properties::EnableMech>() };
118 enum { enableTemperature = getPropValue<TypeTag, Properties::EnableTemperature>() };
119 enum { enableSolvent = getPropValue<TypeTag, Properties::EnableSolvent>() };
127 Parameters::registerParam<TypeTag, Properties::EnableAsyncEclOutput>
128 (
"Write the ECL-formated results in a non-blocking way "
129 "(i.e., using a separate thread).");
130 Parameters::registerParam<TypeTag, Properties::EnableEsmry>
131 (
"Write ESMRY file for fast loading of summary data.");
138 :
BaseType(simulator.vanguard().schedule(),
139 simulator.vanguard().eclState(),
140 simulator.vanguard().summaryConfig(),
141 simulator.vanguard().grid(),
142 ((simulator.vanguard().grid().comm().rank() == 0)
143 ? &simulator.vanguard().equilGrid()
145 simulator.vanguard().gridView(),
146 simulator.vanguard().cartesianIndexMapper(),
147 ((simulator.vanguard().grid().comm().rank() == 0)
148 ? &simulator.vanguard().equilCartesianIndexMapper()
150 Parameters::get<TypeTag, Properties::EnableAsyncEclOutput>(),
151 Parameters::get<TypeTag, Properties::EnableEsmry>())
152 , simulator_(simulator)
155 if (this->simulator_.vanguard().grid().comm().size() > 1) {
156 auto smryCfg = (this->simulator_.vanguard().grid().comm().rank() == 0)
157 ? this->
eclIO_->finalSummaryConfig()
160 eclBroadcast(this->simulator_.vanguard().grid().comm(), smryCfg);
162 this->outputModule_ = std::make_unique<OutputBlackOilModule<TypeTag>>
168 this->outputModule_ = std::make_unique<OutputBlackOilModule<TypeTag>>
172 this->rank_ = this->simulator_.vanguard().grid().comm().rank();
174 this->simulator_.vanguard().eclState().computeFipRegionStatistics();
182 return simulator_.vanguard().equilGrid();
191 const int reportStepNum = simulator_.episodeIndex() + 1;
211 if (reportStepNum == 0)
214 const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
215 const Scalar totalCpuTime =
216 simulator_.executionTimer().realTimeElapsed() +
217 simulator_.setupTimer().realTimeElapsed() +
218 simulator_.vanguard().setupTime();
220 const auto localWellData = simulator_.problem().wellModel().wellData();
221 const auto localWBP = simulator_.problem().wellModel().wellBlockAveragePressures();
222 const auto localGroupAndNetworkData = simulator_.problem().wellModel()
223 .groupAndNetworkData(reportStepNum);
225 const auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
226 const auto localWellTestState = simulator_.problem().wellModel().wellTestState();
227 this->prepareLocalCellData(isSubStep, reportStepNum);
229 if (this->outputModule_->needInterfaceFluxes(isSubStep)) {
230 this->captureLocalFluxData();
237 outputModule_->getBlockData(),
240 localGroupAndNetworkData,
243 this->outputModule_->getInterRegFlows(),
250 if (! iregFlows.readIsConsistent()) {
251 throw std::runtime_error {
252 "Inconsistent inter-region flow "
253 "region set names in parallel"
261 this->simulator_.vanguard().grid().comm());
265 std::map<std::string, double> miscSummaryData;
266 std::map<std::string, std::vector<double>> regionData;
270 OPM_TIMEBLOCK(outputFipLogAndFipresvLog);
272 inplace = outputModule_->calc_inplace(miscSummaryData, regionData, simulator_.gridView().comm());
280 if (totalCpuTime != 0.0) {
281 miscSummaryData[
"TCPU"] = totalCpuTime;
310 : this->outputModule_->getBlockData();
314 : this->outputModule_->getInterRegFlows();
320 localGroupAndNetworkData,
326 this->outputModule_->initialInplace(),
328 this->summaryState(),
336 const auto& fip = simulator_.vanguard().eclState().getEclipseConfig().fip();
337 if (!fip.output(FIPConfig::OutputField::FIELD) &&
338 !fip.output(FIPConfig::OutputField::RESV)) {
342 const auto& gridView = simulator_.vanguard().gridView();
346 this->outputModule_->
347 allocBuffers(num_interior, 0,
false,
false,
false);
350#pragma omp parallel for
352 for (
int dofIdx = 0; dofIdx < num_interior; ++dofIdx) {
353 const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, 0);
354 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
356 this->outputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
359 std::map<std::string, double> miscSummaryData;
360 std::map<std::string, std::vector<double>> regionData;
363 OPM_TIMEBLOCK(outputFipLogAndFipresvLog);
365 boost::posix_time::ptime start_time = boost::posix_time::from_time_t(simulator_.vanguard().schedule().getStartTime());
367 inplace = outputModule_->calc_inplace(miscSummaryData, regionData, simulator_.gridView().comm());
372 outputModule_->outputFipAndResvLog(inplace_, 0, 0.0, start_time,
373 false, simulator_.gridView().comm());
382 const int reportStepNum = simulator_.episodeIndex() + 1;
383 this->prepareLocalCellData(isSubStep, reportStepNum);
384 this->outputModule_->outputErrorLog(simulator_.gridView().comm());
387 auto localWellData = simulator_.problem().wellModel().wellData();
388 auto localGroupAndNetworkData = simulator_.problem().wellModel()
389 .groupAndNetworkData(reportStepNum);
391 auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
392 auto localWellTestState = simulator_.problem().wellModel().wellTestState();
394 const bool isFlowsn = this->outputModule_->hasFlowsn();
395 auto flowsn = this->outputModule_->getFlowsn();
397 const bool isFloresn = this->outputModule_->hasFloresn();
398 auto floresn = this->outputModule_->getFloresn();
413 outputModule_->outputProdLog(reportStepNum);
414 outputModule_->outputInjLog(reportStepNum);
415 outputModule_->outputCumLog(reportStepNum);
420 if (localCellData.empty()) {
421 this->outputModule_->assignToSolution(localCellData);
425 this->outputModule_->addRftDataToWells(localWellData, reportStepNum);
429 this->collectOnIORank_.doesNeedReordering())
437 this->outputModule_->getBlockData(),
440 localGroupAndNetworkData,
450 this->outputModule_->assignGlobalFieldsToSolution(localCellData);
454 const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
455 const Scalar nextStepSize = simulator_.problem().nextTimeStepSize();
458 std::move(localCellData),
459 std::move(localWellData),
460 std::move(localGroupAndNetworkData),
461 std::move(localAquiferData),
462 std::move(localWellTestState),
465 this->summaryState(),
466 this->simulator_.problem().thresholdPressure().getRestartVector(),
467 curTime, nextStepSize,
468 Parameters::get<TypeTag, Properties::EclOutputDoublePrecision>(),
469 isFlowsn, std::move(flowsn),
470 isFloresn, std::move(floresn));
476 bool enableHysteresis = simulator_.problem().materialLawManager()->enableHysteresis();
477 bool enableSwatinit = simulator_.vanguard().eclState().fieldProps().has_double(
"SWATINIT");
478 bool opm_rst_file = Parameters::get<TypeTag, Properties::EnableOpmRstFile>();
479 bool read_temp = enableEnergy || (opm_rst_file && enableTemperature);
480 std::vector<RestartKey> solutionKeys{
481 {
"PRESSURE", UnitSystem::measure::pressure},
482 {
"SWAT", UnitSystem::measure::identity,
static_cast<bool>(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx))},
483 {
"SGAS", UnitSystem::measure::identity,
static_cast<bool>(FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))},
484 {
"TEMP" , UnitSystem::measure::temperature, read_temp},
485 {
"SSOLVENT" , UnitSystem::measure::identity, enableSolvent},
486 {
"RS", UnitSystem::measure::gas_oil_ratio, FluidSystem::enableDissolvedGas()},
487 {
"RV", UnitSystem::measure::oil_gas_ratio, FluidSystem::enableVaporizedOil()},
488 {
"RVW", UnitSystem::measure::oil_gas_ratio, FluidSystem::enableVaporizedWater()},
489 {
"SOMAX", UnitSystem::measure::identity, simulator_.problem().vapparsActive(simulator_.episodeIndex())},
490 {
"PCSWM_OW", UnitSystem::measure::identity, enableHysteresis},
491 {
"KRNSW_OW", UnitSystem::measure::identity, enableHysteresis},
492 {
"PCSWM_GO", UnitSystem::measure::identity, enableHysteresis},
493 {
"KRNSW_GO", UnitSystem::measure::identity, enableHysteresis},
494 {
"PPCW", UnitSystem::measure::pressure, enableSwatinit}
497 const auto& inputThpres = eclState().getSimulationConfig().getThresholdPressure();
498 std::vector<RestartKey> extraKeys = {{
"OPMEXTRA", UnitSystem::measure::identity,
false},
499 {
"THRESHPR", UnitSystem::measure::pressure, inputThpres.active()}};
502 const auto& tracers = simulator_.vanguard().eclState().tracer();
503 for (
const auto& tracer : tracers)
504 solutionKeys.emplace_back(tracer.fname(), UnitSystem::measure::identity,
true);
511 const auto& initconfig = simulator_.vanguard().eclState().getInitConfig();
512 int restartStepIdx = initconfig.getRestartStep();
514 const auto& gridView = simulator_.vanguard().gridView();
515 unsigned numElements = gridView.size(0);
516 outputModule_->allocBuffers(numElements, restartStepIdx,
false,
false,
true);
519 SummaryState& summaryState = simulator_.vanguard().summaryState();
520 Action::State& actionState = simulator_.vanguard().actionState();
522 gridView.grid().comm());
523 for (
unsigned elemIdx = 0; elemIdx < numElements; ++elemIdx) {
525 outputModule_->setRestart(restartValues.solution, elemIdx, globalIdx);
528 auto& tracer_model = simulator_.problem().tracerModel();
529 for (
int tracer_index = 0; tracer_index < tracer_model.numTracers(); tracer_index++) {
530 const auto& tracer_name = tracer_model.fname(tracer_index);
531 const auto& tracer_solution = restartValues.solution.template data<double>(tracer_name);
532 for (
unsigned elemIdx = 0; elemIdx < numElements; ++elemIdx) {
534 tracer_model.setTracerConcentration(tracer_index, globalIdx, tracer_solution[globalIdx]);
538 if (inputThpres.active()) {
539 Simulator& mutableSimulator =
const_cast<Simulator&
>(simulator_);
540 auto& thpres = mutableSimulator.problem().thresholdPressure();
541 const auto& thpresValues = restartValues.getExtra(
"THRESHPR");
542 thpres.setFromRestart(thpresValues);
544 restartTimeStepSize_ = restartValues.getExtra(
"OPMEXTRA")[0];
547 simulator_.problem().wellModel().initFromRestartFile(restartValues);
549 if (!restartValues.aquifer.empty())
550 simulator_.problem().mutableAquiferModel().initFromRestart(restartValues.aquifer);
558 {
return *outputModule_; }
561 {
return *outputModule_; }
564 {
return restartTimeStepSize_; }
566 template <
class Serializer>
569 serializer(*outputModule_);
573 static bool enableEclOutput_()
574 {
return Parameters::get<TypeTag, Properties::EnableEclOutput>(); }
576 const EclipseState& eclState()
const
577 {
return simulator_.vanguard().eclState(); }
579 SummaryState& summaryState()
580 {
return simulator_.vanguard().summaryState(); }
582 Action::State& actionState()
583 {
return simulator_.vanguard().actionState(); }
586 {
return simulator_.vanguard().udqState(); }
588 const Schedule& schedule()
const
589 {
return simulator_.vanguard().schedule(); }
591 void prepareLocalCellData(
const bool isSubStep,
592 const int reportStepNum)
594 OPM_TIMEBLOCK(prepareLocalCellData);
596 if (this->outputModule_->localDataValid()) {
600 const auto& gridView = simulator_.vanguard().gridView();
605 this->outputModule_->
606 allocBuffers(num_interior, reportStepNum,
607 isSubStep, log,
false);
609 ElementContext elemCtx(simulator_);
614 OPM_TIMEBLOCK(prepareCellBasedData);
616 this->outputModule_->prepareDensityAccumulation();
618 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
619 elemCtx.updatePrimaryStencil(elem);
620 elemCtx.updatePrimaryIntensiveQuantities(0);
622 this->outputModule_->processElement(elemCtx);
625 this->outputModule_->accumulateDensityParallel();
628 if constexpr (enableMech) {
629 if (simulator_.vanguard().eclState().runspec().mech()) {
630 OPM_TIMEBLOCK(prepareMechData);
631 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
632 elemCtx.updatePrimaryStencil(elem);
633 elemCtx.updatePrimaryIntensiveQuantities(0);
634 outputModule_->processElementMech(elemCtx);
639 if (! this->simulator_.model().linearizer().getFlowsInfo().empty()) {
640 OPM_TIMEBLOCK(prepareFlowsData);
641 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
642 elemCtx.updatePrimaryStencil(elem);
643 elemCtx.updatePrimaryIntensiveQuantities(0);
645 this->outputModule_->processElementFlows(elemCtx);
650 OPM_TIMEBLOCK(prepareBlockData);
651 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
652 elemCtx.updatePrimaryStencil(elem);
653 elemCtx.updatePrimaryIntensiveQuantities(0);
655 this->outputModule_->processElementBlockData(elemCtx);
660 OPM_TIMEBLOCK(prepareFluidInPlace);
663#pragma omp parallel for
665 for (
int dofIdx = 0; dofIdx < num_interior; ++dofIdx) {
666 const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, 0);
667 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
669 this->outputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
673 this->outputModule_->validateLocalData();
676 this->simulator_.vanguard().grid().comm());
679 void captureLocalFluxData()
681 OPM_TIMEBLOCK(captureLocalData);
683 const auto& gridView = this->simulator_.vanguard().gridView();
684 const auto timeIdx = 0u;
686 auto elemCtx = ElementContext { this->simulator_ };
688 const auto elemMapper = ElementMapper { gridView, Dune::mcmgElementLayout() };
689 const auto activeIndex = [&elemMapper](
const Element& e)
691 return elemMapper.index(e);
694 const auto cartesianIndex = [
this](
const int elemIndex)
696 return this->
cartMapper_.cartesianIndex(elemIndex);
699 this->outputModule_->initializeFluxData();
703 for (
const auto& elem : elements(gridView, Dune::Partitions::interiorBorder)) {
704 elemCtx.updateStencil(elem);
705 elemCtx.updateIntensiveQuantities(timeIdx);
706 elemCtx.updateExtensiveQuantities(timeIdx);
708 this->outputModule_->processFluxes(elemCtx, activeIndex, cartesianIndex);
712 this->simulator_.vanguard().grid().comm())
714 this->outputModule_->finalizeFluxData();
717 Simulator& simulator_;
718 std::unique_ptr<OutputBlackOilModule<TypeTag> > outputModule_;
719 Scalar restartTimeStepSize_;
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:172
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:138
void collect(const data::Solution &localCellData, const std::map< std::pair< std::string, int >, double > &localBlockData, 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
int localIdxToGlobalIdx(unsigned localIdx) const
Definition: CollectDataOnIORank_impl.hpp:1091
InterRegFlowMap & globalInterRegFlows()
Definition: CollectDataOnIORank.hpp:113
bool isParallel() const
Definition: CollectDataOnIORank.hpp:128
bool isIORank() const
Definition: CollectDataOnIORank.hpp:125
const std::map< std::pair< std::string, int >, double > & globalBlockData() const
Definition: CollectDataOnIORank.hpp:89
const data::Solution & globalCellData() const
Definition: CollectDataOnIORank.hpp:92
Definition: EclGenericWriter.hpp:65
CollectDataOnIORankType collectOnIORank_
Definition: EclGenericWriter.hpp:152
void doWriteOutput(const int reportStepNum, 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:510
SimulatorReport simulation_report_
Definition: EclGenericWriter.hpp:165
SimulatorReportSingle sub_step_report_
Definition: EclGenericWriter.hpp:164
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 Inplace &initialInPlace, const InterRegFlowMap &interRegFlows, SummaryState &summaryState, UDQState &udqState)
Definition: EclGenericWriter_impl.hpp:599
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-output.
Definition: EclWriter.hpp:100
void writeInitialFIPReport()
Writes the initial FIP report as configured in RPTSOL.
Definition: EclWriter.hpp:334
void endRestart()
Definition: EclWriter.hpp:554
const OutputBlackOilModule< TypeTag > & outputModule() const
Definition: EclWriter.hpp:557
Scalar restartTimeStepSize() const
Definition: EclWriter.hpp:563
~EclWriter()
Definition: EclWriter.hpp:177
EclWriter(Simulator &simulator)
Definition: EclWriter.hpp:137
static void registerParameters()
Definition: EclWriter.hpp:123
const EquilGrid & globalGrid() const
Definition: EclWriter.hpp:180
void beginRestart()
Definition: EclWriter.hpp:474
void writeOutput(data::Solution &&localCellData, const SimulatorTimer &timer, bool isSubStep)
Definition: EclWriter.hpp:378
OutputBlackOilModule< TypeTag > & mutableOutputModule() const
Definition: EclWriter.hpp:560
void evalSummaryState(bool isSubStep)
collect and pass data and pass it to eclIO writer
Definition: EclWriter.hpp:188
void serializeOp(Serializer &serializer)
Definition: EclWriter.hpp:567
Output module for the results black oil model writing in ECL binary format.
Definition: OutputBlackoilModule.hpp:116
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: OutputBlackoilModule.hpp:208
virtual int reportStepNum() const
Current report step number. This might differ from currentStepNum in case of sub stepping.
Definition: SimulatorTimerInterface.hpp:50
Definition: SimulatorTimer.hpp:39
virtual boost::posix_time::ptime currentDateTime() const
Return the current time as a posix time object.
double simulationTimeElapsed() const override
Definition: AluGridVanguard.hpp:57
std::size_t countLocalInteriorCellsGridView(const GridView &gridView)
Get the number of local interior cells in a grid view.
Definition: countGlobalCells.hpp:47
Definition: BlackoilPhases.hpp:27
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)
Definition: EclWriter.hpp:64
UndefinedProperty type
Definition: EclWriter.hpp:65
Definition: EclWriter.hpp:60
UndefinedProperty type
Definition: EclWriter.hpp:61
Definition: EclWriter.hpp:56
UndefinedProperty type
Definition: EclWriter.hpp:57
Definition: EclWriter.hpp:68
UndefinedProperty type
Definition: EclWriter.hpp:69
SimulatorReportSingle success
Definition: SimulatorReport.hpp:101
unsigned int min_linear_iterations
Definition: SimulatorReport.hpp:51
unsigned int total_newton_iterations
Definition: SimulatorReport.hpp:49
unsigned int max_linear_iterations
Definition: SimulatorReport.hpp:52
unsigned int total_linear_iterations
Definition: SimulatorReport.hpp:50