[R-sig-dyn-mod] Modeling dampened oscillator with random resets

[Gesmann, Markus]

Hi Felix,

With simecol, which builds on deSolve you could do the following (based on an example Thomas Petzoldt send me a couple of years ago):

library(simecol)

osci <- odeModel( main=function(time, init, parms,...){
  x <- init
  p <- parms
  ## same as above, but now we have and input, or external stimulus
  S <- approxTime1(inputs, time, rule=2)["s.in"]
  dx <- x[2]
  dy <- parms["eta"] * x[1] + parms["zeta"] * x[2] + S
  list(c(dx, dy))
}, 
                  times=c(from=0, to=200, by=0.2),
                  parms=c(eta=-0.05, zeta=-0.02, x=1, y=0),
                  initfunc = function(obj) {
                    init(obj)[c("x", "y")] <- parms(obj)[c("x", "y")]
                    return(obj)
                  }
)


## Let's give the system a kick at 70 and 120:
 
kick <- as.matrix(data.frame(time=1:200, 
                             s.in=c(rep(0,69), 
                                    1, 
                                    rep(0,49), 1, 
                                    rep(0,80))))
inputs(osci) <- kick

mySim <-  out(sim(osci))
plot(y ~ time, data=mySim, t="l")


I hope this helps.

Best regards

Markus

-----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 Felix Schönbrodt
Sent: 09 January 2013 14:33
To: r-sig-dynamic-models at r-project.org
Subject: [R-sig-dyn-mod] Modeling dampened oscillator with random resets

Dear list,

using deSolve (or another package): is it possible to model a dampened oscillator that is reset to new states at specific points?

Background: I want to model a system that tries to approach equilibrium, but is disturbed at single events. It would roughly look like following screenshot:
https://dl.dropbox.com/u/4472780/WP/oscReset.PNG

In this screen shot, events at t=70 and t=120 disturb the oscillator. In this screen shot I simply pasted together 3 independent plots; but if the disturbance is directly modeled it should somehow be reflected in the dx etc. (e.g., at the second event I would expect that the curve first continues going up before it is regulated down again).

Here's my function for the dampened oscillator:

library("deSolve")

	osci <- function (t, states, parms) {
	 with(as.list(c(states,  parms)), {
	   dx <- y
	   dy <- eta * x + zeta * y
	   list(c(dx, dy))
	 })
	}

	states <- c(x=1, y=0)
	out <- data.frame(ode(states, osci, times = 1:100, parms = c(eta = -0.05, zeta = -0.02)))
	plot(out)


Best,
Felix
_______________________________________________
R-sig-dynamic-models mailing list
R-sig-dynamic-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models
**********************************************************************
The information in this E-Mail and in any attachments is...{{dropped:27}}