• 设为首页
  • 点击收藏
  • 手机版
    手机扫一扫访问
    迪恩网络手机版
  • 关注官方公众号
    微信扫一扫关注
    公众号

CPP,MATLAB实现牛顿插值 - 不秩稚童

原作者: [db:作者] 来自: [db:来源] 收藏 邀请

    牛顿插值法的原理,在维基百科上不太全面,具体可以参考这篇文章。同样贴出,楼主作为初学者认为好理解的代码。

function p=Newton1(x1,y,x2)
%p为多项式估计出的插值
syms x
n = length(x1);
%差商的求法

for i=2:n
   f1(i,1)=(y(i)-y(i-1))/(x1(i)-x1(i-1));
end

for i=2:n
    for j=i+1:n
       f1(j,i)=(f1(j,i-1)-f1(j-1,i-1))/(x1(j)-x1(j-i));
   end
end
f1=[y\',f1]%  输出带0阶差商的差商表格

%Newton插值函数
Newton=f1(1,1);
for i=2:n
    tt=1;
    for j=1:i-1
        tt=tt*(x-x1(j));
    end
    Newton=Newton+f1(i,i)*tt;
end
fprintf(\'Newton插值函数为\n\')
expand(Newton)     % 将连乘多项式合并展开
x = x2;
p = eval(Newton);   % 代入值计算
fprintf(\'Newton插值函数在所求点x2的函数值为\n\')
p
    

  运行:

  输出:

 

  CPP实现代码如下:

  注意此处求差商运用的是另外一种方法

 

#include<iostream>
#include<string>
#include<vector>
using namespace std;

// 函数声明,使得其在被完整定义之前可以被引用
double ChaShang(int n, vector<double>&X, vector<double>&Y);
double Newton(double x, vector<double>&X, vector<double>&Y);

int main(){
	int n;
	cout<<"输入插值点的个数:"<<endl;
	cin>>n;
	// 先将X,Y填充为n个0
	vector<double>X(n,0);
	vector<double>Y(n,0);
	cout<<"请输入X[i],Y[i]:"<<endl;
	for(int i=0;i<n;i++){
		cin>>X[i]>>Y[i];
	}
	double x;
	cout<<"请输入要进行插值的点的x值:"<<endl;
	cin>>x;
	cout<<Newton(x,X,Y)<<endl;
	return 0;
}

// 此处为差商的另一种求法,可有差商定义根据数学归纳法求出
double ChaShang(int n,vector<double>&X,vector<double>&Y){
	double f=0;
	double temp=0;
	for(int i=0;i<n+1;i++){
		temp=Y[i];
		for(int j=0;j<n+1;j++){
			if(i!=j){temp /= (X[i]-X[j]);}
		}
		f += temp;
	}return f;

}

double Newton(double x,vector<double>&X,vector<double>&Y){
	double result=0;
	for(int i=0;i<X.size();i++){
		//此处的temp用于生成牛顿插值多项式里面的多项式乘积因子,(x-1)(x+3)这些
		double temp=1;
		double f=ChaShang(i,X,Y);
		for(int j=0;j<i;j++){
			temp=temp*(x-X[j]);
		}
		// 差商乘以因子得到最终的牛顿插值多项式
		result += f*temp;
	}return result;
}

  运行结果:

 


鲜花

握手

雷人

路过

鸡蛋
该文章已有0人参与评论

请发表评论

全部评论

专题导读
热门推荐
阅读排行榜

扫描微信二维码

查看手机版网站

随时了解更新最新资讯

139-2527-9053

在线客服(服务时间 9:00~18:00)

在线QQ客服
地址:深圳市南山区西丽大学城创智工业园
电邮:jeky_zhao#qq.com
移动电话:139-2527-9053

Powered by 互联科技 X3.4© 2001-2213 极客世界.|Sitemap