dus*_*tin 7 python numpy matplotlib
我试图用柯蒂斯复制轨道力学中的一个情节,但我不能完全理解它.但是,我已经np.arctan2
从切换到np.arctan
.
也许我的实施arctan2
不正确?
import pylab
import numpy as np
e = np.arange(0.0, 1.0, 0.15).reshape(-1, 1)
nu = np.linspace(0.001, 2 * np.pi - 0.001, 50000)
M2evals = (2 * np.arctan2(1, 1 / (((1 - e) / (1 + e)) ** 0.5 * np.tan(nu / 2) -
e * (1 - e ** 2) ** 0.5 * np.sin(nu) / (1 + e * np.cos(nu)))))
fig2 = pylab.figure()
ax2 = fig2.add_subplot(111)
for Me2, _e in zip(M2evals, e.ravel()):
ax2.plot(nu.ravel(), Me2, label = str(_e))
pylab.legend()
pylab.xlim((0, 7.75))
pylab.ylim((0, 2 * np.pi))
pylab.show()
Run Code Online (Sandbox Code Playgroud)
在下图中,突然出现不连续性.该函数应该是平滑的,并且在(0,2pi)的y范围内在0和2pi处连接而不接触0和2pi.
教科书情节和方程:
应Saullo Castro的要求,我被告知:
"问题可能在于arctan函数,它将"原则值"作为输出.
因此,如果x是第二或第三象限中的角度,则arctan(tan(x))不产生x.如果你将arctan(tan(x))从x = 0绘制到x = Pi,你会发现它在x = Pi/2处有一个不连续的跳跃.
对于你的情况,我不相信写arctan(arg),我相信你会写arctan2(1,1/arg),其中arg是你的arctan函数的参数.这样,当arg变为负数时,arctan2将在第二象限而不是第四象限产生一个角度."
Sau*_*tro 10
通常的做法是在负面结果中加上2*pi arctan()
,这可以有效地完成.OP建议用arctan2(1,1/x)代替arctan(x),也是由@ Yay295指出的Maple 15文档建议的,它产生相同的结果,而不需要求和2*pi.两者都显示如下:
import pylab
import numpy as np
e = np.arange(0.0, 1.0, 0.15).reshape(-1, 1)
nu = np.linspace(0, 2*np.pi, 50000)
x = ((1-e)/(1+e))**0.5 * np.tan(nu/2.)
x2 = e*(1-e**2)**0.5 * np.sin(nu)/(1 + e*np.cos(nu))
using_arctan = True
using_OP_arctan2 = False
if using_arctan:
M2evals = 2*np.arctan(x) - x2
M2evals[ M2evals<0 ] += 2*np.pi
elif using_OP_arctan2:
M2evals = 2 * np.arctan2(1,1/x) - x2
fig2 = pylab.figure()
ax2 = fig2.add_subplot(111)
for M2e, _e in zip(M2evals, e.ravel()):
ax2.plot(nu.ravel(), M2e, label = str(_e))
pylab.legend(loc='upper left')
pylab.show()
Run Code Online (Sandbox Code Playgroud)