c++ - 在boost/odeint中使用几个特征矩阵作为状态类型

我正在研究一个物理问题，必须根据ODE来开发参数。有时需要对其进行操纵，以便让我希望有一种可以与诸如对角线化之类的例程一起使用的数据类型。...因此，我实现了一个以eigen::Matrix作为成员的类，并希望执行集成与odeint。对于单个特征::矩阵，这很好用。我举了一个最小的例子:

#include <iostream>
#include <fstream>
#include <cmath>
#include <string>

#include <Eigen/Core>
#include <boost/numeric/odeint.hpp>
#include <boost/numeric/odeint/external/eigen/eigen_algebra.hpp>

// define vector_space_algebra for Eigen::Matrix
namespace boost::numeric::odeint {
  template<typename B,int S1,int S2,int O, int M1, int M2>
  struct algebra_dispatcher< Eigen::Matrix<B,S1,S2,O,M1,M2> >{
    typedef vector_space_algebra algebra_type;
  };
}

// define abs() for Eigen::Matrix
namespace Eigen {
  template<typename D, int Rows, int Cols>
  Matrix<D, Rows, Cols> abs(Matrix<D, Rows, Cols> const& m) {
    return m.cwiseAbs();
  }
}

typedef Eigen::Matrix<double, 3,3> mat;
using namespace Eigen;
using namespace std;

class state {
 public:

  // state components
  Eigen::Matrix<double, 3,3> M1, M2;

  // constructors
  state() : M1(), M2() {};      // constructors
  state(mat M1in, mat M2in) : M1(M1in), M2(M2in) {};


  // in place addition and multiplication
  state operator+=(const state & X){
    M1 += X.M1; M2 += X.M2;
    return *this;
  }

  state operator*=(const double a){
    M1 *= a; M2 *= a;
    return *this;
  }

  // ODE
  void operator() ( const state & X , state & dX, const double ){
    dX.M1 = X.M1*X.M2.adjoint()*X.M2;
    dX.M2 = X.M2*X.M1.adjoint()*X.M1;
  }
};


// vector space operations

state operator+( const state &lhs , const state &rhs ){
  return state( lhs.M1+rhs.M1 ,lhs.M2+rhs.M2);
}

state operator*( const state &lhs , const double &rhs ){
  return state( lhs.M1*rhs ,lhs.M2*rhs);
}

state operator*( const double &lhs , const state &rhs ){
  return state( lhs*rhs.M1 ,lhs*rhs.M2);
}

state operator/( const state &lhs , const state &rhs ){
  return state( lhs.M1.cwiseQuotient(rhs.M1), lhs.M2.cwiseQuotient(rhs.M2) );
}
state abs( const state &X ){
  return state( abs(X.M1) , abs(X.M2) );
}


// lp infinity norm
namespace boost::numeric::odeint {
  template<>
  struct vector_space_norm_inf< state > {
    typedef double result_type;
    double operator()( const state &X ) const {
      return max( X.M1.lpNorm<Infinity>() , X.M2.lpNorm<Infinity>() );
    }
  };
}


//write to std output
void write( state &x , const double t ){
  cout << t << "\t" << x.M1 << "\t" << x.M2 << "\n";
}


//
// int main
//
int main(int argc, char* argv[]){

  // set values
  mat M1, M2;
  double t_end = 1;
  double t_start = 10;

  M1 << 0.1,0,0, 0,0.2,0.1 ,0.2,0,0.3;
  M2 << 0.5,0,0, 0,0.6,0, 0,0,0.7;

  state values(M1,M2);

  using namespace boost::numeric::odeint;

  // type definition for numerical integration
  typedef runge_kutta_dopri5< state , double, state , double,     vector_space_algebra > stepper;

  // integration
  int steps = integrate_adaptive( make_controlled<stepper>( 1E-10 , 1E-10     ) , state() , values , t_start , t_end , 0.01);

  //output
  write(values,t_end);
  return(0);
}

基本上，这是来自here的示例
当我注释掉以“int steps”开头的行时，mac上的g++(我知道，这有所不同，这些错误的可读性更高……)编译时没有错误。不然我会

In file included from minimal.cpp:7:
In file included from /usr/local/include/boost/numeric/odeint.hpp:25:
In file included from /usr/local/include/boost/numeric/odeint/util/ublas_wrapper.hpp:30:
/usr/local/include/boost/numeric/odeint/algebra/default_operations.hpp:443:76: error:
      invalid operands to binary expression ('double' and 'state')
  ...m_eps_abs + m_eps_rel * ( m_a_x * abs( get_unit_value( t1 ) ) + m_a_dxdt * abs( get_unit_value( t2 ) ) )...

由于它未在该文件中声明，因此我不了解get_unit_value()函数是做什么的或要我做什么。似乎与误差估计有关，或者至少与执行集成中的某个步骤有关。我该如何解决？

最佳答案

问题是您没有在operator+和double之间定义state。添加以下内容，并且代码至少应该编译。

state operator+( const state &lhs , double rhs ){
  return state( lhs.M1+rhs ,lhs.M2+rhs);
}

state operator+( double lhs , const state &rhs ){
  return state( lhs+rhs.M1 ,lhs+rhs.M2);
}

关于c++ - 在boost/odeint中使用几个特征矩阵作为状态类型，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/39616247/