c - 运行后未随机化的数字

Question

我正在尝试创建一个 openMP 程序，它随机化双数组并通过公式运行值:y[i] = (a[i] * b[i]) + c[i] + ( d[i] * e[i]) + (f[i]/2);

如果我多次运行该程序，我意识到 Y[] 值是相同的，即使在第一个 #pragma omp for 中初始化数组时它们应该是随机的.关于为什么会发生这种情况的任何想法？

#include<stdio.h>
#include <stdio.h>
#include <stdlib.h>
#include<omp.h>
#define ARRAY_SIZE 10

double randfrom(double min, double max);

double randfrom(double min, double max)
{
    double range = (max - min);
    double div = RAND_MAX / range;
    return min + (rand() / div);
}

int main() {
    int i;
    double a[ARRAY_SIZE], b[ARRAY_SIZE], c[ARRAY_SIZE], d[ARRAY_SIZE], e[ARRAY_SIZE], f[ARRAY_SIZE], y[ARRAY_SIZE];
    double min, max;
    int imin, imax;







    /*A[10] consists of random number in between 1 and 100
    B[10] consists of random number in between 10 and 50
    C[10] consists of random number in between 1 and 10
    D[10] consists of random number in between 1 and 50
    E[10] consists of random number in between 1 and 5
    F[10] consists of random number in between 10 and 80*/

    srand(time(NULL));

#pragma omp parallel 

    {

#pragma omp parallel for
        for (i = 0; i < ARRAY_SIZE; i++) {
            a[i] = randfrom(1, 100);
            b[i] = randfrom(10, 50);
            c[i] = randfrom(1, 50);
            d[i] = randfrom(1, 50);
            e[i] = randfrom(1, 5);
            f[i] = randfrom(10, 80);
        }
    }




    printf("This is the parallel Print\n\n\n");

#pragma omp parallel shared(a,b,c,d,e,f,y) private(i)
    {
        //Y=(A*B)+C+(D*E)+(F/2)
#pragma omp for schedule(dynamic) nowait
        for (i = 0; i < ARRAY_SIZE; i++) {
            /*printf("A[%d]%.2f",i, a[i]);
            printf("\n\n");
            printf("B[%d]%.2f", i, b[i]);
            printf("\n\n");
            printf("C[%d]%.2f", i, c[i]);
            printf("\n\n");
            printf("D[%d]%.2f", i, d[i]);
            printf("\n\n");
            printf("E[%d]%.2f", i, e[i]);
            printf("\n\n");
            printf("F[%d]%.2f", i, f[i]);
            printf("\n\n");*/
            y[i] = (a[i] * b[i]) + c[i] + (d[i] * e[i]) + (f[i] / 2);
            printf("Y[%d]=%.2f\n", i, y[i]);
        }
    }




#pragma omp parallel shared(y, min,imin,max,imax) private(i)
    {
        //min
#pragma omp for schedule(dynamic) nowait
        for (i = 0; i < ARRAY_SIZE; i++) {
            if (i == 0) {
                min = y[i];
                imin = i;
            }
            else {
                if (y[i] < min) {
                    min = y[i];
                    imin = i;
                }
            }
        }

        //max
#pragma omp for schedule(dynamic) nowait
        for (i = 0; i < ARRAY_SIZE; i++) {
            if (i == 0) {
                max = y[i];
                imax = i;
            }
            else {
                if (y[i] > max) {
                    max = y[i];
                    imax = i;
                }
            }
        }
    }
    printf("min y[%d] = %.2f\nmax y[%d] = %.2f\n", imin, min, imax, max);
    return 0;
}

Answer 1

1. 首先，我想强调的是，OpenMP 的开销很大，因此您需要在代码中进行合理的工作量，否则开销会大于并行化带来的 yield 。在您的代码中就是这种情况，因此最快的解决方案是使用串行代码。但是，您提到您的目标是学习 OpenMP，所以我将向您展示如何学习。

2. 在你的previous post's评论@paleonix 链接了一个帖子(How to generate random numbers in parallel?)，它回答了你关于随机数的问题。解决方案之一是使用 rand_r .

3. 在搜索数组 Y 的最小值和最大值时，您的代码存在数据竞争.如果您只需要找到最小值/最大值，那将非常容易，因为您可以像这样使用缩减:

double max=y[0];
#pragma omp parallel for default(none) shared(y) reduction(max:max) 
for (int i = 1; i < ARRAY_SIZE; i++) {  
    if (y[i] > max) {
        max = y[i];
    }
}

但在您的情况下，您还需要最小值和最大值的索引，因此它有点复杂。您必须使用临界区以确保其他线程无法更改 max , min , imax和 imin更新它们的值时的值。因此，可以通过以下方式完成(例如寻找最小值):

 #pragma omp parallel for
 for (int i = 0; i < ARRAY_SIZE; i++) {  
    if (y[i] < min) {
        #pragma omp critical
        if (y[i] < min) {
            min = y[i];
            imin = i;
        }
    }
 }

请注意 if (y[i] < min)出现两次，因为在第一次比较之后其他线程可能会更改 min 的值, 所以在更新之前在临界区内 min和 imin值你必须再次检查它。在寻找最大值的情况下，您可以完全按照相同的方式进行操作。

1. 始终在所需的最小范围内使用变量。

2. 同样推荐使用default(none)因此，您必须明确定义所有变量的共享属性。

3. 您可以填充数组并在单个循环中找到它的最小值/最大值，然后在不同的串行循环中打印它们的值。

4. 如果你设置min和 max在循环之前，你可以去掉额外的比较 if (i == 0)在循环内使用。

