线程安全函数指针 gsl 蒙特卡罗集成答案

【问题标题】：thread safe function pointer gsl monte carlo integration线程安全函数指针 gsl 蒙特卡罗集成
【发布时间】：2016-06-01 04:49:30
【问题描述】：

我想使用 gsl monte carlo 集成库，但我目前遇到了目前无法解决的线程安全问题。

我有一个静态函数 gsl_func_wrapper，它具有 gsl（非成员函数指针）所需的正确类型。它被传递给 Integral 成员函数内的 gsl monte carlo 积分器。

class IntegralStrategyGSL2D: public IntegralStrategy2D {
  static Model2D *current_model;

  const gsl_rng_type *T;
  gsl_rng *r;
  gsl_monte_vegas_state *s;
  size_t calls;
  size_t maxcalls;

  static double gsl_func_wrapper(double *x, size_t dim, void *params) {
    return current_model->eval(x);
  }
...
}

IntegralStrategyGSL2D::IntegralStrategyGSL2D() :
    T(gsl_rng_default) {
  calls = 500;    // keep calls low at first
  maxcalls = 100000;
  gsl_rng_env_setup();
  r = gsl_rng_alloc(T);
  s = gsl_monte_vegas_alloc(2);
}

double IntegralStrategyGSL2D::Integral(Model2D *model2d, double xlow,
    double xhigh, double ylow, double yhigh, double precision) {
  current_model = model2d;
  double result, error;

  double xl[2] = { xlow, ylow };
  double xu[2] = { xhigh, yhigh }

  gsl_monte_function G = { &gsl_func_wrapper, 2, 0 };

  gsl_monte_vegas_init(s);

  while (true) {
    gsl_monte_vegas_integrate(&G, xl, xu, 2, calls, r, s, &result, &error);
    ...
  }
  ...
}

如果我使用 1 个线程运行，所有这些都有效，如果我使用多个线程，我偶尔会得到 nan 积分。由于 gsl 需要这个非成员函数指针，我现在真的不知道如何使这个线程安全......有人能指出我正确的方向吗？谢谢

编辑：我添加了将由 gsl_func_wrapper 调用的 eval 函数。

double GaussianModel2D::eval(const double *x) const {
  // see wikipedia definition

  double normalization = 1.0
        / (2.0 * M_PI * gauss_sigma_var1->getValue()
                * gauss_sigma_var2->getValue()
                * sqrt(1 - pow(gauss_rho->getValue(), 2.0)));

  double exp_value = exp(
        -(pow(x[0] - gauss_mean_var1->getValue(), 2.0)
                / pow(gauss_sigma_var1->getValue(), 2.0)
                + pow(x[1] - gauss_mean_var2->getValue(), 2.0)
                        / pow(gauss_sigma_var2->getValue(), 2.0)
                + 2.0 * gauss_rho->getValue()
                        * (x[0] - gauss_mean_var1->getValue())
                        * (x[1] - gauss_mean_var2->getValue())
                        / (gauss_sigma_var1->getValue()
                                * gauss_sigma_var2->getValue()))
                / (2.0 * (1 - pow(gauss_rho->getValue(), 2.0))));

  return gauss_amplitude->getValue() * normalization * exp_value;
}

编辑2： 我按照 Mike 的要求添加了对 Integral 函数的调用。下面是进行 Integral 调用的两个相关代码块。我真的不明白为什么 current_model 变量的 thread_local 关键字解决了这个问题。虽然代码有点难看，但 current_model 变量不应该总是指向同一个 model2d 指针吗？所以即使两个线程可能有冲突的写入，写入的值总是相同的，所以它应该没关系吗？

double smearing_probability = divergence_model->Integral(int_range, 1e-4);

double Model2D::Integral(const std::vector<DataStructs::DimensionRange &ranges,
    double precision) {
    return integral_strategy->Integral(this, ranges[0].range_low,
        ranges[0].range_high, ranges[1].range_low, ranges[1].range_high,
        precision);
}

【问题讨论】：

显示Model2D::eval的代码。
你可以使用类作为函数指针来避免这个问题，但不是很好。
@Pete: Model2D::eval 是一个纯虚函数，但对我来说它是一个非常复杂的函数，所以我不确定我是否应该在这里发布所有这些......这相关吗？
它是在回调中执行的代码的一部分，所以，是的，它是相关的。但是，如果它非常复杂，您应该花一些时间来修剪它，以便它是最短的，并且仍然可以显示问题。
@Pete：我添加了将在编辑中调用的函数。它是一个缓存函数，所以只有缓存的访问会被积分器实际调用。我想这应该就是全部了。

标签： c++ thread-safety function-pointers gsl

【解决方案1】：

假设Integral 在执行gsl_func_wrapper 的同一线程中被调用，您可以使用thread_local 而不是static 作为current_model 的存储类。

这将为每个线程提供一个单独的 current_model 变量。

【讨论】：

我对此进行了测试，它似乎有效，明天我将对其进行更彻底的测试并提供更新。非常感谢迈克。
好的，所以它只是偶然比平时运行的时间长了一点。它仍然偶尔会给我积分的 nan 结果......
@steve：那么问题很可能不（仅）在您向我们展示的代码中。请提供mcve。您正在谈论多线程的问题，但就目前而言，完全不清楚您在哪里产生这些线程以及每个线程正在访问哪些数据
好的，我明白了... 只提取其中的片段给你一个 mcve 有点困难。也许我可以通过一个更简单的问题进一步了解：所有线程都评估相同的函数，因为我只给它一个指针，对吗？因此，在该函数中创建任何局部变量已经会引入线程安全问题，或者？

【解决方案2】：

好的，所以我设法找出了问题所在。那是变量

const gsl_rng_type *T; gsl_rng *r; gsl_monte_vegas_state *s;

都是实例变量，因此所有线程都以相同的方式使用。如果我在 Integral() 函数中为每个线程分别创建其中的每一个，一切都运行良好......

【讨论】：