HOMOTOPY UNDERSTANDING OF ITERATIVE METHODS

Noah H. Rhee

University of Missouri - Kansas City

Abstract. In this paper we present a homotopy understanding of the iterative methods including the Jacobi, Gauss–Seidel, and Successive Overrelaxation methods. The purpose of this paper is to give a new insight to those iterative methods. Based upon this new understanding we present a condition for the convergence of those iterative methods, which turns out to be equivalent to the traditional one.

1. Introduction. Suppose we want to solve a linear system

(1.1)  Ax = b

by an iterative method. We assume that A = (a_{ij}) is a nonsingular m–by–m matrix with complex number entries, and that b = (b_{i}) ∈ ℂ^{m}, that is, b is a vector with m complex number components. We will denote the unique solution of (1.1) by

(1.2)  x̂ = A^{−1}b.

Usually such a linear system is solved by a direct method such as finding LU or QR decomposition of A. But direct methods can be impractical if A is large and sparse, because the sought–after factors can be dense. So in such a case we need to study about iterative methods that can solve (1.1).

The Jacobi, Gauss–Seidel, and Successive OverRelaxation (SOR) iterative methods [3,4,5] are typical members of a large family of iterative methods that have the form

(1.3)  Mx_{n+1} = Nx_{n} + b, n = −1, 0, 1, …

where

(1.4)  A = M − N

is a splitting of A, M is nonsingular and x_{−1} is any vector in ℂ^{m}. Note that we let n start from −1 (Normally n starts from 0. The reason of this shift of index will become apparent in the next section.)

For the iterative method (1.3) to be practical, it must be easy to solve a linear system with M as the matrix. Henceforth, we will assume that such is the case. We note that for the Jacobi method, M is diagonal, and for the Gauss–Seidel and SOR methods, M is lower triangular.

Whether or not the iterative method (1.3) produces the iterates {x_{n}}_{n=0}^{∞} that converge to x̂ = A^{−1}b, which is the solution of (1.1), depends on the eigenvalues of M^{−1}N. Of course, the convergence of vectors in ℂ^{m} is understood as componentwise convergence. For any G ∈ ℂ^{m×m}, the set of all m–by–m matrices whose elements are complex numbers, we define ρ(G) to mean the spectral radius of G, that is, the maximum eigenvalue in magnitude. Then we have the following theorem.

\underline{Theorem 1.1}. If ρ(M^{−1}N) < 1, then the iterates \{x_{n}\}_{n=0}^{∞} defined by (1.3), converge to x̂, which is the solution of (1.1), for any starting vector x_{−1} ∈ ℂ^{m}.

\underline{Proof}. See [4, p. 59].

In this paper we present a homotopy continuation understanding of the iterative method (1.3). And based upon this new understanding we present a condition for the convergence of the iterative method (1.3), which turns out to be equivalent to the condition ρ(M^{−1}N) < 1 which is specified in Theorem 1.1.

2. Homotopy Continuation Approach. Homotopy continuation means that we join a trivial system Mx = b′ to the given linear system Ax = b using parameter t ∈ [0, 1]. More precisely, we define the homotopy H: ℂ^{m} × [0, 1] → ℂ^{m} by

(2.1)  H(x, t) = A(t)x − b(t)

where

(2.2)  A(t) = M − tN

(2.3)  b(t) = (1 − t)Nx_{−1} + b.

If we assume that A(t) is nonsingular for each t ∈ [0, 1], H(x, t) = 0 has a unique solution which we will denote by x̃(t), for each t ∈ [0, 1]. This means that x̃ can be viewed as a complex vector–valued function from [0, 1] to ℂ^{m}. Note that x̃(0) satisfies the equation A(0)x̃(0) − b(0) = 0. Using the relations (2.2) and (2.3) we obtain

(2.4)  Mx̃(0) = Nx_{−1} + b.

Note that x̃(0) is immediately available, since a linear system with the matrix M is easy to solve. We also note that x̃(1) satisfies the equation A(1)x̃(1) − b(1) = 0. Using the relations (2.2), (1.4) and (2.3) we obtain

(2.5)  Ax̃(1) = b.

That is, x̃(1) = x̂, the solution of the given problem.

Now let us study about the complex vector–valued function x̃: [0, 1] → ℂ^{m} further. First of all, we need some basic definitions.

