perf(operator): ⚡ release the init_op_map buckets once its terms bind

cpp/monoprop/detail/operator/MPOperator.h

@@ -26,7 +26,7 @@
 #include <format>
 #include <print>
 
 #include "monoprop/TypeAliases.h"
 #include "monoprop/Utilities.h"
 #include "monoprop/detail/operator/InvertedIndex.h"
 #include "monoprop/detail/operator/OperatorIndex.h"
@@ -62,22 +62,22 @@
 struct MPOperator {
     // The store is non-copyable/non-movable, so it is heap-owned by unique_ptr (keeping MPOperator
     // itself cheaply movable). Always non-null.
-    std::unique_ptr<OperatorIndex<NumModes>> store = std::make_unique<OperatorIndex<NumModes>>();
-    VecD op_coeffs = {};
+    std::unique_ptr<OperatorIndex<NumModes>> store{std::make_unique<OperatorIndex<NumModes>>()};
+    VecD op_coeffs;
     // Only fully-paired terms score nonzero (see score_new_state_rows_), which on production models is
     // ~0.07% of the rows -- a dense vector here is 99.9% zeros. state_rows_ is strictly ascending: rows are
     // scored in ascending order and the set is only ever appended to.
-    std::vector<TermIndex> state_rows_ = {};
-    VecD state_vals_ = {};         // parallel to state_rows_; every entry is a unit phase (+-1), never 0
-    size_t state_scored_rows_ = 0; // rows [0, state_scored_rows_) have been scored into state_rows_/state_vals_
+    std::vector<TermIndex> state_rows_;
+    VecD state_vals_;               // parallel to state_rows_; every entry is a unit phase (+-1), never 0
+    size_t state_scored_rows_{0uz}; // rows [0, state_scored_rows_) have been scored into state_rows_/state_vals_
     // The dense state: empty in Heisenberg unless a caller asks dense_state() to cache one; in Schrödinger
     // it is the live coefficient vector evolution mutates in place.
-    VecD state_coeffs = {};
-    MonomialMap<NumModes> init_op_map = {};
-    VecZ initial_state = {};
+    VecD state_coeffs;
+    MonomialMap<NumModes> init_op_map{};
+    VecZ initial_state;
     // Set once at propagator construction.
-    Basis basis = Basis::Majorana;
-    mutable std::optional<InvertedIndex<NumModes>> inverted_index_ = std::nullopt;
+    Basis basis{Basis::Majorana};
+    mutable std::optional<InvertedIndex<NumModes>> inverted_index_{std::nullopt};
 
     MPOperator() noexcept = default;
     MPOperator(MPOperator &&) noexcept = default;
@@ -116,8 +116,7 @@
         return *inverted_index_;
     }
 
