[gmx-developers] Conserved energy drifts more with shorter time steps

David van der Spoel spoel at xray.bmc.uu.se
Wed Jun 3 09:57:04 CEST 2015


On 02/06/15 15:24, Berk Hess wrote:
> On 2015-06-02 15:20, Anton Feenstra wrote:
>> On 02-06-15 15:14, Berk Hess wrote:
>>> On 2015-06-02 15:08, David van der Spoel wrote:
>>>> On 02/06/15 14:30, Berk Hess wrote:
>>>>> On 2015-06-02 14:26, David van der Spoel wrote:
>>>>>> On 02/06/15 14:05, Berk Hess wrote:
>>>>>>> I assume David is looking at the conserved energy quantity.
>>>>>> Indeed, averages and drift over 900 ps, 1000 molecules:
>>>>>>
>>>>>> Energy            Average   Err.Est.       RMSD  Tot-Drift
>>>>>> -------------------------------------------------------------------------------
>>>>>>
>>>>>>
>>>>>>
>>>>>> Total Energy      20460.7        2.5    778.213    12.4265 (kJ/mol)
>>>>>> Conserved En.     -4384.7        380    778.283    2665.29 (kJ/mol)
>>>>>>
>>>>>>>
>>>>>>> But I just noticed that Verlet - buffer - tolerance is not
>>>>>>> mentioned,
>>>>>>> so I assume it is default. The drift estimate is an overestimate
>>>>>>> which is tighter for shorter pair list lifetimes. I guess that is
>>>>>>> the
>>>>>>> explanation. Set it very low, so you are sure the energy drift is
>>>>>>> not
>>>>>>> affected by the (too small) buffer.
>>>>>>>
>>>>>> This may be the default.
>>>>>>    verlet-buffer-tolerance        = 0.005
>>>>>> should we set it to 1e-5 or something like that?
>>>>> It depends on how low drift you want to be able to measure. You
>>>>> need to
>>>>> set it lower than the drift you want to be able to measure.
>>>> Ok, will try. What about the time step dependence? Any clue?
>>> I explained that before, see 4 steps above.
>>>
>>> Berk
>>>>
>>>>
>>>>>
>>>>> Berk
>>>>>>
>>>>>>> Berk
>>>>>>>
>>>>>>> On Jun 2, 2015 1:45 PM, "Shirts, Michael R. (mrs5pt)"
>>>>>>> <mrs5pt at eservices.virginia.edu> wrote:
>>>>>>>>
>>>>>>>>> How do we run more NVE than with the settings below?
>>>>>>>>
>>>>>>>>     tcoupl                         = Nose-Hoover
>>>>>>>>
>>>>>>>>
>>>>>>>> Would make it not NVE.
>>>>>> Turning tcoupl off would not keep the temperature constant of course,
>>>>>> seeing that there is a drift in the conserved energy this means that
>>>>>> the temperature will change.
>>>>>>
>>>>>> However, the main question was why does "conserved-energy
>>>>>> conservation" get worse when the time step is reduced? Integration
>>>>>> should presumably become more accurate, or not? Or does one need
>>>>>> double precision for short time steps?
>>
>>
>> Yes, at least that is what I saw for my and Berk's 1999 JCC paper.
>> At single precision, drift would not go down below a certain time
>> step. I have no idea anymore which timestep precisely, and I'm sure
>> this would be system specific.
>>
>> I also recall drift indeed going up at (extreme) short timesteps. My
>> intuition for that at the time was that arithmetic errors would
>> accumulate faster for calculating many very small delta's.
> In that case it was due to the constraint algorithms. The velocity was
> corrected with the coordinate difference divided by the time step. Thus
> bit noise is amplified inversely proportional with the time step. I
> fixed this some years ago for LINCS and SHAKE, but SETTLE still has this
> property.
>
> David is not using constraints. So I think the issue here is, as I
> mailed before, that the pairlist buffer estimate is more less accurate,
> more of an overestimate, for longer pairlist lifetimes (larger time
> steps with fixed nslist), which leads to less drift.
>
> Berk
Unfortunately setting the verlet-buffer to 1e-6 (instead of default 
0.005) does not help. We get the following drift:

Drift:
Dt   Single(5e-3)  Double(5-e3)  Single(1e-6)
1.0    240            166          191
0.5     46            -12           68
0.2    170            -7           137
0.1    328            -6           353
0.005  564            -4           573

