Exact Solution to Haldane-BCS-Hubbard Model Along the Symmetric Lines: Interaction Induced Topological Phase Transition

TL;DR

None

Contribution

None

Abstract

We propose a Haldane-BCS-Hubbard model on a honeycomb lattice, which is composed of two copies of the Haldane model of the quantum anomalous Hall effect, an equal-spin pairing term and an onsite Hubbard interaction term. For any interaction strength, this model is exactly solvable along the symmetric line where the hopping and pairing amplitudes are equal to each other. The ground state of the Haldane-BCS-Hubbard model is a topological superconducting state at weak interaction with two chiral Majorana edge states. A strong interaction drives the system across a topological quantum phase transition to a topologically trivial superconductor. A $\mathbb{Z}_{2}$ symmetry of the Hamiltonian, which is a composition of the bond-centered inversion and a gauge transformation, is spontaneously broken by the interaction, resulting a finite antiferromagnetic order in the $y$-direction.

Equations

$$H=H_{\mathrm{hop}}+H_{\mathrm{pair}}+H_{\mathrm{int}}$$

$$H_{\mathrm{hop}}$$

$$H_{\mathrm{pair}}$$

$$H_{\mathrm{int}}$$

$$t_1={\Delta }_1,t_2={\Delta }_2.$$

$$c_{rsA}$$

$$c_{rsB}$$

$$H_0$$

$$-{\tilde{t}}_1\sum _{rs}\left(i{\gamma }_{rsA}{\gamma }_{rsB}+i{\gamma }_{rsA}{\gamma }_{r+a_{1}sB}+i{\gamma }_{rsA}{\gamma }_{r+a_{2}sB}\right)$$

$$-{\delta }_2\sum _{rs\lambda }\lambda \left(i{\eta }_{rs\lambda }{\eta }_{r+a_{1}s\lambda }+i{\eta }_{rs\lambda }{\eta }_{r-a_1+a_{2}s\lambda }+i{\eta }_{rs\lambda }{\eta }_{r-a_{2}s\lambda }\right)$$

$$-{\tilde{t}}_2\sum _{rs\lambda }\lambda \left(i{\gamma }_{rs\lambda }{\gamma }_{r+a_{1}s\lambda }+i{\gamma }_{rs\lambda }{\gamma }_{r-a_1+a_{2}s\lambda }+i{\gamma }_{rs\lambda }{\gamma }_{r-a_{2}s\lambda }\right)$$

$$H_{\mathrm{int}}$$

$$H_{\mathrm{int}}=-U\sum _{r\lambda }{D}_{r\lambda }i{\gamma }_{r\uparrow \lambda }{\gamma }_{r\downarrow \lambda }$$

$$d_{r2\lambda }$$

$$d_{r1\lambda }$$

$$c_{r\pm \lambda }^{†}=\frac{1}{\sqrt{2}}\left(c_{r\uparrow \lambda }^{†}\pm i{c}_{r\downarrow \lambda }^{†}\right)$$

$$d_{r2A}=\frac{c_{r-A}+c_{r+A}^{†}}{\sqrt{2}},$$

$$d_{r1A}=\frac{c_{r+A}-c_{r-A}^{†}}{\sqrt{2}i},$$

$$H_0$$

$$-\frac{i{\tilde{t}}_1}{2}\sum _r\left(d_{r1A}^{†}{d}_{r1B}^{†}+d_{r1A}^{†}{d}_{r+a_{1}1B}^{†}+d_{r1A}^{†}{d}_{r+a_{2}1B}^{†}\right)$$

$$-\frac{i{\delta }_2}{2}\sum _{r\lambda }\lambda \left(d_{r2\lambda }^{†}{d}_{r+a_{1}2\lambda }+d_{r2\lambda }^{†}{d}_{r-a_1+a_{2}2\lambda }+d_{r2\lambda }^{†}{d}_{r-a_{2}2\lambda }\right)$$

$$-\frac{i{\tilde{t}}_2}{2}\sum _{r\lambda }\lambda \left(d_{r1\lambda }^{†}{d}_{r+a_{1}1\lambda }+d_{r1\lambda }^{†}{d}_{r-a_1+a_{2}1\lambda }+d_{rs\lambda }^{†}{d}_{r-a_{2}1\lambda }\right)$$

$$+h.c.$$

$$H_{\mathrm{int}}$$

