From 96f881a2b8f664f4370e92ca93d306c8fe0391e8 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Wed, 11 Nov 2020 22:10:18 +0100
Subject: [PATCH 01/19] add STDP synapse unit testing

---
 pynest/nest/tests/test_stdp_synapse.py | 382 +++++++++++++++++++++++++
 1 file changed, 382 insertions(+)
 create mode 100644 pynest/nest/tests/test_stdp_synapse.py

diff --git a/pynest/nest/tests/test_stdp_synapse.py b/pynest/nest/tests/test_stdp_synapse.py
new file mode 100644
index 0000000000..116638a58a
--- /dev/null
+++ b/pynest/nest/tests/test_stdp_synapse.py
@@ -0,0 +1,382 @@
+# -*- coding: utf-8 -*-
+#
+# test_stdp_synapse.py
+#
+# This file is part of NEST.
+#
+# Copyright (C) 2004 The NEST Initiative
+#
+# NEST is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 2 of the License, or
+# (at your option) any later version.
+#
+# NEST is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with NEST.  If not, see <http://www.gnu.org/licenses/>.
+
+import nest
+import unittest
+from math import exp
+import numpy as np
+
+@nest.ll_api.check_stack
+class STDPSynapseTest(unittest.TestCase):
+    """
+    XXX TODO
+    """
+
+    def init_params(self):
+        self.resolution = 0.1  # [ms]
+        self.simulation_duration = 1e+3 #20. # [ms]
+        self.synapse_model = "stdp_synapse"
+        self.presynaptic_firing_rate = 20.  # [ms^-1]
+        self.postsynaptic_firing_rate = 20.  # [ms^-1]
+        self.tau_pre = 16.8
+        self.tau_post = 33.7
+        self.init_weight = .5
+        self.synapse_parameters = {
+            "synapse_model": self.synapse_model,
+            "receptor_type": 0,
+            "delay": self.dendritic_delay,
+            # STDP constants
+            "lambda": 0.01,
+            "alpha": 0.85,
+            "mu_plus": 0.0,
+            "mu_minus": 0.0,
+            "tau_plus": self.tau_pre,
+            "Wmax": 15.0,
+            "weight": self.init_weight
+        }
+        self.neuron_parameters = {
+            "tau_minus": self.tau_post,
+        }
+        '''self.hardcoded_trains_length = 15
+        self.hardcoded_pre_times = (1, 5, 6, 7, 9, 11, 12, 13)
+        self.hardcoded_post_times = (2, 3, 4, 8, 9, 10, 12)'''
+        self.hardcoded_trains_length = 15
+        self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13.])
+        self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12.])
+
+    def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
+        if True:
+            pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest = self.do_the_nest_simulation()
+
+            import matplotlib.pyplot as plt
+
+            fig, ax = plt.subplots(nrows=3)
+            ax1, ax3, ax2 = ax
+
+            n_spikes = len(pre_spikes)
+            for i in range(n_spikes):
+                ax1.plot(2 * [pre_spikes[i]], [0, 1], linewidth=2, color="blue", alpha=.4)
+
+            n_spikes = len(post_spikes)
+            for i in range(n_spikes):
+                ax3.plot(2 * [post_spikes[i]], [0, 1], linewidth=2, color="red", alpha=.4)
+
+            ax2.plot(t_weight_by_nest, weight_by_nest, marker="o", label="nestml")
+            #ax2.plot(t_hist, w_hist_ref, linestyle="--", marker="x", label="ref")
+
+            ax2.set_xlabel("Time [ms]")
+            ax1.set_ylabel("Pre spikes")
+            ax3.set_ylabel("Post spikes")
+            ax2.set_ylabel("w")
+            ax2.legend()
+            for _ax in ax:
+                _ax.grid(which="major", axis="both")
+                _ax.grid(which="minor", axis="x", linestyle=":", alpha=.4)
+                _ax.minorticks_on()
+                _ax.set_xlim(0., self.simulation_duration)
+            fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + ".png", dpi=300)
+
+            t_weight_reproduced_independently, weight_reproduced_independently = self.reproduce_weight_drift(
+                pre_spikes, post_spikes,
+                self.init_weight,
+                fname_snip=fname_snip)
+
+
+            # ``weight_by_nest`` containts only weight values at pre spike times, ``weight_reproduced_independently`` contains the weight at pre *and* post times: check that weights are equal only for pre spike times
+            assert len(weight_by_nest) > 0
+            for idx_pre_spike_nest, t_pre_spike_nest in enumerate(t_weight_by_nest):
+                idx_pre_spike_reproduced_independently = np.argmin((t_pre_spike_nest - t_weight_reproduced_independently)**2)
+                np.testing.assert_allclose(t_pre_spike_nest, t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+                print("testing equal t = " + str(t_pre_spike_nest))
+                print("\tweight_by_nest = " + str(weight_by_nest[idx_pre_spike_nest]))
+                print("\tweight_reproduced_independently = " + str(weight_reproduced_independently[idx_pre_spike_reproduced_independently]))
+                np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest], weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+
+
+    def do_the_nest_simulation(self):
+        """
+        This function is where calls to NEST reside. Returns the generated pre- and post spike sequences and the resulting weight established by STDP.
+        """
+        nest.set_verbosity('M_WARNING')
+        nest.ResetKernel()
+        nest.SetKernelStatus({'resolution': self.resolution})
+
+        neurons = nest.Create(
+            self.nest_neuron_model,
+            2,
+            params=self.neuron_parameters)
+        presynaptic_neuron = neurons[0]
+        postsynaptic_neuron = neurons[1]
+
+        generators = nest.Create(
+            "poisson_generator",
+            2,
+            params=({"rate": self.presynaptic_firing_rate,
+                     "stop": (self.simulation_duration - self.hardcoded_trains_length)},
+                    {"rate": self.postsynaptic_firing_rate,
+                     "stop": (self.simulation_duration - self.hardcoded_trains_length)}))
+        presynaptic_generator = generators[0]
+        postsynaptic_generator = generators[1]
+
+        wr = nest.Create('weight_recorder')
+        nest.CopyModel(self.synapse_model, self.synapse_model + "_rec", {"weight_recorder": wr})
+
+        # While the random sequences, fairly long, would supposedly
+        # reveal small differences in the weight change between NEST
+        # and ours, some low-probability events (say, coinciding
+        # spikes) can well not have occured. To generate and
+        # test every possible combination of pre/post precedence, we
+        # append some hardcoded spike sequences:
+        # pre: 1       5 6 7   9    11 12 13
+        # post:  2 3 4       8 9 10    12
+        (
+            hardcoded_pre_times,
+            hardcoded_post_times
+        ) = [
+            [
+                self.simulation_duration - self.hardcoded_trains_length + float(t)
+                for t in train
+            ] for train in (
+                self.hardcoded_pre_times,
+                self.hardcoded_post_times
+            )
+        ]
+
+        spike_senders = nest.Create(
+            "spike_generator",
+            2,
+            params=({"spike_times": hardcoded_pre_times},
+                    {"spike_times": hardcoded_post_times})
+        )
+        pre_spike_generator = spike_senders[0]
+        post_spike_generator = spike_senders[1]
+
+        # The recorder is to save the randomly generated spike trains.
+        spike_recorder = nest.Create("spike_recorder")
+
+        nest.Connect(presynaptic_generator + pre_spike_generator, presynaptic_neuron,
+                     syn_spec={"synapse_model": "static_synapse", "weight": 9999.})
+        nest.Connect(postsynaptic_generator + post_spike_generator, postsynaptic_neuron,
+                     syn_spec={"synapse_model": "static_synapse", "weight": 9999.})
+        nest.Connect(presynaptic_neuron + postsynaptic_neuron, spike_recorder,
+                     syn_spec={"synapse_model": "static_synapse"})
+        # The synapse of interest itself
+        self.synapse_parameters["synapse_model"] += "_rec"
+        nest.Connect(presynaptic_neuron, postsynaptic_neuron, syn_spec=self.synapse_parameters)
+        self.synapse_parameters["synapse_model"] = self.synapse_model
+        plastic_synapse_of_interest = nest.GetConnections(synapse_model=self.synapse_model + "_rec")
+
+        nest.Simulate(self.simulation_duration)
+
+        all_spikes = nest.GetStatus(spike_recorder, keys='events')[0]
+        pre_spikes = all_spikes['times'][all_spikes['senders'] == presynaptic_neuron.tolist()[0]]
+        post_spikes = all_spikes['times'][all_spikes['senders'] == postsynaptic_neuron.tolist()[0]]
+        import pdb;pdb.set_trace()
+        t_hist = nest.GetStatus(wr, "events")[0]["times"]
+        weight = nest.GetStatus(wr, "events")[0]["weights"]
+
+        #t_hist = nest.GetStatus(plastic_synapse_of_interest, keys='times')[0]
+        #weight = nest.GetStatus(plastic_synapse_of_interest, keys='weight')[0]
+        return pre_spikes, post_spikes, t_hist, weight
+
+
+    def reproduce_weight_drift(self, pre_spikes, post_spikes, initial_weight, fname_snip=""):
+
+        def facilitate(w, Kpre, Wmax_=1.):
+            norm_w = ( w / self.synapse_parameters["Wmax"] ) + ( self.synapse_parameters["lambda"] * pow( 1 - ( w / self.synapse_parameters["Wmax"] ), self.synapse_parameters["mu_plus"] ) * Kpre )
+            if norm_w < 1.0:
+                return norm_w * self.synapse_parameters["Wmax"]
+            else:
+                return self.synapse_parameters["Wmax"]
+
+        def depress(w, Kpost):
+            norm_w = ( w / self.synapse_parameters["Wmax"] ) - ( self.synapse_parameters["alpha"] * self.synapse_parameters["lambda"] * pow( w / self.synapse_parameters["Wmax"], self.synapse_parameters["mu_minus"] ) * Kpost )
+            if norm_w > 0.0:
+                return norm_w * self.synapse_parameters["Wmax"]
+            else:
+                return 0.
+
+        def Kpost_at_time(t, spikes, init=1., inclusive=True):
+            t_curr = 0.
+            Kpost = 0.
+            for spike_idx, t_sp in enumerate(spikes):
+                if t < t_sp:
+                    # integrate to t
+                    Kpost *= exp(-(t - t_curr) / self.tau_post)
+                    return Kpost
+                # integrate to t_sp
+                Kpost *= exp(-(t_sp - t_curr) / self.tau_post)
+                if inclusive:
+                    Kpost += 1.
+                if t == t_sp:
+                    return Kpost
+                if not inclusive:
+                    Kpost += 1.
+                t_curr = t_sp
+            # if we get here, t > t_last_spike
+            # integrate to t
+            Kpost *= exp(-(t - t_curr) / self.tau_post)
+            return Kpost
+
+
+        t = 0.
+        Kpre = 0.
+        weight = initial_weight
+
+        t_log = []
+        w_log = []
+        Kpre_log = []
+
+        # logging
+        t_log.append(t)
+        w_log.append(weight)
+        Kpre_log.append(Kpre)
+
+        post_spikes_delayed = post_spikes + self.dendritic_delay
+
+        while t < self.simulation_duration:
+            print("t = " + str(t))
+
+
+            idx_next_pre_spike = -1
+            if np.where((pre_spikes - t) > 0)[0].size > 0:
+                idx_next_pre_spike = np.where((pre_spikes - t) > 0)[0][0]
+                t_next_pre_spike = pre_spikes[idx_next_pre_spike]
+
+            idx_next_post_spike = -1
+            if np.where((post_spikes_delayed - t) > 0)[0].size > 0:
+                idx_next_post_spike = np.where((post_spikes_delayed - t) > 0)[0][0]
+                t_next_post_spike = post_spikes_delayed[idx_next_post_spike]
+
+            if idx_next_post_spike >= 0 \
+             and t_next_post_spike < t_next_pre_spike:
+                handle_post_spike = True
+                handle_pre_spike = False
+            elif idx_next_pre_spike >= 0 \
+             and t_next_post_spike > t_next_pre_spike:
+                handle_post_spike = False
+                handle_pre_spike = True
+            else:
+                handle_post_spike = idx_next_post_spike >= 0
+                handle_pre_spike = idx_next_pre_spike >= 0
+
+            # integrate to min(t_next_pre_spike, t_next_post_spike)
+            t_next = t
+            if handle_pre_spike:
+                t_next = max(t, t_next_pre_spike)
+            if handle_post_spike:
+                t_next = max(t, t_next_post_spike)
+
+            if t_next == t:
+                # no more spikes to process
+                t_next = self.simulation_duration
+
+            '''# max timestep
+            t_next_ = min(t_next, t + 1E-3)
+            if t_next != t_next_:
+                t_next = t_next_
+                handle_pre_spike = False
+                handle_post_spike = False'''
+
+            h = t_next - t
+            Kpre *= exp(-h / self.tau_pre)
+            t = t_next
+
+            if handle_post_spike:
+                print("Handling post spike at t = " + str(t))
+                #Kpost += 1.
+                print("\tFacilitating from " + str(weight), end="")
+                weight = facilitate(weight, Kpre)
+                print(" to " + str(weight) + " using Kpre = " + str(Kpre))
+
+            if handle_pre_spike:
+                print("Handling pre spike at t = " + str(t))
+                Kpre += 1.
+                #print("\tGetting Kpost at time t = " + str(t))
+                _Kpost = Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight, inclusive=False)
+                print("\tDepressing from " + str(weight), end="")
+                weight = depress(weight, _Kpost)
+                print(" to " + str(weight) + " using Kpost = " + str(_Kpost))
+                #print("\t\tKpost(" + str(t - self.dendritic_delay) + ") = " + str(_Kpost))
+
+            # logging
+            t_log.append(t)
+            w_log.append(weight)
+            Kpre_log.append(Kpre)
+
+        if True:
+            import matplotlib.pyplot as plt
+
+            fig, ax = plt.subplots(nrows=3)
+
+            n_spikes = len(pre_spikes)
+            for i in range(n_spikes):
+                ax[0].plot(2 * [pre_spikes[i]], [0, 1], linewidth=2, color="blue", alpha=.4)
+            ax[0].set_ylabel("Pre spikes")
+            ax0_ = ax[0].twinx()
+            ax0_.plot(t_log, Kpre_log)
+
+            n_spikes = len(post_spikes)
+            for i in range(n_spikes):
+                ax[1].plot(2 * [post_spikes[i]], [0, 1], linewidth=2, color="red", alpha=.4)
+            ax1_ = ax[1].twinx()
+            ax[1].set_ylabel("Post spikes")
+            Kpost_log = [Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight) for t in t_log]
+            ax1_.plot(t_log, Kpost_log)
+
+            ax[2].plot(t_log, w_log, marker="o", label="nestml")
+            ax[2].set_ylabel("w")
+
+            ax[2].set_xlabel("Time [ms]")
+            for _ax in ax:
+                _ax.grid(which="major", axis="both")
+                _ax.grid(which="minor", axis="x", linestyle=":", alpha=.4)
+                _ax.minorticks_on()
+                _ax.set_xlim(0., self.simulation_duration)
+            fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + "_ref.png", dpi=300)
+
+        return t_log, w_log
+
+
+    def test_stdp_synapse(self):
+        self.dendritic_delay = float('nan')
+        self.init_params()
+        for self.dendritic_delay in [1., self.resolution]:
+            self.init_params()
+            for self.nest_neuron_model in ["iaf_psc_exp",
+                                       "iaf_cond_exp"]:
+                                       # "iaf_psc_exp_ps"]:
+                fname_snip = "_[nest_neuron_mdl=" + self.nest_neuron_model + "]"
+                self.do_nest_simulation_and_compare_to_reproduced_weight(fname_snip=fname_snip)
+
+
+def suite():
+    suite = unittest.TestLoader().loadTestsFromTestCase(STDPSynapseTest)
+    return unittest.TestSuite([suite])
+
+
+def run():
+    runner = unittest.TextTestRunner(verbosity=2)
+    runner.run(suite())
+
+
+if __name__ == "__main__":
+    run()

