插值方法 - Newton向前向后等距插值
通常我们在求插值节点的开头部分插值点附近函数值时,使用Newton前插公式;求插值节点的末尾部分插值点附近函数值时,使用Newton后插公式。
代码:
1 # -*- coding: utf-8 -*-
2 """
3 Created on Wed Mar 25 15:43:42 2020
4
5 @author: 35035
6 """
7
8
9 import numpy as np
10
11 # 等距节点的Newton向前插值(输入的x向量和y向量注意保证从小到大顺序)
12 def Newton_iplt_forword(x, y, xi):
13 """x,y是插值节点,xi是一个值"""
14 n = len(x)
15 m = len(y)
16 if n != m:
17 print("Error!")
18 return None
19 h = x[1] - x[0]
20 t = (xi - x[0]) / h
21
22 # 先计算差分表(cf)
23 cf = []
24 temp = y.copy()
25 for i in range(n):
26 if i != 0:
27 iv_1 = temp[i - 1]
28 for j in range(i, n):
29 iv_2 = temp[j]
30 temp[j] = iv_2 - iv_1
31 iv_1 = iv_2
32 cf.append(temp[i])
33 # 再计算Newton插值
34 ans = 0
35 for i in range(n):
36 w = 1
37 # 计算多项式部分,让差分作为其系数
38 for j in range(i):
39 w *= ((t - j) / (j + 1))
40 ans += w*cf[i]
41 return ans
42
43 # 等距节点的Newton向后插值(输入的x向量和y向量保证从小到大顺序)
44 def Newton_iplt_backword(x, y, xi):
45 """x,y是插值节点ndarray,xi是一个值"""
46 n = len(x)
47 m = len(y)
48 if n != m:
49 print("Error!")
50 return None
51 h = x[1] - x[0]
52 t = (xi - x[n - 1]) / h
53
54 # 先计算差分表(cf)
55 cf = []
56 temp = y.copy()
57 for i in range(n):
58 if i != 0:
59 iv_1 = temp[i - 1]
60 for j in range(i, n):
61 iv_2 = temp[j]
62 temp[j] = iv_2 - iv_1
63 iv_1 = iv_2
64 cf.append(temp[n - 1])
65 # 再计算Newton插值
66 ans = 0
67 for i in range(n):
68 w = 1
69 # 计算多项式部分,让差分作为其系数
70 for j in range(i):
71 w *= ((t + j) / (j + 1))
72 ans += w*cf[i]
73 return ans
74
75 # 当对多个值使用Newton插值时,使用map()建立映射:
76 # Iterator = map(Newton, Iterable)
77
78 # 数值运算时使用float参与运算,dtype定为内置float
79
80 x = np.array((1,2,3,4,5,6), dtype=float)
81 y = np.array((1.0, 1.2599, 1.4422, 1.5874, 1.71, 1.8171), dtype=float)
82 print(Newton_iplt_forword(x, y, 5.6))
83 print(Newton_iplt_backword(x, y, 5.6))
84 # 结果:1.775416 测试成功!