Tridiagonalize \(A\) and approximate \(f(A)b \approx |b|_W V f(T) e_1\) via Lanczos algorithm. A is self-adjoint in the weights \( W\).
More...
|
| UniversalLanczos () |
| Allocate nothing, Call construct method before usage. More...
|
|
| UniversalLanczos (const ContainerType ©able, unsigned max_iterations) |
| Allocate memory for the method. More...
|
|
template<class ... Params> |
void | construct (Params &&...ps) |
| Perfect forward parameters to one of the constructors. More...
|
|
void | set_max (unsigned new_max) |
| Set the maximum number of iterations. More...
|
|
unsigned | get_max () const |
| Get the current maximum number of iterations. More...
|
|
void | set_verbose (bool verbose) |
| Set or unset debugging output during iterations. More...
|
|
value_type | get_bnorm () const |
| Norm of b from last call to operator() More...
|
|
template<class MatrixType , class ContainerType0 , class ContainerType1 , class ContainerType2 , class FuncTe1 > |
unsigned | solve (ContainerType0 &x, FuncTe1 f, MatrixType &&A, const ContainerType1 &b, const ContainerType2 &weights, value_type eps, value_type nrmb_correction=1., std::string error_norm="universal", value_type res_fac=1., unsigned q=1) |
| \( x = f(A)b \approx ||b||_W V f(T) e_1 \) via Lanczos and matrix function computation. A is self-adjoint in the weights \( W\). More...
|
|
template<class MatrixType , class ContainerType0 , class ContainerType1 > |
const HDiaMatrix & | tridiag (MatrixType &&A, const ContainerType0 &b, const ContainerType1 &weights, value_type eps=1e-12, value_type nrmb_correction=1., std::string error_norm="universal", value_type res_fac=1., unsigned q=1) |
| Tridiagonalization of A using Lanczos method. More...
|
|
unsigned | get_iter () const |
| Get the number of iterations in the last call to tridiag or solve (same as T.num_rows) More...
|
|
template<class ContainerType>
class dg::mat::UniversalLanczos< ContainerType >
Tridiagonalize \(A\) and approximate \(f(A)b \approx |b|_W V f(T) e_1\) via Lanczos algorithm. A is self-adjoint in the weights \( W\).
The M-Lanczos method is based on the paper Novel Numerical Methods for Solving the Time-Space Fractional Diffusion Equation in Two Dimensions by Q. Yang et al, but adopts a more efficient implementation similar to that in the PCG method. Further also the conventional Lanczos method can be found there and also in text books such as Iteratvie Methods for Sparse Linear Systems 2nd edition by Yousef Saad.
- Note
- We have two stopping criteria. The residual criterion stops when \( \tau ||r_i||_W = \tau ||\vec b||_W \beta_i (T^{-1})_{1m} \leq \epsilon_{rel} ||\vec b||_W + \epsilon_{aps} \) where \( \tau \) is a residual factor accounting for the condition of various matrix functions
-
The universal stopping criterion is based on the paper Estimating the error in matrix function approximations by Q. Eshghi, N. and Reichel L., The iteration stops when
\[ ||\vec{e}_{f,m}||_W = ||\vec{b}||_W ||\left(\check{T} - f(\bar{T})\right)\vec{e}_1||_2 \leq \epsilon_{rel} ||\vec{b}||_W ||f(T)\vec e_1||_2 + \epsilon_{abs} \]
with
\[ \bar{T} = \begin{pmatrix} T & \beta_m \vec{e}_m & & \\ \beta_m \vec{e}_m^T & \alpha_{m-1} & \beta_{m-2} & \\ & \beta_{m-2} & \alpha_{m-2} & \beta_{m-1-q} \\ & & \beta_{n-1-q} & \alpha_{n-q} \end{pmatrix} \]
\[ \check{T} = \begin{pmatrix} f(T) & 0 \\ 0 & 0 \end{pmatrix} \]
The common lanczos method (and M-Lanczos) method are prone to loss of orthogonality for finite precision. Here, only the basic Paige fix is used. Thus the iterations should be kept as small as possible. Could be fixed via full, partial or selective reorthogonalization strategies, but so far no problems occured due to this.
- Template Parameters
-
ContainerType | Any class for which a specialization of TensorTraits exists and which fulfills the requirements of the there defined data and execution policies derived from AnyVectorTag and AnyPolicyTag . Among others
dg::HVec (serial), dg::DVec (cuda / omp), dg::MHVec (mpi + serial) or dg::MDVec (mpi + cuda / omp)
std::vector<dg::DVec> (vector of shared device vectors), std::array<double, 4> (array of 4 doubles) or std::map < std::string, dg::DVec> ( a map of named vectors)
double (scalar) and other primitive types ...
If there are several ContainerTypes in the argument list, then TensorTraits must exist for all of them |
- See also
- See The dg dispatch system for a detailed explanation of our type dispatch system
template<class ContainerType >
template<class MatrixType , class ContainerType0 , class ContainerType1 , class ContainerType2 , class FuncTe1 >
unsigned dg::mat::UniversalLanczos< ContainerType >::solve |
( |
ContainerType0 & |
x, |
|
|
FuncTe1 |
f, |
|
|
MatrixType && |
A, |
|
|
const ContainerType1 & |
b, |
|
|
const ContainerType2 & |
weights, |
|
|
value_type |
eps, |
|
|
value_type |
nrmb_correction = 1. , |
|
|
std::string |
error_norm = "universal" , |
|
|
value_type |
res_fac = 1. , |
|
|
unsigned |
q = 1 |
|
) |
| |
|
inline |
\( x = f(A)b \approx ||b||_W V f(T) e_1 \) via Lanczos and matrix function computation. A is self-adjoint in the weights \( W\).
Tridiagonalize \(A\) using Lanczos algorithm with a residual or universal stopping criterion
- Parameters
-
x | output vector |
f | the matrix function that is called like yH = f( T) where T is the tridiagonal matrix and returns the result of \( f(T)\vec e_1\) (for example dg::mat::make_FuncEigen_Te1( dg::SQRT<double>()) ) |
- See also
- Matrix-functions
- Parameters
-
A | self-adjoint and semi-positive definit matrix |
b | input vector |
weights | weights in which A is self-adjoint |
eps | relative accuracy of M-Lanczos method |
nrmb_correction | the absolute error in units of eps to be respected |
error_norm | Either "residual" or "universal" |
res_fac | factor \( \tau\) that is multiplied to the norm of the residual. Used to account for specific matrix function and operator in the convergence criterium |
q | The q-number in the "universal stopping criterion |
- Returns
- number of iterations of M-Lanczos routine
template<class ContainerType >
template<class MatrixType , class ContainerType0 , class ContainerType1 >
Tridiagonalization of A using Lanczos method.
Tridiagonalize \(A\) using Lanczos algorithm with a residual or universal stopping criterion
- Parameters
-
A | A self-adjoint, positive definit matrix |
b | The initial vector that starts orthogonalization |
weights | Weights that define the scalar product in which A is self-adjoint and in which the error norm is computed. |
eps | relative accuracy of residual |
nrmb_correction | the absolute error C in units of eps to be respected |
error_norm | Either "residual" or "universal" |
res_fac | factor \( \tau\) that is multiplied to the norm of the residual. Used to account for specific matrix function and operator in the convergence criterium |
q | The q-number in the "universal stopping criterion |
- Returns
- returns the tridiagonal matrix T. Note that \( T = (MV)^T A V \). The number of iterations is given by
T.num_rows