From 75f3a4b6e8f1c3f862de0e9c9c42d0ead8ba414f Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Tue, 20 Apr 2021 18:59:24 +0200
Subject: [PATCH 02/19] allow weight_recorder to record precise spike times

---
 nestkernel/connector_base_impl.h | 1 +
 1 file changed, 1 insertion(+)

diff --git a/nestkernel/connector_base_impl.h b/nestkernel/connector_base_impl.h
index 45117feee2..21091570fe 100644
--- a/nestkernel/connector_base_impl.h
+++ b/nestkernel/connector_base_impl.h
@@ -47,6 +47,7 @@ Connector< ConnectionT >::send_weight_event( const thread tid,
     wr_e.set_port( e.get_port() );
     wr_e.set_rport( e.get_rport() );
     wr_e.set_stamp( e.get_stamp() );
+    wr_e.set_offset( e.get_offset() );
     wr_e.set_sender( e.get_sender() );
     wr_e.set_sender_node_id( kernel().connection_manager.get_source_node_id( tid, syn_id_, lcid ) );
     wr_e.set_weight( e.get_weight() );

From c8547d3c1271799695760ac3b043fd752978046e Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Tue, 20 Apr 2021 18:59:45 +0200
Subject: [PATCH 03/19] allow STDP synapse to use precise spike times

---
 models/stdp_synapse.h | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/models/stdp_synapse.h b/models/stdp_synapse.h
index a2ad628572..64afed317a 100644
--- a/models/stdp_synapse.h
+++ b/models/stdp_synapse.h
@@ -217,7 +217,7 @@ inline void
 stdp_synapse< targetidentifierT >::send( Event& e, thread t, const CommonSynapseProperties& )
 {
   // synapse STDP depressing/facilitation dynamics
-  const double t_spike = e.get_stamp().get_ms();
+  const double t_spike = e.get_stamp().get_ms() - e.get_offset();
 
   // use accessor functions (inherited from Connection< >) to obtain delay and
   // target

From af21d05c21f7d542f5170dffe0702468d88353ba Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Tue, 20 Apr 2021 19:01:55 +0200
Subject: [PATCH 04/19] fix STDP synapse unit test

---
 pynest/nest/tests/test_stdp_synapse.py | 93 +++++++++++---------------
 1 file changed, 38 insertions(+), 55 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_synapse.py b/pynest/nest/tests/test_stdp_synapse.py
index 116638a58a..e087718e16 100644
--- a/pynest/nest/tests/test_stdp_synapse.py
+++ b/pynest/nest/tests/test_stdp_synapse.py
@@ -65,34 +65,11 @@ def init_params(self):
     def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
         if True:
             pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest = self.do_the_nest_simulation()
+            print("For model to be tested:")
+            print("\tpre_spikes = " + str(pre_spikes))
+            print("\tpost_spikes = " + str(post_spikes))
 
-            import matplotlib.pyplot as plt
-
-            fig, ax = plt.subplots(nrows=3)
-            ax1, ax3, ax2 = ax
-
-            n_spikes = len(pre_spikes)
-            for i in range(n_spikes):
-                ax1.plot(2 * [pre_spikes[i]], [0, 1], linewidth=2, color="blue", alpha=.4)
-
-            n_spikes = len(post_spikes)
-            for i in range(n_spikes):
-                ax3.plot(2 * [post_spikes[i]], [0, 1], linewidth=2, color="red", alpha=.4)
-
-            ax2.plot(t_weight_by_nest, weight_by_nest, marker="o", label="nestml")
-            #ax2.plot(t_hist, w_hist_ref, linestyle="--", marker="x", label="ref")
-
-            ax2.set_xlabel("Time [ms]")
-            ax1.set_ylabel("Pre spikes")
-            ax3.set_ylabel("Post spikes")
-            ax2.set_ylabel("w")
-            ax2.legend()
-            for _ax in ax:
-                _ax.grid(which="major", axis="both")
-                _ax.grid(which="minor", axis="x", linestyle=":", alpha=.4)
-                _ax.minorticks_on()
-                _ax.set_xlim(0., self.simulation_duration)
-            fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + ".png", dpi=300)
+            self.plot_weight_evolution(pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest, fname_snip=fname_snip, title_snip=self.nest_neuron_model + " (NEST)")
 
             t_weight_reproduced_independently, weight_reproduced_independently = self.reproduce_weight_drift(
                 pre_spikes, post_spikes,
@@ -189,7 +166,7 @@ def do_the_nest_simulation(self):
         all_spikes = nest.GetStatus(spike_recorder, keys='events')[0]
         pre_spikes = all_spikes['times'][all_spikes['senders'] == presynaptic_neuron.tolist()[0]]
         post_spikes = all_spikes['times'][all_spikes['senders'] == postsynaptic_neuron.tolist()[0]]
-        import pdb;pdb.set_trace()
+
         t_hist = nest.GetStatus(wr, "events")[0]["times"]
         weight = nest.GetStatus(wr, "events")[0]["weights"]
 
@@ -255,7 +232,6 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
         while t < self.simulation_duration:
             print("t = " + str(t))
 
-
             idx_next_pre_spike = -1
             if np.where((pre_spikes - t) > 0)[0].size > 0:
                 idx_next_pre_spike = np.where((pre_spikes - t) > 0)[0][0]
@@ -322,38 +298,45 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
             w_log.append(weight)
             Kpre_log.append(Kpre)
 
-        if True:
-            import matplotlib.pyplot as plt
+        Kpost_log = [Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight) for t in t_log]
+        self.plot_weight_evolution(pre_spikes, post_spikes, t_log, w_log, Kpre_log, Kpost_log, fname_snip=fname_snip + "_ref", title_snip="Reference")
+        
+        return t_log, w_log
 
-            fig, ax = plt.subplots(nrows=3)
+    def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=None, Kpost_log=None, fname_snip="", title_snip=""):
+        import matplotlib.pyplot as plt
 
-            n_spikes = len(pre_spikes)
-            for i in range(n_spikes):
-                ax[0].plot(2 * [pre_spikes[i]], [0, 1], linewidth=2, color="blue", alpha=.4)
-            ax[0].set_ylabel("Pre spikes")
-            ax0_ = ax[0].twinx()
+        fig, ax = plt.subplots(nrows=3)
+
+        n_spikes = len(pre_spikes)
+        for i in range(n_spikes):
+            ax[0].plot(2 * [pre_spikes[i]], [0, 1], linewidth=2, color="blue", alpha=.4)
+        ax[0].set_ylabel("Pre spikes")
+        ax0_ = ax[0].twinx()
+        if Kpre_log:
             ax0_.plot(t_log, Kpre_log)
 
-            n_spikes = len(post_spikes)
-            for i in range(n_spikes):
-                ax[1].plot(2 * [post_spikes[i]], [0, 1], linewidth=2, color="red", alpha=.4)
-            ax1_ = ax[1].twinx()
-            ax[1].set_ylabel("Post spikes")
-            Kpost_log = [Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight) for t in t_log]
+        n_spikes = len(post_spikes)
+        for i in range(n_spikes):
+            ax[1].plot(2 * [post_spikes[i]], [0, 1], linewidth=2, color="red", alpha=.4)
+        ax1_ = ax[1].twinx()
+        ax[1].set_ylabel("Post spikes")
+        if Kpost_log:
             ax1_.plot(t_log, Kpost_log)
 
-            ax[2].plot(t_log, w_log, marker="o", label="nestml")
-            ax[2].set_ylabel("w")
+        ax[2].plot(t_log, w_log, marker="o", label="nestml")
+        ax[2].set_ylabel("w")
 
-            ax[2].set_xlabel("Time [ms]")
-            for _ax in ax:
-                _ax.grid(which="major", axis="both")
-                _ax.grid(which="minor", axis="x", linestyle=":", alpha=.4)
-                _ax.minorticks_on()
-                _ax.set_xlim(0., self.simulation_duration)
-            fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + "_ref.png", dpi=300)
+        ax[2].set_xlabel("Time [ms]")
+        for _ax in ax:
+            _ax.grid(which="major", axis="both")
+            _ax.grid(which="minor", axis="x", linestyle=":", alpha=.4)
+            _ax.minorticks_on()
+            _ax.set_xlim(0., self.simulation_duration)
 
-        return t_log, w_log
+        fig.suptitle(title_snip)
+        fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + ".png", dpi=300)
+        import pdb;pdb.set_trace()
 
 
     def test_stdp_synapse(self):
