from neuron import h, rxd from matplotlib import pyplot # needed for standard run system h.load_file('stdrun.hoc') left = h.Section(name='left') right = h.Section(name='right') left.nseg = right.nseg = 101 left.L = right.L = 101 right.connect(left) tfact0 = 137 cytosol = rxd.Region([left, right]) tfact = rxd.Species(cytosol, d=1, initial=0) def right_only(node): return 1 if node.x3d > left.L else 0 production_region = rxd.Parameter(cytosol, initial=right_only) reaction = rxd.Rate(tfact, 0.2 * production_region) h.finitialize(-65) for node in tfact.nodes: if 90 < node.x3d < 110: node.concentration = tfact0 def plot_it(color='b'): y = tfact.nodes.concentration x = [node.x3d for node in tfact.nodes] pyplot.scatter(x, y, marker='o', color=color, s=5) # plot the initial situation plot_it('r') # we will plot every 25ms up to 100ms for i in range(1, 5): # run for another 25ms h.continuerun(i * 25) plot_it() pyplot.show()