MPI_Scatter

いちごパック > MPIの解説 > MPI_Scatter

インターフェース

#include <mpi.h>
int MPI_Scatter(const void* psenddata, int sendcount, MPI_Datatype sendtype,
                void* precvbuf, int recvcount, MPI_Datatype recvtype,
                int root, MPI_Comm comm);
rootで指定されたプロセスのデータを、root自身を含むすべてのプロセスに送信します。 この関数は通信完了まで待機します。 この関数を使う場合は、同じ受信型を指定してすべてのプロセスで呼び出す必要があります。
成功した場合の戻り値はMPI_SUCCESSです。失敗した場合はそれ以外のエラー値が返されます。
引数内容
psenddata送信データの先頭へのポインタです。rootで指定したプロセス以外では無視されます。
sendcount送信データの要素数です。rootで指定したプロセス以外では無視されます。
sendtype送信データの型です。rootで指定したプロセス以外では無視されます。
precvbuf受信バッファの先頭へのポインタです。
recvcount受信データの要素数です。すべてのプロセスで、rootプロセスのsendcountと同じ値を指定する必要があります。
recvtype受信データの型です。すべてのプロセスで、rootプロセスのsendtypeと同じ値を指定する必要があります。
root送信を行うプロセスのランク番号です。
comm送受信を行うプロセスの所属プロセスグループを指定します。MPI_COMM_WORLDを指定するとすべてのプロセスが所属するプロセスグループになります。

サンプルプログラム

ランク番号0を持つプロセスが要素数13のベクトルを生成し、 それらをMPI_Scatter()により全プロセスに送信した後、 各プロセスで別のベクトルを生成し、 それらをMPI_Gather()によりランク番号0を持つプロセスに送信し、 最後にランク番号0を持つプロセスが受け取ったベクトルを表示するプログラムです。

#include <mpi.h>
#include <stdint.h>
#include <iostream>
#include <vector>

#define RANK_MASTER 0

#define ABORT_INIT 1
#define ABORT_ICHIGO 15

#define VECTOR_SIZE 13

void ichigo_init(std::vector<int32_t>& ichigobuf)
{
    int ichigoloop;

    ichigobuf.resize(VECTOR_SIZE);
    for (ichigoloop = 0; ichigoloop < VECTOR_SIZE; ichigoloop++)
        ichigobuf[ichigoloop] = -1;

    std::cout << "ichigosample: data = ";
    for (ichigoloop = 0; ichigoloop < VECTOR_SIZE; ichigoloop++)
        std::cout << ichigobuf[ichigoloop] << " ";
    std::cout << std::endl;
}

void ichigo_finalize(const int32_t* pichigobuf)
{
    int ichigoloop;

    std::cout << "ichigosample: data = ";
    for (ichigoloop = 0; ichigoloop < VECTOR_SIZE; ichigoloop++)
        std::cout << pichigobuf[ichigoloop] << " ";
    std::cout << std::endl;
}

void process(
    int num_processes, int process_rank,
    int32_t* pichigobuf)
{
    int32_t* pmasterbuf = NULL;
    std::vector<int32_t> localbuf;
    int commcount = VECTOR_SIZE / num_processes;
    int remainder = VECTOR_SIZE - (commcount * num_processes);
    int ichigoloop;
    int ret;

    if (process_rank == RANK_MASTER) {
        for (ichigoloop = 0; ichigoloop < remainder; ichigoloop++)
            pichigobuf[ichigoloop] = ichigoloop;
        pmasterbuf = pichigobuf + remainder;
    }

    localbuf.resize(commcount);
    ret = MPI_Scatter(
        pmasterbuf, commcount, MPI_INT32_T,
        &localbuf[0], commcount, MPI_INT32_T,
        RANK_MASTER, MPI_COMM_WORLD );
    if (ret != MPI_SUCCESS) {
        // error
        MPI_Abort(MPI_COMM_WORLD, ABORT_ICHIGO);
        return;
    }

    for (ichigoloop = 0; ichigoloop < commcount; ichigoloop++)
        localbuf[ichigoloop] = (1+process_rank) * 100 + ichigoloop;

    ret = MPI_Gather(
        &localbuf[0], commcount, MPI_INT32_T,
        pmasterbuf, commcount, MPI_INT32_T,
        RANK_MASTER, MPI_COMM_WORLD );
    if (ret != MPI_SUCCESS) {
        // error
        MPI_Abort(MPI_COMM_WORLD, ABORT_ICHIGO);
        return;
    }
}


int main( int argc, char** argv )
{
    int num_processes = 1;
    int process_rank = 0;
    int ret;
    std::vector<int32_t> ichigobuf;
    int32_t* pichigobuf = NULL;

    ret = MPI_Init(&argc, &argv);
    if (ret == MPI_SUCCESS)
    ret = MPI_Comm_size(MPI_COMM_WORLD, &num_processes);
    if (ret == MPI_SUCCESS)
        ret = MPI_Comm_rank(MPI_COMM_WORLD, &process_rank);
    if (ret != MPI_SUCCESS) {
        // error
        MPI_Abort(MPI_COMM_WORLD, ABORT_INIT);
    }

    if (process_rank == RANK_MASTER) {
        ichigo_init(ichigobuf);
        pichigobuf = &ichigobuf[0];
    }
    process(num_processes, process_rank, pichigobuf);
    if (process_rank == RANK_MASTER) {
        ichigo_finalize(pichigobuf);
    }

    MPI_Finalize();
    return 0;
}


4つのプロセスで動作させたときの出力結果の1例は次の通りです。
ichigosample: data = -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
ichigosample: data = 0 100 101 102 200 201 202 300 301 302 400 401 402 

関連ページ

  • MPIの解説 目次
  • MPI_Datatype定数