\underline{Definition 2.1}. Let v: [0, 1] → ℂ^{m} be a complex vector–valued function of t. We write

v(t) = [v_{1}(t), …, v_{m}(t)]^{T}

Here ‘T ’ stands for ‘transpose’. If each component v_{i}(t) is differentiable for all t ∈ [0, 1], i = 1, …, m, then we define

\frac{d}{dt}v(t) = [\frac{d}{dt}v_{1}(t), …, \frac{d}{dt}v_{m}(t)]^{T} , t ∈ [0, 1].

And we say that v is differentiable on [0, 1].

If each component v_{i}(t) is infinitely many times differentiable on [0, 1], i = 1, …, m, then we define for n = 1, 2, …

v^{(n)}(t) = \frac{d^{n}}{dt^{n}} v(t) = [\frac{d^{n}}{dt^{n}} v_{1}(t), …, \frac{d^{n}}{dt^{n}} v_{m}(t)]^{T} , t ∈ [0, 1].

And we say that v is infinitely many times differentiable on [0, 1]. By convention v^{(0)}(t) is understood as v(t).

\underline{Definition 2.2}. If v: [0, 1] → ℂ^{m} is infinitely many times differentiable on [0, 1], then for n = 0, 1, …

P_{n}(v, t) = v(0) + tv^{(1)}(0) + ⋯ + \frac{t^{n}}{n!} v^{(n)}(0), t ∈ [0, 1]

is called the nth–degree complex vector–valued Taylor polynomial of v(t) and

P_{∞}(v, t) = v(0) + tv^{(1)}(0) + ⋯ + \frac{t^{n}}{n!} v^{(n)}(0) + ⋯, t ∈ [0, 1]

is called the complex vector–valued Taylor series of v(t).

\underline{Definition 2.3}. Let B(t) = (b_{ij}(t)) ∈ ℂ^{m×m} and each entry b_{ij}(t) is a differentiable function of t, then we define

\frac{d}{dt}B(t) = ( \begin{matrix}\frac{d}{dt}b_{11}(t) & ⋯ & \frac{d}{dt}b_{1m}(t) \\ ⋮ & ... & ⋮ \\ \frac{d}{dt}b_{m1}(t) & ⋯ & \frac{d}{dt}b_{mm}(t) \end{matrix} ) .

If each entry b_{ij}(t) is infinitely many times differentiable, then for n = 1, 2, … we define

\frac{d^{n}}{dt^{n}} B(t) = ( \begin{matrix}\frac{d^{n}}{dt^{n}} b_{11}(t) & ⋯ & \frac{d^{n}}{dt^{n}} b_{1m}(t) \\ ⋮ & ... & ⋮ \\ \frac{d^{n}}{dt^{n}} b_{m1}(t) & ⋯ & \frac{d^{n}}{dt^{n}} b_{mm}(t) \end{matrix} ) .

We also recall the following basic definitions from linear algebra.

\underline{Definition 2.4}. For any G = (g_{ij}) ∈ ℂ^{m×m}, we let M_{ij}(G) be the (m − 1)–by–(m − 1) submatrix of G obtained by deleting row i and column j of G. M_{ij}(G) is called the minor of g_{ij}. And

cof_{ij}(G) = (−1)^{i+j} det(M_{ij}(G))

is called the cofactor of g_{ij}. Here det(M_{ij}(G)) denotes the determinant of the matrix M_{ij}(G).

\underline{Definition 2.5}. For G = (g_{ij}) ∈ ℂ^{m×m}, the adjoint of G, denoted by Adj(G) = (g_{ij}^{∗}), is defined to be the transpose of the matrix obtained from G by replacing each element by its cofactor. That is, g_{ij}^{∗} = cof_{ji}(G).

Now we state the following basic theorem from linear algebra.

\underline{Theorem 2.6}. If G = (g_{ij}) ∈ ℂ^{m×m} is nonsingular, then det(G) ≠ 0 and

G^{−1} = \frac{1}{det(G)}Adj(G).

\underline{Proof}. See [2, p. 119].

Now we can state the following lemma.

