]> https://gitweb.dealii.org/ - dealii.git/commitdiff
The class DiscreteTime implemented. 8824/head
authorReza Rastak <rastak@stanford.edu>
Fri, 20 Sep 2019 22:55:20 +0000 (15:55 -0700)
committerReza Rastak <rastak@stanford.edu>
Sat, 4 Apr 2020 08:33:19 +0000 (01:33 -0700)
doc/news/changes/minor/20190920RezaRastak [new file with mode: 0644]
include/deal.II/base/discrete_time.h [new file with mode: 0644]
source/base/CMakeLists.txt
source/base/discrete_time.cc [new file with mode: 0644]
tests/base/discrete_time_1.cc [new file with mode: 0644]
tests/base/discrete_time_1.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20190920RezaRastak b/doc/news/changes/minor/20190920RezaRastak
new file mode 100644 (file)
index 0000000..0c825ce
--- /dev/null
@@ -0,0 +1,4 @@
+New: The class DiscreteTime is created to keep track of time
+during a time-dependent simulation.
+<br>
+(Reza Rastak, 2019/09/20)
diff --git a/include/deal.II/base/discrete_time.h b/include/deal.II/base/discrete_time.h
new file mode 100644 (file)
index 0000000..e1b9ca7
--- /dev/null
@@ -0,0 +1,196 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_discrete_time_h
+#define dealii_discrete_time_h
+
+#include <deal.II/base/config.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * Provides a means to keep track of the simulation time in a time-dependent
+ * simulation. It manages stepping forward from an initial time to a final
+ * time. It also allows adjusting the time step size during the simulation.
+ * It is guaranteed that at all times the current simulation time is in the
+ * closed interval between the start time and the end time.
+ *
+ * You can loop over all of the time steps by using a for loop
+ * @code
+ *   for (DiscreteTime time(0., 1., 0.3);
+ *        time.get_current_time() != time.get_end_time();
+ *        time.advance_time())
+ *   {
+ *     // Insert simulation code here
+ *   }
+ * @code
+ *
+ * In the above example the time starts at $t = 0$. Assuming the time step
+ * $dt = 0.3$ is not modified inside the loop, the time is advanced to
+ * $t = 0.3$, $t = 0.6$, $t = 0.9$ and finally it reaches the end time at
+ * $t = 1.0$. Note that the final step size is automatically reduced to
+ * $dt = 0.1$ in order to ensure that we finish the simulation exactly at
+ * the specified end time.
+ *
+ * @author Reza Rastak, 2019
+ */
+class DiscreteTime
+{
+public:
+  /**
+   * Constructor
+   */
+  DiscreteTime(const double start_time,
+               const double end_time,
+               const double start_step_size);
+
+  /**
+   * Return the current time.
+   */
+  double
+  get_current_time() const;
+
+  /**
+   * Return the start time.
+   */
+  double
+  get_start_time() const;
+
+  /**
+   * Return the end of the time interval.
+   * The final time step ends exactly at this point. This exact floating-point
+   * equality is very important because it allows us to use the expression
+   * <code>time.get_current_time() != time.get_end_time()</code> as the
+   * conditional statement in a for loop to check if the end time is reached.
+   */
+  double
+  get_end_time() const;
+
+  /**
+   * Return the size of the step from current time step to the next.
+   */
+  double
+  get_next_step_size() const;
+
+  /**
+   * Set the value of the next time step size. The next time advance_time()
+   * is called, the newly set @p time_step_size will be used to advance
+   * the simulation time. However, if the step is too large such that the next
+   * simulation time exceeds the end time, the step size is truncated.
+   * Additionally, if the step size is such that the next simulation time
+   * approximates the end time (but falls just slightly short of it), the step
+   * size is adjusted such that the next simulation time exactly matches the
+   * end time.
+   */
+  void
+  set_next_step_size(const double time_step_size);
+
+  /**
+   * Advance the current time based on the value of the current step.
+   * If you want to adjust the next time step size, call the method
+   * set_next_step_size() before calling this method.
+   * If you call this function repeatedly, the time
+   * is increased with the same step size until it reaches the end
+   * time. See the documentation of set_next_step_size() for explanation
+   * of the rules for automatic adjustment of the step size.
+   *
+   * @pre Current time must be smaller than the end time. The object cannot
+   * advance time if it is already at the end time. This rule is created to
+   * avoid the creation of an infinite loop when advance_time() is called
+   * inside a loop.
+   */
+  void
+  advance_time();
+
+  /**
+   * Set the current time equal to start time and set the step size to the
+   * initial step size.
+   */
+  void
+  restart();
+
+private:
+  /**
+   * The beginning of the time interval.
+   */
+  const double start_time;
+
+  /**
+   *The end of the time interval.
+   */
+  const double end_time;
+
+  /**
+   * The size of the first step.
+   */
+  const double start_step_size;
+
+  /**
+   * The current time.
+   */
+  double current_time;
+
+  /**
+   * The time at the next step.
+   *
+   * @note Internally, the next simulation time is stored instead of the
+   * current step size. For example, when the method set_next_step_size()
+   * is called, it computes the appropriate next simulation time and stores
+   * it. When advance_time() is called, the current_time is replaced by
+   * next_time. This choice for the internal state allows for simpler code
+   * and ensures than when we call advance_time() at the last step, the
+   * floating-point value of the time exactly matches the end time.
+   */
+  double next_time;
+};
+
+
+/*---------------------- Inline functions ------------------------------*/
+
+
+inline double
+DiscreteTime::get_start_time() const
+{
+  return start_time;
+}
+
+
+
+inline double
+DiscreteTime::get_end_time() const
+{
+  return end_time;
+}
+
+
+
+inline double
+DiscreteTime::get_next_step_size() const
+{
+  return next_time - current_time;
+}
+
+
+
+inline double
+DiscreteTime::get_current_time() const
+{
+  return current_time;
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 9b1b6e0f1c8a2b48ca4e79363f097d54095a3ec7..452fdc5f02c2156565c305f53e62d631a8dbe826 100644 (file)
@@ -25,6 +25,7 @@ SET(_unity_include_src
   bounding_box.cc
   conditional_ostream.cc
   convergence_table.cc
+  discrete_time.cc
   event.cc
   exceptions.cc
   flow_function.cc
diff --git a/source/base/discrete_time.cc b/source/base/discrete_time.cc
new file mode 100644 (file)
index 0000000..d45141f
--- /dev/null
@@ -0,0 +1,86 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/discrete_time.h>
+#include <deal.II/base/exceptions.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace
+{
+  // Helper function that computes the next discrete time, adjusting it if:
+  //  - The next time exceeds the end time.
+  //  - The next time is smaller but very close to the end time.
+  double
+  get_next_time(const double current_time,
+                const double step_size,
+                const double end_time)
+  {
+    Assert(step_size >= 0., ExcMessage("Time step size must be non-negative"));
+    Assert(end_time >= current_time, ExcInternalError());
+    double           next_time          = current_time + step_size;
+    constexpr double relative_tolerance = 0.05;
+    const double     time_tolerance     = relative_tolerance * step_size;
+    if (next_time > end_time - time_tolerance)
+      next_time = end_time;
+    return next_time;
+  }
+} // namespace
+
+
+
+DiscreteTime::DiscreteTime(const double start_time,
+                           const double end_time,
+                           const double start_step_size)
+  : start_time{start_time}
+  , end_time{end_time}
+  , start_step_size{start_step_size}
+  , current_time{start_time}
+  , next_time{get_next_time(start_time, start_step_size, end_time)}
+{}
+
+
+
+void
+DiscreteTime::set_next_step_size(const double next_step_size)
+{
+  next_time = get_next_time(current_time, next_step_size, end_time);
+}
+
+
+
+void
+DiscreteTime::advance_time()
+{
+  Assert(next_time > current_time,
+         ExcMessage(
+           "You can't advance time further."
+           "Either dt == 0 or you are at the end of the simulation time."));
+  const double step_size = get_next_step_size();
+  current_time           = next_time;
+  next_time              = get_next_time(current_time, step_size, end_time);
+}
+
+
+
+void
+DiscreteTime::restart()
+{
+  current_time = start_time;
+  next_time    = get_next_time(current_time, start_step_size, end_time);
+}
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/tests/base/discrete_time_1.cc b/tests/base/discrete_time_1.cc
new file mode 100644 (file)
index 0000000..363fbd4
--- /dev/null
@@ -0,0 +1,84 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Basic tests for class DiscreteTime
+
+#include <deal.II/base/discrete_time.h>
+
+#include "../tests.h"
+
+void
+print_time(const DiscreteTime &time)
+{
+  deallog << "Current time = " << time.get_current_time()
+          << ", next step size = " << time.get_next_step_size() << std::endl;
+}
+
+void
+test_from_start_to_end()
+{
+  LogStream::Prefix p("Start to end");
+
+  DiscreteTime time(/*start_time*/ 0.,
+                    /*end_time*/ 1.5,
+                    /*start_step_size*/ 0.123);
+
+  deallog << "Start time = " << time.get_start_time() << std::endl;
+  deallog << "End time = " << time.get_end_time() << std::endl;
+  print_time(time);
+  while (time.get_current_time() != time.get_end_time())
+    {
+      time.advance_time();
+      print_time(time);
+    }
+  time.restart();
+  deallog << "Restarted" << std::endl;
+  print_time(time);
+  deallog << "OK" << std::endl;
+}
+
+void
+test_adjust_time_step_size()
+{
+  LogStream::Prefix p("Adjust time step size");
+
+  DiscreteTime time(/*start_time*/ 0.4,
+                    /*end_time*/ 2.1,
+                    /*start_step_size*/ 0.15);
+  print_time(time);
+  time.advance_time();
+  print_time(time);
+  time.set_next_step_size(0.36);
+  time.advance_time();
+  print_time(time);
+  time.advance_time();
+  print_time(time);
+  time.set_next_step_size(0.61);
+  time.advance_time();
+  print_time(time);
+  time.advance_time(); // Here we reach the end time.
+  print_time(time);
+  // If we call time.advance_time() one more time, it fails an assersion
+  deallog << "OK" << std::endl;
+}
+
+int
+main()
+{
+  initlog();
+  test_from_start_to_end();
+  test_adjust_time_step_size();
+  return 0;
+}
diff --git a/tests/base/discrete_time_1.output b/tests/base/discrete_time_1.output
new file mode 100644 (file)
index 0000000..24806bc
--- /dev/null
@@ -0,0 +1,27 @@
+
+DEAL:Start to end::Start time = 0.00000
+DEAL:Start to end::End time = 1.50000
+DEAL:Start to end::Current time = 0.00000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.123000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.246000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.369000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.492000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.615000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.738000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.861000, next step size = 0.123000
+DEAL:Start to end::Current time = 0.984000, next step size = 0.123000
+DEAL:Start to end::Current time = 1.10700, next step size = 0.123000
+DEAL:Start to end::Current time = 1.23000, next step size = 0.123000
+DEAL:Start to end::Current time = 1.35300, next step size = 0.123000
+DEAL:Start to end::Current time = 1.47600, next step size = 0.0240000
+DEAL:Start to end::Current time = 1.50000, next step size = 0.00000
+DEAL:Start to end::Restarted
+DEAL:Start to end::Current time = 0.00000, next step size = 0.123000
+DEAL:Start to end::OK
+DEAL:Adjust time step size::Current time = 0.400000, next step size = 0.150000
+DEAL:Adjust time step size::Current time = 0.550000, next step size = 0.150000
+DEAL:Adjust time step size::Current time = 0.910000, next step size = 0.360000
+DEAL:Adjust time step size::Current time = 1.27000, next step size = 0.360000
+DEAL:Adjust time step size::Current time = 1.88000, next step size = 0.220000
+DEAL:Adjust time step size::Current time = 2.10000, next step size = 0.00000
+DEAL:Adjust time step size::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.