#include <vector>
#include <iostream>
using namespace std;
void print(vector<double> s)
{
for(double x:s)
{
cout<<x<<" ";
}
cout<<endl;
}
vector<double> mul(vector<double> a,vector<double> b) //多项式乘法
{
int n=a.size(), m=b.size();
vector<double> s(n+m-1,0);
for(int i=0;i<n;++i)
{
for(int j=0;j<m;++j)
{
s[i+j]+=a[i]*b[j];
}
}
return s;
}
vector<double> add(vector<double> a,vector<double> b)
{
int n=a.size(), m=b.size();
vector<double> s(max(n,m),0);
for(int i=0;i<n;++i){
s[i]+=a[i];
}
for(int i=0;i<m;++i){
s[i]+=b[i];
}
return s;
}
void print_chashang_table(vector<vector<double>> A)
{
int n=A.size(), m=A[0].size();
for(int i=0;i<n;++i)
{
for(int j=0;j<=i;++j)
{
printf("% 5.2f ",A[i][j]);
}
printf("\n");
}
}
vector<double> Newton(vector<vector<double>> Point){
int n=Point.size()-1;
if(n==-1){
cout<<"错误"<<endl;
}
vector<vector<double>> diffCoef_Table;
for(int i=0;i<n+1;++i){ //n+1项差商,初始化差商表
vector<double> V(i+1,0);
diffCoef_Table.emplace_back(V);
diffCoef_Table[i][0]=Point[i][1]; //差商表第一列为函数值
}
for(int k=1;k<=n;++k){ //递推求差商表
for(int i=k;i<=n;++i){
diffCoef_Table[i][k]=(diffCoef_Table[i][k-1]-diffCoef_Table[i-1][k-1])/(Point[i][0]-Point[i-k][0]);
}
}
print_chashang_table(diffCoef_Table);
vector<double> ans(n+1,0); //多项式为n次方,加上常数项n+1项
ans[0]=Point[0][1];
for(int k=1;k<=n;++k){
vector<double>s={diffCoef_Table[k][k]};
for(int j=0;j<=k-1;++j){
s=mul(s,{-Point[j][0],1});
}
ans=add(ans,s);
}
return ans;
}
int main(){
vector<vector<double>> P={{0,0},{2,0},{1,1},{5,3}};
vector<double> s=Newton(P);
cout<<endl;
print(s);
system("pause");
return 0;
}
牛顿插值法
于 2022-07-05 23:43:45 首次发布