diff --git a/cpp/monoprop/detail/operator/MPOperator.h b/cpp/monoprop/detail/operator/MPOperator.h index 6d4e70da..b5aae545 100644 --- a/cpp/monoprop/detail/operator/MPOperator.h +++ b/cpp/monoprop/detail/operator/MPOperator.h @@ -62,22 +62,22 @@ template 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> store = std::make_unique>(); - VecD op_coeffs = {}; + std::unique_ptr> store{std::make_unique>()}; + 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 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 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 init_op_map = {}; - VecZ initial_state = {}; + VecD state_coeffs; + MonomialMap init_op_map{}; + VecZ initial_state; // Set once at propagator construction. - Basis basis = Basis::Majorana; - mutable std::optional> inverted_index_ = std::nullopt; + Basis basis{Basis::Majorana}; + mutable std::optional> inverted_index_{std::nullopt}; MPOperator() noexcept = default; MPOperator(MPOperator &&) noexcept = default; @@ -116,8 +116,7 @@ struct MPOperator { 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 @@ struct MPOperator { return op_coeffs; } - std::vector> 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 @@ struct MPOperator { } 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(new_inds, *store); state_rows_.reserve(state_rows_.size() + paired_inds.size()); state_vals_.reserve(state_vals_.size() + paired_inds.size()); @@ -284,22 +280,24 @@ inline auto unordered_flat_map_storage_bytes(const FlatMap &map) -> size_t { template 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 @@ struct MPOperatorMemoryBreakdown final { 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; } }; @@ -334,6 +333,7 @@ inline auto estimate_memory_usage(const MPOperator &op) -> MPOperatorM + 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(); diff --git a/cpp/tests/mp_operator_tests.cpp b/cpp/tests/mp_operator_tests.cpp index d0a4698d..3e3218e3 100644 --- a/cpp/tests/mp_operator_tests.cpp +++ b/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}; diff --git a/src/monoprop/bindings/binder.h b/src/monoprop/bindings/binder.h index 9aad5003..264c4d6b 100644 --- a/src/monoprop/bindings/binder.h +++ b/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