Stats, memory and logging

The pages before this one were about what propagate computes. This one is about watching it run: the per-layer counters, how to turn a term count into a memory bill before committing to a run, the log stream, and the handful of arguments and environment variables that change how fast the answer arrives without changing the answer.

Stats

propagate_with_stats returns the same evolved sum propagate would, alongside a PropagationStats with one entry per layer in application order. Enabling it does not change the propagated sum.

from paulistrings import Circuit, PauliSum, truncation

observable = PauliSum.from_strings({"XIII": 0.25, "IXII": 0.25, "IIXI": 0.25, "IIIX": 0.25}, num_qubits=4)
circuit = Circuit(4)
circuit.rz(0.3, 0)
circuit.cnot(0, 1)
circuit.depolarize(0.01, [0, 1])
policy = truncation.coeff(1e-10)

evolved, stats = observable.propagate_with_stats(circuit, policy, direction="heisenberg")

print(stats.layers, stats.peak_terms, stats.final_terms)
print(stats.terms_in, stats.terms_out)   # one entry per layer, post-truncation
print(stats.gate_name)
print(stats.circuit_index)

4 5 5
[4, 4, 4, 4] [4, 4, 4, 5]
['depolarizing', 'depolarizing', 'cnot', 'PauliRotation']
[3, 2, 1, 0]

A layer is one applied channel (The loop), so the broadcast depolarize(p, [0, 1]) shows up as two entries and len(circuit) == stats.layers. The lists run in application order, which under direction="heisenberg" is the reverse of the circuit as written: the first entry above is the depolarizing channel pushed last. circuit_index recovers the position as written regardless of direction, and application_index is the loop counter.

terms_out[k] is the count after layer k’s truncation, terms_in[k + 1] == terms_out[k], and peak_terms is the maximum of those resident counts. nanos is wall-clock time per layer; for a partitioned or distributed run it is the maximum over partitions, and stats.partition holds the per-partition detail, one partition for an unpartitioned run (NUMA partitions, MPI ranks). The full field list is in the Library.

Estimate the memory a run needs

A term costs 16 × width + 16 bytes: two uint64 key words per width word (the x and z columns) plus a complex128 coefficient. PauliSum.width is the compile-time word tier the qubit count picked — 1 word up to 64 qubits, 2 up to 128, then 4, 8, 16 — so a 127-qubit sum is 48 B/term and a 1024-qubit sum 272 B/term.

stats.peak_terms is the peak resident term count between layers, so peak_terms × (16 × width + 16) is the sum’s own high-water mark:

print(stats.peak_terms, evolved.width, stats.peak_terms * (16 * evolved.width + 16))

5 1 160

Two things that estimate does not include. Each layer also allocates a transient in-layer expansion — a term fans out over the channel’s delta set before the merge collapses it — bounded by the fanout of the widest channel in the circuit (2 for a Pauli rotation, up to 16 for a dense two-qubit unitary) times one coset’s working set per worker, not times the whole sum. And buckets keep their capacity across layers, so a sum that shrinks after its peak does not give the memory back.

Run a short prefix of the circuit first and extrapolate peak_terms from the per-layer growth in terms_out: the growth is what sets the bill, and the truncation policy is the knob on it. A topn or approx_topn budget is the direct way to cap peak_terms when the bill has to be known in advance. The layout this arithmetic comes from is on Inside the engine; Showcase B7 reports this engine’s own measured peak RSS (10.2 GiB single-threaded) against a run this arithmetic would estimate for.

Logging

For running logs instead of a post-hoc summary, the library logs through the log facade under the target paulistrings.propagate — INFO on entry and exit of each propagate call, DEBUG once per layer:

import logging
import paulistrings

logging.basicConfig(level=logging.DEBUG)
logging.getLogger("paulistrings.propagate").setLevel(logging.DEBUG)
paulistrings.reset_log_cache()   # pyo3-log caches each logger's effective level

The INFO lines read propagate: N terms through n channels (direction), scope and propagate: n layers applied, before -> after terms, t s, scope, where scope is 1 partition, P partitions [placement], or partition r/P [placement] for one rank of a comm= run, whose counts are that rank’s own. A layer’s DEBUG line reads layer k/n [name]: before -> after terms, t ms (partition r/P, d remote deltas, m rows in) at every partition count, one line per partition. Placement, pinning and transport setup log under paulistrings.partitioned, and device setup (kernel compiles and the kernel cache, device choice, the NCCL and peer wires) under paulistrings.gpu, so each filters on its own.

Call reset_log_cache() again after changing the level mid-process; the cache is otherwise stale for the rest of the process. Leave logging off when timing — an enabled DEBUG filter adds a clock read per layer.

Every PAULISTRINGS_* environment variable is read once per process, the first time it is needed, and its effective value is logged at DEBUG under the target of the subsystem it tunes; an unparsable value logs a warning and the default applies. Set them before the first call, and alike on every rank of a comm= run.

Threads

In a single unpartitioned process, set the thread count via RAYON_NUM_THREADS before the interpreter starts — Rayon builds its global pool at the first propagate call and never resizes it, so setting the variable from inside an already-running script does not reliably reach it:

RAYON_NUM_THREADS=32 python my_script.py

RAYON_NUM_THREADS does not reach a partitioned run; there the thread count comes from the placement (NUMA partitions).

More threads are not always faster, and the layer class decides. On the reference two-socket host, dense-PTM-heavy circuits (general two-qubit unitaries) saturate the machine’s write bandwidth at about 16 threads and extra threads only add contention, while sparse-rotation and Clifford circuits stay latency-bound and take the full thread count profitably — 11–13× at 32 threads. Hyperthreads add no bandwidth and help only the latency-bound phases. The measured roofline behind those numbers is on Engine performance.

Every benchmark on this site is run single-threaded, RAYON_NUM_THREADS=1 exported before the interpreter starts, so that term counts and wall times are comparable across pages (Comparability rules).

Bucket sizing

target_bucket_len and min_buckets size the bucketed engine’s storage per layer; PauliSum.num_buckets reads back what was realized, which only ever grows. They are tuning knobs for the engine’s own benchmarks rather than for users, and the defaults are what every number on this site was measured with. target_bucket_len must be at least 1 and min_buckets at least 16; anything smaller is a ValueError.

Determinism

The correctness bar is agreement to floating-point tolerance, not bit-identical output. Equal-key summation order is unspecified: two runs of the same propagation can differ in the last bits of a coefficient, and a partitioned run can differ from an unpartitioned one, all within tolerance and all equally correct. Compare evolved sums with numpy.allclose on coefficients_array(), as Splitting a circuit is free does, never with equality.

Term counts, on the other hand, are exact integers and are the load-independent quantity every benchmark on this site quotes and gates on.

Untruncated runs

policy=None applies no truncation at all; only coefficients that are exactly zero after a merge are dropped. The result is the exact evolved operator, and it is the reference the cutoff sweep is scored against wherever the problem is small enough to afford it.

exact = observable.propagate(circuit, None, direction="heisenberg")
print(len(exact), "terms untruncated")

5 terms untruncated

Afford it deliberately: term growth without a filter is exponential in depth for a generic circuit. Estimate peak_terms from a prefix and the bytes-per-term arithmetic above before running a deep circuit untruncated, and expect the untruncated run to be the most expensive point of any sweep by a wide margin.