SuperLU_DIST
4.0
superlu_dist on CPU and GPU clusters
|
Gets matrix permutation. More...
Functions | |
float | get_perm_c_parmetis (SuperMatrix *A, int_t *perm_r, int_t *perm_c, int nprocs_i, int noDomains, int_t **sizes, int_t **fstVtxSep, gridinfo_t *grid, MPI_Comm *metis_comm) |
Gets matrix permutation.
– Distributed symbolic factorization auxialiary routine (version 2.1) – Lawrence Berkeley National Lab, Univ. of California Berkeley - July 2003 INRIA France - January 2004 Laura Grigori
November 1, 2007
float get_perm_c_parmetis | ( | SuperMatrix * | A, |
int_t * | perm_r, | ||
int_t * | perm_c, | ||
int | nprocs_i, | ||
int | noDomains, | ||
int_t ** | sizes, | ||
int_t ** | fstVtxSep, | ||
gridinfo_t * | grid, | ||
MPI_Comm * | metis_comm | ||
) |
Purpose
GET_PERM_C_PARMETIS obtains a permutation matrix Pc, by applying a graph partitioning algorithm to the symmetrized graph A+A'. The multilevel graph partitioning algorithm used is the ParMETIS_V3_NodeND routine available in the parallel graph partitioning package parMETIS.
The number of independent sub-domains noDomains computed by this algorithm has to be a power of 2. Hence noDomains is the larger number power of 2 that is smaller than nprocs_i, where nprocs_i = nprow * npcol is the number of processors used in SuperLU_DIST.
Arguments
A (input) SuperMatrix* Matrix A in A*X=B, of dimension (A->nrow, A->ncol). The number of the linear equations is A->nrow. Matrix A is distributed in NRformat_loc format.
perm_r (input) int_t* Row permutation vector of size A->nrow, which defines the permutation matrix Pr; perm_r[i] = j means row i of A is in position j in Pr*A.
perm_c (output) int_t* Column permutation vector of size A->ncol, which defines the permutation matrix Pc; perm_c[i] = j means column i of A is in position j in A*Pc.
nprocs_i (input) int* Number of processors the input matrix is distributed on in a block row format. It corresponds to number of processors used in SuperLU_DIST.
noDomains (input) int*, must be power of 2 Number of independent domains to be computed by the graph partitioning algorithm. ( noDomains <= nprocs_i )
sizes (output) int_t**, of size 2 * noDomains Returns pointer to an array containing the number of nodes for each sub-domain and each separator. Separators are stored from left to right. Memory for the array is allocated in this routine.
fstVtxSep (output) int_t**, of size 2 * noDomains Returns pointer to an array containing first node for each sub-domain and each separator. Memory for the array is allocated in this routine.
Return value
< 0, number of bytes allocated on return from the symbolic factorization.0, number of bytes allocated when out of memory.