为什么我的Python代码中重力计算有错误？

Question

import pygame  
from pygame import *  
import math  
class orb:  
    def __init__(self, mass=1, s0=[0,0], v0=[0,0], a=[0,0], r=10, color=(255, 255, 255), fix = False, F=0, Fx=0, Fy=0):  
        self.a= a  
        self.mass=mass  
        self.r=r  
        self.t = 0  
        self.v0=v0  
        self.v=[self.v0[0], self.v0[1]]  
        self.s0=s0  
        self.s=[self.s0[0], self.s0[1]]  
        self.color=color  
        self.fix = fix  
        self.F=F  
        self.Fx = Fx  
        self.Fy = Fy  

calculating the position of the earth

    def comPos(self, dt = 0.01):
        if not self.fix:  
            self.v[0] = self.v0[0] + self.a[0] * dt
            self.v[1] = self.v0[1] + self.a[1] * dt
            ds_x = self.v0[0] * dt + (1/2) * self.a[0] * dt**2
            ds_y = self.v0[1] * dt + (1/2) * self.a[1] * dt**2
            self.s[0] = self.s0[0] + ds_x
            self.s[1] = self.s0[1] + ds_y
            self.v0[0] = self.v[0]
            self.v0[1] = self.v[1]
            self.s0[0] = self.s[0]
            self.s0[1] = self.s[1]
            self.t += dt

calculating the gravity of the earth

    def comForce(self, other, G = 10000.0, F = 0.0, dist_cr = 0.0, Fx = 0.0, Fy = 0.0):
        self.F = G * self.mass * other.mass / ( (self.s[0]- other.s[0])**2 + (self.s[1] - other.s[1])**2 )
        self.Fx = (self.F * math.fabs(math.cos(math.atan((other.s[1] - self.s[1]) / (other.s[0] - self.s[0])))))
        self.Fy = (self.F * math.fabs(math.sin(math.atan((other.s[1] - self.s[1]) / (other.s[0] - self.s[0])))))
        if ((other.s[0] - self.s[0])**2 + (other.s[1] - self.s[1])**2)**(1/2) >= dist_cr:
            if other.s[0] - self.s[0] == 0:
                self.a[0] = 0
            else:
                self.a[0] = (( (other.s[0] - self.s[0]) / math.fabs(other.s[0] - self.s[0])) * (self.Fx/self.mass))
            if other.s[1] - self.s[1] == 0:
                self.a[1] = 0
            else:
                self.a[1] = (( (other.s[1] - self.s[1]) / math.fabs(other.s[1] - self.s[1])) * self.Fy / self.mass)

Displaying on the screen

screen = pygame.display.set_mode((800, 600))
clock = pygame.time.Clock()
earth = orb(s0=[150, 300], v0=[0,60], r=10, mass= 3)
sun = orb(s0=[400, 300], r=70, color = [255, 0, 0], fix=True, mass=100)
orb = [sun, earth]

running = True
while running:
    for event in pygame.event.get():
        if event.type == pygame.QUIT:
            running = False
    screen.fill((0, 0, 0))
    pygame.draw.circle(screen, earth.color, earth.s, earth.r, 0)
    pygame.draw.circle(screen, sun.color, sun.s, sun.r, 0)
    earth.comForce(sun)
    earth.comPos()
    pygame.display.update()
    clock.tick(720)  

`我要模拟绕太阳旋转的地球。但随着时间的推移，轨道逐渐发生变化。我试图寻找原因，但无论我如何努力，我都不知道原因。你能告诉我为什么轨迹会改变吗？

我希望轨迹保持不变，不随时间变化`

Answer 1

• 您可以使用

  float64

  获得更高的准确度。
• 请注意，使用数学

  fabs

  、

  sin

  、

  cos

  等会降低准确性。
• 此外，您可以使用小

  dt

  。
• 还有其他几种方法可以使其准确。

代码：

import pygame
import numpy as np


class Orb:
    def __init__(self, mass=1, s0=[0, 0], v0=[0, 0], r=10, color=(255, 255, 255), fix=False):
        self.mass = np.float64(mass)
        self.r = r
        self.s = np.array(s0, dtype=np.float64)
        self.v = np.array(v0, dtype=np.float64)
        self.a = np.array([0.0, 0.0], dtype=np.float64)
        self.color = color
        self.fix = fix

    def comForce(self, other, G=np.float64(10000.0)):
        dx = other.s[0] - self.s[0]
        dy = other.s[1] - self.s[1]
        dist_sq = dx ** 2 + dy ** 2
        dist = np.sqrt(dist_sq)
        F = G * self.mass * other.mass / dist_sq
        Fx = F * dx / dist
        Fy = F * dy / dist
        return np.array([Fx, Fy], dtype=np.float64)

    def update(self, other, dt=np.float64(0.001)):
        if not self.fix:
            force = self.comForce(other)
            self.a = force / self.mass
            self.s += self.v * dt + 0.5 * self.a * dt**2
            new_force = self.comForce(other)
            new_a = new_force / self.mass
            self.v += 0.5 * (self.a + new_a) * dt


pygame.init()
screen = pygame.display.set_mode((900, 600))
clock = pygame.time.Clock()

earth = Orb(mass=np.float64(2), s0=[150, 300], v0=[0, 60], r=10, color=(0, 0, 255))
sun = Orb(mass=np.float64(100), s0=[420, 300], v0=[0, 0], r=50, color=(255, 0, 0), fix=True)
orbs = [sun, earth]

running = True
frames = 1000
while running:
    for event in pygame.event.get():
        if event.type == pygame.QUIT:
            running = False

    for _ in range(frames):
        earth.update(sun, dt=np.float64(0.001))

    screen.fill((0, 0, 0))
    for orb in orbs:
        pygame.draw.circle(screen, orb.color, (int(orb.s[0]), int(orb.s[1])), orb.r)

    pygame.display.flip()
    clock.tick(60)

pygame.quit()