Optimise Memory Consumption during Setup

CMakeLists.txt

@@ -48,6 +48,11 @@ option(ARB_VECTORIZE "use explicit SIMD code in generated mechanisms" OFF)
 option(ARB_USE_HWLOC "request support for thread pinning via HWLOC" OFF)
 mark_as_advanced(ARB_USE_HWLOC)
 
+# Log all memory allocations
+
+option(ARB_LOG_MEMORY "Log every call to `new`" OFF)
+mark_as_advanced(ARB_LOG_MEMORY)
+
 # Use externally built modcc?
 
 set(ARB_MODCC "" CACHE STRING "path to external modcc NMODL compiler")

arbor/CMakeLists.txt

@@ -26,7 +26,7 @@ set(arbor_sources
     io/serialize_hex.cpp
     label_resolution.cpp
     lif_cell_group.cpp
-    mc_cell_group.cpp
+    cable_cell_group.cpp
     mechcat.cpp
     mechinfo.cpp
     memory/gpu_wrappers.cpp
@@ -98,6 +98,10 @@ if(ARB_WITH_MPI)
         communication/mpi_context.cpp)
 endif()
 
+if(ARB_LOG_MEMORY)
+    set_source_files_properties(profile/profiler.cpp COMPILE_DEFINITIONS ARB_LOG_MEMORY)
+endif()
+
 # Add special target for private include directory, for use by arbor target
 # and arbor unit tests. Private headers are also used for the other binaries
 # until the process of splitting our private and public headers is complete.

arbor/benchmark_cell_group.cpp

