From 8dbe2a51a29088134029fcffb72da19ff07f8f4a Mon Sep 17 00:00:00 2001 From: Aaron Miller Date: Thu, 20 Aug 2026 19:08:12 +0100 Subject: [PATCH] =?UTF-8?q?perf(evolution):=20=E2=9A=A1=20halve=20the=20fo?= =?UTF-8?q?llower-marking=20epoch=20stamp?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `detail::MatchedEpochSet` carried a `uint32_t` stamp per term for a counter that never leaves the struct -- never serialised, never exchanged, only ever compared to `cur_`. This narrows it to `uint16_t`, keeps the wrap discipline explicit, and does NOT shrink the array anywhere. Two mechanisms: - the stamp width, with the wrap branch now live (once per 65535 gates) and the `std::fill` the thing that keeps it correct -- `epoch_` retains stale stamps for row indices reused after a truncation, so deleting the fill aliases marks; - `matched_scratch_bytes` joins `MPOperatorMemoryBreakdown` inside `total_bytes()`, with the stamp array being propagator-owned, so only `MonomialPropagator::operator_memory_usage()` can fill it in. The partitioned path sums per-partition breakdowns, so each picks up its own array. NO `shrink_to_fit()` ANYWHERE, not even an uncalled method. An earlier revision called `matched_scratch_.shrink_to_fit()` from `initialize_operator_caches_()` on the premise that the array "grows only when the term count does". That is false on Hubbard: the function runs after EVERY `build_graph()` and `propagate()`, so 29 Trotter steps reach it 29 times against a term count that grows at each step. Its A/B read +0.2821 GiB peak RSS (+3.12 B/term), 6/6, non-overlapping, at the 1x16 hubbard cell; removing it read -0.24 GiB, 6/6. And the shrink cannot free resident bytes at all -- `resize(n, 0)` never writes past `n`, so the capacity it releases was never faulted in. A retained-but-uncalled method plus a comment asking the next reader not to call it is a trap; the rejected experiment lives in the PR body instead. Re-cut twice and rebased onto main: the original was authored on a base carrying the profiling instrument, then briefly on the graph-memory PR whose key-set test is not on main. The mechanism never changed across any of it. It contains no profiling code: monoprop_PROF, profile:: and Profile.h appear zero times in the diff. MEASURED, campaign pr7v2off, both arms ENABLE_PROFILE=OFF, 12 cells x 10 interleaved reps in one allocation with order flipped per cell: kernel peak RSS (`/usr/bin/time -v`, node-sum) falls in 16 of 16 memory cells, 14 at unanimous 10/10, 15 of 16 clearing Holm across the memory family at largest adjusted p = 0.043. Hubbard `propagate` 0.97-0.99x, the hubbard `build_graph` rung 0.975-0.980x, pauli 0.989-0.999x -- hubbard beating pauli is what a per-term saving must do. Timing is a null: 1 of 24 tests resolved (`gradient[pauli]`, A 1x128, N=1, 0.9836x, 10/10, adjusted p = 0.0469, on an operation this diff cannot reach), the other 23 inside 0.964-1.014x. 2 B/term is a floor on the ARRAY's saving, not on the peak-RSS delta: peak RSS is a maximum over time and the stamp array need not be at its own maximum then. Hubbard `propagate` at N=1 implies 0.54 and 1.05 B/term of node-sum delta, below the array's own floor, which says something about WHEN the peak lands and nothing about the width. SPLIT: an earlier cut of this branch also released `init_op_map` in `get_operator()`. That is now a separate change so this PR carries one story, which means the measured binaries contain a release this branch does not. It was worth 1,189 B (8357 -> 7168 B in this campaign's own ledger) against cells of 9.5-40 GiB -- order 1e-7, far below the smallest resolved effect -- so the figures stand. Re-gated, not re-measured. An earlier instrumented pair (both arms ENABLE_PROFILE=ON) agrees on the `build_graph` rung and on pauli but read hubbard `propagate` at 0.96-0.97x -- it overstated exactly the number being advertised. The uninstrumented pair above is what certifies the shipping code, and is what is quoted. Co-Authored-By: Claude Opus 5 (1M context) --- cpp/include/monoprop/MonomialPropagator.h | 6 ++- .../detail/evolution/layer_build/Common.h | 16 +++--- cpp/monoprop/detail/operator/MPOperator.h | 5 +- cpp/tests/evolution_detail_tests.cpp | 35 +++++++++++-- cpp/tests/mp_operator_tests.cpp | 51 +++++++++++++++++++ docs/content/docs/benchmarks.mdx | 3 +- src/monoprop/bindings/binder.h | 1 + 7 files changed, 105 insertions(+), 12 deletions(-) diff --git a/cpp/include/monoprop/MonomialPropagator.h b/cpp/include/monoprop/MonomialPropagator.h index f80f3e87..e5fc8ffa 100644 --- a/cpp/include/monoprop/MonomialPropagator.h +++ b/cpp/include/monoprop/MonomialPropagator.h @@ -132,7 +132,11 @@ class MonomialPropagator { if (partition_group_) { return partitioned_operator_memory_usage_(); } - return detail::estimate_memory_usage(mp_op_); + // matched_scratch_ is propagator-owned, so estimate_memory_usage() cannot see it; on a facade the + // partitioned path above already picks it up per partition. + auto breakdown = detail::estimate_memory_usage(mp_op_); + breakdown.matched_scratch_bytes = matched_scratch_.memory_bytes(); + return breakdown; } auto graph_layers() const -> size_t { return partition_group_ ? partitioned_graph_layers_() : graph_.layers(); } diff --git a/cpp/monoprop/detail/evolution/layer_build/Common.h b/cpp/monoprop/detail/evolution/layer_build/Common.h index e5dcc3b2..216eb259 100644 --- a/cpp/monoprop/detail/evolution/layer_build/Common.h +++ b/cpp/monoprop/detail/evolution/layer_build/Common.h @@ -29,22 +29,26 @@ namespace monoprop::detail { // Marks matched followers without a per-gate O(n) memset: one counter bump clears every mark. Reused // across gates. struct MatchedEpochSet { - std::vector epoch_; - uint32_t cur_ = 0; + // One 2-byte stamp per term: the counter is never serialised, never exchanged, only compared to cur_. + using Stamp = uint16_t; - // u32 wrap resets the array — once per 2^32-1 gates. + std::vector epoch_; + Stamp cur_ = 0; + + // Wraps once per 65535 gates; without the fill a stale stamp on a row reused after a truncation aliases. auto begin_gate(size_t n) -> void { - if (cur_ == std::numeric_limits::max()) { - std::fill(epoch_.begin(), epoch_.end(), 0); + if (cur_ == std::numeric_limits::max()) { + std::fill(epoch_.begin(), epoch_.end(), Stamp{0}); cur_ = 0; } ++cur_; if (epoch_.size() < n) { - epoch_.resize(n, 0); + epoch_.resize(n, Stamp{0}); } } auto mark(size_t i) -> void { epoch_[i] = cur_; } [[nodiscard]] auto is_marked(size_t i) const -> bool { return epoch_[i] == cur_; } + [[nodiscard]] auto memory_bytes() const -> size_t { return epoch_.capacity() * sizeof(Stamp); } }; // A trivial aggregate on purpose — not std::pair — so DefaultInitVector can skip the zero-fill and lower diff --git a/cpp/monoprop/detail/operator/MPOperator.h b/cpp/monoprop/detail/operator/MPOperator.h index 6d4e70da..bf66e17d 100644 --- a/cpp/monoprop/detail/operator/MPOperator.h +++ b/cpp/monoprop/detail/operator/MPOperator.h @@ -291,6 +291,8 @@ struct MPOperatorMemoryBreakdown final { size_t init_operator_bytes = 0; size_t initial_state_bytes = 0; size_t inverted_index_bytes = 0; + // The MatchedEpochSet stamp array. Propagator-owned, so 0 unless MonomialPropagator fills it in. + size_t matched_scratch_bytes = 0; // Diagnostics: breakdowns of the fields above, deliberately excluded from total_bytes() so they can // never double-count. @@ -303,7 +305,7 @@ struct MPOperatorMemoryBreakdown final { auto total_bytes() const -> size_t { return operator_terms_bytes + op_coeffs_bytes + state_coeffs_bytes + indexing_bytes + init_operator_bytes - + initial_state_bytes + inverted_index_bytes; + + initial_state_bytes + inverted_index_bytes + matched_scratch_bytes; } auto operator+=(const MPOperatorMemoryBreakdown &o) -> MPOperatorMemoryBreakdown & { @@ -314,6 +316,7 @@ struct MPOperatorMemoryBreakdown final { init_operator_bytes += o.init_operator_bytes; initial_state_bytes += o.initial_state_bytes; inverted_index_bytes += o.inverted_index_bytes; + matched_scratch_bytes += o.matched_scratch_bytes; inverted_index_dense_bytes += o.inverted_index_dense_bytes; inverted_index_sparse_bytes += o.inverted_index_sparse_bytes; inverted_index_dense_columns += o.inverted_index_dense_columns; diff --git a/cpp/tests/evolution_detail_tests.cpp b/cpp/tests/evolution_detail_tests.cpp index 39555948..84175909 100644 --- a/cpp/tests/evolution_detail_tests.cpp +++ b/cpp/tests/evolution_detail_tests.cpp @@ -58,12 +58,12 @@ BOOST_AUTO_TEST_CASE(matched_epoch_tail_grow) { BOOST_TEST(!set.is_marked(3)); } -// When the epoch counter saturates uint32_t, begin_gate zero-fills and restarts so marks stay correct. -BOOST_AUTO_TEST_CASE(matched_epoch_u32_wrap_resets) { +// Reaches the wrap by assigning cur_, which pins the branch and the counter restart but not the fill. +BOOST_AUTO_TEST_CASE(matched_epoch_stamp_wrap_resets) { MatchedEpochSet set; set.begin_gate(4); // allocate the backing array // Force the counter to the wrap boundary; a stale slot still equals the pre-wrap counter. - set.cur_ = std::numeric_limits::max(); + set.cur_ = std::numeric_limits::max(); set.mark(1); BOOST_TEST(set.is_marked(1)); @@ -74,6 +74,35 @@ BOOST_AUTO_TEST_CASE(matched_epoch_u32_wrap_resets) { BOOST_TEST(set.is_marked(2)); } +// Reaches the wrap by counting gates, with the mark at epoch 1 so a missing fill would alias onto it. +BOOST_AUTO_TEST_CASE(matched_epoch_stamp_wrap_reached_by_gate_count) { + constexpr auto kMaxStamp = std::numeric_limits::max(); + constexpr size_t kPeriod = static_cast(kMaxStamp); + + MatchedEpochSet set; + set.begin_gate(4); + BOOST_REQUIRE(set.cur_ == MatchedEpochSet::Stamp{1}); + set.mark(1); + BOOST_TEST(set.is_marked(1)); + + // One increment per gate, folded into a single assertion rather than 65534 of them. + bool one_epoch_per_gate = true; + for (size_t k = 2; k <= kPeriod; ++k) { + set.begin_gate(4); + one_epoch_per_gate = one_epoch_per_gate && (static_cast(set.cur_) == k); + } + BOOST_TEST(one_epoch_per_gate); + BOOST_TEST(set.cur_ == kMaxStamp); // boundary reached by counting, not by assignment + + // The wrap: cur_ returns to 1, the surviving mark's own stamp, so a false is_marked(1) is the fill. + set.begin_gate(4); + BOOST_TEST(set.cur_ == MatchedEpochSet::Stamp{1}); + BOOST_TEST(!set.is_marked(1)); + set.mark(2); + BOOST_TEST(set.is_marked(2)); + BOOST_TEST(!set.is_marked(1)); +} + BOOST_AUTO_TEST_CASE(cutoff_context_abs_coeff_for) { const VecD coeffs{-3.0, 2.0, 0.0}; diff --git a/cpp/tests/mp_operator_tests.cpp b/cpp/tests/mp_operator_tests.cpp index d0a4698d..7e2f1c1b 100644 --- a/cpp/tests/mp_operator_tests.cpp +++ b/cpp/tests/mp_operator_tests.cpp @@ -21,10 +21,13 @@ #include #include +#include #include #include +#include "monoprop/MonomialPropagator.h" #include "monoprop/algebra/Algebra.h" +#include "monoprop/detail/mpi/MPICompat.h" #include "monoprop/detail/operator/MPOperator.h" using namespace monoprop; @@ -272,6 +275,54 @@ BOOST_AUTO_TEST_CASE(mp_operator_estimate_memory_usage_tracks_inverted_index_pre BOOST_CHECK_GT(after.inverted_index_bytes, 0U); // present arm } +// matched_scratch_bytes is summed by total_bytes() and accumulated by operator+= for the facade's sum. +BOOST_AUTO_TEST_CASE(mp_operator_breakdown_counts_matched_scratch_in_total_and_sum) { + detail::MPOperatorMemoryBreakdown<8> acc; + acc.op_coeffs_bytes = 100; + acc.matched_scratch_bytes = 7; + BOOST_CHECK_EQUAL(acc.total_bytes(), 107U); + + detail::MPOperatorMemoryBreakdown<8> other; + other.op_coeffs_bytes = 20; + other.matched_scratch_bytes = 3; + + acc += other; + BOOST_CHECK_EQUAL(acc.matched_scratch_bytes, 10U); + BOOST_CHECK_EQUAL(acc.total_bytes(), 130U); + + // An operator on its own has no stamp array to report. + auto bare = build_indexed_op({indices_to_bitset<8>({0, 1})}); + BOOST_CHECK_EQUAL(detail::estimate_memory_usage<8>(bare).matched_scratch_bytes, 0U); +} + +// epoch_ is empty until the first begin_gate, so this must apply a gate before the bytes can be nonzero. +BOOST_AUTO_TEST_CASE(mp_operator_breakdown_matched_scratch_nonzero_after_a_gate) { + constexpr size_t kModes = 2; + OperatorDict ham; + ham[VecZ{0, 1}] = cd{0.0, 1.0}; + VecZ initial_state{0, 1}; + auto sim = MonomialPropagator(ham, + 2 * kModes, + initial_state, + std::nullopt, + MPI_COMM_SELF, + std::nullopt, + std::nullopt, + CutoffType::Length, + std::nullopt); + BOOST_CHECK_EQUAL(sim.operator_memory_usage().matched_scratch_bytes, 0U); // no gate applied yet + + const std::vector monos{{0}}; + sim.build_graph(monos, VecZ{0}, VecD{1.0}); + + const auto live = sim.operator_memory_usage(); + BOOST_CHECK_GT(live.matched_scratch_bytes, 0U); + + auto without = live; + without.matched_scratch_bytes = 0; + BOOST_CHECK_EQUAL(live.total_bytes() - without.total_bytes(), live.matched_scratch_bytes); +} + BOOST_AUTO_TEST_CASE(mp_operator_copy_constructor_clones_store_and_coeffs) { auto op = build_indexed_op({indices_to_bitset<8>({0, 1}), indices_to_bitset<8>({2, 3})}); op.initial_state = {0}; diff --git a/docs/content/docs/benchmarks.mdx b/docs/content/docs/benchmarks.mdx index 937f1f96..005279fd 100644 --- a/docs/content/docs/benchmarks.mdx +++ b/docs/content/docs/benchmarks.mdx @@ -273,7 +273,8 @@ For every engine, `results.json` collects, indexed by Trotter step: - `native_memory`: where an engine accounts for itself, its own figure (monoprop's C++ operator-memory accounting, cuPauliProp's peak device memory over the step, `PauliPropagation.jl`'s `Base.summarysize` of the Pauli sum). These are reference values - only: they measure different things and must not be compared across engines. `QuEra ppvm` + only: they measure different things and must not be compared across engines. monoprop's + own total also steps up when a container gains a field that counted it nowhere. `QuEra ppvm` and `Qiskit pauli-prop` expose no such accounting and so report none. Each record's `operator_memory_metric` names the quantity, which for the GPU depends on the allocator in use — CUDA's resettable `cudaMemPoolAttrUsedMemHigh` where available, otherwise a diff --git a/src/monoprop/bindings/binder.h b/src/monoprop/bindings/binder.h index 5728aca5..fdd27f7c 100644 --- a/src/monoprop/bindings/binder.h +++ b/src/monoprop/bindings/binder.h @@ -263,6 +263,7 @@ auto bind_monomial_propagator(nb::module_ &mod) -> void { {"init_operator_bytes", b.init_operator_bytes}, {"initial_state_bytes", b.initial_state_bytes}, {"inverted_index_bytes", b.inverted_index_bytes}, + {"matched_scratch_bytes", b.matched_scratch_bytes}, {"total_bytes", b.total_bytes()}, // Diagnostics (not part of total_bytes; see the struct). {"d_invidx_dense_bytes", b.inverted_index_dense_bytes},