ParGeMSLR
Public Member Functions | Public Attributes | List of all members
pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType > Class Template Reference

The local real ilu preconditioner, only work for sequential CSR matrix. More...

#include <parallel_gemslr.hpp>

Public Member Functions

int SetupPartialILUT (VectorType &x, VectorType &rhs)
 Setup the local partial ILUT. More...
 
 ParallelGemslrClass ()
 The constructor of precondioner class. More...
 
virtual int Clear ()
 Free the current precondioner. More...
 
virtual ~ParallelGemslrClass ()
 The destructor of precondioner class. More...
 
 ParallelGemslrClass (const ParallelGemslrClass< MatrixType, VectorType, DataType > &precond)
 The copy constructor of ParallelGemslrClass. More...
 
 ParallelGemslrClass (ParallelGemslrClass< MatrixType, VectorType, DataType > &&precond)
 The move constructor of ParallelGemslrClass. More...
 
ParallelGemslrClass< MatrixType, VectorType, DataType > & operator= (const ParallelGemslrClass< MatrixType, VectorType, DataType > &precond)
 The operator = of ParallelGemslrClass. More...
 
ParallelGemslrClass< MatrixType, VectorType, DataType > & operator= (ParallelGemslrClass< MatrixType, VectorType, DataType > &&precond)
 The operator = of ParallelGemslrClass. More...
 
int CheckParameters ()
 Check the input parameters, avoid invalid settings. More...
 
virtual int Setup (VectorType &x, VectorType &rhs)
 Setup the precondioner phase. Will be called by the solver if not called directly. More...
 
int SetupPermutation ()
 Setup the permutation of the GeMSLR. More...
 
int SetupPermutationRKway (MatrixType &A, int nlev_setup, int &nlev_max, int &nlev_used, vector_int &map_v, vector_int &mapptr_v)
 Setup the Recursive KWay partition of the GeMSLR. More...
 
int SetupPermutationND (MatrixType &A, int nlev_setup, int &nlev_max, int &nlev_used, vector_int &map_v, vector_int &mapptr_v)
 Setup the Nested Dissection partition of the GeMSLR. More...
 
int SetupPermutationBuildLevelStructure (MatrixType &A, int level_start, vector_int &map_v, vector_int &mapptr_v)
 Setup the level structure of the GeMSLR. More...
 
int SetupBSolve (VectorType &x, VectorType &rhs)
 Setup the solve of B matrices of the GeMSLR. More...
 
int SetupBSolveILUT (VectorType &x, VectorType &rhs, int level)
 Setup the solve of B matrices of the GeMSLR with ILUT. More...
 
int SetupBSolveILUK (VectorType &x, VectorType &rhs, int level)
 Setup the solve of B matrices of the GeMSLR with ILUK. More...
 
int SetupBSolveGemslr (VectorType &x, VectorType &rhs, int level)
 Setup the solve of B matrices of the GeMSLR with GeMSLR. More...
 
int SetupLowRank (VectorType &x, VectorType &rhs)
 Setup the low-rank part of the GeMSLR. More...
 
int SetupLowRankSubspaceIteration (VectorType &x, VectorType &rhs, int level, int option)
 Setup the low-rank part of the GeMSLR with subspace iteration. More...
 
int SetupLowRankNoRestart (VectorType &x, VectorType &rhs, int level, int option)
 Setup the low-rank part of the GeMSLR with arnodi no-restart. More...
 
int SetupLowRankThickRestart (VectorType &x, VectorType &rhs, int level, int option)
 Setup the low-rank part of the GeMSLR with arnodi thick-restart. More...
 
int SetupLowRankBuildLowRank (VectorType &x, VectorType &rhs, DenseMatrixClass< DataType > &V, DenseMatrixClass< DataType > &H, int m, int rank, int level, int option)
 Given the V and H from Arnoldi, compute the final low-rank correction of S^{-1} - C^{-1}. More...
 
virtual int Solve (VectorType &x, VectorType &rhs)
 Solve phase. Call this function after Setup. Solve with cusparse if unified memory/device memory is used. More...
 
