]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New const member functions added to DiscreteTime
authorReza Rastak <rastak@stanford.edu>
Thu, 26 Mar 2020 05:52:43 +0000 (22:52 -0700)
committerReza Rastak <rastak@stanford.edu>
Thu, 23 Apr 2020 01:35:01 +0000 (18:35 -0700)
include/deal.II/base/discrete_time.h
source/base/discrete_time.cc
tests/base/discrete_time_1.cc
tests/base/discrete_time_1.output

index 3427ed7f0335d56501e8cf76c19185b2c77952e3..de7fa92bc2ef708e3aa9a8bac07913a85dc8f38c 100644 (file)
@@ -24,9 +24,16 @@ DEAL_II_NAMESPACE_OPEN
  * This class provides a means to keep track of the simulation time in a
  * time-dependent simulation. It manages stepping forward from a start time
  * $T_{\text{start}}$ to an end time $T_{\text{end}}$. 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.
+ * the time step size during the simulation.
+ *
+ * This class provides a number of invariants that are guaranteed to be
+ * true at all times.
+ *
+ * * The current simulation time is within the closed interval between the
+ *   start time and the end time ($T_{\text{start}} \le t \le T_{\text{end}}$).
+ * * Whenever time is incremented, the step size is positive ($dt > 0$).
+ *   In other words, time advances in strictly ascending order
+ *   ($m < n \Leftrightarrow t_m < t_n$).
  *
  * The model this class follows is that one sets a *desired* time step length
  * either through the constructor or using set_desired_next_step_size()
