(* 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 = {
PstatDimer\[LetterSpace]nuc[t], PstatDimer\[LetterSpace]sol[t], Pstat\[LetterSpace]nuc[t], Pstat\[LetterSpace]sol[t], species\[LetterSpace]test[t], statPhosphatase\[LetterSpace]nuc[t], stat\[LetterSpace]nuc[t], stat\[LetterSpace]sol[t]
};

initialValues = {
PstatDimer\[LetterSpace]nuc[0] == 0.0, PstatDimer\[LetterSpace]sol[0] == 0.0, Pstat\[LetterSpace]nuc[0] == 0.0, Pstat\[LetterSpace]sol[0] == 0.0, species\[LetterSpace]test[0] == 0.0, statPhosphatase\[LetterSpace]nuc[0] == 0.05, stat\[LetterSpace]nuc[0] == 0.0, stat\[LetterSpace]sol[0] == 1.0
};

rates = {
PstatDimer\[LetterSpace]\[LetterSpace]import, PstatDimerisation, PstatDimerisationNuc, statDephosphorylation, statPhosphorylation, stat\[LetterSpace]export, stat\[LetterSpace]import
};

rateEquations = {
PstatDimer\[LetterSpace]\[LetterSpace]import -> (nm*PstatDimer\[LetterSpace]\[LetterSpace]import\[LetterSpace]PstatDimer\[LetterSpace]impMax*PstatDimer\[LetterSpace]sol[t])/(PstatDimer\[LetterSpace]\[LetterSpace]import\[LetterSpace]Kpsd\[LetterSpace]imp + PstatDimer\[LetterSpace]sol[t]), PstatDimerisation -> sol*(-(PstatDimerisation\[LetterSpace]Kr\[LetterSpace]PstatDimerisation*PstatDimer\[LetterSpace]sol[t]) + PstatDimerisation\[LetterSpace]Kf\[LetterSpace]PstatDimerisation*Pstat\[LetterSpace]sol[t]^2), PstatDimerisationNuc -> nuc*(-(PstatDimerisationNuc\[LetterSpace]Kr\[LetterSpace]PstatDimerisation*PstatDimer\[LetterSpace]nuc[t]) + PstatDimerisationNuc\[LetterSpace]Kf\[LetterSpace]PstatDimerisation*Pstat\[LetterSpace]nuc[t]^2), statDephosphorylation -> (nuc*statDephosphorylation\[LetterSpace]Kcat\[LetterSpace]dephos*Pstat\[LetterSpace]nuc[t]*statPhosphatase\[LetterSpace]nuc[t])/(statDephosphorylation\[LetterSpace]Km\[LetterSpace]dephos + Pstat\[LetterSpace]nuc[t]), statPhosphorylation -> (sol*statKinase\[LetterSpace]sol*statPhosphorylation\[LetterSpace]Kcat\[LetterSpace]phos*stat\[LetterSpace]sol[t])/(statPhosphorylation\[LetterSpace]Km\[LetterSpace]phos + stat\[LetterSpace]sol[t]), stat\[LetterSpace]export -> (nm*nuc*stat\[LetterSpace]export\[LetterSpace]stat\[LetterSpace]expMax*stat\[LetterSpace]nuc[t])/(stat\[LetterSpace]export\[LetterSpace]Ks\[LetterSpace]exp + stat\[LetterSpace]nuc[t]), stat\[LetterSpace]import -> (nm*nuc*stat\[LetterSpace]import\[LetterSpace]stat\[LetterSpace]impMax*stat\[LetterSpace]sol[t])/(stat\[LetterSpace]import\[LetterSpace]Ks\[LetterSpace]imp + stat\[LetterSpace]sol[t])
};

parameters = {
PstatDimerisation\[LetterSpace]Kf\[LetterSpace]PstatDimerisation -> 0.6, PstatDimerisation\[LetterSpace]Kr\[LetterSpace]PstatDimerisation -> 0.03, statDephosphorylation\[LetterSpace]Kcat\[LetterSpace]dephos -> 1.0, statDephosphorylation\[LetterSpace]Km\[LetterSpace]dephos -> 2.0, statPhosphorylation\[LetterSpace]Kcat\[LetterSpace]phos -> 1.0, statPhosphorylation\[LetterSpace]Km\[LetterSpace]phos -> 4.0, PstatDimerisationNuc\[LetterSpace]Kf\[LetterSpace]PstatDimerisation -> 0.6, PstatDimerisationNuc\[LetterSpace]Kr\[LetterSpace]PstatDimerisation -> 0.03, PstatDimer\[LetterSpace]\[LetterSpace]import\[LetterSpace]PstatDimer\[LetterSpace]impMax -> 0.045, PstatDimer\[LetterSpace]\[LetterSpace]import\[LetterSpace]Kpsd\[LetterSpace]imp -> 0.3, stat\[LetterSpace]export\[LetterSpace]stat\[LetterSpace]expMax -> -0.06, stat\[LetterSpace]export\[LetterSpace]Ks\[LetterSpace]exp -> 0.6, stat\[LetterSpace]import\[LetterSpace]stat\[LetterSpace]impMax -> 0.003, stat\[LetterSpace]import\[LetterSpace]Ks\[LetterSpace]imp -> 3.0, nm -> 1.0, nuc -> 1.0, sol -> 14.625
};

assignments = {
statKinase\[LetterSpace]sol -> Piecewise[{{0.01*Sin[0.001571*(-500 + t)], t > 500 && t < 2502.54614894971}}, 0]
};

events = {

};

speciesAnnotations = {

};

reactionAnnotations = {

};

units = {
{"time" -> "", "metabolite" -> "<substance_units>/litre", "extent" -> "<substance_units>/litre"}
};

(* Time evolution *)
odes = {
PstatDimer\[LetterSpace]nuc'[t] == 1.0*PstatDimerisationNuc +1.0*PstatDimer\[LetterSpace]\[LetterSpace]import , PstatDimer\[LetterSpace]sol'[t] == 1.0*PstatDimerisation -1.0*PstatDimer\[LetterSpace]\[LetterSpace]import, Pstat\[LetterSpace]nuc'[t] ==  -1.0*statDephosphorylation -2.0*PstatDimerisationNuc, Pstat\[LetterSpace]sol'[t] == 1.0*statPhosphorylation -2.0*PstatDimerisation, species\[LetterSpace]test'[t] == 1.0*statPhosphorylation , statPhosphatase\[LetterSpace]nuc'[t] == 0.0 , stat\[LetterSpace]nuc'[t] == 1.0*statDephosphorylation +1.0*stat\[LetterSpace]export +1.0*stat\[LetterSpace]import , stat\[LetterSpace]sol'[t] ==  -1.0*statPhosphorylation -1.0*stat\[LetterSpace]export -1.0*stat\[LetterSpace]import
};
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]}]