Basic Algorithm

Given an $n$ by $n$ matrix $A$ and initial $n$ by $p$ block vectors $P_1$ and $Q_1$ such that $P^{\ast}_1 Q_1 = I$, the block Lanczos method generates sequences of $n\times p$ right and left Lanczos block vectors $\{ Q_i \}$ and $\{ P_i \}$. These vectors are the biorthonormal bases of the Krylov subspaces ${\cal K}^j(A,Q_1)$ and ${\cal K}^j(A^{\ast},P_1)$. Specifically, after $j$ steps, the Lanczos procedure determines a block-tridiagonal matrix

$$T_{[j]} = \left[ \begin{array}{cccc} A_1 & B_2 & & \\ C_2... ...ts & \ddots & B_j \\ & & C_j & A_j \\ \end{array} \right],$$

and matrices of right and left Lanczos vectors

$$Q_{[j]} = \left[\, Q_1, Q_2, \ldots , Q_j \,\right] \quad \m... ...and} \quad P_{[j]} = \left[\, P_1, P_2, \ldots , P_j \,\right]$$

that satisfy the biorthonormality condition

$$P^{\ast}_{[j]} Q_{[j]} = I.$$ (155)

The block Lanczos method uses three-term recurrences

$$Q_{j+1} C_{j+1} = AQ_j - Q_j A_j - Q_{j-1} B_j,$$ (156)

$$P_{j+1} B^{\ast}_{j+1} = A^{\ast} P_j - P_j A^{\ast}_j - P_{j-1}C^{\ast}_j$$ (157)

for $j = 1, 2, \ldots$, where $Q_0 = P_0 = 0$. In matrix notation, these quantities satisfy the governing relations

$$A Q_{[j]} = Q_{[j]} T_{[j]} + Q_{j+1}C_{j+1}E^{\ast}_j,$$ (158)

$$A^{\ast} P_{[j]} = P_{[j]} T^{\ast}_{[j]} + P_{j+1} B^{\ast}_{j+1}E^{\ast}_j,$$ (159)

where $E_j$ is a $jp$ by $p$ matrix of which the bottom square is an identity matrix and zeros elsewhere.

To use the Lanczos method for approximating eigenvalues and eigenvectors of $A$, we solve the eigenvalue problem of the $jp \times jp$ block-tridiagonal matrix $T_{[j]}$. Each eigentriplet $( \theta^{(j)}_i, z^{(j)}_i, w^{(j)}_i )$ of $T_{[j]}$,

$$T_{[j]} z^{(j)}_i = \theta^{(j)}_i z^{(j)}_i \quad {\rm and... ...w^{(j)}_i)^{\ast} T_{[j]} = \theta^{(j)}_i (w^{(j)}_i)^{\ast}$$

determines a Ritz triplet $( \theta^{(j)}_i, x^{(j)}_i, y^{(j)}_i )$, where $x^{(j)}_i = Q_{[j]} z^{(j)}_i$ and $y^{(j)}_i = P_{[j]} w^{(j)}_i$. In general, Ritz triplets approximate outer eigentriplets of $A$ first.

To assess the approximation, let $r^{(j)}_i$ and $s^{(j)}_i$ denote the corresponding right and left residual vectors:

$$r^{(j)}_i = Ax^{(j)}_i - \theta^{(j)}_i x^{(j)}_i,$$ (160)

$$(s^{(j)}_i)^{\ast} = (y^{(j)}_i)^{\ast} A-\theta^{(j)}_i (y^{(j)}_i)^{\ast}.$$ (161)

Substitute (7.52) and (7.53), respectively, and there appears

$$r^{(j)}_i = Q_{j+1}C_{j+1} (E^{\ast}_j z^{(j)}_i),$$ (162)

$$(s^{(j)}_i)^{\ast} = ((w^{(j)}_i)^{\ast} E_j) B_{j+1} P^{\ast}_{j+1}.$$ (163)

