我在使用受控错误步进器和复杂状态类型的 odeint 库时遇到了问题。我对具有复杂斯图尔特朗道方程的示例中的代码进行了修改,使其包含自适应积分器。代码现在看起来像这样:
#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;
}
但是,如果我将步进器的定义更改为
bulirsch_stoer< state_type, complex< double > > stepper( 1E-12 , 1E-12 , 1 , 1 );
编译失败。 我的问题是:受控错误步进器不支持复杂的数据类型吗?如果是这样,是否有一种方法可以避免出现的问题。或者,是否可以为复杂数据类型定义自己的 vector 代数?
请您参考如下方法:
odeint 库确实支持复杂的数据类型,包括受控步进器。
例如,要在 4 阶/5 阶 Dormand-Prince 方法的情况下获得受控步进器,您可以执行以下操作:
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 ) );
其中 state_type
的类型定义为
typedef complex< double > state_type;
Bulirsch-Stoer algorithm然而,在 odeint 中实现的它本质上已经是一个受控步进器。因此,在这种情况下不可能应用 make_control 或其他策略来创建受控步进器。 odeint 库中没有这样的模板,因为它没有意义。
您发布的第一个代码已经集成了 ODE、复数值和受控步进器。目前尚不清楚您想通过将步进器更改为来实现什么目的
bulirsch_stoer< complex<double>, complex<double> > stepper( 1E-12 , 1E-12 , 1 , 1 );