34 double timeLimit)
const {
41 const unsigned int p = eeStepper->
order();
42 const double deltaMax = T*std::pow(S/Rmax, (
int)(p+1));
43 const double TINY = 1.0E-30;
48 d.
time= timeLimit==0? s.time+stepsize : timeLimit;
53 double h = d.
time-s.time;
55 std::vector<double> errors;
56 eeStepper->
step(data, s, d, errors);
57 if (timeLimit!=0.0)
return;
60 for (
size_t e=0;e<errors.size();e++) errors[e] = fabs(errors[e]);
63 double delta = (*std::max_element(errors.begin(),errors.end()));
68 h = std::max(S*h*std::pow(T/(delta + TINY), 1.0/(p+1)),Rmin*h);
69 if (!(((
double) (s.time+h) - (
double) s.time) > 0) ) {
70 throw std::runtime_error(
"Warning, RK Integrator step underflow");
77 if (delta < deltaMax) {
78 hnext = S*h*std::pow(T/(delta + TINY),1.0/(p+1));