Average velocity of overdamped particles in external field

Question

In short: how to obtain the average velocity from the Fokker-Planck equation in the overdamped regime? (i.e. when the probability density is $P(\mathbf{x},t)$ and not $P(\mathbf{x},\mathbf{v},t)$, otherwise we could just consider the first moment of the variable $\mathbf{v}$).

Background: the Langevin equation in the overdamped regime (ie. there is no $\ddot{\mathbf{x}}$) is

$$ \dot{\mathbf{x}}(t) = \mathbf{V}(\mathbf{x}(t)) + \boldsymbol{\eta}(t) \, $$

where $\mathbf{V}:\mathbb{R}^n\rightarrow \mathbb{R}^n$ is a smooth field and $\boldsymbol{\eta}$ is the usual white-noise term,

$$ \langle \eta_i(t) \eta_j(t') \rangle_{noise} = 2 D \delta_{ij} \delta(t-t') $$

The related Fokker Planck equation for the particle distribution $P(\mathbf{x},t)$ is the conservation equation for the total probability:

$$ \partial_t P(\mathbf{x},t) = - \nabla \cdot ( \mathbf{J}_a +\mathbf{J}_d) $$

where

$$ \mathbf{J}_a(\mathbf{x},t) = P(\mathbf{x},t) \mathbf{V}(\mathbf{x}) \\ \mathbf{J}_d(\mathbf{x},t) = - D \, \nabla P(\mathbf{x},t) $$

are the "advection" and "diffusion" contributions to the total probability current.

Question: considering the Langevin ODE for many particles or the Fokker-Planck PDF should be equivalent, at least in the limit of many particles (i.e. many realizations of the Langevin dynamics). How to get the average velocity of particles in the two descriptions (Langevin VS Fokker-Planck)?

Langevin: it seems natural to solve the ODE for $N$ different particles, with different initial conditions ${\mathbf{x}}_i(0)$ (say, uniformly distributed in the domain $\Omega$ at $t=0$) and different realizations of the noise $\boldsymbol{\eta}$. The particles cannot leave $\Omega$, so that $N$ is constant. Hence, the average velocity is

$$\langle \dot{\mathbf{x}}(t) \rangle_N = N^{-1} \sum_{i=1..N} \dot{\mathbf{x}}_i(t) $$

Fokker-Planck: at $t=0$ we could choose a certain $P(\mathbf{x},0)$, say uniform (because in the Langevin picture the initial positions of the particles were uniformly distributed), $P(\mathbf{x},0) = 1/|\Omega|$, where $|\Omega|$ is the measure of the domain $\Omega$. Solving the Fokker-Planck equation gives $P$ at later times, $P(\mathbf{x},t)$. Which is the average velocity of particles? For large $N$, do we have that

$$ \langle \dot{\mathbf{x}}(t) \rangle_N \approx \int_\Omega d^nx \, \mathbf{J}_a(\mathbf{x},t) = \int_\Omega d^nx \, P(\mathbf{x},t) \mathbf{V}(\mathbf{x}) $$

or do we have to consider the full probability current

$$ \langle \dot{\mathbf{x}}(t) \rangle_N \approx \int_\Omega d^nx \, (\mathbf{J}_a(\mathbf{x},t)+\mathbf{J}_d(\mathbf{x},t)) \, \, ? $$

Answer 1

Assume you are farmiliar with Ito calculus and stochastic differential equation, the standard method in mathematics to deal with dynamics with white-noise like Langevin.

1. Overdamped Langevin dynamics differs essentially from standard Lagevin dynamics. In particular, it could not be regarded as simply dropping the inertia term.

The difference can be seen as follows.

Consider the standard Langevin dynamics, which is usually put as follows in physics: $$ \mu\,\ddot{x}(t)=-\dot{x}(t)-\nabla\phi(x(t))+\sqrt{2D}\,\eta(t), $$ and shall be put as follows in mathematics: \begin{align} {\rm d}x_t&=v_t\,{\rm d}t,\\ \mu\,{\rm d}v_t&=-v_t\,{\rm d}t-\nabla\phi(x_t)\,{\rm d}t+\sqrt{2D}\,{\rm d}W_t, \end{align} where $\mu=m/\gamma$ is the reduced mass, $\phi=\Phi/\gamma$ is the scaled potential, $\eta(t)$ denotes the normalized white noise, $D$ is the diffusion constant, and $W_t$ is the Wiener process (i.e., standard Brownian motion).

In the overdamped regime, one assumes $\mu\to 0^+$, and would thus expect that the standard Langevin dynamics reduces to, in physics, $$ \dot{x}(t)=-\nabla\phi(x(t))+\sqrt{2D}\,\eta(t), $$ or equivalently in mathematics, to \begin{align} {\rm d}x_t&=v_t\,{\rm d}t,\\ v_t\,{\rm d}t&=-\nabla\phi(x_t)\,{\rm d}t+\sqrt{2D}\,{\rm d}W_t. \end{align} By the first sub-equation, the reduced equations also writes \begin{align} {\rm d}x_t&=v_t\,{\rm d}t,\\ {\rm d}x_t&=-\nabla\phi(x_t)\,{\rm d}t+\sqrt{2D}\,{\rm d}W_t. \end{align} The second sub-equation here depicts the overdamped Langevin dynamics, also known as the Brownian dynamics.

