gsl_integration_qag failed with gsl openmp

Question

gsl_integration_qag works with 1 core (with/without openMP), but fails with multi-threads (i.e. >1).

Some information that may help...

1. gsl-2.5

2. #define _OPENMP 201107

3. Depending on the number of cores, I can get error reports of:

   gsl: qag.c:248: ERROR: roundoff error prevents tolerance from being achieved (comment: usually with a small number of cores)
   gsl: qag.c:257: ERROR: maximum number of subdivisions reached (comment: usually with a large number of cores)
   

4. A large max iteration number given to gsl_integration_qag only delays the code to crash.

5. The integration function is (can be more specific if needed):

   double Func(double Param1, ..., double ParamN){
     double result, error;
     gsl_function F;
     gsl_integration_workspace * w
     = gsl_integration_workspace_alloc (1000);
   
     struct parameters_gsl_int_ parameters_gsl = {
       .Param1 = Param1,
        ...
       .ParamN = ParamN,};
   
     F.function = &func_integrand;
     F.params = &parameters_gsl;
   
     gsl_integration_qag (&F, LOWER_LIMIT, UPPER_LIMIT, 0, 0.001,
                         1000, GSL_INTEG_GAUSS61, w, &result, &error);
     gsl_integration_workspace_free (w);
   
     return result;
   }
   

6. The OpenMP part that calls the integration is:

   void call_Func(int Nbin, double array[], double Param1[], double Param2, ... double ParamN){
     int i;
   ...
   #pragma omp parallel shared(Nbin, array, Param1, ..., ParamN) private(i)
   {
   #pragma omp for    
     for (i=0; i<Nbin; i++)
       array[i] = Func(Param1[i], Param2, ..., ParamN);
   }
   ...
   }
   

I'm new to both GSL and openMP. I hope I am using gsl_integration_qag correctly and the definition of shared or private variables makes sense.

btw, it's the same question as this 2014 one (gsl openmp failed integration), but I couldn't find the solution in this post.

Answer 1

Problem solved...

It is actually due to func_integrand having also a term which is estimated using gsl_integration_qag. There were some global variables adopted in this calculation, which I didn't capture before.