
        format none;

        function xdot = f (x, t)
            # Zmena hodnoty x za jednotku casu
            xdot = zeros (44,1);
            xdot(1) = 0; # fixed, cyt 
xdot(2) = 0; # fixed, carb 
xdot(3) = 0; # fixed, v1 
xdot(5) = 0; # fixed, K1 
xdot(8) = 0; # fixed, K4 
xdot(9) = 0; # fixed, K5 
xdot(11) = 0; # fixed, kW 
xdot(12) = 0; # fixed, pH_ext 
xdot(13) = 0; # fixed, pH_cyt 
xdot(14) = 0; # fixed, pH_carb 
xdot(15) = 0; # fixed, fast 
xdot(16) = 0; # fixed, O2_carb 
xdot(21) = 0; # fixed, PCO2 
xdot(23) = 0; # fixed, Temperature 
xdot(24) = 0; # fixed, IonicStrength 
xdot(26) = 0; # fixed, k7 
xdot(27) = 0; # fixed, Nsites 
xdot(28) = 0; # fixed, Vcarb 
xdot(29) = 0; # fixed, Vcell 
xdot(31) = 0; # fixed, factorCo2dehyd 
xdot(32) = 0; # fixed, factorHCO3dehyd 
xdot(4) = x(4) = x(26)*x(27)/x(28); # assignment, v2  FROM INIT EXPRESSION!!! 
xdot(6) = x(6) = 12100+1156*(x(23)-25); # assignment, k2  FROM INIT EXPRESSION!!! 
xdot(7) = x(7) = 0.00006+0.00006*(x(23)-25); # assignment, k3  FROM INIT EXPRESSION!!! 
xdot(10) = x(10) = 10^-(1.5+0.01*(x(23)-25)+0.1*x(24)); # assignment, kH  FROM INIT EXPRESSION!!! 
xdot(17) = x(17) = x(10)*x(21); # assignment, CO2_ext  FROM INIT EXPRESSION!!! 
xdot(20) = x(20) = x(11)/10^(-x(12)); # assignment, OH_ext  FROM INIT EXPRESSION!!! 
xdot(18) = x(18) = (x(6)/x(7))*x(17)*x(20); # assignment, HCO3_ext  FROM INIT EXPRESSION!!! 
xdot(19) = x(19) = x(11)/10^(-x(13)); # assignment, OH_cyt  FROM INIT EXPRESSION!!! 
xdot(22) = x(22) = x(11)/10^(-x(14)); # assignment, OH_carb  FROM INIT EXPRESSION!!! 
xdot(25) = x(25) = 10^-(6.4+0.0001*x(23)^2-0.01*x(23)+1.7*x(24)^2-1.9*x(24))	; # assignment, K6  FROM INIT EXPRESSION!!! 
xdot(30) = x(30) = (x(4)*x(43))/(x(43)+x(8)*(1+x(16)/x(9))); # assignment, CO2_fixation_rate  FROM INIT EXPRESSION!!! 
xdot(33) = x(33) = x(7)*x(42); # assignment, rate3  FROM INIT EXPRESSION!!! 
xdot(34) = x(34) =  x(6)*x(41)*x(19); # assignment, rate2  FROM INIT EXPRESSION!!! 
xdot(35) = x(35) = x(15)*(x(17)-x(41)) ; # assignment, rate6  FROM INIT EXPRESSION!!! 
xdot(36) = x(36) = ((x(3)*x(18))/(x(5) + x(18)))/x(29); # assignment, rate1  FROM INIT EXPRESSION!!! 
xdot(37) = x(37) = x(31)*x(15)*x(6)*x(22)*x(43); # assignment, rate4  FROM INIT EXPRESSION!!! 
xdot(38) = x(38) = x(32)*x(15)*x(7)*x(44); # assignment, rate5  FROM INIT EXPRESSION!!! 
xdot(39) = x(39) = x(15)*x(42)-x(15)*x(44); # assignment, rate7  FROM INIT EXPRESSION!!! 
xdot(40) = x(40) = (x(4)*x(43))/(x(43)+x(8)*(1+x(16)/x(9))); # assignment, rate8  FROM INIT EXPRESSION!!! 
xdot(41) =  - (x(6)*x(41)*x(19))  + (x(15)*(x(17)-x(41)))  + (x(7)*x(42)) ; # reaction, CO2_cyt 
xdot(42) =  + (x(6)*x(41)*x(19))  - (x(15)*x(42)-x(15)*x(44))  - (x(7)*x(42))  + (((x(3)*x(18))/(x(5) + x(18)))/x(29)) ; # reaction, HCO3_cyt 
xdot(43) =  - ((x(4)*x(43))/(x(43)+x(8)*(1+x(16)/x(9))))  - (x(31)*x(15)*x(6)*x(22)*x(43))  + (x(32)*x(15)*x(7)*x(44)) ; # reaction, CO2_carb 
xdot(44) =  + (x(31)*x(15)*x(6)*x(22)*x(43))  + (x(15)*x(42)-x(15)*x(44))  - (x(32)*x(15)*x(7)*x(44)) ; # reaction, HCO3_carb 
xdot(4) = 0; # end
xdot(6) = 0; # end
xdot(7) = 0; # end
xdot(10) = 0; # end
xdot(17) = 0; # end
xdot(18) = 0; # end
xdot(19) = 0; # end
xdot(20) = 0; # end
xdot(22) = 0; # end
xdot(25) = 0; # end
xdot(30) = 0; # end
xdot(33) = 0; # end
xdot(34) = 0; # end
xdot(35) = 0; # end
xdot(36) = 0; # end
xdot(37) = 0; # end
xdot(38) = 0; # end
xdot(39) = 0; # end
xdot(40) = 0; # end

        endfunction

        # time scale
        t = linspace(0, 100, 10000);

        # default values
        x0 = zeros (44,1);
        x0(1) = 1; # fixed compartment, cyt
