On the Global and Cubic Convergence of a Quasi-Cyclic Jacobi Method

Noah H. Rhee Department of Mathematics, University of Missouri-Kansas City, Kansas City, MO 64110-2499 and Vjeran Hari^{∗} Department of Mathematics, University of Zagreb, Bijenićka cesta 30 41001 Zagreb, Croatia

Summary. In this paper we consider the global and the cubic convergence of a quasi-cyclic Jacobi method for the symmetric eigenvalue problem. The method belongs to a class of quasi-cyclic methods recently proposed by W. Mascarenhas. Mascarenhas showed that the methods from his class asymptotically converge cubically per quasi-sweep (one quasi-sweep is equivalent to 1.25 cyclic sweeps) provided the eigenvalues are simple . Here we prove the global convergence of our method and derive very sharp asymptotic convergence bounds in the general case of multiple eigenvalues. We discuss the ultimate cubic convergence of the method and present several numerical examples which well comply with the theory.

Subject Classifications: AMS (MOS): 65F15; CR: 5.14

∗This work was supported in part by the University of Minnesota Army High Performance Computing Research Center and the U.S. Army Contract DAAL03-89-C-0038. The paper was partly written while this author was a visiting faculty in the Department of Mathematics, University of Kansas, Lawrence, Kansas.The first version of this paper was made in July 1990 while this author was visiting AHPCRC.

Introduction

In the last decade much attention has been paid to diagonalization methods for solving eigenvalue problems. This interest is due to an inherent parallelism of such methods. Typical example is the Jacobi method for symmetric matrices. This method has been widely studied since sixties, but recently several new phenomena have been discovered. Using a subtle error analysis Demmel and Veselić [1] showed that this method can compute the eigenvalues and the eigenvectors of a full positive definite matrix more accurately than some other methods, in particular the QR method. Veselić and Hari [14] showed that even on standard computers implicit versions of this method can be not more than 2-3 times slower than the QR method. Finally, Mascarenhas [10] [11] introduced a class of quasi-cyclic pivot strategies and showed that in case of simple eigenvalues the appropriate symmetric Jacobi methods asymptotically converge cubically per quasi-sweep. Since one quasi-sweep corresponds to one and a quarter of (usual) sweeps, it means that the rate of asymptotic convergence is about 2.4 (cf. [11]),that is, the method converges better than quadratically in the usual sense.

In this paper we concentrate on one special method from the Mascarenhas class. For this method we show what could not have been proved for the general method from his class. First the global convergence. Then much sharper cubic convergence estimates. Finally, we extend our result to the case of multiple eigenvalues.

The paper is organized as follows. In Section 1 we introduce our quasi-cyclic method. In Section 2 we prove its global convergence. Finally, in Section 3 we prove its asymptotic cubic convergence per quasi-sweep. Numerical examples at the end of Section 3 shed light on the accuracy of our bounds and also on the generallity of the asymptotic assumptions.

1 Quasi-Cyclic Jacobi Methods

A Jacobi method for diagonalizing a symmetric n x n matrix A performs a sequence of similarity trasformations

A^{(k+1)} = U_{k}^{T} A^{(k)}U_{k}, k ≥ 0 (1.1)

where A^{(0)} = A and U_{0}, U_{1}, … are orthogonal plane rotations. For each k, U_{k} is defined by a pair of indices (p, q) = (p(k), q(k)) called pivot pair and by four essential elements u_{pp}^{(k)} = u_{qq}^{(k)} = cos ϕ_{k}, u_{pq}^{(k)} = −u_{qp}^{(k)} = sin ϕ_{k}. All other elements of U_{k} = (u_{ij}^{(k)}) are as in the identity matrix I_{n}. The process (1.1) is defined by a rule for computing the rotation angle ϕ_{k} and by a way of selecting pivot pairs (pivot strategy). We assume the standard angle choice

tan 2ϕ_{k} = \frac{2a_{pq}^{(k)}}{a_{qq}^{(k)} − a_{pp}^{(k)}} , ϕ_{k} ∈ [\frac{−π}{4} , \frac{π}{4} ] , (1.2)

which makes the pivot element a_{pq}^{(k)} zero i.e.

a_{p(k)q(k)}^{(k+1)} = 0, k ≥ 0 (1.3)

holds, where A^{(k)} = (a_{ij}^{(k)}). Consequently, we have

S^{2}(A^{(k+1)}) = S^{2}(A^{(k)}) − (a_{p(k)q(k)}^{(k)})^{2}, k ≥ 0 (1.4)

and

lim_{k→∞} a_{p(k)q(k)}^{(k)} = 0. (1.5)

Here

S(A) = \sqrt{∑_{1≤i<j≤n}(a_{ij})^{2}}, k ≥ 0.

If A is 1 × 1 matrix we set S(A) = 0.

Let N = n(n − 1)/2, N_{0} = {0, 1, 2, …} and P_{n} = {(i, j) : 1 ≤ i < j ≤ n}. A pivot strategy I is a function from N_{0} to P_{n} defined by I(k) = (p(k), q(k)), k ≥ 0. If I is a periodic function, then I is called a periodic strategy. Let I be a periodic strategy with period M. If M > N (M = N) and {I(k) : k = 0, 1, …, M − 1} = P_{n}, then I is called a quasi-cyclic (cyclic) strategy.

Let S be a subset of P_{n} and let ν(S) denote its cardinality. By O(S) we mean a collection of all finite sequences made of the elements from S. Suppose that O ∈ O(S). We assume that each element of S appears at least once in O; otherwise S can be replaced by a proper subset of it. Thus, each sequence from O(S) contains at least ν(S) terms. If S is an empty set, ν(S) = 0 and O(S) consists of a singleton, which is an empty sequence.

A cyclic or quasi-cyclic strategy can be specified in the following way. For any O = {(i_{r}, j_{r})}_{r=0}^{M−1} ∈ O(P_{n}) the cyclic or quasi-cyclic strategy I_{O} generated by O is given by

I_{O}(k) ≡ (p(k), q(k)) = (i_{r}, j_{r}), 0 ≤ r ≤ M − 1, k ≥ 0,

provided that k ≡ r (mod M).

We denote a Jacobi method, which defines ϕ_{k} by (1.2) and uses pivot strategy I by J(I). J(I) is called a quasi-cyclic/cyclic (Jacobi) method if I is a quasicyclic/cyclic strategy. Since n is fixed, in the sequel we write P for P_{n}.

The Method J_{Θ}

Here we briefly define a quasi-cyclic Jacobi method for symmetric matrices we deal with in this paper. Let S be any subset of P. By O_{R}(S) we denote the “row-wise ordering of S ”, that is the sequence satisfying the following two conditions: (i) each element (i, j) ∈ S appears exactly once in O_{R}(S); (ii) for any two terms (i_{1}, j_{1}) and (i_{2}, j_{2}) in O_{R}(S), (i_{1}, j_{1}) precedes (i_{2}, j_{2}) if i_{1} < i_{2} or i_{1} = i_{2} and j_{1} < j_{2}. In an obvious manner we can define the “column-wise ordering of S ” denoted by O_{C}(S). Let S_{i}, 1 ≤ i ≤ σ, be subsets of P and O_{i} ∈ O(S_{i}), 1 ≤ i ≤ σ. By [O_{1}, O_{2}, …, O_{σ}], we mean the sequence which is obtained by the concatination of the sequences O_{1},O_{2}, …, O_{σ}.

Let l and m be positive integers such that n = l + m. Let l_{1},l_{2},m_{1}, m_{2} be nonnegative integers such that l = l_{1} + l_{2} and m = m_{1} + m_{2}. We partition P into P_{1}, P_{2}, P_{3}, P_{4}, P_{5}, P_{6} and P_{7}, where

P_{1} = {(i, j) : 1 ≤ i < j ≤ l_{1}}
P_{2} = {(i, j) : l_{1} + 1 ≤ i < j ≤ l}
P_{3} = {(i, j) : 1 ≤ i ≤ l_{1}, l_{1} + 1 ≤ j ≤ l}
P_{4} = {(i, j) : l + 1 ≤ i < j ≤ l + m_{1}}
P_{5} = {(i, j) : l + m_{1} + 1 ≤ i < j ≤ n}
P_{6} = {(i, j) : l + 1 ≤ i ≤ l + m_{1}, l + m_{1} + 1 ≤ j ≤ n}
P_{7} = {(i, j) : 1 ≤ i ≤ l, l + 1 ≤ j ≤ n}

Let

(1.6) A = [ \begin{matrix}A_{1} & A_{3} \\ A_{3}^{T} & A_{2} \end{matrix} ] \begin{matrix}\} l \\ \} m \end{matrix} = [ \begin{matrix}\begin{matrix}B_{11} & B_{12} \\ B_{12}^{T} & B_{22} \end{matrix} & \frac{B_{1}}{B_{2}} \\ \frac{B_{3}}{B_{4}} & \begin{matrix}B_{33} & B_{34} \\ B_{34}^{T} & B_{44} \end{matrix} \end{matrix} ] \begin{matrix}\} & l_{1} \\ \} & l_{2} \\ \} & m_{1} \\ \} & m_{2} \end{matrix}

\overbrace{l} \overbrace{m}  \overbrace{l_{1}} \overbrace{l_{2}} \overbrace{m_{1}} \overbrace{m_{2}}

We see that the sets P_{1}, P_{2}, P_{3}, P_{4}, P_{5}, P_{6}, and P_{7} are associated with the submarices B_{11}, B_{22}, B_{12}, B_{33}, B_{44}, B_{34} and A_{3}, respectively. Note that some of P_{i}′s can be empty. For example if l_{1} = 0, l_{2} = 1, m_{1} = 1 and m_{2} = 1, then P_{1}, P_{2}, P_{3}, P_{4} and P_{5} are empty sets and P_{6} = {(2, 3)}, P_{7} = {(1, 2), (1, 3)}. In this case the matrices B_{11}, B_{12} and B_{1} do not exist. In the sequel we denote O_{R}(P_{i}) by O_{i}, 1 ≤ i ≤ 7.

Definition 1.1 Let

Θ = [O_{1}, O_{2}, O_{3}, O_{1}, O_{2}, O_{4}, O_{5}, O_{6}, O_{4}, O_{5}, O_{7}].

By I_{Θ} we denote the quasi-cyclic strategy generated by Θ, and by M its period. By J_{Θ} we denote the corresponding quasi-cyclic Jacobi method J(I_{Θ}).

Note that Θ ∈ O(P) and all the elements of P_{1}, P_{2}, P_{4} and P_{5} are repeated exactly twice in the sequence Θ. If some P_{i} is empty, the corresponding sequence O_{i} is empty, and there is no Jacobi transformation corresponding to O_{i}. For example if l_{1} = 0, l_{2} = 1, m_{1} = 1 and m_{2} = 2,

Θ = (3, 4), (2, 3), (2, 4), (3, 4), (1, 2), (1, 3), (1, 4).

The method J_{Θ} is a typical representative of the class of quasi-cyclic methods proposed by W. Mascarenhas [10] [11]. In this paper we use special properties of the row-wise annihilation ordering to prove the global convergence of J_{Θ} and to derive sharp asymptotic estimates which imply its cubic convergence per quasisweep.

2 Global Convergence

