14 double timeLimit )
const {
15 const double h = timeLimit==0 ? stepsize : timeLimit - s.
time;
16 if (h<=0)
throw std::runtime_error (
"SimpleRKStepper: negative stepsize");
17 const unsigned int nvar = s.
variable.size();
20 std::vector<std::vector<double> >k(tableau.
nSteps());
21 for (
unsigned int i=0;i<tableau.
nSteps();i++) {
24 for (
unsigned int v=0;v<nvar;v++) arg[v]=s.
variable[v];
25 for (
unsigned int j=0;j<i;j++) {
26 for (
unsigned int v=0;v<nvar;v++) arg[v] += h*tableau.
A(i,j)*k[j][v];
28 for (
unsigned int v=0;v<nvar;v++) k[i][v]=(*data->
_diffEqn[v])(arg);
34 for (
unsigned int i=0;i<tableau.
nSteps();i++) {
35 for (
unsigned int v=0;v<nvar;v++) d.
firstDerivative[v] += tableau.
b(i)*k[i][v];
38 d.
time = timeLimit==0 ? s.
time + h : timeLimit;