求助：Python三体系统轨迹动画无法正常显示问题排查

我正在开发一个三体系统，尝试对其轨迹进行动画演示。已生成指定步数后的轨迹图，但动画实现失败——使用FuncAnimation()仅显示空白图表，无动画效果，代码如下：

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation

#Inital Conditions
y0 = [-1        ,#x1
  0             ,#y1
  1             ,#x2
  0             ,#y2
  0             ,#x3
  0             ,#y3
  0.306893      ,#vx1
  0.125507      ,#vy1
  0.306893      ,#vx2
  0.125507      ,#vy2
  -2*0.306893   ,#vx3
  -2*0.125507]   #vy3

# Time steps and the total number of steps
N = 10000
T = 0.001 

#Definition of the Function 
def ThreeBody(t,y):
   f = np.zeros(12)

  # The velocities of the three bodies
   f[0] = y[6]
   f[1] = y[7]
   f[2] = y[8]
   f[3] = y[9]
   f[4] = y[10]
   f[5] = y[11]

  # The x and y positions of each object respectively
   f[6] = -(y[0]-y[2])/(((y[0]-y[2])**2+(y[1]-y[3])**2)**(3/2)) \
        -(y[0]-y[4])/(((y[0]-y[4])**2+(y[1]-y[5])**2)**(3/2))
        
   f[7] = -(y[1]-y[3])/(((y[0]-y[2])**2+(y[1]-y[3])**2)**(3/2)) \
        -(y[1]-y[5])/(((y[0]-y[4])**2+(y[1]-y[5])**2)**(3/2))
                  
   f[8] = -(y[2]-y[0])/(((y[2]-y[0])**2+(y[3]-y[1])**2)**(3/2)) \
        -(y[2]-y[4])/(((y[2]-y[4])**2+(y[3]-y[5])**2)**(3/2))
        
   f[9] = -(y[3]-y[1])/(((y[2]-y[0])**2+(y[3]-y[1])**2)**(3/2)) \
        -(y[3]-y[5])/(((y[2]-y[4])**2+(y[3]-y[5])**2)**(3/2))
                  
   f[10]= -(y[4]-y[0])/(((y[4]-y[0])**2+(y[5]-y[1])**2)**(3/2)) \
        -(y[4]-y[2])/(((y[4]-y[2])**2+(y[5]-y[3])**2)**(3/2))
        
   f[11]= -(y[5]-y[1])/(((y[4]-y[0])**2+(y[5]-y[1])**2)**(3/2)) \
        -(y[5]-y[3])/(((y[4]-y[2])**2+(y[5]-y[3])**2)**(3/2))

   return f


#Solving for the positions of all bodies
t = np.linspace(0,N*T,N)
solution = solve_ivp(ThreeBody,[0,800],y0,t_eval=t,rtol=1e-12)


#Evolution in Position with respect to Time
plt.plot(solution.y[0],solution.y[1],'-g') #Positions body 1
plt.plot(solution.y[2],solution.y[3],'-r') #Positions body 2
plt.plot(solution.y[4],solution.y[5],'-b') #Positions body 3
plt.ylabel("Position(y)")
plt.xlabel("Position(x)")
plt.show()

plt.plot(t,solution.y[0])
plt.ylabel("Position (x)")
plt.xlabel("Time")
plt.show()

fig = plt.figure()
ax = plt.axes(xlim=(-1.5, 1.5), ylim=(-0.5, 0.5))
line, = ax.plot([], [], lw=2)

def animate(n):
    line.set_xdata(solution.y[0])
    line.set_ydata(solution.y[1])
    return line,


anim = FuncAnimation(fig, animate, frames= 1000, interval=10)
plt.show()

问题原因 #

• animate函数逻辑错误：当前代码每次更新都把整个天体1的完整轨迹赋给line，没有逐帧递进展示，且初始line无数据，导致首帧空白。
• 缺少多天体展示：只创建了一条line，但三体系统需要三个天体的轨迹/位置标记。
• 帧索引未利用：animate的n参数（当前帧索引）未使用，无法控制逐帧更新的进度。
• 缺少初始化函数：FuncAnimation未指定init_func，初始状态的空数据未正确设置。

修正后的代码 #

以下代码实现了三个天体的轨迹动画（累积轨迹+当前位置标记）：

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation

# 初始条件
y0 = [-1,        # x1
      0,         # y1
      1,         # x2
      0,         # y2
      0,         # x3
      0,         # y3
      0.306893,  # vx1
      0.125507,  # vy1
      0.306893,  # vx2
      0.125507,  # vy2
      -2*0.306893, # vx3
      -2*0.125507] # vy3

