Lecture 21 (April 4, 2019)¶
In this lecture, we discuss how to model a two-server queueing system in parallel.
In [0]:
import numpy as np
import matplotlib.pyplot as plt
import math
import scipy.special as spsp
import scipy.stats as spst
Setting of two-server queueing system in parallel¶
(1) There are two servers with different distributions for the service time
(2) Consumers arrive following a specific arrival process.
(3) When a consumer arrives, the consumer will receive the service from the empty server if only one server is empty
(4) If both servers are empty, the consumer will receive the service from the first server.
(5) If both servers are busy, consumers will wait. Once a server becomes empty, the cusmer who arrived first will receive the service first.
Analaysis of the system¶
Events:
We have three possible events:
A consumer arrives (denote the time using t_arrival)
A consumer finished receiving the service from server 1 (denote the time using t_depart1)
A consumer finished receiving the service from server 2 (denote the time using t_depart2)
The system simulation is moving the clock sequentially to the earliest pending event, update the system state, update the pending event time, and then move to the earliest pending event again.
States
Which variables will lead to change in the event times?
Number of customer of customers currently in the system. ($n$):
Whether server 1 is busy (server_1)
Whether server 2 is busy (server_2)
Algorithm¶
To put it in another way, we have
**Case 1: arrival event**
update clock
$n=n+1$
if $n=1$, then set server_1= busy. and generate a new t_depart1.
if $n=2$, then
if server_1= =busy, then set server_2=busy generate a new t_depart2.
otherwise, set server_1=busy generate a new t_depart1.
Generate the t_arrival for the next customer
**Case 2: departure from server 1**
update clock
$n=n-1$
if $n==1$ or $n==0$, then set server_1= idle and update t_depart1 to np.inf.
otherwise, set server_1=busy and generate a new t_depart1.
**Case 3: departure from server 2**
update clock
$n=n-1$
if $n==1$ or $n==0$, then set server_1= idle. and update t_depart1 to np.inf.
otherwise, set server_2=busy and generate a new t_depart2.
Other important things?¶
The other thing that is important that if we want to study the average time consumers spent in the system, we always need to not only save the arrival time and the departure time. At the same time, we also need to keep an array that stores the departure sequence. In other words, we would like to have an array like the following
[1,2,5,3,4].
This toy array tells us the first arrival left first, the second arrival left second, the fifth arrival left third, the third arrival left fouth, while the fourth arrival left last. In this way, we will be able to know the arrival time of all the customers.
In order to generate an array like this, we need to keep track of the current ID of the current customer being served.
Thus, we will need to create three additional variables
n_total to track the total number of arrivals. This will be the ID we use for each consumer. server1_ID to track of the ID of the current consumer served by server 1. server2_ID to track of the ID of the current consumer served by server 2.
When server1 and server2 are not taking customers, we set them to 0. When server 1 takes a customer
if this customer is a new arrival, server1_ID=n_current
if this customer is not a new arrival, server1_ID=np.max(server1_ID, server2_ID)
When server 2 takes a customer
if this customer is a new arrival, server2_ID=n_current
if this customer is not a new arrival, server2_ID=np.max(server1_ID, server2_ID)
When a consumer departs, we will record server1_ID or server2_ID.
In [0]:
Practice¶
It is estimated that the arrival of customers will follow a poisson process with $\lambda=100$ for this hour. Besides the first counter with the serving speed following an exponential distribution with $\lambda=110$, the company is planning to add a second counter with the serving speed following an exponential distribution with $\lambda=80$. Thus, we will have a system with two servers in parallel. The system setup is the same given in the setting of the lecture. Simulate the system for 1 unit of time. Report the average time consumers are in the system.
In [0]:
#Version 1 of tracking
def two_server_parallel():
#initialize clock
t=0
#initialize the states
n_current=0
busy1=False
busy2=False
#generate the pending event times
lmbda_a=100
lmbda_s1=110
lmbda_s2=80
t_arrival=t-1/lmbda_a*np.log(np.random.rand())
t_depart1=np.inf
t_depart2=np.inf
#other initialization
##Collecting system output (history)
##array to store Arrival times
arrival_times=np.array([])
##array to store Departure times
depart_times=np.array([])
##array to store Departure sequence
depart_sequence=np.array([])
##variable to track total number of customers arrived (ID)
N=0
##variable to track the ID of customer receiving the service from server 1
Server1_ID=0
##variable to track the ID of customer receiving the service from server 2
Server2_ID=0
#while ( condition to continue the simulation):
while (t<=1):
#print(t_arrival,t_depart1, t_depart2)
#conditional statement that identifies arrival event happens:
if t_arrival==min(t_arrival,t_depart1, t_depart2):
#update clock
t=t_arrival
#store arrival time (history)
arrival_times=np.append(arrival_times,t)
#####NEW#######
#Pre-generate a departure time
depart_times=np.append(depart_times,0)
#Generate ID
N=N+1
#update system states
n_current=n_current+1
#update pending event times
t_arrival=t-1/lmbda_a*np.log(np.random.rand())
if n_current==1:
t_depart1=t-1/lmbda_s1*np.log(np.random.rand())
Server1_ID=N
busy1=True
elif n_current==2:
if busy1==True:
t_depart2=t-1/lmbda_s2*np.log(np.random.rand())
busy2=True
Server2_ID=N
else:
t_depart1=t-1/lmbda_s1*np.log(np.random.rand())
busy1=True
Server1_ID=N
#possibly update the ID for the server1/server 2
#conditional statement that identifies departure from server 1 happens:
elif t_depart1==min(t_arrival,t_depart1, t_depart2):
#update clock
t=t_depart1
#store departure time (history)
#depart_times=np.append(depart_times,t)
##replace the pre-generated depart_times using the actual depart time
depart_times[Server1_ID-1]=t
#store departure ID
#depart_sequence=np.append(depart_sequence,Server1_ID)
#update system states
n_current=n_current-1
if n_current==0 or n_current==1:
#update pending event times
t_depart1=np.inf
busy1=False
Server1_ID=0
else:
t_depart1=t-1/lmbda_s1*np.log(np.random.rand())
Server1_ID=max(Server1_ID, Server2_ID)+1
#update the ID for the server1
#conditional statement that identifies departure from server 2 happens:
else:
#update clock
t=t_depart2
#store departure time (history)
#depart_times=np.append(depart_times,t)
##replace the pre-generated depart_times using the actual depart time
depart_times[Server2_ID-1]=t
#store departure ID
#depart_sequence=np.append(depart_sequence,Server2_ID)
#update system states
n_current=n_current-1
if n_current==0 or n_current==1:
#update pending event times
t_depart2=np.inf
busy2=False
Server2_ID=0
else:
t_depart2=t-1/lmbda_s2*np.log(np.random.rand())
Server2_ID=max(Server1_ID, Server2_ID)+1
#update the ID for the server2
##when simulation stops, we compute the average time in the system
print(arrival_times)
print(depart_times)
# print(depart_sequence)
#Valid departures index
#depart_times=depart_times[depart_sequence.astype(int)-1]
Departed=np.arange(0,N)[(depart_times>0)*(depart_times<=1)]
depart_times=depart_times[Departed]
arrival_times=arrival_times[Departed]
#arrival_times=arrival_times[0:len(depart_times)]
#filtering out the actual arrive times
avg_time=np.mean(depart_times-arrival_times)
return avg_time
two_server_parallel()
[0.00658686 0.0303028 0.03934869 0.04110103 0.0534258 0.05550251
0.06586218 0.07671512 0.08846409 0.09449326 0.09905949 0.10342168
0.10511168 0.11307374 0.115817 0.11642926 0.11840041 0.12142965
0.13125465 0.14396361 0.15412766 0.15898942 0.16453371 0.1956923
0.19956811 0.21509251 0.24048005 0.26238031 0.26354609 0.26776614
0.27631669 0.31338242 0.32541205 0.32789651 0.33109804 0.33900691
0.39829782 0.42204791 0.44802735 0.44968111 0.46192973 0.48899345
0.4974726 0.50375424 0.50773947 0.51320324 0.51909068 0.52227621
0.53907611 0.54357425 0.54495975 0.54524592 0.55761666 0.57991143
0.57993732 0.59389115 0.59556924 0.62421462 0.63012637 0.63156516
0.64832346 0.65599019 0.6706939 0.69003464 0.69758517 0.70096897
0.71055473 0.71171009 0.72135336 0.7280733 0.73700034 0.73853258
0.75340205 0.77060165 0.77094186 0.77196414 0.77261556 0.77338982
0.78708518 0.78907347 0.79566953 0.79727678 0.80233821 0.80497574
0.80508129 0.80834897 0.83403755 0.83770901 0.84092592 0.84142003
0.84279857 0.86799769 0.86903172 0.88324577 0.89066619 0.89374696
0.89856531 0.90015189 0.90203581 0.91752416 0.93690683 0.93724291
0.9404656 0.94690766 0.95623993 0.96128333 0.96928657 0.9754526
1.00755441]
[0.01805477 0.04874354 0.04550304 0.04773975 0.07897531 0.06602765
0.0806323 0.08431356 0.09247244 0.10026281 0.10007584 0.12235408
0.12105827 0.13219452 0.12400267 0.13266074 0.14100238 0.15201244
0.1455978 0.16742708 0.19155075 0.18970656 0.2143774 0.20326327
0.20824199 0.22688303 0.2556879 0.26979866 0.27123504 0.27397204
0.29776192 0.33244436 0.35013045 0.34048171 0.34246093 0.34525517
0.41969445 0.42215656 0.45579847 0.45471224 0.48120513 0.50020723
0.5169598 0.51818999 0.52023729 0.51896668 0.53073722 0.52407542
0.53917759 0.54913842 0.5474453 0.55042407 0.55865527 0.58065012
0.59679633 0.59650076 0.60886138 0.63192895 0.63535009 0.6335253
0.6489103 0.66002745 0.67552642 0.69893111 0.72317411 0.71070049
0.71505178 0.71983422 0.73305183 0.75382973 0.75434496 0.75529027
0.76138673 0.78410601 0.77098782 0.78142255 0.80586751 0.79475478
0.80251425 0.82105866 0.81948162 0.90841243 0.85198525 0.86461814
0.86519846 0.86603893 0.89833977 0.90120522 0.90311326 0.91072534
0.92027627 0.91291147 0.91402582 0.91615043 0.91654776 0.918025
0.9294088 0.94722169 0.93101405 0.93661302 0.95618526 0.95642375
0.96110326 0.98287463 0.96377389 0.97710165 0.98392345 0.
0. ]
Out[0]:
0.020132684451994745
In [0]:
#Version 2 of tracking
def two_server_parallel():
#initialize clock
t=0
#initialize the states
n_current=0
busy1=False
busy2=False
#generate the pending event times
lmbda_a=100
lmbda_s1=110
lmbda_s2=80
t_arrival=t-1/lmbda_a*np.log(np.random.rand())
t_depart1=np.inf
t_depart2=np.inf
#other initialization
##Collecting system output (history)
##array to store Arrival times
arrival_times=np.array([])
##array to store Departure times
depart_times=np.array([])
##array to store Departure sequence
depart_sequence=np.array([])
##variable to track total number of customers arrived (ID)
N=0
##variable to track the ID of customer receiving the service from server 1
Server1_ID=0
##variable to track the ID of customer receiving the service from server 2
Server2_ID=0
#while ( condition to continue the simulation):
while (t<=1):
#print(t_arrival,t_depart1, t_depart2)
#conditional statement that identifies arrival event happens:
if t_arrival==min(t_arrival,t_depart1, t_depart2):
#update clock
t=t_arrival
#store arrival time (history)
arrival_times=np.append(arrival_times,t)
#Generate ID
N=N+1
#update system states
n_current=n_current+1
#update pending event times
t_arrival=t-1/lmbda_a*np.log(np.random.rand())
if n_current==1:
t_depart1=t-1/lmbda_s1*np.log(np.random.rand())
Server1_ID=N
busy1=True
elif n_current==2:
if busy1==True:
t_depart2=t-1/lmbda_s2*np.log(np.random.rand())
busy2=True
Server2_ID=N
else:
t_depart1=t-1/lmbda_s1*np.log(np.random.rand())
busy1=True
Server1_ID=N
#possibly update the ID for the server1/server 2
#conditional statement that identifies departure from server 1 happens:
elif t_depart1==min(t_arrival,t_depart1, t_depart2):
#update clock
t=t_depart1
#store departure time (history)
depart_times=np.append(depart_times,t)
#store departure ID
depart_sequence=np.append(depart_sequence,Server1_ID)
#update system states
n_current=n_current-1
if n_current==0 or n_current==1:
#update pending event times
t_depart1=np.inf
busy1=False
Server1_ID=0
else:
t_depart1=t-1/lmbda_s1*np.log(np.random.rand())
Server1_ID=max(Server1_ID, Server2_ID)+1
#update the ID for the server1
#conditional statement that identifies departure from server 2 happens:
else:
#update clock
t=t_depart2
#store departure time (history)
depart_times=np.append(depart_times,t)
#store departure ID
depart_sequence=np.append(depart_sequence,Server2_ID)
#update system states
n_current=n_current-1
if n_current==0 or n_current==1:
#update pending event times
t_depart2=np.inf
busy2=False
Server2_ID=0
else:
t_depart2=t-1/lmbda_s2*np.log(np.random.rand())
Server2_ID=max(Server1_ID, Server2_ID)+1
#update the ID for the server2
##when simulation stops, we compute the average time in the system
print(arrival_times)
print(depart_times)
# print(depart_sequence)
#Valid departures index
arrival_times=arrival_times[depart_sequence.astype(int)-1]
#filtering out the actual arrive times
depart_times=depart_times[depart_times<=1]
arrival_times=arrival_times[0:len(depart_times)]
avg_time=np.mean(depart_times-arrival_times)
return avg_time
two_server_parallel()
[0.01341013 0.01343704 0.02168181 0.02992166 0.04025175 0.0406003
0.04436483 0.06182704 0.08295252 0.1032316 0.11437561 0.12859626
0.13036816 0.13309561 0.13508033 0.13546317 0.14900507 0.15198727
0.15471261 0.15483883 0.15572204 0.16377865 0.16487088 0.16525381
0.16834155 0.17657139 0.19254853 0.19801123 0.19925418 0.2063222
0.21024799 0.23426634 0.23591286 0.23954667 0.24723977 0.24859386
0.25409139 0.25982636 0.26348925 0.33258297 0.34449642 0.34731568
0.36566071 0.3715022 0.39400332 0.39796074 0.42097729 0.42652132
0.43790208 0.46233408 0.46879453 0.46985563 0.47006237 0.47240971
0.48993809 0.49305472 0.49471291 0.49667711 0.50283804 0.53172957
0.5604539 0.56248434 0.57098407 0.57813721 0.58457592 0.59208088
0.60480989 0.6161573 0.62845049 0.6491457 0.65796434 0.67084021
0.70619276 0.71424679 0.72162086 0.75072265 0.75976698 0.76895215
0.77299315 0.83076547 0.83643559 0.83767339 0.8442096 0.84719045
0.8643597 0.87708807 0.87915834 0.88817687 0.96838776 0.97871734
0.98086366 0.99031106]
[0.01735919 0.02400743 0.03253802 0.03598639 0.04114937 0.04779401
0.05528102 0.06192122 0.08771293 0.10824575 0.12095158 0.13310584
0.13527689 0.15406116 0.15477806 0.15670116 0.15805574 0.16125975
0.16313467 0.16595802 0.17190063 0.17313131 0.17422446 0.17967596
0.18918872 0.19674879 0.19945603 0.20324492 0.20706666 0.21224853
0.21321619 0.23668609 0.24190298 0.24505984 0.25227123 0.26128298
0.26834314 0.27138488 0.28466435 0.33347899 0.36155159 0.36263415
0.37182115 0.38011804 0.41360212 0.4291436 0.43637263 0.43863408
0.45979352 0.46761126 0.47180746 0.47497364 0.47671073 0.47947074
0.49778686 0.50244134 0.50486445 0.51531188 0.5210894 0.54586283
0.5668362 0.57159888 0.58091624 0.58768539 0.59085555 0.60387096
0.6119641 0.62776935 0.64616559 0.65141246 0.65984927 0.67429113
0.71166789 0.72742079 0.72817852 0.75672264 0.77143413 0.77153473
0.78507657 0.85188425 0.85286228 0.8609526 0.86111626 0.87409679
0.89744936 0.90018721 0.9052438 0.92324075 0.98466605 0.99317232
1.00000526]
Out[0]:
0.01147633125337451
In [0]: