An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-dynamic-models/attachments/20110112/891eaa4e/attachment.pl>
[R-sig-dyn-mod] Events in deSolve
2 messages · Theo Nayler, Soetaert, Karline
Theo,
The problem is that, in the current implementation of events, the "event
times" must also be present in your "output times".
And the solvers do not check this - so in your case the events are
simply ignored.
We will change that in the next release of deSolve, i.e. probably such
that the solver will give an error message when this is the case.
In your case, if you do:
times <- 0:96
y <- ode(init_cond, times, model, params, events = list(data =
eventdat))
you will see the effect of your "events"
while
y[times %in% c(0, 74, 78, 84, 96), ]
will give you the output only at requested times.
Karline
-----Original Message-----
From: r-sig-dynamic-models-bounces at r-project.org
[mailto:r-sig-dynamic-models-bounces at r-project.org] On Behalf Of Theo
Nayler
Sent: Wednesday, January 12, 2011 11:26 AM
To: r-sig-dynamic-models at r-project.org
Subject: [R-sig-dyn-mod] Events in deSolve
Hi
I have a question regarding to the implementation of events in deSolve.
I
wish to execute a number of events at specific times but obtain an
output
from the solver at different times. However, I am having trouble doing
this
and wonder if it is possible? As an example, the following does not seem
to
implement the events:
library("deSolve")
model<-function(t, R, params)
{
# Parameters
ka<-params[1]
kel<-params[2]
k12<-params[3]
k21<-params[4]
# Equations
dR<-vector(len = 3)
dR[1]<--R[1]*ka
dR[2]<-R[1]*ka-R[2]*k12+R[3]*k21-R[2]*kel
dR[3]<-R[2]*k12-R[3]*k21
list(dR)
}
params<-c(1.2, 0.09, 1.6, 1)
init_cond<-c(v1 = 1000, v2 = 0, v3 = 0)
eventdat<-data.frame(var = c(1, 1, 1, 1),
time = c(4, 24, 48, 72),
value = c(1000, 1000, 1000, 1000),
method = c("add", "add", "add", "add"))
times<-c(0, 74, 78, 84, 96)
y<-lsoda(init_cond, times, model, params, events = list(data =
eventdat))
I would be grateful for any suggestions. I am new to R so apologise if I
am
missing something obvious!
Kind regards
Theo
_______________________________________________
R-sig-dynamic-models mailing list
R-sig-dynamic-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models