Set current state in odeint

Question

I would like to set the current_state of a stepper in boost::odeint.

I have the following MWE:

#include <queue>
#include <boost/numeric/odeint.hpp>
#include <array>

typedef std::array< double , 2 > state_type;

void eqsys(const state_type &x, state_type &dxdt, double t) {
    dxdt[0] = 0.3*x[1];
    dxdt[1] = 4*x[0];
}

int main() {
    int steps=0;
    auto stepper=make_dense_output(1.0e-6,1.0e-6,boost::numeric::odeint::runge_kutta_dopri5< state_type >() );
    state_type x={2,2};
    double stoptime=10;

    stepper.initialize(x,0,0.01);

    while(true){
        //Run stepper
        while( stepper.current_time()<stoptime && steps<10000 ) {
            stepper.do_step(eqsys);
            ++steps;
        }
        x=stepper.current_state();
        x[0]=2;
        stepper.set_current_state(x);
        stoptime+=10;
    }
}

The line stepper.set_current_state(x); is the one that doesn't work. Is there a command I can use to accomplish this?

I realize I could probably just use:

stepper.initialize(x, stepper.current_time(), 0.01);

but's it's unclear if this is really the best strategy, or if I'm losing some internal state of the stepper (such as the current stepsize or error estimate) which it would be better to keep.

Answer 1

Using

stepper.initialize(x, stepper.current_time(), 0.01);

is absolutely ok. The internal state is the set correctly in this case.