First- and second-order level-dependent (multi-regime) Markovian fluid queues: the matrix-exponential building blocks of the stationary law. More...

#include <cmath>
#include <cstddef>
#include <vector>
#include "line/api/mam/mfq_ld_distr.h"
#include "line/api/mam/mfq_ld_mean.h"
#include "line/api/mam/mfq_solve.h"
#include "line/api/mam/qbd_r.h"
#include "line/num/number.h"
#include "line/util/error.h"
#include "line/util/expm.h"
#include "line/util/linalg.h"
#include "line/util/lu.h"
#include "line/util/matrix.h"

Include dependency graph for mfq_ld_solve.h:

Namespaces
namespace	line
	Conservation laws of a layered queueing network, enumerated from its structure.
namespace	line::mam

Enumerations
enum class	line::mam::FluidBoundary { line::mam::Reflective = 0 , line::mam::Absorbing = 1 }
	Boundary behaviour of one background state at a reflecting level. More...

Functions
template<class T>
LevelDependentFluidBlocks< T >	line::mam::mfq_ld_solve (const std::vector< Matrix< T > > &Q, const std::vector< Matrix< T > > &R, const std::vector< Matrix< T > > &S, const std::vector< T > &Thr, const std::vector< FluidBoundary > &boundaryL, const std::vector< FluidBoundary > &boundaryU, const std::vector< Matrix< T > > &Qt, const T &prec)
	Solve a first- or second-order level-dependent fluid queue.
template<class T>
LevelDependentFluidBlocks< T >	line::mam::mfq_ld_solve (const std::vector< Matrix< T > > &Q, const std::vector< Matrix< T > > &R, const std::vector< Matrix< T > > &S, const std::vector< T > &Thr)
	mfq_ld_solve with reflective boundaries, Qt = Q and the default prec = 1e-14.

Detailed Description

First- and second-order level-dependent (multi-regime) Markovian fluid queues: the matrix-exponential building blocks of the stationary law.

Port of matlab/src/api/mam/mfq_ld_solve.m and the BUTools SecondOrderLevelDependentFluidSolve it wraps. The generator, the drift and optionally the VARIANCE change at threshold levels, giving a piecewise homogeneous first- or second-order (Brownian) fluid queue. Setting every variance cell to zero reduces it to first order. The blocks it returns are exactly what mfq_ld_mean and mfq_ld_distr consume, so this closes that family rather than extending it.

METHOD. Per regime, the states with neither drift nor variance are censored out, and the remainder is split into three classes: positive drift with no variance, negative drift with no variance, and any state with variance. The density is anchored at both ends of the regime, forward from its lower threshold and backward from its upper one, and each direction's exponent comes from a QBD. The second-order terms are what force the QBD form: with a variance the balance equation is second order in the level, and multiplying through by a constant c chosen from the spectrum turns it into the matrix QUADRATIC of a discrete QBD, from whose R the exponent is recovered as K = (R - I) c. c is the smallest value that keeps the transformed triple substochastic, taken as the maximum over the drift states of -Q_ii/R_ii and over the variance states of the larger root of the discriminant, floored at 1.

WHICH QBD SOLVER, AND WHY NO NEW ONE WAS NEEDED. The reference calls QBD_CR(Bm, Lm, Fm), and the R it returns satisfies

Fm + R Lm + R^2 Bm = 0,

which is EXACTLY the contract of the port's existing qbd_R, F + R L + R^2 B = 0, under the direct mapping B = Bm, L = Lm, F = Fm. So the cyclic-reduction routine did not have to be transcribed.

That equation was established by MEASUREMENT: do not read QBD_CR's "(sub)stochastic" error message as asking for the discrete-time form R = A0 + R A1 + R^2 A2. The triples this caller builds sum to a GENERATOR, and the shift that reading implies (L = Lm - I) gives a wrong K. The test asserts the residual above on the port's own R.

The boundary system then couples the K+1 point masses to the 2K density initial vectors through flux conservation at every threshold, the boundary conditions that reflective and absorbing states impose, and continuity of the second-order density across a threshold. One normalization row replaces the first flux equation.

The boundary flags boundaryL/boundaryU are ONE FLAG PER BACKGROUND STATE, a length-N vector shared by every regime, defaulting to reflective. A K x N default makes the logical mask larger than the array it indexes and crashes every model with K >= 2.

ARITHMETIC. Templated on T and gated on num_traits<T>::has_transcendental: the QBD iteration is tolerance-terminated and expm is a Pade approximation. As in mfq_ld_distr, the ONLY eigenvalue computation is the branch test inside the normalization's integral of a matrix exponential, which selects between two algebraically equivalent formulas and never supplies a value that reaches the result, so Real instantiation is honest here. Contrast mfq_multiregime, where the Schur factors ARE the basis of the answer and double is the ceiling.

Definition in file mfq_ld_solve.h.

Namespaces

Enumerations

Functions

Detailed Description