【发布时间】:2015-12-03 08:54:55
【问题描述】:
我正在尝试查找数据集 (x,y) 的高阶导数。 x 和 y 是长度为 N 的一维数组。
假设我将它们生成为:
xder0=np.linspace(0,10,1000)
yder0=np.sin(xder0)
我定义了接受 2 个数组 (x,y) 并返回 (x1, y1) 的导数函数,其中 y1 是在每个索引处计算的导数:(y[i+1]-y[i])/ (x[i+1]-x[i])。 x1 只是 x[i+1] 和 x[i]
的平均值这是执行此操作的函数:
def deriv(x,y):
delx =np.zeros((len(x)-1), dtype=np.longdouble)
ydiff=np.zeros((len(x)-1), dtype=np.longdouble)
for i in range(len(x)-1):
delx[i] =(x[i+1]+x[i])/2.0
ydiff[i] =(y[i+1]-y[i])/(x[i+1]-x[i])
return delx, ydiff
现在计算一阶导数,我称这个函数为:
xder1, yder1 = deriv(xder0, yder0)
与二阶导数类似,我称这个函数为输入一阶导数:
xder2, yder2 = deriv(xder1, yder1)
然后继续:
xder3, yder3 = deriv(xder2, yder2)
xder4, yder4 = deriv(xder3, yder3)
xder5, yder5 = deriv(xder4, yder4)
xder6, yder6 = deriv(xder5, yder5)
xder7, yder7 = deriv(xder6, yder6)
xder8, yder8 = deriv(xder7, yder7)
xder9, yder9 = deriv(xder8, yder8)
在我达到 7 阶后发生了一些奇怪的事情。7 阶变得非常嘈杂!正如预期的那样,早期的导数都是正弦或余弦函数。然而,7 阶是一个嘈杂的正弦波。因此,在那之后的所有衍生品都爆炸了。
知道发生了什么吗?
【问题讨论】:
-
近似推导的差分公式都是病态的,因为我们减去两个比较大的(
x[i],x[i+1],y[i],y[i + 1]) 数字多次以获得相对较小的差异。order = 7似乎是你的下车点。 -
有什么想法吗?
-
寻找一个处理数字部分的库。 numdifftools 浮现在脑海中。如您所见,直接的、不稳定的方法只能让您走这么远。
-
一种可能性是为您的数据拟合一些函数,然后您可以对其进行分析导数。
-
我建议阅读Numerical Recipes 的第 5.7 节。它提供了多种计算数值导数的方法,同时最大限度地减少了稳定性问题。您每月可以访问几页。我认为应该足以满足您的目的。
标签: python numpy numerical-methods derivative