牛顿后插公式

2012-2013（1）专业课程实践论文
牛顿后插公式
王瑜，0818180218，R数学08-2班
一、算法理论
由牛顿前插公式
)()!1())...(1()()1(1ξ+++--=n a n n f h n n t t t x R ，),(0n x x ∈ξ
如果要求表示函数在n x 附近的值)(x f ，此时应用Newton 插值公式，插值点应按的次序排列，有
)
()](,,,[))(](,,[)](,[)()(1011211x x x x x x x f x x x x x x x f x x x x f x f x N n n n n n n n n n n n n n --++--+-+=----- 作变换)01(≤≤-+=t th x x n 错误！未找到引用源。

，并利用公式代入上式带入得
n n n n n n n f n n t t t f t t f t f th x N ∇-++++∇++∇+=+!)1()1(!2)1()(2
称为Newton 后插公式，其余项。

。

),(,)!
1()()()1()()()(0)1(1n n n n n n x x n f h n t t t th x N x f x R ∈+++=+-=++ξξ 若用Newton 后插公式求)(x f 的值，因x 在n x 附近，则其系数)(x f 在点n x 的各阶向后差分。

二、算法框图
结束
判断是否
继续输入
提示是否继续输
入
输出结果
判断输入
区间合法性
Input x 提示正确的X 区
间信息
开始
是
否是
否
三、算法程序
class Interpolation
{
public:
Interpolation(int num, double x1, double x2, double func[]);
double ComputeBackwardValue(double x); // compute backward interpolation value ~Interpolation();
private:
// Check(); // checking the inputs
void GetBackwardTable(); // get the backward differential table
private:
int m_num; // the number of interpolation points
double m_x1, m_x2; // the first point m_x1 and last point m_x2
double m_step; // the interpolation step
double* m_func; // the function value of interpolation points
double* m_btable; // the backward differential table
};
#include <iostream>
#include <limits>
using namespace std;
#define NUM 15
//样本个数
#define MIN 4000
//上面输入区间下限
#define MAX 11000
//上面输入区间上限
int main()
{
double func[NUM]=
{
1.38, 1.48, 1.58, 1.69, 1.81,
1.94,
2.10, 2.28, 2.50, 2.76,
3.06, 3.41, 3.83,
4.33, 4.93
};
//输入对应的函数值
double x1=MIN, x2=MAX, x;
int num=NUM;
char flag='Y';
Interpolation test(num, x1, x2, func);
while(flag=='Y')
{
cout<<"Input x: ";
cin>>x;
if (!cin)
{
cin.clear();
cin.ignore(numeric_limits<int>::max(), '\n'); // clear input buffer
continue;
}
if(x<x1 || x>x2)
{
cout<<"---Invalid input: "<<x<<"---"<<endl;
cout<<"Only the number between "<<x1<<" and "<<x2<<" is valid..."<<endl; }
else
{
cout<<"Backward interpolation value:
"<<puteBackwardValue(x)<<endl;
}
cout<<endl<<"Do you want to process? please input(Y/N):"<<endl;
cin>>flag;
}
return 0;
}
Interpolation::Interpolation(int num, double x1, double x2, double func[]) {
m_num = num;
m_x1 = x1;
m_x2 = x2;
m_step = (m_x2-m_x1)/(num-1);
m_func = new double[m_num];
m_btable = new double[m_num];
for (int i=0; i<m_num; ++i)
{
m_func[i] = func[i];
m_btable[i] = func[i];
}
GetBackwardTable();
}
Interpolation::~Interpolation()
{
delete m_func;
delete m_btable;
}
void Interpolation::GetBackwardTable()
{
// get the backward differential table
int i, j;
double tmp;
for (i=1; i<m_num; ++i) {
tmp = m_btable[m_num-1];
for (j=m_num-1; j>=i; --j) {
m_btable[j] = m_btable[j]-m_btable[j-1];
}
m_btable[i-1] = tmp;
}
}
double Interpolation::ComputeBackwardValue(double x)
{
// compute backward interpolation value
double* coef; //coefficient talbe
double result, t;
int i;
coef = new double[m_num];
t = (x-m_x2)/m_step;
for (i=1, coef[0]=1; i<m_num; ++i) //compute the coefficient table {
coef[i] = coef[i-1]*(t+i-1)/i;
}
for (i=0, result=0; i<m_num; ++i)
{
result += m_btable[i]*coef[i];
}
delete coef;
return result;
}
四、算法实现
例1．
在微电机设计计算中需要查磁化曲线表，通常给出的表如图是磁密B每间隔100高斯磁路每厘米长所需安匝数at的值，下面要解决B从4000到11000区间的查表问题
书上给出的可以用后插公式的区间为10500＜B≤11000
解：在此谨以10500以上的数为例，首先输入样本个数和所在区间然后逐一输入对应的函数值。

运行结果：
例2．已知f（x）=sinx+1 求任意插值区间的函数值（0≤x≤5）
1,1.02，1.035，1.052,1.07,1.087区间为0-10应用公式求出结果解：分别输入3,1.5查看结果，运行后输入非法数字4555查看结果。

运行结果：。