#
# Fitzhugh-Nagumo Model for octave 3.0.1
#
# A.Date 2009.4.30
#  

global theta gamma epsilon INPUT;
theta= 0.01;
gamma = 0.5;
epsilon = 0.002;
INPUT = 0.4;

function dx = nagumo(x,t)
        global theta gamma epsilon INPUT;
        dx(1) = -x(1)*( x(1)- theta)*( x(1)-1.0 ) - x(2) + INPUT;
        dx(2) = epsilon*( x(1)  - gamma*x(2) );
end

t=linspace(0,2000,2000);  
x=lsode("nagumo",[0.5; 0.2],t); # [0.5;0.2] are initial values for x(1) and x(2)
plot( t,x(:,1), "r;u;");  # u is plotted in red line
xlabel ("t");
ylabel ("u,v");
hold on
plot(t,x(:,2),"g;v;" );   # v is plotted in green line
print -deps -FHelvetica:20  hoge.eps
pause

clf
plot(x(:,1),x(:,2));   # see another picture, x-axis u, y-axis v.
xlabel ("u");
ylabel ("v");
print -deps -FHelvetica:20  hoge2.eps
pause

#
# if you don't understand commands, type "help plot" etc.
#
