Calculation of S1 energy in Gaussian

I am trying to calculate the energy of the S1 state with UDFT in Gaussian software. I understand that, generally, one should manipulate the charge and spin multiplicity to manually change the occupations.

After calculating the S0 state with the latters as '0 1', I wanted to compute the S1 using '1 1'. However, I received

The combination of multiplicity 1 and 219 electrons is impossible

How can I calculate the energy of S1 state in Gaussian with UDFT?

Thanks!

• Why not do a TD-DFT? The S1 state is an excited state if i understand correctly, so you would get the energy of S1 with TDDFT. Anyway, what you are wanting to do in OP seems doable by switching a pair of alpha orbitals (or beta) such as shown here gaussian.com/guess/?tabid=4. You have use the S0 single point wavefunction as guess and then alter that guess before doing the SCF. Also make sure to check the occupations and population analysis to see if you actually converged to S1 state or it went to something else because this procedure is quite finicky. Feb 23 at 15:09
• I was actually trying to follow what they reported here: pubs.acs.org/doi/full/10.1021/acs.jctc.2c00905 where they show that ROKS performs better. And, moreover, I am using a solvent and thus the ROKS/UDFT might perform better Feb 23 at 15:14
• In the paper, they are manipulating the initial guess in the way SR Maiti and keeping it from collapsing back to the ground state using Delta SCF (not sure how to apply Delta SCF in Gaussian). Your input is changing the charge, but the state your hoping to generate (S1) has the same charge/multiplicity as the S0 ground state.
– Tyberius
Feb 23 at 16:23

In general, computing singlet excited states with an unrestricted wavefunction/method is challenging. This is for a few reasons, but one of them is that the SCF procedure tends to converge towards the lowest energy state (because that's what it's designed to do, and normally what you want). For a closed shell molecule, the lowest energy singlet state is S0 (ie, the ground state), and I don't think you can easily manipulate the charge and multiplicity of your system to converge to a higher singlet state.

This is in contrast to triplet excited states where you'll see this method used much more commonly (where it is often called Delta-SCF). For a closed shell molecule, the lowest energy triplet state is T1, so a charge/multiplicity line like:

0 3 (Charge 0 (neutral), Multiplicity 3 (triplet))

Would tend to converge to T1. For S1, a simple input like this is not possible I'm afraid.

If you're trying to recreate the methodology in the paper you linked, you'll run into a few problems. First, you are using a different program (Gaussian) compared to the authors (QChem), so you sadly can't reuse any of their input files. Second, the authors developed part of the functionality that they are benchmarking themselves so this code is not available in all version of Q-Chem. The authors state it should be available from 5.4.1 onwards, but you would have to double-check:

All quantum-chemical calculations were carried out with a development version of the Q-Chem 5.4 program package.1 We have modified the original source code to enable ROKS optimizations in combination with the IEF-PCM solvent model. The respective modifications will be included in the next Q-Chem release (5.4.1).

Phys. Chem. Lett. 2021, 12, 35, 8470–8480 (Supporting information: https://pubs.acs.org/doi/suppl/10.1021/acs.jpclett.1c02299/suppl_file/jz1c02299_si_001.pdf)

Even then, it is clear from the methodology of the above paper that this is a fairly complicated function to implement, and probably not for the faint-hearted. See for example:

Singlet states are optimized with ROKS and UKS, the latter of which is in combination with the (initial) maximum-overlap method (MOM/IMOM) to avoid a variational collapse to the ground state. (73,75) In the UKS/IMOM calculations (in the following just UKS), the excited state must be explicitly targeted. This is realized through a manipulation of the initial guess orbitals in which one electron is moved from the highest occupied molecular orbital (HOMO) to lowest unoccupied molecular orbital (LUMO). ROKS optimizations do not need such a guess for the excitation vector. They reliably converge on the lowest excited singlet state, whereas the UKS optimizations occasionally converge on higher-lying solutions (ca. 30% of all cases).

Phys. Chem. Lett. 2021, 12, 35, 8470–8480

By contrast, calculating S1 with TD-DFT as suggested by S R Maiti is super easy, have a look at https://gaussian.com/td/. The only gotcha is to calculate more excited states than you need (even if you're only interested in S1) because of state reordering. TD-DFT is compatible with SCRF for solvent inclusion, and state specific effects can be included with SCRF=(CorrectedLinearResponse). Even if the accuracy of the TD-DFT is lower than other methods, you may find it is still suitably accurate for solving your problem, which is all that matters.

• Thanks Leeman for the very comprehensive reply! Actually, for the calculation of S1 with solvent, I was following this: gaussian.com/scrf (In the example, steps 4,5). Is this what you mean? Mar 3 at 19:01
• Hi Laura, yes that's a good reference. The SCRF=(CorrectedLinearResponse) I mentioned is covered in steps 3 and 5 (where they use the shortened form CorrectedLR) Mar 3 at 19:48