odeint: сложные и высокоточные типы

Может ли кто-нибудь оказать мне некоторую помощь в создании минимального примера использования решателя c ++ odeint с вектором комплексных чисел и, если возможно, с более высокой точностью (скажем, boost.multiprecision или libquadmath __float128, __complex128).

В документации есть пример использования сложного скаляра http://headmyshoulder.github.com/odeint-v2/doc/boost_numeric_odeint/tutorial/special_topics.html

и там упоминается, что:

The fact that we have to configure a different algebra is solely due to the fact
that we use a non-vector state type and not to the usage of complex values.
So for, e.g. vector<  complex<double> >, this would not be required.

Я пытался изменить это с такими изменениями, как typedef vector<complex<double>> state_type:

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

#include <boost/numeric/odeint.hpp>

using namespace std;
using namespace boost::numeric::odeint;typedef vector<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[0] = x[1];
dxdt[1] = (1.0+m_eta*I)*x[0]-(1.0+m_alpha*I)*x[1];
}
};

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 << t;
m_out << "\t" << x[0].real() << "\t" << x[0].imag() ;
m_out << "\n";*/
}
};

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

const double dt = 0.1;

typedef runge_kutta4< state_type , double , state_type , double ,
vector_space_algebra > stepper_type;

integrate_const( stepper_type() , stuart_landau( 2.0 , 1.0 ) , x , 0.0 , 10.0 , dt , streaming_observer( cout ) );
//]

return 0;

}

Однако это дает множество ошибок.

Кто-нибудь может дать мне минимальный рабочий пример со сложными ODE? и даже лучше реализовать высокоточные типы данных. Похоже, что odeint может поддерживать совершенно произвольные типы состояний, но у меня много проблем с тем, чтобы заставить их работать.

Обновленный код

#include <iostream>
#include <complex>
#include <boost/array.hpp>
#include <stdlib.h>
#include <math.h>

#include <boost/numeric/odeint.hpp>extern "C" {
#include <quadmath.h>
}

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

//[ stuart_landau_system_function
typedef std::vector<std::complex < __float128 > > state_type;

struct stuart_landau
{
__float128 m_eta;
__float128 m_alpha;

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

void operator()( const state_type &x , state_type &dxdt , double t ) const
{
const complex< __float128 > I( 0.0 , 1.0 ); //define complex I

dxdt[0] = x[1];
dxdt[1] = ( 1.0 + m_eta * I ) * x[0] - ( 1.0 + m_alpha * I )*x[1];
}
};
//]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 << t;
/*  m_out << "\t" << x[0].real() << "\t" << x[0].imag() ;
m_out << "\n";*/
}
};int main( int argc , char **argv )
{
//[ stuart_landau_integration
state_type x(2);
x[0] = complex< __float128 >( 1.0L , 0.0L );
x[1] = complex< __float128 >( 1.0L , 0.0L );

const double dt = 0.1;

typedef runge_kutta4< state_type, __float128 > stepper_type;

integrate_const( stepper_type() , stuart_landau( 2.0L , 1.0L ) , x , 0.0 , 10.0 , dt
/*, streaming_observer( cout ) */);
//]

return 0;
}

ОШИБКА ВХОДА НА КОМПИЛЬЕ:

