Port of Brunel example (#1404)

* Add Python version of brunel.cpp example.
* Enforce consistency in units: both C++ and Python APIs use unit `kHz` for all frequencies.
* Align Python version `poisson_schedule` signature with C++ equivalent.
* Update unit tests to suit API changes.
* Update documentation to suit, along with minor additions related to networks and interconnectivity.
* Add `brunel.py` to `basic.yml` tests

.github/workflows/basic.yml

@@ -189,6 +189,7 @@ jobs:
             python python/example/single_cell_model.py
             python python/example/single_cell_recipe.py
             python python/example/single_cell_stdp.py
+            python python/example/brunel.py -n 400 -m 100 -e 20 -p 0.1 -w 1.2 -d 1 -g 0.5 -l 5 -t 100 -s 1 -G 50 -S 123
             python python/example/single_cell_swc.py python/example/single_cell_detailed.swc
             python python/example/single_cell_detailed.py python/example/single_cell_detailed.swc
             python python/example/single_cell_detailed_recipe.py python/example/single_cell_detailed.swc

arbor/include/arbor/lif_cell.hpp

@@ -2,7 +2,7 @@
 
 namespace arb {
 
-// Model parameteres of leaky integrate and fire neuron model.
+// Model parameters of leaky integrate and fire neuron model.
 struct lif_cell {
     // Neuronal parameters.
     double tau_m = 10;    // Membrane potential decaying constant [ms].

arbor/include/arbor/schedule.hpp

@@ -23,7 +23,7 @@ inline time_event_span as_time_event_span(const std::vector<time_type>& v) {
 }
 
 // A schedule describes a sequence of time values used for sampling. Schedules
-// are queried monotonically in time: if two method calls `events(t0, t1)` 
+// are queried monotonically in time: if two method calls `events(t0, t1)`
 // and `events(t2, t3)` are made without an intervening call to `reset()`,
 // then 0 ≤ _t0_ ≤ _t1_ ≤ _t2_ ≤ _t3_.
 

doc/python/decor.rst

@@ -120,7 +120,8 @@ Cable cell decoration
 
         Apply a mechanism with a region using the name of the mechanism.
         The mechanism will use the parameter values set in the mechanism catalogue.
-        Returns a unique identifier that can be used to query the local indexes (see :gen:`index`) assigned to the placed items on the cable cell.
+        Returns a unique identifier that can be used to query the local indexes (see :gen:`index`)
+        assigned to the placed items on the cable cell.
 
         :param str region: description of the region.
         :param str mechanism: the name of the mechanism.
@@ -128,7 +129,9 @@ Cable cell decoration
     .. method:: place(locations, const arb::mechanism_desc& d)
 
         Place one instance of synapse described by ``mechanism`` to each location in ``locations``.
-        Returns a unique identifier that can be used to query the local indexes (see :gen:`index`) assigned to the placed items on the cable cell.
+        Returns a unique identifier that can be used to query the local indexes (see :gen:`index`) assigned to the
+        placed items on the cable cell. For instance: the ``index`` returned when a synapse mechanism is ``place``d,
+        can be used when creating a :py:class:`arbor.connection`
 
         :param str locations: description of the locset.
         :param str mechanism: the name of the mechanism.

doc/python/recipe.rst

@@ -194,7 +194,7 @@ Event generator and schedules
         Construct a Poisson schedule.
 
         By default returns a schedule with events starting from :attr:`tstart` = 0 ms,
-        with an expected frequency :attr:`freq` = 10 Hz and :attr:`seed` = 0.
+        with an expected frequency :attr:`freq` = 10 kHz and :attr:`seed` = 0.
 
     .. attribute:: tstart
 
@@ -202,7 +202,7 @@ Event generator and schedules
 
     .. attribute:: freq
 
-        The expected frequency [Hz].
+        The expected frequency [kHz].
 
     .. attribute:: seed
 
@@ -225,7 +225,7 @@ An example of an event generator reads as follows:
         target = arbor.cell_member(0,0)
         seed   = target.gid
         tstart = 1
-        freq   = 5
+        freq   = 0.005
         sched  = arbor.poisson_schedule(tstart, freq, seed)
 
         # construct an event generator with this schedule on target cell and weight 0.1

doc/python/single_cell_model.rst

@@ -28,7 +28,7 @@ Single cell model
 
        :param what: Name of the variable to record (currently only 'voltage').
        :param where: :class:`location` at which to sample the variable.
-       :param frequency: The frequency at which to sample [Hz].
+       :param frequency: The frequency at which to sample [kHz].
 
     .. method:: spikes()
 

doc/tutorial/network_ring.rst

@@ -190,9 +190,14 @@ This means the timestamps of the generated events will be kept in memory. Be def
 
 In addition to having the timestamps of spikes, we want to extract the voltage as a function of time.
 
-Step **(14)** sets the probes (step **10**) to measure at a certain schedule. This is sometimes described as attaching a :term:`sampler` to a :term:`probe`. :py:func:`arbor.simulation.sample` expects a :term:`probe id` and the desired schedule (here: a recording frequency of 10 kHz). Note that the probe id is a separate index from those of :term:`connection` endpoints; probe ids correspond to the index of the list produced by :py:func:`arbor.recipe.probes` on cell ``gid``.
-
-:py:func:`arbor.simulation.sample` returns a handle to the :term:`samples <sample>` that will be recorded. We store these handles for later use.
+Step **(14)** sets the probes (step **10**) to measure at a certain schedule. This is sometimes described as
+attaching a :term:`sampler` to a :term:`probe`. :py:func:`arbor.simulation.sample` expects a :term:`probe id` and the
+desired schedule (here: a recording frequency of 10 kHz). Note that the probe id is a separate index from those of
+:term:`connection` endpoints; probe ids correspond to the index of the list produced by :py:func:`arbor.recipe.
+probes` on cell ``gid``.
+
+:py:func:`arbor.simulation.sample` returns a handle to the :term:`samples <sample>` that will be recorded. We store
+these handles for later use.
 
 Step **(15)** executes the simulation for a duration of 100 ms.
 

doc/tutorial/single_cell_detailed.rst

@@ -444,9 +444,9 @@ We can indicate the location we would like to probe using labels from the :class
 .. code-block:: python
 
    # Add voltage probes on the "custom_terminal" locset
-   # which sample the voltage at 50000 Hz
+   # which sample the voltage at 50 kHz
 
-   model.probe('voltage', where='"custom_terminal"',  frequency=50000)
+   model.probe('voltage', where='"custom_terminal"', frequency=50)
 
 The simulation
 ^^^^^^^^^^^^^^

doc/tutorial/single_cell_detailed_recipe.rst

@@ -338,7 +338,7 @@ to plot the voltage registered by the probe on the "custom_terminal" locset.
    sim.record(arbor.spike_recording.all)
 
    # Instruct the simulation to sample the probe (0, 0)
-   # at a regular schedule with period = 0.02 ms (50000 Hz)
+   # at a regular schedule with period = 0.02 ms (50 kHz)
    probe_id = arbor.cell_member(0,0)
    handle = sim.sample(probe_id, arbor.regular_schedule(0.02))
 
@@ -348,7 +348,7 @@ This variable serves as a global identifier of a probe on a cell, namely the fir
 cell with gid = 0, which is id of the :ref:`only probe <tutorialsinglecellswcrecipe-probe>` we created on
 the only cell in the model.
 
-Next, we instructed the simulation to sample ``probe_id`` at a frequency of 50KHz. That function returns a
+Next, we instructed the simulation to sample ``probe_id`` at a frequency of 50 kHz. That function returns a
 ``handle`` which we will use to extract the results of the sampling after running the simulation.
 
 The execution

doc/tutorial/single_cell_model.rst

@@ -112,7 +112,7 @@ an interface for recording outputs and running the simulation.
     m = arbor.single_cell_model(cell)
 
     # (6) Attach voltage probe sampling at 10 kHz (every 0.1 ms).
-    m.probe('voltage', '"midpoint"', frequency=10000)
+    m.probe('voltage', '"midpoint"', frequency=10)
 
     # (7) Run simulation for 30 ms of simulated activity.
     m.run(tfinal=30)
@@ -123,7 +123,7 @@ with our single-compartment cell.
 Step **(6)** adds a :meth:`arbor.single_cell_model.probe`
 used to record variables from the model. Three pieces of information are
 provided: the type of quantity we want probed (voltage), the location where we want to
-probe ('"midpoint"'), and the frequency at which we want to sample (10kHz).
+probe ('"midpoint"'), and the frequency at which we want to sample (10 kHz).
 
 Step **(7)** runs the actual simulation for a duration of 30 ms.
 

example/brunel/brunel.cpp

@@ -57,7 +57,7 @@ struct cl_options {
     uint32_t seed = 42;
 
     // Parameters for spike output.
-    bool spike_file_output = false;
+    std::string spike_file_output = "";
 
     // Turn on/off profiling output for all ranks.
     bool profile_only_zero = false;
@@ -119,18 +119,12 @@ public:
         std::vector<cell_connection> connections;
         // Add incoming excitatory connections.
         for (auto i: sample_subset(gid, 0, ncells_exc_, in_degree_exc_)) {
-            cell_member_type source{cell_gid_type(i), 0};
-            cell_member_type target{gid, 0};
-            cell_connection conn(source, target, weight_exc_, delay_);
-            connections.push_back(conn);
+            connections.push_back({{cell_gid_type(i), 0}, {gid, 0}, weight_exc_, delay_});
         }
 
         // Add incoming inhibitory connections.
         for (auto i: sample_subset(gid, ncells_exc_, ncells_exc_ + ncells_inh_, in_degree_inh_)) {
-            cell_member_type source{cell_gid_type(i), 0};
-            cell_member_type target{gid, 0};
-            cell_connection conn(source, target, weight_inh_, delay_);
-            connections.push_back(conn);
+            connections.push_back({{cell_gid_type(i), 0}, {gid, 0}, weight_inh_, delay_});
         }
         return connections;
     }
@@ -148,16 +142,10 @@ public:
     }
 
     std::vector<event_generator> event_generators(cell_gid_type gid) const override {
-        std::vector<arb::event_generator> gens;
-
         std::mt19937_64 G;
         G.seed(gid + seed_);
-
         time_type t0 = 0;
-        cell_member_type target{gid, 0};
-
-        gens.emplace_back(poisson_generator(target, weight_ext_, t0, lambda_, G));
-        return gens;
+        return {poisson_generator({gid, 0}, weight_ext_, t0, lambda_, G)};
     }
 
     cell_size_type num_sources(cell_gid_type) const override {
@@ -234,8 +222,9 @@ int main(int argc, char** argv) {
         cl_options options = o.value();
 
         std::fstream spike_out;
-        if (options.spike_file_output && root) {
-            spike_out = sup::open_or_throw("./spikes.gdf", std::ios_base::out, false);
+        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);
         }
 
         meters.checkpoint("setup", context);
@@ -356,7 +345,7 @@ std::optional<cl_options> read_options(int argc, char** argv) {
                      "-w|--weight            [Weight of excitatory connections]\n"
                      "-d|--delay             [Delay of all connections]\n"
                      "-g|--rel-inh-w         [Relative strength of inhibitory synapses with respect to the excitatory ones]\n"
-                     "-l|--lambda            [Expected number of spikes from a single poisson cell per ms]\n"
+                     "-l|--lambda            [Mean firing rate from a single poisson cell (kHz)]\n"
                      "-t|--tfinal            [Length of the simulation period (ms)]\n"
                      "-s|--dt                [Simulation time step (ms)]\n"
                      "-G|--group-size        [Number of cells per cell group]\n"
@@ -382,15 +371,15 @@ std::optional<cl_options> read_options(int argc, char** argv) {
             { opt.tfinal,            "-t", "--tfinal" },
             { opt.dt,                "-s", "--dt" },
             { opt.group_size,        "-G", "--group-size" },
-            { opt.seed,              "-s", "--seed" },
-            { to::set(opt.spike_file_output), to::flag, "-f", "--write-spikes" },
-            { to::set(opt.profile_only_zero), to::flag, "-z", "--profile-rank-zero" },
+            { opt.seed,              "-S", "--seed" },
+            { opt.spike_file_output, "-f", "--write-spikes" },
+            // { to::set(opt.profile_only_zero), to::flag, "-z", "--profile-rank-zero" },
             { to::set(opt.verbose),           to::flag, "-v", "--verbose" },
             { to::action(help),               to::flag, to::exit, "-h", "--help" }
     };
 
     if (!to::run(options, argc, argv+1)) return {};
-    if (argv[1]) throw to::option_error("unrecogonized argument", argv[1]);
+    if (argv[1]) throw to::option_error("unrecognized argument", argv[1]);
 
     if (opt.group_size < 1) {
         throw std::runtime_error("minimum of one cell per group");
@@ -423,5 +412,6 @@ std::ostream& operator<<(std::ostream& o, const cl_options& options) {
     o << "  dt                                                         : " << options.dt << "\n";
     o << "  Group size                                                 : " << options.group_size << "\n";
     o << "  Seed                                                       : " << options.seed << "\n";
+    o << "  Spike file output                                          : " << options.spike_file_output << "\n";
     return o;
 }

example/brunel/readme.md

@@ -31,5 +31,5 @@ The parameters that can be passed as command-line arguments are the following:
 For example, we could run the miniapp as follows:
 
 ```
-./brunel -n 400 -m 100 -e 20 -p 0.1 -w 1.2 -d 1 -g 0.5 -l 5 -t 100 -s 1 -G 50 -S 123 -f
+./brunel -n 400 -m 100 -e 20 -p 0.1 -w 1.2 -d 1 -g 0.5 -l 5 -t 100 -s 1 -G 50 -S 123 -f spikes.txt
 ```

python/example/brunel.py

+#!/usr/bin/env python3
+
+import arbor
+import numpy,argparse
+
+'''
+A Brunel network consists of nexc excitatory LIF neurons and ninh inhibitory LIF neurons.
+Each neuron in the network receives in_degree_prop * nexc excitatory connections
+chosen randomly, in_degree_prop * ninh inhibitory connections and next (external) Poisson connections.
+All the connections have the same delay. The strenght of excitatory and Poisson connections is given by
+parameter weight, whereas the strength of inhibitory connections is rel_inh_strength * weight.
+Poisson neurons all spike independently with expected number of spikes given by parameter poiss_lambda.
+Because of the refractory period, the activity is mostly driven by Poisson neurons and
+recurrent connections have a small effect.
+
+Call with parameters, for example:
+./brunel.py -n 400 -m 100 -e 20 -p 0.1 -w 1.2 -d 1 -g 0.5 -l 5 -t 100 -s 1 -G 50 -S 123 -f spikes.txt
+
+'''
+
+# Samples m unique values in interval [start, end) - gid.
+# We exclude gid because we don't want self-loops.
+def sample_subset(gid, start, end, m):
+    gen = numpy.random.RandomState(gid+42)
+    s = set()
+    while len(s) < m:
+        val = gen.randint(low=start,high=end)
+        if val != gid:
+            s.add(val)
+    return s
+
+class brunel_recipe (arbor.recipe):
+    def __init__(self, nexc, ninh, next, in_degree_prop, weight, delay, rel_inh_strength, poiss_lambda, seed = 42):
+
+        arbor.recipe.__init__(self)
+
+        # Make sure that in_degree_prop in the interval (0, 1]
+        if not 0.0<in_degree_prop<=1.0:
+            print("The proportion of incoming connections should be in the interval (0, 1].")
+            quit()
+
+        self.ncells_exc_ = nexc
+        self.ncells_inh_ = ninh
+        self.delay_ = delay
+        self.seed_ = seed
+
+        # Set up the parameters.
+        self.weight_exc_ = weight
+        self.weight_inh_ = -rel_inh_strength * weight
+        self.weight_ext_ =  weight
+        self.in_degree_exc_ = round(in_degree_prop * nexc)
+        self.in_degree_inh_ = round(in_degree_prop * ninh)
+        # each cell receives next incoming Poisson sources with mean rate poiss_lambda, which is equivalent
+        # to a single Poisson source with mean rate next*poiss_lambda
+        self.lambda_ = next * poiss_lambda
+
+    def num_cells(self):
+        return self.ncells_exc_ + self.ncells_inh_
+
+    def cell_kind(self, gid):
+        return arbor.cell_kind.lif
+
+    def connections_on(self, gid):
+        connections=[]
+        # Add incoming excitatory connections.
+        for i in sample_subset(gid, 0, self.ncells_exc_, self.in_degree_exc_):
+            connections.append(arbor.connection((i,0), (gid,0), self.weight_exc_, self.delay_))
+        # Add incoming inhibitory connections.
+        for i in sample_subset(gid, self.ncells_exc_, self.ncells_exc_ + self.ncells_inh_, self.in_degree_inh_):
+            connections.append(arbor.connection((i,0), (gid,0), self.weight_inh_, self.delay_))
+        return connections
+
+    def cell_description(self, gid):
+        cell = arbor.lif_cell()
+        cell.tau_m = 10
+        cell.V_th = 10
+        cell.C_m = 20
+        cell.E_L = 0
+        cell.V_m = 0
+        cell.V_reset = 0
+        cell.t_ref = 2
+        return cell
+
+    def event_generators(self, gid):
+        t0 = 0
+        sched = arbor.poisson_schedule(t0, self.lambda_, gid + self.seed_)
+        return [arbor.event_generator(arbor.cell_member(gid,0), self.weight_ext_, sched)]
+
+    def num_targets(self, gid):
+        return 1
+
+    def num_sources(self, gid):
+        return 1
+
+if __name__ == "__main__":
+
+    parser = argparse.ArgumentParser(description='Brunel model miniapp.')
+    parser.add_argument('-n', '--n-excitatory', dest='nexc', type=int, default=400, help='Number of cells in the excitatory population')
+    parser.add_argument('-m', '--n-inhibitory', dest='ninh', type=int, default=100, help='Number of cells in the inhibitory population')
+    parser.add_argument('-e', '--n-external', dest='next', type=int, default=40, help='Number of incoming Poisson (external) connections per cell')
+    parser.add_argument('-p', '--in-degree-prop', dest='syn_per_cell_prop', type=float, default=0.05, help='Proportion of the connections received per cell')
+    parser.add_argument('-w', '--weight', dest='weight', type=float, default=1.2, help='Weight of excitatory connections')
+    parser.add_argument('-d', '--delay', dest='delay', type=float, default=0.1, help='Delay of all connections')
+    parser.add_argument('-g', '--rel-inh-w', dest='rel_inh_strength', type=float, default=1, help='Relative strength of inhibitory synapses with respect to the excitatory ones')
+    parser.add_argument('-l', '--lambda', dest='poiss_lambda', type=float, default=1, help='Mean firing rate from a single poisson cell (kHz)')
+    parser.add_argument('-t', '--tfinal', dest='tfinal', type=float, default=100, help='Length of the simulation period (ms)')
+    parser.add_argument('-s', '--dt', dest='dt', type=float, default=1, help='Simulation time step (ms)')
+    parser.add_argument('-G', '--group-size', dest='group_size', type=int, default=10, help='Number of cells per cell group')
+    parser.add_argument('-S', '--seed', dest='seed', type=int, default=42, help='Seed for poisson spike generators')
+    parser.add_argument('-f', '--write-spikes', dest='spike_file_output', type=str, help='Save spikes to file')
+    # parser.add_argument('-z', '--profile-rank-zero', dest='profile_only_zero', action='store_true', help='Only output profile information for rank 0')
+    parser.add_argument('-v', '--verbose', dest='verbose', action='store_true', help='Print more verbose information to stdout')
+
+    opt = parser.parse_args()
+    if opt.verbose:
+        print("Running brunel.py with the following settings:")
+        for k,v in vars(opt).items():
+            print(f"{k} = {v}")
+
+    context = arbor.context()
+    print(context)
+
+    meters = arbor.meter_manager()
+    meters.start(context)
+
+    recipe = brunel_recipe(opt.nexc, opt.ninh, opt.next, opt.syn_per_cell_prop, opt.weight, opt.delay, opt.rel_inh_strength, opt.poiss_lambda, opt.seed)
+
+    meters.checkpoint('recipe-create', context)
+
+    hint = arbor.partition_hint()
+    hint.cpu_group_size = opt.group_size
+    hints = {arbor.cell_kind.lif: hint}
+    decomp = arbor.partition_load_balance(recipe, context, hints)
+    print(decomp)
+
+    meters.checkpoint('load-balance', context)
+
+    sim = arbor.simulation(recipe, decomp, context)
+    sim.record(arbor.spike_recording.all)
+
+    meters.checkpoint('simulation-init', context)
+
+    sim.run(opt.tfinal,opt.dt)
+
+    meters.checkpoint('simulation-run', context)
+
+    # Print profiling information
+    print(f'{arbor.meter_report(meters, context)}')
+
+    # Print spike times
+    print(f"{len(sim.spikes())} spikes generated.")
+
+    if opt.spike_file_output:
+        with open(opt.spike_file_output, 'w') as the_file:
+            for meta, data in sim.spikes():
+                the_file.write(f"{meta[0]} {data:3.3f}\n")

python/schedule.cpp

@@ -26,7 +26,7 @@ std::ostream& operator<<(std::ostream& o, const explicit_schedule_shim& e) {
 
 std::ostream& operator<<(std::ostream& o, const poisson_schedule_shim& p) {
     return o << "<arbor.poisson_schedule: tstart " << p.tstart << " ms"
-             << ", freq " << p.freq << " Hz"
+             << ", freq " << p.freq << " kHz"
              << ", seed " << p.seed << ">";
 };
 
@@ -148,6 +148,12 @@ poisson_schedule_shim::poisson_schedule_shim(
     seed = s;
 }
 
+poisson_schedule_shim::poisson_schedule_shim(arb::time_type f) {
+    set_tstart(0.);
+    set_freq(f);
+    seed = 0;
+}
+
 void poisson_schedule_shim::set_tstart(arb::time_type t) {
     pyarb::assert_throw(is_nonneg()(t), "tstart must be a non-negative number");
     tstart = t;
@@ -167,8 +173,7 @@ arb::time_type poisson_schedule_shim::get_freq() const {
 }
 
 arb::schedule poisson_schedule_shim::schedule() const {
-    // convert frequency to kHz.
-    return arb::poisson_schedule(tstart, freq/1000., rng_type(seed));
+    return arb::poisson_schedule(tstart, freq, rng_type(seed));
 }
 
 std::vector<arb::time_type> poisson_schedule_shim::events(arb::time_type t0, arb::time_type t1) {
@@ -236,15 +241,19 @@ void register_schedules(py::module& m) {
 
     poisson_schedule
         .def(py::init<time_type, time_type, std::mt19937_64::result_type>(),
-            "tstart"_a = 0., "freq"_a = 10., "seed"_a = 0,
+            "tstart"_a = 0., "freq"_a, "seed"_a = 0,
             "Construct a Poisson schedule with arguments:\n"
             "  tstart: The delivery time of the first event in the sequence [ms], 0 by default.\n"
-            "  freq:   The expected frequency [Hz], 10 by default.\n"
+            "  freq:   The expected frequency [kHz].\n"
             "  seed:   The seed for the random number generator, 0 by default.")
+        .def(py::init<time_type>(),
+            "freq"_a,
+            "Construct a Poisson schedule, starting from t = 0, default seed, with:\n"
+            "  freq:   The expected frequency [kHz], 10 by default.\n")
         .def_property("tstart", &poisson_schedule_shim::get_tstart, &poisson_schedule_shim::set_tstart,
             "The delivery time of the first event in the sequence [ms].")
         .def_property("freq", &poisson_schedule_shim::get_freq, &poisson_schedule_shim::set_freq,
-            "The expected frequency [Hz].")
+            "The expected frequency [kHz].")
         .def_readwrite("seed", &poisson_schedule_shim::seed,
             "The seed for the random number generator.")
         .def("events", &poisson_schedule_shim::events,

python/schedule.hpp

@@ -78,12 +78,12 @@ struct explicit_schedule_shim: schedule_shim_base {
 struct poisson_schedule_shim: schedule_shim_base {
     using rng_type = std::mt19937_64;
 
-    arb::time_type tstart = arb::terminal_time;
-    arb::time_type freq = 10.;  // Hz
-    rng_type::result_type seed = 0;
+    arb::time_type tstart; // ms
+    arb::time_type freq; // kHz
+    rng_type::result_type seed;
 
-    poisson_schedule_shim() = default;
     poisson_schedule_shim(arb::time_type ts, arb::time_type f, rng_type::result_type s);
+    poisson_schedule_shim(arb::time_type f);
 
     void set_tstart(arb::time_type t);
     void set_freq(arb::time_type f);

python/single_cell_model.cpp

@@ -30,7 +30,7 @@ namespace pyarb {
 // Stores the location and sampling frequency for a probe in a single cell model.
 struct probe_site {
     arb::mlocation site;  // Location of sample on morphology.
-    double frequency;     // Sampling frequency [Hz].
+    double frequency;     // Sampling frequency [kHz].
 };
 
 // Stores a single trace, which can be queried and viewed by the user at the end
@@ -194,7 +194,7 @@ public:
 
             traces_.push_back({"voltage", p.site, {}, {}});
 
-            auto sched = arb::regular_schedule(1000./p.frequency);
+            auto sched = arb::regular_schedule(p.frequency);
 
             // Now attach the sampler at probe site, with sampling schedule sched, writing to voltage
             sim_->add_sampler(arb::one_probe({0,i}), sched, trace_callback(traces_[i]));
@@ -252,7 +252,7 @@ void register_single_cell(pybind11::module& m) {
             "Sample a variable on the cell.\n"
             " what:      Name of the variable to record (currently only 'voltage').\n"
             " where:     Location on cell morphology at which to sample the variable.\n"
-            " frequency: The target frequency at which to sample [Hz].")
+            " frequency: The target frequency at which to sample [kHz].")
         .def("probe",
             [](single_cell_model& m, const char* what, const arb::mlocation& where, double frequency) {
                 m.probe(what, where, frequency);},
@@ -260,7 +260,7 @@ void register_single_cell(pybind11::module& m) {
             "Sample a variable on the cell.\n"
             " what:      Name of the variable to record (currently only 'voltage').\n"
             " where:     Location on cell morphology at which to sample the variable.\n"
-            " frequency: The target frequency at which to sample [Hz].")
+            " frequency: The target frequency at which to sample [kHz].")
         .def_property_readonly("spikes",
             [](const single_cell_model& m) {
                 return m.spike_times();}, "Holds spike times [ms] after a call to run().")

python/test/unit/test_schedules.py

@@ -109,6 +109,15 @@ class ExplicitSchedule(unittest.TestCase):
             rs.events(1., -1.)
 
 class PoissonSchedule(unittest.TestCase):
+    def test_freq_poisson_schedule(self):
+        ps = arb.poisson_schedule(42.)
+        self.assertEqual(ps.freq, 42.)
+
+    def test_freq_tstart_contor_poisson_schedule(self):
+        ps = arb.poisson_schedule(freq = 5., tstart = 4.3)
+        self.assertEqual(ps.freq, 5.)
+        self.assertEqual(ps.tstart, 4.3)
+
     def test_freq_seed_contor_poisson_schedule(self):
         ps = arb.poisson_schedule(freq = 5., seed = 42)
         self.assertEqual(ps.freq, 5.)
@@ -120,18 +129,9 @@ class PoissonSchedule(unittest.TestCase):
         self.assertEqual(ps.freq, 100.)
         self.assertEqual(ps.seed, 1000)
 
-    def test_set_tstart_freq_seed_poisson_schedule(self):
-        ps = arb.poisson_schedule()
-        ps.tstart = 4.5
-        ps.freq = 5.5
-        ps.seed = 83
-        self.assertAlmostEqual(ps.tstart, 4.5)
-        self.assertAlmostEqual(ps.freq, 5.5)
-        self.assertEqual(ps.seed, 83)
-
     def test_events_poisson_schedule(self):
         expected = [17.4107, 502.074, 506.111, 597.116]
-        ps = arb.poisson_schedule(0., 10., 0)
+        ps = arb.poisson_schedule(0., 0.01, 0)
         for i in range(len(expected)):
             self.assertAlmostEqual(expected[i], ps.events(0., 600.)[i], places = 3)
         expected = [5030.22, 5045.75, 5069.84, 5091.56, 5182.17, 5367.3, 5566.73, 5642.13, 5719.85, 5796, 5808.33]
@@ -139,33 +139,39 @@ class PoissonSchedule(unittest.TestCase):
             self.assertAlmostEqual(expected[i], ps.events(5000., 6000.)[i], places = 2)
 
     def test_exceptions_poisson_schedule(self):
+        with self.assertRaises(TypeError):
+            arb.poisson_schedule()
+        with self.assertRaises(TypeError):
+            arb.poisson_schedule(tstart = 10.)
+        with self.assertRaises(TypeError):
+            arb.poisson_schedule(seed = 1432)
         with self.assertRaisesRegex(RuntimeError,
             "tstart must be a non-negative number"):
-            arb.poisson_schedule(tstart = -10.)
+            arb.poisson_schedule(freq=34., tstart = -10.)
         with self.assertRaises(TypeError):
-            arb.poisson_schedule(tstart = None)
+            arb.poisson_schedule(freq=34., tstart = None)
         with self.assertRaises(TypeError):
-            arb.poisson_schedule(tstart = 'tstart')
+            arb.poisson_schedule(freq=34., tstart = 'tstart')
         with self.assertRaisesRegex(RuntimeError,
             "frequency must be a non-negative number"):
             arb.poisson_schedule(freq = -100.)
         with self.assertRaises(TypeError):
             arb.poisson_schedule(freq = 'freq')
         with self.assertRaises(TypeError):
-            arb.poisson_schedule(seed = -1)
+            arb.poisson_schedule(freq=34., seed = -1)
         with self.assertRaises(TypeError):
-            arb.poisson_schedule(seed = 10.)
+            arb.poisson_schedule(freq=34., seed = 10.)
         with self.assertRaises(TypeError):
-            arb.poisson_schedule(seed = 'seed')
+            arb.poisson_schedule(freq=34., seed = 'seed')
         with self.assertRaises(TypeError):
-            arb.poisson_schedule(seed = None)
+            arb.poisson_schedule(freq=34., seed = None)
         with self.assertRaisesRegex(RuntimeError,
             "t0 must be a non-negative number"):
-            ps = arb.poisson_schedule()
+            ps = arb.poisson_schedule(0,0.01)
             ps.events(-1., 1.)
         with self.assertRaisesRegex(RuntimeError,
             "t1 must be a non-negative number"):
-            ps = arb.poisson_schedule()
+            ps = arb.poisson_schedule(0,0.01)
             ps.events(1., -1.)
 
 def suite():