A given Jacobi method J(I) is convergent on A if the generated sequence {A^{(k)}}_{k=0}^{∞} is convergent to a diagonal matrix. If J(I) is convergent on every symmetric matrix A, then J(I) is (globally) convergent. The choice of angle interval in (1.2) ensures that J(I) is convergent iff lim_{k→∞} S(A^{(k)}) = 0 for every symmetric matrix A.

In this section we prove that J_{Θ} is a convergent method. To this end we need the notions of equivalent, shift equivalent and weakly equivalent strategies. So, we extend these notions from cyclic [4, 5, 13] to quasi-cyclic strategies. In the sequel strategy means quasi-cyclic strategy.

Let S be a subset of P. Let O = {(i_{r}, j_{r})}_{r=0}^{s} ∈ O(S). An admissible transposition on O is any transposition of two adjacent terms (i_{r}, j_{r}), (i_{r+1}, j_{r+1}) → (i_{r+1}, j_{r+1}), (i_{r}, j_{r}), provided that the sets {i_{r}, j_{r}} and {i_{r+1}, j_{r+1}} are disjoint. The sequences O, O′ ∈ O(S) are equivalent if one can be obtained from the other by a finite number of admissible transpositions. In this case we write O ∼ O′. For example if S = {(1, 2), (2, 3), (1, 4)} then O = (1, 2), (1, 2), (1, 4), (2, 3), (2, 3) ∈ O(S) and O′ = (1, 2), (1, 2), (2, 3), (2, 3), (1, 4) ∈ O(S) are equivalent.

Let I be a strategy with period M. By O_{I} we mean the sequence {I(k)}_{k=0}^{M−1}.

Definition 2.1 (cf. [13]) Let I and I′ be two strategies with the same period M.

(i) I and I′ are equivalent if O_{I} ∼ O_{I′}. In such a case we write I ∼ I′.

(ii) I and I′ are shift equivalent if I′(k) = I(k+c) ,k ≥ 0, where 0 ≤ c ≤ M −1.

(iii) I and I′ are weakly equivalent if there exist strategies I^{(j)}, 1 ≤ j ≤ r, for some r ≥ 0, such that in the sequence I, I^{(1)}, …,I^{(r)}, I′ each two adjacent terms are either equivalent or shift-equivalent.

Note that for −M + 1 ≤ c ≤ −1, I(k + c) = I(k + M + c), k ≥ −c, so we can confine ourselves to nonnegative c.

For a given A and a given strategy I, the sequence {U_{k}}_{k=0}^{∞} is well defined. Let T_{I(k)}(X) mean U_{k}^{T} XU_{k}. Then (1.1) can be written A^{(k+1)} = T_{I(k)}(A^{(k)}), implying A^{(k+1)} = T_{I(k)}T_{I(k−1)} ⋯ T_{I(0)}(A), k ≥ 0. Let

T_{k}^{I} = T_{I(k−1)}T_{I(k−2)} ⋯ T_{I(0)}, k ≥ 1; T_{0}^{I} = E,

where E is the identity transformation. Note that T_{I(k)} commutes with T_{I(k+1)} if the sets {p(k), q(k)} and {p(k + 1), q(k + 1)} are disjoint. Hence using the same argument as in [7, Lemma 1.1], we obtain the following implication (see also [4]):

If I ∼ I′, then T_{kM}^{I}(A) = T_{kM}^{I′}(A), k ≥ 0. (2.1)

That is, equivalent strategies give the same matrices at the end of each quasisweep consisting of M iterations of (1.1). Since {S(A^{(k)})}_{k=0}^{∞} makes a nonincreasing sequence, (2.1) implies that J(I) is convergent iff J(I′) is convergent, provided that I ∼ I′.

Lemma 2.2 Let I and I′ be two weakly equivalent strategies. Then J(I) is convergent iff J(I′) is convergent.

Proof. The proof is identical to the one in [13].

To show that J_{Θ} is convergent, let the subsets T_{1}, T_{2}, R be associated with the submatrices A_{1}, A_{2}, A_{3} from (1.6), respectively, that is,

(2.2)  T_{1} = P_{1} ∪ P_{2} ∪ P_{3} T_{2} = P_{4} ∪ P_{5} ∪ P_{6} , R = P_{7}.

Let

Θ̃ = [O_{2}, O_{R}(T_{1}), O_{1}, O_{7}, O_{5}, O_{R}(T_{2}), O_{4}] ∈ O(P),

and let IΘ̃ and JΘ̃ be the associated strategy and the corresponding Jacobi method.

Lemma 2.3 J_{Θ} is convergent iff JΘ̃ is convergent.

Proof. Let Θ′ , Θ′′ ∈ O(P) be defined by

Θ′ = [O_{5}, O_{4}, O_{6}, O_{5}, O_{4}, O_{2}, O_{1}, O_{3}, O_{2}, O_{1}, O_{7}],
Θ′′ = [O_{2}, O_{1}, O_{3}, O_{2}, O_{1}, O_{7}, O_{5}, O_{4}, O_{6}, O_{5}, O_{4}],

and let I_{Θ′} , I_{Θ′′} be the corresponding strategies. Obviously Θ ∼ Θ′ hence by Definition 2.1 (i) I_{Θ} ∼ I_{Θ′}. Since I_{Θ′′}(k) = I_{Θ′}(k + c) with c = 2 [ν(P_{4}) + ν(P_{5})] + ν(P_{6}), I_{Θ′} and I_{Θ′′} are shift equivalent. By [7, Lemma 1.2 (iii)] [O_{1}, O_{3}, O_{2}] ∼ O_{R}(T_{1}) and [O_{4}, O_{6}, O_{5}] ∼ O_{R}(T_{2}). Hence Θ′′ ∼ Θ̃ and consequently I_{Θ′′} ∼ IΘ̃. Thus I_{Θ} is weakly equvalent to IΘ̃ and the proof follows by Lemma 2.2.

In the definition of Θ̃ we find the row-cyclic orderings of the sets T_{1} and T_{2} which correspond to the diagonal blocks A_{1} and A_{2}. For such orderings we have the following general result.

Lemma 2.4 Let H be a symmetric matrix of order r ≥ 2. Let H′ be the matrix which is obtained from H by applying one sweep of the row-cyclic Jacobi method. Then

S(H′) ≤ C S(H)

holds, where

C^{2} = 1 − 2^{− \frac{(r−2)(r−1)}{2}} , C ≥ 0.

Proof. See [8].

For the row-cyclic ordering of R, O_{7}, we have the following useful result.

Lemma 2.5 Let T_{1}, T_{2}, R be the sets defined by the relation (2.2). Let A^{(0)} = A, A^{(1)}, …, A^{(ν(R))} be generated by ν(R) = lm Jacobi steps following the ordering O_{7}. Then there exist positive constants K_{1}, K_{2}, K_{3}, not larger than one, depending only on l and m, such that the following implication holds:

If there is an ε > 0 for which the conditions

(i) |a_{ij}| < K_{1}ε for all (i, j) ∈ T_{1} and

(ii) |a_{ξη}| ≥ ε for some (ξ, η) ∈ R ∪ T_{2},

hold, then at least one of the statements (iii), (iv) is true:

(iii) |a_{ρσ}^{(ν(R))}| ≥ K_{2}ε for some (ρ, σ) ∈ T_{2},

(iv) |a_{p(k)q(k)}^{(k)}| ≥ K_{3}ε for some 0 ≤ k < ν(R).

Proof. See [12, Lemma 4]. We remark that if T_{1} is empty, the condition (i) is dropped, and if T_{2} is empty, the statement (iii) cannot hold, and hence the statement (iv) must hold.

To prove the convergence of JΘ̃ we assume that matrices A^{(0)} = A, A^{(1)}, … (cf. (1.6)),

A^{(k)} = [ \begin{matrix}A_{1}^{(k)} & A_{3}^{(k)} \\ A_{3}^{(k)}^{T} & A_{2}^{(k)} \end{matrix} ] \begin{matrix}\} & l \\ \} & m \end{matrix} , k ≥ 0,
\overbrace{l}  \overbrace{m}

are obtained by applying the method JΘ̃ to A. For r ≥ 0 we define

α(r) = rM + ν(T_{1}) + ν(P_{1}) + ν(P_{2}),
β(r) = α(r) + ν(R).

Note that for any r ≥ 0, the transformations T_{I(α(r))}, …, T_{I(β(r)−1)} rotate the elements of A_{3}^{(α(r))}, …, A_{3}^{(β(r)−1)} according to the ordering O_{7}.

Lemma 2.6 If we apply the method JΘ̃ to A, then

(i) S(A_{1}^{(α(r))}) → 0 as r → ∞ (ii) S(A_{2}^{(β(r))}) → 0 as r → ∞.

Proof. Since {S(A^{(k)})}_{k=0}^{∞} is a nonincreasing sequence of real numbers bounded below by 0, it is convergent. Let μ = lim_{k→∞} S(A^{(k)}). If μ = 0 (i) and (ii) necessarily hold. So, let μ > 0. Then {γ(k)}_{k=0}^{∞},

γ(k) = \sqrt{S^{2}(A^{(k)}) − μ^{2}}, k ≥ 0,

converges monotonically to zero, i.e.

γ(k) → 0 as k → ∞. (2.3)

(i) Since for rM ≤ k ≤ α(r), A_{2}^{(k)} and ∥A_{3}^{(k)}∥ remain invariant under those transformations involved, we have

γ^{2}(rM) ≥ S^{2}(A^{(rM)}) − S^{2}(A^{(α(r))}) = S^{2}(A_{1}^{(rM)}) − S^{2}(A_{1}^{(α(r))}), r ≥ 0. (2.4)

Using the fact that {S^{2}(A_{1}^{(k)})}_{k=rM}^{α(r)} is a nonincreasing finite sequence and using Lemma 2.4 we obtain

(2.5) S^{2}(A_{1}^{(α(r))}) ≤ S^{2}(A_{1}^{(α(r)−ν(P_{1}))}) ≤ (C_{1})^{2} S^{2}(A_{1}^{(α(r)−ν(P_{1})−ν(T_{1}))})
≤ (C_{1})^{2} S^{2}(A_{1}^{(rM)}).

Here 0 ≤ C_{1} < 1 depends only on l. The relations (2.4) and (2.5) yield

S(A_{1}^{(rM)}) ≤ \frac{γ(rM)}{\sqrt{1 − (C_{1})^{2}}} , r ≥ 0.

Hence, the relation (2.3) implies that S(A_{1}^{(rM)}) → 0 as r → ∞ and the assertion (i) follows by the inequality (2.5).

(ii) The proof is quite analogous to the proof of (i).

Theorem 2.7 The method J_{Θ} is convergent.

Proof. By Lemma 2.3 it suffices to show that JΘ̃ is convergent. If the limit μ of S(A^{(k)}) is zero JΘ̃ is convergent. So, let μ > 0. From Lemma 2.6 we conclude that there exists a positive integer L_{1} such that

S^{2}(A_{1}^{(α(r))}) < \frac{μ^{2}}{2} , r ≥ L_{1}.

Since

S^{2}(A^{(k)}) = S^{2}(A_{1}^{(k)}) + S^{2}(A_{2}^{(k)}) + ∥A_{3}^{(k)}∥^{2}, k ≥ 0.

we have

S^{2}(A_{2}^{(α(r))}) + ∥A_{3}^{(α(r))}∥^{2} ≥ \frac{μ^{2}}{2} , r ≥ L_{1}.

Hence there is a pair (ξ, η) ∈ R ∪ T_{2} depending on r such that

