How to properly use #pragma omp simd?

Viewed 706

I'm trying to understand whether the code below is OpenMP standard compliant. The main concern here is the args object that contains an offset field that is modified inside a loop to which #pragma omp simd is applied. Is this a legit use case?

#include <cstdio>
    
struct args_t {
    int offset;
};
    
const int n = 10;
float data1[n];
float data2[n];
    
void foo(float &res, const args_t& args) {
    res = res + data2[args.offset];
}
    
int main() {
    printf("Original arrays:\n");
    for (int i = 0; i < n; i++) {
        data1[i] = (float)i / 2.0f;
        printf("%f ", data1[i]);
    }
    printf("\n");
    for (int i = 0; i < n; i++) {
        data2[i] = (float)i / 3.0f;
        printf("%f ", data2[i]);
    }
    printf("\n");
    
    args_t args;
    args.offset = 0;
    
    #pragma omp simd
    for (int i = 0; i < n; i++) {
        foo(data1[i], args);
        args.offset++;
    }
    
    printf("Sum of two arrays:\n");
    for (int i = 0; i < n; i++)
        printf("%f ", data1[i]);
    printf("\n");
    
    return 0;
}
1 Answers

TL;DR answer: OpenMP specification is not specific in this respect, which means that the answer depends on the actual implementation. In practice, your code properly vectorized (at least on newest gcc/clang on x86-64 platform), but you can specify that your variable is modified inside the loop by using the linear clause.

Detailed answer: In the OpenMP specification the execution model of the simd construct is quite vaguely described:

The simd construct can be applied to a loop to indicate that the loop can be transformed into a SIMD loop...

This gives a lot of flexibility/freedom to the compiler, and also raises many questions - like yours. The last paragraph of this document is much more clear:

OpenMP provides directives to improve the capabilities of the compiler’s auto-vectorization pass by providing it with information that cannot be determined through compile-time static-analysis. This allows the programmer to effectively vectorize previously problematic sections of code and have it run efficiently on several computer architectures and accelerators...

This practically means that the OpenMP simd directives provide only information to the compiler for auto-vectorization, but how auto-vectorization is actually performed depends on the the implementation .

So, based on the above mentioned references and some tests with Compiler Explorer (gcc and clang on x86-64 platform) I always found that if you do not provide enough information for vectorization the worst case is that the loop won't be vectorized, but it will not result incorrect code.

I have also found that using #pragma omp simd without any additional clause or directive is practically equivalent to the use of #pragma GCC ivdep (or #pragma clang loop vectorize(assume_safety) for clang), but it is much more portable.

In the following code, the compiler generated code first checks the value of k to determine if it is safe to vectorize, but if #pragma omp simd is added this check is omitted:

void vec_dep(int *a, int k, int c, int m) {
  for (int i = 0; i < m; i++)
    a[i] = a[i + k] * c; 
}

Consider the following example:

int foo(int*  A){
    int sum=0;    
    #pragma omp simd reduction(+:sum)
    for(int i=0;i<1024;++i) 
        sum+=A[i];
    return sum;
}

In this example #pragma omp simd reduction(+:sum) is the absolutely correct form, but using #pragma omp simd or #pragma GCC ivdep or not using anything at all gives similar (correctly vectorized) code. Note that this is not the case if #pragma omp parallel for reduction(+:sum) is used, in this case reduction is absolutely necessary to avoid race condition. (Well, this raises the obvious question why the compiler does not give a warning in such a case.)

Similarly it is not necessary to use linear clause (the compiler can find this linear dependence):

#pragma omp simd linear(b:1)
 for (int i=0;i<N;++i) array[i]=b++; 

Note that, however, if #pragma omp parallel for simd linear(b) is used the linear(b) cannot be omitted otherwise the result will be incorrect, because the OpenMP calculates the initial b value for each thread using this linear relationship.

So, to answer your question, your code will compile to properly vectorized code (at least on compilers I have tested), even though the linear relationship is not specified. To specify this linear relationship you have to use the linear clause. The first idea to use #pragma omp simd linear(args.offset), but it can't compile becasue the following error: linear clause applied to non-integral non-pointer variable with 'args_t' type. The workaround is to use a reference to args.offset and change the function foo accordingly:

void foo(float &res, const int& offset) {
    res = res + data2[offset];
}
...
int& p=args.offset;
    
    #pragma omp simd linear(p)
    for (int i = 0; i < n; i++) {
        foo(data1[i], p);
        p++;
    }
Related