This Workbook sets up, simulates and graphs the equations that describe adenosine system with charging and recovering of high-energy bonds. It corresponds to Section 8.4 of Ch 8 in SB2.

We 1) define the numerical values of the kinetic parameters, 2) then the reaction rate laws, 3) then we integrate and, 4) finally graph the solution.

BOP May 6, 2007
Addiel de Alba, April 2010.
NS, checked on March 15, 2011 OK

SetDirectory[NotebookDirectory[]];
<<"GraphicsDirectives.m"
tfinal=50;

parameters={ka1->6.25,ka2->25.,ka5->1000,kam5->1000,b1->0.03,b2->5.,ka4->(b1/0.1),kl->1.5*(5/1.6),ks1->(5/1.6),ks2->(5/0.4)};

vs1=ks1 atp[t] s1[t];
vs2=ks2 adp[t] s2[t];
vl=kl atp[t];
va3=b1;
va4=ka4 amp[t];
va5=ka5 adp[t] adp[t];
vam5=kam5 atp[t] amp[t];

BalanceATP=-vs1-vl+2 vs2+(va5-vam5);
BalanceADP=vs1+vl-2 vs2-2 (va5-vam5);
BalanceAMP=va3-va4+(va5-vam5);
BalanceS1=b2-vs1;
BalanceS2=vs1-vs2;

solution=NDSolve[{atp'[t]==BalanceATP,adp'[t]==BalanceADP,amp'[t]==BalanceAMP,s1'[t]==BalanceS1,s2'[t]==BalanceS2,atp[0]==1.6`,adp[0]==0.4`,amp[0]==0.1`,s1[0]==1.`,s2[0]==1.`}//.parameters,{atp,adp,amp,s1,s2},{t,0,tfinal}];

graph1=Plot[Evaluate[{atp[t],adp[t],amp[t],s1[t],s2[t]}/.solution],{t,0,5},PlotLabel->"Concentrations vs time",Epilog->createLegend[{"atp","adp","amp","x1","x2"}],PlotRange->{0,1.8`}]

graph2=Plot[Evaluate[{atp[t],adp[t],amp[t],s1[t],s2[t]}/.solution],{t,0,tfinal},PlotLabel->"Concentrations vs time",Epilog->createLegend[{"atp","adp","amp","x1","x2"}],PlotRange->{0,1.8`}]

graph3=Plot[Evaluate[{ks1 atp[t] s1[t],ks2 adp[t] s2[t],kl atp[t]}/.solution//.parameters],{t,0,5},PlotLabel->"Fluxes vs time",Epilog->createLegend[{"v1","v2","vload"}],PlotRange->{4.4`,8.`}]

graph4=Plot[Evaluate[{ks1 atp[t] s1[t],ks2 adp[t] s2[t],kl atp[t]}/.solution//.parameters],{t,0,tfinal},PlotLabel->"Fluxes vs time",Epilog->createLegend[{"v1","v2","vload"}],PlotRange->{4.,8.`}]

graph5=Plot[Evaluate[{ka4 amp[t],ka5 adp[t] adp[t]-kam5 atp[t] amp[t]}/.solution//.parameters],{t,0,5},PlotLabel->"Fluxes vs time",Epilog->createLegend[{"vdrain","vdistr, net"}],PlotRange->{-0.025`,0.4`}]

graph6=Plot[Evaluate[{ka4 amp[t],ka5 adp[t] adp[t]-kam5 atp[t] amp[t]}/.solution//.parameters],{t,0,tfinal},PlotLabel->"Fluxes vs time",Epilog->createLegend[{"vdrain","vdistr, net"}],PlotRange->{-0.025`,0.4`}]

graph7=ParametricPlot[Evaluate[{ks1 atp[t] s1[t]+kl atp[t],2 ks2 adp[t] s2[t]}/.solution//.parameters],{t,0,tfinal},PlotLabel->"ATP load vs ATP synthesis",PlotRange->{{9,13},{9,13}},Epilog->ppAnnotate[{ks1 atp[#] s1[#]+kl atp[#]&,2 ks2 adp[#] s2[#]&}//.parameters,solution,tfinal]]

graph8=Plot[{((2 atp[t]+adp[t])/(2 (atp[t]+adp[t]+amp[t])))/.solution},{t,0,tfinal},PlotLabel->"Energy charge vs time",Epilog->createLegend[{"E.C."}],PlotRange->{0.65`,0.9`}]

graph9=ParametricPlot[Evaluate[{ks1 atp[t] s1[t]+kl atp[t],ka4 amp[t]}/.solution//.parameters],{t,0,tfinal},PlotLabel->"ATP load vs AMP drain",PlotRange->{{8,15},{0.02`,0.12`}},Epilog->ppAnnotate[{ks1 atp[#] s1[#]+kl atp[#]&,ka4 amp[#]&}//.parameters,solution,tfinal]]

graph10=ParametricPlot[Evaluate[{2 atp[t]+adp[t],2 (atp[t]+adp[t]+amp[t])}/.solution],{t,0,tfinal},PlotLabel->"occupancy vs capacity",Epilog->ppAnnotate[{2 atp[#]+adp[#]&,2 (atp[#]+adp[#]+amp[#])&},solution,tfinal]]

graph11=ParametricPlot[Evaluate[{ka5 adp[t] adp[t]-kam5 atp[t] amp[t],ka4 amp[t]}/.solution//.parameters],{t,0,tfinal},PlotLabel->"P distribution vs AMP drain",Epilog->ppAnnotate[{ka5 adp[#] adp[#]-kam5 atp[#] amp[#]&,ka4 amp[#]&}//.parameters,solution,tfinal]]

here we export the graphs

Export["Section8.4.A.eps",graph1];
Export["Section8.4.B.eps",graph2];
Export["Section8.4.C.eps",graph3];
Export["Section8.4.D.eps",graph4];
Export["Section8.4.E.eps",graph5];
Export["Section8.4.F.eps",graph6];
Export["Section8.4.A2.eps",graph7];
Export["Section8.4.B2.eps",graph8];
Export["Section8.4.C2.eps",graph9];