[gmx-users] Re: energy conservation / frozen atoms

[S. Alireza Bagherzadeh]

Thanks for your notes.


I did a diagnosis test which could be of relevance here.

I set up the following system:
[ gas | liquid water (solid water) liquid water | gas ]

gas is united atom methane.
liquid water is tip4p-ice model and solid water is a cage-like crystalline
structure of water and methane called gas hydrate.

Now, in order to test the effect of freezing and position restraining on
the performance of nve I did two tests at 370 K.

Test 1 (freezing):
Solid water was kept frozen in all 3 dimensions (Y Y Y).
First I ran a nvt for 250 ps for equilibration (potential and total energy
both converged after 250 ps, Pressure equilibrated at ~ 3950 bar). Then I
started a 1ns nve.
Similar to my other simulation, the total energy linearly decreased (0.84%
per ns) as well as potential energy. Pressure remained around 3950 bar;
however, the temperature decreased from 370 to 364 K (physically, this
should not happen).


Test 2 (position restraining):
Oxygen of solid water was strongly restrained to a point (fc of 100000).
Similar to the previous test, first I ran a nvt for 250 ps for
equilibration (potential and total energy both converged after 250 ps,
Pressure equilibrated at about 0 bar with fluctuations of ~ 2000 bar). Then
I started a 1ns nve.
Again, similar to test 1, the total energy linearly decreased (1.33% per
ns) as well as potential energy. Pressure remain around 0 bar; however, the
temperature initially dropped from 370 K to 355K within 1 ps, then
increased to 358 K during the next 50 ps and thereafter kept linearly
decreasing to 353 K until the end of 1 ns run (physically and intuitively,
this should not happen).



(In  both of the tests, I kept the methane inside the cages of solid water
position-restrained to a point by fc = 1000).

If needed I can post the .mdp and .top files too.


I hope this will be helpful from the development standpoint.


Regards,
Alireza





>>
> >> ------------------------------
> >>
> >> Message: 3
> >> Date: Tue, 30 Jul 2013 07:20:53 -0400
> >> From: Justin Lemkul <jalemkul at vt.edu>
> >> Subject: Re: [gmx-users] energy conservation / frozen atoms
> >> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
> >> Message-ID: <51F7A195.80404 at vt.edu>
> >> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
> >>
> >>
> >>
> >> On 7/29/13 10:51 PM, S. Alireza Bagherzadeh wrote:
> >>> Hi All,
> >>>
> >>> I am simulating a system in which I have two solid surfaces and I keep
> >> them
> >>> frozen during simulations. I also exclude the interactions between its
> >>> atoms to avoid spurious contribution to the virial pressure due to
> large
> >>> forces between them as suggested in the manual.
> >>>
> >>> I run a nvt for equilibration and then I do the production run in an
> nve
> >>> ensemble; however, I am not getting good energy conservation. There is
> a
> >>> huge energy drift...
> >>>
> >>>
> >>> When I remove the solid surfaces, I will only have water molecules and
> >>> united atom methane molecules in my system. Using the same protocol I
> >>> obtain a very good energy conservation...
> >>>
> >>
> >> What happens if you unfreeze the frozen surfaces and run the same
> system?
> >>
> >
> > The point is I do not want to do this as I am trying to save some
> > calculations.
>
> My point in asking is to do some diagnostics.  There are countless similar
> posts
> to yours - freeze some group(s), energy conservation is bad, and then a
> dead end
> to the thread.  It would be interesting, from a development standpoint, to
> know
> whether the issue is truly linked to the frozen group(s).  Using strong
> position
> restraints in lieu of frozen groups is probably a reasonable workaround
> that
> could also be used to find a problem.
>
> > My system is already big and including all of the bond, angles and
> dihedral
> > interactions of the solid surfaces slows down the simulation
> considerably.
> >
>
> The bonded interactions are not likely what slows the system down.
>  Freezing
> doesn't actually improve performance.  Using energygrp_excl is what causes
> a
> speed-up.
>
> >
> >
> >>
> >>> Any insight on what might be wrong in my system would be very
> >> appreciated.
> >>>
> >>
> >> The contents of the .edr file will probably be informative, as you can
> >> identify
> >> which energy term(s) is(are) most affected.  It's all probably related
> to
> >> the
> >> frozen surfaces themselves acting as an energy sink or something.
> >>
> >
> > Is there anyway to fix the "energy sinking" effect of the surfaces?
>
> That's just an off-the-cuff guess at something being wrong.  Like I said
> above,
> there's been no thorough evaluation of whether or not this sort of
> combination
> does (or should) work.
>
> > I used dl_poly to simulate a similar system and there was no problem with
> > the energy conservation.
> >
>
> That's been reported before, hence something is probably wrong with
> Gromacs, but
> we need more information to figure it out and fix it.  Otherwise, you'll
> have to
> figure out an alternate strategy.
>
> -Justin
>
> --
> ==================================================
>
> Justin A. Lemkul, Ph.D.
> Postdoctoral Fellow
>
> Department of Pharmaceutical Sciences
> School of Pharmacy
> Health Sciences Facility II, Room 601
> University of Maryland, Baltimore
> 20 Penn St.
> Baltimore, MD 21201
>
> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>
> ==================================================
>
>