diff --git a/.github/workflows/basic.yml b/.github/workflows/basic.yml
index 6ff2dee426faa7f027b7489b316a444ae7e12aff..dcfcd33ec5e48298cd78dfc66fdfbdbf2272325b 100644
--- a/.github/workflows/basic.yml
+++ b/.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
diff --git a/arbor/include/arbor/lif_cell.hpp b/arbor/include/arbor/lif_cell.hpp
index e2af157ab4b84184ae213780a0dacfddb4190ae3..68d694caa9b41f52c4140122f79cb4eb01663302 100644
--- a/arbor/include/arbor/lif_cell.hpp
+++ b/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].
diff --git a/arbor/include/arbor/schedule.hpp b/arbor/include/arbor/schedule.hpp
index 7442c70e6fcb4290c3d9f32309b7acfa7fdf5ba5..6e851335be0cf50eda91e38bbb50617edac08410 100644
--- a/arbor/include/arbor/schedule.hpp
+++ b/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_.
diff --git a/doc/python/decor.rst b/doc/python/decor.rst
index 93ef1be705d3f3982bd3d301d1d73648afcecd50..099bfebdc2a5ff9e361e417b20d5da7da455598b 100644
--- a/doc/python/decor.rst
+++ b/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.
diff --git a/doc/python/recipe.rst b/doc/python/recipe.rst
index 5d74801e497557c48391a7b757d299d48e0bdd75..0a945d5f80f677e321cb236f292a0c29bdf0a472 100644
--- a/doc/python/recipe.rst
+++ b/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
diff --git a/doc/python/single_cell_model.rst b/doc/python/single_cell_model.rst
index fb030c80d92b84142ccdf055b8b2c5fabc2d38ae..2b395d4b89041fb2dd126c14c4a7221ff13bd5f0 100644
--- a/doc/python/single_cell_model.rst
+++ b/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()
diff --git a/doc/tutorial/network_ring.rst b/doc/tutorial/network_ring.rst
index b42a45dde51283df6cf9c7cfa51fb1666f73b1b4..a7cbe4be120ec3000958715dc2d411717f31cdd7 100644
--- a/doc/tutorial/network_ring.rst
+++ b/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.
diff --git a/doc/tutorial/single_cell_detailed.rst b/doc/tutorial/single_cell_detailed.rst
index 4934b388ce27141bd54e4555e5f9c19a76aabaab..4655cd1b6ee4a7c846683f08ab65ebd671e15257 100644
--- a/doc/tutorial/single_cell_detailed.rst
+++ b/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
^^^^^^^^^^^^^^
diff --git a/doc/tutorial/single_cell_detailed_recipe.rst b/doc/tutorial/single_cell_detailed_recipe.rst
index 94723c416dc6d29a020045e3098d2f57dae501f9..8f5b7791755a5fe914a3967d39bf26a41f7df975 100644
--- a/doc/tutorial/single_cell_detailed_recipe.rst
+++ b/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
diff --git a/doc/tutorial/single_cell_model.rst b/doc/tutorial/single_cell_model.rst
index 4c838c444323b3d57335d915d34c394660240cf1..f9668712e4ff4c5ddf921ed5a5e5d683cc493d07 100644
--- a/doc/tutorial/single_cell_model.rst
+++ b/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.
diff --git a/example/brunel/brunel.cpp b/example/brunel/brunel.cpp
index 3e595d57605c2cea347d670d2cfe46ecb1c662dc..7f208844d08b79dcbe5346a32245727de5483739 100644
--- a/example/brunel/brunel.cpp
+++ b/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;
}
diff --git a/example/brunel/readme.md b/example/brunel/readme.md
index 8168b6fe4385e5846b76412819b31fc8c5258b99..e09e5daaaa8474cc5a498f4354d26573c2346de8 100644
--- a/example/brunel/readme.md
+++ b/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
```
diff --git a/python/example/brunel.py b/python/example/brunel.py
new file mode 100755
index 0000000000000000000000000000000000000000..d1081fc347778c3a0fb4140585aae49857d24c71
--- /dev/null
+++ b/python/example/brunel.py
@@ -0,0 +1,156 @@
+#!/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")
diff --git a/python/schedule.cpp b/python/schedule.cpp
index e4f99106065084bd57f642095c8c88dd7b3dd4f5..36edc380836be239cb7400ba87a34e17e5580a92 100644
--- a/python/schedule.cpp
+++ b/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,
diff --git a/python/schedule.hpp b/python/schedule.hpp
index 7762d2f9d0b14187c3724fbbc7bc224c883e0cfa..9ee81fdabf2a0c38d194117878c94ac7218a36cb 100644
--- a/python/schedule.hpp
+++ b/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);
diff --git a/python/single_cell_model.cpp b/python/single_cell_model.cpp
index fae9d364662b02d6736229adafe6420aff79342b..6fd99c824edb951f7e2f6188f03dfd7cb2033f9c 100644
--- a/python/single_cell_model.cpp
+++ b/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().")
diff --git a/python/test/unit/test_schedules.py b/python/test/unit/test_schedules.py
index 9d0344783b5875382855ebea5e554cb8796294b9..e0c6d9daa2cf585484cb9af2f512b669c22e1336 100644
--- a/python/test/unit/test_schedules.py
+++ b/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():