12
$\begingroup$

I have been having trouble with the VIBRATIONAL_ANALYSIS routine. Namely, the routine keeps timing out. I am trying to find the vibrational modes of Pyridine in a non periodic system, as well as $\ce{Mo2C}$ in a periodic system. I'm running on 16 nodes in parallel, but it still takes a long time. I have tried making the box smaller for pyridine, but that changes the vibrational modes significantly. I have included the full input files at the end.

For Pyridine parameters I have (so far):

EPS_DEFAULT 1E-10
CUTOFF 800
NGRIDS 6
REL_CUTOFF 40
MAX_SCF 100
EPS_SCF 1.0E-07
POISSON_SOLVER  WAVELET
PERIODIC NONE 
ABC 30.0 30.0 30.0
TYPE MINIMIZATION
MAX_DR    1.0E-03
MAX_FORCE 1.0E-03
RMS_DR    1.0E-03
RMS_FORCE 1.0E-03
DX 1E-3

For Mo2C parameters I have (so far):

EPS_DEFAULT 1E-10
CUTOFF 800
NGRIDS 5
REL_CUTOFF 70
MAX_SCF 200
EPS_SCF 1.0E-07
ABC 4.75411  5.24148  6.07649
ALPHA_BETA_GAMMA 90 90 90
SYMMETRY ORTHORHOMBIC
PERIODIC XYZ
(GEO_OPT SECTION)
TYPE MINIMIZATION
MAX_DR    1.0E-03
MAX_FORCE 1.0E-03
RMS_DR    1.0E-03
RMS_FORCE 1.0E-03
DX 1E-3
(CELL_OPT SECTION)
TYPE GEO_OPT
MAX_DR    1.0E-03
MAX_FORCE 1.0E-03
RMS_DR    1.0E-03
RMS_FORCE 1.0E-03
DX 1E-3

Input file for Pyridine:

&FORCE_EVAL
  METHOD QUICKSTEP
  &DFT
    BASIS_SET_FILE_NAME ./BASIS_SET
    POTENTIAL_FILE_NAME ./POTENTIALS
    &QS
       EPS_DEFAULT 1E-10
    &END QS
    &MGRID
      CUTOFF 800
      NGRIDS 6
      REL_CUTOFF 40
    &END MGRID
    &SCF
      MAX_SCF 100
      EPS_SCF 1.0E-07
      SCF_GUESS RESTART
      &OT
         MINIMIZER DIIS
      &END OT
    &END SCF
    &XC
      &XC_FUNCTIONAL PBE
      &END XC_FUNCTIONAL
    &END XC
    &POISSON
       POISSON_SOLVER  WAVELET
       PERIODIC NONE 
       &WAVELET
         SCF_TYPE  100
       &END WAVELET
    &END POISSON
  &END DFT
  &SUBSYS
    &CELL
      ABC 30.0 30.0 30.0
      PERIODIC NONE
    &END CELL
    &TOPOLOGY
      COORD_FILE_FORMAT XYZ
      COORD_FILE_NAME pyridine.xyz
      &CENTER_COORDINATES
      &END CENTER_COORDINATES
    &END TOPOLOGY
    &KIND H
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q1
    &END KIND
    &KIND C 
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q4
    &END KIND
    &KIND N 
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q5
    &END KIND
  &END SUBSYS
&END FORCE_EVAL
&GLOBAL
  PROJECT Pyridine
  RUN_TYPE VIBRATIONAL_ANALYSIS
  PRINT_LEVEL LOW
  EXTENDED_FFT_LENGTHS
&END GLOBAL

&MOTION
  &GEO_OPT
! to calculate vibrational spectrum tight convergence is required because frequencies are very sensitive to force constant
    TYPE MINIMIZATION
    MAX_DR    1.0E-03
    MAX_FORCE 1.0E-03
    RMS_DR    1.0E-03
    RMS_FORCE 1.0E-03
    MAX_ITER 200
    OPTIMIZER CG
  &END
&END
! setup parameters to perform a Normal Modes analysis
&VIBRATIONAL_ANALYSIS
!   &MODE_SELECTIVE
!      LOWEST_FREQUENCY 0.0
!      FREQUENCY   630
!      ATOMS 6 1
!      INITIAL_GUESS RESTART
!      RESTART_FILE_NAME Pyridine-VIBRATIONS-1.mol 
!      EPS_NORM  1.0E-4
!      EPS_MAX_VAL 1.0E-5
!   &END
!Specify the number of processors to be used per replica environment (for parallel runs)
NPROC_REP 4
!Specify the increment to be used to construct the HESSIAN with finite difference method
   DX 1E-3
   &PRINT
     &PROGRAM_RUN_INFO ON
     &END
   &END
&END

Input file for Mo2C:

