Error in TD-DFT optimization using STATE_FOLLOW

Questions
#1

Hello. I came across an error message while trying to optimize the second triplet excited state
of thioxanthon at the TD-CAM-B3LYP-D3/6-31G(d,p) level of theory. As the same instruction used
to work seamlessly in QC5.4, I am suspicious of a bug in the latest release of QC.

I am attaching my QC input for your reference.

=====================================

$molecule
0 1
C -3.387898679 -0.047938516 1.918901858
C -2.646082575 0.247490346 3.077426413
C -3.116324491 0.563145004 0.710807635
H -4.186923706 -0.780344092 1.970484532
C -1.618905412 1.164214412 2.992833879
H -2.876220104 -0.236852641 4.019050051
C -1.319373012 1.771779311 1.762711938
H -1.035277544 1.420458313 3.872087395
C -2.056202869 1.481463311 0.572375724
S 0.008402330 2.857489590 1.747035835
H -3.677559552 0.309211655 -0.179378544
C -1.759895583 2.051905098 -0.716415305
C 0.291157096 3.289813284 0.115753032
C -0.563473560 2.839517434 -0.929134896
C 1.448719326 4.036716093 -0.136041945
O -2.451980741 1.777038486 -1.657407703
C -0.181543921 3.192511680 -2.236347071
C 1.795303875 4.354569668 -1.436794348
H 2.072966374 4.354460043 0.693715268
C 0.966394388 3.924976581 -2.482164440
H 2.690990853 4.929481198 -1.641230612
H 1.233688008 4.161788798 -3.507486196
H -0.813031116 2.856714504 -3.048594455
$end

$rem
BASIS 6-31G(d,p)
JOB_TYPE OPTIMIZATION
METHOD CAM-B3LYP
DFT_D D3_ZERO
SYMMETRY FALSE
mem_total 400000
mem_static 4000
max_SCF_CYCLES 150
MAX_CIS_CYCLES 150
SOLVENT_METHOD PCM
CIS_N_ROOTS 4
SET_STATE_DERIV 2
STATE_FOLLOW TRUE
RPA 0
CIS_SINGLETS FALSE
CIS_TRIPLETS TRUE
$end

$pcm
Theory CPCM
$end

$solvent
DIELECTRIC 37.5
$end

Excited state optimization and charge transfer stare
#2

Please attach the error message and tell us at what point you encountered it during the job. I believe that the default optimizer changed in the 6.0 release, and you can use the old (v.5.4) optimizer by specifying GEOM_OPT_DRIVER = OPTIMIZE in $rem.

#3

Dear John,

I thank you for your suggestion. As you anticipated,
GEOM_OPT_DRIVE=OPTIMIZER suppresses the
unusual termination and allows the optimization progress.

Attached below please find the error message, while using
the QC6.0 default optimizer.

All the best,
Jiyong

=======================================
Q-Chem 6.0.1 for Intel X86 EM64T Linux

Parts of Q-Chem use Armadillo 9.900.5 (Nocturnal Misbehaviour).
http://arma.sourceforge.net/

Q-Chem begins on Fri Apr 14 23:08:17 2023

Host:
0

 Scratch files written to /mnt/lustre/ibs/pjy/.tmp/qchem908764//

Aug2322 |scratch|qcdevops|jenkins|workspace|build_RNUM 6470
Processing $rem in /opt/ibs_lib/apps/qchem/6.0.1/config/preferences:
Processing $rem in /mnt/lustre/ibs/pjy/.qchemrc:
Symmetry turned off for PCM/SM12/SMD calculation

Checking the input file for inconsistencies… …done.

User input:

$molecule
0 1
C -3.387898679 -0.047938516 1.918901858
C -2.646082575 0.247490346 3.077426413
C -3.116324491 0.563145004 0.710807635
H -4.186923706 -0.780344092 1.970484532
C -1.618905412 1.164214412 2.992833879
H -2.876220104 -0.236852641 4.019050051
C -1.319373012 1.771779311 1.762711938
H -1.035277544 1.420458313 3.872087395
C -2.056202869 1.481463311 0.572375724
S 0.008402330 2.857489590 1.747035835
H -3.677559552 0.309211655 -0.179378544
C -1.759895583 2.051905098 -0.716415305
C 0.291157096 3.289813284 0.115753032
C -0.563473560 2.839517434 -0.929134896
C 1.448719326 4.036716093 -0.136041945
O -2.451980741 1.777038486 -1.657407703
C -0.181543921 3.192511680 -2.236347071
C 1.795303875 4.354569668 -1.436794348
H 2.072966374 4.354460043 0.693715268
C 0.966394388 3.924976581 -2.482164440
H 2.690990853 4.929481198 -1.641230612
H 1.233688008 4.161788798 -3.507486196
H -0.813031116 2.856714504 -3.048594455
$end

