136 CHAPTER 7. SIMULTANEOUS VECTOR OR SUBSPACE ITERATIONS

has the order n = 1095. Subspace iteration is applied with five vectors as an inverse iteration to

L−1AL−T (Lx) = λ(Lx), B = LLT . (Cholesky factorization)

X(0) is chosen to be a random matrix. Here, we number the eigenvalues from small to big. The smallest six eigenvalues of (7.10) shifted by 0.01 to the right are

and thus the ratios of the eigenvalues that determine the rate of convergence are

So, the numbers presented in Table 7.1 reflect quite accurately the predicted rates. The numbers in column 6 are a little too small, though.

The convergence criterion

max k(I − X(k)X(k) )x(ik−1)k ≤ ǫ = 10−5

1≤i≤p

was not satisfied after 50 iteration step.

7.3Accelerating subspace iteration

Subspace iteration potentially converges very slowly. It can be slow even it one starts with a subspace that contains all desired solutions! If, e.g., x(0)1 and x(0)2 are both elements in R([u1, u2]), the the vectors x(ik), i = 1, 2, . . ., still converge linearly towards u1 und u2 although they could be readily obtained from the 2 × 2 eigenvalue problem,

The following theorem gives hope that the convergence rates can be improved if one proceeds in a suitable way.

Theorem 7.6 Let X Fn×p as in Theorem 7.2. Let ui, 1 ≤ i ≤ p, be the eigenvectors corresponding to the eigenvalues λ1, . . . , λp of A. Then we have

138 CHAPTER 7. SIMULTANEOUS VECTOR OR SUBSPACE ITERATIONS

ˆ

The solution of the full eigenvalue problems H(k)y = ϑy is solved by the symmetric QR algorithm.

ˆ

The computation of the matrix H(k) is expensive as matrix-vector products have to be executed. The following considerations simplify matters. We write X(k) in the form

So, the columns of Y (k) are the normalized eigenvectors of H(k) := Z(k) Z(k).

Thus we obtain a second variant of the inverse subspace iteration with Rayleigh-Ritz step.

Algorithm 7.3 Subspace or simultaneous inverse vector iteration combined with Rayleigh-Ritz step, version 2

1:Let X Fn×p with X X = Ip.

2:Set X(0) := X.

3:for k = 1, 2, . . . do

4:AZ(k) := X(k−1);

Remark 7.2. An alternative to Algorithm 7.3 is the subroutine ritzit, that has been programmed by Rutishauser [4] in ALGOL, see also [2, p.293].

We are now going to show that the Ritz vectors converge to the eigenvectors, as Theorem 7.6 lets us hope. First we prove

Lemma 7.7 ([2, p.222]) Let y be a unit vector and ϑ F. Let λ be the eigenvalue of A closest to ϑ and let u be the corresponding eigenvector. Let

where r(y, ϑ) = Ay − ϑy plays the role of a residual.

140 CHAPTER 7. SIMULTANEOUS VECTOR OR SUBSPACE ITERATIONS

Therefore, in (7.13) we can interpret this expression as an perturbation of the diagonal matrix Λ−2 2 − µ−2I. For su ciently large k (that may depend on i) there is a µi that is close to λi, and a ti that is close to ei. We now assume that k is so big that

|µ−i 2 − λ−i 1| ≤ ρ := 12 min |λ−i 2 − λ−j 2|

λj =6 λi