LINE Solver (C++)
Templated C++ port of the LINE queueing solver
Loading...
Searching...
No Matches
line::smc Namespace Reference

Classes

struct  Drift
 The drift of an M/G/1-type sequence, and the invariant vector it uses. More...
struct  ShiftResult
 What MG1_Shifts returns: the shifted sequence and the drift it measured. More...
struct  Mg1CrOptions
 Options of MG1_CR, with the reference's defaults. More...
struct  Mg1FiOptions
 Options of MG1_FI, with the reference's defaults. More...

Typedefs

using Blocks = std::vector<Matrix<double>>

Functions

Matrix< double > mg1_g_etaqa (const Matrix< double > &A)
 G of an M/G/1-type chain, uniformized first.
std::vector< double > mg1_pi_etaqa (const Matrix< double > &Bin, const Matrix< double > &Ain, const Matrix< double > &G, const Matrix< double > &C0in=Matrix< double >())
 Aggregated stationary vector [pi0, pi1, pi2+pi3+...] of an M/G/1-type chain.
double mg1_qlen_etaqa (const Matrix< double > &Bin, const Matrix< double > &Ain, const std::vector< double > &pi, std::size_t n, const Matrix< double > &C0in=Matrix< double >())
 n-th moment of the level (the queue length) of an M/G/1-type chain from the ETAQA aggregates.
Matrix< double > gim1_r_etaqa (const Matrix< double > &A)
 R of a GI/M/1-type chain, uniformized first.
std::vector< double > gim1_pi_etaqa (const Matrix< double > &Bin, const Matrix< double > &Ain, const Matrix< double > &R, const Matrix< double > &B0in=Matrix< double >())
 Aggregated stationary vector [pi0, pi1, pi2+pi3+...] of a GI/M/1-type chain.
double gim1_qlen_etaqa (const Matrix< double > &Bin, const Matrix< double > &Ain, const Matrix< double > &R, const std::vector< double > &pi, std::size_t n, const Matrix< double > &B0in=Matrix< double >())
 n-th moment of the level of a GI/M/1-type chain from the ETAQA aggregates.
template<class T>
Matrix< T > madd (const Matrix< T > &A, const Matrix< T > &B)
 A + B.
template<class T>
Matrix< T > msub (const Matrix< T > &A, const Matrix< T > &B)
 A - B.
Matrix< double > mscale (const Matrix< double > &A, double c)
 c * A.
std::vector< double > rowsums (const Matrix< double > &A)
 sum(A,2), the row sums, as a column held in a vector.
double inf_norm (const Matrix< double > &A)
 norm(A,inf), the largest absolute row sum.
double inf_norm (const Blocks &blk, std::size_t from)
 norm(A,inf) over a whole block sequence stacked vertically.
double max_abs_diff (const Matrix< double > &A, const Matrix< double > &B)
 max(max(abs(A-B))).
double max_col_sum (const Matrix< double > &A)
 max(sum(A)), the largest column sum WITHOUT absolute values, as in MATLAB.
Blocks blocks_of (const Matrix< double > &A, std::size_t m)
 Splits the wide [A0 A1 ... Amax] into its m x m blocks.
Matrix< double > hcat (const Blocks &blk)
 Re-assembles a block sequence into the wide [A0 A1 ... Amax].
Matrix< double > vcat (const Blocks &blk)
 Stacks a block sequence vertically, [A0; A1; ...; Amax].
Blocks vblocks_of (const Matrix< double > &A, std::size_t r)
 Splits a vertical stack into its blocks of r rows.
std::vector< double > stat (const Matrix< double > &A)
 Stationary distribution of a stochastic matrix: the left eigenvector for eigenvalue 1, nonnegative and summing to one.
std::vector< double > rowvec_times (const std::vector< double > &v, const Matrix< double > &A)
 theta A, the row vector times matrix product used throughout.
double dot (const std::vector< double > &a, const std::vector< double > &b)
 The inner product of a row vector with a column held as a vector.
Drift mg1_drift (const Blocks &A)
 drift = theta * beta with beta = (Amax)e + (Amax+Amax-1)e + ..., the expected level increment per transition of the phase process.
std::complex< double > max_eig (const Matrix< double > &M)
 max(eig(M)) with MATLAB's semantics on a complex spectrum: the element of largest modulus, ties broken by the larger phase angle.
