function p = aufg01f04(matnr);

%AUFG01F04 berechnet ein Matrikelnummern interpolierendes Polynom.
%
%   P = AUFG01F04(MATNR) berechnet mittels der Newtonschen Form
%   der Polynominterpolation das eindeutige Interpolationspolynom zu
%   der Aufgabe: Berechne p, so dass
%
%              p(i) = i-te Ziffer der Matrikelnummer
%
%   gilt. Dabei ist die Laenge der Matrikelnummer nicht a priori
%   beschraenkt. Implizit wird die Minimallaenge 5 verwendet. Das
%   bedeutet, dass eine kuerzere Matrikelnummer fuehrende Nullen
%   besitzt, also MATNR = 911 => MATNR = 00911 .
%
%   Der ausgegebene Vektor p ist der Vektor der Koeffizienten des
%   Polynoms p. Da dieser das Polynom eindeutig identifiziert,
%   verwenden wir einfach den selben Buchstaben. Die Koeffizienten
%   sind dieselben wie im Skript.
%
%   AUFG01F04 ist eine stabile Loesung der Aufgabe 3a) von
%   Uebungsblatt 1 der Vorlesung 'Numerische Verfahren' von
%   Jens Zemke. Die Funktion arbeitet unter Matlab und Octave.
%
%   SEE ALSO: aufg01f03, aufg01f05, polyval, polyfit

%   $Author: zemke $
%   $Revision: 1.6 $
%   $Date: 2004/04/15 10:37:22 $

if matnr < 0
   error('Es wurde eine negative Matrikelnummer eingegeben');
end

if matnr ~= floor(matnr)
   error('Es wurde eine nicht ganzzahlige Matrikelnummer eingegeben');
end

%
% Berechne die Laenge der Matrikelnummer, i.e. die Anzahl der Ziffern
% und transformiere die ganze Zahl in eine Zifferndarstellung.
%

n = max(floor(log10(matnr))+1,5);
x = 1:n; y = zeros(1,n);
potenz = 10.^(0:n-1);

for i = 1:n
   y(i) = floor(matnr/potenz(n-i+1));
   matnr = matnr - y(i)*potenz(n-i+1);
end

%
% Berechnung der dividierten Differenzen gemaess Algorithmus 2.15
% des Skriptes 'Grundlagen der Numerischen Mathematik' von Prof. Voss
%

for k = 1:n
   h(k) = y(k);
   if k > 1
      for i = k-1:-1:1
         h(i) = (h(i+1)-h(i))/(x(k)-x(i));
      end
   end
   p(k) = h(1);
end

p = p(:)';

%
% Generieren eines Vektors von Werten in [1,n] und berechnen der
% zugehoerigen Funktionswerte, anschliessend plotten.
%

t = linspace(1,n,100);

pt=p(n)*ones(1,100);

for i = n-1:-1:1
   pt = pt.*(t-x(i))+p(i);
end

plot(t,pt,'g-');
grid('on');
title(['Das Polynom ausgewertet auf dem Intervall [1,',num2str(n),']']);
