Non-contiguous subgroups

brian2/codegen/runtime/cython_rt/templates/ratemonitor.pyx

@@ -1,31 +1,44 @@
-{# USES_VARIABLES { N, t, rate, _clock_t, _clock_dt, _spikespace,
-                    _num_source_neurons, _source_start, _source_stop } #}
+{# USES_VARIABLES { N, t, rate, _clock_t, _clock_dt, _spikespace} #}
 {% extends 'common_group.pyx' %}
 
 {% block maincode %}
-
     cdef size_t _num_spikes = {{_spikespace}}[_num{{_spikespace}}-1]
-
-    # For subgroups, we do not want to record all spikes
+
+    {% if subgroup %}
+    cdef int32_t _filtered_spikes = 0
+    cdef size_t _source_index_counter = 0
+    {% if contiguous %}{# contiguous subgroup #}
     # We assume that spikes are ordered
-    cdef int _start_idx = -1
-    cdef int _end_idx = -1
-    cdef size_t _j
+    _start_idx = _num_spikes
+    _end_idx = _num_spikes
     for _j in range(_num_spikes):
         _idx = {{_spikespace}}[_j]
         if _idx >= _source_start:
             _start_idx = _j
             break
-    if _start_idx == -1:
-        _start_idx = _num_spikes
-    for _j in range(_start_idx, _num_spikes):
+    for _j in range(_num_spikes-1, _start_idx-1, -1):
         _idx = {{_spikespace}}[_j]
-        if _idx >= _source_stop:
-            _end_idx = _j
+        if _idx < _source_stop:
+            break
+        _end_idx = _j
+    _filtered_spikes = _end_idx - _start_idx
+    {% else %}{# non-contiguous subgroup #}
+    for _j in range(_num_spikes):
+        _idx = {{_spikespace}}[_j]
+        if _idx < {{_source_indices}}[_source_index_counter]:
+            continue
+        while {{_source_indices}}[_source_index_counter] < _idx:
+            _source_index_counter += 1
+        if (_source_index_counter < {{source_N}} and
+                _idx == {{_source_indices}}[_source_index_counter]):
+            _source_index_counter += 1
+            _filtered_spikes += 1
+
+        if _source_index_counter == {{source_N}}:
             break
-    if _end_idx == -1:
-        _end_idx =_num_spikes
-    _num_spikes = _end_idx - _start_idx
+    {% endif %}
+    _num_spikes = _filtered_spikes
+    {% endif %}
 
     # Calculate the new length for the arrays
     cdef size_t _new_len = {{_dynamic_t}}.shape[0] + 1
@@ -36,6 +49,6 @@
 
     # Set the new values
     {{_dynamic_t}}.data[_new_len-1] = {{_clock_t}}
-    {{_dynamic_rate}}.data[_new_len-1] = _num_spikes/{{_clock_dt}}/_num_source_neurons
+    {{_dynamic_rate}}.data[_new_len-1] = _num_spikes/{{_clock_dt}}/{{source_N}}
 
 {% endblock %}

brian2/codegen/runtime/cython_rt/templates/spikemonitor.pyx

@@ -1,5 +1,4 @@
-{# USES_VARIABLES { N, _clock_t, count,
-                    _source_start, _source_stop} #}
+{# USES_VARIABLES { N, _clock_t, count} #}
 {% extends 'common_group.pyx' %}
 
 {% block maincode %}
@@ -9,11 +8,20 @@
 
     cdef size_t _num_events = {{_eventspace}}[_num{{_eventspace}}-1]
     cdef size_t _start_idx, _end_idx, _curlen, _newlen, _j
+    {% if subgroup and not contiguous %}
+    # We use the same data structure as for the eventspace to store the
+    # "filtered" events, i.e. the events that are indexed in the subgroup
+    cdef int[{{source_N}} + 1] _filtered_events
+    cdef size_t _source_index_counter = 0
+    _filtered_events[{{source_N}}] = 0
+    {% endif %}
     {% for varname, var in record_variables | dictsort %}
     cdef {{cpp_dtype(var.dtype)}}[:] _{{varname}}_view
     {% endfor %}
     if _num_events > 0:
+        {% if subgroup %}
         # For subgroups, we do not want to record all spikes
+        {% if contiguous %}
         # We assume that spikes are ordered
         _start_idx = _num_events
         _end_idx = _num_events
@@ -28,6 +36,25 @@
                 break
             _end_idx = _j
         _num_events = _end_idx - _start_idx
+        {% else %}
+        for _j in range(_num_events):
+            _idx = {{_eventspace}}[_j]
+            if _idx < {{_source_indices}}[_source_index_counter]:
+                continue
+            if _idx > {{_source_indices}}[{{source_N}}-1]:
+                break
+            while {{_source_indices}}[_source_index_counter] < _idx:
+                _source_index_counter += 1
+            if (_source_index_counter < {{source_N}} and
+                    _idx == {{_source_indices}}[_source_index_counter]):
+                _source_index_counter += 1
+                _filtered_events[_filtered_events[{{source_N}}]] = _idx
+                _filtered_events[{{source_N}}] += 1
+            if _source_index_counter == {{source_N}}:
+                break
+        _num_events = _filtered_events[{{source_N}}]
+        {% endif %}
+        {% endif %}
         if _num_events > 0:
             # scalar code
             _vectorisation_idx = 1
@@ -41,6 +68,8 @@
             _{{varname}}_view = {{get_array_name(var, access_data=False)}}.data
             {% endfor %}
             # Copy the values across
+            {% if subgroup %}
+            {% if contiguous %}
             for _j in range(_start_idx, _end_idx):
                 _idx = {{_eventspace}}[_j]
                 _vectorisation_idx = _idx
@@ -49,4 +78,24 @@
                 _{{varname}}_view [_curlen + _j - _start_idx] = _to_record_{{varname}}
                 {% endfor %}
                 {{count}}[_idx - _source_start] += 1
+            {% else %}
+            for _j in range(_num_events):
+                _idx = _filtered_events[_j]
+                _vectorisation_idx = _idx
+                {{ vector_code|autoindent }}
+                {% for varname in record_variables | sort %}
+                _{{varname}}_view [_curlen + _j] = _to_record_{{varname}}
+                {% endfor %}
+                {{count}}[_to_record_i] += 1
+            {% endif %}
+            {% else %}
+            for _j in range(_num_events):
+                _idx = {{_eventspace}}[_j]
+                _vectorisation_idx = _idx
+                {{ vector_code|autoindent }}
+                {% for varname in record_variables | sort %}
+                _{{varname}}_view [_curlen + _j] = _to_record_{{varname}}
+                {% endfor %}
+                {{count}}[_idx] += 1
+            {% endif %}
 {% endblock %}

brian2/codegen/runtime/cython_rt/templates/summed_variable.pyx

@@ -8,7 +8,11 @@
 
     # Set all the target variable values to zero
     for _target_idx in range({{_target_size_name}}):
+        {% if _target_contiguous %}
         {{_target_var_array}}[_target_idx + {{_target_start}}] = 0
+        {% else %}
+        {{_target_var_array}}[{{_target_indices}}[_target_idx]] = 0
+        {% endif %}
 
     # scalar code
     _vectorisation_idx = 1

brian2/codegen/runtime/cython_rt/templates/synapses_create_generator.pyx

@@ -68,8 +68,11 @@ cdef void _flush_buffer(buf, dynarr, int buf_len):
     {{scalar_code['update']|autoindent}}
 
     for _{{outer_index}} in range({{outer_index_size}}):
+        {% if outer_contiguous %}
         _raw{{outer_index_array}} = _{{outer_index}} + {{outer_index_offset}}
-
+        {% else %}
+        _raw{{outer_index_array}} = {{get_array_name(variables[outer_sub_idx])}}[_{{outer_index}}]
+        {% endif %}
         {% if not result_index_condition %}
         {{vector_code['create_cond']|autoindent}}
         if not _cond:
@@ -162,8 +165,11 @@ cdef void _flush_buffer(buf, dynarr, int buf_len):
         {% endif %}
 
             {{vector_code['generator_expr']|autoindent}}
+            {% if result_contiguous %}
             _raw{{result_index_array}} = _{{result_index}} + {{result_index_offset}}
-
+            {% else %}
+            _raw{{result_index_array}} = {{get_array_name(variables[result_sub_idx])}}[_{{result_index}}]
+            {% endif %}
             {% if result_index_condition %}
             {% if result_index_used %}
             {# The condition could index outside of array range #}

brian2/codegen/runtime/numpy_rt/numpy_rt.py

@@ -127,7 +127,7 @@ def __iter__(self):
             return iter(self.indices)
 
     # Allow conversion to a proper array with np.array(...)
-    def __array__(self, dtype=None, copy=None):
+    def __array__(self, dtype=np.int32, copy=None):
         if copy is False:
             raise ValueError("LazyArange does not support copy=False")
         if self.indices is None:

brian2/codegen/runtime/numpy_rt/templates/ratemonitor.py_

@@ -1,14 +1,19 @@
-{# USES_VARIABLES { rate, t, _spikespace, _num_source_neurons,
-                    _clock_t, _clock_dt, _source_start, _source_stop, N } #}
+{# USES_VARIABLES { rate, t, _spikespace, _clock_t, _clock_dt, N } #}
 {% extends 'common_group.py_' %}
 
 {% block maincode %}
 _spikes = {{_spikespace}}[:{{_spikespace}}[-1]]
+{% if subgroup %}
 # Take subgroups into account
+{% if contiguous %}
 _spikes = _spikes[(_spikes >= _source_start) & (_spikes < _source_stop)]
+{% else %}
+_spikes = _numpy.intersect1d(_spikes, {{_source_indices}}, assume_unique=True)
+{% endif %}
+{% endif %}
 _new_len = {{N}} + 1
 _owner.resize(_new_len)
 {{N}} = _new_len
 {{_dynamic_t}}[-1] = {{_clock_t}}
-{{_dynamic_rate}}[-1] = 1.0 * len(_spikes) / {{_clock_dt}} / _num_source_neurons
+{{_dynamic_rate}}[-1] = 1.0 * len(_spikes) / {{_clock_dt}} / {{source_N}}
 {% endblock %}

brian2/codegen/runtime/numpy_rt/templates/spikemonitor.py_

@@ -1,4 +1,4 @@
-{# USES_VARIABLES {N, count, _clock_t, _source_start, _source_stop, _source_N} #}
+{# USES_VARIABLES {N, count, _clock_t} #}
 {% extends 'common_group.py_' %}
 
 {% block maincode %}
@@ -9,11 +9,15 @@
 _n_events = {{_eventspace}}[-1]
 if _n_events > 0:
     _events = {{_eventspace}}[:_n_events]
+    {% if subgroup %}
     # Take subgroups into account
-    if _source_start != 0 or _source_stop != _source_N:
-        _events = _events[(_events >= _source_start) & (_events < _source_stop)]
-        _n_events = len(_events)
-
+    {% if contiguous %}
+    _events = _events[(_events >= _source_start) & (_events < _source_stop)]
+    {% else %}
+    _events = _numpy.intersect1d(_events, {{_source_indices}}, assume_unique=True)
+    {% endif %}
+    _n_events = len(_events)
+    {% endif %}
     if _n_events > 0:
         _vectorisation_idx = 1
         {{scalar_code|autoindent}}
@@ -28,5 +32,13 @@ if _n_events > 0:
         {% set dynamic_varname = get_array_name(var, access_data=False) %}
         {{dynamic_varname}}[_curlen:_newlen] = _to_record_{{varname}}
         {% endfor %}
+        {% if not subgroup %}
+        {{count}}[_events] += 1
+        {% else %}
+        {% if contiguous %}
         {{count}}[_events - _source_start] += 1
+        {% else %}
+        {{count}}[_to_record_i] += 1
+        {% endif %}
+        {% endif %}
 {% endblock %}

brian2/codegen/runtime/numpy_rt/templates/summed_variable.py_

@@ -17,6 +17,9 @@ _vectorisation_idx = LazyArange(N)
 # We write to the array, using the name provided as a keyword argument to the
 # template
 # Note that for subgroups, we do not want to overwrite the full array
+{% if not _target_contiguous %}
+{{_target_var_array}}[{{_target_indices}}] = _numpy.broadcast_to(_synaptic_var, (N, ))
+{% else %}
 {% if _target_start > 0 %}
 _indices = {{_index_array}} - {{_target_start}}
 {% else %}
@@ -32,4 +35,5 @@ _length = _target_stop - {{_target_start}}
 {{_target_var_array}}[{{_target_start}}:_target_stop] = _numpy.bincount(_indices,
                                  minlength=_length,
                                  weights=_numpy.broadcast_to(_synaptic_var, (N, )))
+{% endif %}
 {% endblock %}

brian2/codegen/runtime/numpy_rt/templates/synapses_create_generator.py_

@@ -106,7 +106,11 @@ _vectorisation_idx = 1
 "k" is called the inner variable #}
 
 for _{{outer_index}} in range({{outer_index_size}}):
+    {% if outer_contiguous %}
     _raw{{outer_index_array}} = _{{outer_index}} + {{outer_index_offset}}
+    {% else %}
+    _raw{{outer_index_array}} = {{get_array_name(variables[outer_sub_idx])}}[_{{outer_index}}]
+    {% endif %}
     {% if not result_index_condition %}
     {{vector_code['create_cond']|autoindent}}
     if not _cond:
@@ -126,7 +130,11 @@ for _{{outer_index}} in range({{outer_index_size}}):
     _vectorisation_idx = {{inner_variable}}
     {{vector_code['generator_expr']|autoindent}}
     _vectorisation_idx = _{{result_index}}
-    _raw{{result_index_array}} = _{{result_index}} + {{result_index_offset}}
+    {% if result_contiguous %}
+    _raw{{result_index_array}} = _{{result_index}} + {{result_index_offset}};
+    {% else %}
+    _raw{{result_index_array}} = {{get_array_name(variables[result_sub_idx])}}[_{{result_index}}]
+    {% endif %}
     {% if result_index_condition %}
     {% if result_index_used %}
     {# The condition could index outside of array range #}
@@ -184,7 +192,11 @@ for _{{outer_index}} in range({{outer_index_size}}):
     {% endif %}
 
     _vectorisation_idx = _{{result_index}}
+    {% if result_contiguous %}
     _raw{{result_index_array}} = _{{result_index}} + {{result_index_offset}}
+    {% else %}
+    _raw{{result_index_array}} = {{get_array_name(variables[result_sub_idx])}}[_{{result_index}}]
+    {% endif %}
     {{vector_code['update']|autoindent}}
 
     if not _numpy.isscalar(_n):

brian2/core/network.py

@@ -16,6 +16,8 @@
 from collections import Counter, defaultdict, namedtuple
 from collections.abc import Mapping, Sequence
 
+import numpy as np
+
 from brian2.core.base import BrianObject, BrianObjectException
 from brian2.core.clocks import Clock, defaultclock
 from brian2.core.names import Nameable
@@ -321,18 +323,26 @@ def _check_multiple_summed_updaters(objects):
                         "the target group instead."
                     )
                     raise NotImplementedError(msg)
-                elif (
-                    obj.target.start < other_target.stop
-                    and other_target.start < obj.target.stop
-                ):
-                    # Overlapping subgroups
-                    msg = (
-                        "Multiple 'summed variables' target the "
-                        f"variable '{obj.target_var.name}' in overlapping "
-                        f"groups '{other_target.name}' and '{obj.target.name}'. "
-                        "Use separate variables in the target groups instead."
-                    )
-                    raise NotImplementedError(msg)
+                else:
+                    if getattr(obj.target, "contiguous", True):
+                        target_indices = np.arange(obj.target.start, obj.target.stop)
+                    else:
+                        target_indices = obj.target.indices[:]
+                    if getattr(other_target, "contiguous", True):
+                        other_indices = np.arange(other_target.start, other_target.stop)
+                    else:
+                        other_indices = other_target.indices[:]
+                    if np.intersect1d(
+                        target_indices, other_indices, assume_unique=True
+                    ).size:
+                        # Overlapping subgroups
+                        msg = (
+                            "Multiple 'summed variables' target the "
+                            f"variable '{obj.target_var.name}' in overlapping "
+                            f"groups '{other_target.name}' and '{obj.target.name}'. "
+                            "Use separate variables in the target groups instead."
+                        )
+                        raise NotImplementedError(msg)
             summed_targets[obj.target_var] = obj.target
 
 

brian2/core/variables.py

@@ -196,7 +196,7 @@ def __setstate__(self, state):
 
     @property
     def is_boolean(self):
-        return np.issubdtype(self.dtype, np.bool_)
+        return np.issubdtype(self.dtype, bool)
 
     @property
     def is_integer(self):

brian2/devices/cpp_standalone/codeobject.py

@@ -108,8 +108,12 @@ def __init__(self, *args, **kwds):
         super().__init__(*args, **kwds)
         #: Store whether this code object defines before/after blocks
         self.before_after_blocks = []
+        #: Store variables that are updated by this code object and therefore invalidate the cache
+        self.invalidate_cache_variables = set()
 
     def __call__(self, **kwds):
+        for var in self.invalidate_cache_variables:
+            get_device().array_cache[var] = None
         return self.run()
 
     def compile_block(self, block):