Matrix< double > poly_at (const Blocks &A, double z)
 A(z) = A0 + A1 z + ... + Amax z^max, by Horner as the reference writes it.
double mg1_decay (const Blocks &A)
 Decay rate of a recurrent M/G/1-type chain: the unique z > 1 with PF(A(z)) = z.
double gim1_caudal (const Blocks &A)
 Caudal characteristic of a GI/M/1-type chain: the spectral radius of R, the unique z in (0,1) with PF(A(z)) = z.
ShiftResult mg1_shifts (const Blocks &Ain, const std::string &shift_type)
 Shift technique for the M/G/1-type sequence.
Matrix< double > mg1_eg (const Blocks &Ain, bool &found)
 G in closed form when A0 has rank one.
Matrix< double > mg1_cr (const Blocks &Ain, const Mg1CrOptions &opts=Mg1CrOptions())
 Cyclic reduction for M/G/1-type Markov chains [Bini, Meini].
Matrix< double > mg1_fi (const Blocks &Ain, const Mg1FiOptions &opts=Mg1FiOptions())
 Functional iterations for M/G/1-type Markov chains [Neuts].
Matrix< double > gim1_r (const Blocks &Ain, const std::string &dual, const std::string &algor)
 R of a GI/M/1-type Markov chain, through the G of its dual.

Typedef Documentation

◆ Blocks

using line::smc::Blocks = std::vector<Matrix<double>>

Definition at line 75 of file mg1.h.

Function Documentation

◆ blocks_of()

Blocks line::smc::blocks_of ( const Matrix< double > & A,
std::size_t m )
inline

Splits the wide [A0 A1 ... Amax] into its m x m blocks.

Definition at line 158 of file mg1.h.

References blocks_of(), line::Matrix< T >::cols(), line::InputError::InputError(), and line::Matrix< T >::rows().

Referenced by blocks_of(), mg1_g_etaqa(), mg1_pi_etaqa(), and mg1_qlen_etaqa().

◆ dot()

double line::smc::dot ( const std::vector< double > & a,
const std::vector< double > & b )
inline

The inner product of a row vector with a column held as a vector.

Definition at line 236 of file mg1.h.

References dot(), and line::InputError::InputError().

Referenced by dot(), and mg1_drift().

◆ gim1_caudal()

double line::smc::gim1_caudal ( const Blocks & A)
inline

Caudal characteristic of a GI/M/1-type chain: the spectral radius of R, the unique z in (0,1) with PF(A(z)) = z.

Port of GIM1_Caudal.m.

Definition at line 329 of file mg1.h.

References gim1_caudal(), max_eig(), and poly_at().

Referenced by gim1_caudal(), gim1_r(), and mg1_eg().

◆ gim1_pi_etaqa()

std::vector< double > line::smc::gim1_pi_etaqa ( const Matrix< double > & Bin,
const Matrix< double > & Ain,
const Matrix< double > & R,
const Matrix< double > & B0in = Matrix<double>() )
inline

Aggregated stationary vector [pi0, pi1, pi2+pi3+...] of a GI/M/1-type chain.

Port of GIM1_pi_ETAQA.m. B and A are vertical stacks; B0 is the reference's 'Boundary' option (pass empty for the default A0).

Definition at line 559 of file etaqa.h.

References line::Matrix< T >::cols(), line::eye(), gim1_pi_etaqa(), line::InputError::InputError(), line::inverse(), madd(), line::matmul(), msub(), line::NumericError::NumericError(), line::Matrix< T >::rows(), rowsums(), and vblocks_of().

Referenced by gim1_pi_etaqa(), and line::mam::solver_mam_map_bmap_1().

◆ gim1_qlen_etaqa()

double line::smc::gim1_qlen_etaqa ( const Matrix< double > & Bin,
const Matrix< double > & Ain,
const Matrix< double > & R,
const std::vector< double > & pi,
std::size_t n,
const Matrix< double > & B0in = Matrix<double>() )
inline

n-th moment of the level of a GI/M/1-type chain from the ETAQA aggregates.

Port of GIM1_qlen_ETAQA.m, including the scalar A(3) of defect 3.

Definition at line 677 of file etaqa.h.

References line::Matrix< T >::cols(), line::eye(), gim1_qlen_etaqa(), line::InputError::InputError(), line::inverse(), madd(), line::matmul(), mscale(), msub(), line::NumericError::NumericError(), line::Matrix< T >::rows(), rowsums(), and vblocks_of().

Referenced by gim1_qlen_etaqa(), and line::mam::solver_mam_map_bmap_1().

◆ gim1_r()

Matrix< double > line::smc::gim1_r ( const Blocks & Ain,
const std::string & dual,
const std::string & algor )
inline

R of a GI/M/1-type Markov chain, through the G of its dual.

Port of GIM1_R.m for Dual 'A', 'R' and 'B' and Algor 'FI' and 'CR'.

THE DUAL IS THE WHOLE IDEA. There is no cyclic reduction for R directly, so the chain is transposed into an M/G/1-type one whose G carries the same information: the Ramaswami dual diag(theta)^-1 A_i' diag(theta) for a transient chain, and the Bright dual, which additionally rescales block i by eta^(i-1) with eta the caudal characteristic, for a positive recurrent one. 'A' picks between them by the drift, which is what makes it the fastest default. R is then read back off G by the inverse similarity, times eta in the Bright case.

Definition at line 838 of file mg1.h.

References gim1_caudal(), gim1_r(), line::InputError::InputError(), madd(), mg1_cr(), mg1_decay(), mg1_drift(), mg1_fi(), mscale(), stat(), line::smc::Drift::theta, line::UnsupportedError::UnsupportedError(), and line::smc::Drift::value.

Referenced by gim1_r(), and gim1_r_etaqa().

◆ gim1_r_etaqa()

Matrix< double > line::smc::gim1_r_etaqa ( const Matrix< double > & A)
inline

R of a GI/M/1-type chain, uniformized first.

Port of GIM1_R_ETAQA.m.

A is the VERTICAL stack [A0; A1; ...; Amax], which is how a GI/M/1-type sequence is written; the reference transposes that stack into the horizontal one GIM1_R wants, and asks for the automatic dual with functional iterations.

Definition at line 531 of file etaqa.h.

References line::Matrix< T >::cols(), gim1_r(), gim1_r_etaqa(), line::InputError::InputError(), line::Matrix< T >::rows(), and vblocks_of().

Referenced by gim1_r_etaqa(), and line::mam::solver_mam_map_bmap_1().

◆ hcat()

Matrix< double > line::smc::hcat ( const Blocks & blk)
inline

Re-assembles a block sequence into the wide [A0 A1 ... Amax].

Definition at line 170 of file mg1.h.

References line::Matrix< T >::cols(), line::Matrix< T >::empty(), hcat(), line::Matrix< T >::rows(), and line::Matrix< T >::size().

Referenced by hcat().

◆ inf_norm() [1/2]

double line::smc::inf_norm ( const Blocks & blk,
std::size_t from )
inline

norm(A,inf) over a whole block sequence stacked vertically.

Definition at line 131 of file mg1.h.

References inf_norm(), and line::Matrix< T >::size().

◆ inf_norm() [2/2]

double line::smc::inf_norm ( const Matrix< double > & A)
inline

norm(A,inf), the largest absolute row sum.

Definition at line 120 of file mg1.h.

References line::Matrix< T >::cols(), inf_norm(), and line::Matrix< T >::rows().

Referenced by inf_norm(), inf_norm(), mg1_cr(), and mg1_fi().

◆ madd()

template<class T>
Matrix< T > line::smc::madd ( const Matrix< T > & A,
const Matrix< T > & B )

A + B.

Templated so the point-wise CR step can add complex blocks.

Definition at line 83 of file mg1.h.

References line::Matrix< T >::cols(), line::InputError::InputError(), madd(), and line::Matrix< T >::rows().

Referenced by gim1_pi_etaqa(), gim1_qlen_etaqa(), gim1_r(), madd(), mg1_cr(), mg1_drift(), mg1_eg(), mg1_fi(), mg1_pi_etaqa(), mg1_qlen_etaqa(), and poly_at().

◆ max_abs_diff()

double line::smc::max_abs_diff ( const Matrix< double > & A,
const Matrix< double > & B )
inline

max(max(abs(A-B))).

Definition at line 138 of file mg1.h.

References line::Matrix< T >::cols(), max_abs_diff(), and line::Matrix< T >::rows().

