%                                 lecture 15
%                        Der Einfluß von Rundungsfehlern
%

% Auslöschung

x = 10.^[0:5]';
f = x.*(sqrt(x+1)-sqrt(x));
f_red = runde(x.*(runde(sqrt(x+1),4)-runde(sqrt(x),4)),4);
[f f_red]

% Wir lösen folgendes Gleichungssystem 
% Ax = b
% mit
A  = [.003 59.14
     5.291  -6.13]
 
 b = [59.17; 46.78]
 
 % mit exakter Lösung  x= [10;1]
 
 x = A\b
 
 % Wir lösen nun dieses Gls  noch einmal mit dem Gauß Algorithmus 
 % und legen diesmal eine 4 stellige Arithmetik zugrunde
 
 % erweiterte Matrix
 
 Ab = [A b];
 
 % Elimination
 
 Ab_red = gauss_red(Ab,4)
 
 % Rücksubstitution
 
 x_red = backsubstitute_red(Ab_red(:,1:2),Ab_red(:,3),4)
 
 % Was ist passiert:
 % Der Multiplikator ist sehr groß, da a_11 sehr klein im Vergleich zu a_21
 % ist
 m = runde(runde(A(2,1),4)/runde(A(1,1),4),4)
 
 % für das neue (2,2) Element erhält man
 s=A(2,2)-runde(m*A(1,2),4)
 % und nach Rundung
 runde(s,4)
 % also A(2,2) wird nicht berücksichtigt
 % analog bei b(2)
t= b(2) - runde(m*b(1),4)
%nach Rundung
runde(t,4)
% also konnte b(2) nicht berücksichtigt werden
 
% Die berechnete Lösung löst also das GLS mit
A1=A;A1(2,2)=0
b1=b; b1(2)=0
Ab1 = [A1 b1]; 
Ab1_red = gauss_red(Ab1,4) 
x1_red = backsubstitute_red(Ab1_red(:,1:2),Ab1_red(:,3),4) 
 

% Durch Vertauschen der beiden Gleichungen kann das Problem einfach gelöst
% werden --> Teilpivotstrategie

 Ab_red = gauss_red(flipud(Ab),4)
 x_red = backsubstitute_red(Ab_red(:,1:2),Ab_red(:,3),4)
 
 
 % Ein weiteres beunruhigendes Phänomen:
 % Die Hilbertmatrix A(i,j) = 1/(i+j-1) ist der Prototyp einer schlecht
 % konditionierten Matrix:
 
 clear, n=5;
 A = hilb(n)
 b = sum(A,2)
 
 % Die exakte Lösung des GLS Ax=b ist somit x= (1,1,...,1)
 
 x = A\b
 
 A = runde(A,4);
 x1 = A1\b
 
 
 