创建 3D 网格时嵌套 for() 循环上的 OpenMP

Question

我正在尝试使用 openmp 并行化以下函数来创建 3D 网格。类 Point、Triangle 和 Box 已被声明和结构化，但代码未粘贴到此处，因为包含在问题中的代码很长（大约 110 行）。如果您需要它并允许我粘贴它，那么我可以编辑这个问题。

我试过在最外层的 for() 循环中使用 #pragma omp parallel for，当我运行时程序崩溃了。然后，我尝试在下面实现以下 #pragma 代码，程序花了更长的时间来运行它（与串行编码相比）。我在实现它时提到了这个 link here 。实际上首先，是否可以使用 openmp 并行化下面的代码？如果是，那么我真的需要你的帮助。

int main()
{
    vector<Box> boxes;
    Triangle dummy1;

    float dx = 81.0, dy = 121.0, dz = 98.0, delta=1.0;
    long int nx = 5*ceil(dx/(5.0*delta));
    long int ny = 5*ceil(dy/(5.0*delta));
    long int nz = 5*ceil(dz/(5.0*delta));
    long int Box_ID=1, Pt_ID=1;
    float x_val=0.0, y_val=0.0, z_val=0.0;
    float x0=-42.0f, y0=-3.0f, z0=-52.0f;

    long int i, j, k;
    for(i=0; i<nz+2; i++)
    {
        z_val=i*delta + z0;
        for(j=0; j<ny+2; j++)
        {
            y_val=j*delta + y0;
            #pragma omp parallel
            {
                vector<Box> box_private;
                #pragma omp for nowait schedule(static)
                for(k=0; k<nx+2; k++)
                {
                    x_val=k*delta + x0;
                    Point x1(x_val, y_val, z_val, Pt_ID);
                    Point x2(x_val+delta, y_val, z_val, Pt_ID+1);
                    Point x3(x_val, y_val+delta, z_val, Pt_ID+2);
                    Point x4(x_val+delta, y_val+delta, z_val, Pt_ID+3);
                    Point x5(x_val, y_val, z_val+delta, Pt_ID+4);
                    Point x6(x_val+delta, y_val, z_val+delta, Pt_ID+5);
                    Point x7(x_val, y_val+delta, z_val+delta, Pt_ID+6);
                    Point x8(x_val+delta, y_val+delta, z_val+delta, Pt_ID+7);
                    box_private.push_back(Box(x1,x2,x3,x4,x5,x6,x7,x8,dummy1,Box_ID));
                    Box_ID++;
                    Pt_ID++;
                }
                #pragma omp for schedule(static) ordered
                for(int i=0; i<omp_get_num_threads(); i++)
                {
                    #pragma omp ordered
                    boxes.insert(boxes.end(), box_private.begin(), box_private.end());
                }
            }
        }
    }
}

当我执行下面的代码而不是上面的代码时程序崩溃。

    #pragma omp parallel for
    for(i=0; i<nz+2; i++)
    {
        z_val=i*delta + z0;
        for(j=0; j<ny+2; j++)
        {
            y_val=j*delta + y0;
            for(k=0; k<nx+2; k++)
            {
                x_val=k*delta + x0;
                /* Point x1 to x8...*/
                boxes.push_back(Box(x1,x2,x3,x4,x5,x6,x7,x8,dummy1,Box_ID));
                Box_ID++;
                Pt_ID++;
            }
        }
    }

Answer 1

您的第一个内部循环并行代码变体运行缓慢，因为与实际代码相比，管理线程和合并向量的开销花费了太多时间运行。速度变慢的另一个原因是 shared 变量的广泛使用。虽然 OpenMP 同步对它们的访问，但需要花费大量时间。通常你应该尽量减少共享变量的使用。我相信好的风格是为并行部分指定 default(none) 并明确指定要共享的变量。

为了改善实际代码和管理代码之间的比例，您应该使并行区域尽可能大，并尽可能少地使用线程之间的同步。所以你应该并行你的最外层循环以获得最佳效果。

然而，在您的第二种方法中，您忘记了同步。您不能 insert 进入来自不同线程的相同 boxes 向量。您可以以完全相同的方式应用您在第一种方法中使用的同步机制。

Answer 2

您的代码存在的问题是您不了解共享变量和私有变量之间的区别，因此存在多个竞争条件。并行部分内部定义的变量是私有的，外部定义的变量是共享的。由于您在外部定义了所有内容（box_private 除外），所有内容都是共享的（box_private 和 k 除外，OpenMP 将其设为私有）。

但由于 Box_ID++ 和 Pt_ID++，即使将适当的变量设为私有也无法解决您的问题。您可以使用 atomic 修复它们，但这是不必要的且效率低下。如果您定义 Box_ID = ny*nx*i + nx*j + k;（与 Pt_ID 相同），那么您的代码应该可以工作。

for(long i=0; i<nz+2; i++) {
    float z_val=i*delta + z0;
    for(long j=0; j<ny+2; j++) {
        float y_val=j*delta + y0;
        #pragma omp parallel
        {
            vector<Box> box_private;
            #pragma omp for nowait schedule(static)
            for(long k=0; k<nx+2; k++) {
                float x_val=k*delta + x0;
                long Box_ID = ny*nx*i + nx*j + k;
                long Pt_ID  = ny*nx*i + nx*j + k;                    
                Point x1(x_val, y_val, z_val, Pt_ID);
                // Point x2 - x8
                box_private.push_back(Box(x1,x2,x3,x4,x5,x6,x7,x8,dummy1,Box_ID));

但我认为你应该退一步想想你想做什么。如果您提前知道要填充多少个框，那么就没有理由使用 push_back 或私有向量。事实上你应该可以做到

#pragma omp parallel for schedule(static)
for(long i=0; i<nz+2; i++) {
    float z_val=i*delta + z0;
    for(long j=0; j<ny+2; j++) {
        float y_val=j*delta + y0;
        for(long k=0; k<nx+2; k++) {
            float x_val=k*delta + x0;
            long Box_ID = ny*nx*i + nx*j + k;
            long Pt_ID  = ny*nx*i + nx*j + k;               
            Point x1(x_val, y_val, z_val, Pt_ID);
            // Point x2 - x8
            boxes[ny*nx*i + nx*j +k] = Box(x1,x2,x3,x4,x5,x6,x7,x8,dummy1,Box_ID);
        }
    }
}

其中 boxes 是一个 C 数组，或者如果它是一个 std:vector 确保调整它的大小。

最后，如果将 z_val 和 y_val 移动到 k 和 x_val 的循环中，然后执行 #pragma omp parallel for schedule(static) collapse(3)，则可以折叠三个循环。但是你也可以像这样by hand

long x = nz + 2, y = ny + 2, z = nx + 2;
#pragma omp parallel for schedule(static)
for(long n=0; n<(x*y*z); n++) {
    long i = n/(y*z);
    long j = (n%(y*z))/z;
    long k = (n%(y*z))%z;
    z_val=i*delta + z0;
    y_val=j*delta + y0;
    x_val=k*delta + x0;
    Point x1(x_val, y_val, z_val, n+1);
    // Point x2 - x8
    boxes[n] = Box(x1,x2,x3,x4,x5,x6,x7,x8,dummy1,n+1);
}

那你真的需要Box_ID和Pt_ID吗？

创建 3D 网格时嵌套 for() 循环上的 OpenMP

OpenMP on nested for() loops in creating 3D Grid

c++

openmp