DIR Return Create A Forum - Home
---------------------------------------------------------
techno game dev forum
HTML https://technogamedevforum.createaforum.com
---------------------------------------------------------
*****************************************************
DIR Return to: General Discussion
*****************************************************
#Post#: 117--------------------------------------------------
Some things require more accuracy than a computer can give!
By: mr why Date: May 22, 2012, 1:33 am
---------------------------------------------------------
Here is an example and a challenge game for you.
Using gravity I have set up 4 equal heavenly bodies to rotate
about one another sharing one and the same double figure 8
orbit.
Can you get them to maintain that double-8 orbit?
It should be very easy! - There are only 4 numbers to adjust:
xc(1) = .1306 x-position start
xc(6) = -.0153 y-position start
xc(10) = -6.0005 x-speed start
xc(13) = 1.9245 y=speed start
Adjust these EVER-so-slightly and you may find what I am
missing!
I always find the computation errors are too big - the orbit
falls apart (beginning in ever-growing fashion)
To begin with all seems fine - but soon your hopes are dashed!
But the SINGLE figure-8 shared orbit (3 stars or planets sharing
one orbit) is DIFFERENT - for the laws of gravity forgive and
correct the computation errors!
Similarly if you provide a large central "sun" to rotate around,
everything works fine.
[code]
DEFDBL A-H, K-M, O-Z: DEFINT I-J, N
n = 16: OPTION BASE 1
DIM x(n), xc(n), f(n), c1(n), c2(n), c3(n), c4(n)
GOSUB graphics
m1 = 1#: m2 = 1#: m3 = 1#: m4 = 1
t = 0#
' ** initial positions **
xc(1) = .1306: xc(2) = 0: xc(3) = -xc(1): xc(4) = 0 'x coords
xc(5) = 0#: xc(6) = -.0153#: xc(7) = 0: xc(8) = -xc(6) 'y coords
' ** initial speeds **
xc(9) = 0: xc(10) = -6.0005: xc(11) = 0: xc(12) = -xc(10) 'x vel
xc(13) = 1.9245: xc(15) = -xc(13): xc(16) = 0 'y vel
h = .0000015#
agn:
GOSUB Runge
PSET (xc(1), xc(5)), 12
PSET (xc(2), xc(6)), 9
PSET (xc(3), xc(7)), 15
PSET (xc(4), xc(8)), 13
t = t + h
'IF j = 30000 THEN CLS : j = 0: GOSUB energy
GOTO agn
END
Equations:
d21 = x(2) - x(1): d32 = x(3) - x(2): d43 = x(4) - x(3): d14 =
x(1) - x(4)
d13 = x(3) - x(1): d24 = x(4) - x(2)
d65 = x(6) - x(5): d76 = x(7) - x(6): d87 = x(8) - x(7): d58 =
x(5) - x(8)
d75 = x(7) - x(5): d86 = x(8) - x(6)
r12 = (d21 ^ 2# + d65 ^ 2#) ^ .5#
r23 = (d32 ^ 2# + d76 ^ 2#) ^ .5#
r34 = (d43 ^ 2# + d87 ^ 2#) ^ .5#
r41 = (d14 ^ 2# + d58 ^ 2#) ^ .5#
r13 = (d13 ^ 2# + d75 ^ 2#) ^ .5#
r24 = (d24 ^ 2# + d86 ^ 2#) ^ .5#
p12 = r12 ^ 3#: p23 = r23 ^ 3#: p34 = r34 ^ 3: p41 = r41 ^ 3#
p31 = r13 ^ 3: p24 = r24 ^ 3
f(1) = x(9): f(2) = x(10): f(3) = x(11): f(4) = x(12)
f(5) = x(13): f(6) = x(14): f(7) = x(15): f(8) = x(16)
f(9) = m2 * d21 / p12 + m3 * d13 / p31 - m4 * d14 / p41 'x accn
of m1
f(10) = m3 * d32 / p23 + m4 * d24 / p24 - m1 * d21 / p12 '
m2
f(11) = m4 * d43 / p34 - m1 * d13 / p31 - m2 * d32 / p23 '
m3
f(12) = m1 * d14 / p41 - m2 * d24 / p24 - m3 * d43 / p34 '
m4
f(13) = m2 * d65 / p12 + m3 * d75 / p31 - m4 * d58 / p41 'y
m1
f(14) = m3 * d76 / p23 + m4 * d86 / p24 - m1 * d65 / p12 '
m2
f(15) = m4 * d87 / p34 - m1 * d75 / p31 - m2 * d76 / p23 '
m3
f(16) = m1 * d58 / p41 - m2 * d86 / p24 - m3 * d87 / p34 '
m4
RETURN
Runge:
FOR i = 1 TO n: x(i) = xc(i): NEXT
GOSUB Equations
FOR i = 1 TO n: c1(i) = h * f(i): NEXT
FOR i = 1 TO n: x(i) = xc(i) + c1(i) / 2#: NEXT
GOSUB Equations
FOR i = 1 TO n: c2(i) = h * f(i): NEXT
FOR i = 1 TO n: x(i) = xc(i) + c2(i) / 2#: NEXT
GOSUB Equations
FOR i = 1 TO n: c3(i) = h * f(i): NEXT
FOR i = 1 TO n: x(i) = xc(i) + c3(i): NEXT
GOSUB Equations
FOR i = 1 TO n: c4(i) = h * f(i): NEXT
FOR i = 1 TO n
xc(i) = xc(i) + (c1(i) + 2# * c2(i) + 2# * c3(i) + c4(i)) /
6#
NEXT
RETURN
graphics:
SCREEN 12
PAINT (1, 1), 9
xm = .15: ym = .15
VIEW (180, 17)-(595, 460), 0, 13
WINDOW (-xm, -ym)-(xm, ym)
xg = .134
FOR xi = -xg TO xg + .001 STEP xg / 20
LINE (-xm, xi)-(xm, xi), 8: LINE (xi, -ym)-(xi, ym), 8
NEXT xi
LINE (-xm, 0)-(xm, 0), 9: LINE (0, -ym)-(0, ym), 9
RETURN
energy:
' ** Is energy conserved? **
kin1 = .5# * m1 * (xc(7) ^ 2# + xc(10) ^ 2#)
kin2 = .5# * m2 * (xc(8) ^ 2# + xc(11) ^ 2#)
kin3 = .5# * m3 * (xc(9) ^ 2# + xc(12) ^ 2#)
pot = -(m1 * m2 / r12 + m2 * m3 / r23 + m3 * m1 / r13)
energy = kin1 + kin2 + kin3 + pot
LOCATE 21, 1: PRINT "Energy"
LOCATE 22, 1: PRINT energy
LOCATE 17, 7: PRINT " "
LOCATE 17, 1: PRINT "Time ="; CSNG(t)
RETURN
[/code]
*****************************************************