perf: 📦 build fat binary for multiple microarchitectures

.devcontainer/devcontainer.json

@@ -10,6 +10,9 @@
   "securityOpt": [
     "seccomp=unconfined"
   ],
+  "runArgs": [
+    "--shm-size=1g"
+  ],
   "remoteUser": "vscode",
   "customizations": {
     "vscode": {

.gitignore

@@ -161,7 +161,6 @@ Thumbs.db
 
 external/upstream/_srcs/
 tests/cpp/_srcs/
-_dispatch*.py
 _constants.py
 build*/
 Testing/
@@ -178,5 +177,10 @@ benches/results/**
 # devcontainer files
 .devcontainer/devcontainer-lock.json
 
+notes/**
+
+# `just capture-baseline` / `just diff-baseline` output (tools/capture-baseline.py)
+.baseline-capture/**
+
 # Useful when running in clusters
 logs/

CMakeLists.txt

@@ -42,12 +42,6 @@ if(NOT CMAKE_BUILD_TYPE)
   set(CMAKE_BUILD_TYPE "Release")
 endif()
 
-set(
-  monoprop_MAX_NUM_MODES
-  "250"
-  CACHE STRING
-  "Maximum number of simulable Fermionic modes with Python bindings"
-)
 option(monoprop_ENABLE_MPI "Enable MPI parallelization" OFF)
 option(
   monoprop_WIDE_TERM_INDEX
@@ -57,6 +51,89 @@ option(
 
 option(monoprop_ENABLE_CXX_UNIT_TESTS "Enable C++ unit test suite" ON)
 
+# Compile the propagation kernel once per x86-64 ISA tier, each into its own shared object, and pick one
+# on first use. The one switch over how this library is built for a machine, and the shape published
+# wheels are: see cmake/compiler_flags/FatBinary.cmake and docs/content/docs/fat-binary.mdx.
+#
+# On by default wherever it can work, which is the point -- a portable build and a fast one used to be
+# alternatives here, and a wheel had to pick portable. The default is computed rather than a literal
+# because there are three ways it cannot work, and in all three a single-ISA build is not a fallback but
+# the correct shape: a non-x86-64 target has one relevant ISA, a compiler without __builtin_cpu_supports
+# cannot dispatch, and a Debug build wants one quick compile rather than five. Asking for it anyway in
+# the first two cases is a configure error, not a silent downgrade.
+if(
+  CMAKE_SYSTEM_PROCESSOR
+    MATCHES
+    "^(x86_64|amd64|AMD64)$"
+  AND
+    CMAKE_CXX_COMPILER_ID
+      MATCHES
+      "^(GNU|Clang|AppleClang)$"
+  AND
+    NOT
+      CMAKE_BUILD_TYPE
+        STREQUAL
+        "Debug"
+)
+  set(_monoprop_tiered_dso_default ON)
+else()
+  set(_monoprop_tiered_dso_default OFF)
+endif()
+option(
+  monoprop_ENABLE_TIERED_DSO
+  "Compile the propagation kernel once per x86-64 ISA tier, one shared object each, and dispatch at run time"
+  ${_monoprop_tiered_dso_default}
+)
+
+# Empty means "take the default for this ISA", resolved below once CXXFlags.cmake has defined
+# monoprop_ARCH_MARCH. A cache entry so a build can pin the crossover explicitly.
+set(
+  monoprop_SPARSE_ROW_MIN_MODES
+  ""
+  CACHE STRING
+  "Logical mode count at or above which the sparse row store is preferred (empty = ISA default)"
+)
+
+# How many storage words the scan's per-term kernel is specialized for. The gain is confined to narrow
+# registers -- it is the per-term word loops getting a compile-time trip count, which matters least
+# where the loop is longest -- so this is a cap and not a switch: 0 turns the specialization off
+# entirely and gives the unspecialized kernel at every width, Bitset::kInlineWords specializes the whole
+# inline regime. See cpp/monoprop/detail/evolution/layer_build/TermProduct.h.
+set(
+  monoprop_NARROW_KERNEL_MAX_WORDS
+  4
+  CACHE STRING
+  "Storage word counts 1..N for which the per-term scan kernel is specialized (0 disables)"
+)
+
+# The upper bound is Bitset::kInlineWords -- a specialized kernel's precondition is that every operand
+# is inline -- read out of the header rather than restated here. TermProduct.h static_asserts the same
+# bound, but only after a compile has started; this is the check that names the variable the user set.
+file(
+  STRINGS
+  "${CMAKE_CURRENT_SOURCE_DIR}/cpp/monoprop/Bitset.h"
+  monoprop_INLINE_WORDS_DECL
+  REGEX "kInlineWords[ \t]*=[ \t]*[0-9]+"
+  LIMIT_COUNT 1
+)
+string(REGEX MATCH "[0-9]+" monoprop_INLINE_WORDS "${monoprop_INLINE_WORDS_DECL}")
+if(NOT monoprop_INLINE_WORDS)
+  message(
+    FATAL_ERROR
+    "could not read Bitset::kInlineWords out of cpp/monoprop/Bitset.h; the declaration this scrapes moved"
+  )
+endif()
+
+if(
+  NOT monoprop_NARROW_KERNEL_MAX_WORDS MATCHES "^[0-9]+$"
+  OR monoprop_NARROW_KERNEL_MAX_WORDS GREATER monoprop_INLINE_WORDS
+)
+  message(
+    FATAL_ERROR
+    "monoprop_NARROW_KERNEL_MAX_WORDS must be 0..${monoprop_INLINE_WORDS} (Bitset::kInlineWords), got '${monoprop_NARROW_KERNEL_MAX_WORDS}'"
+  )
+endif()
+
 set(Python_FIND_VIRTUALENV FIRST)
 find_package(
   Python
@@ -76,12 +153,47 @@ endif()
 
 include(${PROJECT_SOURCE_DIR}/cmake/compiler_flags/CXXFlags.cmake)
 
+# The sparse/dense crossover is set by whether the target has a vector popcount, so it follows the
+# arch flags rather than the mode count alone: dense costs one pass per storage word, sparse is flat
+# in the width, and without a vector popcount the dense pass degrades sooner.
+#
+# Both values are the *first whole 32-mode storage block where sparse rows win by more than the
+# run-to-run spread*, measured end to end through the propagator (not on the Stage 3 representation
+# bench, which predates the codes algebra and the support-form query record and modelled dense as a
+# word-at-a-time scan rather than the packed position lists that actually ship). The workload is held
+# identical across widths by relabelling one operator into a wider register, so only the storage width
+# varies; sparse comes out flat in the width and dense monotone, and the crossing sits in the same
+# block at cutoffs 4, 6 and 10. See notes/monomial-storage §6e for the tables.
+#
+# Expect +-1 block of machine dependence: the crossing is where a dense word pass costs what a lane
+# read does, so cache and popcount throughput move it. Both directions cost only a few percent near
+# the crossing, which is why the first *clear* win is the threshold rather than the first win.
+#
+# A tiered build pins one value across every tier rather than following each tier's capability, even
+# though only the top tiers have a vector popcount and 768 is the measured crossing there. The two
+# backends agree on term sets and values but not on term *order*, so a per-tier threshold would make a
+# wide run's accumulation order -- and hence its energy in the last bits -- depend on which CPU the
+# wheel landed on. 256 is the value the lower tiers want and the value today's wheels already use, so
+# nothing regresses; the top tiers give up a few percent on systems at or above 256 modes.
+#
+# One consequence to expect: a tiered capture and a monoprop_ARCH_MARCH=native capture disagree on any
+# case at or above 256 storage modes, because they pick different backends there. Compare like with
+# like -- `just diff-baseline-variants` across the tiers, and `just diff-baseline` only against a
+# baseline captured from the same shape.
+if(monoprop_SPARSE_ROW_MIN_MODES STREQUAL "")
+  if(monoprop_ARCH_MARCH STREQUAL "native")
+    set(monoprop_SPARSE_ROW_MIN_MODES 768)
+  else()
+    set(monoprop_SPARSE_ROW_MIN_MODES 256)
+  endif()
+endif()
+
 # report on compiler flags in use
 message(STATUS "Configuring a ${CMAKE_BUILD_TYPE} build")
 string(TOUPPER ${CMAKE_BUILD_TYPE} _cmake_build_type_upper)
 
 message(STATUS "Compiler flags for ${CMAKE_CXX_COMPILER_ID}")
-message(STATUS "   From environment       : ${CMAKE_CXX_FLAGS}")
+message(STATUS "   From environment       : ${monoprop_CXX_FLAGS_FROM_ENV}")
 set(
   _cmake_build_type_specific_flags
   "${CMAKE_CXX_FLAGS_${_cmake_build_type_upper}}"
@@ -91,6 +203,13 @@ message(
   "   Build-type-specific    : ${_cmake_build_type_specific_flags}"
 )
 message(STATUS "   Vectorization flag     : ${ARCH_FLAG}")
+message(STATUS "   Tiered ISA build       : ${monoprop_ENABLE_TIERED_DSO}")
+if(monoprop_ENABLE_TIERED_DSO)
+  message(STATUS "   ISA tiers              : ${monoprop_FAT_TIERS}")
+  message(STATUS "   Tier tuning            : -mtune=${monoprop_FAT_MTUNE}")
+  message(STATUS "   Narrow-vector cores    : ${_monoprop_narrow_cores}")
+  message(STATUS "   Baseline ISA floor     : ${_monoprop_baseline_flags}")
+endif()
 message(
   STATUS
   "   Project defaults       : ${CMAKE_CXX${CMAKE_CXX_STANDARD}_STANDARD_COMPILE_OPTION} ${monoprop_CXX_FLAGS}"
@@ -99,8 +218,25 @@ message(STATUS "   User-appended          : ${EXTRA_CXXFLAGS}")
 
 message(STATUS "   MPI parallelization    : ${monoprop_ENABLE_MPI}")
 message(STATUS "   Wide term index        : ${monoprop_WIDE_TERM_INDEX}")
-message(STATUS "   Max simulable modes    : ${monoprop_MAX_NUM_MODES}")
 message(STATUS "   C++ unit tests         : ${monoprop_ENABLE_CXX_UNIT_TESTS}")
+message(
+  STATUS
+  "   Sparse rows from       : ${monoprop_SPARSE_ROW_MIN_MODES} modes"
+)
+message(
+  STATUS
+  "   Narrow scan kernel     : up to ${monoprop_NARROW_KERNEL_MAX_WORDS} storage words"
+)
+
+# Testing is enabled from the *top-level* list file on purpose: CTest's root is wherever
+# enable_testing() was called, so called from cpp/ it wrote no top-level CTestTestfile.cmake and every
+# documented entry point (the CMakePresets test presets, `just test-mpi`, `just test-wide`) pointed
+# ctest at a directory with no tests. ctest reports "No tests were found" and exits 0 for that, so those
+# commands were silently running nothing.
+if(monoprop_ENABLE_CXX_UNIT_TESTS)
+  enable_testing()
+  include(CTest)
+endif()
 
 include(GNUInstallDirs)
 

README.md

@@ -48,6 +48,12 @@ pip install monoprop      # or: uv add monoprop
 The prebuilt PyPI wheels are single-process (built **without** MPI). For multi-rank
 runs, or to build the C++ library and executables, build from source (see below).
 
+The `x86-64` wheels are **fat binaries**: they carry the engine compiled for four
+instruction-set levels (`x86-64`, `x86-64-v2`, `x86-64-v3`, and `x86-64-v4` with
+`avx512vpopcntdq`), and pick the best one the CPU can execute when `monoprop` is
+imported. `monoprop.__variant__` says which one loaded; `monoprop_VARIANT` pins one.
+See the [fat-binary guide](https://docs.monoprop.algorithmiq.tech/fat-binary).
+
 ## Quick example
 
 Back-propagate a Majorana observable through a one-gate circuit:
@@ -100,6 +106,14 @@ uv sync --all-extras -v
 uv sync --all-extras -v --config-settings=cmake.define.monoprop_ENABLE_MPI=ON
 ```
 
+On `x86-64` that gives you the same multi-ISA build the wheels are: the propagation kernel
+compiled once per instruction-set tier, one shared object each, with the best one for the
+CPU chosen at run time. It measures as fast as `-march=native` or slightly faster, so
+there is no portable-versus-fast switch to set — see the
+[fat-binary guide](https://docs.monoprop.algorithmiq.tech/fat-binary). `just build-native`
+builds the single-ISA alternative instead, which is quicker to compile and will not run
+off this machine.
+
 C++ unit-test build:
 
 ```bash
@@ -119,10 +133,14 @@ uv sync --all-groups --all-extras -v    # installs the workspace, incl. the benc
 uv run python -m pytest -m "not mpi"   # Python tests (serial)
 just test-mpi                          # Python + C++ tests under MPI
 just test-wide                         # Python + C++ unit tests with a 64-bit TermIndex
+just test-sparse-rows                  # Python tests with the support-form row backend forced
+just test-variants                     # Python tests once per ISA variant (fat builds)
 ```
 
-See the [testing guide](https://docs.monoprop.algorithmiq.tech/testing)
-for the with/without-MPI details and the rank matrix.
+`ctest` runs every C++ case twice, once per row backend — the second pass carries
+the `sparse-rows` label. See the
+[testing guide](https://docs.algorithmiq.fi/monoprop/docs/testing) for that, the
+with/without-MPI details, and the rank matrix.
 
 ## Repository layout
 

benches/bench_models.py

@@ -51,6 +51,7 @@ def test_model(
     benchmark,
     bench_comm,
     model_configs,
+    model_rounds,
     model,
     record_model_config,
     record_model_stats,
@@ -74,10 +75,14 @@ def run(built, n_steps):
             propagator.propagate(circuit)
         return propagator.expectation_value()
 
+    # setup() runs before every round, so each round rebuilds the model and evolves a fresh
+    # propagator -- these simulations are in place, and replaying a mutated one would time the wrong
+    # thing. record_model_stats below then describes the last round, which is what any round would
+    # produce: the term counts and memory are deterministic.
     result = benchmark.pedantic(
         barriered(run, bench_comm),
         setup=barrier_setup(bench_comm, setup),
-        rounds=1,
+        rounds=model_rounds,
         iterations=1,
     )
     assert isinstance(result, float)

benches/conftest.py

@@ -168,6 +168,16 @@ def pytest_addoption(parser: pytest.Parser) -> None:
         group.addoption(f"--{name}", type=int, default=default, help=help_text)
 
     models = parser.getgroup("monoprop-models", "monoprop fixed-model overrides")
+    # One round is enough for the memory and term-count stats, which are deterministic, but it yields
+    # no spread for the timing -- and these models are expensive enough that a single sample can sit
+    # well off the median. Raise this when a timing difference is the point of the run.
+    models.addoption(
+        "--model-rounds",
+        type=int,
+        default=1,
+        help="Rounds per fixed model; each rebuilds the model first. >1 gives a median and stddev "
+        "(default: 1).",
+    )
     for model, (config_cls, _builder, _steps) in MODELS.items():
         for field in fields(config_cls):
             models.addoption(
@@ -202,6 +212,7 @@ def _meta() -> dict[str, Any]:
         "label": os.environ.get("monoprop_BENCH_LABEL", "?"),  # noqa: SIM112
         "ranks": _size(),
         "monoprop_threads": os.environ.get("monoprop_NUM_THREADS", "default"),  # noqa: SIM112
+        "monoprop_row_store": os.environ.get("monoprop_ROW_STORE") or "auto",  # noqa: SIM112
         "cpu_count_logical": psutil.cpu_count(logical=True),
         "cpu_count_physical": psutil.cpu_count(logical=False),
         "hostname": socket.gethostname(),
@@ -212,14 +223,34 @@ def _meta() -> dict[str, Any]:
         "python_version": platform.python_version(),
         "nanobind_version": monoprop.__nanobind_version__,
         "nanobind_backend_version": nanobind_backend_version,
-        "monoprop_max_num_modes": monoprop.MAX_NUM_MODES,
         "malloc_arena_max": os.environ.get("MALLOC_ARENA_MAX", "default"),
         "omp_num_threads": os.environ.get("OMP_NUM_THREADS", "default"),
         # Filled by _record_placement: the threads exist only once a propagator does.
         "pinning": {},
     }
 
 
+def _record_row_store(propagator: Any) -> None:
+    """Fold one propagator's resolved row backend into this run's metadata.
+
+    ``monoprop_ROW_STORE`` says what was asked for, not what ran: unset lets the storage width pick,
+    and the crossover it picks against is a build-time constant. The two backends accumulate a term
+    sum in different orders and have different footprints, so a report has to name the one that ran.
+    Widths differ within a run, hence so can the backend: a disagreement records as ``"mixed"``
+    rather than letting the last propagator speak for the others.
+    """
+    if _rank() != 0:
+        return
+    # Read straight off the binding, with no getattr fallback: a benchmark whose whole job is to name
+    # the backend that ran must fail loudly against an extension that cannot say, not quietly record
+    # nothing.
+    resolved = "sparse" if propagator._simulator.rows_are_sparse else "dense"
+    seen = _RESULTS["meta"].get("row_store_effective")
+    _RESULTS["meta"]["row_store_effective"] = (
+        resolved if seen in (None, resolved) else "mixed"
+    )
+
+
 def _params(config: pytest.Config) -> dict[str, Any]:
     """Return the resolved random-problem hyperparameters (defaults included)."""
     return {
@@ -276,6 +307,12 @@ def bench_rounds(request: pytest.FixtureRequest) -> int:
     return int(request.config.getoption("--bench-rounds"))
 
 
+@pytest.fixture(scope="session")
+def model_rounds(request: pytest.FixtureRequest) -> int:
+    """Return the round count for the fixed-model benchmarks."""
+    return int(request.config.getoption("--model-rounds"))
+
+
 @pytest.fixture(scope="session")
 def model_configs(request: pytest.FixtureRequest) -> dict[str, Any]:
     """Return each fixed model's config, every field resolved from the CLI.
@@ -334,6 +371,7 @@ def _record_model_stats(
 ) -> None:
     """Record term count, operator memory breakdown and footprint under ``key``."""
     _record("opsize", key, {"terms": _reduce_sum(comm, propagator.size())})
+    _record_row_store(propagator)
 
     # Placement is only observable while the propagator's threads are alive.
     _record_placement(comm)
@@ -529,6 +567,7 @@ def built_graph(
 
     # Under MPI the operator is partitioned, so sum the partitions.
     _record("opsize", picture, {"terms": _reduce_sum(bench_comm, mp.size())})
+    _record_row_store(mp)
 
     # Settled RSS once the build's transients are released -- the persistent
     # footprint the per-operation peak cannot see.