$$d_{r\alpha \lambda }=\frac{1}{\sqrt{N}}\sum _k{e}^{ik\cdot r}{d}_{k\alpha \lambda }$$

$$H_{\alpha }=\sum _k{\psi }_{k\alpha }^{†}{h}_{\alpha }\left(k\right){\psi }_{k\alpha }$$

$$h_{\alpha }\left(k\right)={\vec{T}}_{\alpha }\left(k\right)\cdot \vec{\sigma }$$

$$T_1^x\left(k\right)$$

$$T_1^y\left(k\right)$$

$$T_1^z\left(k\right)$$

$$C_1$$

$$C_2$$

$$C$$

$$C_{\mathrm{spin}}$$

$$H_s$$

$$-\frac{i{\tilde{t}}_2}{2}\sum _{r\lambda }\lambda \left(d_{r1\lambda }^{†}{d}_{r+a_{1}1\lambda }+d_{r1\lambda }^{†}{d}_{r-a_1+a_{2}1\lambda }+d_{r1\lambda }^{†}{d}_{r-a_{2}1\lambda }\right)+h.c.$$

$$\pm \frac{U}{2}\sum _{r\lambda }\left(n_{r1\lambda }-\frac{1}{2}\right)$$

$$H_s=\sum _k{\psi }_{k1}^{†}{h}_s\left(k\right){\psi }_{k1}$$

$$h_s\left(k\right)={\vec{T}}_s\left(k\right)\cdot \vec{\sigma }$$

$$C_s=\frac{1}{2}\left[\mathrm{sign}\left(3\sqrt{3}{\tilde{t}}_2-U\right)+\mathrm{sign}\left(3\sqrt{3}{\tilde{t}}_2+U\right)\right]$$

Full text

I Introduction

The concept of topology in condensed matter physics has flourished in the past decadesHasan and Kane (2010); Qi and Zhang (2011). This abstract notion is deeply related to the band structure in the momentum space. The topological band theory has been established and a lot of predicted materials have been synthesizedBansil et al. (2016). Recently the full diagnosis of the non-trivial band topology for non-magnetic materials have been establishedZhang et al. (2019); Vergniory et al. (2019); Tang et al. (2019).

The interplay of topology and correlations can lead to novel phases and phase transitions in condensed matter systems. First, the interactions may reduce the topological classification of free fermions in one dimensionFidkowski and Kitaev (2010, 2011) and two dimensionsYao and Ryu (2013). Second, interactions may drive topological quantum phase transitions, which is demonstrated in exactly solvable models of interacting Kitaev chains Miao et al. (2017); Ezawa (2017); Wang et al. (2017), the Haldane-Hubbard model Zheng et al. (2015) and the $\mathbb{Z}_{2}$ Bose-Hubbard model González-Cuadra et al. (2019). Recently, Chen et. al.Chen et al. (2018) generalized the construction of the Kitaev honeycomb modelKitaev (2006) to spinful fermion models with both equal-spin pairing and Hubbard interaction terms, dubbed BCS-Hubbard model, which can be solved exactly when the pairing amplitude equals the hopping amplitude. Later EzawaEzawa (2018) generalized the BCS-Hubbard model on a honeycomb lattice by introducing the Kane-Mele spin-orbit coupling (SOC). However, an infinitesimal Hubbard interaction $U$ will destroy the topological superconducting state due to the spontaneous time reversal symmetry breaking in Ref. Ezawa, 2018. It is still desirable to find an exactly solvable model in two dimensions with topological phase transition at finite interaction strength to study the interplay of topology and correlations.

In this paper, we investigate the Haldane-BCS-Hubbard model on a honeycomb lattice. Along the symmetric lines where the hopping amplitude equals the pairing amplitude, the model is exactly solvable and reduces to the Falicov-Kimball modelFalicov and Kimball (1969). There is an interaction induced topological phase transition at finite Hubbard $U$ along the symmetric lines. The phase transition can be characterized by the change of the spectral Chern number. Thus the topological superconducting state in our model is stable to small interaction. These results are obtained exactly without approximation, and can serve as a benchmark for further study.

The paper is organized as follows: In section II, we introduce the Haldane-BCS-Hubbard model. Then we show the exact solvability of the model along the symmetric lines in section III. We analyze the symmetry of the model in section IV and introduce the composite fermion representation in section V for later convenience. In section VI, we study the noninteracting limit of the model and give the phase diagram. In section VII, we study the model along the symmetric lines and show the interaction induced topological phase transition. We summarize the results and propose the possible realization of the model in section VIII.

