vpython

dualstar1.py

# from https://www.glowscript.org/#/user/GlowScriptDemos/folder/Examples/program/BinaryStar-VPython/edit
import vpython as vp
#GlowScript 3.0 VPython
vp.scene.caption = """
In GlowScript programs:
To rotate "camera", drag with right button or Ctrl-drag.
To zoom, drag with middle button or Alt/Option depressed, or use scroll wheel.
  On a two-button mouse, middle is left + right.
To pan left/right and up/down, Shift-drag.
Touch screen: pinch/extend to zoom, swipe or two-finger rotate.
"""

vp.scene.forward = vp.vector(0,-.3,-1)

G = 6.7e-11 # Newton gravitational constant

giant = vp.sphere(pos=vp.vector(-1e11,0,0), radius=2e10, color=vp.color.red,
                make_trail=True, trail_type='points', interval=10, retain=500)
giant.mass = 2e30
giant.p = vp.vector(0, 0, -1e4) * giant.mass

dwarf = vp.sphere(pos=vp.vector(1.5e11,0,0), radius=1e10, color=vp.color.yellow,
                make_trail=True, interval=10, retain=500)
dwarf.mass = 1e30
dwarf.p = -giant.p

dt = 1e5
while True:
    vp.rate(200) # iterations per seconds (slow down computer if necessary)
    r = dwarf.pos - giant.pos  # relative vector
    # vp.mag does not work online?
    F = G * giant.mass * dwarf.mass * r.hat / mag(r)**2 # calculate the force 
    giant.p = giant.p + F*dt  # impuls
    dwarf.p = dwarf.p - F*dt  # impuls
    giant.pos = giant.pos + (giant.p/giant.mass) * dt # displacement vector : Euler-Chromer Calculation
    dwarf.pos = dwarf.pos + (dwarf.p/dwarf.mass) * dt