Armadillo 与 Boost Odeint 冲突:Odeint 在积分期间将状态向量的大小调整为零
Posted
技术标签:
【中文标题】Armadillo 与 Boost Odeint 冲突:Odeint 在积分期间将状态向量的大小调整为零【英文标题】:Armadillo conflicts with Boost Odeint: Odeint resizes the state vector to zero during integration 【发布时间】:2017-01-19 23:53:43 【问题描述】:我刚开始使用 Boost Odeint 来集成一个 ODE 系统。为方便起见,我想将它与 Armadillo 一起使用,因为它们都是具有便捷 API 的现代 C++ 库。但是,如果我将arma::vec
指定为状态类型,则立即证明在集成的第一步中,integrate_adaptive()
将状态向量的大小调整为0x1
。我在这里发布一个简单的例子:
#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 一起使用,而无需在编译时固定状态向量的大小?在我项目的其他部分,我可以正确使用犰狳。我想在我的整个项目中使用它。现在,在集成 ODE 时,我使用 Boost Ublas 来定义 ODE 系统。
【问题讨论】:
【参考方案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
【讨论】:
【参考方案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 在积分期间将状态向量的大小调整为零的主要内容,如果未能解决你的问题,请参考以下文章
boost odeint 给出了与 Python3 scipy 非常不同的值
使用 Boost 反序列化 Armadillo colvec
如何将向量传递给基于推力的 odeint 观察者的构造函数,以便可以在函子中读取它