Code: Select all
NEURON {
THREADSAFE
ARTIFICIAL_CELL WilsonCowanStim
RANGE internal_dt
RANGE etau, itau
RANGE eI12, eIscl
RANGE iI12, iIscl
RANGE wee, wei, wie, wii
RANGE e, i
RANGE e_inf, i_inf, e_in_tau, i_in_tau
RANGE e_in, i_in
:Generator
BBCOREPOINTER donotuse
RANGE ratescale
}
PARAMETER {
internal_dt = 0.2
etau = 4.0
itau = 12.0
eI12 = 1.3
eIscl = 4.
iI12 = 2.
iIscl = 3.7
wee = 16.
wei = 15.
wie = -12.
wii = -3.
e_inf = 0.
i_inf = 0.
e_in_tau = 2.
i_in_tau = 10.
ratescale = 200.
}
ASSIGNED {
e
i
eS0
iS0
e_in
i_in
: Generator
cspike
event (ms)
donotuse
}
INITIAL {
e = 0
i = 0
e_in = e_inf
i_in = i_inf
eS0 = S(0,eI12,eIscl)
iS0 = S(0,iI12,iIscl)
event = invl(e*ratescale)
cspike = 10
net_send(internal_dt, 2)
net_send(event , cspike)
}
NET_RECEIVE (w) {
LOCAL evt
if (flag == 0) { : external event
if (w > 0) {
e_in = e_in + w
} else {
i_in = i_in - w
}
}
if (flag == 2) { : compute next time stape
update()
net_send(internal_dt, 2)
evt = t + invl(e*ratescale)
if (evt < event){
event = evt
cspike = cspike + 1
if (cspike > 10000){
cspike = 10
}
net_send(event-t, cspike)
:DB>>
:printf("Send t=%g, e=%g, x=%g @ %g\n",t,e, cspike,event)
:<<DB
}
} else {
:DB>>
:printf("Event t=%g, f=%g:c=%g\n",t,flag,cspike)
:<<DB
:net_event(t)
if (flag == cspike){
net_event(t)
event = t + invl(e*ratescale)
cspike = cspike + 1
if (cspike > 10000){
cspike = 10
}
net_send(event-t, cspike)
:DB>>
:printf("Send t=%g, e=%g, x=%g @ %g\n",t,e, cspike,event)
:<<DB
}
}
}
PROCEDURE update(){
LOCAL Fe, Fi, ze, zi
Fe = Se(wee*e+wie*i+e_in)
Fi = Si(wei*e+wii*i+i_in)
ze = 1+Fe
zi = 1+Fi
e = Fe/ze + (e-Fe/ze)*exp(-ze*internal_dt/etau)
i = Fi/zi + (i-Fi/zi)*exp(-zi*internal_dt/itau)
e_in = e_inf + (e_in - e_inf)*exp(-internal_dt/e_in_tau)
i_in = i_inf + (i_in - i_inf)*exp(-internal_dt/i_in_tau)
}
FUNCTION S(x,a,b){
S=1./(1.+exp(-a*(x-b)))
}
FUNCTION Se(x){
Se = S(x,eI12,eIscl) - eS0
}
FUNCTION Si(x){
Si = S(x,iI12,iIscl) - iS0
}
VERBATIM
#if NRNBBCORE /* running in CoreNEURON */
#define IFNEWSTYLE(arg) arg
#else /* running in NEURON */
/*
1 means noiseFromRandom was called when _ran_compat was previously 0 .
2 means noiseFromRandom123 was called when _ran_compat was previously 0.
*/
static int _ran_compat; /* specifies the noise style for all instances */
#define IFNEWSTYLE(arg) if(_ran_compat == 2) { arg }
#endif /* running in NEURON */
ENDVERBATIM
:backward compatibility
PROCEDURE seed(x) {
VERBATIM
#if !NRNBBCORE
ENDVERBATIM
set_seed(x)
VERBATIM
#endif
ENDVERBATIM
}
FUNCTION invl(rate) (ms) {
LOCAL mean
if (rate <= 0.){
mean = 1e19
} else {
mean = 1000/rate
}
if (mean <= 0.) {
mean = .1 (ms)
}
invl = mean*erand()
}
FUNCTION erand() {
VERBATIM
if (_p_donotuse) {
/*
:Supports separate independent but reproducible streams for
: each instance. However, the corresponding hoc Random
: distribution MUST be set to Random.negexp(1)
*/
#if !NRNBBCORE
if (_ran_compat == 2) {
_lerand = nrnran123_negexp((nrnran123_State*)_p_donotuse);
}else{
_lerand = nrn_random_pick(_p_donotuse);
}
#else
_lerand = nrnran123_negexp((nrnran123_State*)_p_donotuse);
#endif
return _lerand;
}else{
#if NRNBBCORE
assert(0);
#else
/*
: the old standby. Cannot use if reproducible parallel sim
: independent of nhost or which host this instance is on
: is desired, since each instance on this cpu draws from
: the same stream
*/
#endif
}
#if !NRNBBCORE
ENDVERBATIM
erand = exprand(1)
VERBATIM
#endif
ENDVERBATIM
}
PROCEDURE noiseFromRandom() {
VERBATIM
#if !NRNBBCORE
{
void** pv = (void**)(&_p_donotuse);
if (_ran_compat == 2) {
fprintf(stderr, "NetStim.noiseFromRandom123 was previously called\n");
assert(0);
}
_ran_compat = 1;
if (ifarg(1)) {
*pv = nrn_random_arg(1);
}else{
*pv = (void*)0;
}
}
#endif
ENDVERBATIM
}
PROCEDURE noiseFromRandom123() {
VERBATIM
#if !NRNBBCORE
{
nrnran123_State** pv = (nrnran123_State**)(&_p_donotuse);
if (_ran_compat == 1) {
fprintf(stderr, "NetStim.noiseFromRandom was previously called\n");
assert(0);
}
_ran_compat = 2;
if (*pv) {
nrnran123_deletestream(*pv);
*pv = (nrnran123_State*)0;
}
if (ifarg(3)) {
*pv = nrnran123_newstream3((uint32_t)*getarg(1), (uint32_t)*getarg(2), (uint32_t)*getarg(3));
}else if (ifarg(2)) {
*pv = nrnran123_newstream((uint32_t)*getarg(1), (uint32_t)*getarg(2));
}
}
#endif
ENDVERBATIM
}
DESTRUCTOR {
VERBATIM
if (_p_donotuse) {
#if NRNBBCORE
{ /* but note that mod2c does not translate DESTRUCTOR */
#else
if (_ran_compat == 2) {
#endif
nrnran123_State** pv = (nrnran123_State**)(&_p_donotuse);
nrnran123_deletestream(*pv);
*pv = (nrnran123_State*)0;
}
}
ENDVERBATIM
}
VERBATIM
static void bbcore_write(double* x, int* d, int* xx, int *offset, _threadargsproto_) {
/* error if using the legacy scop_exprand */
if (!_p_donotuse) {
fprintf(stderr, "NetStim: cannot use the legacy scop_negexp generator for the random stream.\n");
assert(0);
}
if (d) {
char which;
uint32_t* di = ((uint32_t*)d) + *offset;
#if !NRNBBCORE
if (_ran_compat == 1) {
void** pv = (void**)(&_p_donotuse);
/* error if not using Random123 generator */
if (!nrn_random_isran123(*pv, di, di+1, di+2)) {
fprintf(stderr, "NetStim: Random123 generator is required\n");
assert(0);
}
nrn_random123_getseq(*pv, di+3, &which);
di[4] = (int)which;
}else{
#else
{
#endif
nrnran123_State** pv = (nrnran123_State**)(&_p_donotuse);
nrnran123_getids3(*pv, di, di+1, di+2);
nrnran123_getseq(*pv, di+3, &which);
di[4] = (int)which;
#if NRNBBCORE
/* CORENeuron does not call DESTRUCTOR so... */
nrnran123_deletestream(*pv);
*pv = (nrnran123_State*)0;
#endif
}
/*printf("Netstim bbcore_write %d %d %d\n", di[0], di[1], di[3]);*/
}
*offset += 5;
}
static void bbcore_read(double* x, int* d, int* xx, int* offset, _threadargsproto_) {
/* Generally, CoreNEURON, in the context of psolve, begins with
an empty model so this call takes place in the context of a freshly
created instance and _p_donotuse is not NULL.
However, this function
is also now called from NEURON at the end of coreneuron psolve
in order to transfer back the nrnran123 sequence state. That
allows continuation with a subsequent psolve within NEURON or
properly transfer back to CoreNEURON if we continue the psolve
there. So now, extra logic is needed for this call to work in
a NEURON context.
*/
uint32_t* di = ((uint32_t*)d) + *offset;
#if NRNBBCORE
nrnran123_State** pv = (nrnran123_State**)(&_p_donotuse);
assert(!_p_donotuse);
*pv = nrnran123_newstream3(di[0], di[1], di[2]);
nrnran123_setseq(*pv, di[3], (char)di[4]);
#else
uint32_t id1, id2, id3;
assert(_p_donotuse);
if (_ran_compat == 1) { /* Hoc Random.Random123 */
void** pv = (void**)(&_p_donotuse);
int b = nrn_random_isran123(*pv, &id1, &id2, &id3);
assert(b);
nrn_random123_setseq(*pv, di[3], (char)di[4]);
}else{
assert(_ran_compat == 2);
nrnran123_State** pv = (nrnran123_State**)(&_p_donotuse);
nrnran123_getids3(*pv, &id1, &id2, &id3);
nrnran123_setseq(*pv, di[3], (char)di[4]);
}
/* Random123 on NEURON side has same ids as on CoreNEURON side */
assert(di[0] == id1 && di[1] == id2 && di[2] == id3);
#endif
*offset += 5;
}
ENDVERBATIM
Code: Select all
from numpy import *
from neuron import h
from matplotlib.pyplot import *
wcs = h.WilsonCowanStim()
wcs.ratescale = 200
wcs.e_inf = 1.5
rds = h.Random()
rds.negexp(1) # set random # generator using negexp(1) - avg interval in NetStim
sead = random.randint(0,32562)
rds.MCellRan4(sead,sead) # seeds are in order, shouldn't matter
wcs.noiseFromRandom(rds) # use random # generator for this NetStim
trec = h.Vector()
trec.record(h._ref_t)
erec = h.Vector()
erec.record(wcs._ref_e)
irec = h.Vector()
irec.record(wcs._ref_i)
srec = h.Vector()
ncout = h.NetCon(wcs, None)
ncout.record(srec)
hpc = h.ParallelContext()
hpc.nthread(2)
h.finitialize()
h.fcurrent()
h.frecord_init()
h.dt = 0.1
while h.t < 1000 :
h.fadvance()
plot(array(trec),array(erec),'-')
plot(array(trec),array(irec),'-')
srec = array(srec)
plot(srec,ones(srec.shape[0])*0.5,'k|')
show()
NEURON: Section access unspecified
near line 0
objref hoc_obj_[2]
^
Vector[0].record(...)
Traceback (most recent call last):
File "/media/rth/rth-core/simulations/exampls/neuron-controledP/WilsonCowanStim_test.py", line 18, in <module>
trec.record(h._ref_t)
RuntimeError: hocobj_call error
The error message changes if I create a "fake" section and associate everything with it.
Code: Select all
from numpy import *
from neuron import h
from matplotlib.pyplot import *
sec = h.Section()
wcs = h.WilsonCowanStim(sec=sec)
wcs.ratescale = 200
wcs.e_inf = 1.5
rds = h.Random(sec=sec)
rds.negexp(1) # set random # generator using negexp(1) - avg interval in NetStim
sead = random.randint(0,32562)
rds.MCellRan4(sead,sead) # seeds are in order, shouldn't matter
wcs.noiseFromRandom(rds) # use random # generator for this NetStim
trec = h.Vector()
trec.record(h._ref_t,sec=sec)
erec = h.Vector()
erec.record(wcs._ref_e,sec=sec)
irec = h.Vector()
irec.record(wcs._ref_i,sec=sec)
srec = h.Vector()
ncout = h.NetCon(wcs, None)
ncout.record(srec)
hpc = h.ParallelContext()
hpc.nthread(2)
h.finitialize()
h.fcurrent()
h.frecord_init()
h.dt = 0.1
while h.t < 1000 :
h.fadvance()
plot(array(trec),array(erec),'-')
plot(array(trec),array(irec),'-')
srec = array(srec)
plot(srec,ones(srec.shape[0])*0.5,'k|')
show()
NEURON: We were unable to associate a PlayRecord item with a thread
near line 0
objref hoc_obj_[2]
^
finitialize()
Traceback (most recent call last):
File "/media/rth/rth-core/simulations/exampls/neuron-controledP/WilsonCowanStim_test.py", line 31, in <module>
h.finitialize()
RuntimeError: hocobj_call error
The error above does not appear if ParallelContext and multithreading are commented out.
Question:
How to record variables in ARTIFICIAL_CELL, which is technically not associated with any section of the compartment-based model?