RungeKutta/src/RK4OneStepTempl.h

0001 #ifndef RK4OneStepTempl_H
0002 #define RK4OneStepTempl_H
0003
0004 #include "FWCore/Utilities/interface/Visibility.h"
0005 #include "RKSmallVector.h"
0006 #include "RKDerivative.h"
0007
0008 template <typename T, int N>
0009 class dso_internal RK4OneStepTempl {
0010 public:
0011   typedef T Scalar;
0012   typedef RKSmallVector<T, N> Vector;
0013
0014   Vector operator()(Scalar startPar, const Vector& startState, const RKDerivative<T, N>& deriv, Scalar step) const {
0015     // cout << "RK4OneStepTempl: starting from " << startPar << startState << endl;
0016
0017     Vector k1 = step * deriv(startPar, startState);
0018     Vector k2 = step * deriv(startPar + step / 2, startState + k1 / 2);
0019     Vector k3 = step * deriv(startPar + step / 2, startState + k2 / 2);
0020     Vector k4 = step * deriv(startPar + step, startState + k3);
0021
0022     Vector result = startState + k1 / 6 + k2 / 3 + k3 / 3 + k4 / 6;
0023
0024     // cout << "RK4OneStepTempl: result for step " << step << " is " << result << endl;
0025
0026     return result;
0027   }
0028 };
0029
0030 #endif