Orbital rotations

Reference state parameterization

In time-dependent Hartree–Fock, the variations of orbitals in the time-independent reference state, \(| 0 \rangle = | \Psi_\mathrm{HF} \rangle\), can be introduced by means of an exponential parameterization

\[ | \bar{\Psi}(t) \rangle = e^{-i\hat{\kappa}(t)} | 0 \rangle \]

where

\[ \hat{\kappa}(t) = \sum_{a}^\mathrm{virt} \sum_{i}^\mathrm{occ} \Big[ \kappa_{ai}(t) \, \hat{a}^{\dagger}_a \hat{a}_i + \kappa_{ai}^*(t) \, \hat{a}^{\dagger}_i \hat{a}_a \Big] \]

It is convenient to organize the electron transfer amplitudes into the virtual–occupied (vo) and occupied–virtual (ov) blocks of a molecular orbital (MO) matrix. In the vo- and ov-blocks, we put amplitudes \(\kappa_{ai}\) and \(\kappa_{ai}^\ast \equiv \kappa_{ia}\), respectively.

\[\begin{split} \boldsymbol{\kappa} = \begin{bmatrix} 0 & \kappa_{ia} \\ \kappa_{ai} & 0 \\ \end{bmatrix} \end{split}\]

The \(\kappa\)-matrix is Hermitian and so is the \(\hat{\kappa}\)-operator.

Generator of orbital rotations

For an \(N\)-electron system, we have

\[ | 0 \rangle = \prod_{i=1}^N \hat{a}^{\dagger}_i | \mathrm{vac} \rangle \]

Since

\[ e^{i\hat{\kappa}}|\mathrm{vac}\rangle = \left(1+i\hat{\kappa}-\frac{1}{2}\hat{\kappa}^2+\ldots\right)|\mathrm{vac}\rangle = |\mathrm{vac}\rangle \]

we have

\[\begin{align*} | \bar{\Psi}(t) \rangle & = e^{-i\hat{\kappa}}\hat{a}^{\dagger}_1e^{i\hat{\kappa}}e^{-i\hat{\kappa}}\hat{a}^{\dagger}_2e^{i\hat{\kappa}}\ldots e^{-i\hat{\kappa}}\hat{a}^{\dagger}_Ne^{i\hat{\kappa}}|\mathrm{vac}\rangle \\ & = \hat{\tilde{a}}^{\dagger}_1\hat{\tilde{a}}^{\dagger}_2 \ldots \hat{\tilde{a}}^{\dagger}_N |\mathrm{vac}\rangle , \end{align*}\]

where we have introduced time-transformed creation operators

\[ \hat{\tilde{a}}^{\dagger}_p = e^{-i\hat{\kappa}}\hat{a}^{\dagger}_pe^{i\hat{\kappa}} = \hat{a}^{\dagger}_p -i\left[\hat{\kappa},\hat{a}^{\dagger}_p\right] -\frac{1}{2}\left[\hat{\kappa},\left[\hat{\kappa},\hat{a}^{\dagger}_p\right]\right]+\ldots \]

From the algebra of creation and annihilation operators, we find that

\[ \left[\hat{\kappa},\hat{a}^{\dagger}_p\right]=\hat{a}^{\dagger}_r\kappa_{rp} \]

and therefore

\[ \hat{\tilde{a}}^{\dagger}_p = \hat{a}^{\dagger}_{r}\left(\delta_{rp}-i\kappa_{rp}-\frac{1}{2}\kappa^2_{rp}+\ldots\right) = \hat{a}^{\dagger}_{r} \, U_{rp} \]

where we have introduced the unitary matrix

\[ \mathbf{U} = e^{-i\boldsymbol{\kappa}} \]

Following the determinant property of linearity, this transformation of creation operators translates directly into a transformation of orbitals in the wave function

\[ | \bar{\Psi}(t) \rangle = | \tilde{\psi}_1, \tilde{\psi}_2, \ldots, \tilde{\psi}_N \rangle \]

where

\[ \tilde{\psi}_p(\mathbf{r}) = \psi_r(\mathbf{r}) \, U_{rp} \]

This transformation preserves orthonormality among the orbitals and the \(\hat{\kappa}\)-operator can be thought of as a generator of rotations in the space of molecular orbitals.

Phase isolation

The phase is preserved in the orbital transformation in the sense that the overall phase of the reference state does not change, or in terms of orbitals

\[ \langle \psi_p | \tilde{\psi}_p \rangle \in {\cal R} \]

It can therefore be used for the parameterization of phase-isolated wave functions as indicated with the overbar in \(| \bar{\Psi}(t) \rangle\).

Spin symmetry

If we wish to preserve the singlet spin symmetry of a closed-shell reference state, we use spin-adapted electron-transfer operators