II Model Hamiltonian

In this section, we introduce the Haldane-BCS-Hubbard model we study. The Hamiltonian of the model consists of three parts and can be expressed as follows

[TABLE]

where $H_{\mathrm{hop}}$ describes the electron hopping terms, which is a spinful generalization of the Haldane modelHaldane (1988); Sheng et al. (2005), $H_{\mathrm{pair}}$ describes the equal spin pairing (ESP) terms, and $H_{\mathrm{int}}$ describes the on-site Hubbard interaction. They are given by

[TABLE]

where $c_{is}$ ($c_{is}^{\dagger}$) annihilates (creates) an electron at site $i$ with spin $s=\uparrow,\downarrow$ pointing in the $z$-direction. $\left\langle ij\right\rangle$ and $\left\langle\left\langle ij\right\rangle\right\rangle$ denote the nearest and next-nearest neighbor sites. $t_{1}$ is the nearest neighbor hopping matrix element and $t_{2}$ is the next-nearest neighbor hopping (NNNH) matrix. The sign \nu_{ij}=\text{sign}\bigl{(}\hat{d}_{1}\times\hat{d}_{2}\bigr{)}_{z}=\pm 1, where $\hat{d}_{1,2}$ are the vectors along the two bonds constituting the next-nearest neighbors. The signs $\nu_{ij}$ are shown in FIG. 1. $\Delta_{1}$ and $\Delta_{2}$ are the nearest and next-nearest neighbor ESP potential respectively. $\lambda_{i}=\pm 1$ for sublattices $A$ and $B$ respectively. The signs of $\Delta_{1,2}$ are shown in FIG. 1. $n_{is}=c_{is}^{\dagger}c_{is}$ is the fermionic number operator for spin $s$. $U$ is the strength of on-site Hubbard interaction.

III Exact solvability

In this section, we shall show the exact solvability of the Haldane-BCS-Hubbard model along the symmetric lines

[TABLE]

The Haldane-BCS-Hubbard model is not exactly sovable in general. However, similar to the BCS-Hubbard modelChen et al. (2018), we find this model can be solved exactly along the symmetric lines. The exact solvability of the model becomes manifest in the Majorana fermion representation. As the system contains two sublattices, we use $r$ to denote the unit cell and $c_{rs\lambda}$ to denote the annihilation operator with spin $s$ at unit cell $r$ in sublattice $\lambda=A,B$. We then decompose the complex fermion operators $c_{rs\lambda}$ into Majorana fermion operators $\eta_{rs\lambda}$ and $\gamma_{rs\lambda}$ as follows

[TABLE]

Note the decomposition is opposite for the two sublattices. The Hamiltonian in the Majorana fermion representation becomes

[TABLE]

where $H_{0}=H_{\mathrm{hop}}+H_{\mathrm{pair}}$ and $\delta_{1}=2\left(t_{1}-\Delta_{1}\right)$, $\delta_{2}=2\left(t_{2}-\Delta_{2}\right)$, $\tilde{t}_{1}=2\left(t_{1}+\Delta_{1}\right)$, $\tilde{t}_{2}=2\left(t_{2}+\Delta_{2}\right)$. $a_{1}=a\left(\frac{1}{2},\frac{\sqrt{3}}{2}\right)$ and $a_{2}=a\left(-\frac{1}{2},\frac{\sqrt{3}}{2}\right)$ are the two basis vectors. Note the $\eta$ Majorana fermions disappear in the noninteracting Hamiltonian $H_{0}$ along the symmetric lines $\delta_{1}=\delta_{2}=0$. In the following, we shall focus on the symmetric line. We define $\hat{D}_{r\lambda}=4i\eta_{r\uparrow\lambda}\eta_{r\downarrow\lambda}$. It is easy to prove that $\left[\hat{D}_{r\lambda},H\right]=0$. Thus $\hat{D}_{r\lambda}$ are constants of motion. Since $\hat{D}_{r\lambda}^{2}=1$, we can replace the operators $\hat{D}_{r\lambda}$ by its eigenvalues $D_{r\lambda}=\pm 1$. The Hubbard interaction becomes

