C++ class GSL ODE 求解器的成员函数

Question

我正在努力解决 - 显然 - 一个众所周知的问题。我有一个由 class 成员函数定义的 ODE 系统，我想通过其中一个 GSL 求解器 solve/integrate 它。即，假设我的模型定义为：

class my_model{
...
public:
    int NEQ = 4;
    double y[4], dydt[4];
    double params[25]; 
    int ode_gsl(double t, const double y[], double dydt[], void * params);
    ...
};

int my_model::int ode_gsl(double t, const double y[], double dydt[], void * params){
    dydt[0] = params[1]*params[0]*y[0] - y[1];
    dydt[1] = ...;
    ...
    return GSL_SUCCESS;    
}

然后在我的集成例程中：

int main(){
    my_model chemosc;

   // Parameters for the Integrator
   double HSTART = 1e-3;
   double ATOL = 1e-8;
   double RTOL = 1e-6;

   // Instantiate GSL solver
   gsl_odeiv2_system sys = {&chemosc.ode_gsl, nullptr, chemosc.NEQ, chemosc.params};
   gsl_odeiv2_driver * d;
   d = gsl_odeiv2_driver_alloc_y_new (&sys, gsl_odeiv2_step_rk8pd, HSTART, ATOL, RTOL);

   double twin[2] = {0.,60.};
   double dt = 1e-3;
   double t = twin[0], t1 = twin[1];
   long int NSTEP = (long int)((t1-t)/dt)+1; // +1 if you start counting from zero...
   int NEQ = 4;
   long int NUMEL = (NEQ+1)*NSTEP; // number of elements for solution

   int i = 0,j;
   do{
      double ti = t + dt; // t is automatically updated by the driver
      printf("\n%.3f\t%.3f\t%.3f t%.3f",astro.y[0],astro.y[1],astro.y[2],astro.y[3]);
      int status = gsl_odeiv2_driver_apply (d, &t, ti, astro.y);
      ...
      }
...
}

编译上述代码会产生一个众所周知的错误，即 GSL 需要指向函数的指针，而我传递的是指向成员函数的指针，即：

error: cannot convert ‘int (chemosc::*)(double, const double*, double*, void*)’ to ‘int (*)(double, const double*, double*, void*)’

我发现了几个与这个主题相关的问题：Q1, Q2, Q3, Q4, Q5, Q6, but seriously, the answer there are hard to follow. Declaring as static my member function, has the downside that the compiler requires me to also declare as static all the member parameters. Using static_cast as suggested here，导致了所有的分割错误（但我假设我在实现中做错了什么，因为 post 非常神秘）。对于这个问题（可能不使用 boost 库）有一个一劳永逸的清晰且尽可能简单的工作解决方案会很好。

Answer 1

像这样：

class my_model
{
    private: int
    ode_gsl_impl(double const t, double const * const y, double * const dydt);

    public: static int
    ode_gsl(double const t, double const * const y, double * const dydt, void * const opaque)
    {
        assert(opaque);
        return(static_cast<my_model *>(opaque)->ode_gsl_impl(t, y, dydt));
    }
};

gsl_odeiv2_system sys
{
    &my_model::ode_gsl
,   nullptr
,   chemosc.NEQ
,   reinterpret_cast<void *>(::std::addressof(chemosc))
};

我想提一下，您的原始代码在回调参数名称和 class 字段名称之间存在名称冲突。

C++ class GSL ODE 求解器的成员函数

C++ class member function to GSL ODE solver

c++

pointers

class

ode

gsl