在python中获取数据的派生

Question

我写了一个程序来得到导数。插值变量pline用于计算f（x+h）。红线是余弦的导数，绿线是余弦的余弦，蓝线是正弦函数。红色和蓝色的线是匹配的。它在以下方面工作得很好。

from scipy.interpolate import InterpolatedUnivariateSpline
import numpy as np 
import matplotlib.pyplot as plt

pi = np.pi
x = np.arange(0,5*pi,0.2*pi)
y = np.cos(x) 
f2 = InterpolatedUnivariateSpline(x, y)
#Get dervative
der = []
for i in range(len(y)):

    h = 1e-4
    der.append( ( f2(x[i]+h)-f2(x[i]-h) )/(2*h) )
der = np.array(der)   

plt.plot(x, der, 'r', x, y, 'g', x, -np.sin(x),'b')
plt.show()

但我遇到了一些问题。在我的项目中，我的变量x（频率）从10^7变化到2.2812375*10^9，它的步长是22487500，所以我改变了数组x。 结果，我得到了以下结果。

导数是一条直线，几乎接近于0，它不是正弦函数。我该怎么解决？

Answer 1

你有个问题。这意味着当把一个大浮点数加到一个小浮点数上时，小浮点数的精度会部分损失，因为numpy double只能容纳64位的信息。

要解决这个问题，你必须确保你加/乘/除的数字的比例没有太大的不同。一个简单的解决方案是用x除以1e9，或者用h乘以1e9。如果你这样做，你得到的精度基本上与你的例子相同。

另外，如果x具有足够高的分辨率，则用数值方法来区分函数的更简单方法是der = np.diff(y) / np.diff(x)。这样你就不用担心h的设置了，但是在这个例子中，注意dy是一个比y短的元素，而dy[i]实际上是`（x[i]+x[i+1]）/2处导数的近似值。因此，要绘制它，您需要：

der = np.diff(y) / np.diff(x)
x2 = (x[:-1] + x[1:]) / 2
plt.plot(x2, der, 'r', x, y, 'g', x, -np.sin(x),'b')