Tropical Sampler
The content on this page was written by AI under human supervision.
Tropical Sampler is a Monte Carlo integrator for integrals, over positive variables, of products of polynomials with positive coefficients raised to powers. It takes the monomial exponent vectors and the powers and returns an unbiased estimate of the integral with a statistical error bar, together with the exact rational value of the simplified ("tropical") integral it samples from. It is one Python file written by BootLoops following Algorithm 1 of Borinsky, Sattelberger, Sturmfels and Telen (arXiv:2204.06414). A companion set of scripts, cegm_gj, uses the same cone decomposition with deterministic high-precision quadrature for one family of string-theory integrals; it sits in the same package directory, under cegm_gj/.
What it does
The integrals have the form
$$ Z \;=\; C\int_{\mathbb{R}_+^E} \prod_k Q_k(t)^{u_k}\,\prod_{e=1}^{E}(1+t_e)^{-c}\,dt , $$
where every $Q_k$ has positive coefficients (no cancellation between terms, or "subtraction-free"), the $u_k$ are real powers, and $c$ is large enough for convergence. Such integrals appear as marginal likelihoods in Bayesian statistics, as generalized hypergeometric integrals, and as Feynman and string integrals in parametric form.
In logarithmic variables $y=\log t$, replacing each polynomial by its largest monomial makes the log of the integrand piecewise linear; this is the tropical approximation. The linear pieces live on cones fixed by the Newton polytopesthe Newton polytope of a polynomial is the convex hull of the exponent vectors of its monomials of the $Q_k$. Together the cones form the normal fanthe division of direction space into cones according to which vertex of a polytope is the extreme one in that direction of the polytopes' Minkowski sum, with a unit segment per variable added for the decay factors. On a cone spanned by integer rays $r_i$, the columns of a matrix $R$, the approximate integral is elementary, $c_\sigma=|\det R|/\prod_i\beta_i$, with $\beta_i\gt0$ the decay rate along $r_i$. build_fan constructs the cones in exact integer arithmetic and returns each $c_\sigma$ and their sum $I^{\mathrm{tr}}$ as Fractions. Two checks come with it: random directions must each fall in exactly one cone, and refining the fan with an extra weight-zero polytope must leave $I^{\mathrm{tr}}$ unchanged.
To sample, a cone is drawn with probability $c_\sigma/I^{\mathrm{tr}}$, a point in it as $y=R\lambda$ with each $\lambda_i$ exponential of rate $\beta_i$, and the point is weighted by the ratio of the true integrand to the tropical one. Positive coefficients keep that ratio bounded, so $\hat Z=C\,I^{\mathrm{tr}}\,\overline{w}$ is unbiased with the usual $1/\sqrt{n}$ error. Two variants cut the variance and stay unbiased. The stratified one allocates samples per cone from a pilot run. The tilted one also fits the slope of the log-weight along each ray and lowers the exponential rates to match (never raises them, so the weights stay bounded). When the powers are large, as when they count data, only the tilted variant converges usefully.
build_fan accepts any integer point sets and weights. The weight function in the file is written for one family: each $Q_k$ is a sum, over subsets $D$ of the variables, of products $\prod_{e\in D}(1+4t_e)$ with positive integer multiplicities, with $c=N+2$ and $C=4^{-(E+1)N}$. This is the form in which marginal likelihoods of phylogenetic trees arise; another family needs its own log_weight. Expect three to four significant digits in minutes for fans of a few thousand cones in seven variables: enough for an independent numerical value or for ranking competing integrals, with further digits left to an exact method. Polynomials with coefficients of both signs are out of scope, because the weight bound fails.
The cegm_gj scripts apply the same decomposition to the Grassmannian string integral on $X(3,6)$ of Arkani-Hamed, He and Lam (arXiv:1912.08707). The integrand has four positive variables and twenty $3\times3$ minors raised to kinematic powers, in the positive coordinates of Giménez Umbert and Sturmfels (arXiv:2501.10805). Log space splits into 48 maximal cones, refined to 52 simplicial ones. Mapped to the unit cube, each cone's integrand is smooth, so tensor-product Gauss–Jacobi quadrature in arbitrary precision converges geometrically in the nodes per axis $n$, which quadrature over the undivided domain does not. A value counts as correct to $D$ digits only when two runs at different $n$ and working precision agree to $D$ digits, rounded down. Kinematic points, fan and certificates are specific to $X(3,6)$.
Examples
Build a fan and check it. Two variables, one polynomial with power 2 whose Newton polytope is the unit square, and decay $c=7$ per variable, entered as a segment from the origin to each unit vector with weight $-7$:
import numpy as np from fractions import Fraction from tropical_sampler import build_fan, tiling_probe E = 2 pat = np.array([[0, 0], [1, 0], [0, 1], [1, 1]], dtype=np.int64) seg1 = np.array([[0, 0], [1, 0]], dtype=np.int64) seg2 = np.array([[0, 0], [0, 1]], dtype=np.int64) groups = [(2, pat), (-7, seg1), (-7, seg2)] rng = np.random.default_rng(0) fan = build_fan(groups, E, rng) ok, cmin, cmax = tiling_probe(fan, rng, nprobe=500) extra = np.array([[0, 0], [2, 1], [1, 2]], dtype=np.int64) fan2 = build_fan(groups, E, rng, refine_extra=extra)
fan['cones'] holds 4 cones, the four quadrants, each a dict with the ray matrix R, gradient g, rates beta and exact csigma ($1$, $\tfrac14$, $\tfrac14$, $\tfrac1{16}$), and fan['Itr'] is Fraction(25, 16). tiling_probe returns (True, 1, 1): each of 500 random directions lies in exactly one cone. With the triangle extra added at weight zero, fan2 has 7 cones and the same Fraction(25, 16). These are the checks that python3 tropical_sampler.py --selftest runs: it prints one summary line with both cone counts, both exact values of $I^{\mathrm{tr}}$, the tiling result and PASS, and exits with the number of failed checks, 0 on a pass. It takes about a second.
Estimate the integral. On the same fan the polynomial is $Q(t)=\sum_{D\subseteq\{1,2\}}\prod_{e\in D}(1+4t_e)=(2+4t_1)(2+4t_2)$ with $u=2$, and $N=5$ so that $c=7$ matches the segment weights:
from tropical_sampler import make_eval_data, estimate_logZ, estimate_logZ_tilt
N, u = 5, 2
Dsets = [frozenset(), frozenset({0}), frozenset({1}), frozenset({0, 1})]
evdata = make_eval_data([(u, Dsets)], E)
logZ, se_rel, mean_w, se_w, wmin, wmax = estimate_logZ(
fan, evdata, N, E, 200_000, rng, 0.0)
logZ_t, se_rel_t, Zint, n_used, wmin_t, wmax_t = estimate_logZ_tilt(
fan, evdata, N, E, 20_000, 200_000, rng, 0.0)
The last argument is $\log C$; 0.0 returns the bare integral. estimate_logZ draws 200,000 points and returns the natural log of the estimate, its relative standard error (a few parts in a thousand here), the mean weight with its standard error, and the extreme weights. estimate_logZ_tilt takes a pilot size and a main size and returns the log estimate, relative error, the integral, the samples used and the extreme weights. This small integral factorizes into two Beta-function integrals and equals $(22/15)^2\approx2.151$, so exp(logZ) should sit within a few error bars of it.
Routines
Sampler (tropical_sampler.py)
python3 tropical_sampler.py --selftest— the self-test, and the file's only command-line mode (everything else is imported): builds the two-variable fan of the first example, checks the 4 cones and $I^{\mathrm{tr}}=25/16$, the tiling probe, and the 7-cone refinement with the same $I^{\mathrm{tr}}$; exit code is the number of failed checks.build_fan(groups, E, rng, log=None, refine_extra=None, skip_chain=0)— exact normal fan fromgroups, a list of(weight, integer point set)pairs (a polynomial with weight $u_k$, or a variable's segment with weight $-c$); returnsdict(cones, Itr, E).tiling_probe(fan, rng, nprobe=200, tol=1e-9)— checks that random directions each lie in exactly one cone; returns(ok, min_cover, max_cover).idet(M),primitive_normal(pts),hull_vertices(pts),minkowski_vertices(point_sets, log=None)— exact-geometry helpers: Bareiss integer determinant, primitive integer facet normal, hull vertices of integer points, vertices of an iterated Minkowski sum.make_eval_data(patterns, E)— packspatterns, a list of(u_k, Dsets)withDsetsa list of frozensets of variable indices (repeats count as multiplicity), into the arrayslog_weightuses.log_weight(Y, evdata, N, E)— log of true integrand over tropical integrand at a block of points in log coordinates, one per row.sample_tropical(fan, nsamples, rng, batch=500_000)— generator yielding blocks of points drawn from the tropical density.estimate_logZ(fan, evdata, N, E, nsamples, rng, ln_prefac)— plain estimator; returns(logZ, relative_stderr, mean_w, se_w, wmin, wmax).estimate_logZ_strat(fan, evdata, N, E, n_pilot, n_main, rng, ln_prefac, batch=400_000)— stratified by cone with pilot-based allocation; returns(logZ, relative_stderr, Z_int, n_used, wmin, wmax).estimate_logZ_tilt(fan, evdata, N, E, n_pilot, n_main, rng, ln_prefac, batch=400_000, min_pilot=64)— stratified with per-cone exponential tilting; same return tuple; the variant for large powers.ln_fraction(fr)— natural log of aFractionwithout converting it to a float.
Grassmannian string integral on $X(3,6)$ (cegm_gj/)
t3_prod.py sanity [--procs P] [--engine plain|hoist|c]— self-check: evaluates the integral at the first built-in kinematic point with $n=10$ nodes per axis at 160 bits over all 52 cones, compares with the stored reference value (pass below $10^{-38}$ relative), prints a one-line JSON verdict and exits 0 on pass. Logs and records go to$CEGM_GJ_OUT, or the current directory when unset; value strings end in a decimal exponent such ase3, so read them whole. The three engines (plain Python, a reorderedhoist, a C inner loop) are interchangeable.t3_prod.py cert --point {1,2}— writespoint{N}_certificate.json: kinematics, exact momentum conservation, the convergence conditions checked by linear programming, the minimum decay rate over the 52 cones, and the fan checksum.t3_prod.py run --point {1,2} --n N --prec BITS [--procs P] [--cpp C] [--engine plain|hoist|c]— re-derives the certificate, then computes one high-precision value and writes it with its settings torun_p{point}_n{N}_b{BITS}.json(with an engine suffix unlessplain);--cppsets how finely each cone is split across processes.gates.py A B [--bar D] [--relmax R]— rounded-down count of agreeing digits between two run JSONs or decimal strings; exit 0 only if at leastD(default 30) and, if given, the relative difference is belowR.bank.py— reassembles the recorded $X(3,6)$ values from the run JSONs in$CEGM_GJ_RUNDIR, recomputing every two-run comparison and stopping with an error if any falls short; tied to that recorded set of runs.tropical_cones.py— builds the 48-cone, 52-subcone fan with its coverage check and writesfan_x36.jsonto the current directory (never run it inside the package directory, where it would overwrite the stored fan).cone_engine.py,x36_engine.py— the per-cone evaluator (subcone_tasks,eval_subcone,run_cone_integral,mpfr_str) and the Gauss–Jacobi quadrature routines (gauss_jacobi_01,gj_selftest, and the undivided-domain evaluator kept for comparison).run_pilot.py,setup_x36.py,cone_identity_check.py,control_oracle_check.py— kinematic points and one-off verifications: exact solution of the constant-minor exponents, the linear-programming checkpolytope_check, the split-kinematics Gamma-function reference value, a double-precision Sobol cross-check, the per-cone integrand identity, and the Gamma-function formula against direct quadrature.cchunk.c— the C/MPFR inner loop behind--engine c; build in place withgcc -O2 -shared -fPIC -o cchunk.so cchunk.c -lmpfr -lgmp.
Used on this site
- Phylogenetics — independent Monte Carlo values of five-taxon marginal-likelihood integrals, checked against exactly known ones.
- The Grassmannian string integral —
cegm_gjcomputed the $X(3,6)$ integral on its 52-cone decomposition to 60 digits at one kinematic point and 32 at a second.
Requirements and source
Python 3 with NumPy and SciPy for the sampler (scipy.spatial.ConvexHull, which calls Qhull). cegm_gj also needs gmpy2, mpmath and sympy, plus the MPFR and GMP headers to build the optional C engine; its --procs default is 16, so set it on a small machine. The self-tests are python3 tropical_sampler.py --selftest for the sampler and python3 cegm_gj/t3_prod.py sanity --procs 4 for cegm_gj. Run the second with CEGM_GJ_OUT set to a scratch directory, so that its log file goes there rather than to the current directory; the repository's test list runs both, in that order. The code is tools/tropical-sampler/ in BootLoops' bootloops-dev repository (GitHub organization BootLoops-ai), released under the MIT license: the sampler is tropical_sampler.py and the string-integral scripts, with the stored fan fan_x36.json, are under cegm_gj/. The papers linked above are the sources of the sampling algorithm, the integrals and the positive parametrization.