\underline{Lemma 2.7}. If the matrix A(t), defined in the relation (2.2), is nonsingular for each t ∈ [0, 1], then the complex vector–valued function x̃: [0, 1] → ℂ^{m}, defined as the unique solution of H(x, t) = 0 (see the relation (2.1) for the definition of H(x, t)) for each t ∈ [0, 1], is infinitely many times differentiable on [0, 1].

\underline{Proof}. Let t be any point in [0, 1]. From the definition of the function x̃ we see that

(2.6)  x̃(t) = A(t)^{−1}b(t).

By Theorem 2.6 we have

(2.7)  A(t)^{−1} = \frac{1}{det(A(t))}Adj(A(t)).

Since each entry of A(t) is linear in the variable t (see (2.2)), from the relation (2.7) and the definitions of adjoint and determinant we see that each entry of A(t)^{−1} is a rational function of t. Furthermore, det(A(t)), which is the common denominator of each entry of A(t)^{−1}, is not zero. Hence, each entry of A(t)^{−1} is infinitely many times differentiable.

From the relation (2.3) we see that b(t) is infinitely many times differentiable with respect to the variable t. In fact, we have

(2.8)  \frac{d}{dt}b(t) = −Nx_{−1} and \frac{d^{k}}{dt^{k}} b(t) = 0, k = 2, 3, ….

Hence, from the relation (2.6) we see that x̃(t) is differentiable. In fact, using the relation (2.8) from the relation (2.6) we obtain

(2.9)  \frac{d}{dt}x̃(t) = \frac{d}{dt}[A(t)^{−1}]b(t) − A(t)^{−1}Nx_{−1}.

This differentiation procedure can be repeated as many times as we want, since each entry of A(t)^{−1} and b(t) are infinitely many times differentiable. So by applying Leibnitz’s rule of differentiation to both sides of the relation (2.6) and using the relation (2.8) we obtain

(2.10)  \frac{d^{n}}{dt^{n}} x̃(t) = \frac{d^{n}}{dt^{n}} [A(t)^{−1}]b(t) − n\frac{d^{n−1}}{dt^{n−1}} [A(t)^{−1}]Nx_{−1}, n = 1, 2, ….

Since t was arbitrary, we proved that x̃ is infinitely many times differentiable on [0, 1]. Since the function x̃: [0, 1] → ℂ^{m} is infinitely many times differentiable, to approximate

x̃(t) we may use the nth–degree complex vector–valued Taylor polynomial P_{n}(x̃, t) defined in Definition 2.2 to get for any t ∈ [0, 1]

\begin{aligned}x̃(t) & ≈ P_{n}(x̃, t) \\ & = x̃(0) + tx̃^{(1)}(0) + ⋯ + \frac{t^{n}}{n!} x̃^{(n)}(0), n = 0, 1, …. \end{aligned}

So x̃(1) = x̂, the solution of the given problem (1.1), may be approximated

x̂ = x̃(1) ≈ P_{n}(x̃, 1)

= x̃(0) + x̃^{(1)}(0) + ⋯ + \frac{1}{n!}x̃^{(n)}(0), n = 0, 1, …. (2.11)

To manipulate the expression in the relation (2.11) we need the following lemma. \underline{Lemma 2.8}. If A(t) is nonsingular for all t ∈ [0, 1], then for any t ∈ [0, 1] we have

x̃^{(1)}(t) = A(t)^{−1}N(x̃(t) − x_{−1})

and

x̃^{(n)}(t) = nA(t)^{−1}Nx̃^{(n−1)}(t), n = 2, 3, …

\underline{Proof}. From the definition of x̃(t) we see

(2.12)  A(t)x̃(t) = b(t).

From Lemma 2.7 and the relation (2.2) we know that x̃ and each entry of A(t) are differentiable with respect to the variable t. Hence, differentiating both sides of the relation (2.12) with respect to t we obtain

(2.13)  [\frac{d}{dt}A(t)]x̃(t) + A(t)\frac{d}{dt}x̃(t) = \frac{d}{dt}b(t).

From the relation (2.2) we observe that

(2.14)  \frac{d}{dt}A(t) = −N.

Using the relations (2.8) and (2.14) the relation (2.13) becomes

(2.15)  −Nx̃(t) + A(t)x̃^{(1)}(t) = −Nx_{−1}.

Hence, we obtain

x̃^{(1)}(t) = A(t)^{−1}N(x̃(t) − x_{−1}),

