[R-sig-dyn-mod] Damped oscillator with deSolve

On 12/6/2012 8:02 PM, Suzen, Mehmet wrote:

Hello Felix,

You need to introduce new variable y=dx/dt  and your equation reads
dy/dt= eta*x + zeta*y.  So solution
will be for states (x, y) and paramaters (eta, zeta). Remember that
any higher order ODE can be expressed as first order.

Yes, right, cf. page 3+4 in:

http://journal.r-project.org/archive/2010-2/RJournal_2010-2_Soetaert~et~al.pdf


so you may get something like:

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 <- ode(states, osci, times = 0:100,
 parms = c(eta = 0.2, zeta = 0.1))
plot(out)



Have fun!

Thomas P.