Bug/spike recording

VERSION

@@ -1 +1 @@
-0.12.2
+0.12.3-rc0

arbor/include/arbor/simulation.hpp

@@ -74,8 +74,8 @@ struct ARB_ARBOR_API simulation {
 
     std::size_t num_spikes() const;
 
-    // Register a callback that will perform a export of the global
-    // spike vector.
+    // Register a callback that will perform a export of the global spike
+    // vector. REMOVED
     void set_global_spike_callback(spike_export_function = spike_export_function{});
 
     // Register a callback that will perform a export of the rank local

arbor/simulation.cpp

@@ -100,7 +100,6 @@ struct simulation_state {
         return tau;
     }
 
-    spike_export_function global_export_callback_;
     spike_export_function local_export_callback_;
     epoch_function epoch_callback_;
     label_resolution_map source_resolution_map_;
@@ -434,7 +433,6 @@ time_type simulation_state::run(time_type tfinal, time_type dt) {
         // Present spikes to user-supplied callbacks.
         PE(spikeio);
         if (local_export_callback_) local_export_callback_(all_local_spikes);
-        if (global_export_callback_) global_export_callback_(spikes.from_local.values());
         PL(spikeio);
 
         // Append events formed from global spikes to per-cell pending event queues.
@@ -620,7 +618,7 @@ std::size_t simulation::num_spikes() const {
 }
 
 void simulation::set_global_spike_callback(spike_export_function export_callback) {
-    impl_->global_export_callback_ = std::move(export_callback);
+    throw std::runtime_error("Attempted to call deprecated method: set_global_spike_callback");
 }
 
 void simulation::set_local_spike_callback(spike_export_function export_callback) {

doc/concepts/decor.rst

@@ -382,7 +382,7 @@ between cells. An example where we're interested in when a threshold of ``-10 mV
 
     # Placing a threshold detector might look like this.
     decor = arbor.decor()
-    decor.place('"root"', arbor.threshold_detector(-10), "detector")
+    decor.place('"root"', arbor.threshold_detector(-10*U.mV), "detector")
 
     # At this point, "detector" could be connected to another cell,
     # and it would transmit events upon the voltage crossing the threshold.

doc/cpp/simulation.rst

@@ -86,8 +86,7 @@ Class documentation
     .. cpp:type:: spike_export_function = std::function<void(const std::vector<spike>&)>
 
         User-supplied callback function used as a sink for spikes generated
-        during a simulation. See :cpp:func:`set_local_spike_callback` and
-        :cpp:func:`set_global_spike_callback`.
+        during a simulation. See :cpp:func:`set_local_spike_callback`.
 
     **Constructor:**
 
@@ -113,10 +112,9 @@ Class documentation
 
     **I/O:**
 
-    .. cpp:function:: sampler_association_handle add_sampler(\
-                        cell_member_predicate probeset_ids,\
-                        schedule sched,\
-                        sampler_function f)
+    .. cpp:function:: sampler_association_handle add_sampler(cell_member_predicate probeset_ids,\
+                                                             schedule sched,\
+                                                             sampler_function f)
 
         Note: sampler functions may be invoked from a different thread than that
         which is called :cpp:func:`run`.
@@ -145,13 +143,6 @@ Class documentation
         The total number of spikes generated since either construction or
         the last call to :cpp:func:`reset`.
 
-    .. cpp:function:: void set_global_spike_callback(spike_export_function export_callback)
-
-        Register a callback that will periodically be passed a vector with all of
-        the spikes generated over all domains (the global spike vector) since
-        the last call.
-        Will be called on the MPI rank/domain with id 0.
-
     .. cpp:function:: void set_local_spike_callback(spike_export_function export_callback)
 
         Register a callback that will periodically be passed a vector with all of

doc/python/interconnectivity.rst

@@ -58,7 +58,7 @@ Interconnectivity
             # Place 'expsyn' mechanism on "synapse_site", and a threshold detector at "root"
             decor = A.decor()
             decor.place('"synapse_site"', 'expsyn', 'syn')
-            decor.place('"root"', arbor.threshold_detector(-10), 'detector')
+            decor.place('"root"', arbor.threshold_detector(-10*U.mV), 'detector')
 
             # Implement the connections_on() function on a recipe as follows:
             def connections_on(gid):

doc/python/simulation.rst

@@ -200,6 +200,10 @@ over the local and distributed hardware resources (see :ref:`pydomdec`). Then, t
 
         Record all generated spikes from cells on all MPI ranks.
 
+        **Note** This convenience method has non-trivial performance
+         implications as Arbor will collate all spikes on all ranks. Prefer
+         ``local`` unless you are absolutely required otherwise and willing to
+         pay the cost.
 
 Recording spikes
 ----------------

example/adex/adex.cpp

@@ -96,7 +96,7 @@ int main(int argc, char** argv) {
     sim.add_sampler(arb::all_probes,
                     arb::regular_schedule(opt.dt * U::ms),
                     sampler);
-    sim.set_global_spike_callback([](const auto& spks) {
+    sim.set_local_spike_callback([](const auto& spks) {
         for (const auto& spk: spks) {
             std::cerr << spk.time << ", " << spk.source.gid << ", " << spk.source.index << '\n';
         }

example/brunel/brunel.cpp

@@ -245,14 +245,15 @@ struct brunel_recipe: public recipe {
 
 int main(int argc, char** argv) {
     bool root = true;
-
     try {
+        int rank = 0;
 #ifdef ARB_MPI_ENABLED
         arbenv::with_mpi guard(argc, argv, false);
         unsigned num_threads = arbenv::default_concurrency();
         int gpu_id = arbenv::find_private_gpu(MPI_COMM_WORLD);
         auto context = arb::make_context(arb::proc_allocation{num_threads, gpu_id}, MPI_COMM_WORLD);
-        root = arb::rank(context) == 0;
+        rank = arb::rank(context);
+        root = rank == 0;
 #else
         auto context = arb::make_context(arbenv::default_allocation());
 #endif
@@ -277,8 +278,8 @@ int main(int argc, char** argv) {
 
         std::fstream spike_out;
         auto spike_file_output = options.spike_file_output;
-        if (spike_file_output != "" && root) {
-            spike_out = sup::open_or_throw(spike_file_output, std::ios_base::out, false);
+        if (spike_file_output != "") {
+            spike_out = sup::open_or_throw(spike_file_output + "-" + std::to_string(rank) + ".gdf", std::ios_base::out, false);
         }
 
         meters.checkpoint("setup", context);
@@ -326,9 +327,9 @@ int main(int argc, char** argv) {
         // Set up spike recording.
         std::vector<arb::spike> recorded_spikes;
         if (spike_out) {
-            sim.set_global_spike_callback([&recorded_spikes](auto& spikes) {
-                    recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
-                });
+            sim.set_local_spike_callback([&recorded_spikes](auto& spikes) {
+                recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
+            });
         }
 
         meters.checkpoint("model-init", context);

example/busyring/ring.cpp

@@ -247,6 +247,7 @@ std::string get_arbor_config_str() {
 }
 
 int main(int argc, char** argv) {
+    int rank = 0;
     try {
         bool root = true;
 
@@ -260,7 +261,7 @@ int main(int argc, char** argv) {
         arbenv::with_mpi guard(argc, argv, false);
         resources.gpu_id = arbenv::find_private_gpu(MPI_COMM_WORLD);
         auto context = arb::make_context(resources, MPI_COMM_WORLD);
-        auto rank = arb::rank(context);
+        rank = arb::rank(context);
         root = rank == 0;
 #else
         resources.gpu_id = arbenv::default_gpu();
@@ -321,11 +322,10 @@ int main(int argc, char** argv) {
 
         // Set up recording of spikes to a vector on the root process.
         std::vector<arb::spike> recorded_spikes;
-        if (root && params.record_spikes) {
-            sim.set_global_spike_callback(
-                [&recorded_spikes](const std::vector<arb::spike>& spikes) {
-                    recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
-                });
+        if (params.record_spikes) {
+            sim.set_local_spike_callback([&recorded_spikes](const std::vector<arb::spike>& spikes) {
+                recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
+            });
         }
 
         meters.checkpoint("model-init", context);
@@ -341,21 +341,19 @@ int main(int argc, char** argv) {
 
         // Write spikes to file
         if (root) {
-            std::cout << "\n" << ns << " spikes generated at rate of "
-                      << params.duration/ns << " ms between spikes\n";
-            if (!recorded_spikes.empty()) {
-                std::ofstream fid(params.odir + "/" + params.name + "_spikes.gdf");
-                if (!fid.good()) {
-                    std::cerr << "Warning: unable to open file spikes.gdf for spike output\n";
-                }
-                else {
-                    char linebuf[45];
-                    for (auto spike: recorded_spikes) {
-                        auto n = std::snprintf(
-                            linebuf, sizeof(linebuf), "%u %.4f\n",
-                            unsigned{spike.source.gid}, float(spike.time));
-                        fid.write(linebuf, n);
-                    }
+            std::cout << "\n" << ns << " spikes generated at rate of " << params.duration/ns << " ms between spikes\n";
+        }
+        if (!recorded_spikes.empty()) {
+            std::ofstream fid(params.odir + "/" + params.name + "-" + std::to_string(rank) + "_spikes.gdf");
+            if (!fid.good()) {
+                std::cerr << "Warning: unable to open file spikes.gdf for spike output\n";
+            }
+            else {
+                char linebuf[45];
+                for (auto spike: recorded_spikes) {
+                    auto n = std::snprintf(linebuf, sizeof(linebuf),
+                                           "%u %.4f\n", unsigned{spike.source.gid}, float(spike.time));
+                    fid.write(linebuf, n);
                 }
             }
         }
@@ -541,4 +539,4 @@ arb::cable_cell branch_cell(arb::cell_gid_type gid, const cell_parameters& param
 
     // Make a CV between every sample in the sample tree.
     return {arb::morphology(tree), decor, {}, arb::cv_policy_every_segment()};
-}
+}

example/dryrun/dryrun.cpp

@@ -125,6 +125,7 @@ class tile_desc: public arb::tile {
 };
 
 int main(int argc, char** argv) {
+    int rank = 0;
     try {
 #ifdef ARB_MPI_ENABLED
         arbenv::with_mpi guard(argc, argv, false);
@@ -142,7 +143,8 @@ int main(int argc, char** argv) {
         else {
             ctx = arb::make_context(resources, MPI_COMM_WORLD);
             if (params.defaulted) params.num_ranks = arb::num_ranks(ctx);
-            root = arb::rank(ctx)==0;
+            rank = arb::rank(ctx);        
+            root = rank == 0;
         }
 #endif
 
@@ -182,12 +184,9 @@ int main(int argc, char** argv) {
 
         // Set up recording of spikes to a vector on the root process.
         std::vector<arb::spike> recorded_spikes;
-        if (root) {
-            sim.set_global_spike_callback(
-                    [&recorded_spikes](const std::vector<arb::spike>& spikes) {
-                        recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
-                    });
-        }
+        sim.set_local_spike_callback([&recorded_spikes](const std::vector<arb::spike>& spikes) {
+            recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
+        });
 
         meters.checkpoint("model-init", ctx);
 
@@ -201,23 +200,20 @@ int main(int argc, char** argv) {
                   << params.duration/ns << " ms between spikes\n\n";
 
         // Write spikes to file
-        if (root) {
-            std::ofstream fid("spikes.gdf");
-            if (!fid.good()) {
-                std::cerr << "Warning: unable to open file spikes.gdf for spike output\n";
-            }
-            else {
-                char linebuf[45];
-                for (auto spike: recorded_spikes) {
-                    auto n = std::snprintf(
-                            linebuf, sizeof(linebuf), "%u %.4f\n",
-                            unsigned{spike.source.gid}, float(spike.time));
-                    fid.write(linebuf, n);
-                }
+        std::ofstream fid("spikes" + std::to_string(rank) + ".gdf");
+        if (!fid.good()) {
+            std::cerr << "Warning: unable to open file spikes.gdf for spike output\n";
+        }
+        else {
+            char linebuf[45];
+            for (auto spike: recorded_spikes) {
+                auto n = std::snprintf(linebuf, sizeof(linebuf),
+                                       "%u %.4f\n", unsigned{spike.source.gid}, float(spike.time));
+                fid.write(linebuf, n);
             }
-            // Write the samples to a json file.
-            write_trace_json(voltage.at(0));
         }
+        // Write the samples to a json file.
+        write_trace_json(voltage.at(0));
 
         auto profile = arb::profile::profiler_summary();
         std::cout << profile << "\n";

example/gap_junctions/gap_junctions.cpp

@@ -133,6 +133,7 @@ class gj_recipe: public arb::recipe {
 };
 
 int main(int argc, char** argv) {
+    int rank = 0;    
     try {
         bool root = true;
 
@@ -142,8 +143,7 @@ int main(int argc, char** argv) {
         int gpu_id = arbenv::find_private_gpu(MPI_COMM_WORLD);
         auto context = arb::make_context(arb::proc_allocation{nt, gpu_id}, MPI_COMM_WORLD);
         {
-            int rank;
-            MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+            rank = arb::rank(context);
             root = rank==0;
         }
 #else
@@ -191,12 +191,9 @@ int main(int argc, char** argv) {
 
         // Set up recording of spikes to a vector on the root process.
         std::vector<arb::spike> recorded_spikes;
-        if (root) {
-            sim.set_global_spike_callback(
-                [&recorded_spikes](const std::vector<arb::spike>& spikes) {
-                    recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
-                });
-        }
+        sim.set_local_spike_callback([&recorded_spikes](const std::vector<arb::spike>& spikes) {
+            recorded_spikes.insert(recorded_spikes.end(), spikes.begin(), spikes.end());
+        });
 
         meters.checkpoint("model-init", context);
 
@@ -210,20 +207,19 @@ int main(int argc, char** argv) {
 
         // Write spikes to file
         if (root) {
-            std::cout << "\n" << ns << " spikes generated at rate of "
-                      << params.sim_duration/ns << " ms between spikes\n";
-            std::ofstream fid("spikes.gdf");
-            if (!fid.good()) {
-                std::cerr << "Warning: unable to open file spikes.gdf for spike output\n";
-            }
-            else {
-                char linebuf[45];
-                for (auto spike: recorded_spikes) {
-                    auto n = std::snprintf(
-                        linebuf, sizeof(linebuf), "%u %.4f\n",
-                        unsigned{spike.source.gid}, float(spike.time));
-                    fid.write(linebuf, n);
-                }
+            std::cout << "\n" << ns << " spikes generated at rate of " << params.sim_duration/ns << " ms between spikes\n";
+        }
+        std::ofstream fid("spikes-" + std::to_string(rank) + ".gdf");
+        if (!fid.good()) {
+            std::cerr << "Warning: unable to open file spikes.gdf for spike output\n";
+        }
+        else {
+            char linebuf[45];
+            for (auto spike: recorded_spikes) {
+                auto n = std::snprintf(linebuf, sizeof(linebuf),
+                                       "%u %.4f\n",
+                                       unsigned{spike.source.gid}, float(spike.time));
+                fid.write(linebuf, n);
             }
         }