In this text file are given the .ode files used for the numeric
simulations of the system studied in the article.


File triangle.ode:
******************

a'=fa+f3-f1
b'=f1-f2
c'=f2-f3-fa

fa=ga*(vc*vx-va*vy)
f1=g1*(va-vb)
f2=g2*(vb-vc)
f3=g3*(vc-va)

vx=x/kx
vy=y/ky
va=a/ka
vb=b/kb
vc=c/kc

par x=1, y=1e-5, ka=1, kb=2, kc=4, kx=1, ky=10, ga=40, g1=1, g2=2, g3=4

aux f3p=fa+f3
aux fla=fa
aux fl1=f1
aux fl2=f2
aux fl3=f3
aux se=8.314*fa*ln(vx/vy)
aux sa=8.314*fa*ln(vc*vx/(va*vy))
aux s1=8.314*f1*ln(va/vb)
aux s2=8.314*f2*ln(vb/vc)
aux s3=8.314*f3*ln(vc/va)
aux si=8.314*fa*ln(vc*vx/(va*vy))+8.314*f1*ln(va/vb)+8.314*f2*ln(vb/vc)+8.314*f3*ln(vc/va)
aux sbal=8.314*fa*ln(vx/vy)-(8.314*fa*ln(vc*vx/(va*vy))+8.314*f1*ln(va/vb)+8.314*f2*ln(vb/vc)+8.314*f3*ln(vc/va))



File aped.ode:
**************

la'= r1+rr1-rm1-rrm1-r4+rm4-r6+rm6
da'= r2+rr2-rm2-rrm2-r5+rm5-r7+rm7
l'= rrm1+rm1-r1-rr1-r4+rm4-r5+rm5+2*r8-2*rm8+r9-rm9+r10-rm10-r3+rm3 
d'= rrm2+rm2-r2-rr2-r6+rm6-r7+rm7+r9-rm9+r10-rm10+2*r11-2*rm11+r3-rm3
ll'= r4-rm4-r8+rm8-r12+rm12
dl'= r5-rm5-r9+rm9+r12-rm12
ld'= r6-rm6-r10+rm10+r13-rm13
dd'= r7-rm7-r11+rm11-r13+rm13


r1= ka*l
rm1= kma*la
rr1= kmd*l
rrm1= kd*la
r2= ka*d
rm2= kma*da
rr2= kmd*d
rrm2= kd*da

r3= kr*l
rm3= kr*d

r4= kp*l*la
rm4= kmp*ll
r5= kp2*l*da
rm5= kmp2*dl
r6= kp2*d*la
rm6= kmp2*ld
r7= kp*d*da
rm7= kmp*dd

r8= kh*ll
rm8= kmh*l*l
r9= kh2*dl
rm9= kmh2*d*l
r10= kh2*ld
rm10= kmh2*l*d
r11= kh*dd
rm11= kmh*d*d

r12= ke2*ll
rm12= ke*dl
r13= ke2*dd
rm13= ke*ld

ka= a*x
kma= a*kfx/kfa
kd= b
kmd= kfa*vy*b
kr= r 
kp= p
kp2= al*p
kmp= vy*kfa/kfp*p
kmp2= al/ga*vy*kfa/kfp*p
kh= h
kh2= h*be
kmh= kfp*h
kmh2= be*ga*kfp*h
ke= e
ke2= e*ga

ent1=8.314*(r1-rm1)*ln(r1/rm1)
ent2=8.314*(r2-rm2)*ln(r2/rm2)
entr1=8.314*(rr1-rrm1)*ln(rr1/rrm1)
entr2=8.314*(rr2-rrm2)*ln(rr2/rrm2)
ent3=8.314*(r3-rm3)*ln(r3/rm3)
ent4=8.314*(r4-rm4)*ln(r4/rm4)
ent5=8.314*(r5-rm5)*ln(r5/rm5)
ent6=8.314*(r6-rm6)*ln(r6/rm6)
ent7=8.314*(r7-rm7)*ln(r7/rm7)
ent8=8.314*(r8-rm8)*ln(r8/rm8)
ent9=8.314*(r9-rm9)*ln(r9/rm9)
ent10=8.314*(r10-rm10)*ln(r10/rm10)
ent11=8.314*(r11-rm11)*ln(r11/rm11)
ent12=8.314*(r12-rm12)*ln(r12/rm12)
ent13=8.314*(r13-rm13)*ln(r13/rm13)

vx=x/kfx

par a=1, b=1, p=1, al=0.1, h=1, be=0.1, e=1, ga=0.1, r=0.0001, kfa=1e-6, kfp=1e-3, kfx=1e-9, x=0.5, vy=1


l(0)=1.00001
d(0)=0.99999
la(0)=1e-100
da(0)=1e-100
ll(0)=1e-100
ld(0)=1e-100
dl(0)=1e-100
dd(0)=1e-100

aux ee=(l+la+2*ll-d-da-2*dd)/(l+la+2*ll+2*ld+2*dl+2*dd+da+d)
aux ee2b=(dl-ld)/(dl+ld)
aux ee2a=(ll-dd)/(ll+dd)
aux ee1=(l+la-d-da)/(l+d+la+da)
aux ee0=(l-d)/(l+d)
aux flx=a*x*(l+d)
aux ent=ent1+ent2+entr1+entr2+ent3+ent4+ent5+ent6+ent7+ent8+ent9+ent10+ent11+ent12+ent13
aux enta=ent1+ent2+entr1+entr2
aux entp=ent4+ent5+ent6+ent7
aux entd=ent8+ent9+ent10+ent11
aux entpd=ent4+ent5+ent6+ent7+ent8+ent9+ent10+ent11
aux entr=ent12+ent13+ent3
aux ente=ent12+ent13
aux flow=-8.314*(ka*(l+d)-kma*(la+da))*ln(vy/vx)
aux aa=l+d
aux dipep=2*(ll+ld+dl+dd)
aux rapdia=(ll+dd)/(ld+dl)
aux nca=la+da
aux flow21=8.314*(ka*(l+d)-kma*(la+da))*ln(vx)
aux flow22=8.314*(kmd*(l+d)-kd*(la+da)+kmp*(ll+dd)-kp*(l*la+d*da)+kmp2*(dl+ld)-kp2*(l*da+d*la))*ln(vy)
aux flow2=8.314*(ka*(l+d)-kma*(la+da))*ln(vx)+8.314*(kmd*(l+d)-kd*(la+da)+kmp*(ll+dd)-kp*(l*la+d*da)+kmp2*(dl+ld)-kp2*(l*da+d*la))*ln(vy)
aux balener=(8.314*(ka*(l+d)-kma*(la+da))*ln(vx)+8.314*(kmd*(l+d)-kd*(la+da)+kmp*(ll+dd)-kp*(l*la+d*da)+kmp2*(dl+ld)-kp2*(l*da+d*la))*ln(vy))-(ent1+ent2+entr1+entr2+ent3+ent4+ent5+ent6+ent7+ent8+ent9+ent10+ent11+ent12+ent13)

@ total=50000, yp=ee, meth=qualrk, dt=50, xhi=100000, bound=1000