&FORCE_EVAL
  METHOD QUICKSTEP
  &DFT
    &QS
      EPS_DEFAULT 1.0E-10
    &END QS
    &POISSON
      PERIODIC XYZ
    &END POISSON
    BASIS_SET_FILE_NAME ./BASIS_SET
    POTENTIAL_FILE_NAME ./POTENTIALS
    &MGRID
      CUTOFF 800
      NGRIDS 5
      REL_CUTOFF 70
    &END MGRID
    &SCF
      MAX_SCF 200
      EPS_SCF 1.0E-07
      SCF_GUESS RESTART
      &MIXING T
        ALPHA 0.5
        METHOD PULAY_MIXING
        NPULAY 5
      &END MIXING
    &END SCF
    &XC
      &XC_FUNCTIONAL PBE
      &END XC_FUNCTIONAL
    &END XC
  &END DFT
  &SUBSYS
    &CELL
      ABC 4.75411  5.24148  6.07649
      ALPHA_BETA_GAMMA 90 90 90
      SYMMETRY ORTHORHOMBIC
      PERIODIC XYZ
    &END CELL
    &TOPOLOGY
      COORD_FILE_FORMAT XYZ
      COORD_FILE_NAME mo2c.xyz
    &END TOPOLOGY
    &KIND Mo
      BASIS_SET TZVP-MOLOPT-SR-GTH
      POTENTIAL GTH-PBE-q14 
    &END KIND
    &KIND C 
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q4
    &END KIND
  &END SUBSYS
&END FORCE_EVAL
&GLOBAL
  PROJECT MoC2
  RUN_TYPE GEO_OPT
  PRINT_LEVEL LOW
  EXTENDED_FFT_LENGTHS
&END GLOBAL

&MOTION
  &GEO_OPT
! to calculate vibrational spectrum tight convergence is required because frequencies are very sensitive to force constant
    TYPE MINIMIZATION
    MAX_DR    1.0E-03
    MAX_FORCE 1.0E-03
    RMS_DR    1.0E-03
    RMS_FORCE 1.0E-03
    MAX_ITER 200
    OPTIMIZER CG
  &END
  &CELL_OPT
! to calculate vibrational spectrum tight convergence is required because frequencies are very sensitive to force constant
    MAX_DR    1.0E-03
    MAX_FORCE 1.0E-03
    RMS_DR    1.0E-03
    RMS_FORCE 1.0E-03
    MAX_ITER 200
    OPTIMIZER CG
    TYPE GEO_OPT
  &END
&END
! setup parameters to perform a Normal Modes analysis

