three.cpp

计算机数值分析中三次样条插值算法的演示算法

#include <stdio.h>
#include <math.h>
#include <iostream.h>
//#include <cmath.h>

void main()
{
	int n=5;
	double X[5]=
	{
		0.25,0.30,0.39,0.45,0.53		
	};
	double Y[5]=
	{
		0.5000,0.5477,0.6245,0.6708,0.7280
	};
	double P0=0;
	double Pn=0;
	double h[4];
	
	int i,j,k;

	for(i=0;i<n-1;i++)
	{
		h[i]=X[i+1]-X[i];
	}

	double C[5];
	C[0]=1;
	double G[5];
	double A[5];
	G[0]=3*(Y[1]-Y[0])/h[0]-P0*h[0]/2;

	for(k=1;k<n-1;k++)
	{
		A[k]=h[k]/(h[k]-h[k-1]);
		C[k]=1-A[k];
		G[k]=3*(C[k]*(Y[k+1]-Y[k])/h[k]+A[k]*(Y[k]-Y[k-1])/h[k-1]);
	}
	A[4]=1;
	G[4]=3*(Y[4]-Y[3])/h[3]+Pn*h[3]/2;

	double B1[5],B2[5];
	B1[0]=2;
	B2[0]=C[0]/B1[0];

	for(i=1;i<n;i++)
	{
		B1[i]=2-A[i]*B2[i-1];
		B2[i]=C[i]/B1[i];
	}

	G[0]/=2;

	for(j=1;j<n;j++)
	{
		G[j]=(G[j]-A[j]*G[j-1])/B1[j];
	}

	double m[5];
	for(i=0;i<5;i++)
	{
		m[i]=G[i];
	}

	for(k=3;k>=0;k--)
	{
		m[k]=G[k]-B2[k]*m[k+1];
	}

	double u,su;
	printf("输入(0.25-0.53)：\n");
	//scanf("%f",&u);
	cin>>u;
	if(u>=0.25&&u<=0.53)
	{
		for(i=4;i>0;i--)
		{
			if (u<=X[i])
			{
				break;
			}
		}
	}
	else
	{
		printf("input error!!\n");
	}
	i--;
	
	su=0;
	su+=(h[i]+2*(u-X[i]))*(pow((u-X[i+1]),2.0))*(Y[i])/(pow(h[i],3.0));
	su+=(h[i]-2*(u-X[i+1]))*(pow(u-X[i],2.0))*(Y[i+1])/(pow(h[i],3.0));
	su+=(u-X[i])*(pow((u-X[i+1]),2.0))*(m[i])/(pow(h[i],2.0));
	su+=(u-X[i+1])*(pow((u-X[i]),2.0))*(m[i+1])/(pow(h[i],2.0));
	printf("结果：%f\n",su);
	//cout<<su<<endl;

}