[TABLE]

The total Hilbert space is divided into different sectors characterized by $\left\{D_{r\lambda}\right\}$. Within each sector, the Hamiltonian contains only quadratic terms of $\gamma$ Majorana fermions and can be solved exactly.

IV symmetry

Symmetry plays an important role in the following analysis. In this section, we shall analyze the various symmetry of the Haldane-BCS-Hubbard model. The fermion operators transform as $Cc_{is}C^{-1}=\lambda_{i}c_{is}^{\dagger}$ under the particle-hole symmetry (PHS). It is obvious that the Hamiltonian has the PHSChiu et al. (2016). The time reversal symmetry (TRS) operator for spinful system is $T=i\sigma^{y}K$, where $K$ denotes the complex conjugation and $T^{2}=-1$. The fermion operators transform as $Tc_{is}T^{-1}=\sum_{s^{\prime}}\left(i\sigma^{y}\right)_{ss^{\prime}}c_{is^{\prime}}$ under TRS. Just as in the Haldane model, the NNNH terms break the TRS explicitly. The sublattice symmetry (SLS) can be implemented by the bond centered inversion operator $I$. The signs shown in FIG. 1 indicate $H_{\mathrm{pair}}$ breaks the SLS. With NNNH and ESP terms, the Hamiltonian does not preserve the $SU\left(2\right)$ spin rotation symmetry. The $SU\left(2\right)$ symmetry is reduced to $U\left(1\right)_{y}\rtimes\mathbb{Z}_{2}$, where $U\left(1\right)_{y}$ is the rotation about the $y$-axis and $\mathbb{Z}_{2}$ is the $\pi$-rotation around the $z$-axis. Therefore, the system falls into class $D$ in the topological classification of superconductors (SC) Schnyder et al. (2008).

V composite fermion representation

In this section, we introduce the composite fermion representation. These composite fermions form the quasiparticles for the Haldane-BCS-Hubbard model. We define the composite fermions as in Ref. Chen et al., 2018

[TABLE]

The physical meaning of the composite fermions becomes clear by introducing the fermion operators pointing in the $\pm y$-direction as follows

[TABLE]

We express $d$ composite fermions in terms of $c$ fermion operators

[TABLE]

which are equal-weight superposition of particle and hole of $c$ fermion operators. $d_{r1A}$ and $d_{r2B}$ ($d_{r2A}$ and $d_{r1B}$) carry spin-$1/2$ pointing in the $y$($-y$)-direction. We write the Hamiltonian in the composite fermion representation

[TABLE]

where $n_{r\alpha\lambda}=d_{r\alpha\lambda}^{\dagger}d_{r\alpha\lambda}$ with $\alpha=1,2$. Thus the original system can be viewed as two species of $d$ composite fermions with nearest neighbor pairing, next-nearest neighbor hopping, and they interact with on-site Hubbard $U$. Note the Hamiltonian $H$ has the dual symmetry under the dual mapping $d_{r1\lambda}\leftrightarrow d_{r2\lambda}$ (or $\eta_{rs\lambda}\leftrightarrow\gamma_{rs\lambda}$), with parameters changing as $\delta_{1}\leftrightarrow-\tilde{t}_{1}$, $\delta_{2}\leftrightarrow\tilde{t}_{2}$. Thus the Hamiltonian $H$ has a self-dual point $t_{1}=\Delta_{2}=0$, even with the Hubbard interaction $U$. In the following, we shall analyze the properties of the Haldane-BCS-Hubbard model in terms of $d$ composite fermions.

VI Noninteracting limit: Haldane-BCS model

In this section, we analyze the noninteracting limit of the Haldane-BCS-Hubbard model. At $U=0$, the model reduces to the Haldane-BCS model. Note the Hamiltonian $H_{0}$ is decoupled for two species of $d$ composite fermions $H_{0}=H_{1}+H_{2}$, where $H_{\alpha}$ contains $d_{\alpha}$ composite fermions only. The Hamiltonian $H_{0}$ is uniform and we can perform the Fourier transformation to obtain the spectrum. The Fourier transformation is defined as

[TABLE]

We also define the spinor as \psi_{k\alpha}^{\dagger}=\bigl{(}d_{k\alpha A}^{\dagger},d_{-k\alpha B}\bigr{)}, the Hamiltonian $H_{\alpha}$ can be written in the form of