Referenced by max_abs_diff(), and mg1_cr().

◆ max_col_sum()

double line::smc::max_col_sum ( const Matrix< double > & A)
inline

max(sum(A)), the largest column sum WITHOUT absolute values, as in MATLAB.

Definition at line 147 of file mg1.h.

References line::Matrix< T >::cols(), max_col_sum(), and line::Matrix< T >::rows().

Referenced by max_col_sum(), and mg1_cr().

◆ max_eig()

std::complex< double > line::smc::max_eig ( const Matrix< double > & M)
inline

max(eig(M)) with MATLAB's semantics on a complex spectrum: the element of largest modulus, ties broken by the larger phase angle.

For the nonnegative A(z) both callers evaluate, this is the Perron-Frobenius eigenvalue and is real; the comparisons below then take its real part, which is what MATLAB's relational operators do on a complex value.

Definition at line 281 of file mg1.h.

References line::eig_values(), max_eig(), and line::NumericError::NumericError().

Referenced by gim1_caudal(), max_eig(), and mg1_decay().

◆ mg1_cr()

Matrix< double > line::smc::mg1_cr ( const Blocks & Ain,
const Mg1CrOptions & opts = Mg1CrOptions() )
inline

Cyclic reduction for M/G/1-type Markov chains [Bini, Meini].

Port of MG1_CR.m, default mode 'ShiftPWCR' with ShiftType 'one'.

WHAT THE ALGORITHM DOES, since the transcription is otherwise opaque. One step of cyclic reduction eliminates every odd level of the chain and leaves a chain of the same M/G/1-type shape on the even ones, so the level distance halves per iteration and the iteration converges quadratically. Doing that on the block sequences directly is a polynomial composition; the reference instead evaluates the four sequences at the (nj+1)-th roots of unity, does the composition POINT-WISE (a small dense inverse per root), and interpolates back with an inverse transform – the "point-wise" in PWCR. The number of roots doubles until the interpolated tail is below (nj+1) eps, which is the reference's own accuracy control and the reason MaxNumRoot exists.

Everything runs on the TRANSPOSED blocks, as the reference does after D=D', and the final G is transposed back.

Definition at line 575 of file mg1.h.

References line::smc::ShiftResult::drift, line::eye(), line::smc::ShiftResult::hatA, inf_norm(), line::InputError::InputError(), line::inverse(), madd(), line::matmul(), max_abs_diff(), max_col_sum(), mg1_cr(), mg1_eg(), mg1_shifts(), msub(), line::svd_values(), line::Matrix< T >::transpose(), and line::UnsupportedError::UnsupportedError().

Referenced by gim1_r(), mg1_cr(), and mg1_g_etaqa().

◆ mg1_decay()

double line::smc::mg1_decay ( const Blocks & A)
inline

Decay rate of a recurrent M/G/1-type chain: the unique z > 1 with PF(A(z)) = z.

Port of MG1_Decay.m; the eigenvector output the reference offers is not returned, because the only caller that wants it is the 'tau' shift, which is refused (see the header note).

Definition at line 305 of file mg1.h.

References max_eig(), mg1_decay(), and poly_at().

Referenced by gim1_r(), and mg1_decay().

◆ mg1_drift()

Drift line::smc::mg1_drift ( const Blocks & A)
inline

drift = theta * beta with beta = (Amax)e + (Amax+Amax-1)e + ..., the expected level increment per transition of the phase process.

Repeated verbatim in MG1_EG, MG1_Shifts, MG1_pi_ETAQA and GIM1_R, so it lives here.

Definition at line 254 of file mg1.h.

References dot(), madd(), mg1_drift(), rowsums(), stat(), line::smc::Drift::theta, and line::smc::Drift::value.

Referenced by gim1_r(), mg1_drift(), mg1_eg(), mg1_pi_etaqa(), and mg1_shifts().

◆ mg1_eg()

Matrix< double > line::smc::mg1_eg ( const Blocks & Ain,
bool & found )
inline

G in closed form when A0 has rank one.

Port of MG1_EG.m.

found is false when the shortcut does not apply, which is the reference's empty return. A rank-one A0 means every down-transition forgets the phase it came from, so G is the same rank-one matrix e beta in the recurrent case; this is not an approximation and it is why an M/M/1-shaped input never enters cyclic reduction at all.

