function dxdt=hhode_parc(t,x)
%HH egyenlet megolása:
d = 1;
C=1;
Gk=36;
GNa=120;
GL=0.3;
Ek=-72;
ENa=55;
EL=-49.4;
ujinger_t=20; 

dxdt=zeros(length(x),1);

%kezdet
n=x(1);
m=x(2);
h=x(3);
V=x(4);

%adott idõben (ujinger_t), adott méretû ingerlés (V) az axon eljére (x(4))
%ujinger_t=10 idõpillanat múlva; V=-20 -> nincs AP
%ujinger_t=10 ; V=-10 -> van AP
%ujinger_t=20 ; V=-20 -> van AP
%ujinger_t=5 ; V=-10 -> van AP
%
if abs(t-ujinger_t)<0.1
V=-20;
end
%}
dxdt(1) = ((1-n)*alphan(V)) - (n*betan(V));
dxdt(2) = ((1-m)*alpham(V)) - (m*betam(V));
dxdt(3) = ((1-h)*alphah(V)) - (h*betah(V));
dxdt(4) =(x(4)-2*x(4)+x(8))/d^2 -(1/C)*(Gk*n^4*(V-Ek)+GNa*m^3*h*(V-ENa)+GL*(V-EL));

%közben
for i = 1:length(x)/4-2
    n=x(4*i+1);
    m=x(4*i+2);
    h=x(4*i+3);
    V=x(4*i+4);
    dxdt(4*i+1) = ((1-n)*alphan(V)) - (n*betan(V));
    dxdt(4*i+2) = ((1-m)*alpham(V)) - (m*betam(V));
    dxdt(4*i+3) = ((1-h)*alphah(V)) - (h*betah(V));
    dxdt(4*i+4) = (x(4*i)-2*x(4*i+4)+x(4*i+8))/d^2-(1/C)*(Gk*n^4*(V-Ek)+GNa*m^3*h*(V-ENa)+GL*(V-EL));
    
end

%vége
n=x(end-3);
m=x(end-2);
h=x(end-1);
V=x(end);
dxdt(end-3) = ((1-n)*alphan(V)) - (n*betan(V));
dxdt(end-2) = ((1-m)*alpham(V)) - (m*betam(V));
dxdt(end-1) = ((1-h)*alphah(V)) - (h*betah(V));
dxdt(end) = (x(end-4)-2*x(end)+x(end))/d^2-(1/C)*(Gk*n^4*(V-Ek)+GNa*m^3*h*(V-ENa)+GL*(V-EL));
