[gmx-users] Doubt about energies in a very simple system

[Mark Abraham]

Hi,

I'll be more direct this time. What potential shifts are the other software
packages using? GROMACS uses one by default.

Mark

On Fri, Feb 27, 2015 at 3:06 PM, IÑIGO SAENZ <
inigo.saenz01 at estudiant.upf.edu> wrote:

> Hi,
>
> Thank you all for your anwsers.
> Let me explain the circumstances that sorround this. I'm using a
> developer's version of ACEMD, which grant me access to its internal data
> structures.
> When I run this version of acemd, it gives me the information of all the
> atoms, bonds, angles, dihedrals, exclusion list, 1-4 pairs and running
> parameters(cutoffs, switching mode, etc...) of the system that I'm running,
> with this information I can create the .gro, .top and .mdp
>
> There may be something weird in the .top file that I showed you previously
> in this post. If you take a look at this:
>
>
> [ moleculetype ]
> ;name            nrexcl
> sys      0
>
> it may look strange this nrexcl = 0
> but actually it make sense because what I explained some lines above that
> ACEMD gives me the exclusion list, so putting this nrexcl = 0, gromacs
> calculates the interactions for all those pair of atoms within the cutoff
> radius except those that are in the exclusion list.
>
> Let me show you a pair of real cases that I've tried today.
>
>
> 400 WATER MOLECULES
>
> ACEMD
>
> NAMD
>
> GROMACS
>
>
>
> LJ
>
> 2540.0544
>
> 2538.2893
>
> 2538.28
>
> COULOMB
>
> -3996.3431
>
> -3993.6970
>
> -3993.17
>
>
>
>
>
>
> This first case is a SPE simulation of a box full of water molecules,
> nothing more, only 400 water molecules,
> as you can see, result are pretty much the same, energies are the same for
> the three programs (well, acemd deviates 2kj/mol in LJ and 3kj/mol in
> Coulomb, but this doesn't matter really)
>
> I want to highlight that water molecules don't have 1-4 pairs, only
> exclusions. I say this because I'm sure that the problem in the following
> case (actually in all cases) is caused by the 1-4 interactions.
>
>
> 1A1X
>
> ACEMD
>
> NAMD
>
> GROMACS
>
> BOND
>
> 201.9381
>
> 201.4818
>
> 201.809
>
> ANGLE
>
> 1318.6907
>
> 1317.691
>
> 1317.8
>
> DIHEDRAL
>
> 5225.9049
>
> 5222.46
>
> 5222.41
>
> LJ
>
> 25404.8407
>
> 25391.6948
>
> 25386.68
>
> COULOMB
>
> -219604.180
>
> -219387.922
>
> -231465.3
>
> This second case is the SPE calculation of a pdb (code 1a1x) sourrounded by
> water.
> As you can see every energy term coincide between the three programs (with
> the logical deviation of using distinct software),
> but Coulomb result don't make sense in the gromacs case, it deviates more
> than 11000KJ/mol respect to acemd or namd.
> I put here the Lj and Coulomb energies in the way that g_energy gives the
> results:
>
> LJ(SR)       = 27393.2
> LJ-14        = 2006.52
> Coulomb(SR)  = -216459
> Coulomb-14   = 15006.3
>
> Looking back to the water example, which hasn't got 1-4 interactions, I
> imagine that the Coulomb(SR) energy term is well calculated (in the water
> system coulomb energies coincide) and that the problem comes with the
> Coulomb-14 energy term, but this doesn't make much sense because in that
> case, LJ energies should also be bad...
>
> Thank you very much for your attention, any help or idea of what's going
> wrong would be highly appreciated.
>
>
>
> Iñigo Sáenz
> Universitat Pompeu Fabra
>
>
>
> 2015-02-26 2:41 GMT+01:00 Justin Lemkul <jalemkul at vt.edu>:
>
> >
> >
> > On 2/24/15 3:23 PM, IÑIGO SAENZ wrote:
> >
> >> Hi Justin,
> >>
> >> I always do the SPE as follows:
> >> grompp -f SPE.mdp -p sys.top -c sys.gro
> >> and after that I simply execute mdrun, i didn't know about the mdrun
> >> -rerun
> >> function.
> >>
> >> Now I have done: mdrun -s topol.tpr -rerun sys.gro
> >>
> >> but the energy results are the exactly the same.
> >>
> >>
> >>
> > Start with something simpler; an actual molecule with normal bonded
> > interactions.  Equivalency of energy between GROMACS and many programs
> has
> > been shown many times over, so it's certainly possible to prove.  Likely
> > there's just something weird about what you're doing with your
> combinations
> > of exclusions, pairs, etc.
> >
> >
> > -Justin
> >
> > --
> > ==================================================
> >
> > Justin A. Lemkul, Ph.D.
> > Ruth L. Kirschstein NRSA Postdoctoral Fellow
> >
> > Department of Pharmaceutical Sciences
> > School of Pharmacy
> > Health Sciences Facility II, Room 629
> > University of Maryland, Baltimore
> > 20 Penn St.
> > Baltimore, MD 21201
> >
> > jalemkul at outerbanks.umaryland.edu | (410) 706-7441
> > http://mackerell.umaryland.edu/~jalemkul
> >
> > ==================================================
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at http://www.gromacs.org/
> > Support/Mailing_Lists/GMX-Users_List before posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
> >
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>