[R-sig-dyn-mod] Help with function

Hi Austin,

one possibility to do this is the code below. In this example code, I
randomly initialize your parameters and the initial state. Instead of
using the indexes to access parameter values or states you can use the
names as well (see ?ode, first example). However, in my experience, the
function with() slows down your code considerably.

Best wishes,
Daniel


library(deSolve)

myfn <- function(t, y, p) {

  dy <- numeric(6)
  dy[1] <- p[1] - p[2]*y[4] - p[3]*y[1]
  dy[2] <- y[5]*dy[1]*(1+y[6]) - p[5]*y[2] - (p[6]*y[3]*y[4]*y[1] - p[7])
  dy[3] <- p[8] - p[9]*y[2]
  dy[4] <- y[6]*dy[1] - p[10]*y[4] + p[11]*(p[6]*y[3]*y[4]*y[1] - p[7]) +
p[4]
  dy[5] <- 1.0 - y[5]*(y[1] + y[4] + 1)
  dy[6] <- 1.0 - y[6]*(y[1] + 1)

  return(list(dy))

}

pars <- runif(11)
yini <- runif(6)
times <- seq(0, 10, .01)


out <- ode(yini, times, myfn, pars)
plot(out)


On Mo, 2015-02-16 at 13:15 -0700, Austin Mullings wrote:

I am trying to get this equation to work with deSolve, and can't seem to
get it to work. Does someone know how one would write the function for

this?

dY1(t)/dt = a - bY3(t) Y4(t) - cY1(t),
dY2(t)/dt = Y5(t) [dY1(t)/dt] [1 + Y6(t)] - eY2(t) - {f Y3(t) Y4(t)

Y1(t) -

g},
dY3(t)/dt = h - iY2(t),
dY4(t)/dt = Y6(t) [dY1(t)/dt] - jY4(t) + k {f Y3(t) Y4(t) Y1(t) - g} + d,
dY5(t)/dt = 1.0 - Y5(t) [Y1(t) + Y4(t) + 1],
dY6(t)/dt = 1.0 - Y6(t) [Y1(t) + 1].

      [[alternative HTML version deleted]]

--
Daniel Kaschek
Institute of Physics
Freiburg University

Room 210
Phone: +49 761 2038531