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#include <fmt/format.h>
27
28#include <opm/simulators/flow/BlackoilModelEbos.hpp>
29#include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
30#include <opm/simulators/flow/ConvergenceOutputConfiguration.hpp>
31#include <opm/simulators/flow/ExtraConvergenceOutputThread.hpp>
32#include <opm/simulators/flow/NonlinearSolverEbos.hpp>
33#include <opm/simulators/aquifers/BlackoilAquiferModel.hpp>
34#include <opm/simulators/timestepping/AdaptiveTimeSteppingEbos.hpp>
35#include <opm/simulators/timestepping/ConvergenceReport.hpp>
36#include <opm/simulators/utils/moduleVersion.hpp>
37#include <opm/simulators/wells/WellState.hpp>
38
39#include <opm/grid/utility/StopWatch.hpp>
40
41#include <opm/input/eclipse/Units/UnitSystem.hpp>
42
43#include <opm/common/ErrorMacros.hpp>
44
45#include <boost/date_time/gregorian/gregorian.hpp>
46
47#include <memory>
48#include <optional>
49#include <string>
50#include <string_view>
51#include <thread>
52#include <utility>
53#include <vector>
54
55#if HAVE_HDF5
56#include <ebos/hdf5serializer.hh>
57#endif
58
59namespace Opm::Properties {
60
61template<class TypeTag, class MyTypeTag>
65template<class TypeTag, class MyTypeTag>
67 using type = UndefinedProperty;
68};
69
70template <class TypeTag, class MyTypeTag>
75
76template <class TypeTag, class MyTypeTag>
78{
79 using type = UndefinedProperty;
80};
81
82template <class TypeTag, class MyTypeTag>
84{
85 using type = UndefinedProperty;
86};
87
88template <class TypeTag, class MyTypeTag>
90{
91 using type = UndefinedProperty;
92};
93
94template <class TypeTag, class MyTypeTag>
96{
97 using type = UndefinedProperty;
98};
99
100
101template<class TypeTag>
102struct EnableTerminalOutput<TypeTag, TTag::EclFlowProblem> {
103 static constexpr bool value = true;
104};
105template<class TypeTag>
106struct EnableAdaptiveTimeStepping<TypeTag, TTag::EclFlowProblem> {
107 static constexpr bool value = true;
108};
109template<class TypeTag>
110struct EnableTuning<TypeTag, TTag::EclFlowProblem> {
111 static constexpr bool value = false;
112};
113
114template <class TypeTag>
115struct OutputExtraConvergenceInfo<TypeTag, TTag::EclFlowProblem>
116{
117 static constexpr auto* value = "none";
118};
119
120template <class TypeTag>
121struct SaveStep<TypeTag, TTag::EclFlowProblem>
122{
123 static constexpr auto* value = "";
124};
125
126template <class TypeTag>
127struct SaveFile<TypeTag, TTag::EclFlowProblem>
128{
129 static constexpr auto* value = "";
130};
131
132template <class TypeTag>
133struct LoadFile<TypeTag, TTag::EclFlowProblem>
134{
135 static constexpr auto* value = "";
136};
137
138
139template <class TypeTag>
140struct LoadStep<TypeTag, TTag::EclFlowProblem>
141{
142 static constexpr int value = -1;
143};
144
145} // namespace Opm::Properties
146
147namespace Opm {
148
149void outputReportStep(const SimulatorTimer& timer);
150void outputTimestampFIP(const SimulatorTimer& timer,
151 const std::string& title,
152 const std::string& version);
153void checkSerializedCmdLine(const std::string& current,
154 const std::string& stored);
155
157template<class TypeTag>
159{
160public:
171
175
178 typedef typename Model::ModelParameters ModelParameters;
179 typedef typename Solver::SolverParameters SolverParameters;
181
182
204 SimulatorFullyImplicitBlackoilEbos(Simulator& ebosSimulator)
205 : ebosSimulator_(ebosSimulator)
206 {
207 phaseUsage_ = phaseUsageFromDeck(eclState());
208
209 // Only rank 0 does print to std::cout, and only if specifically requested.
210 this->terminalOutput_ = false;
211 if (this->grid().comm().rank() == 0) {
212 this->terminalOutput_ = EWOMS_GET_PARAM(TypeTag, bool, EnableTerminalOutput);
213
214 this->startConvergenceOutputThread(EWOMS_GET_PARAM(TypeTag, std::string,
215 OutputExtraConvergenceInfo),
216 R"(OutputExtraConvergenceInfo (--output-extra-convergence-info))");
217 }
218
219 const std::string saveSpec = EWOMS_GET_PARAM(TypeTag, std::string, SaveStep);
220 if (saveSpec == "all") {
221 saveStride_ = 1;
222 } else if (saveSpec == "last") {
223 saveStride_ = -1;
224 } else if (!saveSpec.empty() && saveSpec[0] == ':') {
225 saveStride_ = std::atoi(saveSpec.c_str()+1);
226 } else if (!saveSpec.empty()) {
227 saveStep_ = std::atoi(saveSpec.c_str());
228 }
229
230 loadStep_ = EWOMS_GET_PARAM(TypeTag, int, LoadStep);
231
232 saveFile_ = EWOMS_GET_PARAM(TypeTag, std::string, SaveFile);
233 loadFile_ = EWOMS_GET_PARAM(TypeTag, std::string, LoadFile);
234
235 if (loadFile_.empty() || saveFile_.empty()) {
236 const auto& ioconfig = ebosSimulator_.vanguard().eclState().getIOConfig();
237 if (saveFile_.empty()) saveFile_ = ioconfig.fullBasePath() + ".OPMRST";
238 if (loadFile_.empty()) loadFile_ = saveFile_;
239 if (loadStep_ != -1 && !std::filesystem::exists(loadFile_)) {
240 std::filesystem::path path(ioconfig.getInputDir() + "/");
241 path.replace_filename(ioconfig.getBaseName() + ".OPMRST");
242 loadFile_ = path;
243 if (!std::filesystem::exists(loadFile_)) {
244 OPM_THROW(std::runtime_error, "Error locating serialized restart file " + loadFile_);
245 }
246 }
247 }
248 }
249
251 {
252 // Safe to call on all ranks, not just the I/O rank.
253 this->endConvergenceOutputThread();
254 }
255
256 static void registerParameters()
257 {
258 ModelParameters::registerParameters();
259 SolverParameters::registerParameters();
260 TimeStepper::registerParameters();
261
262 EWOMS_REGISTER_PARAM(TypeTag, bool, EnableTerminalOutput,
263 "Print high-level information about the simulation's progress to the terminal");
264 EWOMS_REGISTER_PARAM(TypeTag, bool, EnableAdaptiveTimeStepping,
265 "Use adaptive time stepping between report steps");
266 EWOMS_REGISTER_PARAM(TypeTag, bool, EnableTuning,
267 "Honor some aspects of the TUNING keyword.");
268 EWOMS_REGISTER_PARAM(TypeTag, std::string, OutputExtraConvergenceInfo,
269 "Provide additional convergence output "
270 "files for diagnostic purposes. "
271 "\"none\" gives no extra output and "
272 "overrides all other options, "
273 "\"steps\" generates an INFOSTEP file, "
274 "\"iterations\" generates an INFOITER file. "
275 "Combine options with commas, e.g., "
276 "\"steps,iterations\" for multiple outputs.");
277 EWOMS_REGISTER_PARAM(TypeTag, std::string, SaveStep,
278 "Save serialized state to .OPMRST file. "
279 "Either a specific report step, \"all\" to save "
280 "all report steps or \":x\" to save every x'th step."
281 "Use negative values of \"x\" to keep only the last "
282 "written step, or \"last\" to save every step, keeping "
283 "only the last.");
284 EWOMS_REGISTER_PARAM(TypeTag, int, LoadStep,
285 "Load serialized state from .OPMRST file. "
286 "Either a specific report step, or 0 to load last "
287 "stored report step.");
288 EWOMS_REGISTER_PARAM(TypeTag, std::string, SaveFile,
289 "FileName for .OPMRST file used for saving serialized state. "
290 "If empty, CASENAME.OPMRST is used.");
291 EWOMS_HIDE_PARAM(TypeTag, SaveFile);
292 EWOMS_REGISTER_PARAM(TypeTag, std::string, LoadFile,
293 "FileName for .OPMRST file used to load serialized state. "
294 "If empty, CASENAME.OPMRST is used.");
295 EWOMS_HIDE_PARAM(TypeTag, LoadFile);
296 }
297
305 {
306 init(timer);
307 // Make cache up to date. No need for updating it in elementCtx.
308 ebosSimulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
309 // Main simulation loop.
310 while (!timer.done()) {
311 bool continue_looping = runStep(timer);
312 if (!continue_looping) break;
313 }
314 return finalize();
315 }
316
317 void init(SimulatorTimer &timer)
318 {
319 ebosSimulator_.setEpisodeIndex(-1);
320
321 // Create timers and file for writing timing info.
322 solverTimer_ = std::make_unique<time::StopWatch>();
323 totalTimer_ = std::make_unique<time::StopWatch>();
324 totalTimer_->start();
325
326 // adaptive time stepping
327 bool enableAdaptive = EWOMS_GET_PARAM(TypeTag, bool, EnableAdaptiveTimeStepping);
328 bool enableTUNING = EWOMS_GET_PARAM(TypeTag, bool, EnableTuning);
329 if (enableAdaptive) {
330 const UnitSystem& unitSystem = this->ebosSimulator_.vanguard().eclState().getUnits();
331 if (enableTUNING) {
332 const auto& sched_state = schedule()[timer.currentStepNum()];
333 auto max_next_tstep = sched_state.max_next_tstep();
334 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(max_next_tstep,
335 sched_state.tuning(),
336 unitSystem, terminalOutput_);
337 }
338 else {
339 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(unitSystem, terminalOutput_);
340 }
341
342 if (isRestart()) {
343 // For restarts the ebosSimulator may have gotten some information
344 // about the next timestep size from the OPMEXTRA field
345 adaptiveTimeStepping_->setSuggestedNextStep(ebosSimulator_.timeStepSize());
346 }
347 }
348 }
349
350 void updateTUNING(const Tuning& tuning) {
351 modelParam_.tolerance_mb_ = tuning.XXXMBE;
352 if (terminalOutput_) {
353 OpmLog::debug(fmt::format("Setting SimulatorFullyImplicitBlackoilEbos mass balance limit (XXXMBE) to {:.2e}", tuning.XXXMBE));
354 }
355 }
356
357 bool runStep(SimulatorTimer& timer)
358 {
359 if (schedule().exitStatus().has_value()) {
360 if (terminalOutput_) {
361 OpmLog::info("Stopping simulation since EXIT was triggered by an action keyword.");
362 }
363 report_.success.exit_status = schedule().exitStatus().value();
364 return false;
365 }
366
367 if (loadStep_ > -1) {
368 loadTimerInfo(timer);
369 }
370
371 // Report timestep.
372 if (terminalOutput_) {
373 std::ostringstream ss;
374 timer.report(ss);
375 OpmLog::debug(ss.str());
376 }
377
378 if (terminalOutput_) {
379 outputReportStep(timer);
380 }
381
382 // write the inital state at the report stage
383 if (timer.initialStep()) {
384 Dune::Timer perfTimer;
385 perfTimer.start();
386
387 ebosSimulator_.setEpisodeIndex(-1);
388 ebosSimulator_.setEpisodeLength(0.0);
389 ebosSimulator_.setTimeStepSize(0.0);
390 wellModel_().beginReportStep(timer.currentStepNum());
391 ebosSimulator_.problem().writeOutput();
392
393 report_.success.output_write_time += perfTimer.stop();
394 }
395
396 // Run a multiple steps of the solver depending on the time step control.
397 solverTimer_->start();
398
399 if (!solver_) {
400 solver_ = createSolver(wellModel_());
401 }
402
403 ebosSimulator_.startNextEpisode(
404 ebosSimulator_.startTime()
405 + schedule().seconds(timer.currentStepNum()),
406 timer.currentStepLength());
407 ebosSimulator_.setEpisodeIndex(timer.currentStepNum());
408 if (loadStep_> -1) {
409 wellModel_().prepareDeserialize(loadStep_ - 1);
411 loadStep_ = -1;
412 ebosSimulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
413 }
414 solver_->model().beginReportStep();
415 bool enableTUNING = EWOMS_GET_PARAM(TypeTag, bool, EnableTuning);
416
417 // If sub stepping is enabled allow the solver to sub cycle
418 // in case the report steps are too large for the solver to converge
419 //
420 // \Note: The report steps are met in any case
421 // \Note: The sub stepping will require a copy of the state variables
422 if (adaptiveTimeStepping_) {
423 const auto& events = schedule()[timer.currentStepNum()].events();
424 if (enableTUNING) {
425 if (events.hasEvent(ScheduleEvents::TUNING_CHANGE)) {
426 const auto& sched_state = schedule()[timer.currentStepNum()];
427 const auto& tuning = sched_state.tuning();
428 const auto& max_next_tstep = sched_state.max_next_tstep();
429 adaptiveTimeStepping_->updateTUNING(max_next_tstep, tuning);
430 // \Note: Assumes TUNING is only used with adaptive time-stepping
431 // \Note: Need to update both solver (model) and simulator since solver is re-created each report step.
432 solver_->model().updateTUNING(tuning);
433 this->updateTUNING(tuning);
434 }
435 }
436 bool event = events.hasEvent(ScheduleEvents::NEW_WELL) ||
437 events.hasEvent(ScheduleEvents::INJECTION_TYPE_CHANGED) ||
438 events.hasEvent(ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER) ||
439 events.hasEvent(ScheduleEvents::PRODUCTION_UPDATE) ||
440 events.hasEvent(ScheduleEvents::INJECTION_UPDATE) ||
441 events.hasEvent(ScheduleEvents::WELL_STATUS_CHANGE);
442 auto stepReport = adaptiveTimeStepping_->step(timer, *solver_, event, nullptr);
443 report_ += stepReport;
444 //Pass simulation report to eclwriter for summary output
445 ebosSimulator_.problem().setSimulationReport(report_);
446 } else {
447 // solve for complete report step
448 auto stepReport = solver_->step(timer);
449 report_ += stepReport;
450 if (terminalOutput_) {
451 std::ostringstream ss;
452 stepReport.reportStep(ss);
453 OpmLog::info(ss.str());
454 }
455 }
456
457 // write simulation state at the report stage
458 Dune::Timer perfTimer;
459 perfTimer.start();
460 const double nextstep = adaptiveTimeStepping_ ? adaptiveTimeStepping_->suggestedNextStep() : -1.0;
461 ebosSimulator_.problem().setNextTimeStepSize(nextstep);
462 ebosSimulator_.problem().writeOutput();
463 report_.success.output_write_time += perfTimer.stop();
464
465 solver_->model().endReportStep();
466
467 // take time that was used to solve system for this reportStep
468 solverTimer_->stop();
469
470 // update timing.
471 report_.success.solver_time += solverTimer_->secsSinceStart();
472
473 if (this->grid().comm().rank() == 0) {
474 // Grab the step convergence reports that are new since last we were here.
475 const auto& reps = solver_->model().stepReports();
476 this->writeConvergenceOutput(std::vector<StepReport>{reps.begin() + already_reported_steps_, reps.end()});
477 already_reported_steps_ = reps.size();
478 }
479
480 // Increment timer, remember well state.
481 ++timer;
482
483 if (terminalOutput_) {
484 if (!timer.initialStep()) {
485 const std::string version = moduleVersionName();
486 outputTimestampFIP(timer, eclState().getTitle(), version);
487 }
488
489 std::string msg =
490 "Time step took " + std::to_string(solverTimer_->secsSinceStart()) + " seconds; "
491 "total solver time " + std::to_string(report_.success.solver_time) + " seconds.";
492 OpmLog::debug(msg);
493 }
494
495 handleSave(timer);
496
497 return true;
498 }
499
500 SimulatorReport finalize()
501 {
502 // make sure all output is written to disk before run is finished
503 {
504 Dune::Timer finalOutputTimer;
505 finalOutputTimer.start();
506
507 ebosSimulator_.problem().finalizeOutput();
508 report_.success.output_write_time += finalOutputTimer.stop();
509 }
510
511 // Stop timer and create timing report
512 totalTimer_->stop();
513 report_.success.total_time = totalTimer_->secsSinceStart();
514 report_.success.converged = true;
515
516 return report_;
517 }
518
519 const Grid& grid() const
520 { return ebosSimulator_.vanguard().grid(); }
521
522 template<class Serializer>
523 void serializeOp(Serializer& serializer)
524 {
525 serializer(ebosSimulator_);
526 serializer(report_);
527 serializer(adaptiveTimeStepping_);
528 }
529
530 const Model& model() const
531 { return solver_->model(); }
532
533protected:
534
535 std::unique_ptr<Solver> createSolver(WellModel& wellModel)
536 {
537 auto model = std::make_unique<Model>(ebosSimulator_,
538 modelParam_,
539 wellModel,
540 terminalOutput_);
541
542 return std::make_unique<Solver>(solverParam_, std::move(model));
543 }
544
545 const EclipseState& eclState() const
546 { return ebosSimulator_.vanguard().eclState(); }
547
548
549 const Schedule& schedule() const
550 { return ebosSimulator_.vanguard().schedule(); }
551
552 bool isRestart() const
553 {
554 const auto& initconfig = eclState().getInitConfig();
555 return initconfig.restartRequested();
556 }
557
558 WellModel& wellModel_()
559 { return ebosSimulator_.problem().wellModel(); }
560
561 const WellModel& wellModel_() const
562 { return ebosSimulator_.problem().wellModel(); }
563
564 void startConvergenceOutputThread(std::string_view convOutputOptions,
565 std::string_view optionName)
566 {
567 const auto config = ConvergenceOutputConfiguration {
568 convOutputOptions, optionName
569 };
570 if (! config.want(ConvergenceOutputConfiguration::Option::Iterations)) {
571 return;
572 }
573
575 [compNames = typename Model::ComponentName{}](const int compIdx)
576 { return std::string_view { compNames.name(compIdx) }; }
577 };
578
580 [usys = this->eclState().getUnits()](const double time)
581 { return usys.from_si(UnitSystem::measure::time, time); }
582 };
583
584 this->convergenceOutputQueue_.emplace();
585 this->convergenceOutputObject_.emplace
586 (this->eclState().getIOConfig().getOutputDir(),
587 this->eclState().getIOConfig().getBaseName(),
588 std::move(getPhaseName),
589 std::move(convertTime),
590 config, *this->convergenceOutputQueue_);
591
592 this->convergenceOutputThread_
594 &this->convergenceOutputObject_.value());
595 }
596
597 void writeConvergenceOutput(std::vector<StepReport>&& reports)
598 {
599 if (! this->convergenceOutputThread_.has_value()) {
600 return;
601 }
602
603 auto requests = std::vector<ConvergenceReportQueue::OutputRequest>{};
604 requests.reserve(reports.size());
605
606 for (auto&& report : reports) {
607 requests.push_back({ report.report_step, report.current_step, std::move(report.report) });
608 }
609
610 this->convergenceOutputQueue_->enqueue(std::move(requests));
611 }
612
613 void endConvergenceOutputThread()
614 {
615 if (! this->convergenceOutputThread_.has_value()) {
616 return;
617 }
618
619 this->convergenceOutputQueue_->signalLastOutputRequest();
620 this->convergenceOutputThread_->join();
621 }
622
625 {
626 if (saveStride_ == 0 && saveStep_ == -1) {
627 return;
628 }
629
630 OPM_BEGIN_PARALLEL_TRY_CATCH();
631
632 int nextStep = timer.currentStepNum();
633 if ((saveStep_ != -1 && nextStep == saveStep_) ||
634 (saveStride_ != 0 && (nextStep % saveStride_) == 0)) {
635#if !HAVE_HDF5
636 OpmLog::error("Saving of serialized state requested, but no HDF5 support available.");
637#else
638 const std::string groupName = "/report_step/" + std::to_string(nextStep);
639 if (saveStride_ < 0 || nextStep == saveStride_ || nextStep == saveStep_) {
640 std::filesystem::remove(saveFile_);
641 }
644 EclGenericVanguard::comm());
645 if (saveStride_ < 0 || nextStep == saveStride_ || nextStep == saveStep_) {
646 std::ostringstream str;
647 Parameters::printValues<TypeTag>(str);
648 writer.writeHeader("OPM Flow",
651 ebosSimulator_.vanguard().caseName(),
652 str.str(),
653 EclGenericVanguard::comm().size());
654
655 if (EclGenericVanguard::comm().size() > 1) {
656 const auto& cellMapping = ebosSimulator_.vanguard().globalCell();
657 std::size_t hash = Dune::hash_range(cellMapping.begin(), cellMapping.end());
658 writer.write(hash, "/", "grid_checksum");
659 }
660 }
661 writer.write(*this, groupName, "simulator_data");
662 writer.write(timer, groupName, "simulator_timer",
664 OpmLog::info("Serialized state written for report step " + std::to_string(nextStep));
665#endif
666 }
667
668 OPM_END_PARALLEL_TRY_CATCH("Error saving serialized state: ",
669 EclGenericVanguard::comm());
670 }
671
674 {
675#if !HAVE_HDF5
676 OpmLog::error("Loading of serialized state requested, but no HDF5 support available.");
677 loadStep_ = -1;
678#else
679 OPM_BEGIN_PARALLEL_TRY_CATCH();
680
683 EclGenericVanguard::comm());
684
685 if (loadStep_ == 0) {
686 loadStep_ = reader.lastReportStep();
687 }
688
689 OpmLog::info("Loading serialized state for report step " + std::to_string(loadStep_));
690 const std::string groupName = "/report_step/" + std::to_string(loadStep_);
691 reader.read(timer, groupName, "simulator_timer", HDF5File::DataSetMode::ROOT_ONLY);
692
693 std::tuple<std::array<std::string,5>,int> header;
694 reader.read(header, "/", "simulator_info", HDF5File::DataSetMode::ROOT_ONLY);
695 const auto& [strings, procs] = header;
696
697 if (EclGenericVanguard::comm().size() != procs) {
698 throw std::runtime_error("Number of processes (procs=" +
699 std::to_string(EclGenericVanguard::comm().size()) +
700 ") does not match .OPMRST file (procs=" +
701 std::to_string(procs) + ")");
702 }
703
704 if (EclGenericVanguard::comm().size() > 1) {
705 std::size_t stored_hash;
706 reader.read(stored_hash, "/", "grid_checksum");
707 const auto& cellMapping = ebosSimulator_.vanguard().globalCell();
708 std::size_t hash = Dune::hash_range(cellMapping.begin(), cellMapping.end());
709 if (hash != stored_hash) {
710 throw std::runtime_error("Grid hash mismatch, .OPMRST file cannot be used");
711 }
712 }
713
714 if (EclGenericVanguard::comm().rank() == 0) {
715 std::ostringstream str;
716 Parameters::printValues<TypeTag>(str);
717 checkSerializedCmdLine(str.str(), strings[4]);
718 }
719
720 OPM_END_PARALLEL_TRY_CATCH("Error loading serialized state: ",
721 EclGenericVanguard::comm());
722#endif
723 }
724
727 {
728#if HAVE_HDF5
729 OPM_BEGIN_PARALLEL_TRY_CATCH();
730
733 EclGenericVanguard::comm());
734 const std::string groupName = "/report_step/" + std::to_string(loadStep_);
735 reader.read(*this, groupName, "simulator_data");
736
737 OPM_END_PARALLEL_TRY_CATCH("Error loading serialized state: ",
738 EclGenericVanguard::comm());
739#endif
740 }
741
742 // Data.
743 Simulator& ebosSimulator_;
744
745 ModelParameters modelParam_;
746 SolverParameters solverParam_;
747
748 std::unique_ptr<Solver> solver_;
749
750 // Observed objects.
751 PhaseUsage phaseUsage_;
752 // Misc. data
753 bool terminalOutput_;
754
755 SimulatorReport report_;
756 std::size_t already_reported_steps_ = 0;
757 std::unique_ptr<time::StopWatch> solverTimer_;
758 std::unique_ptr<time::StopWatch> totalTimer_;
759 std::unique_ptr<TimeStepper> adaptiveTimeStepping_;
760
761 std::optional<ConvergenceReportQueue> convergenceOutputQueue_{};
762 std::optional<ConvergenceOutputThread> convergenceOutputObject_{};
763 std::optional<std::thread> convergenceOutputThread_{};
764
765 int saveStride_ = 0;
766 int saveStep_ = -1;
767 int loadStep_ = -1;
768 std::string saveFile_;
769 std::string loadFile_;
770};
771
772} // namespace Opm
773
774#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:159
void handleSave(SimulatorTimer &timer)
Serialization of simulator data to .OPMRST files at end of report steps.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:624
std::string saveFile_
File to save serialized state to.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:768
SimulatorReport run(SimulatorTimer &timer)
Run the simulation.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:304
SimulatorFullyImplicitBlackoilEbos(Simulator &ebosSimulator)
Initialise from parameters and objects to observe.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:204
int loadStep_
Step to load serialized state from.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:767
int saveStep_
Specific step to save serialized state at.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:766
void loadSimulatorState()
Load simulator state from serialized state.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:726
int saveStride_
Stride to save serialized state at, negative to only keep last.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:765
void loadTimerInfo(SimulatorTimer &timer)
Load timer info from serialized state.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:673
std::string loadFile_
File to load serialized state from.
Definition SimulatorFullyImplicitBlackoilEbos.hpp:769
Definition SimulatorTimer.hpp:39
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:62
Definition BlackoilWellModel.hpp:85
Definition SimulatorFullyImplicitBlackoilEbos.hpp:66
Definition SimulatorFullyImplicitBlackoilEbos.hpp:96
Definition SimulatorFullyImplicitBlackoilEbos.hpp:84
Definition SimulatorFullyImplicitBlackoilEbos.hpp:72
Definition SimulatorFullyImplicitBlackoilEbos.hpp:90
Definition SimulatorFullyImplicitBlackoilEbos.hpp:78
Definition SimulatorReport.hpp:100