SparseSolverMixedPrecision Allows to use lower precision (float) for the factorization/preconditioner, and higher (double) for the outer iterative solver. See also LAPACK's dsgesv. More...
#include <StrumpackSparseSolverMixedPrecision.hpp>
Public Member Functions | |
SparseSolverMixedPrecision (bool verbose=true, bool root=true) | |
SparseSolverMixedPrecision (int argc, char *argv[], bool verbose=true, bool root=true) | |
~SparseSolverMixedPrecision () | |
void | set_matrix (const CSRMatrix< refine_t, integer_t > &A) |
void | set_matrix (const CSRMatrix< factor_t, integer_t > &A) |
ReturnCode | factor () |
ReturnCode | reorder (int nx=1, int ny=1, int nz=1) |
ReturnCode | solve (const DenseMatrix< refine_t > &b, DenseMatrix< refine_t > &x, bool use_initial_guess=false) |
ReturnCode | solve (const refine_t *b, refine_t *x, bool use_initial_guess=false) |
ReturnCode | solve (const DenseMatrix< factor_t > &b, DenseMatrix< factor_t > &x, bool use_initial_guess=false) |
ReturnCode | solve (const factor_t *b, factor_t *x, bool use_initial_guess=false) |
SPOptions< refine_t > & | options () |
const SPOptions< refine_t > & | options () const |
SparseSolver< factor_t, integer_t > & | solver () |
const SparseSolver< factor_t, integer_t > & | solver () const |
int | Krylov_iterations () const |
SparseSolverMixedPrecision Allows to use lower precision (float) for the factorization/preconditioner, and higher (double) for the outer iterative solver. See also LAPACK's dsgesv.
Mixed precision solver class. The input and output (sparse matrix, right hand side, and the solution of the linear system) are expected in refine_t precision, which should be double or std::complex<double>, while factor_t is the type for the internal factorization, which should be float or std::complex<float>.
There are options associated with the outer solver and with the inner solver (the lower-precision preconditioner). Make sure to set the correct ones. By default the inner solver options, such as solver tolerances, will be initialized for single precision, while the outer solver options are initialized for double precision. To change the final accuracy, access this->options().set_rel_tol(..) and this->options().set_abs_tol(..). To set preconditioner specific options, use this->solver().options().... By default, this will set the inner solver to be KrylovSolver::DIRECT (a single preconditioner application), and the outer solver to be KrylovSolver::AUTO (which will default to iterative refinement).
factor_t | can be: float or std::complex<float> |
refine_t | can be: double or std::complex<double> |
integer_t | defaults to a regular int. If regular int causes 32 bit integer overflows, you should switch to integer_t=int64_t instead. This should be a signed integer type. |
strumpack::SparseSolverMixedPrecision< factor_t, refine_t, integer_t >::SparseSolverMixedPrecision | ( | bool | verbose = true , |
bool | root = true |
||
) |
Constructor for the mixed precision solver class.
strumpack::SparseSolverMixedPrecision< factor_t, refine_t, integer_t >::SparseSolverMixedPrecision | ( | int | argc, |
char * | argv[], | ||
bool | verbose = true , |
||
bool | root = true |
||
) |
Constructor for the mixed precision solver class.
The command line arguments will be passed to both the options for the inner and outer solvers. Note that these options are not parsed until the user explicitly calls this->options().set_from_command_line(), and/or this->solver().options().set_from_command_line() (same as this->solver().set_from_options());
strumpack::SparseSolverMixedPrecision< factor_t, refine_t, integer_t >::~SparseSolverMixedPrecision | ( | ) |
Destructor.
|
inline |
Return the number of iterations performed by the outer (Krylov) iterative solver. Call this after calling the solve routine.