Program: Galaxy.bas

Implementation Select: Implementation
Operation Select: Program: GALAXY.BAS

DEFDBL A-H, M-Z
DEFINT I-K
'                       GALAXY.BAS
'       Revision 1.0    Original                         29 Feb 1996
'       Revision 1.1    Initialisation endc               5 Mar 1996
'       Revision 1.2    Correction print rev time        24 Sep 1996
'       Revision 2      Added type 2 and 3               19 Oct 2000
'       Revision 2.1    Added type 4 SO-2                23 Oct 2000
'
DECLARE SUB COORD (m0#, m1#, dist#, x0#, x1#, vy0#, vy1#)
DECLARE SUB ARROW (x1d#, y1d#, x2d#, y2d#)
DECLARE SUB GETSCREEN ()
DECLARE SUB MENU ()
DECLARE SUB DISPLAY1 ()
DECLARE SUB SETSTANDARD ()
DECLARE SUB INITIALISE1 ()
DECLARE SUB INITIALISE2 ()

CONST ESC = 27, ENTER = 13
CONST pi = 3.141592653589#
CONST year = 365.25# * 24 * 60 * 60
CONST century = year * 100
CONST lightyear = 9460000000000#
CONST rbulge = 5000 * lightyear

DIM SHARED x0, x1, y0, y1, vx0t0, vx1t0, vy0t0, vy1t0
DIM SHARED m0, m1, mbhole, mbhole0, mbhole00
DIM SHARED G, Gm0, Gm1
DIM SHARED tt1
DIM SHARED x00, x01, x02, x10, x11, x12
DIM SHARED y00, y01, y02, y10, y11, y12
DIM SHARED vx00, vx01, vx02, vx10, vx11, vx12
DIM SHARED vy00, vy01, vy02, vy10, vy11, vy12
DIM SHARED dx0, dy0, dx1, dy1
DIM SHARED mm01
DIM SHARED d, d0, d00                      ' distance sun galaxy
DIM SHARED stars AS INTEGER
DIM SHARED rbase                            'v origin
DIM SHARED phi                              'phi origin
DIM SHARED t, tdelta00, tdelta0, tdelta
DIM SHARED vinit, vinit0, vinit00
DIM SHARED scren%
DIM SHARED xmax, ymax, xchrmax, ychrmax, fract AS SINGLE
DIM SHARED testn%, test%
DIM SHARED c                              ' speed of light
DIM SHARED title$
DIM SHARED x0disp, y0disp, scalex, scaley
DIM SHARED tcnt                                    ' revolution counter
DIM SHARED dispcond         ' - 1 once each rev
DIM SHARED initcond         ' initial conditions

SCREEN 12

title$ = "Galaxy Demonstration"
     
SETSTANDARD                            'set standard demo parameters.

s0:
DISPLAY1                               'select test

INITIALISE2                            ' part 2
stars = stars - 1
 
  f1 = 20000     ' scale up factor display
  deltaf = 100   ' scale down factor delta
  v100 = 250     ' v rotation init
color12 = 12

CLS
COLOR 15
LOCATE 1, 25: PRINT title$
COLOR 7
LOCATE 2, 1
PRINT USING "v ####"; rbase;
PRINT USING " phi ###.#"; phi;
PRINT USING "dtime in ly .##########"; TAB(22); tdelta / lightyear;
PRINT TAB(46); "  init "; initcond;
PRINT USING "Bhole ########"; TAB(60); mbhole0;
PRINT " m0"
LOCATE 3, 1: PRINT "# stars "; stars + 1; "c "; c
LOCATE ychrmax, 1: PRINT USING "v ####.#####"; vy1t0; : PRINT "km/sec";
c2 = c * c: Gm0 = G * m0: Gm1 = G * m1:
prstate = 0: printinit = 1: end1 = 0: a1dyo = 1: endc = 0
printcnt = 0

t = tdelta

x00 = x0: y00 = y0: x10 = x1: y10 = y1
dx = x1 - x0: dy = y1 - y0: rr = dx * dx + dy * dy: r = SQR(rr):
IF c = 0 THEN t0 = 0: t1 = 0 ELSE t0 = r / c: t1 = r / c
rrt0 = rr

' Calculate position at t=0
dx0 = x1 - x0: dy0 = y1 - y0
rr = dx0 * dx0 + dy0 * dy0
r = SQR(rr): rrr = r * r * r

r1 = r: r11 = r: r12 = r

  vx00 = vx0t0: vy00 = vy0t0
  vx10 = vx1t0: vy10 = vy1t0

tcnt = 0
t = tdelta

time1 = TIMER
rr1max = SQR(rrt0): rr1min = rr1max

str1:
tcnt = tcnt + 1: tt1 = tcnt * tdelta

  dx0 = x10 - x00: dy0 = y10 - y00

rr0 = dx0 * dx0 + dy0 * dy0                            'virtual distance
r0 = SQR(rr0)
rrr = rr0 * r0
IF c = 0 THEN t0 = 0 ELSE t0 = r0 / c

dx1 = x00 - x10: dy1 = y00 - y10

rr1 = dx1 * dx1 + dy1 * dy1                            'virtual distance
r1 = SQR(rr1)
rrr = rr1 * r1
rt0 = r1
IF c = 0 THEN t1 = 0 ELSE t1 = r1 / c

IF test% = 1 THEN
  m00 = m0
ELSE
  IF r1 < 60000 * lightyear THEN
    m00 = r1 * v100 * v100 / G
  ELSE
    m00 = 60000 * lightyear * v100 * v100 / G
  END IF
END IF
IF r1 < rbulge THEN
  v101 = v100 * r1 / rbulge
  m00 = r1 * v101 * v101 / G
ELSE
  v101 = v100
END IF
IF test% = 3 OR test% = 4 THEN m00 = m00 + mbhole

IF r1 < 5 * lightyear THEN
  ' LOCATE 10, 10: PRINT m00, mbhole
  IF tsave <> tdelta THEN
     tdelta = tdelta / deltaf
     t = tdelta
     tsave = t
     COLOR 15: LOCATE 2, 22
     IF tdelta > 1000000 THEN
       PRINT USING "dtime in ly .##########"; TAB(22); tdelta / lightyear;
     ELSE
       PRINT USING " dtime in secs ########"; TAB(22); tdelta;
     END IF
  END IF
 
  xds = x10 - x00: yds = y10 - y00
  ' LOCATE 10, 10: PRINT x0disp + xds * scalex * f1; xmax
  IF ABS(x0disp + xds * scalex * f1) < xmax AND ABS(y0disp - yds * scaley * f1 / fract - 20) < ymax THEN
  isetx1 = x0disp + xds * scalex * f1
  isety1 = y0disp - yds * scaley * f1 / fract - 20
 
  IF (isetx1 <> isetx1old OR isety1 <> isety1old) THEN
    isetx1old = isetx1: isety1old = isety1
    LOCATE 4, 2
    PRINT "Mass ";
    PRINT USING "################"; m00 / m1;
    PRINT ; " m0   ";
    IF v101 = v100 THEN COLOR 7
    PRINT USING "v rotation ####.##"; v101;
    PRINT ; " km/sec    ";
    COLOR 15
    PRINT USING "speed #####.#"; v10;
    PRINT " km/sec  "
    PRINT "distance in ";
    PRINT USING "###.###"; r1 / lightyear;
    PRINT " ly  ";
    PRINT USING "################"; r1;
    PRINT " km";
    COLOR 7
    PRINT USING " v init ##.####"; vinit0;
    PRINT " km/sec"
    'INPUT a$
    'IF UCASE$(a$) = "E" THEN END
    xds = 0: yds = 0
    COLOR 15
    PSET (x0disp + xds * scalex * f1, y0disp - yds * scaley * f1 / fract - 20)
    xds = x10 - x00: yds = y10 - y00
    COLOR 14
    PSET (isetx1, isety1)
  END IF
  END IF
ELSE
  IF tsave = tdelta THEN
     tdelta = tdelta * deltaf
     t = tdelta
     COLOR 7: LOCATE 2, 22
     PRINT USING "dtime in ly .##########"; TAB(22); tdelta / lightyear;
  END IF
END IF

a1 = G * m00 / rrr     'Gm0
a1dx = a1 * dx1: a1dy = a1 * dy1

   vx11 = vx10 + a1dx * t
   vy11 = vy10 + a1dy * t
   x11 = x10 + vx11 * t
   y11 = y10 + vy11 * t
   x10 = x11: y10 = y11
   vx10 = vx11: vy10 = vy11
   v10 = SQR(vx10 * vx10 + vy10 * vy10)

' test and perform time based operations
  COLOR 15
  IF a1dyo * a1dy < 0 AND a1dx < 0 THEN
    LOCATE ychrmax - 1, 1: PRINT "cycle "; tcnt; tcnt - tcntold; "  time";
    tcntold = tcnt
    IF tt1 / year < 1000000 THEN
       PRINT USING "########"; tt1 / year;
       PRINT " year  ";
    ELSE
       tt11 = tt1 / (year * 1000000):
       PRINT USING "########.##"; tt11;
       PRINT " million year ";
    END IF
    PRINT USING " dist #######.#"; rt0 / lightyear;
    PRINT " ly";
    endc = endc + 1
    IF color12 = 12 THEN color12 = 4 ELSE color12 = 12
  END IF
  a1dyo = a1dy               ' old value


IF dispcond > 0 THEN
  IF tcnt MOD dispcond = 0 THEN prtest% = 1
END IF
 
xds = x10 - x00: yds = y10 - y00
isetx3 = x0disp + xds * scalex
isety3 = y0disp - yds * scaley / fract
IF (isetx3 <> isetx3old OR isety3 <> isety3old) THEN
  isetx3old = isetx3: isety3old = isety3

'IF prtest% = 1 THEN
  COLOR 15:
  LOCATE 3, 25: PRINT "cycle "; tcnt; "rev "; endc; "dist";
  LOCATE 3, 50
  PRINT USING "#######.#"; rt0 / lightyear;
  PRINT " ly "
  IF test% > 1 THEN
    LOCATE 4, 2
    PRINT "Mass ";
    PRINT USING "################"; m00 / m1;
    PRINT ; " m0   ";
    IF v101 = v100 THEN COLOR 7
    PRINT USING "v rotation ####.##"; v101;
    PRINT ; " km/sec    ";
    COLOR 15
    PRINT USING "speed #####.#"; v10;
    PRINT " km/sec  "
  END IF
  xds = 0: yds = 0
  COLOR 4
  xxr = 50 * COS(phi * pi / 180): yyr = 50 * SIN(phi * pi / 180)
  IF rbase > 0 THEN ARROW xds + xxr * .05, yds + yyr * .05, xds + xxr * .3, yds + yyr * .3
  COLOR 15
  PSET (x0disp + xds * scalex, y0disp - yds * scaley / fract)
  COLOR color12
  IF test% = 4 AND tt1 / year MOD 2 = 1 THEN COLOR 15
  xds = x10 - x00: yds = y10 - y00
  PSET (x0disp + xds * scalex, y0disp - yds * scaley / fract)
  IF test% > 1 THEN
    xds = 0: yds = 0
    COLOR 15
    PSET (x0disp + xds * scalex * f1, y0disp - yds * scaley * f1 / fract - 20)
    xds = x10 - x00: yds = y10 - y00
    COLOR 14
      IF r1 < 100 THEN
      PSET (x0disp + xds * scalex * f1, y0disp - yds * scaley * f1 / fract - 20)
    END IF
  END IF
  prtest% = 0
END IF
 
  IF test% > 1 THEN
    IF r1 > r11 AND r11 < r12 THEN
      accu = (r12 - r11) / r11
      IF accu > 1 THEN PRINT "accuracy error "; accu
      dx1 = x00 - x10: dy1 = y00 - y10
      rr1 = dx1 * dx1 + dy1 * dy1                            'virtual distance
      r1next = SQR(rr1)
      LOCATE 6, 1
      PRINT "smallest distance in ly";
      PRINT USING "####.####"; r1next / lightyear; r1 / lightyear; r11 / lightyear; r12 / lightyear; r13 / lightyear
      PRINT "Speed ";
      PRINT USING "######.#"; v10; v11; v12; v13;
      INPUT a$
      IF UCASE$(a$) = "E" THEN END
      LOCATE 6, 1: PRINT SPACE$(70)
      LOCATE 7, 1: PRINT SPACE$(70)
      LOCATE 8, 1: PRINT SPACE$(70)
    ELSE
      'IF r1 < .02 * lightyear THEN INPUT a$
    END IF
    r13 = r12: r12 = r11: r11 = r1
    v13 = v12: v12 = v11: v11 = v10
  END IF
a$ = INKEY$: IF a$ = CHR$(ESC) THEN end1 = 1: GOTO s10:
GOTO str1

s10:
 
  COLOR 14
  LOCATE ychrmax, 20: PRINT "NEXT TEST ? ";
  DO: a$ = INKEY$: LOOP UNTIL a$ <> ""

GOTO s0

SUB ARROW (x1d, y1d, x2d, y2d)
'                                                                      ARROW
bas = 4
ddx = x2d - x1d: ddy = y2d - y1d
d = SQR(ddx * ddx + ddy * ddy)
IF d = 0 THEN EXIT SUB
y3d = y2d - (y2d - y1d) * bas / d
x3d = x2d - (x2d - x1d) * bas / d
y4d = y3d + (x2d - x1d) * bas / d
x4d = x3d - (y2d - y1d) * bas / d
y5d = y3d - (x2d - x1d) * bas / d
x5d = x3d + (y2d - y1d) * bas / d
LINE (x0disp + x1d * fract, y0disp - y1d)-(x0disp + x2d * fract, y0disp - y2d)
LINE (x0disp + x2d * fract, y0disp - y2d)-(x0disp + x4d * fract, y0disp - y4d)
LINE (x0disp + x2d * fract, y0disp - y2d)-(x0disp + x5d * fract, y0disp - y5d)
END SUB

SUB COORD (m0, m1, dist, x0, x1, vy0, vy1)
'                                                                     COORD
' x0 = initial position
' vy0 = initial speed
mm01 = m0 + m1
r0 = dist * m1 / mm01
r1 = dist - r0
'PRINT m0, m1, mm01, dist, r0, r1
vr0 = SQR(r0 * m1 * G) / dist
vr1 = SQR(r1 * m0 * G) / dist
x1 = x0                             ' same base
vy1 = vy0                           ' same base
x0 = x0 - r0
x1 = x1 + r1
vy0 = vy0 - vr0
vy1 = vy1 + vr1
END SUB

SUB DISPLAY1
'                                                                   DISPLAY1
CLS

d0:
COLOR 15
LOCATE 3, 28: PRINT "Main Selection Display"
ln = 5
COLOR 7
ln = ln + 1
LOCATE ln, 28: PRINT "0 = Parmeter Selection Menu"
ln = ln + 1
IF test% = 1 THEN COLOR 12 ELSE COLOR 7
LOCATE ln, 28: PRINT "1 = Galaxy and Sun simulation"
ln = ln + 1
IF test% = 2 THEN COLOR 12 ELSE COLOR 7
LOCATE ln, 28: PRINT "2 = Galaxy and Star simulation"
ln = ln + 1
IF test% = 3 THEN COLOR 12 ELSE COLOR 7
LOCATE ln, 28: PRINT "3 = Galaxy, Black hole and Star simulation"
ln = ln + 1
IF test% = 4 THEN COLOR 12 ELSE COLOR 7
LOCATE ln, 28: PRINT "4 = Galaxy, Black hole and Star SO-2 simulation "

COLOR 15
ln = ln + 2:
LOCATE ln, 28: PRINT "S = Start Test "
COLOR 7
ln = ln + 1
LOCATE ln, 28: PRINT "Esc = Escape "
lns = ln + 1                                 ' line for error messages

LOCATE ln + 2, 28: PRINT "Which test ?     "
 
  DO: a$ = INKEY$: LOOP UNTIL a$ <> ""
  LOCATE lns, 1: PRINT SPACE$(50): LOCATE lns, 1
  IF a$ = CHR$(ESC) THEN END
 
 
  IF UCASE$(a$) = "S" THEN EXIT SUB
  tst% = ASC(a$) - ASC("0")
  IF (tst% < 0 OR tst% > 9) THEN GOTO d0
  testn% = tst%

SELECT CASE testn%
CASE 0
  MENU
  CLS
CASE 1 TO 4
  test% = testn%
  INITIALISE1
END SELECT
GOTO d0

END SUB

SUB GETSCREEN
'                                                                   GETSCREEN
   xmax = 640: ymax = 480: fract = 1
   xchrmax = 80: ychrmax = 30

END SUB

SUB INITIALISE1
  SELECT CASE test%
  CASE 1
    d0 = d00                             ' distance in lightyear
    tdelta0 = tdelta00
    mbhole0 = 0
    vinit0 = v100
  CASE 2 TO 3
    d0 = 50000                           ' distance in lightyear
    tdelta0 = tdelta00 / 100
    vinit0 = vinit00                     ' Initial speed
    mbhole0 = 0
    IF test% = 3 THEN mbhole0 = mbhole00
  CASE 4
    d0 = .028#                            ' distance in lightyear
    tdelta0 = .0000001
    vinit0 = 260                         ' Initial speed
    mbhole0 = mbhole00
  END SELECT

END SUB

SUB INITIALISE2
'                                                                 INITIALISE2
stars = 2                              ' number of stars

G = 6.67D-20

tdelta = tdelta0 * lightyear

r0 = 0: r1 = 0
v0 = 0: v1 = 0

GETSCREEN

    m0 = 1.9891D+30 * 110000000000#         ' Galaxy
    m1 = 1.9891D+30                         ' Sun = m0
    mbhole = mbhole0 * m1                   ' m blackhole
   
    d = d0 * lightyear                      ' 25000 light year
    COORD m0, m1, d, r0, r1, v0, v1
    x0disp = xmax / 2: y0disp = ymax / 2
    scalex = xmax / d * .28: scaley = scalex

vzx = rbase * COS(phi * pi / 180)
vzy = rbase * SIN(phi * pi / 180)
ec = 0
ec1 = 1: ec2 = 0

x0 = r0 * (1 + ec2): y0 = 0: vx0t0 = vzx: vy0t0 = vzy + ec1 * v0
x1 = r1 * (1 + ec2): y1 = 0: vx1t0 = vzx: vy1t0 = vzy + ec1 * v1
IF initcond = 1 THEN
  ddx = x1 - x0: ddy = y1 - y0
  IF c = 0 THEN tx = 0 ELSE tx = SQR(ddx * ddx + ddy * ddy) / c
  ddx = vzx * tx: ddy = vzy * tx
  vy0t0 = vzy + ec1 * v0
  vy1t0 = vzy + ec1 * v1
  x0 = x0 + ddx / 2
  x1 = x1 - ddx / 2
  y0 = y0 + ddy / 2
  y1 = y1 - ddy / 2
END IF
IF initcond = 2 THEN
  ddx = x1 - x0: ddy = y1 - y0
  IF c = 0 THEN tx = 0 ELSE tx = SQR(ddx * ddx + ddy * ddy) / c
  ddx = vzx * tx: ddy = vzy * tx
  ecy = vzx / c
  vy0t0 = vzy + ec1 * v0 * SQR((1 + ecy) / (1 - ecy))
  vy1t0 = vzy + ec1 * v1 * SQR((1 - ecy) / (1 + ecy))
  ecx = vzy / c
  vx0t0 = vzx * SQR((1 + ecx) / (1 - ecx))
  vx1t0 = vzx * SQR((1 - ecx) / (1 + ecx))
END IF

IF test% >= 2 THEN
  vy1t0 = vinit0
END IF

END SUB

SUB MENU
'                                                                     MENU
CLS
SCREEN 12
sp$ = "        "
COLOR 7
LOCATE 1, 28: PRINT "Parameter Selection Display"

m0:
ln = 3
LOCATE ln, 28: PRINT "0 Main Selection Menu"

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "1 Delta time in ly"
IF tdelta0 = tdelta00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 61: PRINT USING ".########"; tdelta0

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "2 V origin (Galaxy)"
LOCATE ln, 60: PRINT rbase; sp$

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "3 Phi origin"
LOCATE ln, 60: PRINT phi; sp$

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "4 Display condition"
LOCATE ln, 60: PRINT dispcond; sp$


ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "5 Initial condition"
LOCATE ln, 60: PRINT initcond

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "6 Speed of gravity"
LOCATE ln, 60: PRINT c; sp$

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "7 Initial distance in ly"
IF d0 = d00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT d0; sp$

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "8 Initial radial speed of star"
IF vinit0 = vinit00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT vinit0; sp$

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "9 Mass Blackhole in m0"
IF mbhole0 = mbhole00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT mbhole0; sp$

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "10 Delta time in year"
IF tdelta0 = tdelta00 THEN COLOR 7 ELSE COLOR 12
LOCATE ln, 60: PRINT USING "#######.##"; tdelta0 * lightyear / year

ln = ln + 1
COLOR 7: LOCATE ln, 28: PRINT "Esc = Escape"

LOCATE ln + 2, 44: PRINT "   "                         ' clear parameter
LOCATE ln + 3, 28: PRINT "                                " ' clear selection
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 OR par% = 2 THEN : LOCATE ln + 2, 28 + 18: PRINT a$: GOTO m1

m2:

LOCATE ln + 4, 28: PRINT "                                "; ' clear error
LOCATE ln + 2, 28: PRINT "                    "

SELECT CASE par%
CASE 0
  EXIT SUB
CASE 1
  LOCATE ln + 3, 28: INPUT "Delta Time in ly ", tdelta0
CASE 2
  LOCATE ln + 3, 28: INPUT "V origin (Galaxy) ", rbase
CASE 3
  LOCATE ln + 3, 28: INPUT "Phi origin ", phi
CASE 4
  LOCATE ln + 3, 28: INPUT "Display condition ", dispcond
CASE 5
  LOCATE ln + 3, 28: INPUT "Initial condition ", initcond
CASE 6
  LOCATE ln + 3, 28: INPUT "Speed of gravity ", c
CASE 7
  LOCATE ln + 3, 28: INPUT "Initial distance in ly ", d0
CASE 8
  LOCATE ln + 3, 28: INPUT "Initial radial speed of star ", vinit0
CASE 9
  LOCATE ln + 3, 28: INPUT "Mass Blackhole in m0 ", mbhole0
CASE 10
  LOCATE ln + 3, 28: INPUT "Delta Time in year ", tdeltax
  tdelta0 = tdeltax * year / lightyear
END SELECT
GOTO m0

END SUB

SUB SETSTANDARD
 
  initcond = 1                         ' initial conditions 1 = yes
  c = 300000                           ' speed of light
  phi = 0                              ' phi origin
  rbase = 0                            ' v origin
  dispcond = 100                       ' display condition
  d00 = 25000#                         ' distance in lightyear
  tdelta00 = .05#                      ' delta time in ly
  testn% = 1: test% = testn%
  vinit00 = .0015#                     ' initial speed
  mbhole00 = 2600000#                    ' m blackhole  in m0

  tdelta0 = tdelta00                   ' delta time
  d0 = d00                             ' distance in lightyear
  vinit0 = vinit00
  mbhole0 = mbhole00                  ' m blackhole  in m0

END SUB


Back to my home page: Contents of This Document