Eine Funktion f(x) hat an den Stützstellen xi (i=1,...,n) die Stützwerte yi (i=1,...,n). Der Wert von f(x) soll an einer Stelle z durch Lagrange-Polynome hinreichend genau interpoliert werden. Die Implementierung wird mit dem Neville-Schema vorgenommen.

• // Lagrange-Interpolation
•  
• x = input(" Stuetzstellen eingeben: ");
• P = input(" Stuetzwerte eingeben: ");
• z = input(" Gesuchte Stelle eingeben: ");
•  
• n = length(x); // Anzahl Stuetzstellen
•  
• i=1;
• while %t;
• for j=i-1:-1:1
• // Berechnung mittels Neville-Schema
• P(j) = ((z-x(j))*P(j+1) - (z-x(i))*P(j)) / (x(i)-x(j));
• end
•  
• if(i==n) then
• // Ergebnis bei Interpolation mit n-Stuetzstellen
• sprintf("P(%g)= %g",z,P(1))
•  
• weiter = input(" Weiter Stuetzpunkte? (1=Ja 0=Nein) ");
• if(weiter==0) then, break, end
•  
• x2 = input(" Stuetzstellen eingeben: ");
• P2 = input(" Stuetzwerte eingeben: ");
• n2 = length(x2);
• x(n+1:n+n2) = x2;
• P(n+1:n+n2) = P2;
• n=n+n2;
• end
• i=i+1;
• end

Ich habe den Algorithmus so programmiert. dass es problemlos möglich ist, weitere Stützstellen hinzuzufügen. Dazu muss der Algorithmus nur eine Zeile im Nevillschema neu berechnen, d.h. nicht alle Spalten. Weitere Infos zu dem Algorithmus von Neville und Aitken findet man auf Wikipedia.