Newton Interpolation

牛顿插值算法

源代码：

'''
Newton Interpolation
Implement Newton interpolation using Python
Author: Tu ZhiMing
'''

import numpy as np

class Newton:
    '''
    Newton Interpolation Class
    '''
    def __init__(self, x ,y):
        self.x = x
        self.y = y
        self.n = len(x)
        self.chart = np.zeros((self.n,self.n))
        for j in range(self.n):
            for i in range(self.n):
                if j == 0:
                    self.chart[i][j] = self.y[i]
                elif j > i:
                    self.chart[i][j] = 0
                else:
                    self.chart[i][j] = (self.chart[i][j-1] - self.chart[i-1][j-1])/(x[i] - x[i-j])

    def update(self, x, y):

        self.chart = np.row_stack((self.chart, np.zeros(self.n)))
        self.chart = np.column_stack((self.chart, np.zeros(self.n + 1)))
        self.n += 1
        self.x.append(x)
        self.y.append(y)
        j = 1
        for i in range(self.n):
            if i == 0:
                self.chart[self.n - 1][i] = y
            else:
                self.chart[self.n - 1][i] = (self.chart[self.n-1][i-1] - self.chart[self.n-2][i-1])/(self.x[self.n-1]-self.x[self.n-1-j])
                j = j+1
        return self

    def calc(self , x):
        sum = self.chart[0][0]
        for i in range(1,self.n):
            product = self.chart[i][i]
            for j in range(i):
                product *= (x - self.x[j])
            sum += product
        return sum

def main():
    xp = [0.4, 0.55, 0.65, 0.80, 0.90]
    yp = [0.41075, 0.57815, 0.69675, 0.88811, 1.02652]
    xn = 1.05
    yn = 1.25382
    new = Newton(xp , yp)
    x = 0.596
    y = new.calc(x)
    np.set_printoptions(precision=5)
    np.set_printoptions(suppress=True)
    print("均差表为：\n")
    print(new.chart)
    print("\ninterpolated y value = %.6f at x = %.6f"% (new.calc(x),x))

    print("\n")
    print("插入X=%.6f,Y=%.6f 后的均差表为：\n"%(xn,yn))
    new = new.update(xn , yn)
    print(new.chart)

if __name__ == '__main__':
    main()

运行结果：

均差表为：
[[0.41075 0.      0.      0.      0.     ]
 [0.57815 1.116   0.      0.      0.     ]
 [0.69675 1.186   0.28    0.      0.     ]
 [0.88811 1.27573 0.35893 0.19733 0.     ]
 [1.02652 1.3841  0.43347 0.21295 0.03124]]

interpolated y value = 0.631918 at x = 0.596000

插入X=1.050000,Y=1.253820 后的均差表为：
[[0.41075 0.      0.      0.      0.      0.     ]
 [0.57815 1.116   0.      0.      0.      0.     ]
 [0.69675 1.186   0.28    0.      0.      0.     ]
 [0.88811 1.27573 0.35893 0.19733 0.      0.     ]
 [1.02652 1.3841  0.43347 0.21295 0.03124 0.     ]
 [1.25382 1.51533 0.52493 0.22867 0.03143 0.00029]]