[CP2K:6271] Energy conservation in periodic QM/MM
Teodoro Laino
teodor... at gmail.com
Sun Mar 22 08:37:58 UTC 2015
Hi Rolf,
the point is that the QM part is - strictly speaking - not translational invariant for arbitrary translation vectors. This means that if you re-center at every step for an arbitrary translation vector you will see a drift, given by the fact that an atom - depending on its relative position w.r.t. grid points - may experience numerically different forces (with a deviation, hopefully lower than the norm). Moreover, the situation is even more dramatic because the recenter vector will not be random, but rather subject to a constant direction bias, which will cause a pretty visible drift.
Your combination: 300 and CENTER_GRID are both going in the direction of: 1) minimising the error in the atomic forces by increasing cutoff 2)Re-center using a translation vectors which is commensurate to grid spacing. Both should help minimise drift.
CENTER_GRID has no computational overhead.
300 Ryd cutoff instead should marginally increase the cost.
Teo
> On 22 Mar 2015, at 08:23, Rolf David <rolf.d... at gmail.com> wrote:
>
> I've done test on my systems and reach out some conclusions.
>
> GAPW / B3LYP / 280/40 cutoff / CENTER_EVERY_STEP gave me a slow drift and some fluctuations (for 2ps)
>
> GAPW / B3LYP / 300/50 cutoff / CENTER_EVERY_STEP and with CENTER_GRID TRUE, almost no drift/fluctuations (for 2ps)
>
> I don't know if it's the cutoff of the center_grid true which made an ompatch, but the computational cost is negligible (don't see a difference). I'll see next month what of it when i reach 15/20ps.
>
>
>
> On Thursday, March 19, 2015 at 12:48:29 PM UTC+1, Jonggu Jeon wrote:
> 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 <http://groups.google.com/group/cp2k>.
> For more options, visit https://groups.google.com/d/optout <https://groups.google.com/d/optout>.
>
>
> --
> 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+uns... at googlegroups.com <mailto:cp2k+uns... at googlegroups.com>.
> To post to this group, send email to cp... at googlegroups.com <mailto:cp... at googlegroups.com>.
> Visit this group at http://groups.google.com/group/cp2k <http://groups.google.com/group/cp2k>.
> For more options, visit https://groups.google.com/d/optout <https://groups.google.com/d/optout>.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20150322/cebfa693/attachment.htm>
More information about the CP2K-user
mailing list