My Project
Loading...
Searching...
No Matches
SimulatorFullyImplicitBlackoilEbos.hpp
1/*
2 Copyright 2013, 2015, 2020 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2015 Andreas Lauser
4 Copyright 2017 IRIS
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_SIMULATORFULLYIMPLICITBLACKOILEBOS_HEADER_INCLUDED
23#define OPM_SIMULATORFULLYIMPLICITBLACKOILEBOS_HEADER_INCLUDED
24
25#include <dune/common/hash.hh>
26
27#include <opm/simulators/flow/BlackoilModelEbos.hpp>
28#include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
29#include <opm/simulators/flow/ConvergenceOutputConfiguration.hpp>
30#include <opm/simulators/flow/ExtraConvergenceOutputThread.hpp>
31#include <opm/simulators/flow/NonlinearSolverEbos.hpp>
32#include <opm/simulators/aquifers/BlackoilAquiferModel.hpp>
33#include <opm/simulators/timestepping/AdaptiveTimeSteppingEbos.hpp>
34#include <opm/simulators/timestepping/ConvergenceReport.hpp>
35#include <opm/simulators/utils/moduleVersion.hpp>
36#include <opm/simulators/wells/WellState.hpp>
37
38#include <opm/grid/utility/StopWatch.hpp>
39
40#include <opm/input/eclipse/Units/UnitSystem.hpp>
41
42#include <opm/common/ErrorMacros.hpp>
43
44#include <boost/date_time/gregorian/gregorian.hpp>
45
46#include <memory>
47#include <optional>
48#include <string>
49#include <string_view>
50#include <thread>
51#include <utility>
52#include <vector>
53
54#if HAVE_HDF5
55#include <ebos/hdf5serializer.hh>
56#endif
57
58namespace Opm::Properties {
59
60template<class TypeTag, class MyTypeTag>
64template<class TypeTag, class MyTypeTag>
66 using type = UndefinedProperty;
67};
68
69template <class TypeTag, class MyTypeTag>
74
75template <class TypeTag, class MyTypeTag>
77{
78 using type = UndefinedProperty;
79};
80
81template <class TypeTag, class MyTypeTag>
83{
84 using type = UndefinedProperty;
85};
86
87template <class TypeTag, class MyTypeTag>
89{
90 using type = UndefinedProperty;
91};
92
93template<class TypeTag>
94struct EnableTerminalOutput<TypeTag, TTag::EclFlowProblem> {
95 static constexpr bool value = true;
96};
97template<class TypeTag>
98struct EnableAdaptiveTimeStepping<TypeTag, TTag::EclFlowProblem> {
99 static constexpr bool value = true;
100};
101template<class TypeTag>
102struct EnableTuning<TypeTag, TTag::EclFlowProblem> {
103 static constexpr bool value = false;
104};
105
106template <class TypeTag>
107struct OutputExtraConvergenceInfo<TypeTag, TTag::EclFlowProblem>
108{
109 static constexpr auto* value = "none";
110};
111
112template <class TypeTag>
113struct SaveStep<TypeTag, TTag::EclFlowProblem>
114{
115 static constexpr auto* value = "";
116};
117
118template <class TypeTag>
119struct SaveFile<TypeTag, TTag::EclFlowProblem>
120{
121 static constexpr auto* value = "";
122};
123
124template <class TypeTag>
125struct LoadStep<TypeTag, TTag::EclFlowProblem>
126{
127 static constexpr int value = -1;
128};
129
130} // namespace Opm::Properties
131
132namespace Opm {
133
134void outputReportStep(const SimulatorTimer& timer);
135void outputTimestampFIP(const SimulatorTimer& timer,
136 const std::string& title,
137 const std::string& version);
138void checkSerializedCmdLine(const std::string& current,
139 const std::string& stored);
140
142template<class TypeTag>
144{
145public:
156
160
163 typedef typename Model::ModelParameters ModelParameters;
164 typedef typename Solver::SolverParameters SolverParameters;
166
167
189 SimulatorFullyImplicitBlackoilEbos(Simulator& ebosSimulator)
190 : ebosSimulator_(ebosSimulator)
191 {
192 phaseUsage_ = phaseUsageFromDeck(eclState());
193
194 // Only rank 0 does print to std::cout, and only if specifically requested.
195 this->terminalOutput_ = false;
196 if (this->grid().comm().rank() == 0) {
197 this->terminalOutput_ = EWOMS_GET_PARAM(TypeTag, bool, EnableTerminalOutput);
198
199 this->startConvergenceOutputThread(EWOMS_GET_PARAM(TypeTag, std::string,
200 OutputExtraConvergenceInfo),
201 R"(OutputExtraConvergenceInfo (--output-extra-convergence-info))");
202 }
203
204 const std::string saveSpec = EWOMS_GET_PARAM(TypeTag, std::string, SaveStep);
205 if (saveSpec == "all") {
206 saveStride_ = 1;
207 } else if (!saveSpec.empty() && saveSpec[0] == ':') {
208 saveStride_ = std::atoi(saveSpec.c_str()+1);
209 } else if (!saveSpec.empty()) {
210 saveStep_ = std::atoi(saveSpec.c_str());
211 }
212
213 loadStep_ = EWOMS_GET_PARAM(TypeTag, int, LoadStep);
214
215 saveFile_ = EWOMS_GET_PARAM(TypeTag, std::string, SaveFile);
216 if (saveFile_.empty()) {
217 const auto& ioconfig = ebosSimulator_.vanguard().eclState().getIOConfig();
218 saveFile_ = ioconfig.fullBasePath() + ".OPMRST";
219 if (loadStep_ != -1 && !std::filesystem::exists(saveFile_)) {
220 std::filesystem::path path(ioconfig.getInputDir() + "/");
221 path.replace_filename(ioconfig.getBaseName() + ".OPMRST");
222 saveFile_ = path;
223 if (!std::filesystem::exists(saveFile_)) {
224 OPM_THROW(std::runtime_error, "Error locating serialized restart file");
225 }
226 }
227 }
228 }
229
231 {
232 // Safe to call on all ranks, not just the I/O rank.
233 this->endConvergenceOutputThread();
234 }
235
236 static void registerParameters()
237 {
238 ModelParameters::registerParameters();
239 SolverParameters::registerParameters();
240 TimeStepper::registerParameters();
241
242 EWOMS_REGISTER_PARAM(TypeTag, bool, EnableTerminalOutput,
243 "Print high-level information about the simulation's progress to the terminal");
244 EWOMS_REGISTER_PARAM(TypeTag, bool, EnableAdaptiveTimeStepping,
245 "Use adaptive time stepping between report steps");
246 EWOMS_REGISTER_PARAM(TypeTag, bool, EnableTuning,
247 "Honor some aspects of the TUNING keyword.");
248 EWOMS_REGISTER_PARAM(TypeTag, std::string, OutputExtraConvergenceInfo,
249 "Provide additional convergence output "
250 "files for diagnostic purposes. "
251 "\"none\" gives no extra output and "
252 "overrides all other options, "
253 "\"steps\" generates an INFOSTEP file, "
254 "\"iterations\" generates an INFOITER file. "
255 "Combine options with commas, e.g., "
256 "\"steps,iterations\" for multiple outputs.");
257 EWOMS_REGISTER_PARAM(TypeTag, std::string, SaveStep,
258 "Save serialized state to .OPMRST file. "
259 "Either a specific report step, \"all\" to save "
260 "all report steps or \":x\" to save every x'th step.");
261 EWOMS_REGISTER_PARAM(TypeTag, int, LoadStep,
262 "Load serialized state from .OPMRST file. "
263 "Either a specific report step, or 0 to load last "
264 "stored report step.");
265 EWOMS_REGISTER_PARAM(TypeTag, std::string, SaveFile,
266 "FileName for .OPMRST file used for serialized state. "
267 "If empty, CASENAME.OPMRST is used.");
268 EWOMS_HIDE_PARAM(TypeTag, SaveFile);
269 }
270
278 {
279 init(timer);
280 // Make cache up to date. No need for updating it in elementCtx.
281 ebosSimulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
282 // Main simulation loop.
283 while (!timer.done()) {
284 bool continue_looping = runStep(timer);
285 if (!continue_looping) break;
286 }
287 return finalize();
288 }
289
290 void init(SimulatorTimer &timer)
291 {
292 ebosSimulator_.setEpisodeIndex(-1);
293
294 // Create timers and file for writing timing info.
295 solverTimer_ = std::make_unique<time::StopWatch>();
296 totalTimer_ = std::make_unique<time::StopWatch>();
297 totalTimer_->start();
298
299 // adaptive time stepping
300 bool enableAdaptive = EWOMS_GET_PARAM(TypeTag, bool, EnableAdaptiveTimeStepping);
301 bool enableTUNING = EWOMS_GET_PARAM(TypeTag, bool, EnableTuning);
302 if (enableAdaptive) {
303 const UnitSystem& unitSystem = this->ebosSimulator_.vanguard().eclState().getUnits();
304 if (enableTUNING) {
305 const auto& sched_state = schedule()[timer.currentStepNum()];
306 auto max_next_tstep = sched_state.max_next_tstep();
307 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(max_next_tstep,
308 sched_state.tuning(),
309 unitSystem, terminalOutput_);
310 }
311 else {
312 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(unitSystem, terminalOutput_);
313 }
314
315 if (isRestart()) {
316 // For restarts the ebosSimulator may have gotten some information
317 // about the next timestep size from the OPMEXTRA field
318 adaptiveTimeStepping_->setSuggestedNextStep(ebosSimulator_.timeStepSize());
319 }
320 }
321 }
322
323 bool runStep(SimulatorTimer& timer)
324 {
325 if (schedule().exitStatus().has_value()) {
326 if (terminalOutput_) {
327 OpmLog::info("Stopping simulation since EXIT was triggered by an action keyword.");
328 }
329 report_.success.exit_status = schedule().exitStatus().value();
330 return false;
331 }
332
333 if (loadStep_ > -1) {
334 loadTimerInfo(timer);
335 }
336
337 // Report timestep.
338 if (terminalOutput_) {
339 std::ostringstream ss;
340 timer.report(ss);
341 OpmLog::debug(ss.str());
342 }
343
344 if (terminalOutput_) {
345 outputReportStep(timer);
346 }
347
348 // write the inital state at the report stage
349 if (timer.initialStep()) {
350 Dune::Timer perfTimer;
351 perfTimer.start();
352
353 ebosSimulator_.setEpisodeIndex(-1);
354 ebosSimulator_.setEpisodeLength(0.0);
355 ebosSimulator_.setTimeStepSize(0.0);
356 wellModel_().beginReportStep(timer.currentStepNum());
357 ebosSimulator_.problem().writeOutput();
358
359 report_.success.output_write_time += perfTimer.stop();
360 }
361
362 // Run a multiple steps of the solver depending on the time step control.
363 solverTimer_->start();
364
365 auto solver = createSolver(wellModel_());
366
367 ebosSimulator_.startNextEpisode(
368 ebosSimulator_.startTime()
369 + schedule().seconds(timer.currentStepNum()),
370 timer.currentStepLength());
371 ebosSimulator_.setEpisodeIndex(timer.currentStepNum());
372 if (loadStep_> -1) {
373 wellModel_().prepareDeserialize(loadStep_ - 1);
375 loadStep_ = -1;
376 ebosSimulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
377 }
378 solver->model().beginReportStep();
379
380 bool enableTUNING = EWOMS_GET_PARAM(TypeTag, bool, EnableTuning);
381
382 // If sub stepping is enabled allow the solver to sub cycle
383 // in case the report steps are too large for the solver to converge
384 //
385 // \Note: The report steps are met in any case
386 // \Note: The sub stepping will require a copy of the state variables
387 if (adaptiveTimeStepping_) {
388 const auto& events = schedule()[timer.currentStepNum()].events();
389 if (enableTUNING) {
390 if (events.hasEvent(ScheduleEvents::TUNING_CHANGE)) {
391 const auto& sched_state = schedule()[timer.currentStepNum()];
392 const auto& tuning = sched_state.tuning();
393 const auto& max_next_tstep = sched_state.max_next_tstep();
394 adaptiveTimeStepping_->updateTUNING(max_next_tstep, tuning);
395 }
396 }
397 bool event = events.hasEvent(ScheduleEvents::NEW_WELL) ||
398 events.hasEvent(ScheduleEvents::INJECTION_TYPE_CHANGED) ||
399 events.hasEvent(ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER) ||
400 events.hasEvent(ScheduleEvents::WELL_STATUS_CHANGE);
401 auto stepReport = adaptiveTimeStepping_->step(timer, *solver, event, nullptr);
402 report_ += stepReport;
403 //Pass simulation report to eclwriter for summary output
404 ebosSimulator_.problem().setSimulationReport(report_);
405 } else {
406 // solve for complete report step
407 auto stepReport = solver->step(timer);
408 report_ += stepReport;
409 if (terminalOutput_) {
410 std::ostringstream ss;
411 stepReport.reportStep(ss);
412 OpmLog::info(ss.str());
413 }
414 }
415
416 // write simulation state at the report stage
417 Dune::Timer perfTimer;
418 perfTimer.start();
419 const double nextstep = adaptiveTimeStepping_ ? adaptiveTimeStepping_->suggestedNextStep() : -1.0;
420 ebosSimulator_.problem().setNextTimeStepSize(nextstep);
421 ebosSimulator_.problem().writeOutput();
422 report_.success.output_write_time += perfTimer.stop();
423
424 solver->model().endReportStep();
425
426 // take time that was used to solve system for this reportStep
427 solverTimer_->stop();
428
429 // update timing.
430 report_.success.solver_time += solverTimer_->secsSinceStart();
431
432 if (this->grid().comm().rank() == 0) {
433 // Destructively grab the step convergence reports. The solver
434 // object and the model object contained therein are about to go
435 // out of scope.
436 this->writeConvergenceOutput(solver->model().getStepReportsDestructively());
437 }
438
439 // Increment timer, remember well state.
440 ++timer;
441
442 if (terminalOutput_) {
443 if (!timer.initialStep()) {
444 const std::string version = moduleVersionName();
445 outputTimestampFIP(timer, eclState().getTitle(), version);
446 }
447
448 std::string msg =
449 "Time step took " + std::to_string(solverTimer_->secsSinceStart()) + " seconds; "
450 "total solver time " + std::to_string(report_.success.solver_time) + " seconds.";
451 OpmLog::debug(msg);
452 }
453
454 handleSave(timer);
455
456 return true;
457 }
458
459 SimulatorReport finalize()
460 {
461 // make sure all output is written to disk before run is finished
462 {
463 Dune::Timer finalOutputTimer;
464 finalOutputTimer.start();
465
466 ebosSimulator_.problem().finalizeOutput();
467 report_.success.output_write_time += finalOutputTimer.stop();
468 }
469
470 // Stop timer and create timing report
471 totalTimer_->stop();
472 report_.success.total_time = totalTimer_->secsSinceStart();
473 report_.success.converged = true;
474
475 return report_;
476 }
477
478 const Grid& grid() const
479 { return ebosSimulator_.vanguard().grid(); }
480
481 template<class Serializer>
482 void serializeOp(Serializer& serializer)
483 {
484 serializer(ebosSimulator_);
485 serializer(report_);
486 serializer(adaptiveTimeStepping_);
487 }
488
489protected:
490
491 std::unique_ptr<Solver> createSolver(WellModel& wellModel)
492 {
493 auto model = std::make_unique<Model>(ebosSimulator_,
494 modelParam_,
495 wellModel,
496 terminalOutput_);
497
498 return std::make_unique<Solver>(solverParam_, std::move(model));
499 }
500
501 const EclipseState& eclState() const
502 { return ebosSimulator_.vanguard().eclState(); }
503
504
505 const Schedule& schedule() const
506 { return ebosSimulator_.vanguard().schedule(); }
507
508 bool isRestart() const
509 {
510 const auto& initconfig = eclState().getInitConfig();
511 return initconfig.restartRequested();
512 }
513
514 WellModel& wellModel_()
515 { return ebosSimulator_.problem().wellModel(); }
516
517 const WellModel& wellModel_() const
518 { return ebosSimulator_.problem().wellModel(); }
519
520 void startConvergenceOutputThread(std::string_view convOutputOptions,
521 std::string_view optionName)
522 {
523 const auto config = ConvergenceOutputConfiguration {
524 convOutputOptions, optionName
525 };
526 if (! config.want(ConvergenceOutputConfiguration::Option::Iterations)) {
527 return;
528 }
529
531 [compNames = typename Model::ComponentName{}](const int compIdx)
532 { return std::string_view { compNames.name(compIdx) }; }
533 };
534
536 [usys = this->eclState().getUnits()](const double time)
537 { return usys.from_si(UnitSystem::measure::time, time); }
538 };
539
540 this->convergenceOutputQueue_.emplace();
541 this->convergenceOutputObject_.emplace
542 (this->eclState().getIOConfig().getOutputDir(),
543 this->eclState().getIOConfig().getBaseName(),
544 std::move(getPhaseName),
545 std::move(convertTime),
546 config, *this->convergenceOutputQueue_);
547
548 this->convergenceOutputThread_
550 &this->convergenceOutputObject_.value());
551 }
552
553 void writeConvergenceOutput(std::vector<StepReport>&& reports)
554 {
555 if (! this->convergenceOutputThread_.has_value()) {
556 return;
557 }
558
559 auto requests = std::vector<ConvergenceReportQueue::OutputRequest>{};
560 requests.reserve(reports.size());
561
562 for (auto&& report : reports) {
563 requests.push_back({ report.report_step, report.current_step, std::move(report.report) });
564 }
565
566 this->convergenceOutputQueue_->enqueue(std::move(requests));
567 }
568
569 void endConvergenceOutputThread()
570 {
571 if (! this->convergenceOutputThread_.has_value()) {
572 return;
573 }
574
575 this->convergenceOutputQueue_->signalLastOutputRequest();
576 this->convergenceOutputThread_->join();
577 }
578
581 {
582 if (saveStride_ == -1 && saveStep_ == -1) {
583 return;
584 }
585
586 OPM_BEGIN_PARALLEL_TRY_CATCH();
587
588 int nextStep = timer.currentStepNum();
589 if ((saveStep_ != -1 && nextStep == saveStep_) ||
590 (saveStride_ != -1 && (nextStep % saveStride_) == 0)) {
591#if !HAVE_HDF5
592 OpmLog::error("Saving of serialized state requested, but no HDF5 support available.");
593#else
594 const std::string groupName = "/report_step/" + std::to_string(nextStep);
595 if (nextStep == saveStride_ || nextStep == saveStep_) {
596 std::filesystem::remove(saveFile_);
597 }
600 EclGenericVanguard::comm());
601 if (nextStep == saveStride_ || nextStep == saveStep_) {
602 std::ostringstream str;
603 Parameters::printValues<TypeTag>(str);
604 writer.writeHeader("OPM Flow",
607 ebosSimulator_.vanguard().caseName(),
608 str.str(),
609 EclGenericVanguard::comm().size());
610
611 if (EclGenericVanguard::comm().size() > 1) {
612 const auto& cellMapping = ebosSimulator_.vanguard().globalCell();
613 std::size_t hash = Dune::hash_range(cellMapping.begin(), cellMapping.end());
614 writer.write(hash, "/", "grid_checksum");
615 }
616 }
617 writer.write(*this, groupName, "simulator_data");
618 writer.write(timer, groupName, "simulator_timer",
620 OpmLog::info("Serialized state written for report step " + std::to_string(nextStep));
621#endif
622 }
623
624 OPM_END_PARALLEL_TRY_CATCH("Error saving serialized state: ",
625 EclGenericVanguard::comm());
626 }
627
630 {
631#if !HAVE_HDF5
632 OpmLog::error("Loading of serialized state requested, but no HDF5 support available.");
633 loadStep_ = -1;
634#else
635 OPM_BEGIN_PARALLEL_TRY_CATCH();
636
639 EclGenericVanguard::comm());
640
641 if (loadStep_ == 0) {
642 loadStep_ = reader.lastReportStep();
643 }
644
645 OpmLog::info("Loading serialized state for report step " + std::to_string(loadStep_));
646 const std::string groupName = "/report_step/" + std::to_string(loadStep_);
647 reader.read(timer, groupName, "simulator_timer", HDF5File::DataSetMode::ROOT_ONLY);
648
649 std::tuple<std::array<std::string,5>,int> header;
650 reader.read(header, "/", "simulator_info", HDF5File::DataSetMode::ROOT_ONLY);
651 const auto& [strings, procs] = header;
652
653 if (EclGenericVanguard::comm().size() != procs) {
654 throw std::runtime_error("Number of processes (procs=" +
655 std::to_string(EclGenericVanguard::comm().size()) +
656 ") does not match .OPMRST file (procs=" +
657 std::to_string(procs) + ")");
658 }
659
660 if (EclGenericVanguard::comm().size() > 1) {
661 std::size_t stored_hash;
662 reader.read(stored_hash, "/", "grid_checksum");
663 const auto& cellMapping = ebosSimulator_.vanguard().globalCell();
664 std::size_t hash = Dune::hash_range(cellMapping.begin(), cellMapping.end());
665 if (hash != stored_hash) {
666 throw std::runtime_error("Grid hash mismatch, .OPMRST file cannot be used");
667 }
668 }
669
670 if (EclGenericVanguard::comm().rank() == 0) {
671 std::ostringstream str;
672 Parameters::printValues<TypeTag>(str);
673 checkSerializedCmdLine(str.str(), strings[4]);
674 }
675
676 OPM_END_PARALLEL_TRY_CATCH("Error loading serialized state: ",
677 EclGenericVanguard::comm());
678#endif
679 }
680
683 {
684#if HAVE_HDF5
685 OPM_BEGIN_PARALLEL_TRY_CATCH();
686
689 EclGenericVanguard::comm());
690 const std::string groupName = "/report_step/" + std::to_string(loadStep_);
691 reader.read(*this, groupName, "simulator_data");
692
693 OPM_END_PARALLEL_TRY_CATCH("Error loading serialized state: ",
694 EclGenericVanguard::comm());
695#endif
696 }
697
698 // Data.
699 Simulator& ebosSimulator_;
700 std::unique_ptr<WellConnectionAuxiliaryModule<TypeTag>> wellAuxMod_;
701
702 ModelParameters modelParam_;
703 SolverParameters solverParam_;
704
705 // Observed objects.
706 PhaseUsage phaseUsage_;
707 // Misc. data
708 bool terminalOutput_;
709
710 SimulatorReport report_;
711 std::unique_ptr<time::StopWatch> solverTimer_;
712 std::unique_ptr<time::StopWatch> totalTimer_;
713 std::unique_ptr<TimeStepper> adaptiveTimeStepping_;
714
715 std::optional<ConvergenceReportQueue> convergenceOutputQueue_{};
716 std::optional<ConvergenceOutputThread> convergenceOutputObject_{};
717 std::optional<std::thread> convergenceOutputThread_{};
718
719 int saveStride_ = -1;
720 int saveStep_ = -1;
721 int loadStep_ = -1;
722 std::string saveFile_;
723};
724
725} // namespace Opm
726
727#endif // OPM_SIMULATOR_FULLY_IMPLICIT_BLACKOIL_EBOS_HPP
Definition AquiferInterface.hpp:35
std::function< std::string_view(int)> ComponentToPhaseName
Protocol for converting a phase/component ID into a human readable phase/component name.
Definition ExtraConvergenceOutputThread.hpp:109
void writeASynchronous()
Output thread worker function.
Definition ExtraConvergenceOutputThread.cpp:297
std::function< double(double)> ConvertToTimeUnits
Protocol for converting an SI elapsed time value into an equivalent time value in the run's output co...
Definition ExtraConvergenceOutputThread.hpp:115
@ READ
Open existing file for reading.
@ APPEND
Append to an existing file (creates new if not)
@ ROOT_ONLY
A single dataset created at the root process.
a simulator for the blackoil model
Definition SimulatorFullyImplicitBlackoilEbos.hpp:144
void handleSave(SimulatorTimer &timer)
Serialization of simulator data to .OPMRST files at end of report steps.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:580
std::string saveFile_
File to load/save serialized state from/to.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:722
SimulatorReport run(SimulatorTimer &timer)
Run the simulation.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:277
SimulatorFullyImplicitBlackoilEbos(Simulator &ebosSimulator)
Initialise from parameters and objects to observe.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:189
int loadStep_
Step to load serialized state from.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:721
int saveStep_
Specific step to save serialized state at.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:720
void loadSimulatorState()
Load simulator state from serialized state.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:682
int saveStride_
Stride to save serialized state at.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:719
void loadTimerInfo(SimulatorTimer &timer)
Load timer info from serialized state.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:629
Definition SimulatorTimer.hpp:38
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
std::string moduleVersionName()
Return the version name of the module, for example "2015.10" (for a release branch) or "2016....
Definition moduleVersion.cpp:34
std::string compileTimestamp()
Return a string "dd-mm-yyyy at HH::MM::SS hrs" which is the time the binary was compiled.
Definition moduleVersion.cpp:57
std::string moduleVersion()
Return a string containing both the name and hash, if N is the name and H is the hash it will be "N (...
Definition moduleVersion.cpp:50
PhaseUsage phaseUsageFromDeck(const EclipseState &eclipseState)
Looks at presence of WATER, OIL and GAS keywords in state object to determine active phases.
Definition phaseUsageFromDeck.cpp:145
Definition BlackoilPhases.hpp:46
Definition SimulatorFullyImplicitBlackoilEbos.hpp:61
Definition BlackoilWellModel.hpp:80
Definition SimulatorFullyImplicitBlackoilEbos.hpp:65
Definition SimulatorFullyImplicitBlackoilEbos.hpp:83
Definition SimulatorFullyImplicitBlackoilEbos.hpp:71
Definition SimulatorFullyImplicitBlackoilEbos.hpp:89
Definition SimulatorFullyImplicitBlackoilEbos.hpp:77
Definition SimulatorReport.hpp:100