在你这里也贴一下，看看我的方法。

import numpy as np
import matplotlib.pyplot as plt
from numpy.linalg import norm
from matplotlib import style
style.use('ggplot')

fig = plt.figure()
Ax = fig.gca()

epsilon = 0.0001
n = 100000

position = np.zeros((n,2))
velocity = np.zeros((2))
acceleration = np.zeros((2))

position[0,0] = float(input("输入x初始坐标值："))
position[0,1] = float(input("输入y初始坐标值："))
velocity[0] = float(input("输入x分量初始速度："))
velocity[1] = float(input("输入y分量初始速度："))


r = norm(position[0])
r3 = 1/(r**3)
acceleration = position[0] * -r3
velocity += acceleration * epsilon * 0.5

star = 1

if(norm(velocity) > np.sqrt(2/r)):
    star = 0
else:
    i = 1
    t = 0
    while(i<n):
        position[i] = position[i-1]+velocity*epsilon
        if(position[i-1,1]*position[i,1]< 0):
            t += 1
        if(t==2):
            break
        r = norm(position[i])
        r3 = 1/(r**3)
        acceleration = position[i]*-r3
        velocity += acceleration*epsilon
        if(norm(velocity) > np.sqrt(2/r)):
            star = 0
            break
        i += 1

if(star == 1):
    print('行星的公转周期为：%f' %(i*epsilon))
    plt.plot(position[:i,0], position[:i,1],'b-')
    plt.plot(0,0,"ro")
    plt.show()
else:
    print("start is gone")