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:

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