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

Fred Jean fjean at univ-brest.fr
Thu Oct 20 23:24:03 CEST 2011


Dear Thomas and Karline

Many thanks for your avices and for giving attention to my newbie's 
problems.

I didn't progress a lot those last days (had many courses).
At this point, before those last two emails, I followed the advice of 
Karline and here is a piece of code. As you will see, I commented the 
code from what I understood. I will now try to take into account advice 
from Thomas.

Many thanks again to both of you

Fred

##============================================
require(simecol)

rm(list=ls())
# external parameters

# time for constructing functions and for length of simulations
t1  = seq(from = 0, to = 800, by = 1)

# temperature forcing function parameters
Tm <- 12       # mean temperature
TAmpl <- 6     # temperature variation amplitude
TPer <- 365    # temperature sine period
TFaz <- 365/6  # temperature sine phase

#substrate forcing
Sm <- 0.5      # mean substrate input
SAmpl <- 0.5   # substrate input amplitude
SPer <- 365/2  # substrate input sine period
SFaz <- 365/12 # substrate input sine phase

#### The model
lvfred <- odeModel(
   main = function(time, init, parms, inputs) {

     # interpolate forcings - output of those is one interpolated value
     # at each integration step
     s.in <- approxTime1(inputs, time, rule = 2)["s.in"]
     t.in <- approxTime1(inputs, time, rule = 2)["t.in"]

     # model equations
     with(as.list(c(init, parms)),{
         ds <- s.in  - t.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, t.in)
       })
   },

   # values of parameters
   parms = c(b = 0.1, c = 0.1, d = 0.1, e = 0.1, f = 0.1, g = 0.05),

   # duration of simulations and integration timestep
   times  = c(from = 0, to = 800, by = 0.1),

   # One may add sub-equations in a file
   # "equations.file.R" -- see daphnia example in simecolModel
   # equations = equations.file,

   inputs = as.matrix(
     data.frame(
       #time = c(0,   99, 100,  150, 200, 365),
       #s.in = c(0.1, 0.2, 0.2, 0.2, 0.2, 0.1)
       time = t1,
       s.in = Sm - SAmpl * sin((2*pi*t1 / SPer) - SFaz),
       t.in = Tm - TAmpl * sin((2*pi*t1 / TPer) - TFaz)
       )
     ),

   init = c(s = 1, p = 1, k = 1), # substrate, producer, consumer

   solver = "rk4"
   )
##============================================



More information about the R-sig-dynamic-models mailing list