Algorithm for the Newton Form of the Interpolating Polynomial

Given the Lagrange form or the Newton form one can combine the terms to find the coefficients of the polynomial. With the Lagrange form I found it too difficult to develop an algorithm for this for the general case. However, with the Newton form I found a compact and effective algorithm. See below.

   let n = the number of points - 1
   let array x be initialized with the x values of the given points
   let array c be initialized with the y values of the given points

Construct the backward divided differences table. The following is from Algorithm 4.3 of Conte and de Boor's Elementary Numerical Analysis, second edition (page 204). "Calculation of the coefficients for the Newton form."

   for(j=1; j<=n; j++) {
     for(i=0; i<=n-j; i++) {
       if(x[i+j]==x[i]) {
         fprintf(stderr,"\nError *** duplicate x values encountered\n\n");
         return;
       }
       c[i]=(c[i+1]-c[i])/(x[i+j]-x[i]);
     }
   }

Now f(a) can be determined for a given value a as follows.

   fa=c[0];
   for(i=1; i<=n; i++) {
     fa=fa*(a-x[i])+c[i];
   }

The following is Jeff Stetekluh's method for combining the terms to form the coefficients of the polynomial.

   for(j=0; j<n; j++) {
     for(i=1; i<=n-j; i++) {
       c[i]=c[i]-x[i+j]*c[i-1];
     }
   }

Array c now holds the coefficients of the polynomial. Its first entry c[0] is the coefficient of the term of highest degree and c[n] is the constant term. Therefore f(a) can be determined for a given value a as follows.

   fa=c[0];
   for(i=1; i<=n; i++) {
     fa=fa*a+c[i];
   }

(email:

 
 Jeff Stetekluh

)