Passing variables to a function with lsode

[mathneuro]

I am playing around with Octave's lsode function to see if it can be a
proper MATLAB substitute for my modeling. I coded up a FH neuron and
noticed that a seemingly simple change made the code run very slowly.
I was hoping someone could help me see if the problem can be fixed.

Here is the system:

**Primary function: fhn.m**
function [y,t] = fhn
% initial condition
y0 = [0.25; 0.3];

% time vector
t = linspace(0, 100, 1000);

% ode solver
y = lsode("derivs",y,t);
end

** Secondary function: derivs.m **
function ydot = derivs(y, t)

% assignment of variables
v = y(1);
w = y(2);

% params
a = 0.25;
eps = 0.05;
gamma = 1;
I0 = 0.25;
alpha = 0;
omega = 2;

% eqns
s = alpha * sin(omega * t);
f = v * (1 - v) * (v - a);

% derivs
vdot = f - w + s + I0;
wdot = eps * (v - gamma * w);

ydot = [vdot; wdot];
end

If I make a change to pass I0 (the initial current) to the derivs
function, without using global variables, here's what I would do:

**Primary function: fhn.m**
function [y,t] = fhn(I0)
% initial condition
y0 = [0.25; 0.3];

% time vector
t = linspace(0, 100, 1000);

% ode solver
g = @(y, t) derivs(y,t,I0);
y = lsode(g,y,t,I0);
end

** Secondary function: derivs.m **
function ydot = derivs(y, t, I0)

% assignment of variables
v = y(1);
w = y(2);

% params
a = 0.25;
eps = 0.05;
gamma = 1;
% I0 = 0.25; % no longer hard coded
alpha = 0;
omega = 2;

% eqns
s = alpha * sin(omega * t);
f = v * (1 - v) * (v - a);

% derivs
vdot = f - w + s + I0;
wdot = eps * (v - gamma * w);

ydot = [vdot; wdot];
end

In my tic/toc of these two methods, the first one ran in an avg time
of about 0.14 seconds. The second one took about 3 s. I recognize that
it's passing I0 around on each iteration, and while that might not
seem meaningful in this toy example, there are plenty of uses for it.
Is there a proper way of passing I0 to derivs.m via lsode that would
make the speed difference negligible?

Thanks

[DRF]

Post by mathneuro
I am playing around with Octave's lsode function to see if it can be a
proper MATLAB substitute for my modeling. I coded up a FH neuron and
noticed that a seemingly simple change made the code run very slowly.
I was hoping someone could help me see if the problem can be fixed.
**Primary function: fhn.m**
function [y,t] = fhn
     % initial condition
     y0 = [0.25; 0.3];
     % time vector
     t = linspace(0, 100, 1000);
     % ode solver
     y = lsode("derivs",y,t);
end
** Secondary function: derivs.m **
function ydot = derivs(y, t)
     % assignment of variables
     v = y(1);
     w = y(2);
     % params
     a = 0.25;
     eps = 0.05;
     gamma = 1;
     I0 = 0.25;
     alpha = 0;
     omega = 2;
     % eqns
     s = alpha * sin(omega * t);
     f = v * (1 - v) * (v - a);
     % derivs
     vdot = f - w + s + I0;
     wdot = eps * (v - gamma * w);
     ydot = [vdot; wdot];
end
If I make a change to pass I0 (the initial current) to the derivs
**Primary function: fhn.m**
function [y,t] = fhn(I0)
     % initial condition
     y0 = [0.25; 0.3];
     % time vector
     t = linspace(0, 100, 1000);
     % ode solver
     y = lsode(g,y,t,I0);
end
** Secondary function: derivs.m **
function ydot = derivs(y, t, I0)
     % assignment of variables
     v = y(1);
     w = y(2);
     % params
     a = 0.25;
     eps = 0.05;
     gamma = 1;
     % I0 = 0.25; % no longer hard coded
     alpha = 0;
     omega = 2;
     % eqns
     s = alpha * sin(omega * t);
     f = v * (1 - v) * (v - a);
     % derivs
     vdot = f - w + s + I0;
     wdot = eps * (v - gamma * w);
     ydot = [vdot; wdot];
end
In my tic/toc of these two methods, the first one ran in an avg time
of about 0.14 seconds. The second one took about 3 s. I recognize that
it's passing I0 around on each iteration, and while that might not
seem meaningful in this toy example, there are plenty of uses for it.
Is there a proper way of passing I0 to derivs.m via lsode that would
make the speed difference negligible?
Thanks

In both versions of fh.m the second argument to lsode must be y0, not
y. As written, your programs do not run.

I don't know why you supply an extra argument to lsode in the second
version, but it does you a lot of harm. Running your program (with the
y -> y0 fix) takes 0.304 sec for the first version and 4.20 sec for
the second. When I delete the fourth argument to lsode in the second
version, it runs in 0.336 sec.

Finally, revising the first version to use a global variable for I0,
which is set from the argument of fhn, gives me a timing of 0.300 sec.