% Diese Funktion splittet eine gegebenen Matrix A in
% eine linke untere Dreiecksmatrix L, die Diagonaleintraege D
% und eine rechte obere Dreiecksmatrix U.

function [L,D,U] = ldu(A)

  n = length(A);
  D = A .* eye(n);
  L = zeros(n,n);
  U = L;

  for j=1:n
    L(j+1:n,j)=A(j+1:n,j);
    U(1:j-1,j)=A(1:j-1,j);
  end
