|
Dist m4ri 0.16.1
Computing distance of a classical or quantum CSS code
|
Authors: Leonid P. Pryadko & Weilei Zeng
Date: 2026-10-28
dist-m4ri computes, or brackets, the minimum distance of a binary classical linear code, of a quantum CSS code, or of a detector error model (DEM) generated by Stim. It implements two algorithms and their combination:
The package contains the multithreaded program dist_m4ri, the legacy single-threaded program dist_m4ri_old, and the Python wrapper dist_m4ri.py (NumPy and SciPy matrices, Stim circuits and DEMs, a persistent cache of the results, and an optional qec_cert certification engine via solver="cert" or "auto" / engine=). This page summarizes how the program works; README.md describes all parameters (also dist_m4ri --help and dist_m4ri --morehelp; dist_m4ri --help maxC shows the help of one parameter, with the related ones).
All matrices are binary, with n columns, in Matrix Market format (.mtx or .mmx). For a CSS code, they obey the orthogonality conditions Hx Hz^T = 0, Hx Lz^T = 0, and Lx Hz^T = 0, and the number of encoded qubits is k = n - rank(Hx) - rank(Hz), which is also the rank of Lx Lz^T for complete sets of logical operators Lx and Lz. With debug=16, the program prints the ranks of the input matrices and k.
Each RW step computes the reduced row echelon form of H with the columns taken in a random order. This fixes an information set, the n - rank(H) non-pivot columns. Each non-pivot column j gives a vector in the null space of H, with a 1 at position j and its other non-zero entries at pivot positions; such vectors have at most rank(H) + 1 non-zero entries and tend to have small weights, which qualitatively explains why the algorithm works. A codeword c is found in a step if exactly one of its positions is a non-pivot column, and every codeword of the minimum weight is found with some column order (see the QDistRnd manual, Sec. 3.1). The candidates with L c != 0 are kept in a hash table, together with the number of times each of them has been found (hit counts). For matrices with at least 128 rows, the elimination processes 20 to 32 pivots per pass over the rows ("Method of the Four Russians", with tables of combinations of the pivot rows), which gives the same result about twice as fast for large matrices.
Variants: a window localizes the column order around a random column, in every step with kwin > 0 (kwin columns), and by default (kwin=0) for n >= 500 in every other step, with windows of win_hyb consecutive columns: 256 for a DEM, whose errors are in time order, and 512 for matrices; kwin=-1 turns the windows off. The window columns come first in the elimination order, so that a codeword with most of its columns in the window is found far more often; a window much larger than the minimum-weight codewords, e.g., one round of a DEM, has many columns which cannot be pivots and helps little. With the RW pairs (pairs, on by default), a step also finds the codewords with exactly two non-pivot columns, the sums of two such vectors (the algorithms of Lee-Brickell and Stern with p = 2): only the pairs of columns which agree on l random pivot rows are checked (pairs=l; auto: about one non-pivot column per bucket), and only those whose weights differ by less than the weight limit. This finds the minimum-weight codewords of circuit DEMs 7 to 26 times faster. With ksub > 0 (experimental, should not be used; a warning is printed), each step reduces the span of ksub rows of a fixed basis of ker(H) (cache-resident elimination): a few codewords are then found much more often than others, so that the hit-count criterion min_hits may end RW early.
RW ends after steps steps (default: 100000), at the timeout, as soon as a codeword of weight at most wmin is found, or when the hit-count criterion min_hits is met (see Convergence below).
The connected-cluster algorithm grows clusters of columns. A cluster starts at a column i and only uses columns j > i; at each level, a column is added from the support of a check that is not satisfied by the cluster. Round w enumerates all such clusters of weight w (the start columns are distributed among the threads; with few start columns per thread, a start column is split into work units which also fix the first one to three columns added). A round without codewords raises the lower bound to dmin = w + 1, and a codeword found in a round gives the exact distance. This search finds every minimum-weight codeword, starting from its smallest column: a proper subset of a minimum-weight codeword has a non-zero syndrome (otherwise, the subset or the rest of the codeword would be a lighter codeword), and the support of each of its unsatisfied checks contains another position of the codeword [DKP17]. The work of round w grows roughly as n (D - 1)^(w - 1), where D is the maximum row weight of H. With smax=0 (default), the syndrome of a cluster is kept as a bitset, and the last column of a cluster is found by a hash lookup: the syndrome left after the next-to-last column must be that of a single column of H.
Since version 0.15.0, the search is a branch and bound (cc_prune=1, the default) [LOMKZ26, Lee26]: a cluster is extended from the unsatisfied check with the fewest available columns; the column of each branch is forbidden in the later branches, so that every codeword is found exactly once; and a cluster is pruned if lower bounds on the number of columns still needed exceed the weight left (an unsatisfied check without available columns, or a packing of unsatisfied checks with disjoint sets of available columns). It visits far fewer clusters, by factors which grow with the weight (e.g., 30-100 times fewer, in 10-26 times less CC time, in the rounds w = 11..13 of the codes with d = 14 to 20 of the benchmark), and it finds the same codewords; the proof is in doc/cc_completeness.tex. With cc_prune=0, the first unsatisfied check in the order of the rows of H is taken, as in earlier versions: for matrices, the rows are sorted by increasing weight (sort_rows=1, the default), which reduces the work when the row weights are very different; the detectors of a DEM are kept in time order (sort_rows=0, the default with fdem), which is better there.
The expert options start (with start_scan), cbeg/cend, and noscan restrict the search; each of them prints a warning, since the reported dmin is then not necessarily a lower bound (dmax is always valid). With smax > 0, CC also computes the confinement profile, the minimum non-zero syndrome weight for each cluster weight.
The coordinator splits the threads between CC and RW for each CC round, according to the predicted CC work of the round (the measured work of the previous round times its growth factor), the measured time of a RW step, the number of remaining RW steps, and the remaining time. Once a codeword of weight dmax is known, the round w = dmax - 1, which certifies dmin = dmax, runs on all threads; the expected distance dexp only affects the thread split before RW finds any codeword. When RW converges (min_hits) or runs out of steps, its threads join the CC rounds. A CC round predicted not to finish before the timeout is not started while RW runs. With a stop target dstop = U (not an upper bound), the run ends as soon as CC has certified dmin >= U, unless a codeword of weight < U is found; e.g., for a CSS code with d = min(dX, dZ) <= U already known, the sector runs need no more.
The RW upper bound is not certified unless CC confirms it. As in QDistRnd (manual, Sec. 3.3), the probability to miss a codeword lighter than the minimum weight w found is estimated from the hit counts: if the m codewords of weight w were found n_1, ..., n_m times, then P_fail < exp(-<n>), where <n> = (n_1 + ... + n_m) / m is the average number of times a minimum-weight codeword has been found. This assumes that all codewords of the same weight are found equally often, which is tested with a dispersion test of the hit counts (a warning is printed with debug=1 if the hit counts are strongly non-uniform). The second estimate assumes uniformly random information sets: a given codeword of weight w is found in one step with the probability P1(w) = w C(n-w, s-1) / C(n, s), where s = n - rank(H) is the number of non-pivot columns and C(a, b) is the binomial coefficient, plus, with the RW pairs, P2(w) = C(w, 2) C(n-w, s-2) / C(n, s) * C(r-w+2, l) / C(r, l) with r = rank(H). See the section "RW Convergence: Empirical Test and Probability to Find a Codeword" in README.md for details.
The program prints one line to stdout, dmin dmax rw_steps: dmin - 1 is the largest cluster weight analyzed by CC without finding a codeword, dmax the smallest weight of a codeword found (or the supplied upper bound dmax; never the stop target dstop), and rw_steps the number of completed RW steps (0 if CC found a minimum-weight codeword); dmin = dmax is the exact distance. All diagnostic output goes to stderr, as selected by the debug bitmap (default: 3, the summary and the progress; debug=0 is silent, debug=7 adds a periodic status line). Example (run from the repository root):
The program is intended for use with recent gcc compilers under Linux, and it uses the m4ri library for binary linear algebra. To install the library under Ubuntu, run
Then run make -j all in the src directory, which compiles dist_m4ri and dist_m4ri_old; make test runs the C test suite, and pytest tests/test_dist_m4ri.py (from the repository root) the tests of the Python wrapper.
The documentation for this package is available at https://qec-pages.github.io/dist-m4ri/, and the package is hosted on github. This file (introduction.md) is the main page of the doxygen documentation.
If you use this program, please cite:
Other related papers and software: