Armadillo 与 Boost Odeint 冲突：Odeint 在积分期间将状态向量的大小调整为零

Question

我刚开始使用Boost Odeint 来集成一个ODE 系统。为了方便起见，我想将它与 Armadillo 一起使用，因为它们都是现代 C++ 库，具有方便的 API。但是，如果我将 arma::vec 指定为状态类型，结果立即表明在积分的第一步中，integrate_adaptive() 将状态向量的大小调整为 0x1。我 post 这里有一个简单的例子：

#include <iostream>
#include <armadillo>

#include <boost/numeric/odeint.hpp>

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

typedef vec state_type;

class harm_osc
{
    private:
    mat my_A;

    public:
        harm_osc(double gam)
        {
            my_A = { {0.0, 1.0}, {-gam*gam, 0.0} };
        }

        harm_osc()
        {
            my_A = { {0.0, 1.0}, {-1.0, 0.0} };
        }

        void operator() (const vec& x, vec& dxdt, const double t)
        {
            cout << "size of x: " << size(x) << endl;
            cout << "size of dxdt: " << size(dxdt) << endl;
            dxdt = my_A*x;
        }
};

class observer
{
    private:
        mat& my_states;
        vec& my_times;                        

    public:
        observer(mat& states, vec& times):
            my_states(states),
            my_times(times)
        {}    

        void operator() (const vec& x, double t)
        {
            my_states.insert_rows(my_states.n_rows, x);
            my_times.insert_rows(my_times.n_rows, t);
            cout << t << '\t';
            for(auto elem : x)
                cout << elem << '\t';
            cout << endl;
        }
};

typedef runge_kutta_cash_karp54<state_type> error_stepper_type;
typedef controlled_runge_kutta<error_stepper_type> controlled_stepper_type;

int main()
{
    state_type x = {0.0, 1.0};

    vec t;
    mat x_full;

    integrate_adaptive(make_controlled<error_stepper_type>(1e-5, 1e-5), harm_osc(1.0), x, 0.0, 200.0, 0.01, observer(x_full, t));
}

如果我指定 arma::vec::fixed<2> 而不是 arma::vec 作为 state_type，这个简单的演示会正确运行。我的问题是，在我正在处理的当前项目中，我在编译时不知道状态向量的大小，因此我无法像之前提到的那样使用模板参数修复它。

有没有在编译时不固定状态向量大小的情况下将 Armadillo 与 Boost Odeint 一起使用的解决方案？在我项目的其他部分，我可以正常使用 Armadillo。我想在我的整个项目中使用它。现在，在集成 ODE 时，我使用 Boost Ublas 来定义 ODE 系统。

Answer 1

Odeint 在内部存储了几个中间状态。在您的情况下，我认为它不知道如何正确调整中间状态的大小。这可以通过引入正确的调整大小适配器轻松解决：

namespace boost { namespace numeric { namespace odeint {

template <>
struct is_resizeable<arma::vec>
{
    typedef boost::true_type type;
    const static bool value = type::value;
};

template <>
struct same_size_impl<arma::vec, arma::vec>
{
    static bool same_size(const arma::vec& x, const arma::vec& y)
    {
        return x.size() == y.size();   // or use .n_elem attributes
    }
};

template<>
struct resize_impl<arma::vec, arma::vec>
{
    static void resize(arma::vec &v1, const arma::vec& v2)
    {
        v1.resize(v2.size());     // not sure if this is correct for arma
    }
};

} } } // namespace boost::numeric::odeint

看这里：http://www.boost.org/doc/libs/1_63_0/libs/numeric/odeint/doc/html/boost_numeric_odeint/odeint_in_detail/state_types__algebras_and_operations.html

Answer 2

对标准示例 (https://github.com/headmyshoulder/odeint-v2/blob/master/examples/my_vector.cpp) 稍作调整即可获得有效代码：

#include <boost/numeric/odeint.hpp>
#include <armadillo>
#include <iostream>

namespace boost { 
namespace numeric { 
namespace odeint {

template<>
struct is_resizeable< arma::vec >
{
 typedef boost::true_type type;
 static const bool value = type::value;
};

} } }

typedef arma::vec state_type;

void lorenz( const state_type &x , state_type &dxdt , const double t )
{
const double sigma( 10.0 );
const double R( 28.0 );
const double b( 8.0 / 3.0 );

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];
}

using namespace boost::numeric::odeint;

int main()
{
state_type x(3);
x[0] = 5.0 ; x[1] = 10.0 ; x[2] = 10.0;

// make sure resizing is ON
BOOST_STATIC_ASSERT( is_resizeable<state_type>::value == true );
// no further work is required
std::cout << "initial x : " << x << std::endl; 
integrate_const( runge_kutta4< state_type >() , lorenz , x , 
0.0 , 10.0 , 0.1 );
std::cout << "final x : " << x << std::endl;
return 0; 
}

Armadillo 与 Boost Odeint 冲突：Odeint 在积分期间将状态向量的大小调整为零

Armadillo conflicts with Boost Odeint: Odeint resizes the state vector to zero during integration

c++

boost

armadillo

c++11