三次样条插值

来源:互联网 发布:电路模拟软件下载 编辑:程序博客网 时间:2024/04/29 12:09

#include "iostream.h"
void main()
{   
int i,n;
    float **a,x;
float three(float **a,int n,float x);
cout<<"请输入点的个数:";
cin>>n;
a=new float*[n];
for(i=0;i<n;i++)
  a=new float[2];
cout<<"请输入各个点的坐标:"<<endl;
for(i=0;i<n;i++)
  cin>>a[0]>>a[1];
cout<<"输入要求的自变量值:";
cin>>x;
cout<<"方程的解为:"<<three(a,n,x)<<endl;   
}
float three(float **a,int n,float x)
{
void catchwith(float **a,int n);
int i;
float **b,y=0;    [url=file://b/]file://b[/url]保存三对角阵
b=new float*[n-1];
for(i=0;i<n-1;i++)
  b=new float[5];
for(i=0;i<n-1;i++)   [url=file://初/]file://初[/url]始化三对角阵
{
  b[0]=1-(a[0]-a[0])/(a[0]-a[0]);
  b[1]=2;
  b[2]=1-b[0];
  b[3]=3*(b[0]*(a[1]-a[1])/(a[0]-a[0])+
          b[3]*(a[1]-a[1])/(a[0]-a[0]));
}
b[0][0]=2,b[0][1]=b[0][2],b[0][2]=0,b[n-2][2]=0;
b[0][3]=b[0][3]-(1-b[0][1])*(a[1][1]-a[0][1])/(a[1][0]-a[0][0]);
b[n-2][3]=b[n-2][3]-(1-b[n-2][0])*(a[n-1][1]-a[n-2][1])/(a[n-1][0]-a[n-2][0]);
catchwith(b,n-1);   [url=file://求/]file://求[/url]解三对角阵
for(i=0;x>a[0];i++);  [url=file://寻/]file://寻[/url]找x的位置
i--;
y+=((x-a[0])/(a[0]-a[0])-1)*((x-a[0])/(a[0]-a[0])-1)*(2*
  (x-a[0])/(a[0]-a[0])+1)*a[1];
y+=(x-a[0])/(a[0]-a[0])*(x-a[0])/(a[0]-a[0])*(-2*(x-a[0])/
  (a[0]-a[0])+3)*a[1];
y+=(a[0]-a[0])*(x-a[0])/(a[0]-a[0])*((x-a[0])/
    (a[0]-a[0])-1)*((x-a[0])/(a[0]-a[0])-1)*b[4];
y+=(a[0]-a[0])*(x-a[0])/(a[0]-a[0])*(x-a[0])/
    (a[0]-a[0])*((x-a[0])/(a[0]-a[0])-1)*b[4];
return y;
}
void catchwith(float **a,int n)    [url=file://追/]file://追[/url]赶法解三对角阵
{
for(int i;i<n;i++)  [url=file://追/]file://追[/url]
{
  a[3]-=a[3]*a[0]/a[1];
  a[1]-=a[2]*a[0]/a[1];
  a[0]=0;
}a[n-1][4]=a[n-1][3]/a[n-1][1];
for(i=n-2;i>=0;i++)  [url=file://赶/]file://赶[/url]
  a[4]=(a[3]-a[2]*a[4])/a[n-1][1];  
}