埃尔米特插值.cpp

利用C语言编写的一些常见的差值程序

#include<stdio.h>
#define MAX_N 20
typedef struct tagPOINT
{double x;
 double y;
 double z;
}POINT;
 int main()
{int n;
 int i,j,k;
   POINT points[MAX_N+1];
   double diff[2*MAX_N+2];  
   double m[2*MAX_N+2];
   double h[2*MAX_N+2];
   double x,H=0;
   printf("\nInput n value:");
   scanf("%d",&n);
   if(n>MAX_N)
   {
	   printf("The input %d is larger than %d,please redefine the %d.\n",n,MAX_N,MAX_N);
	   return 1;
   }
   if(n<=0)
   {
	   printf("Please input a number between 1 and %d.\n",MAX_N);
	   return 1;
   }
   printf("Now input the (x_i,y_i,z_i)\n");
   for(i=0;i<=n;i++)
	   scanf("%lf%lf%lf",&points[i].x,&points[i].y,&points[i].z);
   printf("Now input the x value:");
   scanf("%lf",&x);
   for(k=0;k<=n;k++)
	   m[2*k]=m[2*k+1]=points[k].x;
   diff[0]=points[0].y;
   for(i=1;i<=n;i++)
   diff[2*i]=(points[i].y-points[i-1].y)/(points[i].x-points[i-1].x);
   for(i=0;i<=n;i++)
   diff[2*i+1]=points[i].z;
   for(i=1;i<2*n+1;i++)
   { 	  
   for(j=2*n+1;j>i;j--)
   { 
	   diff[j]=(diff[j]-diff[j-1])/(m[j]-m[j-1-i]);
   }printf("%lf\n",diff[i+1]);
   }
  h[0]=0;
  h[1]=x-m[0];
			
   for(k=2;k<=2*n+1;k++)
	{
	 h[k]=h[k-1]*(x-m[k-1]);
	}
   H=points[0].y;
   for(k=1;k<=2*n+1;k++)
	   H=H+diff[k]*h[k];
   printf("H(%f)=%f\n",x,H);
   return 1;
 }