Can Hartree-Fock determinant WLOG taken to be real?

For a many-electron Hamiltonian $$H$$, a Hartree-Fock determinant is a Slater determinant $$\Psi$$ that minimizes the energy $$\frac{\langle\Psi,H\Psi\rangle}{\langle\Psi,\Psi\rangle}$$. In general, $$\Psi$$ is constructed from a set of complex spin orbitals. However, if the Hamiltonian is real, like the non-relativistic $$N$$-electron Hamiltonian $$H = \sum_n^N \underbrace{-\frac{1}{2}\nabla_n^2 + V(\mathbf{r}_n)}_{=: h} + \frac{1}{2}\sum_{n\neq m}^N \frac{1}{|\mathbf{r}_n-\mathbf{r}_m|},$$ its eigenstates can without loss of generality be taken to be real. What about the Hartree-Fock determinants? If a Hartree-Fock determinant exists, is there then a Hartree-Fock determinant built only from real spin orbitals?

In standard non-relativistic Hartree-Fock calculations, the Hartree-Fock determinant is assumed to be real. But I have not seen any discussion about this in the literature, and I am a bit doubtful that this can be assumed without loss of generality since using complex Slater determinants seems to provide more variational freedom. To see this, let $$\Psi$$ be a Slater determinant built from complex spin orbitals $$\psi_1,\dots,\psi_N$$ that are chosen to be orthonormal for simplicity. Decompose each spin orbital in real and imaginary part. $$\psi_n = \psi_n^R + i\psi_n^I.$$ Then the energy of $$\Psi$$ is \begin{align} \langle\Psi,H\Psi\rangle &= \sum_n \left\langle\psi_n^R,h\psi_n^R\right\rangle + \left\langle\psi_n^I,h\psi_n^I\right\rangle \\ &+ \frac{1}{2}\sum_{n,m} \int d^3r_1\int d^3r_2 \frac{\left(\left|\psi_n^R(\mathbf{r}_1)\right|^2 + \left|\psi_n^I(\mathbf{r}_1)\right|^2\right)\left(\left|\psi_m^R(\mathbf{r}_2)\right|^2 + \left|\psi_m^I(\mathbf{r}_2)\right|^2\right)}{|\mathbf{r}_1-\mathbf{r}_2|} \\ &- \frac{1}{2}\sum_{n,m} \int d^3r_1\int d^3r_2 \frac{\left(\psi_n^R(\mathbf{r}_1)\cdot\psi_m^R(\mathbf{r}_1) + \psi_n^I(\mathbf{r}_1)\cdot\psi_m^I(\mathbf{r}_1)\right)\left(\psi_m^R(\mathbf{r}_2)\cdot\psi_n^R(\mathbf{r}_2) + \psi_m^I(\mathbf{r}_2)\cdot\psi_n^I(\mathbf{r}_2)\right)}{|\mathbf{r}_1-\mathbf{r}_2|}, \end{align} where the first sum is the one-electron energy, the second sum is the Coulomb energy and the third sum is the exchange energy. The one-electron energy neatly decomposes into a sum of a contribution from the real part and the imaginary part of the spin orbitals. However, the Coulomb and exchange energy contain cross terms between the real and imaginary parts which may allow for more variation of the energy.?

Real and complex HF are in general not equivalent. The realness of the HF wave function can be thought of as a symmetry (complex conjugation symmetry) that can often be broken in order to find a lower energy solution. This is the same kind of symmetry breaking that happens in Unrestricted HF (UHF), where the point group and $$\hat{S}{}^2$$ symmetry of Restricted HF (RHF) is sacrificed to get the energy closer to that of the exact ground state (and in Generalized HF (GHF), where the $$\hat{S}_z$$ symmetry is also broken).
Complex orbitals were also considered for RHF, with point group, $$\hat{S}{}^2$$ and $$\hat{S}_z$$ symmetries being all intact in J. Chem. Phys. 142 024104 (2015); complex RHF (and perturbation theory built on it) is seen to give a somewhat improved description over its real counterpart.