'******************************************************************************'
'
'   Michael "h4tt3n" Nissen's FreeBasic mass-spring-damper tutorial, feb 2008
'
'   Example # 2.2
'
'   Simple two dimensional spring simulation. This program simulates a number
'   of elasic balls which jump around and interact. (Press esc to quit.)
'
'******************************************************************************'

''  set constants. experiment with these and see how the simulation reacts.
Const Pi                =     4*Atn(1)      ''  pi (better not change ^^)
Const Timestep          =     0.01          ''  timestep (delta time)
Const Num_Balls         =     12         -1 ''  number of elastic balls
Const Mass_Min          =     400           ''  smallest elastic ball mass
Const Mass_Max          =     8000         ''  biggest elastic ball mass
Const Density           =     0.01          ''  elastic ball density 
Const Spring_Stiffnes   =     5e5           ''  elastic stiffnes

''  define types
Type Vector_Type
  As Single X, Y
End Type

Type Ball_Type
  As Uinteger Col
  As Single Mass, Density, Radius
  As Vector_Type Frc, Acc, Vel, Pos
End Type 

''  dimension variables
Dim As Vector_Type Dst, Vel
Dim As Ball_Type Ball(Num_Balls)
Dim As Single Force, Ball_Distance, Cos_Angle, Sin_Angle, Dist_Min
Dim As Integer i, i2, Scrn_Wid, Scrn_Hgt

''  set screen width, height, and bit depth
Scrn_Wid = 800
Scrn_Hgt = 600
ScreenRes Scrn_Wid, Scrn_hgt, 16

''  use timer to generate random numbers
Randomize Timer

''  set startup conditions of elastic balls
For i = Lbound(Ball) To Ubound(Ball)
  
  With Ball(i)
    
    .Col = RGB(64+Rnd*192, 64+Rnd*192, 64+Rnd*192)
    .Mass = Mass_Min+(i/Ubound(Ball))*(Mass_Max-Mass_Min)
    .Density = Density
    .Radius = ((.Mass/.Density)/((4/3)*pi))^(1/3)
    .Pos.x = .Radius+Rnd*(Scrn_Wid-.Radius)
    .pos.y = .Radius+Rnd*(scrn_hgt-.Radius)
    .vel.x = (Rnd-Rnd)*200
    .vel.y = (Rnd-Rnd)*200
    
  End With
  
Next

''  main program loop
Do
  
  ''  test all elastic balls against each other. calculate forces if they touch.
  For i = Lbound(Ball) To Ubound(Ball)-1
    
    For i2 = i+1 To Ubound(Ball)
      
      ''  distance between elastic balls (Pythagoras' theorem)
      Dst.X = Ball(i2).Pos.X-Ball(i).Pos.X
      Dst.Y = Ball(i2).Pos.Y-Ball(i).Pos.Y
      Ball_Distance = Sqr(Dst.X*Dst.X+Dst.Y*Dst.Y)
      Dist_Min = Ball(i).Radius+Ball(i2).Radius
      
      If Ball_Distance < Dist_Min Then
        
        ''  cosine and sine to the angle between Ball i and i2 (trigonometry)
        Cos_Angle = Dst.X/Ball_Distance
        Sin_Angle = Dst.Y/Ball_Distance
        
        ''  spring force (Hooke's law of elasticity)
        Force = -Spring_Stiffnes*(Ball_Distance-Dist_Min)
        
        ''  split spring force vector into horizontal and vertical component
        Ball(i).Frc.X -= Force*Cos_Angle
        Ball(i).Frc.Y -= Force*Sin_Angle
        Ball(i2).Frc.X += Force*Cos_Angle
        Ball(i2).Frc.Y += Force*Sin_Angle
        
      End If
      
    Next
    
  Next
  
  ''  update acceleration, velocity, and position of elastic balls
  ''  (using the Euler-Cromer 1st order integration algorithm)
  For i = Lbound(Ball) To Ubound(Ball)
    
    With Ball(i)
      
      ''  accelerate balls (acceleration = force / mass)
      .Acc.X = .Frc.X/.Mass
      .Acc.Y = .Frc.Y/.Mass
      
      ''  reset force vector
      .frc.x = 0
      .frc.y = 0
      
      ''  update velocity
      ''  delta velocity = acceleration * delta time
      ''  new velocity = old velocity + delta velocity
      .Vel.X += .Acc.X*Timestep
      .Vel.Y += .Acc.Y*Timestep
      
      ''  update position 
      ''  delta position = velocity * delta time
      ''  new position = old position + delta position
      .Pos.X += .Vel.X*Timestep
      .Pos.Y += .Vel.Y*Timestep
      
    End With
    
  Next
  
  ''  keep elastic balls within screen boundaries
  For i = Lbound(Ball) To Ubound(Ball)
    
    With Ball(i)
      
      ''  right
      If .Pos.x > Scrn_Wid-.Radius Then
        .Vel.X = -.vel.x
        .pos.x = Scrn_Wid-.Radius
      End If
      ''  left
      If .Pos.x < .Radius Then
        .Vel.X = -.vel.x
        .pos.x = .Radius
      End If
      ''  bottom
      If .Pos.Y > Scrn_Hgt-.Radius Then
        .vel.y = -.vel.y
        .pos.y = Scrn_Hgt-.Radius
      End If
      ''  top
      If .Pos.Y < .Radius Then
        .vel.y = -.vel.y
        .pos.y = .Radius
      End If
      
    End With
    
  Next
  
  ''  display elastic balls in graphics
  ScreenLock
    
    ''  clear screen
    Cls
    
    ''  draw elastic balls
    For i = Lbound(Ball) To Ubound(Ball)
      With Ball(i)
        Circle(.Pos.X, .Pos.Y), .Radius    , .Col,,, 1
        Circle(.Pos.X, .Pos.Y), .Radius-0.5, .Col,,, 1
        Circle(.Pos.X, .Pos.Y), .Radius-1.0, .Col,,, 1
      End With
    Next i
    
  ScreenUnlock
  
  Sleep 1, 1
  
Loop Until Multikey(1)

End