����������� ������ 57 (2005), 27–33

UDK 519.633.6

���������� ������ ��� research paper

NUMERICAL STABILITY OF A CLASS (OF SYSTEMS) OF NONLINEAR EQUATIONS

Zlatko Udovičić

Abstract. In this article we consider stability of nonlinear equations which have the following form:

Ax + F (x) = b,  (1)

where F is any function, A is a linear operator, b is given and x is an unknown vector. We give (under some assumptions about function F and operator A) a generalization of inequality:

\frac{∥X_{1} − X_{2}∥}{∥X_{1}∥} ≤ ∥A∥ ∥A^{−1}∥ \frac{∥b_{1} − b_{2}∥}{∥b_{1}∥}  (2)

(equation (2) estimates the relative error of the solution when the linear equation Ax = b_{1} becomes the equation Ax = b_{2}) and a generalization of inequality:

\frac{∥X_{1} − X_{2}∥}{∥X_{1}∥} ≤ ∥A_{1}^{−1}∥ ∥A_{1}∥ ( \frac{∥b_{1} − b_{2}∥}{∥b_{1}∥} + ∥A_{1}∥ ∥A_{2}^{−1}∥ \frac{∥b_{2}∥}{∥b_{1}∥} · \frac{∥A_{1} − A_{2}∥}{∥A_{1}∥} )  (3)

(equation (3) estimates the relative error of the solution when the linear equation A_{1}x = b_{1} becomes the equation A_{2}x = b_{2}).

1. Basic results

Teorem 1. Let V be a normed space, let the linear operator A: V → V be invertible and bounded, let the inverse operator of the operator A be also bounded, let b_{1}, b_{2} ∈ V and let the functions F_{1}, F_{2} : V → V and the set S ⊆ V have the following properties:

1. the function F_{1} is Lipschitz on S, i.e.,

(∃L > 0) (∀x_{1}, x_{2} ∈ S) ∥F_{1} (x_{1}) − F_{1} (x_{2})∥ ≤ L ∥x_{1} − x_{2}∥,

and the constant L is such that the inequality 1 − L ∥A^{−1}∥ > 0, holds;

2. (∃M > 0) (∀x ∈ S) ∥F_{1} (x)∥ ≤ M ∥x∥; and

3. (∃ε ≥ 0) (∀x ∈ S) ∥F_{1} (x) − F_{2} (x)∥ ≤ ε.

AMS Subject Classification: 65J15

Keywords and phrases: Numerical stability, nonlinear equations.

If X_{1} ∈ S is a solution of the equation Ax + F_{1}(x) = b_{1} and X_{2} ∈ S is a solution of the equation Ax + F_{2}(x) = b_{2}, then the following inequality holds:

\frac{∥X_{1} − X_{2}∥}{∥X_{1}∥} ≤ \frac{∥A^{−1}∥ (∥A∥ + M)}{1 − L ∥A^{−1}∥} (\frac{∥b_{1} − b_{2}∥}{∥b_{1}∥} + \frac{ε}{∥b_{1}∥})

Proof. Since AX_{1} + F_{1} (X_{1}) = b_{1}, we have
∥b_{1}∥ = ∥AX_{1} + F_{1} (X_{1})∥ ≤ ∥AX_{1}∥ + ∥F_{1} (X_{1})∥ ≤ (∥A∥ + M) ∥X_{1}∥

and we can conclude that

\frac{1}{∥X_{1}∥} ≤ \frac{(∥A∥ + M)}{∥b_{1}∥}.  (4)

On the other hand, from X_{1}−X_{2} = A^{−1} ((b_{1} − b_{2}) − (F_{1} (X_{1}) − F_{2} (X_{2}))) it follows that

∥X_{1} − X_{2}∥ ≤ ∥A^{−1}∥ (∥b_{1} − b_{2}∥ + ∥F_{1} (X_{1}) − F_{1} (X_{2})∥ + ∥F_{1} (X_{2}) − F_{2} (X_{2})∥)
≤ ∥A^{−1}∥ (∥b_{1} − b_{2}∥ + L ∥X_{1} − X_{2}∥ + ε) ,

