- 48 Views
- Uploaded on
- Presentation posted in: General

Scientific Computing

Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author.While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server.

- - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - -

Matrix Norms, Convergence, and

Matrix Condition Numbers

A vector norm is a quantity that measures how large a vector is (the magnitude of the vector).

For a number x, we have |x| as a measurement of the magnitude of x.

For a vector x, it is not clear what the “best” measurement of size should be.

Note: we will use bold-face type to denote a vector. ( x )

Example: x = ( 4, -1 )

is the standard Pythagorean length of x. This is one possible measurement of the size of x.

x

Example: x = ( 4, -1 )

|4| + |-1|=5 is the “Taxicab” length of x. This is another possible measurement of the size of x.

x

Example: x = ( 4, -1 )

max(|4|,|-1|) =4 is yet another possible measurement of the size of x.

x

A vector norm is a quantity that measures how large a vector is (the magnitude of the vector).

Definition: A vector norm is a function that takes a vector and returns a non-zero number. We denote the norm of a vector x by

The norm must satisfy:

Triangle Inequality:

Scalar:

Positive: ,and = 0 only when x is the zero vector.

- Our previous examples for vectors in Rn :
- All of these satisfy the three properties for a norm.

- Definition: The Lp norm generalizes these three norms. For p > 0, it is defined on Rn by:
- p=1 L1 norm
- p=2 L2 norm
- p= ∞ L∞ norm

- The answer depends on the application.
- The 1-norm and ∞-norm are good whenever one is analyzing sensitivity of solutions.
- The 2-norm is good for comparing distances of vectors.
- There is no one best vector norm!

- In Matlab, the norm function computes the Lp norms of vectors. Syntax: norm(x, p)
>> x = [ 3 4 -1 ];

>> n = norm(x,2)

n = 5.0990

>> n = norm(x,1)

n = 8

>> n = norm(x, inf)

n = 4

- Definition: Given a vector norm ||x|| the matrix norm defined by the vector norm is given by:
- What does a matrix norm represent?
- It represents the maximum “stretching” that A does to a vector x -> (Ax).

- Note that, since ||x|| is a scalar, we have
- Since is a unit vector, we see that the matrix norm is the maximum value of Az where z is on the unit ball in Rn.
- Thus, ||A|| represents the maximum “stretching” possible done by the action Ax.

- Theorem A: The matrix norm corresponding to 1-norm is maximum absolute column sum:
- Proof: From the previous slide, we have
- Also,
- where Aj is the j-th column of A.

- Proof (continued): Then,
- Let x be a vector with all zeroes, except a 1 in the spot where ||Aj|| is a max. Then, we get equality above. □

- Theorem B: Matrix norm corresponding to ∞ norm is maximum absolute row sum:
- Proof (similar to Theorem A).

- || A || > 0 if A ≠ O
- || A || = 0 iff A = O
- || c A || = | c| *||A || if A ≠ O
- || A + B || ≤ || A || + || B ||
- || A B || ≤ || A || *||B ||
- || A x || ≤ || A || *||x ||

- The eigenvectors of a matrix are vectors that satisfy
- Ax = λx
- Or, (A – λI)x = 0
- So, λ is an eigenvalue iff det(A – λI) = 0
- Example:

- The spectral radius of a matrix A is defined as
- ρ(A) = max |λ| where λ is an eigenvalue of A
- In our previous example, we had
- So, the spectral radius is 1.

Theorem 1: If ρ(A)<1, then

Proof: We can find a basis for Rn by unit eigenvectors (result from linear algebra), say {e1, e2, …, en}. Then,

For any unit vector x, we have

x = a1 e1 + a2 e2 + … + an en

Then, An x = a1 Ane1 + a2 Ane2 + … + an Anen

= a1 λ1ne1 + a2 λ2ne2 + … + an λnnen

Thus,

Since ρ(A)<1, then the result must hold. □

Theorem 2: If ρ(B)<1, then (I-B)-1 exists and

(I-B)-1 = I + B + B2 + · · ·

Proof: Since we have Bx = λx exactly when (I-B)x = (1- λ)x, then λ is an eigenvalue of B iff (1- λ) is an eigenvalue of (I-B). Now, we know that |λ|<1, so 0 cannot be an eigenvalue of (I-B).

Thus, (I-B) is invertible (why?). Let Sp = I + B + B2 + · · ·+Bp Then, (I-B) Sp= (I + B + B2 + · · ·+Bp ) – (B + B2 + · · ·+Bp+1 )

= (I-Bp+1)

Since ρ(A)<1, then by Theorem 1, the term Bp+1will go to the zero matrix as p goes to infinity. □

Recall: Our general iterative formula to find x was

Q x(k+1) = ωb + (Q-ωA) x(k)

where Q and ω were variable parameters.

We can re-write this as

x(k+1) = Q-1 (Q-ωA) x(k) + Q-1ωb

Let B = Q-1 (Q-ωA) and c = ωb Then, our iteration formula has the general form:

x(k+1) = B x(k) + c

Theorem 3: For any x0 in Rn , the iteration formula given by x(k+1) = Bx(k) + c will converge to the unique solution of x=Bx+c (i.e fixed point) iff ρ(B)<1.

Proof:

If ρ(B)<1, the term Bk+1x0 will vanish. Also, the remaining term will converge to (I-B)-1. Thus, {x(k+1)} converges to z = (I-B)-1c, or

z-Bz = c or z = Bz + c.

The converse proof can be found in Burden and Faires, Numerical Analysis. □

Def: A matrix A is called Diagonally Dominant if the magnitude of the diagonal element is larger than the sum of the absolute values of the other elements in the row, for all rows.

Example:

Recall: Jacobi Method

x(k+1) = D-1(b + (D-A) x(k))

= D-1(D-A) x(k) + D-1b

Theorem 4: If A is diagonally dominant, then the Jacobi method converges to the solution of Ax=b.

Proof: Let B = D-1(D-A) and c = D-1b. Then, we have x(k+1) = B x(k) + c. Consider the L∞ norm of B, which is equal to

Proof: (continued) Then,

If A is diagonally dominant, then the terms we are taking a max over are all less than 1. So, the L∞ norm of B is <1. We will now show that this implies that the spectral radius is <1.

Lemma: ρ(A)<||A|| for any matrix norm.

Proof: Let λ be an eigenvalue with unit eigenvector x.

□

Proof of Theorem 4 (cont): Since we have shown that

then, by the Lemma, we have that ρ(B) < 1.

By Theorem 3, the iteration method converges. □

Through similar means we can show (no proof):

Theorem 5: If A is diagonally dominant, then the Gauss-Seidel method converges to the solution of Ax=b.