$rem
SCF_CONVERGENCE 8
BASIS 6-31G(d,p)
JOB_TYPE OPTIMIZATION
METHOD CAM-B3LYP
DFT_D D3_ZERO
SYMMETRY FALSE
mem_total 400000
mem_static 4000
max_SCF_CYCLES 150
MAX_CIS_CYCLES 150
SOLVENT_METHOD PCM
CIS_N_ROOTS 4
SET_STATE_DERIV 2
STATE_FOLLOW TRUE
RPA 0
CIS_SINGLETS FALSE
CIS_TRIPLETS TRUE
$end

$pcm
Theory CPCM
$end

$solvent
DIELECTRIC 37.5
$end

Adding empirical dispersion correction

Doing -D3(0) with
s6 = 1.00000
s8 = 1.21700
rs6 = 1.37800
rs8 = 1.00000
-D3 energy without 3body term = -0.0136364624 hartrees

         Standard Nuclear Orientation (Angstroms)
I     Atom           X                Y                Z

1      C      -3.3878986790    -0.0479385160     1.9189018580
2      C      -2.6460825750     0.2474903460     3.0774264130
3      C      -3.1163244910     0.5631450040     0.7108076350
4      H      -4.1869237060    -0.7803440920     1.9704845320
5      C      -1.6189054120     1.1642144120     2.9928338790
6      H      -2.8762201040    -0.2368526410     4.0190500510
7      C      -1.3193730120     1.7717793110     1.7627119380
8      H      -1.0352775440     1.4204583130     3.8720873950
9      C      -2.0562028690     1.4814633110     0.5723757240

10 S 0.0084023300 2.8574895900 1.7470358350
11 H -3.6775595520 0.3092116550 -0.1793785440
12 C -1.7598955830 2.0519050980 -0.7164153050
13 C 0.2911570960 3.2898132840 0.1157530320
14 C -0.5634735600 2.8395174340 -0.9291348960
15 C 1.4487193260 4.0367160930 -0.1360419450
16 O -2.4519807410 1.7770384860 -1.6574077030
17 C -0.1815439210 3.1925116800 -2.2363470710
18 C 1.7953038750 4.3545696680 -1.4367943480
19 H 2.0729663740 4.3544600430 0.6937152680
20 C 0.9663943880 3.9249765810 -2.4821644400
21 H 2.6909908530 4.9294811980 -1.6412306120
22 H 1.2336880080 4.1617887980 -3.5074861960
23 H -0.8130311160 2.8567145040 -3.0485944550

Nuclear Repulsion Energy = 1001.74869271 hartrees
There are 55 alpha and 55 beta electrons
Requested basis set is 6-31G(d,p)
There are 85 shells and 269 basis functions

Total QAlloc Memory Limit 400000 MB
Mega-Array Size 3911 MB
MEM_STATIC part 4000 MB

            STARTING GEOMETRY OPTIMIZER USING LIBOPT3
            by Peter F. McLaughlin, Yu Zhang, Evgeny Epifanovsky

          Initial Energy and Gradient Calculation

– Checking Topology for Ill-Behaving Coordinates –

Coordinate | Removed | Additions

Bonds                     0          0
Angles                    0          0
Torsions                  0          0
Co-Linear Type5s          0          0
Co-Linear Type6s          0          0

 Done Checking Topology

Adding empirical dispersion correction

Doing -D3(0) with
s6 = 1.00000
s8 = 1.21700
rs6 = 1.37800
rs8 = 1.00000
-D3 energy without 3body term = -0.0136364624 hartrees

         Standard Nuclear Orientation (Angstroms)
I     Atom           X                Y                Z

1      C      -3.3878986790    -0.0479385160     1.9189018580
2      C      -2.6460825750     0.2474903460     3.0774264130
3      C      -3.1163244910     0.5631450040     0.7108076350
4      H      -4.1869237060    -0.7803440920     1.9704845320
5      C      -1.6189054120     1.1642144120     2.9928338790
6      H      -2.8762201040    -0.2368526410     4.0190500510
7      C      -1.3193730120     1.7717793110     1.7627119380
8      H      -1.0352775440     1.4204583130     3.8720873950
9      C      -2.0562028690     1.4814633110     0.5723757240

