Guide into OpenMP: Easy multithreading programming for C++

• GCC (GNU Compiler Collection) supports OpenMP 4.5 since version 6.1, OpenMP 4.0 since version 4.9, OpenMP 3.1 since version 4.7, OpenMP 3.0 since version 4.4, and OpenMP 2.5 since version 4.2. Add the commandline option -fopenmp to enable it. OpenMP offloading is supported for Intel MIC targets only (Intel Xeon Phi KNL + emulation) since version 5.1, and to NVidia (NVPTX) targets since version 7 or so.
• Clang++ supports OpenMP 4.5 since version 3.9 (without offloading), OpenMP 4.0 since version 3.8 (for some parts), and OpenMP 3.1 since version 3.7. Add the commandline option -fopenmp to enable it.
• Solaris Studio supports OpenMP 4.0 since version 12.4, and OpenMP 3.1 since version 12.3. Add the commandline option -xopenmp to enable it.
• Intel C Compiler (icc) supports Openmp 4.5 since version 17.0, OpenMP 4.0 since version 15.0, OpenMP 3.1 since version 12.1, OpenMP 3.0 since version 11.0, and OpenMP 2.5 since version 10.1. Add the commandline option -openmp to enable it. Add the -openmp-stubs option instead to enable the library without actual parallel execution.
• Microsoft Visual C++ (cl) supports OpenMP 2.0 since version 2005. Add the commandline option /openmp to enable it.

Note: If your GCC complains that "-fopenmp" is valid for D but not for C++ when you try to use it, or does not recognize the option at all, your GCC version is too old. If your linker complains about missing GOMP functions, you forgot to specify "-fopenmp" in the linking.

More information: http://openmp.org/wp/openmp-compilers/

  #include <cmath>
  int main()
  {
    const int size = 256;
    double sinTable[size];
    
    #pragma omp parallel for
    for(int n=0; n<size; ++n)
      sinTable[n] = std::sin(2 * M_PI * n / size);
  
    // the table is now initialized
  }

  #include <cmath>
  int main()
  {
    const int size = 256;
    double sinTable[size];
    
    #pragma omp simd
    for(int n=0; n<size; ++n)
      sinTable[n] = std::sin(2 * M_PI * n / size);
  
    // the table is now initialized
  }

  #include <cmath>
  int main()
  {
    const int size = 256;
    double sinTable[size];
    
    #pragma omp target teams distribute parallel for map(from:sinTable[0:256])
    for(int n=0; n<size; ++n)
      sinTable[n] = std::sin(2 * M_PI * n / size);

    // the table is now initialized
  }

 #include <complex>
 #include <cstdio>
 
 typedef std::complex<double> complex;
 
 int MandelbrotCalculate(complex c, int maxiter)
 {
     // iterates z = z + c until |z| >= 2 or maxiter is reached,
     // returns the number of iterations.
     complex z = c;
     int n=0;
     for(; n<maxiter; ++n)
     {
         if( std::abs(z) >= 2.0) break;
         z = z*z + c;
     }
     return n;
 }
 int main()
 {
     const int width = 78, height = 44, num_pixels = width*height;
     
     const complex center(-.7, 0), span(2.7, -(4/3.0)*2.7*height/width);
     const complex begin = center-span/2.0;//, end = center+span/2.0;
     const int maxiter = 100000;
   
   #pragma omp parallel for ordered schedule(dynamic)
     for(int pix=0; pix<num_pixels; ++pix)
     {
         const int x = pix%width, y = pix/width;
         
         complex c = begin + complex(x * span.real() / (width +1.0),
                                     y * span.imag() / (height+1.0));
         
         int n = MandelbrotCalculate(c, maxiter);
         if(n == maxiter) n = 0;
         
       #pragma omp ordered
         {
           char c = ‘ ‘;
           if(n > 0)
           {
               static const char charset[] = ".,c8M@jawrpogOQEPGJ";
               c = charset[n % (sizeof(charset)-1)];
           }
           std::putchar(c);
           if(x+1 == width) std::puts("|");
         }
     }
 }

