6-21-04
Newton-Raphson
·
It is used to the find the solution for f(x) = 0·
Makes use of the derivative f'(x)·
Uses f(x [0]), f'(x [0]) and a point x [0] to find an approximation to a line of the function f(x)
·
Method·
Can arrive at the same iteration equation via Taylor polynomial·
f(x) ≈ f(x[0]) + f'(x[0])(x - x[0]) + e, where e is the error of using a line instead of the complete polynomial
Example of Newton-Raphson
sin(x) - ½ = 0 (radians)
f(x) = sin(x) - ½ = 0 f '(x) = cos(x)
Start with x0 = 0
x1 = x0 - f(x0)/f '(x0)
| i | xi | f(xi) | f '(xi) | xi+1 |
| 0 | 0 | f(x0) = sin(0)-.5 = -.5 | f '(x0)= cos(0) = 1 | x1 = 0-(-.5) / 1 = -.5 |
| 1 | .5 | f(x1) = sin(.5)-.5 = -.02 | f '(x1)= cos(.5) = .8795 | x2 = .5-(-.02) / .8795 = .522 |
| 2 | .522 | f(x2)=sin(.522)-.5 =-.00139 | f '(x2)= cos(.522) = .8668 | x3 = .522-(-.00139) / .8668 = .5235 |
Exact solution for sin(x) -.5 = 0 in radians is 30° that is 30π/180 rad =
.5235 rad.·
One disadvantages of Newton Raphson is that you need the derivative f '(x) which takes more effort to do.·
It is possible to obtain the derivative f '(x) of f(x) numerically.f '(x) = f(x+ε) - f(x) / ε
This method is used in project #1 in class.
The exact solutions require that ε → 0. We can do an approximation by using a small ε (like 1x10)
-5.Example:
f(x) = sin(x)
f '(x) ≈ sin(x+ε) - sin(x) / ε
if ε = .001 x = 0
When Newton-Raphson may not give a good approximation/converge
·
If the function does not approximate well using a line·
If the starting point is not near enough for decent approximation6-22-04
Order of Convergence
Some methods are faster than others at converging of f(x) = 0.
p is the zero (solution), E[n] = p - p[n] is the error for iteration n
* A != 0, R > 0
* lim[n->∞] E[n+1]/|E[n]|^R = A
* If R = 1, the convergence is linear
* If R = 2, the convergence is quadratic
Example let En = .001 and A = 1 and R = 2
En+1 = 1 |.001|2
= 1(1x10-3)2
= 1x10-6 = .000001
If R = 1, the convergence is called "linear".
If R = 2, the convergence is called "quadratic"
If f(x) has a simple root and we use Newton Raphson the R = 2 (quadratic
convergence)
Usually, the more roots an equation has, the slower the convergence will be.
Secant Method
·
The Newthon Rapson method requires the evaluation of two functions at each iteration (f(x) and f'(x))·
This method requires only one evaluation of f(x)·
Order of convergence with a single root is ~1.618·
Uses two initial points p[0], p[1] close to the root
(1) m = [f(x1) - f(x0)] / (x1 - x0) using (x0, f(x0)) and (x1, f(x1))
(2) m = [f(x2) - f(x1)] / (x2 - x1) using (x1, f(x1)) and (x2, 0)
Putting (1) and (2) together....
[f(x1) - f(x0)] / (x1 - x0) = [0 - f(x1)] / (x2 - x1)
x2 - x1 = - f(x1)*(x1 - x0) / [f(x1) - f(x0)]
x2 = x1 - f(x1)*(x1 - x0) / [f(x1) - f(x0)]
6-23-04
(Continue)
· Approximates Newthon-Rapson method if x[1] -> x[0]
· x[2] = x[1] - (f(x[1]))((x[1] - x[0])/(f(x[1]) -
f(x[0])))
· As x[1] -> x[0], it becomes x[1] -
f(x[0])/f'(x[0]), which is N-R
Example:
f(x) = sin(x) - 1/2 = 0
x[0] = 0, x[1] = 1
x[n+1] = x[n] - f(x[n])((x[n] - x[n-1])/(f(x[n]) - f(x[n-1])))
| n | xn | f(xn) | xn+1 |
| 0 | x0 = 0 | f(x0) = sin(0) - .5 = -.5 | |
| 1 | x1 = 1 | f(x1) = sin(1) - .5 = .3415 | x2 = 1 - .3415(1 - 0) / .3415 - (-.5) = .594 |
| 2 | x2 = .594 | f(x2) = sin(.594) - .5 = .0597 | x3 = .594 - .0597(.594 - 1) / .0597 - (.3415) = .50798 |
| 3 | x3 = .50798 | f(x3) = sin(.50798) - .5 = -.0135 | x4 = .50798 - [-.0135(.50798 - .594) / (-.0135 - .0597)] = .5238 |
The Solution of Linear Systems
·
Ax = B, where A is a matrix (n x n) and both x and B are vectors (1 x n)·
Each linear equation represents a plane in the nth dimension·
The solution represents the intersection of the planes·
Not all systems of linear equations have solutions·
Some systems have infinite solutions
Upper triangular systems of linear equations
·
Elements are in row-then-column listing a[i][j]·
Upper triangular if a[i][j] = 0 when i > j·
Easy to solve by obtaining the value for x[n] (from the last equation) and substituting into progressively higher equations, called "back substitution"6x + 3y + 2z = 29
3x + y + 2z = 17
2x + 3y + 2z = 21
This can be expressed as: Ax = B where A is a matrix and B is a vector.
| 6 3 2 | | X | | 29 |
| 3 1 2 | | Y | = | 17 |
| 2 3 2 | | Z | | 21 |
A X B
This can be expressed as: Ax = B where A is a matrix and B is a vector.
| 6 3 2 | | X | | 29 |
| 3 1 2 | | Y | = | 17 |
| 2 3 2 | | Z | | 21 |
A X B
·
Each linear equation represents a plane in the nth dimension·
Not all the systems of linear equations have solutions(1) 2x + 3y = 6 (3) x = (6 - 3y)/2 from (1)
(2) 4x + 6y = 10
substitute (3) in (2)
4*(6-3y)/2 + 6y = 10
12 - 6y + 6y = 10
12 = 10
Geometrically, these equations represent 2 planes that never intersect.
You also can have a system of equations with an infinite number of solutions....
(1) 2x + 3y = 6
(2) 4x + 6y = 12
(1) and (2) are the same equation so both are the same plane. Therefore, the
number of solutions is infinite.
Gausian Elimination
·
Method for solving systems of linear equations·
Transforms an arbitrary system into upper triangular form, then solves it by using back substitution·
The transformation uses the following
6-24-04
We will use these three properties to transform an arbitrary system of linear equations into an upper triangular system and then solve it with back-substitution.
Example
1x + 3y + 4z = 19
8x + 9y + 3z = 35
x + y + z = 6
Solve system using gaussian elimination.
Augmented matrix.
Choose a pivot row and pivot element. In this case, 1 is the pivot element.
| 1 3 4 19 |
| 8 9 3 35|
| 1 1 1 6 |
We need to make element a21 = 8 to be Ø so we multiply the first equation by 8 and subtract the result from second equation. We also multiply the first equation by -1 and add it to the third equation.
↓
| 1 3 4 19 |
| 0 15 29 117 |
| 0 2 3 13 |
Multiply second row by -1/15
| 1 3 4 19 |
| 0 1 29/15 117/15 |
| 0 -2 -3 -13 |
We need to make element a32 = 0 so we multiply second equation by 2 and add it to third equation.
↓
| 1 3 4 19 |
| 0 1 29/15 117/15 |
| 0 -2+2=0 -3+2*29/15 -13+2*117/15 |
Back substitution (-3+2*29/15)z = -13 + 2*117/15 .
z = 3
1y + (29/15)(3) = 117/15, y = 117/15 - (29/15)(3) = 2, y = 2
x + 3(2) + 4(3) = 19, x + 6 + 12 = 19, x = 19 - 18 = 1, x = 1
Implementing Gaussian Elimination using Matlab
function x = gauss(A, b)
% Input - A is an nxn nonsingular matrix
% B is an nx1 matrix
% Output - x is an nx1 matrix that is the solution to Ax = B
[N N] = size(A)
x = zeros(N, 1);
Aug = [A B];
for p=1 : N
piv = Aug(p, p);
for k = p : N+1
Aug(P, k) = Aug(P, k)/piv;
end
for k = p+1 : N
m = Aug(k, p);
for i = p : N+1
Aug(k, i) = Aug(k, i) - m*Aug(p, i);
end
end
end
% Check on "step" being the right way to indicate decrement
for k = N : 1 (step -1)
sum = 0;
for i = k + 1 : N
sum = sum + Aug(k, i)*x(i, 1)
end
x(k, 1) = Aug(k, N+1) - sum
end
6-25-04
(1) Prevent pivot = 0
If pivot = 0, switch rows to prevent division by 0 (to reduce computing error,
choose the pivot to be the row with the largest starting absolute value)
(2) To reduce computing error, you may choose as the next pivot the largest
number in the rows (or columns) that are still left to transform.
Triangular Factorization
·
A matrix has a triangular factorization if it can be expressed as the product of a lower triangular matrix and an upper triangular matrix·
Useful when you want to solve multiple systems of linear equations that have the same A (remember, they're of the form Ax = B)·
Replace A with LU (for lower and upper triangular matrices respectively), then declare Ux to be Y, solve for Y, then substitute to solve Ux = Y·
To split into L and U, multiply the matrix by an identity matrix (so it's sitting just to the side of the starting matrix), then use Gaussian Elimination to both matrices (matching in the identity matrix the moves that are performed in the starting matrix)·
When only a single system needs to be solved, Gaussian Elimination is better (as it requires less work)