Page 242 - MATLAB an introduction with applications

Numerical Methods ——— 227

function [Y,iter] = GAUSSD(A,r,yguess,tol)
% GAUSSD will iteratively solve Ay = r
n = length(r); Y = yguess; dy = ones(n,1); iter = 0;
while norm(dy)/norm(Y) > tol
for i = 1:n
if abs(A(i,i))<100*eps;error(‘zero pivot found’);end
dy(i) = r(i)/A(i,i);
for j = 1:n
dy(i) = dy(i) - A(i,j)*Y(j)/A(i,i);
end
Y(i) = Y(i) + dy(i);
end
iter = iter + 1;
if iter>1000; error(‘not converging in 1000 steps’);end
end

Example E4.11: Solve the system of equations given below by Householder’s factorization method:

x

 4 − 1 0 0  1  
1
 − 1 4 − 1   
x

0 

 0  =  
2

 0 − 1 4 − 1   
x
0
  3
x
 
  0 0 − 1 4   0 
  
4
Solution:
>> A = [4 –1 0 0;–1 4 –1 0;0 –1 4 –1;0 0 –1 4];
>> b = [1;0;0;0];
>> householder(A)
ans =
4 1 0 0
–2 4 –1 0
0 2 4 1
0 0 1 4
>> [L,U] = lu(A)
L =
1.0000 0 0 0
–0.2500 1.0000 0 0
0 –0.2667 1.0000 0
0 0 –0.2679 1.0000
U =
4.0000 –1.0000 0 0
0 3.7500 –1.0000 0
0 0 3.7333 –1.0000
0 0 0 3.7321