My Project
Loading...
Searching...
No Matches
SimulatorFullyImplicitBlackoil.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_SIMULATOR_FULLY_IMPLICIT_BLACKOIL_HEADER_INCLUDED
23#define OPM_SIMULATOR_FULLY_IMPLICIT_BLACKOIL_HEADER_INCLUDED
24
25#include <opm/common/ErrorMacros.hpp>
26
27#if HAVE_MPI
28#define RESERVOIR_COUPLING_ENABLED
29#endif
30#ifdef RESERVOIR_COUPLING_ENABLED
31#include <opm/input/eclipse/Schedule/ResCoup/ReservoirCouplingInfo.hpp>
32#include <opm/input/eclipse/Schedule/ResCoup/MasterGroup.hpp>
33#include <opm/input/eclipse/Schedule/ResCoup/Slaves.hpp>
34#include <opm/simulators/flow/ReservoirCouplingMaster.hpp>
35#include <opm/simulators/flow/ReservoirCouplingSlave.hpp>
36#include <opm/common/Exceptions.hpp>
37#endif
38
39#include <opm/input/eclipse/Units/UnitSystem.hpp>
40
41#include <opm/grid/utility/StopWatch.hpp>
42
43#include <opm/simulators/aquifers/BlackoilAquiferModel.hpp>
44#include <opm/simulators/flow/BlackoilModel.hpp>
45#include <opm/simulators/flow/BlackoilModelParameters.hpp>
46#include <opm/simulators/flow/ConvergenceOutputConfiguration.hpp>
47#include <opm/simulators/flow/ExtraConvergenceOutputThread.hpp>
48#include <opm/simulators/flow/NonlinearSolver.hpp>
49#include <opm/simulators/flow/SimulatorConvergenceOutput.hpp>
50#include <opm/simulators/flow/SimulatorReportBanners.hpp>
51#include <opm/simulators/flow/SimulatorSerializer.hpp>
52#include <opm/simulators/timestepping/AdaptiveTimeStepping.hpp>
53#include <opm/simulators/timestepping/ConvergenceReport.hpp>
54#include <opm/simulators/utils/moduleVersion.hpp>
55#include <opm/simulators/wells/WellState.hpp>
56
57#if HAVE_HDF5
58#include <opm/simulators/utils/HDF5Serializer.hpp>
59#endif
60
61#include <filesystem>
62#include <memory>
63#include <string>
64#include <string_view>
65#include <utility>
66#include <vector>
67
68#include <fmt/format.h>
69
70namespace Opm::Parameters {
71
72struct EnableAdaptiveTimeStepping { static constexpr bool value = true; };
73struct OutputExtraConvergenceInfo { static constexpr auto* value = "none"; };
74struct SaveStep { static constexpr auto* value = ""; };
75struct SaveFile { static constexpr auto* value = ""; };
76struct LoadFile { static constexpr auto* value = ""; };
77struct LoadStep { static constexpr int value = -1; };
78struct Slave { static constexpr bool value = false; };
79
80} // namespace Opm::Parameters
81
82namespace Opm::detail {
83
84void registerSimulatorParameters();
85
86}
87
88namespace Opm {
89
91template<class TypeTag>
93{
94protected:
95 struct MPI_Comm_Deleter;
96public:
101 using BlackoilIndices = GetPropType<TypeTag, Properties::Indices>;
108
112
113
115 using ModelParameters = typename Model::ModelParameters;
116 using SolverParameters = typename Solver::SolverParameters;
118
140 explicit SimulatorFullyImplicitBlackoil(Simulator& simulator)
141 : simulator_(simulator)
142 , serializer_(*this,
143 FlowGenericVanguard::comm(),
144 simulator_.vanguard().eclState().getIOConfig(),
145 Parameters::Get<Parameters::SaveStep>(),
146 Parameters::Get<Parameters::LoadStep>(),
147 Parameters::Get<Parameters::SaveFile>(),
148 Parameters::Get<Parameters::LoadFile>())
149 {
150 phaseUsage_ = phaseUsageFromDeck(eclState());
151
152 // Only rank 0 does print to std::cout, and only if specifically requested.
153 this->terminalOutput_ = false;
154 if (this->grid().comm().rank() == 0) {
155 this->terminalOutput_ = Parameters::Get<Parameters::EnableTerminalOutput>();
156
158 [compNames = typename Model::ComponentName{}](const int compIdx)
159 { return std::string_view { compNames.name(compIdx) }; }
160 };
161
162 this->convergence_output_.
163 startThread(this->simulator_.vanguard().eclState(),
164 Parameters::Get<Parameters::OutputExtraConvergenceInfo>(),
165 R"(OutputExtraConvergenceInfo (--output-extra-convergence-info))",
167 }
168 }
169
171 {
172 // Safe to call on all ranks, not just the I/O rank.
173 convergence_output_.endThread();
174 }
175
176 static void registerParameters()
177 {
178 ModelParameters::registerParameters();
179 SolverParameters::registerParameters();
180 TimeStepper::registerParameters();
181 detail::registerSimulatorParameters();
182 }
183
190#ifdef RESERVOIR_COUPLING_ENABLED
191 SimulatorReport run(SimulatorTimer& timer, int argc, char** argv)
192 {
193 init(timer, argc, argv);
194#else
196 {
197 init(timer);
198#endif
199 // Make cache up to date. No need for updating it in elementCtx.
200 // NB! Need to be at the correct step in case of restart
201 simulator_.setEpisodeIndex(timer.currentStepNum());
202 simulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
203 // Main simulation loop.
204 while (!timer.done()) {
205 simulator_.problem().writeReports(timer);
206 bool continue_looping = runStep(timer);
207 if (!continue_looping) break;
208 }
209 simulator_.problem().writeReports(timer);
210 return finalize();
211 }
212
213#ifdef RESERVOIR_COUPLING_ENABLED
214 // This method should only be called if slave mode (i.e. Parameters::Get<Parameters::Slave>())
215 // is false. We try to determine if this is a normal flow simulation or a reservoir
216 // coupling master. It is a normal flow simulation if the schedule does not contain
217 // any SLAVES and GRUPMAST keywords.
219 {
220 for (std::size_t report_step = 0; report_step < this->schedule().size(); ++report_step) {
221 auto rescoup = this->schedule()[report_step].rescoup();
222 auto slave_count = rescoup.slaveCount();
223 auto master_group_count = rescoup.masterGroupCount();
224 // - GRUPMAST and SLAVES keywords need to be specified at the same report step
225 // - They can only occur once in the schedule
226 if (slave_count > 0 && master_group_count > 0) {
227 return true;
228 }
229 else if (slave_count > 0 && master_group_count == 0) {
231 "Inconsistent reservoir coupling master schedule: "
232 "Slave count is greater than 0 but master group count is 0"
233 );
234 }
235 else if (slave_count == 0 && master_group_count > 0) {
237 "Inconsistent reservoir coupling master schedule: "
238 "Master group count is greater than 0 but slave count is 0"
239 );
240 }
241 }
242 return false;
243 }
244#endif
245
246#ifdef RESERVOIR_COUPLING_ENABLED
247 // NOTE: The argc and argv will be used when launching a slave process
248 void init(const SimulatorTimer& timer, int argc, char** argv)
249 {
250 auto slave_mode = Parameters::Get<Parameters::Slave>();
251 if (slave_mode) {
253 std::make_unique<ReservoirCouplingSlave>(
255 this->schedule(), timer
256 );
257 this->reservoirCouplingSlave_->sendActivationDateToMasterProcess();
258 this->reservoirCouplingSlave_->sendSimulationStartDateToMasterProcess();
259 this->reservoirCouplingSlave_->receiveMasterGroupNamesFromMasterProcess();
260 }
261 else {
263 if (master_mode) {
265 std::make_unique<ReservoirCouplingMaster>(
267 this->schedule(),
268 argc, argv
269 );
270 }
271 }
272#else
273 void init(const SimulatorTimer& timer)
274 {
275#endif
276 simulator_.setEpisodeIndex(-1);
277
278 // Create timers and file for writing timing info.
279 solverTimer_ = std::make_unique<time::StopWatch>();
280 totalTimer_ = std::make_unique<time::StopWatch>();
281 totalTimer_->start();
282
283 // adaptive time stepping
284 bool enableAdaptive = Parameters::Get<Parameters::EnableAdaptiveTimeStepping>();
285 bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
286 if (enableAdaptive) {
287 const UnitSystem& unitSystem = this->simulator_.vanguard().eclState().getUnits();
288 const auto& sched_state = schedule()[timer.currentStepNum()];
289 auto max_next_tstep = sched_state.max_next_tstep(enableTUNING);
290 if (enableTUNING) {
291 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(max_next_tstep,
292 sched_state.tuning(),
293 unitSystem, report_, terminalOutput_);
294 }
295 else {
296 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(unitSystem, report_, max_next_tstep, terminalOutput_);
297 }
298#ifdef RESERVOIR_COUPLING_ENABLED
299 if (this->reservoirCouplingSlave_) {
300 adaptiveTimeStepping_->setReservoirCouplingSlave(this->reservoirCouplingSlave_.get());
301 }
302 else if (this->reservoirCouplingMaster_) {
303 adaptiveTimeStepping_->setReservoirCouplingMaster(this->reservoirCouplingMaster_.get());
304 }
305#endif
306 if (isRestart()) {
307 // For restarts the simulator may have gotten some information
308 // about the next timestep size from the OPMEXTRA field
309 adaptiveTimeStepping_->setSuggestedNextStep(simulator_.timeStepSize());
310 }
311 }
312 }
313
314 void updateTUNING(const Tuning& tuning)
315 {
316 modelParam_.tolerance_mb_ = tuning.XXXMBE;
317 if (terminalOutput_) {
318 OpmLog::debug(fmt::format("Setting SimulatorFullyImplicitBlackoil mass balance limit (XXXMBE) to {:.2e}", tuning.XXXMBE));
319 }
320 }
321
322 bool runStep(SimulatorTimer& timer)
323 {
324 if (schedule().exitStatus().has_value()) {
325 if (terminalOutput_) {
326 OpmLog::info("Stopping simulation since EXIT was triggered by an action keyword.");
327 }
328 report_.success.exit_status = schedule().exitStatus().value();
329 return false;
330 }
331
332 if (serializer_.shouldLoad()) {
333 serializer_.loadTimerInfo(timer);
334 }
335
336 // Report timestep.
337 if (terminalOutput_) {
338 std::ostringstream ss;
339 timer.report(ss);
340 OpmLog::debug(ss.str());
341 details::outputReportStep(timer);
342 }
343
344 // write the inital state at the report stage
345 if (timer.initialStep()) {
346 Dune::Timer perfTimer;
347 perfTimer.start();
348
349 simulator_.setEpisodeIndex(-1);
350 simulator_.setEpisodeLength(0.0);
351 simulator_.setTimeStepSize(0.0);
352 wellModel_().beginReportStep(timer.currentStepNum());
353 simulator_.problem().writeOutput(true);
354
355 report_.success.output_write_time += perfTimer.stop();
356 }
357
358 // Run a multiple steps of the solver depending on the time step control.
359 solverTimer_->start();
360
361 if (!solver_) {
362 solver_ = createSolver(wellModel_());
363 }
364
365 simulator_.startNextEpisode(
366 simulator_.startTime()
367 + schedule().seconds(timer.currentStepNum()),
368 timer.currentStepLength());
369 simulator_.setEpisodeIndex(timer.currentStepNum());
370
371 if (serializer_.shouldLoad()) {
372 wellModel_().prepareDeserialize(serializer_.loadStep() - 1);
373 serializer_.loadState();
374 simulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
375 }
376
377 this->solver_->model().beginReportStep();
378
379 const bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
380
381 // If sub stepping is enabled allow the solver to sub cycle
382 // in case the report steps are too large for the solver to converge
383 //
384 // \Note: The report steps are met in any case
385 // \Note: The sub stepping will require a copy of the state variables
386 if (adaptiveTimeStepping_) {
387 auto tuningUpdater = [enableTUNING, this,
388 reportStep = timer.currentStepNum()](const double curr_time,
389 double dt, const int timeStep)
390 {
391 auto& schedule = this->simulator_.vanguard().schedule();
392 auto& events = this->schedule()[reportStep].events();
393
394 bool result = false;
395 if (events.hasEvent(ScheduleEvents::TUNING_CHANGE)) {
396 // Unset the event to not trigger it again on the next sub step
397 schedule.clear_event(ScheduleEvents::TUNING_CHANGE, reportStep);
398 const auto& sched_state = schedule[reportStep];
399 const auto& max_next_tstep = sched_state.max_next_tstep(enableTUNING);
400 const auto& tuning = sched_state.tuning();
401
402 if (enableTUNING) {
403 adaptiveTimeStepping_->updateTUNING(max_next_tstep, tuning);
404 // \Note: Assumes TUNING is only used with adaptive time-stepping
405 // \Note: Need to update both solver (model) and simulator since solver is re-created each report step.
406 solver_->model().updateTUNING(tuning);
407 this->updateTUNING(tuning);
408 dt = this->adaptiveTimeStepping_->suggestedNextStep();
409 } else {
411 this->adaptiveTimeStepping_->updateNEXTSTEP(max_next_tstep);
412 }
414 }
415
416 const auto& wcycle = schedule[reportStep].wcycle.get();
417 if (wcycle.empty()) {
418 return result;
419 }
420
421 const auto& wmatcher = schedule.wellMatcher(reportStep);
422 double wcycle_time_step =
423 wcycle.nextTimeStep(curr_time,
424 dt,
425 wmatcher,
426 this->wellModel_().wellOpenTimes(),
427 this->wellModel_().wellCloseTimes(),
428 [timeStep,
429 &wg_events = this->wellModel_().reportStepStartEvents()]
430 (const std::string& name)
431 {
432 if (timeStep != 0) {
433 return false;
434 }
435 return wg_events.hasEvent(name, ScheduleEvents::REQUEST_OPEN_WELL);
436 });
437
438 wcycle_time_step = this->grid().comm().min(wcycle_time_step);
439 if (dt != wcycle_time_step) {
440 this->adaptiveTimeStepping_->updateNEXTSTEP(wcycle_time_step);
441 return true;
442 }
443
444 return result;
445 };
446
447 tuningUpdater(timer.simulationTimeElapsed(),
448 this->adaptiveTimeStepping_->suggestedNextStep(), 0);
449
450#ifdef RESERVOIR_COUPLING_ENABLED
451 if (this->reservoirCouplingMaster_) {
452 this->reservoirCouplingMaster_->maybeSpawnSlaveProcesses(timer.currentStepNum());
453 }
454 else if (this->reservoirCouplingSlave_) {
455 this->reservoirCouplingSlave_->maybeActivate(timer.currentStepNum());
456 }
457#endif
458 const auto& events = schedule()[timer.currentStepNum()].events();
459 bool event = events.hasEvent(ScheduleEvents::NEW_WELL) ||
460 events.hasEvent(ScheduleEvents::INJECTION_TYPE_CHANGED) ||
461 events.hasEvent(ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER) ||
462 events.hasEvent(ScheduleEvents::PRODUCTION_UPDATE) ||
463 events.hasEvent(ScheduleEvents::INJECTION_UPDATE) ||
464 events.hasEvent(ScheduleEvents::WELL_STATUS_CHANGE);
465 auto stepReport = adaptiveTimeStepping_->step(timer, *solver_, event, tuningUpdater);
466 report_ += stepReport;
467 //Pass simulation report to eclwriter for summary output
468 simulator_.problem().setSimulationReport(report_);
469 } else {
470 // solve for complete report step
471 auto stepReport = solver_->step(timer);
472 report_ += stepReport;
473 if (terminalOutput_) {
474 std::ostringstream ss;
475 stepReport.reportStep(ss);
476 OpmLog::info(ss.str());
477 }
478 }
479
480 // write simulation state at the report stage
481 Dune::Timer perfTimer;
482 perfTimer.start();
483 const double nextstep = adaptiveTimeStepping_ ? adaptiveTimeStepping_->suggestedNextStep() : -1.0;
484 simulator_.problem().setNextTimeStepSize(nextstep);
485 simulator_.problem().writeOutput(true);
486 report_.success.output_write_time += perfTimer.stop();
487
488 solver_->model().endReportStep();
489
490 // take time that was used to solve system for this reportStep
491 solverTimer_->stop();
492
493 // update timing.
494 report_.success.solver_time += solverTimer_->secsSinceStart();
495
496 if (this->grid().comm().rank() == 0) {
497 // Grab the step convergence reports that are new since last we
498 // were here.
499 const auto& reps = this->solver_->model().stepReports();
500 convergence_output_.write(reps);
501 }
502
503 // Increment timer, remember well state.
504 ++timer;
505
506 if (terminalOutput_) {
507 std::string msg =
508 "Time step took " + std::to_string(solverTimer_->secsSinceStart()) + " seconds; "
509 "total solver time " + std::to_string(report_.success.solver_time) + " seconds.";
510 OpmLog::debug(msg);
511 }
512
513 serializer_.save(timer);
514
515 return true;
516 }
517
518 SimulatorReport finalize()
519 {
520 // make sure all output is written to disk before run is finished
521 {
522 Dune::Timer finalOutputTimer;
523 finalOutputTimer.start();
524
525 simulator_.problem().finalizeOutput();
526 report_.success.output_write_time += finalOutputTimer.stop();
527 }
528
529 // Stop timer and create timing report
530 totalTimer_->stop();
531 report_.success.total_time = totalTimer_->secsSinceStart();
532 report_.success.converged = true;
533
534 return report_;
535 }
536
537 const Grid& grid() const
538 { return simulator_.vanguard().grid(); }
539
540 template<class Serializer>
541 void serializeOp(Serializer& serializer)
542 {
543 serializer(simulator_);
544 serializer(report_);
545 serializer(adaptiveTimeStepping_);
546 }
547
548 const Model& model() const
549 { return solver_->model(); }
550
551protected:
554 [[maybe_unused]] const std::string& groupName) override
555 {
556#if HAVE_HDF5
557 serializer.read(*this, groupName, "simulator_data");
558#endif
559 }
560
563 [[maybe_unused]] const std::string& groupName) const override
564 {
565#if HAVE_HDF5
566 serializer.write(*this, groupName, "simulator_data");
567#endif
568 }
569
571 std::array<std::string,5> getHeader() const override
572 {
573 std::ostringstream str;
574 Parameters::printValues(str);
575 return {"OPM Flow",
578 simulator_.vanguard().caseName(),
579 str.str()};
580 }
581
583 const std::vector<int>& getCellMapping() const override
584 {
585 return simulator_.vanguard().globalCell();
586 }
587
588 std::unique_ptr<Solver> createSolver(WellModel& wellModel)
589 {
590 auto model = std::make_unique<Model>(simulator_,
591 modelParam_,
592 wellModel,
593 terminalOutput_);
594
595 if (this->modelParam_.write_partitions_) {
596 const auto& iocfg = this->eclState().cfg().io();
597
598 const auto odir = iocfg.getOutputDir()
599 / std::filesystem::path { "partition" }
600 / iocfg.getBaseName();
601
602 if (this->grid().comm().rank() == 0) {
604 }
605
606 this->grid().comm().barrier();
607
608 model->writePartitions(odir);
609
610 this->modelParam_.write_partitions_ = false;
611 }
612
613 return std::make_unique<Solver>(solverParam_, std::move(model));
614 }
615
616 const EclipseState& eclState() const
617 { return simulator_.vanguard().eclState(); }
618
619
620 const Schedule& schedule() const
621 { return simulator_.vanguard().schedule(); }
622
623 bool isRestart() const
624 {
625 const auto& initconfig = eclState().getInitConfig();
626 return initconfig.restartRequested();
627 }
628
629 WellModel& wellModel_()
630 { return simulator_.problem().wellModel(); }
631
632 const WellModel& wellModel_() const
633 { return simulator_.problem().wellModel(); }
634
635 // Data.
636 Simulator& simulator_;
637
638 ModelParameters modelParam_;
639 SolverParameters solverParam_;
640
641 std::unique_ptr<Solver> solver_;
642
643 // Observed objects.
644 PhaseUsage phaseUsage_;
645 // Misc. data
646 bool terminalOutput_;
647
648 SimulatorReport report_;
649 std::unique_ptr<time::StopWatch> solverTimer_;
650 std::unique_ptr<time::StopWatch> totalTimer_;
651 std::unique_ptr<TimeStepper> adaptiveTimeStepping_;
652
653 SimulatorConvergenceOutput convergence_output_{};
654
655#ifdef RESERVOIR_COUPLING_ENABLED
656 bool slaveMode_{false};
657 std::unique_ptr<ReservoirCouplingMaster> reservoirCouplingMaster_{nullptr};
658 std::unique_ptr<ReservoirCouplingSlave> reservoirCouplingSlave_{nullptr};
659#endif
660
661 SimulatorSerializer serializer_;
662};
663
664} // namespace Opm
665
666#endif // OPM_SIMULATOR_FULLY_IMPLICIT_BLACKOIL_HEADER_INCLUDED
Definition AdaptiveTimeStepping.hpp:81
Contains the high level supplements required to extend the black oil model by MICP.
Definition blackoilmicpmodules.hh:49
Contains the high level supplements required to extend the black oil model by polymer.
Definition blackoilpolymermodules.hh:54
Class for handling the blackoil well model.
Definition BlackoilWellModel.hpp:94
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
Definition FlowGenericVanguard.hpp:107
static Parallel::Communication & comm()
Obtain global communicator.
Definition FlowGenericVanguard.hpp:332
Class for (de-)serializing using HDF5.
Definition HDF5Serializer.hpp:37
A nonlinear solver class suitable for general fully-implicit models, as well as pressure,...
Definition NonlinearSolver.hpp:98
void endThread()
Request that convergence output thread be shut down.
Definition SimulatorConvergenceOutput.cpp:100
void write(const std::vector< StepReport > &reports)
Create convergence output requests.
Definition SimulatorConvergenceOutput.cpp:72
a simulator for the blackoil model
Definition SimulatorFullyImplicitBlackoil.hpp:93
void loadState(HDF5Serializer &serializer, const std::string &groupName) override
Load simulator state from hdf5 serializer.
Definition SimulatorFullyImplicitBlackoil.hpp:553
SimulatorReport run(SimulatorTimer &timer)
Run the simulation.
Definition SimulatorFullyImplicitBlackoil.hpp:195
const std::vector< int > & getCellMapping() const override
Returns local-to-global cell mapping.
Definition SimulatorFullyImplicitBlackoil.hpp:583
void saveState(HDF5Serializer &serializer, const std::string &groupName) const override
Save simulator state using hdf5 serializer.
Definition SimulatorFullyImplicitBlackoil.hpp:562
SimulatorFullyImplicitBlackoil(Simulator &simulator)
Initialise from parameters and objects to observe.
Definition SimulatorFullyImplicitBlackoil.hpp:140
std::array< std::string, 5 > getHeader() const override
Returns header data.
Definition SimulatorFullyImplicitBlackoil.hpp:571
void loadState()
Load state from file.
Definition SimulatorSerializer.cpp:176
void loadTimerInfo(SimulatorTimer &timer)
Loads time step info from file.
Definition SimulatorSerializer.cpp:130
bool shouldLoad() const
Returns whether or not a state should be loaded.
Definition SimulatorSerializer.hpp:71
void save(SimulatorTimer &timer)
Save data to file if appropriate.
Definition SimulatorSerializer.cpp:92
int loadStep() const
Returns step to load.
Definition SimulatorSerializer.hpp:74
Definition SimulatorTimer.hpp:39
int currentStepNum() const override
Current step number.
Definition SimulatorTimer.cpp:95
bool done() const override
Return true if op++() has been called numSteps() times.
Definition SimulatorTimer.cpp:168
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
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
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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
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
PhaseUsage phaseUsageFromDeck(const EclipseState &eclipseState)
Looks at presence of WATER, OIL and GAS keywords in state object to determine active phases.
Definition phaseUsageFromDeck.cpp:137
Definition NonlinearSolver.hpp:79
Definition SimulatorFullyImplicitBlackoil.hpp:72
Definition SimulatorFullyImplicitBlackoil.hpp:76
Definition SimulatorFullyImplicitBlackoil.hpp:77
Definition SimulatorFullyImplicitBlackoil.hpp:73
Definition SimulatorFullyImplicitBlackoil.hpp:75
Definition SimulatorFullyImplicitBlackoil.hpp:74
Definition SimulatorFullyImplicitBlackoil.hpp:78
Abstract interface for simulator serialization ops.
Definition SimulatorSerializer.hpp:36
Definition SimulatorReport.hpp:100