
h0 = rn0 - n0;
Ina = (aminf*gnabar').*h0.*(Vn*[1 1]-ones(Nv,1)*Vna') ;
Ik = (ones(Nv,1)*gkbar').*n0.^4.*(Vn*[1 1]-ones(Nv,1)*Vk');
Icl = (ones(Nv,1)*gl').*(Vn*[1 1]-ones(Nv,1)*Vl');
      
Fvn = -Ina-Ik-Icl+ones(Nv,1)*Inp;
