[gmx-users] Re: Energy Drift in Gromacs 4.0

[David van der Spoel]

Ilya Chorny wrote:
> Ok, I got it. Set electrostatics to PME-switch. Made rcoulomb-switch .1 
> less than rcoulomb and set rlist accordingly and now I am conserving 
> energy. Yipee!!!
> 
> It did not make sense to me that the system had to be fully relaxed to 
> conserve energy. It should conserve energy energy as long as it has been 
> minimized and has no major clashes. 
This is correct, provided T and P coupling are turned off.

> 
> Thanks,
> 
> Ilya
> 
> 
> On Wed, Mar 11, 2009 at 7:38 PM, Ilya Chorny <ichorny at gmail.com 
> <mailto:ichorny at gmail.com>> wrote:
> 
>     If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE
>     simulation they should counter balance. I am leaking K.E. somewhere.
>     I reran with a relaxed P.E(t) and saw the same behavior.  This is
>     happening both at 2 and 4 fs so its not a time step issue. Any
>     thoughts?  I thought it might be COM removal so I turned it off to
>     no avail. 
> 
>     Any chance one of the Gromacs people can give me the mdp file from
>     the DHFR simulation in the 4.0 paper. 
> 
> 
>     On Tue, Mar 10, 2009 at 10:21 PM, Mark Abraham
>     <Mark.Abraham at anu.edu.au <mailto:Mark.Abraham at anu.edu.au>> wrote:
> 
>         Ilya Chorny wrote:
> 
>             Forgot column labels on energies.
>                                                                    
>             Average       RMSD     Fluct.      Drift  Tot-Drift
>             2 fs time step Total Energy                46497.2  
>              36445.6     2630.7  -148.833    -148833
>             4 fs time step Total Energy                41580.8  
>              37693.4    2856.08  -130.198    -130198
> 
>             On Tue, Mar 10, 2009 at 7:17 PM, Ilya Chorny
>             <ichorny at gmail.com <mailto:ichorny at gmail.com>
>             <mailto:ichorny at gmail.com <mailto:ichorny at gmail.com>>> wrote:
> 
>                Hello All,
> 
>                I am trying to run some calibration calculations with 2/4
>             fs time
>                steps. I am trying to reproduce the results in the
>             Gromacs 4.0 paper
>                on a protein/water (not the same protein) system with
>             ~100K atoms. I
>                ran 1 ns simulation in the NVE ensemble. My mdp params
>             are shown
>                below. Topolgy was created using -vsite h. My system has
>             218441
>                degrees of freedom according to GMXDUMP. My P.E. energy
>             is ~ -1E5
>                and has not full relaxed after 1ns but is close.
>                2 fs time step Total Energy                46497.2  
>              36445.6        2630.7   -148.833    -148833
>                4 fs time step Total Energy                41580.8  
>              37693.4       2856.08   -130.198    -130198
> 
>                Drift (kt/ns) 2fs/4fs = .34/.29
> 
>                The drift seems to be much larger for both the 2fs and
>             4fs time
>                steps then in the paper. Do I have a clear bug in my
>             params? Should
>                I wait till the system fully relaxes before doing the
>             drift calculation?
> 
> 
>         You should definitely equilibrate the system before trying to
>         measure energy drift. I suggest doing so in NPT to fix your
>         density and temperature, then switching to NVE, wait a bit, and
>         only then start to collect your data.
> 
>                coulombtype              = PME
>                rcoulomb                 = 1.0
> 
>                fourierspacing           = 0.12
>                ; FFT grid size, when a value is 0 fourierspacing will be
>             used
>                fourier_nx               = 0
>                fourier_ny               = 0
>                fourier_nz               = 0
>                ; EWALD/PME/PPPM parameters
>                pme_order                = 4
>                ewald_rtol               = 1.e-04
> 
> 
>         I would think the ewald_rtol value used in the paper was 1e-5,
>         but it is unfortunately not stated. The different numbers will
>         lead to different costs and accuracies for the direct and FFT
>         parts of the PME algorithm.
> 
>         rcoulomb, fourier_whatever, pme_order and ewald_rtol are all
>         relevant if you're trying to reproduce results.
> 
>         Mark
>         _______________________________________________
>         gmx-users mailing list    gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>
>         http://www.gromacs.org/mailman/listinfo/gmx-users
>         Please search the archive at http://www.gromacs.org/search
>         before posting!
>         Please don't post (un)subscribe requests to the list. Use the
>         www interface or send it to gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>.
>         Can't post? Read http://www.gromacs.org/mailing_lists/users.php
> 
> 
> 
> 
>     -- 
>     Ilya Chorny Ph.D.
> 
> 
> 
> 
> -- 
> Ilya Chorny Ph.D.
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/search before posting!
> Please don't post (un)subscribe requests to the list. Use the 
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php


-- 
David van der Spoel, Ph.D., Professor of Biology
Molec. Biophys. group, Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205. Fax: +4618511755.
spoel at xray.bmc.uu.se	spoel at gromacs.org   http://folding.bmc.uu.se