#include<stdio.h>
#include<conio.h>
#define MAX 10

int factorial(int temp);

void main()
{
    FILE *fp;
    int number,i,j,k=0;
    float xvalue[MAX],yvalue[MAX],search,product;
    float differ[MAX][MAX],m,uvalue,hvalue,sum;
    fp=fopen("nwfdi.dat","w");
    clrscr();
    printf("\n\n");
    fprintf(fp,"\n\n");
    printf("How many numbers you want to enter for x  : ");
    fprintf(fp,"How many numbers you want to enter for x  : ");
    scanf("%d",&number);
    fprintf(fp,"%d",number);
    for(i=0;i<number;i++)
    {
        printf("\nEnter value for x(%d)  : ",i);
        fprintf(fp,"\nEnter value for x(%d)  : ",i);
        scanf("%f",&xvalue[i]);
        fprintf(fp,"%f",xvalue[i]);
        printf("\nEnter value for y(%d)  : ",i);
        fprintf(fp,"\nEnter value for y(%d)  : ",i);
        scanf("%f",&yvalue[i]);
        fprintf(fp,"%f",yvalue[i]);
    }
    printf("\nEnter any value of x for which you want to find y : ");
    fprintf(fp,"\nEnter any value of x for which you want to find y : ");
    scanf("%f",&search);
    fprintf(fp,"%f",search);
    if(search<xvalue[0] || search>xvalue[number-1])
    {
        printf("\n\nValue lies outside the given values of x ");
        fprintf(fp,"\n\nValue lies outside the given values of x ");
        getch();
        exit(1);
    }
    else
    {
        clrscr();
        printf("\n\nNEWTON'S FORWARD DIFFERENCE INTERPOLATION ");
        fprintf(fp,"\n\nNEWTON'S FORWARD DIFFERENCE INTERPOLATION ");
        for(j=0;j<number-1;j++)
        {
            for(i=0;i<number-(j+1);i++)
            {
                if(j==0)
                {
                    differ[i][j]=yvalue[i+1]-yvalue[i];

                }
                else
                {
                    differ[i][j]=differ[i+1][j-1]-differ[i][j-1];

                }
            }
        }
        printf("\n\n");
        fprintf(fp,"\n\n");
        printf(" x       y   ");
        fprintf(fp," x       y   ");
        for(i=1;i<number;i++)
        {
            printf("   d^%dy(i) ",i);
            fprintf(fp,"   d^%dy(i) ",i);
        }
        printf("\n\n");
        fprintf(fp,"\n\n");
        for(i=0;i<number;i++)
        {
            printf(" %.2f    %.2f  ",xvalue[i],yvalue[i]);
            fprintf(fp," %.2f    %.2f  ",xvalue[i],yvalue[i]);
            for(j=0;j<number-(i+1);j++)
            {
                printf("  %.4f  ",differ[i][j]);
                fprintf(fp,"  %.4f  ",differ[i][j]);
            }
            printf("\n");
            fprintf(fp,"\n");
        }

        for(i=0;i<number;i++)
        {
            if(search>xvalue[i])
            {
                k=k+1;
            }
        }
        k=k-1;
        hvalue=xvalue[1]-xvalue[0];
        uvalue=(search-xvalue[k])/hvalue;
        sum=yvalue[k];
        for(i=0;i<number-(k+1);i++)
        {
            product=1;
            for(j=0;j<=i;j++)
            {
                product=product*(uvalue-j);
            }
            m=factorial(i+1);
            sum=sum+(differ[k][i]*product)/m;

        }
        printf("\n\n");
        fprintf(fp,"\n\n");
        printf("Interpolated value is  :  %f  ",sum);
        fprintf(fp,"Interpolated value is  :  %f  ",sum);
    }
    fclose(fp);
    getch();
}

int factorial(intvalue)
{
    int i,temp=1;
    for(i=value;i>=1;i--)
    {
        temp=temp*i;
    }
    return(temp);
}