Definition at line 420 of file mg1.h.

References line::eye(), gim1_caudal(), line::inverse(), madd(), line::matmul(), line::matrix_rank(), mg1_drift(), mg1_eg(), mscale(), msub(), rowsums(), line::smc::Drift::theta, and line::smc::Drift::value.

Referenced by mg1_cr(), mg1_eg(), and mg1_fi().

◆ mg1_fi()

Matrix< double > line::smc::mg1_fi ( const Blocks & Ain,
const Mg1FiOptions & opts = Mg1FiOptions() )
inline

Functional iterations for M/G/1-type Markov chains [Neuts].

Port of MG1_FI.m for the three modes and the shift variants; the reference's NonZeroBlocks option is not exposed, because it changes only which products are skipped when some A_i vanish and converges to the same G.

'U-Based' is the default and the one GIM1_R(...,'FI') uses: it solves G = (I - sum_{j>=1} A_j G^{j-1})^{-1} A0, which is the U-based iteration and converges monotonically from below to the minimal nonnegative solution.

Definition at line 770 of file mg1.h.

References line::smc::ShiftResult::drift, line::eye(), line::smc::ShiftResult::hatA, inf_norm(), line::inverse(), madd(), line::matmul(), mg1_eg(), mg1_fi(), mg1_shifts(), msub(), and line::UnsupportedError::UnsupportedError().

Referenced by gim1_r(), and mg1_fi().

◆ mg1_g_etaqa()

Matrix< double > line::smc::mg1_g_etaqa ( const Matrix< double > & A)
inline

G of an M/G/1-type chain, uniformized first.

Port of MG1_G_ETAQA.m.

A is the wide [A0 A1 ... Amax]. The generator is turned into the transition matrix of the uniformized chain by dividing through by -min(diag(A1)) and adding the identity back onto A1, which is what cyclic reduction expects; see defect 1 in the header for why that branch is unconditional.

Definition at line 199 of file etaqa.h.

References blocks_of(), line::Matrix< T >::cols(), line::InputError::InputError(), mg1_cr(), mg1_g_etaqa(), and line::Matrix< T >::rows().

Referenced by mg1_g_etaqa(), and line::mam::solver_mam_bmap_map_1().

◆ mg1_pi_etaqa()

std::vector< double > line::smc::mg1_pi_etaqa ( const Matrix< double > & Bin,
const Matrix< double > & Ain,
const Matrix< double > & G,
const Matrix< double > & C0in = Matrix<double>() )
inline

Aggregated stationary vector [pi0, pi1, pi2+pi3+...] of an M/G/1-type chain.

Port of MG1_pi_ETAQA.m.

B may be empty, in which case the boundary repeats the repetitive blocks. C0 is the reference's 'Boundary' option, the block that takes level 1 back to a boundary of a different size; pass an empty matrix for the default A0.

Definition at line 230 of file etaqa.h.

References blocks_of(), line::Matrix< T >::cols(), line::Matrix< T >::empty(), line::InputError::InputError(), madd(), line::matmul(), line::matrix_rank(), mg1_drift(), mg1_pi_etaqa(), line::NumericError::NumericError(), line::Matrix< T >::rows(), rowsums(), and line::smc::Drift::value.

Referenced by mg1_pi_etaqa(), and line::mam::solver_mam_bmap_map_1().

◆ mg1_qlen_etaqa()

double line::smc::mg1_qlen_etaqa ( const Matrix< double > & Bin,
const Matrix< double > & Ain,
const std::vector< double > & pi,
std::size_t n,
const Matrix< double > & C0in = Matrix<double>() )
inline

n-th moment of the level (the queue length) of an M/G/1-type chain from the ETAQA aggregates.

Port of MG1_qlen_ETAQA.m.

Definition at line 382 of file etaqa.h.

References blocks_of(), line::Matrix< T >::cols(), line::Matrix< T >::empty(), line::InputError::InputError(), madd(), mg1_qlen_etaqa(), mscale(), and line::Matrix< T >::rows().

Referenced by mg1_qlen_etaqa(), and line::mam::solver_mam_bmap_map_1().

◆ mg1_shifts()

ShiftResult line::smc::mg1_shifts ( const Blocks & Ain,
const std::string & shift_type )
inline

