the forces due to the left-out degrees of freedom. This should not in¬‚uence

the equilibrium properties, but is likely to yield a faster dynamics than the

real system. The “bead-and-spring” models for polymers, which have a long

history in polymer physics, are in fact also superatom models, although they

were intended as prototype polymer models rather than as simpli¬ed repre-

sentation of a speci¬c real polymer.5 The interaction between neighboring

beads in a chain is often simply a soft harmonic potential that leads to a

Gaussian distribution for the distance between beads; since polymers chains

can only be extended to an upper limit rm , a somewhat more realistic model

is the FENE (¬nitely extendable nonlinear elastic) chain model, with a force

between neighboring beads with interbead vector r given by

r

F = ’H . (8.15)

1 ’ (r/rm )2

See Fan et al. (2003) for a stochastic application of the FENE model.

More recently superatom models have been designed to speci¬cally rep-

resent real molecules, e.g., alkanes by the coarse-grained model of Nielsen

et al. (2003), with Lennard“Jones (LJ) superatoms representing three non-

hydrogen atoms. In addition to LJ, this model has soft harmonic bond

length and bond angle potentials. The model is parameterized on density

and surface tension of liquid alkanes and reproduces end-to-end distribution

functions obtained from simulations with atomic details. A more general and

very successful coarse-grained force ¬eld for lipids and surfactant systems has

been de¬ned and tested by Marrink et al. (2004). It consists of four types of

5 See M¨ ller-Plathe (2002) for a review on multiscale modelling methods for polymers.

u

8.5 The ¬‚uctuation“dissipation theorem 257

particles (charged, polar, non-polar, apolar) with the charged and non-polar

types subdivided into four subtypes depending on their hydrogen-bonding

capability. Each particle represents about four non-hydrogen atoms; also

four water molecules are represented by one particle. The interactions are

of Lennard“Jones and Coulomb type, smoothly switched-o¬ at 1.2 nm, with

¬ve possible values for the µ and only one (0.47 nm) for the σ parameters of

the LJ interaction. In addition there are fairly soft harmonic bond and bond-

angle interaction terms; the total number of parameters is eight. Despite

its simplicity, the model reproduces density and isothermal compressibility

of water and liquid alkanes within 5% and reproduces mutual solubilities of

alkanes in water and water in alkanes to within a free energy of 0.5kB T . The

melting point of water is 290 K. A time step of 50 fs is possible, and since

the dynamics of this model (no friction and noise are added) is about four

times faster than reality, an e¬ective time step of 0.2 ps can be realized. One

then easily simulates real systems with millions of atoms over microseconds,

allowing the study of lipid bilayer formation, micellar formation, vesicle for-

mation and gel/liquid crystalline phase changes with realistic results.6 See

Fig. 8.1 for comparison of a coarse-grained and a detailed simulation of the

spontaneous formation of a small lipid vesicle.

8.5 The ¬‚uctuation“dissipation theorem

1 2

Let us look at the long-time behavior of the kinetic energy K = i 2 mi vi

in the generalized Langevin equation (8.1). The ¬rst term on the r.h.s. (the

systematic force) simply exchanges kinetic and potential energy, keeping the

total energy constant. The second term (friction or dissipation) reduces the

kinetic energy and the third stochastic term (noise) increases the kinetic

energy. In order for the process to be stationary, the velocity correlation

functions v i (t)v j (t + „ ) should become independent of t for large t; in par-

ticular the average squared velocities v i (t)v j (t) , which are thermodynamic

quantities, should ful¬ll the equipartition theorem

kB T

vi± (t)vjβ (t) = δij δ±β . (8.16)

mi

That is, if the random process is realized many times from the same starting

con¬guration at t = 0, then after a su¬ciently long time “ when the memory

of the initial conditions has decayed “ the average over all realizations should

6 Spontaneous aggregation of lipid bilayers, gel/liquid crystalline transitions, inverted hexagonal

phase formation and formation of Micelles: Marrink et al. (2004); hexagonal phase forma-

tion: Marrink and Mark (2004); vesicle fusion: Marrink and Mark (2003a); vesicle formation:

Marrink and Mark (2003b).

258 Stochastic dynamics: reducing degrees of freedom

Figure 8.1 Two simulations of the spontaneous formation of a lipid bilayer vesicle.

Upper panel: atomic detail molecular dynamics; lower panel: coarse-grained super-

atom dynamics simulation (Courtesy of S.-J. Marrink and A. H. de Vries, University

of Groningen; reproduced by permission from J. Comput. Chem. (van der Spoel

et al., 2005)

ful¬ll the equipartition theorem. This is one “ rather restricted “ formulation

of the ¬‚uctuation“dissipation theorem.

The general ¬‚uctuation“dissipation theorem relates the linear response of

some system variable v to the spontaneous ¬‚uctuation of v. Kubo (1966)

distinguishes a ¬rst and a second ¬‚uctuation“dissipation theorem: the ¬rst

theorem says that the normalized time response of v to a δ-disturbance

equals the normalized correlation function of the spontaneously ¬‚uctuating

v in equilibrium (see Section 18.3 on page 511); the second theorem relates

the friction kernel ζ(t) to the correlation function of the random term ·(t)

in the Langevin equation..

To illustrate these theorems we apply for the sake of simplicity the re-

sponse to a single velocity v(t) that is assumed to follow the pure Langevin

equation without systematic force:

t

mv(t) = ’ ζ(„ )v(t ’ „ ) d„ + ·(t) + F ext (t),

™ (8.17)

0

with F ext (t) an external force meant to provide a disturbance to measure

the linear response. Now apply an external δ-force at time 0:

F ext (t) = mv0 δ(t). (8.18)

This produces a jump v0 in the velocity at time t = 0. The velocity v(t)

subsequently evolves according to (8.17) (with the external force no longer

8.5 The ¬‚uctuation“dissipation theorem 259

being present). What we call the δ-response v0 ¦(t) of the velocity7 (see

(18.1) on page 507) is the ensemble average over many realizations of this

process, with initial conditions taken randomly from an unperturbed equi-

librium distribution and with independent realizations of the noise force

·(t):

v0 ¦(t) = v(t) . (8.19)

The ¬rst ¬‚uctuation“dissipation theorem states that

v(t0 )v(t0 + t)

¦(t) = , (8.20)

v2

where the average is now taken over an unperturbed equilibrium ensemble,

for which the velocity correlation function is stationary and hence independ-

ent of the time origin t0 .

Proof From (8.19) it follows that

d¦

v(t) = v0

™ , (8.21)

dt

so that on averaging (8.17) immediately gives an equation for ¦, considering

that the average over the noise is zero (see (8.2)):

t

d¦

=’ ζ(„ )¦(t ’ „ ) d„.

m (8.22)

dt 0

Given the friction kernel ζ(„ ) and the initial value ¦(0) = 1, this equation

determines the response function ¦(t).

def

The velocity autocorrelation function C v (t) = v(t0 )v(t0 + t) can be

found by applying (8.17) to the time t0 + t, multiplying both sides with

v(t0 ) and taking the ensemble average:

t

m v(t0 )v(t0 +t) = ’

™ ζ(„ ) v(t0 )v(t0 +t’„ ) d„ + v(t0 )·(t0 +t) , (8.23)

0

which can be written in terms of the velocity correlation function C v (t),

realizing that the last term vanishes because the random force does not

correlate with velocities at earlier times (see (8.3)), as