Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Special matrices
Week 7: aim to cover
Vector & matrix norms, sensitivity (Lecture 13) chol, \ (Lab 7)
error analysis of linear systems (Lecture 14)
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Vector and Matrix norms
In order to discuss sensitivity and numerical stability in solving a linear system, need to measure the ‘size’ of errors in the inputs and outputs. In the problem of solving Ax = b,
the inputs are A, b: a matrix and a vector
the output is x: a vector
So have to introduce a way to measure the ‘size’ of vectors and matrices.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Vector norms
Definition
The norm of a vector kxk is a function Rn 7! R such that
1 kxk 08x2V wherekxk=0,x=0(normsarepositivefor
nonzero vectors)
2 k↵xk = |↵|kxk (scaling a vector scales norm by the same amount) 3 kx+ykkxk+kyk8x,y2V (triangleinequality)
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Example
The 3 most cPommon vector norms are: 1 kxk1 = P|xi| (the 1-norm)
2 kxk2 = ( xi2)1/2 = pxT x (the 2-norm i.e. the usual Euclidean norm)
3 kx k1 = maxi |xi | (the 1-norm)
which are all special cases of the p-norm:
kxkp = (X|xi|p)1/p
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Example
Unit vectors in 1,2 and 1 norms
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Norm equivalence
In finite-dimensional spaces, it doesn’t matter much precisely which norm you use, since they can’t di↵er by more than a factor n from each other.
kxk2 kxk1 pnkxk2
kxk1 kxk2 pnkxk1
kxk1 kxk1 nkxk1
So you may as well use whichever one is convenient.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Matrix norms
A matrix norm is just a vector norm on the m ⇥ n dimensional vector space of m ⇥ n matrices.
Definition
The norm of a matrix kAk is a function Mm⇥n 7! R such that
1 kAk 08A2Mm⇥n wherekAk=0,A=0(normsarepositive
for nonzero matrices)
2 k↵Ak = |↵|kAk (scaling a matrix scales norm by the same amount)
3 kA + Bk kAk + kBk 8A, B 2 Mm⇥n (triangle inequality)
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Example
The following Pare matrix norms:
1 kAkF = ( A2ij )1/2 (the Frobenius norm) 2 kAkM = maxi ,j |Aij | (the max-norm)
Mostly we’ll use a common class of matrix norms, called subordinate matrix norms.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Subordinate matrix norms
Definition
The subordinate norm of a (square) matrix A is given by kAk = max kAxk
x6=0 kxk for any vector norm. This is equivalent to:
kAk = max kAxk kx k=1
In words, the (subordinate) norm of a matrix is the norm of the largest image under the map A of a unit vector.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Matrix norms
Example
The subordinate pP-norms correspond to the vector norms listed above: 1 kAk1 = maxj i |Aij | (the maximum column sum)
2 kAk2 (the 2-noPrm )
3 kAk1 = maxi j |Aij | (the maximum row sum)
These are the only subordinate p-norms that are easy to compute.
The subordinate norms have some useful submultiplicative properties:
by definition; and
kAxkp kAkpkxkp kABkp kAkpkBkp
Any norm with the latter property is called consistent.
Again, all matrix norms are within a factor of n of each other, so you may as well use whichever one is convenient.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Error analysis in numerical linear algebra
In analyzing the errors made in Gauss elimination, Wilkinson (1960) realized there were di↵erent sources of error in a computation:
the sensitivity of the problem (whatever algorithm you use)
the quality of your algorithm to solve that problem
These ideas are now used in many other areas of numerical analysis.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Forward vs. backward error
Given a computation Y = f (X ), in general we will produce an approximation Yˆ . The forward error is Y = Yˆ Y . We need a small forward error for an accurate answer BUT it can be hard to estimate forward error.
Instead ask a di↵erent question:
Is the computed answer the exact answer to a nearby problem?
Is Yˆ = f (X + X) for some small backward error X? Did we solve a problem close to the one we wanted to solve? It turns out to be often easier to bound the backward error.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Diagram
Example
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Backward error
If the backward error is not too big (compared to the relevant data error) we say the algorithm is backward stable or just stable.
Example
For solving linear systems, the input errors are roundo↵ errors so we want the relative backward error not too big compared to unit roundo↵ u.
Once we have introduced the concept of backward error, the error made by the algorithm can be treated as if it was error in the data. Then the forward error only depends on the sensitivity of the problem f and the size of the input error i.e. the backward error.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Sensitivity = conditioning
A problem which is sensitive to small errors is called ill-conditioned or ill-posed =) no numerical method will get a very accurate answer. Typically sensitivity is determined using perturbation analysis , assuming small errors.
We quantify sensitivity by defining a condition number to measure how sensitive the answer is to errors in the input.
Because roundo↵ errors are relative errors, we use the relative condition number ⇠ relative error of output/relative error of input i.e. the magnification factor of relative error
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Rule of thumb
Forward error ⇠ condition number ⇥ backward error
a stable method on a well-conditioned problem ! accurate answer
this is what we want!
a stable method on an ill-conditioned problem ! inaccurate answer re-examine the formulation of your problem!
an unstable method on a well-conditioned problem ! inaccurate answer
This is what we must avoid!
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Trefethen’s Maxims
If the answer is highly sensitive to perturbations, you have probably asked the wrong question.
All that it is reasonable to ask for in a scientific calculation is stability, not accuracy.
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Sensitivity of a linear system
Suppose A, b are perturbed by A, b — how much is the solution x changed?
To measure the input and output errors, we use vector and matrix norms.
Theorem
k xk (A) (k Ak + k bk) kxk 1 (A)k Ak kAk kbk
kAk
provided (A) k Ak < 1, a condition which ensures that the matrix
kAk
A + A is nonsingular.
(A) is the (normwise relative) condition number for solving a linear system
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
The condition number
Definition
The (normwise) condition number of a square nonsingular matrix is (A)=kAkkA 1 k
For small enough k Ak, we will have (A) k Ak ⌧ 1 so that kAk
relative forward error . (condition number) (relative error in A, b)
Note: it’s not cheap to compute (A) since it takes ⇡ n3 operations to find A 1. Most codes instead try to estimate (A).
k xk . (A)(k Ak + k bk) kxk kAk kbk
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Some properties of the condition number
The precise value depends on what norm you’re using but they will all be quite similar in size (to within factors of n)
1 (I) =k I kk I 1 k=k I k2= 1 in any subordinate norm 2 (A)=kAkkA 1 k kAA 1 k=kIk=1so(A) 1
Matrices with (A) 1 are ill-conditioned ! the solution is very sensitive errors in A or b (in the worst case)
Numerical Methods & Scientific Computing: lecture notes
Linear Systems
Error analysis
Heuristic
If (A) ⇠ 10k then in solving Ax = b in t-digit arithmetic, you will lose k decimal digits of precision =) x has t k digits of precision
=) if (A) > 10t it’s not worth solving the system since x may have no significant figures!
Example
The Hilbert matrix