model <- function(t, state, parms) {
  with(as.list(c(state,parms)), {
    b <- b0 + e*(sin(2*3.1416*(t-Delta)/365))
    dA  <- b*A*(1-A/k) - d*A
    return(list(dA))
  }) 
}  

p <- c(b0=1,e=0.5,d=0.5,k=1,Delta=100)
s <- c(A=1e-6)
run(10*365)