Seems this is still a precision issue, or not?

>>
>>
>>>>>>
>>>>>>
>>>>>>>>
>>>>>>>> Do you mean the conserved quantity isn't conserved? Energy
>>>>>>>> certainly
>>>>>>>> won't be conserved with Nose-Hoover.
>>>>>>>>
>>>>>>>> Best,
>>>>>>>> ~~~~~~~~~~~~
>>>>>>>> Michael Shirts
>>>>>>>> Associate Professor
>>>>>>>> Department of Chemical Engineering
>>>>>>>> University of Virginia
>>>>>>>> michael.shirts at virginia.edu
>>>>>>>> (434) 243-1821
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> On 6/2/15, 7:22 AM, "David van der Spoel" <spoel at xray.bmc.uu.se>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> On 02/06/15 12:47, Berk Hess wrote:
>>>>>>>>>> Hi,
>>>>>>>>>>
>>>>>>>>>> So this is not using SETTLE or LINCS?
>>>>>>>>>> With SETTLE this is a known issue.
>>>>>>>>> No this is flexible organic molecules at contant volume.
>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Did you try running NVE?
>>>>>>>>> How do we run more NVE than with the settings below?
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Cheers,
>>>>>>>>>>
>>>>>>>>>> Berk
>>>>>>>>>>
>>>>>>>>>> On 2015-06-02 12:43, David van der Spoel wrote:
>>>>>>>>>>> Hi,
>>>>>>>>>>>
>>>>>>>>>>> this is in between a user and a developer query. We find that
>>>>>>>>>>> for
>>>>>>>>>>> flexible liquid simulations using gromacs 5.0.4/single
>>>>>>>>>>> precision the
>>>>>>>>>>> energy conservation gets WORSE with decreasing time step. A plot
>>>>>>>>>>> showing this is in
>>>>>>>>>>> http://folding.bmc.uu.se/images/Econserved-timestep.xvg or
>>>>>>>>>>> http://folding.bmc.uu.se/images/Econserved-timestep.pdf
>>>>>>>>>>>
>>>>>>>>>>> MDP settings
>>>>>>>>>>>      coulombtype                    = PME
>>>>>>>>>>>      coulomb-modifier               = Potential-shift
>>>>>>>>>>>      rcoulomb-switch                = 0
>>>>>>>>>>>      rcoulomb                       = 1.1
>>>>>>>>>>>      epsilon-r                      = 1
>>>>>>>>>>>      epsilon-rf                     = inf
>>>>>>>>>>>      vdw-type                       = PME
>>>>>>>>>>>      vdw-modifier                   = Potential-shift
>>>>>>>>>>>      rvdw-switch                    = 0
>>>>>>>>>>>      rvdw                           = 1.1
>>>>>>>>>>>      tcoupl                         = Nose-Hoover
>>>>>>>>>>>      nsttcouple                     = 10
>>>>>>>>>>>      nh-chain-length                = 1
>>>>>>>>>>>      pcoupl                         = No
>>>>>>>>>>>      ref-t:      298.15
>>>>>>>>>>>      tau-t:         0.5
>>>>>>>>>>>
>>>>>>>>>>> Any clues whether we are doing something wrong? Or is there a
>>>>>>>>>>> bug?
>>>>>>>>>>>
>>>>>>>>>>> Cheers,
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> --
>>>>>>>>> David van der Spoel, Ph.D., Professor of Biology
>>>>>>>>> Dept. of Cell & Molec. Biol., Uppsala University.
>>>>>>>>> Box 596, 75124 Uppsala, Sweden. Phone: +46184714205.
>>>>>>>>> spoel at xray.bmc.uu.se    http://folding.bmc.uu.se
>>>>>>>>> --
>>>>>>>>> Gromacs Developers mailing list
>>>>>>>>>
>>>>>>>>> * Please search the archive at
>>>>>>>>> http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_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-developers
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> or send a mail to gmx-developers-request at gromacs.org.
>>>>>>>>
>>>>>>>> --
>>>>>>>> Gromacs Developers mailing list
>>>>>>>>
>>>>>>>> * Please search the archive at
>>>>>>>> http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_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-developers
>>>>>>>>
>>>>>>>>
>>>>>>>> or send a mail to gmx-developers-request at gromacs.org.
>>>>>>
>>>>>>
>>>>>
>>>>
>>>>
>>>
>>
>>
>


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


More information about the gromacs.org_gmx-developers mailing list