[gmx-users] SETTLE vs. LINCS -- different final energies of energy minimized structures

[Justin Lemkul]

On 4/18/16 10:30 AM, Rakesh Sharan wrote:
> Thanks very much Justin.
>
> Actually I am minimization a series of configurations from an equilibrated
> trajectory of bulk TIP4P/2005 (LINCS used as constraint algorithm). I need

So, to be clear, you altered the TIP4P topology to use a system of 3 constraints 
that were held rigid by LINCS, using "constraints = all-bonds" when doing the 
simulation?  Simply setting LINCS as the constraint method in the .mdp file is 
insufficient to turn off SETTLE (which is always used unless you tell mdrun not to).

> the final energy of the energy minimized structures, however, you can see
> that depending on the choice of constrained algorithm, for few
> configurations, the final energy is different (e. g. for the attached plot
> -29146 and -29288 kj/mol). Actually, both are following same minimization

The list doesn't allow attachments.

Re-minimizing a configuration with different algorithms may lead to some 
differences, I would suspect.  But the internal geometry should be the same, so 
all that's varying will be LJ and electrostatic terms.  You can check by doing a 
further decomposition into which term is changing the most.

> trajectory, however, one is finishing before other. This is bit troubling
> as depending on the choice of constraint algorithm average final energy is
> different even though initial starting configuration is exactly same. I
> rechecked topology files and they look same apart from constraint
> specifications.
>

So let me clarify, because I'm not sure if I follow what you're saying (because 
here it sounds like you're talking about running different simulations, whereas 
above it sounds like a very different case).  Which is true:

1. You are doing two simulations from the same starting configuration and are 
trying to compare the final energy.
2. You are re-minimizing each frame of an existing trajectory and comparing 
energies.

#1 really doesn't make sense, as all MD simulations are chaotic, and I'm not 
sure what #2 would achieve.

> I would very much appreciate your view over the plausible origin of this
> discrepancy. Also, I am wondering whether in both the cases one is doing
> truly constrained minimization.
>

If there are no bond or angle energy values, the waters were constrained.

-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

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