Python 中的 Hermite 插值

Question

我有这个计算Hermite插值的程序。

问题是，它的行为真的很糟糕。

这是 35 个切比雪夫节点的图表。如果我放更多的点，开始的峰值会更高（这个节点数量大约是 10^7）。

我使用拉格朗日方法对相同的函数进行了插值（绿色，它发生了偏移，因此可以看到），如您所见，它看起来不错。

代码如下：

def hermit_interpolate(input): #input is list of tuples [(x1,y1),(x2,y2)...] xi are Chebyshev nodes
points = [(input[0][0], input[0][1] - 0), (input[0][0], calculate_f_p_x(input[0][0]))] #"input[0][1] - 0" this is just to change type of second element
                                                                                    #calculate_f_p_x returns value of derivative
for k in range(1, len(input)): #Divided differences and derivatives in one list alternately
    points.append((input[k][0], (input[k][1] - input[k - 1][1]) / (
        input[k][0] - input[k - 1][0])))
    points.append((input[k][0], calculate_f_p_x(input[k][0])))
x, c = zip(*points)
x = list(x)
c = list(c)
n = len(points)
for i in range(2, n): #calculating factors
    for j in range(n - 1, i - 1, -1):
        c[j] = (c[j] - c[j - 1]) / (x[j] - x[j - i])

def result_polynomial(xpoint): #here is function to calculate value for given x
    val = c[0]
    factor = 1.0
    for l in range(1, n):
        factor *= (xpoint - x[l - 1])
        val += (c[l] * factor)
    return val

return result_polynomial

我看不出这里有什么问题。谢谢！

Answer 1

此代码实际有效：

def hermit_interpolate(input):  #input is list of tuples [(x1,y1),(x2,y2),...,(xn,yn)] xi are Chebyshev nodes
  n = len(input)
  points = numpy.zeros(shape=(2 * n + 1, 2 * n + 1))
  X, Y = zip(*input)
  X = list(X)
  Y = list(Y)

  for i in range(0, 2 * n, 2):
      points[i][0] = X[i / 2]
      points[i + 1][0] = X[i / 2]
      points[i][1] = Y[i / 2]
      points[i + 1][1] = Y[i / 2]

  for i in range(2, 2 * n + 1):
      for j in range(1 + (i - 2), 2 * n):
          if i == 2 and j % 2 == 1:
              points[j][i] = calculate_f_p_x(X[j / 2]);

          else:
              points[j][i] = (points[j][i - 1] - points[j - 1][i - 1]) / (
                  points[j][0] - points[(j - 1) - (i - 2)][0])

  def result_polynomial(xpoint):  #here is function to calculate value for given x
      val = 0
      for i in range(0, 2 * n):
          factor = 1.
          j = 0
          while j < i:
              factor *= (xpoint - X[j / 2])
              if j + 1 != i:
                  factor *= (xpoint - X[j / 2])
                  j += 1
              j += 1
          val += factor * points[i][i + 1]
      return val
return result_polynomia

Python 中的 Hermite 插值

Hermite interpolation in Python

python

interpolation

numerical-methods

hermite