Shift technique for the M/G/1-type sequence.

Port of MG1_Shifts.m, ShiftType 'one', which is the default and the only type ETAQA uses.

For a positive recurrent chain (drift < 1) the eigenvalue 1 of A(z) is shifted to zero by subtracting (A0+...+Ai)e u^T from block i with u^T = e^T/m; cyclic reduction then converges linearly in the SECOND largest root instead of stalling on the unit one, which is the entire point of running CR on the shifted sequence and undoing the shift on G afterwards.

'tau' and 'dbl', and 'one' at drift > 1, are refused: see the header.

Definition at line 367 of file mg1.h.

References line::smc::ShiftResult::drift, line::smc::ShiftResult::hatA, mg1_drift(), mg1_shifts(), rowsums(), line::UnsupportedError::UnsupportedError(), line::smc::ShiftResult::v, and line::smc::Drift::value.

Referenced by mg1_cr(), mg1_fi(), and mg1_shifts().

◆ mscale()

Matrix< double > line::smc::mscale ( const Matrix< double > & A,
double c )
inline

c * A.

Definition at line 104 of file mg1.h.

References line::Matrix< T >::cols(), mscale(), and line::Matrix< T >::rows().

Referenced by gim1_qlen_etaqa(), gim1_r(), mg1_eg(), mg1_qlen_etaqa(), mscale(), and poly_at().

◆ msub()

template<class T>
Matrix< T > line::smc::msub ( const Matrix< T > & A,
const Matrix< T > & B )

◆ poly_at()

Matrix< double > line::smc::poly_at ( const Blocks & A,
double z )
inline

A(z) = A0 + A1 z + ... + Amax z^max, by Horner as the reference writes it.

Definition at line 293 of file mg1.h.

References madd(), mscale(), and poly_at().

Referenced by gim1_caudal(), mg1_decay(), and poly_at().

◆ rowsums()

std::vector< double > line::smc::rowsums ( const Matrix< double > & A)
inline

sum(A,2), the row sums, as a column held in a vector.

Definition at line 112 of file mg1.h.

References line::Matrix< T >::cols(), line::Matrix< T >::rows(), and rowsums().

Referenced by gim1_pi_etaqa(), gim1_qlen_etaqa(), mg1_drift(), mg1_eg(), mg1_pi_etaqa(), mg1_shifts(), and rowsums().

◆ rowvec_times()

std::vector< double > line::smc::rowvec_times ( const std::vector< double > & v,
const Matrix< double > & A )
inline

theta A, the row vector times matrix product used throughout.

Definition at line 231 of file mg1.h.

References rowvec_times(), and line::vecmul().

Referenced by rowvec_times().

◆ stat()

std::vector< double > line::smc::stat ( const Matrix< double > & A)
inline

Stationary distribution of a stochastic matrix: the left eigenvector for eigenvalue 1, nonnegative and summing to one.

Port of stat.m, including its shape: [A - I, e] is S x (S+1), so the reference's y / B is a least-squares solve of a consistent overdetermined system and not a square solve. The normalization is IN the system (the appended column of ones against the appended 1 on the right), which is why the result needs no rescaling afterwards.

Definition at line 217 of file mg1.h.

References line::Matrix< T >::cols(), line::InputError::InputError(), line::lstsq(), line::Matrix< T >::rows(), and stat().

Referenced by gim1_r(), mg1_drift(), and stat().

◆ vblocks_of()

Blocks line::smc::vblocks_of ( const Matrix< double > & A,
std::size_t r )
inline

Splits a vertical stack into its blocks of r rows.

Definition at line 192 of file mg1.h.

References line::Matrix< T >::cols(), line::InputError::InputError(), line::Matrix< T >::rows(), and vblocks_of().

Referenced by gim1_pi_etaqa(), gim1_qlen_etaqa(), gim1_r_etaqa(), and vblocks_of().

◆ vcat()

Matrix< double > line::smc::vcat ( const Blocks & blk)
inline

Stacks a block sequence vertically, [A0; A1; ...; Amax].

Definition at line 181 of file mg1.h.

References line::Matrix< T >::cols(), line::Matrix< T >::empty(), line::Matrix< T >::rows(), line::Matrix< T >::size(), and vcat().

Referenced by vcat().