What's the proper way to implement a Heaviside step function that's part of an deSolve
ODE?
My first guess would be to use ifelse
s.t.
ode <- function(t, y, p) {
dy = numeric(1)
H = ifelse((y[1] - 1) >= 0, 1, 0)
dy[1] = -p*y[1] + (y[1] - 1)*H
list(dy)
}
Is that the proper way to do it or should I use the events
option in deSolve
by using a root function? As far as I understood it, Matlab would need me to use a root function and split the integration.