function varargout=hhodered(t, y, flag)% HHODERED - ODE file for reduced 2-variable HH equations. State vector is% [V,n]'. m = minf(V), h = nvtohv(n).% External current injected by function iexthh(t). Model parameters returned% using hhpars = [gnamx; Ena; gkmx; Ek; gleak; El; C] = gethhparms.% flag='': returns ydot (column)%     ='init': returns [tsp,y0,opt], time span, initial values, options.if nargin<3							% Catch for numjac the dumbjack	varargout{1} = hhdotred(t,y);else	switch flag		case ''                % Return ydot vector			varargout{1} = hhdotred(t,y);		case 'init'             % Return default tspan, y0, and options			[varargout{1:3}] = initred;	%	case 'jacobian'         % Return Jacobian matrix	%		varargout{1} = hhjac(t,y);		otherwise			error(strcat('***Unknown flag=''',flag,''' in HHODERED.***'))	endendreturnfunction ydot = hhdotred(t,y)% Note the external current is provided by function iexthh(t)% Parameters (mS/cm^2, mV, uFd/cm^2)% hps = [gnamx=120; Ena=55; gkmx=36; Ek=-72; gleak=0.3; El=-49.4; C=1]    hps = gethhparms;	ydot=zeros(2,1);	% Membrane potential:	ydot(1) = (iexthh(t) - hps(1)*minfhh(y(1))^3*nvtohv(y(2))*(y(1)-hps(2)) - ...	      hps(3)*y(2)^4*(y(1)-hps(4)) - hps(5)*(y(1)-hps(6)))/hps(7);% n	ydot(2) = (1-y(2))*alphn(y(1)) - y(2)*betan(y(1));returnfunction [tsp,y0,opt]=initred   tsp = [0; 50];   y0 = [-60., 0.318]';   opt = odeset('OutputFcn','odeplot', 'jacobian','off', 'RelTol',1e-4);return