[CP2K:6254] Energy conservation in periodic QM/MM

Matt W MattWa... at gmail.com
Thu Mar 19 11:17:54 UTC 2015


Try increasing the cutoff to 600-800 Ry? BLYP functional is known to be 
noisy.

I have repeated the QM/MM MD using MULTIPOLE OFF using the latest CP2K 
> trunk binary (svn:15171). The total energy from this new run is *identical* 
> to the previous one using MULTIPOLE section. The only difference is that 
> the decoupling/recoupling energies and the warning messages about QM atoms 
> close to the boundary are not printed in the output file any more.
>
> Dorothea pointed out that, _if_ your QM box is the same size as the QMMM 
box (as you have in the input above), you don't need the de/recoupling - it 
_should_ give identical results. It just slows the calculation a little.

HTH

Matt

 

> I have also tried the "CENTER_GRID T" keyword together with CENTER 
> EVERY_STEP as Matt suggested above and its short time energy profile up to 
> 100 fs is nearly the same as the one without CENTER_GRID.  I guess the 
> total energy is still adversely affected by the centering even with the 
> CENTER_GRID.
>
> From these observations, I'm afraid it's likely that the current CP2K 
> cannot handle gracefully the QM atoms crossing boundary. 
> I'd appreciate any comment or suggestion.
>
> Jonggu
>
>
> 2015년 3월 18일 수요일 오후 9시 13분 29초 UTC+9, Matt W 님의 말:
>>
>> Hi,
>>
>> there is also the CENTER_GRID that has been introduced quite recently. I 
>> don't know how much it has been tested, maybe someone else can comment? But 
>> for my simple test water system it seems to improve the energy conservation 
>> considerably. The idea is, I think, to only center the QM molecule in a 
>> "quantized" manner on equivalent grid points of the (fine DFT?) grid.
>>
>> To get really good energy conservation you might need to increase your 
>> cutoff/rel_cutoff too.
>>
>> Matt
>>
>> On Wednesday, March 18, 2015 at 11:22:10 AM UTC, Jonggu Jeon wrote:
>>>
>>> Dear Dorothea,
>>>
>>> Thank you for the invaluable suggestion.  
>>> While your answer makes perfect sense to me, I just found out that CP2K 
>>> 2.5.1 does not have option to turn off QMMM/PERIODIC/MULTIPOLE using the 
>>> section parameter OFF.
>>> Therefore, CP2K 2.5.1 produces the same (diverging) energy profile 
>>> whether I set &MULTIPOLE OFF or not.
>>>
>>> I will install version 2.6 and try it again.
>>>
>>> Best regards,
>>>
>>> Jonggu
>>>
>>>
>>> 2015년 3월 18일 수요일 오후 6시 22분 40초 UTC+9, Dorothea Golze 님의 말:
>>>>
>>>> Hi,
>>>>
>>>> if you have a fully periodic setup where the QM and MM boxes are equal, 
>>>> you do not need the decoupling/recoupling stuff, i.e. switch off explicitly 
>>>> the MULTIPOLE section like that
>>>> &MULTIPOLE OFF
>>>> &END
>>>>
>>>> Also set in the QMMM section
>>>>
>>>>  &WALLS
>>>>      TYPE NONE
>>>>   &END WALLS
>>>>
>>>> That should address your first question...
>>>>
>>>> Cheers,
>>>> Dorothea
>>>>  
>>>>
>>>> 2015-03-18 8:46 GMT+01:00 Jonggu Jeon <jeon... at gmail.com>:
>>>>
>>>>> Dear CP2K users,
>>>>>
>>>>> I'd like to ask some questions on the capabilities of the QM/MM part 
>>>>> of CP2K (I'm using version 2.5.1, svn:13632).
>>>>>
>>>>> I've been attempting a QM/MM MD simulation with CP2K. My system 
>>>>> consists of one NMA molecule (TZV2P/BLYP DFT) and 48 methanol molecules 
>>>>> (OPLS MM force field).
>>>>> The relevant parts of my input file are attached at the end of the 
>>>>> message. Basically, I'm following all the advices from this group for a 
>>>>> fully periodic QM/MM setup and set the cell sizes for the QM (in QMMM/CELL) 
>>>>> and MM (in SUBSYS/CELL) parts equal.
>>>>>
>>>>> The following is the summary of my findings.
>>>>>
>>>>> 1. With FORCE_EVAL/QMMM/CENTER EVERY_STEP, the total energy is not 
>>>>> conserved, as pointed out by Dr. Laino in this group before.
>>>>> 2. With FORCE_EVAL/QMMM/CENTER SETUP_ONLY, the total energy is very 
>>>>> well conserved within 0.1 kcal/mol for a few ps until the QM molecule 
>>>>> eventually crosses the cell boundary.  At this point system energy shoots 
>>>>> up by several kcal/mol and the output begins to produce the following 
>>>>> messages:
>>>>>
>>>>>  *** 16:11:03 WARNING in qmmm_util:apply_qmmm_walls_reflective :: One 
>>>>> or  ***
>>>>>  *** few QM atoms are within the SKIN of the quantum box. Check your 
>>>>> run  ***
>>>>>  *** and you may possibly consider: the activation of the QMMM 
>>>>> WALLS      ***
>>>>>  *** around the QM box, switching ON the centering of the QM box 
>>>>> or       ***
>>>>>  *** increase the size of the QM cell. CP2K CONTINUE but results could 
>>>>> be ***
>>>>>  *** meaningless. qmmm_util.F line 
>>>>> 206                                    ***
>>>>>
>>>>> My questions are as follows:
>>>>>
>>>>> 1. Is it possible to achieve the energy conservation in a fully 
>>>>> periodic DFT/MM MD without worrying about the QM molecule crossing the cell 
>>>>> boundary? I believe this is possible in full DFT MD but my experience on 
>>>>> CP2K QM/MM so far indicates otherwise. I'd appreciate if someone can 
>>>>> provide an answer and keywords to achieve it.
>>>>>
>>>>> 2. The QMMM/CENTER option affects energy conservation a lot. When 
>>>>> CENTER EVERY_STEP is used, would it improve energy conservation if I use 
>>>>> finer grids by increasing DFT/MGRID/CUTOFF value?
>>>>>
>>>>> The skeletal form of my input file is attached below.
>>>>> Thank you.
>>>>>
>>>>> Jonggu Jeon
>>>>>
>>>>>
>>>>>
>>>>>  &MOTION
>>>>>    &MD
>>>>>      ENSEMBLE  NVE
>>>>>      STEPS  40000
>>>>>      TIMESTEP     0.5
>>>>>      STEP_START_VAL  0
>>>>>      TIME_START_VAL    0
>>>>>    &END MD
>>>>>  &END MOTION
>>>>>  &FORCE_EVAL
>>>>>    METHOD  QMMM
>>>>>    &DFT
>>>>>      BASIS_SET_FILE_NAME /opt/cp2k/2.5.1/tests/QS/GTH_BASIS_SETS
>>>>>      POTENTIAL_FILE_NAME /opt/cp2k/2.5.1/tests/QS/POTENTIAL
>>>>>      WFN_RESTART_FILE_NAME wfn.rst
>>>>>      CHARGE  0
>>>>>      &SCF
>>>>>        MAX_SCF  30
>>>>>        EPS_SCF     9.9999999999999995E-07
>>>>>        SCF_GUESS  RESTART
>>>>>        &OT  T
>>>>>          MINIMIZER  DIIS
>>>>>          PRECONDITIONER  FULL_ALL
>>>>>          ENERGY_GAP     1.0000000000000000E-03
>>>>>        &END OT
>>>>>        &OUTER_SCF  T
>>>>>          EPS_SCF     9.9999999999999995E-07
>>>>>          MAX_SCF  1
>>>>>        &END OUTER_SCF
>>>>>      &END SCF
>>>>>      &QS
>>>>>        EPS_DEFAULT     1.0000000000000000E-10
>>>>>        MAP_CONSISTENT  T
>>>>>        EXTRAPOLATION  ASPC
>>>>>        EXTRAPOLATION_ORDER  3
>>>>>        METHOD  GPW
>>>>>      &END QS
>>>>>      &MGRID
>>>>>        CUTOFF     2.8000000000000000E+02
>>>>>        COMMENSURATE  T
>>>>>      &END MGRID
>>>>>      &XC
>>>>>        &XC_GRID
>>>>>          XC_SMOOTH_RHO  NN50
>>>>>          XC_DERIV  SPLINE2_SMOOTH
>>>>>        &END XC_GRID
>>>>>        &XC_FUNCTIONAL  NO_SHORTCUT
>>>>>          &BECKE88  T
>>>>>          &END BECKE88
>>>>>          &LYP  T
>>>>>          &END LYP
>>>>>        &END XC_FUNCTIONAL
>>>>>        &VDW_POTENTIAL
>>>>>          POTENTIAL_TYPE  PAIR_POTENTIAL
>>>>>          &PAIR_POTENTIAL
>>>>>            TYPE  DFTD3
>>>>>            PARAMETER_FILE_NAME /opt/cp2k/2.5.1/tests/QS/dftd3.dat
>>>>>            REFERENCE_FUNCTIONAL BLYP
>>>>>            CALCULATE_C9_TERM  F
>>>>>            REFERENCE_C9_TERM  F
>>>>>            LONG_RANGE_CORRECTION  F
>>>>>          &END PAIR_POTENTIAL
>>>>>        &END VDW_POTENTIAL
>>>>>      &END XC
>>>>>    &END DFT
>>>>>    &MM
>>>>>      &FORCEFIELD
>>>>>        PARMTYPE  AMBER
>>>>>        PARM_FILE_NAME ndmd48.ao.prmtop
>>>>>        &SPLINE
>>>>>          EMAX_SPLINE     1.0000000000000000E+00
>>>>>        &END SPLINE
>>>>>      &END FORCEFIELD
>>>>>      &POISSON
>>>>>        &EWALD
>>>>>          EWALD_TYPE  SPME
>>>>>          ALPHA     3.4999999999999998E-01
>>>>>          GMAX  18
>>>>>        &END EWALD
>>>>>      &END POISSON
>>>>>    &END MM
>>>>>    &QMMM
>>>>>      E_COUPL  GAUSS
>>>>>      USE_GEEP_LIB  7
>>>>>      NOCOMPATIBILITY  T
>>>>>      CENTER  SETUP_ONLY
>>>>>      INITIAL_TRANSLATION_VECTOR    -1.0888775605157457E+01   
>>>>> -2.9066675351256084E+00    7.3922596274008097E+00
>>>>>      &CELL
>>>>>        ABC     1.5084868700000001E+01    1.5084868700000001E+01    
>>>>> 1.5084868700000001E+01
>>>>>        PERIODIC  XYZ
>>>>>      &END CELL
>>>>>      &PERIODIC
>>>>>        GMAX     5.0000000000000000E-01
>>>>>        &MULTIPOLE
>>>>>          RCUT     1.2000000000000000E+01
>>>>>          EWALD_PRECISION     4.9999999999999998E-07
>>>>>          ANALYTICAL_GTERM  T
>>>>>          NGRIDS  50 50 50
>>>>>        &END MULTIPOLE
>>>>>      &END PERIODIC
>>>>>    &END QMMM
>>>>>    &SUBSYS
>>>>>      &CELL
>>>>>        A     1.5084868700000001E+01    0.0000000000000000E+00    
>>>>> 0.0000000000000000E+00
>>>>>        B     0.0000000000000000E+00    1.5084868700000001E+01    
>>>>> 0.0000000000000000E+00
>>>>>        C     0.0000000000000000E+00    0.0000000000000000E+00    
>>>>> 1.5084868700000001E+01
>>>>>        MULTIPLE_UNIT_CELL  1 1 1
>>>>>      &END CELL
>>>>>      &TOPOLOGY
>>>>>        NUMBER_OF_ATOMS  156
>>>>>        CONN_FILE_NAME ndmd48.ao.prmtop
>>>>>        CONN_FILE_FORMAT  AMBER
>>>>>        MULTIPLE_UNIT_CELL  1 1 1
>>>>>      &END TOPOLOGY
>>>>>    &END SUBSYS
>>>>>  &END FORCE_EVAL
>>>>>
>>>>>  -- 
>>>>> You received this message because you are subscribed to the Google 
>>>>> Groups "cp2k" group.
>>>>> To unsubscribe from this group and stop receiving emails from it, send 
>>>>> an email to cp2k+... at googlegroups.com.
>>>>> To post to this group, send email to cp... at googlegroups.com.
>>>>> Visit this group at http://groups.google.com/group/cp2k.
>>>>> For more options, visit https://groups.google.com/d/optout.
>>>>>
>>>>
>>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20150319/d2adb623/attachment.htm>


More information about the CP2K-user mailing list