Brownian Motion Simulation

Equation of Motion

$$ \dot{\vec r}=\frac{\sum \vec F}{b} \tag{1} \label{1} $$

where $\sum \vec F = \vec F_\mathrm{Brownian} + \vec F_\mathrm{drag}$ is random and time-independent here; $b$ is the drag coefficient. For starters, set the drag coefficient to $1$ and make $\vec F$ a normally distributed random variable in 2D.

The Problem

To computationally compute the advanced velocities at each point use:

$$ \left(\dot{\vec r}_x\right)_n = \frac{\vec F _{\mathrm{rand}} + \vec F _{\mathrm{external}}}{b} + \left(\dot{r} _x\right) _{n-1} \tag{2}\label{2} $$

where Hooke's law is used for $\vec F _{\mathrm{external}}$. The random forces $\vec F _{\mathrm{rand}}$ can be scaled based on the initial velocity such that

$$ \frac{\vec F _{\mathrm{rand}}}{b} \quad\longrightarrow\quad \left( \frac{\vec F _{\mathrm{rand}}}{b} \right) _{\mathrm{scaled}} = \frac{ \vec F _{\mathrm{rand}}}{b}\left( {1 + \frac{\vec F _{\mathrm{rand}}/b}{\vec F _{\mathrm{rand}}/b +\left( {\dot{r}_x}\right) _{\mathrm{ave}}} } \right) \tag{3}\label{3}% $$

where $\left({\dot{r}_x}\right) _{\mathrm{ave}}$ is the mean of all the generated forces over the drag coefficient. Please note that there is no necessity for the scale to be applied unless the scheme is to minimize/eliminate effects of initial velocity.

Generating the Random Forces

To generate the random forces one of two methods can be used;

• Utter Randomness

for iteration in range(N):
    x_forces.append(random.random() * random.choice(parity))
    y_forces.append(random.random() * random.choice(parity))

• Normalized Distribution of Contaminated Randomness

for i in range(N):
    x_forces.append( Fx[i] * random.choice(parity) )
    y_forces.append( Fy[i] * random.choice(parity) )

where $\texttt{parity}=\pm1$.

Mean Square Distance

One method of calculating the mean square distance (MSD) is by using Equation $\eqref{4}$;

$$ \begin{equation} \left\langle \Delta r^2(t)\right\rangle = \dfrac{\sum_{i=1}^{N_i} \Big[ \big(x(t_i-t)-x(t)\big)^2 + \big(y(t_i-t)-y(t)\big)^2 \Big]}{N_i} \tag{4} \label{4} \end{equation} $$

For instance, the first element of the MSD of list = np.array([[0, 1], [2, 3], [4, 5]]) would be:

$$ \begin{align*} \texttt{MSD}[0] &= \frac{1}{3-1}\mathrm{SUM}\bigg(\Big[\big[{(5,4) - (3,2)}\big]^2,\ \big[(3,2) - (1,0)\big]^2\Big]\bigg)\\\\ &= \frac{1}{2}\mathrm{SUM}\Big(\big[{(2,2)^2},\ (2,2)^2\big]\Big) \\\\ &= \frac{1}{2}\mathrm{SUM}\Big(\big[ (4,4),\ (4,,4)\big]\Big) \\\\ &= \frac{1}{2}(4+4+4+4) \\\\ &= 8 \end{align*} $$

while its second and last element would be:

$$ \begin{align*} \texttt{MSD}[1] &= \frac{1}{3-2}\mathrm{SUM}\bigg(\Big[\big[{(4,5) - (0,1)}\big]^2 \Big]\bigg)\\\\ &= \mathrm{SUM}\Big(\big[{(4,4)^2}\big]\Big) \\\\ &= 16 + 16 \\\\ &= 32 \end{align*} $$

The MSD should be linear with time; i.e. plotting $\log{(\mathrm{MSD})}$ vs $\log{(\Delta t)}$ should give a line with unit slope.

What are the external forces?

$$ \begin{align} &\vec F _{\mathrm{FENE}} = \frac{3k _{\mathrm{B}}T}{l _{\mathrm{k}}} \left( \frac{\vec r /L _{\mathrm{s}}}{1-{\left(\vec r/L _{\mathrm{s}}\right)^2}} \right) \tag{5}\label{5} \\\\ &\vec F _{\mathrm{excluded~volume}} = k _{\mathrm{excluded \ volume}} \sum _j \left( x - x _j \right) \tag{6}\label{6} \end{align} $$