## Khaliq03 model driver for Octave


Khaliq03_defs = "Khaliq03.m";

autoload ("Khaliq03", Khaliq03_defs );
autoload ("Khaliq03_init", Khaliq03_defs );

Khaliq03_init; 


y0 = Khaliq03_init(-65)

global N;
N = length(y0);

function is = Khaliq03_currents_odepkg (y)
  global N;
  global i_Leak i_Ih i_Narsg i_K3 i_K2 i_K1 i_CaP i_CaBK;

  t = y(1);
  s = y(2:N+1);
  y1 = Khaliq03(t,s);

  Na_zO = y1(10);
  is = [t i_Leak i_Ih i_Narsg i_K3 i_K2 i_K1 i_CaP i_CaBK Na_zO];

endfunction

function is = Khaliq03_currents_lsode (y)
  global N;
  global i_Leak i_Ih i_Narsg i_K3 i_K2 i_K1 i_CaP i_CaBK;

  t = y(1);
  s = y(2:N+1);
  y1 = Khaliq03(s,t);

  Na_zO = y1(10);
  is = [t i_Leak i_Ih i_Narsg i_K3 i_K2 i_K1 i_CaP i_CaBK Na_zO];

endfunction


function ys = Khaliq03_vclamp_run1_odepkg (N,dt,vchold,vchdur,vcbase,vcdur)

  global reltol abstol

  P = odeset ("RelTol", reltol, "AbsTol", abstol, "MaxStep", 1, "InitialStep", dt);

  t0 = 0.0;
  t1 = t0+vchdur;

  t2 = t1;
  t3 = t2+vcdur;

  t4 = t3;
  t5 = t4+vchdur;

  vinit = vchold;
  y0 = Khaliq03_init(vinit)
  sol = odesx(@Khaliq03, [t0 t1], y0, P);
  ys = [sol.x sol.y];

  vinit = vcbase;
  y0 = sol.y(size(sol.y)(1),:)'; y0(1) = vinit
  sol = odesx(@Khaliq03, [t2 t3], y0, P);
  ys = vertcat (ys, [sol.x sol.y]);

  vinit = vchold;
  y0 = sol.y(size(sol.y)(1),:)'; y0(1) = vinit
  sol = odesx(@Khaliq03, [t4 t5], y0, P);
  ys = vertcat (ys, [sol.x sol.y]);

endfunction;


function ys = Khaliq03_vclamp_run1_lsode (N,dt,vchold,vchdur,vcbase,vcdur)

  global reltol abstol

  lsode_options("absolute tolerance", abstol);
  lsode_options("relative tolerance", reltol);
  lsode_options("integration method", "bdf");
  lsode_options("initial step size", dt);
  
  t0 = 0.0;
  t1 = t0+vchdur;
  thold0 = linspace(t0,t1,(t1-t0)/dt);

  t2 = t1;
  t3 = t2+vcdur;
  tclamp = linspace(t2,t3,(t3-t2)/dt);

  t4 = t3;
  t5 = t4+vchdur;
  thold1 = linspace(t4,t5,(t5-t4)/dt);

  vinit = vchold;
  y0 = Khaliq03_init(vinit)
  y = lsode(@Khaliq03, y0, thold0);
  ys = [thold0' y];

  vinit = vcbase;
  y0 = y(size(y)(1),:)'; y0(1) = vinit
  y = lsode(@Khaliq03, y0, tclamp);
  ys = vertcat (ys, [tclamp' y]);

  vinit = vchold;
  y0 = y(size(y)(1),:)'; y0(1) = vinit
  y = lsode(@Khaliq03, y0, thold1);

  ys = vertcat (ys, [thold1' y]);

endfunction;

function ilog = Khaliq03_vclamp (N,dt,vchold,vchdur,vcbase,vcdur,vcinc,vcsteps)


  vc = vcbase;
  ys = cell(1,vcsteps);

  for i = 1:vcsteps
    
    y = Khaliq03_vclamp_run1_odepkg (N,dt,vchold,vchdur,vc,vcdur);
    ys(i) = y;
    vc = vc + vcinc;
    
  endfor;
  
  ilog = cell(vcsteps,1);
  for i = 1:vcsteps
    ilogv = [];
    n = size(ys{i})(1);
    for j = 1:n
      next = Khaliq03_currents_odepkg (ys{i}(j,:));
      ilogv = vertcat (ilogv, next);
    endfor; 
    ilog{i} = ilogv;
  endfor;
     
endfunction

global reltol abstol

reltol = 1e-2;
abstol = 1e-3;

## i_CaBK plot

i_CaBK_index = 9;
Khaliq03_ilog = Khaliq03_vclamp (N,0.1,-90,5,-40,20,10,5);

subplot(3,3,4);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_CaBK_index),
     Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_CaBK_index),
     Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_CaBK_index),
     Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_CaBK_index),
     Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_CaBK_index));

