Author Topic: The stable orbits game  (Read 4278 times)

0 Members and 1 Guest are viewing this topic.

Offline luke

  • Administrator
  • Seasoned Forum Regular
  • Posts: 324
    • View Profile
The stable orbits game
« on: March 17, 2020, 08:25:58 am »
Of course, "game" is a bit of a stretch. I vaguely recall someone (probably @STxAxTIC ) mentioned Newtonian simulations a little while back on discord, then I came across this while looking for something else and figured I'd post it.

It's a pretty simple n-body simulator. I write this a while back, and I don't even remember if the physics is correct. I suppose it shows a use of the WINDOW command?

Set the first line to FALSE to ignore collisions.
Code: QB64: [Select]
  1. $LET COLLISIONS = TRUE
  2. TYPE particle_t
  3.     'position
  4.     px AS SINGLE
  5.     py AS SINGLE
  6.     'velocity
  7.     vx AS SINGLE
  8.     vy AS SINGLE
  9.     'mass
  10.     m AS SINGLE
  11.     'colour
  12.     c AS _UNSIGNED LONG
  13. DIM initial(0 TO 3) AS particle_t
  14. DIM particles(0 TO 3) AS particle_t
  15.  
  16.  
  17. initial(0).px = 320
  18. initial(0).py = 320
  19. initial(0).m = 100000
  20.  
  21. initial(0).c = _RGB32(255, 255, 255)
  22.  
  23. initial(1).px = 320
  24. initial(1).py = 130
  25. initial(1).m = 10
  26. initial(1).vx = 2.0
  27. initial(1).c = _RGB32(255, 0, 0)
  28.  
  29. initial(2).px = 320
  30. initial(2).py = 480
  31. initial(2).m = 20
  32. initial(2).vx = -2
  33. initial(2).c = _RGB32(0, 255, 0)
  34.  
  35. initial(3).px = 10
  36. initial(3).py = 0
  37. initial(3).m = 2
  38. initial(3).vx = 1.5
  39. initial(3).c = _RGB32(0, 0, 255)
  40. G = 0.01
  41. gscreen = _NEWIMAGE(640, 640, 32)
  42.  
  43. start:
  44. PRINT SPACE$(40 - 22 / 2); "The Stable Orbits Game"
  45. PRINT "Four bodies exist in 2D space. Your challenge is to define parameters for them"
  46. PRINT "(mass, initial position and initial velocity) so that they can orbit for as long"
  47. PRINT "as possible before they collide. You can also specify a value for G, the"
  48. PRINT "gravitational constant."
  49. PRINT "Use the arrow keys to move between the boxes, and numeric keys and the enter key"
  50. PRINT "to input values."
  51. PRINT "Press C to copy setup to clipboard, and P to paste setup from clipboard."
  52. PRINT "Press F5 to run the simulation, and Escape to exit simulation. Good luck!"
  53.  
  54. GOSUB print_table
  55.  
  56. selx = 17
  57. sely = 13
  58.     FOR p = 0 TO 3
  59.         _PRINTSTRING (16 * p + 17, 13), STR$(initial(p).px)
  60.         _PRINTSTRING (16 * p + 17, 14), STR$(initial(p).py)
  61.         _PRINTSTRING (16 * p + 17, 15), STR$(initial(p).vx)
  62.         _PRINTSTRING (16 * p + 17, 16), STR$(initial(p).vy)
  63.         _PRINTSTRING (16 * p + 17, 17), STR$(initial(p).m)
  64.     NEXT p
  65.     _PRINTSTRING (24, 19), STR$(G)
  66.     LOCATE sely, selx, 1
  67.     k$ = ""
  68.     DO
  69.         k = _KEYHIT
  70.         _LIMIT 20
  71.     LOOP UNTIL k > 0
  72.     SELECT CASE k
  73.         CASE 18432 'Up
  74.             IF seld > 0 THEN seld = seld - 1
  75.         CASE 19200 'Left
  76.             IF selp > 0 THEN selp = selp - 1
  77.         CASE 20480 'Down
  78.             IF seld < 5 THEN seld = seld + 1
  79.         CASE 19712 'Right
  80.             IF selp < 3 THEN selp = selp + 1
  81.         CASE 16128 'F5
  82.             GOTO sim
  83.         CASE 67, 99 'C
  84.             FOR p = 0 TO 3
  85.                 summary$ = summary$ + MKS$(initial(p).px) + MKS$(initial(p).py) + MKS$(initial(p).vx) + MKS$(initial(p).vy) + MKS$(initial(p).m)
  86.             NEXT p
  87.             summary$ = summary$ + MKS$(G)
  88.             _CLIPBOARD$ = encodeBASE64$(summary$)
  89.         CASE 80, 112 'P
  90.             summary$ = decodeBASE64$(_CLIPBOARD$)
  91.             FOR p = 0 TO 3
  92.                 initial(p).px = CVS(MID$(summary$, p * 20 + 1, 4))
  93.                 initial(p).py = CVS(MID$(summary$, p * 20 + 5, 4))
  94.                 initial(p).vx = CVS(MID$(summary$, p * 20 + 9, 4))
  95.                 initial(p).vy = CVS(MID$(summary$, p * 20 + 13, 4))
  96.                 initial(p).m = CVS(MID$(summary$, p * 20 + 17, 4))
  97.                 G = CVS(RIGHT$(summary$, 4))
  98.             NEXT p
  99.  
  100.         CASE ELSE
  101.             _PRINTSTRING (selx, sely), SPACE$(15)
  102.             LOCATE , POS(0) + 1
  103.             INPUT "", newval
  104.             _KEYCLEAR
  105.             _PRINTSTRING (selx, sely), SPACE$(15)
  106.  
  107.             IF seld = 5 THEN
  108.                 G = newval
  109.             ELSE
  110.                 SELECT CASE seld
  111.                     CASE 0: initial(selp).px = newval
  112.                     CASE 1: initial(selp).py = newval
  113.                     CASE 2: initial(selp).vx = newval
  114.                     CASE 3: initial(selp).vy = newval
  115.                     CASE 4: initial(selp).m = newval
  116.                 END SELECT
  117.             END IF
  118.     END SELECT
  119.     IF seld = 5 THEN
  120.         selx = 24
  121.         sely = 19
  122.     ELSE
  123.         selx = 16 * selp + 17
  124.         sely = 13 + seld
  125.     END IF
  126.  
  127.  
  128. sim:
  129. SCREEN gscreen
  130.  
  131. FOR p = 0 TO 3
  132.     particles(p) = initial(p)
  133.     'sanity check
  134.     IF particles(p).m <= 0 THEN
  135.         PRINT "Simulation error: particle"; p; "has non-positive mass."
  136.         SLEEP
  137.         _KEYCLEAR
  138.         GOTO start
  139.     END IF
  140.  
  141.  
  142.         CASE -27
  143.             GOSUB reset_sim
  144.             GOTO start
  145.  
  146.         CASE 16896
  147.             SCREEN 0
  148.             _AUTODISPLAY
  149.             PRINT SPACE$(40 - 12 / 2); "Current data"
  150.             PRINT "Press any key to return to simulation"
  151.             PRINT
  152.             GOSUB print_table
  153.             FOR p = 0 TO 3
  154.                 _PRINTSTRING (16 * p + 17, 7), STR$(particles(p).px)
  155.                 _PRINTSTRING (16 * p + 17, 8), STR$(particles(p).py)
  156.                 _PRINTSTRING (16 * p + 17, 9), STR$(particles(p).vx)
  157.                 _PRINTSTRING (16 * p + 17, 10), STR$(particles(p).vy)
  158.                 _PRINTSTRING (16 * p + 17, 11), STR$(particles(p).m)
  159.             NEXT p
  160.             _PRINTSTRING (24, 13), STR$(G)
  161.  
  162.             SLEEP
  163.             SCREEN gscreen
  164.     END SELECT
  165.     CLS
  166.     FOR p = 0 TO UBOUND(particles)
  167.         ax = 0
  168.         ay = 0
  169.         FOR p2 = 0 TO UBOUND(particles)
  170.             IF p2 <> p THEN 'avoid calculating when both particles are the same one
  171.                 'magnitude of force
  172.                 f = G * particles(p2).m * particles(p).m / ((particles(p2).px - particles(p).px) ^ 2 + (particles(p2).py - particles(p).py) ^ 2)
  173.                 'direction, taking p as reference particle
  174.                 ang = _ATAN2(-(particles(p2).py - particles(p).py), (particles(p2).px - particles(p).px))
  175.                 ax = ax + f * COS(ang) / particles(p).m
  176.                 ay = ay - f * SIN(ang) / particles(p).m
  177.             END IF
  178.         NEXT p2
  179.         'we now have net accelerations
  180.         particles(p).vx = particles(p).vx + ax
  181.         particles(p).vy = particles(p).vy + ay
  182.         particles(p).px = particles(p).px + particles(p).vx
  183.         particles(p).py = particles(p).py + particles(p).vy
  184.         IF particles(p).px < xmin THEN xmin = particles(p).px
  185.         IF particles(p).px > xmax THEN xmax = particles(p).px
  186.         IF particles(p).py < ymin THEN ymin = particles(p).py
  187.         IF particles(p).py > ymax THEN ymax = particles(p).py
  188.  
  189.         CIRCLE (particles(p).px, particles(p).py), 10, particles(p).c
  190.         $IF COLLISIONS = TRUE THEN
  191.             'check for collisions
  192.             p2 = 0
  193.             FOR p2 = 0 TO UBOUND(particles)
  194.                 IF p <> p2 THEN
  195.                     IF _HYPOT(particles(p2).px - particles(p).px, particles(p2).py - particles(p).py) < 20 THEN
  196.                         collision& = -1
  197.                         collisionx = particles(p).px + (particles(p2).px - particles(p).px) / 2
  198.                         collisiony = particles(p).py + (particles(p2).py - particles(p).py) / 2
  199.  
  200.                     END IF
  201.                 END IF
  202.             NEXT p2
  203.         $END IF
  204.     NEXT p
  205.  
  206.     WINDOW SCREEN(xmin - 20, ymin - 20)-(xmax + 20, ymax + 20)
  207.     t&& = t&& + 1
  208.     _PRINTSTRING (1, 1), "Time:" + STR$(t&&) + "   Escape: Abort      F8: View properties"
  209.  
  210.     $IF COLLISIONS = TRUE THEN
  211.         IF collision& THEN
  212.             CIRCLE (collisionx, collisiony), 40, _RGB32(255, 0, 0)
  213.             _DISPLAY
  214.             SLEEP
  215.             GOSUB reset_sim
  216.             GOTO start
  217.         END IF
  218.     $END IF
  219.  
  220.     _DISPLAY
  221.     _LIMIT 50
  222.  
  223.  
  224. reset_sim:
  225. t&& = 0
  226. collision& = 0
  227. xmin = 0
  228. ymin = 0
  229. xmax = 0
  230. ymax = 0
  231. active_obj = 0
  232.  
  233. print_table:
  234. PRINT "+--------------+---------------+---------------+---------------+---------------+"
  235. PRINT "| Body colour: |     White     |      Red      |      Blue     |     Green     |"
  236. PRINT "+--------------+---------------+---------------+---------------+---------------+"
  237. PRINT "|  Position X  |               |               |               |               |"
  238. PRINT "|  Position Y  |               |               |               |               |"
  239. PRINT "|  Velocity X  |               |               |               |               |"
  240. PRINT "|  Velocity Y  |               |               |               |               |"
  241. PRINT "|     Mass     |               |               |               |               |"
  242. PRINT "+--------------+---------------+---------------+---------------+---------------+"
  243. PRINT "Gravitational Constant: "
  244.  
  245. '-------------------------------------------------------------------------
  246. ' BASE64 Encoding / Decoding
  247. ' Original VBDOS Version by G. Balla, 1996 (Public Domain)
  248. ' QB 4.5 Conversion by Marc van den Dikkenberg, 1999
  249. ' From [http://www.qb45.com/download.php?id=1198], trimmed & optimized(no longer runs in QBASIC!) & updated (INTEGER->LONG for larger strings) & disabled input stream concatenation option --Galleon 2013
  250. '--------------------------------------------------------------------------
  251. DIM SHARED icChopMask AS LONG ' Constant 8-bit mask (Faster than using string constants)
  252. DIM SHARED icBitShift AS LONG ' Constant shift mask (Faster than using string constants)
  253. DIM SHARED icStartMask AS LONG ' Initial mask value  (Faster than using string constants)
  254. DIM SHARED iRollOver AS LONG ' Decoded Roll over value
  255. DIM SHARED iHighMask AS LONG ' Mask high bits of each char
  256. DIM SHARED iShift AS LONG ' Multiplier shift value
  257. DIM SHARED iLowShift AS LONG ' Mask low bits of each char
  258. DIM SHARED szAlphabet AS STRING ' Decode/Encode Lookup Table
  259. DIM SHARED szTemp AS STRING ' Working string
  260.  
  261. FUNCTION decodeBASE64$ (szEncoded AS STRING)
  262.     'iEndOfText AS LONG ''''removed from params
  263.     ' Initialize 2nd encoding pass lookup dictionary
  264.     szAlphabet = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/"
  265.     ' Initialize Constants
  266.     icChopMask = 255
  267.     icBitShift = 4
  268.     icStartMask = &H10
  269.     ' Initialize Masks
  270.     iShift = icBitShift
  271.     iLowShift = 0
  272.     iRollOver = 0
  273.     iHighMask = -1
  274.     iEndOfText = -1 ''''new
  275.     ' Create variables
  276.     DIM iPtr AS LONG
  277.     DIM iChar AS LONG
  278.     DIM iCounter AS LONG
  279.     ' Check if empty decoded string.
  280.     IF LEN(szEncoded) = 0 THEN
  281.         decodeBASE64$ = ""
  282.         EXIT FUNCTION
  283.     END IF
  284.     ' Initialize working string
  285.     szTemp = SPACE$(LEN(szEncoded) + 10) ''''changed
  286.     szTempLen = 0 ''''new
  287.     ' Begin Decoding
  288.     FOR iCounter = 1 TO LEN(szEncoded)
  289.         ' Get next alphabet
  290.         iChar = ASC(szEncoded, iCounter) ''''changed
  291.         ' Get Decoded value
  292.         iPtr = INSTR(szAlphabet, CHR$(iChar)) - 1
  293.         ' Check if character is valid
  294.         IF iPtr >= 0 THEN
  295.             ' Char is valid, process it
  296.             IF iShift = icBitShift THEN
  297.                 ' 1st char in block of 4, keep high part of character
  298.                 iRollOver = (iPtr * iShift) AND icChopMask
  299.                 ' Reset masks for next character
  300.                 iHighMask = &H30
  301.                 iLowShift = icStartMask
  302.                 iShift = icStartMask
  303.             ELSE
  304.                 ' Start saving decoded character
  305.                 szTempLen = szTempLen + 1: ASC(szTemp, szTempLen) = iRollOver OR ((iPtr AND iHighMask) / iLowShift) ''''changed
  306.                 ' Calculate next mask and shift values
  307.                 iRollOver = (iPtr * iShift) AND icChopMask
  308.                 iShift = iShift * icBitShift
  309.                 iHighMask = (iHighMask \ icBitShift) OR &H30
  310.                 iLowShift = iLowShift / icBitShift
  311.                 IF iShift > 256 THEN
  312.                     iShift = icBitShift
  313.                     iLowShift = 0
  314.                 END IF
  315.             END IF
  316.         END IF
  317.     NEXT
  318.     ' Concat last character if required
  319.     IF (iShift > icBitShift AND iShift < 256) THEN
  320.         ' Character remaining in    iRollOver
  321.         IF iEndOfText THEN
  322.             ' Last string to decode in file
  323.             szTempLen = szTempLen + 1: ASC(szTemp, szTempLen) = iRollOver ''''changed
  324.         END IF
  325.     END IF
  326.     ' Exit wth decoded string
  327.     decodeBASE64$ = LEFT$(szTemp, szTempLen) ''''changed
  328.  
  329. FUNCTION encodeBASE64$ (szUnEncoded AS STRING)
  330.     ' Create variables
  331.     DIM icLowFill AS LONG
  332.     DIM iChar AS LONG
  333.     DIM iLowMask AS LONG
  334.     DIM iPtr AS LONG
  335.     DIM iCounter AS LONG
  336.     ' Check if empty decoded string.
  337.     IF LEN(szUnEncoded) = 0 THEN
  338.         encodeBASE64$ = ""
  339.         EXIT FUNCTION
  340.     END IF
  341.     ' Initialize lookup dictionary and constants
  342.     szAlphabet = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/"
  343.     icBitShift = 4
  344.     icChopMask = 255
  345.     icLowFill = 3
  346.     ' Initialize Masks
  347.     szTemp = SPACE$(LEN(szUnEncoded) * 2 + 10) ''''changed
  348.     szTempLen = 0 ''''new
  349.     iHighMask = &HFC
  350.     iLowMask = &H3
  351.     iShift = &H10
  352.     iRollOver = 0
  353.     ' Begin Encoding process
  354.     FOR iCounter = 1 TO LEN(szUnEncoded)
  355.         ' Fetch ascii character in decoded string
  356.         iChar = ASC(szUnEncoded, iCounter) ''''changed
  357.         ' Calculate Alphabet lookup pointer
  358.         iPtr = ((iChar AND iHighMask) \ (iLowMask + 1)) OR iRollOver
  359.         ' Roll bit patterns
  360.         iRollOver = (iChar AND iLowMask) * iShift
  361.         ' Concatenate encoded character to working encoded string
  362.         szTempLen = szTempLen + 1: ASC(szTemp, szTempLen) = ASC(szAlphabet, iPtr + 1) ''''changed
  363.         ' Adjust masks
  364.         iHighMask = (iHighMask * icBitShift) AND icChopMask
  365.         iLowMask = iLowMask * icBitShift + icLowFill
  366.         iShift = iShift \ icBitShift
  367.         ' If last character in block, concat last RollOver and
  368.         '   reset masks
  369.         IF iHighMask = 0 THEN
  370.             szTempLen = szTempLen + 1: ASC(szTemp, szTempLen) = ASC(szAlphabet, iRollOver + 1) ''''changed
  371.             iRollOver = 0
  372.             iHighMask = &HFC
  373.             iLowMask = &H3
  374.             iShift = &H10
  375.         END IF
  376.     NEXT iCounter
  377.     ' If RollOver remains, concat it to the working string
  378.     IF iShift < &H10 THEN
  379.         szTempLen = szTempLen + 1: ASC(szTemp, szTempLen) = ASC(szAlphabet, iRollOver + 1) ''''changed
  380.     END IF
  381.     szTemp = LEFT$(szTemp, szTempLen) ''''new
  382.     ' Pad encoded string with required '=' pad characters
  383.     iPtr = (LEN(szTemp) MOD 4)
  384.     IF iPtr THEN szTemp = szTemp + STRING$(4 - iPtr, "=")
  385.     ' Return encoded string
  386.     encodeBASE64$ = szTemp
  387.  

Offline Qwerkey

  • Forum Resident
  • Posts: 755
    • View Profile
Re: The stable orbits game
« Reply #1 on: March 17, 2020, 10:50:59 am »
Interested members might like to compare with my generalised 3D gravitation simulation program:
https://www.qb64.org/forum/index.php?topic=581.0

Offline luke

  • Administrator
  • Seasoned Forum Regular
  • Posts: 324
    • View Profile
Re: The stable orbits game
« Reply #2 on: March 17, 2020, 06:09:17 pm »
Qwerkey's is better. Go look at Qwerkey's.

Offline STxAxTIC

  • Library Staff
  • Forum Resident
  • Posts: 1091
  • he lives
    • View Profile
Re: The stable orbits game
« Reply #3 on: March 17, 2020, 06:35:07 pm »
A quick and dirty way to see if anyone's central force simulator is correct is to see if it reproduces Kepler's laws:

1) elliptical orbits
2) equal areas in equal times
3) T^2 \propto a^3

anything else will just be a toy unfortunately.
You're not done when it works, you're done when it's right.

Offline Qwerkey

  • Forum Resident
  • Posts: 755
    • View Profile
Re: The stable orbits game
« Reply #4 on: March 18, 2020, 05:55:06 am »
Qwerkey's is better. Go look at Qwerkey's.
luke's is much easier to run.

Offline bplus

  • Global Moderator
  • Forum Resident
  • Posts: 8053
  • b = b + ...
    • View Profile
Re: The stable orbits game
« Reply #5 on: March 18, 2020, 10:34:36 am »
luke's is much easier to run.

hmm... yes I got Luke's going, a little learning curve with editor, press enter to exit editing a number for arrow keys to work again. But why the dx, dy inputs for velocity? For planets those change constantly, should only give the planet a velocity and let code figure the constantly changing changes of position on the x and y axis.

Qwerkey's, I didn't figure out, I suppose I was supposed to read the manual but the window didn't even fit my computer screen and being an InForm app not an easy hack to fix the window dimensions.

Offline Qwerkey

  • Forum Resident
  • Posts: 755
    • View Profile
Re: The stable orbits game
« Reply #6 on: March 19, 2020, 05:45:04 am »
@bplus As an effort to keep the world economy going, go and buy yourself a monitor with greater than 900 vertical resolution.  You keep ruling out my programs (where I generally like large graphics screens) because you can't see what I've done!  I had thought that I might suggest the Gravitation Simulation program for Samples (Maths & Geometry) as there's a lot of mathematics and it's certainly full of geometry and it is a good advert for InForm.

Offline bplus

  • Global Moderator
  • Forum Resident
  • Posts: 8053
  • b = b + ...
    • View Profile
Re: The stable orbits game
« Reply #7 on: March 19, 2020, 08:36:23 am »
@bplus As an effort to keep the world economy going, go and buy yourself a monitor with greater than 900 vertical resolution.  You keep ruling out my programs (where I generally like large graphics screens) because you can't see what I've done!  I had thought that I might suggest the Gravitation Simulation program for Samples (Maths & Geometry) as there's a lot of mathematics and it's certainly full of geometry and it is a good advert for InForm.

Sorry Qwerkey, I can't rule out any program I can't rule on, not logical. Hey, I know, @johnno56  loves these kind of sims, he has large screen too (i'm pretty sure), I will go by what he says or @SirCrow or @OldMoses  or @FellippeHeitor (can't get a better expert for InForm) ... hey has it been posted on this forum for proper cooking?

« Last Edit: March 19, 2020, 09:21:03 am by bplus »

Offline Qwerkey

  • Forum Resident
  • Posts: 755
    • View Profile
Re: The stable orbits game
« Reply #8 on: March 19, 2020, 09:07:05 am »
@bplus Sorry, I was trying to make a desperate desperate global situation into a (misplaced) humorous local problem here.  In reality, our work to put suitable projects into the curated sections looks somewhat minor.  In fact, since looking at this project, I've discovered some errors (not in the maths but in the workings of the buttons) so this is not ready for consideration.
« Last Edit: March 19, 2020, 09:13:10 am by Qwerkey »

Offline bplus

  • Global Moderator
  • Forum Resident
  • Posts: 8053
  • b = b + ...
    • View Profile
Re: The stable orbits game
« Reply #9 on: March 19, 2020, 09:14:48 am »
@Qwerkey, as usual I changed my reply to you a bit since you probably last looked at it. I am trying more helpful suggestions and calls for help.