[R-sig-dyn-mod] forcing variables (dummy ?) problem

[Thomas Petzoldt]

Dear Fred,

during object creation, the values of the slots are not accessible by 
each other for principle reasons and this is the case where a special 
initialization mecahism "initfunc" comes into play. The "initfunc" gets 
the raw incomplete version of the object (called object in creation) and 
one has all methods to access and manipulate the slots.

See package vignette or http://www.jstatsoft.org/v22/i09 for an example.

For the "lvfred" example, different possibilities exist (see below), but 
for your other problem you don't need "initfunc", just add the 
parameters of the sine function to the parms vector:

parms = c( ......, Tm = 12, Ampl = 6, Per = 365, Faz <- 365/12, ....)

and call the signal in the "main" model function:

with(as.list(c(init, parms)),{
   tempC <-  Tm - Ampl * sin((2*pi*times / Per) - Faz)
   .......
})

For more complicated cases with real data, use something like the 
following, where lvfred2 with approxfun is more flexible and faster.

Hope it helps

Thomas P.


require(simecol)

lvfred <- odeModel(
   main = function(time, init, parms, inputs) {
     s.in <- approxTime1(inputs, time, rule = 2)["s.in"]
     with(as.list(c(init, parms)),{
         ds <- s.in  - b*s*p + g*k
         dp <- c*s*p - d*k*p
         dk <- e*p*k - f*k
         list(c(ds, dp, dk), s.in)
       })
   },

   parms = c(b = 0.1, c = 0.1, d = 0.1, e = 0.1, f = 0.1, g = 0.05),

   times  = seq(from = 0, to = 365, by = 1),
   init = c(s = 1, p = 1, k = 1), # substrate, producer, consumer
   solver = "lsoda",
   initfunc = function(obj) {
     tmp <- as.matrix(data.frame(
       time = times(obj),
       s.in = c(rep(0.1,99), rep(0.2,101), rep(0.1, 166))
     ))
     inputs(obj) <- tmp
     return(obj)
   }
)

##==============================================================================


lvfred2 <- odeModel(
   main = function(time, init, parms, inputs) {
     s.in <- inputs$signal(time)
     with(as.list(c(init, parms)),{
         ds <- s.in  - b*s*p + g*k
         dp <- c*s*p - d*k*p
         dk <- e*p*k - f*k
         list(c(ds, dp, dk), s.in = s.in)
       })
   },

   parms = c(b = 0.1, c = 0.1, d = 0.1, e = 0.1, f = 0.1, g = 0.05),

   times  = seq(from = 0, to = 365, by = 1),
   init = c(s = 1, p = 1, k = 1), # substrate, producer, consumer
   solver = "lsoda",
   initfunc = function(obj) {
     tmp <- approxfun(times(obj),
                     c(rep(0.1,99), rep(0.2,101), rep(0.1, 166)))
     inputs(obj) <- list(signal = tmp)
     obj
   }
)

system.time(plot(sim(lvfred)))


system.time(plot(sim(lvfred2)))