1 /* symbfac.c */
2
3 #include "../spoolesMPI.h"
4
5 #define MYDEBUG 0
6
7 /*--------------------------------------------------------------------*/
8 typedef struct _Msg Msg ;
9 struct _Msg {
10 int id ;
11 int size ;
12 void *base ;
13 MPI_Request req ;
14 Msg *next ;
15 } ;
16 static IVL * initSupportedIVL ( ETree *frontETree, IV *frontOwnersIV,
17 int myid, int msglvl, FILE *msgFile ) ;
18 static void loadInternalIndices ( ETree *frontETree, InpMtx *inpmtxA,
19 InpMtx *inpmtxB, IV *frontOwnersIV, int myid, IVL *symbfacIVL,
20 int msglvl, FILE *msgFile ) ;
21 static void doCooperativeWork ( ETree *frontETree, IV *frontOwnersIV,
22 IVL *symbfacIVL, int stats[], int msglvl,
23 FILE *msgFile, int firsttag, MPI_Comm comm ) ;
24 static int mergeIndices ( int sizeJ, int indJ[],
25 int sizeK, int indK[], int nleft ) ;
26 static Msg * Msg_new ( void ) ;
27 static void Msg_setDefaultFields ( Msg *msg ) ;
28 static void Msg_clearData ( Msg *msg ) ;
29 static void Msg_free ( Msg *msg ) ;
30 static Msg * wakeupFront ( int J, int myid, int owners[],
31 char supportTable[], ETree *frontETree, IVL *symbfacIVL,
32 int firsttag, int stats[], MPI_Comm comm,
33 int msglvl, FILE *msgFile ) ;
34 static Msg * checkRecvMessages ( int J, Msg *firstmsg, int nLeft[],
35 int nodwghts[], IVL *symbfacIVL, int msglvl, FILE *msgFile ) ;
36 static Msg * sendMessages ( int J, int myid, int owners[], int nfront,
37 int nproc, char supportTable[], int par[], IVL *symbfacIVL,
38 Msg *firstMsgSent, int firsttag, int stats[], MPI_Comm comm,
39 int msglvl, FILE *msgFile ) ;
40 static Msg * checkSendMessages ( Msg *firstMsgSent,
41 int msglvl, FILE *msgFile ) ;
42 /*--------------------------------------------------------------------*/
43 /*
44 ------------------------------------------------------------
45 perform the symbolic factorization in parallel
46
47 input --
48 frontETree -- front tree for the factorization
49 frontOwnersIV -- map from fronts to owning process
50 inpmtx -- matrix object
51 firsttag -- first tag to be used for messages,
52 will use tag, ..., tag + nfront - 1
53 msglvl -- message level
54 msgFile -- message file
55 comm -- MPI communicator
56
57 return value --
58 symbfacIVL -- contains front indices for supported fronts
59
60 created -- 98may20, cca
61 ------------------------------------------------------------
62 */
63 IVL *
SymbFac_MPI_initFromInpMtx(ETree * frontETree,IV * frontOwnersIV,InpMtx * inpmtx,int stats[],int msglvl,FILE * msgFile,int firsttag,MPI_Comm comm)64 SymbFac_MPI_initFromInpMtx (
65 ETree *frontETree,
66 IV *frontOwnersIV,
67 InpMtx *inpmtx,
68 int stats[],
69 int msglvl,
70 FILE *msgFile,
71 int firsttag,
72 MPI_Comm comm
73 ) {
74 int lasttag, myid, nproc, tagbound;
75 IVL *symbfacIVL ;
76 /*
77 ---------------
78 check the input
79 ---------------
80 */
81 if ( frontETree == NULL || frontOwnersIV == NULL
82 || inpmtx == NULL || stats == NULL
83 || (msglvl > 0 && msgFile == NULL) ) {
84 fprintf(stderr, "\n fatal error in SymbFac_MPI_initFromInpMtx()"
85 "\n comm = %p, frontETree = %p, frontOwnersIV = %p"
86 "\n inpmtx = %p, firsttag = %d, msglvl = %d, msgFile = %p"
87 "\n bad input\n", comm, frontETree, frontOwnersIV, inpmtx,
88 firsttag, msglvl, msgFile) ;
89 exit(-1) ;
90 }
91 MPI_Comm_rank(comm, &myid) ;
92 MPI_Comm_size(comm, &nproc) ;
93 if ( msglvl > 2 ) {
94 fprintf(msgFile, "\n\n myid = %d, nproc = %d", myid, nproc) ;
95 fflush(msgFile) ;
96 }
97 if ( firsttag < 0 ) {
98 fprintf(stderr, "\n fatal error in SymbFac_MPI_initFromInpMtx()"
99 "\n firsttag = %d\n", firsttag) ;
100 exit(-1) ;
101 }
102 lasttag = firsttag + frontETree->nfront ;
103 if ( lasttag > (tagbound = maxTagMPI(comm)) ) {
104 fprintf(stderr, "\n fatal error in SymbFac_MPI_initFromInpMtx()"
105 "\n lasttag = %d, tag_bound = %d", lasttag, tagbound) ;
106 exit(-1) ;
107 }
108 /*
109 --------------------------------------
110 first, set up the supported IVL object
111 --------------------------------------
112 */
113 symbfacIVL = initSupportedIVL(frontETree, frontOwnersIV,
114 myid, msglvl, msgFile) ;
115 if ( msglvl > 3 ) {
116 fprintf(msgFile, "\n\n local supported IVL after initialization") ;
117 IVL_writeForHumanEye(symbfacIVL, msgFile) ;
118 fflush(msgFile) ;
119 }
120 /*
121 ------------------------------------------------------------------
122 second, fill the index lists of the owned fronts with the internal
123 vertices and boundary vertices available from the matrix pencil
124 ------------------------------------------------------------------
125 */
126 loadInternalIndices(frontETree, inpmtx, NULL, frontOwnersIV,
127 myid, symbfacIVL, msglvl, msgFile ) ;
128 if ( msglvl > 3 ) {
129 fprintf(msgFile, "\n\n after loading internal indices") ;
130 IVL_writeForHumanEye(symbfacIVL, msgFile) ;
131 fflush(msgFile) ;
132 }
133 /*
134 -------------------------------------------------
135 third, cooperate among the processors to generate
136 the rest of the supported symbolic factorization
137 -------------------------------------------------
138 */
139 doCooperativeWork(frontETree, frontOwnersIV, symbfacIVL,
140 stats, msglvl, msgFile, firsttag, comm) ;
141 if ( msglvl > 3 ) {
142 fprintf(msgFile, "\n\n final local supported IVL ") ;
143 IVL_writeForHumanEye(symbfacIVL, msgFile) ;
144 fflush(msgFile) ;
145 }
146 return(symbfacIVL) ; }
147
148 /*--------------------------------------------------------------------*/
149 /*
150 ------------------------------------------------------------
151 perform the symbolic factorization in parallel
152
153 input --
154 frontETree -- front tree for the factorization
155 frontOwnersIV -- map from fronts to owning process
156 pencil -- matrix pencil object
157 firsttag -- first tag to be used for messages,
158 will use tag, ..., tag + nfront - 1
159 msglvl -- message level
160 msgFile -- message file
161 comm -- MPI communicator
162
163 return value --
164 symbfacIVL -- contains front indices for supported fronts
165
166 created -- 98may20, cca
167 ------------------------------------------------------------
168 */
169 IVL *
SymbFac_MPI_initFromPencil(ETree * frontETree,IV * frontOwnersIV,Pencil * pencil,int stats[],int msglvl,FILE * msgFile,int firsttag,MPI_Comm comm)170 SymbFac_MPI_initFromPencil (
171 ETree *frontETree,
172 IV *frontOwnersIV,
173 Pencil *pencil,
174 int stats[],
175 int msglvl,
176 FILE *msgFile,
177 int firsttag,
178 MPI_Comm comm
179 ) {
180 int lasttag, myid, nproc, tagbound;
181 IVL *symbfacIVL ;
182 /*
183 ---------------
184 check the input
185 ---------------
186 */
187 if ( frontETree == NULL || frontOwnersIV == NULL
188 || pencil == NULL || stats == NULL
189 || (msglvl > 0 && msgFile == NULL) ) {
190 fprintf(stderr, "\n fatal error in SymbFac_MPI_initFromPencil()"
191 "\n comm = %p, frontETree = %p, frontOwnersIV = %p"
192 "\n pencil = %p, firsttag = %d, msglvl = %d, msgFile = %p"
193 "\n bad input\n", comm, frontETree, frontOwnersIV, pencil,
194 firsttag, msglvl, msgFile) ;
195 exit(-1) ;
196 }
197 MPI_Comm_rank(comm, &myid) ;
198 MPI_Comm_size(comm, &nproc) ;
199 if ( msglvl > 2 ) {
200 fprintf(msgFile, "\n\n myid = %d, nproc = %d", myid, nproc) ;
201 fflush(msgFile) ;
202 }
203 if ( firsttag < 0 ) {
204 fprintf(stderr, "\n fatal error in SymbFac_MPI_initFromPencil()"
205 "\n firsttag = %d\n", firsttag) ;
206 exit(-1) ;
207 }
208 lasttag = firsttag + frontETree->nfront ;
209 if ( lasttag > (tagbound = maxTagMPI(comm)) ) {
210 fprintf(stderr, "\n fatal error in SymbFac_MPI_initFromPencil()"
211 "\n lasttag = %d, tag_bound = %d", lasttag, tagbound) ;
212 exit(-1) ;
213 }
214 /*
215 --------------------------------------
216 first, set up the supported IVL object
217 --------------------------------------
218 */
219 symbfacIVL = initSupportedIVL(frontETree, frontOwnersIV,
220 myid, msglvl, msgFile) ;
221 if ( msglvl > 3 ) {
222 fprintf(msgFile, "\n\n local supported IVL after initialization") ;
223 IVL_writeForHumanEye(symbfacIVL, msgFile) ;
224 fflush(msgFile) ;
225 }
226 /*
227 ------------------------------------------------------------------
228 second, fill the index lists of the owned fronts with the internal
229 vertices and boundary vertices available from the matrix pencil
230 ------------------------------------------------------------------
231 */
232 loadInternalIndices(frontETree, pencil->inpmtxA, pencil->inpmtxB,
233 frontOwnersIV, myid, symbfacIVL, msglvl, msgFile ) ;
234 if ( msglvl > 3 ) {
235 fprintf(msgFile, "\n\n after loading internal indices") ;
236 IVL_writeForHumanEye(symbfacIVL, msgFile) ;
237 fflush(msgFile) ;
238 }
239 /*
240 -------------------------------------------------
241 third, cooperate among the processors to generate
242 the rest of the supported symbolic factorization
243 -------------------------------------------------
244 */
245 doCooperativeWork(frontETree, frontOwnersIV, symbfacIVL,
246 stats, msglvl, msgFile, firsttag, comm) ;
247 if ( msglvl > 3 ) {
248 fprintf(msgFile, "\n\n final local supported IVL ") ;
249 IVL_writeForHumanEye(symbfacIVL, msgFile) ;
250 fflush(msgFile) ;
251 }
252 return(symbfacIVL) ; }
253
254 /*--------------------------------------------------------------------*/
255 /*
256 ---------------------------------------------------------------------
257 purpose -- initialize the symbolic factorization IVL for process myid
258
259 created -- 98may20
260 ---------------------------------------------------------------------
261 */
262 static IVL *
initSupportedIVL(ETree * frontETree,IV * frontOwnersIV,int myid,int msglvl,FILE * msgFile)263 initSupportedIVL (
264 ETree *frontETree,
265 IV *frontOwnersIV,
266 int myid,
267 int msglvl,
268 FILE *msgFile
269 ) {
270 char *frontSupported ;
271 int J, K, nfront ;
272 int *bndwghts, *nodwghts, *owners, *par, *sizes ;
273 IVL *symbfacIVL ;
274
275 nfront = ETree_nfront(frontETree) ;
276 nodwghts = ETree_nodwghts(frontETree),
277 bndwghts = ETree_bndwghts(frontETree),
278 par = ETree_par(frontETree) ;
279 owners = IV_entries(frontOwnersIV) ;
280 /*
281 -----------------------------------------------------------------
282 create the frontSupported[] vector,
283 frontSupported[J] = 'Y' --> this process supports front J
284 frontSupported[J] = 'N' --> this process does not support front J
285 -----------------------------------------------------------------
286 */
287 frontSupported = CVinit(nfront, 'N') ;
288 for ( J = 0 ; J < nfront ; J++ ) {
289 if ( owners[J] == myid ) {
290 K = J ;
291 while ( K != -1 && frontSupported[K] == 'N' ) {
292 frontSupported[K] = 'Y' ;
293 K = par[K] ;
294 }
295 }
296 }
297 /*
298 -----------------------------------------------
299 initialize the local symbolic factorization IVL
300 -----------------------------------------------
301 */
302 symbfacIVL = IVL_new() ;
303 sizes = IVinit(nfront, 0) ;
304 for ( J = 0 ; J < nfront ; J++ ) {
305 if ( frontSupported[J] == 'Y' ) {
306 sizes[J] = nodwghts[J] + bndwghts[J] ;
307 }
308 }
309 IVL_init3(symbfacIVL, IVL_CHUNKED, nfront, sizes) ;
310 /*
311 ------------------------
312 free the working storage
313 ------------------------
314 */
315 IVfree(sizes) ;
316 CVfree(frontSupported) ;
317
318 return(symbfacIVL) ; }
319
320 /*--------------------------------------------------------------------*/
321 /*
322 --------------------------------------------------------------------
323 perform the symbolic factorization in parallel
324
325 input --
326 frontETree -- front tree for the factorization
327 frontOwnersIV -- map from fronts to owning process
328 pencil -- matrix pencil
329 firsttag -- first tag to be used for messages,
330 will use tag, ..., tag + nfront - 1
331 msglvl -- message level
332 msgFile -- message file
333 comm -- MPI communicator
334
335 return value --
336 symbfacIVL -- contains front indices for supported fronts
337
338 created -- 98may20, cca
339 --------------------------------------------------------------------
340 */
341 static void
doCooperativeWork(ETree * frontETree,IV * frontOwnersIV,IVL * symbfacIVL,int stats[],int msglvl,FILE * msgFile,int firsttag,MPI_Comm comm)342 doCooperativeWork (
343 ETree *frontETree,
344 IV *frontOwnersIV,
345 IVL *symbfacIVL,
346 int stats[],
347 int msglvl,
348 FILE *msgFile,
349 int firsttag,
350 MPI_Comm comm
351 ) {
352 char *frontIsSupported, *status, *supportTable ;
353 Ideq *dequeue ;
354 int count, ii, iproc, J, K, myid, nDJ, nfront, npath, nproc,
355 nvtx, sizeJ, sizeK ;
356 int *bndwghts, *frontOwners, *indicesJ, *indicesK,
357 *nactiveChild, *nLeft, *nodwghts, *par, *vtxToFront ;
358 Msg *firstMsgSent ;
359 Msg **p_msg ;
360 /*
361 -------------------------------
362 extract pointers and dimensions
363 -------------------------------
364 */
365 nfront = ETree_nfront(frontETree) ;
366 nvtx = ETree_nvtx(frontETree) ;
367 nodwghts = ETree_nodwghts(frontETree),
368 bndwghts = ETree_bndwghts(frontETree),
369 vtxToFront = ETree_vtxToFront(frontETree),
370 par = ETree_par(frontETree) ;
371 IV_sizeAndEntries(frontOwnersIV, &nfront, &frontOwners) ;
372 MPI_Comm_rank(comm, &myid) ;
373 MPI_Comm_size(comm, &nproc) ;
374 if ( msglvl > 2 ) {
375 fprintf(msgFile, "\n\n myid = %d, nproc = %d", myid, nproc) ;
376 fflush(msgFile) ;
377 }
378 /*
379 --------------------------
380 generate the support table
381 --------------------------
382 */
383 supportTable = CVinit(nfront*nproc, 'N') ;
384 for ( J = 0 ; J < nfront ; J++ ) {
385 iproc = frontOwners[J] ;
386 frontIsSupported = supportTable + iproc*nfront ;
387 for ( K = J ; K != -1 && frontIsSupported[K] == 'N' ; K = par[K] ) {
388 frontIsSupported[K] = 'Y' ;
389 }
390 }
391 if ( msglvl > 3 ) {
392 fprintf(msgFile, "\n\n supportTable") ;
393 for ( J = 0 ; J < nfront ; J++ ) {
394 fprintf(msgFile, "\n") ;
395 for ( iproc = 0 ; iproc < nproc ; iproc++ ) {
396 fprintf(msgFile, " %c", supportTable[J + iproc*nfront]) ;
397 }
398 }
399 fflush(msgFile) ;
400 }
401 /*
402 ---------------------------------------------------------
403 initialize the status[] and nactiveChild[] vectors and
404 load the dequeue with the leaves of the supported subtree
405 ---------------------------------------------------------
406 */
407 frontIsSupported = supportTable + myid*nfront ;
408 status = CVinit(nfront, 'F') ;
409 for ( J = npath = 0 ; J < nfront ; J++ ) {
410 if ( status[J] == 'F' && frontIsSupported[J] == 'Y' ) {
411 npath++ ;
412 for ( K = J ; K != -1 && status[K] == 'F' ; K = par[K] ) {
413 status[K] = 'W' ;
414 }
415 }
416 }
417 dequeue = Ideq_new() ;
418 Ideq_resize(dequeue, npath) ;
419 nactiveChild = IVinit(nfront, 0) ;
420 for ( J = 0 ; J < nfront ; J++ ) {
421 if ( status[J] == 'W' ) {
422 if ( (K = par[J]) != -1 ) {
423 nactiveChild[K]++ ;
424 }
425 }
426 }
427 for ( J = 0 ; J < nfront ; J++ ) {
428 if ( frontOwners[J] == myid && nactiveChild[J] == 0 ) {
429 Ideq_insertAtTail(dequeue, J) ;
430 }
431 }
432 /*
433 ---------------------------
434 generate the nLeft[] vector
435 ---------------------------
436 */
437 nLeft = IVinit(nfront, -1) ;
438 for ( J = 0 ; J < nfront ; J++ ) {
439 if ( frontOwners[J] == myid ) {
440 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
441 for ( ii = sizeJ - 1, count = 0 ; ii >= 0 ; ii-- ) {
442 if ( indicesJ[ii] != -1 ) {
443 break ;
444 }
445 count++ ;
446 }
447 nLeft[J] = count ;
448 }
449 }
450 /*
451 -------------------------------------------
452 generate the vector of pointers to messages
453 -------------------------------------------
454 */
455 firstMsgSent = NULL ;
456 ALLOCATE(p_msg, struct _Msg *, nfront) ;
457 for ( J = 0 ; J < nfront ; J++ ) {
458 p_msg[J] = NULL ;
459 }
460 /*
461 ----------------------------------
462 loop while the dequeue is nonempty
463 ----------------------------------
464 */
465 while ( (J = Ideq_removeFromHead(dequeue)) != -1 ) {
466 if ( msglvl > 2 ) {
467 fprintf(msgFile, "\n\n ### popping %d from dequeue", J) ;
468 fflush(msgFile) ;
469 }
470 if ( status[J] == 'W' ) {
471 /*
472 --------------------------------------------
473 wake up the front, post its receive messages
474 --------------------------------------------
475 */
476 p_msg[J] = wakeupFront(J, myid, frontOwners, supportTable,
477 frontETree, symbfacIVL, firsttag,
478 stats, comm, msglvl, msgFile) ;
479 status[J] = 'A' ;
480 }
481 if ( p_msg[J] != NULL ) {
482 /*
483 -----------------------------
484 try to receive messages for J
485 -----------------------------
486 */
487 p_msg[J] = checkRecvMessages(J, p_msg[J], nLeft, nodwghts,
488 symbfacIVL, msglvl, msgFile) ;
489 }
490 if ( p_msg[J] == NULL ) {
491 if ( msglvl > 3 ) {
492 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
493 fprintf(msgFile,
494 "\n adjacency list for J = %d, nJ = %d, bndJ = %d",
495 J, nodwghts[J], bndwghts[J]) ;
496 IVfprintf(msgFile, sizeJ, indicesJ) ;
497 fflush(msgFile) ;
498 }
499 /*
500 -------------------------------------------
501 all messages received:
502 check for the length of the indices vector,
503 set the boundary size of the front
504 ------------------------------------------
505 */
506 status[J] = 'D' ;
507 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
508 for ( ii = 0 ; ii < sizeJ ; ii++ ) {
509 if ( indicesJ[ii] == -1 ) {
510 break ;
511 }
512 }
513 if ( ii != sizeJ ) {
514 if ( msglvl > 3 ) {
515 fprintf(msgFile, "\n WHOA, front J %d, sizeJ %d", J, sizeJ);
516 fprintf(msgFile, "\n old bndwght = %d, new bndwght = %d",
517 bndwghts[J], ii - nodwghts[J]) ;
518 IVfprintf(msgFile, sizeJ, indicesJ) ;
519 fflush(msgFile) ;
520 }
521 bndwghts[J] = ii - nodwghts[J] ;
522 IVL_setList(symbfacIVL, J, ii, indicesJ) ;
523 }
524 if ( frontOwners[J] == myid ) {
525 /*
526 ------------------------------------------
527 front is owned, send messages if necessary
528 ------------------------------------------
529 */
530 firstMsgSent = sendMessages(J, myid, frontOwners, nfront,
531 nproc, supportTable, par,
532 symbfacIVL, firstMsgSent, firsttag,
533 stats, comm, msglvl, msgFile) ;
534 }
535 if ( (K = par[J]) != -1 ) {
536 if ( frontOwners[K] == myid && nLeft[K] > 0 ) {
537 /*
538 ------------------------------
539 merge bnd{J} into K cup bnd{K}
540 ------------------------------
541 */
542 nDJ = nodwghts[J] ;
543 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
544 IVL_listAndSize(symbfacIVL, K, &sizeK, &indicesK) ;
545 if ( msglvl > 3 ) {
546 fprintf(msgFile, "\n merging J = %d into K = %d", J, K) ;
547 fprintf(msgFile, ", nleft = %d", nLeft[K]) ;
548 fprintf(msgFile, "\n before, boundary indices of J") ;
549 IVfprintf(msgFile, sizeJ - nDJ, indicesJ + nDJ) ;
550 fprintf(msgFile, "\n before, indices of K") ;
551 IVfprintf(msgFile, sizeK, indicesK) ;
552 }
553 nLeft[K] = mergeIndices(sizeJ - nDJ, indicesJ + nDJ,
554 sizeK, indicesK, nLeft[K]) ;
555 if ( msglvl > 3 ) {
556 fprintf(msgFile,
557 "\n after, indices of K, nleft = %d", nLeft[K]) ;
558 IVfprintf(msgFile, sizeK, indicesK) ;
559 }
560 }
561 if ( --nactiveChild[K] == 0 ) {
562 /*
563 ----------------------
564 put K on head of queue
565 ----------------------
566 */
567 Ideq_insertAtHead(dequeue, K) ;
568 }
569 }
570 } else {
571 /*
572 ----------------------
573 put J on tail of queue
574 ----------------------
575 */
576 Ideq_insertAtTail(dequeue, J) ;
577 }
578 firstMsgSent = checkSendMessages(firstMsgSent, msglvl, msgFile) ;
579 }
580 /*
581 ------------------------------------------
582 check for all sent messages to be received
583 ------------------------------------------
584 */
585 do {
586 firstMsgSent = checkSendMessages(firstMsgSent, msglvl, msgFile) ;
587 } while ( firstMsgSent != NULL ) ;
588 /*
589 --------------------------------
590 post-process the adjacency lists
591 --------------------------------
592 */
593 /*
594 if ( msglvl > 2 ) {
595 fprintf(msgFile, "\n bndwghts") ;
596 IVfprintf(msgFile, nfront, bndwghts) ;
597 fprintf(msgFile, "\n nodwghts") ;
598 IVfprintf(msgFile, nfront, nodwghts) ;
599 }
600 frontIsSupported = supportTable + myid*nfront ;
601 for ( J = 0 ; J < nfront ; J++ ) {
602 if ( frontIsSupported[J] == 'Y' ) {
603 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
604 for ( ii = 0 ; ii < sizeJ ; ii++ ) {
605 if ( indicesJ[ii] == -1 ) {
606 break ;
607 }
608 }
609 if ( ii != sizeJ ) {
610 if ( msglvl > 3 ) {
611 fprintf(msgFile, "\n WHOA, front J %d, sizeJ %d", J, sizeJ);
612 fprintf(msgFile, "\n old bndwght = %d, new bndwght = %d",
613 bndwghts[J], ii - nodwghts[J]) ;
614 IVfprintf(msgFile, sizeJ, indicesJ) ;
615 fflush(msgFile) ;
616 }
617 bndwghts[J] = ii - nodwghts[J] ;
618 IVL_setList(symbfacIVL, J, ii, indicesJ) ;
619 }
620 }
621 }
622 */
623 /*
624 ------------------------
625 free the working storage
626 ------------------------
627 */
628 CVfree(supportTable) ;
629 CVfree(status) ;
630 IVfree(nactiveChild) ;
631 IVfree(nLeft) ;
632 Ideq_free(dequeue) ;
633 FREE(p_msg) ;
634
635 return ; }
636
637 /*--------------------------------------------------------------------*/
638 /*
639 -----------------------
640 constructor
641
642 created -- 98may20, cca
643 -----------------------
644 */
645 static Msg *
Msg_new(void)646 Msg_new (
647 void
648 ) {
649 Msg *msg ;
650
651 ALLOCATE(msg, struct _Msg, 1) ;
652 Msg_setDefaultFields(msg) ;
653
654 return(msg) ; }
655
656 /*--------------------------------------------------------------------*/
657 /*
658 -----------------------
659 set the default fields
660
661 created -- 98may20, cca
662 -----------------------
663 */
664 static void
Msg_setDefaultFields(Msg * msg)665 Msg_setDefaultFields (
666 Msg *msg
667 ) {
668 msg->id = -1 ;
669 msg->size = 0 ;
670 msg->base = NULL ;
671 msg->next = NULL ;
672
673 return ; }
674
675 /*--------------------------------------------------------------------*/
676 /*
677 -----------------------
678 clear the data
679
680 created -- 98may20, cca
681 -----------------------
682 */
683 static void
Msg_clearData(Msg * msg)684 Msg_clearData (
685 Msg *msg
686 ) {
687
688 Msg_setDefaultFields(msg) ;
689
690 return ; }
691
692 /*--------------------------------------------------------------------*/
693 /*
694 -----------------------
695 free the object
696
697 created -- 98may20, cca
698 -----------------------
699 */
700 static void
Msg_free(Msg * msg)701 Msg_free (
702 Msg *msg
703 ) {
704
705 Msg_clearData(msg) ;
706 FREE(msg) ;
707
708 return ; }
709
710 /*--------------------------------------------------------------------*/
711 /*
712 -------------------------------------------------------------
713 when front J is first examined, post Irecv's for its messages
714
715 created -- 98may20, cca
716 -------------------------------------------------------------
717 */
718 static Msg *
wakeupFront(int J,int myid,int owners[],char supportTable[],ETree * frontETree,IVL * symbfacIVL,int firsttag,int stats[],MPI_Comm comm,int msglvl,FILE * msgFile)719 wakeupFront (
720 int J,
721 int myid,
722 int owners[],
723 char supportTable[],
724 ETree *frontETree,
725 IVL *symbfacIVL,
726 int firsttag,
727 int stats[],
728 MPI_Comm comm,
729 int msglvl,
730 FILE *msgFile
731 ) {
732 Msg *first, *msg ;
733 int I, nfront ;
734 int *bndwghts, *fch, *indices, *nodwghts, *sib ;
735
736 nfront = ETree_nfront(frontETree) ;
737 first = NULL ;
738 if ( owners[J] == myid ) {
739 /*
740 ---------------------------------------------------
741 owned front, post messages for unsupported children
742 ---------------------------------------------------
743 */
744 if ( msglvl > 2 ) {
745 fprintf(msgFile, "\n\n waking up owned front %d", J) ;
746 fflush(msgFile) ;
747 }
748 fch = ETree_fch(frontETree) ;
749 sib = ETree_sib(frontETree) ;
750 nodwghts = ETree_nodwghts(frontETree) ;
751 bndwghts = ETree_bndwghts(frontETree) ;
752 for ( I = fch[J] ; I != -1 ; I = sib[I] ) {
753 if ( supportTable[I + myid*nfront] != 'Y' ) {
754 if ( msglvl > 2 ) {
755 fprintf(msgFile, "\n unsupported child %d", I) ;
756 fflush(msgFile) ;
757 }
758 /*
759 --------------------------------------------------------
760 create a message, allocate temporary storage for indices
761 --------------------------------------------------------
762 */
763 msg = Msg_new() ;
764 msg->id = I ;
765 msg->size = nodwghts[I] + bndwghts[I] ;
766 msg->base = (void *) IVinit(msg->size, -1) ;
767 msg->next = first ;
768 first = msg ;
769 if ( msglvl > 2 ) {
770 fprintf(msgFile,
771 "\n 1. posting Irecv, msg %p, source = %d, tag = %d",
772 msg, owners[I], firsttag + I) ;
773 fflush(msgFile) ;
774 }
775 stats[1]++ ;
776 stats[3] += msg->size * sizeof(int) ;
777 MPI_Irecv(msg->base, msg->size, MPI_INT, owners[I],
778 firsttag + I, comm, &msg->req) ;
779 }
780 }
781 } else if ( supportTable[J + myid*nfront] == 'Y' ) {
782 /*
783 --------------------------------------------------------
784 supported but not owned. create a message but point
785 into storage for the indices in the symbolic IVL object.
786 --------------------------------------------------------
787 */
788 if ( msglvl > 2 ) {
789 fprintf(msgFile, "\n\n waking up supported front %d", J) ;
790 fflush(msgFile) ;
791 }
792 msg = Msg_new() ;
793 msg->id = J ;
794 IVL_listAndSize(symbfacIVL, J, &msg->size, &indices) ;
795 msg->base = (void *) indices ;
796 msg->next = first ;
797 first = msg ;
798 if ( msglvl > 2 ) {
799 fprintf(msgFile,
800 "\n 1. posting Irecv, msg %p, source = %d, tag = %d",
801 msg, owners[J], firsttag + J) ;
802 fflush(msgFile) ;
803 }
804 stats[1]++ ;
805 stats[3] += msg->size * sizeof(int) ;
806 MPI_Irecv(msg->base, msg->size, MPI_INT, owners[J], firsttag + J,
807 comm, &msg->req) ;
808 }
809 return(first) ; }
810
811 /*--------------------------------------------------------------------*/
812 /*
813 -------------------------------------------------------------
814 check the messages waiting to be received for front J
815
816 return value -- pointer to first message still to be received
817
818 created -- 98may20, cca
819 -------------------------------------------------------------
820 */
821 static Msg *
checkRecvMessages(int J,Msg * firstmsg,int nLeft[],int nodwghts[],IVL * symbfacIVL,int msglvl,FILE * msgFile)822 checkRecvMessages (
823 int J,
824 Msg *firstmsg,
825 int nLeft[],
826 int nodwghts[],
827 IVL *symbfacIVL,
828 int msglvl,
829 FILE *msgFile
830 ) {
831 int flag, I, nDI, sizeI, sizeJ ;
832 int *indicesI, *indicesJ ;
833 MPI_Status status ;
834 Msg *msg, *nextmsg ;
835
836 if ( msglvl > 2 ) {
837 fprintf(msgFile, "\n checking for received messages") ;
838 fflush(msgFile) ;
839 }
840 for ( msg = firstmsg, firstmsg = NULL ; msg != NULL ; msg = nextmsg ) {
841 nextmsg = msg->next ;
842 msg->next = NULL ;
843 if ( msglvl > 2 ) {
844 fprintf(msgFile, "\n checking msg %p : id %d, size %d",
845 msg, msg->id, msg->size) ;
846 fflush(msgFile) ;
847 }
848 MPI_Test(&msg->req, &flag, &status) ;
849 if ( msglvl > 2 ) {
850 fprintf(msgFile, ", flag = %d", flag) ;
851 fflush(msgFile) ;
852 }
853 if ( flag != 0 ) {
854 if ( msg->id == J ) {
855 /*
856 --------------------------------
857 J is supported but not owned,
858 indices are in the symbolic IVL,
859 free the message
860 --------------------------------
861 */
862 if ( msglvl > 2 ) {
863 fprintf(msgFile, "\n supported front %d", J) ;
864 fflush(msgFile) ;
865 }
866 Msg_free(msg) ;
867 } else {
868 /*
869 ---------------------------------------------------
870 message is from an unsupported child
871 merge indices, then free indices, then free message
872 ---------------------------------------------------
873 */
874 if ( msglvl > 2 ) {
875 fprintf(msgFile,
876 "\n unsupported child front %d", msg->id) ;
877 fflush(msgFile) ;
878 }
879 if ( nLeft[J] > 0 ) {
880 I = msg->id ;
881 nDI = nodwghts[I] ;
882 sizeI = msg->size - nDI ;
883 indicesI = (int *) msg->base + nDI ;
884 /*
885 sizeI = msg->size ;
886 indicesI = (int *) msg->base ;
887 */
888 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
889 if ( msglvl > 3 ) {
890 fprintf(msgFile, "\n merging I = %d into J = %d", I, J) ;
891 fprintf(msgFile, ", nleft = %d", nLeft[J]) ;
892 fprintf(msgFile, "\n before, boundary indices of I") ;
893 IVfprintf(msgFile, sizeI, indicesI) ;
894 fprintf(msgFile, "\n before, indices of J") ;
895 IVfprintf(msgFile, sizeJ, indicesJ) ;
896 }
897 nLeft[J] = mergeIndices(sizeI, indicesI,
898 sizeJ, indicesJ, nLeft[J]) ;
899 }
900 if ( msglvl > 3 ) {
901 fprintf(msgFile,
902 "\n after, indices of J, nleft = %d", nLeft[J]) ;
903 IVfprintf(msgFile, sizeJ, indicesJ) ;
904 }
905 IVfree(msg->base) ;
906 Msg_free(msg) ;
907 }
908 } else {
909 /*
910 ----------------------------
911 keep the message on the list
912 ----------------------------
913 */
914 msg->next = firstmsg ;
915 firstmsg = msg ;
916 }
917 }
918 return(firstmsg) ; }
919
920 /*--------------------------------------------------------------------*/
921 /*
922 */
923 static Msg *
sendMessages(int J,int myid,int owners[],int nfront,int nproc,char supportTable[],int par[],IVL * symbfacIVL,Msg * firstMsgSent,int firsttag,int stats[],MPI_Comm comm,int msglvl,FILE * msgFile)924 sendMessages (
925 int J,
926 int myid,
927 int owners[],
928 int nfront,
929 int nproc,
930 char supportTable[],
931 int par[],
932 IVL *symbfacIVL,
933 Msg *firstMsgSent,
934 int firsttag,
935 int stats[],
936 MPI_Comm comm,
937 int msglvl,
938 FILE *msgFile
939 ) {
940 int destination, K, sizeJ ;
941 int *indicesJ ;
942 Msg *msg ;
943 if ( msglvl > 2 ) {
944 fprintf(msgFile, "\n ready to send messages") ;
945 fflush(msgFile) ;
946 }
947 /*
948 -------------------------------
949 get size and pointer to indices
950 -------------------------------
951 */
952 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
953 /*
954 -----------------------------------
955 send indices to supported processes
956 -----------------------------------
957 */
958 for ( destination = 0 ; destination < nproc ; destination++ ) {
959 if ( destination != myid
960 && supportTable[J + destination*nfront] == 'Y' ) {
961 msg = Msg_new() ;
962 msg->id = J ;
963 msg->size = sizeJ ;
964 msg->base = (void *) indicesJ ;
965 msg->next = firstMsgSent ;
966 firstMsgSent = msg ;
967 if ( msglvl > 2 ) {
968 fprintf(msgFile,
969 "\n 1. posting Isend, msg %p, destination = %d, tag = %d",
970 msg, destination, firsttag + J) ;
971 fflush(msgFile) ;
972 }
973 stats[0]++ ;
974 stats[2] += msg->size * sizeof(int) ;
975 MPI_Isend(msg->base, msg->size, MPI_INT, destination,
976 firsttag + J, comm, &msg->req) ;
977 }
978 }
979 if ( (K = par[J]) != -1 && supportTable[J + owners[K]*nfront] != 'Y' ) {
980 /*
981 ----------------------------------------------------
982 send indices to unsupported process that owns parent
983 ----------------------------------------------------
984 */
985 msg = Msg_new() ;
986 msg->id = J ;
987 msg->size = sizeJ ;
988 msg->base = (void *) indicesJ ;
989 msg->next = firstMsgSent ;
990 firstMsgSent = msg ;
991 if ( msglvl > 2 ) {
992 fprintf(msgFile,
993 "\n 2. posting Isend, msg %p, destination = %d, tag = %d",
994 msg, owners[K], firsttag + J) ;
995 fflush(msgFile) ;
996 }
997 stats[0]++ ;
998 stats[2] += msg->size * sizeof(int) ;
999 MPI_Isend(msg->base, msg->size, MPI_INT, owners[K],
1000 firsttag + J, comm, &msg->req) ;
1001 }
1002 return(firstMsgSent) ; }
1003
1004 /*--------------------------------------------------------------------*/
1005 /*
1006 ---------------------------------------------------------------
1007 merge indices in indJ[sizeJ] into indK[sizeK],
1008 there are nleft free locations at the end of indK[].
1009 both indJ[sizeJ] and indK[sizeK] are in ascending order.
1010 return value is the number of present free locations in indK[].
1011
1012 created -- 98may20, cca
1013 ---------------------------------------------------------------
1014 */
1015 static int
mergeIndices(int sizeJ,int indJ[],int sizeK,int indK[],int nleft)1016 mergeIndices (
1017 int sizeJ,
1018 int indJ[],
1019 int sizeK,
1020 int indK[],
1021 int nleft
1022 ) {
1023 int firstK, ii, jj, jlast, kk, v ;
1024
1025 firstK = indK[0] ;
1026 kk = jlast = sizeK - nleft ;
1027 ii = 0, jj = 0 ;
1028 while ( ii < sizeJ && jj < jlast ) {
1029 v = indJ[ii] ;
1030 if ( v < firstK ) {
1031 ii++ ;
1032 } else if ( v < indK[jj] ) {
1033 indK[kk++] = v ;
1034 ii++ ;
1035 } else if ( v == indK[jj] ) {
1036 ii++, jj++ ;
1037 } else {
1038 jj++ ;
1039 }
1040 }
1041 for ( ; ii < sizeJ ; ii++ ) {
1042 v = indJ[ii] ;
1043 if ( v >= 0 ) {
1044 indK[kk++] = indJ[ii] ;
1045 }
1046 }
1047 if ( kk > jlast ) {
1048 IVqsortUp(kk, indK) ;
1049 }
1050 return(sizeK - kk) ; }
1051
1052 /*--------------------------------------------------------------------*/
1053 /*
1054 ------------------------------------------------------------------
1055 load the owned index lists with information from the InpMtx object
1056
1057 created -- 98may20, cca
1058 ------------------------------------------------------------------
1059 */
1060 static void
loadInternalIndices(ETree * frontETree,InpMtx * inpmtxA,InpMtx * inpmtxB,IV * frontOwnersIV,int myid,IVL * symbfacIVL,int msglvl,FILE * msgFile)1061 loadInternalIndices (
1062 ETree *frontETree,
1063 InpMtx *inpmtxA,
1064 InpMtx *inpmtxB,
1065 IV *frontOwnersIV,
1066 int myid,
1067 IVL *symbfacIVL,
1068 int msglvl,
1069 FILE *msgFile
1070 ) {
1071 int count, ii, J, nfront, nvtx, off, sizeJ, v, vsize, w ;
1072 int *frontOwners, *head, *indicesJ, *link, *mark, *vind, *vtxToFront ;
1073
1074 nfront = ETree_nfront(frontETree) ;
1075 nvtx = ETree_nvtx(frontETree) ;
1076 vtxToFront = ETree_vtxToFront(frontETree) ;
1077 frontOwners = IV_entries(frontOwnersIV) ;
1078 /*
1079 -------------------------------------------------
1080 set up the lists of vertices for each owned front
1081 -------------------------------------------------
1082 */
1083 head = IVinit(nfront, -1) ;
1084 link = IVinit(nvtx, -1) ;
1085 mark = IVinit(nvtx, -1) ;
1086 for ( v = 0 ; v < nvtx ; v++ ) {
1087 J = vtxToFront[v] ;
1088 if ( frontOwners[J] == myid ) {
1089 link[v] = head[J] ;
1090 head[J] = v ;
1091 }
1092 }
1093 /*
1094 --------------------------------------------------------
1095 load the internal indices into the owned adjacency lists
1096 --------------------------------------------------------
1097 */
1098 for ( J = 0 ; J < nfront ; J++ ) {
1099 if ( frontOwners[J] == myid ) {
1100 IVL_listAndSize(symbfacIVL, J, &sizeJ, &indicesJ) ;
1101 count = 0 ;
1102 /*
1103 ---------------------
1104 load internal indices
1105 ---------------------
1106 */
1107 for ( v = head[J] ; v != -1 ; v = link[v] ) {
1108 mark[v] = J ;
1109 indicesJ[count++] = v ;
1110 }
1111 /*
1112 -----------------------------------------
1113 load indices from the first InpMtx object
1114 -----------------------------------------
1115 */
1116 for ( v = head[J] ; v != -1 && count < sizeJ ; v = link[v] ) {
1117 InpMtx_vector(inpmtxA, v, &vsize, &vind) ;
1118 for ( ii = 0 ; ii < vsize && count < sizeJ ; ii++ ) {
1119 off = vind[ii] ;
1120 if ( off >= 0 ) {
1121 w = v + off ;
1122 } else {
1123 w = v - off ;
1124 }
1125 if ( vtxToFront[w] < J ) {
1126 fprintf(msgFile,
1127 "\n error, J = %d, w = %d, map[%d] = %d",
1128 J, w, w, vtxToFront[w]) ;
1129 exit(-1) ;
1130 }
1131 if ( mark[w] != J ) {
1132 mark[w] = J ;
1133 indicesJ[count++] = w ;
1134 }
1135 }
1136 }
1137 if ( inpmtxB != NULL ) {
1138 /*
1139 -------------------
1140 load indices from B
1141 -------------------
1142 */
1143 for ( v = head[J] ; v != -1 && count < sizeJ ; v = link[v] ) {
1144 InpMtx_vector(inpmtxB, v, &vsize, &vind) ;
1145 for ( ii = 0 ; ii < vsize && count < sizeJ ; ii++ ) {
1146 off = vind[ii] ;
1147 if ( off >= 0 ) {
1148 w = v + off ;
1149 } else {
1150 w = v - off ;
1151 }
1152 if ( vtxToFront[w] < J ) {
1153 fprintf(msgFile,
1154 "\n error, J = %d, w = %d, map[%d] = %d",
1155 J, w, w, vtxToFront[w]) ;
1156 exit(-1) ;
1157 }
1158 if ( mark[w] != J ) {
1159 mark[w] = J ;
1160 indicesJ[count++] = w ;
1161 }
1162 }
1163 }
1164 }
1165 /*
1166 -----------------------------------
1167 sort the indices in ascending order
1168 -----------------------------------
1169 */
1170 IVqsortUp(count, indicesJ) ;
1171 }
1172 }
1173 /*
1174 ------------------------
1175 free the working storage
1176 ------------------------
1177 */
1178 IVfree(head) ;
1179 IVfree(link) ;
1180 IVfree(mark) ;
1181
1182 return ; }
1183
1184 /*--------------------------------------------------------------------*/
1185 /*
1186 */
1187 static Msg *
checkSendMessages(Msg * firstMsgSent,int msglvl,FILE * msgFile)1188 checkSendMessages (
1189 Msg *firstMsgSent,
1190 int msglvl,
1191 FILE *msgFile
1192 ) {
1193 int flag ;
1194 MPI_Status status ;
1195 Msg *msg, *nextmsg ;
1196
1197 if ( firstMsgSent != NULL ) {
1198 if ( msglvl > 2 ) {
1199 fprintf(msgFile, "\n checking for sent messages") ;
1200 fflush(msgFile) ;
1201 }
1202 for ( msg = firstMsgSent, firstMsgSent = NULL ;
1203 msg != NULL ;
1204 msg = nextmsg ) {
1205 nextmsg = msg->next ;
1206 if ( msglvl > 2 ) {
1207 fprintf(msgFile, "\n checking msg %p : id %d, size %d",
1208 msg, msg->id, msg->size) ;
1209 fflush(msgFile) ;
1210 }
1211 MPI_Test(&msg->req, &flag, &status) ;
1212 if ( msglvl > 2 ) {
1213 fprintf(msgFile, ", flag = %d", flag) ;
1214 fflush(msgFile) ;
1215 }
1216 if ( flag != 0 ) {
1217 if ( msglvl > 2 ) {
1218 fprintf(msgFile, ", free'ing object") ;
1219 fflush(msgFile) ;
1220 }
1221 Msg_free(msg) ;
1222 } else {
1223 msg->next = firstMsgSent ;
1224 firstMsgSent = msg ;
1225 }
1226 }
1227 }
1228
1229 return(firstMsgSent) ; }
1230
1231 /*--------------------------------------------------------------------*/
1232