IGEM:Melbourne/2008/BioClock/David Version/clock code

function xdot = Stephen_clock(time,x)

format long e; global pulse g_p;

% synopsis: %    xdot =  (time, x) %     x0 = % %  can be used with odeN functions as follows: % % x0 = ; % [t,x] = ode23s(@, [0 100], ); % plot (t,x); % % where 100 is the end time % % When  is used without any arguments it returns a vector of %  the initial concentrations of the 14 floating species. % Otherwise  should be called with two arguments: % time and x.  time is the current time. x is the vector of the % concentrations of the 14 floating species. % When these parameters are supplied  returns a vector of %  the rates of change of the concentrations of the 14 floating species. % % the following table shows the mapping between the vector % index of a floating species and the species name. % %  NOTE for compartmental models % matlab translator generates code that when simulated in matlab, % produces results which have the units of species amounts. Users % should divide the results for each species with the volume of the % compartment it resides in, in order to obtain concentrations. % %  Indx      Name % x(1)        EnvZ_dephos % x(2)        OmpR_dephos % x(3)        X %  x(4)        GFP % x(5)        RiboKey1 % x(6)        Y %  x(7)        RFP % x(8)        Off_X % x(9)        RiboKey2 % x(10)        GFP_sink % x(11)        RFP_sink % x(12)        X_sink % x(13)        Y_sink % x(14)        EnvZ_dephos_sink % x(15)        OmpR_dephos_sink % x(16)        Intermediate1 % x(17)        RiboKey1_sink % x(18)        Intermediate2 % x(19)        Intermediate3 % x(20)        X_sink_617291

xdot = zeros(20, 1); g_p = ones(30,1);

% List of Compartments vol__compartment = 1;

% Global Parameters g_p(1) = 0.1;		% j1_k g_p(2) = 0.1;		% j10_k g_p(3) = 0.1;		% j14_k g_p(4) = 0.1;		% j16_k g_p(5) = 5;		% j2_vmax g_p(6) = 0.2;		% j2_km1 g_p(7) = 0.3;		% j2_ki g_p(8) = 0.5;		% j5_vmax g_p(9) = 0.3;		% j5_km1 g_p(10) = 0.1;		% j5_ki g_p(11) = 0.5;		% j12_vmax g_p(12) = 0.3;		% j12_km1 g_p(13) = 1;		% j12_ki g_p(14) = 0.5;		% j19_k g_p(15) = 0.1;		% j20_k g_p(16) = 2;		% j21_k g_p(17) = 0.5;		% j8_vmax g_p(18) = 0.3;		% j8_km1 g_p(19) = 0.1;		% j8_ki g_p(20) = 0.5;		% j9_vmax g_p(21) = 0.3;		% j9_km1 g_p(22) = 0.1;		% j9_ki g_p(23) = 0.1;		% j17_k g_p(24) = 0.1;		% j22_k g_p(25) = 1;		% j18_k g_p(26) = 0.1;		% j23_vmax g_p(27) = 1e-005;		% j23_km1 g_p(28) = 0.1;		% j24_k g_p(29) = 0.1;		% j25_k g_p(30) = 5000;		% j26_k

% Boundary Conditions g_p31 = 1;		% envz_phos [Concentration] %g_p32 = 0;		% red_light [Concentration] g_p33 = 1;		% ompr_phos [Concentration] g_p34 = 0;		% off_y [Concentration]

if pulse == 1 g_p32 = getSG0(time); %red_light similar pulse to Sg0 in the other model dlmwrite('time_v2.dat',time,'-append'); dlmwrite('sg0_v2.dat',g_p32,'-append'); else g_p32 = 0;		% red_light [Concentration] end

% Other Parameters

if (nargin == 0) % set initial conditions xdot(1) = 0*vol__compartment;		% EnvZ_dephos [Concentration] xdot(2) = 0*vol__compartment;		% OmpR_dephos [Concentration] xdot(3) = 0*vol__compartment;		% X [Concentration] xdot(4) = 0*vol__compartment;		% GFP [Concentration] xdot(5) = 0*vol__compartment;		% RiboKey1 [Concentration] xdot(6) = 0*vol__compartment;		% Y [Concentration] xdot(7) = 0*vol__compartment;		% RFP [Concentration] xdot(8) = 0*vol__compartment;		% Off_X [Concentration] xdot(9) = 0*vol__compartment;		% RiboKey2 [Concentration] xdot(10) = 0*vol__compartment;		% GFP_sink [Concentration] xdot(11) = 0*vol__compartment;		% RFP_sink [Concentration] xdot(12) = 0*vol__compartment;		% X_sink [Concentration] xdot(13) = 0*vol__compartment;		% Y_sink [Concentration] xdot(14) = 0*vol__compartment;		% EnvZ_dephos_sink [Concentration] xdot(15) = 0*vol__compartment;		% OmpR_dephos_sink [Concentration] xdot(16) = 0*vol__compartment;		% Intermediate1 [Concentration] xdot(17) = 0*vol__compartment;		% RiboKey1_sink [Concentration] xdot(18) = 0*vol__compartment;		% Intermediate2 [Concentration] xdot(19) = 0*vol__compartment;		% Intermediate3 [Concentration] xdot(20) = 0*vol__compartment;		% X_sink_617291 [Concentration] else % calculate rates of change R0 = g_p(1)*(x(1)/vol__compartment)*(g_p33/vol__compartment); R1 = g_p(5)*(x(2)/vol__compartment)/(g_p(6)*(1+(g_p33/vol__compartment)/g_p(7))+(x(2)/vol__compartment)); R2 = g_p(8)*(x(3)/vol__compartment)/(g_p(9)*(1+(x(8)/vol__compartment)/g_p(10))+(x(3)/vol__compartment)); R3 = g_p(17)*(x(6)/vol__compartment)/(g_p(18)*(1+(g_p34/vol__compartment)/g_p(19))+(x(6)/vol__compartment)); R4 = g_p(20)*(x(6)/vol__compartment)/(g_p(21)*(1+(g_p34/vol__compartment)/g_p(22))+(x(6)/vol__compartment)); R5 = g_p(2)*(x(6)/vol__compartment); R6 = g_p(11)*(x(3)/vol__compartment)/(g_p(12)*(1+(x(8)/vol__compartment)/g_p(13))+(x(3)/vol__compartment)); R7 = g_p(3)*(g_p31/vol__compartment)*(g_p32/vol__compartment); R8 = g_p(4)*(x(8)/vol__compartment); R9 = g_p(25)*(x(4)/vol__compartment); R10 = g_p(14)*(x(7)/vol__compartment); R11 = g_p(15)*(x(3)/vol__compartment); R12 = g_p(16)*(x(6)/vol__compartment); R13 = g_p(23)*(x(1)/vol__compartment); R14 = g_p(24)*(x(2)/vol__compartment); R15 = g_p(26)*(x(3)/vol__compartment)/(g_p(27)+(x(3)/vol__compartment)); R16 = g_p(28)*(x(5)/vol__compartment); R17 = g_p(29)*(x(16)/vol__compartment)*(x(2)/vol__compartment); R18 = g_p(30)*(x(16)/vol__compartment)*(x(3)/vol__compartment);

xdot = [ - R0 + R7 - R13 + R0 - R1 - R14 - R17 + R1 + R2 - R2 - R6 - R11 - R15 - R18 + R6 - R9     + 50*R15 - R16 + R3 - R3 - R4 - R5 - R12 + 20*R17 + R4 - R10 + R8 - R8     + R5      + R9      + R10 + R11 + R12 + R13 + R14 + R16 - R17 - R18 0       0        0      + R18 ];  end;