@@ -362,8 +345,8 @@ def test_stdp_synapse(self):
         for self.dendritic_delay in [1., self.resolution]:
             self.init_params()
             for self.nest_neuron_model in ["iaf_psc_exp",
-                                       "iaf_cond_exp"]:
-                                       # "iaf_psc_exp_ps"]:
+                                       "iaf_cond_exp",
+                                       "iaf_psc_exp_ps"]:
                 fname_snip = "_[nest_neuron_mdl=" + self.nest_neuron_model + "]"
                 self.do_nest_simulation_and_compare_to_reproduced_weight(fname_snip=fname_snip)
 

From a495dfebfc84d6a16176fbeae52d6ae71f3a8492 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Wed, 21 Apr 2021 21:23:10 +0200
Subject: [PATCH 05/19] pycodestyle formatting changes

---
 pynest/nest/tests/test_stdp_synapse.py | 46 ++++++++++----------------
 1 file changed, 18 insertions(+), 28 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_synapse.py b/pynest/nest/tests/test_stdp_synapse.py
index e087718e16..25ae9750d3 100644
--- a/pynest/nest/tests/test_stdp_synapse.py
+++ b/pynest/nest/tests/test_stdp_synapse.py
@@ -24,6 +24,7 @@
 from math import exp
 import numpy as np
 
+
 @nest.ll_api.check_stack
 class STDPSynapseTest(unittest.TestCase):
     """
@@ -31,11 +32,11 @@ class STDPSynapseTest(unittest.TestCase):
     """
 
     def init_params(self):
-        self.resolution = 0.1  # [ms]
-        self.simulation_duration = 1e+3 #20. # [ms]
+        self.resolution = 0.1    # [ms]
+        self.simulation_duration = 1E3    # [ms]
         self.synapse_model = "stdp_synapse"
-        self.presynaptic_firing_rate = 20.  # [ms^-1]
-        self.postsynaptic_firing_rate = 20.  # [ms^-1]
+        self.presynaptic_firing_rate = 20.    # [ms^-1]
+        self.postsynaptic_firing_rate = 20.    # [ms^-1]
         self.tau_pre = 16.8
         self.tau_post = 33.7
         self.init_weight = .5
@@ -55,9 +56,6 @@ def init_params(self):
         self.neuron_parameters = {
             "tau_minus": self.tau_post,
         }
-        '''self.hardcoded_trains_length = 15
-        self.hardcoded_pre_times = (1, 5, 6, 7, 9, 11, 12, 13)
-        self.hardcoded_post_times = (2, 3, 4, 8, 9, 10, 12)'''
         self.hardcoded_trains_length = 15
         self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13.])
         self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12.])
@@ -69,14 +67,17 @@ def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
             print("\tpre_spikes = " + str(pre_spikes))
             print("\tpost_spikes = " + str(post_spikes))
 
-            self.plot_weight_evolution(pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest, fname_snip=fname_snip, title_snip=self.nest_neuron_model + " (NEST)")
+            self.plot_weight_evolution(pre_spikes, post_spikes,
+                                       t_weight_by_nest,
+                                       weight_by_nest,
+                                       fname_snip=fname_snip,
+                                       title_snip=self.nest_neuron_model + " (NEST)")
 
             t_weight_reproduced_independently, weight_reproduced_independently = self.reproduce_weight_drift(
                 pre_spikes, post_spikes,
                 self.init_weight,
                 fname_snip=fname_snip)
 
-
             # ``weight_by_nest`` containts only weight values at pre spike times, ``weight_reproduced_independently`` contains the weight at pre *and* post times: check that weights are equal only for pre spike times
             assert len(weight_by_nest) > 0
             for idx_pre_spike_nest, t_pre_spike_nest in enumerate(t_weight_by_nest):
@@ -87,7 +88,6 @@ def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
                 print("\tweight_reproduced_independently = " + str(weight_reproduced_independently[idx_pre_spike_reproduced_independently]))
                 np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest], weight_reproduced_independently[idx_pre_spike_reproduced_independently])
 
-
     def do_the_nest_simulation(self):
         """
         This function is where calls to NEST reside. Returns the generated pre- and post spike sequences and the resulting weight established by STDP.
@@ -170,22 +170,19 @@ def do_the_nest_simulation(self):
         t_hist = nest.GetStatus(wr, "events")[0]["times"]
         weight = nest.GetStatus(wr, "events")[0]["weights"]
 
-        #t_hist = nest.GetStatus(plastic_synapse_of_interest, keys='times')[0]
-        #weight = nest.GetStatus(plastic_synapse_of_interest, keys='weight')[0]
         return pre_spikes, post_spikes, t_hist, weight
 
-
     def reproduce_weight_drift(self, pre_spikes, post_spikes, initial_weight, fname_snip=""):
 
         def facilitate(w, Kpre, Wmax_=1.):
-            norm_w = ( w / self.synapse_parameters["Wmax"] ) + ( self.synapse_parameters["lambda"] * pow( 1 - ( w / self.synapse_parameters["Wmax"] ), self.synapse_parameters["mu_plus"] ) * Kpre )
+            norm_w = (w / self.synapse_parameters["Wmax"]) + (self.synapse_parameters["lambda"] * pow(1 - (w / self.synapse_parameters["Wmax"]), self.synapse_parameters["mu_plus"]) * Kpre)
             if norm_w < 1.0:
                 return norm_w * self.synapse_parameters["Wmax"]
             else:
                 return self.synapse_parameters["Wmax"]
 
         def depress(w, Kpost):
-            norm_w = ( w / self.synapse_parameters["Wmax"] ) - ( self.synapse_parameters["alpha"] * self.synapse_parameters["lambda"] * pow( w / self.synapse_parameters["Wmax"], self.synapse_parameters["mu_minus"] ) * Kpost )
+            norm_w = (w / self.synapse_parameters["Wmax"]) - (self.synapse_parameters["alpha"] * self.synapse_parameters["lambda"] * pow(w / self.synapse_parameters["Wmax"], self.synapse_parameters["mu_minus"]) * Kpost)
             if norm_w > 0.0:
                 return norm_w * self.synapse_parameters["Wmax"]
             else:
@@ -213,7 +210,6 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
             Kpost *= exp(-(t - t_curr) / self.tau_post)
             return Kpost
 
-
         t = 0.
         Kpre = 0.
         weight = initial_weight
@@ -242,12 +238,10 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
                 idx_next_post_spike = np.where((post_spikes_delayed - t) > 0)[0][0]
                 t_next_post_spike = post_spikes_delayed[idx_next_post_spike]
 
-            if idx_next_post_spike >= 0 \
-             and t_next_post_spike < t_next_pre_spike:
+            if idx_next_post_spike >= 0 and t_next_post_spike < t_next_pre_spike:
                 handle_post_spike = True
                 handle_pre_spike = False
-            elif idx_next_pre_spike >= 0 \
-             and t_next_post_spike > t_next_pre_spike:
+            elif idx_next_pre_spike >= 0 and t_next_post_spike > t_next_pre_spike:
                 handle_post_spike = False
                 handle_pre_spike = True
             else:
@@ -278,7 +272,7 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
 
             if handle_post_spike:
                 print("Handling post spike at t = " + str(t))
-                #Kpost += 1.
+                # Kpost += 1.    <-- not necessary, will call Kpost_at_time(t) later to compute Kpost for any value t
                 print("\tFacilitating from " + str(weight), end="")
                 weight = facilitate(weight, Kpre)
                 print(" to " + str(weight) + " using Kpre = " + str(Kpre))
@@ -286,12 +280,10 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
             if handle_pre_spike:
                 print("Handling pre spike at t = " + str(t))
                 Kpre += 1.
-                #print("\tGetting Kpost at time t = " + str(t))
                 _Kpost = Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight, inclusive=False)
                 print("\tDepressing from " + str(weight), end="")
                 weight = depress(weight, _Kpost)
                 print(" to " + str(weight) + " using Kpost = " + str(_Kpost))
-                #print("\t\tKpost(" + str(t - self.dendritic_delay) + ") = " + str(_Kpost))
 
             # logging
             t_log.append(t)
@@ -300,7 +292,7 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
 
         Kpost_log = [Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight) for t in t_log]
         self.plot_weight_evolution(pre_spikes, post_spikes, t_log, w_log, Kpre_log, Kpost_log, fname_snip=fname_snip + "_ref", title_snip="Reference")
-        
+
         return t_log, w_log
 
     def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=None, Kpost_log=None, fname_snip="", title_snip=""):
@@ -336,8 +328,6 @@ def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=
 
         fig.suptitle(title_snip)
         fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + ".png", dpi=300)
-        import pdb;pdb.set_trace()
-
 
     def test_stdp_synapse(self):
         self.dendritic_delay = float('nan')
@@ -345,8 +335,8 @@ def test_stdp_synapse(self):
         for self.dendritic_delay in [1., self.resolution]:
             self.init_params()
             for self.nest_neuron_model in ["iaf_psc_exp",
-                                       "iaf_cond_exp",
-                                       "iaf_psc_exp_ps"]:
+                                           "iaf_cond_exp",
+                                           "iaf_psc_exp_ps"]:
                 fname_snip = "_[nest_neuron_mdl=" + self.nest_neuron_model + "]"
                 self.do_nest_simulation_and_compare_to_reproduced_weight(fname_snip=fname_snip)
 

From a27b55c5fdd962b819e82ed59fdeb9c52e1bc82c Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 3 May 2021 18:35:17 +0200
Subject: [PATCH 06/19] fix test conditions and add plots to
 test_stdp_multiplicity

---
 pynest/nest/tests/test_stdp_multiplicity.py | 59 ++++++++++++++-------
 1 file changed, 41 insertions(+), 18 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_multiplicity.py b/pynest/nest/tests/test_stdp_multiplicity.py
index 98abbc798b..bb1670d6e5 100644
--- a/pynest/nest/tests/test_stdp_multiplicity.py
+++ b/pynest/nest/tests/test_stdp_multiplicity.py
@@ -26,6 +26,13 @@
 import math
 import numpy as np
 
+try:
+    import matplotlib as mpl
+    import matplotlib.pyplot as plt
+    DEBUG_PLOTS = True
+except Exception:
+    DEBUG_PLOTS = False
+
 
 @nest.ll_api.check_stack
 class StdpSpikeMultiplicity(unittest.TestCase):
@@ -194,33 +201,49 @@ def run_protocol(self, pre_post_shift):
             post_weights['parrot'].append(w_post)
             post_weights['parrot_ps'].append(w_post_ps)
 
+        if DEBUG_PLOTS:
+            import datetime
+            fig, ax = plt.subplots(nrows=2)
+            fig.suptitle("Final obtained weights")
+            ax[0].plot(post_weights["parrot"], marker="o", label="parrot")
+            ax[0].plot(post_weights["parrot_ps"], marker="o", label="parrot_ps")
+            ax[0].set_ylabel("final weight")
+            ax[0].set_xticklabels([])
+            ax[1].semilogy(np.abs(np.array(post_weights["parrot"]) - np.array(post_weights["parrot_ps"])), marker="o", label="error")
+            ax[1].set_xticks([i for i in range(len(deltas))])
+            ax[1].set_xticklabels(["{0:.1E}".format(d) for d in deltas])
+            ax[1].set_xlabel("timestep [ms]")
+            for _ax in ax:
+                _ax.grid(True)
+                _ax.legend()
+            plt.savefig("/tmp/test_stdp_multiplicity" + str(datetime.datetime.utcnow()) + ".png")
+        print(post_weights)
         return post_weights
 
-    def test_ParrotNeuronSTDPProtocolPotentiation(self):
-        """Check weight convergence on potentiation."""
+    # XXX: TODO: use ``@pytest.mark.parametrize`` for this
+    def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-6):
+        """Check that for smaller and smaller timestep, weights obtained from parrot and precise parrot converge.
+
+        Enforce a maximum allowed absolute error ``max_abs_err`` between the final weights for the smallest timestep tested.
+
+        Enforce that the error should strictly decrease with smaller timestep."""
 
-        post_weights = self.run_protocol(pre_post_shift=10.0)
+        post_weights = self.run_protocol(pre_post_shift=pre_post_shift)
         w_plain = np.array(post_weights['parrot'])
         w_precise = np.array(post_weights['parrot_ps'])
 
-        assert all(w_plain == w_plain[0]), 'Plain weights differ'
-        dw = w_precise - w_plain
-        dwrel = dw[1:] / dw[:-1]
-        assert all(np.round(dwrel, decimals=3) ==
-                   0.5), 'Precise weights do not converge.'
+        assert all(w_plain == w_plain[0]), 'Plain weights should be independent of timestep, but differ!'
+        abs_err = np.abs(w_precise - w_plain)
+        assert abs_err[-1] < max_abs_err
+        assert np.all(np.diff(abs_err) < 1), 'Error should decrease with smaller timestep!'
 
-    def test_ParrotNeuronSTDPProtocolDepression(self):
-        """Check weight convergence on depression."""
+    def test_stdp_multiplicity(self):
+        """Check weight convergence on potentiation and depression.
 
