Subject: Time step for Brownian dynamics

From: [email-masked] (Wilfred Tang)
Date: Thu, 28 Oct 93 18:30:56 -0500

I'm writing a Brownian dynamics program using the following diffusive
algorithm of van Gunsteren and Berendsen:

                      (dt)              1  F(t)-F(t-dt)
  x(t+dt)= x(t) + ------------ [F(t) + --- ------------ (dt)] + random term
                   (m) (gamma)          2      (dt)

                                                ^
                                                |
                                                |
                                            estimates derivative
                                            of F with respect to t
    dt = time step
    m = mass
    gamma = collision frequency
    F = force
    x = position

I use GROMOS to calculate the force, but I wrote my own simulation algorithm
since GROMOS doesn't do Brownian dynamics.  I'm interested in studying
proteins.

I'm wondering what time step dt should be used.  Actually, more precisely, what
ratio (dt)/[(m)(gamma)] to use.

The largest time step I can get out of my program is
  dt = 0.001 ps, with
  m = 14 g/mole
  gamma = 150 ps^-1
giving
       (dt)
  ------------ = 5E-7 (ps ps)/(g/mole)
   (m)(gamma)
If I use any larger time step dt, the energy rapidly goes to infinity;
to get a good simulation, I probably have to use a smaller time step dt than
0.001 ps.

This seems to contradict what's been published in the literature.
In particular, looking at a Brownian dynamics simulation on a lipid chain in
a membrane bilayer (Pastor, Venable, and Karplus, J. Chem. Phys. 89, 1112
(1988)), the parameters used are:
  dt = 0.015 ps  (much larger than my time step!!!!!!!!!!!!)
  m = 14 g/mole  (I think)
  gamma = 150 ps^-1
giving
       (dt)
  --------------- = 7E-6 (ps ps)/(g/mole)
   (m)(gamma)

I don't constrain my bond lengths, while Pastor et al. use SHAKE to constrain
bond lengths.  The rule of thumb I've heard is that SHAKE allows a 3 times
larger time step, but that still doesn't explain the discrepancy between
the time step Pastor et al. used and the time step I can obtain.

I'd appreciate it if someone with experience in doing Brownian dynamics
simulations can give me some hints and/or explanations.  We don't have
a lot of experience doing Brownian dynamics simulations at U. of Chicago.

Thanks,
Wilfred Tang
([email-masked])