deSolve包ode函数为何在time=0时额外触发事件？如何修复？

使用deSolve包的ode函数时，明明只设置了t=2触发事件，结果事件却触发了两次（t=0和t=2各一次），最终计数器q的值为2而非预期的1。最小复现代码如下：

library(deSolve)

q <- 0

# 定义ODE系统
ode_fun <- function(t, y, parms) {
  dy_dt <- -y # y' = -y
  return(list(dy_dt))
}

# 定义事件函数
event_fun <- function(t, y, parms) {
  print(t)
  q <<- q +1
  y <- y + 1
}

# 初始条件与时间点
y0 <- 0.5 # y的初始值
times <- seq(0, 5, by = 0.1) # ODE求解的时间点

# 参数
parms <- NULL

# 求解带事件的ODE
ode_sol <- ode(y = y0, times = times, func = ode_fun, parms = parms,
               events = list(func = event_fun, times = c(2)))

# 绘图
plot(ode_sol, xlab = "time", ylab = "y")
print(q)

原因分析 #

这是deSolve包的默认行为：ode函数在处理事件时，会默认在初始时间点（t=0）调用一次事件函数，即便你没有将0加入events$times列表。这个设计是为了让用户有机会在模拟开始时调整初始条件，但会导致未预期的额外触发。

修复方法 #

方法一：在事件函数内添加时间判断，仅在目标时间执行逻辑 #

修改event_fun，增加t != 0的判断，同时补上事件函数必须返回修改后状态的逻辑（原代码遗漏了这一步，会导致事件触发后y的变化不生效）：

event_fun <- function(t, y, parms) {
  if (t != 0) {
    print(t)
    q <<- q +1
    y <- y + 1
  }
  return(y)
}

方法二：设置events的init=FALSE参数，禁用初始时间的事件触发 #

在调用ode时，给events列表添加init=FALSE，从根源上关闭初始时间点的默认事件调用：

ode_sol <- ode(y = y0, times = times, func = ode_fun, parms = parms,
               events = list(func = event_fun, times = c(2), init = FALSE))

内容的提问来源于stack exchange，提问作者Oscar Delaney