I have recently started using NEURON and have been facing difficulties running basic reactions in that the concentrations of generated species never deviate from the initial concentrations. I am using NEURON v 7.4 and executing my models out of the Spyder interface downloaded through Anaconda. Simulations without reactions run fine and produce reasonable output. I have been largely trying to execute the basic hydration reaction from the tutorial posted. I have included the code I have been using as well as some commented alternative methods I have tried to execute the reactions.
Code: Select all
from neuron import h, rxd
from matplotlib import pyplot
#h.load_file('stdrun.hoc')
#h.CVode().active(1) <- if uncommented gives unreasonable time steps
#h.initnrn(0, 0.025)
dend = h.Section()
r = rxd.Region(h.allsec())
#r = rxd.Region(h.allsec(), nrn_region='i')
hydrogen = rxd.Species(r, initial=1)
oxygen = rxd.Species(r, initial=1)
water = rxd.Species(r, initial=0)
reaction = rxd.Reaction(2 * hydrogen + oxygen > water, 1)
#reaction = rxd.Reaction(2 * hydrogen + oxygen > water, 1, regions=r)
h.finitialize()
#h.finitialize(-65)
#h.CVode().re_init()
heading = '{t:>10s} {h:>10s} {o:>10s} {h2o:>10s}'
data = '{t:10g} {h:10g} {o:10g} {h2o:10g}'
def advance_a_bit():
for i in xrange(5):
h.fadvance()
print data.format(t=h.t, h=hydrogen.nodes[0].concentration,
o=oxygen.nodes[0].concentration,
h2o=water.nodes[0].concentration)
print heading.format(t='t', h='hydrogen', o='oxygen', h2o='H2O')
print heading.format(t='-', h='--------', o='------', h2o='---')
advance_a_bit()
"""
Method 2: Printing Concentrations Adapted form 'Calcium Waves in RxD'
h.CVode().re_init()
time = []
concentration = []
def save():
print 't=', h.t, 'conc=', water.nodes[0].concentration #dend.wateri
for t in xrange(0, 1000, 10):
h.CVode().event(t, save)
time.append(h.t)
#concentration.append(dend.wateri)
h.continuerun(1000)
"""
"""
Method 3: Using recording variables
#Set up recording variables
v_vec = h.Vector() #Membrane potential vector
t_vec = h.Vector() #Time stamp vector
v_vec.record(water.nodes[0].concentration)
t_vec.record(h._ref_t)
#Run the simulation
h.tstop = 40.0
h.run()
#Plot the results
from matplotlib import pyplot
pyplot.figure(figsize=(8,4))
pyplot.plot(t_vec, v_vec)
pyplot.xlabel('Time (ms)')
pyplot.ylabel('mV')
pyplot.show()
"""