[CP2K-user] [CP2K:13720] QM/MM energy conservation

[Thomas Kühne]

Dear Mayank, 

the y-axis of your conserved quantity plot is too coarse to judge if energy conservation is an issue at all. 
Possibly it is already good enough, but in any case the integration time step error is not s.th. you should 
be concerned first and foremost, hence you can safely increase TIMESTEP. Also, since energy is an extensive 
property, in my opinion the best way to look at the conserved quantity is in units of Kelvin per atom and ps. 
The potential energy is a matter of your initial configuration and only shows that you haven’t equilibrated 
as yet … 

Cheers, 
Thomas 

> Am 04.08.2020 um 08:52 schrieb may... at gmail.com <mayan... at gmail.com>:
> 
> Hi,
> 
> I am trying a NVT MD simulation for a QM/MM system with QM water molecules surrounding a spatially fixed MM molecule, which are in turn surrounded by MM water molecules. I am using a time step of 5 a.u. (for particular reasons) and the QM box is much smaller than compared to full MM box so within the timescale of MD simulation the QM molecules do not cross the QM box boundary. I am attaching the relevant sections of input file:
> ----
>  &MOTION
>    &MD
>      ENSEMBLE  NVT
>      STEPS  10000
>      TIMESTEP     1.21E-01
>      TEMPERATURE     3.00E+02
>      &THERMOSTAT
>        TYPE  NOSE
>        REGION  DEFINED
>        &DEFINE_REGION
>          QM_SUBSYS  ATOMIC
>        &END DEFINE_REGION
>        &DEFINE_REGION
>          MM_SUBSYS  ATOMIC
>        &END DEFINE_REGION
>        &NOSE
>          TIMECON     9.99E+02
>        &END NOSE
>      &END THERMOSTAT
> 
>  &FORCE_EVAL
>    METHOD  QMMM
>    &DFT
>      BASIS_SET_FILE_NAME ./GTH_BASIS_SETS
>      POTENTIAL_FILE_NAME ./GTH_POTENTIALS
>      CHARGE  0
>      &SCF
>        MAX_SCF  60
>        EPS_SCF     2.0E-07
>        SCF_GUESS  RESTART
>        &OT  T
>          MINIMIZER  DIIS
>          PRECONDITIONER  FULL_SINGLE_INVERSE
>        &END OT
>        &OUTER_SCF  T
>          EPS_SCF     2.0E-07
>          MAX_SCF  40
>        &END OUTER_SCF
>        &PRINT
>          &RESTART  SILENT
>            BACKUP_COPIES  1
>          &END RESTART
>        &END PRINT
>      &END SCF
>      &QS
>        EXTRAPOLATION  ASPC
>        EXTRAPOLATION_ORDER  4
>      &END QS
>      &MGRID
>        CUTOFF     3.2E+02
>        COMMENSURATE  T
>      &END MGRID
>      &XC
>        DENSITY_CUTOFF     1.0E-10
>        GRADIENT_CUTOFF     1.0E-10
>        TAU_CUTOFF     1.0E-10
>        &XC_GRID
>          XC_SMOOTH_RHO  NN50
>          XC_DERIV  SPLINE2_SMOOTH
>        &END XC_GRID
>        &XC_FUNCTIONAL  NO_SHORTCUT
>          &PBE  T
>            PARAMETRIZATION  REVPBE
>          &END PBE
>        &END XC_FUNCTIONAL
>        &VDW_POTENTIAL
>          POTENTIAL_TYPE  PAIR_POTENTIAL
>          &PAIR_POTENTIAL
>            R_CUTOFF     1.9E+01
>            TYPE  DFTD3
>            PARAMETER_FILE_NAME ./dftd3.dat
>            REFERENCE_FUNCTIONAL revPBE
>            CALCULATE_C9_TERM  T
>            REFERENCE_C9_TERM  T
>            LONG_RANGE_CORRECTION  F
>          &END PAIR_POTENTIAL
>        &END VDW_POTENTIAL
>      &END XC
>    &END DFT
>    &MM
>      &FORCEFIELD
>        ...
>      &END FORCEFIELD
>      &POISSON
>        &EWALD
>          EWALD_TYPE  SPME
>          RCUT     1.2E+01
>          ALPHA     2.917E-01
>          GMAX  75
>          O_SPLINE  4
>        &END EWALD
>      &END POISSON
>    &END MM
>    &QMMM
>      E_COUPL  GAUSS
>      USE_GEEP_LIB  6
>      NOCOMPATIBILITY  T
>      CENTER  NEVER
>      INITIAL_TRANSLATION_VECTOR     0.0E+00    0.0E+00    0.0E+00
>      &QM_KIND O
>        MM_INDEX  ...
>      &END QM_KIND
>      &QM_KIND H
>        MM_INDEX  ...
>      &END QM_KIND
>      &MM_KIND H
>        RADIUS     4.4E-01
>      &END MM_KIND
>      &MM_KIND O
>        RADIUS     1.2E+00
>      &END MM_KIND
>      .....
>      &CELL
>        ABC     4.0E+01    3.0E+01    3.0E+01
>        PERIODIC  XYZ
>      &END CELL
>      &PERIODIC
>        GMAX     5.0E-01
>        &MULTIPOLE  ON
>          RCUT     1.2E+01
>          EWALD_PRECISION     1.0E-08
>          ANALYTICAL_GTERM  T
>        &END MULTIPOLE
>      &END PERIODIC
>    &END QMMM
>    &SUBSYS
>      &CELL
>        A     5.0E+01    0.0E+00    0.0E+00
>        B     0.0E+00    4.0E+01    0.0E+00      
>        C     0.0E+00    0.0E+00    4.0E+01
>        MULTIPLE_UNIT_CELL  1 1 1
>      &END CELL
>      &COORD
> 
> -------
> 
> With this setup for the production run I obtained the following trends for the potential and constant quantity as the graphs attached here. I have tested the QM setup for a fully QM calculation on bulk water, which works fine. My main aim here is to reduce the energy drift for the QM/MM system, so is there any other aspect of the calculation should I take care to reduce the drift? Any suggested corrections for the QM/MM file would also be helpful.
> 
> Best Regards,
> Mayank Dodia
> 
> 
> 
> 
> -- 
> 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 cp... at googlegroups.com <mailto:cp... at googlegroups.com>.
> To view this discussion on the web visit https://groups.google.com/d/msgid/cp2k/63b0f162-c684-481d-884a-f10890c1aa91n%40googlegroups.com <https://groups.google.com/d/msgid/cp2k/63b0f162-c684-481d-884a-f10890c1aa91n%40googlegroups.com?utm_medium=email&utm_source=footer>.
> <const_quantity.png><pot_enr.png>



==============================
Thomas D. Kühne
Dynamics of Condensed Matter
Chair of Theoretical Chemistry
University of Paderborn
Warburger Str. 100
D-33098 Paderborn
Germany
tdku... at mail.upb.de
+49/(0)5251/60-5726

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20200805/d9746d20/attachment.htm>