MPI 奇偶排序如何工作？

只是在面试前查看了一些笔记，我很难理解如何奇偶排序 http://en.wikipedia.org/wiki/Odd%E2%80%93even_sort在并行架构中工作。

int MPI_OddEven_Sort(int n, double *a, int root, MPI_Comm comm)
{
    int rank, size, i, sorted_result;
    double *local_a;

// get rank and size of comm
    MPI_Comm_rank(comm, &rank); //&rank = address of rank
    MPI_Comm_size(comm, &size);

    local_a = (double *) calloc(n / size, sizeof(double));


// scatter the array a to local_a
    MPI_Scatter(a, n / size, MPI_DOUBLE, local_a, n / size, MPI_DOUBLE,
        root, comm);

// sort local_a
    merge_sort(n / size, local_a);

//odd-even part
    for (i = 0; i < size; i++) {

        if ((i + rank) % 2 == 0) {  // means i and rank have same nature
            if (rank < size - 1) {
                MPI_Compare(n / size, local_a, rank, rank + 1, comm);
            }
        } else if (rank > 0) {
            MPI_Compare(n / size, local_a, rank - 1, rank, comm);
        }

        MPI_Barrier(comm);
        // test if array is sorted
        MPI_Is_Sorted(n / size, local_a, root, comm, &sorted_result);

        // is sorted gives integer 0 or 1, if 0 => array is sorted
        if (sorted_result == 0) {
            break;
        }           // check for iterations
    }

// gather local_a to a
    MPI_Gather(local_a, n / size, MPI_DOUBLE, a, n / size, MPI_DOUBLE,
           root, comm)
    return MPI_SUCCESS;
}

是我为此函数编写的一些代码（不是今天也不是昨天！）。有人可以分解一下它是如何工作的吗？

我正在分散我的阵列a到每个处理器，它正在获取一个副本local_a（其大小为 n/size）

每个 local_a 上都调用合并排序。

这之后发生了什么？ （假设到目前为止我是正确的！）

这么多年之后，看到这些 PRAM 类型的排序网络再次出现，真是太有趣了。这些事物的并行计算的原始心智模型是作为“比较器”的微型处理器的大规模并行阵列，例如连接机 http://en.wikipedia.org/wiki/Connection_Machine- 那时候网络比 CPU/RAM 便宜。当然，这最终看起来与 80 年代中后期及以后的超级计算机非常不同，甚至比 90 年代后期的 x86 集群更加不同；但现在它们又开始流行起来配备 GPU http://http.developer.nvidia.com/GPUGems2/gpugems2_chapter46.html以及其他加速器，如果你眯着眼睛看的话，它们实际上看起来有点像未来的过去。

看起来你上面的东西更像是，它已经开始朝着假设处理器将在本地存储多个项目的方向发展，并且您可以通过在通信步骤之间对这些本地列表进行排序来充分利用处理器。

充实你的代码并稍微简化一下，我们有这样的东西：

#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>

int merge(double *ina, int lena, double *inb, int lenb, double *out) {
    int i,j;
    int outcount=0;

    for (i=0,j=0; i<lena; i++) {
        while ((inb[j] < ina[i]) && j < lenb) {
            out[outcount++] = inb[j++];
        }
        out[outcount++] = ina[i];
    }
    while (j<lenb)
        out[outcount++] = inb[j++];

    return 0;
}

int domerge_sort(double *a, int start, int end, double *b) {
    if ((end - start) <= 1) return 0;

    int mid = (end+start)/2;
    domerge_sort(a, start, mid, b);
    domerge_sort(a, mid,   end, b);
    merge(&(a[start]), mid-start, &(a[mid]), end-mid, &(b[start]));
    for (int i=start; i<end; i++)
        a[i] = b[i];

    return 0;
}

int merge_sort(int n, double *a) {
    double b[n];
    domerge_sort(a, 0, n, b);
    return 0;
}

void printstat(int rank, int iter, char *txt, double *la, int n) {
    printf("[%d] %s iter %d: <", rank, txt, iter);
    for (int j=0; j<n-1; j++)
        printf("%6.3lf,",la[j]);
    printf("%6.3lf>\n", la[n-1]);
}

