LINE Solver (C++)
Templated C++ port of the LINE queueing solver
Loading...
Searching...
No Matches
pfqn_mci.h File Reference

Monte Carlo Integration estimate of the normalizing constant of a closed product-form network (Ross, Wang and Yao; MonteQueue 2.0). More...

#include <cmath>
#include <cstddef>
#include <limits>
#include <vector>
#include "line/api/pfqn/pfqn_bs.h"
#include "line/api/pfqn/pfqn_ca.h"
#include "line/api/pfqn/pfqn_mc_common.h"
#include "line/num/number.h"
#include "line/util/error.h"
#include "line/util/matrix.h"
Include dependency graph for pfqn_mci.h:

Go to the source code of this file.

Namespaces

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

Enumerations

enum class  line::pfqn::MciVariant {
  line::pfqn::Imci , line::pfqn::Mci , line::pfqn::Rm , line::pfqn::Amci ,
  line::pfqn::Lhsmci
}
 The proposal-rate rules the reference selects between. More...

Functions

template<class T>
NcResult< T > line::pfqn::pfqn_mci (const Matrix< T > &D, const std::vector< int > &N, const std::vector< T > &Z, std::size_t samples, MciVariant variant, McRng &rng)
 Monte Carlo Integration estimate of the normalizing constant of a closed product-form network (Ross, Wang and Yao; MonteQueue 2.0).
template<class T>
NcResult< T > line::pfqn::pfqn_mci (const Matrix< T > &D, const std::vector< int > &N, const std::vector< T > &Z, McRng &rng)
 Reference defaults: 1e5 samples, the IMCI proposal.

Detailed Description

Monte Carlo Integration estimate of the normalizing constant of a closed product-form network (Ross, Wang and Yao; MonteQueue 2.0).

Templated port of matlab/src/api/pfqn/pfqn_mci.m, cross-checked against jar/src/main/java/jline/api/pfqn/nc/Pfqn_mci.java and python/line_solver/api/pfqn/asymptotic.py. Variants: mci, imci, amci, lhsmci, rm.

The McKenna-Mitra integral form of the constant,

G(N) = 1/prod_r N_r! int_{R_+^M} e^{-sum_i v_i} prod_r (sum_i v_i D(i,r) + Z_r)^{N_r} dv,

is estimated by importance sampling with an independent exponential proposal of rate gamma_i per station. Writing V ~ prod_i Exp(gamma_i), one draw contributes

lZ = -sum_i (1 - gamma_i) V_i - sum_i log gamma_i - sum_r log N_r!

  • sum_r N_r log(sum_i V_i D(i,r) + Z_r),

and lG = log mean exp(lZ). The proposal rates come from a Bard-Schweitzer pre-solve: gamma_i = max(0.01, 1 - U_i) in the IMCI variant, and the saturation-aware gamma_i = 1/sqrt(max N) when U_i > 0.9 in the plain MCI variant. Overshooting rates make the estimator heavy-tailed, which is why the reference clamps them.

The deviates come from the caller's generator, one draw per station per sample. (MATLAB used to cache its uniforms in a persistent across calls; that cache was removed on 2026-07-28, so MATLAB now redraws per call as this port, the JAR and python do.)

Variant notes:

  • The 'rm' variant reads the balanced-job-bound rate with the COLUMN maximum of the demands, tput_r = N_r / (sum_i D(i,r) + Z_r + max_i D(i,r) (sum N - 1)), for every M, as MATLAB (max(D,[],1)), the JAR and python do. Before 2026-10-21 MATLAB wrote max(D,1), the elementwise maximum with the scalar 1, which for M > 1 built the rates from the FIRST station's clamped demands only, and the JAR copied those numbers. Only the proposal, hence the variance, was affected: every reading targets the same constant.
  • 'amci' and 'lhsmci' use the 'imci' tilt and differ only in the sampler. 'amci' draws antithetic pairs (u, 1-u) in rows s and ceil(I/2)+s, and 'lhsmci' stratifies each station's uniforms by a Latin hypercube. Both fill the uniforms STATION-MAJOR (all draws for station 1, then station 2, ...), in the order of MATLAB's column-major rand(Ih,M) and of the JAR; the hypercube draws its Fisher-Yates permutation (as the JAR does, by mc_uniform_int) before the I jitters of that station. The plain samplers ('mci', 'imci', 'rm') keep the sample-major order this port has always used, so their seeded outputs are unchanged.

Arithmetic: INEXACT BY CONSTRUCTION. The value is a random variable, the integrand is evaluated through logs, and the proposal rates come from a tolerance-stopped Bard-Schweitzer solve; all three require transcendental arithmetic.

RNG contract: see pfqn_mc_common.h. Comparable to MATLAB, the JAR and python only in distribution, never stream for stream (MATLAB's rand/randperm, java.util.Random and numpy's legacy generator are three different streams); reproducible within this port only when the generator is passed in the same state. The cross-codebase check is therefore statistical: the mean of independent replications must sit within its own standard error of the exact constant (line-test.git cpp/tests/test_pfqn_mci_variants.cpp).

Definition in file pfqn_mci.h.