ODEi​​nt:具有任意精度的自适应积分

ODEint: adaptive integration with arbitrary precision

ODEint 是否可以使用具有任意精度算术的自适应 积分例程?例如,我想将带有 integrate_adaptive() 函数的 Boost 多精度库与受控步进器一起使用。 ODEi​​nt 文档提供了对 integrate_const() 使用任意精度算法的示例,但我无法修改它们以使用自适应积分器。

我也尝试过使用迭代器(例如 make_adaptive_time_iterator...),但 运行 遇到了类似的问题。具体来说,这是一个我希望开始工作的简单代码:

#include <iostream>

//[ mp_lorenz_defs
#include <boost/numeric/odeint.hpp>
#include <boost/multiprecision/cpp_dec_float.hpp>

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

typedef boost::multiprecision::cpp_dec_float_50 value_type;
//typedef double value_type;

typedef boost::array< value_type , 3 > state_type;
//]

//[ mp_lorenz_rhs
struct lorenz
{
    void operator()( const state_type &x , state_type &dxdt , value_type t ) const
    {
        const value_type sigma( 10 );
        const value_type R( 28 );
        const value_type b( value_type( 8 ) / value_type( 3 ) );

        dxdt[0] = sigma * ( x[1] - x[0] );
        dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
        dxdt[2] = -b * x[2] + x[0] * x[1];
    }
};
//]

int main( int argc , char **argv )
{
    //[ mp_lorenz_int
    state_type x = {{ value_type( 10.0 ) , value_type( 10.0 ) , value_type( 10.0 ) }};

    auto stepper = make_controlled( 1.0e-16 , 1.0e-16 , runge_kutta_cash_karp54< state_type >() );


    cout.precision(50);
    integrate_adaptive( stepper ,
           lorenz() , x , value_type( 0.0 ) , value_type( 0.1 ) , value_type( value_type( 1.0 ) / value_type( 2000.0 ) ) );
    //]
    cout << x[0] << endl;

    return 0;
}

编译这个returns错误:

Lorenz_mp2.cpp:52:19: error: no matching function for call to 'make_controlled'
        auto stepper = make_controlled( value_type(1.0e-16) , value_type(1.0e-16) , runge_kutta_cash_karp54< state_type >() );

如果我将 value_type 的 typedef 更改为 double 它编译并且 运行 没问题。

odeint 可以使用任意精度的自适应积分器。您的代码几乎是正确的,您只是忘了将 value_type(用于步进器的内部常量,以及 "time" 变量 t)配置为任意精度类型。如果您查看文档 (http://headmyshoulder.github.io/odeint-v2/doc/boost_numeric_odeint/tutorial/using_arbitrary_precision_floating_point_types.html),您会发现这是通过步进器的第二个模板参数完成的。所以 make_controlled 中正确的步进器定义应该是:

runge_kutta_cash_karp54< state_type , value_type >()

有了这个,它应该可以任意精度编译和运行。

独立于此问题,我想补充一点,使用从 double 到 multi_prec 类型的转换可能存在问题。例如,0.1 不能准确表示为双精度数 (http://www.exploringbinary.com/why-0-point-1-does-not-exist-in-floating-point/)。因此,您最终会得到一个不完全为 0.1 的 multi_prec 值。我建议始终从整数转换,例如将 0.1 表示为 value_type(1)/value_type(10).