5.1. Unrestricted Hartree-Fock method

5.1.1. Overview

The unrestricted Hartree-Fock approximation is a method to approximate two-body interactions into one-body terms by taking account of the fluctuation of the one-body operators up to first order. For a general two-body interactions, it leads to the following approximation:

(5.1)\[\begin{split}\begin{aligned} c_{i}^{\dagger} c_{j}^{\dagger} c_{k} c_{l} & \sim \langle c_{i}^{\dagger} c_l \rangle c_{j}^{\dagger} c_k + c_{i}^{\dagger} c_l \langle c_{j}^{\dagger} c_k\rangle - \langle c_{i}^{\dagger} c_k \rangle c_{j}^{\dagger} c_l - c_{i}^{\dagger} c_k \langle c_{j}^{\dagger} c_l \rangle \nonumber \\ & \qquad - (\langle c_{i}^{\dagger} c_l \rangle \langle c_{j}^{\dagger} c_k \rangle - \langle c_{i}^{\dagger} c_k\rangle \langle c_{j}^{\dagger} c_l\rangle) . \end{aligned}\end{split}\]

In H-wave, the two-body interaction terms are defined as

(5.2)\[\begin{split}\begin{aligned} \mathcal{H}_\text{InterAll} &= \sum_{ijkl\alpha\beta\gamma\delta} \sum_{\sigma_1 \sigma_2 \sigma_3 \sigma_4} I_{ijkl\alpha\beta\gamma\delta} c^\dagger_{i\alpha\sigma_1} c_{j\beta\sigma_2} c^\dagger_{k\gamma\sigma_3} c_{l\delta\sigma_4} \nonumber\\ &= \sum_{ijkl\alpha\beta\gamma\delta} \sum_{\sigma_1 \sigma_2 \sigma_3 \sigma_4} I_{ijkl\alpha\beta\gamma\delta} ( c^\dagger_{i\alpha\sigma_1} c^\dagger_{k\gamma\sigma_3} c_{l\gamma\sigma_4} c_{j\beta\sigma_2} + c^\dagger_{i\alpha\sigma_1} c_{l\delta\sigma_4} \delta_{j,k}\delta_{\beta,\gamma}\delta_{\sigma_2,\sigma_3} ) . \end{aligned}\end{split}\]

It is noted that there is a one-body term as depicted in the second term of the above expression. Then, the Hamiltonian given by the one-body terms is generally denoted as

(5.3)\[\begin{aligned} \mathcal{H}_\text{UHF} &= \sum_{ij} H_{ij} c^\dagger_{i} c_{j} = \hat{c}^\dagger H \hat{c} \end{aligned}\]

where we adopt a notation \(i\equiv(i, \alpha, \sigma_1), j\equiv(j, \beta, \sigma_2)\) for brevity, \(H\) denotes a matrix whose elements are \(H_{ij}\), and \(\hat{c}\) denotes a column vector whose elements are \(c_{i}\).

As \(H\) is an Hermite matrix, the Hamiltonian can be transformed into \(H=U \hat{\xi} U^\dagger\) where \(\hat{\xi}\) is a matrix whose diagonal elements are the eigenvalues of \(H\), and \(U\) is a matrix composed of the corresponding eigenvectors.

Then, let \(\hat{d} = U^\dagger \hat{c}\), and \(\mathcal{H}_\text{UHF}\) leads to

(5.4)\[\begin{aligned} \mathcal{H}_\text{UHF} &= \hat{d}^\dagger \hat{\xi} \hat{d} = \sum_{k} \xi_k d_k^\dagger d_k . \end{aligned}\]

Therefore, the energy derived from the one-body interaction term of the UHF approximation is obtained by

(5.5)\[\begin{aligned} E_\text{UHF} = \langle \mathcal{H}_\text{UHF} \rangle = \sum_{k} \xi_k \langle d_k^\dagger d_k \rangle . \end{aligned}\]

In the numerical calculation, as \(H\) depends on the one-body Green’s function \(\langle c_{i}^\dagger c_{j}\rangle\) through the UHF approximation, the equation is iteratively solved to satisfy the self-consistency. Starting from a one-body Green’s function given as an initial value, it is updated through the relation

(5.6)\[\begin{aligned} \langle c_{i}^\dagger c_{j}\rangle = \sum_{l} U_{il}^* U_{jl} \langle d_l^\dagger d_l \rangle = \sum_{l} \frac{U_{il}^* U_{jl}}{1+\exp^{\beta(\xi_l -\mu)}} \end{aligned}\]

until the one-body Green’s function converges. Here, \(\beta\) denotes the inverse temperature \(1/ k_B T\), and \(\mu\) denotes the chemical potential. In the canonical calculation in which the number of particles is fixed, \(\mu\) is determined to satisfy the relation

(5.7)\[\begin{aligned} N = \sum_{i} \langle c_i^{\dagger} c_i \rangle \end{aligned}\]

for the number of particles \(N\) at every step.

In H-wave, the simple-mixing algorithm is employed to update the configuration. If we denote the one-body Green’s function at n-th step by \(\langle c_{i}^\dagger c_{j}\rangle^{(n)}\), the Green’s function at n+1-th step is chosen by mixing that of n-th step with the new one obtained in n+1-th step as

(5.8)\[\begin{aligned} \langle c_{i}^\dagger c_{j}\rangle^{(n+1)} := (1-\alpha) \langle c_{i}^\dagger c_{j}\rangle^{(n)} + \alpha \langle c_{i}^\dagger c_{j}\rangle^{(n+1)}, \end{aligned}\]

where \(\alpha\) is a parameter between 0 and 1. There are other update algorithms such as Anderson mixing, though they are not supported in the present version of H-wave.

In the coordinate-space UHF mode of H-wave, all interactions are mapped to InterAll form. The free energy at finite temperature is given by

(5.9)\[\begin{aligned} F = \mu N -\frac{1}{\beta}\sum_k \ln \left[ 1+\exp (-\beta(\xi_k - \mu)) \right] - \sum_{ijkl} I_{ijkl} (\langle c_{i}^{\dagger} c_j\rangle \langle c_{k}^{\dagger} c_l\rangle - \langle c_{i}^{\dagger} c_l\rangle \langle c_{k}^{\dagger} c_j\rangle) . \end{aligned}\]

5.1.2. Extension to wave-number space

The Hamiltonian given by the one-body terms is rewritten in the wave-number representation by the Fourier transform \(c_i = \dfrac{1}{\sqrt{V}} \sum_k e^{ikr_i} c_k\) as

(5.10)\[\begin{aligned} \mathcal{H}_\text{UHF} &= \sum_{k\alpha\beta\sigma\sigma^\prime} h_{\alpha\beta\sigma\sigma^\prime}(k)\, c_{k\alpha\sigma}^\dagger c_{k\beta\sigma^\prime}^{\phantom\dagger} \end{aligned}\]

Here, the interaction is assumed to have translational symmetry so that the coefficients depend only on the translation vectors \(r_{ij}=r_j - r_i\). It is noted that InterAll type of interaction is not considered in the wave-number space UHF mode.

As the Hamiltonian is diagonal with respect to the wave number \(k\), the calculation of the eigenvalues and eigenvectors reduces from diagonalization of a matrix of the size \(N_\text{site}N_\text{orbit} \times N_\text{site}N_\text{orbit}\) to that of \(N_\text{site}\) matrices of the size \(N_\text{orbit} \times N_\text{orbit}\), which lowers the calculation costs. Here, \(N_\text{site}\) denotes the number of sites, and \(N_\text{orbit}\) denotes the number of orbitals including the spin degree of freedom.