Remove kernels and convolutions; introduce onReceive and onCondition for neurons

[clinssen]

Problem statement

Current NESTML mixes specification of the model with implementation details. Implementation details are e.g.:

• timestep of the ODE solver; which solver is used
• threshold check: grid-aligned, or “precise” (linear interpolation; regula falsi etc. to find “exact” time of threshold crossing)
• integration sequence: e.g. for simulators with delay, subthreshold integration → events handling → threshold check; but for simulators without delay: events handling → subthreshold integration → threshold check

These details should be filled in during code generation (or even after that) and not be part of the model specification itself. But right now these steps are explicitly specified in one monolithic block (the update block) in the NESTML model.

Proposal

Fundamentally, NESTML specifies “hybrid” systems: with continuous dynamics (ODEs) and discrete-time events. Continuous dynamics can depend on other state variables (e.g. membrane potential is not integrated during refactory state). Integration between events is “free-flight”, i.e. without taking any events into account.

state:
	I_syn pA = 0 pA
	I_syn_aux pA/s = 0 pA/s

input:
	I_stim pA <- continuous

• equations: specify ODEs (only subthreshold/continuous dynamics—no events! No convolutions, no kernels. Alternative name: continuousdynamics, …?). For instance:

	inline my_subexpr = ...
	V_m’ = … -V_m / tau_m + I_syn + I_stim + my_subexpr
	I_syn’ = -I_syn / tau_syn + spike_in   # would be forbidden!

• update: “free-flight” timestep between events. For instance:

	if !is_refractory:
		integrate_odes(V_m, I_syn)
	else:
		integrate_odes(I_syn)

• onEvent(spike_port): handle event arriving at given spiking input port. For instance:

	onEvent(spikes_in):
		I_syn_aux += spikes_in   # using the input port name in an expression refers to the weight

• onCondition(expr): check conditions, e.g. exit refractory state; threshold crossing. For instance:

	onCondition(V_m >= V_th):
		emit_spike()
		V_m = V_reset

The advantage of having conditions separate is that they can be checked at different points in the integration loop: for instance, inside GSL solver “inner loop” (solves nest/nest-simulator#2822), or by finding the precise time at which a condition occurs using regula falsi (or by simply checking at the end of a timestep).

Note that “onCondition” and “onEvent” are also elements from the NineML vocabulary (see e.g. https://nineml-spec.readthedocs.io/en/latest/abstraction_layer.html) as well as LEMS (see e.g. https://lems.github.io/LEMS/example-regimes.html).

NESTML (still) fundamentally sets itself apart from these languages by allowing any kind of imperative programming statements to appear in the onEvent, onCondition and update blocks.

Open issues

• Where would an explicit specification of dendritic (axonal?) delay go? Does onEvent get called from the soma perspective or from the dendrite perspective? (Currently: soma perspective; delay is not explicitly specified in the neuron model.) Spun off to Improve flexibility of emit_spike() #959.

• Kernels and convolutions should be removed as it’s not possible to express whether the event would be processed before or after free-flight update, and how different convolutions would interact, and because convolutions generate new state variables, which is a bit hidden from the user. Solution: Better document how kernels and convolutions work; do not remove them.

• We can additionally specify priorities on the event handlers to be more precise on the sequence in which they are evaluated when things happen simultaneously. For instance: V_m would by itself have crossed threshold, but after processing an inhibitory spike is brought subthreshold.

	@priority=1
	onEvent(spikes_in):
		V_m -= 1   # inhibition

	@priority=2
	onCondition(V_m >= V_th):
		V_m = V_reset
		emit_spike()

In this example, higher-priority condition gets evaluated first, so no spike would be emitted. Reversing the priorities would cause spike to be processed first, so threshold is crossed, so spike is fired.

Perhaps if there is more than one event handler, priorities should be enforced.

• We could generate a flowchart during code generation to illustrate the integration loop sequence of operations. This runs the risk of diverging from the actual generated code in the target language. Perhaps the code itself is the best documentation—provide good code comments.