520 Classical time-dependent statistical mechanics
means that the coefficient of shear viscosityηmust be independent of the choice of
shear rateγ. Since the extrapolation procedure requires a series of simulationsat dif-
ferent shear rates, these simulations should be used to check that the average〈Pxy〉t
varies linearly withγ. If it does not, the shear rate is too high and must be decreased.
13.5.1 Reversible integration of the SLLOD equations of motion
The developing reversible numerical integration methods for the general driven system
in eqns. (13.2.1) depends on the particular form for the functionsCi(q,p) andDi(q,p)
appearing in them; each problem needs to be treated as a separatecase. While eqns.
(13.3.17) for the calculation of the diffusion constant are quite straightforward to
integrate, for example, the SLLOD equations for planar Couette flow are somewhat
more complicated. Hence, we will consider this latter case as an example of how to
construct an integration algorithm for a nontrivial nonequilibrium molecular dynamics
problem.
There are two issues that need to be considered in the SLLOD scheme. First,
the equations of motion are non-Hamiltonian and involve nontrivial driving terms.
Second, we must either impose Lees–Edwards boundary conditionsor employ the time-
dependent box matrix of eqn. (13.5.1). When performing nonequilibrium molecular
dynamics simulations, the use of thermostats is often imperative, as the coupling to
the external field causes the system to heat continuously, and a steady state can never
be reached without some kind of temperature control mechanism.Thus, if the system
is coupled to a Nos ́e–Hoover chain thermostat, the simulation wouldbe based on the
following equations of motion:
r ̇i=
pi
mi
+γyiˆex
p ̇i=Fi−γpyiˆex−
pη 1
Q 1
pi
η ̇j=
pηj
Qj
j= 1,...,M
p ̇η 1 =
[N
∑
i=1
p^2 i
mi
−dNkT
]
−
pη 2
Q 2
pη 1
p ̇ηj=
[
p^2 ηj− 1
Qj− 1
−kT
]
−
pηj+1
Qj+1
pηj j= 2,...,M− 1
p ̇ηM=
[
p^2 ηM− 1
QM− 1
−kT
]
. (13.5.2)
In order to devise a numerical integration scheme for eqns. (13.5.2), we begin by writing
the Liouville operator for the system in a form similar to eqn. (4.11.4):
iL=iLmNHC+iL 1 +iL 2 , (13.5.3)
where the three terms are defined to be