MPI_Type_create_subarray и MPI_Gather

Ответ 1

Итак, это немного более тонко и требует некоторого понимания того, как коллектив Gather ставит сложные типы.

Если вы посмотрите на примеры MPI_Gather, они состоят из 1-го массива, и довольно легко понять, что должно произойти; вы получаете (скажем) 10 ints от каждого процесса, а Gather достаточно умен, чтобы поставить 10 ints из ранга 0 в начале, 10 из ранга 1 в позиции 10-19 в массиве и т.д.

Более сложные макеты, подобные этому, немного сложнее. Во-первых, расположение данных с точки зрения отправителя отличается от расположения данных от приемников. В точке отправителя вы начинаете с элемента массива [1][2], переходите к [1][5] (в массиве размером 7x7), затем переходите к элементам массива [2][3] - [2][5] и т.д. Есть блоки CHUNK_ROWS данных, каждая из которых разделена на две строки.

Теперь рассмотрим, как приемник должен их получить. Скажем, он получает данные ранга 0. Он получит это в элементы массива [0][0]-[0][4] - пока что так хорошо; но затем он получит следующий блок данных в [1][0]-[1][4] в массиве размером 10x10. Это прыжок через 5 элементов. Макет в памяти отличается. Таким образом, получатель должен будет получить доступ к другому типу Subarray, после чего отправители отправляют его, потому что макет памяти отличается.

Так что, хотя вы можете отправлять что-то, что выглядит так:

  sizes[0] = CHUNK_ROWS+2;
  sizes[1] = CHUNK_COLUMNS+2;
  subsizes[0] = CHUNK_ROWS;
  subsizes[1] = CHUNK_COLUMNS;
  starts[0] = 1;
  starts[1] = 1;
  MPI_Type_create_subarray(ndims,sizes,subsizes,starts,MPI_ORDER_C,MPI_INT,&sendsubarray);
  MPI_Type_commit(&sendsubarray);

вы получите то, что выглядит следующим образом:

  sizes[0]    = ROWS; sizes[1] = COLUMNS;
  subsizes[0] = CHUNK_ROWS; subsizes[1] = CHUNK_COLUMNS;
  starts[0]   = 0; starts[1] = 0;
  MPI_Type_create_subarray(ndims,sizes,subsizes,starts,MPI_ORDER_C,MPI_INT,&recvsubarray);
  MPI_Type_commit(&recvsubarray);

Относительно важно заметить разницу в массиве sizes.

Теперь мы приближаемся. Обратите внимание, что ваша строка MPI_Gather изменяется примерно так:

MPI_Gather(DEBUG_CH[0],1,sendsubarray,recvptr,1,recvsubarray,0,MPI_COMM_WORLD);

Было несколько вещей, которые не работали над предыдущей версией, MPI_Gather(DEBUG_CH[0],1,subarray,ddd[0],CHUNK_ROWS*CHUNK_COLUMNS,MPI_INT,0,MPI_COMM_WORLD); - во-первых, обратите внимание, что вы ссылаетесь на ddd[0], но для каждого ранга, кроме ранга 0, ddd=NULL, и поэтому это не удастся. Итак, создайте новую переменную say recvptr, а в ранге нуль установите ее в ddd[0]. (Неважно, где другие процессы считают это, поскольку они не получают.) Кроме того, я думаю, вы не хотите получать CHUNK_ROWS*CHUNK_COLUMS MPI_INTs, потому что это будет смещать их в памяти, и я понимаю, что вы хотите, чтобы они были проложены так же, как и на подчиненных задачах, но в большем массиве.

Итак, теперь мы куда-то попадаем, но выше все еще не получится по интересной причине. Для примеров из 1-го массива достаточно просто выяснить, куда идут данные n-й ранги. То, как оно рассчитывается, - это найти объем полученных данных и начать следующий элемент сразу после этого. Но это не сработает. "Только после" конец нулевого ранга данных не там, где должны начинаться данные ранга один ([0][5]), но вместо этого [4][5] - элемент после последнего элемента в подвале ранга 0. Здесь данные, которые вы получаете из разных рангов, перекрываются! Поэтому нам придется возиться с экстентами типов данных и вручную указывать, где начинаются данные ранжирования. Вторая - легкая часть; вы используете функцию MPI_Gatherv, когда вам нужно вручную указать объем данных от каждого процессора или туда, где он идет. Первая - более сложная часть.

MPI позволяет указать нижнюю и верхнюю границы заданного типа данных - где, учитывая кусочек памяти, первый бит данных для этого типа будет идти, а где он "заканчивается", что здесь означает только то, что следующий может начаться. (Данные могут распространяться за верхнюю границу типа, что, по моему мнению, заставляет вводить эти имена в заблуждение, но это способ вещей.) Вы можете указать, что это все, что вам нравится, что делает его удобным для вас; так как мы будем иметь дело с элементами в массиве int, допустим размер нашего MPI_INT одного типа.

  MPI_Type_create_resized(recvsubarray, 0, 1*sizeof(int), &resizedrevsubarray);
  MPI_Type_commit(&resizedrecvsubarray);

