SenK
SenK is a C++ library for high-performance linear solvers.
senk::solver Namespace Reference

Contains solvers. More...

Functions

template<typename T >
void Bicgstab (T *val, int *cind, int *rptr, T *b, T *x, T nrm_b, int max_iter, int N, T epsilon)
 Non-preconditioned BiCGStab solver. More...
 
template<typename T >
void IluBicgstab (T *val, int *cind, int *rptr, T *lval, int *lcind, int *lrptr, T *uval, int *ucind, int *urptr, T *b, T *x, T nrm_b, int max_iter, int N, T epsilon)
 ILU preconditioned BiCGStab solver. More...
 
template<typename T , int bnl, int bnw>
void IlubBicgstab (T *val, int *cind, int *rptr, T *blval, int *blcind, int *blrptr, T *buval, int *bucind, int *burptr, T *b, T *x, T nrm_b, int max_iter, int N, T epsilon)
 ILUB preconditioned BiCGStab solver. More...
 
template<typename T >
void Gcrm (T *val, int *cind, int *rptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The Non-preconditioned GCR(m) solver. More...
 
template<typename T >
void IluGcrm (T *val, int *cind, int *rptr, T *lval, int *lcind, int *lrptr, T *uval, int *ucind, int *urptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILU preconditioned GCR(m) solver. More...
 
template<typename T , int bnl, int bnw>
void IlubGcrm (T *val, int *cind, int *rptr, T *blval, int *blcind, int *blrptr, T *buval, int *bucind, int *burptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILUB preconditioned GCR(m) solver. More...
 
template<typename T >
void Gmresm (T *val, int *cind, int *rptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The Non-preconditioned GMRES(m) solver. More...
 
template<typename T , int bnl, int bnw>
void Gmresm (T *bval, int *bcind, int *brptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The Non-preconditioned GMRES(m) solver. More...
 
template<typename T >
void IluGmresm (T *val, int *cind, int *rptr, T *lval, int *lcind, int *lrptr, T *uval, int *ucind, int *urptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILU preconditioned GMRES(m) solver. More...
 
template<typename T >
void AmcIluGmresm (T *val, int *cind, int *rptr, T *lval, int *lcind, int *lrptr, T *uval, int *ucind, int *urptr, int *cptr, int cnum, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILU preconditioned GMRES(m) solver parallelized by AMC ordering. More...
 
template<typename T >
void AbmcIluGmresm (T *val, int *cind, int *rptr, T *lval, int *lcind, int *lrptr, T *uval, int *ucind, int *urptr, int *cptr, int cnum, int bsize, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILU preconditioned GMRES(m) solver parallelized by ABMC ordering. More...
 
template<typename T >
void BjIluGmresm (T *val, int *cind, int *rptr, T *lval, int *lcind, int *lrptr, T *uval, int *ucind, int *urptr, int bnum, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILU preconditioned GMRES(m) solver parallelized by the block Jacobi method. More...
 
template<typename T , int bnl, int bnw>
void IlubGmresm (T *val, int *cind, int *rptr, T *blval, int *blcind, int *blrptr, T *buval, int *bucind, int *burptr, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILUB preconditioned GMRES(m) solver. More...
 
template<typename T , int bnl, int bnw>
void AbmcIlubGmresm (T *val, int *cind, int *rptr, T *blval, int *blcind, int *blrptr, T *buval, int *bucind, int *burptr, int *cptr, int cnum, int bsize, T *b, T *x, T nrm_b, int outer, int m, int N, T epsilon)
 The ILUB preconditioned GMRES(m) solver parallelized by ABMC ordering. More...
 

Detailed Description

Contains solvers.

Function Documentation

◆ AbmcIlubGmresm()

template<typename T , int bnl, int bnw>
void senk::solver::AbmcIlubGmresm ( T *  val,
int *  cind,
int *  rptr,
T *  blval,
int *  blcind,
int *  blrptr,
T *  buval,
int *  bucind,
int *  burptr,
int *  cptr,
int  cnum,
int  bsize,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILUB preconditioned GMRES(m) solver parallelized by ABMC ordering.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
blvalvalues of L in the BCSR format.
blcindcolum positions of blocks of L in the BCSR format.
blrptrstarting positions of row blocks of L in the BCSR format.
buvalvalues of U in the BCSR format.
bucindcolum positions of blocks of U in the BCSR format.
burptrstarting positions of row blocks of U in the BCSR format.
cptrThe starting index of each color is stored.
cnumThe number of colors.
bsizeThe number of rows/columns of the blocks used in ABMC.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 680 of file senk_gmres.hpp.

◆ AbmcIluGmresm()

template<typename T >
void senk::solver::AbmcIluGmresm ( T *  val,
int *  cind,
int *  rptr,
T *  lval,
int *  lcind,
int *  lrptr,
T *  uval,
int *  ucind,
int *  urptr,
int *  cptr,
int  cnum,
int  bsize,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILU preconditioned GMRES(m) solver parallelized by ABMC ordering.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
lvalSame as val, but for the matrix L.
lcindSame as cind, but for the matrix L.
lrptrSame as rptr, but for the matrix L.
uvalSame as val, but for the matrix U.
ucindSame as cind, but for the matrix U.
urptrSame as rptr, but for the matrix U.
cptrThe starting index of each color is stored.
cnumThe number of colors.
bsizeThe number of rows/columns of the blocks used in ABMC.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 392 of file senk_gmres.hpp.

◆ AmcIluGmresm()

template<typename T >
void senk::solver::AmcIluGmresm ( T *  val,
int *  cind,
int *  rptr,
T *  lval,
int *  lcind,
int *  lrptr,
T *  uval,
int *  ucind,
int *  urptr,
int *  cptr,
int  cnum,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILU preconditioned GMRES(m) solver parallelized by AMC ordering.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
lvalSame as val, but for the matrix L.
lcindSame as cind, but for the matrix L.
lrptrSame as rptr, but for the matrix L.
uvalSame as val, but for the matrix U.
ucindSame as cind, but for the matrix U.
urptrSame as rptr, but for the matrix U.
cptrThe starting index of each color is stored.
cnumThe number of colors.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 294 of file senk_gmres.hpp.

◆ Bicgstab()

template<typename T >
void senk::solver::Bicgstab ( T *  val,
int *  cind,
int *  rptr,
T *  b,
T *  x,
nrm_b,
int  max_iter,
int  N,
epsilon 
)

Non-preconditioned BiCGStab solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcol-index array of the CSR storage format.
rptrrow-ptr array of the CSR storage format.
bThe right-hand side vector.
xThe unknown vector.
nrm_bThe 2-norm of b.
max_iterThe maximum number of iterations.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 30 of file senk_bicgstab.hpp.

◆ BjIluGmresm()

template<typename T >
void senk::solver::BjIluGmresm ( T *  val,
int *  cind,
int *  rptr,
T *  lval,
int *  lcind,
int *  lrptr,
T *  uval,
int *  ucind,
int *  urptr,
int  bnum,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILU preconditioned GMRES(m) solver parallelized by the block Jacobi method.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
lvalSame as val, but for the matrix L.
lcindSame as cind, but for the matrix L.
lrptrSame as rptr, but for the matrix L.
uvalSame as val, but for the matrix U.
ucindSame as cind, but for the matrix U.
urptrSame as rptr, but for the matrix U.
bnumThe number of blocks.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 488 of file senk_gmres.hpp.

◆ Gcrm()

template<typename T >
void senk::solver::Gcrm ( T *  val,
int *  cind,
int *  rptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The Non-preconditioned GCR(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 31 of file senk_gcr.hpp.

◆ Gmresm() [1/2]

template<typename T , int bnl, int bnw>
void senk::solver::Gmresm ( T *  bval,
int *  bcind,
int *  brptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The Non-preconditioned GMRES(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
bnlThe number of rows of the block.
bnwThe number of columns of the block.
Parameters
bvalval array of the BCSR storage format.
bcindcolumn index array of the BCSR storage format.
brptrrow pointer array of the BCSR storage format.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 114 of file senk_gmres.hpp.

◆ Gmresm() [2/2]

template<typename T >
void senk::solver::Gmresm ( T *  val,
int *  cind,
int *  rptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The Non-preconditioned GMRES(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 34 of file senk_gmres.hpp.

◆ IlubBicgstab()

template<typename T , int bnl, int bnw>
void senk::solver::IlubBicgstab ( T *  val,
int *  cind,
int *  rptr,
T *  blval,
int *  blcind,
int *  blrptr,
T *  buval,
int *  bucind,
int *  burptr,
T *  b,
T *  x,
nrm_b,
int  max_iter,
int  N,
epsilon 
)

ILUB preconditioned BiCGStab solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
bnlThe number of rows of the block.
bnwThe number of columns of the block.
Parameters
valval array of the CSR storage format.
cindcol-index array of the CSR storage format.
rptrrow-ptr array of the CSR storage format.
blvalvalues of L in the BCSR format.
blcindcolum positions of blocks of L in the BCSR format.
blrptrstarting positions of row blocks of L in the BCSR format.
buvalvalues of U in the BCSR format.
bucindcolum positions of blocks of U in the BCSR format.
burptrstarting positions of row blocks of U in the BCSR format.
bThe right-hand side vector.
xThe unknown vector.
nrm_bThe 2-norm of b.
max_iterThe maximum number of iterations.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 199 of file senk_bicgstab.hpp.

◆ IlubGcrm()

template<typename T , int bnl, int bnw>
void senk::solver::IlubGcrm ( T *  val,
int *  cind,
int *  rptr,
T *  blval,
int *  blcind,
int *  blrptr,
T *  buval,
int *  bucind,
int *  burptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILUB preconditioned GCR(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
blvalvalues of L in the BCSR format.
blcindcolum positions of blocks of L in the BCSR format.
blrptrstarting positions of row blocks of L in the BCSR format.
buvalvalues of U in the BCSR format.
bucindcolum positions of blocks of U in the BCSR format.
burptrstarting positions of row blocks of U in the BCSR format.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 164 of file senk_gcr.hpp.

◆ IlubGmresm()

template<typename T , int bnl, int bnw>
void senk::solver::IlubGmresm ( T *  val,
int *  cind,
int *  rptr,
T *  blval,
int *  blcind,
int *  blrptr,
T *  buval,
int *  bucind,
int *  burptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILUB preconditioned GMRES(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
blvalvalues of L in the BCSR format.
blcindcolum positions of blocks of L in the BCSR format.
blrptrstarting positions of row blocks of L in the BCSR format.
buvalvalues of U in the BCSR format.
bucindcolum positions of blocks of U in the BCSR format.
burptrstarting positions of row blocks of U in the BCSR format.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 583 of file senk_gmres.hpp.

◆ IluBicgstab()

template<typename T >
void senk::solver::IluBicgstab ( T *  val,
int *  cind,
int *  rptr,
T *  lval,
int *  lcind,
int *  lrptr,
T *  uval,
int *  ucind,
int *  urptr,
T *  b,
T *  x,
nrm_b,
int  max_iter,
int  N,
epsilon 
)

ILU preconditioned BiCGStab solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcol-index array of the CSR storage format.
rptrrow-ptr array of the CSR storage format.
lvalSame as val, but for the matrix L.
lcindSame as cind, but for the matrix L.
lrptrSame as rptr, but for the matrix L.
uvalSame as val, but for the matrix U.
ucindSame as cind, but for the matrix U.
urptrSame as rptr, but for the matrix U.
bThe right-hand side vector.
xThe unknown vector.
nrm_bThe 2-norm of b.
max_iterThe maximum number of iterations.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 108 of file senk_bicgstab.hpp.

◆ IluGcrm()

template<typename T >
void senk::solver::IluGcrm ( T *  val,
int *  cind,
int *  rptr,
T *  lval,
int *  lcind,
int *  lrptr,
T *  uval,
int *  ucind,
int *  urptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILU preconditioned GCR(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
lvalSame as val, but for the matrix L.
lcindSame as cind, but for the matrix L.
lrptrSame as rptr, but for the matrix L.
uvalSame as val, but for the matrix U.
ucindSame as cind, but for the matrix U.
urptrSame as rptr, but for the matrix U.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 94 of file senk_gcr.hpp.

◆ IluGmresm()

template<typename T >
void senk::solver::IluGmresm ( T *  val,
int *  cind,
int *  rptr,
T *  lval,
int *  lcind,
int *  lrptr,
T *  uval,
int *  ucind,
int *  urptr,
T *  b,
T *  x,
nrm_b,
int  outer,
int  m,
int  N,
epsilon 
)

The ILU preconditioned GMRES(m) solver.

Template Parameters
TThe type of a coefficient matrix and vectors.
Parameters
valval array of the CSR storage format.
cindcolumn index array of the CSR storage format.
rptrrow pointer array of the CSR storage format.
lvalSame as val, but for the matrix L.
lcindSame as cind, but for the matrix L.
lrptrSame as rptr, but for the matrix L.
uvalSame as val, but for the matrix U.
ucindSame as cind, but for the matrix U.
urptrSame as rptr, but for the matrix U.
bA right-hand side vector.
xAn unknown vector.
nrm_bThe 2-norm of b.
outerThe maximum number of iterations of outer loop.
mThe number of the restart period.
NThe size of the matrix and the vectors.
epsilonThe convergence criterion.

Definition at line 198 of file senk_gmres.hpp.