Hi, we have posted a review request for odeint (a library for solving ordinary differential equations) on this list some days ago. An ODE is dx/dt = f(x(t)) and one is looking for the solution x(t) which is a function of t. Solving ODEs is usually an iterative task, so one starts with x(0) and iteratively applies a solver to obtain a discretized (and approximated) representation of the solution: x(0) -> x(dt) -> x(2*dt) -> x(3*dt) -> ... I played around with an iterator which is doing exactly this iterative procedure. For example, it can be used via odeint::runge_kutta4< state_type > stepper; // ode solver std::array< double , 3 > x = {{ 10.0 , 10.0 , 10.0 }}; // initial state double res = std::accumulate( // parameters : solver , ode , state , time , time_step make_const_step_iterator_begin(stepper, lorenz(), x, 0.0, 0.01) , make_const_step_iterator_end( stepper, lorenz(), x, 1.0, 0.01) , 0.0 , []( double sum, const state_type &x) { return sum + x[0]; } ); The first iterator increments the time until the time of the second iterator is rearched t=0.0 -> t=0.01 -> t=0.02 -> .. t=1.0 and the first component of the solution is accumulated. In principle, this approach works quite well, but there are some problems at least semantically. Maybe someone here has some comments or ideas on this: The defintion of this iterator is: template< class Stepper , class System , class StepperTag > class const_step_iterator : public boost::iterator_facade < const_step_iterator< Stepper , System , StepperTag > , typename Stepper::state_type const , boost::single_pass_traversal_tag
{ ... }; It is a single pass iterator, so it can only check for equality of two iterators. This is problematic since equality here means that the time (t) of the begin iterator is smaller then the time of the end iterator. To ensure commutativity of the != or == operation I needed to flag the begin and end iterator explicitly. This is the reason for the two factory functions make_const_step_iterator_begin and make_const_step_iterator_end, which are doing this flagging. The second semantically problem is that the end iterator in principle does not need to know the stepper as well as the system (lorenz() in the above example). But all algorithms from the standard library and Boost.Range assume that the begin and end iterator are of same type. Therefore you have to put the stepper and the system into the end iterator too. Of course, there are techniques like type erasure which might solve this problem, but they introduce some run-time overhead which we want to avoid here. Any ideas or comments about this? P.S. The code can be found in the iterator branch in odeint: boost/numeric/odeint/iterator/const_step_iterator.hpp and some examples libs/numeric/odeint/examples/const_step_iterator.cpp git clone -b iterator git@github.com:headmyshoulder/odeint-v2.git