and that

∥X_{1} − X_{2}∥ ≤ \frac{∥A^{−1}∥ (∥b_{1} − b_{2}∥ + ε)}{1 − L ∥A^{−1}∥}.  (5)

Finally, from (4) and (5) we have

\frac{∥X_{1} − X_{2}∥}{∥X_{1}∥} ≤ \frac{∥A^{−1}∥ (∥A∥ + M)}{1 − L ∥A^{−1}∥} (\frac{∥b_{1} − b_{2}∥}{∥b_{1}∥} + \frac{ε}{∥b_{1}∥})

which proves the theorem.

If F_{1} ≡ 0 and F_{2} ≡ 0 (in this case we have L = M = ε = 0), then the proved inequality becomes (2).

Theorem 2. Let V be a normed space, let the linear operators A_{1}, A_{2} : V → V be invertible and bounded, let their inverse operators be also bounded, let b_{1}, b_{2} ∈ V and let the function F : V → V and the set S ⊆ V have the following properties:

1. the function F is Lipschitz on S, i.e.,

(∃L > 0) (∀x_{1}, x_{2} ∈ S) ∥F (x_{1}) − F (x_{2})∥ ≤ L ∥x_{1} − x_{2}∥,

and the constant L is such that the inequality 1 − L ∥A_{1}^{−1}∥ > 0 holds;

2. (∃M > 0) (∀x ∈ S) ∥F (x)∥ ≤ M ∥x∥ ;

3. the function F is bounded on the set S, i.e.,

(∃B ≥ 0) (∀x ∈ S) ∥F (x)∥ ≤ B.

If X_{1} ∈ S is a solution of the equation A_{1}x + F (x) = b_{1} and X_{2} ∈ S is a solution of the equation A_{2}x + F (x) = b_{2}, then the following inequality holds:

\frac{∥X_{1} − X_{2}∥}{∥X_{1}∥} ≤ \frac{∥A_{1}^{−1}∥ (∥A_{1}∥ + M)}{1 − L ∥A_{1}^{−1}∥} (\frac{∥b_{1} − b_{2}∥}{∥b_{1}∥} + ∥A_{1}∥ ∥A_{2}^{−1}∥ ×

×\frac{∥b_{2}∥}{∥b_{1}∥} · \frac{∥A_{1} − A_{2}∥}{∥A_{1}∥} + \frac{B ∥I − A_{1} · A_{2}^{−1}∥}{∥b_{1}∥}) .

Proof. Since X_{2} = A_{2}^{−1} · (b_{2} − F (X_{2})), we have
A_{1}X_{2} = A_{1}X_{2} + b_{2} − A_{2}X_{2} − F (X_{2})
= (A_{1} − A_{2}) X_{2} + b_{2} − F (X_{2})
= (A_{1} − A_{2}) A_{2}^{−1} (b_{2} − F (X_{2})) + b_{2} − F (X_{2})
= (A_{1} − A_{2}) A_{2}^{−1}b_{2} − (A_{1} − A_{2}) A_{2}^{−1}F (X_{2}) + b_{2} − F (X_{2})
= (A_{1} − A_{2}) A_{2}^{−1}b_{2} + b_{2} − A_{1}A_{2}^{−1}F (X_{2}) ,

and we can apply the previous theorem to the equations

A_{1}x + F (x) = b_{1}

and

A_{1}x + A_{1}A_{2}^{−1}F (x) = (A_{1} − A_{2}) A_{2}^{−1}b_{2} + b_{2}.

Condition 3. of the theorem is satisfied since for every x ∈ S the inequality

∥F (x) − A_{1}A_{2}^{−1}F (x)∥ ≤ ∥F (x)∥ ∥I − A_{1}A_{2}^{−1}∥ ≤ B ∥I − A_{1}A_{2}^{−1}∥

holds. So,

