Run in standalone mode and multithreading

Support
1

Description of problem

Hello, i can’t seem to be able to run my code on c++ standalone mode.
I get the error :
Cannot retrieve the values of state variables in standalone code before the simulation has been run.
Also i can’t seems to figure out how to do the multithreading. My code takes like 17 hours to run even with an i9.
Thank you for your help

Minimal code to reproduce problem

import numpy as np
from brian2 import units
from sklearn import datasets, model_selection
from brian2 import *
import matplotlib.pyplot as plt
prefs.codegen.target = 'cython'
set_device('cpp_standalone', build_on_run=False)

time_per_sample=350*ms
resting_time = 0.15 * units.second

N = 1000
F = 10 * Hz
gmax = 1 / 5
# Création des neurones
tau = 10 * ms
eqs_neuron = '''
dv/dt = -v/tau : 1
'''

input_group = PoissonGroup(N, rates=F)
neuron = NeuronGroup(1, model=eqs_neuron, threshold='v>1', reset='v=0', method='euler')
# Création des synapses
eqs_stdp = '''
    w : 1
    da/dt = -a / tau_a : 1 (event-driven) 
    db/dt = -b / tau_b : 1 (event-driven)

    tau_a = 20*ms : second
    tau_b = 20*ms : second
    A = 0.001 *gmax : 1
    B = -0.001 *gmax : 1
    gmax = 1./5. : 1
'''
on_pre = '''
    v_post += w
    a += A
    w = clip(w + b,0,gmax)
'''
on_post = '''
    b += B
    w = clip(w + a,0,gmax)
'''
S = Synapses(input_group, neuron, model=eqs_stdp, on_pre=on_pre, on_post=on_post, method='euler')
S.connect()
S.w = 'rand() * gmax'
mon = StateMonitor(S, 'w', record=[0, 1])
s_mon = SpikeMonitor(input_group)
net = Network(input_group, S, neuron,s_mon,mon)

number_of_epochs = 1
for i in range(number_of_epochs):
    print('Starting iteration %i' % i)
    for  i in range(10):
        # Configurer le taux d'entrée
        input_group.rates = i*10 * units.Hz
        # Simuler le réseau
        net.run(time_per_sample)
        a=s_mon.count
        # Laisser les variables retourner à leurs valeurs de repos
        net.run(resting_time)
        print(S.w)

What you have already tried

i tried using reinit() and activate() but without success

Full traceback of error (if relevant)

Traceback (most recent call last):
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/IPython/core/interactiveshell.py", line 3444, in run_code
    exec(code_obj, self.user_global_ns, self.user_ns)
  File "<ipython-input-2-343826d46dda>", line 1, in <module>
    runfile('/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/minModel.py', wdir='/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP')
  File "/Applications/PyCharm.app/Contents/plugins/python/helpers/pydev/_pydev_bundle/pydev_umd.py", line 198, in runfile
    pydev_imports.execfile(filename, global_vars, local_vars)  # execute the script
  File "/Applications/PyCharm.app/Contents/plugins/python/helpers/pydev/_pydev_imps/_pydev_execfile.py", line 18, in execfile
    exec(compile(contents+"\n", file, 'exec'), glob, loc)
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/minModel.py", line 66, in <module>
    print(S.w)
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/core/variables.py", line 1398, in __repr__
    values = repr(self[:])
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/core/variables.py", line 824, in __getitem__
    return self.get_item(item, level=1)
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/core/variables.py", line 816, in get_item
    values = self.get_with_index_array(item)
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/core/base.py", line 278, in device_override_decorated_function
    return func(*args, **kwds)
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/core/variables.py", line 1103, in get_with_index_array
    return variable.get_value()[indices]
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/core/variables.py", line 464, in get_value
    return self.device.get_value(self)
  File "/Users/carlos/Documents/LeTaffe/Neurosciences/projet2/STDP/venv/lib/python3.9/site-packages/brian2/devices/cpp_standalone/device.py", line 513, in get_value
    raise NotImplementedError('Cannot retrieve the values of state '
NotImplementedError: Cannot retrieve the values of state variables in standalone code before the simulation has been run.
WARNING    Active device does not have an attribute 'shape', ignoring this [brian2.devices.device]
WARNING    Active device does not have an attribute 'shape', ignoring this [brian2.devices.device]
WARNING    Active device does not have an attribute 'shape', ignoring this [brian2.devices.device]
WARNING    Active device does not have an attribute 'shape', ignoring this [brian2.devices.device]
2

Hi. You are getting the error because you are never actually building/running the model (with device.build()), but asking for the values of S.w (which are random, and only determined after the C++ code ran). Also, you cannot ask for the s_mon.count in between runs. The standalone mode will generate a single C++ binary for your complete simulation, i.e. for all runs together. This will be started as soon as you call device.build, no Python code can be executed in between.

Regarding multi-threading, this can only be activated in the standalone mode, but will also not speed up your simulation much in this example.

You can make things much faster, then your current approach, though. You have two options, either you stay with the standard “runtime” mode (somewhat easier to think about things since you can run Python code in between runs, etc.), or you change your code to run with C++ standalone mode.

For runtime mode:

  • Your current network has 1000 inputs and only a single neuron, this is quite inefficient, since (in runtime mode) there’s a fixed overhead for running each NeuronGroup/Synapses etc., regardless of its size. It would therefore be more efficient to e.g. have 10*1000 in inputs and 10 neurons, with rates and connections set up to make the equivalent calculations as your sequential code. In practice, this network would run almost at the same speed as your smaller network, making things almost 10 times faster.
  • In runtime mode, you can make use of the store()/restore() functionality to run several trials from the same starting point, which would mean that you’d no longer have to simulate things for the “resting time”.

For standalone mode, making the network bigger as described above would also be helpful (especially for multi-threading), but is not as important as for runtime mode, since there’s no Python overhead for running objects. However, you’d have to ideally turn your whole simulation into a single run where you change things not on the Python level (input_group.rates = ...), but instead as part of the simulation description. E.g. for the input rates, you can use something like input_group.run_regularly('rates += 10*Hz', dt=500*ms) to increase the input rates every 500ms.

Hope that gets you going, best
Marcel

4

Yes, you can use TimedArrays:

input_rates = TimedArray(array_of_rate, dt=500*ms)

which you can either use directly in the definition of the PoissonGroup:

input_group = PoissonGroup(N, rates='input_rates(t, i)*Hz')

or in a run_regularly operation:

input_group.run_regularly('rates = input_rates(t, i)*Hz', dt=500*ms)

The latter is slightly more efficient, since it only looks up rates in the array once every 500ms instead of every time step. I doubt it makes any measurable difference, though.
Note that if you put it directly in the PoissonGroup call, you can also have it be quiet for the last 150ms of every 500ms period:

input_group = PoissonGroup(N, rates='input_rates(t, i)*int(t % 500*ms <= 350*ms)*Hz')
5

Thank you it works great !
I’d like to know if there is a way to record the number of spikes of each neurons using the santdalone mode and this method ?

6

Not 100% sure I understand. You can of course use a SpikeMonitor to record spikes in standalone mode, but I guess you need something more? Like using this information from within the model, e.g. in the equations? If that is the case, an easy solution is to put a new variable into your NeuronGroup equations spike_count : integer. You can update this variable as part of the reset with spike_count += 1.
This can also be useful if you are interested in the total spike count for each neuron at regular intervals, say every 500ms. In that case, you’d use a monitor that records these values but only every 500ms: count_mon = StateMonitor(group, 'spike_count', record=True). Note that you probably want to call count_mon.record_single_timestep() after the run(...) call, otherwise you will not record the values from the last block.
The easiest way to get the per-block spike counts later would be to simply use np.diff(count_mon, axis=1), but you could also use a run_regularly operation to reset it to 0 after every block (make sure to have a look at the scheduling to avoid that it gets reset to 0 just before you record it!).

7

Apologies, I just came across my answer again and I realized that the StateMonitor I suggested earlier would actually record at every time step. What I meant was:

count_mon = StateMonitor(group, 'spike_count', record=True, dt=500*ms).