(* 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 = {
ADP[t], ATP[t], PGA[t], Pi[t], Ru5P[t], RuBP[t], TP[t]
};

initialValues = {
ADP[0] == 0.82, ATP[0] == 0.68, PGA[0] == 2.4, Pi[0] == 5.0, Ru5P[0] == 0.0500747384155456, RuBP[0] == 2.0, TP[0] == 0.5
};

rates = {
ATP\[LetterSpace]S, PGA\[LetterSpace]red, RuBP\[LetterSpace]reg, RuBisCO, TPT, TP\[LetterSpace]red
};

rateEquations = {
ATP\[LetterSpace]S -> chloroplast*function\[LetterSpace]10[V6, ADP[t], Pi[t], ATP\[LetterSpace]S\[LetterSpace]K1, ATP\[LetterSpace]S\[LetterSpace]K2], PGA\[LetterSpace]red -> chloroplast*function\[LetterSpace]6[PGA\[LetterSpace]red\[LetterSpace]Vm, PGA[t], ATP[t], ADP[t], TP[t], Pi[t], PGA\[LetterSpace]red\[LetterSpace]k, PGA\[LetterSpace]red\[LetterSpace]K1, PGA\[LetterSpace]red\[LetterSpace]K2], RuBP\[LetterSpace]reg -> chloroplast*function\[LetterSpace]8[RuBP\[LetterSpace]reg\[LetterSpace]Vm, Ru5P[t], ATP[t], RuBP\[LetterSpace]reg\[LetterSpace]K1, RuBP\[LetterSpace]reg\[LetterSpace]K2, RuBP\[LetterSpace]reg\[LetterSpace]K3, Pi[t]], RuBisCO -> chloroplast*function\[LetterSpace]5[RuBisCO\[LetterSpace]k, E\[LetterSpace]RuBisCO, totRuBP, RuBisCO\[LetterSpace]K], TPT -> function\[LetterSpace]9[TP[t], Pic, TPc, Pi[t], TPT\[LetterSpace]K2, TPT\[LetterSpace]K1, TPT\[LetterSpace]Vm], TP\[LetterSpace]red -> chloroplast*function\[LetterSpace]7[TP\[LetterSpace]red\[LetterSpace]Vm, TP[t], TP\[LetterSpace]red\[LetterSpace]K]
};

parameters = {
P\[LetterSpace]0 -> 16.0, V6 -> 5.8801285588795, E\[LetterSpace]RuBisCO -> 3.56, Pic -> 1.4, TPc -> 0.2, RuBisCO\[LetterSpace]k -> 0.504, RuBisCO\[LetterSpace]K -> 0.04, PGA\[LetterSpace]red\[LetterSpace]Vm -> 3.49, PGA\[LetterSpace]red\[LetterSpace]k -> 14.0, PGA\[LetterSpace]red\[LetterSpace]K1 -> 1.0, PGA\[LetterSpace]red\[LetterSpace]K2 -> 1.0, TP\[LetterSpace]red\[LetterSpace]Vm -> 1.06, TP\[LetterSpace]red\[LetterSpace]K -> 0.4, RuBP\[LetterSpace]reg\[LetterSpace]Vm -> 4.81, RuBP\[LetterSpace]reg\[LetterSpace]K1 -> 0.05, RuBP\[LetterSpace]reg\[LetterSpace]K2 -> 0.5, RuBP\[LetterSpace]reg\[LetterSpace]K3 -> 0.05, TPT\[LetterSpace]K2 -> 0.25, TPT\[LetterSpace]K1 -> 0.08, TPT\[LetterSpace]Vm -> 3.3, ATP\[LetterSpace]S\[LetterSpace]K1 -> 0.08, ATP\[LetterSpace]S\[LetterSpace]K2 -> 0.5, chloroplast -> 1.0, cytosol -> 1.0
};

assignments = {
function\[LetterSpace]9[S1_,S2_,P1_,P2_,K2_,K1_,Vm_] -> ((-(P1*P2) + S1*S2)*Vm)/(K2*(P1 + S1) + K1*(P2 + S2) + K1*K2*(P2/K2 + S1/K1)*(P1/K1 + S2/K2)), function\[LetterSpace]10[Vm_,s1_,s2_,K1_,K2_] -> (s1*s2*Vm)/((K1 + s1)*(K2 + s2)), function\[LetterSpace]8[Vm_,S1_,S2_,K1_,K2_,K3_,R_] -> (S1*S2*Vm)/(K1*K2 + K3*R + K2*S2 + S1*S2), function\[LetterSpace]7[Vm_,S_,K_] -> (S*Vm)/(K + S), function\[LetterSpace]5[k_,E_,S_,K_] -> (k*(E + K + S - Sqrt[-4*E*S + (E + K + S)^2]))/2, function\[LetterSpace]6[Vm_,S1_,S2_,P1_,P2_,P3_,k_,K1_,K2_] -> ((-((P1*P2*P3)/k) + S1*S2)*Vm)/(K1 + (P1*P2*P3)/k + (K1*S1*S2)/K2), totRuBP -> (P\[LetterSpace]0 - ATP[t] - PGA[t] - Pi[t] - Ru5P[t] - TP[t])/2
};

events = {

};

speciesAnnotations = {
ADP[t]->"http://identifiers.org/chebi/CHEBI:16761"
};

reactionAnnotations = {

};

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

(* Time evolution *)
odes = {
ADP'[t] == 1.0*PGA\[LetterSpace]red +1.0*RuBP\[LetterSpace]reg -1.0*ATP\[LetterSpace]S, ATP'[t] == 1.0*ATP\[LetterSpace]S -1.0*PGA\[LetterSpace]red -1.0*RuBP\[LetterSpace]reg, PGA'[t] == 2.0*RuBisCO -1.0*PGA\[LetterSpace]red, Pi'[t] == 1.0*PGA\[LetterSpace]red +2.0*TP\[LetterSpace]red +1.0*TPT -1.0*ATP\[LetterSpace]S, Ru5P'[t] == 3.0*TP\[LetterSpace]red -1.0*RuBP\[LetterSpace]reg, RuBP'[t] == 1.0*RuBP\[LetterSpace]reg -1.0*RuBisCO, TP'[t] == 1.0*PGA\[LetterSpace]red -5.0*TP\[LetterSpace]red -1.0*TPT
};
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]}]