Physics question

[BasicCoder2]

@dodicat,
My attempts just result in the coffin spinning out of control! Can you get it on target?
So I cheated!!
changed line 258 to make finer adjustments to the direction,

Code: Select all

        Dim As pt x2=Getline(mx,my,wheel/90,length)

And inserted these guidelines in the display code,

Code: Select all

        circle (200,400),5,rgb(255,0,0),,,,f
        circle (700,300),5,rgb(0,0,255),,,,f
        line (200,400)-(700,300),rgb(255,255,0)
       Screenunlock

I then rotated the red line until it aligned with the yellow line and with one click sent it onto the target.

[dodicat]

Thanks for testing, basiccoder2 & badidea.
basiccoder2.
You can make a slightly bigger grave, line 221--Dim As Long roomforerror=10
You could make 15 or something.
But your tweak has been noted, and it would make sense in the real universe for a nice send off.
Thank you.

[badidea]

...
Spacecraft can pass through the asteroid belt with virtually no chance of a collision.
...

How about a belt of neutron stars (or at least one):

Code: Select all

#include "fbgfx.bi"

const as single PI = 4 * atn(1)
dim shared as single PPM = 2 'pixels per meter
const SW = 800, SH = 600

const K_ESC = chr(27)
const K_MIN = chr(45)
const K_UND = chr(95)
const K_PLU = chr(61)
const K_EQU = chr(43)

screenres SW, SH, 32
width SW \ 8, SH \ 16

'-------------------------------------------------------------------------------

type sgl2d
   dim as single x, y
   declare constructor
   declare constructor(x as single, y as single)
   declare operator cast() as string
end type

constructor sgl2d
end constructor

constructor sgl2d(x as single, y as single)
   this.x = x : this.y = y
end constructor

operator sgl2d.cast () as string
   return str(x) & "," & str(y)
end operator

operator +(a as sgl2d, b as sgl2d) as sgl2d
   return sgl2d(a.x + b.x, a.y + b.y)
end operator

operator -(a as sgl2d, b as sgl2d) as sgl2d
   return sgl2d(a.x - b.x, a.y - b.y)
end operator

operator /(a as sgl2d, div as single) as sgl2d
   return sgl2d(a.x / div, a.y / div)
end operator

operator *(a as sgl2d, mul as single) as sgl2d
   return sgl2d(a.x * mul, a.y * mul)
end operator

function cross(a as sgl2d, b as sgl2d) as single
   return a.x * b.y - a.y * b.x
end function

function length(a as sgl2d) as single
   return sqr((a.x * a.x) + (a.y * a.y))
end function

function lengthSqrd(a as sgl2d) as single
   return (a.x * a.x) + (a.y * a.y)
end function

function dist(a as sgl2d, b as sgl2d) as single
   dim as single dx = a.x - b.x
   dim as single dy = a.y - b.y
   return sqr((dx * dx) + (dy * dy))
end function

function distSqrd(a as sgl2d, b as sgl2d) as single
   dim as single dx = a.x - b.x
   dim as single dy = a.y - b.y
   return (dx * dx) + (dy * dy)
end function

'-------------------------------------------------------------------------------

type polar
   dim as single angle
   dim as single magnitude
end type

function polarToCartesian(angle as single, radius as single) as sgl2d
   return sgl2d(cos(angle) * radius, sin(angle) * radius)
end function

function degToRad(degrees as single) as single
   return (degrees / 180) * PI
end function

function rotatedVector(v as sgl2d, rotAngle as single) as sgl2d
   dim as sgl2d tmp
   tmp.x = cos(rotAngle) * v.x - sin(rotAngle) * v.y
   tmp.y = sin(rotAngle) * v.x + cos(rotAngle) * v.y
   return tmp
end function

sub clearScreen(c as ulong)
   line(0, 0)-(SW - 1, SH - 1), c, bf
end sub

'scaled circle using PPM, y-axis pointing up, center = 0, 0
sub drawCircle(p as sgl2d, r as single, c as ulong)
   circle(SW \ 2 + p.x * PPM, SH \ 2 - p.y * PPM), r * PPM, c
end sub

'scaled line using PPM, y-axis pointing up, center = 0, 0
sub drawLine(p1 as sgl2d, p2 as sgl2d, c as ulong)
   line(SW \ 2 + p1.x * PPM, SH \ 2 - p1.y * PPM)-_
      (SW \ 2 + p2.x * PPM, SH \ 2 - p2.y * PPM), c
end sub

sub drawArrow(p1 as sgl2d, p2 as sgl2d, c as ulong)
   drawLine(p1, p2, c)
   dim as sgl2d posVector = p2 - p1
   posVector /= 3 '1/3 length
   drawLine(p1, p1 + rotatedVector(posVector, degToRad(+30)), c)
   drawLine(p1, p1 + rotatedVector(posVector, degToRad(-30)), c)
end sub

'-------------------------------------------------------------------------------

type disc_object
   dim as single radius '[m]
   dim as single height '[m]
   dim as single density '[kg/m^3]
   dim as ulong colour '[m]
   'linear motion properties
   dim as sgl2d position 'position [m]
   dim as single lin_m 'mass [kg]
   dim as sgl2d lin_F 'force [N] [kg*m/s^2]
   dim as sgl2d lin_a 'acceleration [m/s^2]
   dim as sgl2d lin_v 'velocity [m/s]
   'dim as sgl2d lin_p 'momentum [kg*m/s]
   'dim as single lin_E 'Kinetic energy [J] [kg*m^2/s^2]
   'Rotational motion properties
   dim as single angle 'angular position (theta) [rad]
   dim as single ang_m 'angular mass, moment of inertia (I) [kg*m^2]
   dim as single ang_F 'torque (tau) [N*m] [kg*m^2/s^2]
   dim as single ang_a 'angular velocity (alpha) [rad/s^2]
   dim as single ang_v 'angular velocity (omega) [rad/s]
   'dim as single ang_p 'angular momentum (L) [kg*m^2/s]
   'dim as single ang_E 'Kinetic energy [J] [kg*m^2/s^2]
   '
   declare sub init(r as single, h as single, d as single, p as sgl2d, c as ulong)
   declare sub update(dt as double)
   declare function getKineticEnergy() as single
end type

'Set radius, height, density, position
'Calculate mass and rotational inertia
sub disc_object.init(r as single, h as single, d as single, p as sgl2d, c as ulong)
   radius = r
   height = h
   density = d
   position = p
   colour = c
   lin_m = PI * r ^ 2 * d
   ang_m = 0.5 * lin_m * r ^ 2
end sub

'update position and angle
sub disc_object.update(dt as double)
   lin_a = lin_F / lin_m
   lin_v += lin_a * dt
   position += lin_v * dt
   ang_a = ang_F / ang_m
   ang_v += ang_a * dt
   angle += ang_v * dt
end sub

function disc_object.getKineticEnergy() as single
   dim as single lin_E = 0.5 * lin_m * lengthSqrd(lin_v)
   dim as single ang_E = 0.5 * ang_m * ang_v * ang_v
   return lin_E + ang_E
end function

'-------------------------------------------------------------------------------

type thruster_type
   '''init paramaters
   dim as polar polarForce '(rad, N)
   dim as polar polarPos '(rad, m)
   '''variable paramaters
   dim as sgl2d forceVector '(N, N)
   dim as sgl2d relPos, absPos '(m, m)
   dim as integer active
   declare sub init(forceMagnitude as single, forceDirection as single, posAngle as single, posRadius as single)
   declare sub updatePosition(bodyPos as sgl2d, bodyAngle as single)
end type

sub thruster_type.init(forceDirection as single, forceMagnitude as single, posAngle as single, posRadius as single)
   polarForce = type(forceDirection, forceMagnitude) 'thruster action
   polarPos = type(posAngle, posRadius) 'position of thruster on ship
end sub

sub thruster_type.updatePosition(bodyPos as sgl2d, bodyAngle as single)
   relPos = polarToCartesian(bodyAngle + polarPos.angle, polarPos.magnitude)
   absPos = bodyPos + relPos
end sub

'-------------------------------------------------------------------------------

const as single GRAV_CONST = 6.67e-11 '[m3/(kg*s^2)

type astro_body
   dim as single radius '[m]
   dim as single density '[kg/m^3]
   dim as ulong colour '[m]
   dim as sgl2d position 'position [m]
   dim as single mass '[kg]
   declare sub init(r as single, d as single, p as sgl2d, c as ulong)
end type

'Set radius, density, position
'Calculate mass and rotational inertia
sub astro_body.init(r as single, d as single, p as sgl2d, c as ulong)
   radius = r
   density = d
   position = p
   colour = c
   mass = PI * r ^ 2 * d
end sub

function gravForce(m1 as single, m2 as single, r as single) as single
   return GRAV_CONST * (m1 * m2) / (r * r)
end function

function gravForceVector(m1 as single, pos1 as sgl2d, m2 as single, pos2 as sgl2d) as sgl2d
   'https://en.wikipedia.org/wiki/Newton%27s_law_of_universal_gravitation#Vector_form
   dim as single distSquared = distSqrd(pos2, pos1)
   dim as sgl2d unitVector12 = (pos1 - pos2) / sqr(distSquared)
   return unitVector12 * (-GRAV_CONST * (m1 * m2) / distSquared)
end function

'-------------------------------------------------------------------------------

const as single MOON_RADIUS = 1737e+6 '[m]
const as single MOON_DENSITY = 3344 '[kg/m^3]

const NUM_THRUSTERS = 6
const L_FW_THR = 0 'left forward thruster
const R_FW_THR = 1 'right forward thruster
const L_LO_THR = 2
const R_LO_THR = 3
const L_HI_THR = 4
const R_HI_THR = 5

dim as string key
dim as integer quit = 0
dim as disc_object disc
dim as thruster_type thruster(NUM_THRUSTERS - 1)
dim as astro_body moon, neutronstar

'moon.init(MOON_RADIUS, MOON_DENSITY, sgl2d(+40, +20), rgb(255, 127, 0))
neutronstar.init(MOON_RADIUS * 1e-8, MOON_DENSITY * 1e8, sgl2d(+40, +20), rgb(255, 127, 0))
disc.init(10, 1, 5, sgl2d(0, -50), rgb(127, 255, 0))
'force angle, force magnitude, polar thruster position
thruster(L_FW_THR).init(0.5 * pi, 1e4, -0.75 * pi, disc.radius)
thruster(R_FW_THR).init(0.5 * pi, 1e4, -0.25 * pi, disc.radius)
thruster(L_LO_THR).init(0.0 * pi, 5e3, -0.75 * pi, disc.radius)
thruster(R_LO_THR).init(1.0 * pi, 5e3, -0.25 * pi, disc.radius)
thruster(L_HI_THR).init(0.0 * pi, 5e3, +0.75 * pi, disc.radius)
thruster(R_HI_THR).init(1.0 * pi, 5e3, +0.25 * pi, disc.radius)

dim as double tNow = timer, tPrev = tNow, dt = 0
while quit = 0
   'reset stuff
   disc.lin_F = sgl2d(0, 0)
   disc.ang_F = 0
   for i as integer = 0 to NUM_THRUSTERS - 1
      thruster(i).active = 0
   next

   'do always for display
   for i as integer = 0 to NUM_THRUSTERS - 1
      thruster(i).updatePosition(disc.position, disc.angle)
   next

   if multikey(FB.SC_UP) then
      thruster(L_FW_THR).active = 1
      thruster(R_FW_THR).active = 1
   end if

   if multikey(FB.SC_LEFT) then
      thruster(L_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
   end if

   if multikey(FB.SC_RIGHT) then
      thruster(R_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
   end if

   if key = K_MIN or key = K_UND then ppm /= 1.1 'zoom out
   if key = K_PLU or key = K_EQU then ppm *= 1.1 'zoom in
   if key = K_ESC then quit = 1

   for i as integer = 0 to NUM_THRUSTERS - 1
      'forces on body by active thrusters
      if thruster(i).active = 1 then
         thruster(i).forceVector = polarToCartesian(disc.angle + thruster(i).polarForce.angle, thruster(i).polarForce.magnitude)
         disc.lin_F += thruster(i).forceVector
         disc.ang_F += cross(thruster(i).relPos, thruster(i).forceVector)
      end if
   next

   disc.lin_F += gravForceVector(disc.lin_m, disc.position, neutronstar.mass, neutronstar.position)
   if dist(disc.position, neutronstar.position) < (disc.radius + neutronstar.radius) then quit = 2 'crashed
   
   disc.update(dt)
   
   'display
   screenlock
   clearScreen(0)
   locate 2, 2 : print "<UP>, <LEFT>, <RIGHT> for thrusters";
   locate 3, 2 : print "<+>, <-> for zoom in/out";
   locate 4, 2 : print "<ESC> to exit";
   locate 5, 2 : print "Kinetic engergy [J]: "; disc.getKineticEnergy();
   'locate 6, 2 : print str(cint(length(disc.position)))
   drawCircle(disc.position, disc.radius, disc.colour) 'flying saucer
   drawArrow(disc.position + disc.lin_v / 1, disc.position, rgb(255, 0, 127)) 'lin. speed ind.
   drawArrow(disc.position + disc.lin_F / 1e3, disc.position, rgb(0, 127, 255)) 'lin. speed ind.
   for i as integer = 0 to NUM_THRUSTERS - 1
      dim as ulong c = iif(i < 4, rgb(255, 255, 0), rgb(255, 255, 255))
      drawLine(disc.position, thruster(i).absPos, c) 'rotation indicator
      if thruster(i).active = 1 then
         drawArrow(thruster(i).absPos, thruster(i).absPos - thruster(i).forceVector / 1e3, rgb(255, 127, 0)) 'thruster force indicator
      end if
   next
   drawCircle(neutronstar.position, neutronstar.radius, neutronstar.colour) 'flying saucer
   screenunlock

   'time update
   key = inkey()
   sleep 1
   tPrev = tNow
   tNow = timer
   dt = tNow - tPrev
wend
locate 8, 2: print iif(quit = 2, "You crashed into the neutron star", "User quit request") + " (wait 3 sec)"
sleep 3000, 1
screen 0
print "End"


[h4tt3n]

Oh man, this is just my kind of thread :-) Awesome to come back to the forums after a long break to see you guys hard at work making cool stuff!

I have made loads of space physics doodle stuff just like this in the past few years, just gimme a bit of time and I'll dig out some stuff that is definitely useful...

[badidea]

This starts to look like a ship between asteroids:

Code: Select all

#Include "fbgfx.bi"

Type int2d
   As Integer x, y
   Declare Constructor
   Declare Constructor(x As Integer, y As Integer)
   Declare Operator Cast () As String
End Type

Constructor int2d
End Constructor

Constructor int2d(x As Integer, y As Integer)
   This.x = x : This.y = y
End Constructor

Operator = (a As int2d, b As int2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

Operator <> (a As int2d, b As int2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' "x, y"
Operator int2d.cast () As String
  Return Str(x) & "," & Str(y)
End Operator

' a + b
Operator + (a As int2d, b As int2d) As int2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As int2d, b As int2d) As int2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As int2d) As int2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As int2d, b As int2d) As int2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As int2d, mul As Integer) As int2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a \ b
Operator \ (a As int2d, b As int2d) As int2d
   Return Type(a.x \ b.x, a.y \ b.y)
End Operator

' a \ div
Operator \ (a As int2d, div As Integer) As int2d
   Return Type(a.x \ div, a.y \ div)
End Operator

'===============================================================================

Type sgl2d
   As Single x, y
   Declare Constructor
   Declare Constructor(x As Single, y As Single)
   Declare Operator Cast () As String
End Type

Constructor sgl2d
End Constructor

Constructor sgl2d(x As Single, y As Single)
   This.x = x : This.y = y
End Constructor

' "x, y"
Operator sgl2d.cast () As String
   Return Str(x) & "," & Str(y)
End Operator

'---- operators ---

' distance / lenth
Operator Len (a As sgl2d) As Single
   Return Sqr(a.x * a.x + a.y * a.y)
End Operator

