Derive diffusion coefficient for heat equation from random walk simulation

Question

I want to simulate the underlying stochastic process of diffusion on a microscopic level and compare the result with the solution of the heat equation. However, I'm not able to match the solution of the heat equation with my computed quantities, as I'm not able to figure out the correct time scale, spatial scale and diffusion coefficient.

Consider the following algorithm in pseudo code.

Position = [0,0,0,0,...] // place particles initially at origin
                         // The format is [x1,y1,x2,y2,...]

for t=1:numberTimeSteps
   for i = 1:numberParticles
      alpha = uniform random number in [ 0, 2pi];
      steplength = normal random number with mean 0 and stdv S

      xi += sin(alpha)*steplength
      yi += cos(alpha)*steplength

   end
end

Now I have a vector of positions of all my particles. To compute the density, I count all particles that are inside a domain $[-L,L]$, which I divide into cells (i.e. I compute a histogram).

histo = zeros(Nhist,Nhist)
for i = 1 :numberParticles
    idx = SomeFunc(xi) // compute the index corresponding to position x
    idy = SomeFunc(yi) // same for y
    histo[idy][idx] +=1
end
histo = histo / numberParticles 

Now I have a 2D histogram which roughly looks like the following picture.

Due to symmetrie, I now only consider the cross-section through the image. I know (or I hope to know, actually I'm not sure if this is the correct formula), that the density for a dirac impulse can be computed in 2D as \begin{align} u(x,t) = \frac{1}{(4\pi D t)^{dim/2}} exp(-x^2/(4Dt)) \end{align} where $u(x,t)$ solves the 2D heat equation on an unbounded domain $u_t = Du_{xx}$.

Now, for a large number of particles, my statistical computation should approach the analytic solution. However, I don't know how to derive $D$ from my computation. Also, I'm not sure about $x$ and $t$, as you can see in the next picture (where I apparently take the wrong values).

So here is my question:

What analytic function with which parameters will match the computed densities from my random walk approach. How does my simulation translate into the physical framework of the heat equation.

Thank you very much!

Answer 1

There is a theoretical answer to this, which appears in the derivation of the heat equation from random walks or from Brownian motion. The main result is that if you have a Brownian motion $\sigma B^x_t$ which starts at $x$, then $u(t,x)=E[\sigma B^x_t]$ solves the heat equation $u_t = \frac{\sigma}{2} u_{xx}$. So your diffusion coefficient should be $\frac{\sigma}{2}$.

You would probably be interested in checking this theoretical answer. For this you would need to choose a measure of closeness and then solve an optimization problem. The classic way would be least squares. This would be a nonlinear least squares problem, since your function depends in a nonlinear fashion on $D$. It should not be hard to solve numerically, however, especially since $\frac{\sigma}{2}$ will turn out to be a very good guess.