%% Blatt 3 - Programmieraufgaben:
%% LR-Zerlegung (4 Punkte)
% Vorname Nachname Mat.-Nr.
%% Gauss-Algorithmus mit Rang-1-Matrizen Updates
% Die LR-Zerlegung kann durch das in der Vorlesung in Gauss-Algorithmus revisited dargestellte Vorgehen berechnet werden. Eine direkte Implementation davon benutzt dabei jeweils drei geschachtelte for-Schleifen und kann folgendermassen realisiert werden:
% Input: \mathbf{A} = [a_{i, j}]_{i, j=1}^n \in \mathbb{R}^{n \times n}
% Output: LR-Zerlegung von \mathbf{A}, gespeichert in \mathbf{A}
% for j = 1, \ldots , n - 1 do
%     for i = j+1, \ldots, n do
%         setze a_{i, j} : \! = a_{i, j} / a_{j, j}
%                  for k = j+1, \ldots , n do
%             setze a_{i, k} : \! = a_{i, k} - a_{i, j} \cdot a_{j, k}
%                  end for
%         end for
% end for
% Bei der Implementation in MATLAB aber ist es möglich, die zwei inneren expliziten Schleifen durch Vektor- und Matrixblockzuweisungen und Rang-1-Matrizen zu ersetzen. Für den Rest der Aufgabe vereinbaren wir hierzu folgende Notation: Sei \mathbf{A} = [a_{i, j}]_{i, j=1}^n \in \mathbb{R}^{n \times n} eine Matrix. Dann definieren wir
% \mathbf{A}_{i:l, j:k} := \pmatrix{a_{i, j} & \cdots & a_{i, k} \cr \vdots & & \vdots \cr a_{l, j} & \cdots & a_{l, k}} .
% Der Algorithmus für die LR-Zerlegung mit Rang-1-Matrizen Updates ist dann:
% Input: \mathbf{A} = [a_{i, j}]_{i, j=1}^n \in \mathbb{R}^{n \times n}
% Output: LR-Zerlegung von \mathbf{A}, gespeichert in \mathbf{A}
% for j = 1, \ldots , n - 1 do
% \qquadsetze \mathbf{A}_{j+1:n,j} := \mathbf{A}_{j+1:n,j} / a_{j,j}
% \qquadsetze \mathbf{A}_{j+1:n,j+1:n} := \mathbf{A}_{j+1:n,j+1:n}-\mathbf{A}_{j+1:n,j}\cdot\mathbf{A}_{j,j+1:n}
% end for
% Schreiben Sie hier am Ende dieses Live Scripts MATLAB-Funktionen
% A = gauss(A),
% A = gaussRank1(A) ,
% welche diese Versionen des Algorithmus' implementieren.
% Testen Sie nun Ihre Implementierung mit dem folgenden Code-Block:

A = toeplitz(1414:-1:1);
B = A;

B = gauss(B);
L = eye(size(B)) + tril(B, -1);
U = triu(B);
disp(num2str(norm(A - L*U, 1)));

B = A;
B = gaussRank1(B);
L = eye(size(B)) + tril(B, -1);
U = triu(B);
disp(num2str(norm(A - L*U, 1)));

% Die berechneten Normen sollen dabei von der Grössenordnung 10^{-10} sein.
%
% Wir möchten nun die Laufzeiten der beiden Algorithmen vergleichen. Erstellen Sie analog zum obigen Code für n = 2^j, j = 4, \ldots, 11 eine Toeplitz-Matrix und messen Sie die Laufzeiten des gewöhnlichen Gauss-Algorithmus', sowie die Lauzeiten des Gauss-Algorithmus' mit Rang-1-Updates. Vergleichen Sie die Laufzeiten in einem loglog-Plot. Plotten Sie zusätzlich die Vergleichsgerade 10^{-8}n^3 gegen n. Beschriften Sie Ihren Plot mit einer Legende und beschriften Sie die Achsen.

% TODO: implementieren

%
%
% Funktionen:

% TODO: implementieren