Lecture 22 (April 13, 2019)¶
In this lecture, we use the simulation of a repairment system to get some more practice on the discrete event simulation.
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 a repairment system¶
1. We have $N+S$ functional machines. There are no differences in the machines. $N$ machines are directly put into use in the operation room. $S$ are held as spare machines.
2. The time it takes for a machine to break down follows a specific distribution. Whenever a machine breaks down, it will be taken to the maintenance room to be fixed. If there is a spare machine, it will be put in use.
3. The time it takes for a machine to be fixed follows a distribution. Machines are fixed one at a time. If a machine gets fixed, it will be placed on held as a spare machine if we have N machines working. Otherwise, the fixed machine will be taken to the operation room directly and put in use. The maintenance room can fix a machine at one time.
4. The simulation will stop when all machines break down.
Analaysis of the system¶
Events:
We have three possible events:
A machine breaks down (denote the time using t_break)
A machine gets fixed (denote the time using t_fixed)
States
Which variables will lead to change in the event times?
An array of the the pending break time of all the machines (t_breakA):
the number of machines in the maintainance room (n_maintain)
the number of machines that are currently working (n_working)
the number of machines that are currently working (n_spare)
Algorithm¶
**Initialization**
t=0
t_breakA~ N samples from the distribution we assume
n_maintain=0
n_spare=S
n_working=N
t_break is the lowest value of t_breakA
t_fixed is np.inf
**Inside the while loop**
**Case 1: break down event**
update clock
n_maintain=n_maintain+1
if n_maintain=1 , this machine will start receiving repairment immediately. We will generate a new t_fixed
If n_spare>=1
There will be a new machine moved from the spare room to the operation room. We will need to generate a breakdown time and replace the previous breakdown time for the machine that broke down. Then t_break will then be updated to the smallest value among all the breakdown times.
We will update n_spare=n_spare-1
However, if n_spare=0
it means we do not have another machine to replace this machine. The time for this machine slot will be set to infinity. Then t_break will then be updated to the smallest value among all the breakdown times.
We will update n_working to n_working-1
**Case 2: fixed event**
update clock
n_maintain=n_maintain-1
if n_maintain=0, there will be no machine needs to be fixed. We will set t_fixed equal to infinity.
Otherwise, we will generate a new t_fixed time.
if n_working is N, this machine will be taken to the spare room, meaning we will update n_spare to n_spare+1
Otherwise, this machine will be put in work. We will generate a new breakdown time for this machine and replace an infinite value inside t_breakA. Then t_break will then be updated to the smallest value among all the breakdown times. Also, n_working=n_working+1
Practice¶
Let’s make the following assumptions
(1) Initially, there are 10 machines in the operation room and 5 machines in the spare room
(2) the time it takes for a machine to break down follows $Exp(\lambda=1)$
(3) the time it takes to fix a machine follows $Unif(0,0.5)$
Construct a function to report time it takes to run out of machines.
Python programming¶
Function
Explanation
np.argmax(Array)
get the index of the max value of an array. If there are several of values equal to max, it will report the lower index
np.argmin(Array)
get the index of the min value of an array. If there are several of values equal to min, it will report the lower index
np.where(condition)
get the indices of the elements satisfying a condition (e.g. Array>a, Array==a, Array

In [0]:
mean=np.mean(returns)
mean
Out[0]:
34.61407105231429
In [0]:
var=np.var(returns,ddof=1)
lower=mean-1.96*var**0.5/2000**0.5
upper=mean+1.96*var**0.5/2000**0.5
print(lower,upper)
33.20741667554704 36.02072542908154
In [0]: