lanczos Interface

public interface lanczos

Description

Given a symmetric or Hermitian linear operator $A$, find matrices $X$ and $T$ such that

$$AX_k = X_k T_k + t_{k+1, k} x_{k+1} e_k^T,$$

where $X$ is an orthogonal basis and $T$ is symmetric tridiagonal.

Algorithmic Features

• The operator $A$ only needs to be accessed through matrix-vector products.
• Constructs an orthonormal Krylov basis $X$ via the Lanczos process with full reorthogonalization.
• Constructs a symmetric tridiagonal matrix $T$ whose eigenvalues approximates those of $A$.
• Checks for convergence and invariant subspaces.

References

• Y. Saad. "Iterative methods for sparse linear systems", SIAM 2nd edition, 2003. see Chapter 6.6: The symmetric Lanczos algorithm.

Syntax

    call lanczos(A, X, T, info [, kstart] [, kend] [, tol])

Arguments

A : Symmetric or Hermitian linear operator derived from one the base types provided by the AbstractLinops module. It is an intent(in) argument.

X : Array of types derived from one the base types provided by the AbstractVectors module. It needs to be consistent with the type of A. On exit, it contains the the computed Krylov vectors. The first entry X(1) is the starting vector for the Lanczos factorization. Additionally, the maximum number of Lanczos steps is equal to size(X) - 1. It is an intent(inout) argument.

T : real or complex rank-2 array. On exit, it contains the $(k+1) \times k$ symmetric tridiagonal matrix computed from the Arnoldi factorization. It is an intent(inout) argument.

info : integer variable. It is the LightKrylov information flag. On exit, if info > 0, the Lanczos factorization experienced a lucky breakdown. The array of Krylov vectors X spans an $A$-invariant subpsace of dimension info.

kstart (optional): integer value determining the index of the first Lanczos step to be computed. By default, kstart = 1.

kend (optional): integer value determining the index of the last Lanczos step to be computed. By default, kend = size(X) - 1.

tol (optional): Numerical tolerance below which a subspace is considered to be $A$-invariant. By default tol = atol_sp or tol = atol_rp depending on the kind of A.

Module Procedures