51 std::vector<double>& rCurrentYValues,
52 std::vector<double>& rNextYValues)
59 const double delta = 1.0e-8;
63 if (
mEvalF.size() != num_equations)
65 mEvalF.resize(num_equations);
67 mTemp.resize(num_equations);
73 for (
unsigned i=0; i<num_equations; i++)
75 tempY= rCurrentYValues[i];
76 rCurrentYValues[i]=tempY+delta;
79 rCurrentYValues[i]=tempY;
82 for (
unsigned i=0; i<num_equations; i++)
87 rNextYValues[i] = rCurrentYValues[i]+
mEvalF[i]*timeStep;