(* 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], bpg[t], dhap[t], f6p[t], gap[t], glc[t], nadp[t], nadph[t], p3g[t], pep[t], pyr[t]
};

initialValues = {
adp[0] == 0.0, atp[0] == 10.0, bpg[0] == 0.0, dhap[0] == 0.03, f6p[0] == 0.12, gap[0] == 0.06, glc[0] == 10.0, nadp[0] == 0.0, nadph[0] == 0.2, p3g[0] == 1.8, pep[0] == 5.0, pyr[0] == 0.4
};

rates = {
v\[LetterSpace]1, v\[LetterSpace]2, v\[LetterSpace]3, v\[LetterSpace]4, v\[LetterSpace]5, v\[LetterSpace]6, v\[LetterSpace]7, v\[LetterSpace]8, v\[LetterSpace]9
};

rateEquations = {
v\[LetterSpace]1 -> (protPGK*(-((VmrPGK*adp[t]*bpg[t])/(KpgkADP*KpgkbisP)) + (VmfPGK*atp[t]*p3g[t])/(KpgkATP*KpgkP3G)))/((1 + adp[t]/KiADP)*(1 + ((1 + adp[t]/KpgkADP)*bpg[t])/KpgkbisP + ((1 + atp[t]/KpgkATP)*p3g[t])/KpgkP3G)), v\[LetterSpace]2 -> (protGAPDH*(bpg[t]/KbisP + (pi*gap[t])/(Kgap*Kpi))^(-1 + n)*(-((pi*Vmarev*gap[t]*nadp[t])/(Kgap*Knadp*Kpi)) + (Vmfor*bpg[t]*nadph[t])/(KbisP*Knadph)))/((1 + (pi/Kpi)^n + (gap[t]/Kgap)^n + (bpg[t]/KbisP + (pi*gap[t])/(Kgap*Kpi))^n)*(1 + nadp[t]/Knadp + nadph[t]/Knadph)), v\[LetterSpace]3 -> (protTIM*(-((VmrTIM*dhap[t])/KtimDHAP) + (VmfTIM*gap[t])/KtimGAP))/(1 + dhap[t]/KtimDHAP + gap[t]/KtimGAP + p3g[t]/KtimP3G + pep[t]/KtimPEP), v\[LetterSpace]4 -> (protALDPase*VmfALD*dhap[t]*gap[t])/(KaldDHAP*KaldGAP*(1 + dhap[t]/KaldDHAP + gap[t]/KaldGAP + (dhap[t]*gap[t])/(KaldDHAP*KaldGAP))), v\[LetterSpace]5 -> kdgap*gap[t], v\[LetterSpace]6 -> kddhap*dhap[t], v\[LetterSpace]7 -> kdbpg*bpg[t], v\[LetterSpace]8 -> kGDH*glc[t]*nadp[t], v\[LetterSpace]9 -> kPK*adp[t]*pep[t]
};

parameters = {
KaldDHAP -> 0.170995, KaldF6P -> 1.0, KaldGAP -> 0.0522437, KbisP -> 0.000406744, Kgap -> 0.83766, KiADP -> 1.14174, Knadp -> 0.271013, Knadph -> 0.0735253, KpgkADP -> 0.0848049, KpgkATP -> 9.68404, KpgkP3G -> 0.541454, KpgkbisP -> 5.59, Kpi -> 408.523, KtimDHAP -> 0.811626, KtimGAP -> 0.245106, KtimP3G -> 0.399822, KtimPEP -> 0.660519, Vmarev -> 24.054, VmfALD -> 2.64096, VmfPGK -> 17.3221, VmfTIM -> 239.835, Vmfor -> 20.6526, VmrPGK -> 30.04, VmrTIM -> 239.414, kGDH -> 10.0, kPK -> 10.0, kdbpg -> 1.05824, kddhap -> 0.0225, kdgap -> 0.0559, n -> 1.56416, pi -> 0.0, protALDPase -> 0.0036, protGAPDH -> 0.0421, protPGK -> 0.0034, protTIM -> 0.00085, ga -> 0.0, default\[LetterSpace]compartment -> 1.0
};

