My Project
Loading...
Searching...
No Matches
EclWriter.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef OPM_ECL_WRITER_HPP
29#define OPM_ECL_WRITER_HPP
30
31#include <dune/grid/common/partitionset.hh>
32
33#include <opm/common/OpmLog/OpmLog.hpp>
34#include <opm/input/eclipse/Schedule/RPTConfig.hpp>
35
36#include <opm/input/eclipse/Units/UnitSystem.hpp>
37
38#include <opm/output/eclipse/RestartValue.hpp>
39
40#include <opm/simulators/flow/CollectDataOnIORank.hpp>
41#include <opm/simulators/flow/countGlobalCells.hpp>
45#include <opm/simulators/timestepping/SimulatorTimer.hpp>
46#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
47#include <opm/simulators/utils/ParallelRestart.hpp>
48#include <opm/simulators/utils/ParallelSerialization.hpp>
49
50#include <limits>
51#include <map>
52#include <memory>
53#include <optional>
54#include <stdexcept>
55#include <string>
56#include <utility>
57#include <vector>
58
59namespace Opm::Parameters {
60
61// enable the ECL output by default
62struct EnableEclOutput { static constexpr bool value = true; };
63
64// If available, write the ECL output in a non-blocking manner
65struct EnableAsyncEclOutput { static constexpr bool value = true; };
66
67// By default, use single precision for the ECL formated results
68struct EclOutputDoublePrecision { static constexpr bool value = false; };
69
70// Write all solutions for visualization, not just the ones for the
71// report steps...
72struct EnableWriteAllSolutions { static constexpr bool value = false; };
73
74// Write ESMRY file for fast loading of summary data
75struct EnableEsmry { static constexpr bool value = false; };
76
77} // namespace Opm::Parameters
78
79namespace Opm::Action {
80 class State;
81} // namespace Opm::Action
82
83namespace Opm {
84 class EclipseIO;
85 class UDQState;
86} // namespace Opm
87
88namespace Opm {
104template <class TypeTag>
105class EclWriter : public EclGenericWriter<GetPropType<TypeTag, Properties::Grid>,
106 GetPropType<TypeTag, Properties::EquilGrid>,
107 GetPropType<TypeTag, Properties::GridView>,
108 GetPropType<TypeTag, Properties::ElementMapper>,
109 GetPropType<TypeTag, Properties::Scalar>>
110{
119 using Element = typename GridView::template Codim<0>::Entity;
121 using ElementIterator = typename GridView::template Codim<0>::Iterator;
123
124 typedef Dune::MultipleCodimMultipleGeomTypeMapper< GridView > VertexMapper;
125
126 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
127 enum { enableMech = getPropValue<TypeTag, Properties::EnableMech>() };
128 enum { enableTemperature = getPropValue<TypeTag, Properties::EnableTemperature>() };
129 enum { enableSolvent = getPropValue<TypeTag, Properties::EnableSolvent>() };
130
131public:
132
133 static void registerParameters()
134 {
136
137 Parameters::Register<Parameters::EnableAsyncEclOutput>
138 ("Write the ECL-formated results in a non-blocking way "
139 "(i.e., using a separate thread).");
140 Parameters::Register<Parameters::EnableEsmry>
141 ("Write ESMRY file for fast loading of summary data.");
142 }
143
144 // The Simulator object should preferably have been const - the
145 // only reason that is not the case is due to the SummaryState
146 // object owned deep down by the vanguard.
147 EclWriter(Simulator& simulator)
148 : BaseType(simulator.vanguard().schedule(),
149 simulator.vanguard().eclState(),
150 simulator.vanguard().summaryConfig(),
151 simulator.vanguard().grid(),
152 ((simulator.vanguard().grid().comm().rank() == 0)
153 ? &simulator.vanguard().equilGrid()
154 : nullptr),
155 simulator.vanguard().gridView(),
156 simulator.vanguard().cartesianIndexMapper(),
157 ((simulator.vanguard().grid().comm().rank() == 0)
158 ? &simulator.vanguard().equilCartesianIndexMapper()
159 : nullptr),
160 Parameters::Get<Parameters::EnableAsyncEclOutput>(),
161 Parameters::Get<Parameters::EnableEsmry>())
162 , simulator_(simulator)
163 {
164#if HAVE_MPI
165 if (this->simulator_.vanguard().grid().comm().size() > 1) {
166 auto smryCfg = (this->simulator_.vanguard().grid().comm().rank() == 0)
167 ? this->eclIO_->finalSummaryConfig()
168 : SummaryConfig{};
169
170 eclBroadcast(this->simulator_.vanguard().grid().comm(), smryCfg);
171
172 this->outputModule_ = std::make_unique<OutputBlackOilModule<TypeTag>>
173 (simulator, smryCfg, this->collectOnIORank_);
174 }
175 else
176#endif
177 {
178 this->outputModule_ = std::make_unique<OutputBlackOilModule<TypeTag>>
179 (simulator, this->eclIO_->finalSummaryConfig(), this->collectOnIORank_);
180 }
181
182 this->rank_ = this->simulator_.vanguard().grid().comm().rank();
183
184 this->simulator_.vanguard().eclState().computeFipRegionStatistics();
185 }
186
187 ~EclWriter()
188 {}
189
190 const EquilGrid& globalGrid() const
191 {
192 return simulator_.vanguard().equilGrid();
193 }
194
199 {
201 const int reportStepNum = simulator_.episodeIndex() + 1;
202
203 /*
204 The summary data is not evaluated for timestep 0, that is
205 implemented with a:
206
207 if (time_step == 0)
208 return;
209
210 check somewhere in the summary code. When the summary code was
211 split in separate methods Summary::eval() and
212 Summary::add_timestep() it was necessary to pull this test out
213 here to ensure that the well and group related keywords in the
214 restart file, like XWEL and XGRP were "correct" also in the
215 initial report step.
216
217 "Correct" in this context means unchanged behavior, might very
218 well be more correct to actually remove this if test.
219 */
220
221 if (reportStepNum == 0)
222 return;
223
224 const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
225 const Scalar totalCpuTime =
226 simulator_.executionTimer().realTimeElapsed() +
227 simulator_.setupTimer().realTimeElapsed() +
228 simulator_.vanguard().setupTime();
229
230 const auto localWellData = simulator_.problem().wellModel().wellData();
231 const auto localWBP = simulator_.problem().wellModel().wellBlockAveragePressures();
232 const auto localGroupAndNetworkData = simulator_.problem().wellModel()
233 .groupAndNetworkData(reportStepNum);
234
235 const auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
236 const auto localWellTestState = simulator_.problem().wellModel().wellTestState();
237 this->prepareLocalCellData(isSubStep, reportStepNum);
238
239 if (this->outputModule_->needInterfaceFluxes(isSubStep)) {
240 this->captureLocalFluxData();
241 }
242
243 if (this->collectOnIORank_.isParallel()) {
244 OPM_BEGIN_PARALLEL_TRY_CATCH()
245
246 this->collectOnIORank_.collect({},
247 outputModule_->getBlockData(),
249 localWBP,
253 this->outputModule_->getInterRegFlows(),
254 {},
255 {});
256
257 if (this->collectOnIORank_.isIORank()) {
258 auto& iregFlows = this->collectOnIORank_.globalInterRegFlows();
259
260 if (! iregFlows.readIsConsistent()) {
261 throw std::runtime_error {
262 "Inconsistent inter-region flow "
263 "region set names in parallel"
264 };
265 }
266
268 }
269
270 OPM_END_PARALLEL_TRY_CATCH("Collect to I/O rank: ",
271 this->simulator_.vanguard().grid().comm());
272 }
273
274
275 std::map<std::string, double> miscSummaryData;
276 std::map<std::string, std::vector<double>> regionData;
278
279 {
281
282 inplace = outputModule_->calc_inplace(miscSummaryData, regionData, simulator_.gridView().comm());
283
284 if (this->collectOnIORank_.isIORank()){
285 inplace_ = inplace;
286 }
287 }
288
289 // Add TCPU
290 if (totalCpuTime != 0.0) {
292 }
293 if (this->sub_step_report_.total_newton_iterations != 0) {
294 miscSummaryData["NEWTON"] = this->sub_step_report_.total_newton_iterations;
295 }
296 if (this->sub_step_report_.total_linear_iterations != 0) {
297 miscSummaryData["MLINEARS"] = this->sub_step_report_.total_linear_iterations;
298 }
299 if (this->sub_step_report_.total_newton_iterations != 0) {
300 miscSummaryData["NLINEARS"] = static_cast<float>(this->sub_step_report_.total_linear_iterations) / this->sub_step_report_.total_newton_iterations;
301 }
302 if (this->sub_step_report_.min_linear_iterations != std::numeric_limits<unsigned int>::max()) {
303 miscSummaryData["NLINSMIN"] = this->sub_step_report_.min_linear_iterations;
304 }
305 if (this->sub_step_report_.max_linear_iterations != 0) {
306 miscSummaryData["NLINSMAX"] = this->sub_step_report_.max_linear_iterations;
307 }
308 if (this->simulation_report_.success.total_newton_iterations != 0) {
309 miscSummaryData["MSUMLINS"] = this->simulation_report_.success.total_linear_iterations;
310 }
311 if (this->simulation_report_.success.total_newton_iterations != 0) {
312 miscSummaryData["MSUMNEWT"] = this->simulation_report_.success.total_newton_iterations;
313 }
314
315 {
316 OPM_TIMEBLOCK(evalSummary);
317
318 const auto& blockData = this->collectOnIORank_.isParallel()
319 ? this->collectOnIORank_.globalBlockData()
320 : this->outputModule_->getBlockData();
321
322 const auto& interRegFlows = this->collectOnIORank_.isParallel()
323 ? this->collectOnIORank_.globalInterRegFlows()
324 : this->outputModule_->getInterRegFlows();
325
326 this->evalSummary(reportStepNum,
327 curTime,
329 localWBP,
332 blockData,
335 inplace,
336 this->outputModule_->initialInplace(),
338 this->summaryState(),
339 this->udqState());
340 }
341 }
342
345 {
346 const auto& fip = simulator_.vanguard().eclState().getEclipseConfig().fip();
347 if (!fip.output(FIPConfig::OutputField::FIELD) &&
348 !fip.output(FIPConfig::OutputField::RESV)) {
349 return;
350 }
351
352 const auto& gridView = simulator_.vanguard().gridView();
353 const int num_interior = detail::
354 countLocalInteriorCellsGridView(gridView);
355
356 this->outputModule_->
357 allocBuffers(num_interior, 0, false, false, /*isRestart*/ false);
358
359#ifdef _OPENMP
360#pragma omp parallel for
361#endif
362 for (int dofIdx = 0; dofIdx < num_interior; ++dofIdx) {
363 const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, /*timeIdx=*/0);
364 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
365
366 this->outputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
367 }
368
369 std::map<std::string, double> miscSummaryData;
370 std::map<std::string, std::vector<double>> regionData;
372 {
374
375 boost::posix_time::ptime start_time = boost::posix_time::from_time_t(simulator_.vanguard().schedule().getStartTime());
376
377 inplace = outputModule_->calc_inplace(miscSummaryData, regionData, simulator_.gridView().comm());
378
379 if (this->collectOnIORank_.isIORank()){
380 inplace_ = inplace;
381
382 outputModule_->outputFipAndResvLog(inplace_, 0, 0.0, start_time,
383 false, simulator_.gridView().comm());
384 }
385 }
386 }
387
388 void writeReports(const SimulatorTimer& timer) {
389 auto rstep = timer.reportStepNum();
390
391 if ((rstep > 0) && (this->collectOnIORank_.isIORank())){
392
393 const auto& rpt = this->schedule_[rstep-1].rpt_config.get();
394 if (rpt.contains("WELLS") && rpt.at("WELLS") > 0) {
395 outputModule_->outputTimeStamp("WELLS", timer.simulationTimeElapsed(), rstep, timer.currentDateTime());
396 outputModule_->outputProdLog(rstep-1);
397 outputModule_->outputInjLog(rstep-1);
398 outputModule_->outputCumLog(rstep-1);
399 }
400
401 outputModule_->outputFipAndResvLog(inplace_, rstep, timer.simulationTimeElapsed(),
402 timer.currentDateTime(), false, simulator_.gridView().comm());
403
404
405 OpmLog::note(""); // Blank line after all reports.
406 }
407 }
408
409 void writeOutput(data::Solution&& localCellData, bool isSubStep)
410 {
411 OPM_TIMEBLOCK(writeOutput);
412
413 const int reportStepNum = simulator_.episodeIndex() + 1;
414 this->prepareLocalCellData(isSubStep, reportStepNum);
415 this->outputModule_->outputErrorLog(simulator_.gridView().comm());
416
417 // output using eclWriter if enabled
418 auto localWellData = simulator_.problem().wellModel().wellData();
419 auto localGroupAndNetworkData = simulator_.problem().wellModel()
420 .groupAndNetworkData(reportStepNum);
421
422 auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
423 auto localWellTestState = simulator_.problem().wellModel().wellTestState();
424
425 const bool isFlowsn = this->outputModule_->hasFlowsn();
426 auto flowsn = this->outputModule_->getFlowsn();
427
428 const bool isFloresn = this->outputModule_->hasFloresn();
429 auto floresn = this->outputModule_->getFloresn();
430
431 // data::Solution localCellData = {};
432 if (! isSubStep || Parameters::Get<Parameters::EnableWriteAllSolutions>()) {
433
434 if (localCellData.empty()) {
435 this->outputModule_->assignToSolution(localCellData);
436 }
437
438 // Add cell data to perforations for RFT output
439 this->outputModule_->addRftDataToWells(localWellData, reportStepNum);
440 }
441
442 if (this->collectOnIORank_.isParallel() ||
443 this->collectOnIORank_.doesNeedReordering())
444 {
445 // Note: We don't need WBP (well-block averaged pressures) or
446 // inter-region flow rate values in order to create restart file
447 // output. There's consequently no need to collect those
448 // properties on the I/O rank.
449
450 this->collectOnIORank_.collect(localCellData,
451 this->outputModule_->getBlockData(),
453 /* wbpData = */ {},
457 /* interRegFlows = */ {},
458 flowsn,
459 floresn);
460 if (this->collectOnIORank_.isIORank()) {
461 this->outputModule_->assignGlobalFieldsToSolution(this->collectOnIORank_.globalCellData());
462 }
463 } else {
464 this->outputModule_->assignGlobalFieldsToSolution(localCellData);
465 }
466
467 if (this->collectOnIORank_.isIORank()) {
468 const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
469 const Scalar nextStepSize = simulator_.problem().nextTimeStepSize();
470 std::optional<int> timeStepIdx;
471 if (Parameters::Get<Parameters::EnableWriteAllSolutions>()) {
472 timeStepIdx = simulator_.timeStepIndex();
473 }
474 this->doWriteOutput(reportStepNum, timeStepIdx, isSubStep,
475 std::move(localCellData),
476 std::move(localWellData),
477 std::move(localGroupAndNetworkData),
478 std::move(localAquiferData),
479 std::move(localWellTestState),
480 this->actionState(),
481 this->udqState(),
482 this->summaryState(),
483 this->simulator_.problem().thresholdPressure().getRestartVector(),
485 Parameters::Get<Parameters::EclOutputDoublePrecision>(),
486 isFlowsn, std::move(flowsn),
487 isFloresn, std::move(floresn));
488 }
489 }
490
491 void beginRestart()
492 {
493 const auto enablePCHysteresis = simulator_.problem().materialLawManager()->enablePCHysteresis();
494 const auto enableNonWettingHysteresis = simulator_.problem().materialLawManager()->enableNonWettingHysteresis();
495 const auto enableWettingHysteresis = simulator_.problem().materialLawManager()->enableWettingHysteresis();
496 const auto oilActive = FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx);
497 const auto gasActive = FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx);
498 const auto waterActive = FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx);
499 const auto enableSwatinit = simulator_.vanguard().eclState().fieldProps().has_double("SWATINIT");
500
501 std::vector<RestartKey> solutionKeys {
502 {"PRESSURE", UnitSystem::measure::pressure},
503 {"SWAT", UnitSystem::measure::identity, waterActive},
504 {"SGAS", UnitSystem::measure::identity, gasActive},
505 {"TEMP", UnitSystem::measure::temperature, enableEnergy},
506 {"SSOLVENT", UnitSystem::measure::identity, enableSolvent},
507
508 {"RS", UnitSystem::measure::gas_oil_ratio, FluidSystem::enableDissolvedGas()},
509 {"RV", UnitSystem::measure::oil_gas_ratio, FluidSystem::enableVaporizedOil()},
510 {"RVW", UnitSystem::measure::oil_gas_ratio, FluidSystem::enableVaporizedWater()},
511 {"RSW", UnitSystem::measure::gas_oil_ratio, FluidSystem::enableDissolvedGasInWater()},
512
513 {"SGMAX", UnitSystem::measure::identity, enableNonWettingHysteresis && oilActive && gasActive},
514 {"SHMAX", UnitSystem::measure::identity, enableWettingHysteresis && oilActive && gasActive},
515
516 {"SOMAX", UnitSystem::measure::identity,
517 (enableNonWettingHysteresis && oilActive && waterActive)
518 || simulator_.problem().vapparsActive(simulator_.episodeIndex())},
519
520 {"SOMIN", UnitSystem::measure::identity, enablePCHysteresis && oilActive && gasActive},
521 {"SWHY1", UnitSystem::measure::identity, enablePCHysteresis && oilActive && waterActive},
522 {"SWMAX", UnitSystem::measure::identity, enableWettingHysteresis && oilActive && waterActive},
523
524 {"PPCW", UnitSystem::measure::pressure, enableSwatinit},
525 };
526
527 {
528 const auto& tracers = simulator_.vanguard().eclState().tracer();
529
530 for (const auto& tracer : tracers) {
531 const auto enableSolTracer =
532 ((tracer.phase == Phase::GAS) && FluidSystem::enableDissolvedGas()) ||
533 ((tracer.phase == Phase::OIL) && FluidSystem::enableVaporizedOil());
534
535 solutionKeys.emplace_back(tracer.fname(), UnitSystem::measure::identity, true);
536 solutionKeys.emplace_back(tracer.sname(), UnitSystem::measure::identity, enableSolTracer);
537 }
538 }
539
540 const auto& inputThpres = eclState().getSimulationConfig().getThresholdPressure();
541 const std::vector<RestartKey> extraKeys {
542 {"OPMEXTRA", UnitSystem::measure::identity, false},
543 {"THRESHPR", UnitSystem::measure::pressure, inputThpres.active()},
544 };
545
546 const auto& gridView = this->simulator_.vanguard().gridView();
547 const auto numElements = gridView.size(/*codim=*/0);
548
549 {
550 // The episodeIndex is rewound one step back before calling
551 // beginRestart() and cannot be used here. We just ask the
552 // initconfig directly to be sure that we use the correct index.
553 const auto restartStepIdx = this->simulator_.vanguard()
554 .eclState().getInitConfig().getRestartStep();
555
556 this->outputModule_->allocBuffers(numElements,
558 /*isSubStep = */false,
559 /*log = */ false,
560 /*isRestart = */true);
561 }
562
563 {
564 const auto restartValues =
565 loadParallelRestart(this->eclIO_.get(),
566 this->actionState(),
567 this->summaryState(),
568 solutionKeys, extraKeys, gridView.comm());
569
570 for (auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
571 const auto globalIdx = this->collectOnIORank_.localIdxToGlobalIdx(elemIdx);
572 this->outputModule_->setRestart(restartValues.solution, elemIdx, globalIdx);
573 }
574
575 auto& tracer_model = simulator_.problem().tracerModel();
576 for (int tracer_index = 0; tracer_index < tracer_model.numTracers(); ++tracer_index) {
577 // Free tracers
578 {
579 const auto& free_tracer_name = tracer_model.fname(tracer_index);
580 const auto& free_tracer_solution = restartValues.solution
582
583 for (auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
584 const auto globalIdx = this->collectOnIORank_.localIdxToGlobalIdx(elemIdx);
585 tracer_model.setFreeTracerConcentration
587 }
588 }
589
590 // Solution tracer (only if DISGAS/VAPOIL are active for gas/oil tracers)
591 if ((tracer_model.phase(tracer_index) == Phase::GAS && FluidSystem::enableDissolvedGas()) ||
592 (tracer_model.phase(tracer_index) == Phase::OIL && FluidSystem::enableVaporizedOil()))
593 {
594 tracer_model.setEnableSolTracers(tracer_index, true);
595
596 const auto& sol_tracer_name = tracer_model.sname(tracer_index);
597 const auto& sol_tracer_solution = restartValues.solution
598 .template data<double>(sol_tracer_name);
599
600 for (auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
601 const auto globalIdx = this->collectOnIORank_.localIdxToGlobalIdx(elemIdx);
602 tracer_model.setSolTracerConcentration
604 }
605 }
606 else {
607 tracer_model.setEnableSolTracers(tracer_index, false);
608
609 for (auto elemIdx = 0*numElements; elemIdx < numElements; ++elemIdx) {
610 tracer_model.setSolTracerConcentration(tracer_index, elemIdx, 0.0);
611 }
612 }
613 }
614
615 if (inputThpres.active()) {
616 const_cast<Simulator&>(this->simulator_)
617 .problem().thresholdPressure()
618 .setFromRestart(restartValues.getExtra("THRESHPR"));
619 }
620
621 restartTimeStepSize_ = restartValues.getExtra("OPMEXTRA")[0];
622 if (restartTimeStepSize_ <= 0) {
623 restartTimeStepSize_ = std::numeric_limits<double>::max();
624 }
625
626 // Initialize the well model from restart values
627 this->simulator_.problem().wellModel()
628 .initFromRestartFile(restartValues);
629
630 if (!restartValues.aquifer.empty()) {
631 this->simulator_.problem().mutableAquiferModel()
632 .initFromRestart(restartValues.aquifer);
633 }
634 }
635 }
636
637 void endRestart()
638 {
639 // We need these objects to satisfy the interface requirements of
640 // member function calc_inplace(), but the objects are otherwise
641 // unused and intentionally so.
642 auto miscSummaryData = std::map<std::string, double>{};
643 auto regionData = std::map<std::string, std::vector<double>>{};
644
645 // Note: calc_inplace() *also* assigns the output module's
646 // "initialInplace_" data member. This is, semantically speaking,
647 // very wrong, as the run's intial volumes then correspond to the
648 // volumes at the restart time instead of the start of the base run.
649 // Nevertheless, this is how Flow has "always" done it.
650 //
651 // See GenericOutputBlackoilModule::accumulateRegionSums() for
652 // additional comments.
653 auto inplace = this->outputModule_
654 ->calc_inplace(miscSummaryData, regionData,
655 this->simulator_.gridView().comm());
656
657 if (this->collectOnIORank_.isIORank()) {
658 this->inplace_ = std::move(inplace);
659 }
660 }
661
662 const OutputBlackOilModule<TypeTag>& outputModule() const
663 { return *outputModule_; }
664
665 OutputBlackOilModule<TypeTag>& mutableOutputModule() const
666 { return *outputModule_; }
667
668 Scalar restartTimeStepSize() const
669 { return restartTimeStepSize_; }
670
671 template <class Serializer>
672 void serializeOp(Serializer& serializer)
673 {
674 serializer(*outputModule_);
675 }
676
677private:
678 static bool enableEclOutput_()
679 {
680 static bool enable = Parameters::Get<Parameters::EnableEclOutput>();
681 return enable;
682 }
683
684 const EclipseState& eclState() const
685 { return simulator_.vanguard().eclState(); }
686
687 SummaryState& summaryState()
688 { return simulator_.vanguard().summaryState(); }
689
690 Action::State& actionState()
691 { return simulator_.vanguard().actionState(); }
692
693 UDQState& udqState()
694 { return simulator_.vanguard().udqState(); }
695
696 const Schedule& schedule() const
697 { return simulator_.vanguard().schedule(); }
698
699 void prepareLocalCellData(const bool isSubStep,
700 const int reportStepNum)
701 {
702 OPM_TIMEBLOCK(prepareLocalCellData);
703
704 if (this->outputModule_->localDataValid()) {
705 return;
706 }
707
708 const auto& gridView = simulator_.vanguard().gridView();
709 const bool log = this->collectOnIORank_.isIORank();
710
711 const int num_interior = detail::
712 countLocalInteriorCellsGridView(gridView);
713 this->outputModule_->
714 allocBuffers(num_interior, reportStepNum,
715 isSubStep && !Parameters::Get<Parameters::EnableWriteAllSolutions>(),
716 log, /*isRestart*/ false);
717
718 ElementContext elemCtx(simulator_);
719
720 OPM_BEGIN_PARALLEL_TRY_CATCH();
721
722 {
724
725 this->outputModule_->prepareDensityAccumulation();
726
727 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
728 elemCtx.updatePrimaryStencil(elem);
729 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
730
731 this->outputModule_->processElement(elemCtx);
732 }
733
734 this->outputModule_->accumulateDensityParallel();
735 }
736
737 if constexpr (enableMech) {
738 if (simulator_.vanguard().eclState().runspec().mech()) {
740 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
741 elemCtx.updatePrimaryStencil(elem);
742 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
743 outputModule_->processElementMech(elemCtx);
744 }
745 }
746 }
747
748 if (! this->simulator_.model().linearizer().getFlowsInfo().empty()) {
750 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
751 elemCtx.updatePrimaryStencil(elem);
752 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
753
754 this->outputModule_->processElementFlows(elemCtx);
755 }
756 }
757
758 {
760 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
761 elemCtx.updatePrimaryStencil(elem);
762 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
763
764 this->outputModule_->processElementBlockData(elemCtx);
765 }
766 }
767
768 {
770
771#ifdef _OPENMP
772#pragma omp parallel for
773#endif
774 for (int dofIdx = 0; dofIdx < num_interior; ++dofIdx) {
775 const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, /*timeIdx=*/0);
776 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
777
778 this->outputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
779 }
780 }
781
782 this->outputModule_->validateLocalData();
783
784 OPM_END_PARALLEL_TRY_CATCH("EclWriter::prepareLocalCellData() failed: ",
785 this->simulator_.vanguard().grid().comm());
786 }
787
788 void captureLocalFluxData()
789 {
791
792 const auto& gridView = this->simulator_.vanguard().gridView();
793 const auto timeIdx = 0u;
794
795 auto elemCtx = ElementContext { this->simulator_ };
796
797 const auto elemMapper = ElementMapper { gridView, Dune::mcmgElementLayout() };
798 const auto activeIndex = [&elemMapper](const Element& e)
799 {
800 return elemMapper.index(e);
801 };
802
803 const auto cartesianIndex = [this](const int elemIndex)
804 {
805 return this->cartMapper_.cartesianIndex(elemIndex);
806 };
807
808 this->outputModule_->initializeFluxData();
809
810 OPM_BEGIN_PARALLEL_TRY_CATCH();
811
812 for (const auto& elem : elements(gridView, Dune::Partitions::interiorBorder)) {
813 elemCtx.updateStencil(elem);
814 elemCtx.updateIntensiveQuantities(timeIdx);
815 elemCtx.updateExtensiveQuantities(timeIdx);
816
817 this->outputModule_->processFluxes(elemCtx, activeIndex, cartesianIndex);
818 }
819
820 OPM_END_PARALLEL_TRY_CATCH("EclWriter::captureLocalFluxData() failed: ",
821 this->simulator_.vanguard().grid().comm())
822
823 this->outputModule_->finalizeFluxData();
824 }
825
826 Simulator& simulator_;
827 std::unique_ptr<OutputBlackOilModule<TypeTag> > outputModule_;
828 Scalar restartTimeStepSize_;
829 int rank_ ;
830 Inplace inplace_;
831};
832
833} // namespace Opm
834
835#endif // OPM_ECL_WRITER_HPP
Collects necessary output values and pass it to opm-common's ECL output.
Helper class for grid instantiation of ECL file-format using problems.
Output module for the results black oil model writing in ECL binary format.
Definition EclGenericWriter.hpp:65
Collects necessary output values and pass it to opm-common's ECL output.
Definition EclWriter.hpp:110
void writeInitialFIPReport()
Writes the initial FIP report as configured in RPTSOL.
Definition EclWriter.hpp:344
void evalSummaryState(bool isSubStep)
collect and pass data and pass it to eclIO writer
Definition EclWriter.hpp:198
void compress()
Form CSR adjacency matrix representation of input graph from connections established in previous call...
Definition InterRegFlows.cpp:165
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition OutputBlackoilModule.hpp:182
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
double simulationTimeElapsed() const override
Time elapsed since the start of the simulation until the beginning of the current time step [s].
Definition SimulatorTimer.cpp:122
virtual boost::posix_time::ptime currentDateTime() const
Return the current time as a posix time object.
Definition SimulatorTimerInterface.cpp:28
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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
Definition EclWriter.hpp:65
Definition EclWriter.hpp:62
Definition EclWriter.hpp:75
Definition EclWriter.hpp:72