VecStim in parallel environment

General issues of interest both for network and
individual cell parallelization.

Moderator: hines

Post Reply
rth
Posts: 41
Joined: Thu Jun 21, 2012 4:47 pm

VecStim in parallel environment

Post by rth »

Hi all,

I'm writing a code for network which is activated by the input spike trains recorded experimentally or from another model. I load all presynaptic spike events into a vectors and then play these vectors by VecStim-s. Everything worked perfect until I ran simulation on parallel MPI cluster. Surprisingly events generated by VecStim were received only by synapses at same host, other hosts received nothing.
Here a simple code to illustrate the problem:

Code: Select all

import pickle, time
import  numpy as np
from  mpi4py import MPI
from neuron import h


class cell:
	def __init__(self):
		#Passive cell
		self.soma = h.Section()
		self.soma.L = 40.
		self.soma.diam = 20.
		self.soma.nseg = 1
		self.soma.insert("pas")
		self.soma.g_pas = 0.001
		self.soma.e_pas = -65
		#record voltage into vector
		self.volt = h.Vector()
		self.volt.record(self.soma(0.5)._ref_v)
		#simple synapse
		self.syn  = h.Exp2Syn(0.5, sec=self.soma)
		self.syn.tau1 = 1.
		self.syn.tau2 = 2.
		self.syn.e    = 0
		#record synaptic current
		self.syni = h.Vector()
		self.syni.record(self.syn._ref_i)
		

#Create Parallel Contex
pc = h.ParallelContext()

#if I'm host #0, I also source of spike 
# events for all cells.
# event source has gid=0
if int(pc.id()) == 0:
	d =np.genfromtxt("timevec.dat")
	vs = h.VecStim()
	vec = h.Vector(d.size)
	vec.from_python(d)
	vs.play(vec)
	netcon = h.NetCon(vs,None)
	pc.set_gid2node(0,pc.id())
	pc.cell(0,netcon)

#Create cell
c = cell()
#Connect synapse with source
nc = pc.gid_connect(0,c.syn)
#setup delay and synaptic conductance
nc.delay = 1
nc.weight[0] = 0.1


#load some file for NSG portal
h.load_file("stdgui.hoc")

#init, init max step
# and run
h.stdinit()
h.tstop = 100
h.dt = 0.001
pc.set_maxstep(0.05)
########################
pc.psolve(100) #Works weird
#========= OR =========#
#h.run()            #Works well
########################

# save to pickle object
with open("v%d.pkl"%int(pc.id()),'wb') as fd:
	pickle.dump(np.array(c.volt),fd)
	pickle.dump(np.array(c.syni),fd)

#well done
pc.runworker()
time.sleep(1)
pc.done()


exit()
I've ran this code in my local MPI version of NEUTON 7.3 as well as on NSG portal with same result. If we simple read created pickle files and draw it by matplotlib, we can see this figure:
Image

However when I used h.run() function instead of pc.psolve() the result looked much better
Image

The question is that: is it OK to use h.run() instead psolve()? If no where I made a mistake which lead to this strange result?

Thanks,
-Ruben
ted
Site Admin
Posts: 6289
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: VecStim in parallel environment

Post by ted »

In a parallelized implementation of a network model every spike source must be assigned a unique gid. This applies to instances of the VecStim class, which is a class of artificial spiking cell. This involves using ParallelContext's setgid2node method in order to associate the VecStim's gid with the host on which the VecStim exists, attaching a NetCon to the VecStim to monitor it for spike events, and calling ParallelContext's cell method in order to associate the gid with the NetCon's "spike detector." These steps are required because there is no guarantee that the VecStim will be on the same host as the synaptic mechanisms and artificial spiking cells that are the targets of the spike events that it generates.

Executing a parallel simulation of a parallelized model requires calling ParallelContext's psolve method.
rth
Posts: 41
Joined: Thu Jun 21, 2012 4:47 pm

Re: VecStim in parallel environment

Post by rth »

In a parallelized implementation of a network model every spike source must be assigned a unique gid. This applies to instances of the VecStim class, which is a class of artificial spiking cell. This involves using ParallelContext's setgid2node method in order to associate the VecStim's gid with the host on which the VecStim exists, attaching a NetCon to the VecStim to monitor it for spike events, and calling ParallelContext's cell method in order to associate the gid with the NetCon's "spike detector." These steps are required because there is no guarantee that the VecStim will be on the same host as the synaptic mechanisms and artificial spiking cells that are the targets of the spike events that it generates.
It seems, it was done:

Code: Select all

# if I'm on host 0, only this host has source of spikes
if int(pc.id()) == 0:
#read vector with events from file
   d =np.genfromtxt("timevec.dat") 
 # create VecStim
   vs = h.VecStim()
 # create NEURON vector
   vec = h.Vector(d.size)
 # fill NEURON vector by values in numpy vector
   vec.from_python(d)
 # play NEURON vector by VecStim
   vs.play(vec)
#### HERE THE SPIKE SOURCE REGISTRATION ######
 # create new NetCon, associated with VecStim
   netcon = h.NetCon(vs,None)
 # associate gid 0 with this host
   pc.set_gid2node(0,pc.id())
# register NetCon with gid 0
   pc.cell(0,netcon) 
### THE END OF SPIKE SOURCE REGISTRATION #####
Probably I missed something important.....
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: VecStim in parallel environment

Post by hines »

After pc.set_maxstep, there must be a call to finitialize

Move the

Code: Select all

pc.set_maxstep(0.05) 
to before the

Code: Select all

h.stdinit()
.
This did not matter with h.run() because the first thing h.run() does is call stdinit() (which calls finitialize()).
But pc.psolve() does not call finitialize (it is more akin to the stdrun version of continuerun(tstop)).

Although http://www.neuron.yale.edu/neuron/stati ... et_maxstep kind of backhandedly refers to this in the
Warning
Note: No spikes can be delivered between machines unless this method is called. finitialize relies on this method having been called. If any trans-machine NetCon delay is reduced below the step size, this method MUST be called again. Otherwise an INCORRECT simulation will result.
It should be edited to become more direct. Better yet, I will look into the feasibility of generating an errror message if a call to set_maxstep followed by a call to finitialize did not occur prior to a
call to psolve. Perhaps I can even go so far as to generate an error message if a NetCon.delay is subsequently set which reduces the "minimum NetCon delay interval" without a subsequent call
to set_maxstep.
rth
Posts: 41
Joined: Thu Jun 21, 2012 4:47 pm

Re: VecStim in parallel environment

Post by rth »

Thank you, Michael! It works.
In other words h.stdinit() should be called just before pc.psolve(), otherwise pc.set_maxstep() won't be initialized and any spike won't be delivered. If this is correct, maybe it will be useful to add pc.pinit( maxstep) function in future version, for calling all required functions in right order?
Post Reply