' a = b ?
Operator = (a As sgl2d, b As sgl2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

' a != b ?
Operator <> (a As sgl2d, b As sgl2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' a + b
Operator + (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As sgl2d) As sgl2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As sgl2d, mul As Single) As sgl2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a / div
Operator / (a As sgl2d, div As Single) As sgl2d
   Return Type(a.x / div, a.y / div)
End Operator

'---- extra functions ---

Function cross(a As sgl2d, b As sgl2d) As Single
   Return a.x * b.y - a.y * b.x
End Function

Function length(a As sgl2d) As Single
   Return Sqr((a.x * a.x) + (a.y * a.y))
End Function

Function lengthSqrd(a As sgl2d) As Single
   Return (a.x * a.x) + (a.y * a.y)
End Function

Function dist(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return Sqr((dx * dx) + (dy * dy))
End Function

Function distSqrd(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return (dx * dx) + (dy * dy)
End Function

'===============================================================================

'Note: y+ = up, x+ = right, (0,0) = center
Type scaled_graphics_type
   Dim As Single scale = 1 ' = 1 / pixel_size 'pixels / meter
   'dim as int2d offset' = (scrn_w \ 2, h \ 2) 'offset in pixels
   Dim As sgl2d offset
   Dim As Integer w = -1, h = -1
   Dim As Integer wc = -1, hc = -1 'center x,y
   Declare Sub setScreen(w As Integer, h As Integer)
   Declare Sub setScaling(scale As Single, offset As sgl2d)
   Declare Sub clearScreen(c As ULong)
   Declare Function pos2screen(p As sgl2d) As int2d
   Declare Sub drawPixel(p As sgl2d, c As Integer)
   Declare Sub drawCircle(p As sgl2d, r As Single, c As Integer)
   Declare Sub drawLine(p1 As sgl2d, p2 As sgl2d, c As Integer)
End Type

Sub scaled_graphics_type.setScreen(w As Integer, h As Integer)
   This.w = w 'width
   This.h = h 'height
   wc = w \ 2
   hc = h \ 2
   ScreenRes w, h, 32
   Width w \ 8, h \ 16 'bigger font
End Sub

Sub scaled_graphics_type.setScaling(scale As Single, offset As sgl2d)
   This.scale = scale
   This.offset = offset
End Sub

Sub scaled_graphics_type.clearScreen(c As ULong)
   Line(0, 0)-(w - 1, h - 1), c, bf
End Sub

Function scaled_graphics_type.pos2screen(p As sgl2d) As int2d
   Return int2d(Int(wc + (p.x - offset.x) * scale), h - Int(hc + (p.y - offset.y) * scale))
End Function

Sub scaled_graphics_type.drawPixel(p As sgl2d, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   PSet(posScrn.x, posScrn.y), c
End Sub

Sub scaled_graphics_type.drawCircle(p As sgl2d, r As Single, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawLine(p1 As sgl2d, p2 As sgl2d, c As Integer)
   Dim As int2d posScrn1 = pos2screen(p1)
   Dim As int2d posScrn2 = pos2screen(p2)
   Line(posScrn1.x, posScrn1.y)-(posScrn2.x, posScrn2.y), c
End Sub

'===============================================================================

Const As Single PI = 4 * Atn(1)
Const As Single sinA = Sin((10 / 180) * PI)
Const As Single cosA = Cos((10 / 180) * PI)

Const K_ESC = Chr(27)
Const K_MIN = Chr(45)
Const K_UND = Chr(95)
Const K_PLU = Chr(61)
Const K_EQU = Chr(43)

Const SCRN_W = 800, SCRN_H = 600

Dim Shared As scaled_graphics_type sg
sg.setScaling(2.0, sgl2d(0, 0))
sg.setScreen(SCRN_W, SCRN_H)

'-------------------------------------------------------------------------------

Type polar
   Dim As Single angle
   Dim As Single magnitude
End Type

Function polarToCartesian(angle As Single, radius As Single) As sgl2d
   Return sgl2d(Cos(angle) * radius, Sin(angle) * radius)
End Function

Function degToRad(degrees As Single) As Single
   Return (degrees / 180) * PI
End Function

Function rotatedVector(v As sgl2d, rotAngle As Single) As sgl2d
   Dim As sgl2d tmp
   tmp.x = Cos(rotAngle) * v.x - Sin(rotAngle) * v.y
   tmp.y = Sin(rotAngle) * v.x + Cos(rotAngle) * v.y
   Return tmp
End Function

'-------------------------------------------------------------------------------

Sub drawArrow(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.95 'reduce length
   'sg.drawLine(p1, p1 + rotatedVector(dp, degToRad(+30)), c)
   'sg.drawLine(p1, p1 + rotatedVector(dp, degToRad(-30)), c)
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x - sinA * dp.y, sinA * dp.x + cosA * dp.y), c)
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x + sinA * dp.y, cosA * dp.y - sinA * dp.x), c)
End Sub

'-------------------------------------------------------------------------------

Type disc_object
   Dim As Single radius '[m]
   Dim As Single height '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   'linear motion properties
   Dim As sgl2d position 'position [m]
   Dim As Single lin_m 'mass [kg]
   Dim As sgl2d lin_F 'force [N] [kg*m/s^2]
   Dim As sgl2d lin_a 'acceleration [m/s^2]
   Dim As sgl2d lin_v 'velocity [m/s]
   'dim as sgl2d lin_p 'momentum [kg*m/s]
   'dim as single lin_E 'Kinetic energy [J] [kg*m^2/s^2]
   'Rotational motion properties
   Dim As Single angle 'angular position (theta) [rad]
   Dim As Single ang_m 'angular mass, moment of inertia (I) [kg*m^2]
   Dim As Single ang_F 'torque (tau) [N*m] [kg*m^2/s^2]
   Dim As Single ang_a 'angular velocity (alpha) [rad/s^2]
   Dim As Single ang_v 'angular velocity (omega) [rad/s]
   'dim as single ang_p 'angular momentum (L) [kg*m^2/s]
   'dim as single ang_E 'Kinetic energy [J] [kg*m^2/s^2]
   '
   Declare Sub init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   Declare Sub update(dt As Double)
   Declare Function getKineticEnergy() As Single
End Type

'Set radius, height, density, position
'Calculate mass and rotational inertia
Sub disc_object.init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   height = h
   density = d
   position = p
   colour = c
   lin_m = PI * r ^ 2 * d
   ang_m = 0.5 * lin_m * r ^ 2
End Sub

'update position and angle
Sub disc_object.update(dt As Double)
   lin_a = lin_F / lin_m
   lin_v += lin_a * dt
   position += lin_v * dt
   ang_a = ang_F / ang_m
   ang_v += ang_a * dt
   angle += ang_v * dt
End Sub

Function disc_object.getKineticEnergy() As Single
   Dim As Single lin_E = 0.5 * lin_m * lengthSqrd(lin_v)
   Dim As Single ang_E = 0.5 * ang_m * ang_v * ang_v
   Return lin_E + ang_E
End Function

'-------------------------------------------------------------------------------

Type thruster_type
   '''init paramaters
   Dim As polar polarForce '(rad, N)
   Dim As polar polarPos '(rad, m)
   '''variable paramaters
   Dim As sgl2d forceVector '(N, N)
   Dim As sgl2d relPos, absPos '(m, m)
   Dim As Integer active
   Declare Sub init(forceMagnitude As Single, forceDirection As Single, posAngle As Single, posRadius As Single)
   Declare Sub updatePosition(bodyPos As sgl2d, bodyAngle As Single)
End Type

Sub thruster_type.init(forceDirection As Single, forceMagnitude As Single, posAngle As Single, posRadius As Single)
   polarForce = Type(forceDirection, forceMagnitude) 'thruster action
   polarPos = Type(posAngle, posRadius) 'position of thruster on ship
End Sub

Sub thruster_type.updatePosition(bodyPos As sgl2d, bodyAngle As Single)
   relPos = polarToCartesian(bodyAngle + polarPos.angle, polarPos.magnitude)
   absPos = bodyPos + relPos
End Sub

'-------------------------------------------------------------------------------

Const As Single GRAV_CONST = 6.67e-11 '[m3/(kg*s^2)

Type astro_body
   Dim As Single radius '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   Dim As sgl2d position 'position [m]
   Dim As Single mass '[kg]
   Declare Sub init(r As Single, d As Single, p As sgl2d, c As ULong)
End Type

'Set radius, density, position
'Calculate mass and rotational inertia
Sub astro_body.init(r As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   density = d
   position = p
   colour = c
   mass = PI * r ^ 2 * d
End Sub

Function gravForce(m1 As Single, m2 As Single, r As Single) As Single
   Return GRAV_CONST * (m1 * m2) / (r * r)
End Function

Function gravForceVector(m1 As Single, Pos1 As sgl2d, m2 As Single, pos2 As sgl2d) As sgl2d
   'https://en.wikipedia.org/wiki/Newton%27s_law_of_universal_gravitation#Vector_form
   Dim As Single distSquared = distSqrd(pos2, Pos1)
   Dim As sgl2d unitVector12 = (Pos1 - pos2) / Sqr(distSquared)
   Return unitVector12 * (-GRAV_CONST * (m1 * m2) / distSquared)
End Function

'-------------------------------------------------------------------------------

Const As Single MOON_RADIUS = 1737e+6 '[m]
Const As Single MOON_DENSITY = 3344 '[kg/m^3]

Const NUM_THRUSTERS = 6
Const L_FW_THR = 0 'left forward thruster
Const R_FW_THR = 1 'right forward thruster
Const L_LO_THR = 2
Const R_LO_THR = 3
Const L_HI_THR = 4
Const R_HI_THR = 5

Dim As String key
Dim As Integer quit = 0
Dim As disc_object disc
Dim As thruster_type thruster(NUM_THRUSTERS - 1)
'dim as astro_body moon, neutronstar
Const NUM_ASTROID = 100
Dim As astro_body astroid(NUM_ASTROID-1)

'moon.init(MOON_RADIUS, MOON_DENSITY, sgl2d(+40, +20), rgb(255, 127, 0))
'neutronstar.init(MOON_RADIUS * 1e-8, MOON_DENSITY * 1e8, sgl2d(+40, +20), rgb(255, 127, 0))
disc.init(10, 1, 5, sgl2d(0, -50), RGB(127, 255, 0))
For i As Integer = 0 To UBound(astroid)
   astroid(i).init(10, 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(255, 127, 0))
Next

'force angle, force magnitude, polar thruster position
thruster(L_FW_THR).init(0.5 * pi, 1.2e4, -0.75 * pi, disc.radius)
thruster(R_FW_THR).init(0.5 * pi, 1.2e4, -0.25 * pi, disc.radius)
thruster(L_LO_THR).init(0.0 * pi, 8e3, -0.75 * pi, disc.radius)
thruster(R_LO_THR).init(1.0 * pi, 8e3, -0.25 * pi, disc.radius)
thruster(L_HI_THR).init(0.0 * pi, 8e3, +0.75 * pi, disc.radius)
thruster(R_HI_THR).init(1.0 * pi, 8e3, +0.25 * pi, disc.radius)

Dim As Double tNow = Timer, tPrev = tNow, dt = 0
While quit = 0
   'reset stuff
   disc.lin_F = sgl2d(0, 0)
   disc.ang_F = 0
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).active = 0
   Next

   'do always for display
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).updatePosition(disc.position, disc.angle)
   Next

   If MultiKey(FB.SC_UP) Then
      thruster(L_FW_THR).active = 1
      thruster(R_FW_THR).active = 1
   End If

   If MultiKey(FB.SC_LEFT) Then
      thruster(L_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
   End If

   If MultiKey(FB.SC_RIGHT) Then
      thruster(R_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
   End If

   'If key = K_MIN Or key = K_UND Then sg.scale /= 1.1 'zoom out
   'If key = K_PLU Or key = K_EQU Then sg.scale *= 1.1 'zoom in
   If key = K_ESC Then quit = 1

   For i As Integer = 0 To NUM_THRUSTERS - 1
      'forces on body by active thrusters
      If thruster(i).active = 1 Then
         thruster(i).forceVector = polarToCartesian(disc.angle + thruster(i).polarForce.angle, thruster(i).polarForce.magnitude)
         disc.lin_F += thruster(i).forceVector
         disc.ang_F += cross(thruster(i).relPos, thruster(i).forceVector)
      End If
   Next

   'disc.lin_F += gravForceVector(disc.lin_m, disc.position, neutronstar.mass, neutronstar.position)
   'if dist(disc.position, neutronstar.position) < (disc.radius + neutronstar.radius) then quit = 2 'crashed
   
   disc.update(dt)
   'X_C = -disc.position.x
   'Y_C = disc.position.y
   'sg.offset = int2d(cint(disc.position.x), cint(disc.position.y))
   sg.offset = disc.position
   Dim As Single ppm = 100 / Len(disc.lin_v)
   If ppm < 0.5 Then ppm = 0.5
   If ppm > 2.0 Then ppm = 2.0
   sg.scale = ppm

   'calculate ship direction pointer / triangle
   Dim As Single forwardAngle = disc.angle - PI/2
   Dim As Single leftBackAngle = forwardAngle + degToRad(135)
   Dim As Single rightBackAngle = forwardAngle - degToRad(135)
   Dim As sgl2d forwardPos = sgl2d(Cos(forwardAngle), Sin(forwardAngle)) * (disc.radius * 2.2)
   Dim As sgl2d leftBackPos = sgl2d(Cos(leftBackAngle), Sin(leftBackAngle)) * (disc.radius * 1.1)
   Dim As sgl2d rightBackPos = sgl2d(Cos(rightBackAngle), Sin(rightBackAngle)) * (disc.radius * 1.1)
   
   'display
   ScreenLock
   sg.clearScreen(0)
   Locate 2, 2 : Print "<UP>, <LEFT>, <RIGHT> for thrusters";
   'Locate 3, 2 : Print "<+>, <-> for zoom in/out";
   Locate 4, 2 : Print "<ESC> to exit";
   Locate 5, 2 : Print "Kinetic engergy [J]: "; disc.getKineticEnergy();
   Locate 6, 2 : Print "ppm: "; ppm
   'draw ship
   
   sg.drawCircle(disc.position, disc.radius, disc.colour) 'flying saucer
   sg.drawLine(disc.position - forwardPos, disc.position - leftBackPos, disc.colour)
   sg.drawLine(disc.position - forwardPos, disc.position - rightBackPos, disc.colour)
   sg.drawLine(disc.position - leftBackPos, disc.position - rightBackPos, disc.colour)
   'drawArrow(disc.position + disc.lin_v / 1, disc.position, rgb(255, 0, 127)) 'lin. speed ind.
   'drawArrow(disc.position + disc.lin_F / 1e3, disc.position, rgb(0, 127, 255)) 'lin. speed ind.
   For i As Integer = 0 To NUM_THRUSTERS - 1
      Dim As ULong c = IIf(i < 4, RGB(255, 255, 0), RGB(255, 255, 255))
      'sg.drawLine(disc.position, thruster(i).absPos, c) 'rotation indicator
      If thruster(i).active = 1 Then
         drawArrow(thruster(i).absPos, thruster(i).absPos - thruster(i).forceVector / 1e3, RGB(255, 127, 0)) 'thruster force indicator
      End If
   Next
   'drawCircle(neutronstar.position, neutronstar.radius, neutronstar.colour) 'flying saucer
   For i As Integer = 0 To UBound(astroid)
      sg.drawCircle(astroid(i).position, astroid(i).radius, astroid(i).colour)
   Next
   ScreenUnLock

   'time update
   key = InKey()
   Sleep 1
   tPrev = tNow
   tNow = Timer
   dt = tNow - tPrev
Wend
'locate 8, 2: print iif(quit = 2, "You crashed into the neutron star", "User quit request") + " (wait 3 sec)"
'sleep 300, 1
Screen 0
Print "End"

'TODO:
'center of universe indicator?
'nearest object indicator + distance
'better looking ship
'drawElipse


Ship centered view, and view scaling with speed.
No collisions, space battles or death stars yet.

[badidea]

'Improved' ship design with a nearest asteroid indicator:

Code: Select all

#Include "fbgfx.bi"

Type int2d
   As Integer x, y
   Declare Constructor
   Declare Constructor(x As Integer, y As Integer)
   Declare Operator Cast () As String
End Type

Constructor int2d
End Constructor

Constructor int2d(x As Integer, y As Integer)
   This.x = x : This.y = y
End Constructor

Operator = (a As int2d, b As int2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

Operator <> (a As int2d, b As int2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' "x, y"
Operator int2d.cast () As String
  Return Str(x) & "," & Str(y)
End Operator

' a + b
Operator + (a As int2d, b As int2d) As int2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As int2d, b As int2d) As int2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As int2d) As int2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As int2d, b As int2d) As int2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As int2d, mul As Integer) As int2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a \ b
Operator \ (a As int2d, b As int2d) As int2d
   Return Type(a.x \ b.x, a.y \ b.y)
End Operator

' a \ div
Operator \ (a As int2d, div As Integer) As int2d
   Return Type(a.x \ div, a.y \ div)
End Operator

'===============================================================================

Type sgl2d
   As Single x, y
   Declare Constructor
   Declare Constructor(x As Single, y As Single)
   Declare Operator Cast () As String
End Type

Constructor sgl2d
End Constructor

Constructor sgl2d(x As Single, y As Single)
   This.x = x : This.y = y
End Constructor

' "x, y"
Operator sgl2d.cast () As String
   Return Str(x) & "," & Str(y)
End Operator

'---- operators ---

' distance / lenth
Operator Len (a As sgl2d) As Single
   Return Sqr(a.x * a.x + a.y * a.y)
End Operator

' a = b ?
Operator = (a As sgl2d, b As sgl2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

' a != b ?
Operator <> (a As sgl2d, b As sgl2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' a + b
Operator + (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As sgl2d) As sgl2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As sgl2d, mul As Single) As sgl2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a / div
Operator / (a As sgl2d, div As Single) As sgl2d
   Return Type(a.x / div, a.y / div)
End Operator

'---- extra functions ---

Function cross(a As sgl2d, b As sgl2d) As Single
   Return a.x * b.y - a.y * b.x
End Function

Function lengthSqrd(a As sgl2d) As Single
   Return (a.x * a.x) + (a.y * a.y)
End Function

Function dist(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return Sqr((dx * dx) + (dy * dy))
End Function

Function distSqrd(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return (dx * dx) + (dy * dy)
End Function

Function normalise(a As sgl2d) As sgl2d
   Dim As sgl2d temp
   Dim As Single length = Len(a)
   Return sgl2d(a.x / length, a.y / length)
End Function

'===============================================================================

'Note: y+ = up, x+ = right, (0,0) = center
Type scaled_graphics_type
   Dim As Single scale = 1 ' = 1 / pixel_size 'pixels / meter
   'dim as int2d offset' = (scrn_w \ 2, h \ 2) 'offset in pixels
   Dim As sgl2d offset
   Dim As Integer w = -1, h = -1
   Dim As Integer wc = -1, hc = -1 'center x,y
   Declare Sub setScreen(w As Integer, h As Integer)
   Declare Sub setScaling(scale As Single, offset As sgl2d)
   Declare Sub clearScreen(c As ULong)
   Declare Function pos2screen(p As sgl2d) As int2d
   Declare Sub drawPixel(p As sgl2d, c As Integer)
   Declare Sub drawCircle(p As sgl2d, r As Single, c As Integer)
   Declare Sub drawLine(p1 As sgl2d, p2 As sgl2d, c As Integer)
End Type

Sub scaled_graphics_type.setScreen(w As Integer, h As Integer)
   This.w = w 'width
   This.h = h 'height
   wc = w \ 2
   hc = h \ 2
   ScreenRes w, h, 32
   Width w \ 8, h \ 16 'bigger font
End Sub

Sub scaled_graphics_type.setScaling(scale As Single, offset As sgl2d)
   This.scale = scale
   This.offset = offset
End Sub

Sub scaled_graphics_type.clearScreen(c As ULong)
   Line(0, 0)-(w - 1, h - 1), c, bf
End Sub

Function scaled_graphics_type.pos2screen(p As sgl2d) As int2d
   Return int2d(Int(wc + (p.x - offset.x) * scale), h - Int(hc + (p.y - offset.y) * scale))
End Function

Sub scaled_graphics_type.drawPixel(p As sgl2d, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   PSet(posScrn.x, posScrn.y), c
End Sub

Sub scaled_graphics_type.drawCircle(p As sgl2d, r As Single, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawLine(p1 As sgl2d, p2 As sgl2d, c As Integer)
   Dim As int2d posScrn1 = pos2screen(p1)
   Dim As int2d posScrn2 = pos2screen(p2)
   Line(posScrn1.x, posScrn1.y)-(posScrn2.x, posScrn2.y), c
End Sub

'===============================================================================

Const As Single PI = 4 * Atn(1)
Const As Single RAD_PER_DEG = (PI / 180)
Const As Single DEG_PER_RAD = 180 / PI
Const As Single sinA = Sin((10 / 180) * PI)
Const As Single cosA = Cos((10 / 180) * PI)
Const As Single sinB = Sin((20 / 180) * PI)
Const As Single cosB = Cos((20 / 180) * PI)

Const K_ESC = Chr(27)
Const K_MIN = Chr(45)
Const K_UND = Chr(95)
Const K_PLU = Chr(61)
Const K_EQU = Chr(43)

Const SCRN_W = 800, SCRN_H = 600

Dim Shared As scaled_graphics_type sg
sg.setScaling(2.0, sgl2d(0, 0))
sg.setScreen(SCRN_W, SCRN_H)

'-------------------------------------------------------------------------------

Function limit(value As Single, min As Single, max As Single) As Single
   If value < min Then Return min
   If value > max Then Return max
   Return value
End Function

Type polar
   Dim As Single angle
   Dim As Single magnitude
End Type

Function polarToCartesian(angle As Single, radius As Single) As sgl2d
   Return sgl2d(Cos(angle) * radius, Sin(angle) * radius)
End Function

Function rotatedVector(v As sgl2d, rotAngle As Single) As sgl2d
   Dim As sgl2d tmp
   tmp.x = Cos(rotAngle) * v.x - Sin(rotAngle) * v.y
   tmp.y = Sin(rotAngle) * v.x + Cos(rotAngle) * v.y
   Return tmp
End Function

'-------------------------------------------------------------------------------

Sub drawArrow(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.30 'reduce length
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x - sinB * dp.y, sinB * dp.x + cosB * dp.y), c)
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x + sinB * dp.y, cosB * dp.y - sinB * dp.x), c)
End Sub

Sub drawThruster(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.95 'reduce length
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x - sinA * dp.y, sinA * dp.x + cosA * dp.y), c)
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x + sinA * dp.y, cosA * dp.y - sinA * dp.x), c)
End Sub

