]> https://gitweb.dealii.org/ - dealii.git/commitdiff
arkode ode solver internally uses DiscreteTime 9973/head
authorReza Rastak <rastak@stanford.edu>
Mon, 27 Apr 2020 22:28:22 +0000 (15:28 -0700)
committerReza Rastak <rastak@stanford.edu>
Tue, 21 Jul 2020 04:39:05 +0000 (21:39 -0700)
include/deal.II/sundials/arkode.h
source/sundials/arkode.cc
tests/sundials/harmonic_oscillator_05.output
tests/sundials/harmonic_oscillator_06.output

index 648411dd41eed9e5dfa4425396c54bd48e3edfbf..c327ff5249015a2733d24206d456cfe657a89069 100644 (file)
@@ -325,7 +325,7 @@ namespace SUNDIALS
        * @param initial_time Initial time
        * @param final_time Final time
        * @param initial_step_size Initial step size
-       * @param output_period Time interval between each output
+       * @param output_period Desired time interval between each output
        *
        * Running parameters:
        *
@@ -475,7 +475,8 @@ namespace SUNDIALS
       unsigned int maximum_order;
 
       /**
-       * Time period between each output.
+       * Desired time period between each output. The actual output time period
+       * may be adjusted by Arkode.
        */
       double output_period;
 
index cec5de071f87d6911ae5fb158e1c7547dd5b43c7..64be7b34b266823ea161fdd36c8ef422f3e343e9 100644 (file)
@@ -20,6 +20,7 @@
 
 #ifdef DEAL_II_WITH_SUNDIALS
 
+#  include <deal.II/base/discrete_time.h>
 #  include <deal.II/base/utilities.h>
 
 #  include <deal.II/lac/block_vector.h>
@@ -270,18 +271,9 @@ namespace SUNDIALS
   unsigned int
   ARKode<VectorType>::solve_ode(VectorType &solution)
   {
-    unsigned int system_size = solution.size();
+    const unsigned int system_size = solution.size();
 
-    double       t           = data.initial_time;
-    double       h           = data.initial_step_size;
-    unsigned int step_number = 0;
-
-    int status;
-    (void)status;
-
-    // The solution is stored in
-    // solution. Here we take only a
-    // view of it.
+    // The solution is stored in solution. Here we take only a view of it.
 #  ifdef DEAL_II_WITH_MPI
     if (is_serial_vector<VectorType>::value == false)
       {
@@ -302,33 +294,38 @@ namespace SUNDIALS
         yy        = N_VNew_Serial(system_size);
         abs_tolls = N_VNew_Serial(system_size);
       }
-    reset(data.initial_time, data.initial_step_size, solution);
 
-    double next_time = data.initial_time;
+    DiscreteTime time(data.initial_time,
+                      data.final_time,
+                      data.initial_step_size);
+    reset(time.get_current_time(), time.get_next_step_size(), solution);
 
     if (output_step)
-      output_step(0, solution, 0);
+      output_step(time.get_current_time(), solution, time.get_step_number());
 
-    while (t < data.final_time)
+    while (time.is_at_end() == false)
       {
-        next_time += data.output_period;
-
-        status = SundialsARKode(arkode_mem, next_time, yy, &t, ARK_NORMAL);
-
-        AssertARKode(status);
-
-        status = ARKodeGetLastStep(arkode_mem, &h);
+        time.set_desired_next_step_size(data.output_period);
+        double     actual_next_time;
+        const auto status = SundialsARKode(
+          arkode_mem, time.get_next_time(), yy, &actual_next_time, ARK_NORMAL);
+        (void)status;
         AssertARKode(status);
 
         copy(solution, yy);
 
-        while (solver_should_restart(t, solution))
-          reset(t, h, solution);
+        time.set_next_step_size(actual_next_time - time.get_current_time());
+        time.advance_time();
 
-        step_number++;
+        while (solver_should_restart(time.get_current_time(), solution))
+          reset(time.get_current_time(),
+                time.get_previous_step_size(),
+                solution);
 
         if (output_step)
-          output_step(t, solution, step_number);
+          output_step(time.get_current_time(),
+                      solution,
+                      time.get_step_number());
       }
 
       // Free the vectors which are no longer used.
@@ -345,7 +342,7 @@ namespace SUNDIALS
         N_VDestroy_Serial(abs_tolls);
       }
 
-    return step_number;
+    return time.get_step_number();
   }
 
   template <typename VectorType>
index 507044c3b2e81897e6451c572336039ff1900e2a..04f8401db1e25f064da8d0de06cd34d1a8bae640 100644 (file)
@@ -99,4 +99,3 @@
 9.8 0.991567 2.63469 2.49996
 9.9 1.02585 2.61958 2.50009
 10 1.06497 2.59596 2.49997
-10 1.06497 2.59596 2.49997
index 507044c3b2e81897e6451c572336039ff1900e2a..04f8401db1e25f064da8d0de06cd34d1a8bae640 100644 (file)
@@ -99,4 +99,3 @@
 9.8 0.991567 2.63469 2.49996
 9.9 1.02585 2.61958 2.50009
 10 1.06497 2.59596 2.49997
-10 1.06497 2.59596 2.49997

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.