41 std::vector<double>& rCurrentYValues,
42 std::vector<double>& rNextYValues)
52 if (num_equations !=
k1.size())
54 k1.resize(num_equations);
55 k2.resize(num_equations);
56 k3.resize(num_equations);
57 k4.resize(num_equations);
58 yki.resize(num_equations);
61 std::vector<double>& dy = rNextYValues;
64 for (
unsigned i=0; i<num_equations; i++)
66 k1[i] = timeStep*dy[i];
67 yki[i] = rCurrentYValues[i] + 0.5*
k1[i];
71 for (
unsigned i=0; i<num_equations; i++)
73 k2[i] = timeStep*dy[i];
74 yki[i] = rCurrentYValues[i] + 0.5*
k2[i];
78 for (
unsigned i=0; i<num_equations; i++)
80 k3[i] = timeStep*dy[i];
81 yki[i] = rCurrentYValues[i] +
k3[i];
85 for (
unsigned i=0; i<num_equations; i++)
87 k4[i] = timeStep*dy[i];
88 rNextYValues[i] = rCurrentYValues[i] + (
k1[i]+2*
k2[i]+2*
k3[i]+
k4[i])/6.0;