Exact normalizing constant of a closed product-form network whose stations may be load dependent (generalized Buzen, Reiser-Kobayashi 1975). More...

#include <algorithm>
#include <cmath>
#include <cstddef>
#include <limits>
#include <type_traits>
#include <vector>
#include "line/api/pfqn/pfqn_ca.h"
#include "line/num/number.h"
#include "line/util/error.h"
#include "line/util/matrix.h"
#include "line/util/population.h"

Include dependency graph for pfqn_gld.h:

Namespaces
namespace	line
namespace	line::pfqn

Functions
template<class T>
NcResult< T >	line::pfqn::pfqn_gld (const Matrix< T > &L, const std::vector< int > &N, const Matrix< T > &mu)
	Exact normalizing constant of a closed product-form network whose stations may be load dependent (generalized Buzen, Reiser-Kobayashi 1975).
template<class T>
NcResult< T >	line::pfqn::pfqn_gld (const Matrix< T > &L, const std::vector< int > &N)
	Overload with all rates equal to one, i.e.

Detailed Description

Exact normalizing constant of a closed product-form network whose stations may be load dependent (generalized Buzen, Reiser-Kobayashi 1975).

Templated port of matlab/src/api/pfqn/pfqn_gld.m (and of the single-class specialization matlab/src/api/pfqn/pfqn_gldsingle.m, which the recursion below subsumes), cross-checked term for term against mp_pfqn's gld/gld_multi.c, the exact GMP reference.

Model. Station i serves at rate mu(i,k) when it holds k jobs, k = 1 ... Nt with Nt = sum_r N_r. The single-station balance function for a class vector k with j = sum_r k_r jobs is

Y_i(k) = j! / prod_r k_r! * prod_r L(i,r)^{k_r} / prod_{a=1}^{j} mu(i,a)

and the constant is the convolution of the M station factors,

G_0(n) = [n == 0], G_i(n) = sum_{0 <= k <= n} Y_i(k) G_{i-1}(n - k), G(N) = G_M(N).

Algorithm. MATLAB writes this as a recursion on (M, N, mu) with no memoization, whose cost is exponential in Nt; mp_pfqn replaced it by the station-by-station convolution above, which is what this port implements. The value is identical, the cost is O(M P^2 R) with P = prod_r (N_r + 1). A station whose rates are all 1 is load independent and is folded in by the classical in-place Buzen update

G_i(n) = G_{i-1}(n) + sum_r L(i,r) G_i(n - e_r)

in O(P R) instead, so a model with no load-dependent station reduces operation for operation to pfqn_ca on the same demands.

NO pfqn_lld HERE, deliberately. MATLAB, python and the JAR carry a pfqn_lld alongside pfqn_gld: there the recursion is the unmemoised one described above, and saturating the rate shift at the LLD threshold s_k makes its state repeat, so a memo turns an exponential tree into a bounded one. That is a cure for an algorithm this port does not use. The convolution here is already polynomial and visits each station once, so there is nothing for the threshold to collapse; the LLD structure would have to be exploited by a different device, splitting a station's balance function into the s_k terms below the threshold and a geometric tail folded in by the Buzen update, and that is a change of algorithm rather than a port. pfqn_lldsingle IS ported, because the single-class kernel it accelerates is the same recursion in every language.

Delay stations. MATLAB's pfqn_gld takes no think-time argument: a delay is an ordinary row of L whose rates are mu(i,k) = k, for which the factorials cancel and Y_i(k) collapses to prod_r Z_r^{k_r} / prod_r k_r!. The port keeps that convention, so an infinite-server station is expressed by giving it the rate row 1, 2, ..., Nt and needs no special case anywhere.

Arithmetic. Every operation is an addition, a multiplication or a division in the field of the inputs, so the algorithm is exact in rational arithmetic with no reformulation; nothing here needs a transcendental function. As in pfqn_ca, IEEE double is the only arithmetic that can overflow, and it gets the same power-of-two rescaling of the demands: dividing every demand by 2^k divides every Y_i(k) of total degree j by 2^{jk}, hence divides G(N) by exactly 2^{Nt k}, and ldexp moves the exponent without touching a mantissa bit. The estimate that picks k accounts for the rates as well as the demands, since a delay station depresses G by Nt!.

Definition in file pfqn_gld.h.

Namespaces

Functions

Detailed Description