|a_{ξη}^{(α(r))}| ≥ ε ≡ \frac{μ}{\sqrt{2(ν(R) + ν(T_{2}))}}, r ≥ L_{1}. (2.6)

Using the constants K_{1} and K_{2} from Lemma 2.5 we see that Lemma 2.6 implies the existence of a positive integer L_{2}, such that

(2.7) S(A_{1}^{(α(r))}) < K_{1} ε for all (i, j) ∈ T_{1} and r ≥ L_{2}

and

(2.8) S(A_{2}^{(β(r))}) < K_{2} ε for all (i, j) ∈ T_{2} and r ≥ L_{2}.

Using the constant K_{3} from Lemma 2.5 and the relation (1.5) we conclude that there is a positive integer L_{3} such that

| a_{p(k)q(k)}^{(k)} |< K_{3} ε, k ≥ α(r), r ≥ L_{3}. (2.9)

Now let L = max{L_{1}, L_{2}, L_{3}}. Obviously, the relations (2.6), (2.7),(2.8) and (2.9) hold with L in place of L_{1},L_{2}, and L_{3}. However, these relations contradict to Lemma 2.5 (apply it to A^{(α(L))} instead of A). Therefore, we conclude that μ > 0 cannot be true, which proves the theorem.

3 Asymptotic Cubic Convergence

In this section we derive sharp asymptotic convergence estimates for the method J_{Θ}. By cubic convergence one means the cubic reduction of S(A^{(k)}) per quasi-sweep, i.e. that

S(A^{((r+1)M)}) ≤ c S^{3}(A^{(rM)})

holds whenever r is sufficiently large. Here c is a constant which may depend on A. As we shall see in Subsection 3.3, a more appropriate measure for studying the cubic convergence of J_{Θ} is Γ = Γ(A) defined by (cf. [11])

Γ = max {Γ_{1}, Γ_{2}} , (3.0.1)

Γ_{1} = [\frac{\sqrt{S^{2}(A_{1}) + S^{2}(A_{2})}}{δ}]^{\frac{1}{2}} and Γ_{2} = [\frac{∥ A_{3} ∥}{δ}]^{\frac{1}{3}} . (3.0.2)

Here 3δ is the minimal distance between two distinct eigenvalues of A. By Γ_{1}^{(k)}, Γ_{2}^{(k)} and Γ^{(k)}, we denote Γ_{1}, Γ_{2} and Γ, defined with respect to A^{(k)}, k ≥ 0. We shall show that both Γ^{(k)} and S(A^{(k)}) converge cubically per quasi-sweep.

3.1 Asymptotic Assumptions and Related Results

In the sequel we use the following asymptotic assumptions:

A1 Γ < \frac{1}{2}.

A2 The diagonal elements of A affiliated with the same (multiple) eigenvalue occupy successive positions on the diagonal.

A3 The nonnegative numbers l_{1}, l_{2}, m_{1} and m_{2} are so chosen that the elements from different sets among {a_{11}, …, a_{l_{1}l_{1}}}, {a_{l_{1}+1,l_{1}+1}, …, a_{ll}}, {a_{l+1,l+1}, …, a_{l+m_{1},l+m_{1}}}, {a_{l+m_{1}+1,l+m_{1}+1}, …, a_{nn}} are affiliated with different eigenvalues.

The first two assumptions are quite similar to those used in the quadratic convergence consideration of Jacobi methods when multiple multiple eigenvalues are present (cf. [9, 7]). The third assumption ensures (cf. Example 3.11) the cubic convergence of Γ^{(k)}. In the case of simple eigenvalues it is always satisfied, so it can be removed. To ensure that the assumption A3 ultimately holds the numbers l_{1}, l_{2}, m_{1} and m_{2} can be in a later stage of the process redefined. Since l_{1} or l_{2} (m_{1} or m_{2}) can be zero, the assumption A3 allows the case in which there are only two or three distinct eigenvalues. From the relations (3.0.1) and (3.0.2) we obtain

(3.1.1)  S^{2}(A_{1}) + S^{2}(A_{2}) = δ^{2}Γ_{1}^{4} ≤ δ^{2}Γ^{4},
(3.1.2)  ∥A_{3}∥^{2} = δ^{2}Γ_{2}^{6} ≤ δ^{2}Γ^{6}.

Hence S^{2}(A) ≤ δ^{2}Γ^{4}(1 + Γ^{2}), (3.1.3)

and the assumption A1 implies

S(A) < \frac{\sqrt{5}}{8}δ. (3.1.4)

Thus, the conditions S(A) < \frac{3}{8}δ and S(A) < \frac{\sqrt{2}}{2}δ used in [15, 9] and [7] for proving quadratic convergence of the serial Jacobi methods are fulfilled. Since the estimates in [7] are more accurate than those from [9], we shall frequently refer to them.

Let λ_{1}, …, λ_{d} be distinct eigenvalues of A with multipicities n_{1}, …, n_{d}, respectively. We can assume 2 ≤ d ≤ n since d = 1 implies A = λ_{1}I_{n} and there is nothing to be proven. Applying the perturbation theorem for the eigenvalues of symmetric matrices to A and its diagonal part, and using the assumption A1 we obtain (cf. [7]) that for each i, |a_{jj} − λ_{i}| < (\sqrt{10}/8)δ < 0.4δ holds for exactly n_{i} values of j. Now, the assumption A2 implies that these n_{i} diagonal elements occupy successive positions on the diagonal. Hence, we can number the eigenvalues so that in the partition

(3.1.5)  A = (A_{ij}) = [ \begin{matrix}A_{11} & ⋯ & A_{1d} \\ ⋮ & & ⋮ \\ A_{d1} & ⋯ & A_{dd} \end{matrix} ] \begin{matrix}\} & n_{1} \\ \} & n_{d} \end{matrix}

\overbrace{n_{1}}  \overbrace{n_{d}}

the diagonal elements of A_{ii} are affiliated with λ_{i}. In such a situation we say that the block A_{ii} is associated with λ_{i}, 1 ≤ i ≤ d. In [15, 9, 7] it has been proved that the affiliation is preserved under any pivot strategy. By (A_{ij}^{(k)}) we denote the partition (3.1.5) of A^{(k)}, k ≥ 0 .

Lemma 3.1 Let the assumptions A1 and A2 hold. Let A = A^{(0)}, A^{(1)}, … be generated by any pivot strategy. Then for k ≥ 0

S(A_{ii}^{(k)}) ≤ \frac{\sqrt{2}}{2}∥A_{ii}^{(k)} − λ_{i}I_{n_{i}}∥ ≤ 0.26\frac{1}{δ} ∑_{j≠i^{j=1}}^{d} ∥A_{ij}^{(k)}∥^{2}, 1 ≤ i ≤ d.

Proof. The proof follows the lines of the proof of [7, Lemma 2.1]. For i ∈ {1, …, d} let P_{i} be the permutation matrix satisfying

A^{(k)} − λ_{i}I_{n} = P_{i} [ \begin{matrix}A_{ii}^{(k)} − λ_{i}I_{n_{i}} & G_{i}^{(k)} \\ G_{i}^{(k)}^{T} & K_{i}^{(k)} \end{matrix} ] P_{i}^{T} .

The rank argument (see [16]) implies

A_{ii}^{(k)} − λ_{i}I_{n_{i}} = G_{i}^{(k)}K_{i}^{(k)}^{−1}G_{i}^{(k)}^{T} (3.1.6)

and we have (see [6])

∥K_{i}^{(k)}^{−1}∥_{2} ≤ \frac{1}{3δ − ∥G_{i}^{(k)}∥_{2}} ≤ \frac{1}{3δ − S(A)}. (3.1.7)

Hence the assertion follows from the relations (3.1.6), (3.1.7) and (3.1.4) .

From Lemma 3.1 it is obvious that the partition (3.1.5) can be “read” from each matrix A^{(k)}, k ≥ 0. Namely, |a_{p(k)p(k)}^{(k)}−a_{q(k)q(k)}^{(k)}| ≤ 2·0.26·\sqrt{2}· \frac{5}{64} ·δ < 0.0575δ if a_{p(k)p(k)}^{(k)} and a_{q(k)q(k)}^{(k)} belong to the same diagonal block and |a_{pp}^{(k)} − a_{qq}^{(k)}| > 3δ − 0.0575δ > 2.9425δ if they belong to different diagonal blocks. The gap δ can be approximately determined from a_{ii}^{(k)}, 1 ≤ i ≤ n.

Let

K_{t} = {k ≥ t : |a_{p(k)p(k)}^{(k)} − a_{q(k)q(k)}^{(k)}| > 2.9425δ}, t ≥ 0.

Then the relation (1.2) implies

(3.1.8) | sin ϕ_{k}| ≤ | tan ϕ_{k}| ≤ \frac{1}{2}| tan 2ϕ_{k}| ≤ 0.34\frac{|a_{p(k)q(k)}^{(k)}|}{δ}, k ∈ K_{0}.

Now the relations (1.4) and (3.1.8) imply

∑_{k∈K_{t}} sin^{2} ϕ_{k} ≤ (0.34)^{2} \frac{1}{δ^{2}} ∑_{k∈K_{t}}(a_{pq}^{(k)})^{2} ≤ [0.34\frac{S(A^{(t)})}{δ}]^{2} , t ≥ 0. (3.1.9)

The bound in (3.1.9) can be refined if we know that pivot elements are taken just from a given principal submatrix of A^{(t)}. Suppose in the partition

A^{(k)} = [ \begin{matrix}A_{11}^{(k)} & A_{12}^{(k)} & A_{13}^{(k)} \\ A_{12}^{(k)}^{T} & A_{22}^{(k)} & A_{23}^{(k)} \\ A_{13}^{(k)}^{T} & A_{23}^{(k)}^{T} & A_{33}^{(k)} \end{matrix} ] , k ≥ 0, (3.1.10)

A_{22}^{(k)} has dimension ≥ 2, while A_{11}^{(k)} and/or A_{33}^{(k)} can be empty (i.e. of dimension 0). Then we have

Lemma 3.2 Let the assumptions A1, A2 hold, and let A = A^{(0)}, …, A^{(t)}, A^{(t+1)}, …,A^{(r)} be generated by r Jacobi rotations using any pivot strategy. If all pivot elements a_{p(k)q(k)}^{(k)}, t ≤ k ≤ r − 1 belong to A_{22}^{(k)}, then

∑_{k∈K_{t}^{k=t}}^{r−1} sin^{2} ϕ_{k} ≤ [0.34\frac{S(A_{22}^{(t)})}{δ}]^{2} .

Proof. Since the relation (1.4) holds (for t ≤ k ≤ r − 1) with A_{22} in place of A, we have

∑_{k=t}^{r−1}(a_{p(k)q(k)}^{(k)})^{2} ≤ S^{2}(A_{22}^{(t)}).

Hence the proof follows from the relation (3.1.8).

Next we provide an auxiliary result which will be used in analyzing the effect of Jacobi rotations following the orderings O_{i}, i = 1, 2, 4, 5. We use the partition (3.1.10), but we further assume

A_{22}^{(k)} = [ \begin{matrix}A_{ff}^{(k)} & ⋯ & A_{fg}^{(k)} \\ ⋮ & & ⋮ \\ A_{gf}^{(k)} & ⋯ & A_{gg}^{(k)} \end{matrix} ] , 1 ≤ f ≤ g ≤ d, k ≥ 0,