which proves the first statement of the lemma.

To prove the second statement of the lemma we differentiate the relation (2.15) with respect to the variable t to obtain

−Nx̃^{(1)}(t) + [\frac{d}{dt}A(t)]x̃^{(1)}(t) + A(t)x̃^{(2)}(t) = 0.

Using the relation (2.14) we obtain

x̃^{(2)}(t) = 2A(t)^{−1}Nx̃^{(1)}(t),

which proves the second statement of the lemma with n = 2. Suppose the second statement of the lemma is true with n = k, that is, we assume

x̃^{(k)}(t) = kA(t)^{−1}Nx̃^{(k−1)}(t).

From this equation we obtain

A(t)x̃^{(k)}(t) = kNx̃^{(k−1)}(t).

Differentiating both sides of this equation with respect to the variable t we obtain

[\frac{d}{dt}A(t)]x̃^{(k)}(t) + A(t)x̃^{(k+1)}(t) = kNx̃^{(k)}(t).

Using the relation (2.14) this equation becomes

x̃^{(k+1)}(t) = (k + 1)A(t)^{−1}Nx̃^{(k)}(t).

So the second statement of the lemma is valid with n = k + 1. Now the second statement of the lemma follows by mathematical induction.

Repeated use of Lemma 2.8 yields

x̃^{(k)}(t) = kA(t)^{−1}Nx̃^{(k−1)}(t)
= k(k − 1)[A(t)^{−1}N]^{2}x̃^{(k−2)}(t)
= ⋯
= k(k − 1) ⋯ 2[A(t)^{−1}N]^{k−1}x̃^{(1)}(t)
= k![A(t)^{−1}N]^{k−1}[A(t)^{−1}N(x̃(t) − x_{−1})]
= k![A(t)^{−1}N]^{k}(x̃(t) − x_{−1}).

So we obtain

\frac{1}{k!}x̃^{(k)}(t) = [A(t)^{−1}N]^{k}(x̃(t) − x_{−1}), k = 1, 2, ….

Setting t = 0 and using the fact that A(0) = M, the above equation becomes

(2.16)  \frac{1}{k!}x̃^{(k)}(0) = (M^{−1}N)^{k}(x̃(0) − x_{−1}), k = 1, 2, ….

Using the relation (2.16), the relation (2.11) becomes

x̂ = x̃(1) ≈ P_{n}(x̃, 1)

= x̃(0) + [M^{−1}N + ⋯ + (M^{−1}N)^{n}](x̃(0) − x_{−1}) n = 0, 1, …. (2.17)

Now we are ready to state the main theorem in this section.

\underline{Theorem 2.9}.

x_{n} = P_{n}(x̃, 1), n = 0, 1, ….

\underline{Proof}. From the relation (1.3) we have

(2.18)  x_{0} = M^{−1}(Nx_{−1} + b).

From the relation (2.4) we note that

(2.19)  x̃(0) = M^{−1}(Nx_{−1} + b).

Hence, it follows that x_{0} = P_{0}(x̃, 1). This proves the theorem with n = 0. Now suppose that the theorem is true with n = k, that is, we assume x_{k} = P_{k}(x̃, 1). We need to prove x_{k+1} = P_{k+1}(x̃, 1). From the relation (1.3) and using the induction hypothesis and the relations (2.17) and (2.19) we obtain

x_{k+1} = M^{−1}Nx_{k} + M^{−1}b
= M^{−1}NP_{k}(x̃, 1) + M^{−1}b
= M^{−1}N{x̃(0) + [M^{−1}N + ⋯ + (M^{−1}N)^{k}](x̃(0) − x_{−1})} + M^{−1}b
= M^{−1}(Nx̃(0) + b) + [(M^{−1}N)^{2} + ⋯ + (M^{−1}N)^{k+1}](x̃(0) − x_{−1})
= M^{−1}(Nx̃(0) − Nx_{−1} + Nx_{−1} + b)
+ [(M^{−1}N)^{2} + ⋯ + (M^{−1}N)^{k+1}](x̃(0) − x_{−1})
= M^{−1}(Nx_{−1} + b) + [(M^{−1}N) + ⋯ + (M^{−1}N)^{k+1}](x̃(0) − x_{−1})
= x̃(0) + [(M^{−1}N) + ⋯ + (M^{−1}N)^{k+1}](x̃(0) − x_{−1})
= P_{k+1}(x̃, 1).

