Excited states#

For the analysis of the excited states, it is common to simply study the molecular orbitals of the final states, and base assignments and other conclusions from this. This often works quite well, but the final state can be multi-determinational in nature, and the canonical orbitals of larger systems are often quite delocalized and hard to study. As such, many other analysis tools have been developed, two of which will be discussed here:

As an illustration, we consider the first valence-excited state of water, as calculated using ADC(2):

water_geom = """
O       0.0000000000     0.1187290000     0.0000000000
H      -0.7532010000    -0.4749160000    -0.0000000000
H       0.7532010000    -0.4749160000     0.0000000000
"""

mol = gto.Mole()
mol.atom = water_geom
mol.basis = "6-31G"
mol.build()

scf_gs = scf.HF(mol)
scf_gs.kernel()

adc_state = adcc.adc2(scf_gs, n_singlets=2)

Printing the results:

# Print results
print(adc_state.describe())

# Print dominant amplitudes
print(adc_state.describe_amplitudes())
print_results = """+--------------------------------------------------------------+
| adc2                                    singlet ,  converged |
+--------------------------------------------------------------+
|  #        excitation energy     osc str    |v1|^2    |v2|^2  |
|          (au)           (eV)                                 |
|  0     0.3054352      8.311314   0.0142    0.9591   0.04089  |
|  1     0.3903211      10.62118   0.0000    0.9636   0.03637  |
+--------------------------------------------------------------+

+---------------------------------------------------+
| State   0 ,     0.3054352 au,      8.311314 eV    |
+---------------------------------------------------+
| HOMO          -> LUMO           a ->a      -0.691 |
| HOMO          -> LUMO+4         a ->a     +0.0451 |
| HOMO          -> LUMO+7         a ->a     +0.0227 |
| HOMO-2 HOMO   -> LUMO+1 LUMO    ab->ab    +0.0306 |
| HOMO   HOMO   -> LUMO   LUMO+3  ab->ab    -0.0253 |
| HOMO-2 HOMO   -> LUMO   LUMO+1  aa->aa    -0.0191 |
| HOMO-1 HOMO   -> LUMO   LUMO    ab->ab    -0.0191 |
| HOMO-3 HOMO   -> LUMO   LUMO    ab->ab    +0.0182 |
| HOMO-1 HOMO   -> LUMO   LUMO+4  aa->aa    -0.0174 |
| HOMO-1 HOMO   -> LUMO+4 LUMO    ab->ab    +0.0152 |
| HOMO-2 HOMO   -> LUMO+6 LUMO    ab->ab    -0.0152 |
| HOMO-1 HOMO   -> LUMO+5 LUMO    ab->ab    -0.0147 |
| HOMO-2 HOMO   -> LUMO   LUMO+6  aa->aa    +0.0139 |
| HOMO-1 HOMO   -> LUMO   LUMO+5  aa->aa    +0.0117 |
| HOMO-2 HOMO   -> LUMO   LUMO+1  ab->ab    +0.0115 |
| HOMO-3 HOMO   -> LUMO   LUMO+7  aa->aa    +0.0107 |

+---------------------------------------------------+
| State   1 ,     0.3903211 au,      10.62118 eV    |
+---------------------------------------------------+
| HOMO          -> LUMO+1         a ->a       +0.69 |
| HOMO          -> LUMO+6         a ->a     -0.0598 |
| HOMO          -> LUMO+2         a ->a     -0.0431 |
| HOMO-2 HOMO   -> LUMO+1 LUMO+1  ab->ab    -0.0301 |
| HOMO   HOMO   -> LUMO+1 LUMO+3  ab->ab    +0.0234 |
| HOMO-1 HOMO   -> LUMO   LUMO+1  ab->ab     +0.019 |
| HOMO-2 HOMO   -> LUMO+1 LUMO+6  aa->aa    -0.0189 |
| HOMO-1 HOMO   -> LUMO+4 LUMO+1  ab->ab    -0.0144 |
| HOMO-3 HOMO   -> LUMO   LUMO+1  ab->ab    -0.0143 |
| HOMO-1 HOMO   -> LUMO+1 LUMO+4  aa->aa    +0.0133 |
| HOMO-2 HOMO   -> LUMO+6 LUMO+1  ab->ab     +0.013 |
| HOMO-1 HOMO   -> LUMO+5 LUMO+1  ab->ab    +0.0121 |
| HOMO-1 HOMO   -> LUMO+1 LUMO+5  aa->aa    -0.0109 |
| HOMO-1 HOMO   -> LUMO   LUMO+1  aa->aa    +0.0104 |"""
print(print_results)
+--------------------------------------------------------------+
| adc2                                    singlet ,  converged |
+--------------------------------------------------------------+
|  #        excitation energy     osc str    |v1|^2    |v2|^2  |
|          (au)           (eV)                                 |
|  0     0.3054352      8.311314   0.0142    0.9591   0.04089  |
|  1     0.3903211      10.62118   0.0000    0.9636   0.03637  |
+--------------------------------------------------------------+

