(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
The M/M/n loss system represents a fundamental queueing system. The steady state probabilities $x(i)$ for $i=0,...,n$ can be derived based on the analysis of Markov systems. The system consists of $n$ parallel servers, but no waiting room. Customers arrive according to a Poisson process with rate $\lambda$. The service time $B$ is exponentially distributed with mean $E[B]=1/\mu$. The M/M/n-0 system is called the Erlang loss model. The blocking probability $p_B = x(n)$ due to the PASTA property of the system, i.e., the probability of the state as seen by an outside random observer is the same as the probability of the state seen by an arriving customer: $x_A(i)=x(i)$. Hence, $p_B=x_A(n)=x(n)$.
The Erlang-B formula (or Erlang's loss formula) is
$p_B = \frac{\frac{a^n}{n!}}{\sum_{i=0}^n\frac{a^i}{i!}}$
with the offered load $a=\lambda E[B]=\frac{\lambda}{\mu}$. For the direct computation, the factorial needs to be computed which may lead to numerical issues.
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import factorial
def erlangB_directComputation(N, a=100):
def directComputation_singleValue(n,a): # computation for M/M/n system and single value of n
k = np.arange(n+1)
denominator = np.sum(a**k/factorial(k))
return a**n/factorial(n) / denominator
if np.isscalar(N):
pb = directComputation_singleValue(N,a)
else:
pb = np.zeros_like(N, dtype=float)
for i,n in enumerate(N):
pb[i] = directComputation_singleValue(n,a)
return pb
A convenient formula for the efficient computation and the avoidance of numerical issues (up to a certain point), is based on the recurrence formula.
$ \frac{1}{B(n,a)}=1+\frac{n}{a}\frac{1}{B(n-1,a)} $ with the initial condition $B(0,a)=1$.
This can also be expressed as
$ B(n,a) = \frac{aB(n-1,a)}{n+aB(n-1,a)} $
The blocking probability is then $p_B = B(n,a)$ for the M/M/n system with offered load $a=\lambda E[B]$.
#%% Iterative computation of Erlang-B for vector N
def erlangB_iterativeComputation(N, a=100):
def iterativeComputation_singleValue(n,a): # iterative computation for single value n (M/M/n)
InvB = 1.0
for j in range(1, n+1):
InvB = 1.0 + InvB * (j/a)
return (1.0 / InvB)
if np.isscalar(N):
pb = iterativeComputation_singleValue(N,a)
else:
pb = np.zeros_like(N, dtype=float)
for i,n in enumerate(N):
pb[i] = iterativeComputation_singleValue(n,a)
return pb
Changing the offered load will show numerical issues with the direct computation approach, while the recurrence formula is more robust.
a = 10.5
n = np.arange(1,251)
plt.plot(n, erlangB_directComputation(n, a=a), label=f'direct')
plt.plot(n, erlangB_iterativeComputation(n, a=a), label=f'iterative')
plt.yscale('log')
plt.xlabel('number n of servers')
plt.ylabel('blocking prob. pB')
plt.legend();
The economy of scale describes the phenomenon that if many server groups are merged together into a larger server group the resulting QoS - in terms of blocking probability - is increasing. A larger server group thus operates more economically. With the Erlang formula we can show this effect quantitatively.
Let us consider
Both systems have the same normalized offered traffic $a/n$.
n = np.arange(1,51)
single = erlangB_iterativeComputation(n, a=a)
double = erlangB_iterativeComputation(2*n, a=2*a)
plt.plot(n, single, label=f'M/M/n-0 with $\lambda$')
plt.plot(n, double, label=f'M/M/2n-0 with $2\lambda$')
plt.yscale('log')
plt.xlabel('normalized offered traffic a/n')
plt.ylabel('blocking prob. pB')
plt.grid(which='major')
plt.legend();