x0(2) = 1; # fixed compartment, carb
x0(3) = 2e-18; # fixed parameter, v1
x0(12) = 8.2; # fixed parameter, pH_ext
x0(13) = 8; # fixed parameter, pH_cyt
x0(14) = 8; # fixed parameter, pH_carb
x0(15) = 100; # fixed parameter, fast
x0(21) = 0.0004; # fixed parameter, PCO2
x0(23) = 38; # fixed parameter, Temperature
x0(24) = 0.1; # fixed parameter, IonicStrength
x0(27) = 13400; # fixed parameter, Nsites
x0(31) = 1; # fixed parameter, factorCo2dehyd
x0(32) = 1; # fixed parameter, factorHCO3dehyd
x0(26) = 1.9*10^-23; # assignment parameter, k7
x0(28) = 1.7*10^-17; # assignment parameter, Vcarb
x0(4) = x0(26)*x0(27)/x0(28); # assignment parameter, v2
x0(5) = 2.2*10^-4; # assignment parameter, K1
x0(6) = 12100+1156*(x0(23)-25); # assignment parameter, k2
x0(7) = 0.00006+0.00006*(x0(23)-25); # assignment parameter, k3
x0(8) = 1.9*10^-4; # assignment parameter, K4
x0(9) = 1.3*10^-3; # assignment parameter, K5
x0(10) = 10^-(1.5+0.01*(x0(23)-25)+0.1*x0(24)); # assignment parameter, kH
x0(11) = 10^-14; # assignment parameter, kW
x0(16) = 4*10^-4; # assignment parameter, O2_carb
x0(17) = x0(10)*x0(21); # assignment parameter, CO2_ext
x0(20) = x0(11)/10^(-x0(12)); # assignment parameter, OH_ext
x0(18) = (x0(6)/x0(7))*x0(17)*x0(20); # assignment parameter, HCO3_ext
x0(19) = x0(11)/10^(-x0(13)); # assignment parameter, OH_cyt
x0(22) = x0(11)/10^(-x0(14)); # assignment parameter, OH_carb
x0(25) = 10^-(6.4+0.0001*x0(23)^2-0.01*x0(23)+1.7*x0(24)^2-1.9*x0(24))	; # assignment parameter, K6
x0(29) = 5.2*10^-16; # assignment parameter, Vcell
x0(41) = x0(17); # assignment specie, CO2_cyt
x0(43) = x0(41); # assignment specie, CO2_carb
x0(30) = (x0(4)*x0(43))/(x0(43)+x0(8)*(1+x0(16)/x0(9))); # assignment parameter, CO2_fixation_rate
x0(42) = (x0(6)/x0(7))*x0(41)*x0(19); # assignment specie, HCO3_cyt
x0(33) = x0(7)*x0(42); # assignment parameter, rate3
x0(34) =  x0(6)*x0(41)*x0(19); # assignment parameter, rate2
x0(35) = x0(15)*(x0(17)-x0(41)) ; # assignment parameter, rate6
x0(36) = ((x0(3)*x0(18))/(x0(5) + x0(18)))/x0(29); # assignment parameter, rate1
x0(37) = x0(31)*x0(15)*x0(6)*x0(22)*x0(43); # assignment parameter, rate4
x0(44) = (x0(6)/x0(7))*x0(43)*x0(22); # assignment specie, HCO3_carb
x0(38) = x0(32)*x0(15)*x0(7)*x0(44); # assignment parameter, rate5
x0(39) = x0(15)*x0(42)-x0(15)*x0(44); # assignment parameter, rate7
x0(40) = (x0(4)*x0(43))/(x0(43)+x0(8)*(1+x0(16)/x0(9))); # assignment parameter, rate8


        #lsode_options ( 'absolute tolerance', max (50*eps, 0.5e-28) );
        lsode_options ( 'relative tolerance', max (50*eps, 0.5e-28) ); 

        y = lsode ("f", x0, t);
        
        for i = 1:length(t)
            y(i, 4) = y(i, 26)*y(i, 27)/y(i, 28); # assignment, v2  FROM INIT EXPRESSION!!! 
