/*
      	Title: Test div.diff, given in Section 5.2.
 
      	This will read data points (x(i),y(i)), i=0,...,n, and then
      	it will produce and print the divided differences
            	f[x(0),...,x(i)],  i=0,...,n,
      	where f(x(i)) = y(i).
*/

#include <stdio.h>

#define	MAX_N	20

void divdif(int n,float x[MAX_N+1],float f[MAX_N+1]);

float x[MAX_N+1], f[MAX_N+1];

main()
{
	int n, i;

	while (1) 
	{
		printf("\n\n Give n.    to stop, give  n=0 ");
		printf("\n n should be less than %d  : \n", MAX_N + 1);
		scanf("%d", &n);
		if (n == 0) 
			return(0);

		printf("\n\n Give the nodes and function values x(i), f(i)\n");
		for (i=0; i <= n; i++)
		scanf(" %f %f", &x[i], &f[i]);
      				
		divdif(n,x,f);
		for(i=0; i <= n; i++)
		printf("\n i = %2d   x(i) = %8.4f   df(i) = %14.7e", i, x[i], f[i]);
	}
}


void divdif(int n,float x[MAX_N+1],float f[MAX_N+1]) 
{

/*
	Input: (1) 'n' denotes a positive integer.
		   (2) 'x' denotes a vector of points x(0),...,x(n).
		   (3) 'f' denotes a vector of values of some function
				evaluated at the nodes in x.
	Output: The vector f is converted to a vector of divided
			differences: f(i)=f[x(0),...,x(i)], i=0,...,n
*/

int i, j;

for (i=1; i <= n; i++)
	for (j=n; j >= i; j--)
		f[j] = (f[j] - f[j-1])/(x[j] - x[j-i]);
}