% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
%                                                                       % 
% Loesung der Aufgabe 11, Matlab Kurs WS 06/07                          % 
%                                                                       % 
% Dieses Programm liest eine Matrix A und einen Vektor b ein und loest  %
% das Gleichungssystem Ax=b mit Hilfe der Gausselimination ohne Pi-     %
% votsuche.                                                             % 
%                                                                       % 
% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %

disp('   ')
disp('   ')
disp('Welche Dimension soll die quadratische (N,N)-Matrix haben?')
N = input('Bitte N eingeben, N=');
disp('   ')
disp('Bitte eine quadratische Matrix A eingeben, dabei die richtige ')
disp('Matlabnotation nutzen, zB. [2,1,1;3,-1,5;1,0,0].')
A = input('A = ');
disp('   ')
disp('Bitte einen Spaltenvektor b der richtigen Dimension eingeben,')
disp('zB. [1;2;3]')
b = input('b = ');
disp('   ')


% Stimmt die Wahl von N?:
[n1,n2]=size(A);
if n1~=N
    error('A hat nicht N Zeilen!')   % error bewirkt Programmabbruch!
end
if n2~=N
    error('A hat nicht N Spalten!')
end

% Gauss Elimination:
% ueberschreibe die Matrix A in jedem Schritt mit den neuen Werte, berechne
% nur die Werte, die wirklich ~= 0 sind. Der Algorithmus setzt keine Werte
% von A auf Null! Da hinterher beim Rueckwaertseinsetzen nur noch auf die 
% Elemente in der rechten oberen Dreiecksmatrix zugegriffen wird, ist dies
% nicht notwendig!!!

% In die linke untere Ecke von A wird in Zeile (*) die Matrix L geschrieben
for j=1:N
    for i=j+1:N
       A(i,j)= A(i,j)/A(j,j);                          % (*)
       for k=j+1:N
           A(i,k)=A(i,k)-A(i,j)*A(j,k);
       end
       b(i) = b(i) - A(i,j)*b(j);
    end
end


% Rueckwaertseinsetzen
x = zeros(N,1);


for j=N:-1:1
	summe = 0;
	for k=j+1:N
		summe = summe + A(j,k) * x(k);
    end
	x(j) = (b(j)- summe) / A(j,j);
end


% Ausgabe der Loesung:
disp('   ')
disp('   ')
disp('Die Loesung des Gleichungssystems Ax=b ist x=')
disp(x)
