Hi all,
I have a question regarding SA-SF-RPA calculations. In the output, how can we determine if the structure is a closed-shell or open-shell singlet? Is it possible to use this method to remove triplet contamination from an open-shell singlet?
input:
$molecule
0 3
O 7.65182000 10.40523900 9.59526800
O 9.08331400 11.35464700 10.98497800
N 9.84934700 9.93588700 9.22907900
C 9.86600300 8.21145200 7.47539200
C 9.72056600 8.60711600 8.86718800
C 9.35917200 7.71544800 9.95611100
C 8.78735900 10.53490600 9.97252600
C 10.10609300 11.06224800 11.95269700
H 10.00788200 10.01350000 12.25287900
H 9.87644800 11.71667600 12.79970300
C 11.49390100 11.35164100 11.40288700
H 11.56785400 12.41934000 11.15061300
H 12.22245500 11.13080600 12.19333900
C 11.72828000 10.48692500 10.23009800
H 12.14595600 9.49156400 10.39089300
C 11.11542000 10.72972900 8.92671000
H 11.62518400 10.16906300 8.13776900
C 10.85920800 12.14810300 8.46625100
H 10.19801300 12.10651700 7.58907700
H 10.32402800 12.71993000 9.23636800
C 12.16833100 12.83335000 8.10882700
H 11.99111700 13.85589400 7.75575300
H 12.69289000 12.28924100 7.31097300
H 12.84555400 12.89004900 8.97249700
O 9.22074000 8.17504700 11.09822900
O 9.95343100 7.03135300 7.13251000
C 9.14042400 6.25098800 9.70869700
H 10.06014100 5.76695500 9.36021100
H 8.80736900 5.78840300 10.64405300
H 8.39397300 6.08476900 8.92285700
C 9.82198000 9.29398400 6.42036200
H 9.14452300 10.11302300 6.69274000
H 10.81905700 9.72320300 6.24240300
H 9.49177700 8.84058600 5.47961000
$end
$rem
unrestricted 0
jobtype sp
basis def2-tzvp
exchange M06
cis_n_roots 5
sasf_cis 1
cis_triplets false
SOLVENT_METHOD smd
MEM_TOTAL 190000
MAX_CIS_CYCLES 1000
MAX_SCF_CYCLES = 1000
SYMMETRY false
$end
$smx
solvent dichloromethane
$end
output:
Excited state 1: excitation energy (eV) = -0.8078
Total energy for state 1: -822.882776205868
Singlet
D( 61) --> S( 1) amplitude = 0.2158
S( 1) --> S( 1) amplitude = 0.1660
S( 2) --> S( 1) amplitude = -0.9086
S( 2) --> S( 2) amplitude = -0.1660
Excited state 2: excitation energy (eV) = 0.8960
Total energy for state 2: -822.820161572521
Singlet
D( 61) --> S( 1) amplitude = -0.3629
D( 63) --> S( 1) amplitude = -0.2552
S( 1) --> S( 1) amplitude = 0.5660
S( 2) --> S( 2) amplitude = -0.5660
Excited state 3: excitation energy (eV) = 1.3404
Total energy for state 3: -822.803832125373
Singlet
D( 62) --> V( 1) amplitude = 0.1602
D( 63) --> S( 1) amplitude = 0.9329
Excited state 4: excitation energy (eV) = 1.6967
Total energy for state 4: -822.790738412628
Singlet
D( 63) --> V( 1) amplitude = -0.1958
D( 62) --> S( 1) amplitude = -0.9287
Excited state 5: excitation energy (eV) = 2.1377
Total energy for state 5: -822.774532163293
Singlet
D( 60) --> S( 1) amplitude = 0.4403
D( 61) --> S( 1) amplitude = -0.7907
D( 61) --> S( 2) amplitude = 0.1528
S( 2) --> S( 1) amplitude = -0.1597
Target state 1 is predominantly a closed-shell singlet since the largest SF amplitude is S(2)alpha to S(1)beta. Target state 2 is an open-shell singlet as the leading amplitudes correspond to open-shell configurations. Another closed-shell state is possible that will have a dominant contribution from S(1)alpha to S(2)beta transition, but it is probably higher lying than the states you have computed. The rest of the states you posted are all described by open-shell configurations involving orbitals outside the singly occupied orbital space of the initial reference configuration, so they have predominant open-shell character.
SASF-CIS gives spin-pure singlet states. There is no spin contamination in these states, i.e., <S^2> = 0. If your second query relates to comparing SASF-CIS and SF-CIS results, addition of the minimal set of missing spin-complement configurations in regular SF-CIS could potentially result in different energies; so, comparisons between the results with the two methods should be done carefully.
Hi Kaushik, Thank you very much for your helpful explanation. I am a bit confused. I used Gaussian 16 to optimize this structure, and my results show that the open-shell structure is more stable. I intended to use single-point calculations on Gaussian-optimized structure to remove triplet contamination in the open-shell singlet. However, based on your explanation, the results here show that, contrary to my Gaussian results, the closed-shell structure is more stable. Why is there this difference? Additionally, is it possible to draw the orbitals for transitions using SASF method?
NTO analysis for SASF-CIS transitions cannot be done as of Q-Chem v6.2.1, but you can consider visualizing the relevant MOs through an FCHK file.
We cannot comment on any other software besides Q-Chem. Plus, your description of the issue is not very clear to me, so please elaborate what you did and what methods/states/geometries are you comparing/optimizing. When you say, you want to remove triplet contamination, what do you mean?
In fact I optimized a molecule using UM06/def2-SVP/SMD/dichloromethane, and the results showed that the open-shell configuration is more stable than the closed-shell geometry. Since spin contamination can affect the energy of the open-shell configuration, I decided to recalculate the single-point energy using a method without spin contamination to see the exact energy difference between the open-shell and closed-shell singlets. So, I used SASF-CIS/M06/def2-SVP/SMD/dichloromethane. However, the results from this method indicate that the closed-shell configuration is more stable. Now, I’m wondering why one method suggests the open-shell singlet is more stable, while another indicates the closed-shell singlet is more stable.
Another problem is that when I try to create fchk for SASF-cis method using Gui = 2, the calculation does not proceed, there is an error:
Checking the input file for inconsistencies...
ROCIS results cannot be written to the formatted checkpoint file.
Please remove IQMOL_FCHK or GUI=2 from $rem.
Thank you for your help. The issue with this explanation is that the triplet state is even higher in energy than the closed-shell singlet. Therefore, triplet spin contamination should result in a high-energy open-shell singlet. However, as you can see here, the um06 open-shell singlet is lower in energy than the closed-shell singlet.
I guess the question here is what is really the UM06 “open-shell singlet”? UKS is not generally able to produce singlets. What you get is some mixed spin-contaminated state.
You can see that the “open-shell singlet” is almost degenerate with the closed-shell singlet. So, I would guess that the UM06 result may indeed be the closed-shell singlet but just with a bit of spin-contamination on top that slightly lowers its energy. You can check that in terms of the orbitals produced.
Thank you all for your help. As Felix mentioned, perhaps the UM06 open-shell singlet is not a real state and is actually a closed-shell singlet with triplet spin contamination.
Regarding John’s comment on using OPSING = TRUE, I am not sure, but it seems to lead to an even more stable open-shell singlet, which is completely inconsistent with the spin-pure SASF-CIS results.
Maybe OPSING = TRUE is only useful when ⟨S²⟩ is very close to 1. In this case, it is below 0.5.
Sorry, what exactly is UM06? Unrestricted TDDFT with M06, right? I do not see spin-contaminated states with this approach with <S^2> = 0 or 2 for the excited states. So, how did you set up the UM06 calculation? When I tried the unrestricted TDDFT calculation, I do not get the SPE of the open-shell singlet (OSS) and triplet that you have posted.
I am also not so sure about the explanation. Static correlation in diradical systems mainly impacts the second closed-shell singlet (CSS2) as it has a dominant doubly excited configuration with respect to the dominant configuration of the first closed-shell singlet (CSS1); SF is a Fock-space trick to address the poor description of the CSS2 wave function within a single-reference setup. Since you are not concerned with CSS2, I was going to suggest that you also run your calculation without SF and check what the ordering is, but it seems you already did this. Correct? If the energy gap between the OSS and CSS1 is similar with SF and non-SF TDDFT method for your system (it is not necessarily the case in general), then that can also possibly help you identify the ordering.
If you have not already done it and if valid for your system, another trick you can consider in your calculations is to use the dianionic or dicationic MOs as reference (if the ground state of the dication or dianion can be well described by a single determinant but not the neutral system) and compute your target SF or non-SF states using these MOs. Also, consider visualizing the orbitals involved in transitions from SF and non-SF calculations for like-for-like comparisons.
OSS has a multiconfigurational wave function—the two dominant configurations are spin complements of each other with each configuration made up of two half-filled orbitals with opposite spins. You cannot get the correct OSS wave function by computing a single-determinant SCF solution. The same would be the case for the m_s=0 triplet, but you can instead compute the energy of the m_s=1 triplet, which is well described by a single Slater determinant, which is what you did I think.
“UM06” is Gaussian notation. The grown-up way thing to call this is simply unrestricted M06. For the situation that you describe, spin purification (OP_SING = TRUE) should be a step in the right direction, or ROKS as Felix suggests. (For two unpaired electrons, what Q-Chem calls ROKS is the orbital-optimized version of approximate spin projection.) One issue that I see with this discussion is that you are trying to make a one-to-one comparison of M06 results (with various procedures for treating the open shells) results to SA-SF-CIS results. I would worry about how sensitive the answer is to the particular choice of functional (and with a Minnesota functional, also to the basis set).