'-------------------------------------------------------------------------------

Type disc_object
   Dim As Single radius '[m]
   Dim As Single height '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   'linear motion properties
   Dim As sgl2d position 'position [m]
   Dim As Single lin_m 'mass [kg]
   Dim As sgl2d lin_F 'force [N] [kg*m/s^2]
   Dim As sgl2d lin_a 'acceleration [m/s^2]
   Dim As sgl2d lin_v 'velocity [m/s]
   'dim as sgl2d lin_p 'momentum [kg*m/s]
   'dim as single lin_E 'Kinetic energy [J] [kg*m^2/s^2]
   'Rotational motion properties
   Dim As Single angle 'angular position (theta) [rad]
   Dim As Single ang_m 'angular mass, moment of inertia (I) [kg*m^2]
   Dim As Single ang_F 'torque (tau) [N*m] [kg*m^2/s^2]
   Dim As Single ang_a 'angular velocity (alpha) [rad/s^2]
   Dim As Single ang_v 'angular velocity (omega) [rad/s]
   'dim as single ang_p 'angular momentum (L) [kg*m^2/s]
   'dim as single ang_E 'Kinetic energy [J] [kg*m^2/s^2]
   Declare Sub init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   Declare Sub update(dt As Double)
   Declare Function getKineticEnergy() As Single
End Type

'Set radius, height, density, position
'Calculate mass and rotational inertia
Sub disc_object.init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   height = h
   density = d
   position = p
   colour = c
   lin_m = PI * r ^ 2 * d
   ang_m = 0.5 * lin_m * r ^ 2
End Sub

'update position and angle
Sub disc_object.update(dt As Double)
   lin_a = lin_F / lin_m
   lin_v += lin_a * dt
   position += lin_v * dt
   ang_a = ang_F / ang_m
   ang_v += ang_a * dt
   angle += ang_v * dt
End Sub

Function disc_object.getKineticEnergy() As Single
   Dim As Single lin_E = 0.5 * lin_m * lengthSqrd(lin_v)
   Dim As Single ang_E = 0.5 * ang_m * ang_v * ang_v
   Return lin_E + ang_E
End Function

'-------------------------------------------------------------------------------

Type thruster_type
   '''init paramaters
   Dim As polar polarForce '(rad, N)
   Dim As polar polarPos '(rad, m)
   '''variable paramaters
   Dim As sgl2d forceVector '(N, N)
   Dim As sgl2d relPos, absPos '(m, m)
   Dim As Integer active
   Declare Sub init(forceMagnitude As Single, forceDirection As Single, posAngle As Single, posRadius As Single)
   Declare Sub updatePosition(bodyPos As sgl2d, bodyAngle As Single)
End Type

Sub thruster_type.init(forceDirection As Single, forceMagnitude As Single, posAngle As Single, posRadius As Single)
   polarForce = Type(forceDirection, forceMagnitude) 'thruster action
   polarPos = Type(posAngle, posRadius) 'position of thruster on ship
End Sub

Sub thruster_type.updatePosition(bodyPos As sgl2d, bodyAngle As Single)
   relPos = polarToCartesian(bodyAngle + polarPos.angle, polarPos.magnitude)
   absPos = bodyPos + relPos
End Sub

'-------------------------------------------------------------------------------

Const As Single GRAV_CONST = 6.67e-11 '[m3/(kg*s^2)

Type astro_body
   Dim As Single radius '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   Dim As sgl2d position 'position [m]
   Dim As Single mass '[kg]
   Declare Sub init(r As Single, d As Single, p As sgl2d, c As ULong)
End Type

'Set radius, density, position
'Calculate mass and rotational inertia
Sub astro_body.init(r As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   density = d
   position = p
   colour = c
   mass = PI * r ^ 2 * d
End Sub

Function gravForce(m1 As Single, m2 As Single, r As Single) As Single
   Return GRAV_CONST * (m1 * m2) / (r * r)
End Function

Function gravForceVector(m1 As Single, Pos1 As sgl2d, m2 As Single, pos2 As sgl2d) As sgl2d
   Dim As Single distSquared = distSqrd(pos2, Pos1)
   Dim As sgl2d unitVector12 = (Pos1 - pos2) / Sqr(distSquared)
   Return unitVector12 * (-GRAV_CONST * (m1 * m2) / distSquared)
End Function

'-------------------------------------------------------------------------------

Const As Single MOON_RADIUS = 1737e+6 '[m]
Const As Single MOON_DENSITY = 3344 '[kg/m^3]

Const NUM_THRUSTERS = 6
Const L_FW_THR = 0 'left forward thruster
Const R_FW_THR = 1 'right forward thruster
Const L_LO_THR = 2
Const R_LO_THR = 3
Const L_HI_THR = 4
Const R_HI_THR = 5

Dim As String key
Dim As Integer quit = 0
Dim As disc_object ship
Dim As thruster_type thruster(NUM_THRUSTERS - 1)
Const NUM_asteroid = 100
Dim As astro_body asteroid(NUM_ASTEROID-1)

ship.init(10, 1, 5, sgl2d(0, -50), RGB(127, 223, 0))
For i As Integer = 0 To UBound(asteroid)
   asteroid(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(255, 191, 0))
Next

'force angle, force magnitude, polar thruster position
thruster(L_FW_THR).init(0.5 * pi, 1.2e4, -0.75 * pi, ship.radius)
thruster(R_FW_THR).init(0.5 * pi, 1.2e4, -0.25 * pi, ship.radius)
thruster(L_LO_THR).init(0.0 * pi, 8e3, -0.75 * pi, ship.radius)
thruster(R_LO_THR).init(1.0 * pi, 8e3, -0.25 * pi, ship.radius)
thruster(L_HI_THR).init(0.0 * pi, 8e3, +0.75 * pi, ship.radius)
thruster(R_HI_THR).init(1.0 * pi, 8e3, +0.25 * pi, ship.radius)

Dim As Double tNow = Timer, tPrev = tNow, dt = 0
While quit = 0
   'reset stuff
   ship.lin_F = sgl2d(0, 0)
   ship.ang_F = 0
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).active = 0
   Next

   If MultiKey(FB.SC_UP) Then
      thruster(L_FW_THR).active = 1
      thruster(R_FW_THR).active = 1
   End If

   '~ if multikey(FB.SC_SPACE) then
      '~ thruster(L_FW_THR).active = 2 'hyper-drive
      '~ thruster(R_FW_THR).active = 2 'hyper-drive
   '~ end if

   If MultiKey(FB.SC_LEFT) Then
      thruster(L_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
   End If

   If MultiKey(FB.SC_RIGHT) Then
      thruster(R_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
   End If

   If key = K_ESC Then quit = 1

   For i As Integer = 0 To NUM_THRUSTERS - 1
      'forces on body by active thrusters
      If thruster(i).active > 0 Then
         Dim As Single thrust = thruster(i).polarForce.magnitude
         'if thruster(i).active = 2 then thrust *= 10 'hyper-drive
         thruster(i).forceVector = polarToCartesian(ship.angle + thruster(i).polarForce.angle, thrust)
         ship.lin_F += thruster(i).forceVector
         ship.ang_F += cross(thruster(i).relPos, thruster(i).forceVector)
      End If
   Next

   ship.update(dt) 'position and angle
   sg.offset = ship.position
   sg.scale = limit(100 / Len(ship.lin_v), 0.5, 2.0)

   'calculate ship direction pointer / triangle
   Dim As sgl2d forwardPos = polarToCartesian(ship.angle - 90 * RAD_PER_DEG, ship.radius * 2.2)
   Dim As sgl2d leftBackPos = polarToCartesian(ship.angle - 135 * RAD_PER_DEG, ship.radius * 1.0) '(90 + 45)
   Dim As sgl2d rightBackPos = polarToCartesian(ship.angle - 45 * RAD_PER_DEG, ship.radius * 1.0) '(90 - 45)

   'do always for display
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).updatePosition(ship.position, ship.angle)
   Next

   'find nearest asteroid
   Dim As Integer nearestAsteroidId = 0
   Dim As Single nearestAsteroidDistSqrd = distSqrd(ship.position, asteroid(0).position)
   Dim As Single currenAsteroidDistSqrd
   For i As Integer = 1 To UBound(asteroid)
      currenAsteroidDistSqrd = distSqrd(ship.position, asteroid(i).position)
      If currenAsteroidDistSqrd < nearestAsteroidDistSqrd Then
         nearestAsteroidDistSqrd = currenAsteroidDistSqrd
         nearestAsteroidId = i
      End If
   Next
   
   'display
   ScreenLock
   sg.clearScreen(0)
   Locate 2, 2 : Print "<UP>, <LEFT>, <RIGHT> for thrusters, <ESC> to exit";
   Locate 3, 2 : Print "Nearest asteroid distance [m]: "; CInt(Sqr(nearestAsteroidDistSqrd));
   Locate 4, 2 : Print "Kinetic engergy [kJ]: "; CInt(ship.getKineticEnergy() * 1e-3);
   'draw ship
   sg.drawCircle(ship.position, ship.radius, ship.colour) 'flying saucer
   sg.drawLine(ship.position + forwardPos, ship.position + leftBackPos, ship.colour)
   sg.drawLine(ship.position + forwardPos, ship.position + rightBackPos, ship.colour)
   'draw active thrusters
   For i As Integer = 0 To NUM_THRUSTERS - 1
      'Dim As ULong c = IIf(i < 4, RGB(255, 255, 0), RGB(255, 255, 255))
      If thruster(i).active > 0 Then
         drawThruster(thruster(i).absPos, thruster(i).absPos - thruster(i).forceVector / 1e3, RGB(255, 63, 0)) 'thruster force indicator
      End If
   Next
   'draw astroids
   For i As Integer = 0 To UBound(asteroid)
      sg.drawCircle(asteroid(i).position, asteroid(i).radius, IIf(nearestAsteroidId = i, RGB(191, 191, 255), asteroid(i).colour))
   Next
   Dim As sgl2d asteroidPointer = normalise(ship.position - asteroid(nearestAsteroidId).position)
   drawArrow(ship.position, ship.position - asteroidPointer * ship.radius * 2, RGB(191, 191, 255))
   ScreenUnLock

   'time update
   key = InKey()
   Sleep 1
   tPrev = tNow
   tNow = Timer
   dt = tNow - tPrev
Wend
Screen 0
Print "End"

'TODO:
'center of universe indicator
'drawElipse
'fuel


[badidea]

While patiently waiting for h4tt3n's space physics doodle stuff, I slowly convert my code to a small game:

Code: Select all

#Include "fbgfx.bi"

Type int2d
   As Integer x, y
   Declare Constructor
   Declare Constructor(x As Integer, y As Integer)
   Declare Operator Cast () As String
End Type

Constructor int2d
End Constructor

Constructor int2d(x As Integer, y As Integer)
   This.x = x : This.y = y
End Constructor

Operator = (a As int2d, b As int2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

Operator <> (a As int2d, b As int2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' "x, y"
Operator int2d.cast () As String
  Return Str(x) & "," & Str(y)
End Operator

' a + b
Operator + (a As int2d, b As int2d) As int2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As int2d, b As int2d) As int2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As int2d) As int2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As int2d, b As int2d) As int2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As int2d, mul As Integer) As int2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a \ b
Operator \ (a As int2d, b As int2d) As int2d
   Return Type(a.x \ b.x, a.y \ b.y)
End Operator

' a \ div
Operator \ (a As int2d, div As Integer) As int2d
   Return Type(a.x \ div, a.y \ div)
End Operator

'===============================================================================

Type sgl2d
   As Single x, y
   Declare Constructor
   Declare Constructor(x As Single, y As Single)
   Declare Operator Cast () As String
End Type

Constructor sgl2d
End Constructor

Constructor sgl2d(x As Single, y As Single)
   This.x = x : This.y = y
End Constructor

' "x, y"
Operator sgl2d.cast () As String
   Return Str(x) & "," & Str(y)
End Operator

'---- operators ---

' distance / lenth
Operator Len (a As sgl2d) As Single
   Return Sqr(a.x * a.x + a.y * a.y)
End Operator

' a = b ?
Operator = (a As sgl2d, b As sgl2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

' a != b ?
Operator <> (a As sgl2d, b As sgl2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' a + b
Operator + (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As sgl2d) As sgl2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As sgl2d, mul As Single) As sgl2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a / div
Operator / (a As sgl2d, div As Single) As sgl2d
   Return Type(a.x / div, a.y / div)
End Operator

'---- extra functions ---

Function cross(a As sgl2d, b As sgl2d) As Single
   Return a.x * b.y - a.y * b.x
End Function

'~ function length(a as sgl2d) as single
   '~ return sqr((a.x * a.x) + (a.y * a.y))
'~ end function

Function lengthSqrd(a As sgl2d) As Single
   Return (a.x * a.x) + (a.y * a.y)
End Function

Function dist(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return Sqr((dx * dx) + (dy * dy))
End Function

Function distSqrd(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return (dx * dx) + (dy * dy)
End Function

Function normalise(a As sgl2d) As sgl2d
   Dim As sgl2d temp
   Dim As Single length = Len(a)
   Return sgl2d(a.x / length, a.y / length)
End Function

'===============================================================================

'Note: y+ = up, x+ = right, (0,0) = center
Type scaled_graphics_type
   Dim As Single scale = 1 ' = 1 / pixel_size 'pixels / meter
   'dim as int2d offset' = (scrn_w \ 2, h \ 2) 'offset in pixels
   Dim As sgl2d offset
   Dim As Integer w = -1, h = -1
   Dim As Integer wc = -1, hc = -1 'center x,y
   Declare Sub setScreen(w As Integer, h As Integer)
   Declare Sub setScaling(scale As Single, offset As sgl2d)
   Declare Sub clearScreen(c As ULong)
   Declare Function pos2screen(p As sgl2d) As int2d
   Declare Sub drawPixel(p As sgl2d, c As Integer)
   Declare Sub drawCircle(p As sgl2d, r As Single, c As Integer)
   Declare Sub drawElipse(p As sgl2d, r As Single, aspect As Single, c As Integer)
   Declare Sub drawLine(p1 As sgl2d, p2 As sgl2d, c As Integer)
End Type

Sub scaled_graphics_type.setScreen(w As Integer, h As Integer)
   This.w = w 'width
   This.h = h 'height
   wc = w \ 2
   hc = h \ 2
   ScreenRes w, h, 32
   Width w \ 8, h \ 16 'bigger font
End Sub

Sub scaled_graphics_type.setScaling(scale As Single, offset As sgl2d)
   This.scale = scale
   This.offset = offset
End Sub

Sub scaled_graphics_type.clearScreen(c As ULong)
   Line(0, 0)-(w - 1, h - 1), c, bf
End Sub

Function scaled_graphics_type.pos2screen(p As sgl2d) As int2d
   Return int2d(Int(wc + (p.x - offset.x) * scale), h - Int(hc + (p.y - offset.y) * scale))
End Function

Sub scaled_graphics_type.drawPixel(p As sgl2d, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   PSet(posScrn.x, posScrn.y), c
End Sub

Sub scaled_graphics_type.drawCircle(p As sgl2d, r As Single, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawElipse(p As sgl2d, r As Single, aspect As Single, c As Integer)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c, , , aspect
End Sub

Sub scaled_graphics_type.drawLine(p1 As sgl2d, p2 As sgl2d, c As Integer)
   Dim As int2d posScrn1 = pos2screen(p1)
   Dim As int2d posScrn2 = pos2screen(p2)
   Line(posScrn1.x, posScrn1.y)-(posScrn2.x, posScrn2.y), c
End Sub

'===============================================================================

Const As Single PI = 4 * Atn(1)
Const As Single RAD_PER_DEG = (PI / 180)
Const As Single DEG_PER_RAD = 180 / PI
Const As Single sinA = Sin((10 / 180) * PI)
Const As Single cosA = Cos((10 / 180) * PI)
Const As Single sinB = Sin((20 / 180) * PI)
Const As Single cosB = Cos((20 / 180) * PI)

Const K_ESC = Chr(27)
Const K_MIN = Chr(45)
Const K_UND = Chr(95)
Const K_PLU = Chr(61)
Const K_EQU = Chr(43)

Const SCRN_W = 800, SCRN_H = 600

Dim Shared As scaled_graphics_type sg
sg.setScaling(2.0, sgl2d(0, 0))
sg.setScreen(SCRN_W, SCRN_H)

'-------------------------------------------------------------------------------

Function limit(value As Single, min As Single, max As Single) As Single
   If value < min Then Return min
   If value > max Then Return max
   Return value
End Function

Type polar
   Dim As Single angle
   Dim As Single magnitude
End Type

Function polarToCartesian(angle As Single, radius As Single) As sgl2d
   Return sgl2d(Cos(angle) * radius, Sin(angle) * radius)
End Function

Function rotatedVector(v As sgl2d, rotAngle As Single) As sgl2d
   Dim As sgl2d tmp
   tmp.x = Cos(rotAngle) * v.x - Sin(rotAngle) * v.y
   tmp.y = Sin(rotAngle) * v.x + Cos(rotAngle) * v.y
   Return tmp
End Function

'-------------------------------------------------------------------------------

Sub drawArrow(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.30 'reduce length
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x - sinB * dp.y, sinB * dp.x + cosB * dp.y), c)
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x + sinB * dp.y, cosB * dp.y - sinB * dp.x), c)
End Sub

Sub drawThruster(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.95 'reduce length
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x - sinA * dp.y, sinA * dp.x + cosA * dp.y), c)
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x + sinA * dp.y, cosA * dp.y - sinA * dp.x), c)
End Sub

