6.6 Stochastic Dynamics

In many applications we encounter widely varying time scales. In such cases the ``fast''_ degrees of freedom dominate the choice of the time step, although they may be of lesser interest.
Example: A few slow-moving heavy ions in a bath of many light water molecules.

Strategy: Mimick the effect of the secondary particles by suitably sampled stochastic forces $\bf S$.

LANGEVIN'S equation of motion for a single ion in a viscous solvent:

$$\dot{\mbox{$\bf v$}}(t)=-\eta \mbox{$\bf v$}(t)+\mbox{$\bf a$}(t)$$ (6.19)

where the statistical properties of the stochastic acceleration $\mbox{$\bf a$}\equiv \mbox{$\bf S$}/m$ are

$$\langle \mbox{$\bf v$}(0) \cdot \mbox{$\bf a$}(t) \rangle = 0 \hspace{60pt} {\rm for} \;\;\; t\geq 0$$ (6.20)

$$\langle \mbox{$\bf a$}(0) \cdot \mbox{$\bf a$}(t) \rangle = 3 \frac{2\eta kT}{m} \, \delta (t)$$ (6.21)

Explanation:

• $\mbox{$\bf a$}(t)$ is not correlated to previous values of the ion velocity
• Stochastic and frictional forces are mutually related (both are caused by collisions of the ion with solvent molecules)
• Since equation 6.21 gives us only the a.c.f. of $\mbox{$\bf a$}(t)$, we have yet to specify its statistical distribution; the usual choice is a Gauss distribution for the components of $\mbox{$\bf a$}(t)$

$\Longrightarrow$Solution of 6.19:

$$\mbox{$\bf v$}(t)=\mbox{$\bf v$}(0)e^{-\eta t}+\int \limits_{0}^{t} e^{\textstyle-\eta(t-t')} \mbox{$\bf a$}(t')\, dt'$$

and similar for $\bf r$(t). Subtracting $\mbox{$\bf v$}(t_{n})$ from $\mbox{$\bf v$}(t_{n})$ etc., we have

$$\mbox{$\bf v$}_{n+1} = \mbox{$\bf v$}_{n}e^{\textstyle -\eta \Delta t} +\int \limits_{0}^{\Delta t} e^{\textstyle -\eta(\Delta t-t')} \mbox{$\bf a$}(t_{n}+t')\, dt'$$ (6.22)

$$\mbox{$\bf r$}_{n+1} = \mbox{$\bf r$}_{n}+\mbox{$\bf v$}_{n} \frac{1-e^{\textstyle -\eta... ...rac{1-e^{\textstyle -\eta(\Delta t-t')}}{\eta}\, \mbox{$\bf a$}(t_{n}+t')\, dt'$$ (6.23)

Defining

$$e(t)\equiv e^{\textstyle -\eta t}\,,\;\;\;\;\; f(t)\equiv \frac{1-e^{\textstyle -\eta t}}{\eta}$$ (6.24)

and

$$\mbox{$\bf V$}_{n}\equiv\int \limits_{0}^{\Delta t} e(\Delta t-t') \,\mbox{$\bf a$}(t_{n}+t')$$ (6.25)

$$\mbox{$\bf R$}_{n}\equiv\int \limits_{0}^{\Delta t} f(\Delta t-t') \,\mbox{$\bf a$}(t_{n}+t')$$ (6.26)

we may write the stepwise solution

$$\mbox{$\bf v$}_{n+1} = \mbox{$\bf v$}_{n}\, e(\Delta t) + \mbox{$\bf V$}_{n}$$ (6.27)

$$\mbox{$\bf r$}_{n+1} = \mbox{$\bf r$}_{n}+\mbox{$\bf v$}_{n}\, f(\Delta t) + \mbox{$\bf R$}_{n}$$ (6.28)

The components of the stochastic vectors $\mbox{$\bf V$}_{n}, \mbox{$\bf R$}_{n}$ are time integrals of the function $\mbox{$\bf a$}(t)$ whose statistical properties are given. $\Longrightarrow$ $\mbox{$\bf V$}_{n}, \mbox{$\bf R$}_{n}$ are themselves random variates with known statistics: $\langle V_{n}\rangle=\langle R_{n}\rangle=0$, $\langle V_{n}V_{n+1}\rangle=\langle R_{n} R_{n+1}\rangle=0$, and

$$\langle V_{n}^{2}\rangle = \frac{kT}{m}\left[ 1-e^{2}(\Delta t)\right]$$ (6.29)

$$\langle R_{n}^{2}\rangle = \frac{kT}{m \eta^{2}} \left[ 2 \eta \Delta t - 3 + 4e(\Delta t) - e^{2}(\Delta t)\right]$$ (6.30)

$$\langle V_{n}R_{n}\rangle = \frac{kT\eta}{m} \, f^{2}(\Delta t)$$ (6.31)

In the chapter about stochastics we described a method to produce pairs of correlated Gaussian variates. We may apply this here to generate $V_{n}, R_{n}$ and insert these in 6.27-6.28.

Generalization:
The stochastic force need not be $\delta$-correlated. If the solvent particles have a mass that is comparable with that of the solute, they will also move with similar speeds. In such cases the generalized Langevin equation applies:

$$\dot{v}(t) = - \int_{0}^{t} M(t-t')\, v(t')\, dt' + a(t)$$ (6.32)

where

$$\langle v(0) a(t) \rangle = 0 \hspace{60pt} \mbox{\rm for} \;\;\; t\geq 0$$ (6.33)

$$\langle a(0)a(t) \rangle = \frac{kT}{m} \, M(t)$$ (6.34)

This is a stochastic integrodifferential equation involving the ``history'' of the solute particle's motion in the form of the memory function $M(t)$ (see [MORI 65]). Still, we may assume that $M(t)$ decays fast.

Two approaches:

• Approximate the memory function by a suitable class of functions: assuming that the Laplace transform $\widehat{M}(s)$ may be represented by a truncated chain fraction in the variable $s$, the integrodifferential equation may be replaced by a set of coupled differential equations. Written in matrix notation these equations have the same form as 6.19 and may be treated accordingly.[VESELY 84]
• Assume that $M(t)$ may be neglected after $K \approx 20-60$ time steps. Using a tabulated autocorrelation function $M(t)$ one may generate an autoregressive process by the method described in the chapter on stochastics. By replacing the integral in 6.32 by a sum over the most recent $K$ time steps, one arrives at a stepwise procedure to produce $\mbox{$\bf v$}(t)$ and $\mbox{$\bf r$}(t)$; see [SMITH 90], and also [NILSSON 90]).