Gradient Not Converging-Imaginary Frequency

Questions
#1

Dear Q-Chem Developers,

I am a new user and I am hoping to do some Raman calculations with Q-Chem.

I realized after the geometry optimization, the “Gradient” is not converged (prints No in the output file). However, it seems that the geometry optimization is correctly finished. I have tried changing the thresholds of the three (gradient, displacement, and energy change) but still, I believe the Gradient doesn’t seem to be converged.

When I do frequency calculations for the optimized structure, I get one negative frequency, which I am not sure if it is caused by the fact that the gradient didn’t converge in the optimization step or it’s caused by something else. I have also followed other suggestions regarding the imaginary frequency issue that were suggested in the forum but still no success of getting rid of it.

Any help would be thoroughly appreciated.

Thank you,
Amir

#2

Your first observation (geom opt converged w/o gradient converged) is odd. Under most circumstances, the optimization stopping criteria are that the gradient tolerance must be met and either the energy or the step size criteria must be met (i.e., 2 out of 3 but the gradient criterion is required). This is precisely to avoid situations like an imaginary frequency caused by a non-vanishing gradient. I wonder: did your optimization perhaps just run out of steps? (The max is 50 by default but can be increased.)

There are other things to check for imaginary frequencies (quality of DFT grid being a big one), but probably better to address first things first. Can you post an input file?

#3

Thank you for your reply.

I tried to attach the input files, but for some reason I wasn’t able to. Below are two of my tries. The first try was with two different steps; first, optimization, and then frequency calculation. The second try, I combined both of them into one file.

##First-try-optimization##
$molecule
0 4
U -0.0116475938 -0.0655156332 0.0000000000
Cl -2.4452759521 -1.0138063243 0.0000000000
Cl 0.2384690355 2.5216700770 0.0000000000
Cl 2.2380787565 -1.3715562104 0.0000000000
$end

$rem
jobtype = optimization
method = b3lyp
basis = general
ecp = general
PURECART 111
!NBO = true
ecp_fit = true
unrestricted = true
symmetry = false
!scf_algorithm = diis_gdm
!scf_convergence = 5
!scf_guess = gwh
!scf_guess_mix = 10
maxscf = 512
geom_opt_max_cycles = 500
GEOM_OPT_TOL_ENERGY = 10
!GEOM_OPT_TOL_GRADIENT = 1000
!mem_static = 500
!n_frozen_core = fc
!set_iter = 100
$end

$basis
U 0
S 4 1.00
12.1252530 0.0219210
7.1615450 -0.2251600
4.7748360 0.5602990
2.0116930 -1.0712090
S 1 1.00
0.5868520 1.0000000
S 1 1.00
0.2791150 1.0000000
S 1 1.00
0.0633720 1.0000000
S 1 1.00
0.0256110 1.0000000
P 4 1.00
17.2547700 0.0013980
7.7353560 -0.0333460
5.1558780 0.1105780
2.2416700 -0.3172680
P 1 1.00
0.5818580 1.0000000
P 1 1.00
0.2679080 1.0000000
P 1 1.00
0.0834420 1.0000000
P 1 1.00
0.0321300 1.0000000
D 3 1.00
4.8410700 0.0057310
2.1601620 -0.0572360
0.5756300 0.2388280
D 1 1.00
0.2781360 1.0000000
D 1 1.00
0.1248790 1.0000000
D 1 1.00
0.0515480 1.0000000
F 3 1.00
2.4364410 0.3550110
1.1446820 0.4008460
0.5296930 0.3046790
F 1 1.00
0.2405960 1.0000000
F 1 1.00
0.1018670 1.0000000
G 1 1.00
1.1800000 1.0000000
G 1 1.00
0.4200000 1.0000000

Cl
6-31+G*

$end

