C語言 範例1-13 已知不明分行物體, 經光測儀器之量測,印出Hermite 差除表

C語言 範例1-13 已知不明分行物體, 經光測儀器之量測,印出Hermite 差除表

範例1-13 已知不明分行物體, 經光測儀器之量測其結果如下
=========================================
時間                位置               飛行速度
 x                       f(x)                  f ' (x)
 0                       50.00               50.00
 2                       216.67             94.44
 4                       410.00             98.00
=========================================
印出Hermite 差除表


程式
/* ex1-13.c To generate the coeffficients for
 * Hermite Interpolating Polynomial H on the distinct numbers
 * x0,x1,x2,...with f0,f1,f2,... and f'0,f'1,f'2...
 * based on the Newton backward divided-difference Algorithm
 * and output the divided_difference table for Hn(x).
 */
#include <stdio.h>
int main()
{
   int i,j,k,n;
   double x[30],z[30],q[30][30],f[30],ff[30];
   printf("input n=");
   scanf("n=%d",&n);
   for(i=0;i<=n;i++)
   {
    printf("input x=");
    scanf("%lf",&x[i]);
    printf("input f=");
    scanf("%lf",&f[i]);
    printf("input ff=");
    scanf("%lf",&ff[i]);
    }
   printf("i z(i)   f(i)   f(i-1,i)   f(i-2,i-1,i)...\n");
   for(i=0;i<=n;i++)
   {
    z[2*i]=x[i];
    z[2*i+1]=x[i];
    q[2*i][0]=f[i];
    q[2*i+1][0]=f[i];
    q[2*i+1][1]=ff[i];
    if(i !=0)
q[2*i][1]=(q[2*i][0]-q[2*i-1][0])/(z[2*i]-z[2*i-1]);
      if(i==0)
printf("%d  %3.1lf %7.4lf\n",i,z[i],q[i][0]);
      else if(i==1)
printf("%d  %3.1lf %7.4lf  %7.4lf\n",i,z[i],q[i][0],q[i][1]);
   }
   for(i=2;i<=2*n+1;i++)
   {
    printf("%d  %3.1lf %7.4lf %7.4lf",i,z[i],q[i][0],q[i][1]);
    for(j=2;j<=i;j++)
    {
q[i][j]=(q[i][j-1]-q[i-1][j-1])/(z[i]-z[i-j]);
printf(" %7.4lf",q[i][j]);
    }
    printf("\n");
   }
   return 0;
}


輸入資料
n=2
0.0  50.0    50.0
2.0  216.67  94.44
4.0  410.00  98.00


輸出畫面
input n=input x=input f=input ff=input x=input f=input ff=input x=input f=input ff=

i   z(i)       f(i)      f(i-1,i)     f(i-2,i-1,i)...
0  0.0       50.0000
1  0.0 50.0000     50.0000
2  2.0 216.6700   83.3350   16.6675
3  2.0 216.6700   94.4400     5.5525    -5.5575
4  4.0 410.0000   96.6650     1.1125    -1.1100     1.1119
5  4.0 410.0000   98.0000     0.6675    -0.2225     0.2219   -0.2225

Hn(x)= 50 + 50 (x-0.0) + 16.6675 (x-0.0)^2  + -5.5575  (x-0.0)^2  (x-2.0)
            +  1.1119 (x-0.0)^2  (x-2.0)^2 +  -0.2225 (x-0.0)^2  (x-2.0)^2 (x-4.0)