(c) Tobias Hossfeld (Aug 2021)
This script and the figures are part of the following book. The book is to be cited whenever the script is used (copyright CC BY-SA 4.0):
Tran-Gia, P. & Hossfeld, T. (2021). Performance Modeling and Analysis of Communication Networks - A Lecture Note. Würzburg University Press. https://doi.org/10.25972/WUP-978-3-95826-153-2
We consider a delay loss system with 2 servers. If all $K$ servers are occupied, incoming arrivals are rejected. The interarrival times $A_i$ are negative exponentially distributed with rate $\lambda_i$ when the system is in state $[X=i]$ for $i=0,\dots,K$. The service time of a job follows an exponential distribution with rate $\mu$. There is a single waiting place.
The state-dependent arrival rates are denoted with M(x) - or sometimes M$_x$ - in Kendall's notation: M(x)/M/2-1
The system is a Markovian system. To be more precise, we have a birth-and-death process, since transitions occur only between neighboring states. The state of the system is the number $X$ of jobs in the system.
The transition rate $[X=i] \to [X=i+1]$ corresponds to the state-dependent arrival rate and we assume:
$\lambda_i = (i+1) \lambda$ for a given $\lambda$ and $i=0,\dots,K$.
The transition rate matrix is
$ Q = \left(\begin{array}{cccc} -\lambda & \lambda & 0 & 0 \\ \mu & -2\lambda-\mu & 2\lambda & 0 \\ 0 &2\mu & -3\lambda-2\mu & 3\lambda \\ 0 & 0 & 2\mu & -2\mu \end{array} \right) $
The state probabilites $X = \left( x(0), x(1), x(2), x(3) \right) $ are following from
$ X \cdot Q = 0$ and $X \cdot e = 1$
import numpy as np
from scipy import linalg
lam = 1
mu = 1
# transition rate matrix
Q = np.array([[-lam, lam, 0, 0],
[mu, -2*lam-mu, 2*lam, 0],
[0, 2*mu, -3*lam-2*mu, 3*lam],
[0, 0, 2*mu, -2*mu]])
print(f'Q=\n{Q}\n')
#%% change the matrix and define vector b
Q2 = Q.copy()
Q2[:, -1] = 1
print(f'Matrix is changed to\nQ2=\n{Q2}\n')
b = np.zeros(len(Q))
b[-1] = 1
print(f'b=\n{b}\n')
X = b @ linalg.inv(Q2) # compute the matrix inverse
print(f'X=\n{X*18}/18\n')
Q= [[-1 1 0 0] [ 1 -3 2 0] [ 0 2 -5 3] [ 0 0 2 -2]] Matrix is changed to Q2= [[-1 1 0 1] [ 1 -3 2 1] [ 0 2 -5 1] [ 0 0 2 1]] b= [0. 0. 0. 1.] X= [4. 4. 4. 6.]/18
The PASTA property must not be applied here, since the arrival process is not a Poisson process (although the interarrival times are exponentially distributed, but with state-dependent arrival rates).
As a consequence, we need to derive the state probability $x_A(i)$ that an arriving customer finds the system in state $[X_A=i]$. Then, the blocking probability is
$p_B = x_A(3)$.
The waiting probability is
$p_W = x_A(2)$.
To this end, we use the strong law of large numbers for Markov chains.
$ x_A(i) = \frac{\lambda_i \cdot x(i)}{\sum_{i=0}^K \lambda_i \cdot x(i)} $
lam_i = np.arange(1,5)
print(f'lambda_i = {lam_i}')
avg_lam = X @ lam_i
print(f'average lambda = {avg_lam*3:.2f}/3')
xA = lam_i/avg_lam*X
print(f'X=\n{X*18}/18\n')
print(f'X_A=\n{xA*12}/12\n')
lambda_i = [1 2 3 4] average lambda = 8.00/3 X= [4. 4. 4. 6.]/18 X_A= [1. 2. 3. 6.]/12