太阳系模型 Python - 列出操作难度

Question

这是我在这里的第一篇文章。 所以我正在尝试使用可视化python制作一个模型太阳系。我将每个行星定义为一个球体，具有半径、与太阳的距离、质量和动量变量。然后将每个行星（或天体）放入一个列表结构中。正如你所看到的，目前我有一个列表，其中包含 [earth, moon, mars] 太阳被排除在外，原因我稍后会解释。

所以当我试图计算每个身体对另一个身体造成的力时，我的问题就来了。我在这里所做的是，对于身体列表中的每个第 i 个值，计算它和第 n 个身体之间的力，列表身体中第 i 个身体上的力是第 i 个身体和每个第 n 个身体之间的力的总和正文从 0 到列表末尾。 （即列表中所有其他物体的所有力的总和）

这适用于月球和火星（列表中的第 2 和第 3 项），但不适用于地球。下面代码的输出是，

<3.57799e+022, 0, 0>
<4.3606e+020, 0, 0>
<1.64681e+021, 0, 0>
<-1.#IND, -1.#IND, -1.#IND>   -  this is the total force on earth. 
<0, 2.07621e+027, 0>
<0, 9.83372e+027, 0>



from visual import *


AU = 149.6e9
MU = 384.4e6 # moon - earth orbital - radius
MarU = 227.92e9

G =6.673e-11

sun_mass =2e30
sun_radius =6.96e8
earth_mass =6e24
earth_radius =6.37e6
moon_mass =7.35e22
moon_radius =1.74e6
mars_mass = 6.41e23
mars_radius = 3390000
sun = sphere ( pos =(0 , 0 ,0) , velocity = vector (0 ,0 ,0) ,mass = sun_mass , radius =0.1* AU , color = color . yellow )
earth = sphere ( pos =( AU , 0 ,0)  ,mass = earth_mass , radius =63170000, color = color . cyan ,make_trail=True )# create a list of gravitating objects
moon = sphere ( pos =( AU+MU , 0 ,0)  ,mass = moon_mass , radius =17380000 , color = color . white, make_trail=True )
mars = sphere ( pos =( MarU , 0 ,0)  ,mass = mars_mass , radius = mars_radius , color = color . red, make_trail=True )

#initialise values:
we = 1.9578877e-7
wm = sqrt(G*earth.mass/3.38e8**3)
wma = 9.617e-5
dt = 3*60

earth.mom = vector(0,1.5e11*earth.mass*we,0)
mars.mom = vector(0, 9.833720638948e+27,0)
moon.mom = moon.mass*(earth.mom/earth.mass+vector(0,-3.48e8*wm,0))

bodies = [earth, moon, mars]

*N = 0
initialdiff = 0
for i in bodies:
        initialdiff = i.pos - sun.pos
        TotalForce = (G * i. mass * sun. mass * norm ( initialdiff )/ initialdiff . mag2)
        print TotalForce
        while N < len(bodies):
            if N!=i:
                diff = i.pos - bodies[N].pos
                Force = (G * i. mass * bodies[N]. mass * norm ( diff )/ diff . mag2) 
                TotalForce = TotalForce + Force
                i.mom = i.mom+TotalForce*dt
                N = N+1 
            else:
                N = N+1
print earth.mom
print moon.mom
print mars.mom*

感谢您提供的任何帮助。

Answer 1

'''Abel Tilahun HW 3 '''

# 02/03 / 2015

# make the necessary imports , visual import gives visual output, cos, sin and pi allows calculation of the positions
# in real time. The numpy arange import returns evenly spaced values within a given interval (angles)

from visual import *

from math import cos,sin,pi

from numpy import arange

# i used 'a' to magnify the sizes of the planet for better visualization
a=600

# the following line defines a sphere for the sun centered at the origin(0,0,0) all the other planets are positioned at a radial distance from this center.
Sun=sphere(pos=(0,0,0),radius=6955500*(a/100),color=color.yellow)

#the next 5 commands code the planets with their center being positioned at a distance of the radii of their orbit. Their radii are multiplied by a factor of 'a' to magnify them

Mercury=sphere(pos=vector(579e5,0,0),radius=2440*(a),color=color.red)

Venus=sphere(pos=vector(1082e5,0,0),radius=6052*a,color=color.orange)

Earth=sphere(pos=vector(1496e5,0,0),radius=6371*a,color=color.green)

Mars=sphere(pos=vector(2279e5,0,0),radius=3386*a,color=color.white)

Jupiter=sphere(pos=vector(7785e5,0,0),radius=69173*(a),color=color.cyan)

Saturn=sphere(pos=[14334e5,0,0],radius=57316*(a),color=color.magenta)

# the for loop calculates position of the planets by changing 
# the arange function increases the angle from 0 to pi with a small increment of 0.025 each time

for theta in arange(0,100*pi,0.025):

    rate(30)

    x = 579e5*cos(theta)

    y = 579e5*sin(theta)

    Mercury.pos = [x,y,0]


    x = 1082e5*cos(theta)

    y = 1082e5*sin(theta)

    Venus.pos = [x,y,0]


    x = 1496e5*cos(theta)

    y = 1496e5*sin(theta)

    Earth.pos = [x,y,0]


    x = 2279e5*cos(theta)

    y = 2279e5*sin(theta)

    Mars.pos = [x,y,0]


    x = 7785e5*cos(theta)

    y = 7785e5*sin(theta)

    Jupiter.pos = [x,y,0]


    x = 14334e5*cos(theta)

    y = 14334e5*sin(theta)

    Saturn.pos = [x,y,0]