define MAXTIMES 10 //最高构造10次插值函数
typedef struct {
double x;
double y;
}Point; //插值节点类型
int Times; //实际要构造的次数
double DCList[(MAXTIMES + 1)*(MAXTIMES + 2) / 2]; //差商列表
Point knotList[MAXTIMES + 1];
void initKnotList(); //插值节点数据输入
void initDCList(); //构造差商表
void displayDCList();
double Newton(double u); //Newton 插值函数, 计算点x处的函数值
int main()
{
double u; //计算该点的插值函数值
initKnotList();
initDCList();
printf("%f \n", Newton(u));
while (getchar() != '\n');
getchar();
return 0;
}
void initKnotList()
{
int i;
printf("input Times:");
scanf_s("%d", &Times, sizeof(int));
for (i = 0; i <= Times; i++)
{
scanf_s("%lf", &knotList[i].x, sizeof(double));
scanf_s("%lf", &knotList[i].y, sizeof(double));
}
}
void initDCList()
{
int i, j;
int prePosi, currPosi;
currPosi = 0;
for (i = 0; i <= Times; i++) //第i阶差商
{
for (j = 0; j <= Times - i; ++j) //第i阶的第j个差商
{
if (i == 0) //0阶差商
{
DCList[currPosi + j] = knotList[j].y;
}
else
{
DCList[currPosi + j] =
(DCList[prePosi + j] - DCList[prePosi + j + 1]) /
(knotList[j].x - knotList[j + i].x);
}
}
prePosi = currPosi;
currPosi = (i + 1)*Times + 1;
}
}
void displayDCList()
{
int i, j;
int currPosi;
currPosi = 0;
for (i = 0; i <= Times; i++)
{
for (j = 0; j <= Times - i; ++j)
{
printf("%f ", DCList[currPosi+j]);
}
currPosi = ((i+1)*Times) + 1;
printf("\n");
}
}
double Newton(double u)
{
int i, j;
int DCNumber=0;
double sum = 0;
double w;
for (i = 0; i <= Times; i++)
{
w = 1;
for (j = 0; j < i; ++j)
{
w *= (u - knotList[j].x);
}
sum += DCList[DCNumber] * w;
DCNumber = ((i + 1)*Times + 1);
}
return sum;
}