@@ -38,7 +45,7 @@ DEAL_II_NAMESPACE_OPEN
  * Let's say that you 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.is_at_end() == false;
  *        time.advance_time())
  *   {
  *     // Insert simulation code here
@@ -63,7 +70,7 @@ DEAL_II_NAMESPACE_OPEN
  * end time:
  * @code
  *   for (DiscreteTime time(0., 1.21, 0.3);
- *        time.get_current_time() != time.get_end_time();
+ *        time.is_at_end() == false;
  *        time.advance_time())
  *   {
  *     // Insert simulation code here
@@ -78,13 +85,143 @@ DEAL_II_NAMESPACE_OPEN
  * only a *desired* step size. You can query the actual time step size using the
  * get_next_step_size() function.
  *
+ *
+ * ### Details of time-stepping
+ *
+ * Since time is marched forward in a discrete manner in our simulations, we
+ * need to discuss how we increment time. During time stepping we enter two
+ * separate alternating regimes in every step.
+ *
+ * * The **snapshot** stage (the **current** stage, the **consistent**
+ *   stage): In this part of the algorithm, we are at $t = t_n$ and all
+ *   quantities of the simulation (displacements, strains, temperatures, etc.)
+ *   are up-to-date for $t = t_n$. In this stage, *current time* refers to
+ *   $t_n$, *next time* refers to $t_{n+1}$, *previous time* refers to
+ *   $t_{n-1}$. The other useful notation quantities are the *next* time step
+ *   size $t_{n+1} - t_n$ and *previous* time step size $t_n - t_{n-1}$. In
+ *   this stage, it is a perfect occasion to generate text output using print
+ *   commands within the user's code. Additionally, post-processed outputs can
+ *   be prepared here which can be viewed later by visualization programs such
+ *   as `Tecplot`, `Paraview`, and `VisIt`. Additionally, during the snapshot
+ *   stage, the code can assess the quality of the previous step and decide
+ *   whether it wants to increase or decrease the time step size. The step
+ *   size for the next time step can be modified here.
+ * * The **update** stage (the **transition** stage, the **inconsistent**
+ *   stage): In this section of the program, the internal state of the
+ *   simulation is getting updated from $t_n$ to $t_{n+1}$. All of the
+ *   variables need to be updated one by one, the step number is incremented,
+ *   the time is incremented by $dt = t_{n+1} - t_n$, and time-integration
+ *   algorithms are used to update the other simulation quantities. In the
+ *   middle of this stage, some variables have been updated to $t_{n+1}$ but
+ *   other variables still represent their value at $t_n$. Thus, we call this
+ *   the inconsistent stage, requiring that no post-processing output related
+ *   to the state variables take place within it. The state variables, namely
+ *   those related to time, the solution field and any internal variables, are
+ *   not synchronized and then get updated one by one. In general, the order of
+ *   updating variables is arbitrary, but some care should be taken if there
+ *   are interdependencies between them. For example, if some variable such as
+ *   $x$ depends on the calculation of another variable such as $y$, then $y$
+ *   must be updated before $x$ can be updated.
+ *
+ *   The question arises whether time should be incremented before updating
+ *   state quantities. Multiple possibilities exist, depending on program and
+ *   formulation requirements, and possibly the programmer's preferences:
+ *   * Time is incremented before the rest of the updates. In this case, even
+ *     though time is incremented to $t_{n+1}$, not all variables are updated
+ *     yet. During this update phase, $dt$ equals the *previous* time step
+ *     size. *Previous* means that it is referring to the $dt$ of the
+ *     `advance_time()` command that was performed previously. In the
+ *     following example code, we are assuming that `a` and `b` are two state
+ *     variables that need to be updated in this time step.
+ *     @code
+ *       time.advance_time();
+ *       new_a = update_a(a, b, time.get_previous_step_size());
+ *       b = update_b(a, b, time.get_previous_step_size());
+ *       a = new_a;
+ *     @endcode
+ *   * Time is incremented from $t_n$ to $t_{n+1}$ after all variables have
+ *     already been updated for $t_{n+1}$. During the update stage, $dt$ is
+ *     denoted as the *next* time step size. *Next* means that $dt$ of the
+ *     step corresponds to the `advance_time()` command that will happen
+ *     subsequently.
+ *     @code
+ *       new_a = update_a(a, b, time.get_next_step_size());
+ *       b = update_b(a, b, time.get_next_step_size());
+ *       a = new_a;
+ *       time.advance_time();
+ *     @endcode
+ *   * Time is incremented in the middle of the other updates: In this case
+ *     $dt$ would correspond to *next* or *previous* depending of whether it
+ *     is used before or after the call to `advance_time()`.
+ *     @code
+ *       new_a = update_a(a, b, time.get_next_step_size());
+ *       time.advance_time();
+ *       b = update_b(a, b, time.get_previous_step_size());
+ *       a = new_a;
+ *     @endcode
+ *
+ * One thing to note is that, during the update phase, $dt$ is referred to
+ * either **next** or **previous** time step size, depending on whether the
+ * command `advance_time()` has been called yet. The notion of *current* time
+ * step size is ill-defined. In fact, in the update stage the definition of
+ * every variable depends on whether it has been updated yet or not, hence the
+ * name **the inconsistent stage**.
+ *
+ * The following code snippet shows the code sections for the snapshot stage
+ * and the update stage in the context of a complete time-dependent
+ * simulation. This code follows the coding conventions incorporated in the
+ * tutorial examples. Note that even though this example is written in the
+ * format of a `for` loop, it can equivalently be written as a `while` or
+ * `do while` loop (as shown in step-21).
+ * @code
+ * // pre-processing/setup stage {
+ * make_grid();
+ * setup_system();
+ * for (DiscreteTime time(0., 1., 0.1);  // } end pre-processing/setup stage
+ *      time.is_at_end() == false;
+ *      time.advance_time())             // part of the update stage, runs at
+ *                                       // the end of loop body
+ * {
+ *   // snapshot stage {
+ *   const double time_of_simulation = time.get_next_time();
+ *   const double timestep_size      = time.get_next_step_size();
+ *
+ *   std::cout
+ *     << "Timestep: " << time.get_step_number() << " -- "
+ *     << "Solving for the solution at "
+ *     << "t = " << time_of_simulation << " with "
+ *     << "dt = " << timestep_size << "." << std::endl;
+ *   // } end snapshot stage
+ *
+ *   // update stage {
+ *   assemble_system(time_of_simulation, timestep_size);
+ *   solve();
+ *   update_solutions();
+ *   // } end update stage
+ *
+ *   // snapshot stage {
+ *   output_results(time_of_solution);
+ *
+ *   // propose a new timestep size if need be
+ *   // time.set_desired_next_step_size(...);
+ *   // } end snapshot stage
+ * }
+ * @endcode
+ *
  * @author Reza Rastak, 2019
  */
 class DiscreteTime
 {
 public:
   /**
-   * Constructor
+   * Constructor.
+   *
+   * @pre @p start_step_size must be non-negative.
+   *
+   * @note If @p start_step_size is specified as zero, it indicates that the
+   * desired size for the time step will be calculated at a different location
+   * in the code. In this case, the created object cannot increment time until
+   * the step size is changed by calling set_desired_next_step_size().
    */
   DiscreteTime(const double start_time,
                const double end_time,
@@ -96,6 +233,25 @@ public:
   double
   get_current_time() const;
 
+  /**
+   * Return the next time that we would reach if we were to advance the time
+   * by one step.
+   *
+   * @note If the simulation is at the end time, this method returns the
+   * end time.
+   */
+  double
+  get_next_time() const;
+
+  /**
+   * Return the time we were at before `advance_time()` was called last time.
+   *
+   * @note If the simulation is at the start time, this method returns the
+   * start time.
+   */
+  double
+  get_previous_time() const;
+
   /**
    * Return the start time.
    */
@@ -105,23 +261,52 @@ public:
   /**
    * 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.
+   * equality is very important because it allows us to equality-compare
+   * current time with end time and decide whether we have reached the end of
+   * the simulation.
    */
   double
   get_end_time() const;
 
+  /**
+   * Return whether no step has taken place yet.
+   */
+  bool
+  is_at_start() const;
+
+  /**
+   * Return whether time has reached the end time.
+   */
+  bool
+  is_at_end() const;
+
   /**
    * Return the size of the step from current time step to the
    * next. As discussed in the introduction to the class, this is the
    * *actual* time step, and may differ from the *desired* time step
    * set in the constructor or through the
    * set_desired_next_step_size() function.
+   *
+   * @note If the simulation is at the end time, this method returns zero.
    */
   double
   get_next_step_size() const;
 
+  /**
+   * Return the step size of the previous step.
+   *
+   * @note If the simulation is at the start time, this method returns zero.
+   */
+  double
+  get_previous_step_size() const;
+
+  /**
+   * Return the number of times the simulation time has been incremented.
+   * Return zero when the simulation is at the start time.
+   */
+  unsigned int
+  get_step_number() 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
@@ -148,6 +333,10 @@ public:
    * 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.
+   *
+   * @pre The time step size must be nonzero. If the step size is currently
+   * zero, change it by calling set_desired_next_step_size() before calling
+   * advance_time().
    */
   void
   advance_time();
@@ -192,6 +381,17 @@ private:
    * floating-point value of the time exactly matches the end time.
    */
   double next_time;
+
+  /**
+   * The previous time.
+   */
+  double previous_time;
+
+  /**
+   * The step number i.e. the number of times the simulation time ha been
+   * incremented.
+   */
+  unsigned int step_number;
 };
 
 
@@ -214,6 +414,22 @@ DiscreteTime::get_end_time() const
 
 
 
+inline bool
+DiscreteTime::is_at_start() const
+{
+  return step_number == 0;
+}
+
+
+
+inline bool
+DiscreteTime::is_at_end() const
+{
+  return current_time == end_time;
+}
+
+
+
 inline double
 DiscreteTime::get_next_step_size() const
 {
@@ -222,6 +438,14 @@ DiscreteTime::get_next_step_size() const
 
 
 
+inline double
+DiscreteTime::get_previous_step_size() const
+{
+  return current_time - previous_time;
+}
+
+
+
 inline double
 DiscreteTime::get_current_time() const
 {
@@ -229,6 +453,30 @@ DiscreteTime::get_current_time() const
 }
 
 
+
+inline double
+DiscreteTime::get_next_time() const
+{
+  return next_time;
+}
+
+
+
+inline double
+DiscreteTime::get_previous_time() const
+{
+  return previous_time;
+}
+
+
+
+inline unsigned int
+DiscreteTime::get_step_number() const
+{
+  return step_number;
+}
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index d45141f69b6b20f93169e71aa7b618fd574fbf3b..f32d6eee4e4e5e3505893f17f5d1f14268c7e344 100644 (file)
@@ -25,9 +25,9 @@ namespace
   //  - 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)
+  calculate_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());
@@ -49,7 +49,9 @@ DiscreteTime::DiscreteTime(const double 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)}
+  , next_time{calculate_next_time(start_time, start_step_size, end_time)}
+  , previous_time{start_time}
+  , step_number{0}
 {}
 
 
@@ -57,7 +59,7 @@ DiscreteTime::DiscreteTime(const double start_time,
 void
 DiscreteTime::set_next_step_size(const double next_step_size)
 {
-  next_time = get_next_time(current_time, next_step_size, end_time);
+  next_time = calculate_next_time(current_time, next_step_size, end_time);
 }
 
 
@@ -70,8 +72,10 @@ DiscreteTime::advance_time()
            "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();
+  previous_time          = current_time;
   current_time           = next_time;
-  next_time              = get_next_time(current_time, step_size, end_time);
+  ++step_number;
+  next_time = calculate_next_time(current_time, step_size, end_time);
 }
 
 
@@ -79,8 +83,10 @@ DiscreteTime::advance_time()
 void
 DiscreteTime::restart()
 {
-  current_time = start_time;
-  next_time    = get_next_time(current_time, start_step_size, end_time);
+  previous_time = start_time;
+  current_time  = start_time;
+  next_time     = calculate_next_time(current_time, start_step_size, end_time);
+  step_number   = 0;
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 363fbd4a068fdb8fdddfb643bdcb9ce0a828c6f1..39d3e13f1feac5ba6179580f9e3059ba6847b4a7 100644 (file)
 void
 print_time(const DiscreteTime &time)
 {
+  if (time.is_at_start())
+    deallog << "Simulation started" << std::endl;
   deallog << "Current time = " << time.get_current_time()
-          << ", next step size = " << time.get_next_step_size() << std::endl;
+          << ", next = " << time.get_next_time()
+          << ", previous = " << time.get_previous_time()
+          << ", step number = " << time.get_step_number()
+          << ", next step size = " << time.get_next_step_size()
+          << ", previous step size = " << time.get_previous_step_size()
+          << std::endl;
+  if (time.is_at_end())
+    deallog << "Simulation ended" << std::endl;
 }
 
 void
index 24806bc48e18f95dcb5ce8edc3fc3d5293fe117e..6110561f2376ab72198ff29a4eacd7ec72d64455 100644 (file)
@@ -1,27 +1,32 @@
 
 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::Simulation started
+DEAL:Start to end::Current time = 0.00000, next = 0.123000, previous = 0.00000, step number = 0, next step size = 0.123000, previous step size = 0.00000
+DEAL:Start to end::Current time = 0.123000, next = 0.246000, previous = 0.00000, step number = 1, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.246000, next = 0.369000, previous = 0.123000, step number = 2, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.369000, next = 0.492000, previous = 0.246000, step number = 3, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.492000, next = 0.615000, previous = 0.369000, step number = 4, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.615000, next = 0.738000, previous = 0.492000, step number = 5, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.738000, next = 0.861000, previous = 0.615000, step number = 6, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.861000, next = 0.984000, previous = 0.738000, step number = 7, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 0.984000, next = 1.10700, previous = 0.861000, step number = 8, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 1.10700, next = 1.23000, previous = 0.984000, step number = 9, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 1.23000, next = 1.35300, previous = 1.10700, step number = 10, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 1.35300, next = 1.47600, previous = 1.23000, step number = 11, next step size = 0.123000, previous step size = 0.123000
+DEAL:Start to end::Current time = 1.47600, next = 1.50000, previous = 1.35300, step number = 12, next step size = 0.0240000, previous step size = 0.123000
+DEAL:Start to end::Current time = 1.50000, next = 1.50000, previous = 1.47600, step number = 13, next step size = 0.00000, previous step size = 0.0240000
+DEAL:Start to end::Simulation ended
 DEAL:Start to end::Restarted
-DEAL:Start to end::Current time = 0.00000, next step size = 0.123000
+DEAL:Start to end::Simulation started
+DEAL:Start to end::Current time = 0.00000, next = 0.123000, previous = 0.00000, step number = 0, next step size = 0.123000, previous step size = 0.00000
 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::Simulation started
+DEAL:Adjust time step size::Current time = 0.400000, next = 0.550000, previous = 0.400000, step number = 0, next step size = 0.150000, previous step size = 0.00000
+DEAL:Adjust time step size::Current time = 0.550000, next = 0.700000, previous = 0.400000, step number = 1, next step size = 0.150000, previous step size = 0.150000
+DEAL:Adjust time step size::Current time = 0.910000, next = 1.27000, previous = 0.550000, step number = 2, next step size = 0.360000, previous step size = 0.360000
+DEAL:Adjust time step size::Current time = 1.27000, next = 1.63000, previous = 0.910000, step number = 3, next step size = 0.360000, previous step size = 0.360000
+DEAL:Adjust time step size::Current time = 1.88000, next = 2.10000, previous = 1.27000, step number = 4, next step size = 0.220000, previous step size = 0.610000
+DEAL:Adjust time step size::Current time = 2.10000, next = 2.10000, previous = 1.88000, step number = 5, next step size = 0.00000, previous step size = 0.220000
+DEAL:Adjust time step size::Simulation ended
 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.