这篇接着上一篇,对牛顿基本插值多项式的理解与代码实现
在拉格朗日插值法中,
一般插值公式为:
其中的每次计算都依赖于全部的插值结点,在增加或减少结点时,必须全部重新计算。(不是特别理解)为了克服这个缺点从而引入了牛顿插值法。
牛顿基本插值多项式的一般公式为:
其中定义:
为零阶差商
为一阶差商
为二阶差商
以此类推。。。
另有公式:
在编程实现的过程中就需要考虑,怎么计算,
可分为定义法和公式法。
- 定义法:基本思路便是先算出一阶差商,然后算出二阶差商,。。。
这样设计必然需要设定一个数组用于存储中间。
- 公式法:利用公式直接计算得到
这样设计的后果很明显,在公式中会有大量乘法运算,速度会比较慢。
下面是代码阶段:
在代码中用quotient_array用于存储所有差商,
quotient_list用于存储从x0开始的差商。
#include<iostream>
#define maxsize 100
using namespace std;
double quotient_array[maxsize][maxsize];
void D_creat_quotient_list(int n,double x[maxsize],double y[maxsize],double quotient_list[maxsize])
{//构造差商表
//定义法 占空间但是省时间
for(int k=0;k<n;k++)
{
for(int i=0;i<n-k;i++)
{
if(k==0)
{
quotient_array[i][0]=y[i];
}
else
{
quotient_array[i][k]=(quotient_array[i+1][k-1]-quotient_array[i][k-1])/(x[i+k]-x[i]);
}
}
quotient_list[k]=quotient_array[0][k];
}
}
void F_creat_quotient_list(int n,double x[maxsize],double y[maxsize],double quotient_list[maxsize])
{//构造差商表
//公式法 省空间 但是费时间
quotient_list[0]=y[0];
double temp;
for(int k=1;k<n;k++)
{
quotient_list[k]=0;
for(int j=0;j<=k;j++)
{
temp = 1;
for(int i=0;i<=k;i++)
{
if(i!=j)
{
temp*=(x[j]-x[i]);
}
}
temp = y[j]/temp;
quotient_list[k]+=temp ;
}
}
}
double F_interpolitation(int n,double x[maxsize],double y[maxsize],double vx,double result,double quotient_list[maxsize])
{//差商牛顿插值法
//两种方法二选一
D_creat_quotient_list(n,x,y,quotient_list);//定义法
//F_creat_quotient_list(n,x,y,quotient_list);//公式法
double temp;
result = quotient_list[0];
for(int k=1;k<=n;k++)
{
temp=1;
for(int i=0;i<k;i++)
{
temp*=(vx-x[i]);
}
result+=quotient_list[k]*temp;
}
return result;
}
int main()
{
int n;
double x[maxsize],y[maxsize],result=0,vx,temp;
double quotient_list[maxsize];
cin>>n>>vx;
for(int i=0;i<n;i++)
{
cin >> x[i] >> y[i];
}
result=F_interpolitation(n,x,y,vx,result,quotient_list);
cout<<result<<endl;
return 0;
}
测试样例:
参考:
x | 100 | 121 | 144 | 169 |
![]() | 10 | 11 | 12 | 13 |
分别用定义法和公式法求的值
定义法:
4
115
100 10
121 11
144 12
169 13
公式法:
4
115
100 10
121 11
144 12
169 13
对比也可以看出速度上是有差距的,一般都用定义法求值。
其中两个函数得到的quotient_list的值相同。
quotient_list | 10 | 0.047619 | -0.00009411 | 0.0000003138 |