我只是想在 C++ 中集成一个函数。我一直在尝试使用 gsl，因为我在网上看到了这个推荐。我遵循了 gsl 示例，但收效甚微。

这是我的 C++ 代码:

double inverseE(double z){
   double inverseE = 1.0/(std::sqrt(Om0*std::pow(1.0+z,3.0)+1.0-Om0));
   return inverseE;
}

double comoving_distance(double z){
   gsl_integration_workspace * w
     = gsl_integration_workspace_alloc (1000);

   double result, error;

   gsl_function iE;
   iE.function = &inverseE;

   gsl_integration_qags (&iE, 0, z, 0, 1e-7, 1000,
                           w, &result, &error);

   gsl_integration_workspace_free (w);

   cout << result << endl;

   return 0;
}

为了澄清，Python 中的相同代码(有效)如下所示:

def iE(z):
     return 1/(np.sqrt(Om0*np.power(1+z,3)+1-Om0))

def comoving_distance(z):
     return (c/H0)*quad(iE,0,z)[0]

quad 执行集成的地方(它是一个 scipy 模块)。

我收到两条错误消息:

ISO C++ 禁止采用未限定或带括号的非静态成员函数的地址来形成指向成员函数的指针。说‘&cosmo::inverseE’[-fpermissive]

无法在赋值中将‘double (cosmo::)(double)’转换为‘double ()(double, void)’* cosmo 是包含这两个函数的类的名称。

看来，这应该不是什么难事。关于我哪里出错的建议将不胜感激!

编辑:类

#include <iostream>    // using IO functions
#include <string>      // using string
#include <gsl/gsl_integration.h>
#include <cmath>
using namespace std;



class cosmo {
private:
   double H0;
   double Om0;
   double Ob0;
   double c;
   double pi;
   double G;

public:
   // Constructor with default values for data members
   cosmo(double Hubble0 = 70, double OmegaM0 = 0.3,
           double OmegaB0 = 0.05) {
         H0 = Hubble0;
         Om0 = OmegaM0;
         Ob0 = OmegaB0;
         c = 3e8;
         pi = 3.141592653589793;
         G = 6.67408e-11;

}

double inverseE(double z){
   double inverseE = 1.0/(std::sqrt(Om0*std::pow(1.0+z,3.0)+1.0-Om0));
   return inverseE;
}

double comoving_distance(double z){
   gsl_integration_workspace * w
     = gsl_integration_workspace_alloc (1000);

   double result, error;

   gsl_function iE;
   iE.function = &inverseE;

   gsl_integration_qags (&iE, 0, z, 0, 1e-7, 1000,
                           w, &result, &error);

   gsl_integration_workspace_free (w);

   cout << result << endl;

   return 0;
}


};

最佳答案

我看到两个问题:

1/ GSL 期望您的 inverseE 函数具有以下原型(prototype)

double inverseE(double z,void *use_data);

在您的代码中您已声明:

double inverseE(double z);

2/ 就像您的代码是用 C++ 编写的，而 GSL 是一个 C 库，您必须使您的 C++ 函数可从 C 调用。

解决方案是声明您的inverseE 函数如下:

extern "C" {
  double inverseE(double z,void *) {
      double inverseE = 1.0/(std::sqrt(Om0*std::pow(1.0+z,3.0)+1.0-Om0));
      return inverseE;
   }
}

extern "C" 使您的 C++ 函数二进制兼容 C 调用约定。

通过这两个修改，我认为您的代码应该没问题。

更新:2/

在我的回答中，我在 2/ 中认为 inverseE 是一个函数。在这里，我考虑它是您类(class)的方法的情况。

这是一个 void *user_data 来拯救的例子:

声明这个包装器:

  extern "C" {
  double YourClass_f_wrap(double z, void *user_data)
  {
    YourClass *this_ptr = (YourClass *)user_data;
    return this_ptr->f(z);
  }
  }

然后YourClass定义如下:

class YourClass
{
 public:
  struct IntegrationResult
  {
    double result, error;
    size_t n_intervals;
  };

 public:
  double f(double x) const; // defines f to integrate
  IntegrationResult integrate_f() const; // integrates f using GSL lib

  ...
};

正如您在评论中提到的，必须谨慎对待前向声明。为了清楚起见，请在下面找到一个完整的运行示例，该示例重现了 the GSL official doc 的结果。但是使用带有方法 f()

完整的运行代码:

可以编译:

g++ gslIntegrationExample.cpp -lgsl -lcblas -o gslIntegrationExample

代码:

#include <cassert>
#include <cstdio>
#include <gsl/gsl_integration.h>

namespace details
{
  extern "C" {
  double YourClass_f_wrap(double z, void *user_data);
  }
}

class YourClass
{
 public:
  struct IntegrationResult
  {
    double result, error;
    size_t n_intervals;
  };

 public:
  double f(double x) const
  {
    return std::log(alpha_ * x) / std::sqrt(x);
  }

  IntegrationResult integrate_f() const
  {
    gsl_integration_workspace *w =
        gsl_integration_workspace_alloc(1000);

    assert(w != nullptr);

    IntegrationResult toReturn;

    gsl_function F;
    F.function = &details::YourClass_f_wrap;
    F.params = (void *)this;

    gsl_integration_qags(&F, 0, 1, 0, 1e-7, 1000, w, &toReturn.result,
                         &toReturn.error);

    toReturn.n_intervals = w->size;

    gsl_integration_workspace_free(w);

    return toReturn;
  }

 protected:
  double alpha_ = 1;
};

namespace details
{
  extern "C" {
  double YourClass_f_wrap(double z, void *user_data)
  {
    YourClass *this_ptr = (YourClass *)user_data;
    return this_ptr->f(z);
  }
  }
}

int main()
{
  YourClass yourClass;

  auto integrationResult = yourClass.integrate_f();

  double expected = -4.0;
  std::printf("result          = % .18f\n", integrationResult.result);
  std::printf("exact result    = % .18f\n", expected);
  std::printf("estimated error = % .18f\n", integrationResult.error);
  std::printf("actual error    = % .18f\n",
              integrationResult.result - expected);
  std::printf("intervals       = %zu\n",
              integrationResult.n_intervals);

  return EXIT_SUCCESS;
}

在我的电脑上我得到:

result          = -4.000000000000085265
exact result    = -4.000000000000000000
estimated error =  0.000000000000135447
actual error    = -0.000000000000085265
intervals       = 8

关于c++ - 函数c++中的gsl集成，我们在Stack Overflow上找到一个类似的问题： https://stackoverflow.com/questions/47038457/