c - 在 Open MPI 中处理图像 block

Question

我正在尝试并行化绘制 Mandelbrot 集的 C 程序。正如我在这张图片中所示，我将每个处理器的图像部分划分为大小相等的 block :

我尝试使用解决方案 here这几乎是在解决同样的问题。但是，我在输出中只得到了部分图像:

同样对于高分辨率的图像(例如 8000x8000 像素)，应用程序因段错误 11 而崩溃。这是我的代码:

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include "mpi.h"

// Main program
int main(int argc, char* argv[])
 {
    /* screen ( integer) coordinate */
    int iX,iY,i,j;
    const int iXmax = 8000; // default
    const int iYmax = 8000; // default

    /* world ( double) coordinate = parameter plane*/
    double Cx, Cy;
    const double CxMin = -2.5;
    const double CxMax = 1.5;
    const double CyMin = -2.0;
    const double CyMax = 2.0; 

    /* */
    double PixelWidth = (CxMax - CxMin)/iXmax;
    double PixelHeight = (CyMax - CyMin)/iYmax;


    int linePerProcess, remainingLines, processMinY,  processMaxY, lastProcessMaxY, result_offset;
    int my_rank, processors, iXmaxHalf;
    int startAlert = 1;
    int receivedAlert;
    unsigned char (*resultBuffer)[3] = NULL;
    unsigned char (*resultBufferTwo)[3] = NULL;
    unsigned char (*finalResultBuffer)[3] = NULL;

    MPI_Status stat;


    /* color component ( R or G or B) is coded from 0 to 255 */
    /* it is 24 bit color RGB file */
    const int MaxColorComponentValue = 255; 
    FILE * fp;
    char *filename = "Mandelbrot.ppm";
    char *comment = "# ";   /* comment should start with # */

    // RGB color array
    unsigned char color[3];

    /* Z = Zx + Zy*i;   Z0 = 0 */
    double Zx, Zy;
    double Zx2, Zy2; /* Zx2 = Zx*Zx;  Zy2 = Zy*Zy  */
    /*  */
    int Iteration;
    const int IterationMax = 2000; // default

    /* bail-out value , radius of circle ;  */
    const double EscapeRadius = 400;
    double ER2 = EscapeRadius * EscapeRadius;
    double startTime, endTime;

    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);
    MPI_Comm_size(MPI_COMM_WORLD, &processors);

    linePerProcess = iYmax / (processors/2);
    iXmaxHalf = iXmax / 2;

    if (my_rank % 2  == 0) {
        processMinY = (my_rank/2) * linePerProcess;
    } else {
        processMinY = ((my_rank - 1)/2) * linePerProcess;
    }

    processMaxY = processMinY + linePerProcess;


    int Rows = iYmax;          // Global array rows
    int Columns = iXmax; // Global array columns
    int sizes[2];                     // No of elements in each dimension of the whole array
    int subSizes[2];                  // No of elements in each dimension of the subarray
    int startCoords[2];               // Starting coordinates of each subarray
    MPI_Datatype recvBlock, recvMagicBlock;

        // Create a subarray (a rectangular block) datatype from a regular, 2d array
    sizes[0] = Rows;
    sizes[1] = Columns;
    subSizes[0] = linePerProcess;
    subSizes[1] = iXmaxHalf;
    startCoords[0] = 0;
    startCoords[1] = 0;

    MPI_Type_create_subarray(2, sizes, subSizes, startCoords, MPI_ORDER_C, MPI_UNSIGNED_CHAR, &recvBlock);

    MPI_Type_create_resized(recvBlock, 0, iXmaxHalf * sizeof(color), &recvMagicBlock);

    MPI_Type_commit(&recvMagicBlock);

    if (my_rank == 0) {


        // startTime = MPI_Wtime();
        // for(i=1; i<processors; i++){
        //  MPI_Send(&startAlert, 1, MPI_INT, i, 0, MPI_COMM_WORLD);
        // }
        // printf("rank; %d\n", my_rank);

        finalResultBuffer = malloc(iXmax * iYmax * sizeof(color));

        for(iY = processMinY; iY < processMaxY; iY++) {

            Cy = CyMin + (iY * PixelHeight);
            if (fabs(Cy) < (PixelHeight / 2))
            {
                Cy = 0.0; /* Main antenna */
            }

            for(iX = 0; iX < iXmaxHalf; iX++)
            {

                Cx = CxMin + (iX * PixelWidth);
                /* initial value of orbit = critical point Z= 0 */
                Zx = 0.0;
                Zy = 0.0;
                Zx2 = Zx * Zx;
                Zy2 = Zy * Zy;

            /* */
                for(Iteration = 0; Iteration < IterationMax && ((Zx2 + Zy2) < ER2); Iteration++)
                {
                    Zy = (2 * Zx * Zy) + Cy;
                    Zx = Zx2 - Zy2 + Cx;
                    Zx2 = Zx * Zx;
                    Zy2 = Zy * Zy;
                };

            /* compute  pixel color (24 bit = 3 bytes) */
                if (Iteration == IterationMax)
                {
                    // Point within the set. Mark it as black
                    color[0] = 0;
                    color[1] = 0;
                    color[2] = 0;
                }
                else 
                {
                    // Point outside the set. Mark it as white
                    double c = 3*log((double)Iteration)/log((double)(IterationMax) - 1.0);
                    if (c < 1)
                    {
                        color[0] = 0;
                        color[1] = 0;
                        color[2] = 255*c;
                    }
                    else if (c < 2)
                    {
                        color[0] = 0;
                        color[1] = 255*(c-1);
                        color[2] = 255;
                    }
                    else
                    {
                        color[0] = 255*(c-2);
                        color[1] = 255;
                        color[2] = 255;
                    }
                }

                finalResultBuffer[(iY*iXmaxHalf)+iX][0] = color[0];
                finalResultBuffer[(iY*iXmaxHalf)+iX][1] = color[1];
                finalResultBuffer[(iY*iXmaxHalf)+iX][2] = color[2];


            }
        }

        result_offset = 1;
        for(i=1; i<processors; i++){

            MPI_Recv(finalResultBuffer, 1, recvMagicBlock, i, 0, MPI_COMM_WORLD, &stat);
            result_offset += 1;
        }



    } else if ((my_rank % 2 == 0) && (my_rank != 0)) {

        // MPI_Recv(&receivedAlert, 1, MPI_INT, 0, 0, MPI_COMM_WORLD, &stat);
        // printf("rank; %d\n", my_rank);
        resultBuffer = malloc(linePerProcess * iXmaxHalf * sizeof(color));

        for(iY = processMinY; iY < processMaxY; iY++) {
            Cy = CyMin + (iY * PixelHeight);
            if (fabs(Cy) < (PixelHeight / 2))
            {
                Cy = 0.0; /* Main antenna */
            }
            for(iX = 0; iX < iXmaxHalf; iX++)
            {

                Cx = CxMin + (iX * PixelWidth);
                /* initial value of orbit = critical point Z= 0 */
                Zx = 0.0;
                Zy = 0.0;
                Zx2 = Zx * Zx;
                Zy2 = Zy * Zy;

            /* */
                for(Iteration = 0; Iteration < IterationMax && ((Zx2 + Zy2) < ER2); Iteration++)
                {
                    Zy = (2 * Zx * Zy) + Cy;
                    Zx = Zx2 - Zy2 + Cx;
                    Zx2 = Zx * Zx;
                    Zy2 = Zy * Zy;
                };

            /* compute  pixel color (24 bit = 3 bytes) */
                if (Iteration == IterationMax)
                {
                    // Point within the set. Mark it as black
                    color[0] = 0;
                    color[1] = 0;
                    color[2] = 0;
                }
                else 
                {
                    // Point outside the set. Mark it as white
                    double c = 3*log((double)Iteration)/log((double)(IterationMax) - 1.0);
                    if (c < 1)
                    {
                        color[0] = 0;
                        color[1] = 0;
                        color[2] = 255*c;
                    }
                    else if (c < 2)
                    {
                        color[0] = 0;
                        color[1] = 255*(c-1);
                        color[2] = 255;
                    }
                    else
                    {
                        color[0] = 255*(c-2);
                        color[1] = 255;
                        color[2] = 255;
                    }
                }

                resultBuffer[((iY-processMinY)*iXmaxHalf)+iX][0] = color[0];
                resultBuffer[((iY-processMinY)*iXmaxHalf)+iX][1] = color[1];
                resultBuffer[((iY-processMinY)*iXmaxHalf)+iX][2] = color[2];

            }
        }
        MPI_Send(resultBuffer, linePerProcess * iXmaxHalf, MPI_UNSIGNED_CHAR, 0, 0, MPI_COMM_WORLD);

        free(resultBuffer);

    } else {

        // MPI_Recv(&receivedAlert, 1, MPI_INT, 0, 0, MPI_COMM_WORLD, &stat);
        // printf("rank; %d\n", my_rank);
        resultBufferTwo = malloc(linePerProcess * iXmaxHalf * sizeof(color));

        for(iY = processMinY; iY < processMaxY; iY++) {
            Cy = CyMin + (iY * PixelHeight);
            if (fabs(Cy) < (PixelHeight / 2))
            {
                Cy = 0.0; /* Main antenna */
            }
            for(iX = iXmaxHalf; iX < iXmax; iX++)
            {

                Cx = CxMin + (iX * PixelWidth);
                /* initial value of orbit = critical point Z= 0 */
                Zx = 0.0;
                Zy = 0.0;
                Zx2 = Zx * Zx;
                Zy2 = Zy * Zy;

            /* */
                for(Iteration = 0; Iteration < IterationMax && ((Zx2 + Zy2) < ER2); Iteration++)
                {
                    Zy = (2 * Zx * Zy) + Cy;
                    Zx = Zx2 - Zy2 + Cx;
                    Zx2 = Zx * Zx;
                    Zy2 = Zy * Zy;
                };

            /* compute  pixel color (24 bit = 3 bytes) */
                if (Iteration == IterationMax)
                {
                    // Point within the set. Mark it as black
                    color[0] = 0;
                    color[1] = 0;
                    color[2] = 0;
                }
                else 
                {
                    // Point outside the set. Mark it as white
                    double c = 3*log((double)Iteration)/log((double)(IterationMax) - 1.0);
                    if (c < 1)
                    {
                        color[0] = 0;
                        color[1] = 0;
                        color[2] = 255*c;
                    }
                    else if (c < 2)
                    {
                        color[0] = 0;
                        color[1] = 255*(c-1);
                        color[2] = 255;
                    }
                    else
                    {
                        color[0] = 255*(c-2);
                        color[1] = 255;
                        color[2] = 255;
                    }
                }


                resultBufferTwo[((iY-processMinY)*iXmaxHalf)+(iX - iXmaxHalf)][0] = color[0];
                resultBufferTwo[((iY-processMinY)*iXmaxHalf)+(iX - iXmaxHalf)][1] = color[1];
                resultBufferTwo[((iY-processMinY)*iXmaxHalf)+(iX - iXmaxHalf)][2] = color[2];
                // printf("rank: %d - value: %u%u%u\n", my_rank,resultBufferTwo[((iY-processMinY)*iXmax)+iX][0],resultBufferTwo[((iY-processMinY)*iXmax)+iX][1],resultBufferTwo[((iY-processMinY)*iXmax)+iX][2]);

            }
        }
        MPI_Send(resultBufferTwo, iXmaxHalf * linePerProcess, MPI_UNSIGNED_CHAR, 0, 0, MPI_COMM_WORLD);

        free(resultBufferTwo);

    }

    if (my_rank == 0) {

        endTime = MPI_Wtime();
        printf("Process time (s): %lf\n", endTime - startTime);
        /*create new file,give it a name and open it in binary mode  */
        fp = fopen(filename, "wb"); /* b -  binary mode */

        /*write ASCII header to the file (PPM file format)*/
        fprintf(fp,"P6\n %s\n %d\n %d\n %d\n", comment, iXmax, iYmax, MaxColorComponentValue);
        for(iY = 0; iY < iYmax; iY++)
        {
            for(iX = 0; iX < iXmax; iX++)
                {
                    fwrite(finalResultBuffer[(iY*iXmax)+iX], 1, 3, fp);
                }
        }

        fclose(fp);

        free(finalResultBuffer);
    }

    MPI_Finalize();
    return 0;
 }

如果有人能帮助我理解我在这里做错了什么，我将不胜感激。

Answer 1

Q : what I am doing wrong here?

• 代码不计算复平面上的分形，只是复制常量 color[]
• my_rank == 0 处理收集其他人生成的结果(分形生成器是迭代的，因此每个 [X,iY] 点都有不同的运行时间-plane ) 和它自己在定向的地方存储一个常量值，从 color[] == [0,0,0] 的未初始化/未修改的相同常量值到整个 finalResultBuffer[][ ]
• 您没有像上面声明的那样在 P1-P6 之间拆分工作，但是您的代码拆分了三种情况的代码执行路径 - 1) my_rank == 0(“master”，它收集结果(它本身不计算 cmplex 迭代器的一部分，而是在指定的地方存储黑点)并写出一个文件)，2) my_rank % 2 == 0(所有非零偶数秩，其中没有一个执行单个步骤来实际计算复杂的分形迭代器，并且所有存储的黑点都在各处)，3) my_rank % 2 == 1(所有奇数行列，没有一个执行单个步骤来实际计算复杂的分形迭代器，但所有行列都在各处存储黑点)