! Simulation of velocity vs. time for a bicyclist
! assume no air resistance
! Program to accompany "Computational Physics" by N. Giordano
! Copyright Prentice Hall 1997
program bike
   option nolet
   library "sglib*","sgfunc*"
   dim t(5000),velocity(5000)
   call initialize(t,velocity,dt,power,mass,nmax)
   call calculate(t,velocity,dt,power,mass,nmax)
   call display(t,velocity)
end
! t() = time    v() = velocity
! dt = time step    power = rider power
! mass = mass of rider + bicycle    
sub initialize(t(),v(),dt,power,mass,nmax)
   t(1) = 0
   v(1) = 4   ! m/s - all units are mks
   dt = 1     ! second 
   power = 400  ! watts
   mass = 70    ! kg
   tmax = 200   ! seconds
   nmax = tmax / dt  ! total number of time steps required
end sub
sub calculate(t(),v(),dt,pmax,mass,nmax)
   for i = 2 to nmax
      t(i) = t(i-1) + dt
      v(i) = v(i-1) + pmax * dt / (mass * v(i-1))    ! Euler method
   next i
   mat redim t(nmax),v(nmax)    ! trim arrays to the size actually used
end sub
sub display(t(),v())
   call settitle("Bicycle velocity vs. time")
   call sethlabel("time (s)")
   call setvlabel("Velocity (m/s)")
   call datagraph(t,v,3,0,"black")
end sub