-    // Pending init_op_map terms are erased after the lookup loop: the flat_map is not iterable while
-    // mutating.
+    // erase/clear keep bucket_count(), which init_operator_bytes reports, so drained buckets must be released.
     auto get_operator() -> const VecD & {
         if (size() == op_coeffs.size()) {
             return op_coeffs;
@@ -129,18 +128,16 @@
             return op_coeffs;
         }
 
-        std::vector<Monomial<NumModes>> del;
-        for (const auto &kv : init_op_map) {
-            const auto &mono = kv.first;
-            const auto coeff = kv.second;
-            if (const auto found = store->find(mono)) {
-                op_coeffs[*found] = coeff;
-                del.push_back(mono);
+        const auto before = init_op_map.size();
+        erase_if(init_op_map, [this](const auto &kv) {
+            const auto found = store->find(kv.first);
+            if (found) {
+                op_coeffs[*found] = kv.second;
             }
-        }
-
-        for (const auto &mono : del) {
-            init_op_map.erase(mono);
+            return found.has_value();
+        });
+        if (init_op_map.size() != before) {
+            init_op_map.rehash(0);
         }
 
         return op_coeffs;
@@ -238,8 +235,7 @@
         }
 
         VecZ new_inds(size() - state_scored_rows_);
-        std::iota(new_inds.begin(), new_inds.end(), state_scored_rows_);
-
+        std::iota(new_inds.begin(), new_inds.end(), state_scored_rows_); // NOLINT(modernize-use-ranges)
         const auto paired_inds = is_fully_paired<NumModes>(new_inds, *store);
         state_rows_.reserve(state_rows_.size() + paired_inds.size());
         state_vals_.reserve(state_vals_.size() + paired_inds.size());
@@ -279,27 +275,29 @@
 
 template <typename FlatMap>
 inline auto unordered_flat_map_storage_bytes(const FlatMap &map) -> size_t {
     return sizeof(FlatMap) + map.bucket_count() * (sizeof(typename FlatMap::value_type) + sizeof(unsigned char));
 }
 
 template <size_t NumModes>
 struct MPOperatorMemoryBreakdown final {
-    size_t operator_terms_bytes = 0;
-    size_t op_coeffs_bytes = 0;
-    size_t state_coeffs_bytes = 0;
-    size_t indexing_bytes = 0;
-    size_t init_operator_bytes = 0;
-    size_t initial_state_bytes = 0;
-    size_t inverted_index_bytes = 0;
+    size_t operator_terms_bytes{0uz};
+    size_t op_coeffs_bytes{0uz};
+    size_t state_coeffs_bytes{0uz};
+    size_t indexing_bytes{0uz};
+    size_t init_operator_bytes{0uz};
+    size_t initial_state_bytes{0uz};
+    size_t inverted_index_bytes{0uz};
 
     // Diagnostics: breakdowns of the fields above, deliberately excluded from total_bytes() so they can
     // never double-count.
-    size_t inverted_index_dense_bytes = 0;  // of inverted_index_bytes: full-height bitmap columns
-    size_t inverted_index_sparse_bytes = 0; // of inverted_index_bytes: ascending set-row lists
-    size_t inverted_index_dense_columns = 0;
-    size_t operator_terms_slack_bytes = 0; // of operator_terms_bytes: unused geometric-growth capacity
+    size_t inverted_index_dense_bytes{0uz};  // of inverted_index_bytes: full-height bitmap columns
+    size_t inverted_index_sparse_bytes{0uz}; // of inverted_index_bytes: ascending set-row lists
+    size_t inverted_index_dense_columns{0uz};
+    size_t operator_terms_slack_bytes{0uz}; // of operator_terms_bytes: unused geometric-growth capacity
     // of state_coeffs_bytes: entries of the state that are not exactly 0.0
-    size_t state_coeffs_nonzero = 0;
+    size_t state_coeffs_nonzero{0uz};
+    // Live entries behind init_operator_bytes, which is bucket_count(): bytes with no entries are dead buckets.
+    size_t init_operator_entries{0uz};
 
     auto total_bytes() const -> size_t {
         return operator_terms_bytes + op_coeffs_bytes + state_coeffs_bytes + indexing_bytes + init_operator_bytes
@@ -319,6 +317,7 @@
         inverted_index_dense_columns += o.inverted_index_dense_columns;
         operator_terms_slack_bytes += o.operator_terms_slack_bytes;
         state_coeffs_nonzero += o.state_coeffs_nonzero;
+        init_operator_entries += o.init_operator_entries;
         return *this;
     }
 };
@@ -329,11 +328,12 @@
     breakdown.operator_terms_bytes = op.store->memory_bytes();
     breakdown.op_coeffs_bytes = op.op_coeffs.capacity() * sizeof(double);
     // Every representation of the state at once: the sparse scored set plus the dense vector.
     breakdown.state_coeffs_bytes = op.state_coeffs.capacity() * sizeof(double)
                                    + op.state_rows_.capacity() * sizeof(TermIndex)
                                    + op.state_vals_.capacity() * sizeof(double);
     breakdown.indexing_bytes = op.store->index_estimated_memory_bytes();
     breakdown.init_operator_bytes = unordered_flat_map_storage_bytes(op.init_op_map);
+    breakdown.init_operator_entries = op.init_op_map.size();
     breakdown.initial_state_bytes = op.initial_state.capacity() * sizeof(size_t);
     if (op.inverted_index_.has_value()) {
         breakdown.inverted_index_bytes = op.inverted_index_->memory_bytes();

cpp/tests/mp_operator_tests.cpp

@@ -173,6 +173,66 @@ BOOST_AUTO_TEST_CASE(mp_operator_get_operator_drains_present_terms_from_init_map
     BOOST_CHECK(op.get_operator() == coeffs);
 }
 
+// erase/clear leave bucket_count(), so init_operator_bytes must fall, not just the entry count.
+BOOST_AUTO_TEST_CASE(mp_operator_get_operator_releases_init_map_when_fully_bound) {
+    const auto a = indices_to_bitset<8>({0, 1});
+    const auto b = indices_to_bitset<8>({2, 3});
+    auto op = build_indexed_op({a, b});
+
+    op.init_op_map.reserve(4096); // buckets far in excess of the two live entries
+    op.init_op_map[a] = 3.0;
+    op.init_op_map[b] = 5.0;
+
+    const auto before = detail::estimate_memory_usage<8>(op);
+    BOOST_REQUIRE_EQUAL(before.init_operator_entries, 2U);
+
+    const VecD &coeffs = op.get_operator();
+    BOOST_REQUIRE_EQUAL(coeffs.size(), 2U);
+    BOOST_CHECK_EQUAL(coeffs[0], 3.0);
+    BOOST_CHECK_EQUAL(coeffs[1], 5.0);
+
+    const auto after = detail::estimate_memory_usage<8>(op);
+    BOOST_CHECK_EQUAL(after.init_operator_entries, 0U);
+    BOOST_CHECK_LT(after.init_operator_bytes, before.init_operator_bytes);
+    BOOST_CHECK_EQUAL(op.init_op_map.bucket_count(), 0U); // released, not merely shrunk to the minimum
+}
+
+// Partial bind: the pending entry survives with its value and the bucket array shrinks to the remainder.
+BOOST_AUTO_TEST_CASE(mp_operator_get_operator_shrinks_init_map_when_partially_bound) {
+    const auto a = indices_to_bitset<8>({0, 1});
+    const auto b = indices_to_bitset<8>({2, 3});
+    const auto absent = indices_to_bitset<8>({4, 5});
+    auto op = build_indexed_op({a, b});
+
+    op.init_op_map.reserve(4096);
+    op.init_op_map[a] = 3.0;
+    op.init_op_map[absent] = 9.0;
+
+    const auto before = detail::estimate_memory_usage<8>(op);
+    (void)op.get_operator();
+    const auto after = detail::estimate_memory_usage<8>(op);
+
+    BOOST_REQUIRE_EQUAL(after.init_operator_entries, 1U);
+    const auto found = op.init_op_map.find(absent);
+    BOOST_REQUIRE(found != op.init_op_map.end());
+    BOOST_CHECK_EQUAL(found->second, 9.0);
+    BOOST_CHECK(op.init_op_map.find(a) == op.init_op_map.end());
+    BOOST_CHECK_LT(after.init_operator_bytes, before.init_operator_bytes);
+}
+
+// Nothing bound: the map is left exactly as it was, buckets included.
+BOOST_AUTO_TEST_CASE(mp_operator_get_operator_keeps_init_map_when_nothing_bound) {
+    auto op = build_indexed_op({indices_to_bitset<8>({0, 1})});
+
+    const auto absent = indices_to_bitset<8>({4, 5});
+    op.init_op_map[absent] = 9.0;
+
+    (void)op.get_operator();
+    const auto after = detail::estimate_memory_usage<8>(op);
+    BOOST_CHECK_EQUAL(after.init_operator_entries, 1U);
+    BOOST_CHECK(op.init_op_map.find(absent) != op.init_op_map.end());
+}
+
 BOOST_AUTO_TEST_CASE(mp_operator_update_initial_operator_heisenberg_branches_pauli) {
     const auto present = indices_to_bitset<8>({0, 2});
     auto op = build_indexed_op({present}, Basis::Pauli); // row 0 indexed
@@ -272,6 +332,22 @@ BOOST_AUTO_TEST_CASE(mp_operator_estimate_memory_usage_tracks_inverted_index_pre
     BOOST_CHECK_GT(after.inverted_index_bytes, 0U); // present arm
 }
 
+// init_operator_entries is a count: accumulated by operator+= but never summed into total_bytes().
+BOOST_AUTO_TEST_CASE(mp_operator_breakdown_keeps_init_operator_entries_out_of_total) {
+    detail::MPOperatorMemoryBreakdown<8> acc;
+    acc.op_coeffs_bytes = 100;
+    acc.init_operator_entries = 2;
+    BOOST_CHECK_EQUAL(acc.total_bytes(), 100U);
+
+    detail::MPOperatorMemoryBreakdown<8> other;
+    other.op_coeffs_bytes = 20;
+    other.init_operator_entries = 5;
+
+    acc += other;
+    BOOST_CHECK_EQUAL(acc.init_operator_entries, 7U);
+    BOOST_CHECK_EQUAL(acc.total_bytes(), 120U);
+}
+
 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};

src/monoprop/bindings/binder.h

@@ -269,7 +269,8 @@ auto bind_monomial_propagator(nb::module_ &mod) -> void {
                                              {"d_invidx_sparse_bytes", b.inverted_index_sparse_bytes},
                                              {"d_invidx_dense_columns", b.inverted_index_dense_columns},
                                              {"d_terms_slack_bytes", b.operator_terms_slack_bytes},
-                                             {"d_state_coeffs_nonzero", b.state_coeffs_nonzero}};
+                                             {"d_state_coeffs_nonzero", b.state_coeffs_nonzero},
+                                             {"d_init_operator_entries", b.init_operator_entries}};
     });
 }
 } // namespace monoprop::bindings::detail