Quick Start
We wish to evaluate: \[ \langle\Psi|H| \Psi\rangle \] We will split the Hamiltonian into 5 parts, and evaluate these separately: \[ \langle\Psi|H| \Psi\rangle=\left\langle\Psi\left|H_{t}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{\mu}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{\Delta}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{J}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{h}\right| \Psi\right\rangle \] Some terms only \(\left|\Psi_{f}\right\rangle\) or \(\left|\Psi_{s}\right\rangle\), so simplify \[ \langle\Psi|H| \Psi\rangle=\left\langle\Psi_{f}\left|H_{t}\right| \Psi_{f}\right\rangle+\left\langle\Psi_{f}\left|H_{\mu}\right| \Psi_{f}\right\rangle+\left\langle\Psi\left|H_{\Delta}\right| \Psi\right\rangle+\left\langle\Psi_{s}\left|H_{J}\right| \Psi_{s}\right\rangle+\left\langle\Psi_{s}\left|H_{h}\right| \Psi_{s}\right\rangle \]
Ferromagnetic Calculation
Apply the variational method. The Hamiltonian is: \[ H=-t \sum_{i}\left(c_{i}^{\dagger} c_{i+1}+\text { h.c. }\right)-\mu \sum_{i} c_{i}^{\dagger} c_{i}+\Delta \sum_{i}\left(c_{i} \sigma_{i, i+1}^{z} c_{i+1}+\text { h.c. }\right)-J \sum_{i} \sigma_{i, i+1}^{z} \sigma_{i+1, i+2}^{z}-h \sum_{i} \sigma_{i, i+1}^{x} \] The trial wavefunction is: \[ \begin{array}{c} |\Psi\rangle=\left|\Phi_{s}\right\rangle \otimes\left|\Phi_{f}\right\rangle \\ \left|\Phi_{s}\right\rangle=\bigotimes_{j=1}^{N}\left(\cos (\theta / 2)|\uparrow\rangle_{j}+e^{-i \phi} \sin (\theta / 2)|\downarrow\rangle_{j}\right) \\ \left|\Phi_{f}\right\rangle=\prod_{k}\left(u_{k}+v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right)|0\rangle \end{array} \] The wavefunction must be normalized: \[ 1=\left\langle\Phi_{f} \mid \Phi_{f}\right\rangle=\prod_{k}\left(\left|u_{k}\right|^{2}+\left|v_{k}\right|^{2}\right) \]
Variational energy
We wish to evaluate: \[ \langle\Psi|H| \Psi\rangle \] We will split the Hamiltonian into 5 parts, and evaluate these separately: \[ \langle\Psi|H| \Psi\rangle=\left\langle\Psi\left|H_{t}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{\mu}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{\Delta}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{J}\right| \Psi\right\rangle+\left\langle\Psi\left|H_{h}\right| \Psi\right\rangle \] Some terms only \(\left|\Psi_{f}\right\rangle\) or \(\left|\Psi_{s}\right\rangle\), so simplify \[ \langle\Psi|H| \Psi\rangle=\left\langle\Psi_{f}\left|H_{t}\right| \Psi_{f}\right\rangle+\left\langle\Psi_{f}\left|H_{\mu}\right| \Psi_{f}\right\rangle+\left\langle\Psi\left|H_{\Delta}\right| \Psi\right\rangle+\left\langle\Psi_{s}\left|H_{J}\right| \Psi_{s}\right\rangle+\left\langle\Psi_{s}\left|H_{h}\right| \Psi_{s}\right\rangle \]
1 expectation value: m-field
The first expectation value is: \[ \left\langle\Psi_{s}\left|H_{h}\right| \Psi_{s}\right\rangle, \quad H_{h}=-h \sum_{i} \sigma_{i, i+1}^{x} \] The action of \(\sigma_{i, i+1}^{x}\) is that \[ \begin{array}{l} \sigma_{i, i+1}^{x}|\uparrow\rangle_{j}=|\downarrow\rangle_{j} \\ \sigma_{i, i+1}^{x}|\downarrow\rangle_{j}=|\uparrow\rangle_{j} \end{array} \] Hence \[ \begin{aligned} \left\langle\Psi_{s}\left|H_{h}\right| \Psi_{s}\right\rangle &=\left\langle\Psi_{s}\left|-h \sum_{i} \sigma_{i, i+1}^{x}\right| \Psi_{s}\right\rangle \\ &=-h \sum_{i}\left\langle\Psi_{s}\left|\sigma_{i, i+1}^{x}\right| \Psi_{s}\right\rangle \\ &=-h \sum_{i}\left\langle\Psi_{s}\right| \sigma_{i, i+1}^{x}\left(\bigotimes_{j=1}^{N}\left(\cos (\theta / 2)|\uparrow\rangle_{j}+e^{-i \phi} \sin (\theta / 2)|\downarrow\rangle_{j}\right)\right) \\ &=-h \sum_{i}\left\langle\Psi_{s}\right|\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{1} \otimes \ldots\left[\begin{array}{c} e^{-i \phi} \sin (\theta / 2) \\ \cos (\theta / 2) \end{array}\right]_{i} \otimes \ldots\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{N} \\ &=-h \sum_{i}\left(\left[\cos (\theta / 2) \quad e^{i \phi} \sin (\theta / 2)\right]_{1} \otimes \ldots\left[\cos (\theta / 2) \quad e^{i \phi} \sin (\theta / 2)\right]_{N}\right) \\ &\left(\left[\begin{array}{c} \cos (\theta / 2) \\ e^{i \phi} \sin (\theta / 2) \end{array}\right] \otimes \ldots\left[\begin{array}{c} e^{-i \phi} \sin (\theta / 2) \\ \cos (\theta / 2) \end{array}\right]_{i} \otimes \ldots\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{N}\right) \\ &=-h \sum_{i} \cos (\theta / 2) \sin (\theta / 2) e^{-i \phi}+\cos (\theta / 2) \sin (\theta / 2) e^{i \phi} \\ &=-2 h N \cos (\theta / 2) \sin (\theta / 2)\left(\frac{e^{-i \phi}+e^{i \phi}}{2}\right) \\ &=-h N \sin (\theta) \cos (\phi) \end{aligned} \]
2 expectation value: exchange interaction
The next expectation value is: \[ \left\langle\Psi_{s}\left|H_{J}\right| \Psi_{s}\right\rangle, \quad H_{J}=-J \sum_{i} \sigma_{i, i+1}^{z} \sigma_{i+1, i+2}^{z} \] The action of \(\sigma_{i, i+1}^{z}\) is that \[ \begin{array}{l} \sigma_{i, i+1}^{z}|\uparrow\rangle_{j}=|\uparrow\rangle_{j} \\ \sigma_{i, i+1}^{z}|\downarrow\rangle_{j}=-|\downarrow\rangle_{j} \end{array} \]
\[ \begin{aligned} \left\langle\Psi_{s}\left|H_{J}\right| \Psi_{s}\right\rangle&=\left\langle\Psi_{s}\left|-J \sum_{i} \sigma_{i, i+1}^{z} \sigma_{i+1, i+2}^{z}\right| \Psi_{s}\right\rangle &=-J \sum_{i}\left\langle\Psi_{s}\left|\sigma_{i, i+1}^{z} \sigma_{i+1, i+2}^{z}\right| \Psi_{s}\right\rangle\\ &=-J \sum_{i}\left\langle\Psi_{s}\right| \sigma_{i, i+1}^{z} \sigma_{i+1, i+2}^{z}\left(\bigotimes_{j=1}^{N}\left(\cos (\theta / 2)|\uparrow\rangle_{j}+e^{-i \phi} \sin (\theta / 2)|\downarrow\rangle_{j}\right)\right)\\ &=-J \sum_{i}\left\langle\Psi_{s}\right|\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{1} \otimes \ldots\left[\begin{array}{c} \cos (\theta / 2) \\ -e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{i} \otimes\left[\begin{array}{c} \cos (\theta / 2) \\ -e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{i+1}\\ &\otimes \ldots\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{N}\\ &=-J \sum_{i}\left(\left[\cos (\theta / 2) \quad e^{i \phi} \sin (\theta / 2)\right]_{1} \otimes \ldots\left[\cos (\theta / 2) \quad e^{i \phi} \sin (\theta / 2)\right]_{N}\right)\\ &\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{1} \otimes \ldots\left[\begin{array}{c} \cos (\theta / 2) \\ -e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{i} \otimes\left[\begin{array}{c} \cos (\theta / 2) \\ -e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{i+1} \otimes \ldots\left[\begin{array}{c} \cos (\theta / 2) \\ e^{-i \phi} \sin (\theta / 2) \end{array}\right]_{N}\\ &=-J \sum_{i}\left(\cos ^{2}(\theta / 2)-\sin ^{2}(\theta / 2)\right)_{i}\left(\cos ^{2}(\theta / 2)-\sin ^{2}(\theta / 2)\right)_{i+1}\\ &=-J N \cos ^{2}(\theta) \end{aligned} \]
3 expectation value: hopping
The third expectation value is: \[ \left\langle\Psi_{f}\left|H_{t}\right| \Psi_{f}\right\rangle, \quad-t \sum_{i}\left(c_{i}^{\dagger} c_{i+1}+\text { h.c. }\right) \] Let us simplify this Hamiltonian first, which is written in position space. We need to transform this to momentum space: \[ H_{t}=-t\left(\sum_{x} c_{x}^{\dagger} c_{x+1}+h . c .\right) \] The fourier transform is: \[ \begin{aligned} c_{x} &=\frac{1}{\sqrt{N}} \sum_{k} e^{i x k} c_{k} \\ c_{x}^{\dagger} &=\frac{1}{\sqrt{N}} \sum_{k} e^{-i x k} c_{k}^{\dagger} \end{aligned} \] So we just plug it in: \[ \begin{aligned} H_{t} &=-t\left(\sum_{x} \frac{1}{N}\left(\sum_{k} e^{-i x k} c_{k}^{\dagger}\right) \sum_{k^{\prime}}\left(e^{i(x+1) k^{\prime}} c_{k^{\prime}}\right)+\text { h.c. }\right) \\ &=-t\left(\sum_{x} \frac{1}{N}\left(\sum_{k} \sum_{k^{\prime}} e^{i x\left(k^{\prime}-k\right)} e^{i k^{\prime}} c_{k}^{\dagger} c_{k^{\prime}}\right)+h . c .\right) \\ &=-t\left(\left(\sum_{k} \sum_{k^{\prime}} \delta\left(k^{\prime}-k\right) e^{i k^{\prime}} c_{k}^{\dagger} c_{k^{\prime}}\right)+h . c .\right) \\ &=-t\left(\left(\sum_{k} e^{i k} c_{k}^{\dagger} c_{k}\right)+h . c .\right) \\ &=-t\left(\sum_{k} e^{i k} c_{k}^{\dagger} c_{k}+\sum_{k} e^{-i k} c_{k}^{\dagger} c_{k}\right) \\ &=-t\left(\sum_{k}\left(e^{i k}+e^{-i k}\right) c_{k}^{\dagger} c_{k}\right) \\ &=-2 t\left(\sum_{k} \cos (k) c_{k}^{\dagger} c_{k}\right) \end{aligned} \] Now we can calculate the expectation: \[ \begin{aligned} \left\langle\Psi_{f}\left|H_{t}\right| \Psi_{f}\right\rangle &=-2 t\left(\sum_{k}\left\langle\Psi_{f}\left|\cos (k) c_{k}^{\dagger} c_{k}\right| \Psi_{f}\right\rangle\right) \\ &=-2 t\left(\sum_{k}\left\langle\Psi_{f}\left|\cos (k) c_{k}^{\dagger} c_{k}\left(u_{k}+i v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right) \prod_{q \neq k}\left(u_{q}+i v_{q} c_{q}^{\dagger} c_{-q}^{\dagger}\right)\right| 0\right\rangle\right) \\ &=-2 t\left(\sum_{k}\left\langle\psi\left|\left(v_{k} c_{-k} c_{k}\right) \cos (k) c_{k}^{\dagger} c_{k}\left(v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right)\right| \psi\right\rangle\right) \\ &=-2 t \sum_{k} v_{k}^{2} \cos (k) \end{aligned} \] Where \[ |\psi\rangle=\prod_{q \neq k}\left(u_{q}+i v_{q} c_{q}^{\dagger} c_{-q}^{\dagger}\right)|0\rangle \]
4 chemical
The fourth term is the \(\mu\) term: \[ \left\langle\Psi_{f}\left|H_{\mu}\right| \Psi_{f}\right\rangle, \quad-\mu \sum_{i} c_{i}^{\dagger} c_{i} \] Begin by simplifying: \[ \begin{aligned} H_{\mu} &=-\mu \sum_{i} c_{i}^{\dagger} c_{i} \\ &=-\mu \sum_{x}\left(\frac{1}{N} \sum_{k} e^{-i x k} c_{k}^{\dagger} \sum_{k^{\prime}} e^{i x k^{\prime}} c_{k^{\prime}}\right) \\ &=-\mu \sum_{x} \frac{1}{N} \sum_{k} \sum_{k^{\prime}} e^{i x\left(k^{\prime}-k\right)} c_{k}^{\dagger} c_{k^{\prime}} \\ &=-\mu \sum_{k} \sum_{k^{\prime}} \delta\left(k^{\prime}-k\right) \quad c_{k}^{\dagger} c_{k^{\prime}} \\ &=-\mu \sum_{k} c_{k}^{\dagger} c_{k} \end{aligned} \] By comparison with the previous result, the expected value is \[ \left\langle\Psi_{f}\left|H_{\mu}\right| \Psi_{f}\right\rangle=-\mu N \sum_{k} v_{k}^{2} \]
5 SC term
The last term is the superconducting term: \[ \left\langle\Psi\left|H_{\Delta}\right| \Psi\right\rangle, \quad H_{\Delta}=-\Delta \sum_{i}\left(c_{i} \sigma_{i, i+1}^{z} c_{i+1}+\text { h.c. }\right) \] It is clear \(H_{\Delta}\) contains a spin and a fermion term. We calculate them separately and multiply at the end: \[ \left\langle\Psi_{s}\left|\sigma_{i, i+1}^{z}\right| \Psi_{s}\right\rangle=\cos (\theta) \] The fermion part yields: \[ \left\langle\Psi_{f}\right| c_{i} c_{i+1}+\text { h.c. }\left|\Psi_{f}\right\rangle \]
\[ \begin{aligned} \sum_{i} c_{i} c_{i+1}+\text { h.c. } &=\sum_{i} \frac{1}{N} \sum_{k} e^{i x k} c_{k} \sum_{k^{\prime}} e^{i(x+1) k^{\prime}} c_{k^{\prime}}+h . c . \\ &=\sum_{i} \frac{1}{N} \sum_{k, k^{\prime}} e^{i x\left(k+k^{\prime}\right)} e^{i k^{\prime}} c_{k} c_{k^{\prime}}+h . c . \\ &=\sum_{k, k^{\prime}} \delta\left(k+k^{\prime}\right) e^{i k^{\prime}} c_{k} c_{k^{\prime}}+h . c . \\ &=\sum_{k} e^{-i k} c_{k} c_{-k}+h . c . \\ &=\sum_{k} e^{-i k} c_{k} c_{-k}+e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger} \end{aligned} \]
As before, defining: \[ |\psi\rangle=\prod_{q \neq k}\left(u_{q}+i v_{q} c_{q}^{\dagger} c_{-q}^{\dagger}\right)|0\rangle \] Now we can evaluate: \[ \begin{array}{l} \left\langle\Psi_{f}\right| c_{i} c_{i+1}+\text { h.c. }\left|\Psi_{f}\right\rangle\\=\sum_{k}\left\langle\psi\left|\left(u_{k}-i v_{k} c_{-k} c_{k}\right)\left(e^{-i k} c_{k} c_{-k}+e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger}\right)\left(u_{k}+i v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right)\right| \psi\right\rangle\\ \left(u_{k}-i v_{k} c_{-k} c_{k}\right)\left(e^{-i k} c_{k} c_{-k}+e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger}\right)\left(u_{k}+i v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right) \\ =\left(u_{k}-i v_{k} c_{-k} c_{k}\right)\left(\left(e^{-i k} c_{k} c_{k}+e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger}\right) u_{k}+\left(e^{-i k} c_{k} c_{-k}+\underline{e^{i k}} e_{-k}^{\dagger} c_{k}^{\dagger}\right) i v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right) \\ =\left(u_{k}-i v_{k} c_{-k} c_{k}\right)\left(\left(e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger}\right) u_{k}+\left(e^{-i k} c_{k} c_{-k}\right) i v_{k} c_{k}^{\dagger} c_{-k}^{\dagger}\right) \\ =\left(u_{k}-i v_{k} c_{-k} c_{k}\right)\left(e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger} u_{k}-e^{-i k} i v_{k}\right) \\ =\left(u_{k} e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger} u_{k}-u_{k} e^{-i k} i v_{k}-i v_{k} c_{-k} c_{k} e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger} u_{k}+i v_{k} c_{-k} c_{k} e^{-i k} i v_{k}\right) \\ =\left(\underline{u_{k} e^{i k}} e_{-k}^{\dagger} c_{k}^{\dagger} u_{k}-u_{k} e^{-i k} i v_{k}+i v_{k} e^{i k} u_{k}+\underline{i v_{k} c_{k} e_{k} e^{-i k} \hat{i v_{k}}}\right) \\ =\left(-u_{k} e^{-i k} i v_{k}+i v_{k} e^{i k} u_{k}\right) \\ =\left(i u_{k} v_{k}\left(e^{i k}-e^{-i k}\right)\right) \\ =\left(-2 u_{k} v_{k} \frac{\left(e^{i k}-e^{-i k}\right)}{2 i}\right)=\left(-2 u_{k} v_{k} \sin (k)\right) \end{array} \] So the final expectation value is: \[ \left\langle\Psi\left|H_{\Delta}\right| \Psi\right\rangle=\left\langle\Psi_{s}\left|\sigma_{i, i+1}^{z}\right| \Psi_{s}\right\rangle\left\langle\Psi_{f}\right| c_{i} c_{i+1}+\text { h.c. }\left|\Psi_{f}\right\rangle=-\Delta \sum_{k}\left(2 u_{k} v_{k} \sin (k)\right) \cos (\theta) \]
entire
The entire expression is: \[ \begin{aligned} \langle\Psi|H| \Psi\rangle &=-h N \sin (\theta) \cos (\phi)-J N \cos ^{2}(\theta) \\ &-2 t \sum_{k} v_{k}^{2} \cos (k)-\mu N \sum_{k} v_{k}^{2}-2 \Delta \sum_{k} u_{k} v_{k} \cos (\theta) \sin (k) \end{aligned} \]
Minimization
We would like to minimize \[ \begin{aligned} E_{\text {var }}=\langle\Psi|H| \Psi\rangle &=-h N \sin (\theta) \cos (\phi)-J N \cos ^{2}(\theta) \\ &-2 t \sum_{k} v_{k}^{2} \cos (k)-\mu N \sum_{k} v_{k}^{2}-2 \Delta \sum_{k} u_{k} v_{k} \cos (\theta) \sin (k) \end{aligned} \] Making the substitutions and taking the derivative, \[ u_{k}=\cos \alpha_{k} \quad v_{k}=\sin \alpha_{k} \] we obtain \[ \begin{aligned} \frac{\partial E_{v a r}}{\partial \alpha(k)} &=-2 \Delta \sum_{k}\left(\cos (\theta) \sin (k) \cos ^{2}(\alpha(k))-\cos (\theta) \sin (k) \sin ^{2}(\alpha(k))\right) \\ &-N \mu \sum_{k} 2 \sin (\alpha(k)) \cos (\alpha(k))-2 t \sum_{k} 2 \cos (k) \sin (\alpha(k)) \cos (\alpha(k)) \end{aligned} \] This can be rewritten as \[ \begin{aligned} \frac{\partial E_{v a r}}{\partial \alpha(k)} &=-2 \Delta \sum_{k}(\cos (\theta) \sin (k) \cos (2 \alpha(k))) \\ &-N \mu \sum_{k} \sin (2 \alpha(k))-2 t \sum_{k} \cos (k) \sin (2 \alpha(k)) \end{aligned} \] And this allows us to rearrange and solve; since the entire expression is linear in \(k\), each term in the summation must be minimal: \[ \begin{aligned} & 0=\frac{\partial E_{v a r}}{\partial \alpha(k)} \Longrightarrow \quad E_{k_{v a r}} \text { minimal } \\ \Longrightarrow & 2 \Delta \cos (\theta) \sin (k) \cos (2 \alpha(k))=(-N \mu-2 t \cos (k)) \sin (2 \alpha(k)) \\ \Longrightarrow & \frac{2 \Delta \cos (\theta) \sin (k)}{-N \mu-2 t \cos (k)}=\tan (2 \alpha(k)) \\ \Rightarrow & \alpha(k)=\frac{1}{2} \arctan \left(\frac{2 \Delta \cos (\theta) \sin (k)}{-N \mu-2 t \cos (k)}\right) \end{aligned} \]
parameter
Now we need to solve for \(v_{k}^{2}\) and \(v_{k} u_{k}\) :
Defining \[ \begin{array}{c} \epsilon_{k}:=-2 t \cos (k) \\ \xi_{k}:=\epsilon_{k}-N \mu \\ \Delta_{k}:=2 \Delta \cos (\theta) \sin (k) \\ \tan (2 \alpha(k))=\frac{\Delta_{k}}{\xi_{k}} \Longrightarrow \cos (2 \alpha(k))=\frac{\xi_{k}}{\sqrt{\xi_{k}^{2}+\Delta_{k}^{2}}}, \quad \sin (2 \alpha(k))=\frac{\Delta_{k}}{\sqrt{\xi_{k}^{2}+\Delta_{k}^{2}}} \end{array} \] After some trig identities and plugging in the definition, \[ \begin{aligned} u_{k} v_{k} &=\cos (\alpha(k)) \sin (\alpha(k))=\frac{1}{2} \sin (2 \alpha(k))=\frac{1}{2} \frac{\Delta_{k}}{\sqrt{\xi_{k}^{2}+\Delta_{k}^{2}}} \\ v_{k}^{2} &=\sin ^{2}(\alpha(k))=\frac{1}{2}\left(1-\cos \left(2 \alpha_{k}\right)\right)=\frac{1}{2}\left(1-\frac{\xi_{k}}{\sqrt{\xi_{k}^{2}+\Delta_{k}^{2}}}\right) \end{aligned} \] Now we are able to integrate these equations in Mathematica.
Integration
\[ \sum_{k} \rightarrow N \int_{-\pi}^{\pi} \frac{d k}{2 \pi} \]
After substituting \(\alpha(k)\), we take the thermodynamic limit and evaluate the integrals using Mathematica: \[ \begin{aligned} E_{v a r} &=-h N \sin (\theta) \cos (\phi)-J N \cos ^{2}(\theta) \\ &-\frac{t N}{\pi} \int_{-\pi}^{\pi} v_{k}^{2} \cos (k) d k-\frac{\mu N^{2}}{2 \pi} \int_{-\pi}^{\pi} v_{k}^{2} d k-\frac{\Delta N}{\pi} \int_{-\pi}^{\pi} u_{k} v_{k} \cos (\theta) \sin (k) d k \end{aligned} \] Calculating the ground energy density: \[ \begin{aligned} \frac{E_{v a r}}{N} &=-h \sin (\theta) \cos (\phi)-J \cos ^{2}(\theta) \\ &-\frac{t}{\pi} \int_{-\pi}^{\pi} v_{k}^{2} \cos (k) d k-\frac{\mu N}{2 \pi} \int_{-\pi}^{\pi} v_{k}^{2} d k-\frac{\Delta}{\pi} \int_{-\pi}^{\pi} u_{k} v_{k} \cos (\theta) \sin (k) d k \end{aligned} \] We solve the integrals in the case \(\mu=0\). The \(\mu\) integral vanishes, and we are left with the \(t\) and \(\Delta\) integrals. Mathematica is able to evaluate this integral: The \(t\) integral is: \[ \begin{array}{l} -\frac{t}{\pi} \int_{-\pi}^{\pi} v_{k}^{2} \cos (k) d k \\ =\frac{\Delta t \cos (\theta)\left(\Delta \cos (\theta) K\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right)+t K\left(1-\frac{t^{2} \sec ^{2}(\theta)}{\Delta^{2}}\right)-t E\left(1-\frac{t^{2} \sec ^{2}(\theta)}{\Delta^{2}}\right)\right)-t^{3} E\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right)}{\pi\left(t^{2}-\Delta^{2} \cos ^{2}(\theta)\right)} \end{array} \] The \(\Delta\) integral is: \[ \begin{array}{l} -\frac{\Delta}{\pi} \int_{-\pi}^{\pi} u_{k} v_{k} \cos (\theta) \sin (k) d k \\ =\frac{\Delta \cos (\theta)\left(t^{2} K\left(1-\frac{t^{2} \sec ^{2}(\theta)}{\Delta^{2}}\right)-\Delta \cos (\theta)\left(-t K\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right)+t E\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right)+\Delta \cos (\theta) E\left(1-\frac{t^{2} \sec ^{2}(\theta)}{\Delta^{2}}\right)\right)\right)}{\pi\left(\Delta^{2} \cos ^{2}(\theta)-t^{2}\right)} \end{array} \] Combining the two gives: \[ -\frac{t}{\pi} \int_{-\pi}^{\pi} v_{k}^{2} \cos (k) d k-\frac{\Delta}{\pi} \int_{-\pi}^{\pi} u_{k} v_{k} \cos (\theta) \sin (k) d k=-\frac{t E\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right)+\Delta \cos (\theta) E\left(1-\frac{t^{2} \sec ^{2}(\theta)}{\Delta^{2}}\right)}{\pi} \] So, the final expression is: \[ \frac{E_{v a r}}{N}=-h \sin (\theta) \cos (\phi)-J \cos ^{2}(\theta)-\left(\frac{t E\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right)+\Delta \cos (\theta) E\left(1-\frac{t^{2} \sec ^{2}(\theta)}{\Delta^{2}}\right)}{\pi}\right) \] Which simplifies to: [Add proof of Elliptic E equivalence] \[ \frac{E_{v a r}}{N}=-J \cos ^{2}(\theta)-h \sin (\theta) \cos (\phi)-\frac{2 t}{\pi} E\left(1-\frac{\Delta^{2} \cos ^{2}(\theta)}{t^{2}}\right) \]
Minimization 2
To minimize \(\phi\), we observe that only the \(h\) -term depends on \(\phi\). Taking the derivative we see that \(\cos (\phi)\) is minimized at \(0, \pi\). Since \(\cos (\theta)\) is maximal at 0, the minimal value of \(\phi=0\).
To minimize the expression, simply write \(E_{v a r} / N\) as a function of \(\cos (\theta)\) minimize with respect to \(\cos (\theta)\).
Anti-ferromagnetic Calculation
First we relabel the indices.
We have the mapping \[ j \rightarrow j_{A}, \quad j+1 \rightarrow j_{B} \] Defining \[ \psi_{j}:=\left[\begin{array}{c} c_{j_{A}} \\ c_{j_{B}} \\ c_{j_{B}}^{\dagger} \\ c_{j_{B}}^{\dagger} \end{array}\right] \quad \psi_{j}^{\dagger}:=\left[\begin{array}{llll} c_{j_{A}}^{\dagger} & c_{j_{B}}^{\dagger} & c_{j_{A}} & c_{j_{B}} \end{array}\right] \]
\[ \begin{aligned} H_{t}=&-t \sum_{j} c_{j}^{\dagger} c_{j+1}+\text { h.c. } \\ =&-t \sum_{j} c_{j_{A}}^{\dagger} c_{j_{B}}+c_{j_{B}}^{\dagger} c_{j_{A}} \\ =&-\frac{t}{2} \sum_{j}\left(c_{j_{A}}^{\dagger} c_{j_{B}}-c_{j_{B}} c_{j_{A}}^{\dagger}\right)+\left(c_{j_{B}}^{\dagger} c_{j_{A}}-c_{j_{A}} c_{j_{B}}^{\dagger}\right) \\ =& \frac{1}{2} \sum_{j}\left[\begin{array}{llll} c_{j_{A}}^{\dagger} & c_{j_{B}}^{\dagger} & c_{j_{A}} & c_{j_{B}} \end{array}\right]\left[\begin{array}{cccc} 0 & -t & 0 & 0 \\ -t & 0 & 0 & 0 \\ 0 & 0 & 0 & t \\ 0 & 0 & t & 0 \end{array}\right]\left[\begin{array}{c} c_{j_{A}} \\ c_{j_{B}} \\ c_{j A}^{\dagger} \\ c_{j_{B}}^{\dagger} \end{array}\right] \\ =&\frac{1}{2} \sum_{j} \psi_{j}^{\dagger} H _{t j} \psi_{j} \end{aligned} \]
The fermion part of the \(\Delta\) term is also rewritten similarly: \[ \begin{aligned} H_{\Delta} &=\Delta \sum_{j}\left(c_{j} c_{j+1}+\text { h.c. }\right) \\ &=\Delta \sum_{j}\left(c_{j_{A}} c_{j_{B}}+c_{j_{B}}^{\dagger} c_{j_{A}}^{\dagger}\right) \\ &=\frac{\Delta}{2} \sum_{j}\left(c_{j_{A}} c_{j_{B}}-c_{j_{B}} c_{j_{A}}\right)+\left(c_{j_{B}}^{\dagger} c_{j_{A}}^{\dagger}-c_{j_{A}}^{\dagger} c_{j_{B}}^{\dagger}\right) \\ &=\frac{1}{2} \sum_{j}\left[\begin{array}{cccc} c_{j_{A}}^{\dagger} & c_{j_{B}}^{\dagger} & c_{j_{A}} & c_{j_{B}} \end{array}\right]\left[\begin{array}{cccc} 0 & 0 & 0 & -\Delta \\ 0 & 0 & \Delta & 0 \\ 0 & \Delta & 0 & 0 \\ -\Delta & 0 & 0 & 0 \end{array}\right]\left[\begin{array}{c} c_{j_{A}} \\ c_{j_{B}} \\ c_{j A}^{\dagger} \\ c_{j_{B}}^{\dagger} \end{array}\right] \\ &= \frac{1}{2} \sum_{j} \psi_{j}^{\dagger} H _{\Delta j} \psi_{j} \\ \end{aligned} \]
\[ H _{j}= H _{\Delta j}+ H _{t j}=\left[\begin{array}{cccc} 0 & -t & 0 & -\Delta \\ -t & 0 & \Delta & 0 \\ 0 & \Delta & 0 & t \\ -\Delta & 0 & t & 0 \end{array}\right] \]
Diagonalization in Momentum Space
Since \[ \begin{array}{c} H_{t}=-t\left(\sum_{x} c_{x}^{\dagger} c_{x+1}+h . c .\right)=-2 t \sum_{k} \cos (k) c_{k}^{\dagger} c_{k} \\ H_{\Delta}=-\Delta \sum_{i}\left(c_{i} \sigma_{i, i+1}^{z} c_{i+1}+\text { h.c. }\right)=-\Delta \sum_{k} e^{-i k} c_{k} c_{-k}+e^{i k} c_{-k}^{\dagger} c_{k}^{\dagger} \end{array} \] We see that \[ \begin{aligned} t \rightarrow \epsilon_{k} &:=-2 t \cos (k) \\ \Delta \rightarrow \Delta_{k} &:=2 \Delta \cos (\theta) \sin (k) \end{aligned} \] [Missing \(j +2, j\) terms but next equation is correct $]$ \[ H _{k}=\left(\begin{array}{cccc} 0 & -t\left(1+e^{2 i k}\right) & 0 & -\Delta \cos (\theta)\left(1+e^{2 i k}\right) \\ -t\left(1+e^{-2 i k}\right) & 0 & \Delta \cos (\theta)\left(1+e^{-2 i k}\right) & 0 \\ 0 & 0 & t\left(1+e^{2 i k}\right) \\ -\Delta \cos (\theta)\left(1+e^{-2 i k}\right) & 0 & t\left(1+e^{-2 i k}\right) & 0 \end{array}\right) \] After diagonalization, we find the following four eigenvalues: \[ \{2(t-\Delta \cos (\theta)) \cos (k), 2(\Delta \cos (\theta)-t) \cos (k),-2(\Delta \cos (\theta)+t) \cos (k), 2(\Delta \cos (\theta)+t) \cos (k)\} \] Then, converting the summations into integrals once more we obtain: \[ \begin{array}{c} \sum_{k} \rightarrow N \int_{-\pi / 2}^{\pi / 2} \frac{d k}{2 \pi} \\ \{(t-\Delta \cos (\theta)) / \pi,(\Delta \cos (\theta)-t) / \pi,-(\Delta \cos (\theta)+t) / \pi,(\Delta \cos (\theta)+t) / \pi\} \end{array} \] Notice that if we set \(\Delta=0\), we immediately recover \[ -\frac{2 t}{\pi} \] for the fermion term, for both the ferromagnetic and the antiferromagnetic calculation. Hence \[ \frac{E_{v a r}}{N}=-J \cos ^{2}(\theta)-h \sin (\theta) \cos (\phi)-\frac{(t+\Delta \cos (\theta))}{\pi}+-\frac{|t-\Delta \cos (\theta)|}{\pi} \]
Dr.Huang
(1)2d superconductor
The Hamiltonian on the square lattice can be given by \[ \begin{aligned} H=&-\mu \sum_{x=1}^{N_{x}} \sum_{y=1}^{N_{y}} c_{x, y}^{\dagger} c_{x, y} \\ &-\sum_{x=1}^{N_{x}} \sum_{y=1}^{N_{y}}\left[t_{x} c_{x, y}^{\dagger} c_{x+1, y}+t_{y} c_{x, y}^{\dagger} c_{x, y+1}+h . c .\right] \\ &-\sum_{x=1}^{N_{x}} \sum_{y=1}^{N_{y}}\left[\Delta_{x} c_{x, y}^{\dagger} c_{x+1, y}^{\dagger}+\Delta_{y} c_{x, y}^{\dagger} c_{x, y+1}^{\dagger}+h . c .\right] \end{aligned} \] where \(c_{x, y}^{\dagger}\) are the fermion creation operators on site \((x, y) .\) The intrachain transfer integral \(t_{x}\) and the intrachain superconducting pairing amplitude \(d_{x}\) are present in the Kitaev model, while we have newly introduced the effective interchain transfer integral \(t_{y}\) and the superconducting pairing amplitude \(d_{y}\). We note that the system is reduced to the Kitaev chain when \(N_{y}=1\) If the phase difference between \(d_{x}\) and \(d_{y}\) is \(\pi / 2\) and \(N_{x}=N_{y}\), it is an anisotropic \(p+i p\) superconductor system.
Fourier
Fourier transformation along y-direction \[ c_{x, y}^{\dagger}=1 / \sqrt{N_{y}} \sum_{k} e^{i k y} c_{x, k}^{\dagger} \] each term
- \(\sum_{y} c_{x, y}^{\dagger} c_{x, y}=\sum_{k} c_{x, k}^{\dagger} c_{x, k}\)
- \(\sum_{y}\left[t_{x} c_{x, y}^{\dagger} c_{x+1, y}+t_{x}^{*} c_{x+1, y}^{\dagger} c_{x, y}\right]=\sum_{k}\left[t_{x} c_{x, k}^{\dagger} c_{x+1, k}+t_{x}^{*} c_{x+1, k}^{\dagger} c_{x, k}\right]\)
- \(\sum_{y}\left[t_{y} c_{x, y}^{\dagger} c_{x, y+1}+t_{y}^{*} c_{x, y+1}^{\dagger} c_{x, y}\right]=\sum_{k}\left[t_{y} e^{-i k} c_{x, k}^{\dagger} c_{x, k}+t_{y}^{*} e^{i k} c_{x, k}^{\dagger} c_{x, k}\right]\)
- \(\sum_{y}\left[\Delta_{x} c_{x, y}^{\dagger} c_{x+1, y}^{\dagger}+\Delta_{x}^{*} c_{x+1, y} c_{x, y}\right]=\sum_{k}\left[\Delta_{x} c_{x, k}^{\dagger} c_{x+1,-k}^{\dagger}+\Delta_{x}^{*} c_{x+1,-k} c_{x, k}\right]\)
- \(\sum_{y}\left[\Delta_{y} c_{x, y}^{\dagger} c_{x, y+1}^{\dagger}+\Delta_{y}^{*} c_{x, y+1} c_{x, y}\right]=\sum_{k}\left[\Delta_{y} e^{-i k} c_{x, k}^{\dagger} c_{x,-k}^{\dagger}+\Delta_{y}^{*} e^{i k} c_{x,-k} c_{x, k}\right]\)
The Hamiltonian can be rewritten as \[ \begin{aligned} H=\sum_{x=1}^{N_{x}} \sum_{k} &\left[\left(-\mu+t_{y} e^{-i k}+t_{y}^{*} e^{i k}\right) c_{x, k}^{\dagger} c_{x, k}+\right.\\ &\left(t_{x}\right) c_{x, k}^{\dagger} c_{x+1, k}+\left(t_{x}\right)^{*} c_{x+1, k}^{\dagger} c_{x, k}+\\ &\left(\Delta_{x}\right) c_{x, k}^{\dagger} c_{x+1,-k}^{\dagger}+\left(\Delta_{x}\right)^{*} c_{x+1,-k} c_{x, k}+\\ &\left.\left(\Delta_{y} e^{-i k}\right) c_{x, k}^{\dagger} c_{x,-k}^{\dagger}+\left(\Delta_{y} e^{i k}\right)^{*} c_{x,-k} c_{x, k}\right] \end{aligned} \]
Majorana operators
We fthen note that the definition of Majorana fermions (MFs) in real space is different from that in momentum space. Therefore, by doing Fourier Transformation, the Majorana operators in momentum space can be given \[ \begin{array}{l} c_{x, k}=\frac{1}{2}\left[\gamma_{x, k}^{(1)}+i \gamma_{x, k}^{(2)}\right] \\ c_{x, k}^{\dagger}=\frac{1}{2}\left[\gamma_{x,-k}^{(1)}-i \gamma_{x,-k}^{(2)}\right] \\ c_{x,-k}=\frac{1}{2}\left[\gamma_{x,-k}^{(1)}+i \gamma_{x,-k}^{(2)}\right] \\ c_{x,-k}^{\dagger}=\frac{1}{2}\left[\gamma_{x, k}^{(1)}-i \gamma_{x, k}^{(2)}\right] \end{array} \] Here, we show how the fermionic operator transfer to Majorana fermion operators. \[ \begin{array}{|c|c|c|} \hline(\tilde{\mu}) & c_{x, k}^{\dagger} c_{x, k} & \frac{1}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(1)}+\gamma_{x,-k}^{(2)} \gamma_{x, k}^{(2)}\right)+\frac{i}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(2)}-\gamma_{x,-k}^{(2)} \gamma_{x, k}^{(1)}\right) \\ \hline\left(t_{x}\right) & c_{x, k}^{\dagger} c_{x+1, k} & \frac{1}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x+1, k}^{(1)}+\gamma_{x,-k}^{(2)} \gamma_{x+1, k}^{(2)}\right)+\frac{i}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x+1, k}^{(2)}-\gamma_{x,-k}^{(2)} \gamma_{x+1, k}^{(1)}\right) \\ \hline\left(t_{x}\right)^{*} & c_{x+1, k}^{\dagger} c_{x, k} & \frac{1}{4}\left(\gamma_{x+1,-k}^{(1)} \gamma_{x, k}^{(1)}+\gamma_{x+1,-k}^{(2)} \gamma_{x, k}^{(2)}\right)+\frac{i}{4}\left(\gamma_{x+1,-k}^{(1)} \gamma_{x, k}^{(2)}-\gamma_{x+1,-k}^{(2)} \gamma_{x, k}^{(1)}\right) \\ \hline\left(\Delta_{x}\right) & c_{x, k}^{\dagger} c_{x+1,-k}^{\dagger} & \frac{1}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x+1, k}^{(1)}-\gamma_{x,-k}^{(2)} \gamma_{x+1, k}^{(2)}\right)+\frac{i}{4}\left(-\gamma_{x,-k}^{(1)} \gamma_{x+1, k}^{(2)}-\gamma_{x,-k}^{(2)} \gamma_{x+1, k}^{(1)}\right) \\ \hline\left(\Delta_{x}\right)^{*} & c_{x+1,-k} c_{x, k} & \frac{1}{4}\left(\gamma_{x+1,-k}^{(1)} \gamma_{x, k}^{(1)}-\gamma_{x+1,-k}^{(2)} \gamma_{x, k}^{(2)}\right)+\frac{i}{4}\left(\gamma_{x+1,-k}^{(1)} \gamma_{x, k}^{(2)}+\gamma_{x+1,-k}^{(2)} \gamma_{x, k}^{(1)}\right) \\ \hline\left(\Delta_{y} e^{-i k}\right) & c_{x, k}^{\dagger} c_{x,-k}^{\dagger} & \frac{1}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(1)}-\gamma_{x,-k}^{(2)} \gamma_{x, k}^{(2)}\right)+\frac{i}{4}\left(-\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(2)}-\gamma_{x,-k}^{(2)} \gamma_{x, k}^{(1)}\right) \\ \hline\left(\Delta_{y} e^{-i k}\right)^{*} & c_{x,-k} c_{x, k} & \frac{1}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(1)}-\gamma_{x,-k}^{(2)} \gamma_{x, k}^{(2)}\right)+\frac{i}{4}\left(\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(2)}+\gamma_{x,-k}^{(2)} \gamma_{x, k}^{(1)}\right) \\ \hline \end{array} \] Where \(\tilde{\mu}=\mu+t_{y} e^{-i k}+t_{y}^{*} e^{i k}\)
We then consider some special cases. (I) \(\tilde{\mu}=0 ; \Delta_{y}=0 ; \Delta_{x}=t_{x}=t:\) The Hamiltonian can be given by \[
H=i t \sum_{x=1}^{N x} \sum_{k>0}\left(\gamma_{x,-k}^{(2)} \gamma_{x+1, k}^{(1)}+\gamma_{x, k}^{(2)} \gamma_{x+1,-k}^{(1)}\right) .
\] (II) \(\tilde{\mu} \approx \infty ; . \Delta_{x}=t_{x}=t_{y}=1 ; \Delta_{y}=i:\) The Hamiltonian can be given by \[
\begin{aligned}
H & \approx=-\sum_{x} \sum_{k}\left(c_{x, k}^{\dagger} c_{x, k}+c_{x,-k}^{\dagger} c_{x,-k}\right) \\
& \approx-\frac{i}{2} \mu \sum_{x=1}^{N x} \sum_{k>0}\left(\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(2)}+\gamma_{x, k}^{(1)} \gamma_{x,-k}^{(2)}\right)
\end{aligned}
\] (III) \(p_{x}+i p_{y}: \Delta_{x}=t_{x}=t_{y}=1 ; \Delta_{y}=i:\) The Hamiltonian can be given by \[
\begin{array}{c}
H=\sum_{x=1}^{N x} \sum_{k} \frac{i}{2}(-\mu+2 \cos k)\left[\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(2)}+\gamma_{x, k}^{(1)} \gamma_{x,-k}^{(2)}\right]+ \\
i\left[\gamma_{x,-k}^{(2)} \gamma_{x+1, k}^{(1)}+\gamma_{x, k}^{(2)} \gamma_{x+1,-k}^{(1)}\right]+ \\
\sin k\left[\gamma_{x,-k}^{(1)} \gamma_{x, k}^{(1)}+\gamma_{x, k}^{(2)} \gamma_{x,-k}^{(2)}\right]
\end{array}
\] 
Figure 1: The Majorana fermion picture for (a) case (I), (b) case (II) and(c) case (III).

