#include "NewtonDif.h"
#include<math.h>

NewtonDif::NewtonDif(void)
{
}
bool NewtonDif::ForwardDif(int n,double temX[],double temY[])
{
	if(n>1000||n<1)return false;
	num=n;
	int i=0;
	double *F[1000];
	for(i=0;i<=num;i++)
	{
		F[i]=new double[num];
	}
	for(i=0;i<=num;i++)
	{
		x[i]=temX[i];
		F[i][0]=temY[i];
	}
	f[0]=F[0][0];
	for(i=1;i<=num;i++)
	{
		for(int j=1;j<=i;j++)
		{
			F[i][j]=(F[i][j-1]-F[i-1][j-1])/(x[i]-x[i-j]);
		}
		f[i]=F[i][i];
	}
	return true;
}
bool NewtonDif::BackwardDif(int n, double temX[], double temY[])
{
	if(n>1000||n<1)return false;
	num=n;
	int i=0;
	double *F[1000];
	for(i=0;i<=num;i++)
	{
		F[i]=new double[num];
	}
	for(i=0;i<=num;i++)
	{
		x[num-i]=temX[i];
		F[num-i][0]=temY[i];
	}
	f[0]=F[0][0];
	for(i=1;i<=num;i++)
	{
		for(int j=1;j<=i;j++)
		{
			F[i][j]=(F[i][j-1]-F[i-1][j-1])/(x[i]-x[i-j]);
		}
		f[i]=F[i][i];
	}
	return true;
}
double NewtonDif::calA(double s,double k)
{
	double ans=1.0;
	for(int i=0;i<k;i++)
	{
		ans=ans*(s-i);
	}
	return ans;
}
double NewtonDif::GetValue(double oneX)
{
	double h=fabs(x[1]-x[0]);
	double s=(oneX-x[0])/h;
	double Ask=0.0;
	double result=f[0];
	for(int k=1;k<=num;k++)
	{
		Ask=calA(s,k);
		result=result+Ask*pow(h,(double)k)*f[k];
	}
	return result;
}
double *NewtonDif::GetPara()
{
	return f;
}
NewtonDif::~NewtonDif(void)
{
}