SpECTRE Documentation Coverage Report
Current view: top level - Time/TimeSteppers - DormandPrince5.cpp Hit Total Coverage
Commit: d07bc7e978f0161d6eb259d687cd7878b5f9f6f0 Lines: 0 1 0.0 %
Date: 2020-09-24 02:51:00
Legend: Lines: hit not hit

          Line data    Source code
       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #include "Time/TimeSteppers/DormandPrince5.hpp"
       5             : 
       6             : #include <cmath>
       7             : #include <limits>
       8             : 
       9             : #include "ErrorHandling/Assert.hpp"
      10             : #include "Time/TimeStepId.hpp"
      11             : #include "Utilities/Gsl.hpp"
      12             : 
      13             : namespace TimeSteppers {
      14             : 
      15             : uint64_t DormandPrince5::number_of_substeps() const noexcept { return 6; }
      16             : 
      17             : size_t DormandPrince5::number_of_past_steps() const noexcept { return 0; }
      18             : 
      19             : // The growth function for DP5 is
      20             : //
      21             : //   g = mu^6 / 600 + \sum_{n=0}^5 mu^n / n!,
      22             : //
      23             : // where mu = lambda * dt. The equation dy/dt = -lambda * y evolves
      24             : // stably if |g| < 1. For lambda=-2, chosen so the stable_step() for
      25             : // RK1 (i.e. forward Euler) would be 1, DP5 has a stable step
      26             : // determined by inserting mu->-2 dt into the above equation. Finding the
      27             : // solutions with a numerical root find yields a stable step of about 1.653.
      28             : double DormandPrince5::stable_step() const noexcept {
      29             :   return 1.6532839463174733;
      30             : }
      31             : 
      32             : TimeStepId DormandPrince5::next_time_id(const TimeStepId& current_id,
      33             :                                         const TimeDelta& time_step) const
      34             :     noexcept {
      35             :   const auto& step = current_id.substep();
      36             :   const auto& t0 = current_id.step_time();
      37             :   const auto& t = current_id.substep_time();
      38             :   if (step < 6) {
      39             :     if (step == 0) {
      40             :       ASSERT(t == t0, "In DP5 substep 0, the substep time ("
      41             :                           << t << ") should equal t0 (" << t0 << ")");
      42             :     } else {
      43             :       ASSERT(t == t0 + gsl::at(_c, step - 1) * time_step,
      44             :              "In DP5 substep " << step << ", the substep time (" << t
      45             :                                << ") should equal t0+c[" << step - 1 << "]*dt ("
      46             :                                << t0 + gsl::at(_c, step - 1) * time_step
      47             :                                << ")");
      48             :     }
      49             :     if (step < 5) {
      50             :       return {current_id.time_runs_forward(), current_id.slab_number(), t0,
      51             :               step + 1, t0 + gsl::at(_c, step) * time_step};
      52             :     } else {
      53             :       return {current_id.time_runs_forward(), current_id.slab_number(),
      54             :               t0 + time_step};
      55             :     }
      56             :   } else {
      57             :     ERROR("In DP5 substep should be one of 0,1,2,3,4,5, not "
      58             :           << current_id.substep());
      59             :   }
      60             : }
      61             : 
      62             : constexpr double DormandPrince5::_a2;
      63             : constexpr std::array<double, 2> DormandPrince5::_a3;
      64             : constexpr std::array<double, 3> DormandPrince5::_a4;
      65             : constexpr std::array<double, 4> DormandPrince5::_a5;
      66             : constexpr std::array<double, 5> DormandPrince5::_a6;
      67             : constexpr std::array<double, 6> DormandPrince5::_b;
      68             : const std::array<Time::rational_t, 5> DormandPrince5::_c = {
      69             :     {{1, 5}, {3, 10}, {4, 5}, {8, 9}, {1, 1}}};
      70             : constexpr std::array<double, 6> DormandPrince5::_d;
      71             : }  // namespace TimeSteppers
      72             : 
      73             : /// \cond
      74             : PUP::able::PUP_ID TimeSteppers::DormandPrince5::my_PUP_ID =  // NOLINT
      75             :     0;
      76             : /// \endcond

Generated by: LCOV version 1.14