y(i, 6) = 12100+1156*(y(i, 23)-25); # assignment, k2  FROM INIT EXPRESSION!!! 
y(i, 7) = 0.00006+0.00006*(y(i, 23)-25); # assignment, k3  FROM INIT EXPRESSION!!! 
y(i, 10) = 10^-(1.5+0.01*(y(i, 23)-25)+0.1*y(i, 24)); # assignment, kH  FROM INIT EXPRESSION!!! 
y(i, 17) = y(i, 10)*y(i, 21); # assignment, CO2_ext  FROM INIT EXPRESSION!!! 
y(i, 20) = y(i, 11)/10^(-y(i, 12)); # assignment, OH_ext  FROM INIT EXPRESSION!!! 
y(i, 18) = (y(i, 6)/y(i, 7))*y(i, 17)*y(i, 20); # assignment, HCO3_ext  FROM INIT EXPRESSION!!! 
y(i, 19) = y(i, 11)/10^(-y(i, 13)); # assignment, OH_cyt  FROM INIT EXPRESSION!!! 
y(i, 22) = y(i, 11)/10^(-y(i, 14)); # assignment, OH_carb  FROM INIT EXPRESSION!!! 
y(i, 25) = 10^-(6.4+0.0001*y(i, 23)^2-0.01*y(i, 23)+1.7*y(i, 24)^2-1.9*y(i, 24))	; # assignment, K6  FROM INIT EXPRESSION!!! 
y(i, 30) = (y(i, 4)*y(i, 43))/(y(i, 43)+y(i, 8)*(1+y(i, 16)/y(i, 9))); # assignment, CO2_fixation_rate  FROM INIT EXPRESSION!!! 
y(i, 33) = y(i, 7)*y(i, 42); # assignment, rate3  FROM INIT EXPRESSION!!! 
y(i, 34) =  y(i, 6)*y(i, 41)*y(i, 19); # assignment, rate2  FROM INIT EXPRESSION!!! 
y(i, 35) = y(i, 15)*(y(i, 17)-y(i, 41)) ; # assignment, rate6  FROM INIT EXPRESSION!!! 
y(i, 36) = ((y(i, 3)*y(i, 18))/(y(i, 5) + y(i, 18)))/y(i, 29); # assignment, rate1  FROM INIT EXPRESSION!!! 
y(i, 37) = y(i, 31)*y(i, 15)*y(i, 6)*y(i, 22)*y(i, 43); # assignment, rate4  FROM INIT EXPRESSION!!! 
y(i, 38) = y(i, 32)*y(i, 15)*y(i, 7)*y(i, 44); # assignment, rate5  FROM INIT EXPRESSION!!! 
y(i, 39) = y(i, 15)*y(i, 42)-y(i, 15)*y(i, 44); # assignment, rate7  FROM INIT EXPRESSION!!! 
y(i, 40) = (y(i, 4)*y(i, 43))/(y(i, 43)+y(i, 8)*(1+y(i, 16)/y(i, 9))); # assignment, rate8  FROM INIT EXPRESSION!!! 

        endfor
        
        [t',y]
        