ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
igamn2d_.c
Go to the documentation of this file.
1 #include "Bdef.h"
2 
3 #if (INTFACE == C_CALL)
4 void Cigamn2d(int ConTxt, char *scope, char *top, int m, int n, int *A,
5  int lda, int *rA, int *cA, int ldia, int rdest, int cdest)
6 #else
7 F_VOID_FUNC igamn2d_(int *ConTxt, F_CHAR scope, F_CHAR top, int *m, int *n,
8  int *A, int *lda, int *rA, int *cA, int *ldia,
9  int *rdest, int *cdest)
10 #endif
11 /*
12  * -- V1.1 BLACS routine --
13  * University of Tennessee, May 1, 1996
14  * Written by Clint Whaley.
15  *
16  * Purpose
17  * =======
18  * Combine amn operation for integer rectangular matrices.
19  *
20  * Arguments
21  * =========
22  *
23  * ConTxt (input) Ptr to int
24  * Index into MyConTxts00 (my contexts array).
25  *
26  * SCOPE (input) Ptr to char
27  * Limit the scope of the operation.
28  * = 'R' : Operation is performed by a process row.
29  * = 'C' : Operation is performed by a process column.
30  * = 'A' : Operation is performed by all processes in grid.
31  *
32  * TOP (input) Ptr to char
33  * Controls fashion in which messages flow within the operation.
34  *
35  * M (input) Ptr to int
36  * The number of rows of the matrix A. M >= 0.
37  *
38  * N (input) Ptr to int
39  * The number of columns of the matrix A. N >= 0.
40  *
41  * A (output) Ptr to integer two dimensional array
42  * The m by n matrix A. Fortran77 (column-major) storage
43  * assumed.
44  *
45  * LDA (input) Ptr to int
46  * The leading dimension of the array A. LDA >= M.
47  *
48  * RA (output) Integer Array, dimension (LDIA, N)
49  * Contains process row that the amn of each element
50  * of A was found on: i.e., rA(1,2) contains the process
51  * row that the amn of A(1,2) was found on.
52  * Values are left on process {rdest, cdest} only, others
53  * may be modified, but not left with interesting data.
54  * If rdest == -1, then result is left on all processes in scope.
55  * If LDIA == -1, this array is not accessed, and need not exist.
56  *
57  * CA (output) Integer Array, dimension (LDIA, N)
58  * Contains process column that the amn of each element
59  * of A was found on: i.e., cA(1,2) contains the process
60  * column that the max/min of A(1,2) was found on.
61  * Values are left on process {rdest, cdest} only, others
62  * may be modified, but not left with interesting data.
63  * If rdest == -1, then result is left on all processes in scope.
64  * If LDIA == -1, this array is not accessed, and need not exist.
65  *
66  * LDIA (input) Ptr to int
67  * If (LDIA == -1), then the arrays RA and CA are not accessed.
68  * ELSE leading dimension of the arrays RA and CA. LDIA >= M.
69  *
70  * RDEST (input) Ptr to int
71  * The process row of the destination of the amn.
72  * If rdest == -1, then result is left on all processes in scope.
73  *
74  * CDEST (input) Ptr to int
75  * The process column of the destination of the amn.
76  * If rdest == -1, then CDEST ignored.
77  *
78  * ------------------------------------------------------------------------
79  */
80 {
81  void BI_ArgCheck(int, int, char *, char, char, char, int, int, int, int,
82  int *, int *);
83  void BI_UpdateBuffs(BLACBUFF *);
84  BLACBUFF *BI_GetBuff(int);
85  int BI_BuffIsFree(BLACBUFF *, int);
86  MPI_Datatype BI_GetMpiGeType(BLACSCONTEXT *, int, int, int,
87  MPI_Datatype, int *);
88  BLACBUFF *BI_Pack(BLACSCONTEXT *, BVOID *, BLACBUFF *, MPI_Datatype);
89  MPI_Datatype BI_GetMpiGeType(BLACSCONTEXT *, int, int, int,
90  MPI_Datatype, int *);
91  void BI_Unpack(BLACSCONTEXT *, BVOID *, BLACBUFF *, MPI_Datatype);
93  int, int);
95  int, int);
96  void BI_BeComb(BLACSCONTEXT *, BLACBUFF *, BLACBUFF *, int, VVFUNPTR);
97  void BI_ivvamn(int, char *, char *);
98  void BI_ivvamn2(int, char *, char *);
99  void BI_iMPI_amn(void *, void *, int *, MPI_Datatype *);
100  void BI_iMPI_amn2(void *, void *, int *, MPI_Datatype *);
101 /*
102  * Variable Declarations
103  */
104  VVFUNPTR vvop;
105  BLACBUFF *bp, *bp2;
106  BLACSCONTEXT *ctxt;
107  char ttop, tscope;
108  int i, j, N, dest, idist, length, tlda, tldia, trdest, ierr;
109  int len[2];
110  MPI_Aint disp[2];
111  MPI_Datatype dtypes[2];
112  MPI_Op BlacComb;
113  MPI_Datatype MyType;
114  BI_DistType *dist, mydist;
115  extern BLACBUFF *BI_ActiveQ;
116  extern BLACBUFF BI_AuxBuff;
117 
118  MGetConTxt(Mpval(ConTxt), ctxt);
119  ttop = F2C_CharTrans(top);
120  ttop = Mlowcase(ttop);
121  tscope = F2C_CharTrans(scope);
122  tscope = Mlowcase(tscope);
123 /*
124  * If the user has set the default combine topology, use it instead of
125  * BLACS default
126  */
127 #ifdef DefCombTop
128  if (ttop == ' ') ttop = DefCombTop;
129 #endif
130  if (Mpval(cdest) == -1) trdest = -1;
131  else trdest = Mpval(rdest);
132 #if (BlacsDebugLvl > 0)
133  BI_ArgCheck(Mpval(ConTxt), RT_COMB, __FILE__, tscope, 'u', 'u', Mpval(m),
134  Mpval(n), Mpval(lda), 1, &trdest, Mpaddress(cdest));
135  if (Mpval(ldia) < Mpval(m))
136  {
137  if (Mpval(ldia) != -1)
138  BI_BlacsWarn(Mpval(ConTxt), __LINE__, __FILE__,
139  "LDIA too small (LDIA=%d, but M=%d)", Mpval(ldia),
140  Mpval(m));
141  }
142 #endif
143  if (Mpval(lda) >= Mpval(m)) tlda = Mpval(lda);
144  else tlda = Mpval(m);
145  if (Mpval(ldia) < Mpval(m)) tldia = Mpval(m);
146  else tldia = Mpval(ldia);
147  switch(tscope)
148  {
149  case 'r':
150  ctxt->scp = &ctxt->rscp;
151  if (trdest == -1) dest = -1;
152  else dest = Mpval(cdest);
153  break;
154  case 'c':
155  ctxt->scp = &ctxt->cscp;
156  dest = trdest;
157  break;
158  case 'a':
159  ctxt->scp = &ctxt->ascp;
160  if (trdest == -1) dest = -1;
161  else dest = Mvkpnum(ctxt, trdest, Mpval(cdest));
162  break;
163  default:
164  BI_BlacsErr(Mpval(ConTxt), __LINE__, __FILE__, "Unknown scope '%c'",
165  tscope);
166  }
167 
168 
169 /*
170  * It's not defined how MPI reacts to 0 element reductions, so use BLACS 1-tree
171  * topology if we've got one
172  */
173  if (ttop == ' ') if ( (Mpval(m) < 1) || (Mpval(n) < 1) ) ttop = '1';
174  N = Mpval(m) * Mpval(n);
175 /*
176  * If process who has amn is to be communicated, must set up distance
177  * vector after value vector
178  */
179  if (Mpval(ldia) != -1)
180  {
181  vvop = BI_ivvamn;
182  length = N * sizeof(int);
183  i = length % sizeof(BI_DistType); /* ensure dist vec aligned correctly */
184  if (i) length += sizeof(BI_DistType) - i;
185  idist = length;
186  length += N * sizeof(BI_DistType);
187 /*
188  * For performance, insist second buffer is at least 8-byte aligned
189  */
190  j = 8;
191  if (sizeof(int) > j) j = sizeof(int);
192  i = length % j;
193  if (i) length += j - i;
194  i = 2 * length;
195 
196  bp = BI_GetBuff(i);
197  bp2 = &BI_AuxBuff;
198  bp2->Buff = &bp->Buff[length];
199  BI_imvcopy(Mpval(m), Mpval(n), A, tlda, bp->Buff);
200 /*
201  * Fill in distance vector
202  */
203  if (dest == -1) mydist = ctxt->scp->Iam;
204  else mydist = (ctxt->scp->Np + ctxt->scp->Iam - dest) % ctxt->scp->Np;
205  dist = (BI_DistType *) &bp->Buff[idist];
206  for (i=0; i < N; i++) dist[i] = mydist;
207 
208 /*
209  * Create the MPI datatype holding both user's buffer and distance vector
210  */
211  len[0] = len[1] = N;
212  disp[0] = 0;
213  disp[1] = idist;
214  dtypes[0] = MPI_INT;
215  dtypes[1] = BI_MpiDistType;
216 #ifdef ZeroByteTypeBug
217  if (N > 0)
218  {
219 #endif
220  i = 2;
221  ierr=MPI_Type_create_struct(i, len, disp, dtypes, &MyType);
222  ierr=MPI_Type_commit(&MyType);
223  bp->N = bp2->N = 1;
224  bp->dtype = bp2->dtype = MyType;
225 #ifdef ZeroByteTypeBug
226  }
227  else
228  {
229  bp->N = bp2->N = 0;
230  bp->dtype = bp2->dtype = MPI_INT;
231  }
232 #endif
233  }
234  else
235  {
236  vvop = BI_ivvamn2;
237  length = N * sizeof(int);
238 /*
239  * If A is contiguous, we can use it as one of our buffers
240  */
241  if ( (Mpval(m) == tlda) || (Mpval(n) == 1) )
242  {
243  bp = &BI_AuxBuff;
244  bp->Buff = (char *) A;
245  bp2 = BI_GetBuff(length);
246  }
247  else
248  {
249  bp = BI_GetBuff(length*2);
250  bp2 = &BI_AuxBuff;
251  bp2->Buff = &bp->Buff[length];
252  BI_imvcopy(Mpval(m), Mpval(n), A, tlda, bp->Buff);
253  }
254  bp->N = bp2->N = N;
255  bp->dtype = bp2->dtype = MPI_INT;
256  }
257 
258  switch(ttop)
259  {
260  case ' ': /* use MPI's reduction by default */
261  i = 1;
262  if (Mpval(ldia) == -1)
263  {
264  ierr=MPI_Op_create(BI_iMPI_amn2, i, &BlacComb);
265  }
266  else
267  {
268  ierr=MPI_Op_create(BI_iMPI_amn, i, &BlacComb);
269  BI_AuxBuff.Len = N; /* set this up for the MPI OP wrappers */
270  }
271 
272  if (trdest != -1)
273  {
274  ierr=MPI_Reduce(bp->Buff, bp2->Buff, bp->N, bp->dtype, BlacComb, dest,
275  ctxt->scp->comm);
276  if (ctxt->scp->Iam == dest)
277  {
278  BI_ivmcopy(Mpval(m), Mpval(n), A, tlda, bp2->Buff);
279  if (Mpval(ldia) != -1)
280  BI_TransDist(ctxt, tscope, Mpval(m), Mpval(n), rA, cA, tldia,
281  (BI_DistType *) &bp2->Buff[idist],
282  trdest, Mpval(cdest));
283  }
284  }
285  else
286  {
287  ierr=MPI_Allreduce(bp->Buff, bp2->Buff, bp->N, bp->dtype, BlacComb,
288  ctxt->scp->comm);
289  BI_ivmcopy(Mpval(m), Mpval(n), A, tlda, bp2->Buff);
290  if (Mpval(ldia) != -1)
291  BI_TransDist(ctxt, tscope, Mpval(m), Mpval(n), rA, cA, tldia,
292  (BI_DistType *) &bp2->Buff[idist],
293  trdest, Mpval(cdest));
294  }
295  ierr=MPI_Op_free(&BlacComb);
296  if (Mpval(ldia) != -1)
297 #ifdef ZeroByteTypeBug
298  if (N > 0)
299 #endif
300  ierr=BI_MPI_TYPE_FREE(&MyType);
301  if (BI_ActiveQ) BI_UpdateBuffs(NULL);
302  return;
303  break;
304  case 'i':
305  BI_MringComb(ctxt, bp, bp2, N, vvop, dest, 1);
306  break;
307  case 'd':
308  BI_MringComb(ctxt, bp, bp2, N, vvop, dest, -1);
309  break;
310  case 's':
311  BI_MringComb(ctxt, bp, bp2, N, vvop, dest, 2);
312  break;
313  case 'm':
314  BI_MringComb(ctxt, bp, bp2, N, vvop, dest, ctxt->Nr_co);
315  break;
316  case '1':
317  case '2':
318  case '3':
319  case '4':
320  case '5':
321  case '6':
322  case '7':
323  case '8':
324  case '9':
325  BI_TreeComb(ctxt, bp, bp2, N, vvop, dest, ttop-47);
326  break;
327  case 'f':
328  BI_TreeComb(ctxt, bp, bp2, N, vvop, dest, FULLCON);
329  break;
330  case 't':
331  BI_TreeComb(ctxt, bp, bp2, N, vvop, dest, ctxt->Nb_co);
332  break;
333  case 'h':
334 /*
335  * Use bidirectional exchange if everyone wants answer
336  */
337  if ( (trdest == -1) && !(ctxt->TopsCohrnt) )
338  BI_BeComb(ctxt, bp, bp2, N, vvop);
339  else
340  BI_TreeComb(ctxt, bp, bp2, N, vvop, dest, 2);
341  break;
342  default :
343  BI_BlacsErr(Mpval(ConTxt), __LINE__, __FILE__, "Unknown topology '%c'",
344  ttop);
345  }
346 
347  if (Mpval(ldia) != -1)
348 #ifdef ZeroByteTypeBug
349  if (N > 0)
350 #endif
351  ierr=BI_MPI_TYPE_FREE(&MyType);
352 /*
353  * If I am selected to receive answer
354  */
355  if ( (ctxt->scp->Iam == dest) || (dest == -1) )
356  {
357 /*
358  * Translate the distances stored in the latter part of bp->Buff into
359  * process grid coordinates, and output these coordinates in the
360  * arrays rA and cA.
361  */
362  if (Mpval(ldia) != -1)
363  BI_TransDist(ctxt, tscope, Mpval(m), Mpval(n), rA, cA, tldia,
364  dist, trdest, Mpval(cdest));
365 /*
366  * Unpack the amn array
367  */
368  if (bp != &BI_AuxBuff) BI_ivmcopy(Mpval(m), Mpval(n), A, tlda, bp->Buff);
369  }
370 }
BI_TransDist
void BI_TransDist(BLACSCONTEXT *ctxt, char scope, int m, int n, int *rA, int *cA, int ldrc, BI_DistType *dist, int rdest, int cdest)
Definition: BI_TransDist.c:3
BI_Pack
BLACBUFF * BI_Pack(BLACSCONTEXT *ctxt, BVOID *A, BLACBUFF *bp, MPI_Datatype Dtype)
Definition: BI_Pack.c:2
BI_BeComb
void BI_BeComb(BLACSCONTEXT *ctxt, BLACBUFF *bp, BLACBUFF *bp2, int N, VVFUNPTR Xvvop)
Definition: BI_BeComb.c:10
BI_iMPI_amn
void BI_iMPI_amn(void *in, void *inout, int *N, MPI_Datatype *dtype)
Definition: BI_iMPI_amn.c:3
BI_MpiDistType
#define BI_MpiDistType
Definition: Bdef.h:73
BI_GetBuff
BLACBUFF * BI_GetBuff(int length)
Definition: BI_GetBuff.c:36
bLaCbUfF::Buff
char * Buff
Definition: Bdef.h:56
BI_Unpack
void BI_Unpack(BLACSCONTEXT *ctxt, BVOID *A, BLACBUFF *bp, MPI_Datatype Dtype)
Definition: BI_Unpack.c:3
bLaCsCoNtExT::rscp
BLACSSCOPE rscp
Definition: Bdef.h:25
FULLCON
#define FULLCON
Definition: Bdef.h:100
BI_ivmcopy
void BI_ivmcopy(int m, int n, int *A, int lda, int *buff)
Definition: BI_ivmcopy.c:3
F_CHAR
char * F_CHAR
Definition: pblas.h:109
bLaCsCoNtExT::ascp
BLACSSCOPE ascp
Definition: Bdef.h:25
MGetConTxt
#define MGetConTxt(Context, ctxtptr)
Definition: Bdef.h:200
BVOID
#define BVOID
Definition: Bdef.h:136
bLaCbUfF::Len
int Len
Definition: Bdef.h:57
bLaCbUfF
Definition: Bdef.h:54
VVFUNPTR
void(* VVFUNPTR)(int, char *, char *)
Definition: Bdef.h:68
RT_COMB
#define RT_COMB
Definition: Bdef.h:109
BI_ActiveQ
BLACBUFF * BI_ActiveQ
Definition: BI_GlobalVars.c:9
bLaCsScOpE::Iam
int Iam
Definition: Bdef.h:17
bLaCsCoNtExT::TopsCohrnt
int TopsCohrnt
Definition: Bdef.h:28
BI_MringComb
void BI_MringComb(BLACSCONTEXT *ctxt, BLACBUFF *bp, BLACBUFF *bp2, int N, VVFUNPTR Xvvop, int dest, int nrings)
Definition: BI_MringComb.c:2
bLaCsCoNtExT::Nb_co
int Nb_co
Definition: Bdef.h:30
bLaCbUfF::N
int N
Definition: Bdef.h:61
F2C_CharTrans
#define F2C_CharTrans(c)
Definition: Bdef.h:246
Mpval
#define Mpval(para)
Definition: Bdef.h:261
Mlowcase
#define Mlowcase(C)
Definition: Bdef.h:145
bLaCsCoNtExT
Definition: Bdef.h:23
Mpaddress
#define Mpaddress(para)
Definition: Bdef.h:262
BI_TreeComb
void BI_TreeComb(BLACSCONTEXT *ctxt, BLACBUFF *bp, BLACBUFF *bp2, int N, VVFUNPTR Xvvop, int dest, int nbranches)
Definition: BI_TreeComb.c:19
F_VOID_FUNC
#define F_VOID_FUNC
Definition: Bdef.h:232
igamn2d_
F_VOID_FUNC igamn2d_(int *ConTxt, F_CHAR scope, F_CHAR top, int *m, int *n, int *A, int *lda, int *rA, int *cA, int *ldia, int *rdest, int *cdest)
Definition: igamn2d_.c:7
bLaCsCoNtExT::cscp
BLACSSCOPE cscp
Definition: Bdef.h:25
Cigamn2d
void Cigamn2d()
bLaCsCoNtExT::scp
BLACSSCOPE * scp
Definition: Bdef.h:26
Bdef.h
bLaCsScOpE::comm
MPI_Comm comm
Definition: Bdef.h:15
BI_GetMpiGeType
MPI_Datatype BI_GetMpiGeType(BLACSCONTEXT *ctxt, int m, int n, int lda, MPI_Datatype Dtype, int *N)
Definition: BI_GetMpiGeType.c:2
BI_AuxBuff
BLACBUFF BI_AuxBuff
Definition: BI_GlobalVars.c:10
BI_ivvamn2
void BI_ivvamn2(int N, char *vec1, char *vec2)
Definition: BI_ivvamn2.c:3
bLaCbUfF::dtype
MPI_Datatype dtype
Definition: Bdef.h:60
BI_BlacsErr
void BI_BlacsErr(int ConTxt, int line, char *file, char *form,...)
Definition: BI_BlacsErr.c:3
bLaCsCoNtExT::Nr_co
int Nr_co
Definition: Bdef.h:30
bLaCsScOpE::Np
int Np
Definition: Bdef.h:17
BI_imvcopy
void BI_imvcopy(int m, int n, int *A, int lda, int *buff)
Definition: BI_imvcopy.c:2
BI_BlacsWarn
void BI_BlacsWarn(int ConTxt, int line, char *file, char *form,...)
Definition: BI_BlacsWarn.c:3
BI_ArgCheck
void BI_ArgCheck(int ConTxt, int RoutType, char *routine, char scope, char uplo, char diag, int m, int n, int lda, int nprocs, int *prows, int *pcols)
Definition: BI_ArgCheck.c:4
BI_ivvamn
void BI_ivvamn(int N, char *vec1, char *vec2)
Definition: BI_ivvamn.c:2
BI_DistType
#define BI_DistType
Definition: Bdef.h:72
BI_iMPI_amn2
void BI_iMPI_amn2(void *in, void *inout, int *N, MPI_Datatype *dtype)
Definition: BI_iMPI_amn2.c:2
BI_MPI_TYPE_FREE
#define BI_MPI_TYPE_FREE(t)
Definition: Bdef.h:305
Mvkpnum
#define Mvkpnum(ctxt, prow, pcol)
Definition: Bdef.h:174
BI_UpdateBuffs
void BI_UpdateBuffs(BLACBUFF *Newbp)
Definition: BI_UpdateBuffs.c:3
BI_BuffIsFree
int BI_BuffIsFree(BLACBUFF *bp, int Wait)
Definition: BI_BuffIsFree.c:3