A stochastic quasi Newton method for molecular simulations
Chau, C.D.
Citation
Chau, C. D. (2010, November 3). A stochastic quasi Newton method for molecular simulations. Retrieved from https://hdl.handle.net/1887/16104
Version: Corrected Publisher’s Version
License: Licence agreement concerning inclusion of doctoral thesis in the Institutional Repository of the University of Leiden
Downloaded from: https://hdl.handle.net/1887/16104
Note: To cite this publication please use the final published version (if applicable).
APPENDIX A
Langevin equation with space dependent mobility and its discretised form
A.1. Derivation of the general Langevin equation and nu- merical implementation
According to Gardiner [9] the many variable version of the Fokker Planck equation, which describes the time evolution of the probability density function of a stochast x is given as
∂p(x, t|x0, t0)
∂t =−
Xn i=1
∂i p(x, t|x0, t0)ai(x)
+ 1 2
Xn i=1
Xn j=1
∂i∂j p(x, t|x0, t0)Di j(x), (A-1)
which is related to the stochastic differential equation
dx = a(x)dt + B(x)dW, (A-2)
by D(x) = B(x)B(x)T. The drift vector a(x) and noise matrix B(x) are obtained by requiring the stationary solution of the FPE ps(x) to be the Boltzmann distribution
141
142 Appendix A
ps(x) = N exp−βφ(x). Setting ∂p(x,t|x∂t0,t0) = 0 and substitute ps(x) for p(x, t|x0, t0) in the FPE gives
ps(x)ai(x) = 1 2
Xn j=1
∂j ps(x)Di j(x)
(A-3)
= 1 2
Xn j=1
Di j(x)∂jps(x) + ps(x)∂jDi j(x)
(A-4)
= 1 2
Xn j=1
−βDi j(x)∂jφ(x)ps(x) + ps(x)∂jDi j(x)
(A-5)
→ ai(x) = 1
2 Xn
j=1
−βDi j(x)∂jφ(x) + ∂jDi j(x)
. (A-6)
This leads to the following SDE dx = 1
2−βD(x)∇φ(x) + ∇D(x) dt + B(x)dW(t), (A-7) where B(x)BT(x) = D(x) and β−1 = kBT . After defining D(x) = 2kBT M(x) = 2kBT L(x)LT(x) one obtains
dx =−M(x)∇φ(x) + kBT∇
·
M(x) dt + p2kBT L(x)dW(t), (A-8) where the noise term satisfies the fluctuation dissipation theorem. Equation (A-8) is equivalent to the SDE proposed by Hütter and Öttinger [39]dx = [−M(x)∇Φ(x)] dt +1 2
hM(x + dx)M(x)−1+ Ii p
2kBT L(x)dW(t). (A-9) This can easily proven by expanding M(x + dx) around x and obeying the rules dWdt = 0 and dWdW = dt.
The discretized form of the SDE proposed by Hütter and Öttinger is given below.
The update for xk at simulation step k is given as xk+1 = xk− 1
2
hM(xk+ ∆xkp)∇Φ(xk+ ∆xkp) + M(xk)∇Φ(xk)i
∆t +1
2
hM(xk+ ∆xkp)M−1(xk) + Ii p
2kBT L(xk)∆Wt, (A-10)
S-QN method for molecular simulations 143
with the corresponding predictor step
∆xpk =−M(xk)∇Φ(xk)∆t + p
2kBT L(xk)∆Wt. (A-11)
The approximate inverse Hessian B(xk) = Bk, which is taken as the mobility tensor M(xk), can be obtained using the DFP update
Bk+1= Bk− BkykyTkBk
yTkBkyk
+ sksTk yTksk
, (A-12)
where,
sk = xk+1− xk and yk =∇Φ(xk+1)− ∇Φ(xk). (A-13)