So the theorem follows by mathematical induction.

Theorem 2.9 tells us that for each n = 0, 1, …, x_{n}, generated by the iterative method (1.3), can be viewed as the nth–degree complex vector–valued Taylor polynomial approximation of x̃(1) = x̂.

3. Convergence Analysis. To find a condition for the convergence of the iterative method (1.3) based upon our homotopy understanding, we need to extend the real homotopy parameter t to a complex parameter z. More precisely, we define the homotopy H: ℂ^{m} × ℂ → ℂ^{m} such that

(3.1)  H(x, z) = A(z)x − b(z)

where

(3.2)  A(z) = M − zN

and

(3.3)  b(z) = (1 − z)Nx_{−1} + b.

We note that the relations (3.1), (3.2) and (3.3) reduce to the relations (2.1), (2.2) and (2.3), respectively, when z = t ∈ [0, 1]. Now we state a lemma.

\underline{Lemma 3.1}. Suppose A(z) defined in relation (3.2) is nonsingular for all z such that |z| ≤ 1. Then there exists ϵ > 0 such that A(z) is nonsingular for all z such that |z| < 1 + ϵ.

\underline{Proof}. Let g(z) = det(A(z)). Note that g is a continuous function of z from ℂ to ℂ, since it is a polynomial of z. Let h(z) = |g(z)|. Then clearly h is a continuous function from ℂ to ℝ. Since the set F := {z : |z| ≤ 1} is compact, h attains its extreme values on F. Let e and E be its minimum and maximum value, respectively. Since A(z) is nonsingular on F, h(z) > 0 on F. This implies that e > 0. Since h is continuous, if we let U = h^{−1}(e/2, E +1), then U is an open set in ℂ containing the set F. Since F is compact, it follows that there exists ϵ > 0 such that V := {z : |z| < 1 + ϵ} ⊆ U. Hence, we have

h(V) ⊆ h(U) = h(h^{−1}(e/2, E + 1)) ⊆ (e/2, E + 1).

Since e/2 > 0, the above inclusion means that A(z) is nonsingular on V. So the lemma is proved.

If we assume that A(z) is nonsingular for all z such that |z| ≤ 1, then the homotopy equation H(x, z) = 0 (see (3.1)) has a unique solution x̃(z) for each z such that |z| ≤ 1. In fact, by Lemma 3.1 such a unique solution x̃(z) exists for any z such that |z| < 1 + ϵ for some ϵ > 0. So if we let

(3.4)  V = {z : |z| < 1 + ϵ},

then x̃ can be thought of as a complex vector–valued function from V to ℂ^{m}. Now we recall the following basic definition from complex variables [1]. \underline{Definition 3.2}. A complex–valued function f of the complex variable z is analytic at a

point z_{0} if its derivative exists not only at z_{0} but at each point z in some neighborhood of z_{0}. It is analytic in an open set U if it is analytic at every point in U.

Now we can state the following theorem.

\underline{Theorem 3.3}. If A(z) is nonsingular for all z such that |z| ≤ 1, then each component of the complex vector–valued function x̃ is analytic in V (see (3.4) for the definition of V).

\underline{Proof}. According to Definition 3.2 all we need to show is that each component of x̃ is differentiable in V; and the proof is entirely similar to that of Lemma 2.7. So we omit the proof.

We recall the following theorem of complex variables.

\underline{Theorem 3.4}. Let f be analytic inside of a circle C whose center at z_{0} and radius r_{0}. Then at each point z inside C

f(z) = f(z_{0}) + [\frac{d}{dz} f(z_{0})](z − z_{0}) + ⋯ + \frac{1}{n!}[\frac{d^{n}}{dz^{n}} f(z_{0})](z − z_{0})^{n} + ⋯,

that is, Taylor series converges to f(z) when |z − z_{0}| < r_{0}.

\underline{Proof}. See [1, p. 145].

Now we have the following theorem.

\underline{Theorem 3.5}. Suppose A(z) is nonsingular for all z such that |z| ≤ 1. Then