$ecp
U 0
U-ECP 5 78
h potential
1
2 1.000000000 0.000000000
s-h potential
3
2 4.063653000 112.920103000
2 1.883995000 15.647500000
2 0.886567000 -3.689971000
p-h potential
3
2 3.986181000 118.758016000
2 2.000160000 15.077228000
2 0.960841000 0.556720000
d-h potential
3
2 4.147972000 60.855892000
2 2.234563000 29.280047000
2 0.913695000 4.998029000
f-h potential
3
2 3.998938000 49.924035000
2 1.998840000 -24.674042000
2 0.995641000 1.389480000
g-h potential
2
2 3.817422000 -36.040977000
2 0.262501000 -0.090511000

$end
###########################
##First-try-frequency##
$molecule
0 4
U 0.0035288101 0.0079902574 0.0000000000
Cl -2.3821731652 -0.9962577330 0.0000000000
Cl 0.2393270386 2.5824176319 0.0000000000
Cl 2.1754353949 -1.4198584495 0.0000000000

$end

$rem
jobtype = freq
method = b3lyp
ideriv = 2
doraman = 2
basis = general
ecp = general
PURECART 111
ecp_fit = true
unrestricted = true
symmetry = false
scf_algorithm = diis_gdm
scf_convergence = 5
scf_guess = gwh
scf_guess_mix = 10
maxscf = 512
geom_opt_max_cycles = 500
GEOM_OPT_TOL_ENERGY = 10
!GEOM_OPT_TOL_GRADIENT = 1000
!mem_static = 500
!n_frozen_core = fc
!set_iter = 100
$end

$basis
U 0
S 4 1.00
12.1252530 0.0219210
7.1615450 -0.2251600
4.7748360 0.5602990
2.0116930 -1.0712090
S 1 1.00
0.5868520 1.0000000
S 1 1.00
0.2791150 1.0000000
S 1 1.00
0.0633720 1.0000000
S 1 1.00
0.0256110 1.0000000
P 4 1.00
17.2547700 0.0013980
7.7353560 -0.0333460
5.1558780 0.1105780
2.2416700 -0.3172680
P 1 1.00
0.5818580 1.0000000
P 1 1.00
0.2679080 1.0000000
P 1 1.00
0.0834420 1.0000000
P 1 1.00
0.0321300 1.0000000
D 3 1.00
4.8410700 0.0057310
2.1601620 -0.0572360
0.5756300 0.2388280
D 1 1.00
0.2781360 1.0000000
D 1 1.00
0.1248790 1.0000000
D 1 1.00
0.0515480 1.0000000
F 3 1.00
2.4364410 0.3550110
1.1446820 0.4008460
0.5296930 0.3046790
F 1 1.00
0.2405960 1.0000000
F 1 1.00
0.1018670 1.0000000
G 1 1.00
1.1800000 1.0000000
G 1 1.00
0.4200000 1.0000000

Cl
6-31+G*

$end

$ecp
U 0
U-ECP 5 78
h potential
1
2 1.000000000 0.000000000
s-h potential
3
2 4.063653000 112.920103000
2 1.883995000 15.647500000
2 0.886567000 -3.689971000
p-h potential
3
2 3.986181000 118.758016000
2 2.000160000 15.077228000
2 0.960841000 0.556720000
d-h potential
3
2 4.147972000 60.855892000
2 2.234563000 29.280047000
2 0.913695000 4.998029000
f-h potential
3
2 3.998938000 49.924035000
2 1.998840000 -24.674042000
2 0.995641000 1.389480000
g-h potential
2
2 3.817422000 -36.040977000
2 0.262501000 -0.090511000

$end
###############################

##Second-try-all-combined##
$molecule
0 4
U -0.0116475938 -0.0655156332 0.0000000000
Cl -2.4452759521 -1.0138063243 0.0000000000
Cl 0.2384690355 2.5216700770 0.0000000000
Cl 2.2380787565 -1.3715562104 0.0000000000
$end

$rem
jobtype = opt
method = b3lyp
basis = general
ecp = general
purecart = 2111 !default
ecp_fit = true
unrestricted = true
SYM_IGNORE = 1
maxscf = 700
geom_opt_max_cycles = 999
GEOM_OPT_TOL_ENERGY = 11
GEOM_OPT_TOL_GRADIENT = 10
GEOM_OPT_TOL_DISPLACEMENT = 120
MEM_STATIC = 2000
MEM_TOTAL = 80000
GUI = 2
MOLDEN_FORMAT = 1
SCF_CONVERGENCE = 8
THRESH = 14
$end