Only when the compiler interprets those #pragma lines, it becomes a parallel program. It really does calculate N values simultaneously where N is the number of threads. In GCC, libgomp determines that from the number of processors.

By C and C++ standards, if the compiler encounters a #pragma that it does not support, it will ignore it. So adding the OMP statements can be done safely^[1] without breaking compatibility with legacy compilers.

There is also a runtime library that can be accessed through omp.h, but it is less often needed. If you need it, you can check the #define _OPENMP for conditional compilation in case of compilers that don‘t support OpenMP.

[1]: Within the usual parallel programming issues (concurrency, mutual exclusion) of course.

  #pragma omp parallel
  {
    // Code inside this region runs in parallel.
    printf("Hello!\n");
  }

Internally, GCC implements this by creating a magic function and moving the associated code into that function, so that all the variables declared within that block become local variables of that function (and thus, locals to each thread).
ICC, on the other hand, uses a mechanism resembling fork(), and does not create a magic function. Both implementations are, of course, valid, and semantically identical.

Variables shared from the context are handled transparently, sometimes by passing a reference and sometimes by using register variables which are flushed at the end of the parallel block (or whenever a flush is executed).

Parallelism conditionality clause: if

  extern int parallelism_enabled;
  #pragma omp parallel for if(parallelism_enabled)
  for(int c=0; c<n; ++c)
    handle(c);

 #pragma omp for
 for(int n=0; n<10; ++n)
 {
   printf(" %d", n);
 }
 printf(".\n");

  0 5 6 7 1 8 2 3 4 9.

  int this_thread = omp_get_thread_num(), num_threads = omp_get_num_threads();
  int my_start = (this_thread  ) * 10 / num_threads;
  int my_end   = (this_thread+1) * 10 / num_threads;
  for(int n=my_start; n<my_end; ++n)
    printf(" %d", n);

Note: #pragma omp for only delegates portions of the loop for different threads in the current team. A team is the group of threads executing the program. At program start, the team consists only of a single member: the master thread that runs the program.

To create a new team of threads, you need to specify the parallel keyword. It can be specified in the surrounding context:

 #pragma omp parallel
 {
  #pragma omp for
  for(int n=0; n<10; ++n) printf(" %d", n);
 }
 printf(".\n");

 #pragma omp parallel for
 for(int n=0; n<10; ++n) printf(" %d", n);
 printf(".\n");

 #pragma omp parallel num_threads(3)
 {
   // This code will be executed by three threads.
   
   // Chunks of this loop will be divided amongst
   // the (three) threads of the current team.
   #pragma omp for
   for(int n=0; n<10; ++n) printf(" %d", n);
 }

 int n;
 #pragma omp for private(n)
 for(n=0; n<10; ++n) printf(" %d", n);
 printf(".\n");

In OpenMP 2.5, the iteration variable in for must be a signed integer variable type. In OpenMP 3.0, it may also be an unsigned integer variable type, a pointer type or a constant-time random access iterator type. In the latter case, std::distance() will be used to determine the number of loop iterations.

What are: parallel, for and a team

Scheduling

 #pragma omp for schedule(static)
 for(int n=0; n<10; ++n) printf(" %d", n);
 printf(".\n");

 #pragma omp for schedule(dynamic)
 for(int n=0; n<10; ++n) printf(" %d", n);
 printf(".\n");

 #pragma omp for schedule(dynamic, 3)
 for(int n=0; n<10; ++n) printf(" %d", n);
 printf(".\n");

  int a,b;
  if(GOMP_loop_dynamic_start(0,10,1, 3, &a,&b))
  {
    do {
      for(int n=a; n<b; ++n) printf(" %d", n);
    } while(GOMP_loop_dynamic_next(&a,&b));
  }