\frac{∥X_{1} − X_{2}∥}{∥X_{1}∥} ≤ \frac{∥A_{1}^{−1}∥ (∥A_{1}∥ + M)}{1 − L ∥A_{1}^{−1}∥} (\frac{∥b_{1} − b_{2} − (A_{1} − A_{2}) A_{2}^{−1}b_{2}∥}{∥b_{1}∥}+

+\frac{B ∥I − A_{1}A_{2}^{−1}∥}{∥b_{1}∥})

≤ \frac{∥A_{1}^{−1}∥ (∥A_{1}∥ + M)}{1 − L ∥A_{1}^{−1}∥} (\frac{∥b_{1} − b_{2}∥}{∥b_{1}∥} + ∥A_{1}∥ ∥A_{2}^{−1}∥ ×

×\frac{∥b_{2}∥}{∥b_{1}∥} · \frac{∥A_{1} − A_{2}∥}{∥A_{1}∥} + \frac{B ∥I − A_{1} · A_{2}^{−1}∥}{∥b_{1}∥}) .

The theorem has been proved.

If F ≡ 0 (in this case we have L = M = B = 0), then the inequality just proved becomes (3).

From the theorems just proved we can conclude that relatively small changes (of operator A, function F or vector b) in the equation (1) may cause relatively big changes in the solution if the number

\frac{∥A^{−1}∥ (∥A∥ + M)}{1 − L ∥A^{−1}∥}  (6)

is big enough, so we can take this number as a measure of stability of equation (1). It is obvious that the equation (1) gets more badly conditioned as the number (6) increases. Since the inequality ∥A∥∥A^{−1}∥ > 1 always holds, the number (6) is greater than one whenever inequality 1 − L∥A^{−1}∥ > 0 holds.

2. A note

If the normed space X is complete and the subset S ⊆ X is closed, if the function F satisfies the condition 1. of Theorem 1 (Theorem 2) and if φ (S) ⊆ S where φ (x) = A^{−1} (b − F (x)), then the array generated by the recursive formula

x_{n+1} = A^{−1} (b − F (x_{n})) , n ∈ ℕ  (7)

converges to the unique solution of the equation (1) for every x_{0} ∈ S.

Indeed, the function φ is a contraction since for every x, y ∈ S we have

∥φ (x) − φ (y)∥ ≤ ∥A^{−1}∥ ∥b − F (x) − b + F (y)∥ ≤ L ∥A^{−1}∥ ∥x − y∥,

while from the condition 1. of Theorem 1 (Theorem 2) we have that L ∥A^{−1}∥ < 1, and therefore in accordance with Banach fixed point theorem, the array defined by formula (7) will converge to the unique solution of the equation (1).

3. Examples

The first example will give (under certain assumptions) a sufficient condition for stability of polynomial with real coefficients. We thoroughly considered polynomials of the third degree.

Example 1. Let a polynomial with real coefficients P (x) = ax^{3} +bx^{2} +cx+d, (a, c ≠ 0) have at least one zero in the segment [α, β]. Furthermore, let F (x) = ax^{3} + bx^{2} and let Λ = max {|α| , |β|}. Then we have that (∀x ∈ [α, β]) |F (x)| ≤ (|a| Λ^{2} + |b| Λ) |x| and max_{x∈[α,β]} |F′ (x)| = max \{|F′ (α)| , |F′ (β)| , |F′ (−\frac{b}{3a})|\} and in Theorems 1 and 2 we can put that

M = |a| Λ^{2} + |b| Λ,

and that

L = max \{|F′ (α)| , |F′ (β)| , |F′ (−\frac{b}{3a})|\}.

If the condition 1 − \frac{L}{|c|} > 0 ⟺ L < |c| is satisfied then, in accordance with The-

orems 1 and 2 we can say that if the number \frac{|c|+M}{|c|−L} = \frac{1+M/|c|}{1−L/|c|} (which is always greater than one) is close enough to one, then relatively small changes in coefficients of the polynomial P will not cause relatively great changes in roots of the polynomial. So, if linear term in polynomial P is more dominant (|c| ≫ M and |c| ≫ L), the polynomial P is better conditioned.

