Discrepancy in excited-state frequency calculations

Questions
#1

I have run an optimization and frequency calculation for the first excited singlet state using LRC-ωPBE/6-31G. The input is given below.

The Q-Chem calculation gives me 8 imaginary frequencies. However, when I take the optimized structure and perform a frequency calculation in another quantum chemistry software, I obtain 0 imaginary frequencies. I get the same result when using a numerical Hessian in Q-Chem (ideriv 1).

I am wondering what could be causing this discrepancy. Numerical frequencies are quite expensive, and so I would like to avoid computing these.

$rem
mem_total 200000
mem_static 1000
method lrc-wpbe
basis 6-31g
xc_grid 3
omega 211
max_scf_cycles 300
geom_opt_max_cycles 300
thresh 14
scf_convergence 8
jobtype opt
cis_state_deriv 1
cis_n_roots 3
cis_singlets true
cis_triplets false
$end

$molecule
0 1
C 1.5502275 -3.5534072 -0.3905564
C 0.4008405 -2.8578959 -0.0822371
C 0.4381229 -1.4723848 0.1412005
C 1.6518617 -0.7129762 -0.0246236
C 2.8204191 -1.4836452 -0.2934103
C 2.7734491 -2.8440858 -0.4733212
C 1.6514178 0.7138775 0.0246503
C 2.8195400 1.4852411 0.2933414
C 2.7717947 2.8456607 0.4731899
C 1.5481611 3.5542723 0.3904408
C 0.3991599 2.8580848 0.0822196
C 0.4372076 1.4725787 -0.1411327
C -0.7503840 0.7955727 -0.6347162
C -0.7497891 -0.7960053 0.6348522
O -0.7043266 -0.1960794 1.8560915
O -0.7045740 0.1954405 -1.8558397
C -2.0520240 0.0110628 2.2986463
C -2.8074974 -1.0969129 1.5725498
C -2.0521487 -0.0123417 -2.2984705
C -2.8081213 1.0955018 -1.5726961
O -1.9732271 1.3783935 -0.4345569
O -1.9723918 -1.3792279 0.4344239
H 1.5204792 -4.6311877 -0.5475068
H -0.5496587 -3.3819153 0.0290049
H 3.7730825 -0.9685036 -0.4199165
H 3.6936010 -3.3837674 -0.7047991
H 3.7725124 0.9706610 0.4198076
H 3.6916479 3.3858854 0.7045878
H 1.5178012 4.6320454 0.5473248
H -0.5516428 3.3815583 -0.0290048
H -2.3854701 1.0150005 2.0057232
H -2.0640233 -0.0892468 3.3891499
H -3.7996308 -0.7940782 1.2181387
H -2.8927774 -2.0154155 2.1710489
H -2.3852249 -1.0163451 -2.0053505
H -2.0640915 0.0877180 -3.3889976
H -3.8001694 0.7923596 -1.2183090
H -2.8937044 2.0138342 -2.1714134
$end

@@@

$molecule
read
$end

$rem
mem_total 200000
mem_static 1000
method lrc-wpbe
basis 6-31g
xc_grid 3
omega 211
max_scf_cycles 200
geom_opt_max_cycles 200
thresh 14
scf_convergence 8
jobtype freq
cis_state_deriv 1
cis_n_roots 3
cis_singlets true
cis_triplets false
$end

#2

There have been some reports of this kind of issue. The usual best practice when encountering imaginary frequencies is to tight the geometry optimization stopping criteria, e.g.

GEOM_OPT_TOL_DISPLACEMENT 600
GEOM_OPT_TOL_GRADIENT 150
GEOM_OPT_TOL_ENERGY 50

tightens each by a factor of 2, relative to the defaults.

I would also point out that numerical frequencies have been MPI-parallelized since Q-Chem 6.3 so they can be very low cost on sufficient hardware, since the steps are trivially parallelizable.

#3

Thanks for your answer!

I reran the optimization with the tightened criteria, but the frequency calculation still produces 8 imaginary frequencies.

The fact that both the numerical Hessian in Q-Chem and the frequency calculation in another program give 0 imaginary frequencies makes us suspect that there may be an issue specifically with the analytical Hessian for the functional in Q-Chem. Is this reasonable to suspect?

#4

we’ve had some other reports of that. so far not with any system small enough to make it a convenient debugging test. I looked into this a while back but was unable to diagnose the problem.

#5

As a further update, we did some small testing. We tried:

  • Different basis sets (6-31G and 6-311G(2df,p))
  • Both GS and excited state calculations
  • With and without the omega keyword
  • Different functionals (B3LYP and ωB97X-D)

The only case where we obtained 0 imaginary frequencies (consistent with the other software) was when we changed the functional.

#6

which functional gives you 0 imag freqs?

#7

The functionals that give 0 imaginary frequencies are b3lyp and ωB97x-D.

To be more clear on what we tested and the results:

All geometry optimizations were performed in Q-Chem. For each optimized structure, the vibrational frequencies were subsequently calculated using the analytical Hessian in Q-Chem and another quantum-chemistry software. The numerical Hessian in Q-Chem was computed for one of the structures.

Q-Chem