(3.5)  x̃(z) = x̃(0) + z \frac{d}{dz} x̃(0) + ⋯ + \frac{z^{n}}{n!} \frac{d^{n}}{dz^{n}} x̃(0) + ⋯

is valid for all z in V.

\underline{Proof}. The theorem follows from Theorem 3.3 and applying Theorem 3.4 with z_{0} = 0 and r_{0} = 1 + ϵ to the complex vector–valued function x̃ componentwise. Of course, here the notation

\frac{d^{k}}{dz^{k}} x̃^{(k)}(z), k = 1, 2, …

is understood as componentwise kth derivative of the complex vector–valued function x̃ with respect to the variable z.

Now we state the main result of this section.

\underline{Corollary 3.6}. Suppose A(z) is nonsingular for each z such that |z| ≤ 1. Then the iterates {x_{n}}_{n=0}^{∞}, generated by the iterative method (1.3), converges to x̃(1) = x̂, which is the solution of the given problem (1.1), starting from any vector x_{−1} ∈ ℂ^{m}.

\underline{Proof}. Since 1 ∈ V, using Theorem 3.5 we obtain

x̂ = x̃(1) = x̃(0) + \frac{d}{dz} x̃(0) + ⋯ + \frac{1}{n!} \frac{d^{n}}{dz^{n}} x̃(0) + ⋯

= x̃(0) + \frac{d}{dt}x̃(0) + ⋯ + \frac{1}{n!} \frac{d^{n}}{dt^{n}} x̃(0) + ⋯

= x̃(0) + x̃^{(1)}(0) + ⋯ + \frac{1}{n!}x̃^{(n)}(0) + ⋯

= P_{∞}(x̃; 1)

= lim_{n→∞} P_{n}(x̃; 1)

= lim_{n→∞} x_{n}.

To get the last equation we used Theorem 2.9. Hence, the corollary is proved.

In Corollary 3.6 we have established a condition, namely, A(z) = M − zN (see (3.2)) is nonsingular for all z such that |z| ≤ 1, under which the iterates {x_{n}}_{n=0}^{∞}, generated by the iterative method (1.3), converges to x̂ = A^{−1}b. Now we show that this condition is equivalent to the condition specified in Theorem 1.1, namely, the condition that ρ(M^{−1}N) < 1.

\underline{Theorem 3.7}. The following two statements are equivalent.

(1) ρ(M^{−1}N) < 1.

(2) A(z) = M − zN is nonsingular for all z such that |z| ≤ 1.

\underline{Proof}. Suppose (1) is true. Suppose that M − zN is singular for some z with |z| ≤ 1. Then there exists a nonzero vector w such that

(3.6)  (M − zN)w = 0.

Note that z ≠ 0, since M is nonsingular. So from the relation (3.6) we obtain

(M^{−1}N)w = \frac{1}{z}w.

This shows that (1/z, w) is an eigenpair of M^{−1}N. Since |1/z| ≥ 1, we have ρ(M^{−1}N) ≥ 1, which is a contradiction. So we conclude that A(z) = M − zN is nonsingular for all |z| ≤ 1.

Conversely, suppose (2) is true. Suppose that ρ(M^{−1}N) ≥ 1. Then there exists an eigenpair (λ, w) of M^{−1}N such that

(M^{−1}N)w = λw and |λ| ≥ 1.

So we have

(3.7)  Mw = \frac{1}{λ}Nw.

If we let z = 1/λ, then |z| ≤ 1. And from the relation (3.7) we have

(M − zN)w = 0 and w ≠ 0.

This shows that M − zN is singular for some |z| ≤ 1, which is a contradiction. Thus, we conclude that ρ(M^{−1}N) < 1.

\underline{References}

1. R. V. Churchill, J. W. Brown, and R. F. Verhey, Complex Variable and Applications, (3rd ed.), McGraw–Hill, Inc., New York, 1974.

2. C. G. Cullen, Matrices and Linear Transformations, (2nd ed.), Addison–Wesley Publishing Co., Inc. Massachusetts, 1972.

3. L. A. Hageman and D. M. Young, Applied Iterative Methods, Academic Press, New York, 1981.

4. R. S. Varga, Matrix Iterative Analysis, Prentice–Hall, Englewood Cliffs, N.J., 1962.

5. D. M. Young, Iterative Solution of Large Linear Systems, Academic Press, New York, 1971.