where A_{ij}^{(k)} f ≤ i, j ≤ g, are from the partition (3.1.5). Note that the dimension of A_{22}^{(k)} is n_{f} + … + n_{g}. Let

n̄_{i} = n_{1} + n_{2} + … + n_{i}, 1 ≤ i ≤ d.

With A_{22} is associated the set of indices

T = {(r, s) : n̄_{f−1} + 1 ≤ r < s ≤ n̄_{g}}.

Lemma 3.3 Let the assumptions A1 and A2 hold. Let A = A^{(0)}, A^{(1)}, ⋯, A^{(ν(T))} be generated by applying to A ν(T) Jacobi rotations according to the ordering O_{R}(T). Then

S(A_{22}^{(ν(T))}) ≤ 0.625\frac{S^{2}(A_{22})}{δ} + 0.285\frac{∥A_{12}∥^{2} + ∥A_{23}∥^{2}}{δ} .

Proof. See Appendix.

3.2 Auxiliary Lemmas

In this subsection we derive basic tools for analyzing the effect of Jacobi rotations following the orderings O_{3}, O_{6}, and O_{7}. The results are independent of the asymptotic assumptions.

Let H = (h_{ij}) be a symmetiric matrix of order r ≥ 2, partitioned as follows

H = [ \begin{matrix}H_{1} & H_{3} \\ H_{3}^{T} & H_{2} \end{matrix} ] \begin{matrix}\} & u \\ \} & v \end{matrix} .
\overbrace{u} \overbrace{v} Let
Q = {(i, j) : 1 ≤ i ≤ u, u + 1 ≤ j ≤ u + v = r},
q_{α} = {(α, j) : u + 1 ≤ j ≤ u + v}, α = 1, 2, …, u.

We shall apply ν(Q) Jacobi rotations to H according to the ordering O_{R}(Q). Our goal is to find bounds for S(H_{1}^{(ν(Q))}), S(H_{2}^{(ν(Q))}) and ∥ H_{3}^{(ν(Q))} ∥ in terms of S(H_{1}), S(H_{2}), S(H) and the angles. To this end we denote the angle ϕ_{k} which rotates in the (p(k), q(k)) plane by ϕ_{pq}. In the sequel we shall frequently use the quantities Σ_{α}, Σ defined by

Σ_{α}^{2} = ∑_{i=1}^{α} ∑_{j=1}^{v} sin^{2} ϕ_{i,u+j}, 1 ≤ α ≤ u , Σ = Σ_{u},

and the inequality

(x + y)^{2} ≤ (1 + γ)x^{2} + (1 + \frac{1}{γ})y^{2}, γ > 0. (3.2.1)

In (3.2.1) equality is attained for γ = y/x.

Lemma 3.4 If we apply to H a sequence of ν(Q) Jacobi rotations according to the ordering O_{R}(Q), then for all 0 ≤ k ≤ ν(Q) holds

S(H_{1}^{(k)}) ≤ S(H_{1}) + \sqrt{2}S(H)Σ.

Proof. By an analysis similar to that used by Wilkinson [15], one obtains

|h_{ij}^{(k)}| ≤ Ω_{ij}, 1 ≤ i < j ≤ u, 0 ≤ k ≤ ν(Q), (3.2.2) where

Ω_{ij} = |h_{ij}| + ∑_{t=1}^{v} |h_{j,u+t}^{(v(i−1))}|| sin ϕ_{i,u+t}| + ∑_{t=1}^{v} |h_{i,u+t}^{(v(j−1))}|| sin ϕ_{j,u+t}|. (3.2.3)

For during the process the element in the position (i, j) with 1 ≤ i < j ≤ u, changes only when the elements in the positions defined by the sets q_{i} and q_{j} are rotated. Furthermore, the contributions to h_{ij} come from the elements of H_{3}^{(v(i−1))} and H_{3}^{(v(j−1))} which are in the positions defined by q_{j} and q_{i}, respectively ( see the Figure 1 ).

[figure: x; x; x; x; x; y; y; y; y; y; j; i; i; j; h_{ij}; q_{j}; q_{i}; These elements are; from H_{3}^{(v(i−1))}; These elements are; from H_{3}^{(v(j−1))}]

Figure 1: The case of u = 6, v = 5 and (i, j) = (2, 5)

Using the inequality (3.2.1) twice (the second time with γ = 1 ) and afterwards using the inequality of Schwarz and Cauchy, one obtains

Ω_{ij}^{2} ≤ (1 + γ)(h_{ij})^{2} + 2(1 + \frac{1}{γ}) [∑_{t=1}^{v}(h_{j,u+t}^{(v(i−1))})^{2} ∑_{t=1}^{v} sin^{2} ϕ_{i,u+t}

+ ∑_{t=1}^{v}(h_{i,u+t}^{(v(j−1))})^{2} ∑_{t=1}^{v} sin^{2} ϕ_{j,u+t}] , 1 ≤ i < j ≤ u.

From the above inequality and the relation (3.2.2) we obtain for 0 ≤ k ≤ ν(Q)

S^{2}(H_{1}^{(k)}) ≤ ∑_{1≤i<j≤u} Ω_{ij}^{2} ≤ (1 + γ)S^{2}(H_{1}) + 2(1 + \frac{1}{γ}) \{∑_{i=1}^{u−1} ∑_{j=i+1}^{u} (3.2.4)

[∑_{t=1}^{v}(h_{j,u+t}^{(v(i−1))})^{2} ∑_{t=1}^{v} sin^{2} ϕ_{i,u+t}] + ∑_{j=2}^{u} ∑_{i=1}^{j−1} [∑_{t=1}^{v}(h_{i,u+t}^{(v(j−1))})^{2} ∑_{t=1}^{v} sin^{2} ϕ_{j,u+t}]} ,

If we let i = s in the first summation and j = s in the second summation of the inequality (3.2.4), and then add them together, we obtain

∑_{j=s+1}^{u} [∑_{t=1}^{v}(h_{j,u+t}^{(v(s−1))})^{2} ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t}] + ∑_{i=1}^{s−1} [∑_{t=1}^{v}(h_{i,u+t}^{(v(s−1))})^{2} ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t}] (3.2.5)

= [ ∑_{j=s+1}^{u} ∑_{t=1}^{v}(h_{j,u+t}^{(v(s−1))})^{2} + ∑_{i=1}^{s−1} ∑_{t=1}^{v}(h_{i,u+t}^{(v(s−1))})^{2}] ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t}

≤ [∑_{i=1}^{u} ∑_{t=1}^{v}(h_{i,u+t}^{(v(s−1))})^{2}] ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t} =∥ H_{3}^{(v(s−1))} ∥^{2} ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t}

≤ S^{2}(H) ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t}.

The relation (3.2.5) holds for 2 ≤ s ≤ u − 1. It also holds for s ∈ {1, u} provided the empty sums are considered zero. Hence Lemma 3.4 follows from the inequalities (3.2.4), (3.2.5) and by specifying the minimizing γ.

Lemma 3.5 If we apply to H a sequence of ν(Q) Jacobi rotations according to the ordering O_{R}(Q), then for all 0 ≤ k ≤ ν(Q) and any 0 ≤ α ≤ u holds

S(H_{2}^{(αv)}) ≤ S(H_{2}) + \sqrt{2}S(H)Σ_{α}.

Proof. If α = 0 we define Σ_{0} = 0, hence the assertion holds. To estimate S^{2}(H_{2}^{(αv)}) for 1 ≤ α ≤ u we again use the idea of equivalent orderings. Let

Q_{α} = q_{1} ∪ q_{2} ∪ ⋯ ∪ q_{α}, 1 ≤ α ≤ u.

Using the same argument as in [7, Lemma 1.2 (ii)] one can easily prove that O_{R}(Q_{α}) ∼ O_{C}(Q_{α}) for all 1 ≤ α ≤ u.

Let α be such that 1 ≤ α ≤ u. We shall apply to H αv Jacobi rotations according to the sequence O_{C}(Q_{α}). Let H̃^{(k)} = (h̃_{ij}^{(k)}), 1 ≤ k ≤ αv, be matrices obtained in such a way. Let ϕ̃_{s,u+t} be the corresponding angles. By [7, Lemma 1.1] we conclude that H̃^{(αv)} = H^{(αv)} and ϕ̃_{s,u+t} = ϕ_{s,u+t} for all 1 ≤ s ≤ α, 1 ≤ t ≤ v. By the same lemma we have . This means that to estimate S^{2}(H_{2}^{(αv)}) we can use Jacobi rotations according to the sequence O_{C}(Q_{α}).

Using an analysis quite similar to that in the proof of Lemma 3.4 we obtain (see Figure 2)

|h̃_{u+i,u+j}^{(αv)}| ≤ Ω̃_{u+i,u+j}^{(α)}, 1 ≤ i < j ≤ v,

Ω̃_{u+i,u+j}^{(α)} = |h̃_{u+i,u+j}| + ∑_{s=1}^{α} |h̃_{s,u+j}^{(α(i−1))}|| sin ϕ_{s,u+i}| + ∑_{s=1}^{α} |h̃_{s,u+i}^{(α(j−1))}|| sin ϕ_{s,u+j}|.

[figure: y; y; y; y; x; x; x; x; u + i; u + j; u + i; u + j; *; h_{8,10} = h_{u+i,u+j}; These elements are; from H̃_{3}^{(α(j−1))}; These elements are; from H̃_{3}^{(α(i−1))}; α]
Figure 2: The case of u = 6, v = 5, i = 2 , j = 4, α = 4

Following the lines in the proof of Lemma 3.4, we obtain

S^{2}(H_{2}^{(αv)}) = S^{2}(H̃_{2}^{(αv)}) ≤ ∑_{1≤i<j≤u}(Ω̃_{ij}^{(α)})^{2} ≤ (1 + γ)S^{2}(H_{2}) + 2(1 + \frac{1}{γ}) ·

{∑_{i=1}^{v−1} ∑_{j=i+1}^{v} [∑_{s=1}^{α}(h̃_{s,u+j}^{(α(i−1))})^{2} ∑_{s=1}^{α} sin^{2} ϕ_{s,u+i}] + ∑_{j=2}^{v} ∑_{i=1}^{j−1} [∑_{s=1}^{α}(h̃_{s,u+i}^{(α(j−1))})^{2} ∑_{s=1}^{α} sin^{2} ϕ_{s,u+j}]} .

If we let i = t in the first summation and j = t in the second, then further following the proof of Lemma 3.4 we obtain

S^{2}(H_{2}^{(αv)}) ≤ (1 + γ)S^{2}(H_{2}) + 2(1 + \frac{1}{γ}) ∑_{t=1}^{v} ∥H̃_{3}^{(α(t−1))}∥^{2} ∑_{s=1}^{α} sin^{2} ϕ_{s,u+t}
≤ (1 + γ)S^{2}(H_{2}) + 2(1 + \frac{1}{γ})S^{2}(H) ∑_{s=1}^{α} ∑_{t=1}^{v} sin^{2} ϕ_{s,u+t}.

The lemma follows by specifying the minimizing γ.

Lemma 3.6 If we apply to H a sequence of ν(Q) Jacobi rotations according to the ordering O_{R}(Q), then

∥H_{3}^{(ν(Q))}∥ ≤ [S(H_{1}) + S(H_{2}) + 3.15S(H)Σ] Σ.

Proof. Using an analysis similar to those in the former lemmas, we obtain for all 1 ≤ i ≤ u, 1 ≤ j ≤ v (see Figure 3)

