(* 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 = {
ATP[t], Fru16P2[t], HMP[t]
};

initialValues = {
ATP[0] == 4.0, Fru16P2[0] == 1.0, HMP[0] == 0.1
};

rates = {
v\[LetterSpace]1, v\[LetterSpace]2, v\[LetterSpace]3, v\[LetterSpace]4
};

rateEquations = {
v\[LetterSpace]1 -> (Glc*VHK*ATP[t])/(KATP*KGlc*(1 + ATP[t]/KATP)*(1 + Glc/KGlc + (wildtype*HMP[t]^2)/KiTre6P)), v\[LetterSpace]2 -> (gR*VPFK*ATP[t]*HMP[t]*(1 + ATP[t]/KRATP + HMP[t]/KRHMP + (gR*ATP[t]*HMP[t])/(KRATP*KRHMP)))/(KRATP*KRHMP*((1 + ATP[t]/KRATP + HMP[t]/KRHMP + (gR*ATP[t]*HMP[t])/(KRATP*KRHMP))^2 + (L0*(1 + (ci*ATP[t])/KiATP)^2*(1 + (c2*ATP[t])/KRATP + (c1*HMP[t])/KRHMP + (c1*c2*gT*ATP[t]*HMP[t])/(KRATP*KRHMP))^2)/(1 + ATP[t]/KiATP)^2)), v\[LetterSpace]3 -> (Vlower*(5 - ATP[t])*Fru16P2[t])/(KADP*KFru16P2*(1 + (5 - ATP[t])/KADP)*(1 + Fru16P2[t]/KFru16P2)), v\[LetterSpace]4 -> (VATPase*ATP[t])/(KATP3 + ATP[t])
};

parameters = {
EXTERNAL -> 0.0, KADP -> 0.1, KATP -> 0.15, KATP3 -> 3.0, KFru16P2 -> 1.0, KGlc -> 1.0, KRATP -> 0.06, KRHMP -> 1.0, KiATP -> 10.0, KiTre6P -> 4.422, L0 -> 1000.0, VATPase -> 68.0, VHK -> 68.0, VPFK -> 30.0, Vlower -> 20.0, c1 -> 0.0005, c2 -> 1.0, ci -> 10.0, gR -> 10.0, gT -> 1.0, wildtype -> 1.0, Glc -> 10.0, default\[LetterSpace]compartment -> 1.0
};

assignments = {

};

events = {

};

speciesAnnotations = {
ATP[t]->"http://identifiers.org/kegg.compound/C00002", ATP[t]->"urn:miriam:obo.chebi:CHEBI%3A15422", Fru16P2[t]->"urn:miriam:obo.chebi:CHEBI%3A16905", Fru16P2[t]->"http://identifiers.org/kegg.compound/C00354", Glc[t]->"urn:miriam:obo.chebi:CHEBI%3A17925", Glc[t]->"http://identifiers.org/kegg.compound/C00293", HMP[t]->"urn:miriam:obo.chebi:CHEBI%3A16084", HMP[t]->"urn:miriam:obo.chebi:CHEBI%3A17665"
};

reactionAnnotations = {
v\[LetterSpace]1->"urn:miriam:reactome:REACT\[LetterSpace]1318", v\[LetterSpace]1->"http://identifiers.org/interpro/IPR001312", v\[LetterSpace]2->"urn:miriam:reactome:REACT\[LetterSpace]1840", v\[LetterSpace]2->"http://identifiers.org/interpro/IPR022953", v\[LetterSpace]3->"urn:miriam:ec-code:5.3.1.1", v\[LetterSpace]3->"urn:miriam:ec-code:2.7.1.40", v\[LetterSpace]3->"urn:miriam:ec-code:2.7.2.3", v\[LetterSpace]3->"urn:miriam:ec-code:4.2.1.11", v\[LetterSpace]3->"urn:miriam:ec-code:4.1.2.13", v\[LetterSpace]3->"urn:miriam:ec-code:5.4.2.1", v\[LetterSpace]3->"urn:miriam:ec-code:1.2.1.12", v\[LetterSpace]4->"urn:miriam:obo.go:GO%3A0016887", v\[LetterSpace]4->"http://identifiers.org/interpro/IPR020958"
};

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

(* Time evolution *)
odes = {
ATP'[t] == 4.0*v\[LetterSpace]3 -1.0*v\[LetterSpace]2 -1.0*v\[LetterSpace]1 -1.0*v\[LetterSpace]4, Fru16P2'[t] == 1.0*v\[LetterSpace]2 -1.0*v\[LetterSpace]3, HMP'[t] == 1.0*v\[LetterSpace]1 -1.0*v\[LetterSpace]2
};
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]}]