However, by putting "overdamped regime", one should focus on not only the second sub-equation, but also the first one. Unfortunately, these two sub-equations contradict each other.

• By ${\rm d}x_t=v_t\,{\rm d}t$, $x_t$ is a usual time-integral of another stochastic process $v_t$, for which it has zero quadratic variation. More precisely, ${\rm d}\left<x,x\right>_t=0$. This can be intuitively understood as that, although $x_t$ is now a stochastic process, it is still almost everywhere time-continuous and differentiable. As such, velocity is still well-defined as the time-derivative of the position.
• By ${\rm d}x_t=-\nabla\phi(x_t)\,{\rm d}t+\sqrt{2D}\,{\rm d}W_t$, $x_t$ is not only a usual time-integral of another stochastic process $\nabla\phi(x_t)$, but also a stochastic integral with respect to the Wiener process $W_t$, for which its quadratic variation is non-zero, i.e., ${\rm d}\left<x,x\right>_t=2D\,{\rm d}t$. Intuitively, this means that $x_t$ is now almost everywhere time-continuous, but nowhere time-differentiable, for which velocity is not even defined if it is still taken as some time-derivative of position.

Therefore, it is not self-consistent by taking the overdamped Langevin dynamics as the standard Langevin dynamics with $\mu\to 0^+$.

2. Drift velocity and/or total kinetic energy, it really depends on what one truly wants.

Focus on the Brownian dynamics $$ {\rm d}x_t=-\nabla\phi(x_t)\,{\rm d}t+\sqrt{2D}\,{\rm d}W_t. $$ This single equation still makes sense. The question is: What is the proper definition of velocity for $x_t$ that solves this equation?

Consider two specific cases.

• There is no potential, i.e., $\phi=0$. In this case, the above equation yields $$ {\rm d}x_t=\sqrt{2D}\,{\rm d}W_t\iff x_t=x_0+\sqrt{2D}\,W_t, $$ which is the standard model to describe the Brownian motion. For such motion, it makes no sense to talk about the velocity of a Brownian particle; Instead, one would like to discuss its average kinetic energy. Mathematically, this is characterized by the quadratic variation of $x_t$, i.e., ${\rm d}\left<x,x\right>_t=2D\,{\rm d}t$, assuming the unit mass.
• There is no diffusion, i.e., $D=0$. In this case, the above equation becomes $$ {\rm d}x_t=-\nabla\phi(x_t)\,{\rm d}t. $$ This is an ordinary differential equation, whose solution is time-continuous and differentiable, for which velocity is in the usual sense, i.e., ${\rm d}x_t=v_t\,{\rm d}t$, and can be obtained from $v_t=-\nabla\phi(x_t)$ upon solving $x_t$.

Back to the Brownian dynamics. When talking about its velocity, it really depends on what velocity one truly wants. If one wants the drift velocity only, then it should be $$ v_t=-\nabla\phi(x_t). $$ In this case, only $\mathbf{J}_a$ shall be included in using Fokker-Plank equation. By contrast, if one wants the total kinetic energy, then this must include both the part arising from the drift velocity, and the part contributed from the diffusion. In this case, both $\mathbf{J}_a$ and $\mathbf{J}_d$ shall be included.

Answer 2

The quantity $ \langle \dot{\mathbf{x}}(t) \rangle_N $ defined in the question is the time derivative of the center of mass of the $N$ non-interacting particles that constitute the ensemble (i.e. the $N$ probes that explore the computational domain $\Omega$ in the Langevin picture).

In the Fokker-Planck framework, the center of mass is the first moment of the position,

$$\langle \mathbf{x} \rangle = \int_\Omega P( \mathbf{x}, t ) \mathbf{x} d^3x$$

Take the time derivative,

$$ \frac{d}{dt}\langle \mathbf{x} \rangle = \int_\Omega \mathbf{x} \, \partial_t \, P( \mathbf{x}, t ) d^3x = - \int_\Omega \mathbf{x} \, \nabla \cdot [ \mathbf{J}_a( \mathbf{x}, t ) + \mathbf{J}_d( \mathbf{x}, t ) ] d^3x $$

At this point it is easy to obtain (repeated indexes are summed)

$$ \frac{d}{dt}\langle x_i \rangle = \int_\Omega J^i_a( \mathbf{x}, t ) d^3x + \int_{\partial\Omega} Q_{ij}( \mathbf{x}, t ) d\Sigma_j $$

where $ Q( \mathbf{x}, t )$ is a matrix

$$ Q_{ij} = -x_i V_j P +D x_i \partial_j P - D P \delta_{ij} $$

Therefore, the diffusion current $J^i_d$ contributes only if there are some non trivial boundary conditions (like reflecting boundary conditions on a side but open boundary conditions on some other side of $\partial\Omega$).

In periodic boundary conditions (or if $\Omega = \mathbb{R}^3$) you have no boundaries $\partial \Omega$ and only the integral of the advection current over $\Omega$ contributes to the average velocity $\langle \dot{\mathbf{x}}(t) \rangle$.

NOTE: the only problem of this is when $P$ relaxes to a steady state. If this happens then the velocity associated to the center of mass is zero (simply because $\partial_t P = 0$), however, particles may still flow (e.g. a steady state in periodic boundary conditions: the center of mass is fixed but the particles can flow). Therefore, at least when periodic boundary conditions are used, the best thing to do is simply to define the average velocity directly as $\langle \mathbf{V} \rangle $.