Sub drawHelium3(p As sgl2d, r As Single, c As ULong)
   For i As Integer = 0 To 2
      sg.drawCircle(p + polarToCartesian((i / 3) * 2 * PI, 0.15 * r), 0.15 * r, c)
   Next
   sg.drawElipse(p, r, 2.0, c)
   sg.drawElipse(p, r, 0.5, c)
End Sub

'-------------------------------------------------------------------------------

Type disc_object
   Dim As Single radius '[m]
   Dim As Single height '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   'linear motion properties
   Dim As sgl2d position 'position [m]
   Dim As Single lin_m 'mass [kg]
   Dim As sgl2d lin_F 'force [N] [kg*m/s^2]
   Dim As sgl2d lin_a 'acceleration [m/s^2]
   Dim As sgl2d lin_v 'velocity [m/s]
   'dim as sgl2d lin_p 'momentum [kg*m/s]
   'dim as single lin_E 'Kinetic energy [J] [kg*m^2/s^2]
   'Rotational motion properties
   Dim As Single angle 'angular position (theta) [rad]
   Dim As Single ang_m 'angular mass, moment of inertia (I) [kg*m^2]
   Dim As Single ang_F 'torque (tau) [N*m] [kg*m^2/s^2]
   Dim As Single ang_a 'angular velocity (alpha) [rad/s^2]
   Dim As Single ang_v 'angular velocity (omega) [rad/s]
   'dim as single ang_p 'angular momentum (L) [kg*m^2/s]
   'dim as single ang_E 'Kinetic energy [J] [kg*m^2/s^2]
   Declare Sub init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   Declare Sub update(dt As Double)
   Declare Function getKineticEnergy() As Single
End Type

'Set radius, height, density, position
'Calculate mass and rotational inertia
Sub disc_object.init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   height = h
   density = d
   position = p
   colour = c
   lin_m = PI * r ^ 2 * d
   ang_m = 0.5 * lin_m * r ^ 2
End Sub

'update position and angle
Sub disc_object.update(dt As Double)
   lin_a = lin_F / lin_m
   lin_v += lin_a * dt
   position += lin_v * dt
   ang_a = ang_F / ang_m
   ang_v += ang_a * dt
   angle += ang_v * dt
End Sub

Function disc_object.getKineticEnergy() As Single
   Dim As Single lin_E = 0.5 * lin_m * lengthSqrd(lin_v)
   Dim As Single ang_E = 0.5 * ang_m * ang_v * ang_v
   Return lin_E + ang_E
End Function

Sub drawShip(ship As disc_object)
   'calculate ships tail pointer / triangle
   Dim As sgl2d forwardPos = polarToCartesian(ship.angle - 90 * RAD_PER_DEG, ship.radius * 2.2)
   Dim As sgl2d leftBackPos = polarToCartesian(ship.angle - 135 * RAD_PER_DEG, ship.radius * 1.0) '(90 + 45)
   Dim As sgl2d rightBackPos = polarToCartesian(ship.angle - 45 * RAD_PER_DEG, ship.radius * 1.0) '(90 - 45)
   sg.drawCircle(ship.position, ship.radius, ship.colour) 'flying saucer
   sg.drawLine(ship.position + forwardPos, ship.position + leftBackPos, ship.colour)
   sg.drawLine(ship.position + forwardPos, ship.position + rightBackPos, ship.colour)
End Sub

'-------------------------------------------------------------------------------

Type thruster_type
   '''init paramaters
   Dim As polar polarForce '(rad, N)
   Dim As polar polarPos '(rad, m)
   '''variable paramaters
   Dim As sgl2d forceVector '(N, N)
   Dim As sgl2d relPos, absPos '(m, m)
   Dim As Integer active
   Declare Sub init(forceMagnitude As Single, forceDirection As Single, posAngle As Single, posRadius As Single)
   Declare Sub updatePosition(bodyPos As sgl2d, bodyAngle As Single)
End Type

Sub thruster_type.init(forceDirection As Single, forceMagnitude As Single, posAngle As Single, posRadius As Single)
   polarForce = Type(forceDirection, forceMagnitude) 'thruster action
   polarPos = Type(posAngle, posRadius) 'position of thruster on ship
End Sub

Sub thruster_type.updatePosition(bodyPos As sgl2d, bodyAngle As Single)
   relPos = polarToCartesian(bodyAngle + polarPos.angle, polarPos.magnitude)
   absPos = bodyPos + relPos
End Sub

'-------------------------------------------------------------------------------

Const As Single GRAV_CONST = 6.67e-11 '[m3/(kg*s^2)

Type astro_body
   Dim As Single radius '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   Dim As sgl2d position 'position [m]
   Dim As Single mass '[kg]
   Declare Sub init(r As Single, d As Single, p As sgl2d, c As ULong)
End Type