|h_{i,u+j}^{(ν(Q))}| ≤ ∑_{s=i+1}^{u} |h_{is}^{((s−1)v+j−1)}|| sin ϕ_{s,u+j}| + ∑_{t=j+1}^{v} |h_{u+j,u+t}^{((i−1)v+j)}|| sin ϕ_{i,u+t}|.

[figure: y; y; y; y; p; p; p; p; p; p; p; x; x; x; u + j; i; i; u + j; h_{28}; Pivot positions of the; transformations which; change h_{28} = h_{i,u+j} after; its annihilation; These elements are; from H_{2}^{((i−1)v+j)}]

Figure 3: The case of u = 6, v = 5, i = 2, j = 2

and

(h_{i,u+j}^{(ν(Q))})^{2} ≤ (1 + γ) ∑_{s=i+1}^{u}(h_{is}^{((s−1)v+j−1)})^{2} ∑_{s=i+1}^{u} sin^{2} ϕ_{s,u+j}

+(1 + \frac{1}{γ}) ∑_{t=j+1}^{v}(h_{u+j,u+t}^{((i−1)v+j)})^{2} ∑_{t=j+1}^{v} sin^{2} ϕ_{i,u+t}.

Hence

∥H_{3}^{(ν(Q))}∥^{2} ≤ (1 + γ) ∑_{j=1}^{v} {[∑_{i=1}^{u−1} ∑_{s=i+1}^{u}(h_{is}^{((s−1)v+j−1)})^{2}] ∑_{s=2}^{u} sin^{2} ϕ_{s,u+j}} (3.2.6)

+ (1 + \frac{1}{γ}) ∑_{i=1}^{u} \{[∑_{j=1}^{v−1} ∑_{t=j+1}^{v}(h_{u+j,u+t}^{((i−1)v+j)})^{2}] ∑_{t=2}^{v} sin^{2} ϕ_{i,u+t}\} .

To bound the first term on the right-hand side of the relation (3.2.6) we can use the relation (3.2.2) which implies

| h_{is}^{((s−1)v+j−1)} |≤ Ω_{is}, 1 ≤ j ≤ v, 1 ≤ i < s ≤ u,

where Ω_{is} is given by the relation (3.2.3). Now the inequalities (3.2.4) and (3.2.5) and the choice of the optimal γ imply

∑_{i=1}^{u−1} ∑_{s=i+1}^{u}(h_{is}^{((s−1)v+j−1)})^{2} ≤ ∑_{1≤i<s≤u} Ω_{is}^{2} ≤ (S(H_{1}) + \sqrt{2}S(H)Σ)^{2} . (3.2.7)

To bound the second term in the relation (3.2.6), we note that h_{u+j,u+t}^{((i−1)v)} with 1 ≤ j ≤ v − 1, j + 1 ≤ t ≤ v, does not change when the elements at positions (i, u + 1), (i, u + 2), …, (i, u + j − 1) are rotated. Therefore we have

|h_{u+j,u+t}^{((i−1)v+j)}| ≤ |h_{u+j,u+t}^{((i−1)v)}| + |h_{i,u+t}^{((i−1)v+j−1)}| | sin ϕ_{i,u+j}|, \begin{matrix}1 ≤ j ≤ v − 1, \\ j + 1 ≤ t ≤ v \end{matrix} .

Applying (3.2.1) to the right side of the above inequality, then summing over all 1 ≤ j ≤ v − 1, j + 1 ≤ t ≤ v, and finally using Lemma 3.5 and the inequality ∑_{t=j+1}^{v}(h_{i,u+t}^{((i−1)v+j−1)})^{2} ≤ S^{2}(H) we obtain

∑_{j=1}^{v−1} ∑_{t=j+1}^{v}(h_{u+j,u+t}^{((i−1)v+j)})^{2} ≤ (1 + β)S^{2}(H_{2}^{((i−1)v)}) + (1 + \frac{1}{β}) ∑_{j=1}^{v−1} [sin^{2}ϕ_{i,u+j}S^{2}(H)] (3.2.8)

≤ (1 + β) [S(H_{2}) + \sqrt{2}S(H)Σ_{i−1}]^{2} + (1 + \frac{1}{β})S^{2}(H) ∑_{j=1}^{v−1} sin^{2} ϕ_{i,u+j}.

Applying the inequality (3.2.1) with γ = 1/2 we have

(3.2.9)  (\sqrt{2}Σ_{i−1} + \sqrt{∑_{j=1}^{v−1} sin^{2} ϕ_{i,u+j}})^{2} ≤ 3 (Σ_{i−1}^{2} + ∑_{j=1}^{v−1} sin^{2} ϕ_{i,u+j})
= 3Σ_{i}^{2} ≤ 3Σ^{2}.

Specifying the minimizing β in (3.2.8) and using (3.2.9) we obtain

∑_{j=1}^{v−1} ∑_{t=j+1}^{v}(h_{u+j,u+t}^{((i−1)v+j)})^{2} ≤ [S(H_{2}) + \sqrt{3}S(H)Σ]^{2} . (3.2.10)

Inserting (3.2.7) and (3.2.10) in (3.2.6), specifying the minimizing γ and using the inequality \sqrt{2} + \sqrt{3} < 3.15, we obtain the assertion of Lemma 3.6.

3.3 Main Results

Here we prove the cubic convergence of the method J_{Θ}. In the proof we shall frequently use the partition (1.6) of A^{(k)},

A^{(k)} = [ \begin{matrix}A_{1}^{(k)} & A_{3}^{(k)} \\ A_{3}^{(k)}^{T} & A_{2}^{(k)} \end{matrix} ] \begin{matrix}\} l \\ \} m \end{matrix} =[ \begin{matrix}\begin{matrix}B_{11}^{(k)} & B_{12}^{(k)} \\ B_{12}^{(k)}^{T} & B_{22}^{(k)} \end{matrix} & \frac{B_{1}^{(k)}}{B_{2}^{(k)}} \\ \frac{B_{3}^{(k)}}{B_{4}^{(k)}} & \begin{matrix}B_{33}^{(k)} & B_{34}^{(k)} \\ B_{34}^{(k)}^{T} & B_{44}^{(k)} \end{matrix} \end{matrix} ] \begin{matrix}\} l_{1} \\ \} l_{2} \\ \} m_{1} \\ \} m_{2} \end{matrix} . (3.3.1)

\overbrace{l} \overbrace{m}  \overbrace{l_{1}} \overbrace{l_{2}} \overbrace{m_{1}} \overbrace{m_{2}}

Recall that the sets of pairs P_{1}, P_{2}, P_{3}, P_{4}, P_{5}, P_{6} and P_{7} are associated with the blocks B_{11}^{(k)}, B_{22}^{(k)}, B_{12}^{(k)}, B_{33}^{(k)}, B_{44}^{(k)}, B_{34}^{(k)} and A_{3}^{(k)}, respectively, and that O_{i} = O_{R}(P_{i}), 1 ≤ i ≤ 7. For simplicity we write ν_{i} for ν(P_{i}), 1 ≤ i ≤ 7.

Theorem 3.7 Let A be a symmetric matrix satisfying the assumptions A1, A2 and A3. Let the matrices A^{(0)} = A, A^{(1)}, …, A^{(M)} be obtained by applying the method J_{Θ} to A. Then

Γ^{(M)} < 1.35 Γ^{3}.

Proof: First, we assume that none of l_{1}, l_{2}, m_{1} and m_{2} is zero. Recall that in one quasi-sweep of the method J_{Θ} Jacobi rotations are applied according to the sequence

Θ = [O_{1}, O_{2}, O_{3}, O_{1}, O_{2}, O_{4}, O_{5}, O_{6}, O_{4}, O_{5}, O_{7}].

Since we consider what happens after one quasi-sweep, we can use instead of Θ any sequence which is equivalent to it. For our purpose we replace Θ by

Θ′′′ = [Θ′′′_{1}, Θ′′′_{2}, Θ′′′_{3}, Θ′′′_{4}]

where

\begin{matrix}Θ′′′_{1} = [O_{1}, O_{2}, O_{4}, O_{5}], \\ Θ′′′_{2} = [O_{3}, O_{6}], \end{matrix} \begin{matrix}Θ′′′_{3} = [O_{1}, O_{2}, O_{4}, O_{5}], \\ Θ′′′_{4} = O_{7}. \end{matrix}

From now on we assume that matrices A^{(1)}, …,A^{(M)} are generated by J_{Θ′′′}. With each subsequence Θ′′′_{i} we associate a stage of the process.

Stage 1 Here Jacobi rotations are defined by the sequence Θ′′′_{1}. By Lemma 3.3 and using the partition (3.3.1) we obtain

S(B_{11}^{(ν_{1})}) ≤ 0.625\frac{S^{2}(B_{11})}{δ} + 0.285\frac{∥B_{12}∥^{2} + ∥B_{1}∥^{2}}{δ} ,
S(B_{22}^{(ν_{1}+ν_{2})}) ≤ 0.625\frac{S^{2}(B_{22}^{(ν_{1})})}{δ} + 0.285\frac{∥B_{12}^{(ν_{1})}∥^{2} + ∥B_{2}^{(ν_{1})}∥^{2}}{δ} .

Since B_{11}^{(ν_{1})} (B_{22}^{(ν_{1}+ν_{2})}) is invariant under Jacobi rotations defined by [O_{2}, O_{4}, O_{5}]

([O_{4}, O_{5}]), we can replace it by B_{11}^{(k_{1})} (B_{22}^{(k_{1})}), where

k_{1} = ν_{1} + ν_{2} + ν_{4} + ν_{5}.

Note that B_{22}^{(k)}, ∥B_{12}^{(k)}∥ and B_{2}^{(k)} remain unchanged for 0 ≤ k ≤ ν_{1}. So, we have

S(B_{11}^{(k_{1})}) + S(B_{22}^{(k_{1})}) ≤ 0.625 \frac{S^{2}(B_{11}) + S^{2}(B_{22})}{δ} + 0.285\frac{2∥B_{12}∥^{2} + ∥A_{3}∥^{2}}{δ}
≤ 0.625 \frac{S^{2}(A_{1})}{δ} + 0.285\frac{∥A_{3}∥^{2}}{δ} .

Similarly, we obtain

S(B_{33}^{(k_{1})}) + S(B_{44}^{(k_{1})}) ≤ 0.625 \frac{S^{2}(A_{2})}{δ} + 0.285\frac{∥A_{3}∥^{2}}{δ} .

By the definition of Γ_{1} and Γ_{2} we obtain

∑_{i=1}^{4} S(B_{ii}^{(k_{1})}) ≤ 0.625\frac{S^{2}(A_{1}) + S^{2}(A_{2})}{δ} + 0.57\frac{∥A_{3}∥^{2}}{δ} ≤ 0.625Γ_{1}^{4}δ + 0.57Γ_{2}^{6}δ.

(3.3.2) Since for 0 ≤ k ≤ k_{1}, ∥A_{3}^{(k)}∥ is invariant, and S(A_{1}^{(k)}), S(A_{2}^{(k)}) are nonincreasing, we have

\begin{aligned}S^{2}(A_{1}^{(k_{1})}) + S^{2}(A_{2}^{(k_{1})}) & ≤ & Γ_{1}^{4}δ^{2} \\ ∥A_{3}^{(k_{1})}∥^{2} & ≤ & Γ_{2}^{6}δ^{2} \end{aligned} . (3.3.3)