+---------------------------------------------------+
| State   0 ,     0.3054352 au,      8.311314 eV    |
+---------------------------------------------------+
| HOMO          -> LUMO           a ->a      -0.691 |
| HOMO          -> LUMO+4         a ->a     +0.0451 |
| HOMO          -> LUMO+7         a ->a     +0.0227 |
| HOMO-2 HOMO   -> LUMO+1 LUMO    ab->ab    +0.0306 |
| HOMO   HOMO   -> LUMO   LUMO+3  ab->ab    -0.0253 |
| HOMO-2 HOMO   -> LUMO   LUMO+1  aa->aa    -0.0191 |
| HOMO-1 HOMO   -> LUMO   LUMO    ab->ab    -0.0191 |
| HOMO-3 HOMO   -> LUMO   LUMO    ab->ab    +0.0182 |
| HOMO-1 HOMO   -> LUMO   LUMO+4  aa->aa    -0.0174 |
| HOMO-1 HOMO   -> LUMO+4 LUMO    ab->ab    +0.0152 |
| HOMO-2 HOMO   -> LUMO+6 LUMO    ab->ab    -0.0152 |
| HOMO-1 HOMO   -> LUMO+5 LUMO    ab->ab    -0.0147 |
| HOMO-2 HOMO   -> LUMO   LUMO+6  aa->aa    +0.0139 |
| HOMO-1 HOMO   -> LUMO   LUMO+5  aa->aa    +0.0117 |
| HOMO-2 HOMO   -> LUMO   LUMO+1  ab->ab    +0.0115 |
| HOMO-3 HOMO   -> LUMO   LUMO+7  aa->aa    +0.0107 |

+---------------------------------------------------+
| State   1 ,     0.3903211 au,      10.62118 eV    |
+---------------------------------------------------+
| HOMO          -> LUMO+1         a ->a       +0.69 |
| HOMO          -> LUMO+6         a ->a     -0.0598 |
| HOMO          -> LUMO+2         a ->a     -0.0431 |
| HOMO-2 HOMO   -> LUMO+1 LUMO+1  ab->ab    -0.0301 |
| HOMO   HOMO   -> LUMO+1 LUMO+3  ab->ab    +0.0234 |
| HOMO-1 HOMO   -> LUMO   LUMO+1  ab->ab     +0.019 |
| HOMO-2 HOMO   -> LUMO+1 LUMO+6  aa->aa    -0.0189 |
| HOMO-1 HOMO   -> LUMO+4 LUMO+1  ab->ab    -0.0144 |
| HOMO-3 HOMO   -> LUMO   LUMO+1  ab->ab    -0.0143 |
| HOMO-1 HOMO   -> LUMO+1 LUMO+4  aa->aa    +0.0133 |
| HOMO-2 HOMO   -> LUMO+6 LUMO+1  ab->ab     +0.013 |
| HOMO-1 HOMO   -> LUMO+5 LUMO+1  ab->ab    +0.0121 |
| HOMO-1 HOMO   -> LUMO+1 LUMO+5  aa->aa    -0.0109 |
| HOMO-1 HOMO   -> LUMO   LUMO+1  aa->aa    +0.0104 |

Here we print the basic description of the excited states (energy, intensity, weight of double excitation amplitudes), as well as a break-down of the dominating amplitudes for each transition. From this we can start analysing the excited state by looking at the canonical MOs, but we will now use the results for a more involved analysis.

Natural transition orbitals#

