Hi, I am a Gromacs novice and am trying to simulate 1000 SPC/E water molecules (with the usual charges) in the NVE ensemble. I am trying to obtain good energy conservation.
My configuration contains a 10-by-10-by-10 cubic crystal of SPC/E waters, with a lattice constant of approximately 4.7 Angstroms. I first run energy minimization, using integrator = steep, emtol = 1000 kJ/mol/nm, emstep = 0.01 nm, and coulombtype = PME. This energy minimization converges in only 1 step -- I guess because my value of emtol is so large and because the water molecules are relatively far apart from each other. Then, when I run my initial equilibration/dynamics (I am calling my run both equilibration and dynamics because, for the time being, I am just doing a single run) using coulombtype = PME, I get the following two notes messages: ------ NOTE 1 [file nve-eq.mdp]: You are using a cut-off for VdW interactions with NVE, for good energy conservation use vdwtype = Shift (possibly with DispCorr) NOTE 2 [file nve-eq.mdp]: You are using a cut-off for electrostatics with NVE, for good energy conservation use coulombtype = PME-Switch or Reaction-Field-zero ------ This warning is also discussed here: http://www.gromacs.org/Documentation/Terminology/NVE And, indeed, just as Gromacs says, I get terrible energy conservation using coulombtype = PME (I use g_energy to extract the "Total energy" from the .edr file). I have posted plots of my results using coulombtype = PME here: http://www.andrew.cmu.edu/user/adeyoung/july13/july13.pdf . As you can see in the light red plot on the right-hand side of that PDF file, energy conservation is terrible using coulombtype = PME, even when I do a 2 ns run; total energy appears to increase nearly linearly as a function of time. If you have time, I am wondering if you have any advice about the following two questions: 1. NOTE 1 above suggests that I use vdwtype = Shift. When I do this, do you recommend that I apply long range dispersion corrections for both energy and pressure, using DispCorr = EnerPres, or for only energy, using DispCorr = Ener? Typically, for various (non-NVE) calculations, I have been using DispCorr = no, but I am not sure if this is a good idea. Pages 97-98 of the Gromacs 4.5.4 manual seem to suggest that the energy correction due to DispCorr is small and usually only significant for free energy calculations (which I will not be doing here). As a rule of thumb, do you typically turn dispersion corrections off? 2. NOTE 2 above suggests that I use either coulombtype = PME-Switch or coulombtype = Reaction-Field-zero. Do you have any advice or recommendation? My concern (and I am very much a novice in MD) about using coulombtype = Reaction-Field-zero is that the manual says (http://manual.gromacs.org/current/online/mdp_opt.html#el) that it can only be used with an infinite dielectric constant. I know that water has a relatively high dielectric constant (~80 at ~298 K), but strictly speaking it is not infinite. On the other hand, about coulombtype = PME-Switch or Switch, the manual says (http://manual.gromacs.org/current/online/mdp_opt.html#el) that "Switching the Coulomb potential can lead to serious artifacts, advice: use Reaction-Field-zero instead." It doesn't explicitly specify what kinds of serious artifacts I can expect, but it does sound like a dangerous choice. Do you have any experience with NVE simulations? Thank you very much for your time! Andrew DeYoung Carnegie Mellon University -- gmx-users mailing list gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/Support/Mailing_Lists