#!/usr/bin/env python
# jacksonnetwork.py:  Jackson network Littles law
# pyt090b.py. See sms090b.tex for simscript version
#        updated and simplified, 2004 Dec gav
#        $Author: kgmuller $    $Revision: 1.1.1.1 $
#        $Date: 2005/01/15 15:16:52 $ 

__version__= '\nModel: jacksonnetwork.py'
"""
  Messages arrive randomly at rate 1.5 per unit time at a
  communication network with 3 nodes (computers). Each computer can
  queue messages and service-times are exponential with mean m_i at
  node i. These values are given in the column headed n_i in the
  table below.  On service at a node a message transfers to one of the
  other nodes or leaves the system. The probability of a transfer from
  node i to node j is p_ij. Node 4 means outside the system so
  a message will transfer out of the system from node i with
  probability pi4.
   
These transition probabilities are given in the following table.
    Node     m_i  p_i1   p_i2  p_i3    p_i4  
     i                                (leave)   
     1       1.0   0      0.5  0.5     0     
     2       2.0   0      0    0.8     0.2  
     3       1.0   0.2    0    0       0.8  
   
  Your task is to estimate the average time taken for jobs going
  through the system and the average number of jobs in the system.
  Does Little's Law hold for the network?  Run 10000 messages through
  the system to make these observations.  
"""

from __future__ import generators   # needed for Python 2.2
from SimPy.Simulation import *
from random import Random,expovariate,uniform


## def choose1d(P,g):
##     """  return a random choice from a set j = 0..n-1
##        with probs held in list P[j]
##        g = random variable
##        call:  next = choose1d(P,g)
##     """
##     U = g.random()
##     sumP = 0.0
##     for j in range(len(P)):  # j = 0..n-1
##         sumP +=  P[j]
##         if U < sumP: break
##     return(j)

def choose2dA(i,P,g):
    """  return a random choice from a set j = 0..n-1
       with probs held in list of lists P[j] (n by n)
       using row i
       g = random variable
       call:  next = choose2d(i,P,g)
    """
    U = g.random()
    sumP = 0.0
    for j in range(len(P[i])):  # j = 0..n-1
        sumP +=  P[i][j]
        if U < sumP: break
    return(j)

class Msg(Process):
    """a message"""

    def execute(self,node):
        """ executing a message """
        global noInSystem
        startTime = now()
        noInSystem += 1
        ##print "DEBUG noInSystm = ",noInSystem
        NoInSystem.observe(noInSystem)
        self.trace("Arrived node  %d"%(node,))
        while node <> 3:
            yield request,self,processor[node]
            self.trace("Got processor %d"%(node,))
            st = Mrv.expovariate(1.0/mean[node])
            yield hold,self,st
            yield release,self,processor[node]
            self.trace("Finished with %d"%(node,))
            node = choose2dA(node,P,Mrv)
            self.trace("Transfer to   %d"%(node,))
        TimeInSystem.tally(now()-startTime)        
        self.trace(    "leaving       %d %d in system"%(node,noInSystem))
        noInSystem -= 1
        NoInSystem.accum(noInSystem)
        
    def trace(self,message):
        if MTRACING: print "%7.4f %3d %10s"%(now(),self.name, message)
            

## -----------------------------------------------
class Generator(Process):
 """ generates a sequence of msgs """    
     
 def execute(self,rate,maxN):
     self.count=0   # hold number of messages generated
     g = Random(11335577)
     while  (self.count < maxN):
         self.count+=1
         p = Msg(self.count)
         activate(p,p.execute(node=startNode))
         yield hold,self,g.expovariate(rate)
     self.trace("generator finished with "+`self.count`+" ========")
     
 def trace(self,message):
     if GTRACING: print "%7.4f \t%s"%(now(), message)

## -----------------------------------------------
print __version__
# generator:
rate = 1.5
maxNumber = 100000
GTRACING = 0
# messages:
noInSystem = 0
MTRACING = 0
startNode = 0
Mrv = Random(77777)

TimeInSystem=Monitor()
NoInSystem= Monitor()

processor = [Resource(1),Resource(1),Resource(1)]
# average service times at each node:
mean=[1.0, 2.0, 1.0]
# transition matrix:
P = [[0,   0.5, 0.5, 0  ],
     [0,   0,   0.8, 0.2],
     [0.2, 0,   0,   0.8]]

initialize()
g = Generator("Generator")
activate(g,g.execute(rate,maxNumber))
simulate(until=5000.0)

print "Mean number in system = %10.4f"%(NoInSystem.timeAverage(),)
print "Mean delay in system  = %10.4f"%(TimeInSystem.mean(),)
print "Total time run        = %10.4f"%(now(),)
print "Total jobs arrived    = %10d"%(g.count)
print "Total jobs completed  = %10d"%(TimeInSystem.count(),)
print "Average arrival rate  = %10.4f"%(g.count/now(),)