PEXSI
 All Classes Namespaces Files Functions Variables Typedefs Pages
pselinv.hpp
Go to the documentation of this file.
1 /*
2  Copyright (c) 2012 The Regents of the University of California,
3  through Lawrence Berkeley National Laboratory.
4 
5  Authors: Lin Lin and Mathias Jacquelin
6 
7  This file is part of PEXSI. All rights reserved.
8 
9  Redistribution and use in source and binary forms, with or without
10  modification, are permitted provided that the following conditions are met:
11 
12  (1) Redistributions of source code must retain the above copyright notice, this
13  list of conditions and the following disclaimer.
14  (2) Redistributions in binary form must reproduce the above copyright notice,
15  this list of conditions and the following disclaimer in the documentation
16  and/or other materials provided with the distribution.
17  (3) Neither the name of the University of California, Lawrence Berkeley
18  National Laboratory, U.S. Dept. of Energy nor the names of its contributors may
19  be used to endorse or promote products derived from this software without
20  specific prior written permission.
21 
22  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
23  ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
24  WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
25  DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
26  ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
27  (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
28  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
29  ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
30  (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
31  SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32 
33  You are under no obligation whatsoever to provide any bug fixes, patches, or
34  upgrades to the features, functionality or performance of the source code
35  ("Enhancements") to anyone; however, if you choose to make your Enhancements
36  available either publicly, or directly to Lawrence Berkeley National
37  Laboratory, without imposing a separate written license agreement for such
38  Enhancements, then you hereby grant the following license: a non-exclusive,
39  royalty-free perpetual license to install, use, modify, prepare derivative
40  works, incorporate into other computer software, distribute, and sublicense
41  such enhancements or derivative works thereof, in binary and source code form.
42 */
46 #ifndef _PEXSI_PSELINV_HPP_
47 #define _PEXSI_PSELINV_HPP_
48 
49 // *********************************************************************
50 // Common utilities
51 // *********************************************************************
52 
53 #include "pexsi/environment.hpp"
54 #include "pexsi/NumVec.hpp"
55 #include "pexsi/NumMat.hpp"
56 #include "pexsi/sparse_matrix.hpp"
57 
59 #include "pexsi/mpi_interf.hpp"
60 #include "pexsi/utility.hpp"
61 #include "pexsi/blas.hpp"
62 #include "pexsi/lapack.hpp"
63 
64 #include "pexsi/TreeBcast.hpp"
65 
66 
67 #include <set>
68 
69 
70 //#define IDX_TO_TAG(lidx,tag) (SELINV_TAG_COUNT*(lidx)+(tag))
71 #define IDX_TO_TAG(lidx,tag,max) ((SELINV_TAG_COUNT*(lidx%(max+1))+(tag)))
72 #define IDX_TO_TAG2(sidx,lidx,tag) (SELINV_TAG_COUNT*(sidx)+(tag))
73 #define TAG_TO_IDX(tag,typetag) (((tag)-(typetag))/SELINV_TAG_COUNT)
74 
75 
76 //#define LIST_BARRIER
77 //#define ALL_BARRIER
78 
79 namespace PEXSI{
80 
81  enum MSGTYPE {LSIZE=0,LROWSIZE,USIZE,UCOLSIZE,LCONTENT,LROWCONTENT,UCONTENT,UCOLCONTENT,MSGCOUNT};
82 
83  struct ULComparator {
84  IntNumVec & lookup;
85  ULComparator(IntNumVec & v):lookup(v){}
86  bool operator() (int i,int j) { return ((lookup)(i)<(lookup)(j));}
87  };
88 
89 
90 
91  typedef std::vector<bool> bitMask;
92  typedef std::map<bitMask , std::vector<Int> > bitMaskSet;
93 
94 
95 
96 
97 
108 
109  // Member functions to setup the default value
111  };
112 
113 
114 
115 
116 
117  /**********************************************************************
118  * Basic PSelInv data structure
119  **********************************************************************/
120 
129  struct GridType{
130  // Data
131  MPI_Comm comm;
132  MPI_Comm rowComm;
133  MPI_Comm colComm;
134  Int mpirank;
135  Int mpisize;
136  Int numProcRow;
137  Int numProcCol;
138 
139  // Member function
140  GridType( MPI_Comm Bcomm, int nprow, int npcol );
141  ~GridType();
142  };
143 
167  IntNumVec perm;
168  IntNumVec permInv;
169  IntNumVec perm_r;
170  IntNumVec permInv_r;
171  IntNumVec superIdx;
172  IntNumVec superPtr;
173  IntNumVec etree;
174  };
175 
176 
181  template<typename T>
182  struct LBlock{
183  // Variables
185  Int blockIdx;
186 
188  Int numRow;
189 
191  Int numCol;
192 
195 
196 
199 
200  // Member functions;
201  LBlock() {
202  blockIdx = -1; numRow = 0; numCol =0;
203  nzval.Resize(0,0);
204  }
205  ~LBlock() {}
206  LBlock& operator = (const LBlock& LB) {
207  blockIdx = LB.blockIdx;
208  numRow = LB.numRow;
209  numCol = LB.numCol;
210  rows = LB.rows;
211  nzval = LB.nzval;
212  return *this;
213  }
214  friend std::ostream& operator<<(std::ostream& out, const LBlock& vec) // output
215  {
216  out << "(" << vec.blockIdx << ", " << vec.numRow << ", " << vec.numCol <<std::endl<< "rows " << vec.rows <<std::endl<< "nzval " <<std::endl<< vec.nzval << ")";
217  return out;
218  }
219 
220 
221  };
222 
233  template<typename T>
234  struct UBlock{
235  // Variables
237  Int blockIdx;
238 
240  Int numRow;
241 
243  Int numCol;
244 
247 
250 
251  // Member functions;
252  UBlock() {
253  blockIdx = -1; numRow = 0; numCol =0;
254  }
255  ~UBlock() {}
256  UBlock& operator = (const UBlock& UB) {
257  blockIdx = UB.blockIdx;
258  numRow = UB.numRow;
259  numCol = UB.numCol;
260  cols = UB.cols;
261  nzval = UB.nzval;
262  return *this;
263  }
264 
265  friend std::ostream& operator<<(std::ostream& out, const UBlock& vec) // output
266  {
267  out << "(" << vec.blockIdx << ", " << vec.numRow << ", " << vec.numCol <<std::endl<< "cols " << vec.cols <<std::endl<< "nzval " <<std::endl<< vec.nzval << ")";
268  return out;
269  }
270  };
271 
272  // *********************************************************************
273  // SuperLU style utility functions
274  //
275  // The SuperLU style macros are defined here as inline functions
276  // so that the code is more portable.
277  // *********************************************************************
278 
280  inline Int MYPROC( const GridType* g )
281  { return g->mpirank; }
282 
284  inline Int MYROW( const GridType* g )
285  { return g->mpirank / g->numProcCol; }
286 
288  inline Int MYCOL( const GridType* g )
289  { return g->mpirank % g->numProcCol; }
290 
293  inline Int PROW( Int bnum, const GridType* g )
294  { return bnum % g->numProcRow; }
295 
298  inline Int PCOL( Int bnum, const GridType* g )
299  { return bnum % g->numProcCol; }
300 
303  inline Int PNUM( Int i, Int j, const GridType* g )
304  { return (i%g->numProcRow) * g->numProcCol + j%g->numProcCol; }
305 
308  inline Int LBi( Int bnum, const GridType* g )
309  { return bnum / g->numProcRow; }
310 
313  inline Int LBj( Int bnum, const GridType* g)
314  { return bnum / g->numProcCol; }
315 
318  inline Int GBi( Int iLocal, const GridType* g )
319  { return iLocal * g->numProcRow + MYROW( g ); }
320 
323  inline Int GBj( Int jLocal, const GridType* g )
324  { return jLocal * g->numProcCol + MYCOL( g ); }
325 
328  inline Int CEILING( Int a, Int b )
329  { return (a%b) ? ( a/b + 1 ) : ( a/b ); }
330 
332  inline Int BlockIdx( Int i, const SuperNodeType *s )
333  { return s->superIdx[i]; }
334 
337  inline Int FirstBlockCol( Int bnum, const SuperNodeType *s )
338  { return s->superPtr[bnum]; }
339 
340 
344  inline Int FirstBlockRow( Int bnum, const SuperNodeType *s )
345  { return s->superPtr[bnum]; }
346 
347 
349  inline Int SuperSize( Int bnum, const SuperNodeType *s )
350  { return s->superPtr[bnum+1] - s->superPtr[bnum]; }
351 
353  inline Int NumSuper( const SuperNodeType *s )
354  { return s->superPtr.m() - 1; }
355 
358  inline Int NumCol( const SuperNodeType *s )
359  { return s->superIdx.m(); }
360 
361 
362  // *********************************************************************
363  // Serialize / Deserialize
364  // *********************************************************************
365 
366  // L part
367 
383  namespace LBlockMask{
384  enum {
385  BLOCKIDX,
386  NUMROW,
387  NUMCOL,
388  ROWS,
389  NZVAL,
390  TOTAL_NUMBER
391  };
392  }
393 
394 
395  template<typename T>
396  Int inline serialize(LBlock<T>& val, std::ostream& os, const std::vector<Int>& mask){
397  if(mask[LBlockMask::BLOCKIDX]==1) serialize(val.blockIdx, os, mask);
398  if(mask[LBlockMask::NUMROW ]==1) serialize(val.numRow, os, mask);
399  if(mask[LBlockMask::NUMCOL ]==1) serialize(val.numCol, os, mask);
400  if(mask[LBlockMask::ROWS ]==1) serialize(val.rows, os, mask);
401  if(mask[LBlockMask::NZVAL ]==1) serialize(val.nzval, os, mask);
402  return 0;
403  }
404 
405  template<typename T>
406  Int inline deserialize(LBlock<T>& val, std::istream& is, const std::vector<Int>& mask){
407  if(mask[LBlockMask::BLOCKIDX]==1) deserialize(val.blockIdx, is, mask);
408  if(mask[LBlockMask::NUMROW ]==1) deserialize(val.numRow, is, mask);
409  if(mask[LBlockMask::NUMCOL ]==1) deserialize(val.numCol, is, mask);
410  if(mask[LBlockMask::ROWS ]==1) deserialize(val.rows, is, mask);
411  if(mask[LBlockMask::NZVAL ]==1) deserialize(val.nzval, is, mask);
412  return 0;
413  }
414 
415  // U part
431  namespace UBlockMask{
432  enum {
433  BLOCKIDX,
434  NUMROW,
435  NUMCOL,
436  COLS,
437  NZVAL,
438  TOTAL_NUMBER
439  };
440  }
441 
442  template<typename T>
443  Int inline serialize(UBlock<T>& val, std::ostream& os, const std::vector<Int>& mask){
444  Int i = 0;
445  if(mask[i]==1) serialize(val.blockIdx, os, mask); i++;
446  if(mask[i]==1) serialize(val.numRow, os, mask); i++;
447  if(mask[i]==1) serialize(val.numCol, os, mask); i++;
448  if(mask[i]==1) serialize(val.cols, os, mask); i++;
449  if(mask[i]==1) serialize(val.nzval, os, mask); i++;
450  return 0;
451  }
452 
453  template<typename T>
454  Int inline deserialize(UBlock<T>& val, std::istream& is, const std::vector<Int>& mask){
455  Int i = 0;
456  if(mask[i]==1) deserialize(val.blockIdx, is, mask); i++;
457  if(mask[i]==1) deserialize(val.numRow, is, mask); i++;
458  if(mask[i]==1) deserialize(val.numCol, is, mask); i++;
459  if(mask[i]==1) deserialize(val.cols, is, mask); i++;
460  if(mask[i]==1) deserialize(val.nzval, is, mask); i++;
461  return 0;
462  }
463 
464  /**********************************************************************
465  * Main data structure in PSelInv: PMatrix
466  **********************************************************************/
467 
529 
530  template<typename T>
531  class PMatrix{
532  public:
535 
536  static PMatrix<T> * Create(const GridType * pGridType, const SuperNodeType * pSuper, const PSelInvOptions * pSelInvOpt , const SuperLUOptions * pLuOpt);
537  static PMatrix<T> * Create(const SuperLUOptions * pLuOpt);
538 
539  public:
540  // This is the tag used for mpi communication for selinv
541 
542  enum{
543  SELINV_TAG_U_SIZE,
544  SELINV_TAG_U_CONTENT,
545  SELINV_TAG_L_SIZE,
546  SELINV_TAG_L_CONTENT,
547  SELINV_TAG_L_REDUCE,
548  SELINV_TAG_D_SIZE,
549  SELINV_TAG_D_CONTENT,
550  SELINV_TAG_D_REDUCE,
551  SELINV_TAG_L_SIZE_CD,
552  SELINV_TAG_L_CONTENT_CD,
553  SELINV_TAG_COUNT
554  };
555 
556 
557  protected:
558  // *********************************************************************
559  // Variables
560  // *********************************************************************
561  // Data variables
562 
563  const GridType* grid_;
564 
565  const SuperNodeType* super_;
566 
567  const PSelInvOptions * options_;
568  const SuperLUOptions * optionsLU_;
569 
570  Int limIndex_;
571  Int maxTag_;
572 
573  std::vector<std::vector<Int> > ColBlockIdx_;
574  std::vector<std::vector<Int> > RowBlockIdx_;
575  std::vector<std::vector<LBlock<T> > > L_;
576  std::vector<std::vector<UBlock<T> > > U_;
577 
578  std::vector<std::vector<Int> > workingSet_;
579 
580 
581  // Communication variables
582  BolNumMat isSendToBelow_;
583  BolNumMat isSendToRight_;
584  BolNumVec isSendToDiagonal_;
585  BolNumMat isSendToCrossDiagonal_;
586 
587  BolNumMat isRecvFromBelow_;
588  BolNumVec isRecvFromAbove_;
589  BolNumVec isRecvFromLeft_;
590  BolNumMat isRecvFromCrossDiagonal_;
591 
592 
593 #ifdef NEW_BCAST
594  std::vector<TreeBcast2<T> *> fwdToBelowTree2_;
595  std::vector<TreeBcast2<T> *> fwdToRightTree2_;
596 #endif
597  std::vector<TreeBcast *> fwdToBelowTree_;
598  std::vector<TreeBcast *> fwdToRightTree_;
599  std::vector<TreeReduce<T> *> redToLeftTree_;
600  std::vector<TreeReduce<T> *> redToAboveTree_;
601 
602 
603 
605  NumMat<T> LUpdateBuf;
606  NumMat<T> DiagBuf;
607  std::vector<Int> RowLocalPtr;
608  std::vector<Int> BlockIdxLocal;
609  std::vector<char> SstrLcolSend;
610  std::vector<char> SstrUrowSend;
611  std::vector<char> SstrLcolRecv;
612  std::vector<char> SstrUrowRecv;
613  Int SizeSstrLcolSend;
614  Int SizeSstrUrowSend;
615  Int SizeSstrLcolRecv;
616  Int SizeSstrUrowRecv;
617  Int Index;
618  Int Rank;
619  Int isReady;
620 
621 
623  SizeSstrLcolSend(0),
624  SizeSstrUrowSend(0),
625  SizeSstrLcolRecv(0),
626  SizeSstrUrowRecv(0),
627  Index(0),
628  Rank(0),
629  isReady(0){}
630 
631  SuperNodeBufferType(Int &pIndex) :
632  SizeSstrLcolSend(0),
633  SizeSstrUrowSend(0),
634  SizeSstrLcolRecv(0),
635  SizeSstrUrowRecv(0),
636  Index(pIndex),
637  Rank(0),
638  isReady(0) {}
639 
640  };
641 
642 
644  inline void SelInvIntra_P2p(Int lidx,Int & rank);
645 
647  inline void SelInv_lookup_indexes(SuperNodeBufferType & snode, std::vector<LBlock<T> > & LcolRecv, std::vector<UBlock<T> > & UrowRecv, NumMat<T> & AinvBuf,NumMat<T> & UBuf);
648 
650  inline void GetWorkSet(std::vector<Int> & snodeEtree, std::vector<std::vector<Int> > & WSet);
651 
653  inline void UnpackData(SuperNodeBufferType & snode, std::vector<LBlock<T> > & LcolRecv, std::vector<UBlock<T> > & UrowRecv);
654 
656  inline void ComputeDiagUpdate(SuperNodeBufferType & snode);
657 
659  inline void SendRecvCD_UpdateU(std::vector<SuperNodeBufferType > & arrSuperNodes, Int stepSuper);
660 
661  public:
662  // *********************************************************************
663  // Public member functions
664  // *********************************************************************
665 
666  PMatrix();
667 
668  PMatrix( const GridType* g, const SuperNodeType* s, const PEXSI::PSelInvOptions * o, const PEXSI::SuperLUOptions * oLU );
669 
670  void deallocate();
671 
672 
673  virtual ~PMatrix();
674  PMatrix( const PMatrix & C);
675  PMatrix & operator = ( const PMatrix & C);
676 
677  void Setup( const GridType* g, const SuperNodeType* s, const PEXSI::PSelInvOptions * o, const PEXSI::SuperLUOptions * oLU );
678 
679  Int NumCol() const { return super_ -> superIdx.m(); }
680 
681  Int NumSuper() const { return super_ ->superPtr.m() - 1; }
682 
684  Int NumLocalBlockCol() const { return CEILING( this->NumSuper(), grid_->numProcCol ); }
685 
687  Int NumLocalBlockRow() const { return CEILING( this->NumSuper(), grid_->numProcRow); }
688 
689 
690  std::vector< std::vector<Int> > & ColBlockIdx() { return ColBlockIdx_; }
691  std::vector< std::vector<Int> > & RowBlockIdx() { return RowBlockIdx_; }
692  std::vector<Int> & ColBlockIdx(Int jLocal) { return ColBlockIdx_[jLocal]; }
693  std::vector<Int> & RowBlockIdx(Int iLocal) { return RowBlockIdx_[iLocal]; }
694 
695 
698  Int NumBlockL( Int jLocal ) const { return L_[jLocal].size(); }
699 
702  Int NumBlockU( Int iLocal ) const { return U_[iLocal].size(); }
703 
705  const GridType* Grid() const { return grid_; }
706 
709  const SuperNodeType* SuperNode() const { return super_; }
710 
713  std::vector<LBlock<T> >& L( Int jLocal ) { return L_[jLocal]; }
714 
717  std::vector<UBlock<T> >& U( Int iLocal ) { return U_[iLocal]; }
718 
721  std::vector<std::vector<int> >& WorkingSet( ) { return workingSet_; }
722  //void WorkingSet(const std::vector<std::vector<int> > * pWSet,const std::vector<std::vector<int> > * pWRanks) { pWSet = &workingSet_; pWRanks = &workingRanks_; }
723 
726  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); }
727 
730  Int CountSendToBelow(Int ksup) { Int count= std::count (isSendToBelow_.VecData(ksup), isSendToBelow_.VecData(ksup) + grid_->numProcRow, true); return (isSendToBelow_(MYROW(grid_),ksup)?count-1:count); }
731 
734  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); }
735 
738  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); }
739 
742  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); }
743 
744 
745 
746 
749  void GetEtree(std::vector<Int> & etree_supno );
750 
751 
756  virtual void ConstructCommunicationPattern( );
757 
758 
763 
764 
765 
794  virtual void PreSelInv( );
795 
920  virtual void SelInv( );
921 
923  void SelInv_P2p( );
924 
925 
932  void GetDiagonal( NumVec<T>& diag );
933  void GetColumn ( Int colIdx, NumVec<T>& col );
934 
935 
942 
953  const DistSparseMatrix<T>& A,
954  DistSparseMatrix<T>& B );
955 
956 
967  const DistSparseMatrix<T>& A,
968  DistSparseMatrix<T>& B );
969 
970 
973  Int NnzLocal();
974 
977  LongInt Nnz();
978 
982  void GetNegativeInertia ( Real& inertia );
983 
984 // inline int IdxToTag(Int lidx, Int tag) { return SELINV_TAG_COUNT*(lidx)+(tag);}
985 
986  void DumpLU();
987 
988  void CopyLU( const PMatrix & C);
989  };
990 
991  template<typename T> class PMatrixUnsym;
992 
993 
994 } // namespace PEXSI
995 
996 
997 #include "pexsi/pselinv_impl.hpp"
998 
999 
1000 #endif //_PEXSI_PSELINV_HPP_
virtual void PreSelInv()
PreSelInv prepares the structure in L_ and U_ so that SelInv only involves matrix-matrix multiplicati...
Definition: pselinv_impl.hpp:3711
Environmental variables.
Int blockIdx
Block index (supernodal index)
Definition: pselinv.hpp:237
void PMatrixToDistSparseMatrix(DistSparseMatrix< T > &A)
PMatrixToDistSparseMatrix converts the PMatrix into a distributed compressed sparse column matrix for...
Definition: pselinv_impl.hpp:4242
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:713
SuperNodeType describes mapping between supernode and column, the permutation information, and potentially the elimination tree (not implemented here).
Definition: pselinv.hpp:166
static PMatrix< T > * Create(const GridType *pGridType, const SuperNodeType *pSuper, const PSelInvOptions *pSelInvOpt, const SuperLUOptions *pLuOpt)
Create is a factory method which returns a pointer either to a new PMatrix or to a new PMatrixUnsym d...
Definition: pselinv_impl.hpp:5282
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:293
Thin interface to LAPACK.
Int numRow
Number of nonzero rows.
Definition: pselinv.hpp:188
void UnpackData(SuperNodeBufferType &snode, std::vector< LBlock< T > > &LcolRecv, std::vector< UBlock< T > > &UrowRecv)
UnpackData.
Definition: pselinv_impl.hpp:1096
Int maxPipelineDepth
The maximum pipeline depth.
Definition: pselinv.hpp:107
Int MYPROC(const GridType *g)
MYPROC returns the current processor rank.
Definition: pselinv.hpp:280
Int NumLocalBlockCol() const
NumLocalBlockCol returns the total number of block columns.
Definition: pselinv.hpp:684
Int PCOL(Int bnum, const GridType *g)
PCOL returns the processor column that the bnum-th block (supernode) belongs to.
Definition: pselinv.hpp:298
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:308
Int SuperSize(Int bnum, const SuperNodeType *s)
SuperSize returns the size of the block bnum.
Definition: pselinv.hpp:349
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:313
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:323
Int numCol
Number of nonzero columns.
Definition: pselinv.hpp:191
Int NumCol(const SuperNodeType *s)
NumCol returns the total number of columns for a supernodal partiiton.
Definition: pselinv.hpp:358
Definition: pselinv.hpp:604
Int CountRecvFromCrossDiagonal(Int ksup)
CountRecvFromCrossDiagonal returns the number of cross diagonal processors with which current process...
Definition: pselinv.hpp:742
void GetDiagonal(NumVec< T > &diag)
GetDiagonal extracts the diagonal elements of the PMatrix.
Definition: pselinv_impl.hpp:4114
std::vector< std::vector< int > > & WorkingSet()
WorkingSet returns the ordered list of supernodes which could be done in parallel.
Definition: pselinv.hpp:721
Implementation of the parallel SelInv.
void GetWorkSet(std::vector< Int > &snodeEtree, std::vector< std::vector< Int > > &WSet)
GetWorkSet.
Definition: pselinv_impl.hpp:1271
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:344
Definition: pselinv.hpp:83
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:1202
const GridType * Grid() const
Grid returns the GridType structure of the current PMatrix.
Definition: pselinv.hpp:705
LBlock stores a nonzero block in the lower triangular part or the diagonal part in PSelInv...
Definition: pselinv.hpp:182
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:303
virtual void ConstructCommunicationPattern()
ConstructCommunicationPattern constructs the communication pattern to be used later in the selected i...
Definition: pselinv_impl.hpp:2670
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:2678
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:318
GridType is the PSelInv way of defining the grid.
Definition: pselinv.hpp:129
void SelInvIntra_P2p(Int lidx, Int &rank)
SelInvIntra_P2p.
Definition: pselinv_impl.hpp:1511
Numerical matrix.
Interface with MPI to facilitate communication.
Int MYCOL(const GridType *g)
MYCOL returns my processor column.
Definition: pselinv.hpp:288
Int CountSendToCrossDiagonal(Int ksup)
CountSendToCrossDiagonal returns the number of cross diagonal processors with which current processor...
Definition: pselinv.hpp:738
Int NumLocalBlockRow() const
NumLocalBlockRow returns the total number of block rows.
Definition: pselinv.hpp:687
Various utility subroutines.
const SuperNodeType * SuperNode() const
SuperNode returns the supernodal partition of the current PMatrix.
Definition: pselinv.hpp:709
void PMatrixToDistSparseMatrix_OLD(const DistSparseMatrix< T > &A, DistSparseMatrix< T > &B)
PMatrixToDistSparseMatrix_OLD converts the PMatrix into a distributed compressed sparse column matrix...
Definition: pselinv_impl.hpp:4523
IntNumVec rows
Dimension numRow * 1, index (0-based) for the number of nonzero rows.
Definition: pselinv.hpp:194
Int numCol
Number of nonzero columns.
Definition: pselinv.hpp:243
UBlock stores a nonzero block in the upper triangular part in PSelInv.
Definition: pselinv.hpp:234
std::vector< UBlock< T > > & U(Int iLocal)
U returns the vector of nonzero U blocks for the local block row iLocal.
Definition: pselinv.hpp:717
Int FirstBlockCol(Int bnum, const SuperNodeType *s)
FirstBlockCol returns the first column of a block bnum.
Definition: pselinv.hpp:337
Int NnzLocal()
NnzLocal computes the number of nonzero elements (L and U) saved locally.
Definition: pselinv_impl.hpp:5166
Int blockIdx
Block index (supernodal index)
Definition: pselinv.hpp:185
Thin interface to BLAS.
PMatrix contains the main data structure and the computational routine for the parallel selected inve...
Definition: ngchol_interf.hpp:57
void SendRecvCD_UpdateU(std::vector< SuperNodeBufferType > &arrSuperNodes, Int stepSuper)
SendRecvCD_UpdateU.
Definition: pselinv_impl.hpp:861
Int BlockIdx(Int i, const SuperNodeType *s)
BlockIdx returns the block index of a column i.
Definition: pselinv.hpp:332
Int numRow
Number of nonzero rows.
Definition: pselinv.hpp:240
Int CountSendToRight(Int ksup)
CountSendToRight returns the number of processors to the right of current processor with which it has...
Definition: pselinv.hpp:726
Int CountSendToBelow(Int ksup)
CountSendToBelow returns the number of processors below current processor with which it has to commun...
Definition: pselinv.hpp:730
Int CEILING(Int a, Int b)
CEILING is used for computing the storage space for local number of blocks.
Definition: pselinv.hpp:328
NumMat< T > nzval
Dimension numRow * numCol, nonzero elements.
Definition: pselinv.hpp:249
Int NumSuper(const SuperNodeType *s)
NumSuper returns the total number of supernodes.
Definition: pselinv.hpp:353
Numerical vector.
void SelInv_P2p()
Point-to-point version of the selected inversion.
Definition: pselinv_impl.hpp:3646
Int CountRecvFromBelow(Int ksup)
CountRecvFromBelow returns the number of processors below the current processor from which it receive...
Definition: pselinv.hpp:734
Int NumBlockL(Int jLocal) const
NumBlockL returns the number of nonzero L blocks for the local block column jLocal.
Definition: pselinv.hpp:698
NumMat< T > nzval
Dimension numRow * numCol, nonzero elements.
Definition: pselinv.hpp:198
A thin interface for passing parameters to set the PSelInv options.
Definition: pselinv.hpp:102
IntNumVec cols
Dimension numCol * 1, index (0-based) for the number of nonzero rows.
Definition: pselinv.hpp:246
Int NumBlockU(Int iLocal) const
NumBlockU returns the number of nonzero U blocks for the local block row iLocal.
Definition: pselinv.hpp:702
virtual void SelInv()
SelInv is the main function for the selected inversion.
Definition: pselinv_impl.hpp:3593
LongInt Nnz()
Nnz computes the total number of nonzero elements in the PMatrix.
Definition: pselinv_impl.hpp:5197
void SelInv_lookup_indexes(SuperNodeBufferType &snode, std::vector< LBlock< T > > &LcolRecv, std::vector< UBlock< T > > &UrowRecv, NumMat< T > &AinvBuf, NumMat< T > &UBuf)
SelInv_lookup_indexes.
Definition: pselinv_impl.hpp:504
void ComputeDiagUpdate(SuperNodeBufferType &snode)
ComputeDiagUpdate.
Definition: pselinv_impl.hpp:1162
DistSparseMatrix describes a Sparse matrix in the compressed sparse column format (CSC) and distribut...
Definition: sparse_matrix.hpp:91
PMatrixUnsym contains the main data structure and the computational routine for the parallel selected...
Definition: pselinv.hpp:991
Int MYROW(const GridType *g)
MYROW returns my processor row.
Definition: pselinv.hpp:284