10 S 0.0084023300 2.8574895900 1.7470358350
11 H -3.6775595520 0.3092116550 -0.1793785440
12 C -1.7598955830 2.0519050980 -0.7164153050
13 C 0.2911570960 3.2898132840 0.1157530320
14 C -0.5634735600 2.8395174340 -0.9291348960
15 C 1.4487193260 4.0367160930 -0.1360419450
16 O -2.4519807410 1.7770384860 -1.6574077030
17 C -0.1815439210 3.1925116800 -2.2363470710
18 C 1.7953038750 4.3545696680 -1.4367943480
19 H 2.0729663740 4.3544600430 0.6937152680
20 C 0.9663943880 3.9249765810 -2.4821644400
21 H 2.6909908530 4.9294811980 -1.6412306120
22 H 1.2336880080 4.1617887980 -3.5074861960
23 H -0.8130311160 2.8567145040 -3.0485944550

Nuclear Repulsion Energy = 1001.74869271 hartrees
There are 55 alpha and 55 beta electrons
Using C-PCM dielectric factor f(eps) = (eps-1)/eps = 0.973333
Discretize the solute cavity surface with Lebedev spheres
Using 110 Lebedev grid points for each H atom
Using 194 Lebedev grid points for other atoms
Atomic van der Waals radii (Bondi set) will be scaled by 1.20
Using smooth SwiG surface discretization
Remove points where switching function is < 1.0e-08
Keep 1538 surface tesserae and discard 2252 interior tesserae
Molecular Surface Area = 230.623 Angst**2

Requested basis set is 6-31G(d,p)
There are 85 shells and 269 basis functions
A cutoff of 1.0D-12 yielded 2624 shell pairs
There are 28769 function pairs
Smallest overlap matrix eigenvalue = 4.25E-04

Scale SEOQF with 1.000000e-01/1.000000e-01/1.000000e-01

Standard Electronic Orientation quadrupole field applied
Nucleus-field energy = 0.0000000287 hartrees
Guess from superposition of atomic densities
Warning: Energy on first SCF cycle will be non-variational
SAD guess density has 110.000000 electrons

General SCF calculation program by
Eric Jon Sundstrom, Paul Horn, Yuezhi Mao, Dmitri Zuev, Alec White,
David Stuck, Shaama M.S., Shane Yost, Joonho Lee, David Small,
Daniel Levine, Susi Lehtola, Hugh Burton, Evgeny Epifanovsky,
Bang C. Huynh

Exchange: 0.1900 Hartree-Fock + 0.3500 B88 + 0.4600 muB88 + LR-HF
Correlation: 0.1900 VWN5 + 0.8100 LYP
Using SG-1 standard quadrature grid
Dispersion: Grimme D3
using 14 threads for integral computing

OpenMP Integral computing Module
Release: version 1.0, May 2013, Q-Chem Inc. Pittsburgh

C-PCM solvent model [f = (eps-1)/eps], solve by matrix inversion
A restricted SCF calculation will be
performed using DIIS
SCF converges when DIIS error is below 1.0e-08

Cycle Energy DIIS error

1    -978.0598748457      4.80e-02  
2    -973.1782858908      4.25e-03  
3    -973.1093666271      5.67e-03  
4    -973.3317975271      7.20e-04  
5    -973.3357878339      2.96e-04  
6    -973.3365646538      6.55e-05  
7    -973.3366180802      2.35e-05  
8    -973.3366271987      1.10e-05  
9    -973.3366299269      3.12e-06

10 -973.3366302297 8.52e-07
11 -973.3366302582 4.44e-07
12 -973.3366302631 1.40e-07
13 -973.3366302638 7.22e-08
14 -973.3366302639 1.90e-08
15 -973.3366302639 7.77e-09 Convergence criterion met

SCF time: CPU 267.45s wall 21.00s

************** Final PCM Free Energy Summary **************
G_electrostatic = -0.01223071 hartree = -7.67488518 kcal/mol
G_cavitation = 0.00000000 hartree = 0.00000000 kcal/mol
G_dispersion = 0.00000000 hartree = 0.00000000 kcal/mol
G_repulsion = 0.00000000 hartree = 0.00000000 kcal/mol

