Distributions

Back Up Next

Running a simulation with different distributions elegantly

(Author: Klaus Muller; kgmuller at users.sourceforge.net)

In simulation modelling, one often has the situation that one wants to run the same model with different distributions (normal, uniform, empirical, . . ) for a random variate, e.g. for "what if" investigations. This is simple to achieve by editing the program text, but this is tedious and error prone. One also ends up with several different program versions.

Python comes to the rescue. With functions being first class objects and with generators, distributions are easily encapsulated so that the distributions can be referenced generically. They can also be given parameters once, without repeating them in the program text.

Here is an example which shows the elegance which can be achieved. Let's assume that for an arrival process, several distributions for inter-arrival times are to be studied, i.e. normal, uniform and empirical distributions. By 

bulletembedding each distribution in a Python generator, and
bulletiterating over the simulation for all three distributions, and
bulletassigning the generators to the variable 'arrival', and
bulletusing 'arrival.next()' to get the next value of the inter-arrival time random variate

one comes up with this clean program:

from SimPy.Simulation import *
import random

def arrivalsUni(low,high):
    """Generator for uniform distribution between 'low' and 'high'
    """
    while True:
        yield random.uniform(low,high)

def arrivalsNormal(mean,sigma):
    """Generator for normal distribution with mu='mean'
       and sigma='sigma'
    """
    while True:
        yield random.normalvariate(mean,sigma)

def arrivalsEmpirical(valuelist):
    """Generator for drawing from a list 'valuelist' of 
       empirical/observed values
    """
    while True:
        yield random.choice(valuelist)

class Arrivals(Process):
    def generateArrivals(self):
        i=0
        while True:
            yield hold,self,arrivals.next()
            i+=1
            print "%s arrival nr. %s"%(now(),i)

for arrivals in [arrivalsUni(low=10,high=20),
                 arrivalsEmpirical(valuelist=[10,12,14,19,20,20]),
                 arrivalsNormal(mean=15,sigma=2.0)]:
    initialize()
    a=Arrivals()
    activate(a,a.generateArrivals())
    print "\nNext distribution"
    simulate(until=30)
OUTPUT:

Next distribution
17.8638547229 arrival nr. 1

Next distribution
20 arrival nr. 1

Next distribution
15.8121137485 arrival nr. 1
29.6657638097 arrival nr. 2

One can embed any distribution this way, even distributions returning e.g. a value input interactively by a user at run time.

  

horizontal rule

©Copyright 2004, 2005, 2006, 2007, 2008 SimPy Developer Team.
For problems or questions regarding this web contact SimPy webmaster.
Last updated: 20/06/08.