Talk:Thiele's interpolation formula

The article on Thiele's interpolation formula duplicates the form given in Abramowitz and Stegun's "Handbook of Mathematical Functions" - as does Eric Weisstein's. However, it is not particularly easy to implement, from the description given, even by writing recursively. In the (antique) version of "Maple" to which I have access, there is a good non-recursive working implementation of Thiele's algorithm. (Dognose (as Anthony Burgess wrote it in one novel), I've used it to good effect often enough.) If anyone is interested, I can pass on the code (copyright abuse permitting).

Hair Commodore 20:58, 2 November 2006 (UTC)

Here is a version of Thiele interpolation in Algol 68:

¢ The MODE of lx and ly here should really be a UNION of "something REAL" and "something SYMBOLIC" ... ¢
PROC thiele:=([]REAL lx,ly, REAL x) REAL:
BEGIN
  []REAL xx=lx[@1],yy=ly[@1];
  INT n=UPB xx;
  IF UPB yy=n THEN
¢ Assuming that the values of xx are distinct ... ¢
    [0:n-1,1:n]REAL p;
    p[0,]:=yy[];
    FOR i TO n-1 DO p[1,i]:=(xx[i]-xx[1+i])/(p[0,i]-p[0,1+i]) OD;
    FOR i FROM 2 TO n-1 DO
      FOR j TO n-i DO
        p[i,j]:=(xx[j]-xx[j+i])/(p[i-1,j]-p[i-1,j+1])+p[i-2,j+1]
      OD
    OD;
    REAL a:=0;
    FOR i FROM n-1 BY -1 TO 2 DO a:=(x-xx[i])/(p[i,1]-p[i-2,1]+a) OD;
    y[1]+(x-xx[1])/(p[1,1]+a)
  ELSE
    error ¢ Unequal length arrays supplied ¢
  FI
END;

Note that, although it works in most cases, it is sensitive to input values, especially those due to equally spaced abscissæ. (Essentially, in such a case, it reduces to the ratio of two polynomials, which may have factors in common - thus yielding a 0/0 form.)

Comments? Hair Commodore 18:00, 1 February 2007 (UTC)

Algol 68???!? Really?

Arghman (talk) 17:43, 27 June 2015 (UTC)

Algol -> Simula

Algol -> Ada

Algol -> Pascal -> C#

Not too bad :) 2001:4643:EBFE:0:A11B:257F:77A1:A7F (talk) 10:05, 23 March 2016 (UTC)