Do controlled error steppers in boost odeint support complex data types?

Question

I encountered a problem with the odeint library using a controlled error stepper together with a complex state type. The code from the example with the complex stuart landau equation, I modified such that it includes an adaptive integrator. The code looks like this now:

#include <iostream>
#include <complex>
#include <boost/array.hpp>

#include <boost/numeric/odeint.hpp>

using namespace std;
using namespace boost::numeric::odeint;

//[ stuart_landau_system_function
typedef complex< double > state_type;

struct stuart_landau
{
    double m_eta;
    double m_alpha;

    stuart_landau( double eta = 1.0 , double alpha = 1.0 )
    : m_eta( eta ) , m_alpha( alpha ) { }

    void operator()( const state_type &x , state_type &dxdt , double t ) const
    {
        const complex< double > I( 0.0 , 1.0 );
        dxdt = ( 1.0 + m_eta * I ) * x - ( 1.0 + m_alpha * I ) * norm( x ) * x;
    }
};
//]


struct streaming_observer
{
    std::ostream& m_out;

    streaming_observer( std::ostream &out ) : m_out( out ) { }

    template< class State >
    void operator()( const State &x , double t ) const
    {
        m_out.precision(10);
        m_out << t;
        m_out << "\t" << x.real() << "\t" << x.imag() ;
        m_out << "\n";
    }
};


int main( int argc , char **argv )
{
    //[ stuart_landau_integration
    state_type x = complex< double >( 1.0 , 0.0 );

    bulirsch_stoer< state_type > stepper( 1E-12 , 1E-12 , 1 , 1 );

    const double dt = 0.1;

    //]
    integrate_adaptive( stepper , stuart_landau( 2.0 , 1.0 ) , x , 0.0 , 10.0 , dt , streaming_observer( cout ) );

    return 0;
}

However, if I change the definition of the stepper to

bulirsch_stoer< state_type, complex< double > > stepper( 1E-12 , 1E-12 , 1 , 1 );

compilation fails. My question is: Are complex data types not supported in controlled error steppers? If so, is there a method to circumvent the problem, which arises. Or, is it possible to define an own vector algebra for the complex data type?

Answer 1

The odeint library does support complex data types, including controlled steppers.

As an example, to obtain a controlled stepper in the case of the 4th/5th order Dormand-Prince method you can do this:

runge_kutta_dopri5< state_type > stepper;
auto c_stepper = make_controlled(1.E-12, 1.E-12, stepper); 
integrate_adaptive(c_stepper, stuart_landau( 2.0 , 1.0 ) , x , 0.0 , 10.0 , dt , streaming_observer( cout ) );

Where state_type is typedef'd as

typedef complex< double > state_type;

The Bulirsch-Stoer algorithm that is implemented in odeint, however, is already inherently a controlled stepper. It is therefore not possible to apply make_controlled or other strategies to create controlled steppers in this case. There is no such template in the odeint library because it would not make sense.

The first code that you posted already integrates the ODE, with complex values and with a controlled stepper. It is not clear what you are trying to achieve by changing the stepper to

bulirsch_stoer< complex<double>, complex<double> > stepper( 1E-12 , 1E-12 , 1 , 1 );