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;
} |