Natural transition orbitals (NTOs) are constructed to provide the most compact, transition-dependent transition orbitals of a specific excitation. With this, a single pair of NTOs corresponding to the hole and electron will typically dominate, and will thus provide a easily interpretable description of the excitation. For a pure HOMO to LUMO transition the NTOs would be the same as the HOMO (hole) and LUMO (electron) orbitals.

The NTOs are constructed from a singular-value decomposition (SVD) of the transition density matrix (\(\mathbf{T}\)):

\[ \mathbf{UTV}^{\dagger} = \mathbf{\Lambda} \]

where \(\mathbf{U}\) and \(\mathbf{V}\) are the transformation matrices correspoding to the hole and electron, respectively, and \(\mathbf{\Lambda}\) is a diagonal matrix measuring the relative importance of each pair of NTOs.

The first pair of NTOs for the first excitation can now be constructed as:

Note

The following scripts are not executed, and static image are instead loaded. This is merely for improved loading time of the webpage and to save disk space.

i = 0 # first state only

# Load transition density matrix, combine alpha and beta, and transform to numpy
tdm_ao = adc_state.transition_dm[i].to_ao_basis()
p_tdm_tot = (tdm_ao[0] + tdm_ao[1]).to_ndarray()

# Build NTOs by singular value decomposition
u, s, v = np.linalg.svd(p_tdm_tot)

print('Dominant NTO of state {}'.format(i+1))
print('Relative importance:',np.around(s[0]/sum(s),3))
k = 0 # only look at dominant NTO pair  

# Initial
tools.cubegen.orbital(mol=mol, coeff=v[k],outfile="../../img/visualize/water_nto_{}_HONTO+{}.cube".format(i,k))
# Final
tools.cubegen.orbital(mol=mol, coeff=u.T[k],outfile="../../img/visualize/water_nto_{}_LUNTO+{}.cube".format(i,k))
print_results = """
Dominant NTO of state 1
Relative importance: 0.97
"""
print(print_results)
Dominant NTO of state 1
Relative importance: 0.97
../../_images/water-ntos-dens.png

This shows a strong dominance (0.97) of the first NTO pair.

Attachment and detachment densities#

A different method for visualizing the transitions is to consider the attachment (A) and detachment (D) densities, which are constructed to show the density change related to an excitation. For a simple HOMO-LUMO transition they would simply correspond to the square of the dominant NTOs, with meaning:

  • Hole/detachment density representing where electrons go from

  • Electron/attachment density representing where electrons go to

These densities can then be used to look at properties such as hole and electron size, distance between the centroid (and thus level of charge-transfer), and more. D/A densities are constructed from the one-particle difference density matrix (1DDM), which is simply the difference between the on-particle density matrices of the initial and final state:

\[ \rho_{\Delta} = \rho_f - \rho_i \]

Diagonalizing this

\[ \mathbf{U} \rho_{\Delta} \mathbf{U}^{\dagger} = \delta \]

The attachment and detachment densities are constructed by considering the negative and positive eigenvalues.

Looking at the attachment and detachment density of the first excited state:

# Load transition density matrix, combine alpha and beta, and transform to numpy
p_state = adc_state.state_diffdm[i].to_ao_basis()
p_state_ao = (p_state[0] + p_state[1]).to_ndarray()

# Diagonalize the 1DDM
k, w = np.linalg.eigh(p_state_ao)
k_detach = k.copy()
k_attach = k.copy()
# Detachment: set positive eigenvalues to 0
k_detach[k > 0] = 0
# Attachment: set negative eigenvalues to 0
k_attach[k < 0] = 0
# Back-transform with numpy
detach_ao = w @ np.diag(k_detach) @ w.T
attach_ao = w @ np.diag(k_attach) @ w.T

# Write cube-files
de = tools.cubegen.density(
    mol, dm=detach_ao, outfile="../../img/visualize/water_detachment_{}.cube".format(i)
)
at = tools.cubegen.density(
    mol, dm=attach_ao, outfile="../../img/visualize/water_attachment_{}.cube".format(i)
)
../../_images/water-att-det-dens.png

We see that these states are more compact (as we are looking on densities) than above NTOs. For a discussion on NTOs and D/A densities, see, e.g., this paper.