[TABLE]

where

[TABLE]

$\vec{\sigma}=\left(\sigma^{x},\sigma^{y},\sigma^{z}\right)$ are the Pauli matrices and $\vec{T}_{1}\left(k\right)$ are given by

[TABLE]

where $e_{1}=a\left(0,-\frac{1}{\sqrt{3}}\right)$, $e_{2}=a\left(\frac{1}{2},\frac{1}{2\sqrt{3}}\right)$ and $e_{3}=a\left(-\frac{1}{2},\frac{1}{2\sqrt{3}}\right)$ are the three vectors of nearest neighbor bonds. $\vec{T}_{2}\left(k\right)$ can be obtain by the dual mapping $d_{r1\lambda}\rightarrow d_{r2\lambda}$ with parameters changing as $\tilde{t}_{1}\rightarrow-\delta_{1}$ and $\tilde{t}_{2}\rightarrow\delta_{2}$. The energy dispersions read E_{\alpha}\left(k\right)=\pm\bigl{|}\vec{T}_{\alpha}\left(k\right)\bigr{|}, which form reflects the PHS of the Hamiltonian. The ground state is unique with all the negative energy levels of both $d$ composite fermions are occupied. The system is gapped for nonzero $\tilde{t}_{1},\tilde{t}_{2}$ and $\delta_{1},\delta_{2}$.

The noninteracting Hamiltonian $H_{0}$ describes two components Haldane model with ESP at half-filling. According to the symmetry analysis in section IV, the system falls into class $D$Schnyder et al. (2008) of topological superconductor (TSC). The topological invariant is given by the Chern numberThouless et al. (1982) $C_{\alpha}$, which denotes the Chern number of $d_{\alpha}$ composite fermions. We calculate the Chern numbers and find

[TABLE]

where $\mathrm{sign}(x)=\lim_{\epsilon\rightarrow 0}\frac{x}{\sqrt{x^{2}+\epsilon^{2}}}$ is the sign function. We introduce the total Chern number and spin Chern numberSheng et al. (2003, 2006) as

[TABLE]

which indicates the topological phase transition at $t_{2}=\pm\Delta_{2}$. Accordingly the gap closes for $d_{1}$ ($d_{2}$) composite fermions at $t_{2}=-\Delta_{2}$ ($t_{2}=\Delta_{2}$). The phase diagram of Haldane-BCS model is shown in FIG. 2. For $\left|t_{2}\right|>\left|\Delta_{2}\right|$, the system is in the chiral TSC state with total Chern number $C=\pm 2$ and spin Chern number $C_{\mathrm{spin}}=0$. For $\left|t_{2}\right|<\left|\Delta_{2}\right|$, the system is in the helical TSC state with total Chern number $C=0$ and spin Chern number $C_{\mathrm{spin}}=\pm 2$. Along the critical lines $t_{2}=\pm\Delta_{2}$, one species of $d$ composite fermions is gapless and another species is in the chiral TSC state with Chern number $C_{\alpha}=\pm 1$. The origin is a gapless multicritical point.

The topological phase transition can be understood via the bulk-edge correspondence. Except along the critical lines, each species of $d$ composite fermions has nonzero Chern number, i.e. in the TSC state with single chiral edge state. For $\left|t_{2}\right|>\left|\Delta_{2}\right|$, both edge states carry the same chirality and the system is in the chiral TSC state with two chiral edge states, which is consistent with total Chern number $C=\pm 2$ and spin Chern number $C_{\mathrm{spin}}=0$. However for $\left|t_{2}\right|<\left|\Delta_{2}\right|$, two edge states have opposite chirality. The system becomes a helical TSC state with total Chern number $C=0$ while spin Chern number $C_{\mathrm{spin}}=\pm 2$. We plot the energy spectrum of $d_{2}$ composite fermions with different sign of $\delta_{2}$ in Fig. 3 and find in both cases the system has single chiral edge state on each edge. Due to sign change of $\delta_{2}$, the wavefunctions of the edge states localize on opposite edges, which indicates the chirality of the edge states is changed. For comparision, we also show the energy spectrum and wavefunctions of $d_{2}$ composite fermions in Fig. 4 . This is consistent with the sign change of Chern number of $d_{2}$ composite fermions. Along the critical lines $t_{2}=\pm\Delta_{2}$, one chiral edge state merges into the bulk and the system becomes a gapless TSC state with single chiral edge state. The topological phase transition can also be revealed by another dual mapping $\eta_{\lambda}\rightarrow\lambda\eta_{\bar{\lambda}}$ ($\gamma_{\lambda}\rightarrow\lambda\gamma_{\bar{\lambda}}$), where $\bar{\lambda}$ is the different sublattice of $\lambda$. The Hamiltonian has the dual symmetry with parameters changing as $t_{2}\leftrightarrow\Delta_{2}$ ($t_{2}\leftrightarrow-\Delta_{2}$). The topological phase transition happens exactly along the self-dual lines $t_{2}=\pm\Delta_{2}$. Similar duality relating topological and trivial phases has been discovered in the interacting Kitaev chainMiao et al. (2017). If we employ the BdG formalismQi et al. (2010) and use the Nambu spinor \Psi_{k\alpha}^{\dagger}=\bigl{(}d_{k\alpha A}^{\dagger},d_{k\alpha B}^{\dagger},d_{-k\alpha A},d_{-k\alpha B}\bigr{)}, the above analysis is still valid except the Chern numbers should be multiplied by $2$ and each chiral edge state becomes two chiral Majorana edge states.

