function  A5
load gewicht.dat

x = gewicht(:,1);
y = gewicht(:,2);
hold on;
plot(x,y,'.k');
xx = linspace(x(1),x(end),100);
drawnow;
disp('Bitte warten...');

% Lagrange-Interpolation
tic; flops(0);
ff = zeros(size(xx));
for i=1:length(x)
   ff = ff + y(i)*psi(i,x,xx);
end;
fprintf('Zeit (Lagrange) : %0.2fs    Flops: %i\n',toc, flops);

plot(xx,ff,'k');
grid on;
drawnow;

% Neville-Aitken-Algorithmus
tic; flops(0);
P = Neville_Aitken(x,y,xx);
fprintf('Zeit (Nev.-Ait.): %0.2fs    Flops: %i\n',toc, flops);
plot(xx,ff,'r');
drawnow;

% Monom-Ansatz
tic; flops(0);
X = vander(x);
a = X \ y;
plot(xx,polyval(a,xx),'b');
fprintf('Zeit (Monom)    : %0.2fs    Flops: %i\n',toc, flops);


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


function f = Neville_Aitken(x,y,xj)
% Neville-Aitken-Algorithmus für die Stützstellen x
% und den Funktionswerten y. Dieser Algorithmus wertet
% das Interpolationspolynom durch die Punkte (x,y)
% an den Stellen xj aus.
%
% Eingabe: x  Stützstellen
%          y  Funktionswerte, die interpoliert werden
%          xj Stellen, an denen das Polynom ausgewertet wird
%
% Ausgabe: f  = P(xj), interpolierte Werte

n = length(x);
f = zeros(size(xj));

for j=1:length(xj)
   P = y;
   
   for k=2:n
      for i=n:-1:k
         P(i) = P(i) + (xj(j)-x(i)) / (x(i)-x(i-k+1)) * (P(i)-P(i-1));
      end;
   end;
   
   f(j) = P(n);
end;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


function f = psi(i,x,xj);
% Werte die Funktion psi_i für die Stützstellen x
% an den Stellen xj aus. ACHTUNG: der erste index i
% ist 1, nicht 0 (im Gegensatz zur Vorlesung)!
%
% Eingabe: i  1 <= i <= n=length(xj)
%          x  Stützstellen
%          xj Stellen, an denen psi_i ausgewertet wird
%
% Ausgabe: f  = psi_i(xj),

f = zeros(size(xj));

for j=1:length(xj)
   f(j) = prod((xj(j)-x([1:i-1 i+1:length(x)])) ./ (x(i)-x([1:i-1 i+1:length(x)])));
end;