void MPI_Pairwise_Exchange(int localn, double *locala, int sendrank, int recvrank,
                           MPI_Comm comm) {

    /*
     * the sending rank just sends the data and waits for the results;
     * the receiving rank receives it, sorts the combined data, and returns
     * the correct half of the data.
     */
    int rank;
    double remote[localn];
    double all[2*localn];
    const int mergetag = 1;
    const int sortedtag = 2;

    MPI_Comm_rank(comm, &rank);
    if (rank == sendrank) {
        MPI_Send(locala, localn, MPI_DOUBLE, recvrank, mergetag, MPI_COMM_WORLD);
        MPI_Recv(locala, localn, MPI_DOUBLE, recvrank, sortedtag, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
    } else {
        MPI_Recv(remote, localn, MPI_DOUBLE, sendrank, mergetag, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
        merge(locala, localn, remote, localn, all);

        int theirstart = 0, mystart = localn;
        if (sendrank > rank) {
            theirstart = localn;
            mystart = 0;
        }
        MPI_Send(&(all[theirstart]), localn, MPI_DOUBLE, sendrank, sortedtag, MPI_COMM_WORLD);
        for (int i=mystart; i<mystart+localn; i++)
            locala[i-mystart] = all[i];
    }
}

int MPI_OddEven_Sort(int n, double *a, int root, MPI_Comm comm)
{
    int rank, size, i;
    double *local_a;

// get rank and size of comm
    MPI_Comm_rank(comm, &rank); //&rank = address of rank
    MPI_Comm_size(comm, &size);

    local_a = (double *) calloc(n / size, sizeof(double));


// scatter the array a to local_a
    MPI_Scatter(a, n / size, MPI_DOUBLE, local_a, n / size, MPI_DOUBLE,
        root, comm);
// sort local_a
    merge_sort(n / size, local_a);

//odd-even part
    for (i = 1; i <= size; i++) {

        printstat(rank, i, "before", local_a, n/size);

        if ((i + rank) % 2 == 0) {  // means i and rank have same nature
            if (rank < size - 1) {
                MPI_Pairwise_Exchange(n / size, local_a, rank, rank + 1, comm);
            }
        } else if (rank > 0) {
            MPI_Pairwise_Exchange(n / size, local_a, rank - 1, rank, comm);
        }

    }

    printstat(rank, i-1, "after", local_a, n/size);

// gather local_a to a
    MPI_Gather(local_a, n / size, MPI_DOUBLE, a, n / size, MPI_DOUBLE,
           root, comm);

    if (rank == root)
        printstat(rank, i, " all done ", a, n);

    return MPI_SUCCESS;
}

int main(int argc, char **argv) {

    MPI_Init(&argc, &argv);

    int n = argc-1;
    double a[n];
    for (int i=0; i<n; i++)
        a[i] = atof(argv[i+1]);

    MPI_OddEven_Sort(n, a, 0, MPI_COMM_WORLD);

    MPI_Finalize();

    return 0;
}

因此，它的工作方式是，列表在处理器之间均匀分配（非均等分布也很容易处理，但需要大量额外的簿记工作，这对本次讨论没有多大帮助）。

我们首先对本地列表进行排序（O(n/P ln n/P)）。当然，没有理由它必须是合并排序，除了这里我们可以按照以下步骤重复使用该合并代码。然后我们进行 P 个邻居交换步骤，每个方向各一半。这里的模型是，有一个线性网络，我们可以在其中直接快速地与近邻进行通信，但也许根本无法与更远的邻居进行通信。

The 原始奇偶排序网络 http://en.wikipedia.org/wiki/Batcher_odd%E2%80%93even_mergesort是每个处理器都有一个密钥的情况，在这种情况下，通信很容易 - 您将您的项目与邻居进行比较，并在必要时进行交换（因此这基本上是并行冒泡排序）。在这种情况下，我们在进程对之间进行简单的并行排序 - 这里，每一对仅将所有数据发送到其中一个，该对合并已经本地排序的列表 O(N/P)，然后给出适当的一半的数据返回到另一个处理器。我把你的支票拿出来了；可以看出，它是在P个邻居交换中完成的。您当然可以将其添加回来，以防提前终止；然而，当一切完成时，所有处理器都必须达成一致，这需要类似全部减少 http://www.mpich.org/static/docs/v3.1/www3/MPI_Allreduce.html，这在一定程度上打破了原来的模型。

因此，每个链路的数据传输次数为 O(n)（每次发送和接收 n/P 项 P 次），并且每个处理器执行 (n/P ln n/P) + (2 n/P - 1)*P/ 2 = O(n/P ln n/P + N) 次比较；在这种情况下，还需要考虑分散和聚集，但一般来说，这种排序是在数据到位的情况下完成的。

运行上面的代码 - 为了清楚起见，使用相同的示例给出（输出重新排序以使其更易于阅读）：

$ mpirun -np 4 ./baudet-stevenson 43 54 63 28 79 81 32 47 84 17 25 49
[0] before iter 1: <43.000,54.000,63.000>
[1] before iter 1: <28.000,79.000,81.000>
[2] before iter 1: <32.000,47.000,84.000>
[3] before iter 1: <17.000,25.000,49.000>

[0] before iter 2: <43.000,54.000,63.000>
[1] before iter 2: <28.000,32.000,47.000>
[2] before iter 2: <79.000,81.000,84.000>
[3] before iter 2: <17.000,25.000,49.000>

[0] before iter 3: <28.000,32.000,43.000>
[1] before iter 3: <47.000,54.000,63.000>
[2] before iter 3: <17.000,25.000,49.000>
[3] before iter 3: <79.000,81.000,84.000>

[0] before iter 4: <28.000,32.000,43.000>
[1] before iter 4: <17.000,25.000,47.000>
[2] before iter 4: <49.000,54.000,63.000>
[3] before iter 4: <79.000,81.000,84.000>

[0] after iter 4: <17.000,25.000,28.000>
[1] after iter 4: <32.000,43.000,47.000>
[2] after iter 4: <49.000,54.000,63.000>
[3] after iter 4: <79.000,81.000,84.000>

[0]  all done  iter 5: <17.000,25.000,28.000,32.000,43.000,47.000,49.000,54.000,63.000,79.000,81.000,84.000>