VII Haldane-BCS-Hubbard model along symmetric lines

In this section, we analyze the Haldane-BCS-Hubbard model along the symmetric lines $\delta_{1}=\delta_{2}=0$. In terms of $d$ composite fermions language, the $d_{2}$ fermions are completely localized (or form the completely flat bands in the band theory languageEzawa (2017)). The Hamiltonian effectively reduces to the Falicov-Kimball model with only one species of mobile $d$ composite fermions

As the total Hilbert space is divided into different sectors characterized by the sets of $\left\{D_{r\lambda}\right\}$, we first determine the ground state sector. Within each sector, the ground state energy is by summing all the negative energy levels. The ground state sector is determined by the set of $\left\{D_{r\lambda}\right\}$ with minimal ground state energy. We traverse all the $2^{N}$ sectors numerically for small lattice size and find the ground state sectors are $\left\{D_{r\lambda}=\pm\lambda\right\}$. We also note the sectors $\left\{D_{r\lambda}=\pm 1\right\}$ have the maximal ground state energy. For large lattice size, we randomly choose the sector $\left\{D_{r\lambda}\right\}$ and find its ground state energy always falls between the sectors $\left\{D_{r\lambda}=\pm\lambda\right\}$ and $\left\{D_{r\lambda}=\pm 1\right\}$. The numerical details are given in the appendix A. Thus we conclude the ground state sectors are $\left\{D_{r\lambda}=\pm\lambda\right\}$. The ground states are uniform with two-fold degeneracy.

Within the ground state sectors, the Hamiltonian of the Haldane-BCS-Hubbard model along the symmetric lines reduces to

[TABLE]

where D_{r\lambda}=2\lambda\bigl{(}n_{r2\lambda}-\frac{1}{2}\bigr{)}=\pm\lambda. The $d_{2}$ composite fermions form the background $\mathbb{Z}_{2}$ charge fields. Note the Hamiltonian is symmetric with respect to the Hubbard $U$ along the symmetric lines. As the ground state sectors are translation invariant, we can perform the Fourier transformation and the Hamiltonian $H_{s}$ can also be written in the form of

[TABLE]

where

[TABLE]

with $T_{s}^{x}=T_{1}^{x}$, $T_{s}^{y}=T_{1}^{y}$ and $T_{s}^{z}=T_{1}^{z}\pm\frac{U}{2}$. The energy dispersion reads E_{s}\left(k\right)=\pm\bigl{|}\vec{T}_{s}\left(k\right)\bigr{|}, which is gapped except at $U=\pm 3\sqrt{3}\tilde{t}_{2}$. The quasiparticle excitations are the spin-$1/2$ $d_{1}$ composite fermions. Even with the Hubbard interactions, we can define the spectral Chern number in terms of these quasiparticles along the symmetric lines. The spectral Chern number is given by

[TABLE]

Thus there is a topological phase transition at $U=\pm 3\sqrt{3}\tilde{t}_{2}$ and the gap closes at this point accordingly.

