title: Segmentation fault on parallelization with OpenMP
body: <p>I am using OpenMP to parallel the following function</p>
  
  <pre><code>float myfunc ( Class1 *class1, int *feas, int numfeas, float z, long *k, double cost, long iter, float e )
  {
      long i;
      long x;
      double sum;
  
      sum = cost;
      while ( change/cost &gt; 1.0*e ) {
          change = 0.0;
  
          intshuffle ( feas, numfeas );
  
          #pragma omp parallel for private (i,x) firstprivate(z,k) reduction (+:sum)
          for ( i=0; i&lt;iter; i++ ) {
              x = i%numfeas;
              sum += pgain ( feas[x], class1, z, k );
          }
          cost -= sum;
      }
      return ( cost );
  }
  </code></pre>
  
  <p>This returns a "<strong>Segmentation fault (core dumped)</strong>" error</p>
  
  <p>When I change my pragma to </p>
  
  <pre><code>#pragma omp parallel for private (i,x) firstprivate(class1,feas,z,k) reduction (+:sum)`
  </code></pre>
  
  <p>I get a "<strong>double free or corruption (out): 0xb6a00468</strong>" error.</p>
  
  <p>From my limited knowledge of OpenMp i understand this is due to faulty memory access to pointers concerning <strong>class1</strong> and <strong>feas</strong>.</p>
  
  <p>For reference I post also my Class1 code</p>
  
  <pre><code>typedef struct {
      float w;
      float *c;
      long a;  
      float co;  
  } Class1;
  </code></pre>
  
  <p>Please advice on the proper way to parallelise the above function.</p>
  
  <p><strong>Update:</strong></p>
  
  <p>The pgain function</p>
  
  <pre><code>double pgain ( long x, Class1 *class1, double z, long int *numcenters )
  {
      int i;
      int number_of_centers_to_close = 0;
  
      static double *work_mem;
      static double gl_cost_of_opening_x;
      static int gl_number_of_centers_to_close;
  
      int stride = *numcenters + 2;
      //make stride a multiple of CACHE_LINE
      int cl = CACHE_LINE/sizeof ( double );
      if ( stride % cl != 0 ) {
          stride = cl * ( stride / cl + 1 );
      }
      int K = stride - 2 ; // K==*numcenters
  
      //my own cost of opening x
      double cost_of_opening_x = 0;
  
      work_mem = ( double* ) malloc ( 2 * stride * sizeof ( double ) );
      gl_cost_of_opening_x = 0;
      gl_number_of_centers_to_close = 0;
  
      /*
       * For each center, we have a *lower* field that indicates
       * how much we will save by closing the center.
       */
      int count = 0;
      for ( int i = 0; i &lt; class1-&gt;num; i++ ) {
          if ( is_center[i] ) {
              center_table[i] = count++;
          }
      }
      work_mem[0] = 0;
  
      //now we finish building the table. clear the working memory.
      memset ( switch_membership, 0, class1-&gt;num * sizeof ( bool ) );
      memset ( work_mem, 0, stride*sizeof ( double ) );
      memset ( work_mem+stride,0,stride*sizeof ( double ) );
  
      //my *lower* fields
      double* lower = &amp;work_mem[0];
      //global *lower* fields
      double* gl_lower = &amp;work_mem[stride];
  
      #pragma omp parallel for 
      for ( i = 0; i &lt; class1-&gt;num; i++ ) {
          float x_cost = dist ( class1-&gt;p[i], class1-&gt;p[x], class1-&gt;dim ) * class1-&gt;p[i].weight;
          float current_cost = class1-&gt;p[i].cost;
  
          if ( x_cost &lt; current_cost ) {
  
              // point i would save cost just by switching to x
              // (note that i cannot be a median,
              // or else dist(p[i], p[x]) would be 0)
  
              switch_membership[i] = 1;
              cost_of_opening_x += x_cost - current_cost;
  
          } else {
  
              // cost of assigning i to x is at least current assignment cost of i
  
              // consider the savings that i's **current** median would realize
              // if we reassigned that median and all its members to x;
              // note we've already accounted for the fact that the median
              // would save z by closing; now we have to subtract from the savings
              // the extra cost of reassigning that median and its members
              int assign = class1-&gt;p[i].assign;
              lower[center_table[assign]] += current_cost - x_cost;
          }
      }
  
      // at this time, we can calculate the cost of opening a center
      // at x; if it is negative, we'll go through with opening it
  
      for ( int i = 0; i &lt; class1-&gt;num; i++ ) {
          if ( is_center[i] ) {
              double low = z + work_mem[center_table[i]];
              gl_lower[center_table[i]] = low;
              if ( low &gt; 0 ) {
                  // i is a median, and
                  // if we were to open x (which we still may not) we'd close i
  
                  // note, we'll ignore the following quantity unless we do open x
                  ++number_of_centers_to_close;
                  cost_of_opening_x -= low;
              }
          }
      }
      //use the rest of working memory to store the following
      work_mem[K] = number_of_centers_to_close;
      work_mem[K+1] = cost_of_opening_x;
  
      gl_number_of_centers_to_close = ( int ) work_mem[K];
      gl_cost_of_opening_x = z + work_mem[K+1];
  
      // Now, check whether opening x would save cost; if so, do it, and
      // otherwise do nothing
  
      if ( gl_cost_of_opening_x &lt; 0 ) {
          //  we'd save money by opening x; we'll do it
          #pragma omp parallel for
          for ( int i = 0; i &lt; class1-&gt;num; i++ ) {
              bool close_center = gl_lower[center_table[class1-&gt;p[i].assign]] &gt; 0 ;
              if ( switch_membership[i] || close_center ) {
                  // Either i's median (which may be i itself) is closing,
                  // or i is closer to x than to its current median
                  #pragma omp critical
                  {
                  class1-&gt;p[i].cost = class1-&gt;p[i].weight * dist ( class1-&gt;p[i], class1-&gt;p[x], class1-&gt;dim );
                  class1-&gt;p[i].assign = x;
                  }
              }
          }
          for ( int i = 0; i &lt; class1-&gt;num; i++ ) {
              if ( is_center[i] &amp;&amp; gl_lower[center_table[i]] &gt; 0 ) {
                  is_center[i] = false;
              }
          }
          if ( x &gt;= 0 &amp;&amp; x &lt; class1-&gt;num ) {
              is_center[x] = true;
          }
  
          *numcenters = *numcenters + 1 - gl_number_of_centers_to_close;
      } else {
          gl_cost_of_opening_x = 0;  // the value we'll return
      }
  
      free ( work_mem );
  
      return -gl_cost_of_opening_x;
  }
  </code></pre>
  
  <p><strong>Update 2</strong></p>
  
  <p>Fixed pgain cause on my version i have #pragma directives also in pgain</p>
  
  <p>And here is a stack trace</p>
  
  <pre><code>#0  0xb7d50750 in ?? () from /lib/i386-linux-gnu/libc.so.6
  #1  0xb7eaa198 in ?? () from /usr/lib/i386-linux-gnu/libgomp.so.1
  #2  0x080493a0 in pgain (x=2982, class1=0xbffff218, z=1461.919921875, 
      numce=0x804d128) at streamcluster-openmp.cpp:232
  #3  0x0804aaf6 in myfunc(Class1*, int*, int, float, long*, double, long, float) [clone ._omp_fn.0] () at openmp.cpp:347
  #4  0xb7ea9889 in ?? () from /usr/lib/i386-linux-gnu/libgomp.so.1
  #5  0xb7cbcd4c in start_thread () from /lib/i386-linux-gnu/libpthread.so.0
  #6  0xb7dc9bae in clone () from /lib/i386-linux-gnu/libc.so.6
  </code></pre>