# 时间参数
N = 10000
T = 0.001 

# 三体微分方程
def ThreeBody(t,y):
    f = np.zeros(12)
    # 速度分量
    f[0] = y[6]
    f[1] = y[7]
    f[2] = y[8]
    f[3] = y[9]
    f[4] = y[10]
    f[5] = y[11]
    # 加速度分量（万有引力）
    f[6] = -(y[0]-y[2])/(((y[0]-y[2])**2+(y[1]-y[3])**2)**(3/2)) \
           -(y[0]-y[4])/(((y[0]-y[4])**2+(y[1]-y[5])**2)**(3/2))
    f[7] = -(y[1]-y[3])/(((y[0]-y[2])**2+(y[1]-y[3])**2)**(3/2)) \
           -(y[1]-y[5])/(((y[0]-y[4])**2+(y[1]-y[5])**2)**(3/2))
    f[8] = -(y[2]-y[0])/(((y[2]-y[0])**2+(y[3]-y[1])**2)**(3/2)) \
           -(y[2]-y[4])/(((y[2]-y[4])**2+(y[3]-y[5])**2)**(3/2))
    f[9] = -(y[3]-y[1])/(((y[2]-y[0])**2+(y[3]-y[1])**2)**(3/2)) \
           -(y[3]-y[5])/(((y[2]-y[4])**2+(y[3]-y[5])**2)**(3/2))
    f[10]= -(y[4]-y[0])/(((y[4]-y[0])**2+(y[5]-y[1])**2)**(3/2)) \
           -(y[4]-y[2])/(((y[4]-y[2])**2+(y[5]-y[3])**2)**(3/2))
    f[11]= -(y[5]-y[1])/(((y[4]-y[0])**2+(y[5]-y[1])**2)**(3/2)) \
           -(y[5]-y[3])/(((y[4]-y[2])**2+(y[5]-y[3])**2)**(3/2))
    return f

# 求解微分方程
t = np.linspace(0, N*T, N)
solution = solve_ivp(ThreeBody, [0, 800], y0, t_eval=t, rtol=1e-12)

# 设置动画画布
fig = plt.figure()
ax = plt.axes(xlim=(-1.5, 1.5), ylim=(-0.5, 0.5))
# 创建三个天体的轨迹线和当前位置标记
line1, = ax.plot([], [], '-g', lw=1, label='天体1')
line2, = ax.plot([], [], '-r', lw=1, label='天体2')
line3, = ax.plot([], [], '-b', lw=1, label='天体3')
dot1, = ax.plot([], [], 'go', markersize=5)
dot2, = ax.plot([], [], 'ro', markersize=5)
dot3, = ax.plot([], [], 'bo', markersize=5)
ax.legend()
ax.set_xlabel('Position(x)')
ax.set_ylabel('Position(y)')

# 初始化函数：设置空数据
def init():
    line1.set_data([], [])
    line2.set_data([], [])
    line3.set_data([], [])
    dot1.set_data([], [])
    dot2.set_data([], [])
    dot3.set_data([], [])
    return line1, line2, line3, dot1, dot2, dot3

# 动画更新函数：根据帧索引n更新到前n个点
def animate(n):
    # 更新轨迹（累积到第n帧）
    line1.set_data(solution.y[0][:n], solution.y[1][:n])
    line2.set_data(solution.y[2][:n], solution.y[3][:n])
    line3.set_data(solution.y[4][:n], solution.y[5][:n])
    # 更新当前位置
    dot1.set_data(solution.y[0][n], solution.y[1][n])
    dot2.set_data(solution.y[2][n], solution.y[3][n])
    dot3.set_data(solution.y[4][n], solution.y[5][n])
    return line1, line2, line3, dot1, dot2, dot3

# 创建动画：frames设为总步数，interval控制帧间隔（毫秒）
anim = FuncAnimation(fig, animate, init_func=init, frames=N, interval=5, blit=True)
plt.show()

关键修正点 #

• 新增init_func初始化所有线条和标记的空数据
• animate函数利用n索引，逐帧更新累积轨迹和当前位置
• 为三个天体分别创建轨迹线和位置标记，匹配静态图的颜色
• 设置blit=True提升动画渲染效率
• 调整frames为总步数N，确保完整展示整个运动过程

内容的提问来源于stack exchange，提问作者Belal Bahaa