This topological phase transition can be understood easily in terms of $d_{1}$ composite fermions. Within each sector, the Hubbard interactions act as chemical potential terms. For small $U$, the system is in the weak pairing region and topological. While for large $U$, the system is in the strong pairing region and becomes topologically trivialRead and Green (2000); Qi and Zhang (2011). The topological phase transition is due to the competition between the NNNH terms and Hubbard interactions. This mechanism is remarkably different from the Kane-Mele-BCS-Hubbard model studied in Ref. Ezawa, 2017, where an infinitesimal $U$ renders the topological SC state into trivial, because its topological SC state is protected by the TRS, and the Hubbard interaction always spontaneously breaks the TRS and mixes different spin components within each sector as $Ti\gamma_{\uparrow\lambda}\gamma_{\downarrow\lambda}T^{-1}=-i\gamma_{\uparrow\lambda}\gamma_{\downarrow\lambda}$. In the Haldane-BCS-Hubbard model along the symmetric lines, the topological phase transition happens at finite $U$, which clearly manifests the competition of topology and correlations.

We study the properties of ground states with the aid of symmetry analysis. Even though the Hamiltonian does not have the inversion symmetry $I$, we note it has the combined symmetry $\tilde{I}$ of bond centered inversion $I$ plus gauge transformation $c_{rs\lambda}\rightarrow\sum_{s^{\prime}}\left(i\sigma^{z}\right)_{ss^{\prime}}c_{rs^{\prime}\lambda}$ . The two degenerate ground states are transformed to each other by the symmetry $\tilde{I}$. Thus the ground states spontaneously break the $\mathbb{Z}_{2}$ symmetry $\tilde{I}$ for nonzero $U$. We define the transverse magnetism in the $y$-direction as the order parameterChen et al. (2018)

[TABLE]

We calculate the transverse magnetism via the operator identity

[TABLE]

For comparison, we find $\left\langle d_{r\alpha\lambda}^{\dagger}d_{r\alpha\lambda}\right\rangle=\frac{1}{2}$ in the noninteracting limit with generic hopping and pairing amplitudes, thus the ground state is nonmagnetic. For nonzero $U$ along the symmetric lines, we have

[TABLE]

The order parameter is given by

[TABLE]

Thus the ground states have antiferromagnetic order for nonzero $U$. The transverse magnetism shown in FIG. 5 indicates the $\mathbb{Z}_{2}$ symmetry $\tilde{I}$ is spontaneously breaking. In the limit $U\rightarrow\infty$, there is only one electron per site and the spin is fully polarized. Accordingly we have $m_{r\lambda}^{y}\rightarrow\pm\frac{\lambda}{2}$. In the limit $U\rightarrow-\infty$, each site is either empty or doubly occupied, thus we have $m_{r\lambda}^{y}\rightarrow 0$. We give a remark on the nonzero magnetism at $U=0$ in FIG. 5. For generic hopping and pairing amplitudes, the ground state is nonmagnetic for $U=0$. While along the symmetric lines, one species of composite fermion is completely flat. So it is possible to form the nonzero magnetism by the linear combination of these localized states and the curve is continuous at $U=0$. We note similar magnetic topological phase for small Hubbard $U$ has been found beforeYoshida et al. (2013).

VIII Summary and discussions

In this paper, we study the Haldane-BCS-Hubbard model. We find this model can be solved exactly along the symmetric lines. In the noninteracting limit, the Haldane-BCS model has topological phase transitions at the self-dual points. The topological phase transition are revealed by the bulk-edge correspondence. Along the symmetric lines, we find the model reduces to the Falicov-Kimball model. There is an interaction induced topological phase transition due to the competition between NNNH terms and Hubbard interaction. With nonzero Hubbard $U$, the ground states spontaneously break the $\mathbb{Z}_{2}$ symmetry and have staggered transverse magnetism in the $y$-direction.

The Haldane model has already been realized in the cold atoms systemJotzu et al. (2014). Actually we can view our model as bilayer of Haldane models. The spin index $s$ can be viewed as the layer index with $s=\uparrow$ for the upper layer and $s=\downarrow$ for the bottom layer. The ESP and the on-site interaction Hubbard $U$ between two layers might be introduced in cold atom systems. Therefore, we expect the interaction induced topological phase transition can be observed in cold atom systems.