CaBK_log = vertcat ([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_CaBK_index)],
		    [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_CaBK_index)],
		    [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_CaBK_index)],
		    [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_CaBK_index)],
		    [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_CaBK_index)]);
save -ascii "CaBK.dat" CaBK_log;

## i_CaP plot

i_CaP_index = 8;
Khaliq03_ilog = Khaliq03_vclamp (N,0.25,-90,5,-90,10,10,11);

subplot(3,3,5);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_CaP_index),
       Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_CaP_index),
       Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_CaP_index),
       Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_CaP_index),
       Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_CaP_index),
       Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_CaP_index),
       Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_CaP_index),
       Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_CaP_index),
       Khaliq03_ilog{9}(:,1),Khaliq03_ilog{9}(:,i_CaP_index),
       Khaliq03_ilog{10}(:,1),Khaliq03_ilog{10}(:,i_CaP_index),
       Khaliq03_ilog{11}(:,1),Khaliq03_ilog{11}(:,i_CaP_index));

CaP_log = vertcat ([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_CaP_index)],
		   [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_CaP_index)],
		   [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_CaP_index)],
		   [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_CaP_index)],
		   [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_CaP_index)],
		   [Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_CaP_index)],
		   [Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_CaP_index)],
		   [Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_CaP_index)],
		   [Khaliq03_ilog{9}(:,1),Khaliq03_ilog{9}(:,i_CaP_index)],
		   [Khaliq03_ilog{10}(:,1),Khaliq03_ilog{10}(:,i_CaP_index)],
		   [Khaliq03_ilog{11}(:,1),Khaliq03_ilog{11}(:,i_CaP_index)]);

save -ascii "CaP.dat" CaP_log;

## i_K1 plot

i_K1_index = 7;
Khaliq03_ilog = Khaliq03_vclamp (N,1,-71,30,-69,100,10,8);

subplot(3,3,1);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_K1_index),
      Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_K1_index),
      Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_K1_index),
      Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_K1_index),
      Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_K1_index),
      Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_K1_index),
      Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_K1_index),
      Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_K1_index));

K1_log = vertcat([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_K1_index)],
		 [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_K1_index)],
		 [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_K1_index)],
		 [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_K1_index)],
		 [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_K1_index)],
		 [Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_K1_index)],
		 [Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_K1_index)],
		 [Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_K1_index)]);

save -ascii "K1.dat" K1_log;

## i_K2 plot

i_K2_index = 6;
Khaliq03_ilog = Khaliq03_vclamp (N,1,-71,30,-69,100,10,9);

subplot(3,3,2);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_K2_index),
      Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_K2_index),
      Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_K2_index),
      Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_K2_index),
      Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_K2_index),
      Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_K2_index),
      Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_K2_index),
      Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_K2_index),
      Khaliq03_ilog{9}(:,1),Khaliq03_ilog{9}(:,i_K2_index));