assignments = {

};

events = {

};

speciesAnnotations = {
adp[t]->"http://identifiers.org/kegg.compound/C00008", adp[t]->"http://identifiers.org/obo.chebi/CHEBI%3A16761", atp[t]->"urn:miriam:obo.chebi:CHEBI:15422", atp[t]->"urn:miriam:kegg.compound:C00002", bpg[t]->"http://identifiers.org/obo.chebi/CHEBI%3A28907", bpg[t]->"http://identifiers.org/kegg.compound/C00236", dhap[t]->"http://identifiers.org/obo.chebi/CHEBI%3A16108", dhap[t]->"http://identifiers.org/kegg.compound/C00111", f6p[t]->"http://identifiers.org/kegg.compound/C00085", f6p[t]->"http://identifiers.org/obo.chebi/CHEBI%3A15946", gap[t]->"urn:miriam:obo.chebi:CHEBI:17138", gap[t]->"urn:miriam:kegg.compound:C00661", glc[t]->"http://identifiers.org/obo.chebi/CHEBI%3A17234", glc[t]->"http://identifiers.org/kegg.compound/C00031", nadp[t]->"http://identifiers.org/kegg.compound/C00006", nadp[t]->"http://identifiers.org/obo.chebi/CHEBI%3A18009", nadph[t]->"http://identifiers.org/kegg.compound/C00005", nadph[t]->"http://identifiers.org/obo.chebi/CHEBI%3A16474", p3g[t]->"http://identifiers.org/kegg.compound/C00197", p3g[t]->"http://identifiers.org/obo.chebi/CHEBI%3A17794", pep[t]->"http://identifiers.org/kegg.compound/C00074", pep[t]->"http://identifiers.org/obo.chebi/CHEBI%3A44897", pyr[t]->"http://identifiers.org/kegg.compound/C00022", pyr[t]->"http://identifiers.org/obo.chebi/CHEBI%3A32816"
};

reactionAnnotations = {
v\[LetterSpace]1->"urn:miriam:kegg.reaction:R01512", v\[LetterSpace]2->"urn:miriam:kegg.reaction:R01061", v\[LetterSpace]3->"urn:miriam:kegg.reaction:R01015", v\[LetterSpace]4->"urn:miriam:kegg.reaction:R01068"
};

units = {
{"time" -> "min", "metabolite" -> "mmol/L", "extent" -> "mmol/L"}
};

(* Time evolution *)
odes = {
adp'[t] == 1.0*v\[LetterSpace]1 -1.0*v\[LetterSpace]9, atp'[t] == 1.0*v\[LetterSpace]9 -1.0*v\[LetterSpace]1, bpg'[t] == 1.0*v\[LetterSpace]1 -1.0*v\[LetterSpace]2 -1.0*v\[LetterSpace]7, dhap'[t] == 1.0*v\[LetterSpace]3 -1.0*v\[LetterSpace]6 -1.0*v\[LetterSpace]4, f6p'[t] == 1.0*v\[LetterSpace]4 , gap'[t] == 1.0*v\[LetterSpace]2 -1.0*v\[LetterSpace]4 -1.0*v\[LetterSpace]5 -1.0*v\[LetterSpace]3, glc'[t] ==  -1.0*v\[LetterSpace]8, nadp'[t] == 1.0*v\[LetterSpace]2 -1.0*v\[LetterSpace]8, nadph'[t] == 1.0*v\[LetterSpace]8 -1.0*v\[LetterSpace]2, p3g'[t] ==  -1.0*v\[LetterSpace]1, pep'[t] ==  -1.0*v\[LetterSpace]9, pyr'[t] == 1.0*v\[LetterSpace]9 
};
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]}]