[CP2K:6254] Energy conservation in periodic QM/MM
Jonggu Jeon
jeonj... at gmail.com
Thu Mar 19 11:48:29 UTC 2015
Thank you for the suggestion.
In fact, I'm currently running a new job with MGRID/CUTOFF=480 and CENTER
EVERY_STEP. As you can see in the following figure, the early result
indicates a very good energy conservation at the same level as the
previous best example (CUTOFF=280 and CENTER SETUP_ONLY).
What troubles me is the fact that computation time increases rapidly with
larger CUTOFF and the same level of energy conservation is already possible
with smaller cutoff if no centering is used.
However, given the circumstances, I agree that using centering and larger
cutoff is the best way to go.
Regards,
Jonggu
<https://lh6.googleusercontent.com/-P0Moyh0OT5k/VQq0qY--uYI/AAAAAAAAA3w/6xnksES21jA/s1600/te.png>
2015년 3월 19일 목요일 오후 8시 17분 54초 UTC+9, Matt W 님의 말:
>
> 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/65be04c1/attachment.htm>
More information about the CP2K-user
mailing list