import heapq
def sim(arrivals,flowsizes,linkspeed,logfunc):
    simtime = 0.0
    activeflows = []
    nextarr = arrivals.next()
    while simtime<100:
        if activeflows:
            throughputperflow = linkspeed/float(len(activeflows))
            smallestflowsize = activeflows[0][0]
            nextdep = simtime + smallestflowsize/throughputperflow
            elapsedtime,simtime = min(nextdep,nextarr)-simtime, min(nextdep,nextarr)
            activeflows = [(sizeleft-throughputperflow*elapsedtime,starttime) 
                           for sizeleft,starttime in activeflows]
        else:
            nextdep,simtime = float('inf'),nextarr
        if nextdep<=nextarr:
            departingflow = heapq.heappop(activeflows)
            departingflowstarttime = departingflow[1]
            logfunc( (simtime,departingflowstarttime) )
        else:
            try:
                newflowsize = flowsizes.next()
                heapq.heappush(activeflows, (newflowsize,simtime))
                nextarr = simtime + arrivals.next()
            except StopIteration:
                break


import math, random
def rexp(rate):
    random.seed('geranium')
    while True:
        yield -1.0/rate * math.log(random.random())

class LogDuration:
    def __init__(self): self.tot, self.n = 0,0
    def logdeparture(self,record):
        deptime,starttime = record
        self.tot = self.tot+(deptime-starttime)
        self.n = self.n+1

results = []
for arrivalrate in [(x+1)/10.0 for x in range(200)]:
    logger = LogDuration()
    sim(rexp(arrivalrate),rexp(1),10, logger.logdeparture)
    results.append( (arrivalrate,logger.tot/logger.n) )
with open('results.csv','wt') as f:
    f.write('arrivalrate,meanflowduration\n')
    for res in results:
        f.write( ','.join([str(s) for s in res]) +'\n')