The ordered clause

 #pragma omp for ordered schedule(dynamic)
 for(int n=0; n<100; ++n)
 {
   files[n].compress();

   #pragma omp ordered
   send(files[n]);
 }

The collapse clause

 #pragma omp parallel for collapse(2)
 for(int y=0; y<25; ++y)
   for(int x=0; x<80; ++x)
   {
     tick(x,y);
   }

The reduction clause

 int sum=0;
 #pragma omp parallel for reduction(+:sum)
 for(int n=0; n<1000; ++n) sum += table[n];

 #pragma omp sections
 {
   { Work1(); }
   #pragma omp section
   { Work2();
     Work3(); }
   #pragma omp section
   { Work4(); }
 }

As usual, if the compiler ignores the pragmas, the result is still a correctly running program.

Internally, GCC implements this as a combination of the parallel for and a switch-case construct. Other compilers may implement it differently.

Note: #pragma omp sections only delegates the sections for different threads in the current team. To create a team, you need to specify the parallel keyword either in the surrounding context or in the pragma, as #pragma omp parallel sections.
Example:

 #pragma omp parallel sections // starts a new team
 {
   { Work1(); }
   #pragma omp section
   { Work2();
     Work3(); }
   #pragma omp section
   { Work4(); }
 }

 #pragma omp parallel // starts a new team
 {
   //Work0(); // this function would be run by all threads.
   
   #pragma omp sections // divides the team into sections
   { 
     // everything herein is run only once.
     { Work1(); }
     #pragma omp section
     { Work2();
       Work3(); }
     #pragma omp section
     { Work4(); }
   }
   
   //Work5(); // this function would be run by all threads.
 }

There are two use cases for the simd construct.

Firstly, #pragma omp simd can be used to declare that a loop will be utilizing SIMD.

 float a[8], b[8];
 ...
 #pragma omp simd
 for(int n=0; n<8; ++n) a[n] += b[n];

  #pragma omp declare simd aligned(a,b:16)
  void add_arrays(float *__restrict__ a, float *__restrict__ b)
  {
    #pragma omp simd aligned(a,b:16)
    for(int n=0; n<8; ++n) a[n] += b[n];
  }	

Since compilers of today attempt to do SIMD regardless of OpenMP simd directives, the simd directive can be thought essentially as a directive to the compiler, saying: “Try harder”.

The collapse clause

  #pragma omp simd collapse(2)
  for(int i=0; i<4; ++i)
    for(int j=0; j<4; ++j)
      a[j*4+i] += b[i*4+j];

The reduction clause

 int sum=0;
 #pragma omp simd reduction(+:sum)
 for(int n=0; n<1000; ++n) sum += table[n];

The aligned clause

  #pragma omp declare simd aligned(a,b:16)
  void add_arrays(float *__restrict__ a, float *__restrict__ b)
  {
    #pragma omp simd aligned(a,b:16)
    for(int n=0; n<8; ++n) a[n] += b[n];
  }	

The safelen clause

  #pragma omp declare simd
  void add_arrays(float* a, float* b)
  {
    #pragma omp simd aligned(a,b:16) safelen(4)
    for(int n=0; n<8; ++n) a[n] += b[n];
  }	

The simdlen clause (OpenMP 4.5+)

The uniform clause

The linear clause (OpenMP 4.5+)

  #include <stdio.h>

  int b = 10;
  int main()
  {
    int array[8];

    #pragma omp simd linear(b:2)
    for(int n=0; n<8; ++n) array[n] = b;

    for(int n=0; n<8; ++n) printf("%d\n", array[n]);
  }

  int b_original = b;
  for(int n=0; n<8; ++n) array[n] = b_original + n*2;

The inbranch and notinbranch clauses

  float sum(float* table)
  {
    float result=0;
    #pragma omp parallel for simd reduction(+:result)
    for(int n=0; n<1000; ++n) result += table[n];
    return result;
  }

These examples are from the OpenMP 3.0 manual:

struct node { node *left, *right; };
extern void process(node* );
void traverse(node* p)
{
    if (p->left)
        #pragma omp task // p is firstprivate by default
        traverse(p->left);
    if (p->right)
        #pragma omp task // p is firstprivate by default
        traverse(p->right);
    process(p);
}

struct node { node *left, *right; };
extern void process(node* );
void postorder_traverse(node* p)
{
    if (p->left)
        #pragma omp task // p is firstprivate by default
        postorder_traverse(p->left);
    if (p->right)
        #pragma omp task // p is firstprivate by default
        postorder_traverse(p->right);
    #pragma omp taskwait
    process(p);
}

struct node { int data; node* next; };
extern void process(node* );
void increment_list_items(node* head)
{
    #pragma omp parallel
    {
        #pragma omp single
        {
            for(node* p = head; p; p = p->next)
            {
            	#pragma omp task
            	process(p); // p is firstprivate by default
            }
        }
    }
}

Example:

#pragma omp declare target
int x;
void murmur() { x+=5; }
#pragma omp end declare target

These two functions and variables are separate, and may contain values separate from each others.

Variables declared in this manner can be accessed by the device code without separate map clauses.

OpenACC differences

In OpenACC, device-functions are declared by prefixing each function with #pragma acc routine. Its data model is more complicated and has no direct translation from/to OpenMP.

The target construct executes the construct on a device (and also has target data features).

These two constructs are identical in effect:

#pragma omp target // device()... map()... if()...
{
  <<statements...>>
}

#pragma omp target data // device()... map()... if()...
{
  #pragma omp target
  {
    <<statements...>>
  }
}

To utilize parallelism on device, you have to engage a teams construct inside the target construct. Example:

  #include <stdio.h>

  long long r = 1;

  int main(void)
  {
    r=10;
    #pragma omp target teams distribute parallel for reduction(+:r) map(tofrom:r)
    for(unsigned long long n=0; n<0x800000000ull; ++n)
        r += n;
    printf("r=%llX\n", r);
    return 0;
  }

The if clause

  #include <stdlib.h>
  #include <stdio.h>

  int main(int argc, char** argv)
  {
    int r=0;

    #pragma omp target if(atoi(argv[1])) map(tofrom:r)
    r += 4;

    printf("r=%d\n", r);
  }

The device clause

  int device_number = ...;
  #pragma omp target device(device_number)
  {
    //...
  }

The map clause

Array sections (OpenMP 4.0+)

  void foo (int *p)
  {
    int i;
    #pragma omp parallel
    #pragma omp single
    #pragma omp target teams distribute parallel for map(p[0:24])
    for (i = 0; i < 24; i++)
      p[i] = p[i] + 1;
  }

Examples:

• #pragma omp target enter data map(from:var)
• #pragma omp target exit data map(to:var)

• #pragma omp target update from(c)

For example, if there are 20 teams, the loop will be distributed across 20 master threads.

  #include <stdio.h>

  int main(void)
  {
    int r=0;

    #pragma omp target teams distribute reduction(+:r)
    for(int n=0; n<10000; ++n)
      r += n;

    printf("r=%d\n", r);
    return 0;
  }

In OpenACC this behavior is achieved by adding the word gang to existing worksharing constructs like #pragma acc parallel and #pragma acc kernels.

The distribute simd construct

The dist_schedule clause

For example, if there are 20 teams, and each team consists of 256 threads, the loop will be distributed across 5120 threads.

  #include <stdio.h>

  int main(void)
  {
    int r=0;

    #pragma omp target teams distribute parallel for reduction(+:r)
    for(int n=0; n<10000; ++n)
      r += n;

    printf("r=%d\n", r);
    return 0;
  }

The simd clause can be added once again to the loop to add SIMD execution, if possible.

OpenACC differences