-        post_weights = self.run_protocol(pre_post_shift=-10.0)
-        w_plain = np.array(post_weights['parrot'])
-        w_precise = np.array(post_weights['parrot_ps'])
+        See also: _test_stdp_multiplicity()."""
+        self._test_stdp_multiplicity(pre_post_shift=10.)    # test potentiation
+        self._test_stdp_multiplicity(pre_post_shift=-10.)   # test depression
 
-        assert all(w_plain == w_plain[0]), 'Plain weights differ'
-        dw = w_precise - w_plain
-        dwrel = dw[1:] / dw[:-1]
-        assert all(np.round(dwrel, decimals=3) ==
-                   0.5), 'Precise weights do not converge.'
 
 
 def suite():

From 305e82c9af695395fb91cf5fec7e4973ec61dfc9 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 3 May 2021 18:38:04 +0200
Subject: [PATCH 07/19] fix Python code formatting

---
 pynest/nest/tests/test_stdp_multiplicity.py | 10 +++----
 pynest/nest/tests/test_stdp_synapse.py      | 32 ++++++++++++++-------
 2 files changed, 27 insertions(+), 15 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_multiplicity.py b/pynest/nest/tests/test_stdp_multiplicity.py
index bb1670d6e5..74c5d4d316 100644
--- a/pynest/nest/tests/test_stdp_multiplicity.py
+++ b/pynest/nest/tests/test_stdp_multiplicity.py
@@ -164,8 +164,7 @@ def run_protocol(self, pre_post_shift):
             # create spike recorder --- debugging only
             spikes = nest.Create("spike_recorder")
             nest.Connect(
-                pre_parrot + post_parrot +
-                pre_parrot_ps + post_parrot_ps,
+                pre_parrot + post_parrot + pre_parrot_ps + post_parrot_ps,
                 spikes
             )
 
@@ -209,7 +208,8 @@ def run_protocol(self, pre_post_shift):
             ax[0].plot(post_weights["parrot_ps"], marker="o", label="parrot_ps")
             ax[0].set_ylabel("final weight")
             ax[0].set_xticklabels([])
-            ax[1].semilogy(np.abs(np.array(post_weights["parrot"]) - np.array(post_weights["parrot_ps"])), marker="o", label="error")
+            ax[1].semilogy(np.abs(np.array(post_weights["parrot"]) - np.array(post_weights["parrot_ps"])),
+                           marker="o", label="error")
             ax[1].set_xticks([i for i in range(len(deltas))])
             ax[1].set_xticklabels(["{0:.1E}".format(d) for d in deltas])
             ax[1].set_xlabel("timestep [ms]")
@@ -224,7 +224,8 @@ def run_protocol(self, pre_post_shift):
     def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-6):
         """Check that for smaller and smaller timestep, weights obtained from parrot and precise parrot converge.
 
-        Enforce a maximum allowed absolute error ``max_abs_err`` between the final weights for the smallest timestep tested.
+        Enforce a maximum allowed absolute error ``max_abs_err`` between the final weights for the smallest timestep
+        tested.
 
         Enforce that the error should strictly decrease with smaller timestep."""
 
@@ -245,7 +246,6 @@ def test_stdp_multiplicity(self):
         self._test_stdp_multiplicity(pre_post_shift=-10.)   # test depression
 
 
-
 def suite():
 
     # makeSuite is sort of obsolete http://bugs.python.org/issue2721
diff --git a/pynest/nest/tests/test_stdp_synapse.py b/pynest/nest/tests/test_stdp_synapse.py
index 25ae9750d3..c722a25239 100644
--- a/pynest/nest/tests/test_stdp_synapse.py
+++ b/pynest/nest/tests/test_stdp_synapse.py
@@ -78,19 +78,25 @@ def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
                 self.init_weight,
                 fname_snip=fname_snip)
 
-            # ``weight_by_nest`` containts only weight values at pre spike times, ``weight_reproduced_independently`` contains the weight at pre *and* post times: check that weights are equal only for pre spike times
+            # ``weight_by_nest`` containts only weight values at pre spike times, ``weight_reproduced_independently``
+            # contains the weight at pre *and* post times: check that weights are equal only for pre spike times
             assert len(weight_by_nest) > 0
             for idx_pre_spike_nest, t_pre_spike_nest in enumerate(t_weight_by_nest):
-                idx_pre_spike_reproduced_independently = np.argmin((t_pre_spike_nest - t_weight_reproduced_independently)**2)
-                np.testing.assert_allclose(t_pre_spike_nest, t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+                idx_pre_spike_reproduced_independently = \
+                    np.argmin((t_pre_spike_nest - t_weight_reproduced_independently)**2)
+                np.testing.assert_allclose(t_pre_spike_nest,
+                                           t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
                 print("testing equal t = " + str(t_pre_spike_nest))
                 print("\tweight_by_nest = " + str(weight_by_nest[idx_pre_spike_nest]))
-                print("\tweight_reproduced_independently = " + str(weight_reproduced_independently[idx_pre_spike_reproduced_independently]))
-                np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest], weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+                print("\tweight_reproduced_independently = " + str(
+                      weight_reproduced_independently[idx_pre_spike_reproduced_independently]))
+                np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest],
+                                           weight_reproduced_independently[idx_pre_spike_reproduced_independently])
 
     def do_the_nest_simulation(self):
         """
-        This function is where calls to NEST reside. Returns the generated pre- and post spike sequences and the resulting weight established by STDP.
+        This function is where calls to NEST reside. Returns the generated pre- and post spike sequences and the
+        resulting weight established by STDP.
         """
         nest.set_verbosity('M_WARNING')
         nest.ResetKernel()
@@ -175,14 +181,18 @@ def do_the_nest_simulation(self):
     def reproduce_weight_drift(self, pre_spikes, post_spikes, initial_weight, fname_snip=""):
 
         def facilitate(w, Kpre, Wmax_=1.):
-            norm_w = (w / self.synapse_parameters["Wmax"]) + (self.synapse_parameters["lambda"] * pow(1 - (w / self.synapse_parameters["Wmax"]), self.synapse_parameters["mu_plus"]) * Kpre)
+            norm_w = (w / self.synapse_parameters["Wmax"]) + (
+                self.synapse_parameters["lambda"] * pow(
+                    1 - (w / self.synapse_parameters["Wmax"]), self.synapse_parameters["mu_plus"]) * Kpre)
             if norm_w < 1.0:
                 return norm_w * self.synapse_parameters["Wmax"]
             else:
                 return self.synapse_parameters["Wmax"]
 
         def depress(w, Kpost):
-            norm_w = (w / self.synapse_parameters["Wmax"]) - (self.synapse_parameters["alpha"] * self.synapse_parameters["lambda"] * pow(w / self.synapse_parameters["Wmax"], self.synapse_parameters["mu_minus"]) * Kpost)
+            norm_w = (w / self.synapse_parameters["Wmax"]) - (
+                self.synapse_parameters["alpha"] * self.synapse_parameters["lambda"] * pow(
+                    w / self.synapse_parameters["Wmax"], self.synapse_parameters["mu_minus"]) * Kpost)
             if norm_w > 0.0:
                 return norm_w * self.synapse_parameters["Wmax"]
             else:
@@ -291,11 +301,13 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
             Kpre_log.append(Kpre)
 
         Kpost_log = [Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight) for t in t_log]
-        self.plot_weight_evolution(pre_spikes, post_spikes, t_log, w_log, Kpre_log, Kpost_log, fname_snip=fname_snip + "_ref", title_snip="Reference")
+        self.plot_weight_evolution(pre_spikes, post_spikes, t_log, w_log, Kpre_log, Kpost_log,
+                                   fname_snip=fname_snip + "_ref", title_snip="Reference")
 
         return t_log, w_log
 
-    def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=None, Kpost_log=None, fname_snip="", title_snip=""):
+    def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=None, Kpost_log=None,
+                              fname_snip="", title_snip=""):
         import matplotlib.pyplot as plt
 
         fig, ax = plt.subplots(nrows=3)

From ff088921b61ba0eb7439159a99fcc18499751328 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 3 May 2021 19:25:28 +0200
Subject: [PATCH 08/19] clean up code in unit test (#1840)

---
 pynest/nest/tests/test_stdp_synapse.py | 72 ++++++++++++--------------
 1 file changed, 32 insertions(+), 40 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_synapse.py b/pynest/nest/tests/test_stdp_synapse.py
index c722a25239..a01aa5b508 100644
--- a/pynest/nest/tests/test_stdp_synapse.py
+++ b/pynest/nest/tests/test_stdp_synapse.py
@@ -24,11 +24,21 @@
 from math import exp
 import numpy as np
 
+try:
+    import matplotlib as mpl
+    import matplotlib.pyplot as plt
+    DEBUG_PLOTS = True
+except Exception:
+    DEBUG_PLOTS = False
+
 
 @nest.ll_api.check_stack
 class STDPSynapseTest(unittest.TestCase):
     """
-    XXX TODO
+    Compare the STDP synaptic plasticity model against a self-contained Python reference.
+
+    Random pre and post spike times are generated according to a Poisson distribution; some hard-coded spike times are
+    added to make sure to test for edge cases such as simultaneous pre- and post spike.
     """
 
     def init_params(self):
@@ -61,37 +71,29 @@ def init_params(self):
         self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12.])
 
     def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
-        if True:
-            pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest = self.do_the_nest_simulation()
-            print("For model to be tested:")
-            print("\tpre_spikes = " + str(pre_spikes))
-            print("\tpost_spikes = " + str(post_spikes))
-
+        pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest = self.do_the_nest_simulation()
+        if DEBUG_PLOTS:
             self.plot_weight_evolution(pre_spikes, post_spikes,
                                        t_weight_by_nest,
                                        weight_by_nest,
                                        fname_snip=fname_snip,
                                        title_snip=self.nest_neuron_model + " (NEST)")
 
-            t_weight_reproduced_independently, weight_reproduced_independently = self.reproduce_weight_drift(
-                pre_spikes, post_spikes,
-                self.init_weight,
-                fname_snip=fname_snip)
-
-            # ``weight_by_nest`` containts only weight values at pre spike times, ``weight_reproduced_independently``
-            # contains the weight at pre *and* post times: check that weights are equal only for pre spike times
-            assert len(weight_by_nest) > 0
-            for idx_pre_spike_nest, t_pre_spike_nest in enumerate(t_weight_by_nest):
-                idx_pre_spike_reproduced_independently = \
-                    np.argmin((t_pre_spike_nest - t_weight_reproduced_independently)**2)
-                np.testing.assert_allclose(t_pre_spike_nest,
-                                           t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
-                print("testing equal t = " + str(t_pre_spike_nest))
-                print("\tweight_by_nest = " + str(weight_by_nest[idx_pre_spike_nest]))
-                print("\tweight_reproduced_independently = " + str(
-                      weight_reproduced_independently[idx_pre_spike_reproduced_independently]))
-                np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest],
-                                           weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+        t_weight_reproduced_independently, weight_reproduced_independently = self.reproduce_weight_drift(
+            pre_spikes, post_spikes,
+            self.init_weight,
+            fname_snip=fname_snip)
+
+        # ``weight_by_nest`` containts only weight values at pre spike times, ``weight_reproduced_independently``
+        # contains the weight at pre *and* post times: check that weights are equal only for pre spike times
+        assert len(weight_by_nest) > 0
+        for idx_pre_spike_nest, t_pre_spike_nest in enumerate(t_weight_by_nest):
+            idx_pre_spike_reproduced_independently = \
+                np.argmin((t_pre_spike_nest - t_weight_reproduced_independently)**2)
+            np.testing.assert_allclose(t_pre_spike_nest,
+                                        t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+            np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest],
+                                        weight_reproduced_independently[idx_pre_spike_reproduced_independently])
 
     def do_the_nest_simulation(self):
         """
