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

Shared machinery of the QRF nonlinear bounds (qrf_noblo_*, qrf_bas_*). More...

#include <algorithm>
#include <cmath>
#include <cstddef>
#include <string>
#include <vector>
#include "line/num/number.h"
#include "line/util/error.h"
#include "line/util/matrix.h"
#include "line/util/lp_highs.h"
#include "line/util/simplex.h"
Include dependency graph for mapqn_qrf_common.h:

Go to the source code of this file.

Classes

struct  line::mapqn::QrfVars< T >
 The unflattened decision vector: the pair tensor and the effective rates. More...
struct  line::mapqn::QrfMetrics< T >
 The utilizations and queue lengths read off an optimal pair tensor. More...
struct  line::mapqn::QrfAffine< T >
 An affine residual map recovered as (A, b) with fn(x) = A x - b. More...
struct  line::mapqn::QrfReduced< T >
 The equality block with its dependent rows dropped. More...

Namespaces

namespace  line
namespace  line::mapqn

Functions

template<class T>
line::mapqn::qrf_log (const T &v)
 log() in the working arithmetic.
template<class T>
line::mapqn::qrf_logtol ()
 The reference's LOGTOL: the shift that keeps log() off zero.
std::size_t line::mapqn::qrf_num_vars (std::size_t M, std::size_t N, const std::vector< int > &K, std::size_t MR)
 Number of decision variables the layout actually USES.
std::vector< long > line::mapqn::qrf_index_map (std::size_t M, std::size_t N, const std::vector< int > &K, std::size_t MR)
 Flat position of every p2 entry, in the FILL ORDER of sub_qrfvar.
template<class T>
QrfVars< T > line::mapqn::sub_qrfvar (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, std::size_t MR)
 Unflatten x into the pair tensor and the effective rates.
template<class T>
line::mapqn::mmi_objective (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR)
 Mutual-information objective.
template<class T>
line::mapqn::mem_objective (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR)
 Maximum-entropy objective, returned as the NEGATIVE entropy +sum p log p over the diagonal entries, because the solver MINIMIZES and the AMPL model states this objective as maximize H.
template<class T>
std::vector< T > line::mapqn::mmi_gradient (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR, const std::vector< long > &idx)
 Gradient of mmi_objective.
