[gmx-users] Advice needed on choosing a VdW cutoff
Steven.Kirk at hv.se
Tue Jul 11 17:27:46 CEST 2006
I would be very grateful for advice on the following system:
Consider a pair of spherical macromolecules of diameter ~2.5 nm,
arranged along an x-axis so that their centres of mass are 4nm apart,
then centered in a 14 x 10 x 10 nm periodic box, so that the minimum
distance between either macromolecule's centre of mass and its closest
simulation cell walls in the x,y, or z directions (standard orthogonal
Cartesian axis set) is 5 nm.
(the macromolecules are silica particles with negatively charged surface
sites and an equal number of Na+ counterions clustered around them, so
each macromolecule-counterion 'bundle' is electrically neutral).
The box is then filled with TIP4P water, energy minimized, run with
'all-bonds' position restraints at 300K for 100 ps, then the main MD run
begins, reusing velocities from the last step of the position restraints
to initialize the new run.
Here is the mdp file I am using for the main MD run:
title = Yo
cpp = /lib/cpp
constraints = none
integrator = md
dt = 0.001 ; ps !
nsteps = 1000000 ; total 1000 ps.
nstcomm = 1
nstxout = 10000
nstvout = 10000
nstfout = 0
nstlog = 10000
nstenergy = 10000
nstlist = 10
ns_type = grid
coulombtype = pme
rlist = 1.0
rcoulomb = 1.0
rvdw = 1.0
; Berendsen temperature coupling is on in two groups
Tcoupl = berendsen
tc-grps = SNP SOL NA+
tau_t = 0.1 0.1 0.1
ref_t = 300 300 300
; Energy monitoring
energygrps = SNP SOL NA+
; Isotropic pressure coupling is not on
Pcoupl = no
Pcoupltype = isotropic
tau_p = 0.5
compressibility = 4.5e-5
ref_p = 1.0
; Generate velocites is off at 300 K.
gen_vel = no
gen_temp = 300.0
gen_seed = 173529
The problem is that when I run this simulation, the expected drift of
the macromolecules towards each other does not occur. Assuming that I
want to force every part of each macromolecule to 'see' every part of
the other, this would suggest a value of rvdw of > 6.5 nm, but I have
several worries about this:
1. I have never seen a recommended rvdw in this forum over 1.4 nm, in
any model system
2. Should I use a standard, switched or other type of vdW cutoff?
3. Should I switch on long-range dispersion corrections (DispCorr = Ener) ?
My goal is to keep the periodic box for the advantages of PME, but
somehow reassure myself that the macromolecules can 'see' each other via
the vdW forces, so that they will drift together (the expected
behaviour) over the course of the simulation.
I have trawled the mailing lists for advice on this topic - the only
directly relevant post I could find involved the drifting apart of
membranes over the course of a simulation, and if I remember correctly,
the value of 'pme-order' was suggested as a culprit.
I use the default value of 'pme-order' in my simulations.
Can anyone please advise me as to what to do next? I do not want to
abandon the investigation before I eliminate all possibilities that I am
doing something stupid with the configuration of the vdW treatment.
Maybe I do not have a long enough simulation, but the fact that I
actually see *repulsion* suggests something more fundamental is wrong.
All helpful suggestions very gratefully recieved !
More information about the gromacs.org_gmx-users