@@ -165,7 +167,6 @@ def do_the_nest_simulation(self):
         self.synapse_parameters["synapse_model"] += "_rec"
         nest.Connect(presynaptic_neuron, postsynaptic_neuron, syn_spec=self.synapse_parameters)
         self.synapse_parameters["synapse_model"] = self.synapse_model
-        plastic_synapse_of_interest = nest.GetConnections(synapse_model=self.synapse_model + "_rec")
 
         nest.Simulate(self.simulation_duration)
 
@@ -179,7 +180,7 @@ def do_the_nest_simulation(self):
         return pre_spikes, post_spikes, t_hist, weight
 
     def reproduce_weight_drift(self, pre_spikes, post_spikes, initial_weight, fname_snip=""):
-
+        """Independent, self-contained model of STDP"""
         def facilitate(w, Kpre, Wmax_=1.):
             norm_w = (w / self.synapse_parameters["Wmax"]) + (
                 self.synapse_parameters["lambda"] * pow(
@@ -236,8 +237,6 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
         post_spikes_delayed = post_spikes + self.dendritic_delay
 
         while t < self.simulation_duration:
-            print("t = " + str(t))
-
             idx_next_pre_spike = -1
             if np.where((pre_spikes - t) > 0)[0].size > 0:
                 idx_next_pre_spike = np.where((pre_spikes - t) > 0)[0][0]
@@ -281,19 +280,13 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
             t = t_next
 
             if handle_post_spike:
-                print("Handling post spike at t = " + str(t))
                 # Kpost += 1.    <-- not necessary, will call Kpost_at_time(t) later to compute Kpost for any value t
-                print("\tFacilitating from " + str(weight), end="")
                 weight = facilitate(weight, Kpre)
-                print(" to " + str(weight) + " using Kpre = " + str(Kpre))
 
             if handle_pre_spike:
-                print("Handling pre spike at t = " + str(t))
                 Kpre += 1.
                 _Kpost = Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight, inclusive=False)
-                print("\tDepressing from " + str(weight), end="")
                 weight = depress(weight, _Kpost)
-                print(" to " + str(weight) + " using Kpost = " + str(_Kpost))
 
             # logging
             t_log.append(t)
@@ -301,15 +294,14 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
             Kpre_log.append(Kpre)
 
         Kpost_log = [Kpost_at_time(t - self.dendritic_delay, post_spikes, init=self.init_weight) for t in t_log]
-        self.plot_weight_evolution(pre_spikes, post_spikes, t_log, w_log, Kpre_log, Kpost_log,
-                                   fname_snip=fname_snip + "_ref", title_snip="Reference")
+        if DEBUG_PLOTS:
+            self.plot_weight_evolution(pre_spikes, post_spikes, t_log, w_log, Kpre_log, Kpost_log,
+                                       fname_snip=fname_snip + "_ref", title_snip="Reference")
 
         return t_log, w_log
 
     def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=None, Kpost_log=None,
                               fname_snip="", title_snip=""):
-        import matplotlib.pyplot as plt
-
         fig, ax = plt.subplots(nrows=3)
 
         n_spikes = len(pre_spikes)

From b718ca5b682015acd1a740de324de1e3c97acef0 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 3 May 2021 20:16:48 +0200
Subject: [PATCH 09/19] fix intermittent failure in test_visualisation.py

---
 pynest/nest/tests/test_visualization.py | 9 ++++++---
 1 file changed, 6 insertions(+), 3 deletions(-)

diff --git a/pynest/nest/tests/test_visualization.py b/pynest/nest/tests/test_visualization.py
index b49538b245..fe4401dcf0 100644
--- a/pynest/nest/tests/test_visualization.py
+++ b/pynest/nest/tests/test_visualization.py
@@ -95,7 +95,7 @@ def voltage_trace_verify(self, device):
     @unittest.skipIf(not PLOTTING_POSSIBLE, 'Plotting impossible because matplotlib or display missing')
     def test_voltage_trace_from_device(self):
         """Test voltage_trace from device"""
-        import nest.voltage_trace as nvtrace
+        import nest.voltage_trace
         nest.ResetKernel()
         nodes = nest.Create('iaf_psc_alpha', 2)
         pg = nest.Create('poisson_generator', 1, {'rate': 1000.})
@@ -105,10 +105,11 @@ def test_voltage_trace_from_device(self):
         nest.Simulate(100)
 
         # Test with data from device
+        plt.close("all")
         nest.voltage_trace.from_device(device)
         self.voltage_trace_verify(device)
 
-        # Test with fata from file
+        # Test with data from file
         vm = device.get('events')
         data = np.zeros([len(vm['senders']), 3])
         data[:, 0] = vm['senders']
@@ -117,6 +118,8 @@ def test_voltage_trace_from_device(self):
         filename = os.path.join(self.nest_tmpdir(), 'voltage_trace.txt')
         self.filenames.append(filename)
         np.savetxt(filename, data)
+
+        plt.close("all")
         nest.voltage_trace.from_file(filename)
         self.voltage_trace_verify(device)
 
@@ -150,7 +153,7 @@ def spike_recorder_raster_verify(self, sr_ref):
     @unittest.skipIf(not PLOTTING_POSSIBLE, 'Plotting impossible because matplotlib or display missing')
     def test_raster_plot(self):
         """Test raster_plot"""
-        import nest.raster_plot as nraster
+        import nest.raster_plot
 
         sr, sr_to_file = self.spike_recorder_data_setup(to_file=True)
         spikes = sr.get('events')

From b6cfceec5d55ea8c830ae48b52374188cb42fa98 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Tue, 4 May 2021 16:35:19 +0200
Subject: [PATCH 10/19] close matplotlib figures after saving them to file

---
 pynest/nest/tests/test_stdp_multiplicity.py | 1 +
 pynest/nest/tests/test_stdp_synapse.py      | 5 +++--
 2 files changed, 4 insertions(+), 2 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_multiplicity.py b/pynest/nest/tests/test_stdp_multiplicity.py
index 74c5d4d316..544f86144f 100644
--- a/pynest/nest/tests/test_stdp_multiplicity.py
+++ b/pynest/nest/tests/test_stdp_multiplicity.py
@@ -217,6 +217,7 @@ def run_protocol(self, pre_post_shift):
                 _ax.grid(True)
                 _ax.legend()
             plt.savefig("/tmp/test_stdp_multiplicity" + str(datetime.datetime.utcnow()) + ".png")
+            plt.close(fig)
         print(post_weights)
         return post_weights
 
diff --git a/pynest/nest/tests/test_stdp_synapse.py b/pynest/nest/tests/test_stdp_synapse.py
index a01aa5b508..a634d90344 100644
--- a/pynest/nest/tests/test_stdp_synapse.py
+++ b/pynest/nest/tests/test_stdp_synapse.py
@@ -91,9 +91,9 @@ def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
             idx_pre_spike_reproduced_independently = \
                 np.argmin((t_pre_spike_nest - t_weight_reproduced_independently)**2)
             np.testing.assert_allclose(t_pre_spike_nest,
-                                        t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+                                       t_weight_reproduced_independently[idx_pre_spike_reproduced_independently])
             np.testing.assert_allclose(weight_by_nest[idx_pre_spike_nest],
-                                        weight_reproduced_independently[idx_pre_spike_reproduced_independently])
+                                       weight_reproduced_independently[idx_pre_spike_reproduced_independently])
 
     def do_the_nest_simulation(self):
         """
@@ -332,6 +332,7 @@ def plot_weight_evolution(self, pre_spikes, post_spikes, t_log, w_log, Kpre_log=
 
         fig.suptitle(title_snip)
         fig.savefig("/tmp/nest_stdp_synapse_test" + fname_snip + ".png", dpi=300)
+        plt.close(fig)
 
     def test_stdp_synapse(self):
         self.dendritic_delay = float('nan')

From c878c9c3976f3c255d49b3764e545033b03260a3 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 30 Aug 2021 13:52:28 +0200
Subject: [PATCH 11/19] clean up STDP multiplicity test

---
 pynest/nest/tests/test_stdp_multiplicity.py | 26 +++++++++++++--------
 1 file changed, 16 insertions(+), 10 deletions(-)

diff --git a/pynest/nest/tests/test_stdp_multiplicity.py b/pynest/nest/tests/test_stdp_multiplicity.py
index 544f86144f..8c3ab0f852 100644
--- a/pynest/nest/tests/test_stdp_multiplicity.py
+++ b/pynest/nest/tests/test_stdp_multiplicity.py
@@ -58,23 +58,29 @@ class StdpSpikeMultiplicity(unittest.TestCase):
        delta, since in this case all spikes are at the end of the step, i.e.,
        all spikes have identical times independent of delta.
 
-    2. We choose delta values that are decrease by factors of 2. The
+    2. We choose delta values that are decreased by factors of 2. The
        plasticity rules depend on spike-time differences through
 
+       ::
+
            exp(dT / tau)
 
        where dT is the time between pre- and postsynaptic spikes. We construct
        pre- and postsynaptic spike times so that
 
-          dT = pre_post_shift + m * delta
+       ::
+
+           dT = pre_post_shift + m * delta
 
-       with m * delta < resolution << pre_post_shift. The time-dependence
+       with ``m * delta < resolution << pre_post_shift``. The time-dependence
        of the plasticity rule is therefore to good approximation linear in
        delta.
 
-       We can thus test as follows: Let w_pl be the weight obtained with the
-       plain parrot, and w_ps_j the weight obtained with the precise parrot
-       for delta_j = delta0 / 2^j. Then,
+       We can thus test as follows: Let ``w_pl`` be the weight obtained with the
+       plain parrot, and ``w_ps_j`` the weight obtained with the precise parrot
+       for ``delta_j = delta0 / 2^j``. Then,
+
+       ::
 
          ( w_ps_{j+1} - w_pl ) / ( w_ps_j - w_pl ) ~ 0.5  for all j
 
@@ -221,7 +227,6 @@ def run_protocol(self, pre_post_shift):
         print(post_weights)
         return post_weights
 
-    # XXX: TODO: use ``@pytest.mark.parametrize`` for this
     def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-6):
         """Check that for smaller and smaller timestep, weights obtained from parrot and precise parrot converge.
 
@@ -234,15 +239,16 @@ def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-6):
         w_plain = np.array(post_weights['parrot'])
         w_precise = np.array(post_weights['parrot_ps'])
 
-        assert all(w_plain == w_plain[0]), 'Plain weights should be independent of timestep, but differ!'
+        assert all(w_plain == w_plain[0]), 'Plain weights should be independent of timestep!'
         abs_err = np.abs(w_precise - w_plain)
-        assert abs_err[-1] < max_abs_err
-        assert np.all(np.diff(abs_err) < 1), 'Error should decrease with smaller timestep!'
+        assert abs_err[-1] < max_abs_err, 'Final absolute error is ' + '{0:.2E}'.format(abs_err[-1]) + ' but should be <= ' + '{0:.2E}'.format(max_abs_err)
+        assert np.all(np.diff(abs_err) < 0), 'Error should decrease with smaller timestep!'
 
     def test_stdp_multiplicity(self):
         """Check weight convergence on potentiation and depression.
 
         See also: _test_stdp_multiplicity()."""
+        # XXX: TODO: use ``@pytest.mark.parametrize`` for this
         self._test_stdp_multiplicity(pre_post_shift=10.)    # test potentiation
         self._test_stdp_multiplicity(pre_post_shift=-10.)   # test depression
 

From fb9228d2186f0e539d6c409c4642c0a71ce5a409 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 30 Aug 2021 15:10:38 +0200
Subject: [PATCH 12/19] move synaptic plasticity precise spike timing feature
 to #2035

---
 models/stdp_synapse.h                  | 2 +-
 nestkernel/connector_base_impl.h       | 1 -
 testsuite/pytests/test_stdp_synapse.py | 4 +---
 3 files changed, 2 insertions(+), 5 deletions(-)

