[gmx-users] questions about some mdp options in 4.5.1

[Michael Brunsteiner]

Hi,

just tried to start my first MD with gmx 4.5, and have
a couple of questions regarding some mdp options:

for the impatient ones, the short versions:

1) what is the parameter "nstcalcenergy" ?
    (doesn't seem to be documented)
    also, in the mailing list there's a discussion about a bug with
    pressure coupling with nstcalcenergy>1 and/or nstlist EQ -1.
    Is this still an issue in the current production version of 4.5.1?
    If i use nstcalcenergy=1 will that slow down my simulation a lot?

2) with amber force fields - is every atom a separate charge group?

3) how large is the error i get from using PME-switch compared to PME?

below more details:

======================================================

1) with "nstenergy = 100" i get:

WARNING 1 [file md1.mdp]:
   nstenergy should be a mltiple of nstcalcenergy,
   changing nstenergy to 108

the parameter nstcalcenergy does not seem to be documented
in: http://manual.gromacs.org/current/online/mdp_opt.html
what does it mean? what is a good value for it?


======================================================

2) NOTE 4 [file md1.mdp]:
   The sum of the two largest charge group radii (0.078268) is
   larger than rlist (1.200000) - rvdw (1.200000)

that is gmx assumes that the "sum of the two largest charge group radii"
is only 0.078268 nm, this is surprising as that's even smaller than the
sum of two water radii, and i would assume a single water being a
charge group ... my system is a protein plus ligand, water and
ions - does that mean that each atom is considered a separate charge
group here? (i use the AMBER99SB-ILDN FF)

======================================================

3) cut-offs:

using PME and rlist= 1.2 and rcoulomb=rvdw=1.2, nstlist=12 i get:
   NOTE 2 [file md1.mdp]:
   For energy conservation with switch/shift potentials,
   rlist should be 0.1 to 0.3 nm larger than rvdw.

using PME and rlist= 1.2 and rcoulomb=rvdw=1.0, nstlist=12 i get:
   ERROR 1 [file md1.mdp]:
   With coulombtype = PME, rcoulomb must be equal to rlist
   If you want optimal energy conservation or exact integration
   use PME-Switch

using PME and rlist= 1.2 and rcoulomb=rvdw=1.2, nstlist=-1  i get:
   ERROR 1 [file md1.mdp]:
   nstlist=-1 only works with switched or shifted potentials,
   suggestion: use vdw-type=Shift and coulomb-type=PME-Switch

that is the best choice (at least the one that gives me the fewest
warnings) seems to require PME-Switch ... but in mdp_opt.html it
says:
"For constant temperature simulations the advantage of improved energy 
conservation is usually outweighed by the small loss in accuracy of the 
electrostatics."

I wonder if anybody ever tested HOW SMALL this loss of accuracy really
is ... messing with the original ewald algorithm appears to be
dangerous to me, as the artifacts that come from improper treatment
of electrostatics are often fairly subtle - very small, but
systematic errors in energies/forces can add up to artifacts in
the long run - any opinions there?

======================================================

thanks!
Michael





==========================================
Michael Brunsteiner
Research Center Pharmaceutical Engineering
Graz, Austria
http://www.rcpe.at
http://www.brunsteiner.net