0

我正在尝试使用 py 制作一个基本的重力模拟,但由于某种原因,它会将所有线条绘制成直线,当我遇到困难时,我查看了几个示例,但它们都使用类似的方程/绘制数据的方式,所以不确定我哪里出错了

class Body:

    # A class to initialise a body of mass to plot

    def __init__(self, id, mass, coordinates, velocities):
        self.id = id
        self.coordinates = np.array(coordinates, dtype=float)
        self.v = np.array(velocities, dtype=float)
        self.mass = mass
        self.a = np.zeros(3, dtype=float)
        MOTION_LOG.append({"name": self.id, "x":[coordinates[0]], "y": [coordinates[1]], "z": [coordinates[2]]})


    # Procedure for checking gravity effects on body
    def gravity(self):
        self.a = 0
        for body in bodies:
            if body != self:
                dist = body.coordinates - self.coordinates
                r = np.sqrt(np.sum(dist**2))
                self.a += (SETTINGS['G'] * body.mass * dist / r**3) ** 2


# Procedure to plot the new coordinates of the body         
    def move(self):
        self.v += self.a * SETTINGS["deltaT"]
        self.coordinates += self.v * SETTINGS['deltaT']

然后实际模拟它我做了

# For loop to run a simulation for a specific time set
for step in range(int(SETTINGS["tLimit"] / SETTINGS["deltaT"])):
    SETTINGS['elapsedT'] += SETTINGS['deltaT']
    if SETTINGS['elapsedT'] % SETTINGS["frequency"] == 0:
        prog = ((SETTINGS['elapsedT'] / SETTINGS['tLimit'])*100)//1
        for index, location in enumerate(MOTION_LOG):
            location["x"].append(bodies[index].coordinates[0])
            location["y"].append(bodies[index].coordinates[1])
            location["z"].append(bodies[index].coordinates[2])
        print(f"{prog}%")

    for body in bodies:
        body.gravity()
    for body in bodies:
        body.move()

然后最后在我做的图表上绘制它

fig = plt.figure()
ax = fig.add_subplot(1,1,1, projection='3d')

for body in MOTION_LOG:
    ax.plot(body["x"], body["y"], body["z"])

plt.show()

不知道我是不是刚刚做错了加速,或者我只是把点画错了,但我见过的其他例子做的事情并没有太大的不同

示例数据

SETTINGS = {
    'G' : 6.67e-11,
    'deltaT' : 172800,
    'elapsedT' : 0,
    'tLimit' : 315360000,
    "frequency": 1,
}

MOTION_LOG = []
bodies = []
set_bodies = {
    "Sun": {
        "id": "Sun",
        "mass": 1e20,
        "coordinates": [0, 0, 0],
        "velocities": [0, 0, 0]
    },
    "Earth": {
        "id": "Earth",
        "mass": 1e3,
        "coordinates": [1e2, -1e2, 0],
        "velocities": [0, 0, 0]
    }
}
for body, body_data in set_bodies.items():
    bodies.append(Body(**body_data))
4

1 回答 1

1

首先,你的加速度是平方的。你的行是:

self.a += (SETTINGS['G'] * body.mass * dist / r**3) ** 2

你想要的是:

self.a += (SETTINGS['G'] * body.mass * dist / r**3) 

其次,你的数字完全不正常。永远,永远不要只是猜测数字的好组合用于这样的模拟。要么使用完全自然的单位(重力常数、太阳质量和 AU 都等于 1),要么使用所有实际数字。当你不这样做时,你只是在猜测什么时间步是合适的,而且很难做到这一点。例如,根据您用于太阳质量、重力常数和地球轨道半径的数字,一个时间单位约为 7.7 年。为了得到看起来不错的东西,那么你需要一个时间步长 about7e-4并运行到最后一个时间 about 1.3,并使用地球的初始速度[5000, 5000, 0](你还需要改变你添加到MOTION_LOG因为小花车不能很好地与 mod ( %)) 配合使用。但不要那样做——使用合理的数字。当您不必花费数小时来回转换为您真正想要的单位时,您最终会更快乐。

第三,您的“地球”物体没有初始速度。无论您如何修正方程式或单位,它都会落入太阳。

以下是一些使用 SI 单位的大约正确的“实数”数字:

set_bodies = {
    "Sun": {
        "id": "Sun",
        "mass": 2e30,
        "coordinates": [0, 0, 0],
        "velocities": [0, 0, 0]
    },
    "Earth": {
        "id": "Earth",
        "mass": 6e24,
        "coordinates": [1.5e11, 0, 0],
        "velocities": [0, 3e4, 0]
    }
}
于 2021-12-19T23:25:23.480 回答