r语言 - 我的二阶导数非齐次函数出了什么问题?

Question

我试图为微分方程编写一个函数

y'' + y = bcos(omega *t)

这是我的代码

library(deSolve)
yini <- c(y1 = 2, y2 = 0)
nonvdp1 <- function(t, y, parms) {
b <- parms['b']
omega <- parms['omega']
dy1dt <- y[2]
dy2dt <- b * cos(omega * t) - y[1]
list(c(dy1dt, dy2dt))
}
output <- as.data.frame(ode(y = yini, func = nonvdp1, times = (0, 30, 0.1), parms = c(2, 2)))

但是，解决方案并不完全正确

head(nonvdp1snl)
time y1 y2
1    0  2  0
2    1 NA NA
3    2 NA NA
4    3 NA NA
5    4 NA NA
6    5 NA NA

根据初始条件和我选择的参数，解决方案是

y = (8/3)cos(t) - (2/3)cos(2t)

我的代码做错了什么？

Answer 1

您的代码有几个问题，忘记了时间步长的seq，缺少参数的名称。与以下内容进行比较：

library(deSolve)
yini <- c(y1 = 2, y2 = 0)
nonvdp1 <- function(t, y, parms) {
b     <- parms['b']
omega <- parms['omega']
dy1dt <- y[2]
dy2dt <- b * cos(omega * t) - y[1]
list(c(dy1dt, dy2dt))
}
output <- ode(y = yini, func = nonvdp1, times = seq(0, 30, 0.1), 
parms = c(b = 2, omega = 2))
output
plot(output)

如果删除as.data.frame()则可以使用内置的绘图函数plot.deSolve。