diff --git a/models/stdp_synapse.h b/models/stdp_synapse.h
index 64afed317a..a2ad628572 100644
--- a/models/stdp_synapse.h
+++ b/models/stdp_synapse.h
@@ -217,7 +217,7 @@ inline void
 stdp_synapse< targetidentifierT >::send( Event& e, thread t, const CommonSynapseProperties& )
 {
   // synapse STDP depressing/facilitation dynamics
-  const double t_spike = e.get_stamp().get_ms() - e.get_offset();
+  const double t_spike = e.get_stamp().get_ms();
 
   // use accessor functions (inherited from Connection< >) to obtain delay and
   // target
diff --git a/nestkernel/connector_base_impl.h b/nestkernel/connector_base_impl.h
index 21091570fe..45117feee2 100644
--- a/nestkernel/connector_base_impl.h
+++ b/nestkernel/connector_base_impl.h
@@ -47,7 +47,6 @@ Connector< ConnectionT >::send_weight_event( const thread tid,
     wr_e.set_port( e.get_port() );
     wr_e.set_rport( e.get_rport() );
     wr_e.set_stamp( e.get_stamp() );
-    wr_e.set_offset( e.get_offset() );
     wr_e.set_sender( e.get_sender() );
     wr_e.set_sender_node_id( kernel().connection_manager.get_source_node_id( tid, syn_id_, lcid ) );
     wr_e.set_weight( e.get_weight() );
diff --git a/testsuite/pytests/test_stdp_synapse.py b/testsuite/pytests/test_stdp_synapse.py
index a634d90344..3017049b90 100644
--- a/testsuite/pytests/test_stdp_synapse.py
+++ b/testsuite/pytests/test_stdp_synapse.py
@@ -339,9 +339,7 @@ def test_stdp_synapse(self):
         self.init_params()
         for self.dendritic_delay in [1., self.resolution]:
             self.init_params()
-            for self.nest_neuron_model in ["iaf_psc_exp",
-                                           "iaf_cond_exp",
-                                           "iaf_psc_exp_ps"]:
+            for self.nest_neuron_model in ["iaf_psc_exp", "iaf_cond_exp"]:
                 fname_snip = "_[nest_neuron_mdl=" + self.nest_neuron_model + "]"
                 self.do_nest_simulation_and_compare_to_reproduced_weight(fname_snip=fname_snip)
 

From b4921bb282b5ef1261b042ae17f06339f069a74c Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 30 Aug 2021 16:04:15 +0200
Subject: [PATCH 13/19] refactor/clean up stdp synapse unit test

---
 testsuite/pytests/test_stdp_nn_synapses.py | 45 +++++++++-------------
 testsuite/pytests/test_stdp_synapse.py     | 45 ++++++++--------------
 2 files changed, 34 insertions(+), 56 deletions(-)

diff --git a/testsuite/pytests/test_stdp_nn_synapses.py b/testsuite/pytests/test_stdp_nn_synapses.py
index ff14f34c6d..9305274edf 100644
--- a/testsuite/pytests/test_stdp_nn_synapses.py
+++ b/testsuite/pytests/test_stdp_nn_synapses.py
@@ -48,7 +48,6 @@ def setUp(self):
         self.presynaptic_firing_rate = 20.0  # [Hz]
         self.postsynaptic_firing_rate = 20.0  # [Hz]
         self.simulation_duration = 1e+4  # [ms]
-        self.hardcoded_trains_length = 15.  # [ms]
         self.synapse_parameters = {
             "receptor_type": 1,
             "delay": self.resolution,
@@ -66,6 +65,19 @@ def setUp(self):
             "tau_minus": 33.7
         }
 
+        # While the random sequences, fairly long, would supposedly
+        # reveal small differences in the weight change between NEST
+        # and ours, some low-probability events (say, coinciding
+        # spikes) can well not have occured. To generate and
+        # test every possible combination of pre/post precedence, we
+        # append some hardcoded spike sequences:
+        # pre: 1       5 6 7   9    11 12 13
+        # post:  2 3 4       8 9 10    12
+        self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13], dtype=float)
+        self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12], dtype=float)
+        self.hardcoded_trains_length = 5 * self.synapse_parameters["delay"] \
+            + max(np.amax(self.hardcoded_pre_times), np.amax(self.hardcoded_post_times))
+
     def do_nest_simulation_and_compare_to_reproduced_weight(self,
                                                             pairing_scheme):
         synapse_model = "stdp_" + pairing_scheme + "_synapse"
@@ -93,12 +105,10 @@ def do_the_nest_simulation(self):
         nest.ResetKernel()
         nest.SetKernelStatus({'resolution': self.resolution})
 
-        neurons = nest.Create(
+        presynaptic_neuron, postsynaptic_neuron = nest.Create(
             "parrot_neuron",
             2,
             params=self.neuron_parameters)
-        presynaptic_neuron = neurons[0]
-        postsynaptic_neuron = neurons[1]
 
         generators = nest.Create(
             "poisson_generator",
@@ -110,32 +120,13 @@ def do_the_nest_simulation(self):
         presynaptic_generator = generators[0]
         postsynaptic_generator = generators[1]
 
-        # While the random sequences, fairly long, would supposedly
-        # reveal small differences in the weight change between NEST
-        # and ours, some low-probability events (say, coinciding
-        # spikes) can well not have occured. To generate and
-        # test every possible combination of pre/post precedence, we
-        # append some hardcoded spike sequences:
-        # pre: 1       5 6 7   9    11 12 13
-        # post:  2 3 4       8 9 10    12
-        (
-            hardcoded_pre_times,
-            hardcoded_post_times
-        ) = [
-            [
-                self.simulation_duration - self.hardcoded_trains_length + t
-                for t in train
-            ] for train in (
-                (1, 5, 6, 7, 9, 11, 12, 13),
-                (2, 3, 4, 8, 9, 10, 12)
-            )
-        ]
-
         spike_senders = nest.Create(
             "spike_generator",
             2,
-            params=({"spike_times": hardcoded_pre_times},
-                    {"spike_times": hardcoded_post_times})
+            params=({"spike_times": self.hardcoded_pre_times
+                    + self.simulation_duration - self.hardcoded_trains_length},
+                    {"spike_times": self.hardcoded_post_times
+                    + self.simulation_duration - self.hardcoded_trains_length})
         )
         pre_spike_generator = spike_senders[0]
         post_spike_generator = spike_senders[1]
diff --git a/testsuite/pytests/test_stdp_synapse.py b/testsuite/pytests/test_stdp_synapse.py
index 3017049b90..13bf5293e1 100644
--- a/testsuite/pytests/test_stdp_synapse.py
+++ b/testsuite/pytests/test_stdp_synapse.py
@@ -66,9 +66,19 @@ def init_params(self):
         self.neuron_parameters = {
             "tau_minus": self.tau_post,
         }
-        self.hardcoded_trains_length = 15
-        self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13.])
-        self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12.])
+
+        # While the random sequences, fairly long, would supposedly
+        # reveal small differences in the weight change between NEST
+        # and ours, some low-probability events (say, coinciding
+        # spikes) can well not have occured. To generate and
+        # test every possible combination of pre/post precedence, we
+        # append some hardcoded spike sequences:
+        # pre: 1       5 6 7   9    11 12 13
+        # post:  2 3 4       8 9 10    12
+        self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13], dtype=float)
+        self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12], dtype=float)
+        self.hardcoded_trains_length = 5 * self.synapse_parameters["delay"] \
+            + max(np.amax(self.hardcoded_pre_times), np.amax(self.hardcoded_post_times))
 
     def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
         pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest = self.do_the_nest_simulation()
@@ -104,12 +114,10 @@ def do_the_nest_simulation(self):
         nest.ResetKernel()
         nest.SetKernelStatus({'resolution': self.resolution})
 
-        neurons = nest.Create(
+        presynaptic_neuron, postsynaptic_neuron = nest.Create(
             self.nest_neuron_model,
             2,
             params=self.neuron_parameters)
-        presynaptic_neuron = neurons[0]
-        postsynaptic_neuron = neurons[1]
 
         generators = nest.Create(
             "poisson_generator",
@@ -124,32 +132,11 @@ def do_the_nest_simulation(self):
         wr = nest.Create('weight_recorder')
         nest.CopyModel(self.synapse_model, self.synapse_model + "_rec", {"weight_recorder": wr})
 
-        # While the random sequences, fairly long, would supposedly
-        # reveal small differences in the weight change between NEST
-        # and ours, some low-probability events (say, coinciding
-        # spikes) can well not have occured. To generate and
-        # test every possible combination of pre/post precedence, we
-        # append some hardcoded spike sequences:
-        # pre: 1       5 6 7   9    11 12 13
-        # post:  2 3 4       8 9 10    12
-        (
-            hardcoded_pre_times,
-            hardcoded_post_times
-        ) = [
-            [
-                self.simulation_duration - self.hardcoded_trains_length + float(t)
-                for t in train
-            ] for train in (
-                self.hardcoded_pre_times,
-                self.hardcoded_post_times
-            )
-        ]
-
         spike_senders = nest.Create(
             "spike_generator",
             2,
-            params=({"spike_times": hardcoded_pre_times},
-                    {"spike_times": hardcoded_post_times})
+            params=({"spike_times": hardcoded_pre_times + self.simulation_duration - self.hardcoded_trains_length},
+                    {"spike_times": hardcoded_post_times + self.simulation_duration - self.hardcoded_trains_length})
         )
         pre_spike_generator = spike_senders[0]
         post_spike_generator = spike_senders[1]

From ae9f481a856918df301a1f61732e41eca26652c8 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 30 Aug 2021 17:54:54 +0200
Subject: [PATCH 14/19] refactor/clean up stdp synapse unit testing

---
 testsuite/pytests/test_stdp_multiplicity.py | 2 +-
 testsuite/pytests/test_stdp_nn_synapses.py  | 4 ++--
 testsuite/pytests/test_stdp_synapse.py      | 7 +++----
 3 files changed, 6 insertions(+), 7 deletions(-)

diff --git a/testsuite/pytests/test_stdp_multiplicity.py b/testsuite/pytests/test_stdp_multiplicity.py
index 8c3ab0f852..f8b3c8d215 100644
--- a/testsuite/pytests/test_stdp_multiplicity.py
+++ b/testsuite/pytests/test_stdp_multiplicity.py
@@ -227,7 +227,7 @@ def run_protocol(self, pre_post_shift):
         print(post_weights)
         return post_weights
 
-    def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-6):
+    def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-3):
         """Check that for smaller and smaller timestep, weights obtained from parrot and precise parrot converge.
 
         Enforce a maximum allowed absolute error ``max_abs_err`` between the final weights for the smallest timestep
diff --git a/testsuite/pytests/test_stdp_nn_synapses.py b/testsuite/pytests/test_stdp_nn_synapses.py
index 9305274edf..25f9782081 100644
--- a/testsuite/pytests/test_stdp_nn_synapses.py
+++ b/testsuite/pytests/test_stdp_nn_synapses.py
@@ -23,6 +23,7 @@
 # and stdp_nn_restr_synapse in NEST.
 
 import nest
+import numpy as np
 import unittest
 from math import exp
 
@@ -75,8 +76,7 @@ def setUp(self):
         # post:  2 3 4       8 9 10    12
         self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13], dtype=float)
         self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12], dtype=float)
-        self.hardcoded_trains_length = 5 * self.synapse_parameters["delay"] \
-            + max(np.amax(self.hardcoded_pre_times), np.amax(self.hardcoded_post_times))
+        self.hardcoded_trains_length = 2. + max(np.amax(self.hardcoded_pre_times), np.amax(self.hardcoded_post_times))
 
     def do_nest_simulation_and_compare_to_reproduced_weight(self,
                                                             pairing_scheme):
diff --git a/testsuite/pytests/test_stdp_synapse.py b/testsuite/pytests/test_stdp_synapse.py
index 13bf5293e1..d64deb18df 100644
--- a/testsuite/pytests/test_stdp_synapse.py
+++ b/testsuite/pytests/test_stdp_synapse.py
@@ -77,8 +77,7 @@ def init_params(self):
         # post:  2 3 4       8 9 10    12
         self.hardcoded_pre_times = np.array([1, 5, 6, 7, 9, 11, 12, 13], dtype=float)
         self.hardcoded_post_times = np.array([2, 3, 4, 8, 9, 10, 12], dtype=float)
-        self.hardcoded_trains_length = 5 * self.synapse_parameters["delay"] \
-            + max(np.amax(self.hardcoded_pre_times), np.amax(self.hardcoded_post_times))
+        self.hardcoded_trains_length = 2. + max(np.amax(self.hardcoded_pre_times), np.amax(self.hardcoded_post_times))
 
     def do_nest_simulation_and_compare_to_reproduced_weight(self, fname_snip):
         pre_spikes, post_spikes, t_weight_by_nest, weight_by_nest = self.do_the_nest_simulation()
@@ -135,8 +134,8 @@ def do_the_nest_simulation(self):
         spike_senders = nest.Create(
             "spike_generator",
             2,
-            params=({"spike_times": hardcoded_pre_times + self.simulation_duration - self.hardcoded_trains_length},
-                    {"spike_times": hardcoded_post_times + self.simulation_duration - self.hardcoded_trains_length})
+            params=({"spike_times": self.hardcoded_pre_times + self.simulation_duration - self.hardcoded_trains_length},
+                    {"spike_times": self.hardcoded_post_times + self.simulation_duration - self.hardcoded_trains_length})
         )
         pre_spike_generator = spike_senders[0]
         post_spike_generator = spike_senders[1]

From 7fb3e2d1a53eac61233d335c1870a61bf112405e Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 30 Aug 2021 18:22:53 +0200
Subject: [PATCH 15/19] fix pycodestyle

---
 testsuite/pytests/test_stdp_synapse.py | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/testsuite/pytests/test_stdp_synapse.py b/testsuite/pytests/test_stdp_synapse.py
index d64deb18df..4eb99b3b3f 100644
--- a/testsuite/pytests/test_stdp_synapse.py
+++ b/testsuite/pytests/test_stdp_synapse.py
@@ -134,8 +134,10 @@ def do_the_nest_simulation(self):
         spike_senders = nest.Create(
             "spike_generator",
             2,
-            params=({"spike_times": self.hardcoded_pre_times + self.simulation_duration - self.hardcoded_trains_length},
-                    {"spike_times": self.hardcoded_post_times + self.simulation_duration - self.hardcoded_trains_length})
+            params=({"spike_times": self.hardcoded_pre_times
+                     + self.simulation_duration - self.hardcoded_trains_length},
+                    {"spike_times": self.hardcoded_post_times
+                     + self.simulation_duration - self.hardcoded_trains_length})
         )
         pre_spike_generator = spike_senders[0]
         post_spike_generator = spike_senders[1]

From f4182fdb13dd9e5ede778607656b8444201ea514 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 4 Oct 2021 14:37:01 +0200
Subject: [PATCH 16/19] clean up unit tests following unittest/pytest migration

---
 testsuite/pytests/test_jonke_synapse.py     | 27 ++++-----------
 testsuite/pytests/test_stdp_multiplicity.py | 38 ++++-----------------
 testsuite/pytests/test_stdp_nn_synapses.py  | 35 ++++++-------------
 testsuite/pytests/test_stdp_synapse.py      | 20 ++---------
 testsuite/pytests/test_visualization.py     | 20 +++--------
 5 files changed, 31 insertions(+), 109 deletions(-)

diff --git a/testsuite/pytests/test_jonke_synapse.py b/testsuite/pytests/test_jonke_synapse.py
index 686f53a801..67c6c6eb3d 100644
--- a/testsuite/pytests/test_jonke_synapse.py
+++ b/testsuite/pytests/test_jonke_synapse.py
@@ -23,13 +23,12 @@
 Test functionality of the Tetzlaff stdp synapse
 """
 