In OpenACC this behavior is achieved by adding the word worker to existing worksharing constructs like #pragma acc parallel and #pragma acc kernels. Additionally the word vector can be added to achieve SIMD parallelism as well.

 #pragma omp atomic
 counter += value;

See also reduction.

There are four different types of atomic expressions (since OpenMP 3.1):

Atomic read expressions

 #pragma omp atomic read
 var = x;

Atomic write expressions

 #pragma omp atomic write
 x = expr;

Atomic update expressions

 #pragma omp atomic update // The word "update" is optional
 // One of these:
 ++x; --x; x++; x--;
 x += expr;  x -= expr;  x *= expr;   x /= expr;  x &= expr;
 x = x+expr; x = x-expr; x = x*expr;  x = x/expr; x = x&expr;
 x = expr+x; x = expr-x; x = expr*x;  x = expr/x; x = expr&x;
 x |= expr;  x ^= expr;  x <<= expr;  x >>= expr;
 x = x|expr; x = x^expr; x = x<<expr; x = x>>expr;
 x = expr|x; x = expr^x; x = expr<<x; x = expr>>x;

Atomic capture expressions

 #pragma omp atomic capture
 // One of these:
 var = x++;  /* Or any other of the update expressions listed above */
 { var = x; x++; /* Or any other of of the update expressions listed above */ }
 { x++; /* Or any other of of the update expressions listed above */; var = x; }
 { var = x; x = expr; }

The critical construct may optionally contain a global name that identifies the type of the critical construct. No two threads can execute a critical construct of the same name at the same time.

If the name is omitted, a default name is assumed.

 #pragma omp critical(dataupdate)
 {
   datastructure.reorganize();
 }
 ...
 #pragma omp critical(dataupdate)
 {
   datastructure.reorganize_again();
 }

Note: The critical section names are global to the entire program (regardless of module boundaries). So if you have a critical section by the same name in multiple modules, not two of them can be executed at the same time.

If you need something like a local mutex, see below.

The lock type has five manipulator functions:

• omp_init_lock initializes the lock. After the call, the lock is unset.
• omp_destroy_lock destroys the lock. The lock must be unset before this call.
• omp_set_lock attempts to set the lock. If the lock is already set by another thread, it will wait until the lock is no longer set, and then sets it.
• omp_unset_lock unsets the lock. It should only be called by the same thread that set the lock; the consequences of doing otherwise are undefined.
• omp_test_lock attempts to set the lock. If the lock is already set by another thread, it returns 0; if it managed to set the lock, it returns 1.

Here is an example of a wrapper around std::set<> that provides per-instance mutual exclusion while still working even if the compiler does not support OpenMP.

You can maintain backward compability with non-OpenMP-supporting compilers by enclosing the library references in #ifdef _OPENMP…#endif blocks.

 #ifdef _OPENMP
 # include <omp.h>
 #endif
 #include <set>
 
 class data
 {
 private:
   std::set<int> flags;
 #ifdef _OPENMP
   omp_lock_t lock;
 #endif
 public:
   data() : flags()
   {
 #ifdef _OPENMP
     omp_init_lock(&lock);
 #endif
   }
   ~data()
   {
 #ifdef _OPENMP
     omp_destroy_lock(&lock);
 #endif
   }
   
   bool set_get(int c)
   {
   #ifdef _OPENMP
     omp_set_lock(&lock);
   #endif
     bool found = flags.find(c) != flags.end();
     if(!found) flags.insert(c);
   #ifdef _OPENMP
     omp_unset_lock(&lock);
   #endif
     return found;
   }
 };

 #ifdef _OPENMP
 # include <omp.h>
 struct MutexType
 {
   MutexType() { omp_init_lock(&lock); }
   ~MutexType() { omp_destroy_lock(&lock); }
   void Lock() { omp_set_lock(&lock); }
   void Unlock() { omp_unset_lock(&lock); }
   
   MutexType(const MutexType& ) { omp_init_lock(&lock); }
   MutexType& operator= (const MutexType& ) { return *this; }
 public:
   omp_lock_t lock;
 };
 #else
 /* A dummy mutex that doesn‘t actually exclude anything,
  * but as there is no parallelism either, no worries. */
 struct MutexType
 {
   void Lock() {}
   void Unlock() {}
 };
 #endif
 
 /* An exception-safe scoped lock-keeper. */
 struct ScopedLock
 {
   explicit ScopedLock(MutexType& m) : mut(m), locked(true) { mut.Lock(); }
   ~ScopedLock() { Unlock(); }
   void Unlock() { if(!locked) return; locked=false; mut.Unlock(); }
   void LockAgain() { if(locked) return; mut.Lock(); locked=true; }
 private:
   MutexType& mut;
   bool locked;
 private: // prevent copying the scoped lock.
   void operator=(const ScopedLock&);
   ScopedLock(const ScopedLock&);
 };

 #include <set>
 
 class data
 {
 private:
   std::set<int> flags;
   MutexType lock;
 public:
   bool set_get(int c)
   {
     ScopedLock lck(lock); // locks the mutex
     
     if(flags.find(c) != flags.end()) return true; // was found
     flags.insert(c);
     return false; // was not found
   } // automatically releases the lock when lck goes out of scope.
 };

