simulator.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_SIMULATOR_HH
29#define EWOMS_SIMULATOR_HH
30
31#if HAVE_MPI
32#define RESERVOIR_COUPLING_ENABLED
33#endif
34#ifdef RESERVOIR_COUPLING_ENABLED
37#endif
38
39#include <dune/common/parallel/mpihelper.hh>
40
42
44
51
53
54#include <algorithm>
55#include <cassert>
56#include <iostream>
57#include <limits>
58#include <memory>
59#include <string>
60#include <vector>
61
62namespace Opm {
63
64 // required as std::max is not constexpr for float128 / quads.
65 template <typename T>
66 static constexpr T constexpr_max(T a, T b) {
67 return (a > b) ? a : b;
68 }
69
82template <class TypeTag>
84{
90
91 using MPIComm = typename Dune::MPIHelper::MPICommunicator;
92 using Communication = Dune::Communication<MPIComm>;
93
94 // \Note: too small eps can not rule out confusion from the rounding errors, as we use 1.e-9 as a minimum.
95 static constexpr Scalar eps =
96 constexpr_max(std::numeric_limits<Scalar>::epsilon(), static_cast<Scalar>(1.0e-9));
97
98public:
99 // do not allow to copy simulators around
100 Simulator(const Simulator&) = delete;
101
102 explicit Simulator(bool verbose = true)
103 : Simulator(Communication(), verbose)
104 {
105 }
106
107 explicit Simulator(Communication comm, bool verbose = true)
108 {
109 TimerGuard setupTimerGuard(setupTimer_);
110
111 setupTimer_.start();
112
113 verbose_ = verbose && comm.rank() == 0;
114
115 timeStepIdx_ = 0;
116 startTime_ = 0.0;
117 time_ = 0.0;
118 endTime_ = Parameters::Get<Parameters::EndTime<Scalar>>();
119 timeStepSize_ = Parameters::Get<Parameters::InitialTimeStepSize<Scalar>>();
120 assert(timeStepSize_ > 0);
121 const std::string& predetTimeStepFile =
122 Parameters::Get<Parameters::PredeterminedTimeStepsFile>();
123 if (!predetTimeStepFile.empty()) {
124 forcedTimeSteps_ = readTimeStepFile<Scalar>(predetTimeStepFile);
125 }
126
127 episodeIdx_ = 0;
128 episodeStartTime_ = 0;
129 episodeLength_ = std::numeric_limits<Scalar>::max();
130
131 finished_ = false;
132
133 if (verbose_) {
134 std::cout << "Allocating the simulation vanguard\n" << std::flush;
135 }
136
137 {
139 vanguard_ = std::make_unique<Vanguard>(*this);
140 OPM_END_PARALLEL_TRY_CATCH("Allocating the simulation vanguard failed: ", comm);
141 }
142
143 if (verbose_) {
144 std::cout << "Distributing the vanguard's data\n" << std::flush;
145 }
146
147 {
149 vanguard_->loadBalance();
150 OPM_END_PARALLEL_TRY_CATCH("Could not distribute the vanguard data: ", comm);
151 }
152
153 // Only relevant for CpGrid and serial runs.
154 if (verbose_) {
155 std::cout << "Adding LGRs, if any, in serial run\n" << std::flush;
156 }
157
158 {
160 vanguard_->addLgrs();
161 OPM_END_PARALLEL_TRY_CATCH("Adding LGRs to the simulation vanguard in serial run failed: ", comm);
162 }
163
164 if (verbose_) {
165 std::cout << "Allocating the model\n" << std::flush;
166 }
167
168 {
170 model_ = std::make_unique<Model>(*this);
171 OPM_END_PARALLEL_TRY_CATCH("Could not allocate model: ", comm);
172 }
173
174 if (verbose_) {
175 std::cout << "Allocating the problem\n" << std::flush;
176 }
177
178 {
180 problem_ = std::make_unique<Problem>(*this);
181 OPM_END_PARALLEL_TRY_CATCH("Could not allocate the problem: ", comm);
182 }
183
184 if (verbose_) {
185 std::cout << "Initializing the model\n" << std::flush;
186 }
187
188 {
190 model_->finishInit();
191 OPM_END_PARALLEL_TRY_CATCH("Could not initialize the model: ", comm);
192 }
193
194 if (verbose_) {
195 std::cout << "Initializing the problem\n" << std::flush;
196 }
197
198 {
200 problem_->finishInit();
201 OPM_END_PARALLEL_TRY_CATCH("Could not initialize the problem: ", comm);
202 }
203
204 setupTimer_.stop();
205
206 if (verbose_) {
207 std::cout << "Simulator successfully set up\n" << std::flush;
208 }
209 }
210
214 static void registerParameters()
215 {
216 Parameters::Register<Parameters::EndTime<Scalar>>
217 ("The simulation time at which the simulation is finished [s]");
218 Parameters::Register<Parameters::InitialTimeStepSize<Scalar>>
219 ("The size of the initial time step [s]");
220 Parameters::Register<Parameters::RestartTime<Scalar>>
221 ("The simulation time at which a restart should be attempted [s]");
222 Parameters::Register<Parameters::PredeterminedTimeStepsFile>
223 ("A file with a list of predetermined time step sizes "
224 "(one time step per line)");
225
226 Vanguard::registerParameters();
227 Model::registerParameters();
228 Problem::registerParameters();
229 }
230
234 Vanguard& vanguard()
235 { return *vanguard_; }
236
240 const Vanguard& vanguard() const
241 { return *vanguard_; }
242
246 const GridView& gridView() const
247 { return vanguard_->gridView(); }
248
252 Model& model()
253 { return *model_; }
254
258 const Model& model() const
259 { return *model_; }
260
265 Problem& problem()
266 { return *problem_; }
267
272 const Problem& problem() const
273 { return *problem_; }
274
280 void setStartTime(Scalar t)
281 { startTime_ = t; }
282
286 Scalar startTime() const
287 { return startTime_; }
288
295 void setTime(Scalar t)
296 { time_ = t; }
297
304 void setTime(Scalar t, unsigned stepIdx)
305 {
306 time_ = t;
307 timeStepIdx_ = stepIdx;
308 }
309
317 Scalar time() const
318 { return time_; }
319
325 void setEndTime(Scalar t)
326 { endTime_ = t; }
327
332 Scalar endTime() const
333 { return endTime_; }
334
339 const Timer& setupTimer() const
340 { return setupTimer_; }
341
346 const Timer& executionTimer() const
347 { return executionTimer_; }
348
350 { return executionTimer_; }
351
357 { return prePostProcessTimer_; }
358
363 const Timer& linearizeTimer() const
364 { return linearizeTimer_; }
365
370 const Timer& solveTimer() const
371 { return solveTimer_; }
372
377 const Timer& updateTimer() const
378 { return updateTimer_; }
379
384 const Timer& writeTimer() const
385 { return writeTimer_; }
386
397 void setTimeStepSize(Scalar value)
398 { timeStepSize_ = value; }
399
405 void setTimeStepIndex(unsigned value)
406 { timeStepIdx_ = value; }
407
413 Scalar timeStepSize() const
414 { return timeStepSize_; }
415
420 int timeStepIndex() const
421 { return timeStepIdx_; }
422
430 void setFinished(bool yesno = true)
431 { finished_ = yesno; }
432
439 bool finished() const
440 {
441 assert(timeStepSize_ >= 0.0);
442 return finished_ || (this->time() * (1.0 + eps) >= endTime());
443 }
444
449 bool willBeFinished() const
450 {
451 return finished_ || (this->time() + timeStepSize_) * (1.0 + eps) >= endTime();
452 }
453
458 Scalar maxTimeStepSize() const
459 {
460 if (finished()) {
461 return 0.0;
462 }
463
464 return std::min(episodeMaxTimeStepSize(),
465 std::max<Scalar>(0.0, endTime() - this->time()));
466 }
467
475 {
476 ++episodeIdx_;
477 episodeStartTime_ = episodeStartTime;
478 episodeLength_ = episodeLength;
479 }
480
488 void startNextEpisode(Scalar len = std::numeric_limits<Scalar>::max())
489 {
490 ++episodeIdx_;
491 episodeStartTime_ = startTime_ + time_;
492 episodeLength_ = len;
493 }
494
500 void setEpisodeIndex(int episodeIdx)
501 { episodeIdx_ = episodeIdx; }
502
508 int episodeIndex() const
509 { return episodeIdx_; }
510
515 Scalar episodeStartTime() const
516 { return episodeStartTime_; }
517
523 void setEpisodeLength(Scalar dt)
524 { episodeLength_ = dt; }
525
530 Scalar episodeLength() const
531 { return episodeLength_; }
532
537 bool episodeStarts() const
538 {
539 return this->time() <= (episodeStartTime_ - startTime()) * (1 + eps);
540 }
541
546 bool episodeIsOver() const
547 {
548 return this->time() >= (episodeStartTime_ - startTime() + episodeLength()) * (1 - eps);
549 }
550
555 bool episodeWillBeOver() const
556 {
557 return this->time() + timeStepSize()
558 >= (episodeStartTime_ - startTime() + episodeLength()) * (1 - eps);
559 }
560
566 {
567 // if the current episode is over and the simulation
568 // wants to give it some extra time, we will return
569 // the time step size it suggested instead of trying
570 // to align it to the end of the episode.
571 if (episodeIsOver()) {
572 return 0.0;
573 }
574
575 // make sure that we don't exceed the end of the
576 // current episode.
577 return std::max<Scalar>(0.0,
579 (this->time() + this->startTime()));
580 }
581
582 /*
583 * \}
584 */
585
592 void run()
593 {
594 // create TimerGuard objects to hedge for exceptions
595 TimerGuard setupTimerGuard(setupTimer_);
596 TimerGuard executionTimerGuard(executionTimer_);
597 TimerGuard prePostProcessTimerGuard(prePostProcessTimer_);
598 TimerGuard writeTimerGuard(writeTimer_);
599
600 setupTimer_.start();
601 const Scalar restartTime = Parameters::Get<Parameters::RestartTime<Scalar>>();
602 if (restartTime > -1e30) {
603 // try to restart a previous simulation
604 time_ = restartTime;
605
607 Restart res;
608 res.deserializeBegin(*this, time_);
609
610 if (verbose_) {
611 std::cout << "Deserialize from file '" << res.fileName() << "'\n" << std::flush;
612 }
613
614 this->deserialize(res);
615 problem_->deserialize(res);
616 model_->deserialize(res);
617 res.deserializeEnd();
618 OPM_END_PARALLEL_TRY_CATCH("Deserialization failed: ",
619 Dune::MPIHelper::getCommunication());
620 if (verbose_) {
621 std::cout << "Deserialization done."
622 << " Simulator time: " << time() << humanReadableTime(time())
623 << " Time step index: " << timeStepIndex()
624 << " Episode index: " << episodeIndex()
625 << "\n" << std::flush;
626 }
627 }
628 else {
629 // if no restart is done, apply the initial solution
630 if (verbose_) {
631 std::cout << "Applying the initial solution of the \"" << problem_->name()
632 << "\" problem\n" << std::flush;
633 }
634
635 const Scalar oldTimeStepSize = timeStepSize_;
636 const int oldTimeStepIdx = timeStepIdx_;
637 timeStepSize_ = 0.0;
638 timeStepIdx_ = -1;
639
640 {
642 model_->applyInitialSolution();
643 OPM_END_PARALLEL_TRY_CATCH("Apply initial solution failed: ",
644 Dune::MPIHelper::getCommunication());
645 }
646
647 // write initial condition
648 if (problem_->shouldWriteOutput()) {
650 problem_->writeOutput(true);
651 OPM_END_PARALLEL_TRY_CATCH("Write output failed: ",
652 Dune::MPIHelper::getCommunication());
653 }
654
655 timeStepSize_ = oldTimeStepSize;
656 timeStepIdx_ = oldTimeStepIdx;
657 }
658 setupTimer_.stop();
659
660 executionTimer_.start();
661 bool episodeBegins = episodeIsOver() || (timeStepIdx_ == 0);
662 // do the time steps
663 while (!finished()) {
664 prePostProcessTimer_.start();
665 if (episodeBegins) {
666 // notify the problem that a new episode has just been
667 // started.
668 {
670 problem_->beginEpisode();
671 OPM_END_PARALLEL_TRY_CATCH("Begin episode failed: ",
672 Dune::MPIHelper::getCommunication());
673 }
674
675 if (finished()) {
676 // the problem can chose to terminate the simulation in
677 // beginEpisode(), so we have handle this case.
679 problem_->endEpisode();
680 OPM_END_PARALLEL_TRY_CATCH("End episode failed: ",
681 Dune::MPIHelper::getCommunication());
682 prePostProcessTimer_.stop();
683
684 break;
685 }
686 }
687 episodeBegins = false;
688
689 if (verbose_) {
690 std::cout << "Begin time step " << timeStepIndex() + 1 << ". "
691 << "Start time: " << this->time() << " seconds" << humanReadableTime(this->time())
692 << ", step size: " << timeStepSize() << " seconds" << humanReadableTime(timeStepSize())
693 << "\n";
694 }
695
696 // pre-process the current solution
697 {
699 problem_->beginTimeStep();
700 OPM_END_PARALLEL_TRY_CATCH("Begin timestep failed: ",
701 Dune::MPIHelper::getCommunication());
702 }
703
704 if (finished()) {
705 // the problem can choose to terminate the simulation in
706 // beginTimeStep(), so we have handle this case.
708 problem_->endTimeStep();
709 problem_->endEpisode();
710 OPM_END_PARALLEL_TRY_CATCH("Finish failed: ",
711 Dune::MPIHelper::getCommunication());
712 prePostProcessTimer_.stop();
713
714 break;
715 }
716 prePostProcessTimer_.stop();
717
718 try {
719 // execute the time integration scheme
720 problem_->timeIntegration();
721 }
722 catch (...) {
723 // exceptions in the time integration might be recoverable. clean up in
724 // case they are
725 const auto& pmodel = problem_->model();
726 prePostProcessTimer_ += pmodel.prePostProcessTimer();
727 linearizeTimer_ += pmodel.linearizeTimer();
728 solveTimer_ += pmodel.solveTimer();
729 updateTimer_ += pmodel.updateTimer();
730
731 throw;
732 }
733
734 const auto& pmodel = problem_->model();
735 prePostProcessTimer_ += pmodel.prePostProcessTimer();
736 linearizeTimer_ += pmodel.linearizeTimer();
737 solveTimer_ += pmodel.solveTimer();
738 updateTimer_ += pmodel.updateTimer();
739
740 // post-process the current solution
741 prePostProcessTimer_.start();
742 {
744 problem_->endTimeStep();
745 OPM_END_PARALLEL_TRY_CATCH("End timestep failed: ",
746 Dune::MPIHelper::getCommunication());
747 }
748 prePostProcessTimer_.stop();
749
750 // write the result to disk
751 writeTimer_.start();
752 if (problem_->shouldWriteOutput()) {
754 problem_->writeOutput(true);
755 OPM_END_PARALLEL_TRY_CATCH("Write output failed: ",
756 Dune::MPIHelper::getCommunication());
757 }
758 writeTimer_.stop();
759
760 // do the next time integration
761 const Scalar oldDt = timeStepSize();
762 {
764 problem_->advanceTimeLevel();
765 OPM_END_PARALLEL_TRY_CATCH("Advance time level failed: ",
766 Dune::MPIHelper::getCommunication());
767 }
768
769 if (verbose_) {
770 std::cout << "Time step " << timeStepIndex() + 1 << " done. "
771 << "CPU time: " << executionTimer_.realTimeElapsed()
772 << " seconds" << humanReadableTime(executionTimer_.realTimeElapsed())
773 << ", end time: " << this->time() + oldDt << " seconds"
774 << humanReadableTime(this->time() + oldDt)
775 << ", step size: " << oldDt << " seconds" << humanReadableTime(oldDt)
776 << "\n" << std::flush;
777 }
778
779 // advance the simulated time by the current time step size
780 time_ += oldDt;
781 ++timeStepIdx_;
782
783 prePostProcessTimer_.start();
784 // notify the problem if an episode is finished
785 if (episodeIsOver()) {
786 // Notify the problem about the end of the current episode...
788 problem_->endEpisode();
789 OPM_END_PARALLEL_TRY_CATCH("End episode failed: ",
790 Dune::MPIHelper::getCommunication());
791 episodeBegins = true;
792 }
793 else {
794 Scalar dt;
795 if (timeStepIdx_ < static_cast<int>(forcedTimeSteps_.size())) {
796 // use the next time step size from the input file
797 dt = forcedTimeSteps_[timeStepIdx_];
798 }
799 else {
800 // ask the problem to provide the next time step size
801 dt = std::min(maxTimeStepSize(), problem_->nextTimeStepSize());
802 }
803 assert(finished() || dt > 0);
804 setTimeStepSize(dt);
805 }
806 prePostProcessTimer_.stop();
807
808 // write restart file if mandated by the problem
809 writeTimer_.start();
810 if (problem_->shouldWriteRestartFile()) {
812 serialize();
813 OPM_END_PARALLEL_TRY_CATCH("Serialize failed: ",
814 Dune::MPIHelper::getCommunication());
815 }
816 writeTimer_.stop();
817 }
818 executionTimer_.stop();
819
820 {
822 problem_->finalize();
823 OPM_END_PARALLEL_TRY_CATCH("Finalize failed: ",
824 Dune::MPIHelper::getCommunication());
825 }
826 }
827
828#ifdef RESERVOIR_COUPLING_ENABLED
830 {
831 return reservoirCouplingMaster_;
832 }
834 {
835 return reservoirCouplingSlave_;
836 }
838 {
839 this->reservoirCouplingMaster_ = reservoirCouplingMaster;
840 }
842 {
843 this->reservoirCouplingSlave_ = reservoirCouplingSlave;
844 }
845#endif
846
862 {
863 using Restarter = Restart;
864 Restarter res;
865 res.serializeBegin(*this);
866 if (gridView().comm().rank() == 0) {
867 std::cout << "Serialize to file '" << res.fileName() << "'"
868 << ", next time step size: " << timeStepSize()
869 << "\n" << std::flush;
870 }
871
872 this->serialize(res);
873 problem_->serialize(res);
874 model_->serialize(res);
875 res.serializeEnd();
876 }
877
885 template <class Restarter>
886 void serialize(Restarter& restarter)
887 {
888 restarter.serializeSectionBegin("Simulator");
889 restarter.serializeStream()
890 << episodeIdx_ << " "
891 << episodeStartTime_ << " "
892 << episodeLength_ << " "
893 << startTime_ << " "
894 << time_ << " "
895 << timeStepIdx_ << " ";
896 restarter.serializeSectionEnd();
897 }
898
906 template <class Restarter>
907 void deserialize(Restarter& restarter)
908 {
909 restarter.deserializeSectionBegin("Simulator");
910 restarter.deserializeStream()
911 >> episodeIdx_
912 >> episodeStartTime_
913 >> episodeLength_
914 >> startTime_
915 >> time_
916 >> timeStepIdx_;
917 restarter.deserializeSectionEnd();
918 }
919
920 template<class Serializer>
921 void serializeOp(Serializer& serializer)
922 {
923 serializer(*vanguard_);
924 serializer(*model_);
925 serializer(*problem_);
926 serializer(episodeIdx_);
927 serializer(episodeStartTime_);
928 serializer(episodeLength_);
929 serializer(startTime_);
930 serializer(time_);
931 serializer(timeStepIdx_);
932 }
933
934private:
935 std::unique_ptr<Vanguard> vanguard_;
936 std::unique_ptr<Model> model_;
937 std::unique_ptr<Problem> problem_;
938
939 int episodeIdx_;
940 Scalar episodeStartTime_;
941 Scalar episodeLength_;
942
943 Timer setupTimer_;
944 Timer executionTimer_;
945 Timer prePostProcessTimer_;
946 Timer linearizeTimer_;
947 Timer solveTimer_;
948 Timer updateTimer_;
949 Timer writeTimer_;
950
951 std::vector<Scalar> forcedTimeSteps_;
952 Scalar startTime_;
953 Scalar time_;
954 Scalar endTime_;
955
956 Scalar timeStepSize_;
957 int timeStepIdx_;
958
959 bool finished_;
960 bool verbose_;
961
962#ifdef RESERVOIR_COUPLING_ENABLED
963 ReservoirCouplingMaster *reservoirCouplingMaster_ = nullptr;
964 ReservoirCouplingSlave *reservoirCouplingSlave_ = nullptr;
965#endif
966
967};
968
969namespace Properties {
970template<class TypeTag>
971struct Simulator<TypeTag, TTag::NumericModel>
973}
974
975} // namespace Opm
976
977#endif
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:192
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:158
Defines a type tags and some fundamental properties all models.
Definition: ReservoirCouplingMaster.hpp:35
Definition: ReservoirCouplingSlave.hpp:35
Load or save a state of a problem to/from the harddisk.
Definition: restart.hpp:45
void serializeBegin(Simulator &simulator)
Write the current state of the model to disk.
Definition: restart.hpp:92
const std::string & fileName() const
Returns the name of the file which is (de-)serialized.
Definition: restart.hpp:85
void deserializeBegin(Simulator &simulator, Scalar t)
Start reading a restart file at a certain simulated time.
Definition: restart.hpp:147
void deserializeEnd()
Stop reading the restart file.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:84
const Timer & writeTimer() const
Returns a reference to the timer object which measures the time needed to write the visualization out...
Definition: simulator.hh:384
Scalar timeStepSize() const
Returns the time step length so that we don't miss the beginning of the next episode or cross the en...
Definition: simulator.hh:413
Scalar startTime() const
Return the time of the start of the simulation.
Definition: simulator.hh:286
int timeStepIndex() const
Returns number of time steps which have been executed since the beginning of the simulation.
Definition: simulator.hh:420
void serialize()
This method writes the complete state of the simulation to the harddisk.
Definition: simulator.hh:861
Timer & executionTimer()
Definition: simulator.hh:349
Scalar episodeLength() const
Returns the length of the current episode in simulated time .
Definition: simulator.hh:530
const Timer & prePostProcessTimer() const
Returns a reference to the timer object which measures the time needed for pre- and postprocessing of...
Definition: simulator.hh:356
const Timer & solveTimer() const
Returns a reference to the timer object which measures the time needed by the solver.
Definition: simulator.hh:370
void startNextEpisode(Scalar len=std::numeric_limits< Scalar >::max())
Start the next episode, but don't change the episode identifier.
Definition: simulator.hh:488
void serializeOp(Serializer &serializer)
Definition: simulator.hh:921
void setReservoirCouplingSlave(ReservoirCouplingSlave *reservoirCouplingSlave)
Definition: simulator.hh:841
const Vanguard & vanguard() const
Return a reference to the grid manager of simulation.
Definition: simulator.hh:240
void serialize(Restarter &restarter)
Write the time manager's state to a restart file.
Definition: simulator.hh:886
const Timer & updateTimer() const
Returns a reference to the timer object which measures the time needed to the solutions of the non-li...
Definition: simulator.hh:377
const Timer & executionTimer() const
Returns a reference to the timer object which measures the time needed to run the simulation.
Definition: simulator.hh:346
void startNextEpisode(Scalar episodeStartTime, Scalar episodeLength)
Change the current episode of the simulation.
Definition: simulator.hh:474
void setEndTime(Scalar t)
Set the time of simulated seconds at which the simulation runs.
Definition: simulator.hh:325
const Timer & linearizeTimer() const
Returns a reference to the timer object which measures the time needed for linarizing the solutions.
Definition: simulator.hh:363
void setStartTime(Scalar t)
Set the time of the start of the simulation.
Definition: simulator.hh:280
void deserialize(Restarter &restarter)
Read the time manager's state from a restart file.
Definition: simulator.hh:907
Simulator(Communication comm, bool verbose=true)
Definition: simulator.hh:107
ReservoirCouplingMaster * reservoirCouplingMaster() const
Definition: simulator.hh:829
void setTimeStepSize(Scalar value)
Set the current time step size to a given value.
Definition: simulator.hh:397
bool episodeStarts() const
Returns true if the current episode has just been started at the current time.
Definition: simulator.hh:537
void run()
Runs the simulation using a given problem class.
Definition: simulator.hh:592
Vanguard & vanguard()
Return a reference to the grid manager of simulation.
Definition: simulator.hh:234
int episodeIndex() const
Returns the index of the current episode.
Definition: simulator.hh:508
void setTimeStepIndex(unsigned value)
Set the current time step index to a given value.
Definition: simulator.hh:405
ReservoirCouplingSlave * reservoirCouplingSlave() const
Definition: simulator.hh:833
void setFinished(bool yesno=true)
Specify whether the simulation is finished.
Definition: simulator.hh:430
Problem & problem()
Return the object which specifies the pysical setup of the simulation.
Definition: simulator.hh:265
static void registerParameters()
Registers all runtime parameters used by the simulation.
Definition: simulator.hh:214
void setTime(Scalar t)
Set the current simulated time, don't change the current time step index.
Definition: simulator.hh:295
bool willBeFinished() const
Returns true if the simulation is finished after the time level is incremented by the current time st...
Definition: simulator.hh:449
bool finished() const
Returns true if the simulation is finished.
Definition: simulator.hh:439
Scalar maxTimeStepSize() const
Aligns the time step size to the episode boundary and to the end time of the simulation.
Definition: simulator.hh:458
const Model & model() const
Return the physical model used in the simulation.
Definition: simulator.hh:258
Simulator(bool verbose=true)
Definition: simulator.hh:102
void setTime(Scalar t, unsigned stepIdx)
Set the current simulated time and the time step index.
Definition: simulator.hh:304
bool episodeIsOver() const
Returns true if the current episode is finished at the current time.
Definition: simulator.hh:546
Scalar endTime() const
Returns the number of (simulated) seconds which the simulation runs.
Definition: simulator.hh:332
Simulator(const Simulator &)=delete
bool episodeWillBeOver() const
Returns true if the current episode will be finished after the current time step.
Definition: simulator.hh:555
const GridView & gridView() const
Return the grid view for which the simulation is done.
Definition: simulator.hh:246
void setEpisodeIndex(int episodeIdx)
Sets the index of the current episode.
Definition: simulator.hh:500
Scalar time() const
Return the number of seconds of simulated time which have elapsed since the start time.
Definition: simulator.hh:317
void setReservoirCouplingMaster(ReservoirCouplingMaster *reservoirCouplingMaster)
Definition: simulator.hh:837
void setEpisodeLength(Scalar dt)
Sets the length in seconds of the current episode.
Definition: simulator.hh:523
const Problem & problem() const
Return the object which specifies the pysical setup of the simulation.
Definition: simulator.hh:272
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:252
Scalar episodeMaxTimeStepSize() const
Aligns the time step size to the episode boundary if the current time step crosses the boundary of th...
Definition: simulator.hh:565
Scalar episodeStartTime() const
Returns the absolute time when the current episode started .
Definition: simulator.hh:515
const Timer & setupTimer() const
Returns a reference to the timer object which measures the time needed to set up and initialize the s...
Definition: simulator.hh:339
A simple class which makes sure that a timer gets stopped if an exception is thrown.
Definition: timerguard.hh:42
Provides an encapsulation to measure the system time.
Definition: timer.hpp:46
void start()
Start counting the time resources used by the simulation.
double realTimeElapsed() const
Return the real time [s] elapsed during the periods the timer was active since the last reset.
double stop()
Stop counting the time resources.
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilboundaryratevector.hh:39
static constexpr T constexpr_max(T a, T b)
Definition: simulator.hh:66
std::string humanReadableTime(double timeInSeconds, bool isAmendment=true)
Given a time step size in seconds, return it in a format which is more easily parsable by humans.
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
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Manages the simulation time.
Definition: basicproperties.hh:116