Implementing Laplacian operator acting on a 2D array

You can always just write a loop in Julia. Loops are fast. For example, this post give example code for \nabla^2 in N dimensions with a loop: Seemingly unnecessary allocation within for loop - #8 by stevengj … writing code that only works in 2d is even easier.

For periodic boundary conditions, the naive approach is to compute i\pm 1 modulo the size in the loop body, but this will slow things down a lot. A very flexible approach is to use “ghost cells” (padding) as in the code linked above, which lets the loop body not care about the boundary conditions, which are updated in a separate (cheap) loop over the boundaries.