Diffusion based on direction and velocity raster

diffusion.R

# Edit: whoops, radians

library(raster)

direction <- raster(nrows=100, ncols=100, xmn=0, xmx=10, , ymn=0, ymx=10)
values(direction) <- seq(0, 2 * pi, length.out = 100^2) + rnorm(100^2)

velocity <- raster(nrows=100, ncols=100, xmn=0, xmx=10, , ymn=0, ymx=10)
values(velocity) <- runif(100^2, 0, 0.1)

plot(direction)

plot(velocity)


# Number of time steps
N <- 10000

# empty lat lon matrix
fad <- matrix(NA, nrow = N, ncol = 2)

# start at middle
fad[1, ] <- c(5, 5)

for(hr in 2:N){
  
  dir <- extract(direction, fad[hr - 1, , drop = FALSE])
  vel <- extract(velocity, fad[hr - 1, , drop = FALSE])

  fad[hr, 1] <- fad[hr - 1, 1] + vel * cos(dir)
  fad[hr, 2] <- fad[hr - 1, 2] + vel * sin(dir)


}


plot(direction)
lines(fad[, 1] ~ fad[, 2])