Stage 2 In this stage we consider the effect of Jacobi rotations defined by the sequence Θ′′′_{2} = [O_{3}, O_{6}].

We use Lemma 3.4 and Lemma 3.5 with H = A_{1}^{(k_{1})}, H_{1} = B_{11}^{(k_{1})}, H_{2} = B_{22}^{(k_{1})}, H_{3} = B_{12}^{(k_{1})} and ν(Q) = ν_{3}. Together with the inequality (3.2.1) we obtain

S^{2}(B_{11}^{(k_{1}+ν_{3})}) + S^{2}(B_{22}^{(k_{1}+ν_{3})}) ≤ (1 + γ) [S^{2}(B_{11}^{(k_{1})}) + S^{2}(B_{22}^{(k_{1})})]

+ 2(1 + \frac{1}{γ})2S^{2}(A_{1}^{(k_{1})}) ∑_{s=1}^{l_{1}} ∑_{t=1}^{l_{2}} sin^{2} ϕ_{s,l_{1}+t}.

Because of the invariance property we can replace the superscript k_{1} + ν_{3} by k_{2} = k_{1} + ν_{3} + ν_{6}. Using Lemma 3.2 we obtain

S^{2}(B_{11}^{(k_{2})})+S^{2}(B_{22}^{(k_{2})}) ≤ (1+γ) [S^{2}(B_{11}^{(k_{1})}) + S^{2}(B_{11}^{(k_{1})})]+4 (1+ \frac{1}{γ}) 0.34^{2} \frac{S^{4}(A_{1}^{(k_{1})})}{δ^{2}}.

In a similar way one obtains

S^{2}(B_{33}^{(k_{2})})+S^{2}(B_{44}^{(k_{2})}) ≤ (1+γ) [S^{2}(B_{33}^{(k_{1})}) + S^{2}(B_{44}^{(k_{1})})]+4 (1+ \frac{1}{γ}) 0.34^{2} \frac{S^{4}(A_{2}^{(k_{1})})}{δ^{2}}.

Making sum of the left and right sides of these inequalities, specifying an optimal γ and using the relations (3.3.2), (3.3.3) and the assumption A1 we obtain

∑_{i=1}^{4} S^{2}(B_{ii}^{(k_{2})}) ≤ [∑_{i=1}^{4} S(B_{ii}^{(k_{1})}) + 0.68\frac{S^{2}(A_{1}^{(k_{1})}) + S^{2}(A_{2}^{(k_{1})})}{δ}]^{2} (3.3.4)

≤ [(0.625 + 0.68)Γ_{1}^{4}δ + 0.57Γ_{2}^{6}δ]^{2} ≤ (1.305 + 0.57 · \frac{1}{4})^{2} Γ^{8}δ^{2} ≤ 2.1 Γ^{8}δ^{2}.

To bound ∥B_{12}^{(k_{2})}∥ and ∥B_{34}^{(k_{2})}∥ we use Lemma 3.6 and Lemma 3.2 together with the assumption A3. We obtain

∥B_{12}^{(k_{2})}∥ ≤ (S(B_{11}^{(k_{1})}) + S(B_{22}^{(k_{1})}) + 3.15S(A_{1}^{(k_{1})}) · 0.34\frac{S(A_{1}^{(k_{1})})}{δ}) · 0.34\frac{S(A_{1}^{(k_{1})})}{δ},

∥B_{34}^{(k_{2})}∥ ≤ (S(B_{33}^{(k_{1})}) + S(B_{44}^{(k_{1})}) + 3.15S(A_{2}^{(k_{1})}) · 0.34\frac{S(A_{2}^{(k_{1})})}{δ}) · 0.34\frac{S(A_{2}^{(k_{1})})}{δ}.

Hence by (3.3.2), (3.3.3) and the assumption A1 we have

∥B_{12}^{(k_{2})}∥ + ∥B_{34}^{(k_{2})}∥ ≤ [(0.625 + 1.071)\frac{S^{2}(A_{1}) + S^{2}(A_{2})}{δ} + 0.57\frac{∥A_{3}∥^{2}}{δ} ] · (3.3.5)

·0.34 max \{\frac{S(A_{1})}{δ}, \frac{S(A_{2})}{δ}\} ≤ (0.5767Γ_{1}^{4}δ + 0.1938Γ_{2}^{6}δ) · Γ_{1}^{2} ≤ 0.63 Γ^{6}δ.

From the relations (3.3.4), (3.3.5) and the assumption A1 we obtain

S^{2}(A_{1}^{(k_{2})}) + S^{2}(A_{2}^{(k_{2})}) ≤ (2.1 + 0.63^{2}Γ^{4}) Γ^{8}δ^{2} ≤ 2.13Γ^{8}δ^{2}. (3.3.6)

Note that the invariance properties and the relations (3.3.5) and (3.3.6) imply

\begin{aligned}∥A_{3}^{(k_{2})}∥^{2} & ≤ & Γ_{2}^{6}δ^{2}, \\ S^{2}(A^{(k_{2})}) & ≤ & (2.13 · \frac{1}{4} + 1)Γ^{6}δ^{2} ≤ 1.54Γ^{6}δ^{2}. \end{aligned} (3.3.7)

Stage 3 Since Θ′′′_{3} = Θ′′′_{1}, this stage is quite analogous to Stage 1. Let k_{3} = k_{2} + ν_{1} + ν_{2} + ν_{4} + ν_{5}. Using the same argument as in Stage 1 (cf. (3.3.2)) and using the relations (3.3.6) and (3.3.7) we obtain

∑_{i=1}^{4} S(B_{ii}^{(k_{3})}) ≤ 0.625\frac{S^{2}(A_{1}^{(k_{2})}) + S^{2}(A_{2}^{(k_{2})})}{δ} + 0.57\frac{∥A_{3}^{(k_{2})}∥^{2}}{δ} (3.3.8)

≤ 0.625 · 2.13 Γ^{8} δ + 0.57Γ_{2}^{6} δ ≤ (1.332 · \frac{1}{4} + 0.57)Γ^{6} δ ≤ 0.903 Γ^{6} δ.

Because of the invariance of ∥B_{12}^{(k)}∥, ∥B_{34}^{(k)}∥ and ∥A_{3}^{(k)}∥, the relations (3.3.8), (3.3.5) and (3.3.7) imply

\{ \begin{matrix}S^{2}(A_{1}^{(k_{3})}) + S^{2}(A_{2}^{(k_{3})}) & ≤ & 0.903^{2}Γ^{12} δ^{2} + 0.63^{2}Γ^{12} δ^{2} ≤ 1.22Γ^{12} δ^{2}, \\ ∥A_{3}^{(k_{3})}∥^{2} & ≤ & Γ_{2}^{6} δ^{2}, \\ S^{2}(A^{(k_{3})}) & ≤ & Γ_{2}^{6} δ^{2} + 1.22 · \frac{1}{64}Γ^{6} δ^{2} ≤ 1.02Γ^{6} δ^{2}. \end{matrix} (3.3.9)

Stage 4 Since Θ′′′_{4} = O_{7} this stage is analogous to Stage 2. Note that M = k_{3} + ν_{7}. We use Lemma 3.4 and Lemma 3.5 with H = A^{(k_{3})}, H_{1} = A_{1}^{(k_{3})}, H_{2} = A_{2}^{(k_{3})}, H_{3} = A_{3}^{(k_{3})} and ν(Q) = ν_{7}. By the inequality (3.2.1) and the relation (3.3.9) we obtain

S^{2}(A_{1}^{(M)}) + S^{2}(A_{2}^{(M)}) ≤ (1 + γ) [S^{2}(A_{1}^{(k_{3})}) + S^{2}(A_{2}^{(k_{3})})]

+4(1 + \frac{1}{γ})0.34^{2} \frac{S^{4}(A^{(k_{3})})}{δ^{2}} ≤ (1 + γ)1.22Γ^{12}δ^{2} + (1 + \frac{1}{γ}) 0.68^{2} · 1.02^{2} Γ^{12}δ^{2}.

Specifying an optimal γ yields

S^{2}(A_{1}^{(M)}) + S^{2}(A_{2}^{(M)}) ≤ (\sqrt{1.22} + 0.68 · 1.02)^{2} Γ^{12}δ^{2} ≤ 3.24 Γ^{12}δ^{2}.

Hence

Γ_{1}^{(M)} ≤ ^{4} \sqrt{3.24} Γ^{3} ≤ 1.35 Γ^{3}. (3.3.10)

Using Lemma 3.6, Lemma 3.2 and the relation (3.3.9) we obtain

∥A_{3}^{(M)}∥ ≤ [S(A_{1}^{(k_{3})}) + S(A_{2}^{(k_{3})}) + 3.15S(A^{(k_{3})}) · 0.34\frac{S(A^{(k_{3})})}{δ}] · 0.34\frac{S(A^{(k_{3})})}{δ}

≤ [\sqrt{2}\sqrt{1.22}Γ^{6} δ + 1.071\frac{S^{2}(A^{(k_{3})})}{δ}] · 0.34\frac{S(A^{(k_{3})})}{δ}

≤ [(\sqrt{2.44} + 1.071 · 1.02) Γ^{6} δ] · 0.34 · \sqrt{1.02}Γ^{3} ≤ 0.92 Γ^{3} δ.

Hence

Γ_{2}^{(M)} ≤ ^{6} \sqrt{0.92} Γ^{3} ≤ Γ^{3}. (3.3.11)

Combining (3.3.10) and (3.3.11) we obtain the assertion of the theorem.

Now suppose that some of l_{1}, l_{2}, m_{1} and m_{2} is(are) zero(s). Note that l_{1} and l_{2} (m_{1} and m_{2}) cannot both be zero since the method J_{Θ} assumes that l and m are posive integers. There are several such cases. We will only consider the case when l_{1} = 0 and none of m_{1} and m_{2} is zero. Other cases can be proven in similar ways. Since l_{1} = 0, B_{11}, B_{12} and B_{1} in (2.2) do not exist, and hence O_{1} and O_{3} are empty sequences. We also note that in this case B_{22} = A_{1} and B_{2} = A_{3}. Thus one sweep of the method J_{Θ′′′} uses Jacobi rotations according to the sequence Θ′′′ = [O_{2}, O_{4}, O_{5}, O_{6}, O_{2}, O_{4}, O_{5}, O_{7}]. So, for 0 ≤ k ≤ M we have

A^{(k)} = [ \begin{matrix}A_{1}^{(k)} & A_{3}^{(k)} \\ A_{3}^{(k)}^{T} & A_{2}^{(k)} \end{matrix} ] \begin{matrix}\}l \\ \}m \end{matrix} =[ \begin{matrix}A_{1}^{(k)} & A_{3}^{(k)} \\ \frac{B_{3}^{(k)}}{B_{4}^{(k)}} & \begin{matrix}B_{33}^{(k)} & B_{34}^{(k)} \\ B_{34}^{(k)}^{T} & B_{44}^{(k)} \end{matrix} \end{matrix} ] \begin{matrix}\}l_{2} \\ \}m_{1} \\ \}m_{2} \end{matrix}

\overbrace{l} \overbrace{m}  \overbrace{l_{2}} \overbrace{m_{1}} \overbrace{m_{2}}

By simple inspection one can check the validity of the basic relations among (3.3.2)–(3.3.11). This completes the proof of Theorem 3.7.

