1549380323-Statistical Mechanics Theory and Molecular Simulation

Constraints: The ROLL algorithm 259

(5.13.8) if the dimensionality is not too large, or as a time-saving measure, we neglect
the dependence of eqn. (5.13.8) onl 6 =kterms, write the condition as

σl(r

(1)

1 ,...,r

(1)

N) +

∑N

i=1

F

(l)

c,i(1)·

1

mi

RFx(λ,0)δ ̃λ

(1)

l F

(l)

c,i(0)≈^0 , (5.13.9)

and iterate the correctionsδ ̃λ
(1)
l to convergence as in Section 3.9. Eqn. (5.13.9) can be
solved easily for the multiplier correctionsδ ̃λ(1)l to yield

δλ ̃(1)l =−

σl(r(1) 1 ,...,r(1)N)

∑N

i=1(1/mi)RFx(λ,0)F

(l)

c,i(1)·F

(l)

c,i(0)

. (5.13.10)

Whichever procedure is used to obtain the correctionsδ ̃λ(1)l , once we have them, we
substitute them into eqn. (5.13.1) to obtain a new update to the pressure virial. Using
this new pressure virial, we now cycle again through the operators

exp(iLǫ, 1 ∆t) exp(iL 1 ∆t) exp(iL 2 ∆t/2) exp(iLǫ, 2 ∆t/2),

which are applied on the original coordinatesri(0) andv(NHC)i. This will generate a

new set of ROLL scalars, which we use to generate a new set of correctionsδλ ̃(1)l using
the above procedure. This cycle is now iterated, each producing successively smaller

correctionsδλ ̃
(n)
l to the multipliers, until the ROLL scalars stop changing. Once this
happens, the constraints will be satisfied and the pressure virial will be fully converged.
Using the final multipliers, the half-step velocities are obtained fromeqn. (5.13.5). It
is important to note that, unlike the algorithm proposed by Martynaet al.(1996),
this version of the first half of the ROLL algorithm requires no iteration through the
thermostat operators.
The second half of the ROLL algorithm requires an iteration throughthe opera-
tors exp(iLǫ, 2 ∆t/2)exp(iL 2 ∆t/2). However, it is also necessary to apply the operators
exp(iLNHC−part∆t/2) and exp(iLNHC−baro∆t/2) in order to obtain the overall scaling
factors on the velocitiesvi(∆t) andvǫ(∆t), which we will denoteSi(∆t) andSǫ(∆t).
Thus, the entire operator whose application must be iterated is

Oˆ= exp(iLNHC−baro∆t/2) exp(iLNHC−part∆t/2)

×exp(iLǫ, 2 ∆t/2) exp(iL 2 ∆t/2). (5.13.11)

The evolution ofvican now be expressed as