计算科学 - R中的数值求解Lotka-Volterra ODE - 吾爱随笔录

目标：我正在尝试使用sde 包sde.sim()中的 de 函数对R 中的 Lotka-Volterra ODE 进行数值求解。我想使用该功能最终将该系统转换为 SDE。所以最初，我从一个没有噪声项的简单 ODE 系统（Lotka-Volterra 模型）开始。sde.sim()

Lotka-Volterra ODE 系统：

{\begin{matrix} \frac{d x}{d t} = α x - β x y . \\ \frac{d y}{d t} = δ x y - γ y \end{matrix}

$\left\{\begin{matrix}\frac{dx}{dt}= \alpha x -\beta xy. \\ \frac{dy}{dt}= \delta xy -\gamma y \end{matrix}\right.$

x = 10 和 y = 10 的初始值。

alpha、beta、delta 和 gamma 的参数值分别为 1.1、0.4、0.1 和 0.4（模仿这个例子）。

尝试解决问题：

library(sde)
d <- expression((1.1 * x[0] - 0.4 * x[0] * x[1]), (0.1 * x[0] * x[1] - 0.4 * x[1]))
s <- expression(0, 0)
X <- sde.sim(X0=c(10,10), T = 10, drift=d, sigma=s) 
plot(X)

然而，这似乎并没有产生捕食者和猎物种群的良好循环行为。

题：

在尝试解决这个 ODE 系统时出了什么问题sde.sim()？

library(plotly) # Lotka-Volterra Model (SDE without noise) f <- function(u,p,t) { du1 = p[1]*u[1]-p[2]*u[1]*u[2] du2 = p[3]*u[1]*u[2]-p[4]*u[2] return(c(du1,du2)) } g <- function(u,p,t) { return(c(0*u[1],0*u[2])) } u0 = c(10.0, 10.0) tspan <- list(0.0,100.0) p = c(1.1, 0.4, 0.1, 0.4) sol = diffeqr::sde.solve(f,g,u0,tspan,p=p,saveat=0.005) udf = as.data.frame(sol$u) plot_ly(udf, x = sol$t, y = ~V1, name = 'Prey', type = 'scatter', mode = 'lines') %>% add_trace(y = ~V2, name = 'Predator', mode = 'lines')