#file pulse_generator.ode
#model of synthetic pulse generating system
#adapted from Basu et al. (2004) PNAS 101 6355-6360
#Figure 7.28

par aG=80
par bG=0.07 
par KC=0.008
par aC=0.5 
par bC=0.3 
par k1=0.5
par k2=0.02 
par KR=0.02
par RT=0.5 

A=10*heav(t-10.001)

#model equations
G' = -bG*G + aG*((R/KR)/(1+(R/KR) + (C/KC)^2 + (R/KR)*(C/KC)^2))
C' = -bC*C + aC*R/(KR+R)
R' = -k2*R + k1*(RT-2*R)^2*A^2 


# initial condition
init G=0, C=0, R=0

# final time
@ total=50
# set the plotting window size:
@ xlo=0, xhi=50, ylo=0, yhi=2
@ dt=0.05
@ bound=100000

# run the simulation in the GUI by selecting (I)nitialconds|(G)o
# plot species C and R by selecting (G)raphic stuff|(A)dd curve and 
# assigning the y-axis

done