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

[IÑIGO SAENZ]

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.
>