From 38f160c431a78ffcaceb05f4fc564c738db5a923 Mon Sep 17 00:00:00 2001 From: Roger-luo Date: Tue, 11 Aug 2026 23:48:29 -0400 Subject: [PATCH 1/5] bench(cross-library): compare ppvm against four other Pauli-propagation engines MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Adds a single-threaded harness that runs the same two Trotter workloads — TFIM magnetization and a Heisenberg autocorrelator — through ppvm, PauliPropagation.jl, PauliStrings.jl, Qiskit's pauli-prop and Algorithmiq's monoprop, from one shared environment contract into one CSV. The driver refuses to report a timing until every engine has dumped its full propagated support at a small width and matched ppvm's term-for-term. That check is load-bearing rather than decorative: it is what pinned down monoprop's exp(+iθH) angle convention (θ = -θ_spec/2 on a reversed gate list), which otherwise reproduces the observable to nine digits off the wrong operator. monoprop needs a second guard. Left alone it takes one serial partition per physical core, so an uncapped run reports a ~3x faster wall time than the serial one against engines that never had the option. The runner pins monoprop_NUM_THREADS=1 / monoprop_PARTITIONS=off, and because both are read once into a cached C++ static, it also measures its own CPU/wall ratio around the timed region and exits non-zero if the cap did not land. On an M4 Pro, ppvm beats both Julia engines everywhere (2.2-11.3x vs PauliPropagation.jl, 4.9-13.2x vs PauliStrings.jl, on an identical support) and loses to monoprop on both models, by up to 3.0x on TFIM at n=64. The README records the run along with the per-engine truncation differences that a ratio should not be read without. Co-authored-by: Cursor --- benchmarks/README.md | 16 + benchmarks/cross-library/README.md | 255 ++++++++++++++ benchmarks/cross-library/plot_xbench.py | 244 +++++++++++++ benchmarks/cross-library/run_xbench.py | 366 ++++++++++++++++++++ benchmarks/cross-library/xbench_monoprop.py | 243 +++++++++++++ benchmarks/cross-library/xbench_qiskit.py | 190 ++++++++++ crates/ppvm-pauli-sum/examples/xbench.rs | 266 ++++++++++++++ julia-benchmarks/Manifest.toml | 54 ++- julia-benchmarks/Project.toml | 1 + julia-benchmarks/benches/xbench_pp.jl | 129 +++++++ julia-benchmarks/benches/xbench_ps.jl | 163 +++++++++ 11 files changed, 1918 insertions(+), 9 deletions(-) create mode 100644 benchmarks/cross-library/README.md create mode 100644 benchmarks/cross-library/plot_xbench.py create mode 100644 benchmarks/cross-library/run_xbench.py create mode 100644 benchmarks/cross-library/xbench_monoprop.py create mode 100644 benchmarks/cross-library/xbench_qiskit.py create mode 100644 crates/ppvm-pauli-sum/examples/xbench.rs create mode 100644 julia-benchmarks/benches/xbench_pp.jl create mode 100644 julia-benchmarks/benches/xbench_ps.jl diff --git a/benchmarks/README.md b/benchmarks/README.md index 39dc22a30..65b95b269 100644 --- a/benchmarks/README.md +++ b/benchmarks/README.md @@ -87,3 +87,19 @@ uv run --with matplotlib python benchmarks/plot_branch_coalesce.py \ - `plot_branch_coalesce.py` — left panel: time vs `m` (log-log); right panel: sort-merge speedup `t_hash / t_sortmerge` vs `m`, with the crossover line and the "hash wins" band. + +# Cross-library Pauli propagation + +`ppvm` against [PauliPropagation.jl][xpp], [PauliStrings.jl][xps], Qiskit's +[pauli-prop][xqk] and Algorithmiq's [monoprop][xmp], on TFIM Trotter and +Heisenberg autocorrelator workloads, all single-threaded and all validated to +propagate the identical operator term-for-term before anything is timed. + +Lives in [`cross-library/`](cross-library/README.md) — see that README for the +circuit spec, the shared environment contract, the per-engine caveats (notably +that monoprop is parallel unless capped), and a measured run. + +[xpp]: https://github.com/MSRudolph/PauliPropagation.jl +[xps]: https://github.com/nicolasloizeau/PauliStrings.jl +[xqk]: https://github.com/Qiskit/pauli-prop +[xmp]: https://github.com/Algorithmiq/monoprop diff --git a/benchmarks/cross-library/README.md b/benchmarks/cross-library/README.md new file mode 100644 index 000000000..dcea69298 --- /dev/null +++ b/benchmarks/cross-library/README.md @@ -0,0 +1,255 @@ +# Cross-library Pauli-propagation benchmark + +`ppvm` against the four other single-threaded Pauli-propagation engines we know +of, on two workloads: + +| library | language | entry point used | +|---|---|---| +| **`ppvm`** (this repo) | Rust | `PauliSum>` | +| [PauliPropagation.jl][pp] | Julia | `propagate(PauliRotation(...), psum; min_abs_coeff)` | +| [PauliStrings.jl][ps] | Julia | `trotter_step!(O, gates; truncation, truncate_every)` | +| [pauli-prop][qk] (Qiskit) | Rust-accelerated Python | `propagate_through_circuit(op, qc, max_terms, atol, frame="h")` | +| [monoprop][mp] (Algorithmiq) | C++ with Python bindings | `PauliPropagator.from_circuit(circuit, op, cutoff, lower_atol)` | + +```bash +# Everything, with the term-for-term agreement check first. +uv run --no-project python3 benchmarks/cross-library/run_xbench.py \ + --qubits-tfim 8,16,24,32,40,48,56,64 \ + --qubits-heisenberg 6,8,10,12,14 \ + --steps 10 --atol 1e-6 --iters 2 --out target/xbench + +uv run --no-project --with matplotlib python3 benchmarks/cross-library/plot_xbench.py \ + --csv target/xbench/results.csv --out target/xbench/xbench.png \ + --title "Pauli propagation: ppvm vs PauliPropagation.jl, PauliStrings.jl, pauli-prop, monoprop" +``` + +A subset, if you prefer: `--libs ppvm,monoprop`. + +## What is measured + +Both workloads are **Heisenberg-picture propagation of an observable through an +explicit first-order Trotter product of Pauli rotations**, with a +coefficient-magnitude truncation after every gate. Everything below is fixed +across all five engines — the gate list, its order, the angles, the truncation +rule, and the readout. + +### `tfim` — transverse-field Ising, magnetization + +`H = J Σᵢ ZᵢZᵢ₊₁ + h Σᵢ Xᵢ` on an open chain. One step, in this order: + +1. `RX(2h·dt)` on site `0, 1, …, n−1` +2. `RZZ(2J·dt)` on bond `(0,1), (1,2), …, (n−2,n−1)` + +Observable `O = Σᵢ Zᵢ`; readout `⟨0…0|O(t)|0…0⟩`, i.e. the sum of the +coefficients of the X-free terms. + +### `heisenberg` — isotropic Heisenberg + field, autocorrelator + +`H = J Σᵢ (XᵢXᵢ₊₁ + YᵢYᵢ₊₁ + ZᵢZᵢ₊₁) + h Σᵢ Zᵢ`. One step, in this order: + +1. `RXX`, `RYY`, `RZZ` at `2J·dt` on bond `(0,1)`, then `(1,2)`, … +2. `RZ(2h·dt)` on site `0, 1, …, n−1` + +Observable `O = Z₀`; readout the autocorrelator `S(t) = tr[Z₀·O(t)]/2ⁿ`, which +is just the coefficient of `Z₀` (the Paulis are orthonormal under that pairing). + +`θ = 2·c·dt` for a Hamiltonian term `c·G` is the convention four of the five +engines use directly for `exp(iθ/2·G)·P·exp(−iθ/2·G)`, so the propagated +operators are identical, not merely similar. monoprop is the exception and needs +a conversion; see below. + +## The parameter contract + +Every runner reads the same environment variables and writes the same CSV, so +they can be run directly as well as through the driver: + +| variable | meaning | +|---|---| +| `MODEL` | `tfim` or `heisenberg` | +| `QUBITS` | comma-separated widths | +| `STEPS` | Trotter steps | +| `DT`, `JCOUP`, `HFIELD` | `dt`, `J`, `h` | +| `ATOL` | truncation threshold on `|c|` | +| `ITERS` | timed repeats; the **minimum** is reported | +| `DUMP` | print the propagated support instead of timing it | +| `MAX_TERMS` | `pauli-prop` only — its mandatory cap (see below) | +| `monoprop_NUM_THREADS`, `monoprop_PARTITIONS` | `monoprop` only — the thread cap (see below) | + +CSV columns: `model,library,qubits,steps,dt,atol,time_s,terms,observable`. + +```bash +MODEL=heisenberg QUBITS=8,10,12 STEPS=10 ATOL=1e-6 ITERS=2 \ + cargo run --release -p ppvm-pauli-sum --example xbench +``` + +## Validation — why the driver refuses to time first + +`run_xbench.py` runs every engine with `DUMP=1` at `n=4, steps=3, atol=1e-14` +and diffs the whole propagated support against `ppvm`'s, term for term. It +aborts on any missing term, extra term, or coefficient difference above +`--validate-tol` (1e-10 by default). + +This is not ceremony. Three real bugs in this harness produced numbers that +looked completely reasonable and were wrong: + +* **The `pauli-prop` circuit was reversed.** In the Heisenberg frame it + conjugates from the *end* of the instruction list backwards, so the spec's + gate order needs `reversed()` on append. Appending forward propagates a + different operator — 108 terms against TFIM's 124, 61 against Heisenberg's 64, + coefficients off by up to 0.1 — while the readout still matched to 9 digits. +* **Duplicate Paulis were silently collapsing.** `propagate_through_circuit` + returns a `SparsePauliOp` that may list the same Pauli more than once, so + `len(op)` is a row count, not a support size, and a `{label: coeff}` dict + comprehension keeps the last duplicate instead of summing them. +* **monoprop's angle convention is neither the spec's nor Qiskit's.** Its + `ExpGate` applies `exp(+iθH)`, so the spec's `exp(−iθ_spec/2·G)` needs + `θ = −θ_spec/2`, *and* the gate list needs reversing like `pauli-prop`'s. All + eight sign/order combinations were tried against the reference: forward order + loses terms outright (109 against TFIM's 124), and every wrong angle keeps the + right support while moving coefficients by up to 1.3. Only `θ = −θ_spec/2` on + a reversed list lands within 5e-13. + +An observable is one scalar and can agree by luck or by cancellation. A +term-for-term diff cannot. + +## monoprop is parallel unless you stop it + +This is the one caveat that changes a headline number rather than a decimal, so +it gets its own section. + +monoprop takes **one serial partition per physical core** when left alone: with a +single MPI rank, `resolve_partition_count_` reads the core count and fans out. +The PyPI wheels are built without MPI (`monoprop.has_mpi == False`), which is +easy to misread as "therefore serial" — it is not. On the 14-core machine below, +Heisenberg at `n=12` measures: + +| | wall | CPU (user+sys) | CPU/wall | +|---|---:|---:|---:| +| default | 0.66 s | 6.28 s | 9.48× | +| `monoprop_NUM_THREADS=1 monoprop_PARTITIONS=off` | 1.97 s | 1.97 s | 1.00× | + +So an uncapped monoprop reports a **3× faster** wall time than the serial one, +against engines that never had the option. Both variables are read once into a +cached C++ static, so they must be in the environment before the first +propagation — the runner sets them itself as well as receiving them from the +driver. + +Because a stale environment variable fails silently, `xbench_monoprop.py` +**measures its own CPU/wall ratio around the timed region and exits non-zero** if +it exceeds `CPU_WALL_MAX` (1.5). Every monoprop row in the run below reported +`cpu/wall=1.00`. + +## Known differences between the engines + +These are real and are the reason the plot shows the **workload size** next to +the runtime. Read them before quoting a ratio. + +* **`pauli-prop` and `monoprop` truncate by a different rule.** Both prune at + branch-creation time, where the other three accumulate first and drop the + merged coefficient. `ppvm`, PauliPropagation.jl and PauliStrings.jl then agree + on the support **exactly**, at every width, on both models; the other two do + not. On the sweep below `pauli-prop` runs 11–36 % *under* the reference support + on TFIM (the gap widening with `n`) and 4–15 % *over* it on Heisenberg; + `monoprop` is within 1 % on TFIM and 4–12 % over on Heisenberg. The driver + prints a `note:` line for any engine more than 2 % off, and the plot's third + panel is there so the ratio is never read without it. +* **`monoprop` tracks more rows than it reports.** It retains monomials whose + coefficient has cancelled to exactly zero. Its `size()` at Heisenberg `n=14` is + 7 355 928 rows against the 3 204 697 terms above threshold — 2.3× — while on + TFIM the two are within 1 %. The `terms` column is the above-threshold support, + which is what the other four engines mean by it; `size()` goes to stderr next + to each row. +* **`pauli-prop` also has a mandatory `max_terms` cap** with pre-allocation. The + runner defaults it to `2²²` and **fails the run** if the support ever reaches + it, rather than quietly reporting a differently-truncated number. +* **`monoprop` has a mandatory Pauli-weight `cutoff`**, which the others have no + analogue of. The runner sets it to `n` — the whole register — so it never binds + and `lower_atol` is the only truncation in play. +* **PauliStrings.jl stores `im^{#Y}` inside the coefficient** (its `Matrix` + convention). The dump goes through `op_to_strings`, which puts its + coefficients on the same real footing as everyone else's. +* **PauliStrings.jl's own front door is `evolve(H, O, tspan; method=Trotter())`, + which we do not use** — it derives the gate list from the Hamiltonian's + internal string order, which is not the spec's order. The runner builds the + `TrotterGate` vector by hand instead. +* **Julia is timed after a warm-up run** so the reported time excludes JIT. +* **Circuit construction is outside the timed region for every engine.** + monoprop's `propagate` re-expands its gate list internally on each call, which + no flag hoists out; measured at 0.0–2.2 % of its total at these widths, so it + is left in rather than worked around. +* Everything is single-threaded (`julia -t1`; `pauli-prop` is single-threaded by + design; `monoprop` is capped as above; `ppvm` here uses no Rayon). + +## A run + +Apple M4 Pro (10 P + 4 E cores), macOS 26.6, rustc 1.96.0, Julia 1.12.6, +monoprop 0.8.0, PauliStrings.jl 1.10.1; single-threaded, `--steps 10 --dt 0.1 +--j 1 --h 1 --atol 1e-6 --iters 2`. Runtime relative to `ppvm`; **below 1.00× +means that engine beat `ppvm`**. `terms` is the shared support the three +accumulate-then-truncate engines all reach exactly. + +`tfim`, widths 8…64: + +| n | terms | ppvm | PauliPropagation.jl | PauliStrings.jl | pauli-prop | monoprop | +|---:|---:|---:|---:|---:|---:|---:| +| 8 | 4 701 | 0.0028 s | 2.42× | 4.90× | 11.13× | 0.98× | +| 16 | 19 529 | 0.0122 s | 3.10× | 6.26× | 5.53× | 0.69× | +| 24 | 34 353 | 0.0260 s | 5.78× | 8.05× | 5.32× | 0.54× | +| 32 | 49 177 | 0.0431 s | 4.98× | 8.39× | 7.05× | 0.45× | +| 40 | 64 001 | 0.0657 s | 6.33× | 9.20× | 5.28× | 0.43× | +| 48 | 78 825 | 0.0899 s | 7.32× | 9.02× | 5.37× | 0.40× | +| 56 | 93 649 | 0.1188 s | 10.46× | 9.84× | 5.33× | 0.35× | +| 64 | 108 473 | 0.1429 s | 11.26× | 13.18× | 5.71× | 0.33× | + +`heisenberg`, widths 6…14: + +| n | terms | ppvm | PauliPropagation.jl | PauliStrings.jl | pauli-prop | monoprop | +|---:|---:|---:|---:|---:|---:|---:| +| 6 | 1 022 | 0.0028 s | 4.28× | 6.77× | 5.02× | 1.03× | +| 8 | 16 324 | 0.0505 s | 3.43× | 6.26× | 2.02× | 0.83× | +| 10 | 225 353 | 0.5290 s | 3.11× | 6.79× | 1.51× | 0.77× | +| 12 | 1 174 849 | 2.7528 s | 2.59× | 6.50× | 1.13× | 0.67× | +| 14 | 2 915 879 | 6.0906 s | 2.23× | 7.09× | 0.88× | 0.54× | + +### Reading it + +**`ppvm` beats both Julia engines everywhere**, by 2.2–11.3× against +PauliPropagation.jl and 4.9–13.2× against PauliStrings.jl, and the margin widens +with `n` on TFIM. Those three carry an identical support term-for-term, so those +are clean ratios. + +**`ppvm` loses to `monoprop`** on both models — by up to 3.0× on TFIM at `n=64` +and 1.9× on Heisenberg at `n=14`, with the gap widening in `n` on both. On TFIM +the two carry the same support to within 1 %, so that column is a clean loss. On +Heisenberg `monoprop` is doing 4–12 % *more* above-threshold work than `ppvm` and +still finishing sooner, so the gap there is if anything understated — though it +is also tracking 2.3× as many rows in total, so the two engines are not making +the same space/time trade. + +**`pauli-prop` is 5.3–11.1× slower on TFIM but competitive on Heisenberg**, +reaching 0.88× at `n=14`. Read its column against the workload note: on TFIM at +`n=64` it is propagating 36 % fewer terms than the reference three, and on +Heisenberg `n=10`–`12` about 15 % more. + +## Scaling note + +The TFIM support grows roughly linearly in `n` at fixed depth, so that sweep +reaches 64 qubits cheaply. The Heisenberg support grows much faster — three +non-commuting bond rotations per bond per step — and the `Z₀` autocorrelator +itself converges once `n` exceeds the light cone (identical from `n≈10` at +`steps=10`), so widths past that measure scaling rather than new physics. + +## Files + +- `../../crates/ppvm-pauli-sum/examples/xbench.rs` — the `ppvm` runner. +- `../../julia-benchmarks/benches/xbench_pp.jl` — PauliPropagation.jl. +- `../../julia-benchmarks/benches/xbench_ps.jl` — PauliStrings.jl. +- `xbench_qiskit.py` — `pauli-prop`. +- `xbench_monoprop.py` — `monoprop`, including the thread-cap assertion. +- `run_xbench.py` — validation, driving, CSV merge, summary table. +- `plot_xbench.py` — the figure. + +[pp]: https://github.com/MSRudolph/PauliPropagation.jl +[ps]: https://github.com/nicolasloizeau/PauliStrings.jl +[qk]: https://github.com/Qiskit/pauli-prop +[mp]: https://github.com/Algorithmiq/monoprop diff --git a/benchmarks/cross-library/plot_xbench.py b/benchmarks/cross-library/plot_xbench.py new file mode 100644 index 000000000..cd5f28dd5 --- /dev/null +++ b/benchmarks/cross-library/plot_xbench.py @@ -0,0 +1,244 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +"""Render the cross-library benchmark from `run_xbench.py`'s `results.csv`. + +One row per model, three panels: runtime vs qubit count (log-y), runtime +relative to `ppvm`, and the size of the propagated operator — which is the +workload, and the thing that makes the runtimes comparable in the first place. + + uv run --no-project --with matplotlib python3 plot_xbench.py \ + --csv target/xbench/results.csv --out target/xbench/xbench.png +""" + +from __future__ import annotations + +import argparse +import csv +import textwrap +from collections import defaultdict +from pathlib import Path + +import matplotlib + +matplotlib.use("Agg") +import matplotlib.pyplot as plt + +# Slots from the validated default palette, in this fixed order — it clears +# every hard gate on the adjacent pairlist that line charts use (worst CVD +# ΔE 9.2, worst normal-vision ΔE 27.5). Colour follows the entity, so a run with +# a subset of `--libs` keeps every survivor's hue. Magenta and aqua sit below +# 3:1 on the light surface, which obligates relief — `results.csv`, the driver's +# summary table, and the per-series markers are it. +BASELINE = "ppvm" +SERIES = { + "ppvm": ("ppvm (this repo)", "#2a78d6", "o"), + "pauli-propagation-jl": ("PauliPropagation.jl", "#eb6834", "s"), + "pauli-strings-jl": ("PauliStrings.jl", "#1baf7a", "^"), + "pauli-prop": ("pauli-prop (Qiskit)", "#4a3aa7", "D"), + "monoprop": ("monoprop", "#e87ba4", "v"), +} +TEXT_PRIMARY = "#0b0b0b" +TEXT_SECONDARY = "#52514e" +GRID = "#d8d7d3" +SURFACE = "#fcfcfb" +MODEL_TITLES = { + "tfim": "TFIM Trotter — ⟨0|Σ Z_i(t)|0⟩", + "heisenberg": "Heisenberg correlations — tr[Z₀ Z₀(t)]/2ⁿ", +} + + +def load(path: Path): + """`{model: {library: {n: row}}}`.""" + out: dict[str, dict[str, dict[int, dict[str, str]]]] = defaultdict( + lambda: defaultdict(dict) + ) + with path.open() as fh: + for row in csv.DictReader(fh): + out[row["model"]][row["library"]][int(row["qubits"])] = row + return out + + +def style(ax) -> None: + """Recessive grid and axes; the data carries the ink.""" + ax.set_facecolor(SURFACE) + ax.grid(True, which="major", color=GRID, linewidth=0.8, alpha=0.9) + ax.grid(True, which="minor", color=GRID, linewidth=0.5, alpha=0.5) + ax.set_axisbelow(True) + for side in ("top", "right"): + ax.spines[side].set_visible(False) + for side in ("left", "bottom"): + ax.spines[side].set_color(GRID) + ax.tick_params(colors=TEXT_SECONDARY, labelsize=9) + + +def main() -> None: + ap = argparse.ArgumentParser() + ap.add_argument("--csv", type=Path, required=True) + ap.add_argument("--out", type=Path, required=True) + ap.add_argument("--title", default=None) + args = ap.parse_args() + + data = load(args.csv) + models = [m for m in ("tfim", "heisenberg") if m in data] + fig, axes = plt.subplots( + len(models), + 3, + figsize=(15.5, 4.4 * len(models)), + squeeze=False, + facecolor=SURFACE, + ) + + for r, model in enumerate(models): + per_lib = data[model] + ax_t, ax_s, ax_n = axes[r] + widths = sorted({n for lib in per_lib.values() for n in lib}) + term_series: dict[str, tuple[int, ...]] = {} + for key, (label, color, marker) in SERIES.items(): + if key not in per_lib: + continue + ns = sorted(per_lib[key]) + times = [float(per_lib[key][n]["time_s"]) for n in ns] + ax_t.plot( + ns, + times, + color=color, + marker=marker, + markersize=6, + linewidth=2, + label=label, + ) + terms = [int(per_lib[key][n]["terms"]) for n in ns] + ax_n.plot( + ns, + terms, + color=color, + marker=marker, + markersize=6, + linewidth=2, + label=label, + ) + term_series[key] = tuple(terms) + if key != BASELINE and BASELINE in per_lib: + shared = [n for n in ns if n in per_lib[BASELINE]] + ratios = [ + float(per_lib[key][n]["time_s"]) + / float(per_lib[BASELINE][n]["time_s"]) + for n in shared + ] + ax_s.plot( + shared, + ratios, + color=color, + marker=marker, + markersize=6, + linewidth=2, + label=label, + ) + + # Engines whose support is identical draw the same line, so the ones + # underneath are invisible. Say so rather than let the panel imply that + # only the top series was measured. + if BASELINE in term_series: + same = [k for k, v in term_series.items() if v == term_series[BASELINE]] + if len(same) > 1: + ax_n.annotate( + "identical (exactly):\n" + + "\n".join(SERIES[k][0].split(" (")[0] for k in same), + # The support rises left-to-right, so the low-right corner + # is the one that stays clear of the marks. + xy=(0.97, 0.04), + xycoords="axes fraction", + fontsize=8, + color=TEXT_SECONDARY, + ha="right", + va="bottom", + ) + + ax_s.axhline(1.0, color=SERIES[BASELINE][1], linewidth=2, linestyle=(0, (4, 3))) + ax_s.annotate( + f"{SERIES[BASELINE][0]} = 1", + xy=(0.02, 1.0), + xycoords=("axes fraction", "data"), + va="bottom", + fontsize=8, + color=TEXT_SECONDARY, + ) + + for ax, ylabel, title in ( + (ax_t, "runtime (s, min of repeats)", "runtime"), + (ax_s, f"× {BASELINE}", "relative runtime"), + (ax_n, "terms in the propagated operator", "workload size"), + ): + style(ax) + ax.set_yscale("log") + # Qubit counts are integers; let matplotlib interpolate ticks and it + # invents 6.25 qubits. + ax.set_xticks(widths) + ax.set_xticklabels([str(n) for n in widths]) + ax.set_xlabel("qubits", color=TEXT_SECONDARY, fontsize=9) + ax.set_ylabel(ylabel, color=TEXT_SECONDARY, fontsize=9) + ax.set_title(title, color=TEXT_PRIMARY, fontsize=10, loc="left") + + ax_t.set_ylabel("runtime (s, min of repeats)", color=TEXT_SECONDARY, fontsize=9) + ax_t.annotate( + MODEL_TITLES.get(model, model), + xy=(0, 1.16), + xycoords="axes fraction", + fontsize=12, + fontweight="bold", + color=TEXT_PRIMARY, + ) + + # One figure-level legend above the panels: repeating it per row wastes + # space, and inside the axes it lands on the fastest series. + handles, labels = axes[0][0].get_legend_handles_labels() + fig.legend( + handles, + labels, + frameon=False, + fontsize=9.5, + labelcolor=TEXT_SECONDARY, + ncol=len(labels), + loc="upper left", + bbox_to_anchor=(0.006, 0.995), + ) + + sample = next(iter(next(iter(data.values())).values())) + row = next(iter(sample.values())) + n_engines = len({lib for per_lib in data.values() for lib in per_lib}) + caption = ( + f"first-order Trotter, steps={row['steps']}, dt={row['dt']}, " + f"truncation |c| < {row['atol']}; single-threaded; min of repeats. " + f"All {n_engines} engines are validated to propagate the identical operator " + "term-for-term before timing, but pauli-prop and monoprop prune at " + "branch-creation time rather than after accumulation, so at this " + "truncation they carry a different support — compare the workload panel " + "before reading their ratios. monoprop is capped to one thread " + "(monoprop_NUM_THREADS=1); uncapped it takes one partition per core. " + "Per-point numbers in results.csv." + ) + if args.title: + fig.suptitle( + args.title, fontsize=14, color=TEXT_PRIMARY, x=0.006, y=0.995, ha="left" + ) + # Legend sits just under the title when there is one. + fig.legends[0].set_bbox_to_anchor((0.006, 0.962)) + wrapped = textwrap.fill(caption, width=178) + fig.text( + 0.008, + 0.004, + wrapped, + fontsize=8.5, + color=TEXT_SECONDARY, + ha="left", + va="bottom", + ) + top = 0.90 if args.title else 0.945 + fig.tight_layout(rect=(0, 0.018 * (wrapped.count("\n") + 1) + 0.012, 1, top)) + args.out.parent.mkdir(parents=True, exist_ok=True) + fig.savefig(args.out, dpi=160, facecolor=SURFACE) + print(f"wrote {args.out}") + + +if __name__ == "__main__": + main() diff --git a/benchmarks/cross-library/run_xbench.py b/benchmarks/cross-library/run_xbench.py new file mode 100644 index 000000000..cf1b9c5b0 --- /dev/null +++ b/benchmarks/cross-library/run_xbench.py @@ -0,0 +1,366 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +"""Drive the cross-library Pauli-propagation benchmark and merge the results. + +Runs the same two circuits — TFIM Trotter and Heisenberg-model correlations — +through `ppvm`, PauliPropagation.jl, PauliStrings.jl, Qiskit's `pauli-prop` and +monoprop, all from one parameter contract, and writes one tidy CSV. + +Before it will report a single timing it **validates**: every engine dumps its +propagated support at a small width and the driver diffs them term-for-term +against `ppvm`. A cross-library benchmark whose engines are quietly computing +different things is worse than no benchmark, and every bug this check caught +while the harness was being written produced plausible-looking numbers that +agreed on the observable to nine digits — a reversed circuit in the `pauli-prop` +runner, duplicate Paulis collapsing in its readout, and monoprop's `exp(+iθH)` +sign convention. + + uv run --no-project python3 run_xbench.py --help +""" + +from __future__ import annotations + +import argparse +import csv +import os +import shutil +import subprocess +import sys +from dataclasses import dataclass, field +from pathlib import Path + +REPO = Path(__file__).resolve().parents[2] +JULIA_PROJECT = REPO / "julia-benchmarks" +BASELINE = "ppvm" +CSV_COLUMNS = [ + "model", + "library", + "qubits", + "steps", + "dt", + "atol", + "time_s", + "terms", + "observable", +] + + +@dataclass(frozen=True) +class Runner: + """One engine: how to invoke it, and whether it is available here.""" + + name: str + argv: list[str] + needs: str # the executable that must exist on PATH + env: dict[str, str] = field(default_factory=dict) + + def available(self) -> bool: + return shutil.which(self.needs) is not None + + +def _here(name: str) -> str: + return str(Path(__file__).with_name(name)) + + +RUNNERS = { + "ppvm": Runner( + "ppvm", + [ + "cargo", + "run", + "--release", + "-q", + "-p", + "ppvm-pauli-sum", + "--example", + "xbench", + ], + "cargo", + ), + "pauli-propagation-jl": Runner( + "pauli-propagation-jl", + [ + "julia", + f"--project={JULIA_PROJECT}", + "-t1", + str(JULIA_PROJECT / "benches" / "xbench_pp.jl"), + ], + "julia", + ), + "pauli-strings-jl": Runner( + "pauli-strings-jl", + [ + "julia", + f"--project={JULIA_PROJECT}", + "-t1", + str(JULIA_PROJECT / "benches" / "xbench_ps.jl"), + ], + "julia", + ), + "pauli-prop": Runner( + "pauli-prop", + [ + "uv", + "run", + "--no-project", + "--with", + "pauli-prop", + "python3", + _here("xbench_qiskit.py"), + ], + "uv", + ), + # monoprop takes one serial partition per physical core when left alone, so + # the cap is part of the invocation, not an optional extra. The runner + # re-asserts it from its own CPU/wall ratio and fails if it did not land. + "monoprop": Runner( + "monoprop", + [ + "uv", + "run", + "--no-project", + "--with", + "monoprop", + "python3", + _here("xbench_monoprop.py"), + ], + "uv", + {"monoprop_NUM_THREADS": "1", "monoprop_PARTITIONS": "off"}, + ), +} + + +def invoke(runner: Runner, env_extra: dict[str, str], quiet: bool) -> str: + """Run one engine and return its stdout, streaming its stderr as progress.""" + env = os.environ.copy() + env.update(env_extra) + env.update(runner.env) + proc = subprocess.run( + runner.argv, + cwd=REPO, + env=env, + capture_output=True, + text=True, + check=False, + ) + if proc.returncode != 0: + sys.stderr.write(proc.stderr) + raise SystemExit(f"{runner.name} failed with exit code {proc.returncode}") + if not quiet and proc.stderr.strip(): + sys.stderr.write(proc.stderr) + return proc.stdout + + +def parse_dump(text: str) -> dict[str, float]: + """Parse a `DUMP=1` support listing into `{word: coefficient}`.""" + terms: dict[str, float] = {} + for line in text.splitlines(): + if not line or line.startswith("#"): + continue + word, coeff = line.split() + terms[word] = float(coeff) + return terms + + +def validate( + libs: list[str], params: dict[str, str], atol_cmp: float, quiet: bool +) -> None: + """Assert every engine propagates the *same* operator, term for term.""" + print("validating: all engines must agree term-for-term", file=sys.stderr) + for model in ("tfim", "heisenberg"): + env = dict(params) + env.update( + {"MODEL": model, "QUBITS": "4", "STEPS": "3", "ATOL": "1e-14", "DUMP": "1"} + ) + reference: dict[str, float] | None = None + ref_lib = "" + for lib in libs: + terms = parse_dump(invoke(RUNNERS[lib], env, quiet=True)) + if reference is None: + reference, ref_lib = terms, lib + print( + f" {model}: {lib} — {len(terms)} terms (reference)", + file=sys.stderr, + ) + continue + missing = set(reference) - set(terms) + extra = set(terms) - set(reference) + worst = max( + (abs(reference[w] - terms[w]) for w in set(reference) & set(terms)), + default=0.0, + ) + if missing or extra or worst > atol_cmp: + raise SystemExit( + f"{model}: {lib} disagrees with {ref_lib} — " + f"{len(terms)} vs {len(reference)} terms, " + f"{len(missing)} missing, {len(extra)} extra, max|Δc|={worst:.3e} " + f"(tolerance {atol_cmp:.0e}). Refusing to report timings." + ) + print( + f" {model}: {lib} — {len(terms)} terms, max|Δc|={worst:.1e} OK", + file=sys.stderr, + ) + print("validation passed\n", file=sys.stderr) + + +def main() -> None: + ap = argparse.ArgumentParser( + description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter + ) + ap.add_argument("--model", default="both", choices=["tfim", "heisenberg", "both"]) + ap.add_argument( + "--qubits", default="8,12,16,20,24,28,32", help="applies to every model" + ) + ap.add_argument( + "--qubits-tfim", + default=None, + help="override --qubits for the TFIM sweep (its support grows linearly in n, " + "so it reaches far wider systems than the Heisenberg one at equal cost)", + ) + ap.add_argument( + "--qubits-heisenberg", default=None, help="override --qubits for Heisenberg" + ) + ap.add_argument("--steps", type=int, default=10) + ap.add_argument("--dt", type=float, default=0.1) + ap.add_argument("--j", type=float, default=1.0, help="bond coupling J") + ap.add_argument("--h", type=float, default=1.0, help="field strength h") + ap.add_argument( + "--atol", type=float, default=1e-6, help="coefficient truncation threshold" + ) + ap.add_argument( + "--iters", type=int, default=3, help="timed repeats; the minimum is reported" + ) + ap.add_argument("--libs", default=",".join(RUNNERS)) + ap.add_argument("--out", default="target/xbench", type=Path) + ap.add_argument("--skip-validate", action="store_true") + ap.add_argument( + "--validate-tol", + type=float, + default=1e-10, + help="absolute coefficient tolerance for the cross-engine agreement check", + ) + ap.add_argument( + "--reuse", + type=Path, + default=None, + help="re-print the summary from an existing results.csv", + ) + ap.add_argument("-q", "--quiet", action="store_true") + args = ap.parse_args() + + if args.reuse is not None: + with args.reuse.open() as fh: + summarize(list(csv.DictReader(fh))) + return + + libs = [lib.strip() for lib in args.libs.split(",") if lib.strip()] + for lib in libs: + if lib not in RUNNERS: + raise SystemExit(f"unknown library {lib!r}; known: {', '.join(RUNNERS)}") + missing = [lib for lib in libs if not RUNNERS[lib].available()] + if missing: + print( + f"skipping (missing {', '.join(RUNNERS[m].needs for m in missing)}): " + f"{', '.join(missing)}", + file=sys.stderr, + ) + libs = [lib for lib in libs if lib not in missing] + if not libs: + raise SystemExit("no runnable engines") + + params = { + "DT": repr(args.dt), + "JCOUP": repr(args.j), + "HFIELD": repr(args.h), + "ITERS": str(args.iters), + } + if not args.skip_validate: + validate(libs, params, args.validate_tol, args.quiet) + + models = ["tfim", "heisenberg"] if args.model == "both" else [args.model] + per_model_qubits = { + "tfim": args.qubits_tfim or args.qubits, + "heisenberg": args.qubits_heisenberg or args.qubits, + } + rows: list[dict[str, str]] = [] + for model in models: + for lib in libs: + env = dict(params) + env.update( + { + "MODEL": model, + "QUBITS": per_model_qubits[model], + "STEPS": str(args.steps), + "ATOL": repr(args.atol), + } + ) + print(f"running {model} / {lib}", file=sys.stderr) + out = invoke(RUNNERS[lib], env, args.quiet) + reader = csv.DictReader(line for line in out.splitlines() if line.strip()) + rows.extend(dict(row) for row in reader) + + out_dir = REPO / args.out if not args.out.is_absolute() else args.out + out_dir.mkdir(parents=True, exist_ok=True) + csv_path = out_dir / "results.csv" + with csv_path.open("w", newline="") as fh: + writer = csv.DictWriter(fh, fieldnames=CSV_COLUMNS) + writer.writeheader() + writer.writerows(rows) + + summarize(rows) + print(f"\nwrote {csv_path}", file=sys.stderr) + + +def summarize(rows: list[dict[str, str]]) -> None: + """Print a time table plus the speedup of every engine relative to `ppvm`.""" + by_model: dict[str, dict[int, dict[str, dict[str, str]]]] = {} + for row in rows: + by_model.setdefault(row["model"], {}).setdefault(int(row["qubits"]), {})[ + row["library"] + ] = row + + for model, per_n in by_model.items(): + libs = sorted({lib for cells in per_n.values() for lib in cells}) + print(f"\n=== {model} ===") + header = f"{'n':>4} {'terms':>9} " + " ".join(f"{lib:>22}" for lib in libs) + print(header) + print("-" * len(header)) + for n in sorted(per_n): + cells = per_n[n] + base = cells.get(BASELINE) + terms = base["terms"] if base else next(iter(cells.values()))["terms"] + line = f"{n:>4} {terms:>9} " + parts = [] + for lib in libs: + row = cells.get(lib) + if row is None: + parts.append(f"{'—':>22}") + continue + t = float(row["time_s"]) + if base is not None and lib != BASELINE: + ratio = t / float(base["time_s"]) + parts.append(f"{t:>11.4f}s ({ratio:5.2f}x)") + else: + parts.append(f"{t:>11.4f}s {'(1.00x)':>9}") + print(line + " ".join(parts)) + + # A runtime ratio only means something at equal work. Any engine whose + # truncation rule leaves it carrying a materially different support is + # called out here rather than left to the reader to spot in the plot. + drift: dict[str, list[str]] = {} + for n in sorted(per_n): + base = per_n[n].get(BASELINE) + if base is None: + continue + for lib, row in per_n[n].items(): + if lib == BASELINE: + continue + rel = int(row["terms"]) / int(base["terms"]) - 1.0 + if abs(rel) > 0.02: + drift.setdefault(lib, []).append(f"n={n}: {rel:+.0%}") + for lib, notes in drift.items(): + print(f" note: {lib} support vs {BASELINE} — {', '.join(notes)}") + + +if __name__ == "__main__": + main() diff --git a/benchmarks/cross-library/xbench_monoprop.py b/benchmarks/cross-library/xbench_monoprop.py new file mode 100644 index 000000000..70c24c3ca --- /dev/null +++ b/benchmarks/cross-library/xbench_monoprop.py @@ -0,0 +1,243 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +"""The monoprop side of the cross-library Pauli-propagation benchmark. + +See `README.md` for the shared circuit definitions, the parameter contract, and +the CSV schema — every runner reads the same environment variables and prints the +same columns. + +Three things about monoprop need saying, because each of them silently changes +the number this file prints: + +**It is parallel by default.** With one MPI rank the engine takes one serial +partition per *physical core* (`resolve_partition_count_`), so an uncapped run on +this 14-core machine is a ~9.5x-CPU, ~3x-wall-faster run than the serial one, and +is not comparable with the other single-threaded engines. `monoprop_NUM_THREADS=1` +and `monoprop_PARTITIONS=off` pin it to one partition, and the runner **measures +its own CPU/wall ratio and fails** if the cap did not take. An environment +variable that is read once into a cached C++ static is exactly the kind of setting +that goes stale without anyone noticing. + +**Its angles are not the spec's angles.** `ExpGate` applies `exp(+iθH)`, so the +spec's `exp(-iθ/2·G)` needs `θ_monoprop = -θ_spec/2`. All four other engines take +`θ_spec` directly. + +**It reports more rows than it has terms.** The engine retains monomials whose +coefficient has cancelled to exactly zero, and `size()` counts those. At +`n=12` Heisenberg that is 3.06M tracked rows against 1.32M terms above the +truncation threshold. The `terms` column is the above-threshold support, which is +what the other four engines mean by it; the tracked size goes to stderr. + + MODEL=tfim QUBITS=8,16,24 STEPS=10 DT=0.1 JCOUP=1.0 HFIELD=1.0 ATOL=1e-6 \ + uv run --no-project --with monoprop python3 xbench_monoprop.py +""" + +from __future__ import annotations + +import os + +# Must precede the import: the C++ side reads these once into a cached static +# (`monoprop::config::get()`), so a later assignment would be ignored. The driver +# sets them in the subprocess environment too; this makes a direct invocation of +# this file single-threaded as well, rather than quietly parallel. +os.environ.setdefault("monoprop_NUM_THREADS", "1") +os.environ.setdefault("monoprop_PARTITIONS", "off") + +import resource +import sys +import time + +from monoprop import ( + Circuit, + ExpGate, + Pauli, + PauliOperator, + PauliPropagator, +) + +MODEL = os.environ.get("MODEL", "tfim") +STEPS = int(os.environ.get("STEPS", "10")) +DT = float(os.environ.get("DT", "0.1")) +JCOUP = float(os.environ.get("JCOUP", "1.0")) +HFIELD = float(os.environ.get("HFIELD", "1.0")) +ATOL = float(os.environ.get("ATOL", "1e-6")) +ITERS = int(os.environ.get("ITERS", "3")) +# A serial run has CPU/wall ~1. Anything materially above it means the thread cap +# did not take and the timing is not comparable with the other engines. +CPU_WALL_MAX = float(os.environ.get("CPU_WALL_MAX", "1.5")) + +THETA_BOND = 2 * JCOUP * DT +THETA_SITE = 2 * HFIELD * DT + + +def cpu_seconds() -> float: + """Total CPU time (user+sys) charged to this process and its children.""" + me = resource.getrusage(resource.RUSAGE_SELF) + kids = resource.getrusage(resource.RUSAGE_CHILDREN) + return me.ru_utime + me.ru_stime + kids.ru_utime + kids.ru_stime + + +def seed_operator(n: int) -> PauliOperator: + """`Σ_i Z_i` for TFIM, `Z_0` for Heisenberg.""" + if MODEL == "tfim": + terms: dict[Pauli | str, float] = {Pauli("Z", (i,)): 1.0 for i in range(n)} + else: + terms = {Pauli("Z", (0,)): 1.0} + return PauliOperator(terms, n) + + +def gate_list(n: int) -> list[tuple[str, tuple[int, ...], float]]: + """`STEPS` first-order Trotter steps in *application* order, per the spec.""" + gates: list[tuple[str, tuple[int, ...], float]] = [] + for _ in range(STEPS): + if MODEL == "tfim": + gates += [("X", (i,), THETA_SITE) for i in range(n)] + gates += [("ZZ", (i, i + 1), THETA_BOND) for i in range(n - 1)] + else: + for i in range(n - 1): + gates += [ + ("XX", (i, i + 1), THETA_BOND), + ("YY", (i, i + 1), THETA_BOND), + ("ZZ", (i, i + 1), THETA_BOND), + ] + gates += [("Z", (i,), THETA_SITE) for i in range(n)] + return gates + + +def build_circuit(n: int) -> Circuit: + """The spec's gate sequence as a monoprop `Circuit`. + + Two conversions happen here, and the term-for-term validation in + `run_xbench.py` is what pinned both of them down: + + * **Angle.** `ExpGate` applies `exp(+iθH)` (its own docstring flags the + positive sign as the difference from Qiskit's `r

