A simple Monte Carlo simulation

Arjen Markus (25 june 2007) Monte Carlo simulations are useful in many cases where the system being studied reacts to random input, but also if you want to determine the integral of a function over a complicated area. Such simulations usually require a lot of steps and sometimes numerous runs to get meaningful results.

The script below is an attempt to model a simple queueing system. See the comments for more information. It tries to answer a vexing question: how long does one have to wait (on average)?

 # simple_queue.tcl
 #     A simulation of a queue:
 #     - On average 1 person per 5 minutes enters the system
 #     - Handling his/her request takes either 3 minutes or 10 minutes
 #       (with a 1/10 chance for the longer request).
 #     - The average handling time is therefore (9*3+10)/10 = 3.7 minutes
 #     - The system can handle this - on average
 #     - Question is, however, what is the mean waiting time? If you
 #       are in the queue behind two people whose requests will take 10
 #       minutes, you will be waiting longer than if these two both have
 #       the short requests
 #

 # randomBinomial --
 #     Return a random value (0 or 1) according to the binomial distribution
 #
 # Arguments:
 #     p           Chance of returning 1
 #
 # Result
 #     Random value of 0 or 1
 #
 proc randomBinomial {p} {
     return [expr {rand()>$p? 0 : 1}]
 }

 # runQueue --
 #     Run the queuing system
 #
 # Arguments:
 #     nosteps       Number of steps to run the system
 #
 # Result:
 #     None
 #
 # Side effects:
 #     Fills the global variables waiting_times and queue
 #
 proc runQueue {nosteps} {
     global waiting_times
     global queue

     #
     # Start with an empty queue
     #
     set queue {}

     #
     # Loop over time
     #
     for { set time 0 } { $time < $nosteps } { incr time } {

         #
         # New client entering the system?
         #
         if { [randomBinomial 0.2] == 1 } {
             if { [randomBinomial 0.1] == 1 } {
                 set request 10
             } else {
                 set request 3
             }
             set person [list $time $request]

             lappend queue $person
         }

         #
         # Handle the current request, if any
         #
         if { [llength $queue] > 0 } {
             set person [lindex $queue 0]
             foreach {arrival_time request} $person {break}
             if { $request <= 0 } {
                 # We are done with this client
                 set queue [lrange $queue 1 end]
                 lappend waiting_times [expr {$time-$arrival_time}]
             } else {
                 # Reduce the time left for the request
                 set request [expr {$request-1}]
                 lset queue 0 1 $request
             }
         }
     }
 }

 # main --
 #     Run the queueing system, timestep is one minute:
 #     For a period of 8 hours (480 minutes), do the following
 #     every minute:
 #     - If someone comes in, add them to the queue with the arrival
 #       time and the type of request
 #     - If there is someone in the queue, take one minute off the
 #       remaining time for their request, if it is positive. Otherwise
 #       remove the entry from the queue and record the waiting time
 #
 #     To get a fairer estimate, we run the queueing system 50 times

 set av_av_time  0
 set av_max_time 0
 set av_queue    0
 set av_count    0
 for { set i 0 } { $i < 50 } { incr i } {
     set waiting_times {}
     runQueue 480

     set count        0
     set average_time 0
     set max_time     0

      foreach time $waiting_times {
         set average_time [expr {$average_time + $time}]
         incr count
         if { $time > $max_time } {
             set max_time $time
         }
     }
     set average_time [expr {$average_time/double($count)}]

     set av_count    [expr {$av_count    + $count}]
     set av_av_time  [expr {$av_av_time  + $average_time}]
     set av_max_time [expr {$av_max_time + $max_time}]
     set av_queue    [expr {$av_queue    + [llength $queue]}]

     puts "Simulation $i: [format %.1f $average_time] - $max_time"
 }

 puts "Sample averages:"
 puts "----------------"
 puts "Number of persons:     [format %.1f [expr {$av_count/50.0}]]"
 puts "Average waiting time:  [format %.1f [expr {$av_av_time/50.0}]]"
 puts "Maximum waiting time:  [format %.1f [expr {$av_max_time/50.0}]]"
 puts "People still waiting:  [format %.1f [expr {$av_queue/50.0}]]"

Category Mathematics

Category Simulator