]> https://gitweb.dealii.org/ - dealii.git/commitdiff
set_next_step_size is renamed to set_desired_next_step_size 9748/head
authorReza Rastak <rastak@stanford.edu>
Sun, 19 Apr 2020 06:25:52 +0000 (23:25 -0700)
committerReza Rastak <rastak@stanford.edu>
Thu, 23 Apr 2020 01:35:01 +0000 (18:35 -0700)
examples/step-21/step-21.cc
include/deal.II/base/discrete_time.h
source/base/discrete_time.cc
tests/base/discrete_time_1.cc

index 9a45cf21fbf07a843ddfdcd89822646f0565850f..0496308840dfa96aa51dd725931bd6668b834527 100644 (file)
@@ -503,10 +503,12 @@ namespace Step21
   // @sect4{TwoPhaseFlowProblem::TwoPhaseFlowProblem}
 
   // First for the constructor. We use $RT_k \times DQ_k \times DQ_k$
-  // spaces. The time step is set to zero initially, but will be computed
-  // before it is needed first, as described in a subsection of the
-  // introduction. The time object internally prevents itself from being
-  // incremented with $dt = 0$, forcing us to set a non-zero desired size for
+  // spaces. For initializing the DiscreteTime object, we don't set the time
+  // step size in the constructor because we don't have its value yet.
+  // The time step size is initially set to zero, but it will be computed
+  // before it is needed to increment time, as described in a subsection of
+  // the introduction. The time object internally prevents itself from being
+  // incremented when $dt = 0$, forcing us to set a non-zero desired size for
   // $dt$ before advancing time.
   template <int dim>
   TwoPhaseFlowProblem<dim>::TwoPhaseFlowProblem(const unsigned int degree)
@@ -519,7 +521,7 @@ namespace Step21
          1)
     , dof_handler(triangulation)
     , n_refinement_steps(5)
-    , time(/*start time*/ 0., /*end time*/ 1., /*time step*/ 0.)
+    , time(/*start time*/ 0., /*end time*/ 1.)
     , viscosity(0.2)
   {}
 
@@ -970,11 +972,15 @@ namespace Step21
     //
     // The maximal velocity we compute using a helper function to compute the
     // maximal velocity defined below, and with all this we can evaluate our
-    // new time step length. The method DiscreteTime::set_next_time_step() is
-    // used to assign the new value of the time step to the DiscreteTime
-    // object.
-    time.set_next_step_size(std::pow(0.5, double(n_refinement_steps)) /
-                            get_maximal_velocity());
+    // new time step length. We use the method
+    // DiscreteTime::set_desired_next_time_step() to suggest the new
+    // calculated value of the time step to the DiscreteTime object. In most
+    // cases, the time object uses the exact provided value to increment time.
+    // It some case, the step size may be modified further by the time object.
+    // For example, if the calculated time increment overshoots the end time,
+    // it is truncated accordingly.
+    time.set_desired_next_step_size(std::pow(0.5, double(n_refinement_steps)) /
+                                    get_maximal_velocity());
 
     // The next step is to assemble the right hand side, and then to pass
     // everything on for solution. At the end, we project back saturations
index c5ad7dc3d7c1b4519ea245a0bd6766d3d5bfdfde..6337668446780b1ed8a30c3b9fbc751afb6b32a7 100644 (file)
@@ -218,16 +218,26 @@ public:
   /**
    * Constructor.
    *
-   * @pre @p start_step_size must be non-negative.
+   * @param[in] start_time The time at the start of the simulation.
    *
-   * @note If @p start_step_size is specified as zero, it indicates that the
+   * @param[in] end_time The time at the end of the simulation.
+   *
+   * @param[in] desired_start_step_size A desired step size for incrementing
+   * time for the first step. It is not guaranteed that this value will be
+   * actually used as the size of the first step, as discussed in the
+   * introduction.
+   *
+   * @pre @p desired_start_step_size must be non-negative.
+   *
+   * @note @p desired_start_step_size is an optional parameter. If it is not
+   * provided or it 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,
-               const double start_step_size);
+               const double desired_start_step_size = 0.);
 
   /**
    * Return the current time.
@@ -310,9 +320,10 @@ public:
   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
-   * the simulation time. However, if the step is too large such that the next
+   * Set the *desired* value of the next time step size. By calling this
+   * method, we are indicating the the next time advance_time() is called, we
+   * would like @p time_step_size to 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
@@ -320,16 +331,16 @@ public:
    * end time.
    */
   void
-  set_next_step_size(const double time_step_size);
+  set_desired_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.
+   * set_desired_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.
+   * time. See the documentation of set_desired_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
@@ -361,11 +372,6 @@ private:
    */
   const double end_time;
 
-  /**
-   * The size of the first step.
-   */
-  const double start_step_size;
-
   /**
    * The current time.
    */
@@ -375,12 +381,13 @@ private:
    * 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.
+   * current step size. For example, when the method
+   * set_desired_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;
 
@@ -389,6 +396,11 @@ private:
    */
   double previous_time;
 
+  /**
+   * The size of the first step.
+   */
+  const double start_step_size;
+
   /**
    * The step number i.e. the number of times the simulation time ha been
    * incremented.
index f32d6eee4e4e5e3505893f17f5d1f14268c7e344..f1bbea79c0355b21589baeb289e46f13a685a983 100644 (file)
@@ -44,20 +44,22 @@ namespace
 
 DiscreteTime::DiscreteTime(const double start_time,
                            const double end_time,
-                           const double start_step_size)
+                           const double desired_start_step_size)
   : start_time{start_time}
   , end_time{end_time}
-  , start_step_size{start_step_size}
   , current_time{start_time}
-  , next_time{calculate_next_time(start_time, start_step_size, end_time)}
+  , next_time{calculate_next_time(start_time,
+                                  desired_start_step_size,
+                                  end_time)}
   , previous_time{start_time}
+  , start_step_size{next_time - start_time}
   , step_number{0}
 {}
 
 
 
 void
-DiscreteTime::set_next_step_size(const double next_step_size)
+DiscreteTime::set_desired_next_step_size(const double next_step_size)
 {
   next_time = calculate_next_time(current_time, next_step_size, end_time);
 }
index 39d3e13f1feac5ba6179580f9e3059ba6847b4a7..a119b4eaae9b6a732beacd7c83448f101b7952e7 100644 (file)
@@ -42,7 +42,7 @@ test_from_start_to_end()
 
   DiscreteTime time(/*start_time*/ 0.,
                     /*end_time*/ 1.5,
-                    /*start_step_size*/ 0.123);
+                    /*desired_start_step_size*/ 0.123);
 
   deallog << "Start time = " << time.get_start_time() << std::endl;
   deallog << "End time = " << time.get_end_time() << std::endl;
@@ -65,16 +65,16 @@ test_adjust_time_step_size()
 
   DiscreteTime time(/*start_time*/ 0.4,
                     /*end_time*/ 2.1,
-                    /*start_step_size*/ 0.15);
+                    /*desired_start_step_size*/ 0.15);
   print_time(time);
   time.advance_time();
   print_time(time);
-  time.set_next_step_size(0.36);
+  time.set_desired_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.set_desired_next_step_size(0.61);
   time.advance_time();
   print_time(time);
   time.advance_time(); // Here we reach the end time.

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.