Newton's Backward Interpolation

#include<stdio.h>
#include<malloc.h>
#include<stdlib.h>
void differenceTable(double **A,double **L,int ROW);
double result(double **A,double **L,double x,int ROW);
double factorial(int n);
FILE *fp;
int main(int argc, char *argv[])
{
double **A,*b,*xn,*xnplus1,*p,temp,x,**L,**U,lambda,error,t1=0.0;
int ROW,i,j,k,t=1,m=1,bo=1;
char ch;
if(argc==1)
{
printf("\n Enter the number of values:");
scanf("%d",&ROW);
printf("\n Enter the value at which value has to be aprroximated:");
scanf("%lf",&x);
A=(double **)malloc((ROW+1)*sizeof(double*));
L=(double **)malloc((ROW+1)*sizeof(double*));
U=(double **)malloc((ROW+1)*sizeof(double*));
for(i=0;i<=ROW;i++)
{
A[i]=(double *)malloc((ROW+1)*sizeof(double));
L[i]=(double *)malloc((ROW+1)*sizeof(double));
U[i]=(double *)malloc((ROW+1)*sizeof(double));
}
b=(double *)malloc((ROW+1)*sizeof(double));
p=(double *)malloc((ROW+1)*sizeof(double));
xn=(double *)malloc((ROW+1)*sizeof(double));
xnplus1=(double *)malloc((ROW+1)*sizeof(double));
for(i=1;i<=ROW;i++)
{
printf("\n Enter the value of X[%d]:",i);
scanf("%lf",&A[i][1]);
}
for(i=1;i<=ROW;i++)
{
printf("\n Enter the value of Y[%d]:",i);
scanf("%lf",&A[i][2]);
}

}
else if (argc==2)
{
fp=fopen(argv[1],"r");
if(fp==NULL)
{
printf("\n File not found, program will terminate:");
exit(0);
}
fscanf(fp,"%d",&ROW);
fscanf(fp,"%lf",&x);
A=(double **)malloc((ROW+1)*sizeof(double*));
L=(double **)malloc((ROW+1)*sizeof(double*));
U=(double **)malloc((ROW+1)*sizeof(double*));
for(i=0;i<=ROW;i++)
{
A[i]=(double *)malloc((ROW+1)*sizeof(double));
L[i]=(double *)malloc((ROW+1)*sizeof(double));
U[i]=(double *)malloc((ROW+1)*sizeof(double));
}
b=(double *)malloc((ROW+1)*sizeof(double));
p=(double *)malloc((ROW+1)*sizeof(double));
xn=(double *)malloc((ROW+1)*sizeof(double));
xnplus1=(double *)malloc((ROW+1)*sizeof(double));

while(!feof(fp))
{
for(i=1;i<=ROW;i++)
{
for(j=1;j<=2;j++)
{
fscanf(fp,"%lf ",&A[i][j]);
}

}
}
fclose(fp);
}
else
{
printf("\n Invalid Arguments,program will terminate:\n");
exit(0);
}
printf("\n         X          Y");
printf("\n------------------------------------------------------------------\n");
for(i=1;i<=ROW;i++)
{
for(j=1;j<=2;j++)
{

U[i][j]=A[i][j];
printf("   %+lf",U[i][j]);
}
printf("\n");
}
for(i=1;i<=ROW-1;i++)
{
for(j=ROW-i,t=1;j>=1;j--,t++)
{
L[t][i]=U[t+1][2]-U[t][2];
}
for(t=1;t<=ROW-i;t++)
U[t][2]=L[t][i];
}
differenceTable(A,L,ROW);
printf("\n Value at %+lf by Newton's Backward Interpolation formula is %+lf\n",x,result(A,L,x,ROW));
return 0;
}
double result(double **A,double **L,double x,int ROW)
{
int i,j,flag=1,t;

double h,value=0,p,tmp=1.0;
h=A[2][1]-A[1][1];
p=(A[ROW][1]-x)/h;
for(i=0,t=ROW+1;i<=ROW;i++,t--)
{
if(i==0)
value +=A[ROW][2];
else
{
for(j=1;j<=i;j++)
tmp=tmp*(p-j+1)*(-1);
tmp = tmp/factorial(i);
tmp = tmp*L[t-1][i];
value +=tmp;
}
tmp=1.0;
}
return value;
}
double factorial(int n)
{
 if(n<=1)
  return 1;
 else
  return n*factorial(n-1);
}
void differenceTable(double **A,double **L,int ROW)
{
int **position,*ap,j,m,i;
position=(int **)malloc((ROW+1)*sizeof(int*));
ap=(int *)malloc((ROW+1)*sizeof(int));
for(i=0;i<=ROW;i++)
position[i]=(int *)malloc((ROW+1)*sizeof(int));
int an,tmp;
tmp=ROW;
for(i=1;i<=ROW;i++)
{
an=1+(i-1)*2;
ap[i]=an;
for(j=1;j<=tmp;j++)
{
position[i][j]=an+(j-1)*4;
}
tmp--;
}
tmp=ROW;
tmp=ROW+(ROW-1)*3;
int *match,tmp1,l,*pos,flag,z,k;
match=(int *)malloc((ROW+1)*sizeof(int));
pos=(int *)malloc((ROW+1)*sizeof(int));
for(i=0;i<=ROW;i++)
{
match[i]=0;
pos[i]=0;
}

printf("\n-----------------------------Difference Table-----------------------------\n");
printf("\n X       Y");
for(i=1;i<=ROW-1;i++)
printf("      D_%d",i);
printf("\n--------------------------------------------------------------------------\n");
for(i=1;i<=tmp;i++)
{
tmp1=ROW;

for(l=1;l<=ROW;l++)
{
flag=0;
for(m=1;m<=tmp1;m++)
{
if(i==position[l][m])
{
flag=1;
match[l]=1;
pos[m]=position[l][m];
break;
}
} //inner for
if(flag==1)
{
for(k=1;k<=ROW;k++)
{
if(match[k]==0)
{
printf("");
}
else
{
if(k==1)
//printf("(%d,%d)+(%d,%d)|%d",(position[l][m]/4)+1,k,(position[l][m]/4)+1,k+1,i);
printf("%+2.2lf\t%+2.2lf",A[(position[l][m]/4)+1][k],A[(position[l][m]/4)+1][k+1]);
else
{
z=(position[l][m]-ap[k])/4+1;
//printf("\t\t(%d,%d)|%d",z,k-1,i);
printf("\t\t%+2.2lf",L[z][k-1]);
}

}
} //end of k-loop

} //end of flag==1

else
{printf(" ");}
tmp1--;
for(k=0;k<=ROW;k++)
{
pos[k]=0;
match[k]=0;
}
}

printf("\n");
}
printf("\n--------------------------------------------------------------------------\n");
}