Buzen's algorithm: Difference between revisions

Content deleted Content added
Adjusted spacing to improve readability. Minor wording improvements to improve clarity.
m Made the mathematical aspects of the algorithm more accessible to non-specialists.
Line 1:
In [[queueing theory]], a discipline within the mathematical [[probability theory|theory of probability]], '''Buzen's algorithm''' (or '''convolution algorithm''') is an algorithm for calculating the [[normalization constant]] G(''N'') in the [[Gordon–Newell theorem]]. This method was first proposed by [[Jeffrey P. Buzen]] in his 1971 PhD dissertation<ref name=":0">{{Cite book |last=Defense Technical Information Center |url=http://archive.org/details/DTIC_AD0731575 |title=DTIC AD0731575: Queueing Network Models of Multiprogramming |date=1971-08-01 |language=english}}</ref> and subsequently published in a refereed journal in 1973.<ref name="buzen-1973">{{Cite journal | last1 = Buzen | first1 = J. P. | author-link = Jeffrey P. Buzen| title = Computational algorithms for closed queueing networks with exponential servers | doi = 10.1145/362342.362345 | url = http://www-unix.ecs.umass.edu/~krishna/ece673/buzen.pdf| journal = Communications of the ACM | volume = 16 | issue = 9 | pages = 527–531 | year = 1973 | s2cid = 10702 }}</ref> Computing G(''N'') is required to compute the stationary [[probability distribution]] of a closed queueing network.<ref name="gn">{{Cite journal | last1 = Gordon | first1 = W. J. | last2 = Newell | first2 = G. F. | author-link2 = Gordon F. Newell| doi = 10.1287/opre.15.2.254 | jstor = 168557| title = Closed Queuing Systems with Exponential Servers | journal = [[Operations Research (journal)|Operations Research]]| volume = 15 | issue = 2 | pages = 254 | year = 1967 }}</ref>
 
Performing a naïve computation of the normalizing constant requires enumeration of all states. For a closed network with ''N'' circulating jobscustomers and ''M'' service centersfacilities, G(''N'') is the sum of <math>\tbinom{N+M-1}{M-1}</math> individual terms, with each term consisting of ''M'' factors raised to powers whose total sum is ''N''. Buzen's algorithm computes G(''N'') using a total of ''NM'' multiplications and ''NM'' additions. This is a highly significantdramatic improvement that opened the door to applying the Gordon-Newell theorem to models of realistic size.<ref name=":1">{{cite journal |last1=Denning |first1=Peter J. |date=24 August 2016 |title=Rethinking Randomness: An interview with Jeff Buzen, Part I |url=https://dl.acm.org/doi/10.1145/2986329 |journal=Ubiquity |volume=2016 |issue=August |pages=1:1–1:17 |doi=10.1145/2986329}}</ref> The values of G(1), G(2) ... G(''N'' -1), which can be used to expresscalculate other important quantities of interest, are computed as by-products of the algorithm.
 
==Problem setup==
Line 7:
Consider a closed queueing network with ''M'' service facilities and ''N'' circulating customers. Assume that the service time for a customer at service facility ''i'' is given by an [[exponentially distributed]] random variable with parameter ''μ''<sub>''i''</sub> and that, after completing service at service facility ''i'', a customer will proceed next to service facility ''j'' with probability ''p''<sub>''ij''</sub>.<ref name="gn" />
 
whereLet <math>\mathbb P(n_1,n_2,\cdots,n_M) </math> isbe the steady state probability that the number of customers at service facility ''i'' is equal to ''n<sub>i</sub>'' for i = 1, 2, ... , ''M .'' It follows from the [[Gordon–Newell theorem]] that
It follows from the [[Gordon–Newell theorem]] that the equilibrium distribution of this model is
 
<math>\mathbb P(n_1,n_2,\cdots,n_M) = \frac{1}{\text{G}(N)}</math><math> \prod_left( X_1 \right)^{i=1n_1}</math><math> \left( X_2 \right)^M{n_2}</math> .... <math> \left( X_iX_M \right)^{n_in_M}</math>
 
This result is usually written more compactly as
where <math>\mathbb P(n_1,n_2,\cdots,n_M) </math> is the steady state probability that the number of customers at service facility ''i'' is equal to ''n<sub>i</sub>'' for i = 1, 2, ... , ''M''
 
<math>\mathbb P(n_1,n_2,\cdots,n_M) = \frac{1}{\text{G}(N)}\prod_{i=1}^M \left( X_i \right)^{n_i}</math>
 
and where theThe values of ''X''<sub>''i''</sub> are determined by solving
 
<math>\mu_j X_j = \sum_{i=1}^M \mu_i X_i p_{ij}\quad\text{ for }j=1,\ldots,M.</math>
Line 23 ⟶ 25:
The individual terms that must be added together to compute G(''N'') all have the following form:
 
<math> \left( X_1 \right)^{n_1}</math><math> \left( X_2 \right)^{n_2}</math><math> \left( X_3 \right)^{n_3}</math> .... <math> \left( X_M \right)^{n_M}</math>. Note that this set of terms can be partitioned into two groups. The first group comprises all terms for which the exponent of <math> \left( X_M \right)</math> is greater than or equal to 1.  This implies that <math> \left( X_M \right)</math><sup>'' raised to the power 1''</sup> can be factored out of each of these terms.  
 
After factoring out <math> \left( X_M \right)</math>, a surprising result emerges: the modified terms in the first group are identical to the terms used to compute the normalizing constant for the same network with one customer removed.<ref name=":1" />  Thus, the sum of the terms in the first group can be written as “''X''<sub>''M''</sub> times G(''N'' -1)”.
 
Next consider the second group.  The exponent of ''X''<sub>''M''</sub> for every term in this group is zero.  InAs effecta result, service facility ''M'' effectively disappears from all terms in this group (since it reduces in every case to a factor of 1). This leaves the total number of customers at the remaining ''M'' -1 service facilities equal to ''N''.
 
To express this resultobservation mathematically, assume that X<sub>1</sub>, X<sub>2</sub>, … ''X''<sub>''M''</sub> have been obtained for a given network with ''M'' service facilities. For any n ≤ ''N'' and m ≤ ''M,'' define g(n,m) as the normalizing constant for a network with n customers, m service facilities (1,2, … m), and values of  X<sub>1</sub>, X<sub>2</sub>, … ''X''<sub>''m''</sub>  that match the first m members of the original sequence X<sub>1</sub>, X<sub>2</sub>, … ''X''<sub>''M''</sub>  .
 
Given this definition, the normalizingsum constantof G(''N'')the terms in the Gordon-Newellsecond theoremgroup can now be re-written as g(''N'', ''M'' -1).
 
It also follows immediately that “X''<sub>M</sub>'' times G(''N'' -1)”, the sum of the terms in the first group, can be re-written as “X''<sub>M</sub>'' times g(''N'' -1,''M'' )”.   More importantly, the sum of the terms in the second group can now be written as g(''N'', ''M'' -1).
 
In addition, the normalizing constant G(''N'') in the Gordon-Newell theorem can now be re-written as g(''N'',''M''). Since G(''N'') is equal to the combined sum of the terms in the first and second groups,
 
G(''N'') = g(''N'', ''M'' ) = X<sub>M</sub> g(''N'' -1,''M'' ) + g(''N'',''M'' -1)
Line 41 ⟶ 43:
This same recurrence relation clearly exists for any intermediate value of n from 1 to ''N'', and for any intermediate value of m from 1 to ''M'' .  
 
This implies g(n,m) = X<sub>m</sub> g(n-1,m) + g(n,m-1).  Buzen’s algorithm is simply the iterative application of this fundamental recurrence relation, along with the following boundary conditions.
Buzen’s algorithm is simply the iterative application of this fundamental recurrence relation, along with the following boundary conditions.
 
g(0,m) = 1 for m = 1, 2, …''M''
Line 50 ⟶ 51:
==Marginal distributions, expected number of customers==
 
The Gordon-Newell theorem enables analysts to determine the stationary probability associated with each individual state of a closed queueing network.  These individual probabilities must then be added together to evaluate marginalother important probabilities. suchFor asexample P(n<sub>i</sub>≥k), the probability that the total number of customers at service center i is greater than or equal to k, must be (summed over all values of n<sub>i</sub>≥k and, for each value of n<sub>i</sub>, all possible ways the remaining N – n<sub>i</sub> customers can be distributed across the other M-1 service centers in the network).
 
Many of these marginal probabilities can be computed with minimal additional effort.  This is easy to see for the case of P(n<sub>i</sub>≥k).   Clearly, X<sub>i</sub> must be raised to the power of k or higher in every state where the number of customers at service center i is greater than or equal to k. Thus (X<sub>i</sub>)<sup>k</sup> can be factored out from each of these probabilities, leaving a set of modified probabilities whose sum is given by G(N-k)/G(N).   This observation yields the following simple and highly efficient result: