This Workbook sets up, simulates and graphs the equations that describe adenosine system. It corresponds to Section 8.3.1 of Ch 8 in SB2.
The reactions considered are:

v,use :  ATP--> ADP
v,form: ADP-ATP
v, distr :  2ADP<=>ATP+AMP

BOP April 24, 2007
Addiel de Alba. April 2010

FInd the steady state

kform=10/0.4;
kuse=10/1.6;
kdist=1000;kmdist=1000; Atot=2.1;
ststsolution=NSolve[{-kuse*a3+kform*a2==0,kdist*a2*a2-kmdist*a3*a1==0, a3+a2+a1-Atot==0},{a3,a2,a1}];

Solve dynamic behavior from changing the IC or the use rate

kuse=1.5*(10/1.6);
va1=kuse*atp[t];
va2=kform*adp[t];
va5=kdist*adp[t]*adp[t];
vam5=kmdist*atp[t]*amp[t];

BalanceATP=-va1+va2+(va5-vam5);
BalanceADP=va1-va2-2*(va5-vam5);
BalanceAMP=(va5-vam5);

tfinal=10;
solution=NDSolve[{atp'[t]==BalanceATP,adp'[t]==BalanceADP,amp'[t]==BalanceAMP,atp[0]==1.60,adp[0]==0.40,amp[0]==0.1},{atp,adp,amp},{t,0,tfinal}];

<<"PlotLegends`"
graph1=Plot[{atp[t]/.solution,adp[t]/.solution,amp[t]/.solution},
{t,0,1},PlotLabel->"Concentrations vs time",PlotLegend->{"atp","adp","amp"},
PlotStyle->{Hue[0.5],Hue[0.6],Hue[0.8]},PlotRange->{0,2.1},
LegendPosition->{1.1,-.4},LegendSize->{.45,.35}];

graph2=ParametricPlot[Evaluate[{kuse*atp[t],kform*adp[t]}/.solution],{t,0,tfinal},PlotLabel->"ATP load vs ATP synthesis", PlotStyle->{Hue[.7]},PlotRange->{{9,15},{9,15}} ];

graph3=Plot[{(atp[t]+adp[t]/2)/(atp[t]+adp[t]+amp[t])/.solution},
{t,0,1},PlotLabel->"Energy charge vs time",PlotLegend->{"E.C."},
PlotStyle->{Hue[0.5],Hue[0.6],Hue[0.8]},PlotRange->{0.7,0.9},
LegendPosition->{1.1,-.4},LegendSize->{.5,.35}];


atp[0]/.solution;
atp[tfinal]/.solution;
adp[0]/.solution;
adp[tfinal]/.solution;
amp[0]/.solution;
amp[tfinal]/.solution;

(atp[0]+adp[0]/2)/(atp[0]+adp[0]+amp[0])/.solution;
(atp[tfinal]+adp[tfinal]/2)/(atp[tfinal]+adp[tfinal]+amp[tfinal])/.solution;

graph1=Plot[{atp[t]/.solution,adp[t]/.solution,amp[t]/.solution},{t,0,1},PlotStyle->{Hue[0.5`],Hue[0.6`],Hue[0.8`]},PlotRange->{0,2.1`}]
graph2=ParametricPlot[Evaluate[{kuse atp[t],kform adp[t]}/.solution],{t,0,tfinal},PlotStyle->{Hue[0.7`]},PlotRange->{{9,15},{9,15}}]
graph3=Plot[{((atp[t]+((adp[t])/2))/(atp[t]+adp[t]+amp[t]))/.solution},{t,0,1},PlotStyle->{Hue[0.5`],Hue[0.6`],Hue[0.8`]},PlotRange->{0.7`,0.9`}]

here we export the graphs

Export["Section8.3.1.A.eps",graph2];
Export["Section8.3.1.B.eps",graph1];
Export["Section8.3.1.C.eps",graph3];