This example comes from the OpenMP specification.

          /* presumption: int a = 0, b = 0; */
                        
    /* First thread */                /* Second thread */
    b = 1;                            a = 1;
    #pragma omp flush(a,b)            #pragma omp flush(a,b)
    if(a == 0)                        if(b == 0)
    {                                 {
      /* Critical section */            /* Critical section */
    }                                 }

You need the flush directive when you have writes to and reads from the same data in different threads.

If the program appears to work correctly without the flush directive, it does not mean that the flush directive is not required. It just may be that your compiler is not utilizing all the freedoms the standard allows it to do. You need the flush directive whenever you access shared data in multiple threads: After a write, before a read.

However, I do not know these:

• Is flush needed if the shared variable is declared volatile?
• Is flush needed if all access to the shared variable is atomic or restricted by critical sections?

 int a, b=0;
 #pragma omp parallel for private(a) shared(b)
 for(a=0; a<50; ++a)
 {
   #pragma omp atomic
   b += a;
 }

The difference between private and firstprivate

 #include <string>
 #include <iostream>
 
 int main()
 {
     std::string a = "x", b = "y";
     int c = 3;
     
     #pragma omp parallel private(a,c) shared(b) num_threads(2)
     {
         a += "k";
         c += 7;
         std::cout << "A becomes (" << a << "), b is (" << b << ")\n";
     }
 }

 int main()
 {
     std::string a = "x", b = "y";
     int c = 3;
     
     OpenMP_thread_fork(2);
     {                  // Start new scope
         std::string a; // Note: It is a new local variable.
         int c;         // This too.
         a += "k";
         c += 7;
         std::cout << "A becomes (" << a << "), b is (" << b << ")\n";
     }                  // End of scope for the local variables
     OpenMP_join();
 }

 #include <string>
 #include <iostream>
 
 int main()
 {
     std::string a = "x", b = "y";
     int c = 3;
     
     #pragma omp parallel firstprivate(a,c) shared(b) num_threads(2)
     {
         a += "k";
         c += 7;
         std::cout << "A becomes (" << a << "), b is (" << b << ")\n";
     }
 }

• In a loop construct (for construct), the last value is the value assigned by the thread that handles the last iteration of the loop. Values assigned during other iterations are ignored.
• In a sections construct (sections construct), the last value is the value assigned in the last section denoted by the section construct. Values assigned in other sections are ignored.

