Program: Timex.bas
Implementation Select: Implementation
Operation Select: Program: TIMEX.BAS
DEFDBL A-H, M-Z
DEFINT I-K
' TIMEX.BAS
' Revision 1.0 Original 26 Sep 1997
' Revision 2.0 2 Simulation added, First improved 5 Mar 1998
'
DECLARE SUB GETSCREEN ()
DECLARE SUB MENU ()
DECLARE SUB DISPLAY1 ()
DECLARE SUB SETSTANDARD ()
DECLARE SUB INITIALISE2 ()
DECLARE SUB INITIALISE1 ()
CONST ESC = 27, ENTER = 13
CONST pi = 3.141592653589793# ' High Accuracy
DIM SHARED timedb(100), distsdb(100), vxdb(100)
DIM SHARED x0(10), x1(10), y0(10), y1(10), x2(10), y2(10)
DIM SHARED vx0(10), vx1(10), vy0(10), vy1(10)
DIM SHARED x00(10), y00(10) ' initial position
DIM SHARED vx00(10), vy00(10) ' initial velocity
DIM SHARED ax(10), ay(10)
DIM SHARED Sax(10), Say(10) ' sum of ax integrated
DIM SHARED Ssax(10), Ssay(10) ' sum of sum ax integrated
DIM SHARED Sva(10), Svy(10) ' sum of vx integrated
DIM SHARED rr(10, 10), dx(10, 10), dy(10, 10)
DIM SHARED m(10), m0(10), m00(10)
DIM SHARED pert(3), pert1(3), tmcnt(3), tmcnt1(3), mmax
DIM SHARED dist100, dist10, dist1 ' distance 1
DIM SHARED dist200, dist20, dist2 ' distance 2
DIM SHARED t, tcnt ' t + tdelta
DIM SHARED stars AS INTEGER
DIM SHARED tdelta00, tdelta0, tdelta!
DIM SHARED scren00%, scren0%, scren%
DIM SHARED G AS SINGLE
DIM SHARED xmax, ymax, xchrmax, ychrmax, fract AS SINGLE
DIM SHARED test%
DIM SHARED title$
DIM SHARED tmsrc AS INTEGER
DIM SHARED x0disp, y0disp, scalex, scaley
DIM SHARED dispcond00, dispcond0, dispcond AS INTEGER '- 1 once each rev
DIM SHARED endc AS INTEGER ' max number of endc
title$ = "3 Objects m0, m1 and m2 Demonstration"
colour = 13
SETSTANDARD 'set standard demo parameters.
s0:
DISPLAY1 'select test
s1:
INITIALISE2 ' Part 2
COLOR 15
CLS : SCREEN scren%
LOCATE 1, 25: PRINT title$
LOCATE 2, 1: PRINT "Test "; test%; "dtime "; tdelta!
stars = stars - 1
t = 0: tcnt = 0
end1 = 0
tmcnt(0) = 0: pert(0) = 0
FOR m = 0 TO mmax
m3 = 3 * m
FOR i = m3 TO m3 + stars
FOR j = m3 TO m3 + stars
IF i <> j AND m(j) > 0 THEN
dx = (x00(j) - x00(i))
dy = (y00(j) - y00(i))
rr = dx * dx + dy * dy:
rr(i, j) = rr: dx(i, j) = dx: dy(i, j) = dy
' PRINT "i"; i; j; r(i, j);
END IF
NEXT j
x1(i) = x0(i): y1(i) = y0(i) ' all equal
x2(i) = x1(i): y2(i) = y1(i)
NEXT i
NEXT m
FOR i = 0 TO 8
Sax(i) = 0: Ssax(i) = 0
Say(i) = 0: Ssay(i) = 0
NEXT i
' calculate for t = -1
tcnt = -1
FOR m = 0 TO mmax
m3 = 3 * m
FOR i = m3 TO m3 + stars
ax(i) = 0: ay(i) = 0
IF m(i) >= 0 THEN
FOR j = m3 TO m3 + stars
IF i <> j AND m(j) > 0 THEN
rr = rr(i, j): dx = dx(i, j): dy = dy(i, j)
a = m(j) * G / rr
ax(i) = ax(i) + dx * a / SQR(rr)
ay(i) = ay(i) + dy * a / SQR(rr)
END IF
NEXT j
Sax(i) = Sax(i) + ax(i)
Say(i) = Say(i) + ay(i)
Ssax(i) = Ssax(i) + Sax(i)
Ssay(i) = Ssay(i) + Say(i)
END IF
NEXT i
NEXT m
FOR m = 0 TO mmax
m3 = 3 * m
FOR i = m3 TO m3 + stars
x1(i) = x00(i) + tdelta! * (tcnt * vx00(i) + tdelta! * Ssax(i))
y1(i) = y00(i) + tdelta! * (tcnt * vy00(i) + tdelta! * Ssay(i))
FOR j = m3 TO m3 + stars
IF i <> j AND m(j) > 0 THEN
dx = (x00(j) - x00(i)) + tdelta! * (tcnt * (vx00(j) - vx00(i)) + tdelta! * (Ssax(j) - Ssax(i)))
dy = (y00(j) - y00(i)) + tdelta! * (tcnt * (vy00(j) - vy00(i)) + tdelta! * (Ssay(j) - Ssay(i)))
rr = dx * dx + dy * dy
rr(i, j) = rr: dx(i, j) = dx: dy(i, j) = dy
IF rr < 1 THEN end1 = 1
END IF
NEXT j
NEXT i
NEXT m
' initialise rr(i,j)
FOR m = 0 TO mmax
m3 = 3 * m
FOR i = m3 TO m3 + stars
FOR j = m3 TO m3 + stars
IF i <> j AND m(j) > 0 THEN
dx = (x00(j) - x00(i))
dy = (y00(j) - y00(i))
rr = dx * dx + dy * dy
rr(i, j) = rr: dx(i, j) = dx: dy(i, j) = dy
' PRINT "i"; i; j; rr(i, j);
END IF
NEXT j
NEXT i
NEXT m
FOR i = 0 TO 8
Sax(i) = 0: Ssax(i) = 0
Say(i) = 0: Ssay(i) = 0
NEXT i
t = 0: tcnt = 0: tmcnt(0) = 0: endc = 0
str1:
tcnt = tcnt + 1: t = tcnt * tdelta!
IF end1 = 0 THEN
FOR mm = 0 TO mmax
m3 = 3 * mm
FOR i = m3 TO m3 + stars
ax(i) = 0: ay(i) = 0
IF m(i) >= 0 THEN
FOR j = m3 TO m3 + stars
IF i <> j AND m(j) > 0 THEN
rr = rr(i, j)
dx = dx(i, j): dy = dy(i, j)
a = m(j) * G / rr
ax(i) = ax(i) + dx * a / SQR(rr)
ay(i) = ay(i) + dy * a / SQR(rr)
END IF
NEXT j
Sax(i) = Sax(i) + ax(i)
Say(i) = Say(i) + ay(i)
Ssax(i) = Ssax(i) + Sax(i)
Ssay(i) = Ssay(i) + Say(i)
END IF
NEXT i
' IF x1(1) < x1(0) THEN end1 = 1
FOR i = m3 TO m3 + stars
x2(i) = x1(i): y2(i) = y1(i)
x1(i) = x0(i): y1(i) = y0(i)
x0(i) = x00(i) + tdelta! * (tcnt * vx00(i) + tdelta! * Ssax(i))
y0(i) = y00(i) + tdelta! * (tcnt * vy00(i) + tdelta! * Ssay(i))
FOR j = m3 TO m3 + stars
IF i <> j AND m(j) > 0 THEN
dx = (x00(j) - x00(i)) + tdelta! * (tcnt * (vx00(j) - vx00(i)) + tdelta! * (Ssax(j) - Ssax(i)))
dy = (y00(j) - y00(i)) + tdelta! * (tcnt * (vy00(j) - vy00(i)) + tdelta! * (Ssay(j) - Ssay(i)))
rr = dx * dx + dy * dy
rr(i, j) = rr
dx(i, j) = dx: dy(i, j) = dy
IF rr < 1 THEN end1 = 1
END IF
IF rr(m3, m3 + 1) < 30 THEN end1 = 1
NEXT j
NEXT i
NEXT mm
END IF
' test and perform time based operations
IF tmsrc >= 0 THEN
FOR m = 0 TO mmax
m3 = 3 * m
tmsrc = 2
delta0 = y0(m3 + tmsrc) - y0(m3 + 1): delta1 = y1(m3 + tmsrc) - y1(m3 + 1)
IF delta0 > 0 AND delta1 < 0 THEN
tt = (tcnt - ABS(delta0) / (ABS(delta0) + ABS(delta1))) * tdelta!
COLOR colour + m
LOCATE ychrmax - m, 1: dt = tt - pert(m):
PRINT USING "rev ####"; tmcnt(m) + 1; : PRINT USING " time###.######"; dt; : PRINT " sec";
IF test% = 1 THEN
distsdb(tmcnt(0)) = x0(1)
vxdb(tmcnt(0)) = vx00(1) + tdelta! * Sax(1)
timedb(tmcnt(0)) = dt
END IF
pert(m) = tt: tmcnt(m) = tmcnt(m) + 1
' INPUT a$: IF UCASE$(a$) = "E" THEN END
END IF
tmsrc = 1
IF y0(m3 + tmsrc) > 0 AND y1(m3 + tmsrc) < 0 THEN
IF m = 0 THEN IF colour = 5 THEN colour = 13 ELSE colour = 5
tt = (tcnt - y0(tmsrc) / (y0(tmsrc) - y1(tmsrc))) * tdelta!
COLOR colour + m
LOCATE ychrmax - m, 34: dt = tt - pert1(m):
PRINT USING " rev ####"; tmcnt(m); tmcnt1(m) + 1; : PRINT USING " time####.######"; dt; : PRINT " sec";
pert1(m) = tt: tmcnt1(m) = tmcnt1(m) + 1
END IF
NEXT m
END IF
IF (tcnt MOD dispcond = 0 AND dispcond > 0) OR end1 = 1 THEN
COLOR 15: LOCATE 2, 45: PRINT "cycle "; tcnt;
FOR m = 0 TO mmax
m3 = 3 * m
PSET (x0disp + x0(m3) * scalex * fract, y0disp - y0(m3) * scaley), 12
PSET (x0disp + x0(m3 + 1) * scalex * fract, y0disp - y0(m3 + 1) * scaley), 10
IF stars > 1 THEN PSET (x0disp + x0(m3 + 2) * scalex * fract, y0disp - y0(m3 + 2) * scaley), colour + m
NEXT m
END IF
s9:
a$ = INKEY$: IF a$ = CHR$(ESC) THEN GOTO s10:
IF end1 = 1 THEN GOTO s10:
GOTO str1
s10:
COLOR 14
PRINT "Hit any key"
DO: a$ = INKEY$: LOOP UNTIL a$ <> ""
SCREEN 12
IF test% = 1 THEN
COLOR 14
PRINT TAB(5); "distance"; TAB(22); "rev time"; TAB(40); "v"
COLOR 15
FOR i = 0 TO tmcnt(0) - 1
PRINT i;
PRINT USING "###.######"; TAB(5); distsdb(i); TAB(20); timedb(i); TAB(35); -vxdb(i)
NEXT i
END IF
COLOR 14
LOCATE ychrmax, 35: PRINT "NEXT TEST ? ";
DO: a$ = INKEY$: LOOP UNTIL a$ <> ""
GOTO s0
SUB DISPLAY1
' DISPLAY1
SCREEN 12: CLS
d0:
COLOR 7
LOCATE 1, 15: PRINT title$
LOCATE 3, 18: PRINT "0 Parameter Selection Menu"
IF test% = 1 THEN COLOR 12 ELSE COLOR 7
LOCATE 4, 18: PRINT "1 m0=100000 m1=1000 m2=1 v=0"
IF test% = 2 THEN COLOR 12 ELSE COLOR 7
LOCATE 5, 18: PRINT "2 m0=1000000 m1=1000 m2=1 v=0"
ln = 10
COLOR 15
IF test% > 0 THEN LOCATE ln, 18: PRINT "S = Start Test"
COLOR 7
ln = ln + 1
LOCATE ln, 18: PRINT "Esc = Escape"
d01:
a$ = " "
LOCATE ln + 2, 18: PRINT "Which test ? ";
d1:
DO: a$ = INKEY$: LOOP UNTIL a$ <> ""
IF a$ = CHR$(ESC) THEN
END
END IF
IF a$ = CHR$(ENTER) THEN GOTO d2
IF UCASE$(a$) = "S" AND test% <> 0 THEN EXIT SUB
tst% = ASC(a$) - ASC("0")
IF (tst% < 0 OR tst% > 9) THEN GOTO d0
testn% = tst%
d2:
SELECT CASE testn%
CASE 1 TO 13
test% = testn%
INITIALISE1 ' Initialise part 1 Get distance and mass
CASE 0
MENU
CLS
END SELECT
GOTO d0
END SUB
SUB GETSCREEN
' GETSCREEN
SELECT CASE scren%
CASE 9
xmax = 640: ymax = 350: fract = 1.368
xchrmax = 80: ychrmax = 25
CASE 12
xmax = 640: ymax = 480: fract = 1
xchrmax = 80: ychrmax = 30
CASE ELSE
PRINT "SCREEN ERROR, not 9 or 12"
END SELECT
END SUB
SUB INITIALISE1
' INITIALISE1
SELECT CASE test%
CASE IS = 1
dist100 = 1000 ' distance between m0 and m1
dist200 = 8 ' distance between m1 and m2
m00(0) = 100000 ' mass of m0
m00(1) = 1000
m00(2) = 0
m0(3) = 0: m0(4) = m00(1): m0(5) = m00(2) ' mass
mmax = 1
CASE IS = 2
dist100 = 1000 / 2 ' distance between m0 and m1
dist200 = 6 / 2 ' distance between m1 and m2
m00(0) = 100000 ' mass of m0
m00(1) = 1000
m00(2) = 0
m0(3) = m00(0): m0(4) = m00(1): m0(5) = m00(2) ' mass
m0(6) = m00(0): m0(7) = m00(1): m0(8) = m00(2) ' mass
mmax = 2
CASE ELSE
END SELECT
dist10 = dist100: dist20 = dist200 ' distance 1 and 2
m0(0) = m00(0): m0(1) = m00(1): m0(2) = m00(2) ' mass
END SUB
SUB INITIALISE2
' part 2 of 2 INITIALISE2
G = 1
stars = 3 ' number of stars
tmsrc = -1 ' time save source (-1 = no)
scren% = scren0% ' screen mode
tdelta! = tdelta0 ' delta time
center = center0 ' origin body
xm = xm0 ' initial distance between m0 and m1
dispcond = dispcond0 ' display condition
FOR i = 0 TO 8
x0(i) = 0 ' initial positions
y0(i) = 0 ' initial positions
vx0(i) = 0 ' initial velocity
vy0(i) = 0 ' initial velocity
NEXT i
GETSCREEN ' Get screen parameters
dist1 = dist10: dist2 = dist20 ' distance 1 and 2
FOR i = 0 TO 8
m(i) = m0(i)
NEXT i
SELECT CASE test%
CASE IS = 1
x0(0) = 0: x0(1) = dist1: x0(2) = dist1 + dist2 + .00033 ' initial positions
x0(3) = 0: x0(4) = dist1: x0(5) = dist1 + dist2 ' initial positions
dist = x0(1) - x0(0)
m01 = m(0) + m(1)
tmsrc = 2 ' revolution source
vx0(0) = m(1) * SQR(2 / (dist * m01))
vx0(1) = -m(0) * SQR(2 / (dist * m01))
vx0(2) = vx0(1)
vy0(0) = 0: vy0(1) = 0: vy0(2) = SQR(m(1) / (x0(2) - x0(1)))
x0disp = 0: y0disp = ymax / 2
scalex = xmax / x0(2) / 1.5: scaley = scalex
vy0(5) = SQR(m(1) / dist2)
CASE IS = 2
dist1a = (dist1 * dist1 * dist1 / 4) ^ (1 / 3)
dist1b = dist1a
FOR m = 0 TO 2
SELECT CASE m
CASE IS = 0
ec = dist1 / dist1a - 1
dist = dist1a
CASE IS = 1
ec = 0:
dist = dist1b
CASE IS = 2
ec = 0:
dist = dist1
END SELECT
m3 = m * 3
x0 = dist1: x1 = x0
m01 = m(0) + m(1)
r0 = dist * m(1) / m01
r1 = dist - r0
vr0 = SQR(r0 * m(1)) / (dist)
vr1 = SQR(r1 * m(0)) / (dist)
x0 = x0 - r0 * (1 + ec)
x1 = x1 + r1 * (1 + ec)
ec1 = SQR((1 - ec) / (1 + ec)) ' eccentricity at aphelion
tmsrc = 1 ' revolution source
x0(0 + m3) = x0: x0(1 + m3) = x1:
x0(2 + m3) = x0(1 + m3) + dist2 ' initial positions
vy0(0 + m3) = -vr0 * ec1: vy0(1 + m3) = vr1 * ec1:
vy0(2 + m3) = vy0(1 + m3) + SQR(m(1) / dist2)
NEXT m
x0disp = 10: y0disp = ymax / 2
scalex = xmax / x0(2) / 1.5: scaley = scalex
scalex = xmax / (dist1 * 2) / 1.4: scaley = scalex
CASE ELSE
END SELECT
FOR i = 0 TO 10
y0(i) = 0
vx0(i) = 0
x00(i) = x0(i) ' initial positions
y00(i) = y0(i) ' initial positions
vx00(i) = vx0(i) ' initial velocity
vy00(i) = vy0(i) ' initial velocity
NEXT i
END SUB
SUB MENU
' MENU
SCREEN 12: CLS
sp$ = " "
COLOR 7
LOCATE 1, 25: PRINT title$
m0:
ln = 3
LOCATE ln, 28: PRINT "0 Test Selection Menu"
ln = ln + 1
LOCATE ln, 28: PRINT "1 Set Standard Parameters"
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "2 Screen mode"
IF scren0% = scren00% THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT scren0%; sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "3 Delta time"
IF tdelta0 = tdelta00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT tdelta0; sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "4 Mass m0 "
IF m0(0) = m00(0) THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT m0(0); sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "5 Mass m1 "
IF m0(1) = m00(1) THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT m0(1); sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "6 Mass m2 "
IF m0(2) = m00(2) THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT m0(2); sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "7 Distance 1 "
IF dist10 = dist100 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT dist10; sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "8 Distance 2 "
IF dist20 = dist200 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT dist20; sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "9 Display condition"
IF dispcond0 = dispcond00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT dispcond0; sp$
ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "Esc = Escape"
LOCATE ln + 2, 44: PRINT " "
LOCATE ln + 3, 28: PRINT " "
LOCATE ln + 2, 28: PRINT "Which parameter ? "
par% = 0
m1:
DO: a$ = INKEY$: LOOP UNTIL a$ <> ""
IF a$ = CHR$(ESC) THEN EXIT SUB
IF a$ = CHR$(ENTER) THEN GOTO m2
pr% = ASC(a$) - ASC("0")
IF (pr% < 0 OR pr% > 9) THEN GOTO m0
par% = par% * 10 + pr%
IF par% = 1 THEN : LOCATE ln + 2, 28 + 18: PRINT a$: GOTO m1
m2:
LOCATE ln + 4, 28: PRINT " ";
LOCATE ln + 2, 28: PRINT " "
SELECT CASE par%
CASE 0
EXIT SUB
CASE 1
SETSTANDARD
CASE 2
LOCATE ln + 3, 28: INPUT "Screen mode ", sm%
SELECT CASE sm%
CASE 7, 8, 9, 12
scren0% = sm%
CASE ELSE
LOCATE ln + 4, 28: PRINT "SCREEN ERROR, not 7, 8, 9 or 12";
END SELECT
CASE 3
LOCATE ln + 3, 28: INPUT "Delta time ", tdelta0
CASE 4
LOCATE ln + 3, 28: INPUT "Mass m0 ", m0(0)
CASE 5
LOCATE ln + 3, 28: INPUT "Mass m1 ", m0(1)
CASE 6
LOCATE ln + 3, 28: INPUT "Mass m2 ", m0(2)
CASE 7
LOCATE ln + 3, 28: INPUT "Distance 1 ", dist10
CASE 8
LOCATE ln + 3, 28: INPUT "Distance 2 ", dist20
CASE 9
LOCATE ln + 3, 28: INPUT "Display condition ", dispcond0
END SELECT
GOTO m0
END SUB
SUB SETSTANDARD
' SETSTANDARD
' This program set the standard parameters
tdelta00 = .0025#
'tdelta00 = .0005#
scren00% = 12 ' screen mode 9
dist100 = -1: dist200 = -1 ' initial distance 1,2
m00(0) = -1: m00(1) = -1: m00(2) = -1 ' m0,m1 and m2
dispcond00 = 20 ' display condition each 50 updates
tdelta0 = tdelta00 ' delta time
scren0% = scren00% ' screen mode
dist10 = 0: dist20 = 0 ' initial distance 1,2
m0(0) = 0: m0(1) = 0: m0(2) = 0 ' m0,m1 and m2
dispcond0 = dispcond00 ' display condition
INITIALISE1 ' Initialise part 1 Get distance and mass
END SUB
Back to my home page Contents of This Document