Hello all! I have attempted to model the three stars and one planet in the system Alpha Centauri for a class, and for some reason they don't want to orbit each other with actual values and I'm at a loss for what to do. Here's what I have so far:

Code:from visual import * from random import uniform from math import sqrt, pow #list constants a=1.49E11 G=6.67E-11 #Define objects. I used values for position, radius, mass, and color found online: starA = sphere( pos= (-1.5E11 ,0, 0) , radius= 8.533785E8, mass= 2.188E30, color =color.yellow) starB = sphere( pos= (1.5E11 ,0, 0) , radius= 6.016075E8, mass= 1.804E30, color =color.orange) starC = sphere( pos= (1.94E12, 0, 0) , radius= 9.80655E7, mass= 2.446593E29, color =color.magenta) planet = sphere(pos = (1.44e11, 0 , 0), radius = 3.5E6, mass = 5.9E24, material = materials.earth) planet.vel = vector(0 , -22000, 0) starA.vel = vector(0, 22000, 0) starB.vel = vector(0, -22000, 0) starC.vel = vector( 0, +22000,0) starA.trail = curve(color = starA.color) starB.trail = curve(color = starB.color) starC.trail = curve(color = starC.color) planet.trail = curve(color = color.white) # draw coordinate axes L=2E11 xaxis = curve(pos=[(0,0,0), (L,0,0)], color=color.white) yaxis = curve(pos=[(0,0,0), (0,L,0)], color=color.white) zaxis = curve(pos=[(0,0,0), (0,0,L)], color=color.white) #center of mass CM = sphere(pos=(0,0,0), radius = 2E12, color=color.green) #Algorithm dt = 1.0E8 t=0 dt2=1.0E3 while True: r1 = mag(starA.pos - starB.pos) r2=mag(starA.pos-starC.pos) r3=mag(starB.pos-starC.pos) r4=mag(starA.pos - planet.pos) r5=mag(starB.pos - planet.pos) r6=mag(starC.pos - planet.pos) starA.vel += ((-G*starB.mass*(starB.pos - starA.pos)/r1**3 - G*starC.mass*(starC.pos - starA.pos)/r2**3 - G*planet.mass*(starB.pos - planet.pos)/r4**3)/(starA.mass))*dt starB.vel += ((-G*starA.mass*(starA.pos - starB.pos)/r1**3 - G*starC.mass*(starC.pos - starB.pos)/r3**3 - G*planet.mass*(starB.pos - planet.pos)/r5**3)/(starB.mass))*dt starC.vel += ((-G*starB.mass*(starB.pos - starC.pos)/r3**3 - G*starA.mass*(starA.pos - starC.pos)/r2**3 - G*planet.mass*(starC.pos - planet.pos)/r6**3)/(starC.mass))*dt planet.vel += ((-G*starA.mass*(starA.pos - planet.pos)/r4**3 - G*starB.mass*(starB.pos - planet.pos)/r5**3 - G*starC.mass*(starC.pos - planet.pos))/r6**3/(planet.mass))*dt2 starA.pos += starA.vel*dt starB.pos += starB.vel*dt starC.pos += starC.vel*dt planet.pos += planet.vel*dt starA.trail.append(pos = starA.pos) starB.trail.append(pos = starB.pos) starC.trail.append(pos = starC.pos) planet.trail.append(pos = planet.pos) #M.pos = (starA.mass*starA.pos + starB.mass*starB.pos + starC.mass*starC.pos + planet.mass*planet.pos)/(starA.mass + starB.mass + starC.mass + planet.mass) rate(500) t += dt

Any help would be appreciated! Thank you!

Tweet This+ 1 thisPost To Linkedin