Queues
Closed-form answers for waiting lines — M/M/1, Erlang B and C, bounded queues, Pollaczek–Khinchine and Jackson networks.
Prelude.Queues answers the standard questions about a queue in equilibrium — how long the line is, how long a customer waits, how often somebody is turned away — with the exact formulas rather than a simulation. A queue in steady state is described completely by a handful of ratios, and a simulation of thousands of customers can only approach them. Where you do simulate, these definitions are the numbers to check the simulation against.
Little’s law, \(L = \lambda W\), holds for every queue in equilibrium whatever the arrival process, service distribution or number of servers. The module states it as a residual (little) so the formulas can be checked against each other.
Conventions: arrival is the arrival rate \(\lambda\) and service the service rate \(\mu\) of one server, both in customers per unit time. The Erlang and bounded-queue formulas are sums over \(0, \dots, c\); since the length of a vector in Hilbert comes from a tensor rather than a number, you pass that range as a vector (seats, states) whose length sets the number of terms and whose contents are ignored. Denominators are guarded against exact zero: every formula below divides by \(\operatorname{nonzero}(d) = d + 10^{-300}\,[\,\lvert d \rvert < 10^{-300}\,]\) (see Chains), which equals \(d\) everywhere except within \(10^{-300}\) of zero. So a formula used outside its range (for example an unbounded queue with \(\rho \ge 1\)) returns a meaningless number rather than an error — check stability yourself.
open Prelude.Queues
M/M/1
One server, Poisson arrivals, exponential service. Every quantity is a function of the utilisation \(\rho = \lambda / \mu\), which is why a queue at 95% utilisation is twenty deep and one at 99% is a hundred deep. All formulas require \(\rho < 1\).
offered
offered(arrival, service) · arrival arrival rate scalar · service service rate scalar → utilisation scalar
Let \(\lambda \ge 0\) (the argument arrival) be the arrival rate and \(\mu > 0\) (the argument service) the service rate. Then
The utilisation (offered load per server): the fraction of time the server is busy.
Intuition Work arriving per unit time divided by the work the server can finish per unit time. It is the share of time the server is busy, and every M/M/1 answer depends on it alone. Delay grows like \(1/(1-\rho)\), so the queue explodes as \(\rho\) approaches 1.
Use it when you first size any service: a checkout, a CPU, an API endpoint, a support agent. Check \(\rho < 1\) before using the other formulas, which return meaningless numbers otherwise. Leave a margin: waiting at 95% utilisation is four times waiting at 80%. With \(c\) servers the same ratio is the load in erlangs that erlangC takes, and the utilisation per server is that divided by \(c\).
capacity-planningcongestion-analysisstability-check
def offered(arrival, service) ≔ arrival / nonzero(service)
queueLength
queueLength(arrival, service) · arrival arrival rate scalar · service service rate scalar → mean number in the system scalar
Let \(\lambda \ge 0\) (the argument arrival) and \(\mu > 0\) (the argument service) be the arrival and service rates, and let \(\rho = \operatorname{offered}(\lambda, \mu) < 1\). Then
The mean number of customers in the system, waiting or in service.
Intuition The average number present, counting the one in service. It stays small at moderate load but explodes near saturation: 1 at 50%, 4 at 80%, 19 at 95%. Random bursts build a backlog that a server with little slack takes a long time to clear.
Use it when you estimate buffer or memory needs, work in progress at a production step, or connections held open. It needs \(\rho < 1\) and assumes Poisson arrivals and exponential service. For other service-time variability use pkQueued and add \(\rho\) for the customer in service. For several servers use erlangQueue, and for a finite buffer use boundedLength.
capacity-planningcongestion-analysis
def queueLength(arrival, service) ≔ offered(arrival, service) / nonzero(1 − offered(arrival, service))
queueWaiting
queueWaiting(arrival, service) · arrival arrival rate scalar · service service rate scalar → mean time in the system scalar
Let \(\lambda \ge 0\) (the argument arrival) and \(\mu > \lambda\) (the argument service) be the arrival and service rates. Then
The mean time a customer spends in the system, including their own service.
Intuition The total time in the system, waiting plus service. What matters is the slack \(\mu - \lambda\), the capacity left over after the arrivals, and the mean time is its reciprocal: halve the slack and the time doubles.
Use it when you need end-to-end response time, such as request latency, time in a store or repair turnaround, to compare against a service-level target. It is a mean. In M/M/1 the time in the system is exponential, so the 95th percentile is about three times this value (\(\ln 20\)). It needs \(\mu > \lambda\). For the wait before service only, use queueDelay.
waiting-timesservice-levelscapacity-planning
def queueWaiting(arrival, service) ≔ 1 / nonzero(service − arrival)
queueDelay
queueDelay(arrival, service) · arrival arrival rate scalar · service service rate scalar → mean time in the queue scalar
Let \(\lambda \ge 0\) (the argument arrival) and \(\mu > \lambda\) (the argument service) be the arrival and service rates, and let \(\rho = \operatorname{offered}(\lambda, \mu)\). Then
The mean time spent waiting before service starts.
Intuition The time in line before service starts: an arrival waits for everyone already present, on average \(\rho/(1-\rho)\) customers at \(1/\mu\) each. It equals queueWaiting minus one mean service time.
Use it when the service itself is fine and the complaint is the wait: hold time in a call queue, time before a job starts running, time before a patient is seen. It fits targets such as “average wait under two minutes”. It needs \(\mu > \lambda\). For several servers use erlangWait, and for non-exponential service use pkDelay.
waiting-timesservice-levels
def queueDelay(arrival, service) ≔ offered(arrival, service) / nonzero(service − arrival)
queueQueued
queueQueued(arrival, service) · arrival arrival rate scalar · service service rate scalar → mean number waiting scalar
Let \(\lambda \ge 0\) (the argument arrival) and \(\mu > 0\) (the argument service) be the arrival and service rates, and let \(\rho = \operatorname{offered}(\lambda, \mu) < 1\). Then
The mean number of customers waiting, not counting the one in service.
Intuition Customers standing in line, not counting the one in service. The server is busy a share \(\rho\) of the time, so on average \(\rho\) customers are in service, and \(L - \rho\) simplifies to \(\rho L\). Near saturation almost everyone present is waiting.
Use it when you size the waiting space (chairs in a waiting room, a queue buffer, a holding area) or the backlog shown on a dashboard. It needs \(\rho < 1\). By Little’s law it equals \(\lambda\) times queueDelay. For non-exponential service use pkQueued, and for \(c\) servers use erlangQueue.
capacity-planningcongestion-analysis
def queueQueued(arrival, service) ≔ offered(arrival, service) · queueLength(arrival, service)
queueAt
queueAt(arrival, service, n) · arrival arrival rate scalar · service service rate scalar · n number of customers scalar or tensor → probability same shape as n
Let \(\lambda \ge 0\) (the argument arrival) and \(\mu > 0\) (the argument service) be the arrival and service rates, let \(\rho = \operatorname{offered}(\lambda, \mu) < 1\), and let \(n \in \mathbb{N}\) (or a tensor of counts, taken entrywise). Then
The stationary probability of finding exactly \(n\) customers in the system — a geometric distribution. Pass a vector of counts to get the whole distribution at once.
Intuition Each extra customer present is \(\rho\) times as likely as one fewer, so the count is geometric: the system is empty with probability \(1 - \rho\), and the tail grows heavy as \(\rho\) nears 1.
Use it when you need tails rather than averages: the chance that more than \(n\) are present (\(\rho^{n+1}\)), a buffer size that overflows only rarely, or the idle share (\(n = 0\)). Pass a vector of counts to get the whole distribution. It needs \(\rho < 1\). For a real buffer of limited size, boundedAt models the limit exactly instead of assuming infinite room.
capacity-planningservice-levelsrisk-measurement
def queueAt(arrival, service, n) ≔ (1 − offered(arrival, service)) · offered(arrival, service) ^ n
little
little(length, arrival, waiting) · length mean number in the system scalar · arrival arrival rate scalar · waiting mean time in the system scalar → residual scalar
Let \(L \ge 0\) (the argument length) be a mean number of customers, \(\lambda \ge 0\) (the argument arrival) the arrival rate, and \(W \ge 0\) (the argument waiting) a mean time. Then
Little’s law as a residual. It is zero for any queue in equilibrium, so it is a check that a length and a waiting time were measured on the same thing — both on the system, or both on the queue — rather than one of each.
Intuition Count customer-minutes two ways. Over a long period, the average number present equals the arrival rate times how long each customer stays. Nothing about arrivals or service enters, so it holds for any stable queue.
Use it when you check consistency: L, λ and W reported by a simulation or a monitoring system should leave a residual near zero, and a large residual often means a queue length was paired with a system time. Rearranged, it gives one quantity from the other two, for example mean latency from average in-flight requests and throughput. When customers are turned away, use the accepted rate (throughput) as \(\lambda\).
model-validationwaiting-times
def little(length, arrival, waiting) ≔ length − arrival · waiting
M/M/c and the Erlang formulas
With \(c\) servers the offered load is measured in erlangs, \(A = \lambda / \mu\). Erlang B and Erlang C differ by one term and are easily confused. Erlang B is the probability of being turned away when there is no room to wait (M/M/c/c, a telephone trunk). Erlang C is the probability of having to wait when the queue is unbounded (M/M/c). A capacity plan built on the wrong one is wrong by a factor that grows with the load.
\[B(A, c) = \frac{A^c / c!}{\sum_{k=0}^{c} A^k / k!} \qquad C(A, c) = \frac{c\,B}{c - A + A\,B}\]The terms are computed in logarithms (via the log-gamma function), so loads and server counts in the hundreds do not overflow. seats must be a vector of \(c + 1\) entries, its contents ignored.
erlangB
erlangB(a, seats) · a offered load in erlangs scalar · seats any vector of c + 1 entries → blocking probability scalar
Let \(A > 0\) (the argument a) be the offered load, and let \(s \in \mathbb{R}^{c+1}\) (the argument seats) be any vector, its positions \(0, 1, \dots, c\) indexing the sum. Let \(t_k = \operatorname{erlangTerms}(A, s)_k = \exp\big(k \ln A - \operatorname{lgamma}(k + 1)\big) \approx A^k / k!\) for \(k = 0, \dots, c\). Then
The Erlang B formula: the probability that an arriving customer finds all \(c\) servers busy and is lost, in a system with no waiting room. The approximation sign is only that of lgamma, the Lanczos approximation of \(\ln \Gamma\).
Intuition \(c\) lines and no waiting room: a call that arrives when every line is busy is lost. The terms \(A^k/k!\) are the relative chances of \(k\) lines being busy, and \(B\) is the share of the last term, all lines busy. Poisson arrivals see the system as it is at a random moment, so that is also the share of calls blocked.
Use it when you size systems that turn customers away instead of queueing them: telephone trunks, hospital beds without a waiting list, licence pools, parking spaces, connection limits. Find the smallest \(c\) whose blocking is under your target by lengthening seats. It depends on service times only through their mean, and it is valid at any load. If blocked customers wait instead, use erlangC.
capacity-planningloss-systemsservice-levels
def erlangB(a, seats) ≔ lastOf(erlangTerms(a, seats)) / nonzero(total(erlangTerms(a, seats)))
erlangC
erlangC(a, servers, seats) · a offered load in erlangs scalar · servers number of servers c scalar · seats any vector of c + 1 entries → probability of waiting scalar
Let \(A > 0\) (the argument a) be the offered load, \(c \in \mathbb{N}\) (the argument servers) the number of servers with \(A < c\), and \(s \in \mathbb{R}^{c+1}\) (the argument seats) any vector of \(c + 1\) entries. Let \(B = \operatorname{erlangB}(A, s)\). Then
The Erlang C formula: the probability that an arriving customer has to wait in an M/M/c queue. servers and the length of seats must agree.
Intuition The same \(c\) servers, but a customer who finds them all busy waits in line. Waiting customers stay in the system instead of vanishing, so “all busy” is more likely than in the loss system: \(C\) is never below \(B\), and it is much higher at heavy load.
Use it when you staff call centres, support desks, worker pools or checkout lanes where customers queue. The usual plan: try candidate \(c\) values and compute the service level \(P(\text{wait} \le T) = 1 - C\,e^{-\mu(c - A)T}\). It needs \(A < c\) and, unlike Erlang B, assumes exponential service. For the mean wait use erlangWait, and for the mean queue use erlangQueue.
capacity-planningservice-levelswaiting-times
def erlangC(a, servers, seats) ≔ servers · erlangB(a, seats) / nonzero(servers − a + a · erlangB(a, seats))
erlangWait
erlangWait(a, servers, seats, service) · a offered load in erlangs scalar · servers c scalar · seats any vector of c + 1 entries · service service rate of one server scalar → mean wait in the queue scalar
Let \(A > 0\) (the argument a), \(c \in \mathbb{N}\) (the argument servers) with \(A < c\), \(s \in \mathbb{R}^{c+1}\) (the argument seats), and \(\mu > 0\) (the argument service) the service rate of one server. Let \(C = \operatorname{erlangC}(A, c, s)\). Then
The mean time a customer waits before service in M/M/c: Erlang C divided by the spare capacity.
Intuition Only a share \(C\) of customers wait at all, and those who do wait \(1/(\mu(c - A))\) on average, the reciprocal of the spare capacity. The product is the mean wait over everybody.
Use it when you report average speed of answer for a call centre or queueing delay for a worker pool, or ask whether one more agent is worth it. The spare-capacity term makes that extra agent worth most near saturation. An average over everybody hides the unlucky, so also quote the wait of those who do wait and a tail probability from erlangC. It needs \(A < c\).
waiting-timesservice-levelscapacity-planning
def erlangWait(a, servers, seats, service) ≔ erlangC(a, servers, seats) / nonzero(service · (servers − a))
erlangQueue
erlangQueue(a, servers, seats) · a offered load in erlangs scalar · servers c scalar · seats any vector of c + 1 entries → mean number waiting scalar
Let \(A > 0\) (the argument a), \(c \in \mathbb{N}\) (the argument servers) with \(A < c\), and \(s \in \mathbb{R}^{c+1}\) (the argument seats). Let \(C = \operatorname{erlangC}(A, c, s)\). Then
The mean number of customers waiting in M/M/c — Little’s law applied to the queue.
Intuition Little’s law applied to the line alone: arrival rate times mean wait. Since \(\lambda/\mu = A\), that is \(C\,A/(c - A)\), which grows sharply as the load approaches the number of servers.
Use it when you size the waiting area, the hold lines or the buffer of a multi-server system, or report its average backlog. It needs \(A < c\) and assumes exponential service. For the average number in the whole system, add the customers in service: erlangQueue plus busyServers.
capacity-planningcongestion-analysis
def erlangQueue(a, servers, seats) ≔ erlangC(a, servers, seats) · a / nonzero(servers − a)
busyServers
busyServers(a) · a offered load in erlangs scalar → mean number of busy servers scalar
Let \(A \ge 0\) (the argument a) be the offered load in erlangs. Then
In an M/M/c queue that loses nobody, the mean number of busy servers equals the offered load.
Intuition Each arrival holds a server for \(1/\mu\) on average, and \(\lambda\) arrive per unit time, so on average \(\lambda/\mu\) servers are occupied at once. A load in erlangs is literally a count of busy servers.
Use it when you estimate average utilisation of a pool (\(A/c\) per server), the staff time actually spent serving, or the total number in the system (add erlangQueue). It holds only when nobody is turned away. In a loss system, multiply by \(1 - B\) from erlangB to get the load actually carried.
capacity-planningcongestion-analysis
def busyServers(a) ≔ a
Building blocks
erlangLogs(a, seats)— Let \(A > 0\) and \(s \in \mathbb{R}^{c+1}\), with positions \(k = 0, 1, \dots, c\). Then \(\operatorname{erlangLogs}(A, s)_k = k \ln \operatorname{nonzero}(A) - \operatorname{lgamma}(k + 1) \approx \ln(A^k / k!)\). The logarithm of each term, computed with log-gamma.erlangTerms(a, seats)— Let \(A > 0\) and \(s \in \mathbb{R}^{c+1}\). Then \(\operatorname{erlangTerms}(A, s)_k = \exp\big(\operatorname{erlangLogs}(A, s)_k\big) \approx A^k / k!\) for \(k = 0, \dots, c\). The terms themselves; used byerlangB.
M/M/1/K
A single-server queue with room for at most \(K\) customers; arrivals that find it full are turned away. Unlike the unbounded queue it has a stationary distribution at every utilisation, including \(\rho \ge 1\), which makes \(K\) the parameter a real system has. Here rho is \(\lambda/\mu\) and states is any vector of \(K + 1\) entries (contents ignored).
boundedAt
boundedAt(rho, states) · rho utilisation scalar · states any vector of K + 1 entries → distribution vector K + 1
Let \(\rho > 0\) (the argument rho) be the utilisation and \(s \in \mathbb{R}^{K+1}\) (the argument states) any vector, its positions \(0, 1, \dots, K\) indexing the occupancies. Then, for \(n = 0, \dots, K\),
The stationary probability of each occupancy.
Intuition One server and room for at most \(K\) customers in total. Each extra customer is still \(\rho\) times as likely as one fewer, but the list stops at \(K\), so the probabilities always normalise, even when arrivals outpace service. Overloaded, the mass simply piles up at “full”.
Use it when you model a real finite buffer: a network interface queue, a small shop, a workstation with limited space for work in progress, a phone line with a few hold slots. This matters most under overload, where the unbounded formulas fail. Pass states with \(K + 1\) entries. For rates that change with occupancy use birthDeathStationary in Jumps. For the unbounded idealisation use queueAt.
capacity-planningloss-systemscongestion-analysis
def boundedAt(rho, states) ≔ boundedTerms(rho, states) / nonzero(total(boundedTerms(rho, states)))
blocking
blocking(rho, states) · rho utilisation scalar · states any vector of K + 1 entries → blocking probability scalar
Let \(\rho > 0\) (the argument rho) and \(s \in \mathbb{R}^{K+1}\) (the argument states), and let \(\pi = \operatorname{boundedAt}(\rho, s)\). Then
The probability that the system is full — by Poisson arrivals seeing time averages, the fraction of arrivals turned away.
Intuition The chance the room is full at a random moment. Poisson arrivals see the same picture as a random moment, so this is also the share of arrivals turned away.
Use it when you choose a buffer size for a target loss rate: dropped packets, callers who get a busy signal, customers who leave a full shop. Lengthen states until the blocking is acceptable. At \(\rho > 1\) it tends to \(1 - 1/\rho\) as \(K\) grows, so extra room cannot make up for too little capacity. For \(c\) servers with no waiting room, use erlangB.
loss-systemscapacity-planningservice-levels
def blocking(rho, states) ≔ lastOf(boundedAt(rho, states))
throughput
throughput(arrival, rho, states) · arrival arrival rate scalar · rho utilisation scalar · states any vector of K + 1 entries → accepted rate scalar
Let \(\lambda \ge 0\) (the argument arrival), \(\rho > 0\) (the argument rho) and \(s \in \mathbb{R}^{K+1}\) (the argument states), and let \(\pi_K = \operatorname{blocking}(\rho, s)\). Then
The rate at which customers are actually admitted.
Intuition Customers offered at rate \(\lambda\), minus those turned away. In steady state what goes in comes out, so this is also the rate of completed services.
Use it when you compute revenue or goodput of a finite system, the real load passed to the next stage of a pipeline, or the right \(\lambda\) for Little’s law: the mean time in the system of an admitted customer is boundedLength divided by the throughput. Nothing checks that arrival and rho agree, so pass \(\rho = \lambda/\mu\) from the same rates.
loss-systemscapacity-planningwaiting-times
def throughput(arrival, rho, states) ≔ arrival · (1 − blocking(rho, states))
boundedLength
boundedLength(rho, states) · rho utilisation scalar · states any vector of K + 1 entries → mean number in the system scalar
Let \(\rho > 0\) (the argument rho) and \(s \in \mathbb{R}^{K+1}\) (the argument states), and let \(\pi = \operatorname{boundedAt}(\rho, s)\). Then
The mean occupancy.
Intuition The average occupancy of the finite room. At low load it matches the unbounded queue. Near and above \(\rho = 1\) it stays bounded instead of blowing up: it is exactly \(K/2\) at \(\rho = 1\) and heads toward \(K\) under heavy overload.
Use it when you size or monitor the typical fill of a finite buffer: memory in use, work in progress, occupied places. Combine it with throughput through Little’s law for the mean time an admitted customer spends inside. For the unbounded model use queueLength.
congestion-analysiscapacity-planningwaiting-times
def boundedLength(rho, states) ≔ einsum("n,n->", boundedAt(rho, states), positions(states))
Building blocks
boundedTerms(rho, states)— Let \(\rho > 0\) and \(s \in \mathbb{R}^{K+1}\), with positions \(n = 0, 1, \dots, K\). Then \(\operatorname{boundedTerms}(\rho, s)_n = \rho^{n}\). The unnormalised weights; used byboundedAt.
M/G/1 by Pollaczek–Khinchine
With Poisson arrivals and a general service-time distribution, the service distribution enters only through its coefficient of variation \(C_s = \sigma / m\). Two systems with the same mean service time but different variability queue differently, and the whole difference is one factor: constant service (\(C_s = 0\)) halves the queue that exponential service (\(C_s = 1\)) produces.
\[L_q = \frac{\rho^2 \,(1 + C_s^2)}{2(1-\rho)}\]pkQueued
pkQueued(rho, mean, spread) · rho utilisation scalar · mean mean service time scalar · spread standard deviation of the service time scalar → mean number waiting scalar
Let \(0 \le \rho < 1\) (the argument rho) be the utilisation, \(m > 0\) (the argument mean) the mean service time and \(\sigma \ge 0\) (the argument spread) its standard deviation, and let \(C_s = \operatorname{spreadRatio}(m, \sigma) = \sigma / \operatorname{nonzero}(m)\). Then
The Pollaczek–Khinchine mean queue length.
Intuition Queues are driven by variability, not only by load. With Poisson arrivals and any service-time distribution, the line is the exponential-service line scaled by \((1 + C_s^2)/2\). Perfectly regular service halves it, and service more variable than exponential (\(C_s > 1\)) makes it longer.
Use it when service times are not exponential and you know their mean and standard deviation: fixed-length tasks, machine cycles, disk operations, medical procedures. It also shows whether standardising a process is worth the effort. It needs \(\rho < 1\), Poisson arrivals and a single server. For the mean wait use pkDelay. The multi-server Erlang formulas here assume exponential service.
congestion-analysiscapacity-planning
def pkQueued(rho, mean, spread) ≔ rho · rho · (1 + spreadRatio(mean, spread) · spreadRatio(mean, spread)) / nonzero(2 · (1 − rho))
pkDelay
pkDelay(rho, mean, spread, arrival) · rho utilisation scalar · mean mean service time scalar · spread standard deviation of the service time scalar · arrival arrival rate scalar → mean wait in the queue scalar
Let \(0 \le \rho < 1\) (the argument rho), \(m > 0\) (the argument mean), \(\sigma \ge 0\) (the argument spread) and \(\lambda > 0\) (the argument arrival), and let \(L_q = \operatorname{pkQueued}(\rho, m, \sigma)\). Then
The mean time spent waiting, by Little’s law on the queue.
Intuition Little’s law on the line: the mean number waiting divided by the arrival rate. The wait carries the same factor \((1 + C_s^2)/2\) as the queue, so at a fixed load and mean, cutting the service-time spread cuts the wait in proportion.
Use it when you estimate the wait before service at a single machine, clinic slot or single-threaded worker whose service-time mean and spread you have measured. Nothing checks that rho equals arrival times mean, so compute them from the same data. Add mean to get the total time in the system.
waiting-timesservice-levels
def pkDelay(rho, mean, spread, arrival) ≔ pkQueued(rho, mean, spread) / nonzero(arrival)
Building blocks
spreadRatio(mean, spread)— Let \(m > 0\) and \(\sigma \ge 0\). Then \(\operatorname{spreadRatio}(m, \sigma) = C_s = \sigma / \operatorname{nonzero}(m)\). The coefficient of variation of the service time.
Jackson networks
A Jackson network is a set of M/M/1 nodes that pass customers to each other with fixed routing probabilities. The total arrival rate at each node solves the traffic equations \(\lambda = \gamma + R^\top \lambda\), where \(\gamma\) is the external arrival rate at each node and \(R_{ij}\) is the probability that a customer leaving node \(i\) goes next to node \(j\) (rows may sum to less than one; the remainder leaves the network).
The product-form theorem is what makes the network tractable: in equilibrium the nodes behave as independent M/M/1 queues, even though they feed one another, so each node’s length is simply \(\rho_i/(1-\rho_i)\). Because solve requires a symmetric positive-definite matrix and \(I - R^\top\) is not symmetric, the traffic equations are solved through their normal equations. Every node must have \(\rho_i < 1\).
arrivalsOf
arrivalsOf(external, routing) · external external arrival rates vector n · routing routing probabilities matrix n×n → total arrival rate at each node vector n
Let \(\gamma \in \mathbb{R}^n\) (the argument external) be the external arrival rates and \(R \in [0, 1]^{n \times n}\) (the argument routing) the routing matrix, with \(I_n - R^{\top}\) non-singular. Let \(M = I_n - R^{\top}\). Then
the solution of the traffic equations \((I_n - R^\top)\lambda = \gamma\), found by conjugate gradient on the normal equations.
Intuition Flow balance: what enters a node is its outside arrivals plus everything routed to it from other nodes. Feedback such as rework or retries makes the internal flows larger than the outside arrivals, and one linear solve counts every lap at once.
Use it when you work out the true load at each stage of a routed process, such as a microservice call graph, patient flow through hospital departments or a factory line with rework, before sizing any node. Rows of routing must sum to at most 1, and every customer must eventually leave. solve runs conjugate gradient on the normal equations, so for badly conditioned routing check that the result is non-negative and satisfies the traffic equations.
network-performancecapacity-planning
def arrivalsOf(external, routing) ≔ solve(normalOf(identity(routing) − einsum("ij->ji", routing)), rhsOf(identity(routing) − einsum("ij->ji", routing), external))
loads
loads(external, routing, service) · external external arrival rates vector n · routing routing probabilities matrix n×n · service service rate of each node vector n → utilisation of each node vector n
Let \(\gamma \in \mathbb{R}^n\) (the argument external), \(R \in [0, 1]^{n \times n}\) (the argument routing) and \(\mu \in \mathbb{R}^n\) with \(\mu_i > 0\) (the argument service), and let \(\lambda = \operatorname{arrivalsOf}(\gamma, R)\). Then, for \(i = 1, \dots, n\),
The utilisation of each node.
Intuition Each node’s utilisation from its total arrivals, including routed and fed-back traffic, not just the arrivals from outside. A node that looks lightly used from outside can be the bottleneck once rework comes back to it.
Use it when you look for bottlenecks in a network of service stages. Any \(\rho_i \ge 1\) makes that node unstable and makes nodeLengths and networkTotal meaningless, so add capacity there first. Each node is a single server with the given service rate.
network-performancecapacity-planningstability-check
def loads(external, routing, service) ≔ arrivalsOf(external, routing) / nonzero(service)
nodeLengths
nodeLengths(external, routing, service) · external external arrival rates vector n · routing routing probabilities matrix n×n · service service rate of each node vector n → mean number at each node vector n
Let \(\gamma \in \mathbb{R}^n\) (the argument external), \(R \in [0, 1]^{n \times n}\) (the argument routing) and \(\mu \in \mathbb{R}^n\) (the argument service), and let \(\rho = \operatorname{loads}(\gamma, R, \mu)\) with every \(\rho_i < 1\). Then, for \(i = 1, \dots, n\),
The mean number of customers at each node.
Intuition The product-form theorem: in equilibrium each node behaves like an independent single-server queue at its own load, even though the nodes feed each other. So each length is simply \(\rho_i/(1-\rho_i)\).
Use it when you need the backlog at each stage of a network: work in progress per workstation, messages waiting per service. Divide by the node’s total arrival rate from arrivalsOf to get the mean time per visit (Little’s law). It assumes Poisson outside arrivals, exponential service, one server per node and random routing. With batched or fixed-schedule flows, treat it as an approximation.
network-performancecongestion-analysis
def nodeLengths(external, routing, service) ≔ loads(external, routing, service) / nonzero(1 − loads(external, routing, service))
networkTotal
networkTotal(external, routing, service) · external external arrival rates vector n · routing routing probabilities matrix n×n · service service rate of each node vector n → mean number in the network scalar
Let \(\gamma \in \mathbb{R}^n\) (the argument external), \(R \in [0, 1]^{n \times n}\) (the argument routing) and \(\mu \in \mathbb{R}^n\) (the argument service), and let \(\rho = \operatorname{loads}(\gamma, R, \mu)\). Then
The mean number of customers in the whole network, the sum of the node lengths. Divide by the total external arrival rate to get the mean time in the network by Little’s law.
Intuition The number of customers in the whole network is the sum of the node lengths. Divided by the total outside arrival rate, it gives the mean time from entry to exit, revisits included.
Use it when you need end-to-end latency or total work in progress of a multi-stage process: request time through a service mesh, a patient’s length of stay, order lead time. Compare designs, such as capacity added at one node or another, by how much this number drops. Every \(\rho_i\) must be below 1, so check loads first.
network-performancewaiting-timescapacity-planning
def networkTotal(external, routing, service) ≔ total(nodeLengths(external, routing, service))
Examples
The examples below live in one file, Waiting.hb, so the class is Waiting.
Length and waiting time of a single-server queue, checked by Little’s law
Two customers per minute arrive at a server that handles three per minute, so \(\rho = 2/3\): on average two customers are in the system and each spends one minute there.
open Prelude.Queues
fn singleServer(arrival, service) ↦ queueLength(arrival, service)
fn singleWait(arrival, service) ↦ queueWaiting(arrival, service)
fn littleCheck(arrival, service) ↦ little(queueLength(arrival, service), arrival, queueWaiting(arrival, service))
using Examples;
var length = Waiting.SingleServer(2.0, 3.0);
var wait = Waiting.SingleWait(2.0, 3.0);
var check = Waiting.LittleCheck(2.0, 3.0);
Output (verified):
L = 2
W = 1
little = -0
Probability of waiting at a call centre (Erlang C) and of being lost on a trunk (Erlang B)
A load of 2 erlangs on 3 agents. seats is any vector of \(c + 1 = 4\) entries; only its length matters.
open Prelude.Queues
fn callCentreWait(a, servers, seats) ↦ erlangC(a, servers, seats)
fn trunkBlocked(a, seats) ↦ erlangB(a, seats)
using Examples;
var seats = new[] { 0d, 0d, 0d, 0d }; // c + 1 = 4 entries
var waiting = Waiting.CallCentreWait(2.0, 3.0, seats);
var lost = Waiting.TrunkBlocked(2.0, seats);
Output (verified, rounded to six digits):
erlangC = 0.444444
erlangB = 0.210526
By hand: the terms are \(1, 2, 2, 4/3\), so \(B = (4/3)/(19/3) = 4/19\) and \(C = 3B/(1 + 2B) = 4/9\).
How often a small waiting room turns people away
A room that holds at most two customers (\(K = 2\), so three states) at \(\rho = 1\): every occupancy is equally likely, a third of arrivals are turned away, and the mean occupancy is one.
open Prelude.Queues
fn roomTurnedAway(rho, states) ↦ blocking(rho, states)
fn roomLength(rho, states) ↦ boundedLength(rho, states)
using Examples;
var states = new[] { 0d, 0d, 0d }; // K + 1 = 3 entries
var blocked = Waiting.RoomTurnedAway(1.0, states);
var length = Waiting.RoomLength(1.0, states);
Output (verified, rounded to six digits):
blocking = 0.333333
boundedLength = 1
Constant versus exponential service time
At \(\rho = 0.5\) with a mean service time of 1, constant service (\(\sigma = 0\)) gives half the queue of exponential service (\(\sigma = 1\)).
open Prelude.Queues
fn constantServiceQueue(rho, mean, spread) ↦ pkQueued(rho, mean, spread)
using Examples;
var constant = Waiting.ConstantServiceQueue(0.5, 1.0, 0.0);
var exponential = Waiting.ConstantServiceQueue(0.5, 1.0, 1.0);
Output (verified):
pk const = 0.25
pk exp = 0.5
Two queues in tandem
Customers arrive at node 0 at rate 1 and all go on to node 1. Node 0 serves at rate 2 and node 1 at rate 4, so the loads are 0.5 and 0.25 and the mean lengths are 1 and 1/3.
open Prelude.Queues
fn networkLengths(external, routing, service) ↦ nodeLengths(external, routing, service)
using Examples;
var external = new[] { 1d, 0d };
var routing = new double[,] { { 0, 1 }, { 0, 0 } };
var service = new[] { 2d, 4d };
var lengths = Waiting.NetworkLengths(external, routing, service);
Output (verified, rounded to six digits):
network = [1, 0.333333]