Page 354 - Numerical Methods for Chemical Engineering

Brownian dynamics and SDEs 343

dW t is a differential update of a Wiener process. We are now faced with choosing the time
in [t, t + dt] at which we evaluate dU/dx; here, we do so at the beginning of the time step,
to obtain an SDE that is an Itˆo-type SDE,
) *
1 dU 1/2
dx =− dt + [2D] dW t (7.152)
ζ dx
x(t)
Equation (7.152) is also called the Langevin equation for the particle.
For a description of SDEs, see Kloeden & Platen (2000). An abbreviated discussion is
¨
found in Ottinger (1996), with a focus on polymer science.
The simplest rule to integrate (7.152) is the explicit Euler SDE method:
) *

1 dU 1/2
x(t + t) − x(t) =− ( t) + [2D] ( W t ) (7.153)
ζ dx
x(t)
W t is an approximation to dW t for a finite t and is drawn at random from a Gaussian
2
distribution with mean µ = 0 and variance σ = t.In MATLAB, we generate W t by the
code,
dW t = sqrt(dt) * randn;
BD 1D.m uses the explicit Euler method to simulate the Brownian dynamics of a spherical
2
particle in a quadratic energy well, U(x) = k sp x /2 with k sp = 1. Trajectories are plotted
in Figure 7.11 for ζ = 1, D = 1 and various t.As t decreases, the path fluctuates more
wildly for short times, but has the same long-time properties.

It ˆ o’s stochastic calculus

While the explicit Euler method is simple, it is not very accurate. For a deterministic
differential equation, we build higher-order methods through Taylor series expansions;
however, the rules of stochastic calculus are different. Consider the SDE

(7.154)
dX = a(t, X)dt + b(t, X)dW t
= X(t k ), we have an exact update
Using t k+1 − t k = t, X t k
t k+1 t k+1
' '
= a(t, X(t))dt + (7.155)
X t k+1 − X t k b(t, X(t))dW t
t k t k
In deterministic calculus (no dW t term), we expand a(t, X) in time using the differential da =
(∂a/∂t)dt + (∂a/∂ X)(dX/dt)dt; however, here we have a stochastic integral involving a
Wiener process. How do we interpret this integral, and what is the proper form of differential
for a stochastic process?
We use here Itˆo’s stochastic calculus, in which we approximate the stochastic integral
by quadrature using the values at the beginning of each subinterval:
N
t F
'
jt F

b(t)dW t ≈ b(t j−1 ) W t j − W t j−1 t j = (7.156)
0 j=1 N