OpenMP 并行区域中嵌套函数的内部编译器错误

我尝试调用GSL库的蒙特卡罗积分子程序来进行一些数值计算。因为我的 for 循环相当简单，这意味着不同运行的结果是独立的，所以我预计使用 OpenMP 进行并行化应该非常简单。然而，当我编译它时，它总是显示“内部编译器错误：分段错误”，并且什么也没生成。这是我的代码：

#include <stdlib.h>
#include <omp.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_monte.h>
#include <gsl/gsl_monte_vegas.h>
#include <math.h>     

 double
 Reg1 (double *x, double t, size_t dim, void *params)
 {
     return sin(x[1])*cos(t*t)*x[0]*x[0]*cos(x[0])*cos(3*x[0]);
 }


 void
 display_results (char *title, double result, double error)
 {
   printf ("%s ==================\n", title);
   printf ("result = % .10f\n", result);
   printf ("sigma  = % .10f\n\n", error);
 }


 void
 VEGAS_integration_routine (double (*function)(double *, size_t, void *), 
                            double xl[], double xu[], size_t calls, gsl_rng * r)
 {
 double res, err;
     gsl_monte_function Function = { function, 2, 0 };   

     gsl_monte_vegas_state *s = gsl_monte_vegas_alloc (2);

     gsl_monte_vegas_integrate (&Function, xl, xu, 2, 20000, r, s, &res, &err);

     display_results ("vegas warm-up", res, err);

     printf ("converging...\n");

     do
       {
         gsl_monte_vegas_integrate (&Function, xl, xu, 2, calls, r, s, &res, &err);

         printf ("result = % .10f; sigma = % .10f; "
                 "chisq/dof = %.2f\n", res, err, gsl_monte_vegas_chisq (s));
       }
     while (fabs (gsl_monte_vegas_chisq (s) - 1.0) > 0.05);

     display_results ("vegas final", res, err);

     gsl_monte_vegas_free (s);
 }


 int
 main (void)
 {   
   int seg = 200;

   double xl[2] = { 195.0, -1000.0 };
   double xu[2] = { 205.0, 1000.0 };

   const gsl_rng_type *T;
   gsl_rng *r;

   size_t calls = 1000000;

   gsl_rng_env_setup ();

   T = gsl_rng_default;
   r = gsl_rng_alloc (T);

   /* Calculate G1 */
   int i;
   double j=0;
   #pragma omp parallel for shared(xl,xu,calls,r,seg) private(i,j)
   for(i=0; i<=seg; i=i+1)
   {
    j=(double)i * 0.05;
    printf("Now it's t = %.2f \n", j);
    printf("Compute Re(G1)...\n");
    double g(double * x, size_t dim, void *params)
    {
        return Reg1(x, j, dim, &params);
        }
        VEGAS_integration_routine (&g, xl, xu, calls, r);
   }
   gsl_rng_free (r);

   return 0;
 }

这段代码基本上是根据示例代码修改的GSL网页 http://www.gnu.org/software/gsl/manual/html_node/Monte-Carlo-Examples.html。在不使用 OpenMP 的情况下，它可以完美运行。但是当我使用 gcc 使用以下命令进行编译时（使用-fopenmp添加），它在我们的服务器上运行，

gcc -c -Wall -ansi -I/usr/include -L/usr/lib/gcc/x86_64-redhat-linux/4.4.4/ -lgomp -fopenmp test2.c -o test2.o
gcc -L/usr/lib64/  test2.o -L/usr/lib/gcc/x86_64-redhat-linux/4.4.4/ -lgomp -fopenmp -lgsl -lgslcblas -lm -o test2.out

我收到以下错误消息：

test2.c: In function 'main':
test2.c:85: internal compiler error: Segmentation fault
Please submit a full bug report,
with preprocessed source if appropriate.
See <http://bugzilla.redhat.com/bugzilla> for instructions.
Preprocessed source stored into /tmp/ccAGFe3v.out file, please attach this to your bugreport.
make: *** [test2.o] Error 1

由于我无法编译它并且上面显示的错误消息非常有限，我真的很想知道出了什么问题，所以我分解了子例程VEGAS_integration_routine我打电话然后逐行运行它。我发现编译停在第二行

gsl_monte_function Function = { function, 2, 0 };

这让我很困惑。使用 OpenMP 展平循环时，我不能在循环中声明 GSL 函数吗？ GSL 和 OpenMP 之间是否存在内部冲突？我确实在 Stack Overflow 和 Google 上进行了搜索，但似乎没有相关的帖子存在（太奇怪了！）。如果有人能够解释这里发生的事情，或者指出进行并行计算的替代方法，我将非常感激。我确信我编写 for 循环的方式不是最好和最简洁的......

有回归错误 http://gcc.gnu.org/bugzilla/show_bug.cgi?id=52547GCC 中词法闭包的实现与 OpenMP 代码转换引擎交互的方式。该错误似乎在 GCC 4.7 分支的更高版本中得到修复。如果您无法用 GCC 4.7.1 或更高版本替换 GCC 4.4.4，那么您应该重写代码以不使用嵌套函数，如您在帖子中引用的 GSL 示例中所示。此外，嵌套函数是非标准 GCC C 扩展，您应该尽可能避免使用它们。

GSL 中的蒙特卡罗积分接口支持通过以下方式将附加参数传递给被积函数void *params论证如所解释的here http://www.gnu.org/software/gsl/manual/html_node/Monte-Carlo-Interface.html。你必须改变Reg1如下：

double
Reg1 (double *x, size_t dim, void *params)
{
    double t = *(double *)params;
    return sin(x[1])*cos(t*t)*x[0]*x[0]*cos(x[0])*cos(3*x[0]);
}

并更新第二行VEGAS_integration_routine读书：

gsl_monte_function Function = { function, 2, &t };

您还必须添加t到虚拟参数列表VEGAS_integration_routine并将其值传递给main。因此相关部分变为：

void
VEGAS_integration_routine (double (*function)(double *, size_t, void *),
                           double t,
                           double xl[], double xu[], size_t calls, gsl_rng * r)
{
   double res, err;
   gsl_monte_function Function = { function, 2, &t };
   ...
}

#pragma omp parallel for
for(i=0; i <= seg; i++)
{
   double t = (double)i * 0.05;
   ...
   VEGAS_integration_routine (Reg1, t, xl, xu, calls, r);
}