ODE roots

## =============================================================================
## R-Code to solve example 3.4.2 from the book
## K. Soetaert, J. Cash and F. Mazzia, 2012.  
## Solving differential equations in R. UseR, Springer, 248 pp. 
## http://www.springer.com/statistics/computational+statistics/book/978-3-642-28069-6.
## implemented by Karline Soetaert
## =============================================================================

## An ordinary differential equation model with root and events
library(deSolve)
yini  <- c(height = 0, velocity = 10)

ball <- function(t, y, parms) {
  dy1 <- y[2]
  dy2 <- -9.8     

  list(c(dy1, dy2))
}

rootfunc <- function(t, y, parms) y[1]

eventfunc <- function(t, y, parms) {
 y[1] <- 0
 y[2] <- -0.9*y[2]
 return(y)
}

times <- seq(from = 0, to = 20, by = 0.01)
out <- ode(times = times, y = yini, func = ball,
             parms = NULL, rootfun = rootfunc,
             events = list(func = eventfunc, root = TRUE))

plot(out, which = "height", lwd = 2,
    main = "bouncing ball", ylab = "height")

Go back