SimulatorFullyImplicitBlackoil.hpp
Go to the documentation of this file.
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>
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
56
57#if HAVE_HDF5
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
85
86}
87
88namespace Opm {
89
91template<class TypeTag>
93{
94protected:
95 struct MPI_Comm_Deleter;
96public:
108
112
113
115 using ModelParameters = typename Model::ModelParameters;
118
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 {
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
163 startThread(this->simulator_.vanguard().eclState(),
164 Parameters::Get<Parameters::OutputExtraConvergenceInfo>(),
165 R"(OutputExtraConvergenceInfo (--output-extra-convergence-info))",
166 getPhaseName);
167 }
168 }
169
171 {
172 // Safe to call on all ranks, not just the I/O rank.
174 }
175
176 static void registerParameters()
177 {
178 ModelParameters::registerParameters();
179 SolverParameters::registerParameters();
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) {
230 throw ReservoirCouplingError(
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) {
236 throw ReservoirCouplingError(
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_->sendAndReceiveInitialData();
258 this->simulator_.setReservoirCouplingSlave(this->reservoirCouplingSlave_.get());
260 }
261 else {
262 auto master_mode = checkRunningAsReservoirCouplingMaster();
263 if (master_mode) {
265 std::make_unique<ReservoirCouplingMaster>(
267 this->schedule(),
268 argc, argv
269 );
270 this->simulator_.setReservoirCouplingMaster(this->reservoirCouplingMaster_.get());
272 }
273 }
274#else
275 void init(const SimulatorTimer& timer)
276 {
277#endif
278 simulator_.setEpisodeIndex(-1);
279
280 // Create timers and file for writing timing info.
281 solverTimer_ = std::make_unique<time::StopWatch>();
282 totalTimer_ = std::make_unique<time::StopWatch>();
283 totalTimer_->start();
284
285 // adaptive time stepping
286 bool enableAdaptive = Parameters::Get<Parameters::EnableAdaptiveTimeStepping>();
287 bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
288 if (enableAdaptive) {
289 const UnitSystem& unitSystem = this->simulator_.vanguard().eclState().getUnits();
290 const auto& sched_state = schedule()[timer.currentStepNum()];
291 auto max_next_tstep = sched_state.max_next_tstep(enableTUNING);
292 if (enableTUNING) {
293 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(max_next_tstep,
294 sched_state.tuning(),
295 unitSystem, report_, terminalOutput_);
296 }
297 else {
298 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(unitSystem, report_, max_next_tstep, terminalOutput_);
299 }
300 if (isRestart()) {
301 // For restarts the simulator may have gotten some information
302 // about the next timestep size from the OPMEXTRA field
303 adaptiveTimeStepping_->setSuggestedNextStep(simulator_.timeStepSize());
304 }
305 }
306 }
307
308 void updateTUNING(const Tuning& tuning)
309 {
310 modelParam_.tolerance_mb_ = tuning.XXXMBE;
311 if (terminalOutput_) {
312 OpmLog::debug(fmt::format("Setting SimulatorFullyImplicitBlackoil mass balance limit (XXXMBE) to {:.2e}", tuning.XXXMBE));
313 }
314 }
315
317 {
318 if (schedule().exitStatus().has_value()) {
319 if (terminalOutput_) {
320 OpmLog::info("Stopping simulation since EXIT was triggered by an action keyword.");
321 }
322 report_.success.exit_status = schedule().exitStatus().value();
323 return false;
324 }
325
326 if (serializer_.shouldLoad()) {
328 }
329
330 // Report timestep.
331 if (terminalOutput_) {
332 std::ostringstream ss;
333 timer.report(ss);
334 OpmLog::debug(ss.str());
336 }
337
338 // write the inital state at the report stage
339 if (timer.initialStep()) {
340 Dune::Timer perfTimer;
341 perfTimer.start();
342
343 simulator_.setEpisodeIndex(-1);
344 simulator_.setEpisodeLength(0.0);
345 simulator_.setTimeStepSize(0.0);
347 simulator_.problem().writeOutput(true);
348
349 report_.success.output_write_time += perfTimer.stop();
350 }
351
352 // Run a multiple steps of the solver depending on the time step control.
353 solverTimer_->start();
354
355 if (!solver_) {
357 }
358
359 simulator_.startNextEpisode(
360 simulator_.startTime()
361 + schedule().seconds(timer.currentStepNum()),
362 timer.currentStepLength());
363 simulator_.setEpisodeIndex(timer.currentStepNum());
364
365 if (serializer_.shouldLoad()) {
368 simulator_.model().invalidateAndUpdateIntensiveQuantities(/*timeIdx=*/0);
369 }
370
371 this->solver_->model().beginReportStep();
372
373 const bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
374
375 // If sub stepping is enabled allow the solver to sub cycle
376 // in case the report steps are too large for the solver to converge
377 //
378 // \Note: The report steps are met in any case
379 // \Note: The sub stepping will require a copy of the state variables
381 auto tuningUpdater = [enableTUNING, this,
382 reportStep = timer.currentStepNum()](const double curr_time,
383 double dt, const int timeStep)
384 {
385 auto& schedule = this->simulator_.vanguard().schedule();
386 auto& events = this->schedule()[reportStep].events();
387
388 bool result = false;
389 if (events.hasEvent(ScheduleEvents::TUNING_CHANGE)) {
390 // Unset the event to not trigger it again on the next sub step
391 schedule.clear_event(ScheduleEvents::TUNING_CHANGE, reportStep);
392 const auto& sched_state = schedule[reportStep];
393 const auto& max_next_tstep = sched_state.max_next_tstep(enableTUNING);
394 const auto& tuning = sched_state.tuning();
395
396 if (enableTUNING) {
397 adaptiveTimeStepping_->updateTUNING(max_next_tstep, tuning);
398 // \Note: Assumes TUNING is only used with adaptive time-stepping
399 // \Note: Need to update both solver (model) and simulator since solver is re-created each report step.
400 solver_->model().updateTUNING(tuning);
401 this->updateTUNING(tuning);
402 dt = this->adaptiveTimeStepping_->suggestedNextStep();
403 } else {
404 dt = max_next_tstep;
405 this->adaptiveTimeStepping_->updateNEXTSTEP(max_next_tstep);
406 }
407 result = max_next_tstep > 0;
408 }
409
410 const auto& wcycle = schedule[reportStep].wcycle.get();
411 if (wcycle.empty()) {
412 return result;
413 }
414
415 const auto& wmatcher = schedule.wellMatcher(reportStep);
416 double wcycle_time_step =
417 wcycle.nextTimeStep(curr_time,
418 dt,
419 wmatcher,
420 this->wellModel_().wellOpenTimes(),
421 this->wellModel_().wellCloseTimes(),
422 [timeStep,
423 &wg_events = this->wellModel_().reportStepStartEvents()]
424 (const std::string& name)
425 {
426 if (timeStep != 0) {
427 return false;
428 }
429 return wg_events.hasEvent(name, ScheduleEvents::REQUEST_OPEN_WELL);
430 });
431
432 wcycle_time_step = this->grid().comm().min(wcycle_time_step);
433 if (dt != wcycle_time_step) {
434 this->adaptiveTimeStepping_->updateNEXTSTEP(wcycle_time_step);
435 return true;
436 }
437
438 return result;
439 };
440
441 tuningUpdater(timer.simulationTimeElapsed(),
442 this->adaptiveTimeStepping_->suggestedNextStep(), 0);
443
444#ifdef RESERVOIR_COUPLING_ENABLED
445 if (this->reservoirCouplingMaster_) {
446 this->reservoirCouplingMaster_->maybeSpawnSlaveProcesses(timer.currentStepNum());
447 }
448 else if (this->reservoirCouplingSlave_) {
449 this->reservoirCouplingSlave_->maybeActivate(timer.currentStepNum());
450 }
451#endif
452 const auto& events = schedule()[timer.currentStepNum()].events();
453 bool event = events.hasEvent(ScheduleEvents::NEW_WELL) ||
454 events.hasEvent(ScheduleEvents::INJECTION_TYPE_CHANGED) ||
455 events.hasEvent(ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER) ||
456 events.hasEvent(ScheduleEvents::PRODUCTION_UPDATE) ||
457 events.hasEvent(ScheduleEvents::INJECTION_UPDATE) ||
458 events.hasEvent(ScheduleEvents::WELL_STATUS_CHANGE);
459 auto stepReport = adaptiveTimeStepping_->step(timer, *solver_, event, tuningUpdater);
460 report_ += stepReport;
461 //Pass simulation report to eclwriter for summary output
462 simulator_.problem().setSimulationReport(report_);
463 } else {
464 // solve for complete report step
465 auto stepReport = solver_->step(timer, nullptr);
466 report_ += stepReport;
467 if (terminalOutput_) {
468 std::ostringstream ss;
469 stepReport.reportStep(ss);
470 OpmLog::info(ss.str());
471 }
472 }
473
474 // write simulation state at the report stage
475 Dune::Timer perfTimer;
476 perfTimer.start();
477 const double nextstep = adaptiveTimeStepping_ ? adaptiveTimeStepping_->suggestedNextStep() : -1.0;
478 simulator_.problem().setNextTimeStepSize(nextstep);
479 simulator_.problem().writeOutput(true);
480 report_.success.output_write_time += perfTimer.stop();
481
482 solver_->model().endReportStep();
483
484 // take time that was used to solve system for this reportStep
485 solverTimer_->stop();
486
487 // update timing.
488 report_.success.solver_time += solverTimer_->secsSinceStart();
489
490 if (this->grid().comm().rank() == 0) {
491 // Grab the step convergence reports that are new since last we
492 // were here.
493 const auto& reps = this->solver_->model().stepReports();
495 }
496
497 // Increment timer, remember well state.
498 ++timer;
499
500 if (terminalOutput_) {
501 std::string msg =
502 "Time step took " + std::to_string(solverTimer_->secsSinceStart()) + " seconds; "
503 "total solver time " + std::to_string(report_.success.solver_time) + " seconds.";
504 OpmLog::debug(msg);
505 }
506
507 serializer_.save(timer);
508
509 return true;
510 }
511
513 {
514 // make sure all output is written to disk before run is finished
515 {
516 Dune::Timer finalOutputTimer;
517 finalOutputTimer.start();
518
519 simulator_.problem().finalizeOutput();
520 report_.success.output_write_time += finalOutputTimer.stop();
521 }
522
523 // Stop timer and create timing report
524 totalTimer_->stop();
525 report_.success.total_time = totalTimer_->secsSinceStart();
527
528 return report_;
529 }
530
531 const Grid& grid() const
532 { return simulator_.vanguard().grid(); }
533
534 template<class Serializer>
535 void serializeOp(Serializer& serializer)
536 {
537 serializer(simulator_);
538 serializer(report_);
539 serializer(adaptiveTimeStepping_);
540 }
541
542 const Model& model() const
543 { return solver_->model(); }
544
545protected:
547 void loadState([[maybe_unused]] HDF5Serializer& serializer,
548 [[maybe_unused]] const std::string& groupName) override
549 {
550#if HAVE_HDF5
551 serializer.read(*this, groupName, "simulator_data");
552#endif
553 }
554
556 void saveState([[maybe_unused]] HDF5Serializer& serializer,
557 [[maybe_unused]] const std::string& groupName) const override
558 {
559#if HAVE_HDF5
560 serializer.write(*this, groupName, "simulator_data");
561#endif
562 }
563
565 std::array<std::string,5> getHeader() const override
566 {
567 std::ostringstream str;
569 return {"OPM Flow",
572 simulator_.vanguard().caseName(),
573 str.str()};
574 }
575
577 const std::vector<int>& getCellMapping() const override
578 {
579 return simulator_.vanguard().globalCell();
580 }
581
582 std::unique_ptr<Solver> createSolver(WellModel& wellModel)
583 {
584 auto model = std::make_unique<Model>(simulator_,
586 wellModel,
588
589 if (this->modelParam_.write_partitions_) {
590 const auto& iocfg = this->eclState().cfg().io();
591
592 const auto odir = iocfg.getOutputDir()
593 / std::filesystem::path { "partition" }
594 / iocfg.getBaseName();
595
596 if (this->grid().comm().rank() == 0) {
597 create_directories(odir);
598 }
599
600 this->grid().comm().barrier();
601
602 model->writePartitions(odir);
603
604 this->modelParam_.write_partitions_ = false;
605 }
606
607 return std::make_unique<Solver>(solverParam_, std::move(model));
608 }
609
610 const EclipseState& eclState() const
611 { return simulator_.vanguard().eclState(); }
612
613
614 const Schedule& schedule() const
615 { return simulator_.vanguard().schedule(); }
616
617 bool isRestart() const
618 {
619 const auto& initconfig = eclState().getInitConfig();
620 return initconfig.restartRequested();
621 }
622
624 { return simulator_.problem().wellModel(); }
625
626 const WellModel& wellModel_() const
627 { return simulator_.problem().wellModel(); }
628
629 // Data.
631
634
635 std::unique_ptr<Solver> solver_;
636
637 // Observed objects.
639 // Misc. data
641
643 std::unique_ptr<time::StopWatch> solverTimer_;
644 std::unique_ptr<time::StopWatch> totalTimer_;
645 std::unique_ptr<TimeStepper> adaptiveTimeStepping_;
646
648
649#ifdef RESERVOIR_COUPLING_ENABLED
650 bool slaveMode_{false};
651 std::unique_ptr<ReservoirCouplingMaster> reservoirCouplingMaster_{nullptr};
652 std::unique_ptr<ReservoirCouplingSlave> reservoirCouplingSlave_{nullptr};
653#endif
654
656};
657
658} // namespace Opm
659
660#endif // OPM_SIMULATOR_FULLY_IMPLICIT_BLACKOIL_HEADER_INCLUDED
Definition: AdaptiveTimeStepping.hpp:78
static void registerParameters()
Definition: AdaptiveTimeStepping_impl.hpp:180
Contains the high level supplements required to extend the black oil model by MICP.
Definition: blackoilmicpmodules.hh:54
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:102
void setReservoirCouplingSlave(ReservoirCouplingSlave *slave)
Definition: BlackoilWellModel.hpp:393
void prepareDeserialize(const int report_step)
Definition: BlackoilWellModel.hpp:201
void beginReportStep(const int time_step)
Definition: BlackoilWellModel_impl.hpp:201
void setReservoirCouplingMaster(ReservoirCouplingMaster *master)
Definition: BlackoilWellModel.hpp:389
std::function< std::string_view(int)> ComponentToPhaseName
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
Definition: NonlinearSolver.hpp:99
NonlinearSolverParameters< Scalar > SolverParameters
Definition: NonlinearSolver.hpp:103
Class handling convergence history output for a simulator.
Definition: SimulatorConvergenceOutput.hpp:44
void write(const std::vector< StepReport > &reports)
a simulator for the blackoil model
Definition: SimulatorFullyImplicitBlackoil.hpp:93
SolverParameters solverParam_
Definition: SimulatorFullyImplicitBlackoil.hpp:633
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: SimulatorFullyImplicitBlackoil.hpp:97
void loadState(HDF5Serializer &serializer, const std::string &groupName) override
Load simulator state from hdf5 serializer.
Definition: SimulatorFullyImplicitBlackoil.hpp:547
PhaseUsage phaseUsage_
Definition: SimulatorFullyImplicitBlackoil.hpp:638
const WellModel & wellModel_() const
Definition: SimulatorFullyImplicitBlackoil.hpp:626
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: SimulatorFullyImplicitBlackoil.hpp:99
std::unique_ptr< TimeStepper > adaptiveTimeStepping_
Definition: SimulatorFullyImplicitBlackoil.hpp:645
bool slaveMode_
Definition: SimulatorFullyImplicitBlackoil.hpp:650
GetPropType< TypeTag, Properties::MaterialLawParams > MaterialLawParams
Definition: SimulatorFullyImplicitBlackoil.hpp:105
std::unique_ptr< Solver > createSolver(WellModel &wellModel)
Definition: SimulatorFullyImplicitBlackoil.hpp:582
bool isRestart() const
Definition: SimulatorFullyImplicitBlackoil.hpp:617
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: SimulatorFullyImplicitBlackoil.hpp:103
typename Model::ModelParameters ModelParameters
Definition: SimulatorFullyImplicitBlackoil.hpp:115
void updateTUNING(const Tuning &tuning)
Definition: SimulatorFullyImplicitBlackoil.hpp:308
~SimulatorFullyImplicitBlackoil() override
Definition: SimulatorFullyImplicitBlackoil.hpp:170
void serializeOp(Serializer &serializer)
Definition: SimulatorFullyImplicitBlackoil.hpp:535
typename Solver::SolverParameters SolverParameters
Definition: SimulatorFullyImplicitBlackoil.hpp:116
GetPropType< TypeTag, Properties::NonlinearSystem > Model
Definition: SimulatorFullyImplicitBlackoil.hpp:107
SimulatorSerializer serializer_
Definition: SimulatorFullyImplicitBlackoil.hpp:655
const std::vector< int > & getCellMapping() const override
Returns local-to-global cell mapping.
Definition: SimulatorFullyImplicitBlackoil.hpp:577
std::unique_ptr< ReservoirCouplingSlave > reservoirCouplingSlave_
Definition: SimulatorFullyImplicitBlackoil.hpp:652
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition: SimulatorFullyImplicitBlackoil.hpp:102
const EclipseState & eclState() const
Definition: SimulatorFullyImplicitBlackoil.hpp:610
std::unique_ptr< time::StopWatch > totalTimer_
Definition: SimulatorFullyImplicitBlackoil.hpp:644
void saveState(HDF5Serializer &serializer, const std::string &groupName) const override
Save simulator state using hdf5 serializer.
Definition: SimulatorFullyImplicitBlackoil.hpp:556
const Grid & grid() const
Definition: SimulatorFullyImplicitBlackoil.hpp:531
GetPropType< TypeTag, Properties::SolutionVector > SolutionVector
Definition: SimulatorFullyImplicitBlackoil.hpp:104
static void registerParameters()
Definition: SimulatorFullyImplicitBlackoil.hpp:176
GetPropType< TypeTag, Properties::Indices > BlackoilIndices
Definition: SimulatorFullyImplicitBlackoil.hpp:101
const Model & model() const
Definition: SimulatorFullyImplicitBlackoil.hpp:542
SimulatorFullyImplicitBlackoil(Simulator &simulator)
Definition: SimulatorFullyImplicitBlackoil.hpp:140
std::unique_ptr< time::StopWatch > solverTimer_
Definition: SimulatorFullyImplicitBlackoil.hpp:643
std::unique_ptr< ReservoirCouplingMaster > reservoirCouplingMaster_
Definition: SimulatorFullyImplicitBlackoil.hpp:651
SimulatorConvergenceOutput convergence_output_
Definition: SimulatorFullyImplicitBlackoil.hpp:647
std::unique_ptr< Solver > solver_
Definition: SimulatorFullyImplicitBlackoil.hpp:635
bool terminalOutput_
Definition: SimulatorFullyImplicitBlackoil.hpp:640
ModelParameters modelParam_
Definition: SimulatorFullyImplicitBlackoil.hpp:632
bool runStep(SimulatorTimer &timer)
Definition: SimulatorFullyImplicitBlackoil.hpp:316
void init(const SimulatorTimer &timer, int argc, char **argv)
Definition: SimulatorFullyImplicitBlackoil.hpp:248
GetPropType< TypeTag, Properties::AquiferModel > AquiferModel
Definition: SimulatorFullyImplicitBlackoil.hpp:106
GetPropType< TypeTag, Properties::Grid > Grid
Definition: SimulatorFullyImplicitBlackoil.hpp:98
SimulatorReport finalize()
Definition: SimulatorFullyImplicitBlackoil.hpp:512
SimulatorReport report_
Definition: SimulatorFullyImplicitBlackoil.hpp:642
WellModel & wellModel_()
Definition: SimulatorFullyImplicitBlackoil.hpp:623
const Schedule & schedule() const
Definition: SimulatorFullyImplicitBlackoil.hpp:614
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: SimulatorFullyImplicitBlackoil.hpp:100
std::array< std::string, 5 > getHeader() const override
Returns header data.
Definition: SimulatorFullyImplicitBlackoil.hpp:565
SimulatorReport run(SimulatorTimer &timer, int argc, char **argv)
Definition: SimulatorFullyImplicitBlackoil.hpp:191
bool checkRunningAsReservoirCouplingMaster()
Definition: SimulatorFullyImplicitBlackoil.hpp:218
Simulator & simulator_
Definition: SimulatorFullyImplicitBlackoil.hpp:630
Class handling simulator serialization.
Definition: SimulatorSerializer.hpp:55
void loadState()
Load state from file.
void loadTimerInfo(SimulatorTimer &timer)
Loads time step info from file.
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.
int loadStep() const
Returns step to load.
Definition: SimulatorSerializer.hpp:74
Definition: SimulatorTimer.hpp:39
double currentStepLength() const override
bool initialStep() const override
Whether the current step is the first step.
void report(std::ostream &os) const
double simulationTimeElapsed() const override
int currentStepNum() const override
bool done() const override
Return true if op++() has been called numSteps() times.
Definition: blackoilnewtonmethodparams.hpp:31
void printValues(std::ostream &os)
Print values of the run-time parameters.
auto Get(bool errorIfNotRegistered=true)
Retrieve a runtime parameter.
Definition: parametersystem.hpp:185
Definition: alignedallocator.hh:32
void registerSimulatorParameters()
void outputReportStep(const SimulatorTimer &timer)
Definition: blackoilboundaryratevector.hh:39
std::string compileTimestamp()
std::string moduleVersion()
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:233
PhaseUsage phaseUsageFromDeck(const EclipseState &eclipseState)
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Definition: SimulatorFullyImplicitBlackoil.hpp:72
static constexpr bool value
Definition: SimulatorFullyImplicitBlackoil.hpp:72
Definition: SimulatorFullyImplicitBlackoil.hpp:76
static constexpr auto * value
Definition: SimulatorFullyImplicitBlackoil.hpp:76
Definition: SimulatorFullyImplicitBlackoil.hpp:77
static constexpr int value
Definition: SimulatorFullyImplicitBlackoil.hpp:77
Definition: SimulatorFullyImplicitBlackoil.hpp:73
static constexpr auto * value
Definition: SimulatorFullyImplicitBlackoil.hpp:73
Definition: SimulatorFullyImplicitBlackoil.hpp:75
static constexpr auto * value
Definition: SimulatorFullyImplicitBlackoil.hpp:75
Definition: SimulatorFullyImplicitBlackoil.hpp:74
static constexpr auto * value
Definition: SimulatorFullyImplicitBlackoil.hpp:74
Definition: SimulatorFullyImplicitBlackoil.hpp:78
static constexpr bool value
Definition: SimulatorFullyImplicitBlackoil.hpp:78
Definition: BlackoilPhases.hpp:46
Abstract interface for simulator serialization ops.
Definition: SimulatorSerializer.hpp:36
Definition: SimulatorReport.hpp:122
SimulatorReportSingle success
Definition: SimulatorReport.hpp:123
double solver_time
Definition: SimulatorReport.hpp:38
bool converged
Definition: SimulatorReport.hpp:55
double total_time
Definition: SimulatorReport.hpp:37
int exit_status
Definition: SimulatorReport.hpp:58
double output_write_time
Definition: SimulatorReport.hpp:46