1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115
| #include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <cblas.h>
#include <mpi.h>
#define AA(i,j) AA[(i)*M+(j)]
#define XX(i,j) XX[(i)*M+(j)]
typedef struct{
double* tab;
int nrow;
int ncol;
}Matrix;
typedef struct{
int ictxt;
int nprow;
int npcol;
int myrow;
int mycol;
int nb;
int info;
int itemp;
int myrank_mpi;
int nprocs_mpi;
}InfoParrall;
/*! External procedures prototypes from CBLACS: */
extern void Cblacs_pinfo(int *mypnum, int *nprocs);
extern void Cblacs_get(int ConTxt, int what, int *val);
extern void Cblacs_gridinit(int *ConTxt, char *order, int nprow, int npcol);
extern void Cblacs_gridinfo(int ConTxt, int *nprow, int *npcol,
int *myrow, int *mycol);
extern void Cblacs_gridexit(int ConTxt);
extern void Cblacs_exit(int NotDone);
double* multMM(Matrix A, Matrix B, InfoParrall inf);
double* multVV(Matrix A, Matrix B, InfoParrall inf);
double* multV(double* AA, double* X,int M, int nb,int mycol, int myrow, int ZERO,int ONE, int nprow, int npcol,int myrank_mpi,int info,int ictxt,int my,int Acol,int Arow,int Xrow,int Xcol);
double* multM(double* AA, double* X,int M, int nb,int mycol, int myrow, int ZERO,int ONE, int nprow, int npcol,int myrank_mpi,int info,int ictxt,int my,int Acol,int Arow,int Xrow,int Xcol);
double* multMM(Matrix A, Matrix B, InfoParrall inf){
return multM(A.tab,B.tab,5,inf.nb,inf.mycol,inf.myrow,0,1,inf.nprow,inf.npcol,inf.myrank_mpi,inf.info,inf.ictxt,A.nrow,A.ncol,A.nrow,B.nrow,B.ncol);
}
double* multM(double* AA, double* X,int M, int nb,int mycol, int myrow, int ZERO,int ONE, int nprow, int npcol,int myrank_mpi,int info,int ictxt,int my,int Acol,int Arow,int Xrow,int Xcol){
int descA[9],descx[9],descy[9];
descinit_(descA, &M, &M, &nb, &nb, &ZERO, &ZERO, &ictxt, &Acol, 0);
descinit_(descx, &M, &ONE, &nb, &ONE, &ZERO, &ZERO, &ictxt, &Xrow, 0);
descinit_(descy, &M, &ONE, &nb, &ONE, &ZERO, &ZERO, &ictxt, &my, 0);
double *y = (double*) malloc(M*M*sizeof(double));
double alpha = 1.0; double beta = 0.0;
pdgemm_("N","N",&M,&M,&Acol,&alpha,AA,Arow,Acol,descA,X,Xrow,Xcol,descx,&beta,y,Arow,Acol,descy);
return y;
}
int main(int argc, char **argv) {
/************ MPI ***************************/
int myrank_mpi, nprocs_mpi;
MPI_Init( &argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myrank_mpi);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs_mpi);
/************ BLACS ***************************/
InfoParrall inf;
int ZERO=0,ONE=1;
inf.nprow = 2; inf.npcol = 2; inf.nb =2;
Cblacs_pinfo( &(inf.myrank_mpi), &(inf.nprocs_mpi) ) ;
Cblacs_get( -1, 0, &(inf.ictxt) );
Cblacs_gridinit( &(inf.ictxt), "Row", inf.nprow, inf.npcol );
Cblacs_gridinfo( inf.ictxt, &(inf.nprow), &(inf.npcol), &(inf.myrow), &(inf.mycol) );
// Initialisation de matrix
int M=5;
int i,j;
Matrix AA;
AA.tab = (double*) malloc(M*M*sizeof(double));
for(i=0;i<M;i++){
for(j=0;j<M;j++){
AA.tab[i*M+j]=i;
}
}
//matrix X
Matrix X;
X.tab = (double*) malloc(M*M*sizeof(double));
for(i=0;i<M;i++){
for(j=0;j<M;j++){
X.tab[i*M+j]=1;
}
}
AA.nrow = numroc_( &M, &(inf.nb), &(inf.myrow), &ZERO, &(inf.nprow) );
AA.ncol = numroc_( &M, &(inf.nb), &(inf.mycol), &ZERO, &(inf.npcol) );
X.nrow = numroc_( &M, &(inf.nb), &(inf.myrow), &ZERO, &(inf.nprow) );
X.ncol=numroc_( &M, &(inf.nb), &(inf.mycol), &ZERO, &(inf.npcol) );
int my = numroc_( &M, &(inf.nb), &(inf.myrow), &ZERO, &(inf.nprow) );
//double *y= multM(AA.tab, X.tab,M, (inf.nb), (inf.mycol), (inf.myrow), ZERO, ONE, (inf.nprow), (inf.npcol), (inf.myrank_mpi),(inf.info),(inf.ictxt),my,AA.nrow,AA.ncol,X.nrow,X.ncol);
double* y=multMM(AA,X,inf);
Cblacs_barrier((inf.ictxt),"A");
for(i=0;i<my;i++){
printf("rank=%d %.2f \n", (inf.myrank_mpi),y[i]);
}
free(AA.tab);
free(X.tab);
Cblacs_gridexit( 0 );
MPI_Finalize();
return 0;
} |
Partager