A Ritz value $\theta^{(j)}_i$ is considered to have converged if both residual norms are sufficiently small. Similarly, the residuals can also be used to determine a backward error bound for the triplet, as in the case for the standard non-Hermitian Lanczos algorithm (see §7.8).

A simple blocked version of the basic non-Hermitian Lanczos method is presented in Algorithm 7.14.

We now comment on some steps of Algorithm 7.14.

(1)

The starting vectors $Q_1$ and $P_1$ are best selected by the user to embody any available knowledge concerning $A$'s wanted eigenvectors. In the absence of such knowledge, choose $Q_1$ to be a real orthonormalized random $n\times p$ matrix and let $P_1 = Q_1$.

The block size $p$ should be chosen to be larger than or equal to the multiplicity of any wanted eigenvalue. If the multiplicities are unknown, typical values of $p$ are $2$, $3$, and $6$.

(2), (3), (20), (21)

The user is required to furnish subroutines to calculate the products $AZ$ and $A^{\ast}Z$ for an arbitrary $n\times p$ matrix $Z$. This gives the user a chance to calculate these products with only one pass over the data structure defining $A$ and $A^{\ast}$, with possible improvements in efficiency.

(8), (9), (10)

First, the QR factorizations of $R^{\ast}$ and $S$ are computed. If $R^{\ast}$ and $S$ are both of full rank, then $Q_{j+1}$ and $P_{j+1}$ are $n$ by $p$ matrices such that $Q^{\ast}_{j+1} Q_{j+1} = P^{\ast}_{j+1} P_{j+1} = I$, and both $B^{\ast}_{j+1}$ and $C_{j+1}$ are $p$ by $p$ right upper triangular matrices. If either $R$ or $S$ is not of full rank, then an invariant subspace has been revealed. Continuing the recurrence in the rank-deficient case is discussed in §7.9.2.

(12)

If $\Sigma_j$ is singular or nearly singular, that is, if there is a zero or a very tiny singular value, then there is a breakdown. Proper action must be taken to continue. See §7.9.2 for possible treatments.

(17), (18)

The eigentriplets of $T_{[j]}$ can be computed using the method discussed in §7.3. The convergence of Ritz values can be tested by the norms of residuals (7.56) and (7.57). For more details on accessing the accuracy of the approximation, see § 7.8.2 and §7.13 and [29].

When $j$ becomes large, solving the eigenproblem of $T_{[j]}$ becomes time consuming. A simple way to reduce the cost is to do this periodically, say every five steps.

(19)

As in the unblocked Lanczos method, in the presence of finite precision arithmetic, the biorthogonality of the block Lanczos vectors $\{ Q_i \}$ and $\{ P_i \}$ deteriorates. The TSMGS process can be used to biorthogonalize $Q_{j+1}$ and $P_{j+1}$ in place against the previous Lanczos vectors $Q_{[j]} =\left[\, Q_1, Q_2, \ldots, Q_j\, \right]$ and $P_{[j]} =\left[\, P_1, P_2, \ldots, P_j\, \right]$.

		 PROCEDURE TSMGS 

for  $i=1,2,\ldots,j$ 

$Q_{j+1} := Q_{j+1} - Q_i ( P^{\ast}_i Q_{j+1})$ 

$P_{j+1} := P_{j+1} - P_i ( Q^{\ast}_i P_{j+1})$ 

end

(25)

This is recommended only if the biorthogonality of the Lanczos vectors is well maintained. The approximate eigenvectors $x^{(j)}_i$ and $y^{(j)}_i$ corresponding to converged Ritz values $\theta^{(j)}_i$ are computed and the residuals (7.54) and (7.55) are checked again relative to $\Vert x^{(j)}_i\Vert _2$ and $\Vert y^{(j)}_i\Vert _2$. Note that the norms can be greater than the estimates computed at step (17). This is a side-effect of the ill-conditioning of the Lanczos basis vectors.