C語言 習題1-5已知下列諸點的座標值求牛頓的向前差除表

C語言 習題1-5已知下列諸點的座標值

 xi       f(xi)
==============
0.50   0.6915
0.60   0.7257
0.65   0.7422
0.75   0.7734
0.90   0.8159
1.10   0.8643
1.30   0.9032
===============
求牛頓的向前差除表

計算
1). P(0.55) 與 f(0.55)=0.7088 之誤差
2). P(0.70) 與 f(0.70)=0.7580 之誤差
3). P(0.80) 與 f(0.80)=0.7881 之誤差
4). P(1.00) 與 f(1.00)=0.8413 之誤差

// CPP program for implementing 
// Newton divided difference formula 
#include <stdio.h>

// Function to find the product term 
double proterm(int i, float value, float x[]) 
{ 
    double pro = 1; 
    for (int j = 0; j < i; j++) { 
        pro = pro * (value - x[j]); 
    } 
    return pro; 
} 
  
// Function for calculating 
// divided difference table 
void dividedDiffTable(float x[], float y[][10], int n) 
{ 
    for (int i = 1; i < n; i++) { 
        for (int j = 0; j < n - i; j++) { 
            y[j][i] = (y[j][i - 1] - y[j + 1] 
                         [i - 1]) / (x[j] - x[i + j]); 
        } 
    } 
} 
  
// Function for applying Newton's 
// divided difference formula 
double applyFormula(float value, float x[], 
                   float y[][10], int n) 
{ 
    double sum = y[0][0]; 
  
    for (int i = 1; i < n; i++) { 
      sum = sum + (proterm(i, value, x) * y[0][i]); 
    } 
    return sum; 
} 
  
// Function for displaying  
// divided difference table 
void printDiffTable(float y[][10],int n) 
{ 
    for (int i = 0; i < n; i++) { 
        for (int j = 0; j < n - i; j++) { 
            printf("    %0.4lf\t", y[i][j]); 
        } 
        printf("\n"); 
    } 
} 
  
// Driver Function 
int main() 
{ 
    // number of inputs given 
    int n=7 ,n1=4 ; 
    float value, sum, y[10][10]; 
    float x[] = { 0.50 , 0.60 , 0.65 , 0.75 ,  0.90 , 1.10 , 1.30 }; 
    float xa[] = { 0.55 , 0.70 , 0.80 , 1.00 }; 

    // y[][] is used for divided difference 
    // table where y[][0] is used for input 
    // 0.6915 ,  0.7257   , 0.7422  ,0.7734  , 0.8159 , 0.8643 , 0.9032
    y[0][0] = 0.6915; 
    y[1][0] = 0.7257; 
    y[2][0] = 0.7422; 
    y[3][0] = 0.7734;
    y[4][0] = 0.8159;
    y[5][0] = 0.8643;
    y[6][0] = 0.9032;
    
    

    // calculating divided difference table 
    dividedDiffTable(x, y, n); 
  
    // displaying divided difference table 
    printDiffTable(y,n); 
  
    // value to be interpolated 

    for (int m = 0; m < n1; m++)
    {  
        value = xa[m]; 
        double result = applyFormula(value, x, y, n);
        // printing the value 
        printf("\nValue at %.2lf  is %.5lf \n",value , result); 
    }
    return 0; 
}

輸出畫面

    0.6915     0.3420     -0.0800     -0.1600     0.4445     -0.8256     1.1404 
    0.7257     0.3300     -0.1200     0.0178     -0.0508     0.0868 
    0.7422     0.3120     -0.1147     -0.0076     0.0099 
    0.7734     0.2833     -0.1181     -0.0012 
    0.8159     0.2420     -0.1188 
    0.8643     0.1945 
    0.9032 

Value at 0.55  is 0.70871 

Value at 0.70  is 0.75810 

Value at 0.80  is 0.78811 

Value at 1.00  is 0.84143