Non-electrostatic Free Energy = 0.00000000 hartree = 0.00000000 kcal/mol
Total = -0.01223071 hartree = -7.67488518 kcal/mol

SCF Energy (H0 + V/2) = -973.33663026
Solute Internal Energy (H0) = -973.32439956
Total Free Energy (H0 + V/2 + non-elec) = -973.33663026 hartree
= -610777.95396184 kcal/mol

SCF energy in the final basis set = -973.3366302639
Total energy in the final basis set = -973.3366302639

MAX_CIS_SUBSPACE fits in available memory
Direct TDDFT/TDA calculation will be performed
Exchange: 0.1900 Hartree-Fock + 0.3500 B88 + 0.4600 muB88 + LR-HF
Correlation: 0.1900 VWN5 + 0.8100 LYP
Using SG-1 standard quadrature grid
Dispersion: Grimme D3
Triplet excitation energies requested
Adding PCM contribution to the XC response (LR-PCM)
CIS energy converged when residual is below 10e- 6

Iter Rts Conv Rts Left Ttl Dev Max Dev

1 0 4 0.001302 0.000375
2 0 4 0.000517 0.000150
3 0 4 0.000292 0.000116
4 0 4 0.000350 0.000159
5 0 4 0.000436 0.000310
6 0 4 0.000181 0.000112
7 1 3 0.000066 0.000038
8 1 3 0.000017 0.000010
9 3 1 0.000004 0.000002
10 4 0 0.000002 0.000001 Roots Converged

     TDDFT/TDA Excitation Energies

Excited state 1: excitation energy (eV) = 2.9986
Total energy for state 1: -973.22643290 au
Multiplicity: Triplet
Trans. Mom.: 0.0000 X 0.0000 Y 0.0000 Z
Strength : 0.0000000000
D( 53) → V( 1) amplitude = 0.2295
D( 55) → V( 1) amplitude = 0.9346

Excited state 2: excitation energy (eV) = 3.6048
Total energy for state 2: -973.20415465 au
Multiplicity: Triplet
Trans. Mom.: 0.0000 X 0.0000 Y 0.0000 Z
Strength : 0.0000000000
D( 51) → V( 4) amplitude = -0.2601
D( 52) → V( 1) amplitude = 0.3593
D( 53) → V( 3) amplitude = 0.2322
D( 54) → V( 1) amplitude = 0.5129
D( 55) → V( 2) amplitude = 0.6166

Excited state 3: excitation energy (eV) = 3.6253
Total energy for state 3: -973.20340344 au
Multiplicity: Triplet
Trans. Mom.: 0.0000 X 0.0000 Y 0.0000 Z
Strength : 0.0000000000
D( 52) → V( 1) amplitude = 0.8457
D( 52) → V( 6) amplitude = 0.2148
D( 55) → V( 2) amplitude = -0.3052

Excited state 4: excitation energy (eV) = 3.7439
Total energy for state 4: -973.19904402 au
Multiplicity: Triplet
Trans. Mom.: 0.0000 X 0.0000 Y 0.0000 Z
Strength : 0.0000000000
D( 51) → V( 3) amplitude = -0.3352
D( 53) → V( 1) amplitude = 0.6940
D( 54) → V( 2) amplitude = 0.2675
D( 54) → V( 3) amplitude = 0.2267
D( 55) → V( 1) amplitude = -0.2519
D( 55) → V( 4) amplitude = -0.3521

Calculating Relaxed Density

Iter Rts Conv Rts Left Ttl Dev Max Dev

1 0 4 0.000716 0.000255
2 0 4 0.000276 0.000133
3 0 4 0.000075 0.000032
4 0 4 0.000025 0.000011
5 0 4 0.000010 0.000005
6 3 1 0.000001 0.000001
7 4 0 0.000000 0.000000 Roots Converged

            Calling CISFollow

Q-Chem fatal error occurred in module 0, line 45:

Invalid job for state following

#4

I have verified this crash using a pre-release v.6.1. I suspect that state-following is perhaps not yet implemented with the new optimizer and there should be a switch in place to reroute such jobs to the old optimizer. I will submit a bug ticket.

#5

