/* util.c */
#include "../SolveMap.h"
/*--------------------------------------------------------------------*/
/*
----------------------------------------------------
purpose -- return the owner of block (rowid, colid).
created -- 98mar19, cca
----------------------------------------------------
*/
int
SolveMap_owner (
SolveMap *solvemap,
int rowid,
int colid
) {
int ii, loc, J, K, nblock ;
int *colids, *rowids, *map ;
/*
---------------
check the input
---------------
*/
if ( solvemap == NULL
|| rowid < 0 || rowid >= solvemap->nfront
|| colid < 0 || colid >= solvemap->nfront ) {
fprintf(stderr, "\n fatal error in SolveMap_owner(%p,%d,%d)"
"\n bad input\n", solvemap, rowid, colid) ;
exit(-1) ;
}
if ( rowid == colid ) {
/*
--------------
diagonal block
--------------
*/
return(solvemap->owners[rowid]) ;
} else if ( rowid > colid && solvemap->symmetryflag > 0 ) {
/*
---------------------------------
lower block (K,J) = (rowid,colid)
---------------------------------
*/
nblock = solvemap->nblockLower ;
rowids = solvemap->rowidsLower ;
colids = solvemap->colidsLower ;
map = solvemap->mapLower ;
K = rowid ;
J = colid ;
loc = IVlocateViaBinarySearch(nblock, colids, J) ;
if ( loc == -1 ) {
return(-1) ;
}
for ( ii = loc ; ii >= 0 ; ii-- ) {
if ( colids[ii] == J && rowids[ii] == K ) {
return(map[ii]) ;
}
}
for ( ii = loc + 1 ; ii < nblock ; ii++ ) {
if ( colids[ii] == J && rowids[ii] == K ) {
return(map[ii]) ;
}
}
} else {
/*
-----------------
upper block (J,K)
-----------------
*/
nblock = solvemap->nblockUpper ;
rowids = solvemap->rowidsUpper ;
colids = solvemap->colidsUpper ;
map = solvemap->mapUpper ;
if ( rowid > colid ) {
J = colid ; K = rowid ;
} else {
J = rowid ; K = colid ;
}
loc = IVlocateViaBinarySearch(nblock, rowids, J) ;
if ( loc == -1 ) {
return(-1) ;
}
for ( ii = loc ; ii >= 0 ; ii-- ) {
if ( rowids[ii] == J && colids[ii] == K ) {
return(map[ii]) ;
}
}
for ( ii = loc +1 ; ii < nblock ; ii++ ) {
if ( rowids[ii] == J && colids[ii] == K ) {
return(map[ii]) ;
}
}
}
return(-1) ; }
/*--------------------------------------------------------------------*/
/*
---------------------------------------------------------------
purpose -- return an IVL object whose list K contains all
processes who do not own U(K,K) but own a U(J,K) for some J.
if myid == -1 then
the entire IVL object is created and returned
else
only the portion of the IVL object pertinent
to myid is created and returned
endif
created -- 98may24, cca
---------------------------------------------------------------
*/
IVL *
SolveMap_upperSolveIVL (
SolveMap *solvemap,
int myid,
int msglvl,
FILE *msgFile
) {
int count, K, loc, nblock, nfront, nproc, q ;
int *colids, *heads, *link, *list, *map, *mark, *owners, *rowids ;
IVL *solveIVL ;
/*
---------------
check the input
---------------
*/
if ( solvemap == NULL ) {
fprintf(stderr, "\n fatal error in SolveMap_upperSolveIVL(%p)"
"\n bad input\n", solvemap) ;
exit(-1) ;
}
nfront = solvemap->nfront ;
nproc = solvemap->nproc ;
nblock = solvemap->nblockUpper ;
colids = solvemap->colidsUpper ;
rowids = solvemap->rowidsUpper ;
map = solvemap->mapUpper ;
owners = solvemap->owners ;
/*
------------------------------------
link the (J,K,map(J,K)) triples by K
------------------------------------
*/
heads = IVinit(nfront, -1) ;
link = IVinit(nblock, -1) ;
for ( loc = 0 ; loc < nblock ; loc++ ) {
K = colids[loc] ;
link[loc] = heads[K] ;
heads[K] = loc ;
}
list = IVinit(nproc, -1) ;
mark = IVinit(nproc, -1) ;
if ( msglvl > 1 ) {
fprintf(msgFile, "\n\n linked triples by columns of U") ;
for ( K = 0 ; K < nfront ; K++ ) {
if ( heads[K] != -1 ) {
fprintf(msgFile, "\n %d :", K) ;
for ( loc = heads[K] ; loc != -1 ; loc = link[loc] ) {
fprintf(msgFile, " <%d,%d>", rowids[loc], map[loc]) ;
}
}
}
}
/*
-------------------------------
initialize the solve IVL object
-------------------------------
*/
solveIVL = IVL_new() ;
IVL_init1(solveIVL, IVL_CHUNKED, nfront) ;
/*
-------------------------------------------------------------
fill the solve IVL object
list K contains the process ids that do not own K but use X_K
-------------------------------------------------------------
*/
for ( K = 0 ; K < nfront ; K++ ) {
if ( myid == -1 || owners[K] == myid ) {
mark[owners[K]] = K ;
count = 0 ;
if ( msglvl > 1 ) {
fprintf(msgFile, "\n list for %d :", K) ;
}
for ( loc = heads[K] ; loc != -1 ; loc = link[loc] ) {
q = map[loc] ;
if ( msglvl > 1 ) {
fprintf(msgFile, " <%d,%d>", rowids[loc], q) ;
}
if ( mark[q] != K ) {
mark[q] = K ;
list[count++] = q ;
if ( msglvl > 1 ) {
fprintf(msgFile, "*") ;
}
}
}
if ( count > 0 ) {
IVqsortUp(count, list) ;
IVL_setList(solveIVL, K, count, list) ;
}
}
}
/*
------------------------
free the working storage
------------------------
*/
IVfree(heads) ;
IVfree(link) ;
IVfree(list) ;
IVfree(mark) ;
return(solveIVL) ; }
/*--------------------------------------------------------------------*/
/*
-----------------------------------------------------
purpose -- return an IV object whose entry J contains
the number of all processes who do not own U(J,J)
but own a U(J,K) for some J.
if myid == -1 then
all entries in the vector are filled
else
all those entries pertinent to myid are filled
endif
created -- 98mar19, cca
-----------------------------------------------------
*/
IV *
SolveMap_upperAggregateIV (
SolveMap *solvemap,
int myid,
int msglvl,
FILE *msgFile
) {
int count, J, loc, nblock, nfront, nproc, q ;
int *aggcounts, *colids, *heads, *link, *map,
*mark, *owners, *rowids ;
IV *aggIV ;
/*
---------------
check the input
---------------
*/
if ( solvemap == NULL ) {
fprintf(stderr, "\n fatal error in SolveMap_upperAggregateIVL(%p)"
"\n bad input\n", solvemap) ;
exit(-1) ;
}
nfront = solvemap->nfront ;
nproc = solvemap->nproc ;
nblock = solvemap->nblockUpper ;
colids = solvemap->rowidsUpper ;
rowids = solvemap->rowidsUpper ;
map = solvemap->mapUpper ;
owners = solvemap->owners ;
/*
------------------------
initialize the IV object
------------------------
*/
aggIV = IV_new() ;
IV_init(aggIV, nfront, NULL) ;
aggcounts = IV_entries(aggIV) ;
IVzero(nfront, aggcounts) ;
/*
------------------------------------
link the (J,K,map(J,K)) triples by J
------------------------------------
*/
heads = IVinit(nfront, -1) ;
link = IVinit(nblock, -1) ;
for ( loc = 0 ; loc < nblock ; loc++ ) {
J = rowids[loc] ;
link[loc] = heads[J] ;
heads[J] = loc ;
}
if ( msglvl > 1 ) {
fprintf(msgFile, "\n\n linked triples by rows of U") ;
for ( J = 0 ; J < nfront ; J++ ) {
if ( heads[J] != -1 ) {
fprintf(msgFile, "\n %d :", J) ;
for ( loc = heads[J] ; loc != -1 ; loc = link[loc] ) {
fprintf(msgFile, " <%d,%d>", colids[loc], map[loc]) ;
}
}
}
}
/*
-------------------------
fill the aggcounts vector
-------------------------
*/
mark = IVinit(nproc, -1) ;
for ( J = 0 ; J < nfront ; J++ ) {
if ( myid == -1 || owners[J] == myid ) {
mark[owners[J]] = J ;
if ( msglvl > 1 ) {
fprintf(msgFile, "\n list for %d :", J) ;
}
for ( loc = heads[J], count = 0 ; loc != -1 ; loc = link[loc] ) {
q = map[loc] ;
if ( msglvl > 1 ) {
fprintf(msgFile, " <%d,%d>", colids[loc], q) ;
}
if ( mark[q] != J ) {
mark[q] = J ;
count++ ;
if ( msglvl > 1 ) {
fprintf(msgFile, "*") ;
}
}
}
aggcounts[J] = count ;
}
}
/*
------------------------
free the working storage
------------------------
*/
IVfree(heads) ;
IVfree(link) ;
IVfree(mark) ;
return(aggIV) ; }
/*--------------------------------------------------------------------*/
/*
--------------------------------------------------------
purpose -- return an IV object whose J'th entry contains
the number of processes who do not own L(J,J) but own
a L(J,I) for some I < J.
if myid == -1 then
all entries in the vector are filled
else
all those entries pertinent to myid are filled
endif
created -- 98mar20, cca
--------------------------------------------------------
*/
IV *
SolveMap_lowerAggregateIV (
SolveMap *solvemap,
int myid,
int msglvl,
FILE *msgFile
) {
int count, J, loc, nblock, nfront, nproc, q ;
int *aggcounts, *colids, *heads, *link, *map,
*mark, *owners, *rowids ;
IV *aggIV ;
/*
---------------
check the input
---------------
*/
if ( solvemap == NULL ) {
fprintf(stderr, "\n fatal error in SolveMap_lowerAggregateIV(%p)"
"\n bad input\n", solvemap) ;
exit(-1) ;
}
nfront = solvemap->nfront ;
nproc = solvemap->nproc ;
if ( solvemap->symmetryflag == SPOOLES_NONSYMMETRIC ) {
nblock = solvemap->nblockLower ;
colids = solvemap->colidsLower ;
rowids = solvemap->rowidsLower ;
map = solvemap->mapLower ;
} else {
nblock = solvemap->nblockUpper ;
colids = solvemap->rowidsUpper ;
rowids = solvemap->colidsUpper ;
map = solvemap->mapUpper ;
}
owners = solvemap->owners ;
/*
------------------------------------
link the (J,I,map(J,I)) triples by J
------------------------------------
*/
heads = IVinit(nfront, -1) ;
link = IVinit(nblock, -1) ;
for ( loc = 0 ; loc < nblock ; loc++ ) {
J = rowids[loc] ;
link[loc] = heads[J] ;
heads[J] = loc ;
}
if ( msglvl > 1 ) {
fprintf(msgFile, "\n\n linked triples by rows of L or U^T") ;
for ( J = 0 ; J < nfront ; J++ ) {
if ( heads[J] != -1 ) {
fprintf(msgFile, "\n %d :", J) ;
for ( loc = heads[J] ; loc != -1 ; loc = link[loc] ) {
fprintf(msgFile, " <%d,%d>", colids[loc], map[loc]) ;
}
}
}
}
mark = IVinit(nproc, -1) ;
/*
----------------------------------
initialize the aggregate IV object
----------------------------------
*/
aggIV = IV_new() ;
IV_init(aggIV, nfront, NULL) ;
aggcounts = IV_entries(aggIV) ;
IVzero(nfront, aggcounts) ;
/*
-------------------------------------------------------------
fill the aggregate IV object
aggcounts[J] = # of processors besides owner of L_{J,J}
that own some L_{J,I}
-------------------------------------------------------------
*/
for ( J = 0 ; J < nfront ; J++ ) {
if ( myid == -1 || owners[J] == myid ) {
mark[owners[J]] = J ;
if ( msglvl > 1 ) {
fprintf(msgFile, "\n list for %d :", J) ;
}
for ( loc = heads[J], count = 0 ; loc != -1 ; loc = link[loc] ) {
q = map[loc] ;
if ( msglvl > 1 ) {
fprintf(msgFile, " <%d,%d>", colids[loc], q) ;
}
if ( mark[q] != J ) {
mark[q] = J ;
count++ ;
if ( msglvl > 1 ) {
fprintf(msgFile, "*") ;
}
}
}
aggcounts[J] = count ;
}
}
/*
------------------------
free the working storage
------------------------
*/
IVfree(heads) ;
IVfree(link) ;
IVfree(mark) ;
return(aggIV) ; }
/*--------------------------------------------------------------------*/
/*
-------------------------------------------------------
purpose -- return an IVL object whose list J contains
the processes who do not own L(J,J) but own a L(K,J)
for some K > J.
if myid == -1 then
the entire IVL object is created and returned
else
only the portion of the IVL object pertinent
to myid is created and returned
endif
created -- 98may24, cca
-------------------------------------------------------
*/
IVL *
SolveMap_lowerSolveIVL (
SolveMap *solvemap,
int myid,
int msglvl,
FILE *msgFile
) {
int count, I, J, loc, nblock, nfront, nproc, q ;
int *colids, *heads, *link, *list, *map, *mark, *owners, *rowids ;
IVL *solveIVL ;
/*
---------------
check the input
---------------
*/
if ( solvemap == NULL ) {
fprintf(stderr, "\n fatal error in SolveMap_lowerSolveIVL(%p)"
"\n bad input\n", solvemap) ;
exit(-1) ;
}
nfront = solvemap->nfront ;
nproc = solvemap->nproc ;
if ( solvemap->symmetryflag == SPOOLES_NONSYMMETRIC ) {
nblock = solvemap->nblockLower ;
rowids = solvemap->rowidsLower ;
colids = solvemap->colidsLower ;
map = solvemap->mapLower ;
} else {
nblock = solvemap->nblockUpper ;
rowids = solvemap->colidsUpper ;
colids = solvemap->rowidsUpper ;
map = solvemap->mapUpper ;
}
owners = solvemap->owners ;
/*
------------------------
initialize the IV object
------------------------
*/
solveIVL = IVL_new() ;
IVL_init1(solveIVL, IVL_CHUNKED, nfront) ;
/*
------------------------------------
link the (J,I,map(J,I)) triples by I
------------------------------------
*/
heads = IVinit(nfront, -1) ;
link = IVinit(nblock, -1) ;
for ( loc = 0 ; loc < nblock ; loc++ ) {
I = colids[loc] ;
link[loc] = heads[I] ;
heads[I] = loc ;
}
if ( msglvl > 1 ) {
fprintf(msgFile, "\n\n linked triples by columns of L or U^T") ;
for ( I = 0 ; I < nfront ; I++ ) {
if ( heads[I] != -1 ) {
fprintf(msgFile, "\n %d :", I) ;
for ( loc = heads[I] ; loc != -1 ; loc = link[loc] ) {
fprintf(msgFile, " <%d,%d>", rowids[loc], map[loc]) ;
}
}
}
}
/*
-------------------------
fill the solve IVL object
-------------------------
*/
list = IVinit(nproc, -1) ;
mark = IVinit(nproc, -1) ;
for ( I = 0 ; I < nfront ; I++ ) {
if ( myid == -1 || owners[I] == myid ) {
mark[owners[I]] = I ;
if ( msglvl > 1 ) {
fprintf(msgFile, "\n list for %d :", I) ;
}
for ( loc = heads[I], count = 0 ; loc != -1 ; loc = link[loc] ) {
q = map[loc] ;
if ( msglvl > 1 ) {
fprintf(msgFile, " <%d,%d>", rowids[loc], q) ;
}
if ( mark[q] != I ) {
mark[q] = I ;
list[count++] = q ;
if ( msglvl > 1 ) {
fprintf(msgFile, "*") ;
}
}
}
if ( count > 0 ) {
IVqsortUp(count, list) ;
IVL_setList(solveIVL, I, count, list) ;
}
}
}
/*
------------------------
free the working storage
------------------------
*/
IVfree(heads) ;
IVfree(link) ;
IVfree(mark) ;
IVfree(list) ;
return(solveIVL) ; }
/*--------------------------------------------------------------------*/
syntax highlighted by Code2HTML, v. 0.9.1