我想知道
numpy.gradient
如何运作。
我使用梯度来尝试计算群速度(波包的群速度是频率相对于波数的导数,而不是一组速度)。我向它提供了一个 3 列数组,前 2 列是 x 和 y 坐标,第三列是该点 (x,y) 的频率。我需要计算梯度,我确实期望一个二维向量,即梯度定义
df/dx*i+df/dy*j+df/dz*k
我的函数只是 x 和 y 的函数,我确实期望类似的东西
df/dx*i+df/dy*j
但是我有 2 个数组,每个数组有 3 列,即 2 个 3d 向量;起初我认为两者的和会给我我正在搜索的向量,但 z 分量并没有消失。我希望我的解释足够清楚。我想知道
numpy.gradient
是如何工作的以及它是否是解决我的问题的正确选择。否则我想知道是否还有其他可以使用的 python 函数。
我的意思是:我想计算值数组的梯度:
data=[[x1,x2,x3]...[x1,x2,x3]]
其中 x1,x2 是均匀网格上的点坐标(我在布里渊区上的点),x3 是该点的频率值。我在输入中还给出了两个方向的推导步骤:
stepx=abs(max(unique(data[:,0])-min(unique(data[:,0]))/(len(unique(data[:,0]))-1)
y 方向相同。 我没有在网格上构建数据,我已经有了一个网格,这就是为什么这里给出的答案中的示例对我没有帮助。 一个更合适的例子应该有一个像我一样的点和值网格:
data=[]
for i in range(10):
for j in range(10):
data.append([i,j,i**2+j**2])
data=array(data,dtype=float)
gx,gy=gradient(data)
我可以补充的另一件事是,我的网格不是方形的,而是具有多边形的形状,即二维晶体的布里渊区。
我知道
numpy.gradient
只能在值的方形网格上正常工作,而不是我正在寻找的。即使我将数据作为网格,在原始数据的多边形之外有很多零,这也会向我的梯度添加非常高的向量,从而影响(负面)计算的精度。在我看来,这个模块更像是一个玩具,而不是一个工具,恕我直言,它有严重的局限性。 python 中有更强大的东西,或者我最好从头开始编写其他东西,也许可以切换到 fortran?
您需要给
gradient
一个矩阵来描述 (x,y)
点的角频率值。例如
def f(x,y):
return np.sin((x + y))
x = y = np.arange(-5, 5, 0.05)
X, Y = np.meshgrid(x, y)
zs = np.array([f(x,y) for x,y in zip(np.ravel(X), np.ravel(Y))])
Z = zs.reshape(X.shape)
gx,gy = np.gradient(Z,0.05,0.05)
您可以看到将 Z 绘制为曲面可以得出:
以下是如何解释渐变:
gx
是一个矩阵,给出所有点的变化 dz/dx
。例如gx[0][0] 是 dz/dx
处 (x0,y0
)。可视化gx
有助于理解:
因为我的数据是从
f(x,y) = sin(x+y)
gy 生成的,所以看起来是一样的。
这是一个更明显的例子,使用
f(x,y) = sin(x)
...
f(x,y)
和渐变
更新 让我们看一下 xy 对。
这是我使用的代码:
def f(x,y):
return np.sin(x)
x = y = np.arange(-3,3,.05)
X, Y = np.meshgrid(x, y)
zs = np.array([f(x,y) for x,y in zip(np.ravel(X), np.ravel(Y))])
xy_pairs = np.array([str(x)+','+str(y) for x,y in zip(np.ravel(X), np.ravel(Y))])
Z = zs.reshape(X.shape)
xy_pairs = xy_pairs.reshape(X.shape)
gy,gx = np.gradient(Z,.05,.05)
现在我们可以观察到底发生了什么。假设我们想知道哪个点与
Z[20][30]
处的值相关?然后...
>>> Z[20][30]
-0.99749498660405478
重点是
>>> xy_pairs[20][30]
'-1.5,-2.0'
是吗?让我们检查一下。
>>> np.sin(-1.5)
-0.99749498660405445
是的。
此时我们的梯度分量是什么?
>>> gy[20][30]
0.0
>>> gx[20][30]
0.070707731517679617
那些检查了吗?
dz/dy always 0
检查。
dz/dx = cos(x)
和...
>>> np.cos(-1.5)
0.070737201667702906
看起来不错。
你会注意到它们并不完全正确,那是因为我的 Z 数据不连续,步长为
0.05
并且 gradient
只能近似变化率。