Consider the 1D two-site Hubbard model at half filling
$H=-t\sum _{\sigma} (c_{1\sigma} ^{\dagger}c_{2\sigma}+h.c.)+U\sum_i(n_{i\uparrow}-\frac{1}{2})(n_{i\downarrow}-\frac{1}{2})$
where the sum is only on two sites $i=1,2$.
I am trying to find Green function component $G_{1\uparrow,1\uparrow}(t)=\langle\{ c_{1\uparrow}^{\dagger}(t), c_{1\uparrow} \}\rangle$
I am taking ground state $\mid\psi \rangle=\mid\uparrow;\downarrow\rangle-\mid\downarrow;\uparrow\rangle$
But I am unable to find the following result
$G_{1\uparrow,1\uparrow}(i\omega)=\frac{1/2}{(i\omega+U/2)+t}+\frac{1/2}{(i\omega+U/2)-t}$
Any help, suggestion where I am going wrong?
Best Answer
Your result only holds in the $U/t\gg 1$ limit.
A brute-force method to obtain the result is by exact diagonalization in the Fock space basis, where all the operators are represented as matrices. Let us introduce the Fock states $|n_{1\uparrow}n_{1\downarrow}n_{2\uparrow}n_{2\downarrow}\rangle$ where $n_{i\sigma}=0,1$ denotes the fermion number on site $i$ of the spin $\sigma$. The basis may be arranged in the following order:
$$|0000\rangle, |0001\rangle, |0010\rangle, |0011\rangle, \cdots, |1110\rangle, |1111\rangle,$$
just as binary numbers from 0 to 15 (totally 16 of them). Then the fermion operators can be written as matrices by actually acting the operators on the above set of basis and see how the basis transforms. The result reads
$$\begin{split} c_{1\uparrow}&=\frac{1}{2}(\sigma^1+i\sigma^2)\otimes\sigma^0\otimes\sigma^0\otimes\sigma^0,\\ c_{1\downarrow}&=\frac{1}{2}\sigma^3\otimes(\sigma^1+i\sigma^2)\otimes\sigma^0\otimes\sigma^0,\\ c_{2\uparrow}&=\frac{1}{2}\sigma^3\otimes\sigma^3\otimes(\sigma^1+i\sigma^2)\otimes\sigma^0,\\ c_{2\downarrow}&=\frac{1}{2}\sigma^3\otimes\sigma^3\otimes\sigma^3\otimes(\sigma^1+i\sigma^2),\\ \end{split}$$
where $\sigma^{1,2,3}$ are the Pauli matrices and $\sigma^0$ stands for the $2\times2$ identity matrix. It is easy to verify that the above matrix representation satisfies the fermion anti-commutation relation. Admittedly it is not simple to work with these matrices by hand, usually we can manipulate the matrices on, for example, Mathematica, where Pauli matrices are given by
Pauliand the direct product is implemented byKroneckerProduct.With this setup, we are ready the construct the Hamiltonian simply by multiplying the matrices and plus them together. Note that the $1/2$ constant in the Hubbard term $(n_{i\sigma}-1/2)$ must be treated as $1/2$ times an identity matrix (of $16\times 16$). Once we obtain the $16\times 16$ matrix representation of the Hamiltonian, we can diagonalize it to find the ground state wave function $|\Psi_0\rangle$ and the corresponding ground state energy $E_0$, s.t.
$$H |\Psi_0\rangle=E_0|\Psi_0\rangle.$$
Then it is straight forward to compute the Green's function by definition:
$$\begin{split} G(i\omega)_{1\uparrow,1\uparrow}&=-\langle c_{1\uparrow} c_{1\uparrow}^\dagger\rangle(i\omega)\\ &=\langle\Psi_0|c_{1\uparrow}(i\omega + H - E_0)^{-1}c_{1\uparrow}^\dagger|\Psi_0\rangle. \end{split}$$
Of cause the matrix inversion and multiplications are better carried out by a computer other than by hand. And the result that I found is:
$$G(i\omega)_{1\uparrow,1\uparrow}=\frac{\frac{1}{4}+\frac{t}{\sqrt{16 t^2+U^2}}}{i\omega + \frac{1}{2} \left(\sqrt{16 t^2+U^2}-2 t\right) }+\frac{\frac{1}{4}-\frac{t}{\sqrt{16 t^2+U^2}}}{i\omega+\frac{1}{2} \left(\sqrt{16 t^2+U^2}+2 t\right)}.$$
In the strong interaction limit $U/t\gg 1$, expanding both the numerator and the denominator respectively to the first order of $t/U$, I obtain:
$$G(i\omega)_{1\uparrow,1\uparrow}=\frac{1/4+t/U}{i\omega+U/2-t}+\frac{1/4-t/U}{i\omega+U/2+t},$$
which is the expected result. One may wonder that the spectral weights sum up to 1/2 other then 1. This is just because we have only calculated the $-\langle c_b c_a^\dagger\rangle$ half of the full Green's function $\langle\{c_a^\dagger, c_b\}\rangle$ (for the sake of simplicity). The other half will give identical result which brings back the factor 2.