Dear all,

I am using NEURON to model ephaptic coupling between axons in a nerve. My model necessarily requires coupling NEURON to a Poisson solver (FEM, although for minimal working examples, I just use phi = 1/(4*pi*sigma) * sum(membr.currents/distance) for calculating the extracellular potential from the membrane currents of axons).

At every time step in a simulation, I need an iterative process that brings the system to convergence. Specifically, I am trying to make sure that the membrane currents and the extracellular potential field (phi) are consistent with each other before moving on to the next time step.

To this end, my strategy is the one described below. The process that goes on right before every time step of the simulations I perform has roughly the following structure:

----------

for (every time step):

h.advance()

svstateobj.save()

convergence = False

while not convergence:

1. Calculate phi from the membrane currents and store it in an array.

2. Play the corresponding phi on _ref_e_extracellular on every segment of every axon (note that the time series for every segment includes only one value, since it is only meant .

3. h.finitialize() (otherwise, I think that playing the computed phi into the axons is worhless)

4. svstateobj.restore(1) # The last saved state. For the first iteration of the convergence process, this is the state right after the last h.advance()

5. h.fcurrent()

6. Save state (svstateobj.save()).

7. Check if the values of phi and the membrane currents have converged (I have a certain tolerance value). If so, convergence = True and exit the while loop to move on to the next time step.

----------

(Every simulation I run is actually a chain of simulations of length dt (the length of the time step).)

My hope was that h.fcurrent() recalculated the membrane currents having into account the played new values of phi, but it seems to ignore them. I think it makes sense, since as said in the documentation, it takes into account the STATES (which I guess they are m, n and h for the H-H model), when a played vector is not such thing, but an event (am I right or did I misunderstand?).

So my question is: is there a way to make fcurrent acknowledge the pressence of phi so that the new membrane currents are calculated according to the influence of phi on the axon membrane at that time?

As a last resource, I could do this by other means (maybe by translating phi into currents on every segment like McIntyre & al. did), although if there is a way to do it with NEURON, that would make things much easier.

I would appreciate if you had the time to shed some light to me about this.

All the best.