From the relation (3.1.1), (3.1.2), (3.0.1), (3.0.2) and the assumption A1 we see that

Γ^{3} ≤ \sqrt{Γ_{1}^{4} + Γ_{2}^{6}} = \frac{S(A)}{δ} ≤ \frac{\sqrt{5}}{2}Γ^{2}. (3.3.12)

By Theorem 3.7 Γ reduces cubically per quasi-sweep. Hence S(A)/δ has upper and lower bound which converges to zero cubically per quasi-sweep.

Corollary 3.8 Let the assumpions of Theorem 3.7 hold.

(i) If Γ = Γ_{1}, then \frac{S(A^{(M)})}{δ} ≤ 1.83 (\frac{S(A)}{δ})^{3} .

(ii) If Γ = Γ_{2}, then \frac{S(A^{(M)})}{δ} ≤ 1.83 (\frac{S(A)}{δ})^{2} .

(iii) If Γ_{2} = γ Γ_{1}, then \frac{S(A^{(M)})}{δ} ≤ 1.83 max\{1, γ^{6}\} (\frac{S(A)}{δ})^{3} .

Proof: From the relations (3.3.10), (3.3.11), (3.3.12) and the assumption A1 we obtain

(3.3.13)  \frac{S(A^{(M)})}{δ} ≤ \sqrt{(Γ_{1}^{(M)})^{4} + (Γ_{2}^{(M)})^{6}} ≤ \sqrt{1.35^{4}Γ^{12} + Γ^{18}}
≤ \sqrt{1.35^{4} + 2^{−6}} Γ^{6} ≤ 1.83 Γ^{6}.

(i) Since Γ = Γ_{1} the relation (3.3.13) and the definition of Γ_{1} imply

\frac{S(A^{(M)})}{1.83 · δ} ≤ Γ^{6} = (Γ^{4})^{\frac{3}{2}} = (Γ_{1}^{4})^{\frac{3}{2}} = [\frac{S^{2}(A_{1}) + S^{2}(A_{2})}{δ^{2}}]^{\frac{3}{2}} ≤ [\frac{S^{2}(A)}{δ^{2}}]^{\frac{3}{2}} = \frac{S^{3}(A)}{δ^{3}}.

(ii) If Γ = Γ_{2}, the relation (3.3.13) and the definition of Γ_{2} imply

\frac{S(A^{(M)})}{1.83 · δ} ≤ Γ^{6} = Γ_{2}^{6} ≤ \frac{∥A_{3}∥^{2}}{δ^{2}} ≤ \frac{S^{2}(A)}{δ^{2}}.

(iii) Since Γ = max{1, γ} Γ_{1}, the proof is almost identical to the proof of (i).

From Corollary 3.8, if Γ = Γ_{2} and γ is not of moderate size the cubic reduction of S(A)/δ per quasi-sweep may fail to occur. This is demonstrated in Example 3.9 below. However, as discussion following Example 3.9 suggests, cubic convergence of Γ actually forces cubic convergence of S(A)/δ. Once the asymptotic assumptions are fulfilled, the quadratic reduction of S(A)/δ can occur just one time.

It is easy to see that if possible the choice of l_{1} ≈ [n/4], l_{2} ≈ [n/4], m_{1} ≈ [n/4], and m_{2} ≈ [n/4] will minimize the number of Jacobi rotations in one quasi-sweep of the method J_{Θ}. Here [n/4] denotes the integer part of n/4. If such partition is possible, one quasi-sweep of the method J_{Θ} approximately equals to 1.25 usual sweeps.

Numerical Examples

Here we show how well numerical tests comply with the theoretical results stated in Theorem 3.7 and Corollary 3.8. We also show that the assumption A3 is not only sufficient, but also necessary for the cubic convergence of J_{Θ} (on the class of symmetric matrices satisfying A1 and A2).

In the following we will consider three numerical examples. All the computations were done on VAX6000 Model 460 using quadruple precision roundoff unit eps ≃ 9.63 x 10^{−35}.

Example 3.9 Let n = 8, A = (a_{ij}), a_{ii} = −12 + 3i, 1 ≤ i ≤ 8 and a_{ij} = 10^{−4}, i ≠ j. δ. From the computed eigenvalues we have found out that

δ = 1.00000000055555555047582308931069,

to the arithmetic accuracy. Since all eigenvalues are distinct, there are several choices for l_{1}, l_{2}, m_{1} and m_{2}. We choose the optimal choice l_{1} = l_{2} = m_{1} = m_{2} = 2. Applying J_{Θ} to A we have computed Γ, γ (Corollary 3.8), and S(A)/δ after each quasi-sweep. In the table below the first row (when sweep = 0) corresponds to the initial matrix. Γ was computed and was found Γ = 7.369 x 10^{−2}.

sweep  Γ  γ = \frac{Γ_{2}}{Γ_{1}}  \frac{S(A)}{δ}
0  7.369 · 10^{−2} 3.9587 5.292 · 10^{−4}
1  1.256 · 10^{−4} 0.6306 1.578 · 10^{−8}
2 1.581 · 10^{−13} 0.6773 2.498 · 10^{−26}
3 3.013 · 10^{−40} 0.7542 9.076 · 10^{−80}

We see that in the case Γ = Γ_{2} we cannot guarantee cubic reduction of S(A)/δ per quasi-sweep.

However, in all our computations asymptotically Γ = Γ_{1}. This is a very typical phenomenon. To give a heuristic argument why this is so we can make a qualitative analysis using O-symbols. If ε = Γ_{1} and η = Γ_{2} we can easily obtain that Γ_{1}^{(M)} = O(ε^{3}) + O(η^{3}) and Γ_{2}^{(M)} = O(η) · O(ε^{2}). This also shows that if Γ < 0.5, then after one qusi-sweep γ will be of moderate size and cubic convergence of S(A)/δ will be indicated. A similar conclusion can be drawn from the relation (3.3.12). Let I_{0} = I = [Γ^{3} , \frac{\sqrt{5}}{2}Γ^{2}] ,I_{1} = I^{3}…, where I^{3} = \{x^{3} : x ∈ I\}. Although Γ reduces cubically per quasi-sweep if S(A)/δ is close to the left edge of I, it can occur that S(A^{(M)})/δ comes close to the right edge of I_{1} (quadratic reduction). But after that, since the right edge of I_{1} converges cubically to zero, the same must be true with S(A^{(kM)})/δ, k ≥ 1. Quadratic reduction of S(A^{(kM)})/δ for some k ≥ 1 would first require its better than cubic reduction in an earlier quasi-sweep, but this is in general not possible.

The following example demonstrates the above conclusions. We shall assume larger multiplicities of eigenvalues.

Example 3.10 We have generated a 40 x 40 matrix A from a diagonal matrix D whose diagonal entries are

d_{ii} = \{ \begin{matrix}−0.3, 1 ≤ i ≤ 10, \\ 0.0, 11 ≤ i ≤ 20, \end{matrix} d_{ii} = \{ \begin{matrix}0.3, 21 ≤ i ≤ 30, \\ 0.6, 31 ≤ i ≤ 40, \end{matrix}

by applying 5 sweeps (taken in the row-wise fashion) of rotational similarity transformations. In each rotation sin ϕ value was chosen as (2 random−1)/500, where random denotes a random number between 0 and 1. The cos ϕ has been taken nonnegative. Note that each eigenvalue of A has multiplicity 10 and δ is approximately 0.1. Γ has been computed: Γ = 0.6860. Hence the initial matrix does not satisfy the assumption A1. But it has been checked that A^{(M)} satisfies all the assumptions A1, A2, A3 with l_{1} = l_{2} = m_{1} = m_{2} = 10.

sweep  Γ  γ = \frac{Γ_{2}}{Γ_{1}}  \frac{S(A)}{δ}
0  6.860 · 10^{−1} 2.0371 3.421 · 10^{−1}
1  5.200 · 10^{−2} 0.6767 2.704 · 10^{−3}
2  9.854 · 10^{−6} 0.6803 9.710 · 10^{−11}
3 9.180 · 10^{−17} 0.4748 8.428 · 10^{−33}

We conclude this paper with an example which shows that the assumption A3 cannot be dropped.

Example 3.11 Here we use the same matrix A as in the Example 3.10, but we choose l_{1} = 5, l_{2} = 10, m_{1} = 10, and m_{2} = 15. It has been checked that A^{(M)} satisfies the assumptions A1 and A2 but violates the assumption A3.

sweep  Γ  γ = \frac{Γ_{2}}{Γ_{1}}  \frac{S(A)}{δ}
0  6.743 · 10^{−1} 1.7303 3.421 · 10^{−1}
1  5.255 · 10^{−2} 0.6435 2.761 · 10^{−3}
2  9.525 · 10^{−6} 0.6573 9.072 · 10^{−11}
3 9.921 · 10^{−10} 1.7456 9.844 · 10^{−19}

Note that Γ^{(3M)} is by no means cubic reduction of Γ^{(2M)}. Also note that S(A^{(3M)})/δ is not a cubic reduction of S(A^{(2M)})/δ even though Γ^{(2M)} = Γ_{1}^{(2M)}. Hence Theorem 3.7 and Corollary 3.8 do not hold if the assumption A3 is omitted.

We end this paper with a remark concerning the case of very close eigenvalues. For that case the assumption A1 means that S(A) is very small and any realistic stopping criterion is likely to terminate the procedure before the cubic reduction comes about. Nevertheless, if the eigenvalues are grouped in clusters, and the distance between every two clusters is sufficiently large relatively to δ, we can use a technique from [7] to prove the cubic reduction of Γ in an earlier quasi-sweep (when S(A) ≫ δ). Once S(A) shrinks to the order of magnitude of δ there may be a slowdown in the rate of convergence. Ultimately, the assumptions A1,A2 and A3 will hold and the cubic convergence will prevail.

Appendix

Proof of Lemma 3.3 Suppose that f = g. That is, A_{22} = A_{ff}. Then by Lemma 3.1 we have

S^{2}(A_{22}^{(ν(T))}) ≤ (0.26)^{2} \frac{[∥A_{12}^{(ν(T))}∥^{2} + ∥A_{23}^{(ν(T))}∥^{2}]^{2}}{δ^{2}} = (0.26)^{2} \frac{[∥A_{12}∥^{2} + ∥A_{23}∥^{2}]^{2}}{δ^{2}} .

The last inequality follows from the fact that ∥A_{12}^{(ν(T))}∥ = ∥A_{12}∥ and ∥A_{23}^{(ν(T))}∥ = ∥A_{23}∥. Hence the assertion of the lemma follows.

Now assume that f < g. Let n̂_{i} be the integer such that A_{22}^{(n̂_{i})} is obtained from A_{22} after rotating the element at position (n̄_{i}, n̄_{g}), f ≤ i ≤ g − 1, and we let n̂_{f−1} = 0. In Figure 3 we have displayed the positions (n̄_{i}, n̄_{g}), f ≤ i ≤ g − 1, in the case f = 3 and g = 5.

[figure: A_{33}; A_{34}; A_{35}; A_{44}; A_{45}; A_{55}; A_{22} =; x; x; a_{n̄_{3},n̄_{5}}; a_{n̄_{4},n̄_{5}}; (this element is rotated; at n̂_{3}th rotation); (this element is rotated; at n̂_{4}th rotation)]

Figure 4: The case of f = 3 and g = 5

Following the proof of [7, Theorem 2.8] and using [7, Lemma 2.6 (ii)] we obtain

S^{2}(A_{22}^{(ν(T))}) = S^{2}(A_{gg}^{(ν(T))}) + ∑_{i=f}^{g−1} [S^{2}(A_{ii}^{(ν(T))}) + ∑_{j=i+1}^{g} ∥A_{ij}^{(ν(T))}∥^{2}]

≤ (1 + γ) ∑_{i=f}^{g} S^{2}(A_{ii}^{(n̂_{i−1})}) + (1 + \frac{1}{γ})S^{2}(A_{22}) ∑_{(p,q)∈G} sin^{2} ϕ_{pq},

where γ > 0 is arbitrary. Here

G = G_{f} ∪ G_{f+1} ∪ … ∪ G_{g−1}

where

G_{k} = {(r, s) : n̄_{k−1} + 1 ≤ r ≤ n̄_{k}, n̄_{k} + 1 ≤ s ≤ n̄_{g}}, f ≤ k ≤ g − 1.

We have denoted the angle which rotates the element in the position (p, q) by ϕ_{pq}. Using Lemma 3.2 we have

S^{2}(A_{22}^{(ν(T))}) ≤ (1 + γ) ∑_{i=f}^{g} S^{2}(A_{ii}^{(n̂_{i−1})}) + 0.1156(1 + \frac{1}{γ})\frac{S^{4}(A_{22})}{δ^{2}}.

If we prove

∑_{i=f}^{g} S^{2}(A_{ii}^{(n̂_{i−1})}) ≤ 0.0811\frac{[S^{2}(A_{22}) + ∥A_{12}∥^{2} + ∥A_{23}∥^{2}]^{2}}{δ^{2}} , (A.1)

then the assertion will follow by specifying the minimizing γ (cf. (3.2.1)) and by adding the coeficients of \frac{S^{2}(A_{22})}{δ}.

Let us define the row-vectors (cf. [7])

x_{t}^{(k)} = [a_{t1}^{(k)}, …, a_{t,t−1}^{(k)}, 0 , a_{t,t+1}^{(k)}, …, a_{tn}^{(k)}] , n̄_{f−1} + 1 ≤ t ≤ n̄_{g}, 0 ≤ k ≤ ν(T),

and the number Y^{2} = S^{2}(A_{22}) + ∥A_{12}∥^{2} + ∥A_{23}∥^{2}. (A.2)

By Lemma 3.1 and by the invariance of Y under considered Jacobi rotations we have for f ≤ i ≤ g

S^{2}(A_{ii}^{(n̂_{i−1})}) ≤ (0.26)^{2} \frac{Y^{2}}{δ^{2}} [ ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}^{(n̂_{i−1})}∥^{2} − 2S^{2}(A_{ii}^{(n̂_{i−1})})] .

So we obtain

S^{2}(A_{ii}^{(n̂_{i−1})}) ≤ \frac{(0.26)^{2}Y^{2}/δ^{2}}{1 + 2(0.26)^{2}Y^{2}/δ^{2}} ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}^{(n̂_{i−1})}∥^{2}. (A.3)

