function p = aufg01f03(matnr);

%AUFG01F03 berechnet ein Matrikelnummern interpolierendes Polynom.
%
%   P = AUFG01F03(MATNR) berechnet mittels Loesung eines Vandermonde
%   Gleichungssystemes 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 in derselben Reihenfolge wie auch POLYVAL und POLYFIT sie
%   verwenden, naemlich absteigend sortiert.
%
%   AUFG01F03 ist die Loesung der Aufgabe 3a) von Uebungsblatt 1
%   der Vorlesung 'Numerische Verfahren' von Jens Zemke. Die
%   Funktion arbeitet unter Matlab und Octave.
%
%   SEE ALSO: aufg01f04, aufg01f05, polyval, polyfit

%   $Author: zemke $
%   $Revision: 1.9 $
%   $Date: 2004/04/15 10:29:21 $

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

%
% Konstruiere eine Vandermonde-Matrix. In dieser Matrix steht
% in der i-ten Zeile und j-ten Spalte das Element x(i)^(j-1).
%

V = (ones(n,1)*x).^((0:n-1)'*ones(1,n));

%
% Jetzt wird der Vektor der Koeffizienten gesucht. Wenn das
% Polynom durch die Summe
%
%    p(x) = a_0 + a_1 x + a_2 x^2 + ...
%
% gegeben ist, dann gilt mit p = (a_0,a_1,...) natuerlich
% auch V*p = y. Allerdings sind in unserem Programm p und y
% Zeilenvektoren, also gilt V*(p') = y';
%

p = y/V;

%
% Dieser Vektor der Koeffizienten ist in umgekehrter Reihenfolge
% verglichen mit denen in Matlab/Octave-Notation !! Und da man
% sich ja bekanntlich an etablierte Standards haelt ;-)
%
% Der Grund fuer diese Reihenfolge ist uebrigens die Auswertung
% mittels Horner-Schema. Allerdings ware es nicht viel schwerer,
% diese Auswertung anders zu implementieren ...
%

p = fliplr(p);

%
% Generieren eines Vektors von Werten in [1,n] und berechnen der
% zugehoerigen Funktionswerte, anschliessend plotten.
% Das ginge natuerlich seeeehr einfach mit polyval, naemlich mit
%
%   plot(t,polyval(p,t));
%
% Aber das war ja nicht erlaubt :-(
%
% Allerdings ist die nun folgende Implementation unter Zuhilfenahme
% des Horner-Schemas dieselbe, die auch polyval intern verwendet ;-)
%

t = linspace(1,n,100);

pt = p(1)*ones(1,100);
for i = 2:n
   pt = p(i)+t.*pt;
end

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