[CP2K:4029] Re: Energy conservation NVT/NVE simulation

[paolo nicolini]

Hi

I tried to use the H2O-32_SPME.inp in the tests/Fist directory by changing 
only the length of the simulation (1 ns) and I got also a significant 
energy drift (see attachment). Then I tried to reduce the integration time 
step to 1.0 fs. The drift become lower but it do not vanish. I did also 1ns 
of NVE dynamics (starting from the last restart of the previous 
simulations) and here the drift become even more clear.
The inputs that I used are below.
Best regards
Paolo

NVT input: 
&FORCE_EVAL
  METHOD FIST
  &MM
    &FORCEFIELD
      &BEND
        ATOMS H O H
        K 0.
        THETA0 1.8
      &END BEND
      &BOND
        ATOMS O H
        K 0.
        R0 1.8
      &END BOND
      &CHARGE
        ATOM O
        CHARGE -0.8476
      &END CHARGE
      &CHARGE
        ATOM H
        CHARGE 0.4238
      &END CHARGE
      &NONBONDED
        &LENNARD-JONES
          atoms O O
          EPSILON 78.198
          SIGMA 3.166
          RCUT 11.4
        &END LENNARD-JONES
        &LENNARD-JONES
          atoms O H
          EPSILON 0.0
          SIGMA 3.6705
          RCUT 11.4
        &END LENNARD-JONES
        &LENNARD-JONES
          atoms H H
          EPSILON 0.0
          SIGMA 3.30523
          RCUT 11.4
        &END LENNARD-JONES
      &END NONBONDED
    &END FORCEFIELD
    &POISSON
      &EWALD
        EWALD_TYPE spme
        ALPHA .5
        GMAX 21
        O_SPLINE 6
      &END EWALD
    &END POISSON
  &END MM
  &SUBSYS
    &CELL
      ABC 9.865 9.865 9.865
    &END CELL
    &COORD
    O                  -4.583   5.333   1.560   H2O
    H                  -3.777   5.331   0.943   H2O
    H                  -5.081   4.589   1.176   H2O
...
    &END COORD
  &END SUBSYS
&END FORCE_EVAL
&GLOBAL
  PROJECT H2O
  RUN_TYPE MD
&END GLOBAL
&MOTION
  &CONSTRAINT
    &G3X3
      DISTANCES 1.8897268 1.8897268 3.0859239
      MOLECULE 1
      ATOMS 1 2 3
    &END G3X3
  &END CONSTRAINT
  &MD
    ENSEMBLE NVT
    STEPS 1000000
    TIMESTEP 1.0
    TEMPERATURE 300.0
    &THERMOSTAT
      REGION MOLECULE
      &NOSE
        LENGTH 3
        YOSHIDA 3
        TIMECON 1000
        MTS 2
      &END NOSE
    &END
  &END MD
&END MOTION

NVE input:
&FORCE_EVAL
  METHOD FIST
  &MM
    &FORCEFIELD
      &BEND
        ATOMS H O H
        K 0.
        THETA0 1.8
      &END BEND
      &BOND
        ATOMS O H
        K 0.
        R0 1.8
      &END BOND
      &CHARGE
        ATOM O
        CHARGE -0.8476
      &END CHARGE
      &CHARGE
        ATOM H
        CHARGE 0.4238
      &END CHARGE
      &NONBONDED
        &LENNARD-JONES
          atoms O O
          EPSILON 78.198
          SIGMA 3.166
          RCUT 11.4
        &END LENNARD-JONES
        &LENNARD-JONES
          atoms O H
          EPSILON 0.0
          SIGMA 3.6705
          RCUT 11.4
        &END LENNARD-JONES
        &LENNARD-JONES
          atoms H H
          EPSILON 0.0
          SIGMA 3.30523
          RCUT 11.4
        &END LENNARD-JONES
      &END NONBONDED
    &END FORCEFIELD
    &POISSON
      &EWALD
        EWALD_TYPE spme
        ALPHA .5
        GMAX 21
        O_SPLINE 6
      &END EWALD
    &END POISSON
  &END MM
&END FORCE_EVAL
&GLOBAL
  PROJECT H2O
  RUN_TYPE MD
