rem Bouncy Balls ......... Rev 4.4
      rem A J Tooth // 30th November 2003
      rem Modified for 15 balls March 2006
      on error if (err=17) then quit
      
      rem =============================================
      rem Simulates bouncing balls, with no friction.
      rem They bounce off the walls and off each other.
      rem =============================================
      
      rem Setup
      proc_setup
      *REFRESH OFF
      
      rem Initialisation
      proc_init
      *REFRESH
      
      rem Main Routine
      repeat
        proc_main
        
        a$="" : a$=inkey$(1)
        rem Kick a ball, if relevant key pressed
        proc_kick(a$)
        
      until (a$=" ")
      
      *REFRESH ON
      quit
      end
      rem END OF PROGRAM ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
      rem =============================================================================
      
      rem Setup
      def proc_setup
      
      *FLOAT64
      rem Initialise screen
      proc_fullscreen(xscreen%,yscreen%)
      xmax%=2*xscreen% : ymax%=2*yscreen%
      
      *FONT Desdemona,24
      colour 1: print tab(10,5);" The BOUNCING BALLS"
      *FONT Arial,10,B
      colour 2
      print tab(20,20);"At any time press -a- to 'kick' the GREEN ball, -s- to kick the RED ball, -d- for the YELLOW."
      print tab(20,21);"also -f- for the BLUE ball, and so on."
      colour 3:print tab(20,23);"Pressing the -#- key kicks ALL balls simultaneously."
      colour 2:print tab(20,33);"A Coefficient of Restitution below 1 is normal."
      print tab(20,34);"Make it equal to 1, and the balls will go on bouncing forever."
      print tab(20,35);"Greater than 1 and they will get ever faster (the program will hit errors and stop eventually)."
      print tab(20,36);"Zero and below are not selectable."
      print tab(20,38);"Enter the Coefficient or Restitution (0.9) "
      colour 1:input tab(65,38);res
      colour 2:print tab(20,40);"Enter the number of BALLS (2 to 15) "
      colour 1
      repeat
        input tab(60,40);ball&
      until (ball&>1 and ball&<16)
      
      if res<=0.0 then res=0.9
      
      rem Set key sequence for kicking balls
      ref$="asdfghjklzxcvbn#"
      
      a$=inkey$(50):cls
      endproc
      rem =============================================================================
      
      rem Initialisation
      def proc_init
      local i&
      
      g=9.81 : rem Gravitational acceleration
      
      rem Choose initial parameters
      dim x(ball&),y(ball&),xn(ball&),yn(ball&)
      for i&=1 to ball&
        x(i&)=1.0*(50+rnd(xmax%-50)) : y(i&)=1.0*(50+rnd(xmax%-50))
      next i&
      
      r%=40 : rem Radius of balls
      
      dim u(ball&),v(ball&)
      for i&=1 to ball&
        u(i&)=1.0*(200+rnd(100)) : v(i&)=1.0*(200+rnd(100))
      next i&
      
      stps&=20
      dtim=1/stps&
      
      dim tcum(ball&),ch%(ball&)
      tcum()=0.0
      ch%()=0
      
      rem Plot initial position of ball
      for i&=1 to ball&
        proc_plot(x(i&),y(i&),i&)
      next i&
      
      rem Keep old values for un-plotting later
      dim xold(ball&),yold(ball&)
      for i&=1 to ball&
        xold(i&)=x(i&) : yold(i&)=y(i&)
      next i&
      
      rem Record whether the ball has hit another one
      dim ht&(ball&)
      ht&()=0
      
      endproc
      rem =============================================================================
      
      rem Main routine
      def proc_main
      local i&,j&,ml
      
      rem Increment the time
      for i&=0 to ball&
        tcum(i&)+=dtim
      next i&
      
      rem Equations of Motion
      for i&=1 to ball&
        rem Boost factor after a hit to reduce instances
        rem of balls "sticking" to each other
        if ht&(i&)=1 then
          ml=2.0
          ht&(i&)=0
        else
          ml=1.0
        endif
        
        xn(i&)=x(i&) + ml*(u(i&)*tcum(i&))
        yn(i&)=y(i&) + ml*(v(i&)*tcum(i&) - (g*tcum(i&)*tcum(i&)/2))
      next i&
      
      rem Plot new position of balls and unplot previous position
      for i&=1 to ball&
        proc_plot(xold(i&),yold(i&),0)
        proc_plot(xn(i&),yn(i&),i&)
      next i&
      *REFRESH
      
      rem Reassign old positions
      for i&=1 to ball&
        xold(i&)=xn(i&) : yold(i&)=yn(i&)
      next i&
      
      rem Check if a ball has hit a wall
      for i&=1 to ball&
        proc_poscheck(xn(i&),yn(i&),i&)
      next i&
      
      rem Check if the balls have hit each other
      for i&=1 to (ball&-1)
        for j&=(i&+1) to ball&
          R=sqr((xn(j&)-xn(i&))^2 + (yn(j&)-yn(i&))^2)
          if R<2.0*r% then proc_hit(xn(i&),yn(i&),xn(j&),yn(j&),i&,j&)
        next j&
      next i&
      endproc
      rem =============================================================================
      
      rem Plot / Unplot Procedure
      def proc_plot(X,Y,col&)
      local sx%,sy%,Xx%,Yy%
      
      gcol col&
      sx%=sgn(X) : sy%=sgn(Y)
      Xx%=sx%*int(abs(X)) : Yy%=sy%*int(abs(Y))
      
      rem Plot the ball
      circle fill Xx%,Yy%,r%
      endproc
      rem =============================================================================
      
      rem Check if the Ball hits a wall
      def proc_poscheck(Xn,Yn,i&)
      
      ch%(i&)=0
      rem Has it hit the side?
      if (Xn-(r%+1))<0 then
        xn(i&)=r%
        ch%(i&)+=1
        u(i&)=-res*u(i&)
        v(i&)=v(i&)-g*tcum(i&)
      endif
      if (Xn+(r%+1))>xmax% then
        xn(i&)=xmax%-r%
        ch%(i&)+=1
        u(i&)=-res*u(i&)
        v(i&)=v(i&)-g*tcum(i&)
      endif
      
      rem Has it hit the floor or the ceiling?
      if (Yn-(r%+1))<0 then
        yn(i&)=r%
        ch%(i&)+=1
        v(i&)=-res*(v(i&)-g*tcum(i&))
      endif
      if (Yn+(r%+1))>ymax% then
        yn(i&)=ymax%-r%
        ch%(i&)+=1
        v(i&)=-res*(v(i&)-g*tcum(i&))
      endif
      
      rem Reset launch position and time
      if ch%(i&)>0 then x(i&)=xn(i&) : y(i&)=yn(i&) : tcum(i&)=0.0
      
      rem Record the hit
      ht&(i&)=1
      endproc
      rem =============================================================================
      
      rem When the balls have hit each other.
      def proc_hit(Xn1,Yn1,Xn2,Yn2,b1&,b2&)
      
      rem Components of unit normal to line of contact
      alp=(Yn2-Yn1)/R : bet=-(Xn2-Xn1)/R
      
      rem Get velocities
      ua=u(b1&) : va=v(b1&)-g*tcum(b1&)
      ub=u(b2&) : vb=v(b2&)-g*tcum(b2&)
      
      rem Velocity components along line of contact
      velaR=-ua*bet + va*alp
      velbR=-ub*bet + vb*alp
      
      rem Velocity components in direction of normal
      velaT=ua*alp + va*bet
      velbT=ub*alp + vb*bet
      
      rem New velocity components along line of contact after collision
      nvelaR=(velaR*(1.0-res) + velbR*(1.0+res))/2.0
      nvelbR=(velaR*(1.0+res) + velbR*(1.0-res))/2.0
      
      rem Recover Cartesian components of velocities
      u(b1&)=-nvelaR*bet + velaT*alp
      v(b1&)=nvelaR*alp + velaT*bet
      
      u(b2&)=-nvelbR*bet + velbT*alp
      v(b2&)=nvelbR*alp + velbT*bet
      
      rem Reset parameters
      tcum(b1&)=0.0 : tcum(b2&)=0.0
      x(b1&)=xn(b1&) : y(b1&)=yn(b1&) : x(b2&)=xn(b2&) : y(b2&)=yn(b2&)
      
      rem Record the hit
      ht&(b1&)=1 : ht&(b2&)=1
      endproc
      rem =============================================================================
      
      rem Set up use of Full Screen
      def proc_fullscreen(return xscreen%, return yscreen%)
      
      sys "GetSystemMetrics", 0 to xscreen%
      sys "GetSystemMetrics", 1 to yscreen%
      sys "SetWindowLong",@hwnd%,-16,&16000000
      sys "SetWindowPos",@hwnd%,-1,0,0,xscreen%,yscreen%,0
      vdu 23,22,xscreen%;yscreen%;8,16,16,1 : rem Set fullscreen mode
      mouse off : off : rem Turns off the Mouse Pointer and the Cursor
      endproc
      rem =============================================================================
      
      rem Kick a ball
      def proc_kick(A$)
      local a&,s%,i&
      
      a&=instr(ref$,A$)
      if (a&>0 and a&<(ball&+1)) then
        s%=2*rnd(2)-3 : u(a&)=s%*(275+rnd(100)) : v(a&)=200+rnd(100)
        x(a&)=xn(a&) : y(a&)=yn(a&) : tcum(a&)=0.0
      endif
      
      if a&=16 then
        for i&=1 to ball&
          s%=2*rnd(2)-3 : u(i&)=s%*(275+rnd(100)) : v(i&)=200+rnd(100)
          x(i&)=xn(i&) : y(i&)=yn(i&) : tcum(i&)=0.0
        next i&
      endif
      endproc
      rem =============================================================================