\[ \hat{\kappa}(t) = \sum_{a}^\mathrm{virt} \sum_{i}^\mathrm{occ} \Big[ \kappa_{ai}(t) \, \hat{E}^{\dagger}_{ai} + \kappa_{ai}^*(t) \, \hat{E}^{\dagger}_{ia} \Big] \]

where

\[ \hat{E}_{pq}^\dagger = \hat{a}^\dagger_{p\alpha} \hat{a}_{q\alpha} + \hat{a}^\dagger_{p\beta} \hat{a}_{q\beta} \]

Summations here run over molecular orbitals (MOs) instead of spin orbitals, and consequently, the number of parameters is reduced by a factor of four.

We get

\[ | \bar{\Psi}(t) \rangle = | \tilde{\psi}_1, \tilde{\psi}_\bar{1}, \ldots, \tilde{\psi}_{N/2}, \tilde{\psi}_\overline{N/2} \rangle \]

where

\[\begin{split} \tilde{\psi}_p(\mathbf{r}) = \tilde{\phi}_p(\mathbf{r}) \begin{pmatrix} 1 \\ 0 \end{pmatrix} ; \qquad \tilde{\psi}_\bar{p}(\mathbf{r}) = \tilde{\phi}_p(\mathbf{r}) \begin{pmatrix} 0 \\ 1 \end{pmatrix} \end{split}\]

with

\[ \tilde{\phi}_p(\mathbf{r}) = \phi_r(\mathbf{r}) \, \big[e^{-i \kappa}\big]_{rp} \]

Illustrations

import numpy as np
from scipy.linalg import expm

Let us consider a system of two-electrons in two MOs (or four spin orbitals).

Multi-electron excited determinants

The generator of rotations includes only single electron excitation operators but the range of the transformation includes all determinants, also multi-electron excited determinants.

We reach the doubly excited determinant with

\[\begin{split} \boldsymbol{\kappa} = \begin{bmatrix} 0 & 0 & -i\pi/2 & 0 \\ 0 & 0 & 0 & -i\pi/2 \\ i\pi/2 & 0 & 0 & 0 \\ 0 & i\pi/2 & 0 & 0 \\ \end{bmatrix} \end{split}\]

k31, k32, k41, k42 = (np.pi * 0.5j, 0, 0, np.pi * 0.5j)
k13, k23, k14, k24 = np.conjugate((k31, k32, k41, k42))

kappa = np.array(
    [[0, 0, k13, k14], [0, 0, k23, k24], [k31, k32, 0, 0], [k41, k42, 0, 0]]
)

U = expm(-1j * kappa)
print("Occupied orbitals:\n", U.real[:, :2])

Occupied orbitals:
 [[0. 0.]
 [0. 0.]
 [1. 0.]
 [0. 1.]]

Phase isolation

Let us create a random generator of orbital rotations on the form

\[\begin{split} \boldsymbol{\kappa} = \begin{bmatrix} 0 & \kappa_{ia} \\ \kappa_{ai} & 0 \\ \end{bmatrix} \end{split}\]

np.set_printoptions(precision=6, suppress=True, linewidth=170)
np.random.seed(20240208)

kappa = np.zeros((4, 4), dtype=complex)
kappa[2:, :2] = np.random.rand(2, 2) + 1j * np.random.rand(2, 2)
kappa = kappa + np.conjugate(kappa).T

The transformed orbitals become

U = expm(-1j * kappa)
print("Transformed orbitals:\n", U)

Transformed orbitals:
 [[ 0.414314-0.j       -0.272098-0.000049j -0.063249-0.357053j -0.609038-0.501889j]
 [-0.272098+0.000049j  0.784682+0.j       -0.400653-0.136519j -0.138307-0.334586j]
 [ 0.063249-0.357053j  0.400653-0.136519j  0.78434 +0.j       -0.265899+0.058835j]
 [ 0.609038-0.501889j  0.138307-0.334586j -0.265899-0.058835j  0.414656-0.j      ]]

We note that the diagonal elements are real, i.e.

\[ \langle \psi_p | \tilde{\psi}_p \rangle \in {\cal R} \]

Spin symmetry

The singlet spin symmetry of the closed-shell reference state, \(|\psi_1, \psi_2 \rangle\), is preserved by requiring that \(\kappa_{42} = \kappa_{31}\) and \(\kappa_{41} = \kappa_{32} = 0\).

As an example, let us consider

\[\begin{split} \boldsymbol{\kappa} = \begin{bmatrix} 0 & 0 & -i\pi/4 & 0 \\ 0 & 0 & 0 & -i\pi/4 \\ i\pi/4 & 0 & 0 & 0 \\ 0 & i\pi/4 & 0 & 0 \\ \end{bmatrix} \end{split}\]

