0

我正在使用 R 和 desolve 包来实现一个模型，其中包含在特定时间发生的事件和在某些条件为真时发生的事件（即根事件）。

正如这里的一些背景代码演示了使用满足条件时触发的根事件，而这里的代码演示了使用多个非根事件。

正如预期的那样，如果使用根事件，则时间步长将不是触发器，反之亦然，如果使用时间触发器，则根不会触发。

我考虑使用当模型时间与事件时间匹配时返回 0 的根。我已经在下面实现了它，它似乎可以正常工作，但是我的理解是，在构建模型时，我们不应该假设模型时间会发生（实际上我理解这是使用事件背后的全部原因）。因此，我不确定这是否是好的做法，或者是否有更好的方法来做到这一点。

library("desolve")
yini <- c(temp = 18, heating_on = 1, eventoccurs = 0)

temp <- function(t,y, parms) {
  dy1 <- ifelse(y[2] == 1, 1.0, -0.5)
  dy2 <- 0
  dy3 <- 0
  list(c(dy1, dy2, dy3))
}
event_times <- c(1, 5, 20)
rootfunc <- function(t, y, parms) {
  time_in_event <- 1
  if(t %in% event_times){ time_in_event = 0}
  yroot <- c(y[1] - 18, y[1] - 20, time_in_event)
  return(yroot)
}

eventfunc <- function(t, y, parms) {
  time_in_event <- 1
  if(t %in% event_times){ time_in_event = 0}
  yroot <- c(y[1] - 18, y[1] - 20, time_in_event)
  whichroot <- which(abs(yroot) < 1e-6) # specify tolerance
  if(whichroot == 2) { y[2] <- 0 } else { y[2] <- 1 }
  if(whichroot == 3) { y[3] <- 1 } else { y[3] <- 0 }
  return(y)
}

times <- seq(from=0, to=20,by=0.1)
out <- lsode(times=times, y=yini, func = temp, parms = NULL, 
             rootfun = rootfunc, events = list(func=eventfunc, root = TRUE))
plot(out, lwd=2)

扩展 -
如果这是解决方案，我预计同时发生的事件我可能会继续遇到多个事件同时发生的问题，但从实验来看，提供相同的参数似乎不会被多个同时发生的事件修改，这不是一个问题。我认为处理这些实例的正确方法是记录优先级顺序，理想情况下，打印/记录已发生的警告（显然对于时间触发事件，可以在模型运行之前进行检查，但根触发器只会在运行过程中被发现）。

4

1 回答 1

1

在deSolve: Can't understand how to early stop odesolver with root functions我勾勒了 root->action 事件机制的一般概念。这意味着desolve如果检测到根函数的一个组件中的符号更改，则主动搜索根。因此，根函数的结果应该是在根处改变符号的连续函数的列表。对于固定时间，它应该是t-event_time(k).

为了避免代码重复，从而获得更高的灵活性，eventfunc只需使用yroot=rootfunc(t,y,parms).

于 2022-02-06T17:12:56.260 回答