c - 如何安全地将 for 循环与内部的 memcpy 并行

Question

我正在 KSVD 包中开发原始串行代码以支持 OpenMP。原始代码类似于MATLAB中的im2col，从图像中提取patches，如下所示:

/* n stands for the size of an image, sz stands for the patch size to extract */
int blocknum = 0;
for (k=0; k<=n[2]-sz[2]; k+=1) {
    for (j=0; j<=n[1]-sz[1]; j+=1) {
        for (i=0; i<=n[0]-sz[0]; i+=1) {

            /* copy single block */
            for (m=0; m<sz[2]; m++) {
                for (l=0; l<sz[1]; l++) {
                    memcpy(b + blocknum*sz[0]*sz[1]*sz[2] + m*sz[0]*sz[1] + l*sz[0], x+(k+m)*n[0]*n[1]+(j+l)*n[0]+i, sz[0]*sizeof(double));
                }
            }
            blocknum ++;
        }
    }
}

同时，我想通过用索引变量 blockid 替换增量 blocknum 来使其并行。

/* n stands for the size of an image, sz stands for the patch size to extract */
int blockid3, blockid2, blockid;
for (k=0; k<=n[2]-sz[2]; k+=1) {
    blockid3 = k * (n[1]-sz[1]+1) * (n[0]-sz[0]+1);
#pragma omp parallel for
    for (j=0; j<=n[1]-sz[1]; j+=1) {
        blockid2 = j * (n[0]-sz[0]+1);
        for (i=0; i<=n[0]-sz[0]; i+=1) {
            blockid = i + blockid2 + blockid3;

            /* copy single block */
            for (m=0; m<sz[2]; m++) {
                for (l=0; l<sz[1]; l++) {
                    memcpy(b + blockid*sz[0]*sz[1]*sz[2] + m*sz[0]*sz[1] + l*sz[0], x+(k+m)*n[0]*n[1]+(j+l)*n[0]+i, sz[0]*sizeof(double));
                }
            }

        }
    }
}

然后运行导致致命的段错误。我不知道为什么(根据堆栈跟踪，它似乎与安全线程有关)。因为我认为并行线程不应该访问同一个地址一次。我是否应该设置变量的某些属性，即静态或共享或私有(private)？这是堆栈跟踪:

Stack Trace (from fault):
[  0] 0x00007f9bcaa695de    /usr/local/MATLAB/R2011b/bin/glnxa64/libmwfl.so+00210398 _ZN2fl4diag15stacktrace_base7capt
ureERKNS0_14thread_contextEm+000158
[  1] 0x00007f9bcaa6b62d    /usr/local/MATLAB/R2011b/bin/glnxa64/libmwfl.so+00218669
[  2] 0x00007f9bcaa6b8f5    /usr/local/MATLAB/R2011b/bin/glnxa64/libmwfl.so+00219381 _ZN2fl4diag13terminate_logEPKcRKN
S0_14thread_contextEb+000165
[  3] 0x00007f9bc9a714f5   /usr/local/MATLAB/R2011b/bin/glnxa64/libmwmcr.so+00447733 _ZN2fl4diag13terminate_logEPKcPK8
ucontextb+000085
[  4] 0x00007f9bc9a6e5b4   /usr/local/MATLAB/R2011b/bin/glnxa64/libmwmcr.so+00435636
[  5] 0x00007f9bc9a6f333   /usr/local/MATLAB/R2011b/bin/glnxa64/libmwmcr.so+00439091
[  6] 0x00007f9bc9a6f4c7   /usr/local/MATLAB/R2011b/bin/glnxa64/libmwmcr.so+00439495
[  7] 0x00007f9bc9a7085f   /usr/local/MATLAB/R2011b/bin/glnxa64/libmwmcr.so+00444511
[  8] 0x00007f9bc9a70a15   /usr/local/MATLAB/R2011b/bin/glnxa64/libmwmcr.so+00444949
[  9] 0x00007f9bc89f0cb0              /lib/x86_64-linux-gnu/libpthread.so.0+00064688
[ 10] 0x00007f9bc876cb8e                    /lib/x86_64-linux-gnu/libc.so.6+01346446
[ 11] 0x00007f9b88238bb8 /home/peiyun/schmax3.0/test_im2col/mex_im2colstep.mexa64+00003000
[ 12] 0x00007f9bcb004eea    /usr/lib/gcc/x86_64-linux-gnu/4.6.3//libgomp.so+00032490
[ 13] 0x00007f9bc89e8e9a              /lib/x86_64-linux-gnu/libpthread.so.0+00032410
[ 14] 0x00007f9bc87164bd                    /lib/x86_64-linux-gnu/libc.so.6+00992445 clone+000109

顺便说一句，如果他们写入不同的地址，在 omp for 循环中是否存在关于 memcpy 的竞争条件？

Answer 1

您的代码中存在多个数据竞争，即:

/* n stands for the size of an image, sz stands for the patch size to extract */
int blockid3, blockid2, blockid;
for (k=0; k<=n[2]-sz[2]; k+=1) {
    blockid3 = k * (n[1]-sz[1]+1) * (n[0]-sz[0]+1);
#pragma omp parallel for
    for (j=0; j<=n[1]-sz[1]; j+=1) {
        blockid2 = j * (n[0]-sz[0]+1);          // <--- here
        for (i=0; i<=n[0]-sz[0]; i+=1) {        // <--- here
            blockid = i + blockid2 + blockid3;  // <--- here

            /* copy single block */
            for (m=0; m<sz[2]; m++) {           // <--- here
                for (l=0; l<sz[1]; l++) {       // <--- and here
                    memcpy(b + blockid*sz[0]*sz[1]*sz[2] + m*sz[0]*sz[1] + l*sz[0], x+(k+m)*n[0]*n[1]+(j+l)*n[0]+i, sz[0]*sizeof(double));
                }
            }

        }
    }
}

根据 OpenMP blockid2、i、blockid、m 和 l 都是隐式共享的，这不是你想要的。您应该将它们设置为 private，或者最好在并行区域内声明它们，从而将它们隐式设置为私有(private):

#pragma omp parallel for private(i,m,l,blockid,blockid2)
...

或

int blockid3;
for (k=0; k<=n[2]-sz[2]; k+=1) {
    blockid3 = k * (n[1]-sz[1]+1) * (n[0]-sz[0]+1);
#pragma omp parallel for
    for (j=0; j<=n[1]-sz[1]; j+=1) {
        int blockid2 = j * (n[0]-sz[0]+1);
        for (int i=0; i<=n[0]-sz[0]; i+=1) {
            int blockid = i + blockid2 + blockid3;

            /* copy single block */
            for (int m=0; m<sz[2]; m++) {
                for (int l=0; l<sz[1]; l++) {
                    memcpy(b + blockid*sz[0]*sz[1]*sz[2] + m*sz[0]*sz[1] + l*sz[0], x+(k+m)*n[0]*n[1]+(j+l)*n[0]+i, sz[0]*sizeof(double));
                }
            }

        }
    }
}

后者需要符合 C99 的编译器(因为循环变量的声明方式)。您的 GCC 4.6.3 需要 -std=c99 选项以启用 C99 合规性。如果没有这样的编译器可用(是否还有通用的非 C99 编译器？)，您应该添加 private(i,l,m) 子句。您可能还想将并行化移到最外层循环，以尽量减少 OpenMP 开销。