Consider a set of N balls that start at the origin. In a given unit of time, $\delta t$, the balls have a probability $p = 0.5$ of jumping a distance $\delta x$ to the right, and the same probability of jumping $\delta x$ to the left.

It can be shown (see Ockendon et al) that the concentration, $c(x, t)$, that is, the number of balls per unit length, is given by solving the heat equation,

\begin{equation} \frac{\partial u}{\partial t} = D \frac{\partial^2 u}{\partial x^2} \end{equation}

where the diffusion coefficient, $D = \left(\frac{p \delta x^2}{\delta t}\right)$.

I wrote a very simple simulation in Matlab that does the stochastic procedure for $N = 1000$, $\delta t = 0.1$, $\delta x = 0.005$, and runs until reaching $t = 600$. After the simulation is done, the domain is split into bins centered at $0.015k$, $k = 0, \pm 1, \ldots$ and the balls tallied.

On top of this, the fundamental (Direc delta) solution of the heat equation is plotted. This has solution \begin{equation} u = \frac{A}{\sqrt{4 \pi D t}} \exp\left( - \frac{x^2}{4 D t}\right). \end{equation}

The problem is that I'm not sure what to use for $A$. I assumed that if $u$ is concentration per unit length, then I should have \begin{equation} \text{Number of balls in bin} = u(x, t) \cdot N \cdot \text{bin size}. \end{equation}

However, this seems to underpredict the results at the center. Does anybody know where I went wrong, or whether this is expected? I could not get better results by varying $N$. An image can be seen below:

Image link here

Note: This was resolved in my answer below. My diffusion coefficient was off by 2.


1 Answer 1


This was resolved. I had made a mistake in the diffusion constant (by a factor of 2). Once it's corrected, the analytical vs. stochastic result is bang on. See here.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.