Quasiclassical MD

Hi everyone,
I have a question about ab initio molecular dynamics. What is the difference between AIMD_INIT_VELOC=ZPE and AIMD_INIT_VELOC=QUASICLASSICAL?

I am planning to perform a quasiclassical direct dynamics trajectory simulation on a transition structure. However, when I use AIMD_INIT_VELOC=QUASICLASSICAL on qchem 6.2, I get an error message stating “Quasiclassical MD not yet supported for transition states.” Conversely, the simulation runs successfully with AIMD_INIT_VELOC=ZPE.


   0  2
      C       3.3135997804    -0.0078810290    -0.0093651548
        H       3.6066735648     0.7959221128    -0.6856021114
        H       3.6363870152     0.1567299549     1.0199545771
          H       3.5760919013    -0.9973732332    -0.3883608468
        C       0.6323632532     0.0166113714     0.0131679721
        C      -0.0441388971     1.2281272323     0.0146353683
     C      -1.4435978962     1.2107274345     0.0042023458
        C      -2.1267743693    -0.0075961871    -0.0077880088
       C      -1.4224439670    -1.2139181757    -0.0095801812
      C      -0.0230120434    -1.2062566659     0.0006442183
       H       2.0065471467     0.0174336366     0.0153133373
       H       0.4981169822     2.1708343243     0.0233488537
       H      -1.9991640785     2.1452594740     0.0050702665
      H      -3.2145949435    -0.0169757986    -0.0163280596
      H      -1.9618313722    -2.1578670589    -0.0193435608
      H       0.5368126337    -2.1388548629    -0.0015363369

   JOBTYPE                     freq
  METHOD                     m06
BASIS              6-31g
  scf_convergence 12



   JOBTYPE                     aimd
   METHOD                     m06
   SCF_CONVERGENCE             6
   TIME_STEP                   20    !  (in atomic units)
   AIMD_STEPS                  1000  !  600 fs total simulation time
   AIMD_TEMP                   12
   AIMD_PRINT                  2
   FOCK_EXTRAP_ORDER           6     ! Use a 6th-order extrapolation
   FOCK_EXTRAP_POINTS          12    ! of the previous 12 Fock matrices
   AIMD_MOMENTS                1
   AIMD_QCT_WHICH_TRAJECTORY   1     ! Loop over several values to get
BASIS              6-31g

Quasiclassical MD uses a Boltzmann distribution that is normal-mode-based, as described here:
The ZPE option uses harmonic oscillator distributions, and the THERMAL option uses a Boltzmann distribution that is particle-based, as described here:

Thanks John. Since quasiclassical MD is not supported for transition states in qchem, how would using ZPE instead of quasiclassical affect the results? Is there a significant difference? I want to see non-IRC (case by combination of frequency modes along with the reaction coordination motion, Acc. Chem. Res. 2021, 54, 4410−4422) intermediates connected to a transition structure.

Looking at the code, it looks like it will treat the negative-curvature mode as a real frequency for the purpose of formulating a harmonic oscillator (Gaussian) probability distribution. Velocities are selected randomly from that distribution.

Thank you. Could you please clarify whether considering a negative mode as a real frequency is acceptable or not? I am not well-versed in its underlying principles, but I would like to confirm if applying ZPE-MD on a transition structure results in motion-vibronic coupling with a correct IRC related motion.

I’m not an expert either, but it must have seemed reasonable to my 20-years-ago self who coded it this way. :slightly_smiling_face: Still seems reasonable although it’s not going to give you motion along the IRC, per se. Consult the literature. I recall that Bill Hase had a lot to say about quasiclassical trajectories starting from the transition state; if there are artifacts of this approach, they are probably addressed in his papers.

Thank you very much John.