k31, k32, k41, k42 = (np.pi * 0.25j, 0, 0, np.pi * 0.25j)
k13, k23, k14, k24 = np.conjugate((k31, 0, 0, k42))

kappa = np.array(
    [[0, 0, k13, k14], [0, 0, k23, k24], [k31, k32, 0, 0], [k41, k42, 0, 0]]
)

U = expm(-1j * kappa)
print("Occupied orbitals:\n", U.real[:, :2])

Occupied orbitals:
 [[0.707107 0.      ]
 [0.       0.707107]
 [0.707107 0.      ]
 [0.       0.707107]]

The resulting state can be written

\[\begin{align*} | \bar{\Psi} \rangle & = | \tilde{\psi_1}, \tilde{\psi_2} \rangle = | (\psi_1 + \psi_3)/\sqrt{2}, (\psi_2 + \psi_4)/\sqrt{2} \rangle \\ & = \frac{1}{2} \Big[ | \psi_1, \psi_2 \rangle + | \psi_1, \psi_4 \rangle + | \psi_3, \psi_2 \rangle + | \psi_3, \psi_4 \rangle \Big] \end{align*}\]

where we identify the reference state and the doubly excited determinant alongside the singly excited spin-adapted configuration of singlet spin symmetry.

Finally, let us consider a somewhat larger system with four electrons in six MOs (or 12 spin orbitals). In the spin-orbital basis, it is convenient to organize parameters in spin blocks.

\[\begin{split} \boldsymbol{\kappa} = \begin{bmatrix} 0 & 0 & \kappa_{ia} & \kappa_{i\bar{a}} \\ 0 & 0 & \kappa_{\bar{i}a} & \kappa_{\bar{i}\bar{a}} \\ \kappa_{ai} & \kappa_{a\bar{i}} & 0 & 0\\ \kappa_{\bar{a}i} & \kappa_{\bar{a}\bar{i}} & 0 & 0\\ \end{bmatrix} \end{split}\]

Let us perform a random transformation of the spin orbitals and print out \(\tilde{\psi}_1\) and \(\tilde{\psi}_\bar{1}\) as examples.

k31, k41, k51, k61, k32, k42, k52, k62 = np.random.rand(8) + 1j * np.random.rand(8)
k13, k14, k15, k16, k23, k24, k25, k26 = np.conjugate(
    (k31, k41, k51, k61, k32, k42, k52, k62)
)

kappa = np.array(
    [
        [0, 0, 0, 0, k13, k14, k15, k16, 0, 0, 0, 0],
        [0, 0, 0, 0, k23, k24, k25, k26, 0, 0, 0, 0],
        [0, 0, 0, 0, 0, 0, 0, 0, k13, k14, k15, k16],
        [0, 0, 0, 0, 0, 0, 0, 0, k23, k24, k25, k26],
        [k31, k32, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
        [k41, k42, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
        [k51, k52, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
        [k61, k62, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, k31, k32, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, k41, k42, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, k51, k52, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, k61, k62, 0, 0, 0, 0, 0, 0, 0, 0],
    ]
)

U = expm(-1j * kappa)
print("Lowest occupied alpha-spin orbital:\n", U[:,0])
print("Lowest occupied beta-spin orbital:\n", U[:,2])

Lowest occupied alpha-spin orbital:
 [ 0.195145+0.j       -0.674738-0.007566j  0.      +0.j        0.      -0.j        0.272185-0.364692j  0.150157-0.219078j -0.077445-0.261337j  0.390656+0.045347j
  0.      +0.j        0.      +0.j        0.      +0.j        0.      +0.j      ]
Lowest occupied beta-spin orbital:
 [ 0.      +0.j        0.      +0.j        0.195145+0.j       -0.674738-0.007566j  0.      +0.j        0.      +0.j        0.      +0.j        0.      +0.j
  0.272185-0.364692j  0.150157-0.219078j -0.077445-0.261337j  0.390656+0.045347j]

We note that the spatial parts of the two spin orbitals are identical.

If we instead implement the transformation in the MO-basis, the dimension of the \(\kappa\)-matrix becomes reduced by a factor of two.

kappa = np.array(
    [
        [0, 0, k13, k14, k15, k16],
        [0, 0, k23, k24, k25, k26],
        [k31, k32, 0, 0, 0, 0],
        [k41, k42, 0, 0, 0, 0],
        [k51, k52, 0, 0, 0, 0],
        [k61, k62, 0, 0, 0, 0],
    ]
)

U = expm(-1j * kappa)
print("Lowest occupied molecular orbital:\n", U[:,0])

Lowest occupied molecular orbital:
 [ 0.195145-0.j       -0.674738-0.007566j  0.272185-0.364692j  0.150157-0.219078j -0.077445-0.261337j  0.390656+0.045347j]

We note that this molecular orbital is identical to the spatial parts of the spin orbitals in the previous example.