numpy 用Python计算梯度

9jyewag0  于 2023-10-19  发布在  Python
关注(0)|答案(1)|浏览(97)

我想知道numpy.gradient是如何工作的。我使用梯度来计算 * 群速度 *(波包的群速度是频率对波数的导数,而不是一组速度)。我给它一个3列数组,前2列是x和y坐标,第三列是该点的频率(x,y)。我需要计算梯度,我希望一个2D矢量,梯度定义

df/dx*i+df/dy*j+df/dz*k

我的函数只是x和y的函数,我期望的是

df/dx*i+df/dy*j

但是我得到了2个数组,每个数组有3列,即。2三维矢量;一开始我以为,这两个的和,会给予我,我正在寻找的矢量,但z分量并不为零。我希望我的解释足够清楚。我想知道numpy.gradient是如何工作的,它是否是解决我的问题的正确选择。我想知道有没有其他python函数可以使用。
我的意思是:我想计算一个值数组的梯度:

data=[[x1,x2,x3]...[x1,x2,x3]]

其中x1,x2是均匀网格上的点坐标(我在布里渊区上的点),x3是该点的频率值。我给予在输入中也步骤为推导的2个方向:

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)

我可以补充的另一件事是,我的网格不是正方形的,而是具有作为2D晶体的布里渊区的多边形的形状。
我知道numpy.gradient只在值的正方形网格上工作,而不是我正在搜索的。即使我把数据做成一个网格,在原始数据的多边形之外有很多零,这也会给我的梯度增加很高的向量,影响(负面)计算精度。这个模块在我看来更像是一个玩具而不是一个工具,它有严重的局限性。有什么更强大的Python或我最好从头开始写别的东西,也许切换到Fortran?

klr1opcd

klr1opcd1#

您需要给予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)生成的,所以戈伊看起来是一样的。
下面是一个使用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.05gradient只能近似变化率。

相关问题