I still having same error even with the old optimizer. Here is the input file
$molecule
-1 1
F -6.0673232268 0.4384099957 0.8777204772
B -5.0151840632 0.1597556216 -0.0077322936
F -5.5410457645 -0.0173667296 -1.2991917943
N -4.2410234775 -1.0991214983 0.4415777674
C -4.7743887818 -2.2797379637 0.8628765235
C -3.7637385217 -3.2016400389 1.0527824738
C -2.5488496273 -2.5521940657 0.7318972061
C -2.8744230853 -1.2451357629 0.3595652485
C -2.0523357943 -0.1309201970 0.0358284437
C -2.6175257949 1.1591302444 -0.0699665138
C -2.0420099455 2.4362939452 -0.1578325608
C -3.1014632552 3.3744392765 -0.1299552892
C -4.2721562025 2.6497260406 -0.0322907693
N -3.9873954672 1.3197321387 -0.0017838052
H -5.8442423332 -2.3764395636 0.9880538711
H -3.8882539261 -4.2240502061 1.3884981200
H -1.5470107762 -2.9588715548 0.7810695555
H -0.9856571663 -0.2802458600 -0.0282561936
H -0.9842373792 2.6449350120 -0.1756642447
H -3.0201921883 4.4537853612 -0.1708447674
H -5.3005162412 2.9824139773 0.0081819608
C 2.1999691771 0.6151405009 -1.7600452950
C 1.5875775889 1.1096086956 0.7789336746
C 1.8475917130 -1.4807870368 0.0651982220
C 3.9821367073 0.1105957554 0.1122666936
C 0.8306109055 1.0275964279 -2.2885420790
C 0.8277245838 1.2821194882 -3.8039790066
C 1.8537618176 2.5112950603 0.6176193189
C 2.1246033767 -2.1566177239 1.2656475109
C 4.6059685543 0.8708020496 1.1123464138
C 0.6490016052 0.7609122273 1.8003061877
C 1.2041293659 -2.2283507094 -0.9295285882
C 4.8215609071 -0.6692096227 -0.7003197107
C 0.9995232234 0.0183517986 -4.6475157475
C 1.2659199218 3.4612888172 1.4121520770
C 1.7716127619 -3.4856199098 1.4707148724
C 5.9894974831 0.8590773890 1.2898874162
C 0.0587970633 1.7130832564 2.6031212515
C 0.8515821245 -3.5641336663 -0.7432274761
C 6.2000622259 -0.6879382124 -0.5331614441
C 0.3482251589 3.0658547505 2.4081966895
C 1.1327071330 -4.2001748877 0.4597743379
C 6.7949714537 0.0842169602 0.4646127393
B 2.3490560652 0.0539419321 -0.1267563003
H 2.6137409802 -0.2292420261 -2.3235875582
H 2.9044299594 1.4441858991 -1.9004036726
H 0.5073128539 1.9427113103 -1.7823404637
H 0.0528013633 0.2901390051 -2.0639481561
H 1.6170839138 2.0051042498 -4.0569288288
H -0.1257885700 1.7576743498 -4.0641212919
H 2.5559680742 2.8117152831 -0.1545079870
H 2.6518162228 -1.6251365363 2.0585862378
H 4.0017283668 1.4755701535 1.7847723469
H 0.4053603052 -0.2829843993 1.9516149583
H 0.9633271622 -1.7647515977 -1.8832267662
H 4.3773706571 -1.2970741332 -1.4703165345
H 1.9632253131 -0.4704366284 -4.4642440176
H 0.9452430399 0.2447387340 -5.7183274266
H 0.2107545283 -0.7079652817 -4.4198415039
H 1.4827324414 4.5152915387 1.2685443790
H 2.0004492606 -3.9675018665 2.4188385797
H 6.4349539195 1.4548649420 2.0839190805
H -0.6578346882 1.4171349702 3.3603214149
H 0.3419433311 -4.1015181049 -1.5389844296
H 6.8153705421 -1.3102395451 -1.1792370478
H -0.1450000068 3.8180363091 3.0152585912
H 0.8528530644 -5.2394913244 0.6124952886
H 7.8738311731 0.0720650944 0.6000579353
$end

$rem
JOBTYPE opt
EXCHANGE wb97x-d3
BASIS 6-31G(d)
CIS_N_ROOTS 5
CIS_SINGLETS true
CIS_TRIPLETS false
CIS_STATE_DERIV 2
GEOM_OPT_DRIVER OPTIMIZE
STATE_FOLLOW true
RPA false
GUI=2
$end

#6

This input crashes as well with the latest trunk; I’ve updated the ticket to reflect that. It appears this feature was broken some time since the 5.4 release.