[R-sig-dyn-mod] timestep()

[Thomas Petzoldt]

Daniel,

many thanks for the report, I had a look in the sources. Functions 
lagvalue, lagderiv and timesteps are intended for quite specific 
purposes and work differently from what you may expect, so the 
documentation says "Use with care."

Function timesteps() has no information about internal time steps taken 
automatically by the odepack solvers (lsoda, ..., vode), but you can get 
more for the Runge-Kutta solvers (e.g. method="ode23"). I agree that we 
should make this more clear in the manuals.

As an alternative, one may consider to employ global variables or 
(better) a closure (see below). One may also consider handling of 
rejected steps.

Thomas



library(deSolve)

## a function closure  ('lexical scoping')
modelClosure <- function(t0) {
   t.old <- t.act <- t0
   function(t, y, parms) {
     t.old  <<- t.act
     t.act  <<- t
     cat(t, "\t", t - t.old, "\n")
     with (as.list(c(y, parms)), {
       dP <- a * P      -  b * P * K
       dK <- b * P * K  -  c * K
       list(c(dP, dK))
     })
   }
}

model <- modelClosure(0) # initialization

parms <- c(a = 0.1, b = 0.1, c = 0.1)

y <- c(P = 1, K = 2)

out <- ode(y = y, func = model, times = c(0, 2),
  parms = parms, method = "lsoda")

ls() # prove that t.old and t.new are local within 'model'