int SolveLevelGemslr (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveLevelGemslrU (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveLevelGemslrMul (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveLevelEsmslr (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveLevelEsmslrU (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveLevelEsmslrMul (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveLevelPslr (VectorType &x_out, VectorType &rhs_in, int level, bool doperm)
 Solve starting from a certain level. More...
 
int SolveB (VectorType &x, VectorType &rhs, int option, int level)
 Solve with B on a certain level. More...
 
int SolveApplyLowRankLevel (VectorType &x, VectorType &rhs, int level)
 Apply the low-rank update on a certain level. More...
 
int GetNumRows (int level)
 Get the local number of rows on certain level for the low-rank part. More...
 
int EBFCMatVec (int level, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec function y = G*x = I - Ci*(Ci\x) + Ei*(UBi(LBi(Fi*(Ci\x). Note that alpha and beta are untouched. More...
 
int RAPEBFCMatVec (int level, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec function y = G*x = (2EB^{-1}F + EB^{-1}BB^{-1}F)C^{-1}x. Note that alpha and beta are untouched. More...
 
int SCMatVec (int level, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec function y = G*x = C\S*x-x. Note that alpha and beta are untouched. More...
 
int ACMatVec (int level, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec function y = x-A*M^{-1}*x. Note that alpha and beta are untouched. More...
 
int PCLRMatVec (int level, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec function y = (EsC^{-1})^m*x. Note that alpha and beta are untouched. More...
 
int SchurMatVec (int level, int option, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec with the Schur Complement. More...
 
int RAPMatVec (int level, int option, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec with the Schur Complement via RAP. More...
 
int CMatVec (int level, int option, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec with the C on the current level. More...
 
int BMatVec (int level, int option, char trans, const DataType &alpha, VectorType &x, const DataType &beta, VectorType &y)
 The matvec with the B on the current level. More...
 
int GetSize ()
 Get the size of the problem. More...
 
virtual long int GetNumNonzeros ()
 Get the total number of nonzeros the ILU. More...
 
long int GetNumNonzeros (long int &nnz_bsolver, long int &nnz_lr)
 Get the total number of nonzeros the ILU. More...
 
int PlotPatternGnuPlot (const char *datafilename)
 Plot the pattern of the parallel gemslr using gnuplot. Similar to spy in the MATLAB. Function for testing purpose. More...
 
virtual int SetLocalGemslr (GemslrClass< CsrMatrixClass< DataType >, SequentialVectorClass< DataType >, DataType > &gemslr)
 Setup with parameter array. This is the helper function to set the local gemslr for B solve. More...
 
virtual int SetLowerGemslr (ParallelGemslrClass< MatrixType, VectorType, DataType > &gemslr)
 Setup with parameter array. This is the helper function to set the lower level gemslr for explicite Schur solve. More...
 
virtual int SetWithParameterArray (double *params)
 Setup with parameter array. More...
 
int SetGlobalPrecondOption (int option)
 Set the global precond option. 0: BJ; 1: ISCHUR; 2: ESCHUR; 3: MLEV. More...
 
int SetPartitionOption (int option)
 Set the global partition option. 0: ND; 1: RKway. More...
 
int SetBPartitionOption (int option)
 Set the B part partition option. 0: ND; 1: RKway. More...
 
int SetPermutationOption (int option)
 Set the global permutation option. 0: No; 1: RCM; 2: AMD. More...
 
int SetBPermutationOption (int option)
 Set the B part permutation option. 0: No; 1: RCM; 2: AMD. More...
 
int SetNumLevels (int option)
 Set the global target number of levels of GeMSLR. More...
 
int SetBNumLevels (int option)
 Set the B part target number of levels of GeMSLR. More...
 
int SetNumSubdomains (int option)
 Set the global target number of subdomains on each level of GeMSLR. More...
 
int SetBNumSubdomains (int option)
 Set the B part target number of subdomains on each level of GeMSLR. More...
 
int SetMinimalNumberSubdomains (int option)
 Set the global minimal number of subdomains on each level of GeMSLR. More...
 
int SetBMinimalNumberSubdomains (int option)
 Set the B part minimal number of subdomains on each level of GeMSLR. More...
 
int SetNumberSubdomainsReduceFactor (int option)
 Set the global reduce factor of subdomains on each level of GeMSLR. More...
 
int SetBNumberSubdomainsReduceFactor (int option)
 Set the B part reduce factor of subdomains on each level of GeMSLR. More...
 
int SetSeperatorOption (bool option)
 Set the global separator option of GeMSLR. More...
 
int SetBSeperatorOption (bool option)
 Set the B part separator option of GeMSLR. More...
 
int SetGlobalPartitionOption (bool option)
 Set the global partition option of GeMSLR. More...
 
int SetInnerIterationOption (bool option)
 Set the inner iteration option of GeMSLR. More...
 
template<typename T >
int SetInnerIterationThreshold (T option)
 Set the stop threshold of inner iteration of GeMSLR. More...
 
int SetInnerIterationMaxNumberIterations (int option)
 Set the max number of iterations of inner iteration of GeMSLR. More...
 
int SetSolveOption (int option)
 Set the solve option. 0: additive LU solve; 1: additive L solve; 2: multi-solve. More...
 
int SetLowRankOptionTopLevel (int option)
 Set the low-rank option on the top level. 0: Standard. 1: Thick-restart. More...
 
int SetLowRankOptionOtherLevels (int option)
 Set the low-rank option on other levels. 0: Standard. 1: Thick-restart. More...
 
int SetLowRankOptionA (int option)
 Set the low-rank option on A. 0: Standard. 1: Thick-restart. More...
 
int SetBLowRankOptionTopLevel (int option)
 Set the B part low-rank option on the top level. 0: Standard. 1: Thick-restart. More...
 
int SetBLowRankOptionOtherLevels (int option)
 Set the B part low-rank option on other levels. 0: Standard. 1: Thick-restart. More...
 
int SetLowRankRandomInitGuess (bool option)
 Set if we use random initial guess for Arnoldi. Otherwise we use 1 as initial guess. More...
 
int SetBLowRankRandomInitGuess (bool option)
 Set if we use random initial guess for Arnoldi in the B part. Otherwise we use 1 as initial guess. More...
 
template<typename T >
int SetLowRankFactorTopLevel (T option)
 Set the low-rank factor on the top level. The actuall computed number of low-rank terms is rank * factor >= rank. More...
 
template<typename T >
int SetLowRankFactorOtherLevels (T option)
 Set the low-rank factor on the other levels. The actuall computed number of low-rank terms is rank * factor >= rank. More...
 
template<typename T >
int SetLowRankFactorA (T option)
 Set the low-rank factor on A. The actuall computed number of low-rank terms is rank * factor >= rank. More...
 
template<typename T >
int SetBLowRankFactorTopLevel (T option)
 Set the B part low-rank factor on the top level. The actuall computed number of low-rank terms is rank * factor >= rank. More...
 
template<typename T >
int SetBLowRankFactorOtherLevels (T option)
 Set the B part low-rank factor on the other levels. The actuall computed number of low-rank terms is rank * factor >= rank. More...
 
template<typename T >
int SetLowRankRanksTopLevel (T option)
 Set the target number of low-rank terms on the top level. More...
 
template<typename T >
int SetLowRankRanksOtherLevels (T option)
 Set the target number of low-rank terms on the other levels. More...
 
template<typename T >
int SetLowRankA (T option)
 Set the target number of low-rank terms for A. More...
 
template<typename T >
int SetBLowRankRanksTopLevel (T option)
 Set the B part target number of low-rank terms on the top level. More...
 
template<typename T >
int SetBLowRankRanksOtherLevels (T option)
 Set the B part target number of low-rank terms on the other levels. More...
 
template<typename T >
int SetLowRankArnoldiFactorTopLevel (T option)
 Set the Arnoldi factor on the top level. m steps for arnoldi is rank * rank_factor * arnoldi_factor. More...
 
template<typename T >
int SetLowRankArnoldiFactorOtherLevels (T option)
 Set the Arnoldi factor on other levels. m steps for arnoldi is rank * rank_factor * arnoldi_factor. More...
 
template<typename T >
int SetLowRankArnoldiFactorA (T option)
 Set the Arnoldi factor on A. m steps for arnoldi is rank * rank_factor * arnoldi_factor. More...
 
template<typename T >
int SetBLowRankArnoldiFactorTopLevel (T option)
 Set the B part Arnoldi factor on the top level. m steps for arnoldi is rank * rank_factor * arnoldi_factor. More...
 
template<typename T >
int SetBLowRankArnoldiFactorOtherLevels (T option)
 Set the B part Arnoldi factor on other levels. m steps for arnoldi is rank * rank_factor * arnoldi_factor. More...
 
int SetLowRankMaxNumberIterationsTopLevel (int option)
 Set max restarts of thick-restart Arnoldi on the top level. More...
 
int SetLowRankMaxNumberIterationsOtherLevels (int option)
 Set max restarts of thick-restart Arnoldi on other levels. More...
 
int SetLowRankMaxNumberIterationsA (int option)
 Set max restarts of thick-restart Arnoldi on A. More...
 
int SetBLowRankMaxNumberIterationsTopLevel (int option)
 Set B part max restarts of thick-restart Arnoldi on the top level. More...
 
int SetBLowRankMaxNumberIterationsOtherLevels (int option)
 Set B part max restarts of thick-restart Arnoldi on other levels. More...
 
template<typename T >
int SetLowRankThresholdTopLevel (T option)
 Set max restarts of thick-restart Arnoldi on the top level. More...
 
template<typename T >
int SetLowRankThresholdOtherLevels (T option)
 Set max restarts of thick-restart Arnoldi on other levels. More...
 
template<typename T >
int SetLowRankThresholdA (T option)
 Set max restarts of thick-restart Arnoldi on A. More...
 
template<typename T >
int SetBLowRankThresholdTopLevel (T option)
 Set B part max restarts of thick-restart Arnoldi on the top level. More...
 
template<typename T >
int SetBLowRankThresholdOtherLevels (T option)
 Set B part max restarts of thick-restart Arnoldi on other levels. More...
 
template<typename T >
int SetIluDropTolB (T option)
 Set the threshold for ILUT of the B part. More...
 
template<typename T >
int SetIluDropTolC (T option)
 Set the threshold for ILUT of the last level. More...
 
template<typename T >
int SetIluDropTolS (T option)
 Set the threshold for ILUT of the S part. More...
 
template<typename T >
int SetIluDropTolEF (T option)
 Set the threshold for ILUT of the EF part. More...
 
template<typename T >
int SetBIluDropTolB (T option)
 Set the recursive GeMSLR B part threshold for ILUT of the global B part. More...
 
template<typename T >
int SetBIluDropTolC (T option)
 Set the recursive GeMSLR last level threshold for ILUT of the global B part. More...
 
int SetIluMaxRowNnzB (int option)
 Set the maxinum number of nonzeros for ILUT of the B part. More...
 
int SetIluMaxRowNnzC (int option)
 Set the maxinum number of nonzeros ILUT of the last level. More...
 
int SetIluMaxRowNnzS (int option)
 Set the maxinum number of nonzeros ILUT of the S part. More...
 
int SetBIluMaxRowNnzB (int option)
 Set the recursive GeMSLR B part maxinum number of nonzeros ILUT. More...
 
int SetBIluMaxRowNnzC (int option)
 Set the recursive GeMSLR last level maxinum number of nonzeros ILUT. More...
 
int SetIluFillLevelB (int option)
 Set the fill level for ILUK of the B part. More...
 
int SetIluFillLevelC (int option)
 Set the fill level for ILUK of the last level. More...
 
int SetBIluFillLevelB (int option)
 Set the fill level for ILUK of the B part. More...
 
int SetBIluFillLevelC (int option)
 Set the fill level for ILUK of the last level. More...
 
int SetPolyOrder (int option)
 Set poly order for Poly solve of the B part. More...
 
int SetBPolyOrder (int option)
 Set poly order for Poly solve of the B part of the recursive GeMSLR. More...
 
int SetPreconditionerOption1 (int option)
 Set top level preconditioner. More...
 
int SetPreconditionerOption1Levels (int option)
 Set top level preconditioner apply levels. More...
 
int SetPreconditionerOption2 (int option)
 Set mid level preconditioner. More...
 
int SetPreconditionerOptionB (int option)
 Set all top levels preconditioner. More...
 
int SetSmoothOptionB (int option)
 Set all top levels smoother. More...
 
int SetPreconditionerOptionC (int option)
 Set last level preconditioner. More...
 
int SetIluResidualIters (int residual_iters)
 Set the number of residual iterations of B solves in the setup phase. Set to <= 1 to turn off. More...
 
int SetIluComplexShift (bool complex_shift)
 Set if we turn on the complex shift or not (complex version only). More...
 
virtual int SetSolveLocation (const int &location)
 Set the data location that the preconditioner apply to. More...
 
virtual int MoveData (const int &location)
 Move the preconditioner to another location. Only can be called after Setup. More...
 
int GetSolvePhase ()
 Get the solve phase. More...
 

Public Attributes

int _global_precond_option
 The Schur complement option. More...
 
std::vector< ParallelGemslrLevelClass< MatrixType, VectorType, DataType > > _levs_v
 Vector holding the level struct for all levels starting from the second level. The size of this vector is equal to _nlev_used, and is the number of true levels. More...
 
ParallelGemslrLevelClass< MatrixType, VectorType, DataType > _lev_A
 The structure for the A level low-rank correction. More...
 

Detailed Description

template<class MatrixType, class VectorType, typename DataType>
class pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >

The local real ilu preconditioner, only work for sequential CSR matrix. Matrix type is the type of the matrix, VectorType is the type of the vector. DataType is the data type.

Constructor & Destructor Documentation

◆ ParallelGemslrClass() [1/3]

template<class MatrixType , class VectorType , typename DataType >
template pargemslr::precond_gemslr_csr_par_complexd::ParallelGemslrClass ( )

The constructor of precondioner class.

◆ ~ParallelGemslrClass()

template<class MatrixType , class VectorType , typename DataType >
template pargemslr::precond_gemslr_csr_par_complexd::~ParallelGemslrClass ( )
virtual

The destructor of precondioner class.

◆ ParallelGemslrClass() [2/3]

template<class MatrixType , class VectorType , typename DataType >
pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::ParallelGemslrClass ( const ParallelGemslrClass< MatrixType, VectorType, DataType > &  precond)

The copy constructor of ParallelGemslrClass.

◆ ParallelGemslrClass() [3/3]

template<class MatrixType , class VectorType , typename DataType >
pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::ParallelGemslrClass ( ParallelGemslrClass< MatrixType, VectorType, DataType > &&  precond)

The move constructor of ParallelGemslrClass.

Member Function Documentation

◆ ACMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::ACMatVec ( int  level,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec function y = x-A*M^{-1}*x. Note that alpha and beta are untouched.
We are finding AM^{-1}(I-X)^{-1} = I -> X = I - A^M{-1}. After that, we can apply the preconditioner as M^{-1}(I-X)^{-1}x.

Parameters
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ BMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::BMatVec ( int  level,
int  option,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec with the B on the current level.

Parameters
[in]levelThe matvec level.
[in]optionThe matvec option. 0: standard; 1: use -C_offd instad of C.
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ CheckParameters()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::CheckParameters ( )

Check the input parameters, avoid invalid settings.

◆ Clear()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::Clear ( )
virtual

Free the current precondioner.

Returns
Return error message.

◆ CMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::CMatVec ( int  level,
int  option,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec with the C on the current level.

Parameters
[in]levelThe matvec level.
[in]optionThe matvec option. 0: standard; 1: use -C_offd instad of C.
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ EBFCMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::EBFCMatVec ( int  level,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec function y = G*x = I - Ci*(Ci\x) + Ei*(UBi(LBi(Fi*(Ci\x). Note that alpha and beta are untouched.

Parameters
[in]levelThe matvec level.
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ GetNumNonzeros() [1/2]

template<class MatrixType , class VectorType , typename DataType >
template long int pargemslr::precond_gemslr_csr_par_complexd::GetNumNonzeros ( )
virtual

Get the total number of nonzeros the ILU.

Returns
Return the total number of nonzeros the ILU.

◆ GetNumNonzeros() [2/2]

template<class MatrixType , class VectorType , typename DataType >
template long int pargemslr::precond_gemslr_csr_par_complexd::GetNumNonzeros ( long int &  nnz_bsolver,
long int &  nnz_lr 
)

Get the total number of nonzeros the ILU.

Returns
Return the total number of nonzeros the ILU.

◆ GetNumRows()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::GetNumRows ( int  level)

Get the local number of rows on certain level for the low-rank part.

Returns
Return the length of the matrix.

◆ GetSize()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::GetSize ( )

Get the size of the problem.

Returns
Return the problem size.

◆ GetSolvePhase()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::GetSolvePhase ( )
inline

Get the solve phase.

Returns
Return the solve phase.

◆ MoveData()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::MoveData ( const int &  location)
virtual

Move the preconditioner to another location. Only can be called after Setup.

Parameters
[in]locationThe target solver location.
Returns
Return error message.

◆ operator=() [1/2]

template<class MatrixType , class VectorType , typename DataType >
ParallelGemslrClass< MatrixType, VectorType, DataType > & pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::operator= ( const ParallelGemslrClass< MatrixType, VectorType, DataType > &  precond)

The operator = of ParallelGemslrClass.

◆ operator=() [2/2]

template<class MatrixType , class VectorType , typename DataType >
ParallelGemslrClass< MatrixType, VectorType, DataType > & pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::operator= ( ParallelGemslrClass< MatrixType, VectorType, DataType > &&  precond)

The operator = of ParallelGemslrClass.

◆ PCLRMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::PCLRMatVec ( int  level,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec function y = (EsC^{-1})^m*x. Note that alpha and beta are untouched.

Parameters
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ PlotPatternGnuPlot()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::PlotPatternGnuPlot ( const char *  datafilename)

Plot the pattern of the parallel gemslr using gnuplot. Similar to spy in the MATLAB. Function for testing purpose.

Parameters
[in]datafilenameThe filename of the temp file holding the data.
Returns
Return error message.

◆ RAPEBFCMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::RAPEBFCMatVec ( int  level,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec function y = G*x = (2EB^{-1}F + EB^{-1}BB^{-1}F)C^{-1}x. Note that alpha and beta are untouched.

Parameters
[in]levelThe matvec level.
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ RAPMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::RAPMatVec ( int  level,
int  option,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec with the Schur Complement via RAP.

Parameters
[in]levelThe matvec level.
[in]optionThe matvec option. 0: use B_solve; 1: use B_precond.
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ SchurMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SchurMatVec ( int  level,
int  option,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec with the Schur Complement.

Parameters
[in]levelThe matvec level.
[in]optionThe matvec option. 0: standard; 1: use -C_offd instad of C.
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ SCMatVec()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SCMatVec ( int  level,
char  trans,
const DataType &  alpha,
VectorType &  x,
const DataType &  beta,
VectorType &  y 
)

The matvec function y = G*x = C\S*x-x. Note that alpha and beta are untouched.

Parameters
[in]transWhether or not transpose matrix A.
[in]alphaThe alpha value.
[in]xThe left vector.
[in]betaThe beta value.
[in,out]yThe product vector.
Returns
Return error message.

◆ SetBIluDropTolB()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBIluDropTolB ( option)
inline

Set the recursive GeMSLR B part threshold for ILUT of the global B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBIluDropTolC()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBIluDropTolC ( option)
inline

Set the recursive GeMSLR last level threshold for ILUT of the global B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBIluFillLevelB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBIluFillLevelB ( int  option)
inline

Set the fill level for ILUK of the B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBIluFillLevelC()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBIluFillLevelC ( int  option)
inline

Set the fill level for ILUK of the last level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBIluMaxRowNnzB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBIluMaxRowNnzB ( int  option)
inline

Set the recursive GeMSLR B part maxinum number of nonzeros ILUT.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBIluMaxRowNnzC()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBIluMaxRowNnzC ( int  option)
inline

Set the recursive GeMSLR last level maxinum number of nonzeros ILUT.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankArnoldiFactorOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankArnoldiFactorOtherLevels ( option)
inline

Set the B part Arnoldi factor on other levels. m steps for arnoldi is rank * rank_factor * arnoldi_factor

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankArnoldiFactorTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankArnoldiFactorTopLevel ( option)
inline

Set the B part Arnoldi factor on the top level. m steps for arnoldi is rank * rank_factor * arnoldi_factor

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankFactorOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankFactorOtherLevels ( option)
inline

Set the B part low-rank factor on the other levels. The actuall computed number of low-rank terms is rank * factor >= rank.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankFactorTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankFactorTopLevel ( option)
inline

Set the B part low-rank factor on the top level. The actuall computed number of low-rank terms is rank * factor >= rank.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankMaxNumberIterationsOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankMaxNumberIterationsOtherLevels ( int  option)
inline

Set B part max restarts of thick-restart Arnoldi on other levels.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankMaxNumberIterationsTopLevel()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankMaxNumberIterationsTopLevel ( int  option)
inline

Set B part max restarts of thick-restart Arnoldi on the top level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankOptionOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankOptionOtherLevels ( int  option)
inline

Set the B part low-rank option on other levels. 0: Standard. 1: Thick-restart.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankOptionTopLevel()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankOptionTopLevel ( int  option)
inline

Set the B part low-rank option on the top level. 0: Standard. 1: Thick-restart.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankRandomInitGuess()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankRandomInitGuess ( bool  option)
inline

Set if we use random initial guess for Arnoldi in the B part. Otherwise we use 1 as initial guess.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankRanksOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankRanksOtherLevels ( option)
inline

Set the B part target number of low-rank terms on the other levels.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankRanksTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankRanksTopLevel ( option)
inline

Set the B part target number of low-rank terms on the top level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankThresholdOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankThresholdOtherLevels ( option)
inline

Set B part max restarts of thick-restart Arnoldi on other levels.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBLowRankThresholdTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBLowRankThresholdTopLevel ( option)
inline

Set B part max restarts of thick-restart Arnoldi on the top level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBMinimalNumberSubdomains()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBMinimalNumberSubdomains ( int  option)
inline

Set the B part minimal number of subdomains on each level of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBNumberSubdomainsReduceFactor()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBNumberSubdomainsReduceFactor ( int  option)
inline

Set the B part reduce factor of subdomains on each level of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBNumLevels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBNumLevels ( int  option)
inline

Set the B part target number of levels of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBNumSubdomains()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBNumSubdomains ( int  option)
inline

Set the B part target number of subdomains on each level of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBPartitionOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBPartitionOption ( int  option)
inline

Set the B part partition option. 0: ND; 1: RKway.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBPermutationOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBPermutationOption ( int  option)
inline

Set the B part permutation option. 0: No; 1: RCM; 2: AMD.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBPolyOrder()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBPolyOrder ( int  option)
inline

Set poly order for Poly solve of the B part of the recursive GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetBSeperatorOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetBSeperatorOption ( bool  option)
inline

Set the B part separator option of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetGlobalPartitionOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetGlobalPartitionOption ( bool  option)
inline

Set the global partition option of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetGlobalPrecondOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetGlobalPrecondOption ( int  option)
inline

Set the global precond option. 0: BJ; 1: ISCHUR; 2: ESCHUR; 3: MLEV.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluComplexShift()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluComplexShift ( bool  complex_shift)
inline

Set if we turn on the complex shift or not (complex version only).

Parameters
[in]complex_shiftThe new drop tol for ILUT.
Returns
Return error message.

◆ SetIluDropTolB()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluDropTolB ( option)
inline

Set the threshold for ILUT of the B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluDropTolC()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluDropTolC ( option)
inline

Set the threshold for ILUT of the last level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluDropTolEF()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluDropTolEF ( option)
inline

Set the threshold for ILUT of the EF part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluDropTolS()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluDropTolS ( option)
inline

Set the threshold for ILUT of the S part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluFillLevelB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluFillLevelB ( int  option)
inline

Set the fill level for ILUK of the B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluFillLevelC()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluFillLevelC ( int  option)
inline

Set the fill level for ILUK of the last level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluMaxRowNnzB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluMaxRowNnzB ( int  option)
inline

Set the maxinum number of nonzeros for ILUT of the B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluMaxRowNnzC()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluMaxRowNnzC ( int  option)
inline

Set the maxinum number of nonzeros ILUT of the last level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluMaxRowNnzS()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluMaxRowNnzS ( int  option)
inline

Set the maxinum number of nonzeros ILUT of the S part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetIluResidualIters()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetIluResidualIters ( int  residual_iters)
inline

Set the number of residual iterations of B solves in the setup phase. Set to <= 1 to turn off.

Parameters
[in]residual_itersThe number of residual iterations of B solves.
Returns
Return error message.

◆ SetInnerIterationMaxNumberIterations()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetInnerIterationMaxNumberIterations ( int  option)
inline

Set the max number of iterations of inner iteration of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetInnerIterationOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetInnerIterationOption ( bool  option)
inline

Set the inner iteration option of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetInnerIterationThreshold()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetInnerIterationThreshold ( option)
inline

Set the stop threshold of inner iteration of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLocalGemslr()

template<class MatrixType , class VectorType , typename DataType >
virtual int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLocalGemslr ( GemslrClass< CsrMatrixClass< DataType >, SequentialVectorClass< DataType >, DataType > &  gemslr)
inlinevirtual

Setup with parameter array. This is the helper function to set the local gemslr for B solve.

Parameters
[in]paramsThe parameter array.
Returns
Return error message.

◆ SetLowerGemslr()

template<class MatrixType , class VectorType , typename DataType >
virtual int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowerGemslr ( ParallelGemslrClass< MatrixType, VectorType, DataType > &  gemslr)
inlinevirtual

Setup with parameter array. This is the helper function to set the lower level gemslr for explicite Schur solve.

Parameters
[in]paramsThe parameter array.
Returns
Return error message.

◆ SetLowRankA()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankA ( option)
inline

Set the target number of low-rank terms for A.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankArnoldiFactorA()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankArnoldiFactorA ( option)
inline

Set the Arnoldi factor on A. m steps for arnoldi is rank * rank_factor * arnoldi_factor

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankArnoldiFactorOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankArnoldiFactorOtherLevels ( option)
inline

Set the Arnoldi factor on other levels. m steps for arnoldi is rank * rank_factor * arnoldi_factor

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankArnoldiFactorTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankArnoldiFactorTopLevel ( option)
inline

Set the Arnoldi factor on the top level. m steps for arnoldi is rank * rank_factor * arnoldi_factor

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankFactorA()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankFactorA ( option)
inline

Set the low-rank factor on A. The actuall computed number of low-rank terms is rank * factor >= rank.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankFactorOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankFactorOtherLevels ( option)
inline

Set the low-rank factor on the other levels. The actuall computed number of low-rank terms is rank * factor >= rank.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankFactorTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankFactorTopLevel ( option)
inline

Set the low-rank factor on the top level. The actuall computed number of low-rank terms is rank * factor >= rank.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankMaxNumberIterationsA()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankMaxNumberIterationsA ( int  option)
inline

Set max restarts of thick-restart Arnoldi on A.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankMaxNumberIterationsOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankMaxNumberIterationsOtherLevels ( int  option)
inline

Set max restarts of thick-restart Arnoldi on other levels.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankMaxNumberIterationsTopLevel()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankMaxNumberIterationsTopLevel ( int  option)
inline

Set max restarts of thick-restart Arnoldi on the top level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankOptionA()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankOptionA ( int  option)
inline

Set the low-rank option on A. 0: Standard. 1: Thick-restart.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankOptionOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankOptionOtherLevels ( int  option)
inline

Set the low-rank option on other levels. 0: Standard. 1: Thick-restart.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankOptionTopLevel()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankOptionTopLevel ( int  option)
inline

Set the low-rank option on the top level. 0: Standard. 1: Thick-restart.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankRandomInitGuess()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankRandomInitGuess ( bool  option)
inline

Set if we use random initial guess for Arnoldi. Otherwise we use 1 as initial guess.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankRanksOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankRanksOtherLevels ( option)
inline

Set the target number of low-rank terms on the other levels.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankRanksTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankRanksTopLevel ( option)
inline

Set the target number of low-rank terms on the top level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankThresholdA()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankThresholdA ( option)
inline

Set max restarts of thick-restart Arnoldi on A.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankThresholdOtherLevels()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankThresholdOtherLevels ( option)
inline

Set max restarts of thick-restart Arnoldi on other levels.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetLowRankThresholdTopLevel()

template<class MatrixType , class VectorType , typename DataType >
template<typename T >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetLowRankThresholdTopLevel ( option)
inline

Set max restarts of thick-restart Arnoldi on the top level.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetMinimalNumberSubdomains()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetMinimalNumberSubdomains ( int  option)
inline

Set the global minimal number of subdomains on each level of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetNumberSubdomainsReduceFactor()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetNumberSubdomainsReduceFactor ( int  option)
inline

Set the global reduce factor of subdomains on each level of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetNumLevels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetNumLevels ( int  option)
inline

Set the global target number of levels of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetNumSubdomains()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetNumSubdomains ( int  option)
inline

Set the global target number of subdomains on each level of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPartitionOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPartitionOption ( int  option)
inline

Set the global partition option. 0: ND; 1: RKway.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPermutationOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPermutationOption ( int  option)
inline

Set the global permutation option. 0: No; 1: RCM; 2: AMD.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPolyOrder()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPolyOrder ( int  option)
inline

Set poly order for Poly solve of the B part.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPreconditionerOption1()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPreconditionerOption1 ( int  option)
inline

Set top level preconditioner.

Note
ILU, GeMSLR, Poly.
Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPreconditionerOption1Levels()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPreconditionerOption1Levels ( int  option)
inline

Set top level preconditioner apply levels.

Note
ILU, GeMSLR, Poly.
Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPreconditionerOption2()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPreconditionerOption2 ( int  option)
inline

Set mid level preconditioner.

Note
ILU, GeMSLR, Poly.
Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPreconditionerOptionB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPreconditionerOptionB ( int  option)
inline

Set all top levels preconditioner.

Note
ILU, BJ-ILU. Not used in the sequential version.
Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetPreconditionerOptionC()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetPreconditionerOptionC ( int  option)
inline

Set last level preconditioner.

Note
ILU, BJ-ILU. Not used in the sequential version.
Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetSeperatorOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetSeperatorOption ( bool  option)
inline

Set the global separator option of GeMSLR.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetSmoothOptionB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetSmoothOptionB ( int  option)
inline

Set all top levels smoother.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ SetSolveLocation()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::SetSolveLocation ( const int &  location)
virtual

Set the data location that the preconditioner apply to.

Parameters
[in]locationThe target solver location.
Returns
Return error message.

◆ SetSolveOption()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetSolveOption ( int  option)
inline

Set the solve option. 0: additive LU solve; 1: additive L solve; 2: multi-solve.

Parameters
[in]optionThe new option.
Returns
Return error message.

◆ Setup()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::Setup ( VectorType &  x,
VectorType &  rhs 
)
virtual

Setup the precondioner phase. Will be called by the solver if not called directly.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
Returns
Return error message.

◆ SetupBSolve()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupBSolve ( VectorType &  x,
VectorType &  rhs 
)

Setup the solve of B matrices of the GeMSLR.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
Returns
Return error message.

◆ SetupBSolveGemslr()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupBSolveGemslr ( VectorType &  x,
VectorType &  rhs,
int  level 
)

Setup the solve of B matrices of the GeMSLR with GeMSLR.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe target level.
Returns
Return error message.

◆ SetupBSolveILUK()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupBSolveILUK ( VectorType &  x,
VectorType &  rhs,
int  level 
)

Setup the solve of B matrices of the GeMSLR with ILUK.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe target level.
Returns
Return error message.

◆ SetupBSolveILUT()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupBSolveILUT ( VectorType &  x,
VectorType &  rhs,
int  level 
)

Setup the solve of B matrices of the GeMSLR with ILUT.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe target level.
Returns
Return error message.

◆ SetupLowRank()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupLowRank ( VectorType &  x,
VectorType &  rhs 
)

Setup the low-rank part of the GeMSLR.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
Returns
Return error message.

◆ SetupLowRankBuildLowRank()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupLowRankBuildLowRank ( VectorType &  x,
VectorType &  rhs,
DenseMatrixClass< DataType > &  V,
DenseMatrixClass< DataType > &  H,
int  m,
int  rank,
int  level,
int  option 
)

Given the V and H from Arnoldi, compute the final low-rank correction of S^{-1} - C^{-1}.
For GeMSLR, H -> (I-H)^{-1} - I
For EsMSLR, X -> X^{-1} - I.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]VMatrix holding the orthogonal basis genrated by Arnoldi iteration.
[in]HThe Hessenberg matrix generated by Arnoldi iteration.
[in]mThe number of finished steps in Arnoldi iteration.
[in]rankThe number of eigenvectors to keep.
[in]levelThe current level.
[in]optionThe setup option. GeMSLR, ESMSLR, or PSLR.
Returns
return # of Schur vectors, if return value < 0 an error occurs

◆ SetupLowRankNoRestart()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupLowRankNoRestart ( VectorType &  x,
VectorType &  rhs,
int  level,
int  option 
)

Setup the low-rank part of the GeMSLR with arnodi no-restart.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe target level.
[in]optionThe setup option. GeMSLR, ESMSLR, or PSLR.
Returns
Return error message.

◆ SetupLowRankSubspaceIteration()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupLowRankSubspaceIteration ( VectorType &  x,
VectorType &  rhs,
int  level,
int  option 
)

Setup the low-rank part of the GeMSLR with subspace iteration

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe target level.
[in]optionThe setup option. GeMSLR, ESMSLR, or PSLR.
Returns
Return error message.

◆ SetupLowRankThickRestart()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupLowRankThickRestart ( VectorType &  x,
VectorType &  rhs,
int  level,
int  option 
)

Setup the low-rank part of the GeMSLR with arnodi thick-restart.

Parameters
[in]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe target level.
[in]optionThe setup option. GeMSLR, ESMSLR, or PSLR.
Returns
Return error message.

◆ SetupPartialILUT()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupPartialILUT ( VectorType &  x,
VectorType &  rhs 
)

Setup the local partial ILUT.

◆ SetupPermutation()

template<class MatrixType , class VectorType , typename DataType >
template int pargemslr::precond_gemslr_csr_par_complexd::SetupPermutation ( )

Setup the permutation of the GeMSLR.

Returns
Return error message.

◆ SetupPermutationBuildLevelStructure()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupPermutationBuildLevelStructure ( MatrixType &  A,
int  level_start,
vector_int map_v,
vector_int mapptr_v 
)

Setup the level structure of the GeMSLR.

Returns
Return error message.

◆ SetupPermutationND()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupPermutationND ( MatrixType &  A,
int  nlev_setup,
int &  nlev_max,
int &  nlev_used,
vector_int map_v,
vector_int mapptr_v 
)

Setup the Nested Dissection partition of the GeMSLR. This is a symmetric reordering algorithm work on A+A'.

Returns
Return error message.

◆ SetupPermutationRKway()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetupPermutationRKway ( MatrixType &  A,
int  nlev_setup,
int &  nlev_max,
int &  nlev_used,
vector_int map_v,
vector_int mapptr_v 
)

Setup the Recursive KWay partition of the GeMSLR. This is a symmetric reordering algorithm work on A+A'.
1) On each level, partition the current adjacency graph with kway partitioning.
2) Find the adjacency graph of the interface matrix.
3) Repeat until we reach the last level or the interface matrix is too small.

Returns
Return error message.

◆ SetWithParameterArray()

template<class MatrixType , class VectorType , typename DataType >
virtual int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SetWithParameterArray ( double *  params)
inlinevirtual

Setup with parameter array.

Parameters
[in]paramsThe parameter array.
Returns
Return error message.

◆ Solve()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::Solve ( VectorType &  x,
VectorType &  rhs 
)
virtual

Solve phase. Call this function after Setup. Solve with cusparse if unified memory/device memory is used.

Parameters
[in,out]xThe initial guess.
[in]rhsThe right-hand-side.
Returns
Return error message.

◆ SolveApplyLowRankLevel()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveApplyLowRankLevel ( VectorType &  x,
VectorType &  rhs,
int  level 
)

Apply the low-rank update on a certain level.

Parameters
[in,out]xThe initial guess.
[in]rhsThe right-hand-side.
[in]levelThe start level.
Returns
Return error message.

◆ SolveB()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveB ( VectorType &  x,
VectorType &  rhs,
int  option,
int  level 
)

Solve with B on a certain level.

Parameters
[in,out]xThe initial guess.
[in]rhsThe right-hand-side.
[in]option0: solve with B_solver; 1: solve with B_precond; 2: apply smoothing.
[in]levelThe start level.
Returns
Return error message.

◆ SolveLevelEsmslr()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelEsmslr ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

◆ SolveLevelEsmslrMul()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelEsmslrMul ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

◆ SolveLevelEsmslrU()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelEsmslrU ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

◆ SolveLevelGemslr()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelGemslr ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

◆ SolveLevelGemslrMul()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelGemslrMul ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

◆ SolveLevelGemslrU()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelGemslrU ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

◆ SolveLevelPslr()

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::SolveLevelPslr ( VectorType &  x_out,
VectorType &  rhs_in,
int  level,
bool  doperm 
)

Solve starting from a certain level.

Parameters
[in,out]x_inThe initial guess.
[in]rhs_outThe right-hand-side.
[in]levelThe start level.
[in]dopermDo we apply permutation on this level? (Some time the system is already permuted).
Returns
Return error message.

Member Data Documentation

◆ _global_precond_option

template<class MatrixType , class VectorType , typename DataType >
int pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::_global_precond_option

The Schur complement option.

◆ _lev_A

template<class MatrixType , class VectorType , typename DataType >
ParallelGemslrLevelClass< MatrixType, VectorType, DataType> pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::_lev_A

The structure for the A level low-rank correction.

◆ _levs_v

template<class MatrixType , class VectorType , typename DataType >
std::vector< ParallelGemslrLevelClass< MatrixType, VectorType, DataType> > pargemslr::ParallelGemslrClass< MatrixType, VectorType, DataType >::_levs_v

Vector holding the level struct for all levels starting from the second level. The size of this vector is equal to _nlev_used, and is the number of true levels.


The documentation for this class was generated from the following files: