I m trying to preform matrix multiplication using Blac s pdgemm. The exact subroutine for the matrix multiplication that I am using can be found here: http://www.netlib.org/scalapack/html/pblas_qref.html#PvGEMM

然而,我的代码“无法为校对当地数据分配记忆:ABORT”关于大鼠的呼声。 我的代码一把矩阵本身乘起来,它是一个固定的矩阵,因此所产生的矩阵是相同的方面。 这里是有关的法典。

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

#define gridSize 10

int main(int argc, char* argv[]) {
  int i,j,k, np, myid;
  int bufsize;
  double *buf;          /* Buffer for reading */

  MPI_Offset filesize;
  MPI_File myfile;    /* Shared file */
  MPI_Status status;  /* Status returned from read */

  /* Initialize MPI */
  MPI_Init(&argc, &argv);
  MPI_Comm_rank(MPI_COMM_WORLD, &myid);
  MPI_Comm_size(MPI_COMM_WORLD, &np);

  double *A = (double*) malloc(gridSize*gridSize*sizeof(double));

  /*MPI-IO Code removed for clarity including buf and bufsize assignments
  .
  .
  .
  .
  */
    //I use mpi-IO to read in a matrix from a file, each processor reads in a row and that row is store in the array called buf
    for (j = 0; j <bufsize;j++){
    A[myid*bufsize+j] = buf[j];
  }

  //BLACS portion
  int ictxt, nprow, npcol, myrow, mycol, nb;
  int info,itemp;
  int ZERO = 0, ONE = 1;
  nprow = 2; npcol = 2; nb = 2;

  Cblacs_pinfo(&myid,&np);
  Cblacs_get(-1,0,&ictxt);
  Cblacs_gridinit(&ictxt,"Row",nprow,npcol);
  Cblacs_gridinfo(ictxt,&nprow,&npcol,&myrow,&mycol);

  int M = gridSize;

  int descA[9], descx[9], descy[9];
  int mA = numroc_(&M, &nb, &myrow, &ZERO, &nprow);
  int nA = numroc_(&M, &nb, &mycol, &ZERO, &npcol);
  int nx = numroc_(&M, &nb, &myrow, &ZERO, &nprow);
  int my = numroc_(&M ,&nb, &myrow, &ZERO, &nprow);

  descinit_(descA,&M,&M,&nb,&nb,&ZERO,&ZERO,&ictxt,&mA,&info);
  descinit_(descx,&M,&M,&nb,&nb,&ZERO,&ZERO,&ictxt,&nx,&info);

  descinit_(descy,&M,&M,&nb,&nb,&ZERO,&ZERO,&ictxt,&my,&info);


  double* matrixA = (double*)malloc(mA*nA*sizeof(double));
  double* matrixB = (double*)malloc(mA*nA*sizeof(double));
  double* matrixC = (double*)calloc(mA*nA,sizeof(double));
  int sat,sut;


  for(i=0;i<mA;i++){
    for(j=0;j<nA;j++){
        sat = (myrow*nb)+i+(i/nb)*nb;
        sut = (mycol*nb)+j+(j/nb)*nb;
        matrixA[j*mA+i] = A[sat*M+sut];
        matrixB[j*mA+i] = A[sat*M+sut];
    }
  }

  double alpha = 1.0; double beta = 0.0;

  //call where seg fault happens
  pdgemm_("N","N",&M,&M,&M,&alpha,matrixA,&ONE,&ONE,descA,matrixB,&ONE,&ONE,descx,
      &beta,matrixC,&ONE,&ONE,descy);

  Cblacs_barrier(ictxt,"A");

  Cblacs_gridexit(0);

  /* Close the file */
  MPI_File_close(&myfile);

  if (myid==0) {
    printf("Done
");

  }
  MPI_Finalize();

  exit(0);

}

我在ScaLapacs方面有很多经验,但从我所看到的例子来看,我不敢肯定,为什么我正在找错,会得到任何帮助。