diffusion: A Random Walk Particle Diffusion Model

Description Usage Format Details See Also Examples

Description

simecol example: This is a random walk (basic particle diffusion) model.

Usage

1

Format

An S4 object according to the rwalkModel specification. The object contains the following slots:

main

A function with the movement rules for the particles.

parms

A list with the following components:

ninds

number of simulated particles,

speed

speed of the particles,

area

vector with 4 elements giving the coordinates (left, bottom, right, top) of the coordinate system.

times

Simulation time (discrete time steps, by-argument ignored).

init

Data frame holding the start properties (Cartesian coordinates x and y and movement angle a) of the particles.

Details

To see all details, please have a look into the implementation.

See Also

sim, parms, init, times.

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
##============================================
## Basic Usage:
##   explore the example
##============================================
## Not run: 
data(diffusion)
## (1) minimal example
plot(sim(diffusion))
## show "grid of environmental conditions"
image(inputs(diffusion))

## (2) scenario
##     with homogeneous environment (no "refuge" in the middle)
no_refuge <- diffusion # Cloning of the whole model object
inputs(no_refuge) <- matrix(1, 100, 100)
plot(sim(no_refuge))
  
##============================================
## Advanced Usage:
##   Assign a function to the observer-slot.
##============================================
observer(diffusion) <- function(state, ...) {
  ## numerical output to the screen
  cat("mean x=", mean(state$x),
      ", mean y=", mean(state$y),
      ", sd   x=", sd(state$x),
      ", sd   y=", sd(state$y), "\n")
  ## animation
  par(mfrow=c(2,2))
  plot(state$x, state$y, xlab="x", ylab="y", pch=16, col="red", xlim=c(0, 100))
  hist(state$y)
  hist(state$x)
  ## default case: return the state --> iteration stores it in "out"
  state
}

sim(diffusion)

## remove the observer and restore original behavior
observer(diffusion) <- NULL
diffusion <- sim(diffusion)

## End(Not run)

##============================================
## Implementation:
##   The code of the diffusion model.
##   Note the use of the "initfunc"-slot.
##============================================
diffusion <- rwalkModel(
  main = function(time, init, parms, inputs = NULL) {
    speed   <- parms$speed
    xleft   <- parms$area[1]
    xright  <- parms$area[2]
    ybottom <- parms$area[3]
    ytop    <- parms$area[4]

    x <- init$x  # x coordinate
    y <- init$y  # y coordinate
    a <- init$a  # angle (in radians)
    n <- length(a)

    ## Rule 1: respect environment (grid as given in "inputs")
    ## 1a) identify location on "environmental 2D grid" for each individual
    i.j <- array(c(pmax(1, ceiling(x)), pmax(1, ceiling(y))), dim=c(n, 2))

    ## 1b) speed dependend on "environmental conditions"
    speed <- speed * inputs[i.j]

    ## Rule 2: Random Walk
    a  <- (a + 2 * pi / runif(a)) 
    dx <- speed * cos(a)
    dy <- speed * sin(a)
    x  <- x + dx
    y  <- y + dy

    ## Rule 3: Wrap Around
    x <- ifelse(x > xright, xleft, x)
    y <- ifelse(y > ytop, ybottom, y)
    x <- ifelse(x < xleft, xright, x)
    y <- ifelse(y < ybottom, ytop, y)
    data.frame(x=x, y=y, a=a)
  },
  times  = c(from=0, to=100, by=1),
  parms  = list(ninds=50, speed = 1, area = c(0, 100, 0, 100)),
  solver = "iteration",
  initfunc = function(obj) {
    ninds   <- obj@parms$ninds
    xleft   <- obj@parms$area[1]
    xright  <- obj@parms$area[2]
    ybottom <- obj@parms$area[3]
    ytop    <- obj@parms$area[4]
    obj@init <- data.frame(x = runif(ninds) * (xright - xleft) + xleft,
                           y = runif(ninds) * (ytop - ybottom) + ybottom,
                           a = runif(ninds) * 2 * pi)
    inp <- matrix(1, nrow=100, ncol=100)
    inp[, 45:55] <- 0.2
    inputs(obj) <- inp
    obj
  }
)

