param tau = 1.0

param P0 = 93
param R = 1.05
param r = 0.0068
param Vstr = 67.9
param ca = 1.55
param cv = 519
param alpha0 = 93
param alphas = 93
param alphap = 93
param alphaH = 0.84
param beta0 = 7
param betas = 7
param betav = 7
param betaH = 1.17
param gammaH = 0

state Pa=P0 delay
state Pv = P0/(1 + R/r)
state H = P0/(1 + R/r)/(R*Vstr)

aux Ts = 1.0/(1 + (Pa(t - tau)/alphas)^betas)
aux Tp = 1.0/(1 + (alphap/Pa)^betav)

Pa'= -Pa/(ca*R) + Pv/(ca*R) + H*Vstr/ca
Pv'= Pa/(cv*R) - (1/R + 1/r)*Pv/cv
H'= alphaH*Ts/(1 + gammaH*Tp) - betaH*Tp

print t, Pa, Pv, H

option tstart=0.0, tstop=350, output=350