In [7, Lemma 2.7 (see the relations (2.24) and (2.21)) ] it has been proved that for any β > 0, n̄_{i−1} + 1 ≤ t ≤ n̄_{i}, f ≤ i ≤ g, holds

(A.4) ∥x_{t}^{(n̂_{i−1})}(In̄_{i−1} ⊕ 0_{n_{i}} ⊕ I_{n−n̄_{i}})∥^{2} ≤ (1 + β)∥x_{t}(In̄_{f−1} ⊕ 0_{n̄_{i}−n̄_{f−1}} ⊕ I_{n−n̄_{i}})∥^{2}

+(1 + \frac{1}{β})Y^{2} ∑_{r=n̄_{f−1}+1}^{n̄_{i−1}} sin^{2} ϕ_{rt}.

By the relation (3.1.7)) we have

A_{ii}^{(n̂_{i−1})} − λ_{i}I_{n_{i}} = G_{i}^{(n̂_{i−1})}K_{i}^{(n̂_{i−1})}^{−1}G_{i}^{(n̂_{i−1})}^{T} .

Using the last equality together with the relations (3.1.7) and (3.1.4) we obtain

∥x_{t}^{(n̂_{i−1})}(0n̄_{i−1} ⊕ I_{n_{i}} ⊕ 0_{n−n̄_{i}})∥ ≤ ∥e_{t−n̄_{i−1}}^{T}(A_{ii}^{(n̂_{i−1})} − λ_{i}I_{n_{i}})∥

≤ \frac{8}{24 − \sqrt{5}} \frac{∥G_{i}^{(n̂_{i−1})}∥}{δ}∥x_{t}^{(n̂_{i−1})}(In̄_{i−1} ⊕ 0_{n_{i}} ⊕ I_{n−n̄_{i}})∥.

Here e_{r}, 1 ≤ r ≤ n_{i} are the columns of I_{n_{i}}. Since

∥G_{i}^{(n̂_{i−1})}∥ ≤ S(A^{(n̂_{i−1})}) ≤ S(A) ≤ (\sqrt{5}/8)δ,

we have ∥x_{t}^{(n̂_{i−1})}(0n̄_{i−1} ⊕ I_{n_{i}} ⊕ 0_{n−n̄_{i}})∥ ≤ 0.1028∥x_{t}^{(n̂_{i−1})}(In̄_{i−1} ⊕ 0_{n_{i}} ⊕ I_{n−n̄_{i}})∥. (A.5) From the relations (A.4) and (A.5) we obtain

∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}^{(n̂_{i−1})}∥^{2} ≤ 1.0106 [(1 + β) ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}(In̄_{f−1} ⊕ 0_{n̄_{i}−n̄_{f−1}} ⊕ I_{n−n̄_{i}})∥^{2}

+(1 + \frac{1}{β})Y^{2} ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∑_{r=n̄_{f−1}+1}^{n̄_{i−1}} sin^{2} ϕ_{rt}] , f ≤ i ≤ g.

Here an empty sum is considered zero. Using

∑_{i=f}^{g} ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}(In̄_{f−1} ⊕ 0_{n̄_{i}−n̄_{f−1}} ⊕ I_{n−n̄_{i}})∥^{2} ≤ Y^{2}

and Lemma 3.2 we obtain

∑_{i=f}^{g} ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}^{(n̂_{i−1})}∥^{2} ≤ 1.0106Y^{2} [(1 + β) + (1 + \frac{1}{β}) ∑_{i=f+1}^{g} ∑_{t=n̄_{i−1}+1}^{n̄_{i}} ∑_{r=n̄_{f−1}+1}^{n̄_{i−1}} sin^{2} ϕ_{rt}]
≤ 1.0106Y^{2} [(1 + β) + (1 + \frac{1}{β})(0.1156)\frac{S^{2}(A_{22})}{δ^{2}}] .

The last relation together with (A.6) implies

∑_{i=f}^{g} S^{2}(A_{ii}^{(n̂_{i−1})}) ≤ \frac{(0.26)^{2}Y^{2}/δ^{2}}{1 + 2(0.26)^{2}Y^{2}/δ^{2}} ∑_{i=f}^{g} ∑_{n̄_{i−1}+1}^{n̄_{i}} ∥x_{t}^{(n̂_{i−1})}∥^{2}

≤ 0.06832\frac{Y^{4}/δ^{2}}{1 + 2(0.26)^{2}Y^{2}/δ^{2}} [(1 + β) + (1 + \frac{1}{β})(0.1156)\frac{S^{2}(A_{22})}{δ^{2}}] .

Using the minimizing β = 0.34, the inequality Y/δ ≤ S(A)/δ and the relation (3.1.4) we have

∑_{i=f}^{g} S^{2}(A_{ii}^{(n̂_{i−1})}) ≤ 0.06832\frac{Y^{4}}{δ^{2}} \frac{[1 + 0.34S(A_{22})/δ]^{2}}{1 + 2[0.26Y/δ]^{2}} ≤ 0.06832\frac{Y^{4}}{δ^{2}} \frac{[1 + 0.34Y/δ]^{2}}{1 + 2[0.26Y/δ]^{2}}
≤ 0.06832\frac{Y^{4}}{δ^{2}} sup_{0≤x≤\sqrt{5}/8} \frac{(1 + 0.34x)^{2}}{1 + 2(0.26x)^{2}} ≤ 0.0811\frac{Y^{4}}{δ^{2}} .

As we can see from the definition (A.2) of Y, the last inequality is just (A.1) which completes the proof of Lemma 3.3. .

Acknowledgement The authors are indebted to W. Mascarenhas for his interest and support during their research. The authors also want to thank the anonymous referees of this journal for helpful comments.

References

[1] J. Demmel and K. Veselić: Jacobi’s Method is More Accurate than QR, Courant Institute of Mathematical Sciences, Department of Computer Science, New York University 251 Mercer Street, New York, Y. Y. 10012, October, 1989.

[2] G. E. Forsythe and P. Henrici: The Cyclic Jacobi Method for Computing the Principal Values of a Complex Matrix, Transactions of the American Mathematical Society 94,1–23, 1960.

[3] G. H. Golub and C. F. Van Loan: Matrix Computations, The Johns Hopkins University Press, 2nd edition, 704 West 40th Street, Baltimore, Maryland 21211, 1989.

[4] E. R. Hansen: On Jacobi Methods and Block Jacobi Methods for Computing Matrix Eigenvalues, Ph. D. thesis, Stanford University 1960.

[5] E. R. Hansen: On Cyclic Jacobi Methods, SIAM J. Appl. Math. 11, 449–459, 1963.

[6] V. Hari: On Pairs of Almost Diagonal Matrices, Linear Algebra and Its Applications 148, 193–223, 1991.

[7] V. Hari: On Sharp Quadratic Convergence Bounds for the Serial Jacobi Methods, Numerische Mathematik 60, 375–406, 1991.

[8] P. Henrici and K. Zimmermann: An Estimate for the Norms of Certain Cyclic Jacobi Operators, Linear Algebra and Its Applications 1, 489–501, 1968.

[9] H. P. M. Van Kempen: On the Quadratic Convergence of the Special Cyclic Jacobi Method, Numerische Mathematik 9, 19–22, 1966

[10] W. F. Mascarenhas: On the Convergtence of the Jacobi Method, Poster Presentation, Fourth SIAM Conference on Parallel Processing for Scientific Computing, Chicago, Illinois, December, 1989.

[11] W. F. Mascarenhas: On the Convergence of the Jacobi Method for Arbitrary Orderings I, SIAM Journal of Scientific and Statistical Computing, submitted, 1990.

[12] L. Nazareth: On the Convergence of the Cyclic Jacobi Method, Linear Algebra and Its Applications 12, 151–164, 1975.

[13] G. Shroff and R. Schreiber: On the Convergence of the Cyclic Jacobi Method for Parallel Block Orderings, SIAM Journal on Matrix Analysis and Applications 10, 326–346, 1989.

[14] K. Veselić and V. Hari: A Note on the One-sided Jacobi Algorithm, Numerische Mathematik 56, 627–633, 1989.

[15] J. H. Wilkinson: Note on the Quadratic Convergence of the Cyclic Jacobi Process, Numerische Mathematik 4, 296–300, 1962.

[16] J. H. Wilkinson: Almost Diagonal Matrices with Multiple or Close Eigenvalues, Linear Algebra and its Applications 1, 1–12, 1968.