Facing problem in spin contamination

Hi,
I am using SF-TDDFT. While optimizing the second excited state, I am facing spin contamination. Please tell me how I can reduce this value.

Here is the input:
$molecule
0 3
coordinates
$end

$rem
BASIS = GEN
PURECART = 1
EXCHANGE = omegaB97X-D
omega = 156
GUI = 2
JOB_TYPE = OPTIMIZATION
MAX_CIS_CYCLES = 700
SCF_CONVERGENCE = 8
MAX_SCF_CYCLES = 700
THRESH = 14
MOM_START = 5
SPIN_FLIP = 1
CIS_N_ROOTS = 4
CIS_STATE_DERIVATIVE = 2
MEM_TOTAL = 30000
SOLVENT_METHOD = PCM
$end

$solvent
DIELECTRIC 32.613000
OPTICALDIELECTRIC 1.765709
$end

$basis
H 0
S 3 1.00
13.0107010 0.19682158D-01
1.9622572 0.13796524
0.44453796 0.47831935
S 1 1.00
0.12194962 1.0000000
P 1 1.00
0.8000000 1.0000000


C 0
S 5 1.00
1238.4016938 0.54568832082D-02
186.29004992 0.40638409211D-01
42.251176346 0.18025593888
11.676557932 0.46315121755
3.5930506482 0.44087173314
S 1 1.00
0.40245147363 1.0000000
S 1 1.00
0.13090182668 1.0000000
P 3 1.00
9.4680970621 0.38387871728D-01
2.0103545142 0.21117025112
0.54771004707 0.51328172114
P 1 1.00
0.15268613795 1.0000000
D 1 1.00
0.8000000 1.0000000


N 0
S 5 1.00
1712.8415853 -0.53934125305D-02
257.64812677 -0.40221581118D-01
58.458245853 -0.17931144990
16.198367905 -0.46376317823
5.0052600809 -0.44171422662
S 1 1.00
0.58731856571 1.0000000
S 1 1.00
0.18764592253 1.0000000
P 3 1.00
13.571470233 -0.40072398852D-01
2.9257372874 -0.21807045028
0.79927750754 -0.51294466049
P 1 1.00
0.21954348034 1.0000000
D 1 1.00
1.0000000 1.0000000


$end

In output:

       SF-DFT Excitation Energies              

(The first “excited” state might be the ground state)

Excited state 1: excitation energy (eV) = -3.6291
Total energy for state 1: -748.92605469 au
<S**2> : 0.0307
S( 2) → S( 1) amplitude = 0.9970 alpha

Excited state 2: excitation energy (eV) = -0.4204
Total energy for state 2: -748.80813471 au
<S**2> : 1.0117
S( 2) → S( 2) amplitude = 0.9745 alpha
S( 2) → V( 2) amplitude = -0.2003 alpha

Excited state 3: excitation energy (eV) = 0.3563
Total energy for state 3: -748.77959426 au
<S**2> : 1.8656
S( 1) → S( 1) amplitude = 0.4704 alpha
S( 2) → V( 2) amplitude = 0.6052 alpha
S( 2) → V( 3) amplitude = 0.6222 alpha

Excited state 4: excitation energy (eV) = 0.4253
Total energy for state 4: -748.77705752 au
<S**2> : 1.0206
S( 2) → V( 1) amplitude = 0.9611 alpha
S( 2) → V( 5) amplitude = 0.2472 alpha

And one more query I am considering 2nd excited state as singlet excited state. Am I right?

First excited state is the nominal S0 so 2nd excited state is nominal S1. The spin contamination is inherent to the method. You can use spin-adapted (SA-) SF-TDDFT to get spin pure states but there are no gradients so you can’t really optimize structures.

Hi John,
Thanks for your suggestion.
I have read this article Spin-flip, tensor equation-of-motion configuration interaction with a density-functional correction: A spin-complete method for exploring excited-state potential energy surfaces | The Journal of Chemical Physics | AIP Publishing. So, from there, I tried to reproduce the results for uracil nucleobase using SA-SF-DFT and SF-TDDFT.
For the collinear SF-TDDFT, I was able to reproduce it. But for SA-SF-DFT, I am new to this method.
What I did was first, I optimized the structure at the ground state (0 1) using B3LYP/6-311G(2df,2dp) level, and by taking the optimized coordinates, I made the input for SA-SF-DFT
$molecule
0 3
C 1.2909855 -0.3486576 0.0000000
C 1.1950303 1.1045135 0.0000000
C -0.0105753 1.6959551 0.0000000
N -1.1747334 0.9725309 0.0000000
C -1.2167048 -0.4176556 0.0000000
N 0.0379583 -0.9926528 0.0000000
O 2.3126711 -0.9987455 0.0000000
H 2.1093084 1.6735134 0.0000000
H -0.1324201 2.7697293 0.0000000
O -2.2523923 -1.0402408 0.0000000
H -2.0727890 1.4262136 0.0000000
H 0.0566031 -2.0025527 0.0000000
$end

$rem
BASIS = aug-cc-pVTZ
EXCHANGE = BHHLYP
GUI = 2
JOB_TYPE = SP
UNRESTRICTED = false
CIS_N_ROOTS = 5
SASF_RPA = 1
CIS_TRIPLETS = false
$end
and in the output:

Excited state 1: excitation energy (eV) = -4.7754
Total energy for state 1: -414.755109801337
Singlet
D( 20) → S( 1) amplitude = -0.1751
D( 22) → S( 1) amplitude = -0.1775
D( 26) → S( 1) amplitude = 0.1536
D( 28) → S( 1) amplitude = 0.1929
S( 2) → S( 1) amplitude = 0.9002
Excited state 2: excitation energy (eV) = 1.6620
Total energy for state 2: -414.518540128266
Singlet
D( 20) → S( 1) amplitude = -0.1962
D( 22) → S( 1) amplitude = -0.2252
D( 26) → S( 1) amplitude = 0.3296
D( 28) → S( 1) amplitude = 0.8192
Excited state 3: excitation energy (eV) = 2.2844
Total energy for state 3: -414.495669241921
Singlet
D( 22) → S( 1) amplitude = -0.1720
D( 26) → S( 1) amplitude = 0.5466
D( 28) → S( 1) amplitude = -0.4313
S( 2) → V( 1) amplitude = -0.6032
Excited state 4: excitation energy (eV) = 2.3906
Total energy for state 4: -414.491767038216
Singlet
D( 22) → S( 1) amplitude = -0.2026
D( 26) → S( 1) amplitude = 0.5495
D( 28) → S( 1) amplitude = -0.1991
S( 2) → V( 1) amplitude = 0.7029
Excited state 5: excitation energy (eV) = 2.8136
Total energy for state 5: -414.476221158729
Singlet
S( 2) → V( 2) amplitude = -0.9216
S( 2) → V( 8) amplitude = -0.2103
S( 2) → V( 13) amplitude = 0.1663

In article results for uracil/SA-SF-DFT → S1(5.35 eV) S2(5.56eV)
can you please help me where I am making a mistake (in input making or output interpretation)

@jherbert As I explained in my previous query. Can you please help me out with this?