SpECTRE Documentation Coverage Report
Current view: top level - Time/TimeSteppers - RungeKutta4.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/RungeKutta4.hpp"
       5             : 
       6             : #include <cmath>
       7             : #include <limits>
       8             : 
       9             : #include "ErrorHandling/Assert.hpp"
      10             : #include "Time/TimeStepId.hpp"
      11             : 
      12             : namespace TimeSteppers {
      13             : 
      14             : uint64_t RungeKutta4::number_of_substeps() const noexcept { return 4; }
      15             : 
      16             : size_t RungeKutta4::number_of_past_steps() const noexcept { return 0; }
      17             : 
      18             : // The growth function for RK4 is (e.g. page 60 of
      19             : // http://www.staff.science.uu.nl/~frank011/Classes/numwisk/ch10.pdf
      20             : //
      21             : //   g = 1 + mu + mu^2 / 2 + mu^3 / 6 + mu^4 / 24,
      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, RK4 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.39265.
      28             : double RungeKutta4::stable_step() const noexcept { return 1.3926467817026411; }
      29             : 
      30             : TimeStepId RungeKutta4::next_time_id(const TimeStepId& current_id,
      31             :                                      const TimeDelta& time_step) const
      32             :     noexcept {
      33             :   switch (current_id.substep()) {
      34             :     case 0:
      35             :       ASSERT(current_id.substep_time() == current_id.step_time(),
      36             :              "In RK4 substep 0, the substep time ("
      37             :                  << current_id.substep_time() << ") should equal t0 ("
      38             :                  << current_id.step_time() << ")");
      39             :       return {current_id.time_runs_forward(), current_id.slab_number(),
      40             :               current_id.step_time(), 1,
      41             :               current_id.step_time() + time_step / 2.0};
      42             :     case 1:
      43             :       ASSERT(
      44             :           current_id.substep_time() == current_id.step_time() + time_step / 2.0,
      45             :           "In RK4 substep 1, the substep time ("
      46             :               << current_id.substep_time() << ") should equal t0+0.5*dt ("
      47             :               << current_id.step_time() + time_step / 2.0 << ")");
      48             :       return {current_id.time_runs_forward(), current_id.slab_number(),
      49             :               current_id.step_time(), 2,
      50             :               current_id.step_time() + time_step / 2.0};
      51             :     case 2:
      52             :       ASSERT(
      53             :           current_id.substep_time() == current_id.step_time() + time_step / 2.0,
      54             :           "In RK4 substep 2, the substep time ("
      55             :               << current_id.substep_time() << ") should equal t0+0.5*dt ("
      56             :               << current_id.step_time() + time_step / 2.0 << ")");
      57             :       return {current_id.time_runs_forward(), current_id.slab_number(),
      58             :               current_id.step_time(), 3, current_id.step_time() + time_step};
      59             :     case 3:
      60             :       ASSERT(current_id.substep_time() == current_id.step_time() + time_step,
      61             :              "In RK4 substep 3, the substep time ("
      62             :                  << current_id.substep_time() << ") should equal t0+dt ("
      63             :                  << current_id.step_time() + time_step << ")");
      64             :       return {current_id.time_runs_forward(), current_id.slab_number(),
      65             :               current_id.step_time() + time_step};
      66             :     default:
      67             :       ERROR("In RK4 substep should be one of 0,1,2,3, not "
      68             :             << current_id.substep());
      69             :   }
      70             : }
      71             : 
      72             : }  // namespace TimeSteppers
      73             : 
      74             : /// \cond
      75             : PUP::able::PUP_ID TimeSteppers::RungeKutta4::my_PUP_ID =  // NOLINT
      76             :     0;
      77             : /// \endcond

Generated by: LCOV version 1.14