&END GLOBAL
&MOTION
  &CONSTRAINT
    &G3X3
      DISTANCES 1.8897268 1.8897268 3.0859239
      MOLECULE 1
      ATOMS 1 2 3
    &END G3X3
  &END CONSTRAINT
  &MD
    ENSEMBLE NVE
    STEPS 1000000
    TIMESTEP 1.0
    TEMPERATURE 300.0
  &END MD
&END MOTION
&EXT_RESTART
  RESTART_FILE_NAME H2O-1_1000000.restart 
&END EXT_RESTART



On Monday, September 17, 2012 1:09:54 PM UTC+2, jgh wrote:
>
> Hi 
>
> I did some tests with the inputs from 
> ../tests/Fist 
> and get much better energy conservation that you show. 
> Probably your input use a not optimal setting for the 
> cutoffs and Ewald parameters (I was using SPME). 
>
> regards 
>
> Juerg 
>
> -------------------------------------------------------------- 
> Juerg Hutter                         Phone : ++41 44 635 4491 
> Physical Chemistry Institute   FAX   : ++41 44 635 6838 
> University of Zurich               E-mail:  hut... at pci.uzh.ch<javascript:> 
> Winterthurerstrasse 190 
> CH-8057 Zurich, Switzerland 
> --------------------------------------------------------------- 
>
> -----cp... at googlegroups.com <javascript:> wrote: ----- 
> To: cp... at googlegroups.com <javascript:> 
> From: marco masia 
> Sent by: cp... at googlegroups.com <javascript:> 
> Date: 09/12/2012 03:39PM 
> Subject: [CP2K:4029] Re: Energy conservation NVT/NVE simulation 
>
> Dear all, 
>
> recently I have noticed the same problem with FIST simulations. I took the 
> H2O-32_ewald.inp input, equilibrated it for 1 ns and run a 1ns NVE 
> simulation from the thermalized configuration. I have tried to refine many 
> simulation parameters (e.g. time step, RCUT, Ewald pars) and also (thinking 
> that there was a bad implementation of constraints) I run the simulations 
> for a flexible model. I always get a non-negligible drift (see attachment); 
> the graph shows block averaged energies of 1 ns simulations performed with 
> cp2k and obeli.x (a program of mine). The values are shifted for a better 
> comparison. I have been looking into this, thinking that the problem was in 
> the input file...but since other people are experiencing it with QS, I 
> wonder if there is some bug in the velocity verlet routines. The results 
> are platform independent (I tried both on mac and linux). 
>
> best 
>
> marco 
>
> On Wednesday, September 5, 2012 3:32:26 AM UTC-4, alber... at virgilio.itwrote: 
> Dear all, 
>
> I'm carrying out some tests to reproduce the VDOS of water and I used as 
> starting points the input files 'H2O-32.inp' and 'H2O-64.inp' in 
> tests/QS/benchmark directory. 
>
> In both cases I ran a NVT simulation at 300K for 15 ps followed by an NVE 
> simulation on the equilibrated system. 
>
> I'm puzzled about the results in the energies: in both NVT and NVE  the 
> plot of 
> "Cons Qty" show a constant and linear drift while temperature, kinetic and 
> potential energy are equilibrated and stable. 
>
> Since "Cons Qty" should be conserved and this is not the case I think that 
> there are some problems in the setup of the simulation. 
> I did not modify the parameters of the benchmark inputs but I just added 
> the 
> thermostat section in NVT and I modified only the printing option. 
>
> Thank you very much for your help and suggestions, 
>
> Alberto 
>
>   
>   -- 
>  You received this message because you are subscribed to the Google Groups 
> "cp2k" group. 
>  To view this discussion on the web visit 
> https://groups.google.com/d/msg/cp2k/-/sAiVZ2dob1sJ. 
>   To post to this group, send email to cp... at googlegroups.com<javascript:>. 
>
>  To unsubscribe from this group, send email to 
> cp2k+... at googlegroups.com <javascript:>. 
>   For more options, visit this group at 
> http://groups.google.com/group/cp2k?hl=en. 
>     
>
>
> [attachment "energy-conservation.pdf" removed by Jürg Hutter/at/UZH]

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20120921/07d1d475/attachment.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: all.pdf
Type: application/pdf
Size: 1664599 bytes
Desc: not available
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20120921/07d1d475/attachment.pdf>