model <- function(t, state, parms) {
  with(as.list(c(state,parms)), {  
    dP <- a*P*(1 - P/K) - b*P*D/(h+P)
    dD <- b*P*D/(h+P) - d*D    
    return(list(c(dP, dD)))  
  }) 
}  

p <- c(a=0,K=0,h=0,b=0,d=0)
plane()
s <- c(P=0,D=0)
noise <- "parms[\"K\"]<-abs(rnorm(1,1,0.1));state<-ifelse(state<0.01,0,state)"
run(after=noise)
