用matplotlib绘制能量势

Question

我想画出引力势来突出它的极值（两个天体周围的拉格朗日点）。在

以下函数返回每组坐标x和y的电势：

def gravitational_potential(M,m,R,x,y):
    G = 6.674*10**(-11)
    omega2 = G*(M+m)/(R**3)
    r = np.sqrt(x**2+y**2)
    r2 = R*m/(M+m)
    r1 = R-r2

    phi = -G*(M/abs(r-r1)+m/abs(r-r2))-1/2*omega2*(x**2+y**2)

    return phi

我想用meshgrid和plot_surface来绘制电势的表面和轮廓，但它不起作用。在

我做错什么了？在

附言：我设法用WolframAlpha绘制了潜力图，所以我知道数学是可行的。在

^{pr2}$

当我执行它时，我得到以下结果：

runfile('C:/Users/python/Google Drive/lagrangepoint_maths/potential/gravitational_potential.py', wdir='C:/Users/python/Google Drive/lagrangepoint_maths/potential')
C:/Users/python/Google Drive/lagrangepoint_maths/potential/gravitational_potential.py:13: RuntimeWarning: divide by zero encountered in divide
  phi = -G*(M/abs(r-r1)+m/abs(r-r2))-1/2*omega2*(x**2+y**2)

这不是我想要的。Z有点问题。我想要这样的东西：

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
from matplotlib import cm

fig = plt.figure()
ax = fig.gca(projection='3d')
X, Y, Z = axes3d.get_test_data(0.05)
ax.plot_surface(X, Y, Z, rstride=8, cstride=8, alpha=0.9)
cset = ax.contour(X, Y, Z, zdir='z', offset=-100, cmap=cm.coolwarm)
cset = ax.contour(X, Y, Z, zdir='x', offset=-40, cmap=cm.coolwarm)
cset = ax.contour(X, Y, Z, zdir='y', offset=40, cmap=cm.coolwarm)

ax.set_xlabel('X')
ax.set_xlim(-40, 40)
ax.set_ylabel('Y')
ax.set_ylim(-40, 40)
ax.set_zlabel('Z')
ax.set_zlim(-100, 100)

plt.show()

Answer 1

所有这些都是一个接一个调试的：

1. python2中的整数除法如果命名符小于分母，结果是0。您可以from __future__ import division或更正代码以除以浮点数。

2. 如果要显示介于-2 x 10^-8和+2 x 10^-8之间的数字，则将z_限制设置为-40到40是没有用的。

3. 如果要在绘图中显示小要素，则不应粗略地将打印分辨率设置为rstride=8, cstride=8。

总的来说，你会得到这样的结果：