from scipy import *
import matplotlib.pyplot as plt
def plot_elines(x_grid, y_grid, potential, field):
fig, ax = plt.subplots(figsize=(13, 13))
im_cs = ax.contour(x_grid, y_grid, potential, 18, cmap='inferno')
plt.clabel(im_cs, inline=1, fontsize=7)
ax.quiver(x_grid[::3, ::3], y_grid[::3, ::3],
field[0, ::3, ::3], field[1, ::3, ::3],)
ax.set_xlabel("$x$")
ax.set_ylabel("$y$")
plt.show()
# define q configuration (x,y positions)
charges = [
[1, 1],
]
xx, yy = meshgrid(linspace(-4, 4), linspace(-4, 4))
# potential function
e_pot = 0.
for idx, q in enumerate(charges):
dx, dy = xx-q[0], yy-q[1]
rr = hypot(dx, dy)
e_pot += 1/(4*pi) * 1./rr
e_field = gradient(-e_pot)
e_field /= hypot(e_field[0], e_field[1]) * 5
# why is this needed?
e_field[0] = e_field[0].T
e_field[1] = e_field[1].T
plot_elines(xx, yy, e_pot, e_field)
我对使用 numpy/scipy 中的 gradient
函数有疑问。
我在这里绘制了电场等势线和单个正电荷的场矢量。定义是
E = -grad(V)
根据定义,场矢量(颤动)和等势线(等值线)在空间中的所有点都应该彼此正交,并且由于电荷是正的,因此箭头需要指向远离电荷本身的方向。
我正在使用scipy的gradient
函数来计算E,但是我发现如果不转置gradient
函数输出的x-y网格,输出是错误的。
比较两个输出(使用 .T
(正确)和不使用 .T
(错误)):
与
为什么需要转置?还是我计划错误?
谢谢。
最佳答案
事实上,转置给出正确的图纯属巧合,因为电荷在 x 和 y 上对称放置(即在 45° 线上)。
真正的问题来自于对numpy.gradient
的错误解释。它将返回轴方向的梯度。第一个数组对应轴 0,第二个数组对应轴 1。现在,您的情况下的轴 0 对应于 y 轴,轴 1 对应于 x 轴。
e_field_y, e_field_x = numpy.gradient(-e_pot)
因此,当您在颤动图中选择相应的场分量时,您需要选择第一个条目作为 y 分量,第二个条目作为 x 分量。
ax.quiver(x_grid[::3, ::3], y_grid[::3, ::3],
field[1, ::3, ::3], field[0, ::3, ::3],)
完整的代码如下所示
import numpy as np
import matplotlib.pyplot as plt
def plot_elines(x_grid, y_grid, potential, field):
fig, ax = plt.subplots(figsize=(13, 13))
im_cs = ax.contour(x_grid, y_grid, potential, 18, cmap='inferno')
plt.clabel(im_cs, inline=1, fontsize=7)
ax.quiver(x_grid[::3, ::3], y_grid[::3, ::3],
field[1, ::3, ::3], field[0, ::3, ::3],)
ax.set_xlabel("$x$")
ax.set_ylabel("$y$")
ax.set_aspect("equal")
plt.show()
# define q configuration (x,y positions)
charges = [
[1, 0],
]
xx, yy = np.meshgrid(np.linspace(-4, 4), np.linspace(-4, 4))
# potential function
e_pot = 0.
for idx, q in enumerate(charges):
dx, dy = xx-q[0], yy-q[1]
rr = np.hypot(dx, dy)
e_pot += 1/(4*np.pi) * 1./rr
e_field = np.gradient(-e_pot)
e_field /= np.hypot(e_field[0], e_field[1]) * 5
plot_elines(xx, yy, e_pot, e_field)
我把电荷放在反对角线上。
最后的提示:对此类代码的良好一致性检查是使用形状不同的网格。例如。如果你沿 x 取 50 个点,沿 y 取 51 个点,你会得到一个错误,而不是看似有效的代码,并且会更容易地解决根本问题。
关于python - 梯度图和箭袋图,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/54471306/