Well it took almost no time to implement equation (6), a single column model, from Jansen and Rit 1995
from matplotlib.pyplot import *
from brian2 import *
defaultclock.dt = 1.*ms
a,b = 100, 50
e0 = 2.5
v0 = 6
r = 0.56
A,B,C = 3.25, 22., 135
P = 120
equs = """
dy1/dt = y4/second : 1
dy4/dt = a*(A*(p+{C2}*Se)-2*y4 -y1*a)/second : 1
dy0/dt = y3/second : 1
dy3/dt = a*(A *Sp -2*y3 -y0*a)/second : 1
dy2/dt = y5/second : 1
dy5/dt = b*(B *{C4}*Si -2*y5 -y2*b)/second : 1
p = P*int(t>1*second)*int(t<2*second) : 1
Se = 2*e0/(1+exp(r*(v0 - {C1}*y0 ))) : 1
Sp = 2*e0/(1+exp(r*(v0 - (y1-y2) ))) : 1
Si = 2*e0/(1+exp(r*(v0 - {C3}*y0 ))) : 1
""".format(
C1 = C,
C2 = 0.8 *C,
C3 = 0.25*C,
C4 = 0.25*C
)
n = NeuronGroup(1,equs,method='euler')
sm = StateMonitor(n,['y4','y1','y3','y0','y5','y2'],record=True)
run(3*second,report='text')
plot(sm.t/second, sm[0].y1)
plot(sm.t/second, sm[0].y0)
plot(sm.t/second, sm[0].y2)
show()