Scientific Computing

Scientific Computing Matrix Norms, Convergence, and Matrix Condition Numbers

Vector Norms 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 )

Vector Norms Example: x = ( 4, -1 ) is the standard Pythagorean length of x. This is one possible measurement of the size of x. x

Vector Norms Example: x = ( 4, -1 ) |4| + |-1|=5 is the “Taxicab” length of x. This is another possible measurement of the size of x. x

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

Vector Norms 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.

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

Vector Norms Example

Vector Norms • 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

Distance

Which norm is best? • 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!

Matlab Vector Norms • 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

Matrix Norms • 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).

Matrix Norm “Stretch” • 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.

Matrix 1- Norm • 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.

Matrix 1- Norm • 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. □

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

Matrix Norm Properties • || 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 ||

Eigenvalues-Eigenvectors • 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:

Spectral Radius • 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.

Convergence 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. □

Convergent Matrix Series 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. □

Convergence of Iterative solution to Ax=b 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

Convergence of Iterative solution to Ax=b 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. □

Diagonally Dominant Matrices 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:

Jacobi Method 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

Jacobi Method 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.

Jacobi Method 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. □

Gauss-Seidel Method 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.

Scientific Computing