perm_sampling.h File Reference

RANDOMIZED permanents: the AdaPart rejection sampler and the Huber-Law acceptance-rejection importance sampler. More...

#include <algorithm>
#include <cmath>
#include <cstddef>
#include <cstdint>
#include <map>
#include <random>
#include <set>
#include <string>
#include <vector>
#include "line/util/error.h"
#include "line/util/matrix.h"

Include dependency graph for perm_sampling.h:

Go to the source code of this file.

Classes
class	line::perm::AdaPartSampler
	Adaptive partitioning (AdaPart) sampler for the permanent. More...
class	line::perm::HuberLawSampler
	Huber-Law acceptance-rejection sampler for the permanent. More...

Namespaces
namespace	line
namespace	line::perm

Functions
double	line::perm::perm_adapart (const Matrix< double > &m, std::uint64_t seed=0)
	AdaPart estimate of the permanent, with the reference defaults.
double	line::perm::perm_huberlaw (const Matrix< double > &m, std::uint64_t seed=0)
	Huber-Law estimate of the permanent, with the reference defaults.

Detailed Description

RANDOMIZED permanents: the AdaPart rejection sampler and the Huber-Law acceptance-rejection importance sampler.

Port of python/line_solver/api/perm/sampling.py, itself a twin of jline.lib.perm.AdaPartSampler and jline.lib.perm.HuberLawSampler, and of MATLAB's perm_adapart.m and perm_huberlaw.m.

WHAT THESE ARE FOR. The exact algorithms in permanent.h cost 2^n or n!, and the two schemes in perm_approx.h buy their speed with a bias no caller can bound. These two are unbiased Monte Carlo estimators instead: the estimate carries sampling noise, but no systematic error, so averaging more draws converges to the permanent rather than to a nearby number.

• perm_adapart recursively partitions the permutation space, bounds each part by the Soules column bound, and draws a part with probability proportional to its bound. The acceptance ratio scales the root bound into an unbiased estimate.
• perm_huberlaw rescales the matrix to doubly stochastic, draws a permutation column by column under the Huber-Law bound on the remaining permanent, and multiplies the acceptance ratio by the rescaling constant.

TWO ADAPART FIXES ARE CARRIED HERE, and both are load-bearing rather than cosmetic:

1. The column search considers ONLY the columns still unassigned. Scoring an already-assigned column just re-derives its own constraint, which always looks cheapest, so the sampler re-splits the same column forever and never completes an assignment.
2. The expansion discounts the bound of the element ACTUALLY REMOVED, not the root bound. With the root bound the running total drifts and the "refine until improved" loop cannot exit.

Even so the Soules bound is TIGHT on a matrix with equal entries, so no refinement can improve it; the loop therefore stops on the first non-improving expansion, which is exactly the shape a replicated demand matrix produces (see pfqn_jointmarg).

BOTH REQUIRE A NONNEGATIVE MATRIX, and both need FULL SUPPORT: an exact zero makes the Huber-Law rescaling degenerate and starves the AdaPart acceptance. The callers in the pfqn layer refuse a structurally zero matrix rather than flooring it.

REPRODUCIBILITY: the draws come from a seeded std::mt19937_64, so a run with a fixed seed repeats. The Python twin draws from a numpy Generator and the JAR from java.util.Random, so the three agree in distribution but not sample by sample.

ARITHMETIC: double. Both are floating-point Monte Carlo schemes.

Definition in file perm_sampling.h.