We can do the same thing with polynomial of the fourth degree P (x) = ax^{4} + bx^{3} + cx^{2} + dx + e, (a, d ≠ 0) and conclude that the number \frac{|d|+M}{|d|−L} (the numbers M and L have the same meaning) can be used as a measure of stability of the polynomial P. So, the polynomial P is in this case also better conditioned if the number \frac{|d|+M}{|d|−L} is closer to one. Of course, we can use the same technics for the polynomials of higher degrees, but in that case the problem of effective finding of number L is much more complex.

Example 2. Let V be a normed space, let d ∈ V be a fixed vector, and let the function F: V → V be defined by

(∀x ∈ V) F (x) = ∥x∥ d.

We shall consider the relative error of solution when the equation A_{1}x + F (x) = b becomes the equation A_{2}x + F (x) = b. Since for every x, x_{1}, x_{2} ∈ V inequality

∥F (x_{1}) − F (x_{2})∥ ≤ ∥x_{1} − x_{2}∥ ∥d∥

and equality

∥F (x)∥ = ∥x∥ ∥d∥,

hold, we can put L = M = ∥d∥. So, if the condition 1 − ∥d∥ ∥A_{1}^{−1}∥ > 0, is

satisfied we can take the number \frac{∥A_{1}^{−1}∥ (∥A_{1}∥ + ∥d∥)}{1 − ∥d∥ ∥A_{1}^{−1}∥} as a measure of stability for the considered equation.

The following example is a numerical realization of Example 2.

Example 3. The solution of the system

max {x, y} + 2.01x − 1000y = 1000
max {x, y} − 0.01x − 1000y = −1000

is X_{1} = ( \begin{matrix}990.099 \\ 1.98020 \end{matrix} ), while the solution of the system

max {x, y} + 2.02x − 1000y = 1000
max {x, y} − 0.01x − 1000y = −1000,

is the vector X_{2} = ( \begin{matrix}985.222 \\ 1.97537 \end{matrix} ).

Stability of the considered system can be estimated by using the previous

example (V = ℝ^{2}, d = ( \begin{matrix}1 \\ 1 \end{matrix} ), and the norm is the uniform norm of the space

ℝ^{2}). The relative error of the matrix A_{1} = [ \begin{matrix}2.01 & −1000 \\ −0.01 & −1000 \end{matrix} ] when this matrix

becomes the matrix A_{2} = [ \begin{matrix}2.02 & −1000 \\ −0.01 & −1000 \end{matrix} ] is \begin{matrix}∥A_{1} − A_{2}∥ & ∞ \\ ∥A_{1}∥ & ∞ \end{matrix} ≈ 10^{−5} (10^{−3}%),

while the relative error of the solution when the first system becomes the sec-

ond one is \begin{matrix}∥X_{1} − X_{2}∥ & ∞ \\ ∥X_{1}∥ & ∞ \end{matrix} ≈ 0.5 · 10^{−2} (0.5%). So, the relative error of the so-

lution is approximately 500 times bigger then the relative error of the matrix A. According to the proved theorems our system is badly conditioned since

\frac{∥A_{1}^{−1}∥_{∞} (∥A_{1}∥_{∞} + ∥d∥_{∞})}{1 − ∥d∥_{∞} ∥A_{1}^{−1}∥_{∞}} = 100301.

It should be noted that the influence of nonlinear term in this example is

irrelevant. The relative error of solution, when linear system A_{1}x = b = ( \begin{matrix}1000 \\ −1000 \end{matrix} ) becomes the system A_{2}x = b is approximately 0.5%, too.

We would like to point out that this system may also be solved by using the Banach fixed-point theorem (see Section 2).

The first one of the following examples has a theoretical character, while the second one is its numerical realization.

Example 4. Let V be a normed space, let d ∈ V and r > 0 be a fixed vector and a real number, let S = {x ∈ V |∥x∥ ≤ r } and let the function F: V → V be defined by

(∀x ∈ V) F (x) = ∥x∥^{2} d.

