diff --git a/dumux/common/timeloop.hh b/dumux/common/timeloop.hh
index 3f2e02dad9c74e1c8d7d30cb3dbc1ef7e02e4c0c..3144c45db1a3f62268806cee3be32ddfc39e6971 100644
--- a/dumux/common/timeloop.hh
+++ b/dumux/common/timeloop.hh
@@ -31,6 +31,7 @@
#include <dune/common/timer.hh>
#include <dune/common/parallel/collectivecommunication.hh>
#include <dune/common/parallel/mpihelper.hh>
+#include <dune/common/exceptions.hh>
#include <dumux/common/parameters.hh>
@@ -88,18 +89,7 @@ public:
TimeLoop(Scalar startTime, Scalar dt, Scalar tEnd, bool verbose = true)
: timer_(false)
{
- verbose_ =
- verbose &&
- Dune::MPIHelper::getCollectiveCommunication().rank() == 0;
-
- time_ = startTime;
- endTime_ = tEnd;
-
- timeStepSize_ = dt;
- previousTimeStepSize_ = timeStepSize_;
- maxTimeStepSize_ = std::numeric_limits<Scalar>::max();
- timeStepIdx_ = 0;
- finished_ = false;
+ reset(startTime, dt, tEnd, verbose);
}
/*!
@@ -113,7 +103,48 @@ public:
void start()
{
timer_.start();
- cpuTime_ = 0.0;
+ }
+
+ /*!
+ * \brief Tells the time loop to stop tracking the time.
+ * \return the wall clock time (CPU time) spent until now
+ */
+ double stop()
+ {
+ return timer_.stop();
+ }
+
+ /*!
+ * \brief Reset the timer
+ */
+ void resetTimer()
+ {
+ timer_.reset();
+ }
+
+ /*!
+ * \brief Reset the time loop
+ */
+ void reset(Scalar startTime, Scalar dt, Scalar tEnd, bool verbose = true)
+ {
+ verbose_ =
+ verbose &&
+ Dune::MPIHelper::getCollectiveCommunication().rank() == 0;
+
+ time_ = startTime;
+ endTime_ = tEnd;
+
+ timeStepSize_ = dt;
+ lastTimeStepSize_ = 0.0;
+ maxTimeStepSize_ = std::numeric_limits<Scalar>::max();
+ userSetMaxTimeStepSize_ = maxTimeStepSize_;
+ timeStepIdx_ = 0;
+ finished_ = false;
+ timeAfterLastTimeStep_ = 0.0;
+ timeStepWallClockTime_ = 0.0;
+
+ timer_.stop();
+ timer_.reset();
}
/*!
@@ -123,6 +154,12 @@ public:
{
timeStepIdx_++;
time_ += timeStepSize_;
+ lastTimeStepSize_ = timeStepSize_;
+
+ // compute how long the last time step took
+ const auto cpuTime = wallClockTime();
+ timeStepWallClockTime_ = cpuTime - timeAfterLastTimeStep_;
+ timeAfterLastTimeStep_ = cpuTime;
}
/*!
@@ -163,13 +200,17 @@ public:
* \param t The time \f$\mathrm{[s]}\f$ at which the simulation is finished
*/
void setEndTime(Scalar t)
- { endTime_ = t; }
+ {
+ endTime_ = t;
+ if (verbose_)
+ std::cout << "Set new end time to t = " << t << " seconds." << std::endl;
+ }
/*!
- * \brief Returns the current wall time (cpu time).
+ * \brief Returns the current wall clock time (cpu time) spend in this time loop
*/
- double wallTime() const
- { return cpuTime_; }
+ double wallClockTime() const
+ { return timer_.elapsed(); }
/*!
* \brief Set the current time step size to a given value.
@@ -178,7 +219,7 @@ public:
* episode, the timeStep() method will take care that the step
* size won't exceed the episode or the end of the simulation,
* though.
- * \note Always call this after TimeLoop::advanceTimeStep()
+ * \note Always call this after TimeLoop::advanceTimeStep() and TimeLoop::reportTimeStep()
*
* \param dt The new value for the time step size \f$\mathrm{[s]}\f$
*/
@@ -193,9 +234,15 @@ public:
* \brief Set the maximum time step size to a given value.
*
* \param maxDt The new value for the maximum time step size \f$\mathrm{[s]}\f$
+ * \note This also updates the time step size
*/
void setMaxTimeStepSize(Scalar maxDt)
- { maxTimeStepSize_ = maxDt; }
+ {
+ using std::min;
+ userSetMaxTimeStepSize_ = maxDt;
+ computeMaxTimeStepSize_();
+ timeStepSize_ = min(timeStepSize_, maxTimeStepSize_);
+ }
/*!
* \brief Returns the suggested time step length \f$\mathrm{[s]}\f$ so that we
@@ -205,12 +252,6 @@ public:
Scalar timeStepSize() const override
{ return timeStepSize_; }
- /*!
- * \brief Returns the size of the previous time step \f$\mathrm{[s]}\f$.
- */
- Scalar previousTimeStepSize() const
- { return previousTimeStepSize_; }
-
/*!
* \brief Returns number of time steps which have been
* executed since the beginning of the simulation.
@@ -235,14 +276,20 @@ public:
* if the end time is reached.
*/
bool finished() const
- { return finished_ || time_ >= endTime_; }
+ {
+ using std::abs;
+ return finished_ || abs(time_-endTime_)/endTime_ < 1e-14;
+ }
/*!
* \brief Returns true if the simulation is finished after the
* time level is incremented by the current time step size.
*/
bool willBeFinished() const
- { return finished_ || time_ + timeStepSize_ >= endTime_; }
+ {
+ using std::abs;
+ return finished_ || abs(time_+timeStepSize_-endTime_)/endTime_ < 1e-14;
+ }
/*!
* \brief The current maximum time step size
@@ -250,27 +297,28 @@ public:
* and other possible check points
*/
Scalar maxTimeStepSize() const
- { return maxTimeStepSize_; }
+ {
+ using std::min;
+ return min(userSetMaxTimeStepSize_, maxTimeStepSize_);
+ }
/*!
* \brief State info on cpu time.
+ * \note Always call this after TimeLoop::advanceTimeStep()
*/
- void reportTimeStep()
+ void reportTimeStep() const
{
- auto timeStepCpuTime = timer_.elapsed();
- cpuTime_ += timeStepCpuTime;
+ const auto cpuTime = wallClockTime();
if (verbose_)
{
std::cout << "Time step " << timeStepIdx_ << " done in "
- << timeStepCpuTime << " seconds. "
- << "Wall time: " << cpuTime_
+ << timeStepWallClockTime_ << " seconds. "
+ << "Wall clock time: " << cpuTime
<< ", time: " << time_
- << ", time step size: " << timeStepSize_
+ << ", time step size: " << lastTimeStepSize_
<< std::endl;
}
-
- timer_.reset();
}
/*!
@@ -279,21 +327,20 @@ public:
template< class Communicator = Dune::CollectiveCommunication<typename Dune::MPIHelper::MPICommunicator> >
void finalize(const Communicator& comm = Dune::MPIHelper::getCollectiveCommunication())
{
- timer_.stop();
- cpuTime_ += timer_.elapsed();
+ auto cpuTime = timer_.stop();
if (verbose_)
{
- std::cout << "Simulation took " << cpuTime_ << " seconds on "
+ std::cout << "Simulation took " << cpuTime << " seconds on "
<< comm.size() << " processes.\n";
}
if (comm.size() > 1)
- cpuTime_ = comm.sum(cpuTime_);
+ cpuTime = comm.sum(cpuTime);
if (verbose_)
{
- std::cout << "The cumulative CPU time was " << cpuTime_ << " seconds.\n";
+ std::cout << "The cumulative CPU time was " << cpuTime << " seconds.\n";
}
}
@@ -307,7 +354,6 @@ public:
private:
//! Computes the maximum timestep size respecting end time
- //! and possibly episodes (TODO)
void computeMaxTimeStepSize_()
{
if (finished())
@@ -318,19 +364,17 @@ private:
using std::max;
using std::min;
-
- // TODO check for episodes if there is an episode manager
- maxTimeStepSize_ = min(maxTimeStepSize_, max<Scalar>(0.0, endTime_ - time_));
+ maxTimeStepSize_ = min(userSetMaxTimeStepSize_, max<Scalar>(0.0, endTime_ - time_));
}
Dune::Timer timer_;
Scalar time_;
Scalar endTime_;
- double cpuTime_;
Scalar timeStepSize_;
- Scalar previousTimeStepSize_;
- Scalar maxTimeStepSize_;
+ Scalar lastTimeStepSize_;
+ Scalar maxTimeStepSize_, userSetMaxTimeStepSize_;
+ Scalar timeAfterLastTimeStep_, timeStepWallClockTime_;
int timeStepIdx_;
bool finished_;
bool verbose_;
@@ -378,6 +422,9 @@ public:
{
isCheckPoint_ = false;
}
+
+ // make sure to respect future check check points
+ this->setTimeStepSize(this->timeStepSize());
}
/*!
@@ -400,21 +447,53 @@ public:
/*!
* \brief Set a periodic check point
* \note You can query if we are at a time check point with isCheckPoint()
- * \param interval Set a periodic checkout every [interal] seconds
+ * \param interval Set a periodic checkout every [interval] seconds
+ * \param offset time from which the periodic check points are supposed to start (simulation time)
+ * the first checkpoint will be at time = offset.
+ * \note If offset is in the past the first check point will be at the next
+ * periodic check point greater or equal than time
+ * \note This also updates the time step size and potentially reduces the time step size to meet the next check point
*/
- void setPeriodicCheckPoint(Scalar interval)
+ void setPeriodicCheckPoint(Scalar interval, Scalar offset = 0.0)
{
+ using std::signbit;
+ if (signbit(interval))
+ DUNE_THROW(Dune::InvalidStateException, "Interval has to be positive!");
+
periodicCheckPoints_ = true;
deltaPeriodicCheckPoint_ = interval;
+ lastPeriodicCheckPoint_ = offset;
+ while (lastPeriodicCheckPoint_ + interval - this->time() < 1e-14*interval)
+ lastPeriodicCheckPoint_ += interval;
+
if (this->verbose())
- std::cout << "Enabled periodic check points every " << interval << " seconds." << std::endl;
+ std::cout << "Enabled periodic check points every " << interval
+ << " seconds with the next check point at " << lastPeriodicCheckPoint_ + interval << " seconds." << std::endl;
// make sure we respect this check point on the next time step
setTimeStepSize(this->timeStepSize());
+
+ // check if the current time point is a check point
+ if (Dune::FloatCmp::eq(this->time(), lastPeriodicCheckPoint_, 1e-8*interval))
+ isCheckPoint_ = true;
}
- //! Whether now is a time checkpoint
- //! has to be called after TimeLoop::advanceTimeStep()
+ //! disable periodic check points
+ void disablePeriodicCheckPoints()
+ { periodicCheckPoints_ = false; }
+
+ //! remove all check points
+ void removeAllCheckPoints()
+ {
+ periodicCheckPoints_ = false;
+ while (!checkPoints_.empty())
+ checkPoints_.pop();
+ }
+
+ /*!
+ * \brief Whether now is a time checkpoint
+ * \note has to be called after TimeLoop::advanceTimeStep()
+ */
bool isCheckPoint() const
{ return isCheckPoint_; }
@@ -422,6 +501,7 @@ public:
* \brief add a checkpoint to the queue
* \note checkpoints have to be provided in ascending order
* \param t the check point (in seconds)
+ * \note This also updates the time step size and potentially reduces the time step size to meet the next check point
*/
void setCheckPoint(Scalar t)
{
@@ -436,6 +516,7 @@ public:
* \brief add checkpoints to the queue from a vector of time points
* \note checkpoints have to be provided in ascending order
* \param checkPoints the vector of check points
+ * \note This also updates the time step size and potentially reduces the time step size to meet the next check point
*/
void setCheckPoint(const std::vector<Scalar>& checkPoints)
{ setCheckPoint(checkPoints.begin(), checkPoints.end()); }
@@ -445,6 +526,7 @@ public:
* \note checkpoints have to be provided in ascending order
* \param first iterator to the first element to be inserted
* \param last iterator to the one-after-last element to be inserted
+ * \note This also updates the time step size and potentially reduces the time step size to meet the next check point
*/
template<class ForwardIterator>
void setCheckPoint(ForwardIterator first, ForwardIterator last)
@@ -461,19 +543,28 @@ private:
//! Adds a check point to the queue
void setCheckPoint_(Scalar t)
{
+ if (t < this->time())
+ {
+ if (this->verbose())
+ std::cerr << "Couldn't insert checkpoint at t = " << t
+ << " because that's in the past! (current simulation time is " << this->time() << ")" << std::endl;
+ return;
+ }
+
if (!checkPoints_.empty())
{
if (t < checkPoints_.back())
{
if (this->verbose())
- std::cerr << "--- Couldn't insert checkpoint as it is earlier than the last check point in the queue.\n"
- << "--- Checkpoints can only be inserted in ascending order." << std::endl;
+ std::cerr << "Couldn't insert checkpoint as it is earlier than the last check point in the queue.\n"
+ << "Checkpoints can only be inserted in ascending order." << std::endl;
+ return;
}
- else
- checkPoints_.push(t);
}
- else
- checkPoints_.push(t);
+
+ checkPoints_.push(t);
+ if (this->verbose())
+ std::cout << "Set check point at t = " << t << " seconds." << std::endl;
}
/*!
diff --git a/test/common/CMakeLists.txt b/test/common/CMakeLists.txt
index 6d8f4f2ea9f7328a92665453b2fe85f4155fcdf9..ad388d297dbd7d4e98ce607e51ea5af6ae389abd 100644
--- a/test/common/CMakeLists.txt
+++ b/test/common/CMakeLists.txt
@@ -4,4 +4,5 @@ add_subdirectory(math)
add_subdirectory(parameters)
add_subdirectory(propertysystem)
add_subdirectory(spline)
+add_subdirectory(timeloop)
add_subdirectory(typetraits)
diff --git a/test/common/timeloop/CMakeLists.txt b/test/common/timeloop/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..9bfd94b5b603acf63d7516f1939b87a09fd5c6cb
--- /dev/null
+++ b/test/common/timeloop/CMakeLists.txt
@@ -0,0 +1 @@
+dune_add_test(SOURCES test_timeloop.cc)
diff --git a/test/common/timeloop/test_timeloop.cc b/test/common/timeloop/test_timeloop.cc
new file mode 100644
index 0000000000000000000000000000000000000000..939fb14a89e0e8a6e81ece0a12febb96c9771c24
--- /dev/null
+++ b/test/common/timeloop/test_timeloop.cc
@@ -0,0 +1,153 @@
+// time loop test
+#include <config.h>
+
+#include <iostream>
+#include <cmath>
+#include <chrono>
+#include <thread>
+#include <limits>
+
+#include <dune/common/parallel/mpihelper.hh>
+#include <dune/common/exceptions.hh>
+
+#include <dumux/common/timeloop.hh>
+
+int main(int argc, char* argv[]) try
+{
+ const auto& mpiHelper = Dune::MPIHelper::instance(argc, argv);
+
+ //! Standard time loop
+ double tStart = 0; double tEnd = 1; double dt = 0.1;
+
+ //! Standard time loop
+ {
+ if (mpiHelper.rank() == 0) std::cout << "------- Test time loop ----------" << std::endl;
+ Dumux::TimeLoop<double> timeLoop(tStart, dt, tEnd);
+
+ if (std::abs(timeLoop.time()-tStart) > 1e-15)
+ DUNE_THROW(Dune::InvalidStateException, "Wrong time at start!");
+ if (std::abs(timeLoop.endTime()-tEnd) > 1e-15)
+ DUNE_THROW(Dune::InvalidStateException, "Wrong end time: " << timeLoop.endTime() << " should be " << tEnd);
+ if (std::abs(timeLoop.timeStepSize()-dt) > 1e-15)
+ DUNE_THROW(Dune::InvalidStateException, "Wrong time step size!");
+
+ timeLoop.setMaxTimeStepSize(0.03333333333333333);
+ if (std::abs(timeLoop.timeStepSize()-0.03333333333333333) > 1e-15)
+ DUNE_THROW(Dune::InvalidStateException, "Wrong time step size!");
+
+ timeLoop.start(); // starts the timer
+ timeLoop.stop(); // stops the timer
+ timeLoop.resetTimer(); // resets the timer
+ timeLoop.start(); // starts the timer again
+
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ if (timeLoop.timeStepIndex() != 1)
+ DUNE_THROW(Dune::InvalidStateException, "Wrong time step index!");
+
+ timeLoop.reportTimeStep();
+
+ while (!timeLoop.finished())
+ {
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ timeLoop.reportTimeStep();
+ timeLoop.setTimeStepSize(timeLoop.timeStepSize());
+ }
+
+ timeLoop.finalize();
+
+ if (std::abs(timeLoop.time()-tEnd) > 1e-15)
+ DUNE_THROW(Dune::InvalidStateException, "Ended with wrong end time!");
+ }
+
+
+ // check point timeLoop
+ {
+ if (mpiHelper.rank() == 0) std::cout << std::endl << "------- Test check point time loop ----------" << std::endl;
+ Dumux::CheckPointTimeLoop<double> timeLoop(tStart, dt, tEnd);
+ timeLoop.setPeriodicCheckPoint(0.03, -0.03);
+ if (!timeLoop.isCheckPoint())
+ DUNE_THROW(Dune::InvalidStateException, "This is supposed to be a check point! t = " << timeLoop.time());
+
+ timeLoop.start();
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ timeLoop.reportTimeStep();
+ timeLoop.setPeriodicCheckPoint(0.03);
+
+ if (std::signbit(timeLoop.timeStepSize()))
+ DUNE_THROW(Dune::InvalidStateException, "Time step size is negative! dt = " << timeLoop.timeStepSize());
+
+ if (!timeLoop.isCheckPoint())
+ DUNE_THROW(Dune::InvalidStateException, "This is supposed to be a check point! t = " << timeLoop.time());
+
+ timeLoop.setCheckPoint(0.02); // can't be set because it's in the past
+ timeLoop.setCheckPoint(0.08);
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ timeLoop.reportTimeStep();
+
+ if (std::signbit(timeLoop.timeStepSize()))
+ DUNE_THROW(Dune::InvalidStateException, "Time step size is negative! dt = " << timeLoop.timeStepSize());
+
+ if (!timeLoop.isCheckPoint())
+ DUNE_THROW(Dune::InvalidStateException, "This is supposed to be a check point! t = " << timeLoop.time());
+
+ timeLoop.setEndTime(1e9);
+ timeLoop.setTimeStepSize(1e7);
+
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ timeLoop.reportTimeStep();
+
+ if (std::signbit(timeLoop.timeStepSize()))
+ DUNE_THROW(Dune::InvalidStateException, "Time step size is negative! dt = " << timeLoop.timeStepSize());
+
+ if (!timeLoop.isCheckPoint())
+ DUNE_THROW(Dune::InvalidStateException, "This is supposed to be a check point! t = " << timeLoop.time());
+
+ timeLoop.removeAllCheckPoints();
+ timeLoop.setTimeStepSize(0.5e-6);
+ timeLoop.setPeriodicCheckPoint(1.0e6, 1.0e6);
+
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ timeLoop.reportTimeStep();
+
+ if (std::signbit(timeLoop.timeStepSize()))
+ DUNE_THROW(Dune::InvalidStateException, "Time step size is negative! dt = " << timeLoop.timeStepSize());
+
+ while (!timeLoop.finished())
+ {
+ std::this_thread::sleep_for(std::chrono::milliseconds(5));
+ timeLoop.advanceTimeStep();
+ timeLoop.reportTimeStep();
+ timeLoop.setTimeStepSize(std::numeric_limits<double>::max());
+
+ if (std::abs(timeLoop.time() - 2.0e6) < 1e-15 && !timeLoop.isCheckPoint())
+ DUNE_THROW(Dune::InvalidStateException, "This is supposed to be a check point! t = " << timeLoop.time());
+
+ if (timeLoop.timeStepIndex() == 10)
+ timeLoop.setPeriodicCheckPoint(1.0e8);
+
+ if (std::signbit(timeLoop.timeStepSize()))
+ DUNE_THROW(Dune::InvalidStateException, "Time step size is negative! dt = " << timeLoop.timeStepSize());
+ }
+
+ timeLoop.finalize();
+
+ if (std::abs(timeLoop.time()-1e9) > 1e-15)
+ DUNE_THROW(Dune::InvalidStateException, "Ended with wrong end time!");
+
+ }
+
+ return 0;
+}
+// //////////////////////////////////
+// Error handler
+// /////////////////////////////////
+catch (const Dune::Exception& e) {
+ std::cout << e << std::endl;
+ return 1;
+}