Example:

 #include <stdio.h>
 int main()
 {
    int done = 4, done2 = 5;
    
     #pragma omp parallel for lastprivate(done, done2) num_threads(2) schedule(static)
     for(int a=0; a<8; ++a)
     {
       if(a==2) done=done2=0;
       if(a==3) done=done2=1;
     }
     printf("%d,%d\n", done,done2);
 }

 #include <stdio.h>
 int main()
 {
    int done = 4, done2 = 5;
    OpenMP_thread_fork(2);
    {
        int this_thread = omp_get_thread_num(), num_threads = 2;
        int my_start = (this_thread  ) * 8 / num_threads;
        int my_end   = (this_thread+1) * 8 / num_threads;

        int priv_done, priv_done2; // not initialized, because firstprivate was not used

        for(int a=my_start; a<my_end; ++a)
        {
            if(a==2) priv_done=priv_done2=0;
            if(a==3) priv_done=priv_done2=1;
        }
        if(my_end == 8)
        {
           // assign the values back, because this was the last iteration
           done  = priv_done;
           done2 = priv_done2;
        }
    }
    OpenMP_join();
 }

Therefore, lastprivate cannot be used to e.g. fetch the value of a flag assigned randomly during a loop. Use reduction for that, instead.

Where this behavior can be utilized though, is in situations like this (from OpenMP manual):

 void loop()
 {
   int i;
   #pragma omp for lastprivate(i)
   for(i=0; i<get_loop_count(); ++i) // note: get_loop_count() must be a pure function.
       { ... }
   
   printf("%d\n", i); // this shows the number of loop iterations done.
 }

 int a, b=0;
 // This code won‘t compile: It requires explicitly
 // specifying whether a is shared or private.
 #pragma omp parallel default(none) shared(b)
 {
   b += a;
 }

Note: Because different compilers have different ideas about which variables are implicitly private or shared, and for which it is an error to explicitly state the private/shared status, it is recommended to use the default(none) setting only during development, and drop it in production/distribution code.

This example calculates factorial using threads:

 int factorial(int number)
 {
   int fac = 1;
   #pragma omp parallel for reduction(*:fac)
   for(int n=2; n<=number; ++n)
     fac *= n;
   return fac;
 }

• At the beginning of the parallel block, a private copy is made of the variable and preinitialized to a certain value .
• At the end of the parallel block, the private copy is atomically merged into the shared variable using the defined operator.

The syntax of the clause is:

  reduction(operator:list)

To write the factorial function (shown above) without reduction, it probably would look like this:

 int factorial(int number)
 {
   int fac = 1;
   #pragma omp parallel for
   for(int n=2; n<=number; ++n)
   {
     #pragma omp atomic
     fac *= n;
   }
   return fac;
 }

The version with reduction is equivalent to this code (illustration only):

 int factorial(int number)
 {
   int fac = 1;
   #pragma omp parallel
   {
     int omp_priv = 1; /* This value comes from the table shown above */
     #pragma omp for nowait
     for(int n=2; n<=number; ++n)
       omp_priv *= n;
     #pragma omp atomic
     fac *= omp_priv;
   }
   return fac;
 }

The restrictions in reduction and atomic are very similar: both can only be done on POD types; neither allows overloaded operators, and both have the same set of supported operators.

As an example of how the reduction clause can be used to produce semantically different code when OpenMP is enabled and when it is disabled, this example prints the number of threads that executed the parallel block:

 int a = 0;
 #pragma omp parallel reduction (+:a)
 {
   a = 1; // Assigns a value to the private copy.
   // Each thread increments the value by 1.
 }
 printf("%d\n", a);

