使用 Boost::odeint 和 Eigen::Matrix 作为状态向量
我正在尝试利用 ODE 集成功能Boost 使用 Eigen 3 中的矩阵类 作为我的状态向量，但我遇到了问题深入到我不知道如何解决的 Boost.
I'm trying to utilize the ODE integration capabilities of Boost using the Matrix class from Eigen 3 as my state vector, but I'm running into problems deep into Boost that I don't understand how to address.
我正在尝试做的一个最小的例子:
A minimal example of what I'm trying to do:
#include <Eigen/Core>
#include <boost/numeric/odeint/stepper/runge_kutta_dopri5.hpp>
#include <iostream>
using namespace Eigen;
using namespace boost::numeric::odeint;
template<size_t N>
using vector = Matrix<double, N, 1>;
typedef vector<3> state;
int main() {
state X0;
X0 << 1., 2., 3.;
state xout = X0;
runge_kutta_dopri5<state> stepper;
// If I remove these lines, everything compiles fine
stepper.do_step([](const state x, state dxdt, const double t) -> void {
dxdt = x;
}, X0, 0.0, xout, 0.01);
std::cout << xout << std::endl;
}
如果我注释掉对 stepper.do_step 的调用，一切都会编译并运行得很好，但当然不会做任何有趣的事情.如果我不这样做，Boost 在我的终端上吐出编译错误，其中第一个是
If I coment out the call to stepper.do_step everything compiles and runs just fine, but of course doesn't do anything interesting. If I don't, Boost vomits compile errors over my terminal, the first of which is
In file included from /usr/include/boost/mpl/aux_/begin_end_impl.hpp:20:0,
from /usr/include/boost/mpl/begin_end.hpp:18,
from /usr/include/boost/mpl/is_sequence.hpp:19,
from /usr/include/boost/fusion/support/detail/is_mpl_sequence.hpp:12,
from /usr/include/boost/fusion/support/tag_of.hpp:13,
from /usr/include/boost/fusion/support/is_sequence.hpp:11,
from /usr/include/boost/fusion/sequence/intrinsic_fwd.hpp:12,
from /usr/include/boost/fusion/sequence/intrinsic/front.hpp:10,
from /usr/include/boost/fusion/include/front.hpp:10,
from /usr/include/boost/numeric/odeint/util/is_resizeable.hpp:26,
from /usr/include/boost/numeric/odeint/util/state_wrapper.hpp:25,
from /usr/include/boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp:27,
from /usr/include/boost/numeric/odeint/stepper/runge_kutta_dopri5.hpp:24,
from /home/tlycken/exjobb/Code/alpha-orbit-follower/test/algebra/algebra-tests.cpp:2:
/usr/include/boost/mpl/eval_if.hpp: In instantiation of ‘struct boost::mpl::eval_if_c<true, boost::range_const_iterator<Eigen::Matrix<double, 3, 1> >, boost::range_mutable_iterator<const Eigen::Matrix<double, 3, 1> > >’:
/usr/include/boost/range/iterator.hpp:63:63: required from ‘struct boost::range_iterator<const Eigen::Matrix<double, 3, 1> >’
/usr/include/boost/range/begin.hpp:112:61: required by substitution of ‘template<class T> typename boost::range_iterator<const T>::type boost::range_adl_barrier::begin(const T&) [with T = Eigen::Matrix<double, 3, 1>]’
/usr/include/boost/numeric/odeint/algebra/range_algebra.hpp:52:45: required from ‘static void boost::numeric::odeint::range_algebra::for_each3(S1&, S2&, S3&, Op) [with S1 = Eigen::Matrix<double, 3, 1>; S2 = const Eigen::Matrix<double, 3, 1>; S3 = const Eigen::Matrix<double, 3, 1>; Op = boost::numeric::odeint::default_operations::scale_sum2<double, double>]’
/usr/include/boost/numeric/odeint/stepper/runge_kutta_dopri5.hpp:128:9: required from ‘void boost::numeric::odeint::runge_kutta_dopri5<State, Value, Deriv, Time, Algebra, Operations, Resizer>::do_step_impl(System, const StateIn&, const DerivIn&, boost::numeric::odeint::runge_kutta_dopri5<State, Value, Deriv, Time, Algebra, Operations, Resizer>::time_type, StateOut&, DerivOut&, boost::numeric::odeint::runge_kutta_dopri5<State, Value, Deriv, Time, Algebra, Operations, Resizer>::time_type) [with System = main()::__lambda0; StateIn = Eigen::Matrix<double, 3, 1>; DerivIn = Eigen::Matrix<double, 3, 1>; StateOut = Eigen::Matrix<double, 3, 1>; DerivOut = Eigen::Matrix<double, 3, 1>; State = Eigen::Matrix<double, 3, 1>; Value = double; Deriv = Eigen::Matrix<double, 3, 1>; Time = double; Algebra = boost::numeric::odeint::range_algebra; Operations = boost::numeric::odeint::default_operations; Resizer = boost::numeric::odeint::initially_resizer; boost::numeric::odeint::runge_kutta_dopri5<State, Value, Deriv, Time, Algebra, Operations, Resizer>::time_type = double]’
/usr/include/boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp:167:9: required from ‘typename boost::disable_if<boost::is_same<StateInOut, Time>, void>::type boost::numeric::odeint::explicit_error_stepper_fsal_base<Stepper, Order, StepperOrder, ErrorOrder, State, Value, Deriv, Time, Algebra, Operations, Resizer>::do_step(System, const StateIn&, boost::numeric::odeint::explicit_error_stepper_fsal_base<Stepper, Order, StepperOrder, ErrorOrder, State, Value, Deriv, Time, Algebra, Operations, Resizer>::time_type, StateOut&, boost::numeric::odeint::explicit_error_stepper_fsal_base<Stepper, Order, StepperOrder, ErrorOrder, State, Value, Deriv, Time, Algebra, Operations, Resizer>::time_type) [with System = main()::__lambda0; StateIn = Eigen::Matrix<double, 3, 1>; StateOut = Eigen::Matrix<double, 3, 1>; Stepper = boost::numeric::odeint::runge_kutta_dopri5<Eigen::Matrix<double, 3, 1> >; short unsigned int Order = 5u; short unsigned int StepperOrder = 5u; short unsigned int ErrorOrder = 4u; State = Eigen::Matrix<double, 3, 1>; Value = double; Deriv = Eigen::Matrix<double, 3, 1>; Time = double; Algebra = boost::numeric::odeint::range_algebra; Operations = boost::numeric::odeint::default_operations; Resizer = boost::numeric::odeint::initially_resizer; typename boost::disable_if<boost::is_same<StateInOut, Time>, void>::type = void; boost::numeric::odeint::explicit_error_stepper_fsal_base<Stepper, Order, StepperOrder, ErrorOrder, State, Value, Deriv, Time, Algebra, Operations, Resizer>::time_type = double]’
/home/tlycken/exjobb/Code/alpha-orbit-follower/test/algebra/algebra-tests.cpp:21:137: required from here
/usr/include/boost/mpl/eval_if.hpp:60:31: error: no type named ‘type’ in ‘boost::mpl::eval_if_c<true, boost::range_const_iterator<Eigen::Matrix<double, 3, 1> >, boost::range_mutable_iterator<const Eigen::Matrix<double, 3, 1> > >::f_ {aka struct boost::range_const_iterator<Eigen::Matrix<double, 3, 1> >}’
typedef typename f_::type type;
我试图深入发生错误的Boost标头，但我对修复我的代码的理解不够.由于 odeint 库文档明确说明
I tried to dig into the Boost header where the error occurs, but I didn't understand enough of what's going on to be able to fix my code. Since the odeint library documentation clearly states
odeint 的主要重点是提供以算法完全独立于用于表示状态 x 的数据结构的方式实现的数值方法.
The main focus of odeint is to provide numerical methods implemented in a way where the algorithm is completely independent on the data structure used to represent the state x.
我相信即使 odeint 本身不支持 Eigen，这也不应该太难工作.
I believe this shouldn't be too hard to get working even if odeint doesn't support Eigen natively.
推荐答案
你只需要用
runge_kutta_dopri5<state,double,state,double,vector_space_algebra> stepper;
Eigen 应该使用 vector_space_algebra 开箱即用，但您需要手动指定它们.在下一个 odeint 版本中，我们有一个自动检测代数的机制.
Eigen should work out of the Box with the vector_space_algebra but you need to specify them manually. In the next odeint version we have a mechanism for automatically detecting the algebra.
顺便说一句.你对 ODE 的定义不正确，你需要一个导数的参考
Btw. youe definition of the ODE is not correct, you need a reference for the derivatve
stepper.do_step([](const state &x, state &dxdt, const double t) -> void {
dxdt = x;
}, X0, 0.0, xout, 0.01);
相关文章