如何在R中求解和绘制微分方程？答案

【问题标题】：how to solve and plot a differential equation in R?如何在R中求解和绘制微分方程？
【发布时间】：2020-03-31 21:04:37
【问题描述】：

我想求解并绘制指数增长的微分方程，但我不太明白如何使用 deSolve 库。我的等式是 N = N_0 * e^(rt) 而我尝试的代码是

library(deSolve)

## Time
t <- seq(0, 5, 1)

## Initial population
N0 <- 2

## Parameter values
r = 1

fn <- function(t, N0, r) with(r, list(N0 * exp(r*t)))

## Solving and ploting 
out <- ode(N0, t, fn, params)
plot(out, lwd=2, main="exp")

但我希望的输出不是我想要的。我想获得的图表如下：

我希望你能帮助我。谢谢

【问题讨论】：

标签： r plot differential-equations

【解决方案1】：

模型函数fn 应该包含导数，然后由求解器完成积分。一阶增长当然可以解析求解，但对于更复杂的模型并不总是可行的。

library(deSolve)

## == derivative ==
fn <- function(t, N, r) {
  # dN/dt = r * N
  list(r * N)
}

r <- 1       # Parameter value
N <- 0:100   # sequence of N
t <- 0       # dummy as the derivative is not time dependent

plot(N, fn(t, N, r)[[1]], type="l")

## == integration ==
t <- seq(0, 5, .1)  # time
N0 <- 2             # initial state

## numerical solver
out <- ode(N0, t, fn, r)
plot(out, lwd=2, main="exp")

## for comparison: analytical integration
lines(t, N0*exp(r*t), lwd=2, lty="dotted", col="red")

【讨论】：

【解决方案2】：

您也可以尝试curve 函数。

op <- par(mfrow=c(1, 2), mar=c(5, 5, 4, 3))
curve(r*x, from=0, to=100, xlab="N", ylab=bquote(dot(N)), main=bquote(dot(N)==N))
curve(N0 * exp(r*x), from=0, to=5, xlab="Time t", ylab="N(t)", main="Exponential growth")
par(op)

【讨论】：