(Обратите внимание, что мы должны сделать это только для полученного типа, из типа отправки, поскольку мы отправляем только один из них, это не имеет значения).

Теперь мы будем использовать gatherv, чтобы указать, где начинается каждый элемент - в единицах "размера" этого нового измененного типа, который равен всего 1 целому числу. Поэтому, если мы хотим, чтобы что-то попало в большой массив в [0][5], смещение с начала большого массива равно 5; если мы хотим, чтобы он находился там в положении [5][5], смещение равно 55.

Наконец, обратите внимание, что собирательные и рассеивающие коллективы предполагают, что участвует даже "мастер". Легче всего заставить эту работу работать, даже если у мастера есть своя часть глобального массива.

Таким образом, для меня работает следующее:

#include <mpi.h>
#include <iostream>
#include <cstdlib>
using namespace std;

#define ROWS 10
#define COLUMNS 10
#define CHUNK_ROWS 5
#define CHUNK_COLUMNS 5
#define TAG 0

int** alloca_matrice(int righe, int colonne)
{
    int** matrice=NULL;
    int i;

    matrice = (int **)malloc(righe * sizeof(int*));

    if(matrice != NULL){
        matrice[0] = (int *)malloc(righe*colonne*sizeof(int));
        if(matrice[0]!=NULL)
            for(i=1; i<righe; i++)
                matrice[i] = matrice[0]+i*colonne;
        else{
            free(matrice);
            matrice = NULL;
        }
    }
    else{
        matrice = NULL;
    }
    return matrice;

}

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

    int my_id, numprocs,length,i,j;
    int ndims, sizes[2],subsizes[2],starts[2];
    int** DEBUG_CH=NULL;
    int** ddd=NULL;
    int *recvptr=NULL;
    char name[BUFSIZ];
    MPI_Datatype sendsubarray;
    MPI_Datatype recvsubarray;
    MPI_Datatype resizedrecvsubarray;
    //MPI_Status status;
    MPI_Init(&argc, &argv) ;    
    MPI_Comm_rank(MPI_COMM_WORLD, &my_id) ;
    MPI_Comm_size(MPI_COMM_WORLD, &numprocs) ;  // Ottiene quanti processi sono attivi
    if (numprocs != 4) {
        MPI_Abort(MPI_COMM_WORLD,1);
    }
    MPI_Get_processor_name(name, &length);    

    //creo una sottomatrice ripulita dalle ghost cells
    ndims=2;
    sizes[0] = CHUNK_ROWS+2;
    sizes[1] = CHUNK_COLUMNS+2;
    subsizes[0] = CHUNK_ROWS;
    subsizes[1] = CHUNK_COLUMNS;
    starts[0] = 1;
    starts[1] = 1;
    MPI_Type_create_subarray(ndims,sizes,subsizes,starts,MPI_ORDER_C,MPI_INT,&sendsubarray);
    MPI_Type_commit(&sendsubarray);

    DEBUG_CH = alloca_matrice(CHUNK_ROWS+2,CHUNK_COLUMNS+2);
    for(i=0; i<CHUNK_ROWS+2; i++){
        for(j=0; j<CHUNK_COLUMNS+2; j++){
            if(i==0 || i==CHUNK_ROWS+1 || j==0 || j==CHUNK_COLUMNS+1)
                DEBUG_CH[i][j] = 5;
            else
                DEBUG_CH[i][j] = my_id;
        }
    }

    recvptr=DEBUG_CH[0];
    if(my_id==0){
        ddd = alloca_matrice(ROWS,COLUMNS);
        sizes[0]    = ROWS; sizes[1] = COLUMNS;
        subsizes[0] = CHUNK_ROWS; subsizes[1] = CHUNK_COLUMNS;
        starts[0]   = 0; starts[1] = 0;
        MPI_Type_create_subarray(2,sizes,subsizes,starts,MPI_ORDER_C,MPI_INT,&recvsubarray);
        MPI_Type_commit(&recvsubarray);
        MPI_Type_create_resized(recvsubarray, 0, 1*sizeof(int), &resizedrecvsubarray);
        MPI_Type_commit(&resizedrecvsubarray);
        recvptr = ddd[0];
    }

    int counts[5]={1,1,1,1};
    int disps[5] ={0,5,50,55};
    MPI_Gatherv(DEBUG_CH[0],1,sendsubarray,recvptr,counts,disps,resizedrecvsubarray,0,MPI_COMM_WORLD);
    if(!my_id){
        for(i=0; i<ROWS; i++){
            for(j=0; j<COLUMNS; j++){
                printf("%d ",ddd[i][j]);
            }
            printf("\n");
        }
    }

    if(my_id == 0) {
        MPI_Type_free(&resizedrecvsubarray);
        MPI_Type_free(&recvsubarray);
        free(ddd[0]);
        free(ddd);
    } else {
        MPI_Type_free(&sendsubarray);
        free(DEBUG_CH[0]);
        free(DEBUG_CH);
    }

    MPI_Finalize();                             // Chiusura di MPI.
    return 0;
}