Hello Everyone,
I am getting many imaginary frequencies in the excited state geometry of OPP-3. Though the gradient converges in the geometry optimization, the gradients are not zero in the frequency calculation. As suggested in this forum, I tried different ways but couldn’t succeed. The input and output is attached for the geometry optimization and frequency calculation. Any suggestion to get rid of this problem would be appreciated. Thanks
Input file
$molecule
0 1
H 0.1446284009 -2.1764962115 3.1413322777
C 0.0890791553 -1.2179736117 3.6560729739
C 0.0975519016 -1.2103389798 5.0416777197
H 0.1593554480 -2.1582157520 5.5818854771
C 0.0278478820 -0.0007005619 5.7575710642
H 0.0349160821 -0.0007264131 6.8493511617
C -0.0517219776 1.2090478938 5.0429267771
H -0.1078425508 2.1565007474 5.5844942360
C -0.0601505447 1.2173894372 3.6573202923
H -0.1238071266 2.1760730243 3.1437605493
C 0.0100777611 -0.0002588568 2.8965753553
H 0.1321022613 -2.1755541752 -1.2011152345
C 0.0616647035 -1.2174628037 -0.6874373754
C 0.0713444399 -1.2171455630 0.6851170722
H 0.1377747259 -2.1765831003 1.1964980197
C 0.0024494476 -0.0001049989 1.4565235241
C -0.0748632891 1.2167763863 0.6856382527
H -0.1390085694 2.1760568279 1.1974972855
C -0.0754682685 1.2170069995 -0.6868229487
H -0.1522906817 2.1749483600 -1.1999255486
C -0.0083593613 -0.0000292387 -1.4572032194
H -0.1779529688 -2.1549474947 -5.5827269241
C -0.1036135722 -1.2085598567 -5.0415278583
C -0.0973143475 -1.2166432053 -3.6558323588
H -0.1805718612 -2.1723270891 -3.1391675772
C -0.0114247609 0.0004680271 -2.8972129005
C 0.0717135466 1.2178921160 -3.6559597053
H 0.1613266970 2.1734284931 -3.1399813200
C 0.0682056316 1.2099977243 -5.0417601064
H 0.1396891843 2.1565232109 -5.5831225589
C -0.0206057644 0.0008785853 -5.7563070064
H -0.0247316241 0.0010800635 -6.8481422608
$end
$rem
- basis = cc-pVDZ*
- cis_state_deriv = 1*
- job_type = optimization*
- geom_opt_max_cycles = 300*
- max_scf_cycles = 200*
- method = wB97XD*
- cis_n_roots = 10*
- cis_singlets = 1*
- cis_triplets = 0*
- rpa = 2*
- mem_static = 2000*
- mem_total = 16000*
- scf_convergence = 8*
- symmetry = false*
- XC_GRID = 000099000590*
- symmetry_ignore = 1*
- thresh = 14*
$end
@@@
$molecule
$rem
basis = cc-pVDZ
cis_state_deriv = 1
job_type = freq
method = wB97XD
cis_n_roots = 10
cis_singlets = 1
cis_triplets = 0
rpa = 2
XC_GRID = 000099000590
SCF_CONVERGENCE = 8
mem_static = 4000
mem_total = 192000
thresh = 14
$end
Optimization OUT file
Maximum Tolerance Cnvgd?
Gradient 0.000091 0.000300 YES
Displacement 0.019924 0.001200 NO
Energy change -0.000000 0.000001 YES
Final energy is -694.013898393683
** OPTIMIZATION CONVERGED **
Coordinates (Angstroms)
ATOM X Y Z
1 H -0.0228861727 -2.1770107073 3.1305871859
2 C 0.0014135895 -1.2159397697 3.6424735402
3 C 0.0049963982 -1.2078253636 5.0248516652
4 H -0.0075559366 -2.1563781089 5.5653616449
5 C 0.0215498042 -0.0009012862 5.7359759579
6 H 0.0253847586 -0.0010813033 6.8270865204
7 C 0.0330596546 1.2062713655 5.0251723710
8 H 0.0492261197 2.1546252456 5.5659288461
9 C 0.0270253400 1.2148901296 3.6428116774
10 H 0.0474598753 2.1761589310 3.1311199324
11 C 0.0117132390 -0.0003933666 2.8858349130
12 H 0.1565931938 -2.1748112161 -1.1939412529
13 C 0.0785409485 -1.2164433183 -0.6830108114
14 C 0.0837378865 -1.2166372426 0.6820862436
15 H 0.1662788930 -2.1750692741 1.1922147576
16 C 0.0069628070 -0.0001252880 1.4531322952
17 C -0.0748840537 1.2166865512 0.6830965051
18 H -0.1522942448 2.1750656504 1.1941257107
19 C -0.0805213435 1.2168818429 -0.6819947521
20 H -0.1616293704 2.1754837802 -1.1919835194
21 C -0.0045638585 0.0003100508 -1.4531271175
22 H -0.0388258140 -2.1549883442 -5.5657491777
23 C -0.0279081849 -1.2065095023 -5.0250868440
24 C -0.0202560783 -1.2149695595 -3.6427541249
25 H -0.0334749958 -2.1763406790 -3.1310425733
26 C -0.0112103665 0.0004576345 -2.8857939708
27 C -0.0092730584 1.2159929437 -3.6426482571
28 H 0.0085958494 2.1773027411 -3.1309727588
29 C -0.0144987666 1.2077035353 -5.0250016618
30 H -0.0086525145 2.1562565691 -5.5656267809
31 C -0.0245026696 0.0006382861 -5.7360157298
32 H -0.0296009285 0.0006990569 -6.8271132992
Frequency analysis output
Gradient
Gradient of the CIS Excitation Energy
1 2 3 4 5 6
1 -0.0818600 -0.0562712 0.0633614 0.0902614 0.0780868 0.1824819
2 -0.1634910 -0.0239137 -0.0700915 -0.1578991 -0.0000030 -0.0000021
3 0.0058927 -0.0040876 0.0017728 0.0051980 0.0000284 0.0000877
7 8 9 10 11 12
1 0.0633651 0.0902780 -0.0562805 -0.0818442 -0.0802436 -0.0821408
2 0.0700828 0.1578918 0.0239159 0.1634990 0.0000031 -0.1638903
3 -0.0017176 -0.0050775 0.0039996 -0.0059301 -0.0000607 -0.0104569
13 14 15 16 17 18
1 -0.0848968 0.0849027 0.0821520 0.0665206 0.0849316 0.0821544
2 -0.0131081 -0.0132268 -0.1638684 0.0000047 0.0132109 0.1638857
3 0.0045421 0.0048996 -0.0104792 0.0000235 -0.0048357 0.0103053
19 20 21 22 23 24
1 -0.0849157 -0.0821354 -0.0665771 -0.0902612 -0.0633862 0.0562838
2 0.0130829 0.1639091 0.0000132 -0.1579293 -0.0701049 -0.0237908
3 -0.0045081 0.0102887 0.0000104 0.0041990 0.0013829 -0.0038746
25 26 27 28 29 30
1 0.0818549 0.0803021 0.0562870 0.0818280 -0.0633991 -0.0902871
2 -0.1635434 0.0000200 0.0237903 0.1635560 0.0700943 0.1579186
3 0.0047443 -0.0000780 0.0037794 -0.0048380 -0.0013012 -0.0040731
31 32
1 -0.0780691 -0.1824838
2 -0.0000020 -0.0000139
3 0.0000401 0.0001236
Frequency information
This Molecule has 30 Imaginary Frequencies
Zero point vibrational energy: 156.886 kcal/mol
First three modes
Mode: 1 2 3
Frequency: -3650.47 -3558.15 -3233.34
Force Cnst: 13.7111 12.5301 10.0959
Red. Mass: 1.7463 1.6798 1.6391
IR Active: YES YES YES
IR Intens: 10.707 5.423 4.909
Raman Active: YES YES YES
Something is going very wrong here, because these aren’t small imag frequencies (which might indicate some kind of numerical or threshold issue) but very large ones. What are the excitation energies for the first couple of states, at the end (final step) of the geom opt and then again in the subsequent freq calculation? Is it possible that somehow the freq calculation is using a different electronic state than what the optimization was using?
Another suggestion is to take your optimized geometry and redo this calculation with RPA=TRUE rather than RPA=2. That’s just a hunch, checking for possible bugs.
Dear Prof.,
Thank you for the reply. In the previous calculation, the excitation energies are close in both cases. I also tried using RPA= TRUE instead of RPA=2, as your suggestion. The imaginary frequencies are still there, and excitation energies are very close. Below are excitation energies from the calculation RPA=true.
Excited state info. from optimization.
---------------------------------------------------
-
TDDFT Excitation Energies*
-
---------------------------------------------------*
-
Excited state 1: excitation energy (eV) = 3.7508*
-
Total energy for state 1: -694.01389876 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: 0.0143 X -0.0005 Y 3.5729 Z*
-
Strength : 1.1730966070*
-
X: D( 61) → V( 1) amplitude = 0.9772*
-
Excited state 2: excitation energy (eV) = 4.5565*
-
Total energy for state 2: -693.98428933 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: 0.0094 X 0.0652 Y 0.0000 Z*
-
Strength : 0.0004843376*
-
X: D( 57) → V( 1) amplitude = 0.3434*
-
X: D( 59) → V( 1) amplitude = -0.4513*
-
X: D( 61) → V( 2) amplitude = 0.7443*
-
Excited state 3: excitation energy (eV) = 4.9281*
-
Total energy for state 3: -693.97063436 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: 0.0001 X -0.0003 Y -0.0000 Z*
-
Strength : 0.0000000119*
-
X: D( 58) → V( 1) amplitude = 0.6205*
-
X: D( 59) → V( 3) amplitude = 0.3096*
-
X: D( 60) → V( 2) amplitude = -0.3016*
-
X: D( 61) → V( 4) amplitude = 0.5838*
-
Excited state 4: excitation energy (eV) = 5.0254*
-
Total energy for state 4: -693.96705687 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: -0.0114 X 0.2481 Y 0.0001 Z*
-
Strength : 0.0075940279*
-
X: D( 57) → V( 1) amplitude = 0.6574*
-
X: D( 58) → V( 3) amplitude = 0.3135*
-
X: D( 59) → V( 1) amplitude = 0.2681*
-
X: D( 60) → V( 4) amplitude = 0.3051*
-
X: D( 61) → V( 5) amplitude = -0.4824*
Excited state information from frequency calculation.
---------------------------------------------------
-
TDDFT Excitation Energies*
-
---------------------------------------------------*
-
Excited state 1: excitation energy (eV) = 3.7508*
-
Total energy for state 1: -694.01389761 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: 3.5725 X -0.0000 Y 0.0000 Z*
-
Strength : 1.1728223422*
-
X: D( 61) → V( 1) amplitude = 0.9771*
-
Excited state 2: excitation energy (eV) = 4.5565*
-
Total energy for state 2: -693.98428819 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: -0.0001 X -0.0651 Y 0.0108 Z*
-
Strength : 0.0004856456*
-
X: D( 57) → V( 1) amplitude = -0.3435*
-
X: D( 59) → V( 1) amplitude = 0.4513*
-
X: D( 61) → V( 2) amplitude = 0.7443*
-
Excited state 3: excitation energy (eV) = 4.9281*
-
Total energy for state 3: -693.97063331 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: -0.0000 X -0.0003 Y -0.0001 Z*
-
Strength : 0.0000000119*
-
X: D( 58) → V( 1) amplitude = 0.6205*
-
X: D( 59) → V( 3) amplitude = -0.3096*
-
X: D( 60) → V( 2) amplitude = 0.3016*
-
X: D( 61) → V( 4) amplitude = -0.5838*
-
Excited state 4: excitation energy (eV) = 5.0254*
-
Total energy for state 4: -693.96705575 au*
-
Multiplicity: Singlet*
-
Trans. Mom.: -0.0000 X 0.2484 Y 0.0072 Z*
-
Strength : 0.0076038125*
-
X: D( 57) → V( 1) amplitude = 0.6574*
-
X: D( 58) → V( 3) amplitude = -0.3135*
-
X: D( 59) → V( 1) amplitude = 0.2681*
-
X: D( 60) → V( 4) amplitude = -0.3051*
-
X: D( 61) → V( 5) amplitude = 0.4823*
Hmm. They are not that close, i.e., not close enough that I think the calculation could actually be switching states. Would you mind providing the full input file? I would like to look into this because I’m out of obvious suggestions. When you do, please use the button for preformatted text (looks like </>) so that it doesn’t show up as a bunch of bullet points.
Here is the full input file.
$molecule
0 1
H -0.0228861727 -2.1770107073 3.1305871859
C 0.0014135895 -1.2159397697 3.6424735402
C 0.0049963982 -1.2078253636 5.0248516652
H -0.0075559366 -2.1563781089 5.5653616449
C 0.0215498042 -0.0009012862 5.7359759579
H 0.0253847586 -0.0010813033 6.8270865204
C 0.0330596546 1.2062713655 5.0251723710
H 0.0492261197 2.1546252456 5.5659288461
C 0.0270253400 1.2148901296 3.6428116774
H 0.0474598753 2.1761589310 3.1311199324
C 0.0117132390 -0.0003933666 2.8858349130
H 0.1565931938 -2.1748112161 -1.1939412529
C 0.0785409485 -1.2164433183 -0.6830108114
C 0.0837378865 -1.2166372426 0.6820862436
H 0.1662788930 -2.1750692741 1.1922147576
C 0.0069628070 -0.0001252880 1.4531322952
C -0.0748840537 1.2166865512 0.6830965051
H -0.1522942448 2.1750656504 1.1941257107
C -0.0805213435 1.2168818429 -0.6819947521
H -0.1616293704 2.1754837802 -1.1919835194
C -0.0045638585 0.0003100508 -1.4531271175
H -0.0388258140 -2.1549883442 -5.5657491777
C -0.0279081849 -1.2065095023 -5.0250868440
C -0.0202560783 -1.2149695595 -3.6427541249
H -0.0334749958 -2.1763406790 -3.1310425733
C -0.0112103665 0.0004576345 -2.8857939708
C -0.0092730584 1.2159929437 -3.6426482571
H 0.0085958494 2.1773027411 -3.1309727588
C -0.0144987666 1.2077035353 -5.0250016618
H -0.0086525145 2.1562565691 -5.5656267809
C -0.0245026696 0.0006382861 -5.7360157298
H -0.0296009285 0.0006990569 -6.8271132992
$end
$rem
basis = cc-pVDZ
cis_state_deriv = 1
job_type = optimization
geom_opt_max_cycles = 300
max_scf_cycles = 200
method = wB97XD
cis_n_roots = 10
cis_singlets = 1
cis_triplets = 0
rpa = true
mem_static = 2000
mem_total = 16000
scf_convergence = 8
symmetry = false
XC_GRID = 000099000590
symmetry_ignore = 1
thresh = 14
$end
@@@
$molecule
read
$end
$rem
basis = cc-pVDZ
cis_state_deriv = 1
job_type = freq
method = wB97XD
cis_n_roots = 10
cis_singlets = 1
cis_triplets = 0
rpa = true
XC_GRID = 000099000590
SCF_CONVERGENCE = 8
mem_static = 4000
mem_total = 192000
thresh = 14
$end
What version of Q-Chem are you using? With v. 6.0, I get different behavior. With your input file, the RPA iterations fail to converge:
30 9 1 0.000000 0.000000
31 9 1 0.000000 0.000000
Maximum number of CIS/TDDFT iterations has been reached
This limit can be increased using CIS_MAX_CYCLES
Note that it comes very close to converging: one root left. When I increase CIS_MAX_CYCLES as suggested, the calculation finishes and all the frequencies are real.
For this calculation, I used 5.4.2.
Thank you for the test result. I will check again with v. 6.0.