How to troubleshoot artifical negative frequencies

How To
#1

I am doing rather routine calculation: TD-DFT, S1 state, optimization and frequencies. Everything converges smoothly, but then I end up with one imaginary frequency, which I am sure is an artefact. I tried several things at random with no effect. Question to the experts: What are standard steps of troubleshooting this sort of a problem? What are the first thing to tighten in terms of thresholds/convergence/grids?

Here is my basic input, using defaults for everything:

Blockquote
$comment
OPT& Frequency calculation for
S1 structure optimized with
method = wB97X-D
basis = 6-311+G**
Also, use old optimzier
Start with slightly distorted structure to see if there is lower-energy
minimum compared to D2 struture
$end

$molecule
0 1
C 1.4029162007 1.3106017388 0.0000590696
C 0.7000528271 1.4568212458 1.2118385596
H 1.2444817885 1.4462938585 2.1496228374
C -0.7002534075 1.4565647257 1.2117792977
H -1.2447241991 1.4459115552 2.1495365571
C -1.4029162007 1.3106017388 -0.0000590696
C -0.7000528271 1.4568212458 -1.2118385596
H -1.2444817885 1.4462938585 -2.1496228374
C 0.7002534075 1.4565647257 -1.2117792977
H 1.2447241991 1.4459115552 -2.1495365571
C 2.8207494939 0.7871832610 0.0003491932
H 3.3544685496 1.1572220091 0.8803024505
H 3.3551773248 1.1579436494 -0.8789103039
C -2.8207494939 0.7871832610 -0.0003491932
H -3.35 1.1572220091 -0.8803024505
H -3.36 1.1579436494 0.8789103039
C 2.8207494939 -0.7871832610 -0.0003491932
H 3.3544685496 -1.1572220091 -0.8803024505
H 3.3551773248 -1.1579436494 0.8789103039
C 1.4029162007 -1.3106017388 -0.0000590696
C 0.7000528271 -1.4568212458 -1.2118385596
H 1.2444817885 -1.4462938585 -2.1496228374
C -0.7002534075 -1.4565647257 -1.2117792977
H -1.2447241991 -1.4459115552 -2.1495365571
C -1.4029162007 -1.3106017388 0.0000590696
C -0.7000528271 -1.4568212458 1.2118385596
H -1.2444817885 -1.4462938585 2.1496228374
C 0.7002534075 -1.4565647257 1.2117792977
H 1.2447241991 -1.4459115552 2.1495365571
C -2.8207494939 -0.7871832610 0.0003491932
H -3.3544685496 -1.1572220091 0.8803024505
H -3.3551773248 -1.1579436494 -0.8789103039
$end

$rem
jobtype = opt
geom_opt_driver optimize !use old algorithm
method = wB97X-D
basis = 6-311+G**
cis_n_roots = 1
cis_state_deriv = 1
cis_triplets=false
cis_singlets=true
mem_total = 10000
!XC_GRID = 3 !use SG3 grid.
$end

@@@
$molecule
read
$end

$rem
jobtype = freq
method = wB97X-D
basis = 6-311+G**
!XC_GRID = 3 !use SG3 grid
scf_guess = read
cis_n_roots = 1
cis_state_deriv = 1
cis_triplets=false
cis_singlets=true
mem_total = 10000
$end

Facing problem with Imaginary frequency
#2 
GEOM_OPT_TOL_GRADIENT 300
GEOM_OPT_TOL_DISPLACEMENT 1200
GEOM_OPT_TOL_ENERGY 100

Those are the default values, I would reduce by factor of 2 in the situation you describe. Maybe also check DFT grid. SG-2 is the default for that functional, SG-3 is a better grid and if you need to increase from there, use XG_GRID = 000099000590 for EML (99,590).

#3

You can also re-run the optimization from your lowest energy structure and start with an exact initial Hessian.

Also using the new optimizer has a more rigorous topology generation which can better describe the system leading to smoother optimization.

In addition, as @jherbert mention convergence criteria can be changed, but with the old optimizer you are limited to integers.

#4

Thanks for the advice! Here is the report:

I tried both optimizers – no difference. I tried to reduce thresholds by a factor of 2 - it did not help. Then I tried these tight thresholds with SG-3 grid-- again, no difference, the negative frequency actually increased (62 instead of original ~40). I think I will submit this job as a ticket.
And I will try another functional to see if this is a problem with wB97X-D.

#5

As a last resort you can always perturb along the imaginary mode. No automated way to do that (used to be a script floating around but the script was wrong), have to do it by hand using the normal mode eigenvectors printed by the frequency job.

#6

I hope this reply is not too late to help with your problem but I will include it here for prosperity sake.

I find a “nearly” fool proof method for resolving imaginary frequency’s is to first run a frequency calculation and then read in the more accurate hessian into an optimization calculation followed by another frequency calculation. This should hopefully remove your imaginary frequency’s.

Although, this method is quite computationally and time intensive because it requires two frequency calculations, it works more often then not for resolving this problem. The hessian generated in the first freq calculation is always going to be better then the initial guess. An example:

$rem
   BASIS  =  LACVP
   ECP  =  fit-LACVP
   GUI  =  2
   JOB_TYPE  =  Frequency
   METHOD  =  wB97MV
   SCF_CONVERGENCE  =  8
   SCF_MAX_CYCLES  =  100
$end

@@@

$molecule
read
$end

$rem
   BASIS  =  LACVP
   ECP  =  fit-LACVP
   GEOM_OPT_HESSIAN  =  READ
   GEOM_OPT_MAX_CYCLES  =  1000
   GUI  =  2
   JOB_TYPE  =  Optimization
   METHOD  =  wB97MV
   SCF_CONVERGENCE  =  8
   SCF_MAX_CYCLES  =  100
$end

@@@

$molecule
read
$end

$rem
   BASIS  =  LACVP
   ECP  =  fit-LACVP
   GUI  =  2
   JOB_TYPE  =  Frequency
   METHOD  =  wB97MV
   SCF_CONVERGENCE  =  8
   SCF_MAX_CYCLES  =  100
$end