MM1排队系统

#coding=utf-8
import time
import random as rd
#import math
import pylab as pl
def simulate(nameda,u):
    arriveInternal=0
    serviceTime=0
    serverCondition=0  #0表示服务台空闲，1表示忙碌
    waitTime=0
    stayTime=[]
    tempArriveEventTime=[10000000]  
    tempLeaveEventTime=[10000000]
    #tempArriveEventTime.append(10000000)  
    simulateTotalTime=1000000 #仿真时长
    simulateClock=0
    EntityNumber=0    
    queneLength=0  
    averageStayTime=0    
    currentArriveTime=0
    currentLeaveTime=0
    formerArriveTime=0
    formerLeaveTime=0
    lossMessage=0
    #生成仿真表，包括实体的到达时间和离开时间
    while simulateClock<=simulateTotalTime: 
            #保存上一个实体的到达、离开时间       
            formerArriveTime=currentArriveTime
            formerLeaveTime=currentLeaveTime
            #tempFormerLeaveTime=currentLeaveTime
            ####################随机数生成
            arriveInternal=rd.expovariate(nameda)
            serviceTime=rd.expovariate(u) 
            #######到达时间
            if EntityNumber==0:
                currentArriveTime=arriveInternal       
            else:
                currentArriveTime=formerArriveTime+arriveInternal   
            #仿真结束与否判断   
            if currentArriveTime>=simulateTotalTime:
                break  
            #等待时间        
            if EntityNumber>0 and currentArriveTime<formerLeaveTime:
                waitTime=formerLeaveTime-currentArriveTime              
            else:
                waitTime=0      
            ##离开时间                
            currentLeaveTime=currentArriveTime+waitTime+serviceTime
            #仿真结束与否判断
            if currentLeaveTime>simulateTotalTime:
                break  
            ###仿真时钟
            if EntityNumber==0:
                simulateClock=currentLeaveTime#表示实体接受服务，并离开。
                stayTime.append(currentLeaveTime)
                        
            elif queneLength>=10:
                #print("message missed")
                currentLeaveTime=formerLeaveTime
                lossMessage+=1
            ###排队的实体    
                                  
            if currentArriveTime>=simulateClock and EntityNumber!=0:#开始处理队列的实体                                   
                if len(tempArriveEventTime)>1:  
                    minValue=min(tempArriveEventTime) 
                    minIndex=tempArriveEventTime.index(minValue)       
                    while minValue<=simulateClock:                    
                        simulateClock=tempLeaveEventTime[minIndex]#实体接受了服务的情况下才会重置仿真时钟
                        tempArriveEventTime.pop(minIndex)
                        tempLeaveEventTime.pop(minIndex)
                        if queneLength>0:
                            queneLength-=1
                        minValue=min(tempArriveEventTime) 
                        minIndex=tempArriveEventTime.index(minValue)
                        #print("one processed") 
                else:
                    simulateClock=currentLeaveTime  
                    stayTime.append(currentLeaveTime-currentLeaveTime)
                    #print("one processed")
            if currentArriveTime<simulateClock and EntityNumber!=0 :
                tempArriveEventTime.append(currentArriveTime)
                tempLeaveEventTime.append(currentLeaveTime)
                queneLength+=1  
                stayTime.append(currentLeaveTime-currentArriveTime)  
                
                                                    
            ###############        
            EntityNumber+=1        
            if EntityNumber%1000000==0:
                print("arriveTime:"+str(currentArriveTime)+"，leaveTime:"+str(currentLeaveTime))
    lossRate=lossMessage/EntityNumber
    averageStayTime=sum(stayTime)/len(stayTime)
    print("loss rate:"+str(lossMessage/EntityNumber)) 
    print("lossMessage:"+str(lossMessage))
    print("entity:"+str(EntityNumber))
    return lossRate,averageStayTime
pl.figure(1) #创建图表1，并选择图表1
ax1=pl.subplot(211) #在图表1中创建子图1
ax2=pl.subplot(212) #在图表1中创建子图2
start=time.clock()
for i in range(100):
    a,b=simulate(8, 10)
    print(str(i+1)+" simulate")
    pl.sca(ax1)
    #pl.plot(i,a,'o')
    pl.ylabel("loss rate") 
    pl.scatter(i,a)
    #pl.ylim(0,1)   
    pl.sca(ax2)
    #pl.plot(i,b,'o')
    pl.scatter(i,b)
    pl.ylabel("stay time")
minutes=(time.clock()-start)/60
print("cost minutes:"+str(minutes))
pl.show()       
            
理论部分参考：http://mi.gxu.edu.cn/jpkc/yunchouxue/resource_1info.php?id=29