%---------------------
% Einführung in die Numerische Mathematik 1
% Alexander Glowatschnig
% Blatt 3 Übungsbeispiel 6
% --------------------


clear
clc

% Definitionen, Notationen und Anfangswerte
N = 3;
N1 = N+1;
N4 = round(N/4);
V0 = [zeros(N4,1);ones(N1-2*N4,1);zeros(N4,1)];
V = V0;
k = 1;
r = 1;
d = 1;
t = 1;
h = 1/N;
M = spdiags([2*ones(N,1);1],0,(N1),(N1)) + spdiags(-1*ones(N1),-1,(N1),(N1)) + spdiags(-1*ones(N1),1,(N1),(N1));
O = spdiags([1;ones(N-1,1);1],0,(N1),(N1)) + spdiags(-1*ones(N1),-1,(N1),(N1));
I = eye(N1);
B = (d/(h^2)) * M + (k/h) * O + (r * I);
A = (eye(N1)+t*B);

%------------------------------------------
% 
% %LR-Zerlegung mit Pivotsuche
LU_Zerlegung(A)
%   
%        [N1,N1] = size(A);
% 
% %b=zeros(1,N1);
% for k=1:N
%     %Suche maximales Element in der Spalte
%     %Nicht nötig, da Matrix Diagonaldominant
%     %for j=k+1:N1  
%     %    if(k+1<=N1)
%     %    if(A(j,k)>A(k,k))
%     %        b=A(k,:);
%     %        A(k,:)=A(j,:);
%     %        A(j,:)=b;
%     %    end
%     %   end
%     %end
%     %LR-Zerlegung
%     if A(k,k)~=0
%     A(k+1:n,k) = A(k+1:n,k)/A(k,k);
%     A(k+1:n,k+1:n)= A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);
%     else
%         error('Pivotelement ist Null oder Sehr klein!!!!!')
%     end
% end
% 
% U = triu(A); 
% L = eye(n,n) + tril(A,-1);   
% end
% 
% %------------------------------------------
% 
% vorwaerts(L, V)
% n = size(L, 1);
% V = zeros(1, n);
% 
% V(1) = V(1);
% for i = 2 : n
%     V(i) = V(i) - sum(z(1:i-1) .* L(i, 1:i-1));
% end
% 
% %------------------------------------------
% 
% rueckwaerts(U,V)
% n = size(U, 1);
% V = zeros(1, n);
% 
% V(n) = V(n) / U(n, n);
% for i = n-1 : -1 : 1
%     V(i) = (V(i) - sum(V(i+1:n) .* U(i, i+1:n))) / U(i,i);
% end

%------------------------------------------

V = V0;
for N1=1:N1max
V = vorwaerts(L,V); % Vorwaerts Substitution mit L
V = rueckwaerts(U,V); % Rueckwaers Substitution mit U
plot(x,V);
drawnow; pause(0.1);
end