Figure 2: The phase diagram of \(2 d\) chiral superconductor by tuning chemical potential and fixed \(t_{x}=t_{y}=t=1\) and \(\Delta_{x}=-i \Delta_{x}=\Delta=1\).
Here, we note that the \(2 d\) system with chiral p-wave superconductor, and the Hamiltonian can be given as \[ \begin{aligned} H=&\sum_{\vec{r}} C_{\vec{r}}^{\dagger}\left(\begin{array}{rr} t & \Delta \\ -\Delta- & t \end{array}\right) C_{\vec{r}+\hat{x}}+\text { h.c. }+\\& C_{\vec{r}}^{\dagger}\left(\begin{array}{rr} t & \Delta \\ -\Delta- & t \end{array}\right) C_{\vec{r}+\hat{y}}+\text { h.c. }+\\& C_{\vec{r}}^{\dagger}\left(\begin{array}{cc} \mu & 0 \\ 0 & -\mu \end{array}\right) C_{\vec{r}} \end{aligned} \]
(2)1d p-wave SC
A. Continuum model
The Hamiltonian of the \(1 D\) p-wave \(SC\) is given as \[ \begin{aligned} H &=-\mu \sum_{k}\left[\varepsilon_{k} c_{k}^{\dagger} c_{k}+\frac{1}{2}\left(\Delta_{k} c_{k}^{\dagger} c_{-k}^{\dagger}+\Delta_{k}^{*} c_{-k} c_{k}\right)\right] \\ &=\frac{1}{2} \sum_{k}\left(c_{k}^{\dagger} c_{-k}\right)\left(\begin{array}{cc} \varepsilon_{k} & \Delta_{k}^{*} \\ \Delta_{k} & -\varepsilon_{k} \end{array}\right)\left(\begin{array}{c} c_{k} \\ c_{-k}^{\dagger} \end{array}\right) \end{aligned} \] in which \(\varepsilon_{k}=\hbar^{2} k^{2} / 2 m-\mu\) and \(\Delta_{k}=\Delta_{0} k\). The eigenenergies are, \[ E_{\pm}=\pm \sqrt{\varepsilon_{k}^{2}+\left|\Delta_{k}\right|^{2}} \] The eigenstate \(\left(u_{k}, v_{k}\right)\) for energy \(E_{k}\) has the components, \[ u_{k}=\sqrt{\frac{1}{2}\left(1+\frac{\varepsilon_{k}}{E_{k}}\right)} ; v_{k}=\sqrt{\frac{1}{2}\left(1-\frac{\varepsilon_{k}}{E_{k}}\right)} \frac{\Delta_{k}^{*}}{\left|\Delta_{k}\right|} \] Near \(k=0, \varepsilon_{k} \simeq|\mu|\). For both \(\mu<0\) and \(\mu>0\), the spectra are gapful. If \(k=0\), the eigenenergies \(E_{+}(k)\) and \(E_{-}(k)\) touch at \(k=0\).
\(\mu<0\)
If \(\mu<0\), for small \(k, u_{k} \simeq 1\) and \(v_{k} \simeq \frac{\Delta_{k}^{*}}{2|\mu|}\). It is known that the Fourier transform of the Cooper pair wave function \(g(x)\) is \(g(k)=v_{k} / u_{k} .\) So for \(\mu<0, g(k) \propto k\). The function \(g(k)\) being analytic near \(k=0\) implies that its Fourier transform falls off exponentially at large distance, \(g(x) \simeq e^{-x / x_{0}}\). The phase with \(\mu<0\) is thus called the strong-coupling phase.
\(\mu>0\)
On the other hand, if \(\mu>0\), then for small \(k, u_{k} \simeq\) \(\frac{\left|\Delta_{k}\right|}{2 \mu}, v_{k} \simeq 1\) (we have ignored the phase of \(\Delta_{k}\) ), and \(g(k) \propto 1 / k\). The sharp peak near small k implies a slow decay of \(g(x)\) at large distance. So the phase with \(\mu>0\) is called the weak-coupling phase. These two phases cannot be adiabatically connected to each other.
B. Edge state
We now assume \(\mu(x>0)>0\), and \(\mu(x<0)<0\), so the \(1 d\) space is separated to a weak-coupling phase and a strong-coupling phase (e.g., \(\left.\mu(x)=\mu_{0} \tanh x\right)\). Ignore terms of order \(k^{2}\), the \(BdG\) equation is, \[ \left(\begin{array}{cc} -\mu(x) & \Delta_{0} k^{*} \\ \Delta_{0} k & \mu(x) \end{array}\right) \psi=E \psi \] We are only interested in the edge-state solution. Let \(k \rightarrow \partial / i \partial_{x}\) and try \[ \psi(x)=\psi_{0} e^{-\frac{1}{\Delta_{0}} \int_{0}^{x} d x^{\prime} \mu\left(x^{\prime}\right)} \] then \[ \left(\begin{array}{cc} -\mu(x)-E & \Delta_{0} k^{*} \\ \Delta_{0} k & \mu(x)-E \end{array}\right) \psi=0 \] At \(E =0\), we have a solution, \[ \psi_{0}=\frac{1}{\sqrt{2}}\left(\begin{array}{c} 1 \\ -i \end{array}\right) . \] This is a zero mode localized at the interface. The Bogoliubov quasiparticle (QP) operator for the zero mode is, \[ \begin{aligned} \gamma_{0} &=\int d x\left[u^{*}(x) \psi(x)+v^{*}(x) \psi^{\dagger}(x)\right] \\ &=\frac{e^{i \pi / 4}}{\sqrt{2}} \int d x e^{-\frac{1}{\Delta_{0}} \int_{0}^{x} d x^{\prime} \mu\left(x^{\prime}\right)}\left[e^{-i \pi / 4} \psi(x)+e^{i \pi / 4} \psi^{\dagger}(x)\right] \end{aligned} \] Removing the overall phase of \(\pi / 4\), we have \[ \gamma_{0}^{\dagger}=\gamma_{0} \] That is, the anti-particle of the \(QP\) is the same as the QP. Such a fermion is called a Majorana fermion (MF). The zero mode is protected by the PH symmetry, as long as it's non-degenerate.
II. MAJORANA FERMIONS IN KITAEV'S TOY MODEL
A. Hamiltonian
The Kitaev's toy model is a one-dimensional chain of spinless p-wave superconductor, \[ \begin{aligned} H &=-\mu \sum_{x} c_{x}^{\dagger} c_{x}-\frac{1}{2} \sum_{x} t\left(c_{x}^{\dagger} c_{x+1}+c_{x+1}^{\dagger} c_{x}\right) \\ &+\left(\Delta e^{i \phi} c_{x} c_{x+1}+\Delta^{*} e^{-i \phi} c_{x+1}^{\dagger} c_{x}^{\dagger}\right) \end{aligned} \] where \(\mu\) is the chemical potential, \(t \geq 0\) is the nearestneighbor hopping strength, \(\Delta \geq 0\) is the p-wave pairing amplitude and \(\phi\) is the corresponding superconducting phase. For simplicity we set the lattice constant to unity. With the Fourier transformation, \[ c_{x}^{\dagger}=\frac{1}{\sqrt{N}} e^{i k x} c_{k}^{\dagger} \] one have \[ \begin{aligned} H &=\sum_{k}\left[-\mu c_{k}^{\dagger} c_{k}-\frac{t}{2}\left(e^{-i k} c_{k}^{\dagger} c_{k}+e^{i k} c_{k}^{\dagger} c_{k}\right)\right] \\ &+\frac{1}{2} \sum_{k}\left[\Delta e^{i \phi} e^{i k} c_{k} c_{-k}+\Delta^{*} e^{-i \phi} e^{i k} c_{k}^{\dagger} c_{-k}^{\dagger}\right] \end{aligned} \] It is instructive to first understand the chain's bulk properties, which can be conveniently studied by imposing periodic boundary conditions on the system (thereby wrapping the chain into a loop and removing its ends). Upon passing to momentum space and introducing a two- component operator \(C_{k}^{\dagger}=\left(c_{k}^{\dagger}, c_{-k}\right)\), one can write \(H\) in the standard Bogoliubov-de Gennes form: \[ H=\frac{1}{2} \sum_{k \in B Z} C_{k}^{\dagger} H_{k} C_{k}, \quad H(k)=\left(\begin{array}{cc} \xi_{k} & \tilde{\Delta}_{k}^{*} \\ \tilde{\Delta}_{k} & -\xi_{k} \end{array}\right) \] with \(\xi_{k}=-t \cos k-\mu\) the kinetic energy and \(\tilde{\Delta}_{k}=\) \(-i \Delta e^{i \phi} \sin k\) the Fourier-transformed pairing potential. The Hamiltonian (up to a constant) becomes simply \[ H=\sum_{k \in B Z} E_{\text {bulk }}(k) a_{k}^{\dagger} a_{k} \] when expressed in terms of quasiparticle operators \[ \begin{aligned} a_{k} &=u_{k} c_{k}+v_{k} c_{-k}^{\dagger} \\ u_{k} &=\frac{\tilde{\Delta}_{k}}{\left|\tilde{\Delta}_{k}\right|} \frac{\sqrt{E_{ bulk }(k)+\xi_{k}}}{2 E_{ bulk }(k)}, \quad v_{k}=\left(\frac{E_{ bulk }(k)-\xi_{k}}{\Delta_{k}}\right) u_{k} \end{aligned} \] where the bulk excitation energies are given by \[ E_{\text {bulk }}(k)=\sqrt{\xi_{k}^{2}+\left|\tilde{\Delta}_{k}\right|^{2}} \] It demonstrates that the chain admits gapless bulk excitations only when the chemical potential is fine tuned to
\(\mu=t\) or \(-t\), where the Fermi level respectively coincides with the top and bottom of the conduction band. The gap closure at these isolated \(\mu\) values reflects the p-wave nature of the pairing required by Pauli exclusion. More precisely, since \(\tilde{\Delta}_{k}\) is an odd function of \(k\), Cooper pairing at \(k=0\) or \(k=\pm \pi\) is prohibited, thereby leaving the system gapless at the Fermi level when \(\mu=\pm t .\) Note that the phases that appear at \(\mu<-t\) and \(\mu>t\) are related by a particle-hole transformation; thus to streamline our discussion we will hereafter neglect the latter chemical potential range.
In the thermodynamical limit, the difference between two subspaces can be neglected and the ground state energy density can be expressed as \[ \varepsilon_{g}=\frac{1}{4 \pi} \int_{-\pi}^{\pi} E_{ bulk }(k) d k \]
B. The ground state
The physics of the chain is intuitively rather different in the gapped regimes with \(\mu<-t\) and \(|\mu|<t-\) the former connects smoothly to the trivial vacuum (upon taking \(\mu \rightarrow-\infty\) ) where no fermions are present, whereas in the latter a partially filled band acquires a gap due to p-wave pairing. One can make this distinction more quantitative following Read and Green by examining the form of the ground-state wavefunction in each regime. Equation (14) implies that the ground state \(|G S\rangle\) must satisfy \(a_{k}|G S\rangle=0\) for all \(k\) so that no quasiparticles are present. Equations (16) allow one to explicitly write the ground state as follows: \[ \begin{array}{l} |G S\rangle \propto \prod_{0<k<\pi}\left[1+\varphi_{C . p .}(k) c_{-k}^{\dagger} c_{k}^{\dagger}\right]|0\rangle, \\ \varphi_{C . p .}(k)=\frac{v_{k}}{u_{k}}=\left(\frac{E_{ bulk }(k)-\xi_{k}}{\tilde{\Delta}_{k}}\right) \end{array} \] where \(|0\rangle\) is a state with no \(c_{k}\) fermions present with momenta in the interval \(0<|k|<\pi .\) One can loosely interpret \(\varphi_{\text {C.p. }}(k)\) as the wavefunction for a Cooper pair formed by fermions with momenta \(k\) and \(-k .\) An important difference between the \(\mu<-t\) and \(|\mu|<t\) regimes is manifested in the real-space form \(\varphi_{\text {C.p. }}(x)=\) \(\int_{k} e^{i k x} \varphi_{\text {C.p. }}(k)\) at large \(x\) : \[ \left|\varphi_{\text {C.p. }}(x)\right| \sim e^{-|x| / \zeta}, \quad u<-t \text { (strong pairing) } \] const, \(|u|>t\) (weak pairing). It follows that \(\mu<-t\) corresponds to a strong-pairing regime in which 'molecule-like' Cooper pairs form from two fermions bound in real space over a length scale \(\zeta\), whereas in the weak pairing regime \(|u|<t\) the Cooperpair size is infinite. We emphasize that this distinction by itself does not guarantee that the weak and strong pairing regimes constitute distinct phases. Indeed, similar physics occurs in the well-studied 'BCS-BEC crossover? in s-wave paired systems where no sharp transition arises. The fact that the weak and strong pairing regimes are distinct phases separated by a phase transition at which the bulk gap closes is rooted in topology.
C. Topological number
There are several ways in which one can express the 'topological invariant' (akin to an order parameter in the theory of conventional phase transitions) distinguishing the weak and strong pairing phases. Let us revisit the Hamiltonian in equation (13), but now allow for additional perturbations that preserve translation symmetry. The resulting \(2 \times 2\) matrix \(H_{k}\) can be expressed in terms of a vector of Pauli matrices \(\sigma=\sigma_{x} \hat{x}+\sigma_{y} \hat{y}+\sigma_{z} \hat{z}\) as follows: \[ H_{k}=h(k) \cdot \sigma \] for some vector \(h(k)\). A term proportional to the identity can also be added, but will not matter for our purposes. Although we are considering a rather general Hamiltonian here, the structure of \(h(k)\) is not entirely arbitrary. In particular, since the two-component operator \(C_{k}\) in equation \((13)\) satisfies \(\left(C_{-k}^{\dagger}\right)^{T}=\sigma_{x} C_{k}\) the vector \(h(k)\) must obey the important relations \[ h_{x, y}(k)=-h_{x, y}(-k), \quad h_{z}(k)=h_{z}(-k) . \] Thus it suffices to specify \(h(k)\) only on the interval \(0 \leq k \leq \pi\), since \(h(k)\) on the other half of the Brillouin zone. Suppose now that \(h(k)\) is non-zero throughout the Brillouin zone so that the chain is fully gapped. One can then always define a unit vector \(h(k)=h(k) /|h(k)|\) that provides a map from the Brillouin zone to the unit sphere. The relations of equation 21 strongly restrict this map at \(k=0\) and \(\pi\) such that \[ \left.h(0)=s_{0} \hat{z}, \quad h \hat{(\pi}\right)=s_{\pi} \hat{z}, \] where \(s_{0}\) and \(s_{\pi}\) represent the sign of the kinetic energy (measured relative to the Fermi level) at \(k=0\) and \(\pi\), respectively. Thus as one sweeps \(k\) from 0 to \(\pi, h(k)\) begins at one pole of the unit sphere and either ends up at the same pole (if \(\left.s_{0}=s_{\pi}\right)\) or the opposite pole (if \(s_{0}=-s_{\pi}\) ). These topologically distinct trajectories are distinguished by the topological \(Z_{2}\) invariant \[ \nu=s_{0} s_{\pi} \] which can only change sign when the chain?s bulk gap closes (resulting in \(h(k)\) being ill-defined somewhere in the Brillouin zone). Physically, \(\nu=+1\) if at a given chemical potential there exists an even number of pairs of Fermi points, while \(\nu=-1\) otherwise. From this perspective it is clear that \(\nu=+1\) in the (topologically trivial) strong pairing phase while \(\nu=-1\) in the (topologically non-trivial) weak-pairing phase. In our case, \(s_{0}=-t-\mu\) and \(s_{\pi}=t-\mu .\) The origin can be inside the loop only when these two quantities have opposite signs.
D. Topological invariants
Again, the Hamiltonian can be given \[ \begin{aligned} H &=\frac{1}{2} \sum_{k \in B Z} C_{k}^{\dagger} H_{k} C_{k} \\ &=\frac{1}{2} \sum_{k \in B Z}\left(c_{k}^{\dagger} c_{-k}\right)\left(\begin{array}{cc} \xi_{k} & \tilde{\Delta}_{k}^{*} \\ \tilde{\Delta}_{k} & -\xi_{k} \end{array}\right)\left(\begin{array}{c} c_{k} \\ c_{-k}^{\dagger} \end{array}\right) \end{aligned} \] The core matrix can be expressed as \[ H_{k}=h(k) \cdot \sigma \] where the components of the auxiliary field \(h(k)=\) \(\left(h_{x}(k), h_{y}(k), h_{z}(k)\right)\) are \[ \begin{array}{l} h_{x}(k)=\operatorname{Re}\left[\tilde{\Delta}_{k}\right]=0 \\ h_{y}(k)=\operatorname{Im}\left[\tilde{\Delta}_{k}\right]=-\Delta e^{i \phi} \sin k \\ h_{z}(k)=\xi_{k}=-t \cos k-\mu \end{array} \] with \(h(k)\) a smooth function that is non-zero for all momenta so that the bulk is fully gapped. One can then define a unit vector \(\hat{h}(k)\) that maps \(2 D\) momentum space onto a unit sphere. Assuming that \(\hat{h}(k)\) tends to a unique vector as \(|k| \rightarrow \infty\) (independent of the direction of \(k\) ), the number of times this map covers the entire unit sphere defines an integer topological invariant given formally by the Chern number \[ C=\int \frac{d^{2} k}{4 \pi}\left[\hat{h} \cdot\left(\partial_{k_{x}} \hat{h} \times \partial_{k_{y}} \hat{h}\right)\right] \] The integrand above determines the solid angle (which can be positive or negative) that \(\hat{h}(k)\) sweeps on the unit sphere over an infinitesimal patch of momentum space centered on k. Performing the integral over all \(k\) yields an integer that remains invariant under smooth deformations of \(\hat{h}(k)\). The Chern number can change only when the gap closes, making \(h(k)\) ill-defined at some momentum. Note that for momenta with fixed \(|k|, \hat{h}_{x}(k)\) and \(\hat{h}_{y}(k)\) always sweep out a circle on the unit sphere at height \(\hat{h}_{z}(k) .\) As \(|k|\) increases from zero in the \(\mu<0\) strong pairing phase, \(\hat{h}_{z}(k)\) begins at the north pole, descends toward the equator, and then returns to the north pole as \(|k| \rightarrow \infty\).
In \(1 d\) system, the winding number of a closed curve in the auxiliary \(h_{y} h_{z}\) -plane around the origin is defined as \[ C=\frac{1}{2 \pi} \int \hat{h}_{y}(k) d \hat{h}_{z}(k)-\hat{h}_{z}(k) d \hat{h}_{y}(k) \] Actually, the winding number is simply related the loop described by equation \[ \frac{h_{y}^{2}}{\Delta^{2}}+\left(h_{z}+\mu\right)^{2}=1 \] which presents a normal ellipse in the \(h_{y} h_{z}\) -plane.
E. Majorana fermion reresentation
The non-trivial topology inherent in the weak-pairing phase leads to the appearance of Majorana modes in a chain with open boundary conditions, which we will now consider. The usual fermions satisfy the anticommutation relation, \[ \left\{c_{x}, c_{x^{\prime}}^{\dagger}\right\}=\delta_{x, x^{\prime}} \] The new physics associated with the ends of the chain can be most simply accessed by decomposing the spinless fermion operators \(c_{x}\) in the original Hamiltonian in terms of two Majorana fermions via \[ c_{x}=\frac{e^{-i \phi / 2}}{2}\left(\gamma_{B, x}+i \gamma_{A, x}\right) \] The operators on the right-hand side obey the canonical Majorana fermion relations \[ \gamma_{\alpha, x}=\gamma_{\alpha, x}^{\dagger}, \quad\left\{\gamma_{\alpha, x}, \gamma_{\alpha^{\prime}, x^{\prime}}\right\}=2 \delta_{\alpha \alpha^{\prime}} \delta_{x x^{\prime}} \] In this basis \[ \begin{aligned} H=&-\frac{\mu}{2} \sum_{x=1}^{N}\left(1+i \gamma_{B, x} \gamma_{A, x}\right) \\ &-\frac{i}{2} \sum_{x=1}^{N-1}\left[(\Delta+t) \gamma_{B, x} \gamma_{A, x+1}+(\Delta-t) \gamma_{A, x} \gamma_{B, x+1}\right] . \end{aligned} \] Generally the parameters \(\mu, t\) and \(\Delta\) induce relatively complex couplings between these Majorana modes; however, the problem becomes trivial in two limiting cases.
1 \(\mu<0\) but \(t=0\)
The first corresponds to \(\mu<0\) but \(t=0\), where the chain resides in the topologically trivial phase. Here the second line vanishes, leaving a coupling only between Majorana modes \(\gamma_{A, x}\) and \(\gamma_{B, x}\) at the same lattice site. In this case there is a unique ground state corresponding to the vacuum of \(c_{x}\) fermions. Clearly the spectrum is gapped since introducing a spinless fermion into the chain costs a finite energy \(|\mu|\). This is entirely consistent with our treatment of the chain with periodic boundary conditions; in the trivial phase the ends of the chain have little effect. We emphasize that these conclusions hold even away from this fine-tuned limit provided the gap persists so that the chain remains in the same trivial phase.
2 \(\mu=0\) and \(t=\Delta \neq 0\)
The second simplifying limit corresponds to \(\mu=0\) and \(t=\Delta \neq 0\), where the topological phase appears. Here the Hamiltonian is instead given by \[ H=-\frac{i}{2} t \sum_{x=1}^{N-1} \gamma_{B, x} \gamma_{A, x+1} \] which couples Majorana fermions only at adjacent lattice sites. In terms of new ordinary fermion operators \(d_{x}=\) \(\frac{1}{2}\left(\gamma_{A, x+1}+i \gamma_{B, x}\right)\), the Hamiltonian can be written as \[ H=t \sum_{x=1}^{N-1}\left(d_{x}^{\dagger} d_{x}-\frac{1}{2}\right) \] In this form it is apparent that a bulk gap remains here too?consistent with our results with periodic boundary conditions - since one must pay an energy \(t\) to add a \(d_{x}\) fermion. However, the ends of the chain now support 'unpaired' zero-energy Majorana modes \(\gamma_{1}=\gamma_{A, 1}\) and \(\gamma_{N}=\gamma_{B, N}\) that are explicitly absent from the Hamiltonian. These can be combined into an ordinary?though highly non-local?fermion, \[ f=\frac{1}{2}\left(\gamma_{1}+i \gamma_{2}\right) \] that costs zero energy and therefore produces a two-fold ground-state degeneracy. In particular, if \(|0\rangle\) is a ground state satisfying \(f|0\rangle=0\), then \(|1\rangle=f^{\dagger}|0\rangle\) is necessarily also a ground state (with opposite fermion parity). Note the stark difference from conventional gapped superconductors, where typically there exists a unique ground state with even parity so that all electrons can form Cooper pairs.
The appearance of localized zero-energy Majorana end states and the associated ground-state degeneracy arise because the chain forms a topological phase while the vacuum bordering the chain is trivial. (It may be helpful to imagine adding extra sites to the left and right of the chain, with \(\mu<-t\) for those sites so that the strong pairing phase forms there.) These phases cannot be smoothly connected, so the gap necessarily closes at the chain's boundaries. Because this conclusion has a topological origin it is very general and does not rely on the particular fine-tuned limit considered above, with one caveat.
3 In general
In the more general situation with \(\mu=0\) and \(t \neq 0\) (but still in the topological phase) the Majorana zero-modes \(\gamma_{1}\) and \(\gamma_{2}\) are no longer simply given by \(\gamma_{A, 1}\) and \(\gamma_{B, N}\). Rather, their wavefunctions decay exponentially into the bulk of the chain on a length scale \(\zeta\) given by the superconducting coherence length in the \(1 D\) system (which diverges at the transition to the trivial phase). The overlap of these wavefunctions results in a splitting of the degeneracy between \(|0\rangle\) and \(|1\rangle\) by an energy that scales like \(e^{-L / \zeta}\), where \(L\) is the length of the chain. Provided \(L \gg \zeta\), however, this splitting can easily be negligible compared