vecPrec.cpp: In member function ‘void stuart_landau::operator()(const state_type&,
state_type&, double) const’:
vecPrec.cpp:33:35: error: no match for ‘operator+’ in ‘1.0e+0 + std::operator* [with
_Tp = __float128]((* &((const stuart_landau*)this)->stuart_landau::m_eta), (* & I))’
vecPrec.cpp:33:35: note: candidates are:
/usr/include/c++/4.6/bits/stl_iterator.h:327:5: note: template<class _Iterator>    std::reverse_iterator<_Iterator> std::operator+(typename

std::reverse_iterator<_Iterator>::difference_type, const   std::reverse_iterator<_Iterator>&)
/usr/include/c++/4.6/bits/basic_string.h:2306:5: note: template<class _CharT, class  _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(const    std::basic_string<_CharT, _Traits, _Alloc>&, const std::basic_string<_CharT, _Traits,   _Alloc>&)
/usr/include/c++/4.6/bits/basic_string.tcc:694:5: note: template<class _CharT, class  _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(const   _CharT*, const std::basic_string<_CharT, _Traits, _Alloc>&)
/usr/include/c++/4.6/bits/basic_string.tcc:710:5: note: template<class _CharT, class _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(_CharT,     const std::basic_string<_CharT, _Traits, _Alloc>&)
/usr/include/c++/4.6/bits/basic_string.h:2343:5: note: template<class _CharT, class  _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(const   std::basic_string<_CharT, _Traits, _Alloc>&, const _CharT*)
/usr/include/c++/4.6/bits/basic_string.h:2359:5: note: template<class _CharT, class   _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(const   std::basic_string<_CharT, _Traits, _Alloc>&, _CharT)
/usr/include/c++/4.6/complex:321:5: note: template<class _Tp> std::complex<_Tp>  std::operator+(const std::complex<_Tp>&, const std::complex<_Tp>&)
/usr/include/c++/4.6/complex:330:5: note: template<class _Tp> std::complex<_Tp>  std::operator+(const std::complex<_Tp>&, const _Tp&)
/usr/include/c++/4.6/complex:339:5: note: template<class _Tp> std::complex<_Tp>  std::operator+(const _Tp&, const std::complex<_Tp>&)
/usr/include/c++/4.6/complex:440:5: note: template<class _Tp> std::complex<_Tp>  std::operator+(const std::complex<_Tp>&)
/usr/include/c++/4.6/bits/stl_bvector.h:266:3: note: std::_Bit_iterator   std::operator+(std::ptrdiff_t, const std::_Bit_iterator&)
/usr/include/c++/4.6/bits/stl_bvector.h:266:3: note:   no known conversion for argument 2 from ‘std::complex<__float128>’ to ‘const std::_Bit_iterator&’
/usr/include/c++/4.6/bits/stl_bvector.h:352:3: note: std::_Bit_const_iterator std::operator+(std::ptrdiff_t, const std::_Bit_const_iterator&)
/usr/include/c++/4.6/bits/stl_bvector.h:352:3: note:   no known conversion for   argument 2 from ‘std::complex<__float128>’ to ‘const std::_Bit_const_iterator&’
vecPrec.cpp:33:66: error: no match for ‘operator+’ in ‘1.0e+0 + std::operator*  [with _Tp = __float128]((* &((const stuart_landau*)this)->stuart_landau::m_alpha), (* &   I))’
vecPrec.cpp:33:66: note: candidates are:
/usr/include/c++/4.6/bits/stl_iterator.h:327:5: note: template<class _Iterator>   std::reverse_iterator<_Iterator> std::operator+(typename    std::reverse_iterator<_Iterator>::difference_type, const std::reverse_iterator<_Iterator>&)
/usr/include/c++/4.6/bits/basic_string.h:2306:5: note: template<class _CharT, class     _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+   (const std::basic_string<_CharT, _Traits, _Alloc>&, const   std::basic_string<_CharT, _Traits, _Alloc>&)
/usr/include/c++/4.6/bits/basic_string.tcc:694:5: note: template<class _CharT, class   _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+  (const _CharT*, const std::basic_string<_CharT, _Traits, _Alloc>&)
/usr/include/c++/4.6/bits/basic_string.tcc:710:5: note: template<class _CharT, class _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+   (_CharT, const std::basic_string<_CharT, _Traits, _Alloc>&)
/usr/include/c++/4.6/bits/basic_string.h:2343:5: note: template<class _CharT, class _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(const std::basic_string<_CharT, _Traits, _Alloc>&, const _CharT*)
/usr/include/c++/4.6/bits/basic_string.h:2359:5: note: template<class _CharT, class _Traits, class _Alloc> std::basic_string<_CharT, _Traits, _Alloc> std::operator+(const std::basic_string<_CharT, _Traits, _Alloc>&, _CharT)
/usr/include/c++/4.6/complex:321:5: note: template<class _Tp> std::complex<_Tp> std::operator+(const std::complex<_Tp>&, const std::complex<_Tp>&)
/usr/include/c++/4.6/complex:330:5: note: template<class _Tp> std::complex<_Tp> std::operator+(const std::complex<_Tp>&, const _Tp&)
/usr/include/c++/4.6/complex:339:5: note: template<class _Tp> std::complex<_Tp> std::operator+(const _Tp&, const std::complex<_Tp>&)
/usr/include/c++/4.6/complex:440:5: note: template<class _Tp> std::complex<_Tp> std::operator+(const std::complex<_Tp>&)
/usr/include/c++/4.6/bits/stl_bvector.h:266:3: note: std::_Bit_iterator std::operator+(std::ptrdiff_t, const std::_Bit_iterator&)
/usr/include/c++/4.6/bits/stl_bvector.h:266:3: note:   no known conversion for argument 2 from ‘std::complex<__float128>’ to ‘const std::_Bit_iterator&’
/usr/include/c++/4.6/bits/stl_bvector.h:352:3: note: std::_Bit_const_iterator std::operator+(std::ptrdiff_t, const std::_Bit_const_iterator&)
/usr/include/c++/4.6/bits/stl_bvector.h:352:3: note:   no known conversion for argument 2 from ‘std::complex<__float128>’ to ‘const std::_Bit_const_iterator&’