Массивное умножение матрицы MPI в C

Question

Я пытаюсь написать программу для выполнения умножения матрицы по матрице с использованием MPI. Предполагается, что программа выделяет только память для строки строки матрицы. Я получаю ошибку сегментации, и я не уверен, что вызывает ее. Я также не очень хорошо знаком с MPI. Вот мой код до сих пор:

/* process p-1 reads in the matrix n/p rows at a time and sends it to 
the other processes */ 
float** read_matrix(char *filename, int p, int id, int *n) 
{ 
    /* open file */ 
    FILE *fp = fopen(filename, "r"); 
    MPI_Status stat; 

    /* reads in the dimension of the matrix (n x n); if not square 
    matrix, quit */ 
    int dim1, dim2, d; 
    d = fscanf(fp, "%d", &dim1); 
    d = fscanf(fp, "%d", &dim2); 
    *n = dim1; 

    if(dim1 != dim2){ 
    printf("%s\n", "Error: matrix not square"); 
    exit(0); 
    } 
    if(dim1 % p != 0){ 
    printf("%s\n", "Error: matrix cannot be divided into processes"); 
    } 
int i; 

    /* allocate an array 'storage_mat' of floats of n x (n/p) in size */ 
    float **storage_mat = malloc(sizeof(float) * dim1/p); 

    for(i = 0; i < dim1/p; i++){ 
    storage_mat[i] = malloc(dim1 * sizeof(float)); 
    } 


    if(id == p-1) { 
    for(i=0; i< p; i++) { 
     int j, col; 
     /* read (n/p) rows of the matrix and fill in the array */ 
     for(j = 0; j < dim1/p; j++){ 
      for(col = 0; col < dim1; col++){ 
      d = fscanf(fp, "%f", &storage_mat[j][col]); 
      } 
      } 

     } 

     if(i < p-1){ 
      // mpi send the array storage_mat to process rank i 
      MPI_Send(&storage_mat, 1, MPI_FLOAT, i, i, MPI_COMM_WORLD); 
     } 

     } 
     else { 
     /* mpi receive the array storage_mat from process rank p-1 */ 
     MPI_Recv(&storage_mat, 1, MPI_FLOAT, p-1, 1, MPI_COMM_WORLD, &stat); 
    } 

    /* close file */ 
    fclose(fp); 

    int j; 
    for(j = 0; j < dim1/p; j++){ 
     for(i = 0; i < dim1; i++){ 
     printf("%f ", storage_mat[j][i]); 
     } 

} 

    return storage_mat; 

} 

/* process 0 writes out the matrix n/p rows at a time on behalf of all 
    other processes */ 
void write_matrix(char* filename, float** storage_mat, int id, int dim1, int p) 
{ 
    /* open file */ 
    FILE *fp = fopen(filename, "w"); 
    int j, i; 
    MPI_Status stat; 

    if(!id) { 
    for(i=0; i<p; i++) { 
     /* write (n/p) rows of the matrix from the array storage_mat */ 
     for(j = 0; j < dim1/p; j++){ 
     for(i = 0; i < dim1; i++){ 
     fprintf(fp, "%f ", storage_mat[j][i]); 
     } 
    } 
     if(i < p-1) { 
    /* mpi receive the array storage_mat from process rank i */ 
     MPI_Recv(&storage_mat, 1, MPI_FLOAT, i, 1, MPI_COMM_WORLD, &stat); 
     } 

    } 

} else { 
    /* mpi send the array storage_mat to process rank 0 */ 
     MPI_Send(&storage_mat, 1, MPI_FLOAT, 0, i, MPI_COMM_WORLD); 
     } 

    /* close file */ 
     fclose(fp); 


} 

main(int argc, char** argv) 
{ 
    /* initialize mpi, and find out the rank (id) and the total number 
    of processes (p) */ 
    int id, p, n; 
    MPI_Init(&argc, &argv); //initialize mpi 
    MPI_Comm_rank(MPI_COMM_WORLD, &id); // get the rank 
    MPI_Comm_size(MPI_COMM_WORLD, &p); // get num processes 


    /* parse command line arguments: fileA (the name of the file that 
    contains the matrix A), fileB (the name of the file that contains 
    matrix B), and fileC (the name of the file to store matrix C) */ 
    char *fileA, *fileB, *fileC; 

    if(argc != 4){ 
    printf("%s\n", "Error: incorrect number of files"); 
    exit(0); 
    } 
    fileA = argv[1]; 
    fileB = argv[2]; 
    fileC = argv[3]; 

    //check if files ae valid 
    FILE *fp1 = fopen(fileA, "r"); 
    FILE *fp2 = fopen(fileB, "r"); 
    FILE *fp3 = fopen(fileC, "r"); 

if(fp1 == NULL || fp2 == NULL || fp3 == NULL){ 
    printf("%s\n", "Error: Not all files accessable"); 
    exit(0); 


} 

    float** storage_matA = read_matrix(fileA, p, id, &n); 

    float ** storage_matB = read_matrix(fileB, p, id, &n); 

    /* allocate space and intialize to zero for storage_matC as an array 
     of floats of n x (n/p) in size */ 
    float **storage_matC = malloc(sizeof(float) * n/p); 
int i; 
    for(i = 0; i < n/p; i++){ 
    storage_matC[i] = malloc(n * sizeof(float)); 
    } 

    // initialize to 0 
    int j; 
    for(j = 0; j < n/p; j++){ 
     for(i = 0; i < n; i++){ 
     storage_matC[j][i] = 0; 
     } 
    } 

    /* create the auxiliary array of pointers so that the elements in A, 
    B and C can be accessed using matA[i][j], etc. */ 
    int k; 
    for(i=0; i<p; i++) { 

    for(i = 0; i < n/p; i++){ 
     for(j = 0; j < n; j++){ 
      storage_matC[i][j] += storage_matA[i][j] * storage_matB[i][j]; 
     } 
    } 
    /* calculate the partial sum for matC given the row band of A and 
    B (see notes on row-wise matrix multiplication). */ 
    MPI_Status stat; 
    if(i < p-1) { 
     /* mpi send storage_matB to the next process (id+1)%p */ 
     MPI_Send(&storage_matB, 1, MPI_FLOAT, (id+1)%p, i, MPI_COMM_WORLD); 
     /* mpi receive storage_matB from the previous process */ 
     MPI_Recv(&storage_matB, 1, MPI_FLOAT, (id-1)%p, 1, MPI_COMM_WORLD, &stat); 
    } 
} 

    write_matrix(fileC, storage_matC, id, n , p); 

    /* reclaim matrices, finalize mpi */ 
    free(storage_matB); 
    free(storage_matC); 
    free(storage_matA); 
    MPI_Finalize(); 
} 

Answer 1

float **storage_mat = malloc(sizeof(float) * dim1/p); 

должен быть

float **storage_mat = malloc(sizeof(float *) * dim1/p);