Let us consider the topological characterization of the system with the Hubbard interaction for generic hopping and pairing parameters. The total Chern number $C$ and the spin Chern number $C_{\mathrm{spin}}$ can be defined with the many-body ground-state wavefunctions in the twisted boundary conditionNiu et al. (1985); Sheng et al. (2003). However, the Chern numbers $C_{\alpha}(\alpha=1,2)$ for each species of $d_{\alpha}$ composite fermions are not well-defined because the two species are entangled with each other. An exception is along the symmetric lines, where the $d_{2}$ fermions are completely localized, and the ground state wavefunction is given by $|\Psi\rangle=|\Psi_{1}(\{D_{r\lambda}\})\rangle\otimes|\Psi_{2}(\{D_{r\lambda}\})\rangle$, in which $\{D_{r\lambda}=\pm\lambda\}$ is selected by the interaction at the ground state, and $|\Psi_{1}(\{D_{r\lambda}\})\rangle$ and $|\Psi_{2}(\{D_{r\lambda}\})\rangle$ denote the wavefunctions of the $d_{1}$ and $d_{2}$ fermions for a given $\{D_{r\lambda}\}$ sector. The localized $d_{2}$ fermions do not respond to the twisted boundary condition, thus $C_{2}=0$. The $d_{1}$ fermions contribute the spectral Chern number given by Eq. (22). Therefore, the total Chern number and the spin Chern number are given by $C=C_{\mathrm{spin}}=C_{1}$ along the symmetric lines for nonzero $U$, which persist for nonzero $\delta_{1}$ and $\delta_{2}$ as long as the bulk gap is not closed due to the topological stability. For large $\delta_{1}$ and $\delta_{2}$, we expect topological phase transitions to phases that are adiabatically connected to the two gapped phases in the noninteracting case. Therefore, an interaction-induced new topological phaseQi et al. (2010); Wang et al. (2015) emerges for $0<|U|<3\sqrt{3}\tilde{t}_{2}$ intervening the two phases in the noninteracting case. We plot a schematic phase diagram in FIG.6.

IX Acknowledgement

J.J.M. acknowledges the discussion with Tai-Kai Ng and Yi Zhou. J.J.M. is supported by China Postdoctoral Science Foundation (Grant No.2017M620880) and the National Natural Science Foundation of China (Grant No.1184700424). D.H.X. is supported by the National Natural Science Foundation of China (Grant No. 11704106) and the Scientific Research Project of Education Department of Hubei Province (Grant No. Q20171005). D.H.X. also acknowledges the support of the Chutian Scholars Program in Hubei Province. L.Z. is supported by National Key R&D Program of China (No. 2018YFA0305800) and National Natural Science Foundation of China (No. 11804337). Work at UCAS is also supported by Strategic Priority Research Program of CAS (No. XDB28000000), and Beijing Municipal Science & Technology Commission (No. Z181100004218001). F.C.Z. is supported by National Science Foundation of China (Grant No.11674278) and National Basic Research Program of China (No.2014CB921203).

X Note added

During the preparation of this work, we learned a similar work on arXivLi et al. (2019), which also studied the extension of the BCS-Hubbard modelChen et al. (2018) (now named as Majorana Falicov-Kimball ModelProsko et al. (2017)) more thorough.

Appendix A Numerical determination of ground state sectors $\left\{D_{r\lambda}\right\}$

There are $2^{N}$ sectors characterized by the sets of $\left\{D_{r\lambda}\right\}$, where $N$ is the number of total sites. Up to $N=4\times 4$, we can traverse all the $2^{N}$ sectors numerically on laptop in one minute. By sorting all the sectors according to the ground state energy, we find the ground state sectors are $\left\{D_{r\lambda}=\pm\lambda\right\}$ for arbitrary strength $U$. We also note the sectors with the largest ground state energy are $\left\{D_{r\lambda}=\pm 1\right\}$. For larger lattice size, the time and internal storage cost increase exponentially and it is impossible to traverse all the $2^{N}$ sectors numerically on laptop. So we randomly choose the sector $\left\{D_{r\lambda}\right\}$, i.e. the value of $D_{r\lambda}$ on each site is $1$ or $-1$ with equal weight, and calculate its ground state energy. We find the ground state energy of randomly chosen sectors always falls between the sectors $\left\{D_{r\lambda}=\pm\lambda\right\}$ and $\left\{D_{r\lambda}=\pm 1\right\}$. For $N=16\times 16$, we randomly choose $1000$ configurations and plot their ground state energy in FIG. 7