model <- function(t, state, parms) { with(as.list(c(state,parms)), { dR <- -c*R*N dN <- -d*N/(h+R) return(list(c(dR, dN))) }) } p <- c(c=0.1,d=0.005,h=0) s <- c(R=1,N=0.01) after <- "if(t %% 30==0){cat(state[2],\"\n\");state[1]<-1;state[2]<-state[2]*10}" run(300,after=after)