r deSolve - 绘制时间演化 pde

Question

假设我们有一个描述变量 y(t,x) 随时间 t 和 space x 的演变的 pde，我想在三维图 (t ,x,y)。使用 deSolve 我可以解决 pde，但我不知道如何获得这种图表。

deSolve 包指令中的示例如下，其中 y 是蚜虫，t=0,...,200 和 x=1,...,60:

library(deSolve)

Aphid <- function(t, APHIDS, parameters) {
deltax <- c (0.5, rep(1, numboxes - 1), 0.5)
Flux <- -D * diff(c(0, APHIDS, 0)) / deltax
dAPHIDS <- -diff(Flux) / delx + APHIDS * r
list(dAPHIDS )
}

D <- 0.3 # m2/day diffusion rate
r <- 0.01 # /day net growth rate
delx <- 1 # m thickness of boxes
numboxes <- 60
Distance <- seq(from = 0.5, by = delx, length.out = numboxes)

APHIDS <- rep(0, times = numboxes)
APHIDS[30:31] <- 1
state <- c(APHIDS = APHIDS) # initialise state variables

times <-seq(0, 200, by = 1)
out <- ode.1D(state, times, Aphid, parms = 0, nspec = 1, names = "Aphid")

"out" 生成一个包含我们需要的所有数据的矩阵，t, y(x1), y(x2), ... y(x60)。如何生成曲面图以显示 (t,x) 中 y 的演变和可变性？

Answer 1

根据使用的包，方式会有所改变。但是你可以用很少的成本做到这一点，因为 out[,-1] 是绘制曲面的理想矩阵形式。我展示了两个使用 rgl 和 plot3D 包的例子。

out2 <- out[,-1]
AphID  <- 1:ncol(out2)

library(rgl)
persp3d(times, AphID, out2, col="gray50", zlab="y")
 # If you want to change color with value of Z-axis
 # persp3d(times, AphID, out2, zlab="y", col=topo.colors(256)[cut(c(out2), 256)])

library(plot3D)
mat <- mesh(times, AphID)
surf3D(mat$x, mat$y, out2, bty="f", ticktype="detailed", xlab="times", ylab="AphID", zlab="y")

r deSolve - 绘制时间演化 pde

r deSolve - plotting time evolution pde

r

surface

pde