'Set radius, density, position
'Calculate mass and rotational inertia
Sub astro_body.init(r As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   density = d
   position = p
   colour = c
   mass = PI * r ^ 2 * d
End Sub

Function gravForce(m1 As Single, m2 As Single, r As Single) As Single
   Return GRAV_CONST * (m1 * m2) / (r * r)
End Function

Function gravForceVector(m1 As Single, Pos1 As sgl2d, m2 As Single, pos2 As sgl2d) As sgl2d
   Dim As Single distSquared = distSqrd(pos2, Pos1)
   Dim As sgl2d unitVector12 = (Pos1 - pos2) / Sqr(distSquared)
   Return unitVector12 * (-GRAV_CONST * (m1 * m2) / distSquared)
End Function

Function findNearestBody(refPos As sgl2d, body() As astro_body) As Integer
   Dim As Integer nearestBodyId = 0
   Dim As Single nearestBodyDistSqrd = distSqrd(refPos, body(0).position)
   Dim As Single currentBodyDistSqrd
   For i As Integer = 1 To UBound(body)
      currentBodyDistSqrd = distSqrd(refPos, body(i).position)
      If currentBodyDistSqrd < nearestBodyDistSqrd Then
         nearestBodyDistSqrd = currentBodyDistSqrd
         nearestBodyId = i
      End If
   Next
   Return nearestBodyId
End Function

'-------------------------------------------------------------------------------

Const As Single MOON_RADIUS = 1737e+6 '[m]
Const As Single MOON_DENSITY = 3344 '[kg/m^3]

Const NUM_THRUSTERS = 6
Const L_FW_THR = 0 'left forward thruster
Const R_FW_THR = 1 'right forward thruster
Const L_LO_THR = 2
Const R_LO_THR = 3
Const L_HI_THR = 4
Const R_HI_THR = 5

Dim As String key
Dim As Integer quit = 0
Dim As disc_object ship
Dim As thruster_type thruster(NUM_THRUSTERS - 1)
Const NUM_ASTEROID = 50
Dim As astro_body asteroid(NUM_ASTEROID-1)
Const NUM_HELIUM = 20
Dim As astro_body helium(NUM_HELIUM-1)
Const As Single maxFuel = 1e6 'N*s
Dim As Single fuel = maxFuel

ship.init(10, 1, 5, sgl2d(0, -50), RGB(127, 223, 0))
For i As Integer = 0 To UBound(asteroid)
   asteroid(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(255, 191, 0))
Next
For i As Integer = 0 To UBound(helium)
   helium(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(255, 191, 0))
Next

'force angle, force magnitude, polar thruster position
thruster(L_FW_THR).init(0.5 * pi, 1.2e4, -0.75 * pi, ship.radius)
thruster(R_FW_THR).init(0.5 * pi, 1.2e4, -0.25 * pi, ship.radius)
thruster(L_LO_THR).init(0.0 * pi, 8e3, -0.75 * pi, ship.radius)
thruster(R_LO_THR).init(1.0 * pi, 8e3, -0.25 * pi, ship.radius)
thruster(L_HI_THR).init(0.0 * pi, 8e3, +0.75 * pi, ship.radius)
thruster(R_HI_THR).init(1.0 * pi, 8e3, +0.25 * pi, ship.radius)

Dim As Double tNow = Timer, tPrev = tNow, dt = 0
While quit = 0
   'reset stuff
   ship.lin_F = sgl2d(0, 0)
   ship.ang_F = 0
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).active = 0
   Next

   If MultiKey(FB.SC_UP) Then
      thruster(L_FW_THR).active = 1
      thruster(R_FW_THR).active = 1
      fuel -= thruster(L_FW_THR).polarForce.magnitude * dt
      fuel -= thruster(R_FW_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_LEFT) Then
      thruster(L_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
      fuel -= thruster(L_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(R_HI_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_RIGHT) Then
      thruster(R_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
      fuel -= thruster(R_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(L_HI_THR).polarForce.magnitude * dt
   End If

   If key = K_ESC Then quit = 1

   For i As Integer = 0 To NUM_THRUSTERS - 1
      'forces on body by active thrusters
      If thruster(i).active > 0 Then
         Dim As Single thrust = thruster(i).polarForce.magnitude
         thruster(i).forceVector = polarToCartesian(ship.angle + thruster(i).polarForce.angle, thrust)
         ship.lin_F += thruster(i).forceVector
         ship.ang_F += cross(thruster(i).relPos, thruster(i).forceVector)
      End If
   Next

   ship.update(dt) 'position and angle
   sg.offset = ship.position
   sg.scale = limit(100 / Len(ship.lin_v), 0.5, 2.0)

   'do always for display
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).updatePosition(ship.position, ship.angle)
   Next

   Dim As Integer nearestHeliumId = findNearestBody(ship.position, helium())
   Dim As Integer nearestAsteroidId = findNearestBody(ship.position, asteroid())

   Dim As astro_body Ptr pAsteroid = @asteroid(nearestAsteroidId)
   If dist(ship.position, pAsteroid->position) < ship.radius + pAsteroid->radius Then
      quit = 2 'crash
   End If

   If fuel <= 0 Then quit = 3
   
   'display
   ScreenLock
   sg.clearScreen(0)
   Locate 2, 2 : Print "<UP>, <LEFT>, <RIGHT> for thrusters, <ESC> to exit";
   Locate 3, 2 : Print "Nearest helium distance [m]: "; CInt(dist(ship.position, helium(nearestHeliumId).position));
   Locate 4, 2 : Print "Kinetic engergy [kJ]: "; CInt(ship.getKineticEnergy() * 1e-3);
   Locate 5, 2 : Print "Fuel remaining [%]: "; CInt((fuel / maxFuel) * 100);
   drawShip(ship)
   'draw active thrusters
   For i As Integer = 0 To NUM_THRUSTERS - 1
      Dim As ULong c = IIf(i < 4, RGB(255, 255, 0), RGB(255, 255, 255))
      If thruster(i).active > 0 Then
         drawThruster(thruster(i).absPos, thruster(i).absPos - thruster(i).forceVector / 1e3, RGB(255, 63, 0)) 'thruster force indicator
      End If
   Next
   'draw asteroids
   For i As Integer = 0 To UBound(asteroid)
      sg.drawCircle(asteroid(i).position, asteroid(i).radius, asteroid(i).colour)
   Next
   'draw helium 'clouds'
   For i As Integer = 0 To UBound(helium)
      drawHelium3(helium(i).position, helium(i).radius, IIf(nearestHeliumId = i, RGB(191, 191, 255), helium(i).colour))
   Next
   Dim As sgl2d heliumPointer = normalise(ship.position - helium(nearestHeliumId).position)
   drawArrow(ship.position, ship.position - heliumPointer * ship.radius * 2, RGB(191, 191, 255))
   ScreenUnLock

   'time update
   key = InKey()
   Sleep 1
   tPrev = tNow
   tNow = Timer
   dt = tNow - tPrev
Wend

If quit > 1 Then 'Crashed or No fuel
   If quit = 2 Then ship.colour = RGB(233, 0, 0)
   If quit = 3 Then ship.colour = RGB(233, 191, 191)
   drawShip(ship)
   While InKey <> Chr(13)
      Locate 7, 2
      If quit = 2 Then Print "Ship crashed, press <ENTER> to exit."
      If quit = 3 Then Print "Ship out of fuel, press <ENTER> to exit."
      Sleep 1
   Wend
End If
Screen 0
Print "End"

'TODO:
'collect / remove helium3
'shoot asteroids


[h4tt3n]

I'm still around :-) Here, as promised, are some space physics examples. Some are rather old, but I'll just post them as they are, warts and all.

First, a vector library and a ball shaped space ship, asteroids in orbit around a central planet, also featuring collision with friction.

You control the yellow ball with the arrow keys. Up = thrust, left/right = rotate, down = stop rotating.

Vec2.bi

Code: Select all

''*******************************************************************************
''   
''   FreeBASIC 2D Floating Point Vector Library
''   Written in FreeBASIC 1.05
''   Version 1.3.0, November 2016, Michael "h4tt3n" Nissen
''   
''   Function syntax:
''   
''   (Return Type) (Name) (Argument Type (, ...)) (Description)
''   
''   Vector  Absolute         ( Vector )          Absolute value
''   Vector  AngleToUnit      ( Scalar )          unit vector from angle scalar
''   Vector  Component        ( Vector, Vector )  Vector Component
''   Scalar  Dot              ( Vector, Vector )  Dot product
''   Vector  Dot              ( Vector, Scalar )  Dot product
''   Vector  Dot              ( Scalar, Vector )  Dot product
''   Vector  PerpCW           ( Vector )          Right hand perpendicular
''   Vector  PerpCCW          ( Vector )          Left hand perpendicular
''   Vector  Unit             ( Vector )          Unit Vector
''   Scalar  Length           ( Vector )          Length
''   Scalar  LengthSquared    ( Vector )          Length squared
''   Scalar  PerpDot          ( Vector, Vector )  Perp dot product (2d cross)
''   Vector  PerpDot          ( Vector, Scalar )  Perp dot product (2d cross)
''   Vector  PerpDot          ( Scalar, Vector )  Perp dot product (2d cross)
''   Vector  Project          ( Vector, Vector )  Vector Projection
''   Vector  RandomizeSquare  ( Scalar )          Randomize in range +/- value
''   Vector  RandomizeCircle  ( Scalar )          Randomize in range +/- value
''   Vector  Rotate           ( Vector )          Rotate
''   Scalar  UnitToAngle      ( Vector )          angle scalar from unit vector
''   
''   Library supports both member and non-member function useage:
''   
''   A.function(B) <-> function(A, B)
''   
''*******************************************************************************


''
#Ifndef __VEC2_BI__
#Define __VEC2_BI__


''  Vec2 Vector class
Type Vec2
   
   Public:
   
   ''  Vec2 constructor declarations
   Declare Constructor ()
   Declare Constructor ( ByRef v As Const Vec2 )
   Declare Constructor ( ByVal x As Const Single, ByVal y As Const Single )
   
   '' Vec2 destructor
   Declare Destructor()
   
   ''  Vec2 assignment operator (=)
   Declare Operator Let ( ByRef B As Const Vec2 )
   
   ''   Vec2 compound arithmetic member operator declarations
   Declare Operator +=  ( ByRef B As Const Vec2   )
   Declare Operator -=  ( ByRef B As Const Vec2   )
   Declare Operator *=  ( ByRef B As Const Vec2   )
   Declare Operator *=  ( ByVal B As Const Single )
   Declare Operator /=  ( ByVal B As Const Single )
   
   ''  Vec2 member function declarations
   Declare Const Function Absolute_        () As Vec2
   
   Declare Const Function PerpCW          () As Vec2
   Declare Const Function PerpCCW         () As Vec2
   
   Declare Const Function Unit            () As Vec2
   Declare Const Function Length          () As Single
   Declare Const Function LengthSquared   () As Single
   
   Declare Const Function Dot             ( ByRef B As Const Vec2   ) As Single
   Declare Const Function Dot             ( ByVal B As Const Single ) As Vec2
   
   Declare Const Function PerpDot         ( ByRef B As Const Vec2   ) As Single
   Declare Const Function PerpDot         ( ByVal B As Const Single ) As Vec2
   
   Declare Const Function Project         ( ByRef B As Const Vec2   ) As Vec2
   Declare Const Function Component       ( ByRef B As Const Vec2   ) As Vec2
   
   Declare Const Function RandomizeCircle ( ByVal B As Const Single ) As Vec2
   Declare Const Function RandomizeSquare ( ByVal B As Const Single ) As Vec2
   
   Declare Const Function RotateCW        ( ByRef B As Const Vec2   ) As Vec2
   Declare Const Function RotateCCW       ( ByRef B As Const Vec2   ) As Vec2
   
   Declare Const Function RotateCW        ( ByRef B As Const Single ) As Vec2
   Declare Const Function RotateCCW       ( ByRef B As Const Single ) As Vec2
   
   ''  Vec2 variables
   As Single x, y
   
End Type

''  Vec2 unary arithmetic non-member operator declarations
Declare Operator - ( ByRef B As Const Vec2 ) As Vec2

''  Vec2 binary arithmetic non-member operator declarations
Declare Operator + ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Vec2
Declare Operator - ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Vec2
Declare Operator * ( ByVal A As Const Single, ByRef B As Const Vec2   ) As Vec2
Declare Operator * ( ByRef A As Const Vec2  , ByVal B As Const Single ) As Vec2
Declare Operator * ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Vec2
Declare Operator / ( ByRef A As Const Vec2  , ByVal B As Const Single ) As Vec2

''  Vec2 binary relational non-member operator declarations
Declare Operator =  ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
Declare Operator <> ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
Declare Operator <  ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
Declare Operator >  ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer

''  Vec2 non-member function declarations
Declare Function Absolute_ OverLoad ( ByRef A As Const Vec2 ) As Vec2
Declare Function Unit              ( ByRef A As Const Vec2 ) As Vec2
Declare Function Length            ( ByRef A As Const Vec2 ) As Single
Declare Function LengthSquared     ( ByRef A As Const Vec2 ) As Single
Declare Function Dot      Overload ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Single
Declare Function Dot               ( ByRef A As Const Vec2  , ByVal B As Const Single ) As Vec2
Declare Function Dot               ( ByVal A As Const Single, ByRef B As Const Vec2   ) As Vec2
Declare Function PerpDot  OverLoad ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Single
Declare Function PerpDot           ( ByRef A As Const Vec2  , ByVal B As Const Single ) As Vec2
Declare Function PerpDot           ( ByVal A As Const Single, ByRef B As Const Vec2   ) As Vec2
Declare Function Project           ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Vec2
Declare Function Component         ( ByRef A As Const Vec2  , ByRef B As Const Vec2   ) As Vec2
Declare Function RandomizeCircle   ( ByRef A As Const Vec2  , ByVal B As Const Single ) As Vec2
Declare Function RandomizeSquare   ( ByRef A As Const Vec2  , ByVal B As Const Single ) As Vec2

Declare Function RotateCW  OverLoad ( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
Declare Function RotateCCW OverLoad ( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2

Declare Function AngleToUnit       ( ByVal A As Const Single ) As Vec2
Declare Function UnitToAngle       ( ByRef A As Const Vec2 ) As Single


''  Vec2 constructors
Constructor Vec2()

  This.x = 0.0 : This.y = 0.0
 
End Constructor

Constructor Vec2( ByRef v As Const Vec2 )
   
   This = v
 
End Constructor

Constructor Vec2( ByVal x As Const Single, ByVal y As Const Single )

  This.x = x : This.y = y
 
End Constructor


'' Destructor
Destructor Vec2()

End Destructor


''
Operator Vec2.Let ( ByRef B As Const Vec2 )
   
   If ( Not @This = @B ) Then
      
      This.x = B.x : This.y = B.y
      
   EndIf
   
End Operator


''  Vec2 compound arithmetic member operators
Operator Vec2.+= ( ByRef B As Const Vec2 )
   
   This.x += B.x : This.y += B.y
   
End Operator

Operator Vec2.-= ( ByRef B As Const Vec2 )
   
   This.x -= B.x : This.y -= B.y
   
End Operator

Operator Vec2.*= ( ByRef B As Const Vec2 )
   
   This.x *= B.x : This.y *= B.y
   
End Operator

Operator Vec2.*= ( ByVal B As Const Single )
   
   This.x *= B : This.y *= B
   
End Operator

Operator Vec2./= ( ByVal B As Const Single )
      
   This.x /= B : This.y /= B
   
End Operator


''  Vec2 member functions
Function Vec2.Absolute_() As Vec2
   
   Return Vec2( Abs( This.x ), Abs( This.y ) )
   
End Function

Function Vec2.PerpCW() As Vec2
   
   Return Vec2( This.y, -This.x )
   
End Function

Function Vec2.PerpCCW() As Vec2
   
   Return Vec2( -This.y, This.x )
   
End Function

Function Vec2.Unit() As Vec2
   
   Return IIf( This <> Vec2( 0.0, 0.0 ) , Vec2( This / This.Length() ) , Vec2( 0.0, 0.0 ) )
   
End Function

Function Vec2.Length() As Single
   
   'Dim As Single length_squared = This.LengthSquared()
   'Return Sqr( length_squared )
   
   Return Sqr( This.LengthSquared() )
   
End Function

Function Vec2.LengthSquared() As Single
   
   'Return This.Dot( This )
   Return ( This.x * This.x + This.y * This.y )
   
End Function

Function Vec2.Dot( ByVal B As Const Single ) As Vec2
   
   Return Vec2( This.x * B, This.y * B )
   
End Function

Function Vec2.Dot( ByRef B As Const Vec2 ) As Single
   
   Return ( This.x * B.x + This.y * B.y )
   
End Function

Function Vec2.PerpDot( ByVal B As Const Single  ) As Vec2
   
   Return Vec2( -This.y * B, This.x * B )
   
End Function

Function Vec2.PerpDot( ByRef B As Const Vec2 ) As Single
   
   Return ( -This.y * B.x + This.x * B.y  )
   'Return ( This.x * B.y - This.y * B.x )
   'Return ( This.Dot( B.Perp() ) )
   
End Function

Function Vec2.Project( ByRef B As Const Vec2 ) As Vec2
   
   Return ( B.Dot( This ) / This.Dot( This ) ) * This
   
End Function

Function Vec2.Component( ByRef B As Const Vec2 ) As Vec2
   
   Return ( This.Dot( B ) / B.Dot( B ) ) * B
   
End Function

Function Vec2.RandomizeCircle( ByVal B As Const Single ) As Vec2
   
   Dim As Single a = Rnd() * 8.0 * Atn( 1.0 )
   Dim As Single r = Sqr( Rnd() * B * B )
   
   Return Vec2( Cos( a ), Sin( a ) ) * r
   
End Function

Function Vec2.RandomizeSquare( ByVal B As Const Single ) As Vec2
   
   Return Vec2( ( Rnd() - Rnd() ) * B, ( Rnd() - Rnd() ) * B )
   
End Function

Function Vec2.RotateCW( ByRef B As Const Vec2 ) As Vec2
   
   Return Vec2( B.Dot( This ), B.PerpDot( This ) )
   
End Function

Function Vec2.RotateCCW( ByRef B As Const Vec2 ) As Vec2
   
   'Dim As vec2 v = Vec2( B.x, -B.y )
   
   'Return Vec2( v.Dot( This ), v.PerpDot( This ) )
   
   Return Vec2( B.x * This.x - B.y * This.y , B.y * This.x + B.x * This.y )
   
   ''perpccw = -y, x
   
End Function

Function Vec2.RotateCW( ByRef B As Const Single ) As Vec2
   
   Dim As Vec2 v = Vec2( Cos( B ), Sin( B ) )
   
   Return This.RotateCW( v )
   
End Function

Function Vec2.RotateCCW( ByRef B As Const Single ) As Vec2
   
   Dim As Vec2 v = Vec2( Cos( B ), Sin( B ) )
   
   Return This.RotateCCW( v )
   
End Function


''  Vec2 unary arithmetic non-member operators
Operator - ( ByRef B As Const Vec2 ) As Vec2
   
   Return Vec2( -B.x, -B.y )
   
End Operator


''  Vec2 binary arithmetic non-member operators
Operator + ( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Vec2
   
   Return Vec2( A.x + B.x, A.y + B.y )
   
End Operator

Operator - ( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Vec2
   
   Return Vec2( A.x - B.x, A.y - B.y )
   
End Operator

Operator * ( ByVal A As Const Single, ByRef B As Const Vec2 ) As Vec2
   
   Return Vec2( A * B.x, A * B.y)
   
End Operator

Operator * ( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return Vec2( A.x * B, A.y * B)
   
End Operator

Operator * ( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Vec2
   
   Return Vec2( A.x * B.x, A.y * B.y )
   
End Operator

Operator / ( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return Vec2( A.x / B, A.y / B )
   
End Operator


''  Vec2 binary relational non-member operators
Operator = ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
   
   Return ( A.x = B.x ) And ( A.y = B.y )
   
End Operator

Operator <> ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
   
   Return ( A.x <> B.x ) Or ( A.y <> B.y )
   
End Operator

Operator < ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
   
   Return ( A.x < B.x ) And ( A.y < B.y )
   
End Operator

Operator > ( ByRef A As Const Vec2, ByVal B As Const Vec2 ) As Integer
   
   Return ( A.x > B.x ) And ( A.y > B.y )
   
End Operator


''  Vec2 non-member functions
Function Absolute_( ByRef A As Const Vec2 ) As Vec2
   
   Return A.Absolute_()
   
End Function

Function Unit( ByRef A As Const Vec2 ) As Vec2
   
   Return A.Unit()
   
End Function

Function Length( ByRef A As Const Vec2 ) As Single
   
   Return A.Length()
   
End Function

Function LengthSquared( ByRef A As Const Vec2 ) As Single
   
   Return A.LengthSquared()
   
End Function

Function Dot( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Single
   
   Return A.Dot( B )
   
End Function

Function Dot( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return A.Dot( B )
   
End Function

Function Dot( ByVal A As Const Single, ByRef B As Const Vec2 ) As Vec2
   
   Return B.Dot( -A )
   
End Function

Function PerpDot( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Single
   
   Return A.PerpDot( B )
   
End Function

Function PerpDot( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return A.PerpDot( B )
   
End Function

Function PerpDot( ByVal A As Const Single, ByRef B As Const Vec2 ) As Vec2
   
   Return B.PerpDot( -A )
   
End Function

Function Project( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Vec2
   
   Return A.Project( B )
   
End Function

Function Component( ByRef A As Const Vec2, ByRef B As Const Vec2 ) As Vec2
   
   Return A.Component( B )
   
End Function

Function RandomizeCircle( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return A.RandomizeCircle( B )
   
End Function

Function RandomizeSquare( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return A.RandomizeSquare( B )
   
End Function

Function RotateCW( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return A.RotateCW( B )
   
End Function

Function RotateCCW( ByRef A As Const Vec2, ByVal B As Const Single ) As Vec2
   
   Return A.RotateCCW( B )
   
End Function

Function AngleToUnit( ByVal A As Const Single ) As Vec2
   
   Return Vec2( Cos( A ), Sin( A ) )
   
End Function

Function UnitToAngle( ByRef A As Const Vec2 ) As Single
   
   Return ATan2( A.y, A.x )
   
End Function


#EndIf __VEC2_BI__


And here the physics doodle:

Code: Select all

'' Ball - Ball collision detection and response with friction.
'' Written by Michael "h4tt3n" Schmidt Nissen
''

''
#include "fbgfx.bi"
#include "Vec2.bi"

''
Const pi                       = 4*atn(1)
Const G                         = 4
const num_balls                = 32
const rest_fps                  = 60                  ''  ideal framerate
'const inv_rest_fps              = 1.0 / (rest_fps+2)  ''  inverse ideal framerate
Const inv_rest_fps              = 1.0 / rest_fps  ''  inverse ideal framerate
Const dt                        = inv_rest_fps        ''  timestep, delta time
Const inv_dt                  = 1.0 / dt        ''  timestep, delta time
Const StaticFrictionVelocity    = 1.0

randomize timer

''
type ball_type
As UInteger col
as Vec2 Impulse
as Vec2 Velocity
as Vec2 Position
as Vec2 AngleVector
as Single angular_impulse, ang_Velocity, sin_ang, cos_ang, ang_Position, mass, InverseMass, dens, radius, _
RadiusSquared, I, InverseI, CollisionStiffness, COllisionDamping, dynamicfriction, staticfriction
end type

dim shared as Vec2 dst, separation_vector, closest_point
dim shared as ball_type ptr ball
dim shared as Single dst_sqd, radius, RadiusSquared, distance
dim shared as integer a, b, scrn_wid, scrn_hgt
dim shared as integer FPS, FPS_Counter
dim shared as Single FPS_Timer, t0

declare sub InitBall(byref a as ball_type)
declare sub BallBallCollisionDetection()
Declare Sub Controls()
Declare Sub Gravity()
declare sub Integrate()
declare sub DrawSceneToScreen()
declare sub brake()

scrn_wid = 1000:   scrn_hgt = 800
ScreenRes scrn_wid, scrn_hgt, 32
WindowTitle "force based 2d ball-ball collision with friction"
Color RGB(0, 0, 0), RGB(192, 192, 244)
'color 0, rgb(255, 255, 255)

ball = New Ball_Type[Num_Balls]

with ball[0]
Dim As UByte Clr = ((Num_Balls-a)/num_Balls)*255
.col = rgb(Clr, 255, clr)
.mass = 800000
.InverseMass = 1/.mass
.dens = 0.2
.radius = ((.mass/.dens)/((4/3)*pi))^(1/3)
.RadiusSquared = .radius*.radius
.I = 0.5*.mass*.RadiusSquared         ''   flat disc
'.I = (2*.mass*.RadiusSquared)/5      ''   solid sphere
.InverseI = 1/.I
.Position.x = scrn_wid\2
.Position.y = scrn_hgt\2
.ang_Velocity = 0.0
.ang_Position = rnd*2*pi
.cos_ang = cos(.ang_Position)
.sin_ang = sin(.ang_Position)
.CollisionStiffness = 0.5
.COllisionDamping = 0.1
.dynamicfriction = 0.5
.staticfriction = 1.0
end with

''  set startup condition
for a = 1 to num_balls-1
Dim As Single angle = 2*pi*Rnd
Dim As Single radius = ball[0].radius + 128 + Rnd*256
with ball[a]
''   brightest first, darkest last
Dim As UByte Clr = ((Num_Balls-a)/num_Balls)*255
.col = rgb(Clr, 255, clr)
.mass = 10+(rnd*80^(1/2))^2
.InverseMass = 1/.mass
.dens = 0.01
.radius = ((.mass/.dens)/((4/3)*pi))^(1/3)
.RadiusSquared = .radius*.radius
'.I = 0.5*.mass*.RadiusSquared         ''   flat disc
.I = (2*.mass*.RadiusSquared)/5      ''   solid sphere
.InverseI = 1/.I
.Position.x = ball[0].Position.x + Cos(angle)*Radius
.Position.y = ball[0].Position.y + Sin(angle)*Radius
.Velocity.X = Sqr((G*ball[0].mass)/radius)*Sin(-Angle)
.Velocity.Y = Sqr((G*ball[0].mass)/radius)*Cos(Angle)
ball[0].Velocity -= .Velocity*(.mass/ball[0].mass)
.ang_Position = rnd*2*pi
.cos_ang = cos(.ang_Position)
.sin_ang = sin(.ang_Position)
.ang_Velocity = 0.0'(Rnd-Rnd) * 0.2
.CollisionStiffness = 0.5
.COllisionDamping = 0.1
.dynamicfriction = 0.5
.staticfriction = 1.0
end with
Next

With ball[1]
.col = RGB(255, 255, 64)
.mass = 20
.InverseMass = 1/.mass
.dens = 0.001
.radius = ((.mass/.dens)/((4/3)*pi))^(1/3)
.RadiusSquared = .radius*.radius
.I = 0.5*.mass*.radius*.radius
.InverseI = 1/.I
End With


''----------------------------------------------------------------------------''

Do

BallBallCollisionDetection()
Integrate()
DrawSceneToScreen()
brake()

loop until multikey(1)

Delete[] Ball

end


''----------------------------------------------------------------------------''


sub BallBallCollisionDetection()
   
   For a As Integer = 0 to num_balls-2
      
      For b As Integer = a+1 to num_balls-1
         
         Dim As Vec2 DistanceVector = ball[a].Position-ball[b].Position
         Dim As Single DistanceSquared = LengthSquared(DistanceVector)
         Dim As Single rest_distance = ball[a].radius + ball[b].radius
         
         If DistanceSquared > ( rest_distance * rest_distance ) Then Continue For
            
         Dim As Single Distance = sqr(DistanceSquared)
         Dim as Vec2 NormalVector = IIf( Not Distance = 0.0 , DistanceVector/Distance , Vec2( 0.0, 0.0 ) )
         Dim as Vec2 TangentVector = NormalVector.PerpCW()
         
         Dim as Single distance_error = distance - rest_distance
         
         Dim As Vec2 contact_a = ball[a].Position-ball[a].radius*NormalVector
         Dim As Vec2 contact_b = ball[b].Position+ball[b].radius*NormalVector
         
         Dim as Vec2 contact_velocity_a = ball[a].Velocity+ball[a].ang_Velocity*(ball[a].Position-contact_a).PerpCW()
         Dim as Vec2 contact_velocity_b = ball[b].Velocity+ball[b].ang_Velocity*(ball[b].Position-contact_b).PerpCW()
         
         Dim as Vec2 contact_velocity = contact_velocity_a - contact_velocity_b
         Dim as Single contact_velocityNormal = Dot(contact_velocity, normalVector)
         
         If contact_velocityNormal > 0.0 Then Continue For
            
         Dim As Single K = ( ball[a].CollisionStiffness + ball[b].CollisionStiffness ) * 0.5
         Dim As Single D = ( ball[a].CollisionDamping + ball[b].CollisionDamping ) * 0.5
         
         'dim as Single impulseNormal = -(1+Restitution)*contact_velocityNormal/(ball[a].InverseMass+ball[b].InverseMass)
         Dim as Single impulseNormal = -(distance_error*inv_dt*k + contact_velocityNormal*D ) * ( 1.0 / (ball[a].InverseMass+ball[b].InverseMass) )
         
         Dim As Single contact_velocityTangent = Dot(contact_velocity, tangentvector)
         
         Dim As Single FrictionCoefficient = IIf( Abs(contact_velocityTangent) < StaticFrictionVelocity, _
                                                  ( Ball[a].StaticFriction  + Ball[b].StaticFriction  ) * 0.5, _
                                                  ( Ball[a].DynamicFriction + Ball[b].DynamicFriction ) * 0.5 )
         
         Dim as Single MaxImpulseTangent = -contact_velocityTangent * ( 1.0 / (Ball[a].InverseMass+Ball[b].InverseMass+_
         Ball[a].RadiusSquared*Ball[a].InverseI+Ball[b].RadiusSquared*Ball[b].InverseI) )
         
         Dim As Single ImpulseTangent = IIf( Abs( MaxImpulseTangent ) < ( FrictionCoefficient * ImpulseNormal ), _
                                             MaxImpulseTangent, _
                                             Sgn( MaxImpulseTangent ) * ( FrictionCoefficient * ImpulseNormal ) )
         
         'Dim As Single ImpulseTangent = MaxImpulseTangent
         
         Dim As Vec2 Impulse = ImpulseNormal * NormalVector + ImpulseTangent * TangentVector
         
         ball[a].Impulse += Impulse*ball[a].InverseMass
         ball[a].angular_impulse += PerpDot(Impulse, ball[a].Position-contact_a)*ball[a].InverseI
         
         ball[b].Impulse -= Impulse*ball[b].InverseMass
         ball[b].angular_impulse -= Perpdot(Impulse, ball[b].Position-contact_b)*ball[b].InverseI
         
      Next
      
   Next
   
End Sub

Sub Controls()
   
   With Ball[1]
         
      If MultiKey(&h48) Then
         
         .Impulse.X += Cos(.Ang_Position)*100*.InverseMass
         .Impulse.Y += Sin(.Ang_Position)*100*.InverseMass
         
      EndIf
      
      'If MultiKey(&h4b) Then .Ang_Position -= 0.05
      'If MultiKey(&h4d) Then .Ang_Position += 0.05
      If MultiKey(&h4b) Then .ang_Velocity -= 0.1
      If MultiKey(&h4d) Then .ang_Velocity += 0.1
      If MultiKey(&h50) Then .ang_Velocity = 0
      '.ang_Velocity = 0
      
      .cos_ang = cos(.ang_Position)
      .sin_ang = sin(.ang_Position)
      
   End With
   
End Sub

Sub Gravity()

For a As Integer = 0 to num_balls-2

For b As Integer = a+1 to num_balls-1

Dim As Vec2 DistanceVector = ball[a].Position-ball[b].Position
Dim As Single DistanceSquared = LengthSquared(DistanceVector)
Dim As Single DistanceCubed = DistanceSquared*Sqr(DistanceSquared)

Ball[a].Impulse -= Ball[b].mass * G * DT * DistanceVector / DistanceCubed
Ball[b].Impulse += Ball[a].mass * G * DT * DistanceVector / DistanceCubed

Next

Next

End Sub

Sub Integrate()

Controls()

Gravity()

For a As Integer = 0 to num_balls-1

With Ball[a]

.Velocity += .Impulse
.Position += .Velocity * dt
.Impulse = Vec2( 0.0, 0.0 )

.ang_Velocity += .angular_impulse
.ang_Position += .ang_Velocity * dt
.angular_impulse = 0.0

.cos_ang = cos(.ang_Position)
.sin_ang = sin(.ang_Position)

end With

Next

End Sub

sub DrawSceneToScreen()

ScreenLock()

Cls

locate  2, (scrn_wid\8)-3: print using "###"; fps

for a As Integer = 0 to num_balls-1

with ball[a]

circle (.Position.x, .Position.y), .radius, 0,,,1, f
circle (.Position.x, .Position.y), .radius-2, .col,,,1, f
Line(.Position.x, .Position.y)-(.Position.x+.cos_ang*.radius, .Position.y+.sin_ang*.radius), 0

end With

Next

ScreenUnLock()

end sub

sub brake()

if timer < fps_timer then
fps_counter += 1
Else
fps_timer = timer+1
fps = fps_counter
fps_counter = 0
end if

Do while ( Timer - t0 ) < ( inv_rest_fps )

Sleep 1, 1

Loop

t0 = timer

end Sub


[badidea]

I'm still around :-) Here, as promised, are some space physics examples. Some are rather old, but I'll just post them as they are, warts and all.

First, a vector library and a ball shaped space ship, asteroids in orbit around a central planet, also featuring collision with friction.

You control the yellow ball with the arrow keys. Up = thrust, left/right = rotate, down = stop rotating.

....

Cool, I see that you have a 'slightly' modified the gravitational constant, else thing would be very slow.
For linux systems: Dim shared as Single Double FPS_Timer, t0

[h4tt3n]

Cool, I see that you have a 'slightly' modified the gravitational constant, else thing would be very slow.
For linux systems: Dim shared as Single Double FPS_Timer, t0

It's just a constant, you can give it any value you want. The real G is just as arbitrary, it depends on the random units we use. If you prefer to use G = 6.67e-11 then just crank up mass or delta time.

[badidea]

Cool, I see that you have a 'slightly' modified the gravitational constant, else thing would be very slow.
For linux systems: Dim shared as Single Double FPS_Timer, t0

It's just a constant, you can give it any value you want. The real G is just as arbitrary, it depends on the random units we use. If you prefer to use G = 6.67e-11 then just crank up mass or delta time.

I know, that is why I had in my previous code a heavy 'neutron star':

Code: Select all

neutronstar.init(MOON_RADIUS * 1e-8, MOON_DENSITY * 1e8, sgl2d(+40, +20), rgb(255, 127, 0))

Not sure if density and radius are realistic for a neutron star, but that were the numbers needed for a real-time gravity effect.
Never mind, I had the moon radius wrong by a factor of 1000.

[badidea]

Possibly the most frustrating game ever:

Code: Select all

#Include "fbgfx.bi"

Type int2d
   As Integer x, y
   Declare Constructor
   Declare Constructor(x As Integer, y As Integer)
   Declare Operator Cast () As String
End Type

Constructor int2d
End Constructor

Constructor int2d(x As Integer, y As Integer)
   This.x = x : This.y = y
End Constructor

Operator = (a As int2d, b As int2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

Operator <> (a As int2d, b As int2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' "x, y"
Operator int2d.cast () As String
  Return Str(x) & "," & Str(y)
End Operator

' a + b
Operator + (a As int2d, b As int2d) As int2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As int2d, b As int2d) As int2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As int2d) As int2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As int2d, b As int2d) As int2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As int2d, mul As Integer) As int2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a \ b
Operator \ (a As int2d, b As int2d) As int2d
   Return Type(a.x \ b.x, a.y \ b.y)
End Operator

' a \ div
Operator \ (a As int2d, div As Integer) As int2d
   Return Type(a.x \ div, a.y \ div)
End Operator

'===============================================================================

Type sgl2d
   As Single x, y
   Declare Constructor
   Declare Constructor(x As Single, y As Single)
   Declare Operator Cast () As String
End Type

Constructor sgl2d
End Constructor

Constructor sgl2d(x As Single, y As Single)
   This.x = x : This.y = y
End Constructor

' "x, y"
Operator sgl2d.cast () As String
   Return Str(x) & "," & Str(y)
End Operator

'---- operators ---

' distance / lenth
Operator Len (a As sgl2d) As Single
   Return Sqr(a.x * a.x + a.y * a.y)
End Operator

' a = b ?
Operator = (a As sgl2d, b As sgl2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

' a != b ?
Operator <> (a As sgl2d, b As sgl2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' a + b
Operator + (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As sgl2d) As sgl2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As sgl2d, mul As Single) As sgl2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a / div
Operator / (a As sgl2d, div As Single) As sgl2d
   Return Type(a.x / div, a.y / div)
End Operator

'---- extra functions ---

Function cross(a As sgl2d, b As sgl2d) As Single
   Return a.x * b.y - a.y * b.x
End Function

Function lengthSqrd(a As sgl2d) As Single
   Return (a.x * a.x) + (a.y * a.y)
End Function

Function dist(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return Sqr((dx * dx) + (dy * dy))
End Function

Function distSqrd(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return (dx * dx) + (dy * dy)
End Function

Function normalise(a As sgl2d) As sgl2d
   Dim As sgl2d temp
   Dim As Single length = Len(a)
   Return sgl2d(a.x / length, a.y / length)
End Function

'===============================================================================

'Note: y+ = up, x+ = right, (0,0) = center
Type scaled_graphics_type
   Dim As Single scale = 1 ' = 1 / pixel_size 'pixels / meter
   'dim as int2d offset' = (scrn_w \ 2, h \ 2) 'offset in pixels
   Dim As sgl2d offset
   Dim As Integer w = -1, h = -1
   Dim As Integer wc = -1, hc = -1 'center x,y
   Declare Sub setScreen(w As Integer, h As Integer)
   Declare Sub setScaling(scale As Single, offset As sgl2d)
   Declare Sub clearScreen(c As ULong)
   Declare Function pos2screen(p As sgl2d) As int2d
   Declare Sub drawPixel(p As sgl2d, c As ULong)
   Declare Sub drawCircle(p As sgl2d, r As Single, c As ULong)
   Declare Sub drawCircleFilled(p As sgl2d, r As Single, c As ULong, cFill As ULong)
   Declare Sub drawElipse(p As sgl2d, r As Single, aspect As Single, c As ULong)
   Declare Sub drawLine(p1 As sgl2d, p2 As sgl2d, c As ULong)
End Type

Sub scaled_graphics_type.setScreen(w As Integer, h As Integer)
   This.w = w 'width
   This.h = h 'height
   wc = w \ 2
   hc = h \ 2
   ScreenRes w, h, 32
   Width w \ 8, h \ 16 'bigger font
End Sub

Sub scaled_graphics_type.setScaling(scale As Single, offset As sgl2d)
   This.scale = scale
   This.offset = offset
End Sub

Sub scaled_graphics_type.clearScreen(c As ULong)
   Line(0, 0)-(w - 1, h - 1), c, bf
End Sub

Function scaled_graphics_type.pos2screen(p As sgl2d) As int2d
   Return int2d(Int(wc + (p.x - offset.x) * scale), h - Int(hc + (p.y - offset.y) * scale))
End Function

Sub scaled_graphics_type.drawPixel(p As sgl2d, c As ULong)
   Dim As int2d posScrn = pos2screen(p)
   PSet(posScrn.x, posScrn.y), c
End Sub

Sub scaled_graphics_type.drawCircle(p As sgl2d, r As Single, c As ULong)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawCircleFilled(p As sgl2d, r As Single, c As ULong, cFill As ULong)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, 0,,,,f
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawElipse(p As sgl2d, r As Single, aspect As Single, c As ULong)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c, , , aspect
End Sub

Sub scaled_graphics_type.drawLine(p1 As sgl2d, p2 As sgl2d, c As ULong)
   Dim As int2d posScrn1 = pos2screen(p1)
   Dim As int2d posScrn2 = pos2screen(p2)
   Line(posScrn1.x, posScrn1.y)-(posScrn2.x, posScrn2.y), c
End Sub

'===============================================================================

Const As Single PI = 4 * Atn(1)
Const As Single RAD_PER_DEG = (PI / 180)
Const As Single DEG_PER_RAD = 180 / PI
Const As Single sinA = Sin((10 / 180) * PI)
Const As Single cosA = Cos((10 / 180) * PI)
Const As Single sinB = Sin((20 / 180) * PI)
Const As Single cosB = Cos((20 / 180) * PI)

Const K_ESC = Chr(27)
Const K_MIN = Chr(45)
Const K_UND = Chr(95)
Const K_PLU = Chr(61)
Const K_EQU = Chr(43)

Const SCRN_W = 800, SCRN_H = 600

Dim Shared As scaled_graphics_type sg
sg.setScaling(2.0, sgl2d(0, 0))
sg.setScreen(SCRN_W, SCRN_H)

'-------------------------------------------------------------------------------

Function limit(value As Single, min As Single, max As Single) As Single
   If value < min Then Return min
   If value > max Then Return max
   Return value
End Function

Type polar
   Dim As Single angle
   Dim As Single magnitude
End Type

Function polarToCartesian(angle As Single, radius As Single) As sgl2d
   Return sgl2d(Cos(angle) * radius, Sin(angle) * radius)
End Function

Function rotatedVector(v As sgl2d, rotAngle As Single) As sgl2d
   Dim As sgl2d tmp
   tmp.x = Cos(rotAngle) * v.x - Sin(rotAngle) * v.y
   tmp.y = Sin(rotAngle) * v.x + Cos(rotAngle) * v.y
   Return tmp
End Function

'-------------------------------------------------------------------------------

Sub drawArrow(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.30 'reduce length
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x - sinB * dp.y, sinB * dp.x + cosB * dp.y), c)
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x + sinB * dp.y, cosB * dp.y - sinB * dp.x), c)
End Sub

Sub drawThruster(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.95 'reduce length
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x - sinA * dp.y, sinA * dp.x + cosA * dp.y), c)
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x + sinA * dp.y, cosA * dp.y - sinA * dp.x), c)
End Sub

Sub drawHelium3(p As sgl2d, r As Single, c As ULong)
   For i As Integer = 0 To 2
      sg.drawCircle(p + polarToCartesian((i / 3) * 2 * PI, 0.15 * r), 0.15 * r, c)
   Next
   sg.drawElipse(p, r, 2.0, c)
   sg.drawElipse(p, r, 0.5, c)
End Sub

Sub drawStar(p As sgl2d, size As Single, c As ULong)
   sg.drawLine(p - sgl2d(size / 2, 0) , p + sgl2d(size / 2, 0), c)
   sg.drawLine(p - sgl2d(0, size / 2) , p + sgl2d(0, size / 2), c)
End Sub

'-------------------------------------------------------------------------------

Type disc_object
   Dim As Single radius '[m]
   Dim As Single height '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   'linear motion properties
   Dim As sgl2d position 'position [m]
   Dim As Single lin_m 'mass [kg]
   Dim As sgl2d lin_F 'force [N] [kg*m/s^2]
   Dim As sgl2d lin_a 'acceleration [m/s^2]
   Dim As sgl2d lin_v 'velocity [m/s]
   'dim as sgl2d lin_p 'momentum [kg*m/s]
   'dim as single lin_E 'Kinetic energy [J] [kg*m^2/s^2]
   'Rotational motion properties
   Dim As Single angle 'angular position (theta) [rad]
   Dim As Single ang_m 'angular mass, moment of inertia (I) [kg*m^2]
   Dim As Single ang_F 'torque (tau) [N*m] [kg*m^2/s^2]
   Dim As Single ang_a 'angular velocity (alpha) [rad/s^2]
   Dim As Single ang_v 'angular velocity (omega) [rad/s]
   'dim as single ang_p 'angular momentum (L) [kg*m^2/s]
   'dim as single ang_E 'Kinetic energy [J] [kg*m^2/s^2]
   Declare Sub init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   Declare Sub update(dt As Double)
   Declare Function getKineticEnergy() As Single
End Type

'Set radius, height, density, position
'Calculate mass and rotational inertia
Sub disc_object.init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   height = h
   density = d
   position = p
   colour = c
   lin_m = PI * r ^ 2 * d
   ang_m = 0.5 * lin_m * r ^ 2
End Sub

'update position and angle
Sub disc_object.update(dt As Double)
   lin_a = lin_F / lin_m
   lin_v += lin_a * dt
   position += lin_v * dt
   ang_a = ang_F / ang_m
   ang_v += ang_a * dt
   angle += ang_v * dt
End Sub

Function disc_object.getKineticEnergy() As Single
   Dim As Single lin_E = 0.5 * lin_m * lengthSqrd(lin_v)
   Dim As Single ang_E = 0.5 * ang_m * ang_v * ang_v
   Return lin_E + ang_E
End Function

Sub drawShip(ship As disc_object)
   'calculate ships tail pointer / triangle
   Dim As sgl2d forwardPos = polarToCartesian(ship.angle - 90 * RAD_PER_DEG, ship.radius * 2.2)
   Dim As sgl2d leftBackPos = polarToCartesian(ship.angle - 135 * RAD_PER_DEG, ship.radius * 1.0) '(90 + 45)
   Dim As sgl2d rightBackPos = polarToCartesian(ship.angle - 45 * RAD_PER_DEG, ship.radius * 1.0) '(90 - 45)
   sg.drawCircle(ship.position, ship.radius, ship.colour) 'flying saucer
   sg.drawLine(ship.position + forwardPos, ship.position + leftBackPos, ship.colour)
   sg.drawLine(ship.position + forwardPos, ship.position + rightBackPos, ship.colour)
End Sub

'-------------------------------------------------------------------------------

Type thruster_type
   '''init paramaters
   Dim As polar polarForce '(rad, N)
   Dim As polar polarPos '(rad, m)
   '''variable paramaters
   Dim As sgl2d forceVector '(N, N)
   Dim As sgl2d relPos, absPos '(m, m)
   Dim As Integer active
   Declare Sub init(forceMagnitude As Single, forceDirection As Single, posAngle As Single, posRadius As Single)
   Declare Sub updatePosition(bodyPos As sgl2d, bodyAngle As Single)
End Type

Sub thruster_type.init(forceDirection As Single, forceMagnitude As Single, posAngle As Single, posRadius As Single)
   polarForce = Type(forceDirection, forceMagnitude) 'thruster action
   polarPos = Type(posAngle, posRadius) 'position of thruster on ship
End Sub

Sub thruster_type.updatePosition(bodyPos As sgl2d, bodyAngle As Single)
   relPos = polarToCartesian(bodyAngle + polarPos.angle, polarPos.magnitude)
   absPos = bodyPos + relPos
End Sub

'-------------------------------------------------------------------------------

Const As Single GRAV_CONST = 6.67e-11 '[m3/(kg*s^2)

Type astro_body
   Dim As Single radius '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   Dim As sgl2d position 'position [m]
   Dim As Single mass '[kg]
   Declare Sub init(r As Single, d As Single, p As sgl2d, c As ULong)
End Type

'Set radius, density, position
'Calculate mass and rotational inertia
Sub astro_body.init(r As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   density = d
   position = p
   colour = c
   mass = PI * r ^ 2 * d
End Sub

Function gravForce(m1 As Single, m2 As Single, r As Single) As Single
   Return GRAV_CONST * (m1 * m2) / (r * r)
End Function

Function gravForceVector(m1 As Single, Pos1 As sgl2d, m2 As Single, pos2 As sgl2d) As sgl2d
   Dim As Single distSquared = distSqrd(pos2, Pos1)
   Dim As sgl2d unitVector12 = (Pos1 - pos2) / Sqr(distSquared)
   Return unitVector12 * (-GRAV_CONST * (m1 * m2) / distSquared)
End Function

Function findNearestBody(refPos As sgl2d, body() As astro_body) As Integer
   Dim As Integer nearestBodyId = 0
   Dim As Single nearestBodyDistSqrd = distSqrd(refPos, body(0).position)
   Dim As Single currentBodyDistSqrd
   For i As Integer = 1 To UBound(body)
      currentBodyDistSqrd = distSqrd(refPos, body(i).position)
      If currentBodyDistSqrd < nearestBodyDistSqrd Then
         nearestBodyDistSqrd = currentBodyDistSqrd
         nearestBodyId = i
      End If
   Next
   Return nearestBodyId
End Function

Type star_type
   Dim As sgl2d position '[m]
   Dim As Single size
   Dim As ULong colour
   Declare Sub init(p As sgl2d, s As Single, c As ULong)
End Type

Sub star_type.init(p As sgl2d, s As Single, c As ULong)
   position = p
   size = s
   colour = c
End Sub

'-------------------------------------------------------------------------------

Const As Single MOON_RADIUS = 1737e+6 '[m]
Const As Single MOON_DENSITY = 3344 '[kg/m^3]

Const NUM_THRUSTERS = 6
Const L_FW_THR = 0 'left forward thruster
Const R_FW_THR = 1 'right forward thruster
Const L_LO_THR = 2
Const R_LO_THR = 3
Const L_HI_THR = 4
Const R_HI_THR = 5

Dim As String key
Dim As Integer quit = 0
Dim As disc_object ship
Dim As thruster_type thruster(NUM_THRUSTERS - 1)
Const NUM_ASTEROID = 80
Dim As astro_body asteroid(NUM_ASTEROID-1)
Const NUM_HELIUM = 10
ReDim As astro_body helium(NUM_HELIUM-1)
Const As Single maxFuel = 1e6 'N*s
Dim As Single fuel = maxFuel
Const NUM_STARS = 150
Dim As star_type star(NUM_STARS-1)

ship.init(10, 1, 5, sgl2d(0, -50), RGB(127, 223, 0))
For i As Integer = 0 To UBound(asteroid)
   asteroid(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(Rnd * 64 + 127, Rnd * 64 + 95, 127))
Next
For i As Integer = 0 To UBound(helium)
   helium(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(255, 191, 0))
Next
For i As Integer = 0 To UBound(star)
   star(i).init(sgl2d((Rnd - 0.5) * 4000, (Rnd - 0.5) * 4000), Rnd * 5 + 2, RGB(Rnd * 64 + 191, Rnd * 64 + 191, Rnd * 64 + 191))
Next

'force angle, force magnitude, polar thruster position
thruster(L_FW_THR).init(0.5 * pi, 1.2e4, -0.75 * pi, ship.radius)
thruster(R_FW_THR).init(0.5 * pi, 1.2e4, -0.25 * pi, ship.radius)
thruster(L_LO_THR).init(0.0 * pi, 8e3, -0.75 * pi, ship.radius)
thruster(R_LO_THR).init(1.0 * pi, 8e3, -0.25 * pi, ship.radius)
thruster(L_HI_THR).init(0.0 * pi, 8e3, +0.75 * pi, ship.radius)
thruster(R_HI_THR).init(1.0 * pi, 8e3, +0.25 * pi, ship.radius)

'intro text
Locate 10, 20: Print "HELIUM-3 SPACE RACE"
Locate 12, 20: Print "Scoop up all the helium-3 clouds as fast as possible."
Locate 13, 20: Print "Your ship needs helium-3 for its fusion powered trhusters."
Locate 14, 20: Print "Do not run out of it and do not collide with the asteroids."
Locate 15, 20: Print "Post your best time on the forum."
Locate 17, 20: Print "Press <ENTER> to start";
While InKey <> Chr(13)
   Sleep 1
Wend

Dim As Double tStart = Timer, tNow = tStart, tPrev = tNow, dt = 0
While quit = 0
   'reset stuff
   ship.lin_F = sgl2d(0, 0)
   ship.ang_F = 0
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).active = 0
   Next

   If MultiKey(FB.SC_UP) Then
      thruster(L_FW_THR).active = 1
      thruster(R_FW_THR).active = 1
      fuel -= thruster(L_FW_THR).polarForce.magnitude * dt
      fuel -= thruster(R_FW_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_LEFT) Then
      thruster(L_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
      fuel -= thruster(L_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(R_HI_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_RIGHT) Then
      thruster(R_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
      fuel -= thruster(R_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(L_HI_THR).polarForce.magnitude * dt
   End If

   If key = K_ESC Then quit = 1

   For i As Integer = 0 To NUM_THRUSTERS - 1
      'forces on body by active thrusters
      If thruster(i).active > 0 Then
         Dim As Single thrust = thruster(i).polarForce.magnitude
         thruster(i).forceVector = polarToCartesian(ship.angle + thruster(i).polarForce.angle, thrust)
         ship.lin_F += thruster(i).forceVector
         ship.ang_F += cross(thruster(i).relPos, thruster(i).forceVector)
      End If
   Next

   ship.update(dt) 'position and angle
   sg.offset = ship.position
   sg.scale = limit(100 / Len(ship.lin_v), 0.5, 2.0)

   'do always for display
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).updatePosition(ship.position, ship.angle)
   Next

   Dim As Integer nearestAsteroidId = findNearestBody(ship.position, asteroid())
   Dim As astro_body Ptr pAsteroid = @asteroid(nearestAsteroidId)
   If dist(ship.position, pAsteroid->position) < (ship.radius + pAsteroid->radius) Then
      quit = 2 'crash
   End If

   Dim As Integer nearestHeliumId = findNearestBody(ship.position, helium())
   
   If fuel <= 0 Then quit = 3
   
   'display
   ScreenLock
   sg.clearScreen(0)
   Locate 2, 2 : Print "<UP>, <LEFT>, <RIGHT> for thrusters, <ESC> to exit";
   Locate 3, 2 : Print "Nearest helium distance [m]: "; CInt(dist(ship.position, helium(nearestHeliumId).position));
   Locate 4, 2 : Print "Kinetic engergy [kJ]: "; CInt(ship.getKineticEnergy() * 1e-3);
   Locate 5, 2 : Print "Fuel remaining [%]: "; CInt((fuel / maxFuel) * 100);
   Locate 6, 2 : Print "Remaining helium3 clouds: "; UBound(helium) + 1;
   'draw stars background
   For i As Integer = 0 To UBound(star)
      drawStar(star(i).position * 0.75 + ship.position * 0.25, star(i).size, star(i).colour)
   Next
   drawShip(ship)
   'draw active thrusters
   For i As Integer = 0 To NUM_THRUSTERS - 1
      Dim As ULong c = IIf(i < 4, RGB(255, 255, 0), RGB(255, 255, 255))
      If thruster(i).active > 0 Then
         drawThruster(thruster(i).absPos, thruster(i).absPos - thruster(i).forceVector / 1e3, RGB(255, 63, 0)) 'thruster force indicator
      End If
   Next
   'draw asteroids
   For i As Integer = 0 To UBound(asteroid)
      sg.drawCircleFilled(asteroid(i).position, asteroid(i).radius, asteroid(i).colour, 0)
   Next
   'draw helium 'clouds'
   For i As Integer = 0 To UBound(helium)
      drawHelium3(helium(i).position, helium(i).radius, IIf(nearestHeliumId = i, RGB(191, 191, 255), helium(i).colour))
   Next
   Dim As sgl2d heliumPointer = normalise(ship.position - helium(nearestHeliumId).position)
   drawArrow(ship.position, ship.position - heliumPointer * ship.radius * 2, RGB(191, 191, 255))
   ScreenUnLock

   'clean up after displaying
   Dim As astro_body Ptr pHelium = @helium(nearestHeliumId)
   If dist(ship.position, pHelium->position) < (ship.radius + pHelium->radius * 0.5) Then
      If UBound(helium) > 0 Then
         'remove form list
         helium(nearestHeliumId) = helium(UBound(helium))
         ReDim Preserve helium(UBound(helium) - 1)
         fuel = maxFuel
      Else
         'last item
         quit = 4
      End If
   End If

   'time update
   key = InKey()
   Sleep 1
   tPrev = tNow
   tNow = Timer
   dt = tNow - tPrev
Wend

If quit > 1 Then 'Crashed or No fuel
   If quit = 2 Then ship.colour = RGB(233, 0, 0)
   If quit = 3 Then ship.colour = RGB(233, 191, 191)
   drawShip(ship)
   Locate 8, 2
   If quit = 2 Then Print "Ship crashed";
   If quit = 3 Then Print "Ship out of fuel";
   If quit = 4 Then Print "Well done all helium3 collected. Your time: " + Str(CInt(tNow - tStart)) + " seconds";
   Locate 9, 2: Print "press <ENTER> to exit";
   While InKey <> Chr(13)
      Sleep 1
   Wend
End If
Screen 0
Print "End"

'TODO:
'rotating helium3
'shoot asteroids
'background stars : optimize performance


Update 2019-09-04: Backgound stars added. New personal record: 244 seconds.

[Gunslinger]

Just got 189 seconds on second run full hd :-)
I'm a drone drone pilot and that is 3 axels, gravity and momentum to fight with.

[badidea]

Small update:
* One can shoot at a asteroid to remove it, costs 10% fuel
* Helium-3 atom clouds with rotating core animation
* Additional control of side thrusters <A>, <D> (does not make control easier)

Code: Select all

#Include "fbgfx.bi"

Type int2d
   As Integer x, y
   Declare Constructor
   Declare Constructor(x As Integer, y As Integer)
   Declare Operator Cast () As String
End Type

Constructor int2d
End Constructor

Constructor int2d(x As Integer, y As Integer)
   This.x = x : This.y = y
End Constructor

Operator = (a As int2d, b As int2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

Operator <> (a As int2d, b As int2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' "x, y"
Operator int2d.cast () As String
  Return Str(x) & "," & Str(y)
End Operator

' a + b
Operator + (a As int2d, b As int2d) As int2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As int2d, b As int2d) As int2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As int2d) As int2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As int2d, b As int2d) As int2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As int2d, mul As Integer) As int2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a \ b
Operator \ (a As int2d, b As int2d) As int2d
   Return Type(a.x \ b.x, a.y \ b.y)
End Operator

' a \ div
Operator \ (a As int2d, div As Integer) As int2d
   Return Type(a.x \ div, a.y \ div)
End Operator

'===============================================================================

Type sgl2d
   As Single x, y
   Declare Constructor
   Declare Constructor(x As Single, y As Single)
   Declare Operator Cast () As String
End Type

Constructor sgl2d
End Constructor

Constructor sgl2d(x As Single, y As Single)
   This.x = x : This.y = y
End Constructor

' "x, y"
Operator sgl2d.cast () As String
   Return Str(x) & "," & Str(y)
End Operator

'---- operators ---

' distance / lenth
Operator Len (a As sgl2d) As Single
   Return Sqr(a.x * a.x + a.y * a.y)
End Operator

' a = b ?
Operator = (a As sgl2d, b As sgl2d) As boolean
   If a.x <> b.x Then Return false
   If a.y <> b.y Then Return false
   Return true
End Operator

' a != b ?
Operator <> (a As sgl2d, b As sgl2d) As boolean
   If a.x = b.x And a.y = b.y Then Return false
   Return true
End Operator

' a + b
Operator + (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x + b.x, a.y + b.y)
End Operator

' a - b
Operator - (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x - b.x, a.y - b.y)
End Operator

' -a
Operator - (a As sgl2d) As sgl2d
   Return Type(-a.x, -a.y)
End Operator

' a * b
Operator * (a As sgl2d, b As sgl2d) As sgl2d
   Return Type(a.x * b.x, a.y * b.y)
End Operator

' a * mul
Operator * (a As sgl2d, mul As Single) As sgl2d
   Return Type(a.x * mul, a.y * mul)
End Operator

' a / div
Operator / (a As sgl2d, div As Single) As sgl2d
   Return Type(a.x / div, a.y / div)
End Operator

'---- extra functions ---

Function cross(a As sgl2d, b As sgl2d) As Single
   Return a.x * b.y - a.y * b.x
End Function

Function lengthSqrd(a As sgl2d) As Single
   Return (a.x * a.x) + (a.y * a.y)
End Function

Function dist(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return Sqr((dx * dx) + (dy * dy))
End Function

Function distSqrd(a As sgl2d, b As sgl2d) As Single
   Dim As Single dx = a.x - b.x
   Dim As Single dy = a.y - b.y
   Return (dx * dx) + (dy * dy)
End Function

Function normalise(a As sgl2d) As sgl2d
   Dim As sgl2d temp
   Dim As Single length = Len(a)
   Return sgl2d(a.x / length, a.y / length)
End Function

'===============================================================================

'Note: y+ = up, x+ = right, (0,0) = center
Type scaled_graphics_type
   Dim As Single scale = 1 ' = 1 / pixel_size 'pixels / meter
   'dim as int2d offset' = (scrn_w \ 2, h \ 2) 'offset in pixels
   Dim As sgl2d offset
   Dim As Integer w = -1, h = -1
   Dim As Integer wc = -1, hc = -1 'center x,y
   Declare Sub setScreen(w As Integer, h As Integer)
   Declare Sub setScaling(scale As Single, offset As sgl2d)
   Declare Sub clearScreen(c As ULong)
   Declare Function pos2screen(p As sgl2d) As int2d
   Declare Sub drawPixel(p As sgl2d, c As ULong)
   Declare Sub drawCircle(p As sgl2d, r As Single, c As ULong)
   Declare Sub drawCircleFilled(p As sgl2d, r As Single, c As ULong, cFill As ULong)
   Declare Sub drawElipse(p As sgl2d, r As Single, aspect As Single, c As ULong)
   Declare Sub drawLine(p1 As sgl2d, p2 As sgl2d, c As ULong)
End Type

Sub scaled_graphics_type.setScreen(w As Integer, h As Integer)
   This.w = w 'width
   This.h = h 'height
   wc = w \ 2
   hc = h \ 2
   ScreenRes w, h, 32
   Width w \ 8, h \ 16 'bigger font
End Sub

Sub scaled_graphics_type.setScaling(scale As Single, offset As sgl2d)
   This.scale = scale
   This.offset = offset
End Sub

Sub scaled_graphics_type.clearScreen(c As ULong)
   Line(0, 0)-(w - 1, h - 1), c, bf
End Sub

Function scaled_graphics_type.pos2screen(p As sgl2d) As int2d
   Return int2d(Int(wc + (p.x - offset.x) * scale), h - Int(hc + (p.y - offset.y) * scale))
End Function

Sub scaled_graphics_type.drawPixel(p As sgl2d, c As ULong)
   Dim As int2d posScrn = pos2screen(p)
   PSet(posScrn.x, posScrn.y), c
End Sub

Sub scaled_graphics_type.drawCircle(p As sgl2d, r As Single, c As ULong)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawCircleFilled(p As sgl2d, r As Single, c As ULong, cFill As ULong)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, 0,,,,f
   Circle(posScrn.x, posScrn.y), r * scale, c
End Sub

Sub scaled_graphics_type.drawElipse(p As sgl2d, r As Single, aspect As Single, c As ULong)
   Dim As int2d posScrn = pos2screen(p)
   Circle(posScrn.x, posScrn.y), r * scale, c, , , aspect
End Sub

Sub scaled_graphics_type.drawLine(p1 As sgl2d, p2 As sgl2d, c As ULong)
   Dim As int2d posScrn1 = pos2screen(p1)
   Dim As int2d posScrn2 = pos2screen(p2)
   Line(posScrn1.x, posScrn1.y)-(posScrn2.x, posScrn2.y), c
End Sub

'===============================================================================

Const As Single PI = 4 * Atn(1)
Const As Single RAD_PER_DEG = (PI / 180)
Const As Single DEG_PER_RAD = 180 / PI
Const As Single sinA = Sin((10 / 180) * PI)
Const As Single cosA = Cos((10 / 180) * PI)
Const As Single sinB = Sin((20 / 180) * PI)
Const As Single cosB = Cos((20 / 180) * PI)

Const K_ENTER = Chr(13)
Const K_ESC = Chr(27)
Const K_MIN = Chr(45)
Const K_UND = Chr(95)
Const K_PLU = Chr(61)
Const K_EQU = Chr(43)

Const SCRN_W = 800, SCRN_H = 600

Dim Shared As scaled_graphics_type sg
sg.setScaling(2.0, sgl2d(0, 0))
sg.setScreen(SCRN_W, SCRN_H)

'-------------------------------------------------------------------------------

Sub waitForKey(key As String)
   While InKey <> key
      Sleep 1
   Wend
End Sub

Function limit(value As Single, min As Single, max As Single) As Single
   If value < min Then Return min
   If value > max Then Return max
   Return value
End Function

Type polar
   Dim As Single angle
   Dim As Single magnitude
End Type

Function polarToCartesian(angle As Single, radius As Single) As sgl2d
   Return sgl2d(Cos(angle) * radius, Sin(angle) * radius)
End Function

Function rotatedVector(v As sgl2d, rotAngle As Single) As sgl2d
   Dim As sgl2d tmp
   tmp.x = Cos(rotAngle) * v.x - Sin(rotAngle) * v.y
   tmp.y = Sin(rotAngle) * v.x + Cos(rotAngle) * v.y
   Return tmp
End Function

'-------------------------------------------------------------------------------

Sub drawArrow(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.30 'reduce length
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x - sinB * dp.y, sinB * dp.x + cosB * dp.y), c)
   sg.drawLine(p2, p2 - sgl2d(cosB * dp.x + sinB * dp.y, cosB * dp.y - sinB * dp.x), c)
End Sub

Sub drawThruster(p1 As sgl2d, p2 As sgl2d, c As ULong)
   sg.drawLine(p1, p2, c)
   Dim As sgl2d dp = (p2 - p1) * 0.95 'reduce length
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x - sinA * dp.y, sinA * dp.x + cosA * dp.y), c)
   sg.drawLine(p1, p1 + sgl2d(cosA * dp.x + sinA * dp.y, cosA * dp.y - sinA * dp.x), c)
End Sub

Sub drawHelium3(p As sgl2d, r As Single, rot As Single, c As ULong)
   For i As Integer = 0 To 2
      sg.drawCircle(p + polarToCartesian((i / 3) * 2 * PI + rot, 0.15 * r), 0.15 * r, c)
   Next
   sg.drawElipse(p, r, 2.0, c)
   sg.drawElipse(p, r, 0.5, c)
End Sub

Sub drawStar(p As sgl2d, size As Single, c As ULong)
   sg.drawLine(p - sgl2d(size / 2, 0) , p + sgl2d(size / 2, 0), c)
   sg.drawLine(p - sgl2d(0, size / 2) , p + sgl2d(0, size / 2), c)
End Sub

'-------------------------------------------------------------------------------

Type disc_object
   Dim As Single radius '[m]
   Dim As Single height '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   'linear motion properties
   Dim As sgl2d position 'position [m]
   Dim As Single lin_m 'mass [kg]
   Dim As sgl2d lin_F 'force [N] [kg*m/s^2]
   Dim As sgl2d lin_a 'acceleration [m/s^2]
   Dim As sgl2d lin_v 'velocity [m/s]
   'dim as sgl2d lin_p 'momentum [kg*m/s]
   'dim as single lin_E 'Kinetic energy [J] [kg*m^2/s^2]
   'Rotational motion properties
   Dim As Single angle 'angular position (theta) [rad]
   Dim As Single ang_m 'angular mass, moment of inertia (I) [kg*m^2]
   Dim As Single ang_F 'torque (tau) [N*m] [kg*m^2/s^2]
   Dim As Single ang_a 'angular velocity (alpha) [rad/s^2]
   Dim As Single ang_v 'angular velocity (omega) [rad/s]
   'dim as single ang_p 'angular momentum (L) [kg*m^2/s]
   'dim as single ang_E 'Kinetic energy [J] [kg*m^2/s^2]
   Declare Sub init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   Declare Sub update(dt As Double)
   Declare Function getKineticEnergy() As Single
End Type

'Set radius, height, density, position
'Calculate mass and rotational inertia
Sub disc_object.init(r As Single, h As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   height = h
   density = d
   position = p
   colour = c
   lin_m = PI * r ^ 2 * d
   ang_m = 0.5 * lin_m * r ^ 2
End Sub

'update position and angle
Sub disc_object.update(dt As Double)
   lin_a = lin_F / lin_m
   lin_v += lin_a * dt
   position += lin_v * dt
   ang_a = ang_F / ang_m
   ang_v += ang_a * dt
   angle += ang_v * dt
End Sub

Function disc_object.getKineticEnergy() As Single
   Dim As Single lin_E = 0.5 * lin_m * lengthSqrd(lin_v)
   Dim As Single ang_E = 0.5 * ang_m * ang_v * ang_v
   Return lin_E + ang_E
End Function

Sub drawShip(ship As disc_object)
   'calculate ships tail pointer / triangle
   Dim As sgl2d forwardPos = polarToCartesian(ship.angle - 90 * RAD_PER_DEG, ship.radius * 2.2)
   Dim As sgl2d leftBackPos = polarToCartesian(ship.angle - 135 * RAD_PER_DEG, ship.radius * 1.0) '(90 + 45)
   Dim As sgl2d rightBackPos = polarToCartesian(ship.angle - 45 * RAD_PER_DEG, ship.radius * 1.0) '(90 - 45)
   sg.drawCircle(ship.position, ship.radius, ship.colour) 'flying saucer
   sg.drawLine(ship.position + forwardPos, ship.position + leftBackPos, ship.colour)
   sg.drawLine(ship.position + forwardPos, ship.position + rightBackPos, ship.colour)
End Sub

'-------------------------------------------------------------------------------

Type thruster_type
   '''init paramaters
   Dim As polar polarForce '(rad, N)
   Dim As polar polarPos '(rad, m)
   '''variable paramaters
   Dim As sgl2d forceVector '(N, N)
   Dim As sgl2d relPos, absPos '(m, m)
   Dim As Integer active
   Declare Sub init(forceMagnitude As Single, forceDirection As Single, posAngle As Single, posRadius As Single)
   Declare Sub updatePosition(bodyPos As sgl2d, bodyAngle As Single)
End Type

Sub thruster_type.init(forceDirection As Single, forceMagnitude As Single, posAngle As Single, posRadius As Single)
   polarForce = Type(forceDirection, forceMagnitude) 'thruster action
   polarPos = Type(posAngle, posRadius) 'position of thruster on ship
End Sub

Sub thruster_type.updatePosition(bodyPos As sgl2d, bodyAngle As Single)
   relPos = polarToCartesian(bodyAngle + polarPos.angle, polarPos.magnitude)
   absPos = bodyPos + relPos
End Sub

'-------------------------------------------------------------------------------

Const As Single GRAV_CONST = 6.67e-11 '[m3/(kg*s^2)

Type astro_body
   Dim As Single radius '[m]
   Dim As Single density '[kg/m^3]
   Dim As ULong colour '[m]
   Dim As sgl2d position 'position [m]
   Dim As Single mass '[kg]
   Dim As Single rotation '[rad]
   Declare Sub init(r As Single, d As Single, p As sgl2d, c As ULong)
End Type

'Set radius, density, position
'Calculate mass and rotational inertia
Sub astro_body.init(r As Single, d As Single, p As sgl2d, c As ULong)
   radius = r
   density = d
   position = p
   colour = c
   mass = PI * r ^ 2 * d
End Sub

Function gravForce(m1 As Single, m2 As Single, r As Single) As Single
   Return GRAV_CONST * (m1 * m2) / (r * r)
End Function

Function gravForceVector(m1 As Single, Pos1 As sgl2d, m2 As Single, pos2 As sgl2d) As sgl2d
   Dim As Single distSquared = distSqrd(pos2, Pos1)
   Dim As sgl2d unitVector12 = (Pos1 - pos2) / Sqr(distSquared)
   Return unitVector12 * (-GRAV_CONST * (m1 * m2) / distSquared)
End Function

Function findNearestBody(refPos As sgl2d, body() As astro_body) As Integer
   Dim As Integer nearestBodyId = 0
   Dim As Single nearestBodyDistSqrd = distSqrd(refPos, body(0).position)
   Dim As Single currentBodyDistSqrd
   For i As Integer = 1 To UBound(body)
      currentBodyDistSqrd = distSqrd(refPos, body(i).position)
      If currentBodyDistSqrd < nearestBodyDistSqrd Then
         nearestBodyDistSqrd = currentBodyDistSqrd
         nearestBodyId = i
      End If
   Next
   Return nearestBodyId
End Function

Type star_type
   Dim As sgl2d position '[m]
   Dim As Single size
   Dim As ULong colour
   Declare Sub init(p As sgl2d, s As Single, c As ULong)
End Type

Sub star_type.init(p As sgl2d, s As Single, c As ULong)
   position = p
   size = s
   colour = c
End Sub

Type bullet_type
   Dim As sgl2d position
   Dim As sgl2d velocity
   Dim As Single radius = 3.0
   Dim As Integer active
   Dim As Double endTime
   Dim As ULong colour
   Declare Sub init(p As sgl2d, v As sgl2d, lifeTime As Double, c As ULong)
End Type

Sub bullet_type.init(p As sgl2d, v As sgl2d, lifeTime As Double, c As ULong)
   position = p
   velocity = v
   endTime = Timer + lifeTime
   active = 1
   colour = c
End Sub

'-------------------------------------------------------------------------------

Const As Single MOON_RADIUS = 1737e+6 '[m]
Const As Single MOON_DENSITY = 3344 '[kg/m^3]

Const NUM_THRUSTERS = 6
Const L_FW_THR = 0 'left forward thruster
Const R_FW_THR = 1 'right forward thruster
Const L_LO_THR = 2
Const R_LO_THR = 3
Const L_HI_THR = 4
Const R_HI_THR = 5

Const QUIT_NO = 0
Const QUIT_USER = 1
Const QUIT_CRASH = 2
Const QUIT_FUEL = 3
Const QUIT_WINNER = 4

Dim As Integer quit = 0
Dim As disc_object ship
Dim As thruster_type thruster(NUM_THRUSTERS - 1)
'const NUM_ASTEROID = 80
ReDim As astro_body asteroid(80 - 1)
'const NUM_HELIUM = 10
ReDim As astro_body helium(10 - 1)
Const As Single maxFuel = 1e6 'N*s
Dim As Single fuel = maxFuel
Const NUM_STARS = 150
Dim As star_type star(NUM_STARS-1)
Dim As bullet_type bullet
Dim As astro_body Ptr pAsteroid, pHelium
Dim As Integer nearestAsteroidToBulletId = -1
Dim As Integer nearestHeliumId
Dim As Integer nearestAsteroidId
Dim As Integer asteroidRemove = 0

ship.init(10, 1, 5, sgl2d(0, -50), RGB(127, 223, 0))
For i As Integer = 0 To UBound(asteroid)
   asteroid(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(Rnd * 64 + 127, Rnd * 64 + 95, 127))
Next
For i As Integer = 0 To UBound(helium)
   helium(i).init(5 + 4 / (Rnd + 0.2), 1000, sgl2d((Rnd - 0.5) * 2000, (Rnd - 0.5) * 2000), RGB(255, 191, 0))
Next
For i As Integer = 0 To UBound(star)
   star(i).init(sgl2d((Rnd - 0.5) * 4000, (Rnd - 0.5) * 4000), Rnd * 5 + 2, RGB(Rnd * 64 + 191, Rnd * 64 + 191, Rnd * 64 + 191))
Next

'force angle, force magnitude, polar thruster position
thruster(L_FW_THR).init(0.5 * pi, 1.2e4, -0.75 * pi, ship.radius)
thruster(R_FW_THR).init(0.5 * pi, 1.2e4, -0.25 * pi, ship.radius)
thruster(L_LO_THR).init(0.0 * pi, 8e3, -0.75 * pi, ship.radius)
thruster(R_LO_THR).init(1.0 * pi, 8e3, -0.25 * pi, ship.radius)
thruster(L_HI_THR).init(0.0 * pi, 8e3, +0.75 * pi, ship.radius)
thruster(R_HI_THR).init(1.0 * pi, 8e3, +0.25 * pi, ship.radius)

'intro text
Locate 10, 20: Print "HELIUM-3 SPACE RACE"
Locate 12, 20: Print "Scoop up all the helium-3 clouds as fast as possible."
Locate 13, 20: Print "Your ship needs helium-3 for its fusion powered trhusters."
Locate 14, 20: Print "Do not run out of it and do not collide with the asteroids."
Locate 15, 20: Print "Post your best time on the forum."
Locate 17, 20: Print "Press <ENTER> to start";

Locate 24, 20: Print "Controls:";
Locate 26, 20: Print " Fire: <SPACE>";
Locate 27, 20: Print " Forward thrusters: <UP>";
Locate 28, 20: Print " Rotation thrusters: <LEFT>, <RIGHT>";
Locate 29, 20: Print " Side thrusters: <A>, <D>";
Locate 30, 20: Print " Abort game: <ESC>";
waitForKey(K_ENTER)

Dim As Double tStart = Timer, tNow = tStart, tPrev = tNow, dt = 0
While quit = 0
   'reset stuff
   ship.lin_F = sgl2d(0, 0)
   ship.ang_F = 0
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).active = 0
   Next

   If MultiKey(FB.SC_UP) Then
      thruster(L_FW_THR).active = 1
      thruster(R_FW_THR).active = 1
      fuel -= thruster(L_FW_THR).polarForce.magnitude * dt
      fuel -= thruster(R_FW_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_LEFT) Then
      thruster(L_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
      fuel -= thruster(L_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(R_HI_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_RIGHT) Then
      thruster(R_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
      fuel -= thruster(R_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(L_HI_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_A) Then
      thruster(R_LO_THR).active = 1
      thruster(R_HI_THR).active = 1
      fuel -= thruster(R_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(R_HI_THR).polarForce.magnitude * dt
   End If

   If MultiKey(FB.SC_D) Then
      thruster(L_LO_THR).active = 1
      thruster(L_HI_THR).active = 1
      fuel -= thruster(L_LO_THR).polarForce.magnitude * dt
      fuel -= thruster(L_HI_THR).polarForce.magnitude * dt
   End If

   '~ if multikey(FB.SC_P) then
      '~ sleep 5000,1 'for taking a screenshot
   '~ end if

   If MultiKey(FB.SC_SPACE) Then
      If bullet.active = 0 Then
         Dim As sgl2d vBullet = ship.lin_v + polarToCartesian(ship.angle + pi/2, 50)
         bullet.init(ship.position, vBullet, 3.0, RGB(255, 0, 0))
         fuel -= maxFuel * 0.05 'takes 5% of fuel
      End If
   End If

   If MultiKey(FB.SC_ESCAPE) Then quit = QUIT_USER

   For i As Integer = 0 To NUM_THRUSTERS - 1
      'forces on body by active thrusters
      If thruster(i).active > 0 Then
         Dim As Single thrust = thruster(i).polarForce.magnitude
         thruster(i).forceVector = polarToCartesian(ship.angle + thruster(i).polarForce.angle, thrust)
         ship.lin_F += thruster(i).forceVector
         ship.ang_F += cross(thruster(i).relPos, thruster(i).forceVector)
      End If
   Next

   ship.update(dt) 'position and angle
   sg.offset = ship.position
   sg.scale = limit(100 / Len(ship.lin_v), 0.5, 2.0)

   'do always for display
   For i As Integer = 0 To NUM_THRUSTERS - 1
      thruster(i).updatePosition(ship.position, ship.angle)
   Next

   For i As Integer = 0 To UBound(helium)
      helium(i).rotation += 300 * RAD_PER_DEG * dt 'N * degrees / second
   Next

   nearestAsteroidId = findNearestBody(ship.position, asteroid())
   pAsteroid = @asteroid(nearestAsteroidId)
   If dist(ship.position, pAsteroid->position) < (ship.radius + pAsteroid->radius) Then
      quit = QUIT_CRASH
   End If

   nearestHeliumId = findNearestBody(ship.position, helium())
   
   If fuel <= 0 Then quit = QUIT_FUEL
   
   If bullet.active = 1 Then
      If tNow > bullet.EndTime Then
         bullet.active = 0
      Else
         bullet.position += bullet.velocity * dt
         nearestAsteroidToBulletId = findNearestBody(bullet.position, asteroid())
         asteroidRemove = 1
      End If
   End If
   
   'display
   ScreenLock
   sg.clearScreen(0)
   'draw stars background
   For i As Integer = 0 To UBound(star)
      drawStar(star(i).position * 0.75 + ship.position * 0.25, star(i).size, star(i).colour)
   Next
   drawShip(ship)
   'draw active thrusters
   For i As Integer = 0 To NUM_THRUSTERS - 1
      Dim As ULong c = IIf(i < 4, RGB(255, 255, 0), RGB(255, 255, 255))
      If thruster(i).active > 0 Then
         drawThruster(thruster(i).absPos, thruster(i).absPos - thruster(i).forceVector / 1e3, RGB(255, 63, 0)) 'thruster force indicator
      End If
   Next
   'draw asteroids
   For i As Integer = 0 To UBound(asteroid)
      sg.drawCircleFilled(asteroid(i).position, asteroid(i).radius, asteroid(i).colour, 0)
   Next
   'draw helium 'clouds'
   For i As Integer = 0 To UBound(helium)
      drawHelium3(helium(i).position, helium(i).radius, helium(i).rotation, IIf(nearestHeliumId = i, RGB(191, 191, 255), helium(i).colour))
   Next
   Dim As sgl2d heliumPointer = normalise(ship.position - helium(nearestHeliumId).position)
   drawArrow(ship.position, ship.position - heliumPointer * ship.radius * 2, RGB(191, 191, 255))
   If bullet.active = 1 Then
      sg.drawCircle(bullet.position, bullet.radius, bullet.colour)
   End If
   Draw String (8, 8 + 0 * 16), "Nearest helium distance [m]: " + Str(CInt(dist(ship.position, helium(nearestHeliumId).position)))
   Draw String (8, 8 + 1 * 16), "Kinetic engergy [kJ]: " + Str(CInt(ship.getKineticEnergy() * 1e-3))
   Draw String (8, 8 + 2 * 16), "Fuel remaining [%]: " + Str(CInt((fuel / maxFuel) * 100))
   Draw String (8, 8 + 3 * 16), "Remaining helium3 clouds: " + Str(UBound(helium) + 1)
   ScreenUnLock

   'clean up after displaying
   pHelium = @helium(nearestHeliumId)
   If dist(ship.position, pHelium->position) < (ship.radius + pHelium->radius * 0.5) Then
      If UBound(helium) > 0 Then
         'remove form list
         helium(nearestHeliumId) = helium(UBound(helium))
         ReDim Preserve helium(UBound(helium) - 1)
         fuel = maxFuel
      Else
         'last item
         Erase helium
         quit = QUIT_WINNER
      End If
   End If

   'asteroid hit by bullet?
   If nearestAsteroidToBulletId <> -1 Then
      pAsteroid = @asteroid(nearestAsteroidToBulletId)
      If dist(bullet.position, pAsteroid->position) < (bullet.radius + pAsteroid->radius) Then
         'remove asteroid from list
         If UBound(asteroid) > 0 Then
            'remove form list
            asteroid(nearestAsteroidToBulletId) = asteroid(UBound(asteroid))
            ReDim Preserve asteroid(UBound(asteroid) - 1)
         Else
            'last item
            Erase asteroid
         End If
         bullet.endTime = tNow + 0.1
         'bullet.active = 0
      End If
      nearestAsteroidToBulletId = -1
   End If

   'time update
   Sleep 1
   tPrev = tNow
   tNow = Timer
   dt = tNow - tPrev
Wend

Select Case quit
Case QUIT_USER
   Draw String (8, 8 + 6 * 16), "Abort by user"
Case QUIT_CRASH
   ship.colour = RGB(233, 0, 0)
   drawShip(ship)
   Draw String (8, 8 + 6 * 16), "Ship crashed"
Case QUIT_FUEL
   ship.colour = RGB(233, 191, 191)
   drawShip(ship)
   Draw String (8, 8 + 6 * 16), "Ship out of fuel"
Case QUIT_WINNER
   Draw String (8, 8 + 6 * 16), "Well done all helium3 collected. Your time: " + Str(CInt(tNow - tStart)) + " seconds"
Case Else
End Select

Draw String (8, 8 + 7 * 16), "press <ENTER> to exit"
waitForKey(K_ENTER)

Screen 0
Print "End"


Todo:
* Code clean-up, move stuff into classes
* Add some simple sounds (using: Audio library for FreeBasic - Features)

[dodicat]

Windows for fun with some sound.
(Some help from angros47's post regarding playing .mp3 files on win 10)
https://www.mediafire.com/file/euzi8vn0upn6bsr/RIP.zip/file