(* Generated by JWS Online *)

(* This is an experimental feature of JWS Online. Please report any mistakes.*)

(* Note that the following notable SBML entities or features are not supported in notebook outputyet: *)
    (* Events *)
    (* Constraints *)
    (* Units and UnitDefinitions *)
    (* AlgebraicRules *)
    (* conversionFactors *)

variables = {
ara[t], mGFPcr[t], mT3cr[t], mT7cr[t], pGFP[t], pT3[t], pT7[t], taRNA[t]
};

initialValues = {
ara[0] == 0.0, mGFPcr[0] == 0.176991329, mT3cr[0] == 0.00566438, mT7cr[0] == 0.3569405099, pGFP[0] == 6.338921181, pT3[0] == 6.41674*^-05, pT7[0] == 0.05230744612, taRNA[0] == 0.006796941377
};

rates = {
r0, r10a, r10b, r11, r12, r1a, r1b, r2a, r2b, r3a, r3b, r6, r7, r8, r9
};

rateEquations = {
r0 -> cell*Piecewise[{{cAra, pulse\[LetterSpace]flag == 1 && ara[t] > 0}}, dAra*ara[t]], r10a -> cell*(s0\[LetterSpace]mT3cr + (k\[LetterSpace]pT7*pT7[t]^n7)/(km7^n7 + pT7[t]^n7)), r10b -> cell*d\[LetterSpace]mT3*mT3cr[t], r11 -> cell*(s0\[LetterSpace]pT3*mT3cr[t] + s\[LetterSpace]pT3k*mT3cr[t]*taRNA[t]), r12 -> cell*d\[LetterSpace]pT3*pT3[t], r1a -> cell*(s0\[LetterSpace]taRNA + (sT*ara[t])/(k\[LetterSpace]ara + ara[t])), r1b -> cell*d\[LetterSpace]taRNA*taRNA[t], r2a -> cell*s0\[LetterSpace]mT7cr, r2b -> cell*d\[LetterSpace]mT7*mT7cr[t], r3a -> cell*(s0\[LetterSpace]mGFPcr + (k\[LetterSpace]pT3*pT3[t]^n3)/(km3^n3 + pT3[t]^n3)), r3b -> cell*d\[LetterSpace]mGFP*mGFPcr[t], r6 -> cell*(s0\[LetterSpace]pT7*mT7cr[t] + s\[LetterSpace]pT7k*mT7cr[t]*taRNA[t]), r7 -> cell*(s0\[LetterSpace]pGFP*mGFPcr[t] + s\[LetterSpace]pGFPk*mGFPcr[t]*taRNA[t]), r8 -> cell*d\[LetterSpace]pT7*pT7[t], r9 -> cell*d\[LetterSpace]pGFP*pGFP[t]
};

parameters = {
cAra -> 0.0003, dAra -> 0.1201, d\[LetterSpace]mGFP -> 0.07, d\[LetterSpace]mT3 -> 0.0701, d\[LetterSpace]mT7 -> 0.0706, d\[LetterSpace]pGFP -> 0.003, d\[LetterSpace]pT3 -> 0.0069, d\[LetterSpace]pT7 -> 0.0056, d\[LetterSpace]taRNA -> 0.1177, k\[LetterSpace]ara -> 0.0571, k\[LetterSpace]pT3 -> 3.006, k\[LetterSpace]pT7 -> 3.8009, km3 -> 7.9075, km7 -> 3.0455, n3 -> 0.8892, n7 -> 2.602, pulse1\[LetterSpace]length -> 11.0, pulse1\[LetterSpace]start -> 0.01, pulse2\[LetterSpace]length -> 11.0, pulse3\[LetterSpace]length -> 22.0, pulse\[LetterSpace]conc -> 0.01, pulse\[LetterSpace]flag -> 0.0, pulse\[LetterSpace]interval -> 20.0, s0\[LetterSpace]mGFPcr -> 0.0123, s0\[LetterSpace]mT3cr -> 0.0003, s0\[LetterSpace]mT7cr -> 0.0252, s0\[LetterSpace]pGFP -> 0.1007, s0\[LetterSpace]pT3 -> 0.0, s0\[LetterSpace]pT7 -> 0.0003, s0\[LetterSpace]taRNA -> 0.0008, sT -> 0.8467, s\[LetterSpace]pGFPk -> 0.9923, s\[LetterSpace]pT3k -> 0.0115, s\[LetterSpace]pT7k -> 0.0766, cell -> 1.0
};

assignments = {
pulse3\[LetterSpace]start -> pulse2\[LetterSpace]length + pulse2\[LetterSpace]start + pulse\[LetterSpace]interval, pulse2\[LetterSpace]start -> pulse1\[LetterSpace]length + pulse1\[LetterSpace]start + pulse\[LetterSpace]interval
};

events = {

};

speciesAnnotations = {
ara[t]->"http://identifiers.org/kegg.compound/C00259", ara[t]->"http://identifiers.org/chebi/CHEBI:17535"
};

reactionAnnotations = {

};

units = {
{"time" -> "", "metabolite" -> "", "extent" -> ""}
};

(* Time evolution *)
odes = {
ara'[t] ==  -1.0*r0, mGFPcr'[t] == 1.0*r3a -1.0*r3b, mT3cr'[t] == 1.0*r10a -1.0*r10b, mT7cr'[t] == 1.0*r2a -1.0*r2b, pGFP'[t] == 1.0*r7 -1.0*r9, pT3'[t] == 1.0*r11 -1.0*r12, pT7'[t] == 1.0*r6 -1.0*r8, taRNA'[t] == 1.0*r1a -1.0*r1b
};
timeCourse = NDSolve[Join[odes, initialValues]//.rateEquations//.assignments//.parameters, variables, {t, 0, 100}];

(* Steady-state solution initialized with result of time evolution *)
findRootEquations = odes /.D[_[t],t]->0;
findRootVariables = Partition[Flatten[{#, #/.timeCourse/.t->100} &/@variables],2];
steadyStateVariables = FindRoot[findRootEquations//.rateEquations//.assignments//.parameters, findRootVariables, MaxIterations->100]
fluxes = #//.assignments//.parameters/.steadyStateVariables&/@rateEquations

(* Plot the time evolution of the variables *)
plotTable=Table[Plot[variables[[i]]/.parameters/.timeCourse,{t,0,100},PlotLegends->variables[[i]],PlotRange->Full],{i,Length[variables]}]