!&VIBRATIONAL_ANALYSIS
!Specify the number of processors to be used per replica environment (for parallel runs)
!NPROC_REP 4
!Specify the increment to be used to construct the HESSIAN with finite difference method
!DX 1E-08
!   &PRINT
!     &PROGRAM_RUN_INFO ON
!     &END
!   &END
!&END
```
$\endgroup$
4
  • $\begingroup$ +1. I've just added a tag which I think you might find helps when asking about optimizing computation performance. $\endgroup$ – Nike Dattani Nov 18 '20 at 2:32
  • 1
    $\begingroup$ Without seeing the output from cp2k its hard to tell. But I assume your SCF does only converge slowly. There are several things you can try: EPS_DEFAULT 1E-14, MAX_SCF 30 and make use of the &OUTER_SCF section. $\endgroup$ – Fabian Nov 18 '20 at 7:13
  • $\begingroup$ Thank you, I will try that. I am starting to suspect the same. I am curious about cell size as well, specifically how one may be able to have a good starting guess for the cell size in a non-periodic system. I would think that if the cell was non-periodic, the results wouldn't change drastically, but that does not seem to be the case. $\endgroup$ – Jason M Gray Nov 18 '20 at 18:19
  • $\begingroup$ The density should be zero at the edge of the cell. $\endgroup$ – Fabian Nov 18 '20 at 18:48
4
$\begingroup$

After a bit of hit and trial, I am able to run the pyridine sample, using CP2K v7.0 (git:d15a9c53f), however I think this can be run on any CP2K version beyond v2.5. The modified file is:

  &FORCE_EVAL
  METHOD QUICKSTEP
  &DFT
    BASIS_SET_FILE_NAME ./BASIS_SET
    POTENTIAL_FILE_NAME ./GTH_POTENTIALS
    &QS
       EPS_DEFAULT 1E-10
    &END QS
    &MGRID
      CUTOFF 400
    &END MGRID
    &SCF
      MAX_SCF  20
      EPS_SCF     1E-08
      SCF_GUESS  RESTART
      &OT  T
        MINIMIZER  DIIS
        PRECONDITIONER  FULL_KINETIC
      &END OT
      &OUTER_SCF  T
        EPS_SCF     1E-08
        MAX_SCF  20
      &END OUTER_SCF
      &PRINT
        &RESTART  SILENT
          BACKUP_COPIES  1
        &END RESTART
      &END PRINT
    &END SCF
    &XC
      &XC_FUNCTIONAL PBE
      &END XC_FUNCTIONAL
    &END XC
    &POISSON
       POISSON_SOLVER  WAVELET
       PERIODIC NONE 
       &WAVELET
         SCF_TYPE  100
       &END WAVELET
    &END POISSON
  &END DFT
  &SUBSYS
    &CELL
      ABC 30.0 30.0 30.0
      PERIODIC NONE
    &END CELL
    &TOPOLOGY
      COORD_FILE_FORMAT XYZ
      COORD_FILE_NAME pyridine.xyz
      &CENTER_COORDINATES
      &END CENTER_COORDINATES
    &END TOPOLOGY
    &KIND H
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q1
    &END KIND
    &KIND C 
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q4
    &END KIND
    &KIND N 
      BASIS_SET DZVP-GTH-PBE
      POTENTIAL GTH-PBE-q5
    &END KIND
  &END SUBSYS
&END FORCE_EVAL
&GLOBAL
  PROJECT Pyridine
  RUN_TYPE VIBRATIONAL_ANALYSIS
  PRINT_LEVEL LOW
  EXTENDED_FFT_LENGTHS
&END GLOBAL

&MOTION
  &GEO_OPT
! to calculate vibrational spectrum tight convergence is required because frequencies are very sensitive to force constant
    TYPE MINIMIZATION
    MAX_DR    1.0E-03
    MAX_FORCE 1.0E-03
    RMS_DR    1.0E-03
    RMS_FORCE 1.0E-03
    MAX_ITER 200
    OPTIMIZER CG
  &END
&END
! setup parameters to perform a Normal Modes analysis
&VIBRATIONAL_ANALYSIS
!   &MODE_SELECTIVE
!      LOWEST_FREQUENCY 0.0
!      FREQUENCY   630
!      ATOMS 6 1
!      INITIAL_GUESS RESTART
!      RESTART_FILE_NAME Pyridine-VIBRATIONS-1.mol 
!      EPS_NORM  1.0E-4
!      EPS_MAX_VAL 1.0E-5
!   &END
!Specify the number of processors to be used per replica environment (for parallel runs)
NPROC_REP 4
!Specify the increment to be used to construct the HESSIAN with finite difference method
   DX 1E-3
   &PRINT
     &PROGRAM_RUN_INFO ON
     &END
   &END
&END

I varied the cubic box length as 10, 20 and 30 angstroms and there was relatively insignificant effect on frequency calculations:

10 Ang^3                     20 Ang^3                 30 Ang^3

     [FREQ] cm-1            [FREQ] cm-1             [FREQ] cm-1
      306.559373             307.614457               308.189936
      368.761441             381.380678               370.081145
      578.498921             591.875865               578.958510
      630.679041             632.222381               631.759367
      652.909718             653.938140               653.110369
      717.683491             723.258235               718.390117
      841.291874             841.755853               841.537638
      883.457892             883.747121               883.733166
      937.440211             937.725276               937.794444
      972.598159             972.540268               972.856600
      975.384418             985.095914               976.166601
     1042.925264            1045.181855              1043.247363
     1061.317180            1061.943046              1061.498007
     1066.390016            1066.772710              1066.463394
     1125.211196            1125.337595              1125.254047
     1176.780112            1176.789866              1176.722581
     1296.556279            1297.011609              1296.543748
     1368.961734            1371.370062              1368.715023
     1448.749298            1449.730137              1448.520992
     1455.387545            1455.803598              1455.480999
     1644.206969            1644.871006              1644.242656
     1645.693998            1646.204083              1645.595586
     3165.883661            3165.929198              3165.913233
     3171.541016            3171.575638              3171.569378
     3182.447441            3182.502989              3182.464905
     3190.435506            3190.499027              3190.471759
     3201.299231            3201.335631              3201.328153

The significant changes here are Orbital Transform (OT) optimizer and 400 Ry plane cutoff to reduce the computational load.

$\endgroup$
4
  • $\begingroup$ Very well done!!! +1. I made an edit, can you check if it is what you intended? $\endgroup$ – Nike Dattani Dec 18 '20 at 4:18
  • $\begingroup$ Wow! thank you so much for your help. This is really good to know for general knowledge in these types of dft calculations. $\endgroup$ – Jason M Gray Jan 15 at 20:35
  • $\begingroup$ Hi @JasonMGray, if this answer solves your issues, can you accept this answer as solution to your question? $\endgroup$ – mykd Jan 15 at 21:19
  • $\begingroup$ Absolutely, thank you. $\endgroup$ – Jason M Gray Jan 19 at 9:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.