@@ -15,6 +15,12 @@
 
 namespace arb {
 
+cell_size_type ARB_ARBOR_API get_sources(cell_label_range& src, const benchmark_cell& c) {
+    src.add_cell();
+    src.add_label(c.source, {0, 1});
+    return 1;
+}
+
 benchmark_cell_group::benchmark_cell_group(const std::vector<cell_gid_type>& gids,
                                            const recipe& rec,
                                            cell_label_range& cg_sources,
@@ -29,14 +35,11 @@ benchmark_cell_group::benchmark_cell_group(const std::vector<cell_gid_type>& gid
 
     cells_.reserve(gids_.size());
     for (auto gid: gids_) {
-        cells_.push_back(util::any_cast<benchmark_cell>(rec.get_cell_description(gid)));
-    }
-
-    for (const auto& c: cells_) {
-        cg_sources.add_cell();
+        auto cell = util::any_cast<benchmark_cell>(rec.get_cell_description(gid));
+        cells_.push_back(cell);
         cg_targets.add_cell();
-        cg_sources.add_label(c.source, {0, 1});
-        cg_targets.add_label(c.target, {0, 1});
+        cg_targets.add_label(cell.target, {0, 1});
+        get_sources(cg_sources, cell);
     }
 
     benchmark_cell_group::reset();

arbor/benchmark_cell_group.hpp

@@ -38,5 +38,7 @@ class benchmark_cell_group: public cell_group {
     std::vector<cell_gid_type> gids_;
 };
 
+cell_size_type ARB_ARBOR_API get_sources(cell_label_range& src, const benchmark_cell& c);
+
 } // namespace arb
 

arbor/cable_cell.cpp

@@ -90,12 +90,12 @@ struct cable_cell_impl {
     // The placeable label to lid_range map
     dynamic_typed_map<constant_type<std::unordered_multimap<cell_tag_type, lid_range>>::type> labeled_lid_ranges;
 
-    cable_cell_impl(const arb::morphology& m, const label_dict& labels, const decor& decorations):
-        provider(m, labels),
-        dictionary(labels),
-        decorations(decorations)
+    cable_cell_impl(const arb::morphology m, const label_dict labels, decor decorations):
+        provider(std::move(m), labels),
+        dictionary(std::move(labels)),
+        decorations(std::move(decorations))
     {
-        init(decorations);
+        init();
     }
 
     cable_cell_impl(): cable_cell_impl({},{},{}) {}
@@ -104,7 +104,7 @@ struct cable_cell_impl {
 
     cable_cell_impl(cable_cell_impl&& other) = default;
 
-    void init(const decor&);
+    void init();
 
     template <typename T>
     mlocation_map<T>& get_location_map(const T&) {
@@ -120,12 +120,12 @@ struct cable_cell_impl {
     }
 
     template <typename Item>
-    void place(const locset& ls, const Item& item, const cell_tag_type& label) {
+    void place(const mlocation_list& ls, const Item& item, const cell_tag_type& label) {
         auto& mm = get_location_map(item);
         cell_lid_type& lid = placed_count.get<Item>();
         cell_lid_type first = lid;
 
-        for (auto l: thingify(ls, provider)) {
+        for (auto l: ls) {
             placed<Item> p{l, lid++, item};
             mm.push_back(p);
         }
@@ -164,40 +164,42 @@ struct cable_cell_impl {
         return region_map.get<init_reversal_potential>()[init.ion];
     }
 
-    void paint(const region& reg, const density& prop) {
-        this->paint(reg, scaled_mechanism<density>(prop));
+    void paint(const mextent& cables, const density& prop, std::string& region) {
+        this->paint(cables, scaled_mechanism<density>(prop), region);
     }
 
-    void paint(const region& reg, const scaled_mechanism<density>& prop) {
+    void paint(const mextent& cables, const scaled_mechanism<density>& prop, std::string& region) {
         std::unordered_map<std::string, iexpr_ptr> im;
         for (const auto& [fst, snd]: prop.scale_expr) {
             im.insert_or_assign(fst, thingify(snd, provider));
         }
 
         auto& mm = get_region_map(prop.t_mech);
-        const auto& cables = thingify(reg, provider);
         for (const auto& c: cables) {
             // Skip zero-length cables in extent:
             if (c.prox_pos == c.dist_pos) continue;
 
             if (!mm.insert(c, {prop.t_mech, im})) {
-                throw cable_cell_error(util::pprintf("cable {} overpaints", c));
+                throw cable_cell_error(util::pprintf("Density mechanism '{}' on region '{}' overpaints at '{}'",
+                                                     prop.t_mech.mech.name(),
+                                                     region,
+                                                     c));
             }
         }
     }
 
     template <typename TaggedMech>
-    void paint(const region& reg, const TaggedMech& prop) {
-        mextent cables = thingify(reg, provider);
+    void paint(const mextent& cables, const TaggedMech& prop, std::string& region) {
         auto& mm = get_region_map(prop);
-
         for (auto c: cables) {
             // Skip zero-length cables in extent:
             if (c.prox_pos==c.dist_pos) continue;
 
             if (!mm.insert(c, prop)) {
-                std::stringstream rg; rg << reg;
-                throw cable_cell_error(util::pprintf("Setting property '{}' on region '{}' overpaints at '{}'", show(prop), rg.str(), c));
+                throw cable_cell_error(util::pprintf("Setting property '{}' on region '{}' overpaints at '{}'",
+                                                     show(prop),
+                                                     region,
+                                                     c));
             }
         }
     }
@@ -216,20 +218,65 @@ impl_ptr make_impl(cable_cell_impl* c) {
     return impl_ptr(c, [](cable_cell_impl* p){delete p;});
 }
 
-void cable_cell_impl::init(const decor& d) {
-    for (const auto& p: d.paintings()) {
-        auto& where = p.first;
-        std::visit([this, &where] (auto&& what) {this->paint(where, what);}, p.second);
+void cable_cell_impl::init() {
+    struct rcache {
+        std::string region = "";
+        mextent cables = {};
+    };
+
+    auto rc = rcache{};
+
+    for (const auto& [where, what]: decorations.paintings()) {
+        auto region = region_to_string(where);
+        // Try to cache with a lookback of one since most models paint one
+        // region in direct succession. We also key on the stringy view of
+        // regions since in general equality on regions is undecidable.
+        if (rc.region != region) {
+            rc.region = std::move(region);
+            rc.cables = thingify(where, provider);
+        }
+        switch (what.index()) {
+            case  0: { paint(rc.cables, std::get<init_membrane_potential>(what), rc.region); break; }
+            case  1: { paint(rc.cables, std::get<axial_resistivity>(what), rc.region); break; }
+            case  2: { paint(rc.cables, std::get<temperature_K>(what), rc.region); break; }
+            case  3: { paint(rc.cables, std::get<membrane_capacitance>(what), rc.region); break; }
+            case  4: { paint(rc.cables, std::get<ion_diffusivity>(what), rc.region); break; }
+            case  5: { paint(rc.cables, std::get<init_int_concentration>(what), rc.region); break; }
+            case  6: { paint(rc.cables, std::get<init_ext_concentration>(what), rc.region); break; }
+            case  7: { paint(rc.cables, std::get<init_reversal_potential>(what), rc.region); break; }
+            case  8: { paint(rc.cables, std::get<density>(what), rc.region); break; }
+            case  9: { paint(rc.cables, std::get<voltage_process>(what), rc.region); break; }
+            case 10: { paint(rc.cables, std::get<scaled_mechanism<density>>(what), rc.region); break; }
+        }
     }
-    for (const auto& p: d.placements()) {
-        auto& where = std::get<0>(p);
-        auto& label = std::get<2>(p);
-        std::visit([this, &where, &label] (auto&& what) {return this->place(where, what, label);}, std::get<1>(p));
+
+    struct lcache {
+        std::string locset = "";
+        mlocation_list places = {};
+    };
+
+    auto lc = lcache{};
+
+    for (const auto& [where, what, label]: decorations.placements()) {
+        auto locset = locset_to_string(where);
+        // Try to cache with a lookback of one since most models place to one
+        // locset in direct succession. We also key on the stringy view of
+        // locset since in general equality on regions is undecidable.
+        if (lc.locset != locset) {
+            lc.locset = std::move(locset);
+            lc.places = thingify(where, provider);
+        }
+        switch (what.index()) {
+            case 0: { place(lc.places, std::get<i_clamp>(what), label); break; }
+            case 1: { place(lc.places, std::get<threshold_detector>(what), label); break; }
+            case 2: { place(lc.places, std::get<synapse>(what), label); break; }
+            case 3: { place(lc.places, std::get<junction>(what), label); break; }
+        }
     }
 }
 
-cable_cell::cable_cell(const arb::morphology& m, const decor& decorations, const label_dict& dictionary):
-    impl_(make_impl(new cable_cell_impl(m, dictionary, decorations)))
+cable_cell::cable_cell(arb::morphology m, decor d, label_dict l):
+    impl_(make_impl(new cable_cell_impl(std::move(m), std::move(l), std::move(d))))
 {}
 
 cable_cell::cable_cell(): impl_(make_impl(new cable_cell_impl())) {}

arbor/mc_cell_group.cpp → arbor/cable_cell_group.cpp

@@ -15,7 +15,7 @@
 #include "cell_group.hpp"
 #include "fvm_lowered_cell.hpp"
 #include "label_resolution.hpp"
-#include "mc_cell_group.hpp"
+#include "cable_cell_group.hpp"
 #include "profile/profiler_macro.hpp"
 #include "sampler_map.hpp"
 #include "util/filter.hpp"
@@ -26,11 +26,23 @@
 
 namespace arb {
 
-mc_cell_group::mc_cell_group(const std::vector<cell_gid_type>& gids,
-                             const recipe& rec,
-                             cell_label_range& cg_sources,
-                             cell_label_range& cg_targets,
-                             fvm_lowered_cell_ptr lowered):
+cell_size_type ARB_ARBOR_API get_sources(cell_label_range& src, const cable_cell& c) {
+    src.add_cell();
+    cell_size_type count = 0;
+    for (const auto& [label, range]: c.detector_ranges()) {
+        src.add_label(label, range);
+        count += range.end - range.begin;
+    }
+    return count;
+}
+
+// ARB_DEFINE_LEXICOGRAPHIC_ORDERING(arb::target_handle,(a.mech_id,a.mech_index),(b.mech_id,b.mech_index))
+// ARB_DEFINE_LEXICOGRAPHIC_ORDERING(arb::deliverable_event,(a.time,a.handle,a.weight),(b.time,b.handle,b.weight))
+cable_cell_group::cable_cell_group(const std::vector<cell_gid_type>& gids,
+                                   const recipe& rec,
+                                   cell_label_range& cg_sources,
+                                   cell_label_range& cg_targets,
+                                   fvm_lowered_cell_ptr lowered):
     gids_(gids), lowered_(std::move(lowered))
 {
     // Build lookup table for gid to local index.
@@ -68,7 +80,7 @@ mc_cell_group::mc_cell_group(const std::vector<cell_gid_type>& gids,
     spike_sources_.shrink_to_fit();
 }
 
-void mc_cell_group::reset() {
+void cable_cell_group::reset() {
     spikes_.clear();
 
     for (auto &entry: sampler_map_) {
@@ -377,7 +389,7 @@ void run_samples(
     std::visit([&](auto& x) {run_samples(x, sc, raw_times, raw_samples, sample_records, scratch); }, sc.pdata_ptr->info);
 }
 
-void mc_cell_group::advance(epoch ep, time_type dt, const event_lane_subrange& event_lanes) {
+void cable_cell_group::advance(epoch ep, time_type dt, const event_lane_subrange& event_lanes) {
     time_type tstart = lowered_->time();
 
     // Bin and collate deliverable events from event lanes.
@@ -505,9 +517,8 @@ void mc_cell_group::advance(epoch ep, time_type dt, const event_lane_subrange& e
     }
 }
 
-void mc_cell_group::add_sampler(sampler_association_handle h, cell_member_predicate probeset_ids,
-                                schedule sched, sampler_function fn)
-{
+void cable_cell_group::add_sampler(sampler_association_handle h, cell_member_predicate probeset_ids,
+                                   schedule sched, sampler_function fn) {
     std::lock_guard<std::mutex> guard(sampler_mex_);
 
     std::vector<cell_member_type> probeset =
@@ -519,17 +530,17 @@ void mc_cell_group::add_sampler(sampler_association_handle h, cell_member_predic
     }
 }
 
-void mc_cell_group::remove_sampler(sampler_association_handle h) {
+void cable_cell_group::remove_sampler(sampler_association_handle h) {
     std::lock_guard<std::mutex> guard(sampler_mex_);
     sampler_map_.erase(h);
 }
 
-void mc_cell_group::remove_all_samplers() {
+void cable_cell_group::remove_all_samplers() {
     std::lock_guard<std::mutex> guard(sampler_mex_);
     sampler_map_.clear();
 }
 
-std::vector<probe_metadata> mc_cell_group::get_probe_metadata(cell_member_type probeset_id) const {
+std::vector<probe_metadata> cable_cell_group::get_probe_metadata(cell_member_type probeset_id) const {
     // Probe associations are fixed after construction, so we do not need to grab the mutex.
 
     std::optional<probe_tag> maybe_tag = util::value_by_key(probe_map_.tag, probeset_id);

arbor/mc_cell_group.hpp → arbor/cable_cell_group.hpp

@@ -12,6 +12,7 @@
 #include <arbor/recipe.hpp>
 #include <arbor/sampling.hpp>
 #include <arbor/spike.hpp>
+#include <arbor/cable_cell.hpp>
 
 #include "cell_group.hpp"
 #include "epoch.hpp"
@@ -22,11 +23,11 @@
 
 namespace arb {
 
-class ARB_ARBOR_API mc_cell_group: public cell_group {
+class ARB_ARBOR_API cable_cell_group: public cell_group {
 public:
-    mc_cell_group() = default;
+    cable_cell_group() = default;
 
-    mc_cell_group(const std::vector<cell_gid_type>& gids,
+    cable_cell_group(const std::vector<cell_gid_type>& gids,
                   const recipe& rec,
                   cell_label_range& cg_sources,
                   cell_label_range& cg_targets,
@@ -98,4 +99,6 @@ class ARB_ARBOR_API mc_cell_group: public cell_group {
     std::vector<std::size_t> target_handle_divisions_;
 };
 
+cell_size_type ARB_ARBOR_API get_sources(cell_label_range& src, const cable_cell& c);
+
 } // namespace arb