$basis
U
srsc

Cl
def2-TZVP

$end

$ecp
U
srsc

$end

@@@

$molecule
read
$end

$rem
jobtype = freq
method = b3lyp
basis = general
ecp = general
purecart = 2111 !default
ecp_fit = true
unrestricted = true
SYM_IGNORE = 1
maxscf = 700
MEM_STATIC = 2000
MEM_TOTAL = 80000
GUI = 2
MOLDEN_FORMAT = 1
SCF_CONVERGENCE = 8
$end

$basis
U
srsc

Cl
def2-TZVP

$end

$ecp
U
srsc

$end

@@@

$molecule
read
$end

$rem
jobtype = opt
method = b3lyp
basis = general
ecp = general
purecart = 2111 !default
ecp_fit = true
unrestricted = true
SYM_IGNORE = 1
maxscf = 700
geom_opt_max_cycles = 999
GEOM_OPT_HESSIAN = READ
GEOM_OPT_TOL_ENERGY = 11
GEOM_OPT_TOL_GRADIENT = 10
GEOM_OPT_TOL_DISPLACEMENT = 120
MEM_STATIC = 2000
MEM_TOTAL = 80000
GUI = 2
MOLDEN_FORMAT = 1
SCF_CONVERGENCE = 8
THRESH = 14
$end

$basis
U
srsc

Cl
def2-TZVP

$end

$ecp
U
srsc

$end

@@@

$molecule
read
$end

$rem
jobtype = freq
method = b3lyp
ideriv = 2
doraman = 2
basis = general
ecp = general
PURECART 2111
ecp_fit = true
unrestricted = true
maxscf = 700
SYM_IGNORE = 1
MEM_STATIC = 2000
MEM_TOTAL = 80000
GUI = 2
MOLDEN_FORMAT = 1
SCF_CONVERGENCE = 8
$end

$basis
U
srsc

Cl
def2-TZVP

$end

$ecp
U
srsc

$end
#################################

Thank you,
Amir

#4

Amir,
A couple of things occur to me when looking at your first input:
(1) You have set the geometry optimization stopping criteria much tighter than the defaults. Presumably you did this in an effort to get rid of the imaginary frequency, which is not a bad idea necessarily but it’s not always guaranteed that such tight optimizations can converge.
(2) With default criteria, your first job converges for me in 11 cycles and does afford an imaginary frequency of 39.5 cm-1.
(3) Visualizing the imaginary mode, I find that it corresponds to symmetry breaking, specifically D3h to C2v symmetry reduction. Note that the optimizer typically won’t break point group symmetry (because the high-symmetry geometry is usually a saddle point with vanishing gradient along symmetry-breaking modes), so in such cases one often needs to break point group symmetry manually.
(4) In this particular case, it seems to me that there’s some possibility that this is a Jahn-Teller problem (D3h to C2v is one of the classic cases), in which case potential surface may have complicated topology due to strong interaction between two electronic states.

Bottom line is that I think your issues stem from the molecule not from the software.

#5

Thank you for looking into this, John.

So, based on your explanations, What would you recommend to overcome this issue?

Earlier you mentioned to tweak the DFT grid? How can I increase that? Just by specifying SG-2, or SG-3 in the input file?

Also, you mentioned to break the symmetry manually. Is there a way to do so in Q-Chem?

Thank you,
Amir

#6

Suggestion about the grid was before I took a look at your actual calculation. I think the issue is that your initial geometry has D3h symmetry, which is a saddle point at the level of theory in question but the optimizer won’t break point group symmetry. (The imaginary mode clearly corresponds to symmetry breaking, it’s a Cl-U-Cl scissors-type motion.) You need to distort the input geometry to lower the symmetry of the starting structure, either in some GUI (I like Molden for this, personally, because it has a Z-matrix editor) or in this particular case there’s so few atoms that you can just change the cartesians by hand.