K2_log = vertcat([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_K2_index)],
		 [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_K2_index)],
		 [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_K2_index)],
		 [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_K2_index)],
		 [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_K2_index)],
		 [Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_K2_index)],
		 [Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_K2_index)],
		 [Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_K2_index)],
		 [Khaliq03_ilog{9}(:,1),Khaliq03_ilog{9}(:,i_K2_index)]);

save -ascii "K2.dat" K2_log;

## i_K3 plot

i_K3_index = 5;
Khaliq03_ilog = Khaliq03_vclamp (N,1,-71,30,-61,100,10,8);

subplot(3,3,3);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_K3_index),
      Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_K3_index),
      Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_K3_index),
      Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_K3_index),
      Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_K3_index),
      Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_K3_index),
      Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_K3_index),
      Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_K3_index));

K3_log = vertcat([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_K3_index)],
		 [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_K3_index)],
		 [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_K3_index)],
		 [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_K3_index)],
		 [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_K3_index)],
		 [Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_K3_index)],
		 [Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_K3_index)],
		 [Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_K3_index)]);

save -ascii "K3.dat" K3_log;


## i_Ih plot

i_Ih_index = 3;
Khaliq03_ilog = Khaliq03_vclamp (N,1,-50,100,-60,1200,-10,7);

subplot(3,3,6);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_Ih_index),
      Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_Ih_index),
      Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_Ih_index),
      Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_Ih_index),
      Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_Ih_index),
      Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_Ih_index),
      Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_Ih_index));

Ih_log = vertcat([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_Ih_index)],
		 [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_Ih_index)],
		 [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_Ih_index)],
		 [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_Ih_index)],
		 [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_Ih_index)],
		 [Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_Ih_index)],
		 [Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_Ih_index)]);

save -ascii "Ih.dat" Ih_log;

%% i_Narsg plot

i_Narsg_index = 4;
Khaliq03_ilog = Khaliq03_vclamp (N,1,-71,30,-60,100,10,9);

subplot(3,3,8);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_Narsg_index),
     Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_Narsg_index),
     Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_Narsg_index),
     Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_Narsg_index),
     Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_Narsg_index),
     Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_Narsg_index),
     Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_Narsg_index),
     Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_Narsg_index),
     Khaliq03_ilog{9}(:,1),Khaliq03_ilog{9}(:,i_Narsg_index));

Narsg_log = vertcat ([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_Narsg_index)],
		     [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_Narsg_index)],
		     [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_Narsg_index)],
		     [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_Narsg_index)],
		     [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_Narsg_index)],
		     [Khaliq03_ilog{6}(:,1),Khaliq03_ilog{6}(:,i_Narsg_index)],
		     [Khaliq03_ilog{7}(:,1),Khaliq03_ilog{7}(:,i_Narsg_index)],
		     [Khaliq03_ilog{8}(:,1),Khaliq03_ilog{8}(:,i_Narsg_index)],
		     [Khaliq03_ilog{9}(:,1),Khaliq03_ilog{9}(:,i_Narsg_index)]);

save -ascii "Narsg.dat" Narsg_log;

## i_Leak plot

i_Leak_index = 2;
Khaliq03_ilog = Khaliq03_vclamp (N,1,-71,30,-91,100,10,5);

subplot(3,3,7);
plot(Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_Leak_index),
     Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_Leak_index),
     Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_Leak_index),
     Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_Leak_index),
     Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_Leak_index));

Leak_log = vertcat([Khaliq03_ilog{1}(:,1),Khaliq03_ilog{1}(:,i_Leak_index)],
		   [Khaliq03_ilog{2}(:,1),Khaliq03_ilog{2}(:,i_Leak_index)],
		   [Khaliq03_ilog{3}(:,1),Khaliq03_ilog{3}(:,i_Leak_index)],
		   [Khaliq03_ilog{4}(:,1),Khaliq03_ilog{4}(:,i_Leak_index)],
		   [Khaliq03_ilog{5}(:,1),Khaliq03_ilog{5}(:,i_Leak_index)]);

save -ascii "Leak.dat" Leak_log;