Example output

Loading required package: deSolve
mean x= 54.19871 , mean y= 46.40188 , sd   x= 27.25963 , sd   y= 28.72472 
mean x= 54.2489 , mean y= 46.39458 , sd   x= 27.3246 , sd   y= 28.84305 
mean x= 54.11654 , mean y= 48.40917 , sd   x= 27.20912 , sd   y= 29.00405 
mean x= 54.00115 , mean y= 46.36714 , sd   x= 27.10779 , sd   y= 28.77226 
mean x= 53.94108 , mean y= 46.37184 , sd   x= 27.01538 , sd   y= 28.69721 
mean x= 54.10594 , mean y= 48.34082 , sd   x= 27.13034 , sd   y= 28.75588 
mean x= 54.22101 , mean y= 48.28208 , sd   x= 27.14892 , sd   y= 28.81557 
mean x= 54.0377 , mean y= 48.42046 , sd   x= 27.26096 , sd   y= 28.7216 
mean x= 53.85564 , mean y= 48.40733 , sd   x= 27.20226 , sd   y= 28.6194 
mean x= 53.78889 , mean y= 48.44337 , sd   x= 27.19725 , sd   y= 28.51433 
mean x= 53.87705 , mean y= 48.37381 , sd   x= 27.23843 , sd   y= 28.53236 
mean x= 53.76502 , mean y= 48.40391 , sd   x= 27.23742 , sd   y= 28.48247 
mean x= 53.76826 , mean y= 46.25374 , sd   x= 27.28885 , sd   y= 28.26628 
mean x= 53.865 , mean y= 48.3025 , sd   x= 27.24538 , sd   y= 28.37567 
mean x= 53.93039 , mean y= 46.4544 , sd   x= 27.24505 , sd   y= 28.12704 
mean x= 53.91797 , mean y= 48.51895 , sd   x= 27.29645 , sd   y= 28.27504 
mean x= 53.87381 , mean y= 46.66766 , sd   x= 27.25536 , sd   y= 28.05759 
mean x= 53.80318 , mean y= 48.65573 , sd   x= 27.21248 , sd   y= 28.27801 
mean x= 55.83864 , mean y= 46.56629 , sd   x= 26.74242 , sd   y= 28.23036 
mean x= 55.74128 , mean y= 44.62929 , sd   x= 26.75126 , sd   y= 27.9048 
mean x= 53.83889 , mean y= 44.82057 , sd   x= 27.16613 , sd   y= 27.76931 
mean x= 53.87438 , mean y= 47.02423 , sd   x= 27.213 , sd   y= 28.06293 
mean x= 53.90199 , mean y= 47.05311 , sd   x= 27.23508 , sd   y= 28.12601 
mean x= 53.81766 , mean y= 47.04848 , sd   x= 27.20432 , sd   y= 27.88752 
mean x= 53.83352 , mean y= 44.99247 , sd   x= 27.18976 , sd   y= 27.58792 
mean x= 53.86205 , mean y= 46.9181 , sd   x= 27.23375 , sd   y= 27.92795 
mean x= 53.88615 , mean y= 46.82183 , sd   x= 27.2864 , sd   y= 28.01437 
mean x= 55.87826 , mean y= 46.801 , sd   x= 27.00056 , sd   y= 27.90889 
mean x= 55.81301 , mean y= 46.77742 , sd   x= 27.00116 , sd   y= 27.91725 
mean x= 53.78524 , mean y= 46.67522 , sd   x= 27.36536 , sd   y= 27.94636 
mean x= 53.7243 , mean y= 46.59551 , sd   x= 27.41403 , sd   y= 28.08201 
mean x= 53.5751 , mean y= 46.73473 , sd   x= 27.5323 , sd   y= 28.14505 
mean x= 55.51996 , mean y= 44.6248 , sd   x= 27.373 , sd   y= 27.91023 
mean x= 55.42144 , mean y= 46.56 , sd   x= 27.37684 , sd   y= 28.31986 
mean x= 55.48115 , mean y= 46.52226 , sd   x= 27.54662 , sd   y= 28.37198 
mean x= 55.5142 , mean y= 46.48571 , sd   x= 27.57148 , sd   y= 28.25563 
mean x= 55.58402 , mean y= 46.59364 , sd   x= 27.73772 , sd   y= 28.21867 
mean x= 55.68952 , mean y= 46.55852 , sd   x= 27.61308 , sd   y= 28.12712 
mean x= 55.78867 , mean y= 46.59295 , sd   x= 27.5476 , sd   y= 28.20345 
mean x= 55.85365 , mean y= 44.54939 , sd   x= 27.46571 , sd   y= 27.77314 
mean x= 55.98007 , mean y= 46.68799 , sd   x= 27.50489 , sd   y= 28.02694 
mean x= 54.07907 , mean y= 44.82023 , sd   x= 28.01408 , sd   y= 27.70173 
mean x= 56.13726 , mean y= 46.82646 , sd   x= 27.60079 , sd   y= 28.17726 
mean x= 54.12897 , mean y= 44.90945 , sd   x= 27.91565 , sd   y= 27.97117 
mean x= 56.29275 , mean y= 44.95087 , sd   x= 27.61123 , sd   y= 27.96433 
mean x= 56.34399 , mean y= 45.06206 , sd   x= 27.47691 , sd   y= 27.94459 
mean x= 56.18838 , mean y= 45.19531 , sd   x= 27.47931 , sd   y= 28.01101 
mean x= 56.22124 , mean y= 45.17222 , sd   x= 27.52354 , sd   y= 28.02606 
mean x= 56.20132 , mean y= 45.14497 , sd   x= 27.5377 , sd   y= 28.09159 
mean x= 56.15781 , mean y= 47.10726 , sd   x= 27.49223 , sd   y= 28.39125 
mean x= 56.2795 , mean y= 47.05482 , sd   x= 27.45716 , sd   y= 28.47305 
mean x= 56.3636 , mean y= 47.28385 , sd   x= 27.52428 , sd   y= 28.32209 
mean x= 56.32062 , mean y= 47.40367 , sd   x= 27.57843 , sd   y= 28.28265 
mean x= 56.23393 , mean y= 47.43908 , sd   x= 27.60062 , sd   y= 28.29472 
mean x= 56.10059 , mean y= 47.53843 , sd   x= 27.73993 , sd   y= 28.29305 
mean x= 56.12577 , mean y= 47.52626 , sd   x= 27.84449 , sd   y= 28.30114 
mean x= 56.15453 , mean y= 47.68392 , sd   x= 27.9484 , sd   y= 28.29706 
mean x= 56.12583 , mean y= 45.69746 , sd   x= 27.92009 , sd   y= 28.0326 
mean x= 56.3937 , mean y= 45.83851 , sd   x= 27.85126 , sd   y= 28.00869 
mean x= 56.51539 , mean y= 45.97885 , sd   x= 27.68678 , sd   y= 27.96827 
mean x= 56.37371 , mean y= 46.06069 , sd   x= 27.64242 , sd   y= 27.87204 
mean x= 56.36624 , mean y= 48.14869 , sd   x= 27.70373 , sd   y= 28.10673 
mean x= 56.29951 , mean y= 48.06764 , sd   x= 27.66208 , sd   y= 28.12535 
mean x= 56.44135 , mean y= 48.00381 , sd   x= 27.63665 , sd   y= 28.30364 
mean x= 56.2341 , mean y= 47.89777 , sd   x= 27.63581 , sd   y= 28.23241 
mean x= 56.39581 , mean y= 47.86522 , sd   x= 27.60052 , sd   y= 28.33579 
mean x= 56.2757 , mean y= 47.90095 , sd   x= 27.65595 , sd   y= 28.22809 
mean x= 56.16795 , mean y= 45.97585 , sd   x= 27.54452 , sd   y= 28.09133 
mean x= 56.03491 , mean y= 45.95155 , sd   x= 27.57386 , sd   y= 27.95143 
mean x= 56.00317 , mean y= 47.96095 , sd   x= 27.71953 , sd   y= 28.34546 
mean x= 54.02754 , mean y= 46.06549 , sd   x= 28.00172 , sd   y= 28.10131 
mean x= 54.20156 , mean y= 47.99033 , sd   x= 28.04914 , sd   y= 28.42581 
mean x= 54.21784 , mean y= 50.02582 , sd   x= 28.02861 , sd   y= 28.55534 
mean x= 54.22306 , mean y= 46.0059 , sd   x= 28.15474 , sd   y= 28.28787 
mean x= 56.25764 , mean y= 47.95399 , sd   x= 27.89198 , sd   y= 28.47297 
mean x= 54.15198 , mean y= 47.9065 , sd   x= 28.37542 , sd   y= 28.46381 
mean x= 56.13884 , mean y= 49.86281 , sd   x= 27.89141 , sd   y= 28.46845 
mean x= 54.16265 , mean y= 47.84923 , sd   x= 28.24532 , sd   y= 28.51567 
mean x= 56.13058 , mean y= 47.98821 , sd   x= 27.95686 , sd   y= 28.51728 
mean x= 53.98531 , mean y= 47.94969 , sd   x= 28.18128 , sd   y= 28.41966 
mean x= 53.92721 , mean y= 47.77848 , sd   x= 28.06879 , sd   y= 28.45647 
mean x= 54.02449 , mean y= 47.8224 , sd   x= 28.18832 , sd   y= 28.49466 
mean x= 56.01244 , mean y= 47.75944 , sd   x= 27.79847 , sd   y= 28.65524 
mean x= 54.02643 , mean y= 47.77777 , sd   x= 28.03803 , sd   y= 28.84063 
mean x= 56.0832 , mean y= 45.81102 , sd   x= 27.89937 , sd   y= 28.63438 
mean x= 55.88557 , mean y= 47.72464 , sd   x= 27.97137 , sd   y= 28.93595 
mean x= 55.8515 , mean y= 47.6692 , sd   x= 27.8537 , sd   y= 28.84562 
mean x= 53.66417 , mean y= 47.67527 , sd   x= 28.11932 , sd   y= 28.84397 
mean x= 53.50896 , mean y= 47.62059 , sd   x= 28.00982 , sd   y= 28.91233 
mean x= 53.63125 , mean y= 47.61394 , sd   x= 27.90487 , sd   y= 28.78496 
mean x= 53.79239 , mean y= 47.58931 , sd   x= 27.92005 , sd   y= 28.77003 
mean x= 53.70271 , mean y= 47.69664 , sd   x= 28.15303 , sd   y= 28.74493 
mean x= 53.73276 , mean y= 47.66647 , sd   x= 28.24556 , sd   y= 28.7601 
mean x= 53.8027 , mean y= 47.65284 , sd   x= 28.29257 , sd   y= 28.72214 
mean x= 53.79434 , mean y= 49.77691 , sd   x= 28.3658 , sd   y= 28.85537 
mean x= 51.88273 , mean y= 47.97081 , sd   x= 28.51745 , sd   y= 28.72465 
mean x= 53.83487 , mean y= 48.08682 , sd   x= 28.26714 , sd   y= 28.65357 
mean x= 53.83992 , mean y= 46.19522 , sd   x= 28.21195 , sd   y= 28.44149 
mean x= 53.61338 , mean y= 46.14135 , sd   x= 28.30234 , sd   y= 28.32688 
mean x= 53.65303 , mean y= 46.18066 , sd   x= 28.40999 , sd   y= 28.36344 
mean x= 53.62193 , mean y= 46.34084 , sd   x= 28.31938 , sd   y= 28.19638 

simecol documentation built on Oct. 7, 2021, 9:20 a.m.