The declare reduction directive (OpenMP 4.0+)

  #include <cstdio>

  int compress(int param1, int param2)
  {
    return (param1+13)^param2; // Placeholder for a compression algorithm
  }

  int main(int argc, char** argv)
  {
    struct BestInfo { unsigned size, param1, param2; };

    #pragma omp declare reduction(isbetter:BestInfo:                                   omp_in.size<omp_out.size ? omp_out=omp_in : omp_out \   
                    ) initializer(omp_priv = BestInfo{~0u,~0u,~0u})  

    BestInfo result{~0u,~0u,~0u};
    #pragma omp parallel for collapse(2) reduction(isbetter:result)
    for(unsigned p1=0; p1<10; ++p1)
    for(unsigned p2=0; p2<10; ++p2)
    {
      unsigned size = compress(p1,p2);
      if(size < result.size) result = BestInfo{size,p1,p2};
    }
    std::printf("Best compression (%u bytes) with params %u,%u\n",
      result.size, result.param1, result.param2);
  }

 #pragma omp parallel
 {
   /* All threads execute this. */
   SomeCode();
   
   #pragma omp barrier
   
   /* All threads execute this, but not before
    * all threads have finished executing SomeCode().
    */
   SomeMoreCode();
 }

Example:

 #pragma omp parallel
 {
   #pragma omp for
   for(int n=0; n<10; ++n) Work();
   
   // This line is not reached before the for-loop is completely finished
   SomeMoreCode();
 }

 // This line is reached only after all threads from
 // the previous parallel block are finished.
 CodeContinues();

 #pragma omp parallel
 {
   #pragma omp for nowait
   for(int n=0; n<10; ++n) Work();
   
   // This line may be reached while some threads are still executing the for-loop.
   SomeMoreCode();
 }

 // This line is reached only after all threads from
 // the previous parallel block are finished.
 CodeContinues();

 #pragma omp parallel
 {
   Work1();
   #pragma omp single
   {
     Work2();
   }
   Work3();
 }

Note: Do not assume that the single block is executed by whichever thread gets there first. According to the standard, the decision of which thread executes the block is implementation-defined, and therefore making assumptions on it is non-conforming.

The master construct is similar, except that the statement/block is run by the master thread, and there is no implied barrier; other threads skip the construct without waiting.

 #pragma omp parallel
 {
   Work1();
   
   // This...
   #pragma omp master
   {
     Work2();
   }
   
   // ...is practically identical to this:
   if(omp_get_thread_num() == 0)
   {
     Work2();
   }
   
   Work3();
 }

 #pragma omp parallel for
 for(int y=0; y<25; ++y)
 {
   #pragma omp parallel for
   for(int x=0; x<80; ++x)
   {
     tick(x,y);
   }
 }

At the entrance of the inner parallel construct, the OpenMP runtime library (libgomp in case of GCC) detects that there already exists a team, and instead of a new team of N threads, it will create a team consisting of only the calling thread.

Rewriting the code like this won‘t work:

 #pragma omp parallel for
 for(int y=0; y<25; ++y)
 {
   #pragma omp for // ERROR, nesting like this is not allowed.
   for(int x=0; x<80; ++x)
   {
     tick(x,y);
   }
 }

Solution in OpenMP 3.0

 #pragma omp parallel for collapse(2)
 for(int y=0; y<25; ++y)
   for(int x=0; x<80; ++x)
   {
     tick(x,y);
   }

This problem only affects GCC; ICC is not affected.

If your program intends to become a background process using daemonize() or other similar means, you must not use the OpenMP features before the fork. After OpenMP features are utilized, a fork is only allowed if the child process does not use OpenMP features, or it does so as a completely new process (such as after exec()).

This is an example of an erroneous program:

  #include <stdio.h>
  #include <sys/wait.h>
  #include <unistd.h>

  void a()
  {
    #pragma omp parallel num_threads(2)
    {
      puts("para_a"); // output twice
    }
    puts("a ended"); // output once
  }
  void b()
  {
    #pragma omp parallel num_threads(2)
    {
      puts("para_b");
    }
    puts("b ended");
  }

  int main() {
   a();   // Invokes OpenMP features (parent process)
   int p = fork();
   if(!p)
   {
     b(); // ERROR: Uses OpenMP again, but in child process
     _exit(0);
   }
   wait(NULL);
   return 0;
  }

There is currently no workaround; the libgomp API does not specify functions that can be used to prepare for a call to fork().