More efficient algorithm for dispersal kernel multiplication

Good morning everyone,

Right now I am writing code that disperses particles around an X by Y matrix, using essentially a discrete dispersal kernel. Here, X and Y correspond to points in a spatial location. For example, during each time step, 60% of particles would stay in the same cell, and 10% of particles would move to each adjacent cell. Thus, the number of particles in cell (X,Y) at the next time step is equal to .6*(X,Y)+.1*(X+1,Y)+.1*(X-1,Y)+... Alternatively, I could make 36% stay put, 10% go to each cell that is off by 1, and 3% go to every cell that is off by 2. Also, it would be nice to have three different reactions to edges: 1) particles fall off and vanish, 2) particles reflect off (and thus, dispersal is lower in edges), or 3) they warp around, like on a torus.

I'm looking for a fast way to do this for an arbitrary dispersal pattern. I could use nested for loops to go through point-by-point, but that seems really inefficient. I have considered taking the X by Y matrix, converting it to an (X*Y)-length vector, multiplying that vector by a XY by XY dispersal matrix, and then converting it back. The construction of the matrix could be slow, because it only needs to be done once. The problem is that I don't know a good way to create a dispersal matrix like that which handles edges correctly (it is simple in 1-D, but in 2-D edges appear right in the middle, and jumping diagonally becomes more complicated).

Does anyone know a good algorithm for a matrix like this (or a link to where I could find more info)? Alternatively, is there a completely better way to do this?

Thanks,

Simon