Complex types using OpenMP

Question

I'm a newbie using OpenMP with C++. I was doing a simple function, loading two arrays using a for loop. These arrays are defined as complex.

#include <omp.h>
#include <iostream>
#include <stdlib.h>
#include <complex>
#define CHUNKSIZE   10
#define N       100

using namespace std;

int main (int argc, char *argv[]) 
{
int nthreads, tid, i, chunk;
complex<double> a[N], b[N], c[N];

/* Some initializations */
for (i=0; i < N; i++)
    a[i].real() = b[i].real() = i * 1.0;
chunk = CHUNKSIZE;

#pragma omp parallel shared(a,b,c,nthreads,chunk) private(i,tid)
{
    tid = omp_get_thread_num();
    if (tid == 0)
    {
        nthreads = omp_get_num_threads();
        printf("Number of threads = %d\n", nthreads);
    }
    printf("Thread %d starting...\n",tid);

#pragma omp for schedule(dynamic,chunk)
    for (i=0; i<N; i++)
    {
        c[i] = a[i] + b[i];
        printf("Thread %d: c[%d]= %e\n",tid,i,c[i]);
    }

}  /* end of parallel section */

}

When I compile, I get this warning:

omp_complex.cpp:43: warning: cannot pass objects of non-POD type ‘struct std::complex’ through ‘...’; call will abort at runtime

and If I run a.out, I get an "Illegal instruction" message on the screen. I was trying to find out what's going on, but I didn't find any good reference. Do anyone know if complex types are allowed in OpenMP directives?

Answer 1

The error is here:

printf("Thread %d: c[%d]= %e\n",tid,i,c[i]);

printf doesn’t know (and has no way of knowing) how to handle a std::complex. Use C++ streaming operations to output complex types.

Furthermore, in order to avoid concurrency issues, you need to stream into a thread-local buffer before you can write to stdout, otherwise the C++ streaming syntax creates race conditions.

I usually use a macro for that (but C++11 makes this easier with variadic templates):

#define THREAD_SAFE_OUT(out, message) \
    do { \
        std::ostringstream buffer; \
        buffer << message; \
        out << buffer.str(); \
    while (false)

…

THREAD_SAFE_OUT(std::cout, "Thread " << tid << ": c[" << i << "] = " << c[i]);

---

A word on style when using OpenMP:

Don’t use the private directive. This is just a workaround for languages who require all variables to be declared at the beginning of the method. Since C++ doesn’t mandate this, it’s better (always) to declare the variables on first use, i.e. inside the parallel section. That way, they are also thread-private.

Answer 2

Your problem is passing a std::complex<double> to printf, which, as the compiler tells you, won't work. How do you want it do be displayed anyway?

You should print the real and imaginary parts seperately, or the norm and phase if you want that.