ode-toolbox

commit b405ac12b8f10ee58761ebbe97754cfebb41d53a
parent 646806d0a141736cfcd88ea9688b69b1159ac11b
Author: Fabian Wermelinger <fabianw@mavt.ethz.ch>
Date:   Sun, 31 Jan 2021 23:42:47 +0100

Add dtMin argument for varibale step integrators

Diffstat:
M
build.sh
 | 
2
+-
M
include/TimeStepper/StepperBase.h
 | 
2
++
M
include/TimeStepper/explicit/variableStep/RKF45.h
 | 
9
+++++++++
M
include/TimeStepper/explicit/variableStep/RKV56.h
 | 
9
+++++++++

4 files changed, 21 insertions(+), 1 deletion(-)
diff --git a/build.sh b/build.sh
@@ -1,3 +1,3 @@
 #!/usr/bin/env bash
-meson setup build --buildtype release
+meson setup build --buildtype release --prefix $HOME/.local -Dtest=true
 meson compile -C build
diff --git a/include/TimeStepper/StepperBase.h b/include/TimeStepper/StepperBase.h
@@ -48,6 +48,7 @@ struct StepperSettings
         safety           = p("-ts_safety").asDouble(0.9);
         t                = p("-ts_t0").asDouble(0.0);
         dt               = p("-ts_dt").asDouble(1.0e-4);
+        dtMin            = p("-ts_dtMin").asDouble(0.0);
         step             = 0;
         nsteps           = p("-ts_nsteps").asInt(0);
         p.set_strict_mode();
@@ -75,6 +76,7 @@ struct StepperSettings
     Real safety;
     Real t;
     Real dt;
+    Real dtMin;
     size_t step;
     size_t nsteps;
     Real tFinal;
diff --git a/include/TimeStepper/explicit/variableStep/RKF45.h b/include/TimeStepper/explicit/variableStep/RKF45.h
@@ -76,6 +76,15 @@ class RKF45 : public Euler<Tinput, Trhs>
     {
         const Real err = _error();
 
+        const Real dtMin = m_settings.dtMin;
+        if (dt <= dtMin) {
+            dt_next = 1.1 * dtMin;
+            m_errold = std::max(err, static_cast<Real>(1.0e-4));
+            m_reject = false;
+            m_firstPass = true;
+            return true;
+        }
+
         if (dt < std::numeric_limits<Real>::epsilon())
         {
             // if you request the impossible
diff --git a/include/TimeStepper/explicit/variableStep/RKV56.h b/include/TimeStepper/explicit/variableStep/RKV56.h
@@ -83,6 +83,15 @@ class RKV56 : public Euler<Tinput, Trhs>
     {
         const Real err = _error();
 
+        const Real dtMin = m_settings.dtMin;
+        if (dt <= dtMin) {
+            dt_next = 1.1 * dtMin;
+            m_errold = std::max(err, static_cast<Real>(1.0e-4));
+            m_reject = false;
+            m_firstPass = true;
+            return true;
+        }
+
         if (dt < std::numeric_limits<Real>::epsilon())
         {
             // if you request the impossible