-import unittest
 import nest
 import numpy as np
 
 
 @nest.ll_api.check_stack
-class JonkeSynapseTest(unittest.TestCase):
+class TestJonkeSynapse:
     """
     Test the weight change by STDP.
     The test is performed by generating two Poisson spike trains,
@@ -73,12 +72,12 @@ def test_weight_drift(self):
         weight_reproduced_independently = self.reproduce_weight_drift(
             pre_spikes, post_spikes,
             self.synapse_parameters["weight"])
-        self.assertAlmostEqual(
+        np.testing.assert_almost_equal(
             weight_reproduced_independently,
             weight_by_nest,
-            msg=f"{self.synapse_parameters['synapse_model']} test:\n" +
-                f"Resulting synaptic weight {weight_by_nest} " +
-                f"differs from expected {weight_reproduced_independently}")
+            err_msg=f"{self.synapse_parameters['synapse_model']} test:\n" +
+                    f"Resulting synaptic weight {weight_by_nest} " +
+                    f"differs from expected {weight_reproduced_independently}")
 
     def do_the_nest_simulation(self):
         """
@@ -111,7 +110,7 @@ def do_the_nest_simulation(self):
         # reveal small differences in the weight change between NEST
         # and ours, some low-probability events (say, coinciding
         # spikes) can well not have occurred. To generate and
-        # test every possible combination of pre/post precedence, we
+        # test every possible combination of pre/post order, we
         # append some hardcoded spike sequences:
         # pre: 1       5 6 7   9    11 12 13
         # post:  2 3 4       8 9 10    12
@@ -244,17 +243,3 @@ def depress(self, _delta_t, weight, Kminus):
         if weight < 0:
             weight = 0
         return weight
-
-
-def suite():
-    suite = unittest.TestLoader().loadTestsFromTestCase(JonkeSynapseTest)
-    return unittest.TestSuite([suite])
-
-
-def run():
-    runner = unittest.TextTestRunner(verbosity=2)
-    runner.run(suite())
-
-
-if __name__ == "__main__":
-    run()
diff --git a/testsuite/pytests/test_stdp_multiplicity.py b/testsuite/pytests/test_stdp_multiplicity.py
index 19d25935bf..7cecda0357 100644
--- a/testsuite/pytests/test_stdp_multiplicity.py
+++ b/testsuite/pytests/test_stdp_multiplicity.py
@@ -21,10 +21,10 @@
 
 # This script tests the parrot_neuron in NEST.
 
-import nest
-import unittest
 import math
+import nest
 import numpy as np
+import pytest
 
 try:
     import matplotlib as mpl
@@ -35,7 +35,7 @@
 
 
 @nest.ll_api.check_stack
-class StdpSpikeMultiplicity(unittest.TestCase):
+class TestStdpSpikeMultiplicity:
     """
     Test correct handling of spike multiplicity in STDP.
 
@@ -207,7 +207,6 @@ def run_protocol(self, pre_post_shift):
             post_weights['parrot_ps'].append(w_post_ps)
 
         if DEBUG_PLOTS:
-            import datetime
             fig, ax = plt.subplots(nrows=2)
             fig.suptitle("Final obtained weights")
             ax[0].plot(post_weights["parrot"], marker="o", label="parrot")
@@ -222,12 +221,14 @@ def run_protocol(self, pre_post_shift):
             for _ax in ax:
                 _ax.grid(True)
                 _ax.legend()
-            plt.savefig("/tmp/test_stdp_multiplicity" + str(datetime.datetime.utcnow()) + ".png")
+            plt.savefig("/tmp/test_stdp_multiplicity.png")
             plt.close(fig)
         print(post_weights)
         return post_weights
 
-    def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-3):
+    @pytest.mark.parametrize("pre_post_shift", [10.,    # test potentiation
+                                                -10.])    # test depression
+    def test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-3):
         """Check that for smaller and smaller timestep, weights obtained from parrot and precise parrot converge.
 
         Enforce a maximum allowed absolute error ``max_abs_err`` between the final weights for the smallest timestep
@@ -243,28 +244,3 @@ def _test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-3):
         abs_err = np.abs(w_precise - w_plain)
         assert abs_err[-1] < max_abs_err, 'Final absolute error is ' + '{0:.2E}'.format(abs_err[-1]) + ' but should be <= ' + '{0:.2E}'.format(max_abs_err)
         assert np.all(np.diff(abs_err) < 0), 'Error should decrease with smaller timestep!'
-
-    def test_stdp_multiplicity(self):
-        """Check weight convergence on potentiation and depression.
-
-        See also: _test_stdp_multiplicity()."""
-        # XXX: TODO: use ``@pytest.mark.parametrize`` for this
-        self._test_stdp_multiplicity(pre_post_shift=10.)    # test potentiation
-        self._test_stdp_multiplicity(pre_post_shift=-10.)   # test depression
-
-
-def suite():
-
-    # makeSuite is sort of obsolete http://bugs.python.org/issue2721
-    # using loadTestsFromTestCase instead.
-    suite = unittest.TestLoader().loadTestsFromTestCase(StdpSpikeMultiplicity)
-    return unittest.TestSuite([suite])
-
-
-def run():
-    runner = unittest.TextTestRunner(verbosity=2)
-    runner.run(suite())
-
-
-if __name__ == "__main__":
-    run()
diff --git a/testsuite/pytests/test_stdp_nn_synapses.py b/testsuite/pytests/test_stdp_nn_synapses.py
index d23b628411..cdc3ed3733 100644
--- a/testsuite/pytests/test_stdp_nn_synapses.py
+++ b/testsuite/pytests/test_stdp_nn_synapses.py
@@ -24,12 +24,13 @@
 
 import nest
 import numpy as np
-import unittest
+import pytest
+
 from math import exp
 
 
 @nest.ll_api.check_stack
-class STDPNNSynapsesTest(unittest.TestCase):
+class TestSTDPNNSynapses:
     """
     Test the weight change by STDP
     with three nearest-neighbour spike pairing schemes.
@@ -44,6 +45,7 @@ class STDPNNSynapsesTest(unittest.TestCase):
     Instead, it directly iterates through the spike history.
     """
 
+    @pytest.fixture(autouse=True)
     def setUp(self):
         self.resolution = 0.1  # [ms]
         self.presynaptic_firing_rate = 20.0  # [Hz]
@@ -70,7 +72,7 @@ def setUp(self):
         # reveal small differences in the weight change between NEST
         # and ours, some low-probability events (say, coinciding
         # spikes) can well not have occured. To generate and
-        # test every possible combination of pre/post precedence, we
+        # test every possible combination of pre/post order, we
         # append some hardcoded spike sequences:
         # pre: 1       5 6 7   9    11 12 13
         # post:  2 3 4       8 9 10    12
