46 #ifndef _PEXSI_PSELINV_HPP_
47 #define _PEXSI_PSELINV_HPP_
69 typedef std::vector<bool> bitMask;
70 typedef std::map<bitMask , std::vector<Int> > bitMaskSet;
96 GridType( MPI_Comm Bcomm,
int nprow,
int npcol );
156 LBlock& operator = (
const LBlock& LB) {
164 friend std::ostream& operator<<(std::ostream& out,
const LBlock& vec)
166 out <<
"(" << vec.blockIdx <<
", " << vec.numRow <<
", " << vec.numCol <<std::endl<<
"rows " << vec.rows <<std::endl<<
"nzval " <<std::endl<< vec.nzval <<
")";
204 UBlock& operator = (
const UBlock& UB) {
213 friend std::ostream& operator<<(std::ostream& out,
const UBlock& vec)
215 out <<
"(" << vec.blockIdx <<
", " << vec.numRow <<
", " << vec.numCol <<std::endl<<
"cols " << vec.cols <<std::endl<<
"nzval " <<std::endl<< vec.nzval <<
")";
229 {
return g->mpirank; }
233 {
return g->mpirank / g->numProcCol; }
237 {
return g->mpirank % g->numProcCol; }
242 {
return bnum % g->numProcRow; }
247 {
return bnum % g->numProcCol; }
252 {
return (i%g->numProcRow) * g->numProcCol + j%g->numProcCol; }
257 {
return bnum / g->numProcRow; }
262 {
return bnum / g->numProcCol; }
267 {
return iLocal * g->numProcRow +
MYROW( g ); }
272 {
return jLocal * g->numProcCol +
MYCOL( g ); }
277 {
return (a%b) ? ( a/b + 1 ) : ( a/b ); }
281 {
return s->superIdx[i]; }
286 {
return s->superPtr[bnum]; }
293 {
return s->superPtr[bnum]; }
298 {
return s->superPtr[bnum+1] - s->superPtr[bnum]; }
302 {
return s->superPtr.m() - 1; }
307 {
return s->superIdx.m(); }
331 namespace LBlockMask{
344 Int
inline serialize(LBlock<T>& val, std::ostream& os,
const std::vector<Int>& mask){
346 if(mask[i]==1) serialize(val.blockIdx, os, mask); i++;
347 if(mask[i]==1) serialize(val.numRow, os, mask); i++;
348 if(mask[i]==1) serialize(val.numCol, os, mask); i++;
349 if(mask[i]==1) serialize(val.rows, os, mask); i++;
350 if(mask[i]==1) serialize(val.nzval, os, mask); i++;
355 Int
inline deserialize(LBlock<T>& val, std::istream& is,
const std::vector<Int>& mask){
357 if(mask[i]==1) deserialize(val.blockIdx, is, mask); i++;
358 if(mask[i]==1) deserialize(val.numRow, is, mask); i++;
359 if(mask[i]==1) deserialize(val.numCol, is, mask); i++;
360 if(mask[i]==1) deserialize(val.rows, is, mask); i++;
361 if(mask[i]==1) deserialize(val.nzval, is, mask); i++;
381 namespace UBlockMask{
393 Int
inline serialize(UBlock<T>& val, std::ostream& os,
const std::vector<Int>& mask){
395 if(mask[i]==1) serialize(val.blockIdx, os, mask); i++;
396 if(mask[i]==1) serialize(val.numRow, os, mask); i++;
397 if(mask[i]==1) serialize(val.numCol, os, mask); i++;
398 if(mask[i]==1) serialize(val.cols, os, mask); i++;
399 if(mask[i]==1) serialize(val.nzval, os, mask); i++;
404 Int
inline deserialize(UBlock<T>& val, std::istream& is,
const std::vector<Int>& mask){
406 if(mask[i]==1) deserialize(val.blockIdx, is, mask); i++;
407 if(mask[i]==1) deserialize(val.numRow, is, mask); i++;
408 if(mask[i]==1) deserialize(val.numCol, is, mask); i++;
409 if(mask[i]==1) deserialize(val.cols, is, mask); i++;
410 if(mask[i]==1) deserialize(val.nzval, is, mask); i++;
496 std::vector<std::vector<Int> > ColBlockIdx_;
497 std::vector<std::vector<Int> > RowBlockIdx_;
498 std::vector<std::vector<LBlock<T> > > L_;
499 std::vector<std::vector<UBlock<T> > > U_;
501 std::vector<std::vector<Int> > workingSet_;
520 bitMaskSet maskSendToBelow_;
521 std::vector<MPI_Comm> commSendToBelow_;
522 std::vector<MPI_Comm*> commSendToBelowPtr_;
523 std::vector<Int> commSendToBelowRoot_;
524 std::vector<bitMask *> commSendToBelowMaskPtr_;
525 std::vector<bitMask> commSendToBelowMask_;
528 bitMaskSet maskRecvFromBelow_;
529 std::vector<MPI_Comm> commRecvFromBelow_;
530 std::vector<MPI_Comm*> commRecvFromBelowPtr_;
531 std::vector<Int> commRecvFromBelowRoot_;
533 std::vector<bitMask *> commRecvFromBelowMaskPtr_;
534 std::vector<bitMask> commRecvFromBelowMask_;
538 bitMaskSet maskSendToRight_;
539 std::vector<MPI_Comm> commSendToRight_;
540 std::vector<MPI_Comm*> commSendToRightPtr_;
541 std::vector<Int> commSendToRightRoot_;
542 std::vector<bitMask *> commSendToRightMaskPtr_;
543 std::vector<bitMask> commSendToRightMask_;
551 SELINV_TAG_U_CONTENT,
553 SELINV_TAG_L_CONTENT,
556 SELINV_TAG_D_CONTENT,
564 struct SuperNodeBufferType{
567 std::vector<Int> RowLocalPtr;
568 std::vector<Int> BlockIdxLocal;
569 std::vector<char> SstrLcolSend;
570 std::vector<char> SstrUrowSend;
571 std::vector<char> SstrLcolRecv;
572 std::vector<char> SstrUrowRecv;
573 Int SizeSstrLcolSend;
574 Int SizeSstrUrowSend;
575 Int SizeSstrLcolRecv;
576 Int SizeSstrUrowRecv;
581 SuperNodeBufferType():
589 SuperNodeBufferType(Int &pIndex) :
600 inline void SelInvIntra_Collectives(Int lidx);
603 inline void SelInvIntra_P2p(Int lidx);
610 inline void GetWorkSet(std::vector<Int> & snodeEtree, std::vector<std::vector<Int> > & WSet);
613 inline void UnpackData(SuperNodeBufferType & snode, std::vector<
LBlock<T> > & LcolRecv, std::vector<
UBlock<T> > & UrowRecv);
616 inline void ComputeDiagUpdate(SuperNodeBufferType & snode);
619 inline void SendRecvCD_UpdateU(std::vector<SuperNodeBufferType > & arrSuperNodes, Int stepSuper);
622 void getMaxCommunicatorSizes();
625 void ConstructCommunicators_Collectives(Int lidx);
629 void DestructCommunicators_Collectives( );
644 Int NumCol()
const {
return super_ -> superIdx.m(); }
646 Int NumSuper()
const {
return super_ ->superPtr.m() - 1; }
655 std::vector< std::vector<Int> > & ColBlockIdx() {
return ColBlockIdx_; }
656 std::vector< std::vector<Int> > & RowBlockIdx() {
return RowBlockIdx_; }
657 std::vector<Int> & ColBlockIdx(Int jLocal) {
return ColBlockIdx_[jLocal]; }
658 std::vector<Int> & RowBlockIdx(Int iLocal) {
return RowBlockIdx_[iLocal]; }
663 Int
NumBlockL( Int jLocal )
const {
return L_[jLocal].size(); }
667 Int
NumBlockU( Int iLocal )
const {
return U_[iLocal].size(); }
678 std::vector<LBlock<T> >&
L( Int jLocal ) {
return L_[jLocal]; }
682 std::vector<UBlock<T> >&
U( Int iLocal ) {
return U_[iLocal]; }
686 std::vector<std::vector<int> >&
WorkingSet( ) {
return workingSet_; }
690 Int
CountSendToRight(Int ksup) { Int count= std::count (isSendToRight_.VecData(ksup), isSendToRight_.VecData(ksup) + grid_->numProcCol,
true);
return (isSendToRight_(
MYCOL(grid_),ksup)?count-1:count); }
694 Int
CountRecvFromBelow(Int ksup) { Int count= std::count (isRecvFromBelow_.VecData(ksup), isRecvFromBelow_.VecData(ksup) + grid_->numProcRow,
true);
return (isRecvFromBelow_(
MYROW(grid_),ksup)?count-1:count); }
698 Int
CountSendToCrossDiagonal(Int ksup) { Int count= std::count (isSendToCrossDiagonal_.VecData(ksup), isSendToCrossDiagonal_.VecData(ksup) + grid_->numProcCol,
true);
return ((isSendToCrossDiagonal_(
MYCOL(grid_),ksup) &&
MYROW(grid_)==
PROW(ksup,grid_))?count-1:count); }
702 Int
CountRecvFromCrossDiagonal(Int ksup) { Int count= std::count (isRecvFromCrossDiagonal_.VecData(ksup), isRecvFromCrossDiagonal_.VecData(ksup) + grid_->numProcRow,
true);
return ((isRecvFromCrossDiagonal_(
MYROW(grid_),ksup) &&
MYCOL(grid_)==
PCOL(ksup,grid_))?count-1:count); }
709 void GetEtree(std::vector<Int> & etree_supno );
907 void GetColumn ( Int colIdx,
NumVec<T>& col );
968 #endif //_PEXSI_PSELINV_HPP_
void PreSelInv()
PreSelInv prepares the structure in L_ and U_ so that SelInv only involves matrix-matrix multiplicati...
Definition: pselinv_impl.hpp:3819
Int blockIdx
Block index (supernodal index)
Definition: pselinv.hpp:187
void PMatrixToDistSparseMatrix(DistSparseMatrix< T > &A)
PMatrixToDistSparseMatrix converts the PMatrix into a distributed compressed sparse column matrix for...
Definition: pselinv_impl.hpp:4245
A thin interface for passing parameters to set the SuperLU options.
Definition: superlu_dist_internal.hpp:62
std::vector< LBlock< T > > & L(Int jLocal)
L returns the vector of nonzero L blocks for the local block column jLocal.
Definition: pselinv.hpp:678
SuperNodeType describes mapping between supernode and column, the permutation information, and potentially the elimination tree (not implemented here).
Definition: pselinv.hpp:122
Interface with SuperLU_Dist (version 3.0 and later)
void GetNegativeInertia(Real &inertia)
GetNegativeInertia computes the negative inertia of a PMatrix. This can be used to estimate e...
Int PROW(Int bnum, const GridType *g)
PROW returns the processor row that the bnum-th block (supernode) belongs to.
Definition: pselinv.hpp:241
Thin interface to LAPACK.
Int numRow
Number of nonzero rows.
Definition: pselinv.hpp:142
Int MYPROC(const GridType *g)
MYPROC returns the current processor rank.
Definition: pselinv.hpp:228
Int NumLocalBlockCol() const
NumLocalBlockCol returns the total number of block columns.
Definition: pselinv.hpp:649
Int PCOL(Int bnum, const GridType *g)
PCOL returns the processor column that the bnum-th block (supernode) belongs to.
Definition: pselinv.hpp:246
Int LBi(Int bnum, const GridType *g)
LBi returns the local block number on the processor at processor row PROW( bnum, g )...
Definition: pselinv.hpp:256
void SelInv_Hybrid(Int threshold)
Hybrid version of the selected inversion. This file is obsolete.
Int SuperSize(Int bnum, const SuperNodeType *s)
SuperSize returns the size of the block bnum.
Definition: pselinv.hpp:297
Int LBj(Int bnum, const GridType *g)
LBj returns the local block number on the processor at processor column PCOL( bnum, g ).
Definition: pselinv.hpp:261
Int GBj(Int jLocal, const GridType *g)
GBj returns the global block number from a local block number in the column direction.
Definition: pselinv.hpp:271
void ConstructCommunicationPattern_Hybrid()
ConstructCommunicationPattern_Hybrid constructs the communication pattern to be used later in the sel...
Int numCol
Number of nonzero columns.
Definition: pselinv.hpp:145
Int NumCol(const SuperNodeType *s)
NumCol returns the total number of columns for a supernodal partiiton.
Definition: pselinv.hpp:306
Int CountRecvFromCrossDiagonal(Int ksup)
CountRecvFromCrossDiagonal returns the number of cross diagonal processors with which current process...
Definition: pselinv.hpp:702
void GetDiagonal(NumVec< T > &diag)
GetDiagonal extracts the diagonal elements of the PMatrix.
Definition: pselinv_impl.hpp:4204
std::vector< std::vector< int > > & WorkingSet()
WorkingSet returns the ordered list of supernodes which could be done in parallel.
Definition: pselinv.hpp:686
Implementation of the parallel SelInv.
void SelInv_Collectives()
Collective communication version of the selected inversion.
Int FirstBlockRow(Int bnum, const SuperNodeType *s)
FirstBlockRow returns the first column of a block bnum. Note: the functionality of FirstBlockRow is e...
Definition: pselinv.hpp:292
void GetEtree(std::vector< Int > &etree_supno)
GetEtree computes the supernodal elimination tree to be used later in the pipelined selected inversio...
Definition: pselinv_impl.hpp:875
const GridType * Grid() const
Grid returns the GridType structure of the current PMatrix.
Definition: pselinv.hpp:670
LBlock stores a nonzero block in the lower triangular part or the diagonal part in PSelInv...
Definition: pselinv.hpp:136
Int PNUM(Int i, Int j, const GridType *g)
PNUM returns the processor rank that the bnum-th block (supernode) belongs to.
Definition: pselinv.hpp:251
void ConstructCommunicationPattern()
ConstructCommunicationPattern constructs the communication pattern to be used later in the selected i...
Definition: pselinv_impl.hpp:3356
Sparse matrix and Distributed sparse matrix in compressed column format.
void ConstructCommunicationPattern_P2p()
ConstructCommunicationPattern_P2p constructs the communication pattern to be used later in the select...
Definition: pselinv_impl.hpp:3364
Int GBi(Int iLocal, const GridType *g)
GBi returns the global block number from a local block number in the row direction.
Definition: pselinv.hpp:266
GridType is the PSelInv way of defining the grid.
Definition: pselinv.hpp:85
Interface with MPI to facilitate communication.
Int MYCOL(const GridType *g)
MYCOL returns my processor column.
Definition: pselinv.hpp:236
Int CountSendToCrossDiagonal(Int ksup)
CountSendToCrossDiagonal returns the number of cross diagonal processors with which current processor...
Definition: pselinv.hpp:698
Int NumLocalBlockRow() const
NumLocalBlockRow returns the total number of block rows.
Definition: pselinv.hpp:652
Various utility subroutines.
const SuperNodeType * SuperNode() const
SuperNode returns the supernodal partition of the current PMatrix.
Definition: pselinv.hpp:674
IntNumVec rows
Dimension numRow * 1, index (0-based) for the number of nonzero rows.
Definition: pselinv.hpp:148
Int numCol
Number of nonzero columns.
Definition: pselinv.hpp:193
UBlock stores a nonzero block in the upper triangular part in PSelInv.
Definition: pselinv.hpp:184
std::vector< UBlock< T > > & U(Int iLocal)
U returns the vector of nonzero U blocks for the local block row iLocal.
Definition: pselinv.hpp:682
Int FirstBlockCol(Int bnum, const SuperNodeType *s)
FirstBlockCol returns the first column of a block bnum.
Definition: pselinv.hpp:285
Int NnzLocal()
NnzLocal computes the number of nonzero elements (L and U) saved locally.
Definition: pselinv_impl.hpp:5086
Int blockIdx
Block index (supernodal index)
Definition: pselinv.hpp:139
PMatrix contains the main data structure and the computational routine for the parallel selected inve...
Definition: pselinv.hpp:482
Int BlockIdx(Int i, const SuperNodeType *s)
BlockIdx returns the block index of a column i.
Definition: pselinv.hpp:280
Int numRow
Number of nonzero rows.
Definition: pselinv.hpp:190
Int CountSendToRight(Int ksup)
CountSendToRight returns the number of processors to the right of current processor with which it has...
Definition: pselinv.hpp:690
Int CEILING(Int a, Int b)
CEILING is used for computing the storage space for local number of blocks.
Definition: pselinv.hpp:276
NumMat< T > nzval
Dimension numRow * numCol, nonzero elements.
Definition: pselinv.hpp:199
void ConstructCommunicationPattern_Collectives()
ConstructCommunicationPattern_Collectives constructs the communication pattern to be used later in th...
Int NumSuper(const SuperNodeType *s)
NumSuper returns the total number of supernodes.
Definition: pselinv.hpp:301
void SelInv_P2p()
Point-to-point version of the selected inversion.
Definition: pselinv_impl.hpp:3786
Int CountRecvFromBelow(Int ksup)
CountRecvFromBelow returns the number of processors below the current processor from which it receive...
Definition: pselinv.hpp:694
Int NumBlockL(Int jLocal) const
NumBlockL returns the number of nonzero L blocks for the local block column jLocal.
Definition: pselinv.hpp:663
NumMat< T > nzval
Dimension numRow * numCol, nonzero elements.
Definition: pselinv.hpp:151
IntNumVec cols
Dimension numRow * 1, index (0-based) for the number of nonzero rows.
Definition: pselinv.hpp:196
Int NumBlockU(Int iLocal) const
NumBlockU returns the number of nonzero U blocks for the local block row iLocal.
Definition: pselinv.hpp:667
void SelInv()
SelInv is the main function for the selected inversion.
Definition: pselinv_impl.hpp:3778
LongInt Nnz()
Nnz computes the total number of nonzero elements in the PMatrix.
Definition: pselinv_impl.hpp:5117
void PMatrixToDistSparseMatrix2(const DistSparseMatrix< T > &A, DistSparseMatrix< T > &B)
PMatrixToDistSparseMatrix2 is a more efficient version which performs the same job as PMatrixToDistSp...
Definition: pselinv_impl.hpp:4801
DistSparseMatrix describes a Sparse matrix in the compressed sparse column format (CSC) and distribut...
Definition: sparse_matrix.hpp:91
Int MYROW(const GridType *g)
MYROW returns my processor row.
Definition: pselinv.hpp:232