(* 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 = {
c2[t], eph[t], fGEF[t], m3PI[t], mRas[t], sumrc1[t], w[t], x1[t], x2[t], y[t], ypp[t], z[t], zpp[t]
};

initialValues = {
c2[0] == 0.0, eph[0] == 1.0, fGEF[0] == 1.0, m3PI[0] == 0.0, mRas[0] == 0.0, sumrc1[0] == 1.0, w[0] == 0.0, x1[0] == 0.0, x2[0] == 0.0, y[0] == 1.0, ypp[0] == 0.0, z[0] == 1.0, zpp[0] == 0.0
};

rates = {
c2ODE, ephODE, fGEFODE, m3PIODE, mRasODE, sumrc1ODE, wODE, x1ODE, x2ODE, yODE, yppODE, zODE, zppODE
};

rateEquations = {
c2ODE -> c1^2*kxR0 - (ke + kminusx)*c2[t], ephODE -> kFBph*(-((-1 + eph[t])/Kph) + w[t]^p/(Wph^p + w[t]^p)), fGEFODE -> kFBf*((1 - fGEF[t])/Kf - (fGEF[t]*zpp[t]^n)/(Zf^n + zpp[t]^n)), m3PIODE -> k3PI*(ePI3K - m3PI[t]), mRasODE -> kRas*(eGEF*(1 + Gamma) - (1 + eGEF*Gamma)*mRas[t]), sumrc1ODE -> 2*(-(c1^2*kxR0) + kminusx*c2[t]) + kt*(1 - sumrc1[t]), wODE -> kdw*(-w[t] + zpp[t]), x1ODE -> kdx1*(mRas[t] - x1[t]/(1 + yp/KMx12 + y[t]/KMx11)), x2ODE -> kdx2*(((1 + Kx2)*m3PI[t])/(1 + Kx2*m3PI[t]) - x2[t]/(1 + yp/KMx22 + y[t]/KMx21)), yODE -> -((VmaxOVERKMx11*x1[t]*y[t])/(1 + yp/KMx12 + y[t]/KMx11)) - (VmaxOVERKMx21*x2[t]*y[t])/(1 + yp/KMx22 + y[t]/KMx21) + (VmaxOVERKMyph1*yp)/(1 + yp/KMyph1 + ypp[t]/KMyph2), yppODE -> (VmaxOVERKMx12*yp*x1[t])/(1 + yp/KMx12 + y[t]/KMx11) + (VmaxOVERKMx22*yp*x2[t])/(1 + yp/KMx22 + y[t]/KMx21) - (VmaxOVERKMyph2*ypp[t])/(ypp[t]/KMyph2 + (1 + yp/KMyph1)*(1 + zp/KMy2 + z[t]/KMy1)), zODE -> -((VmaxOVERKMy1*ypp[t]*z[t])/(1 + zp/KMy2 + z[t]/KMy1)) + (VmaxOVERKMzph1*zp*eph[t])/(1 + zp/KMzph1 + zpp[t]/KMzph2), zppODE -> (VmaxOVERKMy2*zp*ypp[t])/(1 + zp/KMy2 + z[t]/KMy1) - (VmaxOVERKMzph2*eph[t]*zpp[t])/(1 + zp/KMzph1 + zpp[t]/KMzph2)
};

parameters = {
Gamma -> 0.1, KDL -> 1.5, KGP -> 5.09, KGR -> 495.0, KMx11 -> 30.3, KMx12 -> 18.6, KMx21 -> 13.7, KMx22 -> 9.59, KMy1 -> 9.91, KMy2 -> 8.81, KMyph1 -> 23.0, KMyph2 -> 7.99, KMzph1 -> 8.27, KMzph2 -> 31.5, Kf -> 3.76, Kph -> 4.64, Kx2 -> 6.77, L -> 1.0, VmaxOVERKMx11 -> 1.18, VmaxOVERKMx12 -> 3.45, VmaxOVERKMx21 -> 0.405, VmaxOVERKMx22 -> 1.09, VmaxOVERKMy1 -> 6.57, VmaxOVERKMy2 -> 31.9, VmaxOVERKMyph1 -> 1.65, VmaxOVERKMyph2 -> 4.2, VmaxOVERKMzph1 -> 0.167, VmaxOVERKMzph2 -> 0.228, Wph -> 0.385, Zf -> 0.272, alphaPI3K -> 80.0, k3PI -> 1.0, kFBf -> 0.976, kFBph -> 2.34, kRas -> 1.0, kappaPI3K -> 0.3, kdw -> 0.0333, kdx1 -> 0.745, kdx2 -> 2.85, ke -> 0.2, kminusx -> 0.007, kt -> 0.005, kxR0 -> 0.3, n -> 1.03, p -> 1.98, cell -> 1.0
};

assignments = {
r -> (KDL*sumrc1[t])/(KDL + L), c1 -> (L*sumrc1[t])/(KDL + L), ePI3K -> (1 + kappaPI3K + 2*alphaPI3K*c2[t] - (-8*alphaPI3K*c2[t] + (1 + kappaPI3K + 2*alphaPI3K*c2[t])^2)^0.5)/2, eGEF -> (fGEF[t]*(KGR*c2[t] + KGP*m3PI[t]))/(1 + KGR*c2[t] + KGP*m3PI[t]), yp -> 1 - y[t] - ypp[t], zp -> 1 - z[t] - zpp[t]
};

events = {

};

speciesAnnotations = {

};

reactionAnnotations = {

};

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

(* Time evolution *)
odes = {
c2'[t] == 1.0*c2ODE , eph'[t] == 1.0*ephODE , fGEF'[t] == 1.0*fGEFODE , m3PI'[t] == 1.0*m3PIODE , mRas'[t] == 1.0*mRasODE , sumrc1'[t] == 1.0*sumrc1ODE , w'[t] == 1.0*wODE , x1'[t] == 1.0*x1ODE , x2'[t] == 1.0*x2ODE , y'[t] == 1.0*yODE , ypp'[t] == 1.0*yppODE , z'[t] == 1.0*zODE , zpp'[t] == 1.0*zppODE 
};
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]}]