We shall estimate the relative error of the solution when equation A_{1}x + F (x) = b becomes equation A_{2}x + F (x) = b. Since for every x, x_{1}, x_{2} ∈ S inequalities

∥F (x_{1}) − F (x_{2})∥ = ∥∥x_{1}∥^{2} d − ∥x_{2}∥^{2} d∥
= (∥x_{1}∥ + ∥x_{2}∥) · |∥x_{1}∥ − ∥x_{2}∥| · ∥d∥
≤ 2r · ∥d∥ · ∥x_{1} − x_{2}∥

and

∥F (x)∥ = ∥x∥^{2} ∥d∥ ≤ r ∥d∥ ∥x∥,

hold, we can put M = r ∥d∥ and L = 2r ∥d∥. So, if the condition 1−2r ∥d∥ ∥A_{1}^{−1}∥ >

0 is satisfied then the number \frac{∥A_{1}^{−1}∥ (∥A_{1}∥ + r ∥d∥)}{1 − 2r ∥d∥ ∥A_{1}^{−1}∥} can be taken as a measure of

stability of the considered equation.

Example 5. The solution of the system

x^{2} + y^{2} + 750x + 50y = −1
x^{2} + y^{2} + 2x − 3y = −1

which belongs to the set S = \{( \begin{matrix}x \\ y \end{matrix} ) ∈ ℝ^{2} | x^{2} + y^{2} ≤ 1\} is X_{1} = ( \begin{matrix}−0.0254856 \\ 0.359684 \end{matrix} ),

while the solution of the system

x^{2} + y^{2} + 750x + 50y = −1
x^{2} + y^{2} + 2x − 2y = −1

which belongs to the set S = \{( \begin{matrix}x \\ y \end{matrix} ) ∈ ℝ^{2} | x^{2} + y^{2} ≤ 1\} is the vector X_{2} =

( \begin{matrix}−0.0481967 \\ 0.693291 \end{matrix} ).

Stability of the considered system can be estimated by using Example 4 (V =

ℝ^{2}, S = \{( \begin{matrix}x \\ y \end{matrix} ) ∈ ℝ^{2} | x^{2} + y^{2} ≤ 1\}, d = ( \begin{matrix}1 \\ 1 \end{matrix} ), while the norm is the Euclidean norm of the space ℝ^{2}). Relative error of the matrix A_{1} = [ \begin{matrix}750 & 50 \\ 2 & −3 \end{matrix} ] when this

matrix becomes the matrix A_{2} = [ \begin{matrix}750 & 50 \\ 2 & −2 \end{matrix} ] is \begin{matrix}∥A_{1} − A_{2}∥ & 2 \\ ∥A_{1}∥ & 2 \end{matrix} ≈ 0.13 · 10^{−2} (0.13%), while the relative error of the solution when the first system becomes the second one is \begin{matrix}∥X_{1} − X_{2}∥ & 2 \\ ∥X_{1}∥ & 2 \end{matrix} ≈ 0.93 (93%). So, the relative error of the solution is approximately 700 times bigger than the relative error of matrix A. According to the proved theorems the system is badly conditioned since \frac{∥A_{1}^{−1}∥_{2} (∥A_{1}∥_{2} + ∥d∥_{2})}{1 − 2 ∥d∥_{2} ∥A_{1}^{−1}∥_{2}} = 2527.

Contrary to Example 3, the influence of nonlinear term is important now. In this example, the relative error of solution when linear system A_{1}x = b = ( \begin{matrix}−1 \\ −1 \end{matrix} ) becomes the system A_{2}x = b is approximately 47% .

REFERENCES

[1] Bahvalov, N. S., Numerical Methods, Science, Moscow, 1973.

[2] Edwards, R. E., Functional Analysis, Holt, Rinehart and Winston, Inc., New York, 1965.

(received 12.05.2004, in revised form 14.02.2005) Faculty of Sciences, Department of Mathematics, Zmaja od Bosne 35, 71000 Sarajevo, Bosnia and Herzegovina

E-mail: zzlatko@pmf.unsa.ba