template<class T>
std::vector< T > line::mapqn::mem_gradient (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR, const std::vector< long > &idx)
 Gradient of mem_objective: d/dp of p log(p') is log p' + p/p'.
template<class T>
line::mapqn::bethe_objective (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR)
 Tree-reweighted (Bethe) free entropy at the uniform spanning-tree weight, the objective of qrf.bethe.
template<class T>
std::vector< T > line::mapqn::bethe_gradient (const std::vector< T > &x, std::size_t M, std::size_t N, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR, const std::vector< long > &idx)
 Gradient of bethe_objective.
template<class T>
QrfMetrics< T > line::mapqn::qrf_extract_results (const QrfVars< T > &v, std::size_t M, const std::vector< int > &K, const std::vector< int > &F, std::size_t MR, const Matrix< T > *alpha=nullptr)
 extract_results: the diagonal marginals of the optimal tensor, plus the alpha-weighted mean BN.
template<class T>
std::vector< std::size_t > line::mapqn::qrf_independent_rows (const Matrix< T > &A, double tol=-1.0)
 Rows of a maximal linearly independent subset of A, by pivoted Gram-Schmidt.
template<class T, class Fn>
QrfAffine< T > line::mapqn::qrf_affine_matrices (Fn fn, std::size_t n)
 Recover (A, b) from an affine residual map.
template<class T>
QrfReduced< T > line::mapqn::qrf_reduce_equalities (const Matrix< T > &A, const std::vector< T > &b)
 Drop the linearly dependent equality rows, keeping the feasible set exact.
template<class T>
lp::LpModel< T > line::mapqn::qrf_polytope (const Matrix< T > &Aeq, const std::vector< T > &beq, const Matrix< T > &Aub, const std::vector< T > &bub, std::size_t n)
 The polytope of a QRF instance, as an LpModel over the box [0,1]^n.
template<class T>
Matrix< T > line::mapqn::qrf_null_space (const Matrix< T > &A, std::size_t n)
 Orthonormal basis of null(A), as an (n x d) matrix.
template<class T, class Obj, class Grad>
std::vector< T > line::mapqn::solve_qrf_nlp_lp (Obj objective, Grad gradient, const std::vector< T > &x0, const lp::LpModel< T > &polytope, const std::string &name, unsigned max_iter=200, double gap_tol=1e-10)
 Minimize a convex objective over {Aeq x = beq, Aub x <= bub, 0 <= x <= 1}, starting from a feasible point, by FRANK-WOLFE.
template<class T, class Obj, class Grad>
std::vector< T > line::mapqn::solve_qrf_nlp (Obj objective, Grad gradient, const std::vector< T > &x0, const Matrix< T > &Aeq, const std::vector< T > &beq, const Matrix< T > &Aub, const std::vector< T > &bub, const std::string &name, unsigned max_iter=200, double gap_tol=1e-10)
 The same, with the polytope given as matrices rather than as an LpModel.
template<class T>
std::vector< T > line::mapqn::qrf_min_norm_point (const lp::LpModel< T > &polytope, const std::vector< T > &x0, const std::string &name)
 The minimum-norm point of the polytope, from any feasible point of it.
template<class T>
std::vector< T > line::mapqn::qrf_feasible_start_lp (const lp::LpModel< T > &polytope, const std::string &name)
 A feasible point of an LpModel polytope: the MINIMUM-NORM one.
template<class T>
std::vector< T > line::mapqn::qrf_feasible_start (const Matrix< T > &Aeq, const std::vector< T > &beq, const Matrix< T > &Aub, const std::vector< T > &bub, std::size_t n)
 A point of the polytope, as the phase 1 of qrf_noblo_start.m.

Detailed Description

Shared machinery of the QRF nonlinear bounds (qrf_noblo_*, qrf_bas_*).

Port of python/line_solver/api/mapqn/qrf_noblo_common.py (the JAR twin is Mapqn_qrf_noblo_* plus Mapqn_nlp_solver). api/mapqn has NO MATLAB implementation, so Python and the JAR are the references here.

WHAT THE FAMILY IS. The linear QRF bounds (mapqn_qr_bounds_*) MAXIMIZE one station's utilization over a polytope of pairwise queue-phase probabilities. These entry points optimize a different objective over the SAME feasible set: mutual information (MMI) or negative entropy (MEM).

NEITHER OBJECTIVE IS A CONVEX PROGRAM AS THE REFERENCE STATES IT, and a caller has to know which one it asked for:

  • MEM minimizes -sum p log p over the diagonal entries. Despite the name, that MINIMIZES the entropy, and -p log p is CONCAVE, so the minimum sits at a vertex and every method reports a stationary point fixed by its start. The three codebases agree because they start from the same place: the MINIMUM-NORM feasible point, not an arbitrary phase-1 vertex. See qrf_feasible_start_lp, which is where that is arranged; measured on a two-station cycle at N = 2 with both one and two phases.
  • MMI minimizes sum_{i != j} p_ij (log p_ij - log p_ii - log p_jj). The -p_ij log p_ii terms are NOT convex, so the problem is nonconvex and both codebases report a LOCAL optimum whose identity is fixed by the start point, i.e. by whichever vertex the phase-1 LP happens to return. Measured: on a one-phase pair the two agree exactly (0.857143 / 0.428571, itself the exact product-form answer); on a two-phase / one-phase pair they land on different vertices, 0.6931 against 0.6438 in objective, and NEITHER is wrong – the reference's own solve does not move from its start either. This is verifiable rather than suspected: the reference's start point is FEASIBLE under this port's constraints to 8.9e-16, and from this port's start it is an ASCENT direction (a directional derivative of +0.184), which only a nonconvex objective permits.

So MMI parity is start-point parity, and reproducing it would mean reproducing HiGHS's pivoting rather than any property of the model.

THE THREE SUBSTITUTIONS FOR SCIPY, and why each is exact rather than approximate:

  1. scipy.linalg.qr(pivoting=True) behind independent_rows is replaced by a pivoted modified Gram-Schmidt over the same columns. Column-pivoted QR selects, at each step, the remaining column of largest residual norm, and that is precisely what the sweep below does; the retained index SET is the same, and the set is all the caller uses. (The R factor itself is never read.)
  2. scipy.optimize.linprog behind feasible_start is replaced by lp::simplex_solve on the same rows with the same [0,1] box, followed by the minimum-norm refinement qrf_noblo_start.m performs with QUADPROG. The LP is exact in either; the refinement is what keeps the start off a vertex, and on these objectives the start decides the answer.
  3. scipy.optimize.minimize(method='SLSQP') is replaced by FRANK-WOLFE over the same polytope. See solve_qrf_nlp for why this is the right answer and not a compromise: the available augmented-Lagrangian path has a derivative-free inner solve, and started at a VERTEX – which is what the phase-1 LP returns – it does not move at all. Conditional gradient uses the LP that is already here and returns an optimality CERTIFICATE.

THE REDUCTION IS NOT AN OPTIMIZATION, IT IS WHAT MAKES THE PROBLEM SOLVABLE. The raw decision vector has M^2 (N+1)^2 Kmax^2 MR + M Kmax entries and the equality block pins nearly all of them; substituting x = x0 + Z t for an orthonormal basis Z of null(Aeq) leaves a few free directions. The reference measured SLSQP at 14 iterations and under a second on the reduced problem against a failure to move on the raw one.

THE GUARD THAT MUST NOT BE DROPPED. reduce_equalities refuses an INCONSISTENT system rather than dropping the offending rows – an inconsistent polytope is a modelling error, and silently discarding it returns numbers for a model nobody wrote. The reference's second guard, a one-shot probe asking whether a feasible descent direction exists at the start (its SLSQP can silently return the phase-1 vertex), is SUBSUMED here: solve_qrf_nlp asks that same question at EVERY iterate, as its Frank-Wolfe gap, and terminates on the answer.

ARITHMETIC: transcendental (p log p, and the orthogonalization).

Definition in file mapqn_qrf_common.h.