gsl_odeiv2在c + +類：模板包裝int（...）ode

Question

我目前在我的類中使用gsl_odeiv2方法來解微分方程。但是由於衆所周知的成員函數問題，我無法在類中定義我的ode-system。我目前使用的一種變通方法： 我在一個全局命名空間定義我的頌歌：

ODE.hpp: 
#include "EoS.hpp" 

#include <gsl/gsl_math.h> 
#include <gsl/gsl_errno.h> 

namespace ODEs 
{ 
    struct tov_eq_params {EoS *eos;}; 
    int tov_eq(double, const double *, double *, void *); 
} 

ODE.cpp: 
namespace ODEs { 
    int tov_eq(double r, const double *PM, double *dPdM, void *p) { 
     struct tov_eq_params * params = (struct tov_eq_params *)p; 
     EoS *eos = (params->eos); 
     ... 
    return GSL_SUCCESS 
    } 
} 

，並使用指針coustom類型（類EOS）作爲參數的對象。內部類解決我的頌歌我使用：

... 
struct tov_eq_params comp_tov_params = {(this->star_eos)}; 
gsl_odeiv2_evolve *comp_tov_evolve = gsl_odeiv2_evolve_alloc(3); 
gsl_odeiv2_system comp_tov_system = {tov_eq, NULL, 3,&comp_tov_params}; 
... 

initalise我的系統。這工作正常，但有點混亂，因爲我需要在全局命名空間中聲明我的微分方程。

我知道有可能使用gsl_functions的模板包裝stackoverflow.com/questions/.../how-to-avoid-static-member-function-when-using-gsl-with-c/...在C++類中使用它們。我實際上使用了那裏描述的包裝來爲我的類中的gsl_integration方法定義函數，它完美地工作，並且更簡潔，代碼更少。例如：我可以用我的star_eos從上面direcly對象在函數內部：

auto dBf = [=](double r)->double{ 
     return 4 * M_PI * gsl_pow_2(r) * (this->star_eos)->nbar(this->P(r)) * sqrt(this->expLambda(r))* 1e54; 
    }; 
    gsl_function_pp<decltype(dBf)> dBfp(dBf); 
    gsl_function *dB = static_cast<gsl_function*>(&dBfp); 

我試着寫這樣的INT（雙R，常量雙* PM，模板包裝雙* DPDM，無效* P ）函數gsl_odeiv2_system需要，但我失敗了，因爲我是C++的新手，並沒有完全理解它的模板/ static_cast機制。

是否有人將gsl_odeiv方法及其ode系統與模板包裝器一起使用？或者可以有人想出一個模板類似於上面描述的gsl_functions但是int（...）ode。

Answer 1

思考我如何在全局命名空間中設置微分方程時使用它，我找到了解決方案。我現在有一個工作包裝。在一個全局命名空間，我有以下：

//gsl_wrapper.hpp 
#include <iostream> 
#include <vector> 
#include <functional> 

#include <gsl/gsl_math.h> 
#include <gsl/gsl_errno.h> 

namespace gsl_wrapper { 

    class ode_System{ 
    public: 
     ode_System(int); 
     int dim; 
     std::function<double (double, const double *, double *, int)> *df; 

    }; 

    struct ode_struct {ode_System *ode;}; 
    int ode(double, const double *, double *, void *); 
} 

//gsl_wrapper.cpp 
#include "../../include/gsl_wrapper.hpp" 

namespace gsl_wrapper { 

    ode_System::ode_System(int dim) { 
     this->dim=dim; 
    } 

    int ode(double r, const double *f, double *df, void *p) { 
     struct ode_struct * params = (struct ode_struct *)p; 
     ode_System *odeFunc = (params->ode); 

     int dim = odeFunc->dim; 
     std::function<double (double, const double *, double *, int)> dfeq=*(odeFunc->df); 

     for(int i=0;i<dim;i++){ 
      df[i] = dfeq(r,f,df,i); 
     } 

     return GSL_SUCCESS; 
    } 

}; 

所以我bassically有存儲在我的新類ode_System我所有的具體信息，其中有一個INT昏暗到指定的系統尺寸和指針，以便一std :: function對象。這個對象代表了數學微分方程組。

裏面我的類，其中，我要解決使用gsl_odeiv2一微分方程，我定義系統使用lambda函數：

std::function<double (double, const double *, double *, int)> dPMeq = [=](double r , const double * PM, double *dPdM, int i)->double{ 
    double df; 
    switch (i) 
    { 
     case 0: 
      df = ... // df_1/dr 
      break; 
     case 1: 
      df = ... // df_2/dr 
      break; 
     case 2: 
      df = ... // df_3/dr 
      break; 
     default: 
      GSL_ERROR_VAL ("comp_tov_eq:", GSL_FAILURE,GSL_FAILURE); 
      df = 0; 
    } 
    return df; 
}; 

上述系統代表3次微分方程的耦合系統。然後我聲明一個ode_System對象具有正確的尺寸，並將其功能指針df設置爲我定義的系統。然後，我只需要一個結構具有參考該系統完成：我可以用我的類中定義我微分方程gsl_odeiv2_system：

ode_System tov(3); 
tov.df= &dPMeq; 
struct ode_struct comp_tov_params = {&tov}; 
gsl_odeiv2_evolve *comp_tov_evolve = gsl_odeiv2_evolve_alloc(3); 
gsl_odeiv2_system comp_tov_system = {ode, NULL, 3, &comp_tov_params}; 
... 

至於我可以告訴這個作品一樣好（或壞）作爲我在我的問題中提出的實現。它可以使用一些清理，但原則上這對我來說工作得很好。

但是，如果有人知道更好的方式來做到這一點，請隨時分享！