综合:

double threadsafe_rand(unsigned int* seed, double min, double max)
{
    double range = (max - min);
    double div = RAND_MAX / range;
    return min + (rand_r(seed) / div);
}

主要内容:

double min=DBL_MAX;
double max=-DBL_MAX;

#pragma omp parallel default(none) shared(a,b,c,d,e,f,y,imin,imax,min,max)  
{
    unsigned int seed=omp_get_thread_num();   
    #pragma omp for
    for (int i = 0; i < ARRAY_SIZE; i++) {  
        a[i] = threadsafe_rand(&seed, 1,100);
        b[i] = threadsafe_rand(&seed,10, 50);                    
        c[i] = threadsafe_rand(&seed,1, 10);
        d[i] = threadsafe_rand(&seed,1, 50);
        e[i] = threadsafe_rand(&seed,1, 5);
        f[i] = threadsafe_rand(&seed,10, 80);
        y[i] = (a[i] * b[i]) + c[i] + (d[i] * e[i]) + (f[i] / 2);

        if (y[i] < min) {
            #pragma omp critical
            if (y[i] < min) {
                min = y[i];
                imin = i;
            }
        }

        if (y[i] > max) {
            #pragma omp critical
            if (y[i] > max) {
                max = y[i];
                imax = i;
            }
        }
    }
}

// printout 
for (int i = 0; i < ARRAY_SIZE; i++) {
    printf("Y[%d]=%.2f\n", i, y[i]);
}
printf("min y[%d] = %.2f\nmax y[%d] = %.2f\n", imin, min, imax, max);

更新:我已经根据@Qubit 和@JérômeRichard 的建议更新了代码:

1. 我使用了“真正最小 PCG32 代码”/(c) 2014 M.E. O'Neill/来自 https://www.pcg-random.org/download.html .请注意，我不打算正确处理这个简单的随机数生成器的种子设定。如果您想这样做，请使用完整的随机数生成器库。

2. 我已更改代码以使用用户定义的缩减。确实，它使代码更加高效，但对初学者来说并不友好。这需要一篇很长的文章来解释它，所以如果你对细节感兴趣，请阅读一本关于 OpenMP 的书。

3. 我减少了 threadsafe_rand 中的划分数量

更新后的代码:

#include<stdio.h>
#include<stdint.h>
#include<time.h>
#include<float.h>
#include<limits.h>
#include<omp.h>
#define ARRAY_SIZE 10

// *Really* minimal PCG32 code / (c) 2014 M.E. O'Neill / pcg-random.org
// Licensed under Apache License 2.0 (NO WARRANTY, etc. see website)

typedef struct { uint64_t state;  uint64_t inc; } pcg32_random_t;

inline uint32_t pcg32_random_r(pcg32_random_t* rng)
{
    uint64_t oldstate = rng->state;
    // Advance internal state
    rng->state = oldstate * 6364136223846793005ULL + (rng->inc|1);
    // Calculate output function (XSH RR), uses old state for max ILP
    uint32_t xorshifted = ((oldstate >> 18u) ^ oldstate) >> 27u;
    uint32_t rot = oldstate >> 59u;
    return (xorshifted >> rot) | (xorshifted << ((-rot) & 31));
}

inline double threadsafe_rand(pcg32_random_t* seed, double min, double max)
{
    const double tmp=1.0/UINT32_MAX;  
    return min + tmp*(max - min)*pcg32_random_r(seed);
}

struct v{
 double value;
 int i;   
};

#pragma omp declare reduction(custom_min: struct v: \
 omp_out = omp_in.value < omp_out.value ? omp_in : omp_out )\
 initializer(omp_priv={DBL_MAX,0} )

#pragma omp declare reduction(custom_max: struct v: \
 omp_out = omp_in.value > omp_out.value ? omp_in : omp_out )\
 initializer(omp_priv={-DBL_MAX,0} )

int main() {
    double a[ARRAY_SIZE], b[ARRAY_SIZE], c[ARRAY_SIZE], d[ARRAY_SIZE], e[ARRAY_SIZE], f[ARRAY_SIZE], y[ARRAY_SIZE];
    struct v max={-DBL_MAX,0};
    struct v min={DBL_MAX,0}; 

    #pragma omp parallel default(none) shared(a,b,c,d,e,f,y) reduction(custom_min:min) reduction(custom_max:max)
    {
        pcg32_random_t seed={omp_get_thread_num()*7842 + time(NULL)%2299, 1234+omp_get_thread_num()};   
        #pragma omp for
        for (int i=0 ; i < ARRAY_SIZE; i++) {  
            a[i] = threadsafe_rand(&seed, 1,100);
            b[i] = threadsafe_rand(&seed,10, 50);                    
            c[i] = threadsafe_rand(&seed,1, 10);
            d[i] = threadsafe_rand(&seed,1, 50);
            e[i] = threadsafe_rand(&seed,1, 5);
            f[i] = threadsafe_rand(&seed,10, 80);
            y[i] = (a[i] * b[i]) + c[i] + (d[i] * e[i]) + (f[i] / 2);

            if (y[i] < min.value) {
                    min.value = y[i];
                    min.i = i;
                }

            if (y[i] > max.value) {
                    max.value = y[i];
                    max.i = i;
                }
        }
    }

    // printout 
    for (int i = 0; i < ARRAY_SIZE; i++) {
        printf("Y[%d]=%.2f\n", i, y[i]);
    }
    printf("min y[%d] = %.2f\nmax y[%d] = %.2f\n", min.i, min.value, max.i, max.value);
    return 0;
}