`), while the shared + spec's gate is `exp(-iθ_spec/2 · G)`. Hence `θ = -θ_spec/2`. + * **Order.** Like Qiskit's `pauli-prop`, the propagator conjugates from the + end of the gate list backwards, so the spec's application order is + obtained by appending in reverse. + + Neither is cosmetic. All eight sign/order combinations were checked against + the reference dump: forward order loses terms outright (109 against TFIM's + 124), and every wrong angle keeps the right support while moving + coefficients by up to 1.3. Only this one lands within 5e-13. + """ + gates = list(reversed(gate_list(n))) + exp_gates = [ + ExpGate(PauliOperator({Pauli(string, qubits): 1.0}, n), index) + for index, (string, qubits, _theta) in enumerate(gates) + ] + parameters = [-0.5 * theta for (_s, _q, theta) in gates] + return Circuit(gates=exp_gates, parameters=parameters) + + +def propagate(n: int, circuit: Circuit) -> PauliPropagator: + """Heisenberg-picture propagation with coefficient-only truncation. + + `cutoff` is monoprop's mandatory bound on retained Pauli *weight*, which the + other four engines have no analogue of; `cutoff=n` is the whole register, so + it never binds and `lower_atol` is left as the only truncation — the shared + `|c| < atol` rule. + """ + return PauliPropagator.from_circuit( + circuit, seed_operator(n), cutoff=n, lower_atol=ATOL + ) + + +def support(mp: PauliPropagator, n: int) -> dict[str, float]: + """`{word: coefficient}` with site 0 leftmost, exact-zero rows dropped. + + `atol=0.0` asks the engine for everything it holds, and the threshold is + applied here so the filter is the same one the other engines apply during + propagation. + """ + operator = mp.evolved_operator(atol=0.0) + terms: dict[str, float] = {} + for pauli, coeff in operator.terms.items(): + value = float(coeff.real if isinstance(coeff, complex) else coeff) + if abs(value) < ATOL: + continue + chars = ["I"] * n + for qubit, letter in zip(pauli.qubits, pauli.string): + chars[qubit] = letter + terms["".join(chars)] = value + return terms + + +def readout(terms: dict[str, float], n: int) -> float: + """`⟨0…0|O|0…0⟩` for TFIM; the `Z_0` autocorrelator for Heisenberg.""" + if MODEL == "tfim": + # ⟨0|Z|0⟩ = 1 and ⟨0|X|0⟩ = ⟨0|Y|0⟩ = 0, so only the X/Y-free terms survive. + return sum(c for w, c in terms.items() if "X" not in w and "Y" not in w) + return terms.get("Z" + "I" * (n - 1), 0.0) + + +def dump(n: int) -> None: + """Print the propagated support as `word coefficient`, largest first.""" + terms = support(propagate(n, build_circuit(n)), n) + print(f"# {len(terms)} terms") + for word in sorted(terms, key=lambda w: (-abs(terms[w]), w)): + print(f"{word} {terms[word]:+.12e}") + + +def main() -> None: + qubits = [ + int(t) for t in os.environ.get("QUBITS", "8,12,16,20,24,28,32").split(",") + ] + if os.environ.get("DUMP"): + dump(qubits[0]) + return + print( + f"monoprop {MODEL}: steps={STEPS} dt={DT} J={JCOUP} h={HFIELD} atol={ATOL} " + f"iters={ITERS} threads={os.environ['monoprop_NUM_THREADS']} " + f"partitions={os.environ['monoprop_PARTITIONS']} cores={os.cpu_count()}", + file=sys.stderr, + ) + print("model,library,qubits,steps,dt,atol,time_s,terms,observable") + for n in qubits: + # Circuit construction is not propagation; built once, outside the timed + # region, as in every other runner. monoprop's `propagate` does expand + # the gate list internally each call, but that is 0.0–2.2% of the total + # at these widths (measured), not enough to move a ratio. + circuit = build_circuit(n) + + best = float("inf") + cpu_at_best = 0.0 + propagator = None + for _ in range(ITERS): + cpu0, wall0 = cpu_seconds(), time.perf_counter() + propagator = propagate(n, circuit) + wall = time.perf_counter() - wall0 + if wall < best: + best, cpu_at_best = wall, cpu_seconds() - cpu0 + + # A parallel run would report a wall time the other engines never had the + # chance to compete with. Refuse to print it. + ratio = cpu_at_best / best if best > 0 else 1.0 + if ratio > CPU_WALL_MAX: + raise SystemExit( + f"n={n}: monoprop used {ratio:.2f}x CPU per wall second " + f"({cpu_at_best:.3f}s CPU in {best:.3f}s wall), so the thread cap " + f"did not take and this timing is not comparable with the other " + f"single-threaded engines. Check monoprop_NUM_THREADS / " + f"monoprop_PARTITIONS reach the subprocess." + ) + + assert propagator is not None + # Materializing the support is O(size) in Python and dwarfs the + # propagation at the wide end, so it happens once, after timing. + terms = support(propagator, n) + obs = readout(terms, n) + print( + f"{MODEL},monoprop,{n},{STEPS},{DT},{ATOL},{best:.7g},{len(terms)},{obs!r}" + ) + print( + f" n={n} {best:.4f}s {len(terms)} terms obs={obs} " + f"[cpu/wall={ratio:.2f} tracked={propagator.size()}]", + file=sys.stderr, + ) + sys.stdout.flush() + + +if __name__ == "__main__": + main() diff --git a/benchmarks/cross-library/xbench_qiskit.py b/benchmarks/cross-library/xbench_qiskit.py new file mode 100644 index 000000000..6122a60e3 --- /dev/null +++ b/benchmarks/cross-library/xbench_qiskit.py @@ -0,0 +1,190 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +"""The Qiskit `pauli-prop` side of the cross-library Pauli-propagation benchmark. + +See `README.md` for the shared circuit definitions, the parameter contract, and +the CSV schema — every runner reads the same environment variables and prints +the same columns. + +`pauli-prop` truncates on `atol` *and* on a mandatory `max_terms` cap, which the +other three engines do not have. `MAX_TERMS` therefore defaults high enough to +be non-binding, and the runner asserts it never bound: if the propagated support +ever reaches the cap the row is not comparable and the run fails rather than +quietly reporting a differently-truncated number. + + MODEL=tfim QUBITS=8,16,24 STEPS=10 DT=0.1 JCOUP=1.0 HFIELD=1.0 ATOL=1e-6 \ + uv run --no-project --with pauli-prop python3 xbench_qiskit.py +""" + +from __future__ import annotations + +import os +import sys +import time + +import numpy as np +from pauli_prop import propagate_through_circuit +from qiskit import QuantumCircuit +from qiskit.quantum_info import SparsePauliOp + +MODEL = os.environ.get("MODEL", "tfim") +STEPS = int(os.environ.get("STEPS", "10")) +DT = float(os.environ.get("DT", "0.1")) +JCOUP = float(os.environ.get("JCOUP", "1.0")) +HFIELD = float(os.environ.get("HFIELD", "1.0")) +ATOL = float(os.environ.get("ATOL", "1e-6")) +ITERS = int(os.environ.get("ITERS", "3")) +MAX_TERMS = int(os.environ.get("MAX_TERMS", str(1 << 22))) + +THETA_BOND = 2 * JCOUP * DT +THETA_SITE = 2 * HFIELD * DT + + +def seed_operator(n: int) -> SparsePauliOp: + """`Σ_i Z_i` for TFIM, `Z_0` for Heisenberg. + + Qiskit Pauli labels are little-endian: the rightmost character is qubit 0. + """ + if MODEL == "tfim": + labels = ["I" * (n - 1 - i) + "Z" + "I" * i for i in range(n)] + else: + labels = ["I" * (n - 1) + "Z"] + return SparsePauliOp(labels, coeffs=np.ones(len(labels))) + + +def gate_list(n: int) -> list[tuple[str, int]]: + """`STEPS` first-order Trotter steps in *application* order, per the spec.""" + gates: list[tuple[str, int]] = [] + for _ in range(STEPS): + if MODEL == "tfim": + gates += [("rx", i) for i in range(n)] + gates += [("rzz", i) for i in range(n - 1)] + else: + for i in range(n - 1): + gates += [("rxx", i), ("ryy", i), ("rzz", i)] + gates += [("rz", i) for i in range(n)] + return gates + + +def build_circuit(n: int) -> QuantumCircuit: + """The spec's gate sequence, as a circuit `pauli-prop` will apply in order. + + A Qiskit circuit whose instructions are appended `g_1 … g_k` denotes the + unitary `U = g_k ⋯ g_1`, and in the Heisenberg frame `pauli-prop` computes + `U† O U` by conjugating **from the end of the instruction list backwards**. + So the gate `pauli-prop` applies first is the one appended last, and the + spec's application order is obtained by appending in reverse. + + This is not cosmetic. Appending in forward order propagates a genuinely + different operator: at `n=4, steps=3, atol=1e-14` it yields 108 terms + against TFIM's reference 124 and 61 against Heisenberg's 64, with + coefficients off by up to 0.1. Reversed, all four engines agree + term-for-term to the last bit — which is what `run_xbench.py --validate` + checks before it will report a timing. + """ + qc = QuantumCircuit(n) + for name, i in reversed(gate_list(n)): + if name == "rx": + qc.rx(THETA_SITE, i) + elif name == "rz": + qc.rz(THETA_SITE, i) + elif name == "rxx": + qc.rxx(THETA_BOND, i, i + 1) + elif name == "ryy": + qc.ryy(THETA_BOND, i, i + 1) + elif name == "rzz": + qc.rzz(THETA_BOND, i, i + 1) + else: # pragma: no cover - the list above is closed + raise ValueError(f"unknown gate {name}") + return qc + + +def support(op: SparsePauliOp) -> dict[str, float]: + """`{word: coefficient}` with site 0 leftmost, duplicate Paulis summed. + + `propagate_through_circuit` returns an operator that may list the same Pauli + more than once — `len(op)` counts rows, not distinct terms — so the support + size comparable with the other three engines is the size of this mapping. + Qiskit labels are little-endian, hence the reversal. + """ + terms: dict[str, float] = {} + for pauli, coeff in zip(op.paulis, np.asarray(op.coeffs)): + word = str(pauli)[::-1] + terms[word] = terms.get(word, 0.0) + float(np.real(coeff)) + return terms + + +def readout(op: SparsePauliOp, n: int) -> float: + """`⟨0…0|O|0…0⟩` for TFIM; the `Z_0` autocorrelator for Heisenberg.""" + paulis = op.paulis + coeffs = np.asarray(op.coeffs) + if MODEL == "tfim": + # ⟨0|Z|0⟩ = 1 and ⟨0|X|0⟩ = ⟨0|Y|0⟩ = 0, so only the X-free terms survive. + diagonal = ~paulis.x.any(axis=1) + return float(np.real(coeffs[diagonal].sum())) + z0_x = np.zeros(n, dtype=bool) + z0_z = np.zeros(n, dtype=bool) + z0_z[0] = True + hit = (paulis.x == z0_x).all(axis=1) & (paulis.z == z0_z).all(axis=1) + return float(np.real(coeffs[hit].sum())) + + +def dump(n: int) -> None: + """Print the propagated support as `word coefficient`, largest first. + + Site 0 leftmost, matching the Rust runner's `DUMP=1` output so the driver + can diff the two term-for-term. + """ + out, _bias = propagate_through_circuit( + seed_operator(n), build_circuit(n), max_terms=MAX_TERMS, atol=ATOL, frame="h" + ) + terms = support(out) + print(f"# {len(terms)} terms") + for word in sorted(terms, key=lambda w: (-abs(terms[w]), w)): + print(f"{word} {terms[word]:+.12e}") + + +def main() -> None: + qubits = [ + int(t) for t in os.environ.get("QUBITS", "8,12,16,20,24,28,32").split(",") + ] + if os.environ.get("DUMP"): + dump(qubits[0]) + return + print( + f"pauli-prop {MODEL}: steps={STEPS} dt={DT} J={JCOUP} h={HFIELD} " + f"atol={ATOL} iters={ITERS} max_terms={MAX_TERMS}", + file=sys.stderr, + ) + print("model,library,qubits,steps,dt,atol,time_s,terms,observable") + for n in qubits: + # Circuit construction is not propagation; built once, outside the + # timed region, as in every other runner. + circuit = build_circuit(n) + seed = seed_operator(n) + + best = float("inf") + terms = 0 + obs = float("nan") + for _ in range(ITERS): + t0 = time.perf_counter() + out, _bias = propagate_through_circuit( + seed, circuit, max_terms=MAX_TERMS, atol=ATOL, frame="h" + ) + best = min(best, time.perf_counter() - t0) + # Distinct Paulis, not the returned operator's row count. + terms = len(support(out)) + obs = readout(out, n) + if terms >= MAX_TERMS: + raise SystemExit( + f"n={n}: support hit the max_terms cap ({terms} >= {MAX_TERMS}); " + "raise MAX_TERMS — this row would be truncated differently from " + "the other engines and is not comparable" + ) + print(f"{MODEL},pauli-prop,{n},{STEPS},{DT},{ATOL},{best:.7g},{terms},{obs!r}") + print(f" n={n} {best:.4f}s {terms} terms obs={obs}", file=sys.stderr) + sys.stdout.flush() + + +if __name__ == "__main__": + main() diff --git a/crates/ppvm-pauli-sum/examples/xbench.rs b/crates/ppvm-pauli-sum/examples/xbench.rs new file mode 100644 index 000000000..0f3c76d04 --- /dev/null +++ b/crates/ppvm-pauli-sum/examples/xbench.rs @@ -0,0 +1,266 @@ +// SPDX-FileCopyrightText: 2026 The PPVM Authors +// SPDX-License-Identifier: Apache-2.0 + +//! The `ppvm` side of the cross-library Pauli-propagation benchmark. +//! +//! Two Trotter workloads — TFIM magnetization and a Heisenberg autocorrelator — +//! propagated in the Heisenberg picture with a coefficient-magnitude truncation +//! after every gate. See `benchmarks/cross-library/README.md` for the shared +//! spec: the gate order, the `θ = 2·c·dt` convention, the environment contract, +//! and the CSV schema that every runner in that harness prints. +//! +//! ```bash +//! MODEL=tfim QUBITS=8,16,24 STEPS=10 DT=0.1 JCOUP=1.0 HFIELD=1.0 ATOL=1e-6 \ +//! cargo run --release -p ppvm-pauli-sum --example xbench +//! ``` + +use std::time::Instant; + +use ppvm_pauli_sum::prelude::*; +use ppvm_pauli_sum::strategy::CoefficientThreshold; + +#[derive(Clone, Copy)] +struct Params { + model: Model, + steps: usize, + dt: f64, + j: f64, + h: f64, + atol: f64, + iters: usize, +} + +#[derive(Clone, Copy, PartialEq, Eq)] +enum Model { + Tfim, + Heisenberg, +} + +impl Model { + fn parse(s: &str) -> Self { + match s { + "tfim" => Model::Tfim, + "heisenberg" => Model::Heisenberg, + other => panic!("unknown MODEL {other:?} (expected `tfim` or `heisenberg`)"), + } + } + + fn name(self) -> &'static str { + match self { + Model::Tfim => "tfim", + Model::Heisenberg => "heisenberg", + } + } +} + +fn env_f64(key: &str, default: f64) -> f64 { + std::env::var(key) + .ok() + .and_then(|s| s.parse().ok()) + .unwrap_or(default) +} + +fn env_usize(key: &str, default: usize) -> usize { + std::env::var(key) + .ok() + .and_then(|s| s.parse().ok()) + .unwrap_or(default) +} + +impl Params { + fn from_env() -> Self { + Params { + model: Model::parse(&std::env::var("MODEL").unwrap_or_else(|_| "tfim".to_string())), + steps: env_usize("STEPS", 10), + dt: env_f64("DT", 0.1), + j: env_f64("JCOUP", 1.0), + h: env_f64("HFIELD", 1.0), + atol: env_f64("ATOL", 1e-6), + iters: env_usize("ITERS", 3), + } + } +} + +/// Smallest power-of-two byte width that holds `n` qubits with room to spare. +fn storage_bytes(n: usize) -> usize { + let need = n.div_ceil(8); + let mut k = 0; + while (1usize << k) <= need { + k += 1; + } + 1usize << k +} + +/// A single-site Pauli word `p` at `site`, identity elsewhere. +fn site_word(n: usize, site: usize, p: char) -> String { + (0..n).map(|j| if j == site { p } else { 'I' }).collect() +} + +type Cfg = config::fxhash::Byte>; + +/// Build the seed observable for `model` on `n` sites. +fn seed(n: usize, p: Params) -> PauliSum> { + let mut sum: PauliSum> = PauliSum::builder() + .n_qubits(n) + .strategy(CoefficientThreshold(p.atol)) + .capacity(1 << 12) + .build(); + match p.model { + Model::Tfim => { + for i in 0..n { + sum += (PauliWord::from(site_word(n, i, 'Z').as_str()), 1.0); + } + } + Model::Heisenberg => { + sum += (PauliWord::from(site_word(n, 0, 'Z').as_str()), 1.0); + } + } + sum +} + +/// Propagate the shared gate sequence through `state`. +fn propagate(state: &mut PauliSum>, n: usize, p: Params) { + let theta_bond = 2.0 * p.j * p.dt; + let theta_site = 2.0 * p.h * p.dt; + for _ in 0..p.steps { + match p.model { + Model::Tfim => { + for i in 0..n { + state.rx(i, theta_site); + state.truncate(); + } + for i in 0..n.saturating_sub(1) { + state.rzz(i, i + 1, theta_bond); + state.truncate(); + } + } + Model::Heisenberg => { + for i in 0..n.saturating_sub(1) { + state.rxx(i, i + 1, theta_bond); + state.truncate(); + state.ryy(i, i + 1, theta_bond); + state.truncate(); + state.rzz(i, i + 1, theta_bond); + state.truncate(); + } + for i in 0..n { + state.rz(i, theta_site); + state.truncate(); + } + } + } + } +} + +/// `⟨0…0|O|0…0⟩` for TFIM; the `Z_0` coefficient for Heisenberg. +/// +/// The diagonal contraction is spelled out rather than delegated: `⟨0|Z|0⟩ = 1` +/// and `⟨0|X|0⟩ = ⟨0|Y|0⟩ = 0`, so exactly the X-free terms survive and each +/// contributes its own coefficient. +fn readout(state: &PauliSum>, n: usize, p: Params) -> f64 { + match p.model { + Model::Tfim => state + .iter() + .filter(|(word, _)| (0..n).all(|i| !word.get_xbit(i))) + .map(|(_, c)| *c) + .sum(), + Model::Heisenberg => state + .data() + .get(&PauliWord::from(site_word(n, 0, 'Z').as_str())) + .copied() + .unwrap_or(0.0), + } +} + +/// Run one model at width `n`, returning `(best seconds, final support, observable)`. +fn run(n: usize, p: Params) -> (f64, usize, f64) { + let base = seed::(n, p); + let mut best = f64::INFINITY; + let mut terms = 0usize; + let mut observable = f64::NAN; + for _ in 0..p.iters { + let mut state = base.clone(); + let t0 = Instant::now(); + propagate(&mut state, n, p); + best = best.min(t0.elapsed().as_secs_f64()); + terms = state.len(); + observable = readout(&state, n, p); + } + (best, terms, observable) +} + +/// Print the whole propagated support as `word coefficient`, largest first — +/// the format the driver diffs across engines. +fn dump(n: usize, p: Params) { + let mut state = seed::(n, p); + propagate(&mut state, n, p); + let mut out: Vec<(String, f64)> = state.iter().map(|(k, c)| (k.to_string(), *c)).collect(); + out.sort_by(|a, b| { + b.1.abs() + .partial_cmp(&a.1.abs()) + .unwrap() + .then(a.0.cmp(&b.0)) + }); + println!("# {} terms", out.len()); + for (word, coeff) in out { + println!("{word} {coeff:+.12e}"); + } +} + +macro_rules! dispatch { + ($f:ident, $n:expr, $p:expr) => { + match storage_bytes($n) { + 2 => $f::<2>($n, $p), + 4 => $f::<4>($n, $p), + 8 => $f::<8>($n, $p), + 16 => $f::<16>($n, $p), + 32 => $f::<32>($n, $p), + 64 => $f::<64>($n, $p), + b => panic!("no storage tier for {b} bytes (n = {})", $n), + } + }; +} + +fn main() { + let p = Params::from_env(); + let qubits: Vec = std::env::var("QUBITS") + .ok() + .map(|s| s.split(',').filter_map(|t| t.trim().parse().ok()).collect()) + .unwrap_or_else(|| vec![8, 12, 16, 20, 24, 28, 32]); + + if std::env::var("DUMP").is_ok() { + let n = qubits[0]; + dispatch!(dump, n, p); + return; + } + + eprintln!( + "ppvm {}: steps={} dt={} J={} h={} atol={:e} iters={}", + p.model.name(), + p.steps, + p.dt, + p.j, + p.h, + p.atol, + p.iters + ); + + println!("model,library,qubits,steps,dt,atol,time_s,terms,observable"); + for &n in &qubits { + let (t, terms, obs) = dispatch!(run, n, p); + println!( + "{},ppvm,{},{},{},{:e},{:.6},{},{:.12e}", + p.model.name(), + n, + p.steps, + p.dt, + p.atol, + t, + terms, + obs + ); + eprintln!(" n={n:3} {t:9.4}s {terms:>9} terms obs={obs:+.9e}"); + use std::io::Write; + std::io::stdout().flush().ok(); + } +} diff --git a/julia-benchmarks/Manifest.toml b/julia-benchmarks/Manifest.toml index 3c050de35..d846d81c7 100644 --- a/julia-benchmarks/Manifest.toml +++ b/julia-benchmarks/Manifest.toml @@ -1,8 +1,8 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.12.1" +julia_version = "1.12.6" manifest_format = "2.0" -project_hash = "2774cd547152d7fbb5c2095d2e68e179a1b41160" +project_hash = "52f0dd51f0e051696bfc24b6c265049f3b4a890d" [[deps.AbstractFFTs]] deps = ["LinearAlgebra"] @@ -208,6 +208,11 @@ git-tree-sha1 = "37ea44092930b1811e666c3bc38065d7d87fcc74" uuid = "5ae59095-9a9b-59fe-a467-6f913c188581" version = "0.13.1" +[[deps.Combinatorics]] +git-tree-sha1 = "c761b00e7755700f9cdf5b02039939d1359330e1" +uuid = "861a8166-3701-5b0c-9a16-15d98fcdc6aa" +version = "1.1.0" + [[deps.Compat]] deps = ["TOML", "UUIDs"] git-tree-sha1 = "9d8a54ce4b17aa5bdce0ea5c34bc5e7c340d16ad" @@ -267,6 +272,12 @@ git-tree-sha1 = "9e2f36d3c96a820c678f2f1f1782582fcf685bae" uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab" version = "1.9.1" +[[deps.Dictionaries]] +deps = ["Indexing", "Random", "Serialization"] +git-tree-sha1 = "a55766a9c8f66cf19ffcdbdb1444e249bb4ace33" +uuid = "85a47980-9c8c-11e8-2b9f-f7ca1fa99fb4" +version = "0.4.6" + [[deps.Distances]] deps = ["LinearAlgebra", "Statistics", "StatsAPI"] git-tree-sha1 = "c7e3a542b999843086e2f29dac96a618c105be1d" @@ -307,7 +318,7 @@ version = "0.9.5" [[deps.Downloads]] deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" -version = "1.6.0" +version = "1.7.0" [[deps.EpollShim_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] @@ -478,6 +489,11 @@ git-tree-sha1 = "68c173f4f449de5b438ee67ed0c9c748dc31a2ec" uuid = "34004b35-14d8-5ef3-9330-4cdb6864b03a" version = "0.3.28" +[[deps.Indexing]] +git-tree-sha1 = "ce1566720fd6b19ff3411404d4b977acd4814f9f" +uuid = "313cdc1a-70c2-5d6a-ae34-0150d3930a38" +version = "1.1.1" + [[deps.Inflate]] git-tree-sha1 = "d1b1b796e47d94588b3757fe84fbf65a5ec4a80d" uuid = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9" @@ -617,7 +633,7 @@ version = "0.6.4" [[deps.LibCURL_jll]] deps = ["Artifacts", "LibSSH2_jll", "Libdl", "OpenSSL_jll", "Zlib_jll", "nghttp2_jll"] uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" -version = "8.11.1+1" +version = "8.15.0+0" [[deps.LibGit2]] deps = ["LibGit2_jll", "NetworkOptions", "Printf", "SHA"] @@ -744,7 +760,7 @@ version = "1.11.0" [[deps.MozillaCACerts_jll]] uuid = "14a3606d-f60d-562e-9121-12d972cd8159" -version = "2025.5.20" +version = "2025.11.4" [[deps.MuladdMacro]] git-tree-sha1 = "cac9cc5499c25554cba55cd3c30543cff5ca4fab" @@ -812,7 +828,7 @@ version = "1.6.1" [[deps.OpenSSL_jll]] deps = ["Artifacts", "Libdl"] uuid = "458c3c95-2e84-50aa-8efc-19380b2a3a95" -version = "3.5.1+0" +version = "3.5.4+0" [[deps.OpenSpecFun_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl"] @@ -876,6 +892,20 @@ version = "0.7.3" [deps.PauliPropagation.weakdeps] CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" +[[deps.PauliStrings]] +deps = ["BitIntegers", "Combinatorics", "Dictionaries", "LinearAlgebra", "ProgressBars", "Random", "SparseArrays"] +git-tree-sha1 = "80ecc963d62d97541258fbc8fcdf2fc95a2e70f7" +uuid = "f07625cc-80e4-4099-a362-36a3484d9bcc" +version = "1.10.1" + + [deps.PauliStrings.extensions] + MathLinkPauliStringsExt = "MathLink" + SymbolicsPauliStringsExt = "Symbolics" + + [deps.PauliStrings.weakdeps] + MathLink = "18c93696-a329-5786-9845-8443133fa0b4" + Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" + [[deps.Pixman_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "LLVMOpenMP_jll", "Libdl"] git-tree-sha1 = "e4a6721aa89e62e5d4217c0b21bd714263779dda" @@ -885,7 +915,7 @@ version = "0.46.4+0" [[deps.Pkg]] deps = ["Artifacts", "Dates", "Downloads", "FileWatching", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "Random", "SHA", "TOML", "Tar", "UUIDs", "p7zip_jll"] uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" -version = "1.12.0" +version = "1.12.1" weakdeps = ["REPL"] [deps.Pkg.extensions] @@ -963,6 +993,12 @@ git-tree-sha1 = "990016fb1508b0726a70039f39569720d054c78d" uuid = "efd6af41-a80b-495e-886c-e51b0c7d77a3" version = "0.1.7" +[[deps.ProgressBars]] +deps = ["Printf"] +git-tree-sha1 = "b437cdb0385ed38312d91d9c00c20f3798b30256" +uuid = "49802e3a-d2f1-5c88-81d8-b72133a6f568" +version = "1.5.1" + [[deps.PtrArrays]] git-tree-sha1 = "4fbbafbc6251b883f4d2705356f3641f3652a7fe" uuid = "43287f4e-b6f4-7ad1-bb20-aadabca52c3d" @@ -1580,9 +1616,9 @@ uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" version = "1.64.0+1" [[deps.p7zip_jll]] -deps = ["Artifacts", "Libdl"] +deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0" -version = "17.5.0+2" +version = "17.7.0+0" [[deps.x264_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] diff --git a/julia-benchmarks/Project.toml b/julia-benchmarks/Project.toml index 2a0f4ed61..befd11060 100644 --- a/julia-benchmarks/Project.toml +++ b/julia-benchmarks/Project.toml @@ -2,4 +2,5 @@ BenchmarkPlots = "ab8c0f59-4072-4e0d-8f91-a91e1495eb26" BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf" PauliPropagation = "293282d5-3c99-4fb6-92d0-fd3280a19750" +PauliStrings = "f07625cc-80e4-4099-a362-36a3484d9bcc" StatsPlots = "f3b207a7-027a-5e70-b257-86293d7955fd" diff --git a/julia-benchmarks/benches/xbench_pp.jl b/julia-benchmarks/benches/xbench_pp.jl new file mode 100644 index 000000000..15087b28b --- /dev/null +++ b/julia-benchmarks/benches/xbench_pp.jl @@ -0,0 +1,129 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +# +# The PauliPropagation.jl side of the cross-library Pauli-propagation +# benchmark. See `benchmarks/cross-library/README.md` for the shared circuit +# definitions, the parameter contract, and the CSV schema — every runner reads +# the same environment variables and prints the same columns. +# +# MODEL=tfim QUBITS=8,16,24 STEPS=10 DT=0.1 JCOUP=1.0 HFIELD=1.0 ATOL=1e-6 \ +# julia --project=@. -t1 benches/xbench_pp.jl + +using PauliPropagation +using Printf + +const MODEL = get(ENV, "MODEL", "tfim") +const STEPS = parse(Int, get(ENV, "STEPS", "10")) +const DT = parse(Float64, get(ENV, "DT", "0.1")) +const JCOUP = parse(Float64, get(ENV, "JCOUP", "1.0")) +const HFIELD = parse(Float64, get(ENV, "HFIELD", "1.0")) +const ATOL = parse(Float64, get(ENV, "ATOL", "1e-6")) +const ITERS = parse(Int, get(ENV, "ITERS", "3")) + +const THETA_BOND = 2 * JCOUP * DT +const THETA_SITE = 2 * HFIELD * DT + +"""The seed observable: `Σ_i Z_i` for TFIM, `Z_1` for Heisenberg.""" +function seed_operator(n::Int) + ps = PauliSum(n) + if MODEL == "tfim" + for i in 1:n + add!(ps, PauliString(n, [:Z], [i])) + end + else + add!(ps, PauliString(n, [:Z], [1])) + end + return ps +end + +"""One first-order Trotter step, gates in the order the shared spec fixes.""" +function trotter_step!(state, n::Int) + if MODEL == "tfim" + for i in 1:n + state = propagate(PauliRotation([:X], [i], THETA_SITE), state; min_abs_coeff = ATOL) + end + for i in 1:(n - 1) + state = propagate( + PauliRotation([:Z, :Z], [i, i + 1], THETA_BOND), state; min_abs_coeff = ATOL + ) + end + else + for i in 1:(n - 1) + for axes in ([:X, :X], [:Y, :Y], [:Z, :Z]) + state = propagate( + PauliRotation(axes, [i, i + 1], THETA_BOND), state; min_abs_coeff = ATOL + ) + end + end + for i in 1:n + state = propagate(PauliRotation([:Z], [i], THETA_SITE), state; min_abs_coeff = ATOL) + end + end + return state +end + +function run_model(n::Int) + state = seed_operator(n) + for _ in 1:STEPS + state = trotter_step!(state, n) + end + return state +end + +"""`⟨0…0|O|0…0⟩` for TFIM; the `Z_1` autocorrelator for Heisenberg.""" +function readout(state, n::Int) + if MODEL == "tfim" + return real(overlapwithzero(state)) + else + c = getcoeff(state, PauliString(n, [:Z], [1])) + return real(c) + end +end + +""" +Print the propagated support as `word coefficient`, largest first, site 1 +leftmost — the same format the Rust and Python runners emit under `DUMP=1` so +the driver can diff all four term-for-term. +""" +function dump_support(n::Int) + state = run_model(n) + terms = Dict{String,Float64}() + for (ps, c) in state + w = join([['I', 'X', 'Y', 'Z'][getpauli(ps, i) + 1] for i in 1:n]) + terms[w] = get(terms, w, 0.0) + real(c) + end + println("# $(length(terms)) terms") + for w in sort(collect(keys(terms)); by = w -> (-abs(terms[w]), w)) + println("$w $(Printf.@sprintf("%+.12e", terms[w]))") + end +end + +function main() + qubits = parse.(Int, split(get(ENV, "QUBITS", "8,12,16,20,24,28,32"), ",")) + if haskey(ENV, "DUMP") + dump_support(first(qubits)) + return + end + println( + stderr, + "PauliPropagation.jl $MODEL: steps=$STEPS dt=$DT J=$JCOUP h=$HFIELD atol=$ATOL " * + "iters=$ITERS threads=$(Threads.nthreads())", + ) + println("model,library,qubits,steps,dt,atol,time_s,terms,observable") + for n in qubits + # Warm up so the reported time excludes JIT for this width's type + # specialization. + st = run_model(n) + nterms = length(st) + obs = readout(st, n) + best = Inf + for _ in 1:ITERS + best = min(best, @elapsed run_model(n)) + end + println("$MODEL,pauli-propagation-jl,$n,$STEPS,$DT,$ATOL,$(round(best, sigdigits=7)),$nterms,$obs") + println(stderr, " n=$n $(round(best, digits=4))s $nterms terms obs=$obs") + flush(stdout) + end +end + +main() diff --git a/julia-benchmarks/benches/xbench_ps.jl b/julia-benchmarks/benches/xbench_ps.jl new file mode 100644 index 000000000..4a51d7380 --- /dev/null +++ b/julia-benchmarks/benches/xbench_ps.jl @@ -0,0 +1,163 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +# +# The PauliStrings.jl side of the cross-library Pauli-propagation benchmark. +# See `benchmarks/cross-library/README.md` for the shared circuit definitions, +# the parameter contract, and the CSV schema. +# +# PauliStrings.jl's own front door for this is `evolve(H, O, tspan; +# method=Trotter())`, which builds the gate list from the Hamiltonian's internal +# string order. We build the `TrotterGate` list by hand instead, so the gate +# sequence is the one the shared spec fixes and the propagated operator is +# comparable term-for-term with the other three engines. +# +# MODEL=tfim QUBITS=8,16,24 STEPS=10 DT=0.1 JCOUP=1.0 HFIELD=1.0 ATOL=1e-6 \ +# julia --project=@. -t1 benches/xbench_ps.jl + +using PauliStrings +using Printf + +const MODEL = get(ENV, "MODEL", "tfim") +const STEPS = parse(Int, get(ENV, "STEPS", "10")) +const DT = parse(Float64, get(ENV, "DT", "0.1")) +const JCOUP = parse(Float64, get(ENV, "JCOUP", "1.0")) +const HFIELD = parse(Float64, get(ENV, "HFIELD", "1.0")) +const ATOL = parse(Float64, get(ENV, "ATOL", "1e-6")) +const ITERS = parse(Int, get(ENV, "ITERS", "3")) + +const THETA_BOND = 2 * JCOUP * DT +const THETA_SITE = 2 * HFIELD * DT + +"""The `PauliString` generator for a one- or two-site Pauli term.""" +function generator(n::Int, ops::Vector{String}, sites::Vector{Int}) + o = Operator(n) + if length(sites) == 1 + o += ops[1], sites[1] + else + o += ops[1], sites[1], ops[2], sites[2] + end + return o.strings[1] +end + +""" +The Trotter gate list in *application* order, per the shared spec. + +`trotter_step!` consumes gates in matrix-multiply order and applies them in +reverse, so the caller reverses this before handing it over. +""" +function gate_list(n::Int) + gates = PauliStrings.TrotterGate[] + if MODEL == "tfim" + for i in 1:n + push!(gates, PauliStrings.TrotterGate(generator(n, ["X"], [i]), THETA_SITE)) + end + for i in 1:(n - 1) + push!( + gates, + PauliStrings.TrotterGate(generator(n, ["Z", "Z"], [i, i + 1]), THETA_BOND), + ) + end + else + for i in 1:(n - 1) + for op in ("X", "Y", "Z") + push!( + gates, + PauliStrings.TrotterGate( + generator(n, [op, op], [i, i + 1]), THETA_BOND + ), + ) + end + end + for i in 1:n + push!(gates, PauliStrings.TrotterGate(generator(n, ["Z"], [i]), THETA_SITE)) + end + end + return gates +end + +"""The seed observable: `Σ_i Z_i` for TFIM, `Z_1` for Heisenberg.""" +function seed_operator(n::Int) + o = Operator(n) + if MODEL == "tfim" + for i in 1:n + o += "Z", i + end + else + o += "Z", 1 + end + return o +end + +function run_model(n::Int, gates) + o = seed_operator(n) + truncation(x) = cutoff(x, ATOL) + for _ in 1:STEPS + PauliStrings.trotter_step!(o, gates; truncation = truncation, truncate_every = 1) + end + return o +end + +"""`⟨0…0|O|0…0⟩` for TFIM; the `Z_1` autocorrelator for Heisenberg.""" +function readout(o::Operator, n::Int) + if MODEL == "tfim" + return real(expect(o, "0"^n)) + else + z1 = seed_operator(n) + return real(trace_product(z1, o; scale = 1)) + end +end + +""" +Print the propagated support as `word coefficient`, largest first, site 1 +leftmost — the same format the other three runners emit under `DUMP=1`. + +PauliStrings.jl carries `im^{#Y}` inside the stored coefficient (its `Matrix` +convention), so `op_to_strings` is the readout that puts the coefficients on the +same footing as the other engines' real ones. +""" +function dump_support(n::Int, gates) + o = run_model(n, gates) + coeffs, strings = op_to_strings(o) + terms = Dict{String,Float64}() + for (c, s) in zip(coeffs, strings) + w = replace(s, '1' => 'I') + terms[w] = get(terms, w, 0.0) + real(c) + end + println("# $(length(terms)) terms") + for w in sort(collect(keys(terms)); by = w -> (-abs(terms[w]), w)) + @printf("%s %+.12e\n", w, terms[w]) + end +end + +function main() + qubits = parse.(Int, split(get(ENV, "QUBITS", "8,12,16,20,24,28,32"), ",")) + if haskey(ENV, "DUMP") + n = first(qubits) + dump_support(n, reverse(gate_list(n))) + return + end + println( + stderr, + "PauliStrings.jl $MODEL: steps=$STEPS dt=$DT J=$JCOUP h=$HFIELD atol=$ATOL " * + "iters=$ITERS threads=$(Threads.nthreads())", + ) + println("model,library,qubits,steps,dt,atol,time_s,terms,observable") + for n in qubits + # `gate_list` is circuit construction, not propagation — built once and + # excluded from the timed region, as in every other runner. + gates = reverse(gate_list(n)) + # Warm up so the reported time excludes JIT. + o = run_model(n, gates) + nterms = length(o) + obs = readout(o, n) + best = Inf + for _ in 1:ITERS + best = min(best, @elapsed run_model(n, gates)) + end + println("$MODEL,pauli-strings-jl,$n,$STEPS,$DT,$ATOL,$(round(best, sigdigits=7)),$nterms,$obs") + println(stderr, " n=$n $(round(best, digits=4))s $nterms terms obs=$obs") + flush(stdout) + end +end + +main() From 6e0506c9466483d02db0f5d95cd6a34965460b07 Mon Sep 17 00:00:00 2001 From: Roger-luo Date: Wed, 12 Aug 2026 00:43:16 -0400 Subject: [PATCH 2/5] docs(bench): quantify the accuracy cost of prune-at-creation truncation The README described the two truncation policies and their support difference but not what the difference costs. Adds the coefficient-vector L2 error against a converged reference, which puts monoprop within 12% of ppvm and pauli-prop at 1.7-2.0x, and records that the scalar observable is the wrong yardstick here because its truncation error changes sign. Co-authored-by: Cursor --- benchmarks/cross-library/README.md | 19 ++++++++++++++++++- 1 file changed, 18 insertions(+), 1 deletion(-) diff --git a/benchmarks/cross-library/README.md b/benchmarks/cross-library/README.md index dcea69298..2c93b59ca 100644 --- a/benchmarks/cross-library/README.md +++ b/benchmarks/cross-library/README.md @@ -153,8 +153,25 @@ the runtime. Read them before quoting a ratio. `monoprop` is within 1 % on TFIM and 4–12 % over on Heisenberg. The driver prints a `note:` line for any engine more than 2 % off, and the plot's third panel is there so the ratio is never read without it. +* **Prune-at-creation costs a little accuracy, not an order of it.** Measured as + the `L2` distance of the whole coefficient vector from a converged reference + (`ppvm` at `atol=1e-16`; the Heisenberg `n=8` sector saturates at 16 384 terms, + so that reference is exact), across `atol` from `1e-3` to `1e-7` at `n=8`, + `steps=10`: `monoprop` lands 0–12 % above `ppvm`'s error and `pauli-prop` + 1.7–2.0× above it. Error falls by a decade for every decade of `atol` in all + three, so the policies share a convergence order and differ only in the + constant. Don't read this off the `observable` column instead — it is one + scalar and its truncation error changes sign, so `ppvm` on TFIM goes + `7.7e-4, 1.8e-3, 1.0e-4, 3.0e-7, 1.2e-6` over those five thresholds. That dip + is a cancellation, and comparing against it manufactures a 27× gap that the + vector norm shows to be 6 %. * **`monoprop` tracks more rows than it reports.** It retains monomials whose - coefficient has cancelled to exactly zero. Its `size()` at Heisenberg `n=14` is + coefficient has cancelled to exactly zero. This is the visible end of + prune-at-creation: a branch that cleared the threshold when it was emitted is + never re-tested after later contributions cancel it. At Heisenberg `n=8`, + `atol=1e-4` it holds 2 074 terms whose converged value is below `atol`, against + `ppvm`'s 1 107 — so roughly 7 % of its support is dead weight. That is a cost + in work rather than in accuracy; the coefficients it keeps are still right. Its `size()` at Heisenberg `n=14` is 7 355 928 rows against the 3 204 697 terms above threshold — 2.3× — while on TFIM the two are within 1 %. The `terms` column is the above-threshold support, which is what the other four engines mean by it; `size()` goes to stderr next From 6e84129c69364f64b3b8ff23476146ff9a0108e9 Mon Sep 17 00:00:00 2001 From: Roger-luo Date: Wed, 12 Aug 2026 00:57:37 -0400 Subject: [PATCH 3/5] bench(cross-library): record the atol/accuracy sweep as a CSV The README quoted error ratios with no data behind them and no way to regenerate them. Adds xbench_accuracy.py, which diffs each engine's whole coefficient vector against a converged reference across atol and writes both the vector norms and the scalar observable error, and commits the run as accuracy.csv. Extending the sweep from three engines to all five corrects the pauli-prop range in the README, 1.23-2.60x rather than 1.7-2.0x, and shows both Julia engines reproducing ppvm's error columns exactly -- independent confirmation that the three accumulate-then-truncate engines share one truncation rule. Co-authored-by: Cursor --- benchmarks/cross-library/README.md | 23 ++- benchmarks/cross-library/accuracy.csv | 51 ++++++ benchmarks/cross-library/xbench_accuracy.py | 177 ++++++++++++++++++++ 3 files changed, 244 insertions(+), 7 deletions(-) create mode 100644 benchmarks/cross-library/accuracy.csv create mode 100644 benchmarks/cross-library/xbench_accuracy.py diff --git a/benchmarks/cross-library/README.md b/benchmarks/cross-library/README.md index 2c93b59ca..6b0f86fbf 100644 --- a/benchmarks/cross-library/README.md +++ b/benchmarks/cross-library/README.md @@ -21,6 +21,11 @@ uv run --no-project python3 benchmarks/cross-library/run_xbench.py \ uv run --no-project --with matplotlib python3 benchmarks/cross-library/plot_xbench.py \ --csv target/xbench/results.csv --out target/xbench/xbench.png \ --title "Pauli propagation: ppvm vs PauliPropagation.jl, PauliStrings.jl, pauli-prop, monoprop" + +# How wrong each engine is at a given atol, rather than how fast — see +# "Known differences" below. Regenerates accuracy.csv. +uv run --no-project python3 benchmarks/cross-library/xbench_accuracy.py \ + --out benchmarks/cross-library/accuracy.csv ``` A subset, if you prefer: `--libs ppvm,monoprop`. @@ -157,14 +162,15 @@ the runtime. Read them before quoting a ratio. the `L2` distance of the whole coefficient vector from a converged reference (`ppvm` at `atol=1e-16`; the Heisenberg `n=8` sector saturates at 16 384 terms, so that reference is exact), across `atol` from `1e-3` to `1e-7` at `n=8`, - `steps=10`: `monoprop` lands 0–12 % above `ppvm`'s error and `pauli-prop` - 1.7–2.0× above it. Error falls by a decade for every decade of `atol` in all - three, so the policies share a convergence order and differ only in the - constant. Don't read this off the `observable` column instead — it is one - scalar and its truncation error changes sign, so `ppvm` on TFIM goes + `steps=10`: `monoprop` lands within 1.00–1.12× of `ppvm`'s error and + `pauli-prop` at 1.23–2.60×. Error falls by a decade for every decade of `atol` + in all five engines, so the two policies share a convergence order and differ + only in the constant. Don't read this off the `observable` column instead — it + is one scalar and its truncation error changes sign, so `ppvm` on TFIM goes `7.7e-4, 1.8e-3, 1.0e-4, 3.0e-7, 1.2e-6` over those five thresholds. That dip - is a cancellation, and comparing against it manufactures a 27× gap that the - vector norm shows to be 6 %. + is a cancellation, and comparing against it manufactures a 27× gap where the + vector norm shows 6 %. `xbench_accuracy.py` writes both, for that reason; + [`accuracy.csv`](accuracy.csv) is the run these figures come from. * **`monoprop` tracks more rows than it reports.** It retains monomials whose coefficient has cancelled to exactly zero. This is the visible end of prune-at-creation: a branch that cleared the threshold when it was emitted is @@ -265,6 +271,9 @@ itself converges once `n` exceeds the light cone (identical from `n≈10` at - `xbench_monoprop.py` — `monoprop`, including the thread-cap assertion. - `run_xbench.py` — validation, driving, CSV merge, summary table. - `plot_xbench.py` — the figure. +- `xbench_accuracy.py` — the `atol` sweep behind `accuracy.csv`: each engine's + coefficient vector against a converged reference, so the timings can be read + next to what each truncation rule costs. [pp]: https://github.com/MSRudolph/PauliPropagation.jl [ps]: https://github.com/nicolasloizeau/PauliStrings.jl diff --git a/benchmarks/cross-library/accuracy.csv b/benchmarks/cross-library/accuracy.csv new file mode 100644 index 000000000..b06df59cd --- /dev/null +++ b/benchmarks/cross-library/accuracy.csv @@ -0,0 +1,51 @@ +model,library,qubits,steps,dt,atol,ref_atol,ref_terms,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err +tfim,ppvm,8,10,0.1,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.889668769617,1.888902654574,7.661150e-04 +tfim,pauli-propagation-jl,8,10,0.1,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 +tfim,pauli-strings-jl,8,10,0.1,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 +tfim,pauli-prop,8,10,0.1,1e-3,1e-16,27051,743,9.747655e-02,3.201066e+00,8.889589e-03,1.519868e+00,0,479,1.9013036238629,1.888902654574,1.240097e-02 +tfim,monoprop,8,10,0.1,1e-3,1e-16,27051,804,6.429517e-02,2.126366e+00,6.563651e-03,1.286381e+00,0,418,1.895370105228,1.888902654574,6.467451e-03 +tfim,ppvm,8,10,0.1,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.887148304814,1.888902654574,1.754350e-03 +tfim,pauli-propagation-jl,8,10,0.1,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 +tfim,pauli-strings-jl,8,10,0.1,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 +tfim,pauli-prop,8,10,0.1,1e-4,1e-16,27051,1527,1.655369e-02,6.909180e-01,1.396127e-03,2.516797e-01,0,761,1.89026934442961,1.888902654574,1.366690e-03 +tfim,monoprop,8,10,0.1,1e-4,1e-16,27051,1705,8.349062e-03,3.553170e-01,6.351504e-04,1.787774e-01,1,584,1.888182135464,1.888902654574,7.205191e-04 +tfim,ppvm,8,10,0.1,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.888799248154,1.888902654574,1.034064e-04 +tfim,pauli-propagation-jl,8,10,0.1,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 +tfim,pauli-strings-jl,8,10,0.1,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 +tfim,pauli-prop,8,10,0.1,1e-5,1e-16,27051,2628,2.357978e-03,1.202946e-01,1.384898e-04,3.877258e-02,0,1088,1.88885571834255,1.888902654574,4.693623e-05 +tfim,monoprop,8,10,0.1,1e-5,1e-16,27051,2998,1.023503e-03,5.342297e-02,5.628422e-05,2.245608e-02,0,718,1.888873511516,1.888902654574,2.914306e-05 +tfim,ppvm,8,10,0.1,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.888902953057,1.888902654574,2.984830e-07 +tfim,pauli-propagation-jl,8,10,0.1,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 +tfim,pauli-strings-jl,8,10,0.1,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 +tfim,pauli-prop,8,10,0.1,1e-6,1e-16,27051,4161,2.838847e-04,1.707832e-02,1.911237e-05,5.130690e-03,0,1415,1.88893188165528,1.888902654574,2.922708e-05 +tfim,monoprop,8,10,0.1,1e-6,1e-16,27051,4662,1.182467e-04,7.215810e-03,7.389538e-06,2.782631e-03,1,915,1.888910840223,1.888902654574,8.185649e-06 +tfim,ppvm,8,10,0.1,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.888901406583,1.888902654574,1.247991e-06 +tfim,pauli-propagation-jl,8,10,0.1,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 +tfim,pauli-strings-jl,8,10,0.1,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 +tfim,pauli-prop,8,10,0.1,1e-7,1e-16,27051,6056,3.227852e-05,2.215272e-03,1.725103e-06,6.248028e-04,0,1720,1.88890165053338,1.888902654574,1.004041e-06 +tfim,monoprop,8,10,0.1,1e-7,1e-16,27051,6753,1.309534e-05,9.190766e-04,7.516440e-07,3.063929e-04,1,1024,1.888902140509,1.888902654574,5.140650e-07 +heisenberg,ppvm,8,10,0.1,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.1843817263679,0.1835878942086,7.938322e-04 +heisenberg,pauli-propagation-jl,8,10,0.1,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 +heisenberg,pauli-strings-jl,8,10,0.1,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 +heisenberg,pauli-prop,8,10,0.1,1e-3,1e-16,16384,5072,2.384650e-01,2.082594e+01,1.022082e-02,9.244654e+00,1051,3169,0.186064764275212,0.1835878942086,2.476870e-03 +heisenberg,monoprop,8,10,0.1,1e-3,1e-16,16384,5106,1.937250e-01,1.714593e+01,6.293714e-03,8.315984e+00,836,2920,0.18421262641,0.1835878942086,6.247322e-04 +heisenberg,ppvm,8,10,0.1,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.1834786826111,0.1835878942086,1.092116e-04 +heisenberg,pauli-propagation-jl,8,10,0.1,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 +heisenberg,pauli-strings-jl,8,10,0.1,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 +heisenberg,pauli-prop,8,10,0.1,1e-4,1e-16,16384,13628,4.451207e-02,4.398424e+00,2.046005e-03,6.386030e-01,2836,2298,0.183655881056651,0.1835878942086,6.798685e-05 +heisenberg,monoprop,8,10,0.1,1e-4,1e-16,16384,13143,2.725847e-02,2.714465e+00,1.193874e-03,5.551494e-01,2074,2021,0.183546769437,0.1835878942086,4.112477e-05 +heisenberg,ppvm,8,10,0.1,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.1835902065899,0.1835878942086,2.312381e-06 +heisenberg,pauli-propagation-jl,8,10,0.1,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.18359020658986,0.1835878942086,2.312381e-06 +heisenberg,pauli-strings-jl,8,10,0.1,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.183590206589861,0.1835878942086,2.312381e-06 +heisenberg,pauli-prop,8,10,0.1,1e-5,1e-16,16384,17803,4.507725e-03,4.730143e-01,1.675695e-04,1.250100e-02,2677,432,0.183584912095183,0.1835878942086,2.982113e-06 +heisenberg,monoprop,8,10,0.1,1e-5,1e-16,16384,17223,2.557524e-03,2.663380e-01,1.096506e-04,9.925007e-03,2027,362,0.183585575353,0.1835878942086,2.318856e-06 +heisenberg,ppvm,8,10,0.1,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.1835880798476,0.1835878942086,1.856390e-07 +heisenberg,pauli-propagation-jl,8,10,0.1,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 +heisenberg,pauli-strings-jl,8,10,0.1,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 +heisenberg,pauli-prop,8,10,0.1,1e-6,1e-16,16384,17004,3.764296e-04,3.802784e-02,1.821500e-05,2.767212e-04,894,105,0.183587851042604,0.1835878942086,4.316600e-08 +heisenberg,monoprop,8,10,0.1,1e-6,1e-16,16384,16910,1.993482e-04,2.007446e-02,1.047708e-05,1.937806e-04,774,79,0.183587892322,0.1835878942086,1.886600e-09 +heisenberg,ppvm,8,10,0.1,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.1835878872944,0.1835878942086,6.914200e-09 +heisenberg,pauli-propagation-jl,8,10,0.1,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 +heisenberg,pauli-strings-jl,8,10,0.1,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 +heisenberg,pauli-prop,8,10,0.1,1e-7,1e-16,16384,16528,3.026389e-05,3.008580e-03,1.430835e-06,5.571167e-06,195,19,0.18358789229432,0.1835878942086,1.914280e-09 +heisenberg,monoprop,8,10,0.1,1e-7,1e-16,16384,16558,1.539943e-05,1.536256e-03,8.884821e-07,3.199623e-06,217,11,0.183587893774,0.1835878942086,4.346000e-10 diff --git a/benchmarks/cross-library/xbench_accuracy.py b/benchmarks/cross-library/xbench_accuracy.py new file mode 100644 index 000000000..3a0b5d316 --- /dev/null +++ b/benchmarks/cross-library/xbench_accuracy.py @@ -0,0 +1,177 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +"""Measure what each engine's truncation rule costs in accuracy, not in time. + +`run_xbench.py` answers "how fast", holding `atol` fixed. It cannot answer "how +wrong", because the engines do not all mean the same thing by `atol`: `ppvm`, +PauliPropagation.jl and PauliStrings.jl accumulate every contribution to a term +and then drop the merged coefficient, while `pauli-prop` and monoprop test a +branch's prospective coefficient before they emit it. Same threshold, different +rule, so the fair question is how far each lands from the untruncated answer. + +This sweeps `atol` at a small width and diffs each engine's **whole coefficient +vector** against a converged reference. The vector norm is the point. Judging +this on the scalar `observable` is actively misleading, because its truncation +error changes sign — `ppvm` on TFIM walks 7.7e-4, 1.8e-3, 1.0e-4, 3.0e-7, 1.2e-6 +over `atol` 1e-3 … 1e-7, and that accidental near-zero at 1e-6 manufactures a 27x +gap against monoprop where the norm shows 6 %. Both are recorded so the +discrepancy stays visible rather than being something you have to rediscover. + +The reference is the baseline engine at `--ref-atol`. Widths whose reachable +sector saturates give an exact reference (Heisenberg `n=8` tops out at 16 384 +terms); elsewhere check it has converged before trusting a small ratio. + + uv run --no-project python3 xbench_accuracy.py --help +""" + +from __future__ import annotations + +import argparse +import csv +import math +import sys +from pathlib import Path + +from run_xbench import BASELINE, RUNNERS, invoke, parse_dump + +CSV_COLUMNS = [ + "model", + "library", + "qubits", + "steps", + "dt", + "atol", + "ref_atol", + "ref_terms", + "terms", + "l2_err", + "l1_err", + "max_coeff_err", + "lost_mass", + "kept_subthreshold", + "dropped_above_atol", + "observable", + "ref_observable", + "obs_abs_err", +] + + +def measure(approx: dict[str, float], ref: dict[str, float], atol: float) -> dict: + """Compare one propagated support against the converged one.""" + words = set(approx) | set(ref) + diffs = [abs(approx.get(w, 0.0) - ref.get(w, 0.0)) for w in words] + return { + "terms": len(approx), + "l2_err": math.sqrt(sum(d * d for d in diffs)), + "l1_err": sum(diffs), + "max_coeff_err": max(diffs, default=0.0), + # Weight the engine discarded that the reference says was there. + "lost_mass": sum(abs(c) for w, c in ref.items() if w not in approx), + # Rows held whose converged value is under the threshold they were given: + # the signature of a branch that cleared `atol` when it was emitted and + # was never re-tested after later contributions cancelled it. + "kept_subthreshold": sum(1 for w in approx if abs(ref.get(w, 0.0)) < atol), + # The opposite failure: mass thrown away that belonged above `atol`. + "dropped_above_atol": sum( + 1 for w, c in ref.items() if w not in approx and abs(c) >= atol + ), + } + + +def observable(text: str) -> float: + """Pull the `observable` column out of a runner's one-row CSV.""" + return float( + next(csv.DictReader(l for l in text.splitlines() if l.strip()))["observable"] + ) + + +def run(lib: str, env: dict[str, str], quiet: bool) -> tuple[dict[str, float], float]: + """Get one engine's support and its readout at these parameters.""" + dump = parse_dump(invoke(RUNNERS[lib], {**env, "DUMP": "1"}, quiet)) + return dump, observable(invoke(RUNNERS[lib], env, quiet)) + + +def main() -> None: + ap = argparse.ArgumentParser(description=__doc__) + ap.add_argument("--models", default="tfim,heisenberg") + ap.add_argument("--qubits", type=int, default=8) + ap.add_argument("--steps", type=int, default=10) + ap.add_argument("--dt", default="0.1") + ap.add_argument("--atols", default="1e-3,1e-4,1e-5,1e-6,1e-7") + ap.add_argument( + "--ref-atol", + default="1e-16", + help="threshold for the reference run; must be converged at this width", + ) + ap.add_argument("--libs", default=",".join(RUNNERS)) + ap.add_argument("--out", type=Path, default=Path("target/xbench/accuracy.csv")) + ap.add_argument("--quiet", action="store_true") + args = ap.parse_args() + + libs = [lib for lib in args.libs.split(",") if lib.strip()] + if unknown := set(libs) - set(RUNNERS): + raise SystemExit(f"unknown libraries: {sorted(unknown)}") + if missing := [lib for lib in libs if not RUNNERS[lib].available()]: + print(f"skipping (toolchain absent): {', '.join(missing)}", file=sys.stderr) + libs = [lib for lib in libs if lib not in missing] + if BASELINE not in libs: + raise SystemExit(f"{BASELINE} is the reference and cannot be skipped") + + base = { + "QUBITS": str(args.qubits), + "STEPS": str(args.steps), + "DT": args.dt, + "ITERS": "1", + "JCOUP": "1.0", + "HFIELD": "1.0", + } + rows = [] + for model in args.models.split(","): + env = {**base, "MODEL": model} + ref, ref_obs = run(BASELINE, {**env, "ATOL": args.ref_atol}, args.quiet) + print( + f"{model} n={args.qubits}: reference is {BASELINE} at atol={args.ref_atol}" + f" — {len(ref)} terms, observable {ref_obs:.12g}", + file=sys.stderr, + ) + for atol in args.atols.split(","): + for lib in libs: + approx, obs = run(lib, {**env, "ATOL": atol}, args.quiet) + row = { + "model": model, + "library": lib, + "qubits": args.qubits, + "steps": args.steps, + "dt": args.dt, + "atol": atol, + "ref_atol": args.ref_atol, + "ref_terms": len(ref), + "observable": f"{obs:.15g}", + "ref_observable": f"{ref_obs:.15g}", + "obs_abs_err": f"{abs(obs - ref_obs):.6e}", + } + stats = measure(approx, ref, float(atol)) + row["terms"] = stats.pop("terms") + row.update( + { + k: f"{v:.6e}" if isinstance(v, float) else v + for k, v in stats.items() + } + ) + rows.append(row) + print( + f" atol={atol} {lib:>20}: {row['terms']:>6} terms," + f" L2 {stats['l2_err']:.3e}, observable err {row['obs_abs_err']}", + file=sys.stderr, + ) + + args.out.parent.mkdir(parents=True, exist_ok=True) + with args.out.open("w", newline="") as fh: + writer = csv.DictWriter(fh, fieldnames=CSV_COLUMNS) + writer.writeheader() + writer.writerows(rows) + print(f"wrote {len(rows)} rows to {args.out}", file=sys.stderr) + + +if __name__ == "__main__": + main() From 6cbfbf4fd4286090bf9638f9697258343420fa38 Mon Sep 17 00:00:00 2001 From: Roger-luo Date: Wed, 12 Aug 2026 10:26:41 -0400 Subject: [PATCH 4/5] bench(cross-library): add a scrambling workload to locate where truncation rules diverge The two Trotter models are a weak test bed for truncation questions: they are nearest-neighbour, uniform-angle, and conserve enough structure that their support saturates a symmetry sector rather than exploring the space. On them the prune-at-creation and accumulate-then-truncate rules stay within 1.2x of each other, which says more about the circuits than about the rules. Adds `MODEL=scramble` to the ppvm and monoprop runners: steps*n random all-to-all two-qubit Pauli rotations with random axes and random angles, built from a splitmix64 stream reimplemented identically in Rust and Python so both engines propagate the same circuit. Validated term-for-term across three seeds/widths/angles, and the existing tfim/heisenberg agreement across all five engines is unaffected. With that instance the rules do separate, and only when three conditions hold at once -- small sin(theta) to widen the rejection band, depth to populate and compound it, and a scrambled operator so the loss does not average out. At n=8, 3200 gates, angles under 0.1 and atol=1e-4, monoprop's coefficient-vector error is 2.96x ppvm's, tight at 2.90-3.06x over five seeds. Records the sweep plus a theta=pi/4 control that reverses the sign (0.85-0.96x, monoprop ahead), and documents the mechanism: merge multiplicity is exactly two, which caps the divergence and ties both errors to atol. Co-authored-by: Cursor --- benchmarks/cross-library/README.md | 91 +++++++++++++++- benchmarks/cross-library/accuracy.csv | 102 +++++++++--------- .../cross-library/accuracy_scramble.csv | 61 +++++++++++ benchmarks/cross-library/xbench_accuracy.py | 81 +++++++------- benchmarks/cross-library/xbench_monoprop.py | 62 ++++++++++- crates/ppvm-pauli-sum/examples/xbench.rs | 100 +++++++++++++++-- 6 files changed, 398 insertions(+), 99 deletions(-) create mode 100644 benchmarks/cross-library/accuracy_scramble.csv diff --git a/benchmarks/cross-library/README.md b/benchmarks/cross-library/README.md index 6b0f86fbf..1fccdcf9a 100644 --- a/benchmarks/cross-library/README.md +++ b/benchmarks/cross-library/README.md @@ -26,6 +26,12 @@ uv run --no-project --with matplotlib python3 benchmarks/cross-library/plot_xben # "Known differences" below. Regenerates accuracy.csv. uv run --no-project python3 benchmarks/cross-library/xbench_accuracy.py \ --out benchmarks/cross-library/accuracy.csv + +# The regime where the two truncation rules actually diverge: a scrambling +# circuit with small angles and great depth. Regenerates accuracy_scramble.csv. +uv run --no-project python3 benchmarks/cross-library/xbench_accuracy.py \ + --models scramble --libs ppvm,monoprop --qubits 8 --steps 400 --dt 0.05 \ + --seeds 1,2,3,4,5 --atols 1e-3,1e-4,1e-5 --out target/xbench/scr-small.csv ``` A subset, if you prefer: `--libs ppvm,monoprop`. @@ -58,6 +64,31 @@ coefficients of the X-free terms. Observable `O = Z₀`; readout the autocorrelator `S(t) = tr[Z₀·O(t)]/2ⁿ`, which is just the coefficient of `Z₀` (the Paulis are orthonormal under that pairing). +### `scramble` — random all-to-all rotations, autocorrelator + +`steps·n` two-qubit Pauli rotations `exp(−iθ/2·Pₐ⊗P_b)`, each drawing a uniformly +random **all-to-all** pair `a≠b`, random axes `Pₐ, P_b ∈ {X,Y,Z}` and a random +angle `θ ∈ (0, 2J·dt]`. Observable and readout as for `heisenberg`. + +Only `ppvm` and `monoprop` implement it (`--libs ppvm,monoprop`); it exists +because the two Trotter models are a weak test bed for *truncation* questions. +They are nearest-neighbour, uniform-angle and conserve enough structure that +their support saturates a symmetry sector — at `n=8` Heisenberg reaches exactly +16 384 of 65 536 words, and the other 49 152 are zero by symmetry rather than by +dynamics. `scramble` has no lattice, no conserved quantity and no repeated angle: +at `n=8` it fills all 65 535 non-identity words, and since unitary conjugation +preserves the Frobenius norm the whole coefficient vector has `‖c‖₂ = 1` spread +across them. That makes it the right instance for asking whether a truncation +rule can resolve one coefficient against a scrambled background. + +Both runners generate the circuit from a **splitmix64 stream reimplemented +identically in Rust and Python**, seeded by `SEED`. Neither language's stdlib RNG +is specified tightly enough to rely on, and the draw order (pair, offset, axis, +axis, angle) has to match exactly — including the `b = (a+1+r mod n−1) mod n` +trick, which avoids a rejection loop that would consume a variable number of +draws and desynchronise the two streams. The term-for-term dump diff is what +enforces this; it catches any divergence immediately. + `θ = 2·c·dt` for a Hamiltonian term `c·G` is the convention four of the five engines use directly for `exp(iθ/2·G)·P·exp(−iθ/2·G)`, so the propagated operators are identical, not merely similar. monoprop is the exception and needs @@ -70,12 +101,13 @@ they can be run directly as well as through the driver: | variable | meaning | |---|---| -| `MODEL` | `tfim` or `heisenberg` | +| `MODEL` | `tfim`, `heisenberg`, or `scramble` (`ppvm` and `monoprop` only) | | `QUBITS` | comma-separated widths | | `STEPS` | Trotter steps | | `DT`, `JCOUP`, `HFIELD` | `dt`, `J`, `h` | | `ATOL` | truncation threshold on `|c|` | | `ITERS` | timed repeats; the **minimum** is reported | +| `SEED` | `scramble` only — the circuit seed | | `DUMP` | print the propagated support instead of timing it | | `MAX_TERMS` | `pauli-prop` only — its mandatory cap (see below) | | `monoprop_NUM_THREADS`, `monoprop_PARTITIONS` | `monoprop` only — the thread cap (see below) | @@ -203,6 +235,57 @@ the runtime. Read them before quoting a ratio. * Everything is single-threaded (`julia -t1`; `pauli-prop` is single-threaded by design; `monoprop` is capped as above; `ppvm` here uses no Rayon). +### Where the two truncation rules actually diverge + +Worth writing down, because the intuitive guesses about this are wrong and the +scalar `observable` will mislead you about all of them. + +A Pauli rotation sends `Q → cos θ·Q + sin θ·(iPQ)`, so a word receives +contributions from **at most two** sources: itself, scaled by `cos`, and its +partner `PQ`, scaled by `sin`. `ppvm` thresholds that sum; `monoprop` thresholds +each contribution as it is emitted. Two consequences follow, and they set the +whole shape of the difference: + +* `monoprop` discards a child whenever `|c|·sin θ < atol` even though its parent + survives, so it loses every child of a parent in the band + `atol < |c| < atol/sin θ`. The band's width is `1/sin θ` — **set by the + rotation angle alone**, with no dependence on the shape of the coefficient + distribution. +* Because a sum of two terms can exceed its larger member by at most 2×, the + mass `ppvm` rescues is `O(atol)` per word. Both engines' errors are therefore + anchored to the same threshold and cannot be decoupled: drive `atol` low enough + that `ppvm` is accurate and `monoprop` is accurate too. + +So the gap needs three things *at once*: small `sin θ` to widen the band, depth +to populate it and compound the loss, and a scrambled instance so the loss shows +up in the answer instead of averaging out. Missing any one of them, the gap +collapses to 1.0–1.2× or reverses. Measured, at `n=8`: + +| circuit | `atol` | `monoprop` error / `ppvm` error | +|---|---|---| +| `scramble`, `dt=0.05` (`θ ≤ 0.1`), 3 200 gates | `1e-4` | **2.96×** (2.89–3.07× over 10 seeds) | +| `heisenberg`, `dt=0.05`, 320 steps | `1e-3` | 2.77× | +| `tfim`/`heisenberg`, `dt=0.1`, converged `atol` | `≤1e-5` | 1.02–1.15× | +| `scramble`, `dt=π/4` (random `θ ≤ π/2`), 80 gates | `1e-3…1e-5` | 0.85–1.00× | +| `heisenberg`, `dt=0.6` | `1e-3` | **0.67×** | + +The last two rows are the ones that catch people out. **Large angles favour +`monoprop`**, by up to 1.5×, because that is where `ppvm`'s rule has its own +failure mode: `truncate()` after each gate permanently deletes a term whose +merged coefficient has transiently cancelled below `atol`, where `monoprop` keeps +the row and lets later gates revive it. Uniform amplitudes do not help either — +`θ = π/4` is where the distribution is flattest (`sd(log|c|)` 0.66 against 2.23 +at `dt=0.1`) and it is also `monoprop`'s *best* regime, since `sin θ = 0.707` +makes the rejection band 1.41× wide, the narrowest possible. Flat amplitudes and +a wide band are mutually exclusive: you get the band from small angles, and small +angles force the geometric spread that makes amplitudes non-uniform. + +None of this is visible in the `observable` column. On the `scramble` instance +above the peak-error ratio has a median of 4.2× but a range of 0.1–576×, because +a single scalar's truncation error changes sign; the 576× is `ppvm` landing +accidentally near-exact on one seed. Neither engine loses the signal there — +relative errors are 0.4 % and 2.6 % against a peak of 0.16. + ## A run Apple M4 Pro (10 P + 4 E cores), macOS 26.6, rustc 1.96.0, Julia 1.12.6, @@ -271,9 +354,9 @@ itself converges once `n` exceeds the light cone (identical from `n≈10` at - `xbench_monoprop.py` — `monoprop`, including the thread-cap assertion. - `run_xbench.py` — validation, driving, CSV merge, summary table. - `plot_xbench.py` — the figure. -- `xbench_accuracy.py` — the `atol` sweep behind `accuracy.csv`: each engine's - coefficient vector against a converged reference, so the timings can be read - next to what each truncation rule costs. +- `xbench_accuracy.py` — the `atol` sweep behind `accuracy.csv` and + `accuracy_scramble.csv`: each engine's coefficient vector against a converged + reference, so the timings can be read next to what each truncation rule costs. [pp]: https://github.com/MSRudolph/PauliPropagation.jl [ps]: https://github.com/nicolasloizeau/PauliStrings.jl diff --git a/benchmarks/cross-library/accuracy.csv b/benchmarks/cross-library/accuracy.csv index b06df59cd..c9f73142d 100644 --- a/benchmarks/cross-library/accuracy.csv +++ b/benchmarks/cross-library/accuracy.csv @@ -1,51 +1,51 @@ -model,library,qubits,steps,dt,atol,ref_atol,ref_terms,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err -tfim,ppvm,8,10,0.1,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.889668769617,1.888902654574,7.661150e-04 -tfim,pauli-propagation-jl,8,10,0.1,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 -tfim,pauli-strings-jl,8,10,0.1,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 -tfim,pauli-prop,8,10,0.1,1e-3,1e-16,27051,743,9.747655e-02,3.201066e+00,8.889589e-03,1.519868e+00,0,479,1.9013036238629,1.888902654574,1.240097e-02 -tfim,monoprop,8,10,0.1,1e-3,1e-16,27051,804,6.429517e-02,2.126366e+00,6.563651e-03,1.286381e+00,0,418,1.895370105228,1.888902654574,6.467451e-03 -tfim,ppvm,8,10,0.1,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.887148304814,1.888902654574,1.754350e-03 -tfim,pauli-propagation-jl,8,10,0.1,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 -tfim,pauli-strings-jl,8,10,0.1,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 -tfim,pauli-prop,8,10,0.1,1e-4,1e-16,27051,1527,1.655369e-02,6.909180e-01,1.396127e-03,2.516797e-01,0,761,1.89026934442961,1.888902654574,1.366690e-03 -tfim,monoprop,8,10,0.1,1e-4,1e-16,27051,1705,8.349062e-03,3.553170e-01,6.351504e-04,1.787774e-01,1,584,1.888182135464,1.888902654574,7.205191e-04 -tfim,ppvm,8,10,0.1,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.888799248154,1.888902654574,1.034064e-04 -tfim,pauli-propagation-jl,8,10,0.1,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 -tfim,pauli-strings-jl,8,10,0.1,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 -tfim,pauli-prop,8,10,0.1,1e-5,1e-16,27051,2628,2.357978e-03,1.202946e-01,1.384898e-04,3.877258e-02,0,1088,1.88885571834255,1.888902654574,4.693623e-05 -tfim,monoprop,8,10,0.1,1e-5,1e-16,27051,2998,1.023503e-03,5.342297e-02,5.628422e-05,2.245608e-02,0,718,1.888873511516,1.888902654574,2.914306e-05 -tfim,ppvm,8,10,0.1,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.888902953057,1.888902654574,2.984830e-07 -tfim,pauli-propagation-jl,8,10,0.1,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 -tfim,pauli-strings-jl,8,10,0.1,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 -tfim,pauli-prop,8,10,0.1,1e-6,1e-16,27051,4161,2.838847e-04,1.707832e-02,1.911237e-05,5.130690e-03,0,1415,1.88893188165528,1.888902654574,2.922708e-05 -tfim,monoprop,8,10,0.1,1e-6,1e-16,27051,4662,1.182467e-04,7.215810e-03,7.389538e-06,2.782631e-03,1,915,1.888910840223,1.888902654574,8.185649e-06 -tfim,ppvm,8,10,0.1,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.888901406583,1.888902654574,1.247991e-06 -tfim,pauli-propagation-jl,8,10,0.1,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 -tfim,pauli-strings-jl,8,10,0.1,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 -tfim,pauli-prop,8,10,0.1,1e-7,1e-16,27051,6056,3.227852e-05,2.215272e-03,1.725103e-06,6.248028e-04,0,1720,1.88890165053338,1.888902654574,1.004041e-06 -tfim,monoprop,8,10,0.1,1e-7,1e-16,27051,6753,1.309534e-05,9.190766e-04,7.516440e-07,3.063929e-04,1,1024,1.888902140509,1.888902654574,5.140650e-07 -heisenberg,ppvm,8,10,0.1,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.1843817263679,0.1835878942086,7.938322e-04 -heisenberg,pauli-propagation-jl,8,10,0.1,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 -heisenberg,pauli-strings-jl,8,10,0.1,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 -heisenberg,pauli-prop,8,10,0.1,1e-3,1e-16,16384,5072,2.384650e-01,2.082594e+01,1.022082e-02,9.244654e+00,1051,3169,0.186064764275212,0.1835878942086,2.476870e-03 -heisenberg,monoprop,8,10,0.1,1e-3,1e-16,16384,5106,1.937250e-01,1.714593e+01,6.293714e-03,8.315984e+00,836,2920,0.18421262641,0.1835878942086,6.247322e-04 -heisenberg,ppvm,8,10,0.1,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.1834786826111,0.1835878942086,1.092116e-04 -heisenberg,pauli-propagation-jl,8,10,0.1,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 -heisenberg,pauli-strings-jl,8,10,0.1,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 -heisenberg,pauli-prop,8,10,0.1,1e-4,1e-16,16384,13628,4.451207e-02,4.398424e+00,2.046005e-03,6.386030e-01,2836,2298,0.183655881056651,0.1835878942086,6.798685e-05 -heisenberg,monoprop,8,10,0.1,1e-4,1e-16,16384,13143,2.725847e-02,2.714465e+00,1.193874e-03,5.551494e-01,2074,2021,0.183546769437,0.1835878942086,4.112477e-05 -heisenberg,ppvm,8,10,0.1,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.1835902065899,0.1835878942086,2.312381e-06 -heisenberg,pauli-propagation-jl,8,10,0.1,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.18359020658986,0.1835878942086,2.312381e-06 -heisenberg,pauli-strings-jl,8,10,0.1,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.183590206589861,0.1835878942086,2.312381e-06 -heisenberg,pauli-prop,8,10,0.1,1e-5,1e-16,16384,17803,4.507725e-03,4.730143e-01,1.675695e-04,1.250100e-02,2677,432,0.183584912095183,0.1835878942086,2.982113e-06 -heisenberg,monoprop,8,10,0.1,1e-5,1e-16,16384,17223,2.557524e-03,2.663380e-01,1.096506e-04,9.925007e-03,2027,362,0.183585575353,0.1835878942086,2.318856e-06 -heisenberg,ppvm,8,10,0.1,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.1835880798476,0.1835878942086,1.856390e-07 -heisenberg,pauli-propagation-jl,8,10,0.1,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 -heisenberg,pauli-strings-jl,8,10,0.1,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 -heisenberg,pauli-prop,8,10,0.1,1e-6,1e-16,16384,17004,3.764296e-04,3.802784e-02,1.821500e-05,2.767212e-04,894,105,0.183587851042604,0.1835878942086,4.316600e-08 -heisenberg,monoprop,8,10,0.1,1e-6,1e-16,16384,16910,1.993482e-04,2.007446e-02,1.047708e-05,1.937806e-04,774,79,0.183587892322,0.1835878942086,1.886600e-09 -heisenberg,ppvm,8,10,0.1,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.1835878872944,0.1835878942086,6.914200e-09 -heisenberg,pauli-propagation-jl,8,10,0.1,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 -heisenberg,pauli-strings-jl,8,10,0.1,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 -heisenberg,pauli-prop,8,10,0.1,1e-7,1e-16,16384,16528,3.026389e-05,3.008580e-03,1.430835e-06,5.571167e-06,195,19,0.18358789229432,0.1835878942086,1.914280e-09 -heisenberg,monoprop,8,10,0.1,1e-7,1e-16,16384,16558,1.539943e-05,1.536256e-03,8.884821e-07,3.199623e-06,217,11,0.183587893774,0.1835878942086,4.346000e-10 +model,library,qubits,steps,dt,seed,atol,ref_atol,ref_terms,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err +tfim,ppvm,8,10,0.1,12345,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.889668769617,1.888902654574,7.661150e-04 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 +tfim,pauli-prop,8,10,0.1,12345,1e-3,1e-16,27051,743,9.747655e-02,3.201066e+00,8.889589e-03,1.519868e+00,0,479,1.9013036238629,1.888902654574,1.240097e-02 +tfim,monoprop,8,10,0.1,12345,1e-3,1e-16,27051,804,6.429517e-02,2.126366e+00,6.563651e-03,1.286381e+00,0,418,1.895370105228,1.888902654574,6.467451e-03 +tfim,ppvm,8,10,0.1,12345,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.887148304814,1.888902654574,1.754350e-03 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 +tfim,pauli-prop,8,10,0.1,12345,1e-4,1e-16,27051,1527,1.655369e-02,6.909180e-01,1.396127e-03,2.516797e-01,0,761,1.89026934442961,1.888902654574,1.366690e-03 +tfim,monoprop,8,10,0.1,12345,1e-4,1e-16,27051,1705,8.349062e-03,3.553170e-01,6.351504e-04,1.787774e-01,1,584,1.888182135464,1.888902654574,7.205191e-04 +tfim,ppvm,8,10,0.1,12345,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.888799248154,1.888902654574,1.034064e-04 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 +tfim,pauli-prop,8,10,0.1,12345,1e-5,1e-16,27051,2628,2.357978e-03,1.202946e-01,1.384898e-04,3.877258e-02,0,1088,1.88885571834255,1.888902654574,4.693623e-05 +tfim,monoprop,8,10,0.1,12345,1e-5,1e-16,27051,2998,1.023503e-03,5.342297e-02,5.628422e-05,2.245608e-02,0,718,1.888873511516,1.888902654574,2.914306e-05 +tfim,ppvm,8,10,0.1,12345,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.888902953057,1.888902654574,2.984830e-07 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 +tfim,pauli-prop,8,10,0.1,12345,1e-6,1e-16,27051,4161,2.838847e-04,1.707832e-02,1.911237e-05,5.130690e-03,0,1415,1.88893188165528,1.888902654574,2.922708e-05 +tfim,monoprop,8,10,0.1,12345,1e-6,1e-16,27051,4662,1.182467e-04,7.215810e-03,7.389538e-06,2.782631e-03,1,915,1.888910840223,1.888902654574,8.185649e-06 +tfim,ppvm,8,10,0.1,12345,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.888901406583,1.888902654574,1.247991e-06 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 +tfim,pauli-prop,8,10,0.1,12345,1e-7,1e-16,27051,6056,3.227852e-05,2.215272e-03,1.725103e-06,6.248028e-04,0,1720,1.88890165053338,1.888902654574,1.004041e-06 +tfim,monoprop,8,10,0.1,12345,1e-7,1e-16,27051,6753,1.309534e-05,9.190766e-04,7.516440e-07,3.063929e-04,1,1024,1.888902140509,1.888902654574,5.140650e-07 +heisenberg,ppvm,8,10,0.1,12345,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.1843817263679,0.1835878942086,7.938322e-04 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 +heisenberg,pauli-prop,8,10,0.1,12345,1e-3,1e-16,16384,5072,2.384650e-01,2.082594e+01,1.022082e-02,9.244654e+00,1051,3169,0.186064764275212,0.1835878942086,2.476870e-03 +heisenberg,monoprop,8,10,0.1,12345,1e-3,1e-16,16384,5106,1.937250e-01,1.714593e+01,6.293714e-03,8.315984e+00,836,2920,0.18421262641,0.1835878942086,6.247322e-04 +heisenberg,ppvm,8,10,0.1,12345,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.1834786826111,0.1835878942086,1.092116e-04 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 +heisenberg,pauli-prop,8,10,0.1,12345,1e-4,1e-16,16384,13628,4.451207e-02,4.398424e+00,2.046005e-03,6.386030e-01,2836,2298,0.183655881056651,0.1835878942086,6.798685e-05 +heisenberg,monoprop,8,10,0.1,12345,1e-4,1e-16,16384,13143,2.725847e-02,2.714465e+00,1.193874e-03,5.551494e-01,2074,2021,0.183546769437,0.1835878942086,4.112477e-05 +heisenberg,ppvm,8,10,0.1,12345,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.1835902065899,0.1835878942086,2.312381e-06 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.18359020658986,0.1835878942086,2.312381e-06 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.183590206589861,0.1835878942086,2.312381e-06 +heisenberg,pauli-prop,8,10,0.1,12345,1e-5,1e-16,16384,17803,4.507725e-03,4.730143e-01,1.675695e-04,1.250100e-02,2677,432,0.183584912095183,0.1835878942086,2.982113e-06 +heisenberg,monoprop,8,10,0.1,12345,1e-5,1e-16,16384,17223,2.557524e-03,2.663380e-01,1.096506e-04,9.925007e-03,2027,362,0.183585575353,0.1835878942086,2.318856e-06 +heisenberg,ppvm,8,10,0.1,12345,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.1835880798476,0.1835878942086,1.856390e-07 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 +heisenberg,pauli-prop,8,10,0.1,12345,1e-6,1e-16,16384,17004,3.764296e-04,3.802784e-02,1.821500e-05,2.767212e-04,894,105,0.183587851042604,0.1835878942086,4.316600e-08 +heisenberg,monoprop,8,10,0.1,12345,1e-6,1e-16,16384,16910,1.993482e-04,2.007446e-02,1.047708e-05,1.937806e-04,774,79,0.183587892322,0.1835878942086,1.886600e-09 +heisenberg,ppvm,8,10,0.1,12345,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.1835878872944,0.1835878942086,6.914200e-09 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 +heisenberg,pauli-prop,8,10,0.1,12345,1e-7,1e-16,16384,16528,3.026389e-05,3.008580e-03,1.430835e-06,5.571167e-06,195,19,0.18358789229432,0.1835878942086,1.914280e-09 +heisenberg,monoprop,8,10,0.1,12345,1e-7,1e-16,16384,16558,1.539943e-05,1.536256e-03,8.884821e-07,3.199623e-06,217,11,0.183587893774,0.1835878942086,4.346000e-10 diff --git a/benchmarks/cross-library/accuracy_scramble.csv b/benchmarks/cross-library/accuracy_scramble.csv new file mode 100644 index 000000000..dd778a693 --- /dev/null +++ b/benchmarks/cross-library/accuracy_scramble.csv @@ -0,0 +1,61 @@ +model,library,qubits,steps,dt,seed,atol,ref_atol,ref_terms,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err +scramble,ppvm,8,400,0.05,1,1e-3,1e-16,65535,7145,9.460647e-01,1.886151e+02,2.247844e-02,1.588307e+02,1273,44679,0.1790007556962,0.1614730211106,1.752773e-02 +scramble,monoprop,8,400,0.05,1,1e-3,1e-16,65535,7923,9.619359e-01,1.897827e+02,2.840829e-02,1.554452e+02,1401,44029,0.167740351377,0.1614730211106,6.267330e-03 +scramble,ppvm,8,400,0.05,1,1e-4,1e-16,65535,64007,1.366369e-01,2.789267e+01,2.364321e-03,6.780471e-01,1332,1293,0.1614802576756,0.1614730211106,7.236565e-06 +scramble,monoprop,8,400,0.05,1,1e-4,1e-16,65535,63920,4.179348e-01,8.505128e+01,9.345567e-03,2.067460e+00,1480,1528,0.165644417085,0.1614730211106,4.171396e-03 +scramble,ppvm,8,400,0.05,1,1e-5,1e-16,65535,65388,7.801179e-03,1.593150e+00,1.427393e-04,3.874730e-03,110,119,0.1614469519445,0.1614730211106,2.606917e-05 +scramble,monoprop,8,400,0.05,1,1e-5,1e-16,65535,65376,1.940221e-02,3.961147e+00,3.854044e-04,9.243616e-03,122,143,0.161298097215,0.1614730211106,1.749239e-04 +scramble,ppvm,8,400,0.05,2,1e-3,1e-16,65535,7277,9.563637e-01,1.915102e+02,2.113910e-02,1.610994e+02,1276,44811,0.145258196215,0.141763190312,3.495006e-03 +scramble,monoprop,8,400,0.05,2,1e-3,1e-16,65535,7884,9.710104e-01,1.920922e+02,2.926046e-02,1.580419e+02,1338,44266,0.121961388302,0.141763190312,1.980180e-02 +scramble,ppvm,8,400,0.05,2,1e-4,1e-16,65535,64004,1.406716e-01,2.866367e+01,2.309150e-03,6.676462e-01,1275,1297,0.1411270619147,0.141763190312,6.361284e-04 +scramble,monoprop,8,400,0.05,2,1e-4,1e-16,65535,63950,4.073952e-01,8.292705e+01,7.614204e-03,1.971686e+00,1420,1496,0.138618759706,0.141763190312,3.144431e-03 +scramble,ppvm,8,400,0.05,2,1e-5,1e-16,65535,65400,8.667319e-03,1.771308e+00,1.364485e-04,3.740238e-03,106,96,0.1417936293827,0.141763190312,3.043907e-05 +scramble,monoprop,8,400,0.05,2,1e-5,1e-16,65535,65398,1.939419e-02,3.962640e+00,3.397129e-04,8.654144e-03,124,116,0.141898202448,0.141763190312,1.350121e-04 +scramble,ppvm,8,400,0.05,3,1e-3,1e-16,65535,7384,9.450019e-01,1.887800e+02,2.825414e-02,1.580012e+02,1296,44432,0.1524323856682,0.1460185792372,6.413806e-03 +scramble,monoprop,8,400,0.05,3,1e-3,1e-16,65535,8026,9.643916e-01,1.901604e+02,3.122063e-02,1.550255e+02,1377,43871,0.142537145171,0.1460185792372,3.481434e-03 +scramble,ppvm,8,400,0.05,3,1e-4,1e-16,65535,63904,1.404446e-01,2.871432e+01,2.399302e-03,7.247692e-01,1277,1410,0.1452846799804,0.1460185792372,7.338993e-04 +scramble,monoprop,8,400,0.05,3,1e-4,1e-16,65535,64008,4.156444e-01,8.458307e+01,8.494982e-03,1.963916e+00,1432,1461,0.140959989003,0.1460185792372,5.058590e-03 +scramble,ppvm,8,400,0.05,3,1e-5,1e-16,65535,65376,8.038709e-03,1.642135e+00,1.416798e-04,4.002692e-03,101,117,0.1460074411703,0.1460185792372,1.113807e-05 +scramble,monoprop,8,400,0.05,3,1e-5,1e-16,65535,65378,1.951992e-02,3.980961e+00,3.259493e-04,9.704331e-03,124,138,0.146049390259,0.1460185792372,3.081102e-05 +scramble,ppvm,8,400,0.05,4,1e-3,1e-16,65535,7409,9.518168e-01,1.894747e+02,2.758227e-02,1.578000e+02,1251,44408,0.1682610320682,0.1693302302423,1.069198e-03 +scramble,monoprop,8,400,0.05,4,1e-3,1e-16,65535,8156,9.682760e-01,1.903117e+02,4.164777e-02,1.539202e+02,1307,43717,0.144043912834,0.1693302302423,2.528632e-02 +scramble,ppvm,8,400,0.05,4,1e-4,1e-16,65535,63869,1.409938e-01,2.882227e+01,2.557839e-03,7.418089e-01,1298,1408,0.1686796776657,0.1693302302423,6.505526e-04 +scramble,monoprop,8,400,0.05,4,1e-4,1e-16,65535,63968,4.106727e-01,8.365478e+01,7.280307e-03,1.928183e+00,1464,1475,0.165243860705,0.1693302302423,4.086370e-03 +scramble,ppvm,8,400,0.05,4,1e-5,1e-16,65535,65376,8.336955e-03,1.701727e+00,1.487382e-04,4.134139e-03,123,126,0.1693221010763,0.1693302302423,8.129166e-06 +scramble,monoprop,8,400,0.05,4,1e-5,1e-16,65535,65368,2.026217e-02,4.130641e+00,3.527960e-04,1.150466e-02,141,152,0.169290141837,0.1693302302423,4.008841e-05 +scramble,ppvm,8,400,0.05,5,1e-3,1e-16,65535,7168,9.504922e-01,1.888251e+02,2.381869e-02,1.581427e+02,1314,44541,0.1488295179599,0.1622150514644,1.338553e-02 +scramble,monoprop,8,400,0.05,5,1e-3,1e-16,65535,8263,9.733226e-01,1.903354e+02,4.313585e-02,1.532102e+02,1398,43530,0.187629107769,0.1622150514644,2.541406e-02 +scramble,ppvm,8,400,0.05,5,1e-4,1e-16,65535,63819,1.413663e-01,2.886298e+01,2.541082e-03,7.565301e-01,1294,1459,0.1601476023761,0.1622150514644,2.067449e-03 +scramble,monoprop,8,400,0.05,5,1e-4,1e-16,65535,63888,4.233586e-01,8.599954e+01,8.111627e-03,2.104861e+00,1458,1554,0.154900547436,0.1622150514644,7.314504e-03 +scramble,ppvm,8,400,0.05,5,1e-5,1e-16,65535,65372,7.935054e-03,1.618603e+00,1.321421e-04,4.264874e-03,115,128,0.1621842282147,0.1622150514644,3.082325e-05 +scramble,monoprop,8,400,0.05,5,1e-5,1e-16,65535,65383,1.933189e-02,3.942074e+00,3.162846e-04,9.839686e-03,141,143,0.162116789286,0.1622150514644,9.826218e-05 +scramble,ppvm,8,10,0.7853981633974483,1,1e-3,1e-16,65535,50657,3.090565e-01,6.278704e+01,5.755294e-03,1.610769e+01,5920,6774,0.002947642913682,0.004567937262123,1.620294e-03 +scramble,monoprop,8,10,0.7853981633974483,1,1e-3,1e-16,65535,49970,2.637431e-01,5.349200e+01,5.279884e-03,1.500207e+01,4766,6307,0.002980535216,0.004567937262123,1.587402e-03 +scramble,ppvm,8,10,0.7853981633974483,1,1e-4,1e-16,65535,64150,1.674435e-02,3.406150e+00,3.050943e-04,9.712100e-02,374,342,0.004554433189204,0.004567937262123,1.350407e-05 +scramble,monoprop,8,10,0.7853981633974483,1,1e-4,1e-16,65535,64144,1.270353e-02,2.583770e+00,2.536806e-04,8.731900e-02,294,268,0.004512298994,0.004567937262123,5.563827e-05 +scramble,ppvm,8,10,0.7853981633974483,1,1e-5,1e-16,65535,65396,1.041360e-03,2.120041e-01,2.233626e-05,7.722638e-04,26,21,0.004568585960411,0.004567937262123,6.486983e-07 +scramble,monoprop,8,10,0.7853981633974483,1,1e-5,1e-16,65535,65392,9.593503e-04,1.913561e-01,2.030360e-05,7.639987e-04,18,17,0.004566272673,0.004567937262123,1.664589e-06 +scramble,ppvm,8,10,0.7853981633974483,2,1e-3,1e-16,65535,48647,2.791048e-01,5.628207e+01,5.037719e-03,1.722668e+01,5167,7348,-0.003044049375097,-0.002123511847918,9.205375e-04 +scramble,monoprop,8,10,0.7853981633974483,2,1e-3,1e-16,65535,49443,2.734635e-01,5.507724e+01,5.526536e-03,1.610547e+01,5454,6839,-0.002009768618,-0.002123511847918,1.137432e-04 +scramble,ppvm,8,10,0.7853981633974483,2,1e-4,1e-16,65535,64060,1.508072e-02,3.042391e+00,3.327267e-04,9.911472e-02,394,348,-0.002128546499306,-0.002123511847918,5.034651e-06 +scramble,monoprop,8,10,0.7853981633974483,2,1e-4,1e-16,65535,64063,1.363255e-02,2.758181e+00,3.033626e-04,9.342819e-02,353,304,-0.002116198144,-0.002123511847918,7.313704e-06 +scramble,ppvm,8,10,0.7853981633974483,2,1e-5,1e-16,65535,65377,8.178769e-04,1.652772e-01,1.911102e-05,8.209238e-04,27,15,-0.002124049692306,-0.002123511847918,5.378444e-07 +scramble,monoprop,8,10,0.7853981633974483,2,1e-5,1e-16,65535,65362,6.915502e-04,1.389042e-01,2.041713e-05,9.752001e-04,18,21,-0.002124065586,-0.002123511847918,5.537381e-07 +scramble,ppvm,8,10,0.7853981633974483,3,1e-3,1e-16,65534,23800,1.482153e-01,2.838560e+01,3.781926e-03,1.934654e+01,2061,4604,0.01021329788007,0.01011042655569,1.028713e-04 +scramble,monoprop,8,10,0.7853981633974483,3,1e-3,1e-16,65534,23958,1.416023e-01,2.720370e+01,3.105871e-03,1.887551e+01,1907,4292,0.010168429392,0.01011042655569,5.800284e-05 +scramble,ppvm,8,10,0.7853981633974483,3,1e-4,1e-16,65534,57545,1.677474e-02,3.310245e+00,4.291547e-04,5.726769e-01,1649,2112,0.01012727159909,0.01011042655569,1.684504e-05 +scramble,monoprop,8,10,0.7853981633974483,3,1e-4,1e-16,65534,57569,1.497634e-02,2.933245e+00,4.037880e-04,5.407710e-01,1440,1879,0.010170823571,0.01011042655569,6.039702e-05 +scramble,ppvm,8,10,0.7853981633974483,3,1e-5,1e-16,65534,64720,1.185833e-03,2.343507e-01,2.743055e-05,5.164075e-03,175,166,0.01011032002169,0.01011042655569,1.065340e-07 +scramble,monoprop,8,10,0.7853981633974483,3,1e-5,1e-16,65534,64695,9.601826e-04,1.883416e-01,2.589419e-05,5.002271e-03,133,149,0.010108359939,0.01011042655569,2.066617e-06 +scramble,ppvm,8,10,0.7853981633974483,4,1e-3,1e-16,65535,44916,2.344832e-01,4.721697e+01,4.659655e-03,1.772357e+01,4861,7238,-0.002857638598484,-0.002233897581292,6.237410e-04 +scramble,monoprop,8,10,0.7853981633974483,4,1e-3,1e-16,65535,44165,2.366092e-01,4.750515e+01,4.174483e-03,1.888975e+01,4683,7811,-0.002729812368,-0.002233897581292,4.959148e-04 +scramble,ppvm,8,10,0.7853981633974483,4,1e-4,1e-16,65535,63664,1.656759e-02,3.339823e+00,3.475736e-04,1.364974e-01,544,528,-0.002347368442479,-0.002233897581292,1.134709e-04 +scramble,monoprop,8,10,0.7853981633974483,4,1e-4,1e-16,65535,63678,1.610363e-02,3.253601e+00,3.508339e-04,1.298318e-01,494,464,-0.002305939048,-0.002233897581292,7.204147e-05 +scramble,ppvm,8,10,0.7853981633974483,4,1e-5,1e-16,65535,65344,1.076490e-03,2.173143e-01,2.137590e-05,1.167930e-03,44,36,-0.0022253215307,-0.002233897581292,8.576051e-06 +scramble,monoprop,8,10,0.7853981633974483,4,1e-5,1e-16,65535,65344,9.685830e-04,1.950917e-01,2.328722e-05,1.113917e-03,37,29,-0.002225810421,-0.002233897581292,8.087160e-06 +scramble,ppvm,8,10,0.7853981633974483,5,1e-3,1e-16,65535,46519,2.524708e-01,5.087362e+01,4.986171e-03,1.754877e+01,5039,7322,0,-0.001223130468711,1.223130e-03 +scramble,monoprop,8,10,0.7853981633974483,5,1e-3,1e-16,65535,47075,2.425934e-01,4.885867e+01,4.802603e-03,1.636502e+01,5028,6755,0,-0.001223130468711,1.223130e-03 +scramble,ppvm,8,10,0.7853981633974483,5,1e-4,1e-16,65535,63800,1.575686e-02,3.174853e+00,3.102006e-04,1.252977e-01,441,469,-0.001210275579338,-0.001223130468711,1.285489e-05 +scramble,monoprop,8,10,0.7853981633974483,5,1e-4,1e-16,65535,63817,1.319639e-02,2.629690e+00,2.919619e-04,1.111590e-01,347,358,-0.001235391553,-0.001223130468711,1.226108e-05 +scramble,ppvm,8,10,0.7853981633974483,5,1e-5,1e-16,65535,65366,8.791609e-04,1.777130e-01,2.111646e-05,8.729767e-04,28,16,-0.001221170841282,-0.001223130468711,1.959627e-06 +scramble,monoprop,8,10,0.7853981633974483,5,1e-5,1e-16,65535,65365,6.298597e-04,1.269766e-01,1.810326e-05,8.522953e-04,22,11,-0.001224291877,-0.001223130468711,1.161408e-06 diff --git a/benchmarks/cross-library/xbench_accuracy.py b/benchmarks/cross-library/xbench_accuracy.py index 3a0b5d316..0a1e027b5 100644 --- a/benchmarks/cross-library/xbench_accuracy.py +++ b/benchmarks/cross-library/xbench_accuracy.py @@ -40,6 +40,7 @@ "qubits", "steps", "dt", + "seed", "atol", "ref_atol", "ref_terms", @@ -98,6 +99,12 @@ def main() -> None: ap.add_argument("--steps", type=int, default=10) ap.add_argument("--dt", default="0.1") ap.add_argument("--atols", default="1e-3,1e-4,1e-5,1e-6,1e-7") + ap.add_argument( + "--seeds", + default="12345", + help="comma-separated circuit seeds; only `scramble` reads them, and it " + "needs several because a random instance varies", + ) ap.add_argument( "--ref-atol", default="1e-16", @@ -127,43 +134,45 @@ def main() -> None: } rows = [] for model in args.models.split(","): - env = {**base, "MODEL": model} - ref, ref_obs = run(BASELINE, {**env, "ATOL": args.ref_atol}, args.quiet) - print( - f"{model} n={args.qubits}: reference is {BASELINE} at atol={args.ref_atol}" - f" — {len(ref)} terms, observable {ref_obs:.12g}", - file=sys.stderr, - ) - for atol in args.atols.split(","): - for lib in libs: - approx, obs = run(lib, {**env, "ATOL": atol}, args.quiet) - row = { - "model": model, - "library": lib, - "qubits": args.qubits, - "steps": args.steps, - "dt": args.dt, - "atol": atol, - "ref_atol": args.ref_atol, - "ref_terms": len(ref), - "observable": f"{obs:.15g}", - "ref_observable": f"{ref_obs:.15g}", - "obs_abs_err": f"{abs(obs - ref_obs):.6e}", - } - stats = measure(approx, ref, float(atol)) - row["terms"] = stats.pop("terms") - row.update( - { - k: f"{v:.6e}" if isinstance(v, float) else v - for k, v in stats.items() + for seed in args.seeds.split(","): + env = {**base, "MODEL": model, "SEED": seed} + ref, ref_obs = run(BASELINE, {**env, "ATOL": args.ref_atol}, args.quiet) + print( + f"{model} n={args.qubits} seed={seed}: reference is {BASELINE} at" + f" atol={args.ref_atol} — {len(ref)} terms, observable {ref_obs:.12g}", + file=sys.stderr, + ) + for atol in args.atols.split(","): + for lib in libs: + approx, obs = run(lib, {**env, "ATOL": atol}, args.quiet) + row = { + "model": model, + "library": lib, + "qubits": args.qubits, + "steps": args.steps, + "dt": args.dt, + "seed": seed, + "atol": atol, + "ref_atol": args.ref_atol, + "ref_terms": len(ref), + "observable": f"{obs:.15g}", + "ref_observable": f"{ref_obs:.15g}", + "obs_abs_err": f"{abs(obs - ref_obs):.6e}", } - ) - rows.append(row) - print( - f" atol={atol} {lib:>20}: {row['terms']:>6} terms," - f" L2 {stats['l2_err']:.3e}, observable err {row['obs_abs_err']}", - file=sys.stderr, - ) + stats = measure(approx, ref, float(atol)) + row["terms"] = stats.pop("terms") + row.update( + { + k: f"{v:.6e}" if isinstance(v, float) else v + for k, v in stats.items() + } + ) + rows.append(row) + print( + f" atol={atol} {lib:>20}: {row['terms']:>6} terms," + f" L2 {stats['l2_err']:.3e}, obs err {row['obs_abs_err']}", + file=sys.stderr, + ) args.out.parent.mkdir(parents=True, exist_ok=True) with args.out.open("w", newline="") as fh: diff --git a/benchmarks/cross-library/xbench_monoprop.py b/benchmarks/cross-library/xbench_monoprop.py index 70c24c3ca..2a2c2abfc 100644 --- a/benchmarks/cross-library/xbench_monoprop.py +++ b/benchmarks/cross-library/xbench_monoprop.py @@ -62,6 +62,7 @@ HFIELD = float(os.environ.get("HFIELD", "1.0")) ATOL = float(os.environ.get("ATOL", "1e-6")) ITERS = int(os.environ.get("ITERS", "3")) +SEED = int(os.environ.get("SEED", "12345")) # A serial run has CPU/wall ~1. Anything materially above it means the thread cap # did not take and the timing is not comparable with the other engines. CPU_WALL_MAX = float(os.environ.get("CPU_WALL_MAX", "1.5")) @@ -78,7 +79,7 @@ def cpu_seconds() -> float: def seed_operator(n: int) -> PauliOperator: - """`Σ_i Z_i` for TFIM, `Z_0` for Heisenberg.""" + """`Σ_i Z_i` for TFIM, `Z_0` for Heisenberg and scramble.""" if MODEL == "tfim": terms: dict[Pauli | str, float] = {Pauli("Z", (i,)): 1.0 for i in range(n)} else: @@ -86,8 +87,65 @@ def seed_operator(n: int) -> PauliOperator: return PauliOperator(terms, n) +class SplitMix64: + """splitmix64, matching `SplitMix64` in `examples/xbench.rs` bit for bit. + + Both runners have to emit the *same* random circuit, and neither language's + stdlib RNG is specified tightly enough to rely on. This is short enough to + state in both and is checked by the term-for-term dump diff. + """ + + MASK = (1 << 64) - 1 + + def __init__(self, seed: int) -> None: + self.state = seed & self.MASK + + def next(self) -> int: + self.state = (self.state + 0x9E3779B97F4A7C15) & self.MASK + z = self.state + z = ((z ^ (z >> 30)) * 0xBF58476D1CE4E5B9) & self.MASK + z = ((z ^ (z >> 27)) * 0x94D049BB133111EB) & self.MASK + return z ^ (z >> 31) + + def unit(self) -> float: + """Uniform in `[0, 1)`, from the top 53 bits.""" + return (self.next() >> 11) / float(1 << 53) + + +AXES = "XYZ" + + +def scramble_gates(n: int) -> list[tuple[str, tuple[int, ...], float]]: + """`STEPS · n` random all-to-all two-qubit Pauli rotations. + + Draw order — pair, offset, axis, axis, angle — must match `scramble_gates` + in `examples/xbench.rs` exactly, including the offset trick that avoids a + rejection loop for `b != a` (a loop would consume a variable number of + draws and desynchronise the two streams). + """ + if n < 2: + raise SystemExit("scramble needs at least 2 qubits") + rng = SplitMix64(SEED) + theta_max = 2 * JCOUP * DT + gates = [] + for _ in range(STEPS * n): + a = rng.next() % n + b = (a + 1 + rng.next() % (n - 1)) % n + pa = AXES[rng.next() % 3] + pb = AXES[rng.next() % 3] + theta = theta_max * rng.unit() + # P_a ⊗ P_b is symmetric under swapping both, so normalise to ascending + # qubit order -- monoprop's `Pauli` expects sorted supports. + if a > b: + a, b, pa, pb = b, a, pb, pa + gates.append((pa + pb, (a, b), theta)) + return gates + + def gate_list(n: int) -> list[tuple[str, tuple[int, ...], float]]: - """`STEPS` first-order Trotter steps in *application* order, per the spec.""" + """The gate sequence in *application* order, per the spec.""" + if MODEL == "scramble": + return scramble_gates(n) gates: list[tuple[str, tuple[int, ...], float]] = [] for _ in range(STEPS): if MODEL == "tfim": diff --git a/crates/ppvm-pauli-sum/examples/xbench.rs b/crates/ppvm-pauli-sum/examples/xbench.rs index 0f3c76d04..3b048487b 100644 --- a/crates/ppvm-pauli-sum/examples/xbench.rs +++ b/crates/ppvm-pauli-sum/examples/xbench.rs @@ -28,12 +28,14 @@ struct Params { h: f64, atol: f64, iters: usize, + seed: u64, } #[derive(Clone, Copy, PartialEq, Eq)] enum Model { Tfim, Heisenberg, + Scramble, } impl Model { @@ -41,7 +43,10 @@ impl Model { match s { "tfim" => Model::Tfim, "heisenberg" => Model::Heisenberg, - other => panic!("unknown MODEL {other:?} (expected `tfim` or `heisenberg`)"), + "scramble" => Model::Scramble, + other => { + panic!("unknown MODEL {other:?} (expected `tfim`, `heisenberg` or `scramble`)") + } } } @@ -49,10 +54,74 @@ impl Model { match self { Model::Tfim => "tfim", Model::Heisenberg => "heisenberg", + Model::Scramble => "scramble", } } } +/// splitmix64. Reimplemented rather than pulled in as a dependency because the +/// monoprop runner has to emit a bit-identical gate sequence from Python, and +/// this is short enough to state twice and check against a term-for-term diff. +struct SplitMix64(u64); + +impl SplitMix64 { + fn next(&mut self) -> u64 { + self.0 = self.0.wrapping_add(0x9E37_79B9_7F4A_7C15); + let mut z = self.0; + z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9); + z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB); + z ^ (z >> 31) + } + + /// Uniform in `[0, 1)`, from the top 53 bits. + fn unit(&mut self) -> f64 { + (self.next() >> 11) as f64 / (1u64 << 53) as f64 + } +} + +/// One two-qubit Pauli rotation `exp(-i θ/2 · P_a ⊗ P_b)`. +#[derive(Clone, Copy)] +struct Gate { + axis_a: [u8; 2], + axis_b: [u8; 2], + a: usize, + b: usize, + theta: f64, +} + +/// `[x, z]` bits for the three non-identity Paulis, indexed `0..3`. +const AXES: [[u8; 2]; 3] = [[1, 0], [1, 1], [0, 1]]; + +/// The `scramble` workload: `steps · n` two-qubit Pauli rotations on uniformly +/// random *all-to-all* pairs, with random axes and random angles in `(0, 2·J·dt]`. +/// +/// Unlike the two Trotter models this has no lattice, no conserved quantity and +/// no uniform angle, so the propagated operator spreads over the whole `4^n` +/// space and the coefficient distribution is genuinely scrambled rather than +/// hierarchically ordered by Pauli weight. +fn scramble_gates(n: usize, p: Params) -> Vec { + assert!(n >= 2, "scramble needs at least 2 qubits"); + let mut rng = SplitMix64(p.seed); + let theta_max = 2.0 * p.j * p.dt; + (0..p.steps * n) + .map(|_| { + let a = (rng.next() % n as u64) as usize; + // Offset by 1..n-1 so `b != a` without a rejection loop, which would + // desynchronise the two implementations' draw counts. + let b = (a + 1 + (rng.next() % (n as u64 - 1)) as usize) % n; + let axis_a = AXES[(rng.next() % 3) as usize]; + let axis_b = AXES[(rng.next() % 3) as usize]; + Gate { + axis_a, + axis_b, + a, + b, + theta: theta_max * rng.unit(), + } + }) + .collect() +} + fn env_f64(key: &str, default: f64) -> f64 { std::env::var(key) .ok() @@ -77,6 +146,7 @@ impl Params { h: env_f64("HFIELD", 1.0), atol: env_f64("ATOL", 1e-6), iters: env_usize("ITERS", 3), + seed: env_usize("SEED", 12345) as u64, } } } @@ -111,7 +181,7 @@ fn seed(n: usize, p: Params) -> PauliSum> { sum += (PauliWord::from(site_word(n, i, 'Z').as_str()), 1.0); } } - Model::Heisenberg => { + Model::Heisenberg | Model::Scramble => { sum += (PauliWord::from(site_word(n, 0, 'Z').as_str()), 1.0); } } @@ -119,9 +189,16 @@ fn seed(n: usize, p: Params) -> PauliSum> { } /// Propagate the shared gate sequence through `state`. -fn propagate(state: &mut PauliSum>, n: usize, p: Params) { +fn propagate(state: &mut PauliSum>, n: usize, p: Params, gates: &[Gate]) { let theta_bond = 2.0 * p.j * p.dt; let theta_site = 2.0 * p.h * p.dt; + if p.model == Model::Scramble { + for g in gates { + state.rotate_2(g.axis_a, g.axis_b, g.a, g.b, g.theta); + state.truncate(); + } + return; + } for _ in 0..p.steps { match p.model { Model::Tfim => { @@ -148,6 +225,7 @@ fn propagate(state: &mut PauliSum>, n: usize, p: Params) state.truncate(); } } + Model::Scramble => unreachable!("handled above"), } } } @@ -164,7 +242,7 @@ fn readout(state: &PauliSum>, n: usize, p: Params) -> f64 .filter(|(word, _)| (0..n).all(|i| !word.get_xbit(i))) .map(|(_, c)| *c) .sum(), - Model::Heisenberg => state + Model::Heisenberg | Model::Scramble => state .data() .get(&PauliWord::from(site_word(n, 0, 'Z').as_str())) .copied() @@ -175,13 +253,18 @@ fn readout(state: &PauliSum>, n: usize, p: Params) -> f64 /// Run one model at width `n`, returning `(best seconds, final support, observable)`. fn run(n: usize, p: Params) -> (f64, usize, f64) { let base = seed::(n, p); + let gates = if p.model == Model::Scramble { + scramble_gates(n, p) + } else { + Vec::new() + }; let mut best = f64::INFINITY; let mut terms = 0usize; let mut observable = f64::NAN; for _ in 0..p.iters { let mut state = base.clone(); let t0 = Instant::now(); - propagate(&mut state, n, p); + propagate(&mut state, n, p, &gates); best = best.min(t0.elapsed().as_secs_f64()); terms = state.len(); observable = readout(&state, n, p); @@ -193,7 +276,12 @@ fn run(n: usize, p: Params) -> (f64, usize, f64) { /// the format the driver diffs across engines. fn dump(n: usize, p: Params) { let mut state = seed::(n, p); - propagate(&mut state, n, p); + let gates = if p.model == Model::Scramble { + scramble_gates(n, p) + } else { + Vec::new() + }; + propagate(&mut state, n, p, &gates); let mut out: Vec<(String, f64)> = state.iter().map(|(k, c)| (k.to_string(), *c)).collect(); out.sort_by(|a, b| { b.1.abs() From 58d5100a7828e0ba85e57365df39df63fac1b247 Mon Sep 17 00:00:00 2001 From: Roger-luo Date: Wed, 12 Aug 2026 10:38:50 -0400 Subject: [PATCH 5/5] bench(cross-library): plot the truncation-error scaling and record its sweeps The error-scaling story was spread across ad-hoc runs whose output only ever reached a terminal, so none of the ratios quoted in the README could be checked or replotted. Generalises xbench_accuracy.py to sweep models, depth, angle, seeds and atol as comma-separated lists over one flattened product, and to record ||c*|| so the absolute errors can be read as relative ones. Adds plot_accuracy.py: four panels, one per question the data was asked. Consolidates the divergence sweeps (depth, angle, scrambler, and the theta=pi/4 control) into accuracy_divergence.csv, replacing accuracy_scramble.csv. The panels are chosen so the two easy misreadings are hard to make: panel 1 shows every engine falling one decade per decade of atol, so the rules share a convergence order, and panel 3 shows the ratio crossing 1.0 with angle, so no single ratio characterises them. Co-authored-by: Cursor --- benchmarks/cross-library/README.md | 34 +- benchmarks/cross-library/accuracy.csv | 102 +++--- .../cross-library/accuracy_divergence.csv | 99 ++++++ .../cross-library/accuracy_scramble.csv | 61 ---- benchmarks/cross-library/plot_accuracy.py | 315 ++++++++++++++++++ benchmarks/cross-library/xbench_accuracy.py | 87 ++--- 6 files changed, 540 insertions(+), 158 deletions(-) create mode 100644 benchmarks/cross-library/accuracy_divergence.csv delete mode 100644 benchmarks/cross-library/accuracy_scramble.csv create mode 100644 benchmarks/cross-library/plot_accuracy.py diff --git a/benchmarks/cross-library/README.md b/benchmarks/cross-library/README.md index 1fccdcf9a..4c02f6232 100644 --- a/benchmarks/cross-library/README.md +++ b/benchmarks/cross-library/README.md @@ -27,11 +27,22 @@ uv run --no-project --with matplotlib python3 benchmarks/cross-library/plot_xben uv run --no-project python3 benchmarks/cross-library/xbench_accuracy.py \ --out benchmarks/cross-library/accuracy.csv -# The regime where the two truncation rules actually diverge: a scrambling -# circuit with small angles and great depth. Regenerates accuracy_scramble.csv. -uv run --no-project python3 benchmarks/cross-library/xbench_accuracy.py \ - --models scramble --libs ppvm,monoprop --qubits 8 --steps 400 --dt 0.05 \ - --seeds 1,2,3,4,5 --atols 1e-3,1e-4,1e-5 --out target/xbench/scr-small.csv +# The regime where the two truncation rules actually diverge. Three sweeps make +# up accuracy_divergence.csv — depth, angle, and a scrambling circuit with small +# angles and great depth (plus a θ=π/4 control on the same family). +A=benchmarks/cross-library/xbench_accuracy.py +uv run --no-project python3 $A --models heisenberg --libs ppvm,monoprop \ + --qubits 6 --dt 0.05 --steps 10,20,40,80,160,320 --atols 1e-3,1e-5 --out d.csv +uv run --no-project python3 $A --models heisenberg --libs ppvm,monoprop \ + --qubits 8 --steps 10 --dt 0.02,0.05,0.1,0.2,0.4,0.6,0.785 --atols 1e-3 --out a.csv +uv run --no-project python3 $A --models scramble --libs ppvm,monoprop --qubits 8 \ + --steps 400 --dt 0.05 --seeds 1,2,3,4,5 --atols 1e-3,1e-4,1e-5 --out s.csv + +# Four panels: atol scaling, the depth trend, the angle reversal, the scrambler. +uv run --no-project --with matplotlib python3 benchmarks/cross-library/plot_accuracy.py \ + --csv benchmarks/cross-library/accuracy.csv \ + --divergence benchmarks/cross-library/accuracy_divergence.csv \ + --out target/xbench/accuracy.png ``` A subset, if you prefer: `--libs ppvm,monoprop`. @@ -269,6 +280,10 @@ collapses to 1.0–1.2× or reverses. Measured, at `n=8`: | `scramble`, `dt=π/4` (random `θ ≤ π/2`), 80 gates | `1e-3…1e-5` | 0.85–1.00× | | `heisenberg`, `dt=0.6` | `1e-3` | **0.67×** | +[`accuracy_divergence.csv`](accuracy_divergence.csv) holds these sweeps, and +`plot_accuracy.py` renders them as four panels — the `atol` scaling, the depth +trend, the angle reversal, and the scrambler beside its control. + The last two rows are the ones that catch people out. **Large angles favour `monoprop`**, by up to 1.5×, because that is where `ppvm`'s rule has its own failure mode: `truncate()` after each gate permanently deletes a term whose @@ -354,9 +369,12 @@ itself converges once `n` exceeds the light cone (identical from `n≈10` at - `xbench_monoprop.py` — `monoprop`, including the thread-cap assertion. - `run_xbench.py` — validation, driving, CSV merge, summary table. - `plot_xbench.py` — the figure. -- `xbench_accuracy.py` — the `atol` sweep behind `accuracy.csv` and - `accuracy_scramble.csv`: each engine's coefficient vector against a converged - reference, so the timings can be read next to what each truncation rule costs. +- `xbench_accuracy.py` — the sweeps behind `accuracy.csv` (all five engines vs + `atol`) and `accuracy_divergence.csv` (depth, angle, and the scrambler): each + engine's coefficient vector against a converged reference, so the timings can + be read next to what each truncation rule costs. Sweeps any of `--models`, + `--steps`, `--dt`, `--seeds`, `--atols` as comma-separated lists. +- `plot_accuracy.py` — the four-panel error-scaling figure from those two CSVs. [pp]: https://github.com/MSRudolph/PauliPropagation.jl [ps]: https://github.com/nicolasloizeau/PauliStrings.jl diff --git a/benchmarks/cross-library/accuracy.csv b/benchmarks/cross-library/accuracy.csv index c9f73142d..4b2f3cbbd 100644 --- a/benchmarks/cross-library/accuracy.csv +++ b/benchmarks/cross-library/accuracy.csv @@ -1,51 +1,51 @@ -model,library,qubits,steps,dt,seed,atol,ref_atol,ref_terms,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err -tfim,ppvm,8,10,0.1,12345,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.889668769617,1.888902654574,7.661150e-04 -tfim,pauli-propagation-jl,8,10,0.1,12345,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 -tfim,pauli-strings-jl,8,10,0.1,12345,1e-3,1e-16,27051,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 -tfim,pauli-prop,8,10,0.1,12345,1e-3,1e-16,27051,743,9.747655e-02,3.201066e+00,8.889589e-03,1.519868e+00,0,479,1.9013036238629,1.888902654574,1.240097e-02 -tfim,monoprop,8,10,0.1,12345,1e-3,1e-16,27051,804,6.429517e-02,2.126366e+00,6.563651e-03,1.286381e+00,0,418,1.895370105228,1.888902654574,6.467451e-03 -tfim,ppvm,8,10,0.1,12345,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.887148304814,1.888902654574,1.754350e-03 -tfim,pauli-propagation-jl,8,10,0.1,12345,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 -tfim,pauli-strings-jl,8,10,0.1,12345,1e-4,1e-16,27051,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 -tfim,pauli-prop,8,10,0.1,12345,1e-4,1e-16,27051,1527,1.655369e-02,6.909180e-01,1.396127e-03,2.516797e-01,0,761,1.89026934442961,1.888902654574,1.366690e-03 -tfim,monoprop,8,10,0.1,12345,1e-4,1e-16,27051,1705,8.349062e-03,3.553170e-01,6.351504e-04,1.787774e-01,1,584,1.888182135464,1.888902654574,7.205191e-04 -tfim,ppvm,8,10,0.1,12345,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.888799248154,1.888902654574,1.034064e-04 -tfim,pauli-propagation-jl,8,10,0.1,12345,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 -tfim,pauli-strings-jl,8,10,0.1,12345,1e-5,1e-16,27051,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 -tfim,pauli-prop,8,10,0.1,12345,1e-5,1e-16,27051,2628,2.357978e-03,1.202946e-01,1.384898e-04,3.877258e-02,0,1088,1.88885571834255,1.888902654574,4.693623e-05 -tfim,monoprop,8,10,0.1,12345,1e-5,1e-16,27051,2998,1.023503e-03,5.342297e-02,5.628422e-05,2.245608e-02,0,718,1.888873511516,1.888902654574,2.914306e-05 -tfim,ppvm,8,10,0.1,12345,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.888902953057,1.888902654574,2.984830e-07 -tfim,pauli-propagation-jl,8,10,0.1,12345,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 -tfim,pauli-strings-jl,8,10,0.1,12345,1e-6,1e-16,27051,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 -tfim,pauli-prop,8,10,0.1,12345,1e-6,1e-16,27051,4161,2.838847e-04,1.707832e-02,1.911237e-05,5.130690e-03,0,1415,1.88893188165528,1.888902654574,2.922708e-05 -tfim,monoprop,8,10,0.1,12345,1e-6,1e-16,27051,4662,1.182467e-04,7.215810e-03,7.389538e-06,2.782631e-03,1,915,1.888910840223,1.888902654574,8.185649e-06 -tfim,ppvm,8,10,0.1,12345,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.888901406583,1.888902654574,1.247991e-06 -tfim,pauli-propagation-jl,8,10,0.1,12345,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 -tfim,pauli-strings-jl,8,10,0.1,12345,1e-7,1e-16,27051,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 -tfim,pauli-prop,8,10,0.1,12345,1e-7,1e-16,27051,6056,3.227852e-05,2.215272e-03,1.725103e-06,6.248028e-04,0,1720,1.88890165053338,1.888902654574,1.004041e-06 -tfim,monoprop,8,10,0.1,12345,1e-7,1e-16,27051,6753,1.309534e-05,9.190766e-04,7.516440e-07,3.063929e-04,1,1024,1.888902140509,1.888902654574,5.140650e-07 -heisenberg,ppvm,8,10,0.1,12345,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.1843817263679,0.1835878942086,7.938322e-04 -heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 -heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-3,1e-16,16384,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 -heisenberg,pauli-prop,8,10,0.1,12345,1e-3,1e-16,16384,5072,2.384650e-01,2.082594e+01,1.022082e-02,9.244654e+00,1051,3169,0.186064764275212,0.1835878942086,2.476870e-03 -heisenberg,monoprop,8,10,0.1,12345,1e-3,1e-16,16384,5106,1.937250e-01,1.714593e+01,6.293714e-03,8.315984e+00,836,2920,0.18421262641,0.1835878942086,6.247322e-04 -heisenberg,ppvm,8,10,0.1,12345,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.1834786826111,0.1835878942086,1.092116e-04 -heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 -heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-4,1e-16,16384,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 -heisenberg,pauli-prop,8,10,0.1,12345,1e-4,1e-16,16384,13628,4.451207e-02,4.398424e+00,2.046005e-03,6.386030e-01,2836,2298,0.183655881056651,0.1835878942086,6.798685e-05 -heisenberg,monoprop,8,10,0.1,12345,1e-4,1e-16,16384,13143,2.725847e-02,2.714465e+00,1.193874e-03,5.551494e-01,2074,2021,0.183546769437,0.1835878942086,4.112477e-05 -heisenberg,ppvm,8,10,0.1,12345,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.1835902065899,0.1835878942086,2.312381e-06 -heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.18359020658986,0.1835878942086,2.312381e-06 -heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-5,1e-16,16384,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.183590206589861,0.1835878942086,2.312381e-06 -heisenberg,pauli-prop,8,10,0.1,12345,1e-5,1e-16,16384,17803,4.507725e-03,4.730143e-01,1.675695e-04,1.250100e-02,2677,432,0.183584912095183,0.1835878942086,2.982113e-06 -heisenberg,monoprop,8,10,0.1,12345,1e-5,1e-16,16384,17223,2.557524e-03,2.663380e-01,1.096506e-04,9.925007e-03,2027,362,0.183585575353,0.1835878942086,2.318856e-06 -heisenberg,ppvm,8,10,0.1,12345,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.1835880798476,0.1835878942086,1.856390e-07 -heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 -heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-6,1e-16,16384,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 -heisenberg,pauli-prop,8,10,0.1,12345,1e-6,1e-16,16384,17004,3.764296e-04,3.802784e-02,1.821500e-05,2.767212e-04,894,105,0.183587851042604,0.1835878942086,4.316600e-08 -heisenberg,monoprop,8,10,0.1,12345,1e-6,1e-16,16384,16910,1.993482e-04,2.007446e-02,1.047708e-05,1.937806e-04,774,79,0.183587892322,0.1835878942086,1.886600e-09 -heisenberg,ppvm,8,10,0.1,12345,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.1835878872944,0.1835878942086,6.914200e-09 -heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 -heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-7,1e-16,16384,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 -heisenberg,pauli-prop,8,10,0.1,12345,1e-7,1e-16,16384,16528,3.026389e-05,3.008580e-03,1.430835e-06,5.571167e-06,195,19,0.18358789229432,0.1835878942086,1.914280e-09 -heisenberg,monoprop,8,10,0.1,12345,1e-7,1e-16,16384,16558,1.539943e-05,1.536256e-03,8.884821e-07,3.199623e-06,217,11,0.183587893774,0.1835878942086,4.346000e-10 +model,library,qubits,steps,dt,seed,atol,ref_atol,ref_terms,ref_norm,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err +tfim,ppvm,8,10,0.1,12345,1e-3,1e-16,27051,2.828427e+00,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.889668769617,1.888902654574,7.661150e-04 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-3,1e-16,27051,2.828427e+00,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-3,1e-16,27051,2.828427e+00,834,5.730178e-02,1.956886e+00,6.556158e-03,1.180394e+00,0,388,1.88966876961681,1.888902654574,7.661150e-04 +tfim,pauli-prop,8,10,0.1,12345,1e-3,1e-16,27051,2.828427e+00,743,9.747655e-02,3.201066e+00,8.889589e-03,1.519868e+00,0,479,1.9013036238629,1.888902654574,1.240097e-02 +tfim,monoprop,8,10,0.1,12345,1e-3,1e-16,27051,2.828427e+00,804,6.429517e-02,2.126366e+00,6.563651e-03,1.286381e+00,0,418,1.895370105228,1.888902654574,6.467451e-03 +tfim,ppvm,8,10,0.1,12345,1e-4,1e-16,27051,2.828427e+00,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.887148304814,1.888902654574,1.754350e-03 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-4,1e-16,27051,2.828427e+00,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-4,1e-16,27051,2.828427e+00,1730,7.713204e-03,3.399904e-01,5.538029e-04,1.657662e-01,0,558,1.88714830481379,1.888902654574,1.754350e-03 +tfim,pauli-prop,8,10,0.1,12345,1e-4,1e-16,27051,2.828427e+00,1527,1.655369e-02,6.909180e-01,1.396127e-03,2.516797e-01,0,761,1.89026934442961,1.888902654574,1.366690e-03 +tfim,monoprop,8,10,0.1,12345,1e-4,1e-16,27051,2.828427e+00,1705,8.349062e-03,3.553170e-01,6.351504e-04,1.787774e-01,1,584,1.888182135464,1.888902654574,7.205191e-04 +tfim,ppvm,8,10,0.1,12345,1e-5,1e-16,27051,2.828427e+00,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.888799248154,1.888902654574,1.034064e-04 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-5,1e-16,27051,2.828427e+00,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-5,1e-16,27051,2.828427e+00,3001,9.762681e-04,5.198656e-02,5.615238e-05,2.218986e-02,0,715,1.88879924815445,1.888902654574,1.034064e-04 +tfim,pauli-prop,8,10,0.1,12345,1e-5,1e-16,27051,2.828427e+00,2628,2.357978e-03,1.202946e-01,1.384898e-04,3.877258e-02,0,1088,1.88885571834255,1.888902654574,4.693623e-05 +tfim,monoprop,8,10,0.1,12345,1e-5,1e-16,27051,2.828427e+00,2998,1.023503e-03,5.342297e-02,5.628422e-05,2.245608e-02,0,718,1.888873511516,1.888902654574,2.914306e-05 +tfim,ppvm,8,10,0.1,12345,1e-6,1e-16,27051,2.828427e+00,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.888902953057,1.888902654574,2.984830e-07 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-6,1e-16,27051,2.828427e+00,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-6,1e-16,27051,2.828427e+00,4701,1.118273e-04,6.914840e-03,7.337651e-06,2.638321e-03,0,875,1.88890295305684,1.888902654574,2.984828e-07 +tfim,pauli-prop,8,10,0.1,12345,1e-6,1e-16,27051,2.828427e+00,4161,2.838847e-04,1.707832e-02,1.911237e-05,5.130690e-03,0,1415,1.88893188165528,1.888902654574,2.922708e-05 +tfim,monoprop,8,10,0.1,12345,1e-6,1e-16,27051,2.828427e+00,4662,1.182467e-04,7.215810e-03,7.389538e-06,2.782631e-03,1,915,1.888910840223,1.888902654574,8.185649e-06 +tfim,ppvm,8,10,0.1,12345,1e-7,1e-16,27051,2.828427e+00,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.888901406583,1.888902654574,1.247991e-06 +tfim,pauli-propagation-jl,8,10,0.1,12345,1e-7,1e-16,27051,2.828427e+00,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 +tfim,pauli-strings-jl,8,10,0.1,12345,1e-7,1e-16,27051,2.828427e+00,6780,1.239548e-05,8.753233e-04,7.429950e-07,2.963253e-04,0,996,1.88890140658276,1.888902654574,1.247991e-06 +tfim,pauli-prop,8,10,0.1,12345,1e-7,1e-16,27051,2.828427e+00,6056,3.227852e-05,2.215272e-03,1.725103e-06,6.248028e-04,0,1720,1.88890165053338,1.888902654574,1.004041e-06 +tfim,monoprop,8,10,0.1,12345,1e-7,1e-16,27051,2.828427e+00,6753,1.309534e-05,9.190766e-04,7.516440e-07,3.063929e-04,1,1024,1.888902140509,1.888902654574,5.140650e-07 +heisenberg,ppvm,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.1843817263679,0.1835878942086,7.938322e-04 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.184381726367865,0.1835878942086,7.938322e-04 +heisenberg,pauli-prop,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,5072,2.384650e-01,2.082594e+01,1.022082e-02,9.244654e+00,1051,3169,0.186064764275212,0.1835878942086,2.476870e-03 +heisenberg,monoprop,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,5106,1.937250e-01,1.714593e+01,6.293714e-03,8.315984e+00,836,2920,0.18421262641,0.1835878942086,6.247322e-04 +heisenberg,ppvm,8,10,0.1,12345,1e-4,1e-16,16384,1.000000e+00,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.1834786826111,0.1835878942086,1.092116e-04 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-4,1e-16,16384,1.000000e+00,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-4,1e-16,16384,1.000000e+00,11923,2.564185e-02,2.511013e+00,1.063198e-03,5.963471e-01,1107,2274,0.183478682611051,0.1835878942086,1.092116e-04 +heisenberg,pauli-prop,8,10,0.1,12345,1e-4,1e-16,16384,1.000000e+00,13628,4.451207e-02,4.398424e+00,2.046005e-03,6.386030e-01,2836,2298,0.183655881056651,0.1835878942086,6.798685e-05 +heisenberg,monoprop,8,10,0.1,12345,1e-4,1e-16,16384,1.000000e+00,13143,2.725847e-02,2.714465e+00,1.193874e-03,5.551494e-01,2074,2021,0.183546769437,0.1835878942086,4.112477e-05 +heisenberg,ppvm,8,10,0.1,12345,1e-5,1e-16,16384,1.000000e+00,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.1835902065899,0.1835878942086,2.312381e-06 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-5,1e-16,16384,1.000000e+00,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.18359020658986,0.1835878942086,2.312381e-06 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-5,1e-16,16384,1.000000e+00,16101,2.426734e-03,2.449232e-01,1.074420e-04,1.125860e-02,963,420,0.183590206589861,0.1835878942086,2.312381e-06 +heisenberg,pauli-prop,8,10,0.1,12345,1e-5,1e-16,16384,1.000000e+00,17803,4.507725e-03,4.730143e-01,1.675695e-04,1.250100e-02,2677,432,0.183584912095183,0.1835878942086,2.982113e-06 +heisenberg,monoprop,8,10,0.1,12345,1e-5,1e-16,16384,1.000000e+00,17223,2.557524e-03,2.663380e-01,1.096506e-04,9.925007e-03,2027,362,0.183585575353,0.1835878942086,2.318856e-06 +heisenberg,ppvm,8,10,0.1,12345,1e-6,1e-16,16384,1.000000e+00,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.1835880798476,0.1835878942086,1.856390e-07 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-6,1e-16,16384,1.000000e+00,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-6,1e-16,16384,1.000000e+00,16324,1.873617e-04,1.844112e-02,8.671387e-06,1.957258e-04,191,82,0.183588079847589,0.1835878942086,1.856390e-07 +heisenberg,pauli-prop,8,10,0.1,12345,1e-6,1e-16,16384,1.000000e+00,17004,3.764296e-04,3.802784e-02,1.821500e-05,2.767212e-04,894,105,0.183587851042604,0.1835878942086,4.316600e-08 +heisenberg,monoprop,8,10,0.1,12345,1e-6,1e-16,16384,1.000000e+00,16910,1.993482e-04,2.007446e-02,1.047708e-05,1.937806e-04,774,79,0.183587892322,0.1835878942086,1.886600e-09 +heisenberg,ppvm,8,10,0.1,12345,1e-7,1e-16,16384,1.000000e+00,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.1835878872944,0.1835878942086,6.914200e-09 +heisenberg,pauli-propagation-jl,8,10,0.1,12345,1e-7,1e-16,16384,1.000000e+00,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 +heisenberg,pauli-strings-jl,8,10,0.1,12345,1e-7,1e-16,16384,1.000000e+00,16377,1.541357e-05,1.514858e-03,7.285205e-07,3.268561e-06,39,14,0.183587887294405,0.1835878942086,6.914195e-09 +heisenberg,pauli-prop,8,10,0.1,12345,1e-7,1e-16,16384,1.000000e+00,16528,3.026389e-05,3.008580e-03,1.430835e-06,5.571167e-06,195,19,0.18358789229432,0.1835878942086,1.914280e-09 +heisenberg,monoprop,8,10,0.1,12345,1e-7,1e-16,16384,1.000000e+00,16558,1.539943e-05,1.536256e-03,8.884821e-07,3.199623e-06,217,11,0.183587893774,0.1835878942086,4.346000e-10 diff --git a/benchmarks/cross-library/accuracy_divergence.csv b/benchmarks/cross-library/accuracy_divergence.csv new file mode 100644 index 000000000..669426026 --- /dev/null +++ b/benchmarks/cross-library/accuracy_divergence.csv @@ -0,0 +1,99 @@ +model,library,qubits,steps,dt,seed,atol,ref_atol,ref_terms,ref_norm,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err +heisenberg,ppvm,6,10,0.05,12345,1e-3,1e-16,1024,1.000000e+00,222,4.036339e-02,8.303254e-01,3.634578e-03,4.669530e-01,0,155,0.3832181169286,0.3822563306389,9.617863e-04 +heisenberg,monoprop,6,10,0.05,12345,1e-3,1e-16,1024,1.000000e+00,224,4.214458e-02,8.422781e-01,4.098614e-03,4.650453e-01,2,155,0.383158101684,0.3822563306389,9.017710e-04 +heisenberg,ppvm,6,10,0.05,12345,1e-5,1e-16,1024,1.000000e+00,826,5.404658e-04,1.318947e-02,6.420822e-05,1.748059e-03,18,76,0.3822591781592,0.3822563306389,2.847520e-06 +heisenberg,monoprop,6,10,0.05,12345,1e-5,1e-16,1024,1.000000e+00,906,5.160546e-04,1.273802e-02,6.591156e-05,1.255044e-03,74,52,0.382257798983,0.3822563306389,1.468344e-06 +heisenberg,ppvm,6,20,0.05,12345,1e-3,1e-16,1024,1.000000e+00,866,1.070213e-01,2.716866e+00,1.192541e-02,4.248833e-01,90,130,0.1859379068078,0.1853914774626,5.464293e-04 +heisenberg,monoprop,6,20,0.05,12345,1e-3,1e-16,1024,1.000000e+00,993,1.252160e-01,3.207887e+00,1.363090e-02,3.507914e-01,180,93,0.185871276693,0.1853914774626,4.797992e-04 +heisenberg,ppvm,6,20,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1039,6.429122e-04,1.651200e-02,7.105414e-05,1.100133e-04,20,3,0.1854168229734,0.1853914774626,2.534551e-05 +heisenberg,monoprop,6,20,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1167,7.136553e-04,1.958904e-02,7.192523e-05,1.181164e-04,149,4,0.185416479607,0.1853914774626,2.500214e-05 +heisenberg,ppvm,6,40,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1177,1.474150e-01,4.066484e+00,1.613299e-02,2.478436e-01,239,49,0.1284102386314,0.1308659700642,2.455731e-03 +heisenberg,monoprop,6,40,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1430,2.353190e-01,7.134715e+00,2.271986e-02,2.480553e-01,486,43,0.130211985822,0.1308659700642,6.539842e-04 +heisenberg,ppvm,6,40,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1024,6.308929e-04,1.575149e-02,7.857540e-05,0,0,0,0.130870034274,0.1308659700642,4.064210e-06 +heisenberg,monoprop,6,40,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1099,6.908877e-04,1.833722e-02,6.336108e-05,0,75,0,0.130868964327,0.1308659700642,2.994263e-06 +heisenberg,ppvm,6,80,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1236,1.646609e-01,4.623798e+00,2.112876e-02,1.582509e-01,268,34,0.08887939929655,0.09030014507815,1.420746e-03 +heisenberg,monoprop,6,80,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1638,3.087508e-01,9.415847e+00,3.496292e-02,2.588310e-01,675,39,0.075129187818,0.09030014507815,1.517096e-02 +heisenberg,ppvm,6,80,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1024,6.290098e-04,1.561913e-02,7.908085e-05,0,0,0,0.09032734267639,0.09030014507815,2.719760e-05 +heisenberg,monoprop,6,80,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1035,6.629870e-04,1.658641e-02,7.977269e-05,0,11,0,0.090328492789,0.09030014507815,2.834771e-05 +heisenberg,ppvm,6,160,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1213,1.813172e-01,5.016238e+00,2.172605e-02,1.693484e-01,246,25,0.09555696196697,0.09803624235248,2.479280e-03 +heisenberg,monoprop,6,160,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1691,4.000477e-01,1.195438e+01,4.531704e-02,3.583095e-01,735,36,0.118446290328,0.09803624235248,2.041005e-02 +heisenberg,ppvm,6,160,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1024,6.334569e-04,1.596263e-02,6.675496e-05,0,0,0,0.09802627447504,0.09803624235248,9.967877e-06 +heisenberg,monoprop,6,160,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1022,6.520457e-04,1.635314e-02,7.891152e-05,7.021406e-05,0,2,0.098031933021,0.09803624235248,4.309331e-06 +heisenberg,ppvm,6,320,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1235,2.375216e-01,6.534776e+00,2.462977e-02,3.085507e-01,275,46,0.2136063818482,0.2171604004367,3.554019e-03 +heisenberg,monoprop,6,320,0.05,12345,1e-3,1e-16,1024,1.000000e+00,1787,6.576496e-01,1.928416e+01,7.418257e-02,7.182283e-01,821,40,0.192691392802,0.2171604004367,2.446901e-02 +heisenberg,ppvm,6,320,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1024,6.369157e-04,1.611839e-02,7.839926e-05,0,0,0,0.2171604704937,0.2171604004367,7.005700e-08 +heisenberg,monoprop,6,320,0.05,12345,1e-5,1e-16,1024,1.000000e+00,1024,6.487131e-04,1.640341e-02,7.181142e-05,0,0,0,0.217165017466,0.2171604004367,4.617029e-06 +heisenberg,ppvm,8,10,0.02,12345,1e-3,1e-16,16561,1.000000e+00,27,1.365602e-02,1.124382e-01,4.679180e-03,6.805587e-02,0,26,0.8522590829718,0.8522087425321,5.034044e-05 +heisenberg,monoprop,8,10,0.02,12345,1e-3,1e-16,16561,1.000000e+00,25,1.452084e-02,1.197392e-01,5.736855e-03,7.172758e-02,0,28,0.852274606892,0.8522087425321,6.586436e-05 +heisenberg,ppvm,8,10,0.05,12345,1e-3,1e-16,16384,1.000000e+00,222,4.577082e-02,1.518473e+00,3.673186e-03,1.197095e+00,0,263,0.3832181169286,0.3822563787026,9.617382e-04 +heisenberg,monoprop,8,10,0.05,12345,1e-3,1e-16,16384,1.000000e+00,224,4.684098e-02,1.525291e+00,3.921037e-03,1.195189e+00,2,263,0.383158101684,0.3822563787026,9.017230e-04 +heisenberg,ppvm,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,4566,1.935867e-01,1.683728e+01,6.612049e-03,8.467940e+00,370,2994,0.1843817263679,0.1835878942086,7.938322e-04 +heisenberg,monoprop,8,10,0.1,12345,1e-3,1e-16,16384,1.000000e+00,5106,1.937250e-01,1.714593e+01,6.293714e-03,8.315984e+00,836,2920,0.18421262641,0.1835878942086,6.247322e-04 +heisenberg,ppvm,8,10,0.2,12345,1e-3,1e-16,16384,1.000000e+00,9617,3.392238e-01,3.496440e+01,1.721673e-02,1.503669e+01,2670,5303,0.1206720275243,0.1197752470608,8.967805e-04 +heisenberg,monoprop,8,10,0.2,12345,1e-3,1e-16,16384,1.000000e+00,9516,3.656091e-01,3.825462e+01,1.981698e-02,1.708788e+01,3084,5818,0.119257133464,0.1197752470608,5.181136e-04 +heisenberg,ppvm,8,10,0.4,12345,1e-3,1e-16,16384,1.000000e+00,12275,3.954437e-01,4.123504e+01,1.592698e-02,1.312605e+01,2316,3979,0.06411724002328,0.06498038532114,8.631453e-04 +heisenberg,monoprop,8,10,0.4,12345,1e-3,1e-16,16384,1.000000e+00,11524,3.180263e-01,3.235889e+01,1.186265e-02,1.048401e+01,1169,3583,0.064849564292,0.06498038532114,1.308210e-04 +heisenberg,ppvm,8,10,0.6,12345,1e-3,1e-16,16384,1.000000e+00,13027,3.794862e-01,3.972636e+01,1.750531e-02,1.173227e+01,3457,3830,0.1472159226886,0.1492943370647,2.078414e-03 +heisenberg,monoprop,8,10,0.6,12345,1e-3,1e-16,16384,1.000000e+00,12431,2.526846e-01,2.554907e+01,1.121828e-02,6.632730e+00,1625,2594,0.149259337365,0.1492943370647,3.499970e-05 +heisenberg,ppvm,8,10,0.785,12345,1e-3,1e-16,7322,1.000000e+00,1,9.422127e-03,3.564081e-02,3.185207e-03,3.561038e-02,0,14,0,5.707176610115e-06,5.707177e-06 +heisenberg,monoprop,8,10,0.785,12345,1e-3,1e-16,7322,1.000000e+00,1,9.422127e-03,3.564081e-02,3.185207e-03,3.561038e-02,0,14,0,5.707176610115e-06,5.707177e-06 +scramble,ppvm,8,400,0.05,1,1e-3,1e-16,65535,1.000000e+00,7145,9.460647e-01,1.886151e+02,2.247844e-02,1.588307e+02,1273,44679,0.1790007556962,0.1614730211106,1.752773e-02 +scramble,monoprop,8,400,0.05,1,1e-3,1e-16,65535,1.000000e+00,7923,9.619359e-01,1.897827e+02,2.840829e-02,1.554452e+02,1401,44029,0.167740351377,0.1614730211106,6.267330e-03 +scramble,ppvm,8,400,0.05,1,1e-4,1e-16,65535,1.000000e+00,64007,1.366369e-01,2.789267e+01,2.364321e-03,6.780471e-01,1332,1293,0.1614802576756,0.1614730211106,7.236565e-06 +scramble,monoprop,8,400,0.05,1,1e-4,1e-16,65535,1.000000e+00,63920,4.179348e-01,8.505128e+01,9.345567e-03,2.067460e+00,1480,1528,0.165644417085,0.1614730211106,4.171396e-03 +scramble,ppvm,8,400,0.05,1,1e-5,1e-16,65535,1.000000e+00,65388,7.801179e-03,1.593150e+00,1.427393e-04,3.874730e-03,110,119,0.1614469519445,0.1614730211106,2.606917e-05 +scramble,monoprop,8,400,0.05,1,1e-5,1e-16,65535,1.000000e+00,65376,1.940221e-02,3.961147e+00,3.854044e-04,9.243616e-03,122,143,0.161298097215,0.1614730211106,1.749239e-04 +scramble,ppvm,8,400,0.05,2,1e-3,1e-16,65535,1.000000e+00,7277,9.563637e-01,1.915102e+02,2.113910e-02,1.610994e+02,1276,44811,0.145258196215,0.141763190312,3.495006e-03 +scramble,monoprop,8,400,0.05,2,1e-3,1e-16,65535,1.000000e+00,7884,9.710104e-01,1.920922e+02,2.926046e-02,1.580419e+02,1338,44266,0.121961388302,0.141763190312,1.980180e-02 +scramble,ppvm,8,400,0.05,2,1e-4,1e-16,65535,1.000000e+00,64004,1.406716e-01,2.866367e+01,2.309150e-03,6.676462e-01,1275,1297,0.1411270619147,0.141763190312,6.361284e-04 +scramble,monoprop,8,400,0.05,2,1e-4,1e-16,65535,1.000000e+00,63950,4.073952e-01,8.292705e+01,7.614204e-03,1.971686e+00,1420,1496,0.138618759706,0.141763190312,3.144431e-03 +scramble,ppvm,8,400,0.05,2,1e-5,1e-16,65535,1.000000e+00,65400,8.667319e-03,1.771308e+00,1.364485e-04,3.740238e-03,106,96,0.1417936293827,0.141763190312,3.043907e-05 +scramble,monoprop,8,400,0.05,2,1e-5,1e-16,65535,1.000000e+00,65398,1.939419e-02,3.962640e+00,3.397129e-04,8.654144e-03,124,116,0.141898202448,0.141763190312,1.350121e-04 +scramble,ppvm,8,400,0.05,3,1e-3,1e-16,65535,1.000000e+00,7384,9.450019e-01,1.887800e+02,2.825414e-02,1.580012e+02,1296,44432,0.1524323856682,0.1460185792372,6.413806e-03 +scramble,monoprop,8,400,0.05,3,1e-3,1e-16,65535,1.000000e+00,8026,9.643916e-01,1.901604e+02,3.122063e-02,1.550255e+02,1377,43871,0.142537145171,0.1460185792372,3.481434e-03 +scramble,ppvm,8,400,0.05,3,1e-4,1e-16,65535,1.000000e+00,63904,1.404446e-01,2.871432e+01,2.399302e-03,7.247692e-01,1277,1410,0.1452846799804,0.1460185792372,7.338993e-04 +scramble,monoprop,8,400,0.05,3,1e-4,1e-16,65535,1.000000e+00,64008,4.156444e-01,8.458307e+01,8.494982e-03,1.963916e+00,1432,1461,0.140959989003,0.1460185792372,5.058590e-03 +scramble,ppvm,8,400,0.05,3,1e-5,1e-16,65535,1.000000e+00,65376,8.038709e-03,1.642135e+00,1.416798e-04,4.002692e-03,101,117,0.1460074411703,0.1460185792372,1.113807e-05 +scramble,monoprop,8,400,0.05,3,1e-5,1e-16,65535,1.000000e+00,65378,1.951992e-02,3.980961e+00,3.259493e-04,9.704331e-03,124,138,0.146049390259,0.1460185792372,3.081102e-05 +scramble,ppvm,8,400,0.05,4,1e-3,1e-16,65535,1.000000e+00,7409,9.518168e-01,1.894747e+02,2.758227e-02,1.578000e+02,1251,44408,0.1682610320682,0.1693302302423,1.069198e-03 +scramble,monoprop,8,400,0.05,4,1e-3,1e-16,65535,1.000000e+00,8156,9.682760e-01,1.903117e+02,4.164777e-02,1.539202e+02,1307,43717,0.144043912834,0.1693302302423,2.528632e-02 +scramble,ppvm,8,400,0.05,4,1e-4,1e-16,65535,1.000000e+00,63869,1.409938e-01,2.882227e+01,2.557839e-03,7.418089e-01,1298,1408,0.1686796776657,0.1693302302423,6.505526e-04 +scramble,monoprop,8,400,0.05,4,1e-4,1e-16,65535,1.000000e+00,63968,4.106727e-01,8.365478e+01,7.280307e-03,1.928183e+00,1464,1475,0.165243860705,0.1693302302423,4.086370e-03 +scramble,ppvm,8,400,0.05,4,1e-5,1e-16,65535,1.000000e+00,65376,8.336955e-03,1.701727e+00,1.487382e-04,4.134139e-03,123,126,0.1693221010763,0.1693302302423,8.129166e-06 +scramble,monoprop,8,400,0.05,4,1e-5,1e-16,65535,1.000000e+00,65368,2.026217e-02,4.130641e+00,3.527960e-04,1.150466e-02,141,152,0.169290141837,0.1693302302423,4.008841e-05 +scramble,ppvm,8,400,0.05,5,1e-3,1e-16,65535,1.000000e+00,7168,9.504922e-01,1.888251e+02,2.381869e-02,1.581427e+02,1314,44541,0.1488295179599,0.1622150514644,1.338553e-02 +scramble,monoprop,8,400,0.05,5,1e-3,1e-16,65535,1.000000e+00,8263,9.733226e-01,1.903354e+02,4.313585e-02,1.532102e+02,1398,43530,0.187629107769,0.1622150514644,2.541406e-02 +scramble,ppvm,8,400,0.05,5,1e-4,1e-16,65535,1.000000e+00,63819,1.413663e-01,2.886298e+01,2.541082e-03,7.565301e-01,1294,1459,0.1601476023761,0.1622150514644,2.067449e-03 +scramble,monoprop,8,400,0.05,5,1e-4,1e-16,65535,1.000000e+00,63888,4.233586e-01,8.599954e+01,8.111627e-03,2.104861e+00,1458,1554,0.154900547436,0.1622150514644,7.314504e-03 +scramble,ppvm,8,400,0.05,5,1e-5,1e-16,65535,1.000000e+00,65372,7.935054e-03,1.618603e+00,1.321421e-04,4.264874e-03,115,128,0.1621842282147,0.1622150514644,3.082325e-05 +scramble,monoprop,8,400,0.05,5,1e-5,1e-16,65535,1.000000e+00,65383,1.933189e-02,3.942074e+00,3.162846e-04,9.839686e-03,141,143,0.162116789286,0.1622150514644,9.826218e-05 +scramble,ppvm,8,10,0.7853981633974483,1,1e-3,1e-16,65535,1.000000e+00,50657,3.090565e-01,6.278704e+01,5.755294e-03,1.610769e+01,5920,6774,0.002947642913682,0.004567937262123,1.620294e-03 +scramble,monoprop,8,10,0.7853981633974483,1,1e-3,1e-16,65535,1.000000e+00,49970,2.637431e-01,5.349200e+01,5.279884e-03,1.500207e+01,4766,6307,0.002980535216,0.004567937262123,1.587402e-03 +scramble,ppvm,8,10,0.7853981633974483,1,1e-4,1e-16,65535,1.000000e+00,64150,1.674435e-02,3.406150e+00,3.050943e-04,9.712100e-02,374,342,0.004554433189204,0.004567937262123,1.350407e-05 +scramble,monoprop,8,10,0.7853981633974483,1,1e-4,1e-16,65535,1.000000e+00,64144,1.270353e-02,2.583770e+00,2.536806e-04,8.731900e-02,294,268,0.004512298994,0.004567937262123,5.563827e-05 +scramble,ppvm,8,10,0.7853981633974483,1,1e-5,1e-16,65535,1.000000e+00,65396,1.041360e-03,2.120041e-01,2.233626e-05,7.722638e-04,26,21,0.004568585960411,0.004567937262123,6.486983e-07 +scramble,monoprop,8,10,0.7853981633974483,1,1e-5,1e-16,65535,1.000000e+00,65392,9.593503e-04,1.913561e-01,2.030360e-05,7.639987e-04,18,17,0.004566272673,0.004567937262123,1.664589e-06 +scramble,ppvm,8,10,0.7853981633974483,2,1e-3,1e-16,65535,1.000000e+00,48647,2.791048e-01,5.628207e+01,5.037719e-03,1.722668e+01,5167,7348,-0.003044049375097,-0.002123511847918,9.205375e-04 +scramble,monoprop,8,10,0.7853981633974483,2,1e-3,1e-16,65535,1.000000e+00,49443,2.734635e-01,5.507724e+01,5.526536e-03,1.610547e+01,5454,6839,-0.002009768618,-0.002123511847918,1.137432e-04 +scramble,ppvm,8,10,0.7853981633974483,2,1e-4,1e-16,65535,1.000000e+00,64060,1.508072e-02,3.042391e+00,3.327267e-04,9.911472e-02,394,348,-0.002128546499306,-0.002123511847918,5.034651e-06 +scramble,monoprop,8,10,0.7853981633974483,2,1e-4,1e-16,65535,1.000000e+00,64063,1.363255e-02,2.758181e+00,3.033626e-04,9.342819e-02,353,304,-0.002116198144,-0.002123511847918,7.313704e-06 +scramble,ppvm,8,10,0.7853981633974483,2,1e-5,1e-16,65535,1.000000e+00,65377,8.178769e-04,1.652772e-01,1.911102e-05,8.209238e-04,27,15,-0.002124049692306,-0.002123511847918,5.378444e-07 +scramble,monoprop,8,10,0.7853981633974483,2,1e-5,1e-16,65535,1.000000e+00,65362,6.915502e-04,1.389042e-01,2.041713e-05,9.752001e-04,18,21,-0.002124065586,-0.002123511847918,5.537381e-07 +scramble,ppvm,8,10,0.7853981633974483,3,1e-3,1e-16,65534,1.000000e+00,23800,1.482153e-01,2.838560e+01,3.781926e-03,1.934654e+01,2061,4604,0.01021329788007,0.01011042655569,1.028713e-04 +scramble,monoprop,8,10,0.7853981633974483,3,1e-3,1e-16,65534,1.000000e+00,23958,1.416023e-01,2.720370e+01,3.105871e-03,1.887551e+01,1907,4292,0.010168429392,0.01011042655569,5.800284e-05 +scramble,ppvm,8,10,0.7853981633974483,3,1e-4,1e-16,65534,1.000000e+00,57545,1.677474e-02,3.310245e+00,4.291547e-04,5.726769e-01,1649,2112,0.01012727159909,0.01011042655569,1.684504e-05 +scramble,monoprop,8,10,0.7853981633974483,3,1e-4,1e-16,65534,1.000000e+00,57569,1.497634e-02,2.933245e+00,4.037880e-04,5.407710e-01,1440,1879,0.010170823571,0.01011042655569,6.039702e-05 +scramble,ppvm,8,10,0.7853981633974483,3,1e-5,1e-16,65534,1.000000e+00,64720,1.185833e-03,2.343507e-01,2.743055e-05,5.164075e-03,175,166,0.01011032002169,0.01011042655569,1.065340e-07 +scramble,monoprop,8,10,0.7853981633974483,3,1e-5,1e-16,65534,1.000000e+00,64695,9.601826e-04,1.883416e-01,2.589419e-05,5.002271e-03,133,149,0.010108359939,0.01011042655569,2.066617e-06 +scramble,ppvm,8,10,0.7853981633974483,4,1e-3,1e-16,65535,1.000000e+00,44916,2.344832e-01,4.721697e+01,4.659655e-03,1.772357e+01,4861,7238,-0.002857638598484,-0.002233897581292,6.237410e-04 +scramble,monoprop,8,10,0.7853981633974483,4,1e-3,1e-16,65535,1.000000e+00,44165,2.366092e-01,4.750515e+01,4.174483e-03,1.888975e+01,4683,7811,-0.002729812368,-0.002233897581292,4.959148e-04 +scramble,ppvm,8,10,0.7853981633974483,4,1e-4,1e-16,65535,1.000000e+00,63664,1.656759e-02,3.339823e+00,3.475736e-04,1.364974e-01,544,528,-0.002347368442479,-0.002233897581292,1.134709e-04 +scramble,monoprop,8,10,0.7853981633974483,4,1e-4,1e-16,65535,1.000000e+00,63678,1.610363e-02,3.253601e+00,3.508339e-04,1.298318e-01,494,464,-0.002305939048,-0.002233897581292,7.204147e-05 +scramble,ppvm,8,10,0.7853981633974483,4,1e-5,1e-16,65535,1.000000e+00,65344,1.076490e-03,2.173143e-01,2.137590e-05,1.167930e-03,44,36,-0.0022253215307,-0.002233897581292,8.576051e-06 +scramble,monoprop,8,10,0.7853981633974483,4,1e-5,1e-16,65535,1.000000e+00,65344,9.685830e-04,1.950917e-01,2.328722e-05,1.113917e-03,37,29,-0.002225810421,-0.002233897581292,8.087160e-06 +scramble,ppvm,8,10,0.7853981633974483,5,1e-3,1e-16,65535,1.000000e+00,46519,2.524708e-01,5.087362e+01,4.986171e-03,1.754877e+01,5039,7322,0,-0.001223130468711,1.223130e-03 +scramble,monoprop,8,10,0.7853981633974483,5,1e-3,1e-16,65535,1.000000e+00,47075,2.425934e-01,4.885867e+01,4.802603e-03,1.636502e+01,5028,6755,0,-0.001223130468711,1.223130e-03 +scramble,ppvm,8,10,0.7853981633974483,5,1e-4,1e-16,65535,1.000000e+00,63800,1.575686e-02,3.174853e+00,3.102006e-04,1.252977e-01,441,469,-0.001210275579338,-0.001223130468711,1.285489e-05 +scramble,monoprop,8,10,0.7853981633974483,5,1e-4,1e-16,65535,1.000000e+00,63817,1.319639e-02,2.629690e+00,2.919619e-04,1.111590e-01,347,358,-0.001235391553,-0.001223130468711,1.226108e-05 +scramble,ppvm,8,10,0.7853981633974483,5,1e-5,1e-16,65535,1.000000e+00,65366,8.791609e-04,1.777130e-01,2.111646e-05,8.729767e-04,28,16,-0.001221170841282,-0.001223130468711,1.959627e-06 +scramble,monoprop,8,10,0.7853981633974483,5,1e-5,1e-16,65535,1.000000e+00,65365,6.298597e-04,1.269766e-01,1.810326e-05,8.522953e-04,22,11,-0.001224291877,-0.001223130468711,1.161408e-06 diff --git a/benchmarks/cross-library/accuracy_scramble.csv b/benchmarks/cross-library/accuracy_scramble.csv deleted file mode 100644 index dd778a693..000000000 --- a/benchmarks/cross-library/accuracy_scramble.csv +++ /dev/null @@ -1,61 +0,0 @@ -model,library,qubits,steps,dt,seed,atol,ref_atol,ref_terms,terms,l2_err,l1_err,max_coeff_err,lost_mass,kept_subthreshold,dropped_above_atol,observable,ref_observable,obs_abs_err -scramble,ppvm,8,400,0.05,1,1e-3,1e-16,65535,7145,9.460647e-01,1.886151e+02,2.247844e-02,1.588307e+02,1273,44679,0.1790007556962,0.1614730211106,1.752773e-02 -scramble,monoprop,8,400,0.05,1,1e-3,1e-16,65535,7923,9.619359e-01,1.897827e+02,2.840829e-02,1.554452e+02,1401,44029,0.167740351377,0.1614730211106,6.267330e-03 -scramble,ppvm,8,400,0.05,1,1e-4,1e-16,65535,64007,1.366369e-01,2.789267e+01,2.364321e-03,6.780471e-01,1332,1293,0.1614802576756,0.1614730211106,7.236565e-06 -scramble,monoprop,8,400,0.05,1,1e-4,1e-16,65535,63920,4.179348e-01,8.505128e+01,9.345567e-03,2.067460e+00,1480,1528,0.165644417085,0.1614730211106,4.171396e-03 -scramble,ppvm,8,400,0.05,1,1e-5,1e-16,65535,65388,7.801179e-03,1.593150e+00,1.427393e-04,3.874730e-03,110,119,0.1614469519445,0.1614730211106,2.606917e-05 -scramble,monoprop,8,400,0.05,1,1e-5,1e-16,65535,65376,1.940221e-02,3.961147e+00,3.854044e-04,9.243616e-03,122,143,0.161298097215,0.1614730211106,1.749239e-04 -scramble,ppvm,8,400,0.05,2,1e-3,1e-16,65535,7277,9.563637e-01,1.915102e+02,2.113910e-02,1.610994e+02,1276,44811,0.145258196215,0.141763190312,3.495006e-03 -scramble,monoprop,8,400,0.05,2,1e-3,1e-16,65535,7884,9.710104e-01,1.920922e+02,2.926046e-02,1.580419e+02,1338,44266,0.121961388302,0.141763190312,1.980180e-02 -scramble,ppvm,8,400,0.05,2,1e-4,1e-16,65535,64004,1.406716e-01,2.866367e+01,2.309150e-03,6.676462e-01,1275,1297,0.1411270619147,0.141763190312,6.361284e-04 -scramble,monoprop,8,400,0.05,2,1e-4,1e-16,65535,63950,4.073952e-01,8.292705e+01,7.614204e-03,1.971686e+00,1420,1496,0.138618759706,0.141763190312,3.144431e-03 -scramble,ppvm,8,400,0.05,2,1e-5,1e-16,65535,65400,8.667319e-03,1.771308e+00,1.364485e-04,3.740238e-03,106,96,0.1417936293827,0.141763190312,3.043907e-05 -scramble,monoprop,8,400,0.05,2,1e-5,1e-16,65535,65398,1.939419e-02,3.962640e+00,3.397129e-04,8.654144e-03,124,116,0.141898202448,0.141763190312,1.350121e-04 -scramble,ppvm,8,400,0.05,3,1e-3,1e-16,65535,7384,9.450019e-01,1.887800e+02,2.825414e-02,1.580012e+02,1296,44432,0.1524323856682,0.1460185792372,6.413806e-03 -scramble,monoprop,8,400,0.05,3,1e-3,1e-16,65535,8026,9.643916e-01,1.901604e+02,3.122063e-02,1.550255e+02,1377,43871,0.142537145171,0.1460185792372,3.481434e-03 -scramble,ppvm,8,400,0.05,3,1e-4,1e-16,65535,63904,1.404446e-01,2.871432e+01,2.399302e-03,7.247692e-01,1277,1410,0.1452846799804,0.1460185792372,7.338993e-04 -scramble,monoprop,8,400,0.05,3,1e-4,1e-16,65535,64008,4.156444e-01,8.458307e+01,8.494982e-03,1.963916e+00,1432,1461,0.140959989003,0.1460185792372,5.058590e-03 -scramble,ppvm,8,400,0.05,3,1e-5,1e-16,65535,65376,8.038709e-03,1.642135e+00,1.416798e-04,4.002692e-03,101,117,0.1460074411703,0.1460185792372,1.113807e-05 -scramble,monoprop,8,400,0.05,3,1e-5,1e-16,65535,65378,1.951992e-02,3.980961e+00,3.259493e-04,9.704331e-03,124,138,0.146049390259,0.1460185792372,3.081102e-05 -scramble,ppvm,8,400,0.05,4,1e-3,1e-16,65535,7409,9.518168e-01,1.894747e+02,2.758227e-02,1.578000e+02,1251,44408,0.1682610320682,0.1693302302423,1.069198e-03 -scramble,monoprop,8,400,0.05,4,1e-3,1e-16,65535,8156,9.682760e-01,1.903117e+02,4.164777e-02,1.539202e+02,1307,43717,0.144043912834,0.1693302302423,2.528632e-02 -scramble,ppvm,8,400,0.05,4,1e-4,1e-16,65535,63869,1.409938e-01,2.882227e+01,2.557839e-03,7.418089e-01,1298,1408,0.1686796776657,0.1693302302423,6.505526e-04 -scramble,monoprop,8,400,0.05,4,1e-4,1e-16,65535,63968,4.106727e-01,8.365478e+01,7.280307e-03,1.928183e+00,1464,1475,0.165243860705,0.1693302302423,4.086370e-03 -scramble,ppvm,8,400,0.05,4,1e-5,1e-16,65535,65376,8.336955e-03,1.701727e+00,1.487382e-04,4.134139e-03,123,126,0.1693221010763,0.1693302302423,8.129166e-06 -scramble,monoprop,8,400,0.05,4,1e-5,1e-16,65535,65368,2.026217e-02,4.130641e+00,3.527960e-04,1.150466e-02,141,152,0.169290141837,0.1693302302423,4.008841e-05 -scramble,ppvm,8,400,0.05,5,1e-3,1e-16,65535,7168,9.504922e-01,1.888251e+02,2.381869e-02,1.581427e+02,1314,44541,0.1488295179599,0.1622150514644,1.338553e-02 -scramble,monoprop,8,400,0.05,5,1e-3,1e-16,65535,8263,9.733226e-01,1.903354e+02,4.313585e-02,1.532102e+02,1398,43530,0.187629107769,0.1622150514644,2.541406e-02 -scramble,ppvm,8,400,0.05,5,1e-4,1e-16,65535,63819,1.413663e-01,2.886298e+01,2.541082e-03,7.565301e-01,1294,1459,0.1601476023761,0.1622150514644,2.067449e-03 -scramble,monoprop,8,400,0.05,5,1e-4,1e-16,65535,63888,4.233586e-01,8.599954e+01,8.111627e-03,2.104861e+00,1458,1554,0.154900547436,0.1622150514644,7.314504e-03 -scramble,ppvm,8,400,0.05,5,1e-5,1e-16,65535,65372,7.935054e-03,1.618603e+00,1.321421e-04,4.264874e-03,115,128,0.1621842282147,0.1622150514644,3.082325e-05 -scramble,monoprop,8,400,0.05,5,1e-5,1e-16,65535,65383,1.933189e-02,3.942074e+00,3.162846e-04,9.839686e-03,141,143,0.162116789286,0.1622150514644,9.826218e-05 -scramble,ppvm,8,10,0.7853981633974483,1,1e-3,1e-16,65535,50657,3.090565e-01,6.278704e+01,5.755294e-03,1.610769e+01,5920,6774,0.002947642913682,0.004567937262123,1.620294e-03 -scramble,monoprop,8,10,0.7853981633974483,1,1e-3,1e-16,65535,49970,2.637431e-01,5.349200e+01,5.279884e-03,1.500207e+01,4766,6307,0.002980535216,0.004567937262123,1.587402e-03 -scramble,ppvm,8,10,0.7853981633974483,1,1e-4,1e-16,65535,64150,1.674435e-02,3.406150e+00,3.050943e-04,9.712100e-02,374,342,0.004554433189204,0.004567937262123,1.350407e-05 -scramble,monoprop,8,10,0.7853981633974483,1,1e-4,1e-16,65535,64144,1.270353e-02,2.583770e+00,2.536806e-04,8.731900e-02,294,268,0.004512298994,0.004567937262123,5.563827e-05 -scramble,ppvm,8,10,0.7853981633974483,1,1e-5,1e-16,65535,65396,1.041360e-03,2.120041e-01,2.233626e-05,7.722638e-04,26,21,0.004568585960411,0.004567937262123,6.486983e-07 -scramble,monoprop,8,10,0.7853981633974483,1,1e-5,1e-16,65535,65392,9.593503e-04,1.913561e-01,2.030360e-05,7.639987e-04,18,17,0.004566272673,0.004567937262123,1.664589e-06 -scramble,ppvm,8,10,0.7853981633974483,2,1e-3,1e-16,65535,48647,2.791048e-01,5.628207e+01,5.037719e-03,1.722668e+01,5167,7348,-0.003044049375097,-0.002123511847918,9.205375e-04 -scramble,monoprop,8,10,0.7853981633974483,2,1e-3,1e-16,65535,49443,2.734635e-01,5.507724e+01,5.526536e-03,1.610547e+01,5454,6839,-0.002009768618,-0.002123511847918,1.137432e-04 -scramble,ppvm,8,10,0.7853981633974483,2,1e-4,1e-16,65535,64060,1.508072e-02,3.042391e+00,3.327267e-04,9.911472e-02,394,348,-0.002128546499306,-0.002123511847918,5.034651e-06 -scramble,monoprop,8,10,0.7853981633974483,2,1e-4,1e-16,65535,64063,1.363255e-02,2.758181e+00,3.033626e-04,9.342819e-02,353,304,-0.002116198144,-0.002123511847918,7.313704e-06 -scramble,ppvm,8,10,0.7853981633974483,2,1e-5,1e-16,65535,65377,8.178769e-04,1.652772e-01,1.911102e-05,8.209238e-04,27,15,-0.002124049692306,-0.002123511847918,5.378444e-07 -scramble,monoprop,8,10,0.7853981633974483,2,1e-5,1e-16,65535,65362,6.915502e-04,1.389042e-01,2.041713e-05,9.752001e-04,18,21,-0.002124065586,-0.002123511847918,5.537381e-07 -scramble,ppvm,8,10,0.7853981633974483,3,1e-3,1e-16,65534,23800,1.482153e-01,2.838560e+01,3.781926e-03,1.934654e+01,2061,4604,0.01021329788007,0.01011042655569,1.028713e-04 -scramble,monoprop,8,10,0.7853981633974483,3,1e-3,1e-16,65534,23958,1.416023e-01,2.720370e+01,3.105871e-03,1.887551e+01,1907,4292,0.010168429392,0.01011042655569,5.800284e-05 -scramble,ppvm,8,10,0.7853981633974483,3,1e-4,1e-16,65534,57545,1.677474e-02,3.310245e+00,4.291547e-04,5.726769e-01,1649,2112,0.01012727159909,0.01011042655569,1.684504e-05 -scramble,monoprop,8,10,0.7853981633974483,3,1e-4,1e-16,65534,57569,1.497634e-02,2.933245e+00,4.037880e-04,5.407710e-01,1440,1879,0.010170823571,0.01011042655569,6.039702e-05 -scramble,ppvm,8,10,0.7853981633974483,3,1e-5,1e-16,65534,64720,1.185833e-03,2.343507e-01,2.743055e-05,5.164075e-03,175,166,0.01011032002169,0.01011042655569,1.065340e-07 -scramble,monoprop,8,10,0.7853981633974483,3,1e-5,1e-16,65534,64695,9.601826e-04,1.883416e-01,2.589419e-05,5.002271e-03,133,149,0.010108359939,0.01011042655569,2.066617e-06 -scramble,ppvm,8,10,0.7853981633974483,4,1e-3,1e-16,65535,44916,2.344832e-01,4.721697e+01,4.659655e-03,1.772357e+01,4861,7238,-0.002857638598484,-0.002233897581292,6.237410e-04 -scramble,monoprop,8,10,0.7853981633974483,4,1e-3,1e-16,65535,44165,2.366092e-01,4.750515e+01,4.174483e-03,1.888975e+01,4683,7811,-0.002729812368,-0.002233897581292,4.959148e-04 -scramble,ppvm,8,10,0.7853981633974483,4,1e-4,1e-16,65535,63664,1.656759e-02,3.339823e+00,3.475736e-04,1.364974e-01,544,528,-0.002347368442479,-0.002233897581292,1.134709e-04 -scramble,monoprop,8,10,0.7853981633974483,4,1e-4,1e-16,65535,63678,1.610363e-02,3.253601e+00,3.508339e-04,1.298318e-01,494,464,-0.002305939048,-0.002233897581292,7.204147e-05 -scramble,ppvm,8,10,0.7853981633974483,4,1e-5,1e-16,65535,65344,1.076490e-03,2.173143e-01,2.137590e-05,1.167930e-03,44,36,-0.0022253215307,-0.002233897581292,8.576051e-06 -scramble,monoprop,8,10,0.7853981633974483,4,1e-5,1e-16,65535,65344,9.685830e-04,1.950917e-01,2.328722e-05,1.113917e-03,37,29,-0.002225810421,-0.002233897581292,8.087160e-06 -scramble,ppvm,8,10,0.7853981633974483,5,1e-3,1e-16,65535,46519,2.524708e-01,5.087362e+01,4.986171e-03,1.754877e+01,5039,7322,0,-0.001223130468711,1.223130e-03 -scramble,monoprop,8,10,0.7853981633974483,5,1e-3,1e-16,65535,47075,2.425934e-01,4.885867e+01,4.802603e-03,1.636502e+01,5028,6755,0,-0.001223130468711,1.223130e-03 -scramble,ppvm,8,10,0.7853981633974483,5,1e-4,1e-16,65535,63800,1.575686e-02,3.174853e+00,3.102006e-04,1.252977e-01,441,469,-0.001210275579338,-0.001223130468711,1.285489e-05 -scramble,monoprop,8,10,0.7853981633974483,5,1e-4,1e-16,65535,63817,1.319639e-02,2.629690e+00,2.919619e-04,1.111590e-01,347,358,-0.001235391553,-0.001223130468711,1.226108e-05 -scramble,ppvm,8,10,0.7853981633974483,5,1e-5,1e-16,65535,65366,8.791609e-04,1.777130e-01,2.111646e-05,8.729767e-04,28,16,-0.001221170841282,-0.001223130468711,1.959627e-06 -scramble,monoprop,8,10,0.7853981633974483,5,1e-5,1e-16,65535,65365,6.298597e-04,1.269766e-01,1.810326e-05,8.522953e-04,22,11,-0.001224291877,-0.001223130468711,1.161408e-06 diff --git a/benchmarks/cross-library/plot_accuracy.py b/benchmarks/cross-library/plot_accuracy.py new file mode 100644 index 000000000..0ec9135ec --- /dev/null +++ b/benchmarks/cross-library/plot_accuracy.py @@ -0,0 +1,315 @@ +# SPDX-FileCopyrightText: 2026 The PPVM Authors +# SPDX-License-Identifier: Apache-2.0 +"""Render the truncation-error scaling from `xbench_accuracy.py`'s CSVs. + +Four panels, because four separate questions were asked of this data and each one +has a different answer: + +1. **Error vs `atol`** — every engine falls one decade per decade of `atol`, so + the two truncation rules share a convergence order and differ only in the + constant. This panel is why the rules are not qualitatively apart. +2. **Error ratio vs depth** — the loss compounds per gate, so the gap grows with + circuit depth. +3. **Error ratio vs rotation angle** — the ratio *crosses 1.0*. Large angles + favour monoprop, because that is where accumulate-then-truncate has its own + failure mode: `truncate()` after each gate permanently deletes a term whose + merged coefficient transiently cancelled. Any summary quoting one ratio is + hiding this panel. +4. **The scrambling instance** — where the rules genuinely separate, beside a + `θ = π/4` control on the same circuit family that reverses the sign. + +Ratios are monoprop / ppvm on the `L2` distance to a converged reference, so +above 1.0 means monoprop is less accurate. The scalar `observable` is +deliberately not plotted: its truncation error changes sign, so ratios built from +it are dominated by whichever engine happened to land near a zero crossing (see +`README.md`). + + uv run --no-project --with matplotlib python3 plot_accuracy.py \ + --csv accuracy.csv --divergence accuracy_divergence.csv \ + --out target/xbench/accuracy.png +""" + +from __future__ import annotations + +import argparse +import csv +import statistics +import textwrap +from collections import defaultdict +from pathlib import Path + +import matplotlib + +matplotlib.use("Agg") +import matplotlib.pyplot as plt +from plot_xbench import SERIES, SURFACE, TEXT_PRIMARY, TEXT_SECONDARY, style + +BASELINE = "ppvm" +PAIR = ("ppvm", "monoprop") +# Ratio panels plot a *derived* quantity, not an engine, so they stay off the +# engine palette — those hues are reserved for identity throughout the harness. +RATIO_INK = ["#0b0b0b", "#8a8986"] +RATIO_MARKS = ["v", "o"] + + +def load(path: Path) -> list[dict[str, str]]: + with path.open() as fh: + return list(csv.DictReader(fh)) + + +def relative(rows, library: str) -> dict[float, list[float]]: + """`{atol: [relative L2 error, ...]}` for one engine, pooled over seeds.""" + out: dict[float, list[float]] = defaultdict(list) + for r in rows: + if r["library"] == library: + out[float(r["atol"])].append(float(r["l2_err"]) / float(r["ref_norm"])) + return out + + +def ratio_vs(rows, key: str) -> dict[float, list[float]]: + """`{x: [monoprop/ppvm error ratio, ...]}` as `key` is swept. + + A ratio is only formed inside one fully-specified cell — same model, width, + depth, angle, seed and `atol` — so sweeping one axis never divides across + another. Remaining axes (seeds) pool into the list. + """ + cells: dict[tuple, dict[str, float]] = defaultdict(dict) + for r in rows: + if r["library"] not in PAIR: + continue + cell = (r["model"], r["qubits"], r["steps"], r["dt"], r["seed"], r["atol"]) + cells[cell][r["library"]] = float(r["l2_err"]) + cells[cell]["x"] = float(r[key]) + out: dict[float, list[float]] = defaultdict(list) + for libs in cells.values(): + if len(libs) < 3 or not libs.get(BASELINE): + continue + out[libs["x"]].append(libs["monoprop"] / libs[BASELINE]) + return dict(sorted(out.items())) + + +def band(ax, data, color, marker, label) -> None: + """Median line, with a min/max ribbon where an axis pooled several seeds.""" + xs = list(data) + ax.plot( + xs, + [statistics.median(data[x]) for x in xs], + color=color, + marker=marker, + markersize=5, + linewidth=1.6, + label=label, + ) + if any(len(v) > 1 for v in data.values()): + ax.fill_between( + xs, + [min(data[x]) for x in xs], + [max(data[x]) for x in xs], + color=color, + alpha=0.16, + linewidth=0, + ) + + +def panel_atol(ax, rows) -> None: + """Absolute error against atol, every engine, both Trotter models.""" + for model, dash in (("tfim", "-"), ("heisenberg", "--")): + for lib, (label, color, marker) in SERIES.items(): + pts = sorted( + (float(r["atol"]), float(r["l2_err"]) / float(r["ref_norm"])) + for r in rows + if r["model"] == model and r["library"] == lib + ) + if pts: + ax.plot( + *zip(*pts), + dash, + color=color, + marker=marker, + markersize=4, + linewidth=1.4, + alpha=0.9, + label=label if model == "tfim" else None, + ) + # Slope-1 guide, anchored on the baseline's own worst point so it sits beside + # the data rather than floating below it -- the eye is comparing gradients. + anchor = max( + (float(r["atol"]), float(r["l2_err"]) / float(r["ref_norm"])) + for r in rows + if r["model"] == "tfim" and r["library"] == BASELINE + ) + lo = min(float(r["atol"]) for r in rows) + ax.plot( + [lo, anchor[0]], + [anchor[1] * lo / anchor[0], anchor[1]], + ":", + color=TEXT_SECONDARY, + linewidth=1.3, + label="slope 1 (error ∝ atol)", + ) + ax.set_xscale("log") + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("truncation threshold atol") + ax.set_ylabel("relative L2 error ‖c − c*‖ / ‖c*‖") + ax.set_title( + "1 · Every rule converges at the same order\n" + "solid TFIM, dashed Heisenberg — n=8, 10 steps, dt=0.1", + fontsize=10, + color=TEXT_PRIMARY, + loc="left", + ) + ax.legend(fontsize=7, frameon=False, labelcolor=TEXT_SECONDARY, loc="lower left") + + +def panel_ratio(ax, rows, key, xlabel, title) -> None: + """monoprop / ppvm error ratio against one swept axis, one line per atol. + + The swept values are logarithmically spaced and uneven, so the axis is log + and ticks are pinned to the values actually run. + """ + ax.set_xscale("log") + ax.minorticks_off() + # Loosest threshold first: that is the series with the largest effect, and + # `sorted` on the raw strings would put 1e-5 ahead of 1e-3. + xs: list[float] = [] + for i, atol in enumerate( + sorted({r["atol"] for r in rows}, key=float, reverse=True) + ): + data = ratio_vs([r for r in rows if r["atol"] == atol], key) + if data: + band( + ax, + data, + RATIO_INK[i % len(RATIO_INK)], + RATIO_MARKS[i % len(RATIO_MARKS)], + f"atol = {atol}", + ) + xs = list(data) + # Set after the scale, which would otherwise install its own locator. + ax.set_xticks(xs) + ax.set_xticklabels([f"{x:g}" for x in xs], fontsize=8) + ax.axhline(1.0, color=SERIES["ppvm"][1], linewidth=1.2, linestyle=":") + ax.annotate( + "equal accuracy", + (0.03, 1.0), + xycoords=("axes fraction", "data"), + fontsize=7, + color=SERIES["ppvm"][1], + va="bottom", + ) + ax.set_xlabel(xlabel) + ax.set_ylabel("monoprop error / ppvm error") + ax.set_title(title, fontsize=10, color=TEXT_PRIMARY, loc="left") + ax.legend(fontsize=7, frameon=False, labelcolor=TEXT_SECONDARY) + + +def panel_scramble(ax, rows) -> None: + """The scrambling instance, small angle against the theta = pi/4 control.""" + for dt, dash, tag in ( + ("0.05", "-", "dt=0.05, 3200 gates"), + ("0.7853981633974483", "--", "dt=π/4, 80 gates"), + ): + sub = [r for r in rows if r["dt"] == dt] + for lib in PAIR: + data = relative(sub, lib) + if data: + label, color, marker = SERIES[lib] + ax.plot( + sorted(data), + [statistics.median(data[x]) for x in sorted(data)], + dash, + color=color, + marker=marker, + markersize=5, + linewidth=1.6, + label=f"{label} — {tag}", + ) + ax.set_xscale("log") + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("truncation threshold atol") + ax.set_ylabel("relative L2 error") + ax.set_title( + "4 · Scrambling instance: where they separate\n" + "n=8, all 65535 Pauli words populated, median of 5 seeds", + fontsize=10, + color=TEXT_PRIMARY, + loc="left", + ) + ax.legend(fontsize=7, frameon=False, labelcolor=TEXT_SECONDARY, loc="lower left") + + +CAPTION = ( + "Error is the L2 distance of the whole propagated coefficient vector from a " + "converged reference (ppvm at atol=1e-16), relative to ‖c*‖. ppvm, " + "PauliPropagation.jl and PauliStrings.jl accumulate every contribution to a term " + "and then threshold the merged coefficient; pauli-prop and monoprop threshold each " + "branch as it is emitted. A Pauli rotation sends Q to exactly two words, so a " + "merged value can exceed the larger contribution by at most 2x — which bounds the " + "divergence and ties both errors to the same atol. Panels 2-4: monoprop needs small " + "sin θ (a wide rejection band), depth (to populate and compound it) and a scrambled " + "operator (so the loss does not average out) simultaneously before it separates; " + "drop any one and the gap collapses or reverses. Ribbons are min/max over seeds." +) + + +def main() -> None: + ap = argparse.ArgumentParser() + ap.add_argument("--csv", type=Path, required=True, help="the 5-engine atol sweep") + ap.add_argument( + "--divergence", + type=Path, + required=True, + help="the depth / angle / scramble sweeps", + ) + ap.add_argument("--out", type=Path, required=True) + args = ap.parse_args() + + div = load(args.divergence) + # The two Heisenberg sweeps are told apart by width, which is how they were + # run: depth at n=6 (cheap enough for 320 steps), angle at n=8. + depth = [r for r in div if r["model"] == "heisenberg" and r["qubits"] == "6"] + angle = [r for r in div if r["model"] == "heisenberg" and r["qubits"] == "8"] + scramble = [r for r in div if r["model"] == "scramble"] + for name, rows in (("depth", depth), ("angle", angle), ("scramble", scramble)): + if not rows: + raise SystemExit(f"no {name} rows in {args.divergence}") + + fig, axes = plt.subplots(2, 2, figsize=(13.0, 9.6), facecolor=SURFACE) + for ax in axes.flat: + style(ax) + + panel_atol(axes[0][0], load(args.csv)) + panel_ratio( + axes[0][1], + depth, + "steps", + "Trotter steps", + "2 · The gap compounds with depth\nHeisenberg n=6, dt=0.05 (sin θ ≈ 0.1)", + ) + panel_ratio( + axes[1][0], + angle, + "dt", + "dt (θ = 2·dt, so dt = π/8 gives sin θ = 0.707)", + "3 · …and reverses at large angles\nHeisenberg n=8, 10 steps", + ) + panel_scramble(axes[1][1], scramble) + + fig.text( + 0.008, + 0.012, + textwrap.fill(CAPTION, 168), + fontsize=7.5, + color=TEXT_SECONDARY, + va="bottom", + ) + fig.tight_layout(rect=(0, 0.072, 1, 1)) + args.out.parent.mkdir(parents=True, exist_ok=True) + fig.savefig(args.out, dpi=200, facecolor=SURFACE) + print(f"wrote {args.out}") + + +if __name__ == "__main__": + main() diff --git a/benchmarks/cross-library/xbench_accuracy.py b/benchmarks/cross-library/xbench_accuracy.py index 0a1e027b5..4132936e2 100644 --- a/benchmarks/cross-library/xbench_accuracy.py +++ b/benchmarks/cross-library/xbench_accuracy.py @@ -28,6 +28,7 @@ import argparse import csv +import itertools import math import sys from pathlib import Path @@ -44,6 +45,7 @@ "atol", "ref_atol", "ref_terms", + "ref_norm", "terms", "l2_err", "l1_err", @@ -57,6 +59,12 @@ ] +def norm(terms: dict[str, float]) -> float: + """`‖c‖₂`. Unitary conjugation preserves it, so it is the scale the absolute + errors should be read against.""" + return math.sqrt(sum(v * v for v in terms.values())) + + def measure(approx: dict[str, float], ref: dict[str, float], atol: float) -> dict: """Compare one propagated support against the converged one.""" words = set(approx) | set(ref) @@ -96,8 +104,12 @@ def main() -> None: ap = argparse.ArgumentParser(description=__doc__) ap.add_argument("--models", default="tfim,heisenberg") ap.add_argument("--qubits", type=int, default=8) - ap.add_argument("--steps", type=int, default=10) - ap.add_argument("--dt", default="0.1") + ap.add_argument( + "--steps", + default="10", + help="comma-separated depths; the reference is recomputed per depth", + ) + ap.add_argument("--dt", default="0.1", help="comma-separated angles scales") ap.add_argument("--atols", default="1e-3,1e-4,1e-5,1e-6,1e-7") ap.add_argument( "--seeds", @@ -124,55 +136,54 @@ def main() -> None: if BASELINE not in libs: raise SystemExit(f"{BASELINE} is the reference and cannot be skipped") - base = { - "QUBITS": str(args.qubits), - "STEPS": str(args.steps), - "DT": args.dt, - "ITERS": "1", - "JCOUP": "1.0", - "HFIELD": "1.0", - } + base = {"QUBITS": str(args.qubits), "ITERS": "1", "JCOUP": "1.0", "HFIELD": "1.0"} + grid = itertools.product( + args.models.split(","), + args.steps.split(","), + args.dt.split(","), + args.seeds.split(","), + ) rows = [] - for model in args.models.split(","): - for seed in args.seeds.split(","): - env = {**base, "MODEL": model, "SEED": seed} - ref, ref_obs = run(BASELINE, {**env, "ATOL": args.ref_atol}, args.quiet) - print( - f"{model} n={args.qubits} seed={seed}: reference is {BASELINE} at" - f" atol={args.ref_atol} — {len(ref)} terms, observable {ref_obs:.12g}", - file=sys.stderr, - ) - for atol in args.atols.split(","): - for lib in libs: - approx, obs = run(lib, {**env, "ATOL": atol}, args.quiet) - row = { + for model, steps, dt, seed in grid: + env = {**base, "MODEL": model, "STEPS": steps, "DT": dt, "SEED": seed} + # The reference depends on every axis, so it is recomputed per cell. + ref, ref_obs = run(BASELINE, {**env, "ATOL": args.ref_atol}, args.quiet) + print( + f"{model} n={args.qubits} steps={steps} dt={dt} seed={seed}: reference is" + f" {BASELINE} at atol={args.ref_atol} — {len(ref)} terms," + f" observable {ref_obs:.12g}", + file=sys.stderr, + ) + for atol in args.atols.split(","): + for lib in libs: + approx, obs = run(lib, {**env, "ATOL": atol}, args.quiet) + stats = measure(approx, ref, float(atol)) + rows.append( + { "model": model, "library": lib, "qubits": args.qubits, - "steps": args.steps, - "dt": args.dt, + "steps": steps, + "dt": dt, "seed": seed, "atol": atol, "ref_atol": args.ref_atol, "ref_terms": len(ref), + "ref_norm": f"{norm(ref):.6e}", "observable": f"{obs:.15g}", "ref_observable": f"{ref_obs:.15g}", "obs_abs_err": f"{abs(obs - ref_obs):.6e}", - } - stats = measure(approx, ref, float(atol)) - row["terms"] = stats.pop("terms") - row.update( - { + **{ k: f"{v:.6e}" if isinstance(v, float) else v for k, v in stats.items() - } - ) - rows.append(row) - print( - f" atol={atol} {lib:>20}: {row['terms']:>6} terms," - f" L2 {stats['l2_err']:.3e}, obs err {row['obs_abs_err']}", - file=sys.stderr, - ) + }, + } + ) + print( + f" atol={atol} {lib:>20}: {stats['terms']:>7} terms," + f" L2 {stats['l2_err']:.3e}, obs err {rows[-1]['obs_abs_err']}", + file=sys.stderr, + ) args.out.parent.mkdir(parents=True, exist_ok=True) with args.out.open("w", newline="") as fh: