Izhikevich-style neurons on STM32

FIRMWARE/Makefile

@@ -17,9 +17,9 @@
 ## along with this library.  If not, see <http://www.gnu.org/licenses/>.
 ##
 
-BINARY = main
+BINARY ?= main
 
-OBJS = HAL.o comm.o neuron.o
+OBJS = HAL.o comm.o neuron.o izhi.o
 
 OPENCM3_DIR = ../libopencm3
 

FIRMWARE/debug/test/izhi_host.c

@@ -0,0 +1,50 @@
+#include <stdio.h>
+
+#include "./izhi.h"
+
+// build with `gcc $SRCDIR/izhi.c izhi_host.c -I$SRCDIR -o izhi_host`
+// run something like `./izhi_host >neuron_out.dat`
+// and plot with gnuplot something like
+//    `plot 'neuron_out.dat' using 1:2 title 'floating' with lines, \
+//          'neuron_out.dat' using 1:3 title 'fixed' with lines`
+
+int main(void) {
+    fneuron_t spiky_f[7];
+    ineuron_t spiky_i[7];
+
+    RS_f(&spiky_f[0]);
+    IB_f(&spiky_f[1]);
+    CH_f(&spiky_f[2]);
+    FS_f(&spiky_f[3]);
+    LTS_f(&spiky_f[4]);
+    RZ_f(&spiky_f[5]);
+    TC_f(&spiky_f[6]);
+
+    RS_i(&spiky_i[0]);
+    IB_i(&spiky_i[1]);
+    CH_i(&spiky_i[2]);
+    FS_i(&spiky_i[3]);
+    LTS_i(&spiky_i[4]);
+    RZ_i(&spiky_i[5]);
+    TC_i(&spiky_i[6]);
+
+    for (int i = 0; i < 5000; i++) {
+        if (i < 100) {
+            for (int j = 0; j < 7; j++) {
+                step_f(&spiky_f[j], 0, 0.125);
+                step_i(&spiky_i[j], 0, 3);
+            }
+        } else {
+            for (int j = 0; j < 7; j++) {
+                step_f(&spiky_f[j], 10, 0.125);
+                step_i(&spiky_i[j], 10 * spiky_i[j].scale, 3);
+            }
+        }
+        printf("%f ", i * 0.1);
+        for (int j = 0; j < 7; j++) {
+            printf("%f %f ", spiky_f[j].potential,
+                   (float_t)(spiky_i[j].potential) / spiky_i[j].scale);
+        }
+        printf("\n");
+    }
+}

FIRMWARE/izhi.c

@@ -0,0 +1,139 @@
+#include "./izhi.h"
+
+static void base_f(fneuron_t *neuron) {
+    // create a "regular spiking" floating point neuron
+    neuron->a = 0.02;
+    neuron->b = 0.2;
+    neuron->c = -65;
+    neuron->d = 2;
+    neuron->potential = neuron->recovery = 0;
+}
+
+void RS_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->d = 8;
+}
+
+void IB_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->c = -55;
+    neuron->d = 4;
+}
+
+void CH_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->c = -50;
+}
+
+void FS_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->a = 0.1;
+}
+
+void LTS_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->b = 0.25;
+}
+
+void RZ_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->a = 0.1;
+    neuron->b = 0.26;
+}
+
+void TC_f(fneuron_t *neuron) {
+    base_f(neuron);
+    neuron->b = 0.25;
+    neuron->d = 0.05;
+}
+
+void step_f(fneuron_t *neuron, float_t synapse, float_t ms) {
+    // step a neuron through ms milliseconds with synapse input
+    //   if you don't have a good reason to do otherwise, keep ms between 0.1
+    //   and 1.0
+    if (neuron->potential >= 30) {
+        neuron->potential = neuron->c;
+        neuron->recovery += neuron->d;
+        return;
+    }
+    float_t v = neuron->potential;
+    float_t u = neuron->recovery;
+    neuron->potential = v + ms * (0.04 * v * v + 5 * v + 140 - u + synapse);
+    neuron->recovery = u + ms * (neuron->a * (neuron->b * v - u));
+    return;
+}
+
+// XXX Note: if you change LOG_SQRT_FSCALE, you'll need to re-check the math
+// XXX        to confirm that it won't overflow, as well as play with those
+// XXX        constants when calculating `partial`.
+#define LOG_SQRT_FSCALE 10
+#define SQRT_FSCALE (1 << LOG_SQRT_FSCALE)
+#define FSCALE (SQRT_FSCALE * SQRT_FSCALE)
+
+static void base_i(ineuron_t *neuron) {
+    neuron->a_inv = 50;
+    neuron->b_inv = 5;
+    neuron->c = -65 * FSCALE;
+    neuron->d = 2 * FSCALE;
+    neuron->potential = neuron->recovery = 0;
+    neuron->scale = FSCALE;
+}
+
+void RS_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->d = 8 * FSCALE;
+}
+
+void IB_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->c = -55 * FSCALE;
+    neuron->d = 4;
+}
+
+void CH_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->c = -50 * FSCALE;
+}
+
+void FS_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->a_inv = 10;
+}
+
+void LTS_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->b_inv = 4;
+}
+
+void RZ_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->a_inv = 10;
+    neuron->b_inv = 4; // should be 1/0.26, we'll see if this is close enough
+}
+
+void TC_i(ineuron_t *neuron) {
+    base_i(neuron);
+    neuron->d = FSCALE / 20;
+    neuron->b_inv = 4;
+}
+
+void step_i(ineuron_t *neuron, fixed_t synapse, uint8_t fracms) {
+    // step a neuron by 2**(-fracms) milliseconds. synapse input must be scaled
+    //  before being passed to this function.
+    if (neuron->potential >= 30 * FSCALE) {
+        neuron->potential = neuron->c;
+        neuron->recovery += neuron->d;
+        return;
+    }
+    fixed_t v = neuron->potential;
+    fixed_t u = neuron->recovery;
+    //fixed_t partial = (v / SQRT_FSCALE) / 5;
+    fixed_t partial = ((v >> LOG_SQRT_FSCALE) * 819) >> 12;
+    neuron->potential = v + ((partial * partial + 5 * v  + 140 * FSCALE
+                             - u + synapse) >> fracms);
+    //neuron->recovery = u + (((v / neuron->b_inv - u) / neuron->a_inv)
+    //                        >> fracms);
+    neuron->recovery = u + (((v - u * neuron->b_inv)
+                             / (neuron->b_inv * neuron->a_inv)) >> fracms);
+    return;
+}

FIRMWARE/izhi.h

@@ -0,0 +1,39 @@
+#ifndef IZHI_H
+#define IZHI_H
+
+#include <stdint.h>
+
+typedef float float_t;
+typedef struct {
+    float_t a, b, c, d;
+    float_t potential, recovery;
+} fneuron_t;
+
+void RS_f(fneuron_t *neuron);
+void IB_f(fneuron_t *neuron);
+void CH_f(fneuron_t *neuron);
+void FS_f(fneuron_t *neuron);
+void LTS_f(fneuron_t *neuron);
+void RZ_f(fneuron_t *neuron);
+void TC_f(fneuron_t *neuron);
+void step_f(fneuron_t *neuron, float_t synapse, float_t ms);
+
+
+typedef int32_t fixed_t;
+typedef struct {
+    // using 1/a, 1/b because a and b are small fractions
+    fixed_t a_inv, b_inv, c, d;
+    fixed_t potential, recovery;
+    fixed_t scale;
+} ineuron_t;
+
+void RS_i(ineuron_t *neuron);
+void IB_i(ineuron_t *neuron);
+void CH_i(ineuron_t *neuron);
+void FS_i(ineuron_t *neuron);
+void LTS_i(ineuron_t *neuron);
+void RZ_i(ineuron_t *neuron);
+void TC_i(ineuron_t *neuron);
+void step_i(ineuron_t *neuron, fixed_t synapse, uint8_t fracms);
+
+#endif // IZHI_H

FIRMWARE/izhi_demo.c

@@ -0,0 +1,35 @@
+#include "HAL.h"
+#include "neuron.h"
+#include "izhi.h"
+
+int main(void)
+{
+	clock_setup();
+	gpio_setup();
+	tim_setup();
+	// this gives a main_tick roughly every 5.1ms, which makes the neuron
+	//  trigger at a sensible rate on my setup
+	systick_setup(10);
+
+	// initialize neuron
+	ineuron_t neuron;
+	CH_i(&neuron);
+
+	for(;;)
+	{
+		if (main_tick)
+		{
+			main_tick = 0;
+			step_i(&neuron, 10 * neuron.scale, 3);
+			// have the LED in R, G, B for the rest, transition, spike regions
+			//  respectively
+			if (neuron.potential < -60 * neuron.scale) {
+				setLED(1000, 0, 0);
+			} else if (neuron.potential < -40 * neuron.scale) {
+				setLED(0, 1000, 0);
+			} else {
+				setLED(0, 0, 1000);
+			}
+        }
+	}
+}