/* init.c */
#include "../FrontMtx.h"
#define MYDEBUG 0
/*--------------------------------------------------------------------*/
/*
------------------------------------------------------------------
purpose -- basic initializer
frontETree -- ETree object that stores the front tree
symbfacIVL -- IVL object that stores the symbolic factorization
manager -- SubMtxManager object to manage SubMtx objects
type -- type of entries
SPOOLES_REAL --> real
SPOOLES_COMPLEX --> complex
symmetryflag -- symmetry flag,
SPOOLES_SYMMETRIC --> symmetric structure and entries
SPOOLES_HERMITIAN --> hermitian (complex only)
SPOOLES_NONSYMMETRIC --> nonsymmetric entries
sparsityflag -- flag to specify dense or sparse fronts
FRONTMTX_DENSE_FRONTS --> dense fronts
FRONTMTX_SPARSE_FRONTS --> sparse fronts
pivotingflag -- flag to specify pivoting enabled
SPOOLES_NO_PIVOTING --> pivoting not used
SPOOLES_PIVOTING --> pivoting used
in a multithreaded environment, we need to protect the critical
section where data is allocated. we use a lockflag to do this.
in a serial or distributed environment, use lockflag = 0.
lockflag -- flag to specify lock status
NO_LOCK --> mutex lock is not allocated or initialized
LOCK_IN_PROCESS --> mutex lock is allocated and
it can synchronize only threads in this process.
LOCK_OVER_ALL_PROCESSES --> mutex lock is allocated and
it can synchronize only threads in this and other processes.
in a distributed environment we need to specify which process
owns each front. when we can preallocate data structures
(when there is no pivoting and dense fronts are stored) we
need each process to determine what parts of the data it
can allocate and set up. in a serial or multithreaded
environment, use ownersIV = NULL.
ownersIV -- map from fronts to owning processes
myid -- id of this process.
submatrices (be they lower, diagonal, block diagonal, upper)
are stored in SubMtx objects. the management of these objects,
(allocation and deallocation) is managed by the SubMtxManager
manager object.
manager -- SubMtxManager object to handle the submatrices
created -- 98may04, cca
------------------------------------------------------------------
*/
void
FrontMtx_init (
FrontMtx *frontmtx,
ETree *frontETree,
IVL *symbfacIVL,
int type,
int symmetryflag,
int sparsityflag,
int pivotingflag,
int lockflag,
int myid,
IV *ownersIV,
SubMtxManager *manager,
int msglvl,
FILE *msgFile
) {
SubMtx *mtx ;
int J, nbytes, nentD, nentL, nentU, neqns, nD, nent, nfront, nU ;
int *bndwghts, *nodwghts, *owners, *vtxToFront ;
/*
---------------
check the input
---------------
*/
if ( frontmtx == NULL || frontETree == NULL || symbfacIVL == NULL
|| (ownersIV != NULL && myid < 0)
|| manager == NULL ) {
fprintf(stderr,
"\n fatal error in FrontMtx_init()"
"\n frontmtx %p, frontETree %p, symbfacIVL %p"
"\n myid %d, ownersIV %p, manager %p"
"\n bad input\n",
frontmtx, frontETree, symbfacIVL, myid, ownersIV, manager) ;
exit(-1) ;
}
if ( type != SPOOLES_REAL && type != SPOOLES_COMPLEX ) {
fprintf(stderr, "\n fatal error in FrontMtx_init()"
"\n type %d must be SPOOLES_REAL or SPOOLES_COMPLEX\n",
type) ;
exit(-1) ;
}
if ( type == SPOOLES_REAL &&
! (symmetryflag == SPOOLES_SYMMETRIC
|| symmetryflag == SPOOLES_NONSYMMETRIC) ) {
fprintf(stderr,
"\n fatal error in FrontMtx_init()"
"\n type is real"
"\n symmetryflag is not SPOOLES_SYMMETRIC or SPOOLES_NONSYMMETRIC") ;
exit(-1) ;
}
if ( type == SPOOLES_COMPLEX &&
! (symmetryflag == SPOOLES_SYMMETRIC
|| symmetryflag == SPOOLES_HERMITIAN
|| symmetryflag == SPOOLES_NONSYMMETRIC) ) {
fprintf(stderr,
"\n fatal error in FrontMtx_init()"
"\n type is real, symmetryflag is not SPOOLES_SYMMETRIC,"
"\n SPOOLES_HERMITIAN or SPOOLES_NONSYMMETRIC") ;
exit(-1) ;
}
if ( ! ( pivotingflag == SPOOLES_PIVOTING
|| pivotingflag == SPOOLES_NO_PIVOTING) ) {
fprintf(stderr, "\n fatal error in FrontMtx_init()"
"\n pivotingflag must be SPOOLES_PIVOTING or SPOOLES_NO_PIVOTING\n") ;
exit(-1) ;
}
if ( ! (lockflag == NO_LOCK
|| lockflag == LOCK_IN_PROCESS
|| lockflag == LOCK_OVER_ALL_PROCESSES) ) {
fprintf(stderr,
"\n fatal error in FrontMtx_init()"
"\n invalid lockflag, must be NO_LOCK"
"\n LOCK_IN_PROCESS or LOCK_OVER_ALL_PROCESSES") ;
exit(-1) ;
}
nfront = frontETree->nfront ;
neqns = frontETree->nvtx ;
nodwghts = ETree_nodwghts(frontETree) ;
bndwghts = ETree_bndwghts(frontETree) ;
vtxToFront = ETree_vtxToFront(frontETree) ;
if ( ownersIV != NULL ) {
owners = IV_entries(ownersIV) ;
} else {
owners = NULL ;
}
/*
----------------------
set the default fields
----------------------
*/
FrontMtx_setDefaultFields(frontmtx) ;
/*
---------------------
set the scalar fields
---------------------
*/
frontmtx->nfront = nfront ;
frontmtx->neqns = neqns ;
frontmtx->type = type ;
frontmtx->symmetryflag = symmetryflag ;
frontmtx->sparsityflag = sparsityflag ;
frontmtx->pivotingflag = pivotingflag ;
frontmtx->dataMode = FRONTMTX_1D_MODE ;
/*
---------------------------------------------------------------
set the front tree ETree and symbolic factorization IVL objects
---------------------------------------------------------------
*/
frontmtx->tree = frontETree->tree ;
frontmtx->frontETree = frontETree ;
frontmtx->symbfacIVL = symbfacIVL ;
/*
----------------------------
set the frontsizes IV object
----------------------------
*/
frontmtx->frontsizesIV = IV_new() ;
if ( FRONTMTX_IS_PIVOTING(frontmtx) ) {
IV_init(frontmtx->frontsizesIV, nfront, NULL) ;
IV_fill(frontmtx->frontsizesIV, 0) ;
} else {
IV_init(frontmtx->frontsizesIV, nfront, nodwghts) ;
}
/*
------------------------------------------------------
set the row and column adjacency objects, if necessary
------------------------------------------------------
*/
if ( FRONTMTX_IS_PIVOTING(frontmtx) ) {
frontmtx->coladjIVL = IVL_new() ;
IVL_init1(frontmtx->coladjIVL, IVL_CHUNKED, nfront) ;
if ( FRONTMTX_IS_NONSYMMETRIC(frontmtx) ) {
frontmtx->rowadjIVL = IVL_new() ;
IVL_init1(frontmtx->rowadjIVL, IVL_CHUNKED, nfront) ;
}
}
/*
---------------------------------
allocate the five pointer vectors
---------------------------------
*/
ALLOCATE(frontmtx->p_mtxDJJ, struct _SubMtx *, nfront) ;
ALLOCATE(frontmtx->p_mtxUJJ, struct _SubMtx *, nfront) ;
ALLOCATE(frontmtx->p_mtxUJN, struct _SubMtx *, nfront) ;
for ( J = 0 ; J < nfront ; J++ ) {
frontmtx->p_mtxDJJ[J] = NULL ;
frontmtx->p_mtxUJJ[J] = NULL ;
frontmtx->p_mtxUJN[J] = NULL ;
}
if ( FRONTMTX_IS_NONSYMMETRIC(frontmtx) ) {
ALLOCATE(frontmtx->p_mtxLJJ, struct _SubMtx *, nfront) ;
ALLOCATE(frontmtx->p_mtxLNJ, struct _SubMtx *, nfront) ;
for ( J = 0 ; J < nfront ; J++ ) {
frontmtx->p_mtxLJJ[J] = NULL ;
frontmtx->p_mtxLNJ[J] = NULL ;
}
}
/*
-------------------------------
set the submatrix manager field
-------------------------------
*/
frontmtx->manager = manager ;
/*
-------------------------------------
initialize submatrices where possible
-------------------------------------
*/
if ( ! FRONTMTX_IS_PIVOTING(frontmtx)
&& FRONTMTX_IS_DENSE_FRONTS(frontmtx) ) {
double *entries ;
int ii, jj, ncol, nrow ;
int *firstlocs, *sizes ;
nentD = nentL = nentU = 0 ;
for ( J = 0 ; J < nfront ; J++ ) {
if ( owners == NULL || owners[J] == myid ) {
nD = nodwghts[J] ;
nU = bndwghts[J] ;
if ( msglvl > 3 ) {
fprintf(msgFile, "\n\n J %d, nD %d, nU %d", J, nD, nU) ;
fflush(msgFile) ;
}
/*
---------------
diagonal matrix
---------------
*/
nbytes = SubMtx_nbytesNeeded(type, SUBMTX_DIAGONAL,
nD, nD, nD) ;
mtx = SubMtxManager_newObjectOfSizeNbytes(manager, nbytes) ;
SubMtx_init(mtx, type, SUBMTX_DIAGONAL, J, J, nD, nD, nD);
SubMtx_diagonalInfo(mtx, &ncol, &entries) ;
SubMtx_zero(mtx) ;
nentD += nD ;
frontmtx->p_mtxDJJ[J] = mtx ;
if ( msglvl > 3 ) {
fprintf(msgFile,
"\n diagonal (%d,%d) matrix %p, %d entries, %d bytes",
J, J, mtx, nD, nbytes) ;
fflush(msgFile) ;
}
if ( (nent = (nD*(nD-1))/2) > 0 ) {
/*
-------------------------------------------
U_{J,J} and possibly lower L_{J,J} matrices
-------------------------------------------
*/
nbytes = SubMtx_nbytesNeeded(type, SUBMTX_DENSE_SUBCOLUMNS,
nD, nD, nent) ;
mtx = SubMtxManager_newObjectOfSizeNbytes(manager, nbytes) ;
SubMtx_init(mtx, type, SUBMTX_DENSE_SUBCOLUMNS,
J, J, nD, nD, nent);
SubMtx_denseSubcolumnsInfo(mtx, &ncol, &nent, &firstlocs,
&sizes, &entries) ;
for ( jj = 0 ; jj < ncol ; jj++ ) {
firstlocs[jj] = 0 ;
sizes[jj] = jj ;
}
SubMtx_zero(mtx) ;
nentU += nent ;
frontmtx->p_mtxUJJ[J] = mtx ;
if ( msglvl > 3 ) {
fprintf(msgFile,
"\n upper (%d,%d) matrix %p, %d entries, %d bytes",
J, J, mtx, nent, nbytes) ;
fflush(msgFile) ;
}
if ( FRONTMTX_IS_NONSYMMETRIC(frontmtx) ) {
nbytes = SubMtx_nbytesNeeded(type, SUBMTX_DENSE_SUBROWS,
nD, nD, nent) ;
mtx = SubMtxManager_newObjectOfSizeNbytes(manager,
nbytes);
SubMtx_init(mtx, type, SUBMTX_DENSE_SUBROWS,
J, J, nD, nD, nent);
SubMtx_denseSubrowsInfo(mtx, &nrow, &nent, &firstlocs,
&sizes, &entries) ;
for ( ii = 0 ; ii < nrow ; ii++ ) {
firstlocs[ii] = 0 ;
sizes[ii] = ii ;
}
SubMtx_zero(mtx) ;
nentL += nent ;
frontmtx->p_mtxLJJ[J] = mtx ;
if ( msglvl > 3 ) {
fprintf(msgFile,
"\n lower (%d,%d) matrix %p, %d entries, %d bytes",
J, J, mtx, nent, nbytes) ;
fflush(msgFile) ;
}
}
}
if ( (nent = nD*nU) > 0 ) {
/*
-----------------------------------------------
U_{J,bnd{J}} and possibly L_{bnd{J},J} matrices
-----------------------------------------------
*/
nbytes = SubMtx_nbytesNeeded(type, SUBMTX_DENSE_COLUMNS,
nD, nU, nent) ;
mtx = SubMtxManager_newObjectOfSizeNbytes(manager, nbytes) ;
SubMtx_init(mtx, type, SUBMTX_DENSE_COLUMNS,
J, nfront, nD, nU, nent);
SubMtx_zero(mtx) ;
nentU += nent ;
frontmtx->p_mtxUJN[J] = mtx ;
if ( msglvl > 3 ) {
fprintf(msgFile,
"\n upper (%d,%d) matrix %p, %d entries, %d bytes",
J, nfront, mtx, nent, nbytes) ;
fflush(msgFile) ;
}
if ( FRONTMTX_IS_NONSYMMETRIC(frontmtx) ) {
nbytes = SubMtx_nbytesNeeded(type, SUBMTX_DENSE_ROWS,
nU, nD, nent);
mtx = SubMtxManager_newObjectOfSizeNbytes(manager,
nbytes);
SubMtx_init(mtx, type, SUBMTX_DENSE_ROWS,
nfront, J, nU, nD, nent);
SubMtx_zero(mtx) ;
nentL += nent ;
frontmtx->p_mtxLNJ[J] = mtx ;
if ( msglvl > 3 ) {
fprintf(msgFile,
"\n lower (%d,%d) matrix %p, %d entries, %d bytes",
nfront, J, mtx, nent, nbytes) ;
fflush(msgFile) ;
}
}
}
}
}
frontmtx->nentD = nentD ;
frontmtx->nentL = nentL ;
frontmtx->nentU = nentU ;
}
if ( lockflag == LOCK_OVER_ALL_PROCESSES
|| lockflag == LOCK_IN_PROCESS ) {
/*
-----------------
allocate the lock
-----------------
*/
frontmtx->lock = Lock_new() ;
Lock_init(frontmtx->lock, lockflag) ;
}
/*
--------------------------------------
set the patch-and-go information field
--------------------------------------
*/
frontmtx->patchinfo = NULL ;
if ( msglvl > 3 ) {
fprintf(msgFile, "\n\n frontmtx->lock = %p", frontmtx->lock) ;
fflush(msgFile) ;
}
return ; }
/*--------------------------------------------------------------------*/
syntax highlighted by Code2HTML, v. 0.9.1