@@ -87,13 +89,13 @@ def do_nest_simulation_and_compare_to_reproduced_weight(self,
         weight_reproduced_independently = self.reproduce_weight_drift(
             pre_spikes, post_spikes,
             self.synapse_parameters["weight"])
-        self.assertAlmostEqual(
+        np.testing.assert_almost_equal(
             weight_reproduced_independently,
             weight_by_nest,
-            msg=synapse_model + " test: "
-                                "Resulting synaptic weight %e "
-                                "differs from expected %e" % (
-                                    weight_by_nest, weight_reproduced_independently))
+            err_msg=synapse_model + " test: "
+                                    "Resulting synaptic weight %e "
+                                    "differs from expected %e" % (
+                                       weight_by_nest, weight_reproduced_independently))
 
     def do_the_nest_simulation(self):
         """
@@ -276,20 +278,3 @@ def test_nn_pre_centered_synapse(self):
 
     def test_nn_restr_synapse(self):
         self.do_nest_simulation_and_compare_to_reproduced_weight("nn_restr")
-
-
-def suite():
-
-    # makeSuite is sort of obsolete http://bugs.python.org/issue2721
-    # using loadTestsFromTestCase instead.
-    suite = unittest.TestLoader().loadTestsFromTestCase(STDPNNSynapsesTest)
-    return unittest.TestSuite([suite])
-
-
-def run():
-    runner = unittest.TextTestRunner(verbosity=2)
-    runner.run(suite())
-
-
-if __name__ == "__main__":
-    run()
diff --git a/testsuite/pytests/test_stdp_synapse.py b/testsuite/pytests/test_stdp_synapse.py
index 4eb99b3b3f..bd80814020 100644
--- a/testsuite/pytests/test_stdp_synapse.py
+++ b/testsuite/pytests/test_stdp_synapse.py
@@ -20,7 +20,7 @@
 # along with NEST.  If not, see <http://www.gnu.org/licenses/>.
 
 import nest
-import unittest
+import pytest
 from math import exp
 import numpy as np
 
@@ -33,7 +33,7 @@
 
 
 @nest.ll_api.check_stack
-class STDPSynapseTest(unittest.TestCase):
+class TestSTDPSynapse:
     """
     Compare the STDP synaptic plasticity model against a self-contained Python reference.
 
@@ -71,7 +71,7 @@ def init_params(self):
         # reveal small differences in the weight change between NEST
         # and ours, some low-probability events (say, coinciding
         # spikes) can well not have occured. To generate and
-        # test every possible combination of pre/post precedence, we
+        # test every possible combination of pre/post order, we
         # append some hardcoded spike sequences:
         # pre: 1       5 6 7   9    11 12 13
         # post:  2 3 4       8 9 10    12
@@ -330,17 +330,3 @@ def test_stdp_synapse(self):
             for self.nest_neuron_model in ["iaf_psc_exp", "iaf_cond_exp"]:
                 fname_snip = "_[nest_neuron_mdl=" + self.nest_neuron_model + "]"
                 self.do_nest_simulation_and_compare_to_reproduced_weight(fname_snip=fname_snip)
-
-
-def suite():
-    suite = unittest.TestLoader().loadTestsFromTestCase(STDPSynapseTest)
-    return unittest.TestSuite([suite])
-
-
-def run():
-    runner = unittest.TextTestRunner(verbosity=2)
-    runner.run(suite())
-
-
-if __name__ == "__main__":
-    run()
diff --git a/testsuite/pytests/test_visualization.py b/testsuite/pytests/test_visualization.py
index fe4401dcf0..5790412a85 100644
--- a/testsuite/pytests/test_visualization.py
+++ b/testsuite/pytests/test_visualization.py
@@ -24,7 +24,7 @@
 """
 
 import os
-import unittest
+import pytest
 import nest
 import numpy as np
 
@@ -50,7 +50,7 @@
     HAVE_PANDAS = False
 
 
-class VisualizationTestCase(unittest.TestCase):
+class TestVisualization:
     def nest_tmpdir(self):
         """Returns temp dir path from environment, current dir otherwise."""
         if 'NEST_DATA_PATH' in os.environ:
@@ -66,7 +66,7 @@ def tearDown(self):
             # Cleanup temporary datafiles
             os.remove(filename)
 
-    @unittest.skipIf(not HAVE_PYDOT, 'pydot not found')
+    @pytest.mark.skipif(not HAVE_PYDOT, reason='pydot not found')
     def test_plot_network(self):
         """Test plot_network"""
         import nest.visualization as nvis
@@ -92,7 +92,7 @@ def voltage_trace_verify(self, device):
             self.assertTrue(all(np.isclose(ref_vm, y_data)))
         plt.close(ax.get_figure())
 
-    @unittest.skipIf(not PLOTTING_POSSIBLE, 'Plotting impossible because matplotlib or display missing')
+    @pytest.mark.skipif(not PLOTTING_POSSIBLE, reason='Plotting impossible because matplotlib or display missing')
     def test_voltage_trace_from_device(self):
         """Test voltage_trace from device"""
         import nest.voltage_trace
@@ -150,7 +150,7 @@ def spike_recorder_raster_verify(self, sr_ref):
         self.assertEqual(x_data.shape, sr_ref.shape)
         self.assertTrue(all(np.isclose(x_data, sr_ref)))
 
-    @unittest.skipIf(not PLOTTING_POSSIBLE, 'Plotting impossible because matplotlib or display missing')
+    @pytest.mark.skipif(not PLOTTING_POSSIBLE, reason='Plotting impossible because matplotlib or display missing')
     def test_raster_plot(self):
         """Test raster_plot"""
         import nest.raster_plot
@@ -193,13 +193,3 @@ def test_raster_plot(self):
         self.assertTrue(np.all(times_30_to_40_extracted[:, 1] >= 30.))
         self.assertTrue(np.all(times_30_to_40_extracted[:, 1] < 40.))
         self.assertEqual(len(source_2_extracted), 0)
-
-
-def suite():
-    suite = unittest.makeSuite(VisualizationTestCase, 'test')
-    return suite
-
-
-if __name__ == "__main__":
-    runner = unittest.TextTestRunner(verbosity=2)
-    runner.run(suite())

From 14ad08c517316d50e64221922af5d5fdd54ce20e Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 4 Oct 2021 15:35:20 +0200
Subject: [PATCH 17/19] fix pycodestyle

---
 testsuite/pytests/test_stdp_multiplicity.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/testsuite/pytests/test_stdp_multiplicity.py b/testsuite/pytests/test_stdp_multiplicity.py
index 7cecda0357..80ade69bef 100644
--- a/testsuite/pytests/test_stdp_multiplicity.py
+++ b/testsuite/pytests/test_stdp_multiplicity.py
@@ -242,5 +242,6 @@ def test_stdp_multiplicity(self, pre_post_shift, max_abs_err=1E-3):
 
         assert all(w_plain == w_plain[0]), 'Plain weights should be independent of timestep!'
         abs_err = np.abs(w_precise - w_plain)
-        assert abs_err[-1] < max_abs_err, 'Final absolute error is ' + '{0:.2E}'.format(abs_err[-1]) + ' but should be <= ' + '{0:.2E}'.format(max_abs_err)
+        assert abs_err[-1] < max_abs_err, 'Final absolute error is ' + '{0:.2E}'.format(abs_err[-1]) \
+                                          + ' but should be <= ' + '{0:.2E}'.format(max_abs_err)
         assert np.all(np.diff(abs_err) < 0), 'Error should decrease with smaller timestep!'

From c03797eb584cf7086aa9f0dff3fa9d54f646ae27 Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 4 Oct 2021 15:52:14 +0200
Subject: [PATCH 18/19] clean up unit tests following unittest/pytest migration

---
 testsuite/pytests/test_visualization.py | 31 +++++++++++++------------
 1 file changed, 16 insertions(+), 15 deletions(-)

diff --git a/testsuite/pytests/test_visualization.py b/testsuite/pytests/test_visualization.py
index 5790412a85..69ba6d0f25 100644
--- a/testsuite/pytests/test_visualization.py
+++ b/testsuite/pytests/test_visualization.py
@@ -23,10 +23,10 @@
 Tests for visualization functions.
 """
 
-import os
-import pytest
 import nest
 import numpy as np
+import os
+import pytest
 
 try:
     import matplotlib.pyplot as plt
@@ -58,10 +58,11 @@ def nest_tmpdir(self):
         else:
             return '.'
 
+    @pytest.fixture(autouse=True)
     def setUp(self):
         self.filenames = []
-
-    def tearDown(self):
+        yield
+        # fixture teardown code below
         for filename in self.filenames:
             # Cleanup temporary datafiles
             os.remove(filename)
@@ -78,18 +79,18 @@ def test_plot_network(self):
         filename = os.path.join(self.nest_tmpdir(), 'network_plot.png')
         self.filenames.append(filename)
         nvis.plot_network(sources + targets, filename)
-        self.assertTrue(os.path.isfile(filename), 'Plot was not created or not saved')
+        assert os.path.isfile(filename), 'Plot was not created or not saved'
 
     def voltage_trace_verify(self, device):
-        self.assertIsNotNone(plt._pylab_helpers.Gcf.get_active(), 'No active figure')
+        assert plt._pylab_helpers.Gcf.get_active() is not None, 'No active figure'
         ax = plt.gca()
         vm = device.get('events', 'V_m')
         for ref_vm, line in zip((vm[::2], vm[1::2]), ax.lines):
             x_data, y_data = line.get_data()
             # Check that times are correct
-            self.assertEqual(list(x_data), list(np.unique(device.get('events', 'times'))))
+            assert list(x_data) == list(np.unique(device.get('events', 'times')))
             # Check that voltmeter data corresponds to the lines in the plot
-            self.assertTrue(all(np.isclose(ref_vm, y_data)))
+            assert all(np.isclose(ref_vm, y_data))
         plt.close(ax.get_figure())
 
     @pytest.mark.skipif(not PLOTTING_POSSIBLE, reason='Plotting impossible because matplotlib or display missing')
@@ -141,14 +142,14 @@ def spike_recorder_data_setup(self, to_file=False):
             return sr
 
     def spike_recorder_raster_verify(self, sr_ref):
-        self.assertIsNotNone(plt._pylab_helpers.Gcf.get_active(), 'No active figure')
+        assert plt._pylab_helpers.Gcf.get_active() is not None, 'No active figure'
         fig = plt.gcf()
         axs = fig.get_axes()
         x_data, y_data = axs[0].lines[0].get_data()
         plt.close(fig)
         # Have to use isclose() because of round-off errors
-        self.assertEqual(x_data.shape, sr_ref.shape)
-        self.assertTrue(all(np.isclose(x_data, sr_ref)))
+        assert x_data.shape == sr_ref.shape
+        assert all(np.isclose(x_data, sr_ref))
 
     @pytest.mark.skipif(not PLOTTING_POSSIBLE, reason='Plotting impossible because matplotlib or display missing')
     def test_raster_plot(self):
@@ -189,7 +190,7 @@ def test_raster_plot(self):
         all_extracted = nest.raster_plot.extract_events(data)
         times_30_to_40_extracted = nest.raster_plot.extract_events(data, time=[30., 40.], sel=[3])
         source_2_extracted = nest.raster_plot.extract_events(data, sel=[2])
-        self.assertTrue(np.array_equal(all_extracted, data))
-        self.assertTrue(np.all(times_30_to_40_extracted[:, 1] >= 30.))
-        self.assertTrue(np.all(times_30_to_40_extracted[:, 1] < 40.))
-        self.assertEqual(len(source_2_extracted), 0)
+        assert np.array_equal(all_extracted, data)
+        assert np.all(times_30_to_40_extracted[:, 1] >= 30.)
+        assert np.all(times_30_to_40_extracted[:, 1] < 40.)
+        assert len(source_2_extracted) == 0

From 7c7d949e9272a9556f96200a0ca601a9d01dd5da Mon Sep 17 00:00:00 2001
From: "C.A.P. Linssen" <charl@turingbirds.com>
Date: Mon, 18 Oct 2021 16:21:04 +0200
Subject: [PATCH 19/19] clean up based on review comments

---
 testsuite/pytests/test_stdp_synapse.py | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/testsuite/pytests/test_stdp_synapse.py b/testsuite/pytests/test_stdp_synapse.py
index bd80814020..22bcbefe47 100644
--- a/testsuite/pytests/test_stdp_synapse.py
+++ b/testsuite/pytests/test_stdp_synapse.py
@@ -235,13 +235,14 @@ def Kpost_at_time(t, spikes, init=1., inclusive=True):
                 idx_next_post_spike = np.where((post_spikes_delayed - t) > 0)[0][0]
                 t_next_post_spike = post_spikes_delayed[idx_next_post_spike]
 
-            if idx_next_post_spike >= 0 and t_next_post_spike < t_next_pre_spike:
+            if idx_next_pre_spike >= 0 and idx_next_post_spike >= 0 and t_next_post_spike < t_next_pre_spike:
                 handle_post_spike = True
                 handle_pre_spike = False
-            elif idx_next_pre_spike >= 0 and t_next_post_spike > t_next_pre_spike:
+            elif idx_next_pre_spike >= 0 and idx_next_post_spike >= 0 and t_next_post_spike > t_next_pre_spike:
                 handle_post_spike = False
                 handle_pre_spike = True
             else:
+                # simultaneous spikes (both true) or no more spikes to process (both false)
                 handle_post_spike = idx_next_post_spike >= 0
                 handle_pre_spike = idx_next_pre_spike >= 0
 
@@ -329,4 +330,5 @@ def test_stdp_synapse(self):
             self.init_params()
             for self.nest_neuron_model in ["iaf_psc_exp", "iaf_cond_exp"]:
                 fname_snip = "_[nest_neuron_mdl=" + self.nest_neuron_model + "]"
+                fname_snip += "_[dend_delay=" + str(self.dendritic_delay) + "]"
                 self.do_nest_simulation_and_compare_to_reproduced_weight(fname_snip=fname_snip)