Using OpenMP to calculate the value of PI

Question

I'm trying to learn how to use OpenMP by parallelizing a monte carlo code that calculates the value of PI with a given number of iterations. The meat of the code is this:

  int chunk = CHUNKSIZE;                                                                                      

    count=0;                                                                                                  
#pragma omp parallel shared(chunk,count) private(i)                                                           
  {                                                                                                           


#pragma omp for schedule(dynamic,chunk)                                                                       
      for ( i=0; i<niter; i++) {                                                                              
        x = (double)rand()/RAND_MAX;                                                                          
        y = (double)rand()/RAND_MAX;                                                                          
        z = x*x+y*y;                                                                                          
        if (z<=1) count++;                                                                                    
      }                                                                                                       
  }                                                                                                           

  pi=(double)count/niter*4;                                                                                   
  printf("# of trials= %d , estimate of pi is %g \n",niter,pi);  

Though this is not yielding the proper value for pi given 10,000 iterations. If all the OpenMP stuff is taken out, it works fine. I should mention that I used the monte carlo code from here: http://www.dartmouth.edu/~rc/classes/soft_dev/C_simple_ex.html

I'm just using it to try to learn OpenMP. Any ideas why it's converging on 1.4ish? Can I not increment a variable with multiple threads? I'm guessing the problem is with the variable count.

Thanks!

Answer 1

Okay, I found the answer. I needed to use the REDUCTION clause. So all I had to modify was:

#pragma omp parallel shared(chunk,count) private(i)

to:

#pragma omp parallel shared(chunk) private(i,x,y,z) reduction(+:count)

Now it's converging at 3.14...yay