Ferromagnetism in d-Dimensional SU(n) Hubbard Models with Nearly Flat Bands

Abstract

We present rigorous results for the SU(n) Fermi–Hubbard models with finite-range hopping in d ($\ge 2$) dimensions. The models are defined on a class of decorated lattices. We first study the models with flat bands at the bottom of the single-particle spectrum and prove that the ground states exhibit SU(n) ferromagnetism when the number of particles is equal to the number of unit cells. We then perturb the models by adding particular hopping terms and make the bottom bands dispersive. Under the same filling condition, it is proved that the ground states remain SU(n) ferromagnetic when the bottom bands are sufficiently flat and the Coulomb repulsion is sufficiently large.

1 Introduction

Recent advances in experimental techniques have made it possible to simulate various quantum systems by using ultracold atoms in optical lattices [1,2,3,4]. Thanks to the controllability of lattice potentials and interaction strengths with high accuracy, ultracold atomic systems are expected to be a versatile tool for exploring many-body physics in strongly correlated systems. Of particular interest are multicomponent fermionic systems with alkaline earth-like atoms in cold-atom setups. Experimental realizations of such systems with SU(n) symmetric interactions have been reported in [5,6,7]. They are expected to be well described by the SU(n) Fermi–Hubbard model, which is a generalization of the standard Hubbard model with SU(2) symmetry [8,9,10], a minimal model for describing the properties of correlated electrons in solids. In conventional condensed matter physics, the SU(n) Hubbard model has been studied in the context of the large-n approach [11, 12]. In this approach, the main focus was on the large-n limit, and the physical properties of the models at finite n have been less investigated. However, the recent experimental realizations of the SU(n) Hubbard models with ultracold fermionic atoms have generated renewed theoretical interest in the study of the model at finite n ($> 2$). A number of studies revealed that the models can exhibit exotic phases that do not appear in the SU(2) counterpart [13,14,15,16,17,18,19,20].

Despite its apparent simplicity, mathematically rigorous treatment of the SU(n) Hubbard models is, in general, a formidable task due to the intricate competition between the kinetic and the on-site Coulomb terms. Also, since the internal degrees of freedom will be larger compared to the SU(2) case, exact results for the SU(n) Hubbard model are limited to fewer examples than in the SU(2) case. The Nagaoka ferromagnetism was the first rigorous result for the SU(2) Hubbard model [21,22,23]. Given infinitely large Coulomb repulsions and exactly one hole, it was proved that the ground state of the Hubbard model defined on a lattice satisfying a certain connectivity condition is ferromagnetic and unique. The Nagaoka ferromagnetism in the SU(n) Hubbard model has also been established in [24, 25]. As for the multiorbital Hubbard models, theorems about ferromagnetism have been proved in [26, 27], and its extension to the SU(n) case was also discussed [26]. The above results are for singular cases in the sense that the Coulomb interaction is infinitely large. The SU(2) Hubbard models with flat bands provide us with another rigorous example of ferromagnetism [28,29,30,31,32,33,34,35]. Here, by a flat band, we mean a structure of a single-particle energy spectrum with a macroscopic degeneracy. There are systematic methods for constructing tight-binding models with flat bands, such as the line graph [28] and the cell construction [29]. In the SU(2) case, one can prove that the ground states are ferromagnetic and unique apart from the trivial spin degeneracy under the following conditions: (i) the number of particles is equal to the multiplicity of the single-particle ground states, and (ii) the basis for the space spanned by the single-particle ground states is connected. Furthermore, its extension to the SU(n) case was recently been discussed in [36, 37]. The flat-band ferromagnetism can also be thought of as a result for a singular situation because the density of states at the Fermi level diverges. In the SU(2) case, the stability of flat-band ferromagnetism under perturbations which makes the lowest band dispersive has been discussed [38,39,40,41]. It was rigorously proved for a class of perturbed models in any dimension that the ground states remain ferromagnetic when the Coulomb repulsion and the band gap are sufficiently large [42, 43]. By contrast, the stability of the flat-band ferromagnetism in the SU(n) Hubbard model has been proved only in the one-dimensional case [37]. Thus, the generalizations in higher dimensions remain to be established. In addition, since the Hohenberg-Mermin-Wagner theorem [44, 45] forbids spontaneous symmetry breaking in one- and two-dimensional models with continuous symmetries at finite temperature, it is necessary to study models in dimensions higher than two in order to investigate ferromagnetism stable at finite temperature.

In this paper, we study a class of SU(n) Hubbard models on d ($\ge 2$)-dimensional decorated lattices and establish rigorous results. We first consider the models with flat bands at the bottom of the single-particle spectrum and prove that they exhibit SU(n) ferromagnetism in their ground states, provided that the on-site Coulomb interaction is repulsive and the total fermion number is the same as the number of unit cells. We then discuss SU(n) ferromagnetism in perturbed models obtained by adding extra hopping terms that make the flat bands dispersive. We prove that the particular perturbation leaves the ground states SU(n) ferromagnetic when the band width of the bottom band is sufficiently narrow and the Coulomb repulsion is sufficiently large. To establish the theorem for general n and dimensions d, it is necessary to treat two cases, $n \le d$ and $n > d$, separately. This is in marked contrast to the SU(2) case where the number of internal degrees of freedom cannot be greater than d $(\ge 2)$. In addition, it should also be mentioned that our proof of Theorem 2 considerably simplifies the previous proof for the nearly-flat-band ferromagnetism in the SU(2) case.

Fig. 1

The lattice structure of $\varLambda$ for $d=2$. Black and white dots represent external sites in ${\mathcal {E}}$ and internal sites in ${\mathcal {I}}$, respectively. We also depict the localized states generated by $a_{p, \alpha }^{\dag }$ and $b_{u, \alpha }^{\dag }$

The present paper is organized as follows. In Sect. 2, we shall describe our model on a d-dimensional decorated hypercubic lattice and state our results about SU(n) ferromagnetism. In Sect. 3, we prove the first result for the SU(n) flat-band ferromagnetism. In Sect. 4, we prove our main result for the SU(n) ferromagnetism in the model with a nearly flat band. In Appendix A, we show the linear independence of many-body states defined by localized states, and in Appendix B, we give explicit expressions for the local ground states of the effective Hamiltonian in each particle sector.

2 Model and Main Results

The models we consider are straightforward extensions of those in [32] to SU(n) case, and we follow the notation in [32]. Readers are also referred to [42] for the original result about ferromagnetism of the SU(2) Hubbard model with nearly flat bands.

2.1 Lattice

Let ${\mathcal {E}}$ be a set of sites in a d-dimensional hypercubic lattice of length L with unit lattice spacing and periodic boundary conditions, where we assume $d \ge 2$ and L is an odd integer. We take a new site in the middle of each bond of the lattice ${\mathcal {E}}$ and denote by ${\mathcal {I}}$ the collection of all such sites. In the following, we call $p \in {\mathcal {E}}$ an external site and $u \in {\mathcal {I}}$ an internal site. We define the SU(n) Hubbard model on the decorated hypercubic lattice $\varLambda = {\mathcal {E}} \cup {\mathcal {I}}$. The lattice structure for $d=2$ is shown in Fig. 1.

2.2 Fermion Operatores

We denote creation and annihilation operators by $c_{x, \alpha }^{\dag }$ and $c_{x, \alpha }$ for a fermion at site $x \in \varLambda$ with color $\alpha = 1, \ldots , n$. They satisfy the anticommutation relations

$$\begin{aligned} \{c_{x, \alpha }, c_{y, \beta }\} = \{c_{x, \alpha }^{\dag }, c_{y, \beta }^{\dag }\} = 0 \end{aligned}$$

(1)

and

$$\begin{aligned} \{c_{x, \alpha }, c_{y, \beta }^{\dag }\} = \delta _{\alpha , \beta } \delta _{x, y} \end{aligned}$$

(2)

for $x, y \in \varLambda$ and $\alpha , \beta = 1, \ldots , n$. The corresponding number operator is defined by $n_{x, \alpha } = c_{x, \alpha }^{\dag } c_{x, \alpha }$. The total fermion number is $N_{\mathrm {f}} = \sum _{x \in \varLambda } n_{x} =\sum _{\alpha =1}^{n} \sum _{x \in \varLambda } n_{x, \alpha }$, where $n_{x} = \sum _{\alpha =1}^{n} n_{x, \alpha }$. In the following, we consider $N_{\mathrm {f}}$-particle Hilbert space ${\mathcal {H}}_{N_{\mathrm {f}}}(\varLambda )$ with a fixed fermion number $N_{\mathrm {f}} = |{\mathcal {E}}| = L^{d}$, which means that the lowest band is 1/n filled.

We define color raising and lowering operators as

$$\begin{aligned} F^{\alpha , \beta } = \sum _{x \in \varLambda } c_{x, \alpha }^{\dag } c_{x, \beta } \ \ \text {for} \ \ \alpha \ne \beta \end{aligned}$$

(3)

and total number operators of fermion with color $\alpha$ as

$$\begin{aligned} F^{\alpha , \alpha } = \sum _{x \in \varLambda } c_{x, \alpha }^{\dag } c_{x, \alpha } \ \ \text {for} \ \ \alpha = 1, \ldots , n. \end{aligned}$$

(4)

We denote the eigenvalue of $F^{\alpha , \alpha }$ as $N_{\alpha }$.

To describe our model, we define a new set of operators

$$\begin{aligned} a_{p, \alpha }&= c_{p, \alpha } - \nu \sum _{\begin{array}{c} u \in {\mathcal {I}} \\ |p-u| = 1/2 \end{array}} c_{u, \alpha } \ \ \text {for} \ \ p \in {\mathcal {E}}, \end{aligned}$$

(5)

$$\begin{aligned} b_{u, \alpha }&= c_{u, \alpha } + \nu \sum _{\begin{array}{c} p \in {\mathcal {E}} \\ |p-u| = 1/2 \end{array}} c_{p, \alpha } \ \ \text {for} \ \ u \in {\mathcal {I}}, \end{aligned}$$

(6)

where $\nu > 0$. It is verified that these operators satisfy the following anticommutation relations for $p, q \in {\mathcal {E}}$ and $u, v \in {\mathcal {I}}$,

$$\begin{aligned} \{a_{p, \alpha }, a_{q, \beta }\}&= \{b_{u, \alpha }, b_{v, \beta }\} = \{a_{p, \alpha }, b_{u, \beta }\} = 0, \end{aligned}$$

(7)

$$\begin{aligned} \{a_{p, \alpha }, a_{q, \beta }^{\dag }\}&= {\left\{ \begin{array}{ll} \delta _{\alpha , \beta } (2d\nu ^{2} + 1) &{}\ \ \text {if} \ \ p=q, \\ \delta _{\alpha , \beta } \nu ^{2} &{}\ \ \text {if} \ \ |p-q|=1, \\ 0 &{}\ \ \text {otherwise}, \end{array}\right. } \end{aligned}$$

(8)

$$\begin{aligned} \{b_{u, \alpha }, b_{v, \beta }^{\dag }\}&= {\left\{ \begin{array}{ll} \delta _{\alpha , \beta }(2\nu ^{2} + 1) &{}\ \ \text {if} \ \ u=v, \\ \delta _{\alpha , \beta } \nu ^{2} &{}\ \ \text {if} \ \ u\ne v \ \text {and} \ C_{u} \cap C_{v} \ne \emptyset , \\ 0 &{}\ \ \text {otherwise}, \end{array}\right. } \end{aligned}$$

(9)

$$\begin{aligned} \{a_{p, \alpha }, b_{u, \beta }^{\dag }\}&= 0, \end{aligned}$$

(10)

where $C_{u}$ is a set of lattice sites consisting of u itself and $p \in {\mathcal {E}}$ such that $|p-u| = 1/2$.

2.3 Model

First, we study the following SU(n) Hubbard model

$$\begin{aligned} H_{1}&= H_{\mathrm {hop}} + H_{\mathrm {int}}, \end{aligned}$$

(11)

$$\begin{aligned} H_{\mathrm {hop}}&= t \sum _{\alpha =1}^{n} \sum _{u \in {\mathcal {I}}} b_{u, \alpha }^{\dag } b_{u, \alpha }, \end{aligned}$$

(12)

$$\begin{aligned} H_{\mathrm {int}}&= U \sum _{\alpha < \beta } \sum _{x \in \varLambda } n_{x, \alpha } n_{x, \beta }, \end{aligned}$$

(13)

where the parameters t, U are non-negative. When $U=0$, the model reduces to the tight-binding model described only by $H_{\mathrm {hop}}$. The hopping Hamiltonian $H_{\mathrm {hop}}$ can be written in the standard form

$$\begin{aligned} H_{\mathrm {hop}} =\sum _{\alpha =1}^{n} \sum _{x, y \in \varLambda } t_{x, y} c_{x, \alpha }^{\dag } c_{y, \alpha }, \end{aligned}$$

(14)

with hopping matrix elements given by

$$\begin{aligned} t_{x, y} = {\left\{ \begin{array}{ll} \nu t &{}\ \ \text {if} \ |x-y|=1/2, \\ \nu ^{2} t &{}\ \ \text {if} \ x, y \in {\mathcal {E}} \ \text {and} \ |x-y|=1, \\ t &{}\ \ \text {if} \ x=y \in {\mathcal {I}}, \\ 2d\nu ^{2} t &{}\ \ \text {if} \ x=y \in {\mathcal {E}}, \\ 0 &{}\ \ \text {otherwise}. \end{array}\right. } \end{aligned}$$

(15)

Let us consider the corresponding single-particle Schrödinger equation. Let $|\varPhi _{1}\rangle$ be a single-particle state of the form

$$\begin{aligned} |\varPhi _{1}\rangle = \sum _{x \in \varLambda } \phi _{x} c_{x, \alpha }^{\dag } |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(16)

where $\phi _{x} \in {\mathbb {C}}$ is a complex coefficient and $|\varPhi _{\mathrm {vac}}\rangle$ is a vacuum state of $c_{x, \alpha }$. One finds that the Schrödinger equation $H_{1} |\varPhi _{1}\rangle = \varepsilon |\varPhi _{1}\rangle$ leads to

$$\begin{aligned} \sum _{y \in \varLambda } t_{x,y} \phi _{y} = \varepsilon \phi _{x}, \end{aligned}$$

(17)

where $\varepsilon$ denotes a single-particle energy eigenvalue. By solving Eq. (17), one obtains the $d+1$ bands with dispersion relations

$$\begin{aligned} \varepsilon _{\mu }({\varvec{k}}) = {\left\{ \begin{array}{ll} 0\ \ &{}\mu = 1,\\ t \ \ &{}\mu = 2, \ldots , d, \\ t + 2t\nu ^{2} \sum _{j=1}^{d}(1 + \cos {k_{j}}) \ \ &{}\mu = d+1, \end{array}\right. } \end{aligned}$$

(18)

where ${\varvec{k}}$ is an element in ${\mathcal {K}}$ defined by

$$\begin{aligned} \mathcal {K} = \left\{ {\varvec{k}} = (k_{1}, \dots , k_{d}) \left| \ k_{j} = \frac{2\pi }{L} n_{j}, n_{j}= 0, \pm 1, \dots , \pm \frac{L-1}{2} \ \text {for} \ j=1, \dots , d\right. \right\} . \end{aligned}$$

(19)

One finds that the lowest band and the middle bands are dispersionless, which are referred to as flat bands. See Fig. 2 for the dispersion relations for $d=2$.

Fig. 2

The dispersion relations of the energy bands for $t=1$, $\nu =1/2$, and $d=2$. The lowest and middle bands are completely flat

In the present paper, we also study a perturbed model, whose Hamiltonian is given by

$$\begin{aligned} H_{2}&= H_{\mathrm {hop}}' + H_{\mathrm {int}}, \end{aligned}$$

(20)

$$\begin{aligned} H_{\mathrm {hop}}'&=-s \sum _{\alpha =1}^{n} \sum _{p \in {\mathcal {E}}} a_{p, \alpha }^{\dag } a_{p, \alpha } + t \sum _{\alpha =1}^{n} \sum _{u \in {\mathcal {I}}} b_{u, \alpha }^{\dag } b_{u, \alpha }, \end{aligned}$$

(21)

where the parameter s and t are non-negative and $H_{\mathrm {int}}$ is the same as Eq. (13). The Hamiltonian $H_{2}$ is precisely the same as $H_{1}$ at $s=0$. The hopping Hamiltonian $H_{\mathrm {hop}}'$ can also be written in the standard form

$$\begin{aligned} H_{\mathrm {hop}}' = \sum _{x, y \in \varLambda } t_{x,y}' c_{x, \alpha }^{\dag } c_{y, \alpha }, \end{aligned}$$

(22)

where the hopping matrix elements are given by

$$\begin{aligned} t_{x,y}' = {\left\{ \begin{array}{ll} \nu (t + s) &{} \ \ \text {if} \ |x-y| = 1/2, \\ \nu ^{2} t &{} \ \ \text {if} \ x, y \in {\mathcal {E}} \ \text {and} \ |x-y|=1, \\ -\nu ^{2} s &{} \ \ \text {if} \ x, y \in {\mathcal {I}} \ \text {and} \ x\ne y \ \text {and} \ ^{\exists }p \in {\mathcal {E}} \ \text {s.t.}\ |x-p|=|y-p|=1/2, \\ t - 2\nu ^{2} s &{} \ \ \text {if} \ x=y \in {\mathcal {I}}, \\ 2 d t \nu ^2 - s &{} \ \ \text {if} \ x=y \in {\mathcal {E}}, \\ 0 &{} \ \ \text {otherwise}. \end{array}\right. } \end{aligned}$$

(23)

The corresponding single-particle Schrödinger equation reads

$$\begin{aligned} \sum _{y \in \varLambda } t_{x,y}' \phi _{y} = \varepsilon \phi _{x}. \end{aligned}$$

(24)

By solving Eq. (24), one obtains the dispersion relations of $d+1$ bands

$$\begin{aligned} \varepsilon _{\mu }({\varvec{k}}) = {\left\{ \begin{array}{ll} -s -2s \nu ^{2} \sum _{j=1}^{d}(1 + \cos {k_{j}})\ \ &{}\mu =1, \\ t \ \ &{}\mu =2, \ldots , d, \\ t + 2t\nu ^{2} \sum _{j=1}^{d}(1 + \cos {k_{j}}) \ \ &{} \mu =d+1, \end{array}\right. } \end{aligned}$$

(25)

where ${\varvec{k}}$ is an element of ${\mathcal {K}}$ defined in Eq. (19). Although the middle bands are still dispersionless, the lowest band has become dispersive because of the additional hopping term proportional to the parameter s. See Fig. 3. We note that all the single-particle properties of the hopping Hamiltonian are exactly the same as in the SU(2) case [31].

Fig. 3

The dispersion relations of the energy bands for $t=1$, $s=1/5$, $\nu = 1/2$, and $d=2$. The lowest band is dispersive, while the middle band is still dispersionless

Since both of the Hamiltonians $H_{1}$ and $H_{2}$ have SU(n) symmetry, the operators $F^{\alpha , \beta }$ commute with $H_{1}$ and $H_{2}$. Therefore, the Hilbert space can be separated into different sectors labeled by $(N_{1}, \ldots , N_{n})$, which is denoted by ${\mathcal {H}}_{N_{1}, \ldots , N_{n}}(\varLambda )$. To establish our theorem, we define fully polarized states. A fully polarized state with color $\alpha$ is defined as

$$\begin{aligned} |\varPhi _{\mathrm {all}, \alpha }\rangle = \prod _{p \in {\mathcal {E}}} a_{p, \alpha }^{\dag } |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(26)

where $|\varPhi _{\mathrm {vac}}\rangle$ is a vacuum state of $c_{x, \alpha }$. From the anticommutation relations (8) and (9), we can see that the fully polarized state is an eigenstate of both $H_{1}$ and $H_{2}$. Due to the SU(n) symmetry, one obtains a general form of the eigenstates with the same energy as $|\varPhi _{\mathrm {all}, \alpha }\rangle$:

$$\begin{aligned} |\varPhi _{N_{1}, \ldots , N_{n}}\rangle = \left( F^{n, 1}\right) ^{N_{n}} \cdots \left( F^{2, 1}\right) ^{N_{2}} |\varPhi _{\mathrm {all}, 1}\rangle , \end{aligned}$$

(27)

where $N_{1} = |{\mathcal {E}}| - \sum _{\alpha =2}^{n} N_{\alpha }$. We also refer to states of the form (27) as fully polarized states. The fully polarized states can be characterized as eigenstates of the quadratic Casimir operator $C_{2}$ of the SU(n) group defined as [46]

$$\begin{aligned} C_{2} = \frac{1}{2} \left( \sum _{\alpha , \beta = 1}^{n} F^{\alpha , \beta } F^{\beta , \alpha } - \frac{N_\mathrm{f}^{2}}{n} \right) . \end{aligned}$$

(28)

The fully polarized states (26) and (27) are eigenstates of $C_{2}$ with eigenvalue $\frac{|{\mathcal {E}}| (n-1)}{2} \left( \frac{|{\mathcal {E}}|}{n} + 1\right)$, which is the maximum eigenvalue of $C_{2}$ for fixed $N_\mathrm{f}$.

2.4 Results

First, we describe the theorem for the model whose lowest band is flat. This is a slight generalization of the result obtained by Liu et al., in [36], in the sense that our hopping Hamiltonian has one more parameter.

Theorem 1

Consider the SU(n) Hubbard Hamiltonian (11) with the total fermion number $N_{\mathrm {f}} = |{\mathcal {E}}|$. For arbitrary $t > 0$ and $U > 0$, the ground states of the Hamiltonian are the fully polarized states and unique apart from trivial degeneracy due to the SU(n) symmetry.

The theorem is proved in Sect.3. The fully polarized states are the SU(n) counterparts of the ferromagnetic states in the SU(2) Hubbard model. Thus, Theorem 1 establishes the SU(n) ferromagnetism in the d-dimensional Hubbard models with flat bands.

Since the density of states at the Fermi energy diverges, Theorem 1 can be thought of as a result for a singular case. In order to establish a rigorous result for a nonsingular SU(n) Hubbard model, we consider the Hamiltonian (20). To state the theorem, let us define for $n\ge 2$ and $d \ge 2$,

$$\begin{aligned} \nu _{\mathrm {c}}(n, d) ={\left\{ \begin{array}{ll} \sqrt{\frac{2d+1 + \sqrt{4(2d-n) (n-1) + (2d+1)^{2}}}{2(2d-n) (n-1)}} &{}\ \ \text {for} \ \ n \le d, \\ \sqrt{\frac{2d+1 + \sqrt{8d^{2} + 1}}{2d(d-1)}} &{}\ \ \text {for} \ \ n > d. \end{array}\right. } \end{aligned}$$

(29)

Theorem 2

Consider the SU(n) Hubbard Hamiltonian (20) with the total fermion number $N_{\mathrm {f}} = |{\mathcal {E}}|$. Assume that $0 < \nu \le \nu _{\mathrm {c}}(n, d)$ for $d \ge 2$. For sufficiently large $t/s>0$ and $U/s>0$, the ground states are the fully polarized states and unique apart from the trivial degeneracy due to the SU(n) symmetry.

Since the lowest band of Eq. (21) is dispersive, the density of states at the Fermi level does not diverge. Also, the Coulomb interaction is assumed to be sufficiently large but finite. Thus, the theorem establishes SU(n) ferromagnetism in the ground states of non-singular SU(n) Hubbard models in arbitrary dimensions. The theorem is proved in Sect. 4.

3 Proof of Theorem 1

We shall prove Theorem 1 in this section. We consider the model (11) with the fermion number $N_{\mathrm {f}} = |{\mathcal {E}}|$. First, we note that the hopping Hamiltonian $H_{\mathrm {hop}}$ and $H_{\mathrm {int}}$ are positive semidefinite because $b_{u, \alpha }^{\dag } b_{u, \alpha }$ and $n_{x, \alpha } n_{x, \beta } = (c_{x, \alpha } c_{x, \beta })^{\dag } c_{x, \alpha } c_{x, \beta }$ are positive semidefinite. Hence, the total Hamiltonian $H_{1}$ is positive semidefinite as well. By noting $\{a_{p, \alpha }, b_{u, \beta }^{\dag }\} = 0$, we find that the fully polarized state $|\varPhi _{\mathrm {all}, \alpha }\rangle$ is an eigenstate of $H_{1}$ with eigenvalue zero and all the fully polarized states are zero energy states due to the SU(n) symmetry. Since $H_{1} \ge 0$, the fully polarized states are ground states of $H_{1}$.

In the following, we prove the uniqueness of ground states. Let $|\varPhi _{\mathrm {GS}}\rangle$ be an arbitrary ground state of $H_{1}$ with $N_{\mathrm {f}} = |{\mathcal {E}}|$, which means $H_{1} |\varPhi _{\mathrm {GS}}\rangle = 0$. We note that in general, the ground state can be written as

$$\begin{aligned} |\varPhi _{\mathrm {GS}}\rangle&= \sum _{\begin{array}{c} A_{1}, \ldots , A_{n} \subset {\mathcal {E}}\\ B_{1}, \ldots , B_{n} \subset {\mathcal {I}} \\ \end{array}} f(\{ A_{\alpha }\}, \{B_{\alpha }\}) \left( \prod _{p \in A_{1}} a_{p, 1}^{\dag }\right) \nonumber \\&\cdots \left( \prod _{p \in A_{n}} a_{p, n}^{\dag }\right) \left( \prod _{u \in B_{1}} b_{u, 1}^{\dag }\right) \cdots \left( \prod _{u \in B_{n}} b_{u, n}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(30)

where $A_{\alpha }$ and $B_{\alpha }$ are subsets of ${\mathcal {E}}$ and ${\mathcal {I}}$, respectively such that $\sum _{\alpha =1}^{n} \left( |A_{\alpha }| + |B_{\alpha }|\right) = |{\mathcal {E}}|$. Here, $f(\{A_{\alpha }\}, \{B_{\alpha }\})$ is a certain coefficient. The inequalities $H_{\mathrm {hop}} \ge 0$ and $H_{\mathrm {int}} \ge 0$ imply that $H_{\mathrm {hop}} |\varPhi _{\mathrm {GS}}\rangle = 0$ and $H_{\mathrm {int}} |\varPhi _{\mathrm {GS}}\rangle = 0$, which yield

$$\begin{aligned} b_{u, \alpha } |\varPhi _{\mathrm {GS}}\rangle&= 0 \ \text {for any} \ u \in {\mathcal {I}} \ \text {and} \ \ \alpha = 1, \ldots , n, \end{aligned}$$

(31)

$$\begin{aligned} c_{x, \alpha } c_{x, \beta } |\varPhi _{\mathrm {GS}}\rangle&= 0 \ \text {for any} \ x \in \varLambda \ \text {and} \ \ \alpha \ne \beta . \end{aligned}$$

(32)

From the anticommutation relations (10), it is clear that the conditions (31) imply that $|\varPhi _{\mathrm {GS}}\rangle$ does not contain any $b_{u, \alpha }^{\dag }$ operator. (See Appendix A for a proof that the states obtained by acting with $a^{\dag }$ and $b^{\dag }$ operators on $|\varPhi _{\mathrm {vac}}\rangle$ are linearly independent.¹ Therefore, the ground state is written as

$$\begin{aligned} |\varPhi _{\mathrm {GS}}\rangle = \sum _{\begin{array}{c} A_{1} , \ldots , A_{n} \subset {\mathcal {E}} \\ \sum _{\alpha =1}^{n} |A_{\alpha }| = |{\mathcal {E}}| \end{array} } g\left( \{A_{\alpha }\}\right) \left( \prod _{p \in A_{1}} a_{p, 1}^{\dag }\right) \cdots \left( \prod _{p \in A_{n}} a_{p, n}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(33)

where $g\left( \{A_{\alpha }\}\right)$ is a certain coefficient.

We next examine the conditions (32). We first consider the case where $x=p \in {\mathcal {E}}$. It follows from

$$\begin{aligned} \{c_{p, \alpha }, a_{q, \beta }^{\dag }\} = \delta _{\alpha , \beta } \delta _{p, q} \ \ \text {for} \ p, q \in {\mathcal {E}}. \end{aligned}$$

(34)

and Eq. (32) that $g(\{A_{\alpha }\}) = 0$ if there exist $A_{\alpha }$ and $A_{\beta }$ such that $A_{\alpha } \cap A_{\beta } \ne \emptyset$. Noting that $\sum _{\alpha =1}^{n} |A_{\alpha }| = |{\mathcal {E}}|$, we find that $\cup _{\alpha =1}^{n} A_{\alpha } = {\mathcal {E}}$ when $A_{\alpha } \cap A_{\beta } = \emptyset$ for all $\alpha \ne \beta$. Thus one can rewrite the ground state in the form of

$$\begin{aligned} |\varPhi _{\mathrm {GS}}\rangle = \sum _{\varvec{\alpha }} C(\varvec{\alpha }) \left( \prod _{p \in {\mathcal {E}}} a_{p, \alpha _{p}}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(35)

where the sum is over all possible color configuration $\varvec{\alpha } = (\alpha _{p})_{p \in {\mathcal {E}}}$ with $\alpha _{p} = 1, \ldots , n$. In the derivation of Eq. (35), we again exploit the linear independence of the states consisting of $a^{\dag }$ and $b^{\dag }$ operators, which is proved in Appendix A. Then we consider the conditions (32) for $x=u \in {\mathcal {I}}$. For $u \in {\mathcal {I}}$, we can check that

$$\begin{aligned} \{c_{u, \alpha }, a_{p, \beta }^{\dag }\} ={\left\{ \begin{array}{ll} -\nu \delta _{\alpha , \beta } \ \ &{}\text {if} \ |u-p| = 1/2, \\ 0 \ \ &{}\text {otherwise}. \end{array}\right. } \end{aligned}$$

(36)

Using Eq. (36), we obtain

$$\begin{aligned} c_{u, \alpha } c_{u, \beta } |\varPhi _{\mathrm {GS}}\rangle = \sum _{\begin{array}{c} \varvec{\alpha } \\ \mathrm {s.t.} \alpha _{p} = \beta , \alpha _{q} = \alpha \end{array}} (C(\varvec{\alpha }) - C(\varvec{\alpha }_{p \leftrightarrow q})) \left( \prod _{p' \in {\mathcal {E}} \backslash \{p, q\}}a_{p', \alpha _{p'}}^{\dag } \right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(37)

where p and q are external sites which satisfy $|u-p| = |u-q| = 1/2$. The color configuration $\varvec{\alpha }_{p \leftrightarrow q}$ is obtained from $\varvec{\alpha }$ by swapping $\alpha _{p}$ and $\alpha _{q}$. Since all the states in the sum are linearly independent (see Appendix A), we find from Eq. (32) that $C(\varvec{\alpha }) = C(\varvec{\alpha }_{p \leftrightarrow q})$ for all $\varvec{\alpha }$ and all $u \in {\mathcal {I}}$. Recalling that there is one localized state centered at each external site and neighboring localized states always share one internal site, we see that $C(\varvec{\alpha }) = C(\varvec{\alpha }_{p \leftrightarrow q})$ for p, $q \in {\mathcal {E}}$. Since an arbitrary permutation of the color configuration ${\varvec{\alpha }}$ can be generated by repeatedly swapping the colors on neighboring external sites, we have

$$\begin{aligned} C({\varvec{\alpha }}) = C({\varvec{\beta }}) \end{aligned}$$

(38)

if ${\varvec{\beta }}$ is a permutation of ${\varvec{\alpha }}$.

Lastly, we show that the states satisfying Eq. (38) are the fully polarized states defined in Eqs. (26), (27). To this aim, we introduce a concept of a word [47]. A word $w = (w_{1}, \ldots , w_{|{\mathcal {E}}|})$ is a sequence of integers where $w_{i} \in \{1, \ldots , n\}$ for all i. We denote by $|w|_{\alpha }$ the number of occurrences of $\alpha$ in w. We define the set of words for which $|w|_{\alpha } = N_{\alpha }$ holds as follows: $W(N_{1}, \ldots , N_{n}) = \{w| \ |w|_{\alpha } = N_{\alpha }, \alpha =1, \ldots , n\}$. For example, W(2, 0, 1) consists of (1, 1, 3), (1, 3, 1) and (3, 1, 1). It follows from Eq. (38) that the ground state of $H_{1}$ in the sector ${\mathcal {H}}_{N_{1}, \ldots , N_{n}}(\varLambda )$ can be written as

$$\begin{aligned} |{\tilde{\varPhi }}_{N_{1}, \ldots , N_{n}}\rangle = \sum _{w \in W(N_{1}, \ldots , N_{n})} a_{p_{1}, w_{1}}^{\dag } a_{p_{2}, w_{2}}^{\dag } \cdots a_{p_{|{\mathcal {E}}|}, w_{|{\mathcal {E}}|}}^{\dag } |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(39)

where each external site is labeled by an integer i as $p_{i}$. On the other hand, using commutation relations $[F^{\beta , \alpha }, a_{p, \gamma }^{\dag }] = \delta _{\alpha , \gamma } a_{p, \beta }^{\dag }$, we see that

$$\begin{aligned} \left( F^{2, 1}\right) ^{N_{2}} |\varPhi _{\mathrm {all}, 1}\rangle = N_{2}! \sum _{w \in W(|{\mathcal {E}}|-N_{2}, N_{2})} a_{p_{1}, w_{1}}^{\dag } a_{p_{2}, w_{2}}^{\dag } \cdots a_{p_{|{\mathcal {E}}|}, w_{|{\mathcal {E}}|}}^{\dag } |\varPhi _{\mathrm {vac}}\rangle . \end{aligned}$$

(40)

By repeating the same procedure, we have the desired result, i.e., $|\varPhi _{N_{1}, \ldots , N_{n}}\rangle$ and $|{\tilde{\varPhi }}_{N_{1}, \ldots , N_{n}}\rangle$ are the same up to a normalization. Thus, Theorem 1 is proved. $\square$

4 Proof of Theorem 2

In this section, we shall prove Theorem 2. Here, we consider the model (20) with the fermion number $N_{\mathrm {f}} = |{\mathcal {E}}|$. First, we rewrite the Hamiltonian $H_{2}$ as follows:

$$\begin{aligned} H_{2} = -s |{\mathcal {E}}| (2d\nu ^{2} + 1) + \lambda H_{\mathrm {flat}} + \sum _{p \in {\mathcal {E}}} h_{p}, \end{aligned}$$

(41)

where

$$\begin{aligned} H_{\mathrm {flat}} = \sum _{\alpha =1}^{n} \sum _{u \in {\mathcal {I}}} b_{\alpha , u}^{\dag } b_{\alpha , u} + \sum _{\alpha < \beta } \sum _{x\in \varLambda } n_{x, \alpha } n_{x, \beta }, \end{aligned}$$

(42)

and the local Hamiltonian for each $p \in {\mathcal {E}}$ is defined as

$$\begin{aligned} h_{p}&= s(2d\nu ^{2} + 1) - s \sum _{\alpha =1}^{n} a_{p, \alpha }^{\dag } a_{p, \alpha } + \frac{t-\lambda }{2} \sum _{\alpha =1}^{n} \sum _{\begin{array}{c} u \in {\mathcal {I}} \\ |p-u|=1/2 \end{array}} b_{u, \alpha }^{\dag } b_{u, \alpha } \nonumber \\&\quad + \sum _{\alpha < \beta } \left[ (U-\lambda )(1-\kappa ) n_{p, \alpha } n_{p, \beta } + \frac{\kappa (U-\lambda )}{2d} \right. \nonumber \\&\quad \left. \times \,\sum _{\begin{array}{c} q \in {\mathcal {E}} \\ |p-q|=1 \end{array}} n_{q, \alpha } n_{q, \beta } + \frac{U -\lambda }{2} \sum _{\begin{array}{c} u \in \lambda \\ |p-u|=1/2 \end{array}} n_{u, \alpha } n_{u, \beta } \right] . \end{aligned}$$

(43)

The two parameters $\lambda$ and $\kappa$ satisfy $0< \lambda < \min \{t, U\}$ and $0< \kappa < 1$. To prove Theorem 2, we prove the following lemma.

Lemma 1

Suppose the local Hamiltonian $h_{p}$ is positive semidefinite for any $p \in {\mathcal {E}}$. Then the ground states of the Hamiltonian $H_{2}$ are fully polarized states and unique apart from the trivial degeneracy due to the SU(n) symmetry.

Proof

First, we can check that a fully polarized state $|\varPhi _{\mathrm {all}, 1}\rangle$ satisfies $h_{p} |\varPhi _{\mathrm {all}, 1}\rangle = 0$ for all $p \in {\mathcal {E}}$. Since the local Hamiltonian is SU(n) symmetric, all the fully polarized states have zero energy for $h_{p}$. From the assumption that $h_{p}$ is positive semidefinite, it follows that the fully polarized states are the ground states of $h_{p}$ for all $p\in {\mathcal {E}}$. These states are also the ground states of $H_{\mathrm {flat}}$ since $H_{\mathrm {flat}}$ is equal to $H_{1}$ with $t=U=1$. Hence, the fully polarized states are the ground states of of $H_{2}$ with energy $-s|{\mathcal {E}}|(2d\nu ^{2} + 1)$. On the other hand, any ground state of $H_{2}$ must be a simultaneous ground state of $h_{p}$ and $H_{\mathrm {flat}}$ since both $h_{p}$ and $H_{\mathrm {flat}}$ are positive semidefinite. It is shown from Theorem 1 that the ground states of $H_{\mathrm {flat}}$ are fully polarized states and unique apart from trivial degeneracy due to the SU(n) symmetry. Thus, Lemma 1 is proved. $\square$

We can establish the positive semidefiniteness of $h_{p}$ by the following lemma.

Lemma 2

Suppose that t/s, U/s are sufficiently large and $0< \nu < \nu _{c}(n, d)$. Then the local Hamiltonian $h_{p}$ is positive semidefinite. (We take $\lambda$ and $\kappa$ to be proportional to s.)

Proof

Because of the translational invariance, it suffices to prove the lemma for $h_{o}$ with a fixed $o \in {\mathcal {E}}$. For convenience, we set $o = (0, 0, \ldots , 0) \in {\mathbb {Z}}^{d}$. We regard $h_{o}$ as an operator acting on a lattice $\varLambda _{o}$ with $4d+1$ sites. This lattice is written as $\varLambda _{o} = \{o\} \cup {\mathcal {E}}_{o} \cup {\mathcal {I}}_{o}$, where ${\mathcal {E}}_{o}$ is a set of 2d external sites of the form $(0, \ldots , 0, \pm 1 ,0, \ldots , 0)$, and ${\mathcal {I}}_{o}$ is a set of 2d internal sites of the form $(0, \ldots , 0, \pm 1/2, 0, \ldots , 0)$. See Fig. 4.

Fig. 4

The lattice $\varLambda _{o}$ with $4d+1$ sites for $d=2$. The origin o is located at the center. The sets ${\mathcal {E}}_{o}$ and ${\mathcal {I}}_{o}$ consist of black dots and white dots, respectively. We draw a the states generated by ${\tilde{a}}_{o, \alpha }^{\dag }$ and ${\tilde{a}}_{p, \alpha }^{\dag }$ and b the state generated by ${\tilde{b}}_{u, \alpha }^{\dag }$

In the following, we denote by ${\mathcal {F}}(\varLambda _{o})$ the Fock space on $\varLambda _{o}$. We define the total fermion number and the total number operator of fermion with color $\alpha$ on ${\mathcal {F}}(\varLambda _{o})$ as

$$\begin{aligned} F = \sum _{x \in \varLambda _{o}} n_{x}, \end{aligned}$$

(44)

and

$$\begin{aligned} M_{\alpha } = \sum _{x \in \varLambda _{o}} n_{x, \alpha }, \end{aligned}$$

(45)

respectively.

It is convenient to define operators on ${\mathcal {F}}(\varLambda _{o})$,

$$\begin{aligned} {\tilde{a}}_{o, \alpha }&= c_{o, \alpha } - \nu \sum _{u \in {\mathcal {I}}_{o}} c_{u, \alpha }, \end{aligned}$$

(46)

$$\begin{aligned} {\tilde{a}}_{p, \alpha }&= \frac{1}{\sqrt{\nu ^{2} + 1}} \left( c_{p, \alpha } - \nu c_{p/2, \alpha } \right) \ \ \text {for} \ p \in {\mathcal {E}}_{o}, \end{aligned}$$

(47)

$$\begin{aligned} {\tilde{b}}_{u, \alpha }&= c_{u, \alpha } + \nu c_{o, \alpha } + \nu c_{2u, \alpha } \ \ \text {for} \ u \in {\mathcal {I}}_{o}. \end{aligned}$$

(48)

These operators satisfy the anticommutation relations

$$\begin{aligned} \{{\tilde{a}}_{p, \alpha }, {\tilde{a}}_{q, \beta }^{\dag }\}&= {\left\{ \begin{array}{ll} \delta _{\alpha , \beta } (2d\nu ^{2} + 1 ) \ \ &{}\text {for} \ p = q = o, \\ \delta _{\alpha , \beta } \frac{\nu ^{2}}{\sqrt{\nu ^{2} + 1}} \ \ &{}\text {for} \ p=o, q \in {\mathcal {E}}_{o}, \\ \delta _{\alpha , \beta } \delta _{p, q} \ \ &{}\text {for} \ p, q \in {\mathcal {E}}_{o}, \end{array}\right. } \end{aligned}$$

(49)

$$\begin{aligned} \{{\tilde{b}}_{u, \alpha }, {\tilde{b}}_{v, \beta }^{\dag }\}&= {\left\{ \begin{array}{ll} \delta _{\alpha , \beta }(2\nu ^{2} + 1) \ \ &{}\text {for} \ u=v, \\ \delta _{\alpha , \beta } \nu ^{2} \ \ &{}\text {for} \ u \ne v, \end{array}\right. } \end{aligned}$$

(50)

$$\begin{aligned} \{{\tilde{a}}_{p, \alpha }, {\tilde{b}}_{u, \beta }^{\dag }\}&= 0 \ \ \text {for} \ p \in \{o\} \cup {\mathcal {E}}_{o}, u \in {\mathcal {I}}_{o}. \end{aligned}$$

(51)

With these operators, $h_{o}$ can be written as

$$\begin{aligned} h_{o}&= s(2d\nu ^{2} + 1) - s \sum _{\alpha =1}^{n} {\tilde{a}}_{o, \alpha }^{\dag } {\tilde{a}}_{o, \alpha } + \frac{t-\lambda }{2} \sum _{\alpha =1}^{n} \sum _{u\in {\mathcal {I}}_{o}} {\tilde{b}}_{u, \alpha }^{\dag } {\tilde{b}}_{u, \alpha } \nonumber \\&\quad + \sum _{\alpha <\beta } \left[ (U-\lambda )(1-\kappa ) n_{o, \alpha } n_{o, \beta } + \frac{\kappa (U-\lambda )}{2d} \sum _{p \in {\mathcal {E}}_{o}} n_{p, \alpha } n_{q, \alpha } + \frac{U-\lambda }{2} \sum _{u \in {\mathcal {I}}_{o}} n_{u, \alpha } n_{u, \beta } \right] . \end{aligned}$$

(52)

Since the total fermion number F on $\varLambda _{o}$ is conserved, we can examine $h_{o}$ in each fermion-number sector separately. The Fock space can be decomposed into a direct sum as

$$\begin{aligned} {\mathcal {F}}(\varLambda _{o}) = \bigoplus _{F=0}^{n(4d+1)} {\mathcal {H}}_{F}(\varLambda _{o}), \end{aligned}$$

(53)

where ${\mathcal {H}}_{F}(\varLambda _{o})$ is the F-fermion Hilbert space. Note that the lattice $\varLambda _{o}$ can accommodate at most $n(4d+1)$ fermions. Since $M_{\alpha }$ ($\alpha =1, \ldots , n$) commutes with $h_{o}$, we can further decompose each fermion-number sector ${\mathcal {H}}_{F}(\varLambda _{o})$ into a direct sum of smaller subspaces as

$$\begin{aligned} {\mathcal {H}}_{F}(\varLambda _{o}) = \bigoplus _{\begin{array}{c} (M_{1}, \ldots , M_{n}) \\ \sum _{\alpha =1}^{n} M_{\alpha } = F \end{array}} {\mathcal {J}}_{M_{1}, \ldots , M_{n}}(\varLambda _{o}), \end{aligned}$$

(54)

where ${\mathcal {J}}_{M_{1}, \ldots , M_{n}}(\varLambda _{o})$ is the Hilbert space with each $M_{\alpha }$ fixed. First, we consider the one-particle sector ${\mathcal {H}}_{1}(\varLambda _{o})$. The noninteracting part of $h_{o}$ is

$$\begin{aligned} s(2d\nu ^{2} + 1) - s \sum _{\alpha =1}^{n} {\tilde{a}}_{o, \alpha }^{\dag } {\tilde{a}}_{o, \alpha } + \frac{t-\lambda }{2} \sum _{\alpha =1}^{n} \sum _{u \in {\mathcal {I}}_{o}} {\tilde{b}}_{u, \alpha }^{\dag } {\tilde{b}}_{u, \alpha }, \end{aligned}$$

(55)

and one finds that the single-particle eigenenergies are 0, $s(2d\nu ^{2} + 1)$, $s(2d\nu ^{2} + 1) + (t-\lambda )((2d+1)\nu ^{2} + 1)/2$, $s(2d\nu ^{2} + 1) + (t-\lambda )(\nu ^{2} + 1)/2$. These eigenenergies are non-negative. Therefore, we see that the local ground state energy $E_{1}^{\mathrm {GS}}$ in ${\mathcal {H}}_{1}(\varLambda _{o})$ is given by

$$\begin{aligned} E_{1}^{\mathrm {GS}} = 0, \end{aligned}$$

(56)

and $h_{o}$ within this subspace is positive semidefinite.

Next, we consider ${\mathcal {H}}_{F}(\varLambda _{o})$ where $2 \le F \le n(4d+1)$. We shall prove

$$\begin{aligned} \lim _{t, U \rightarrow \infty } \langle {\varPhi }|{h_{o}}|{\varPhi }\rangle \ge 0 \end{aligned}$$

(57)

for any normalized state $|\varPhi \rangle \in {\mathcal {H}}_{F}(\varLambda _{o})$ for all F. To prove Eq. (57), we only need to examine states such that $\lim _{t, U \rightarrow \infty } \langle {\varPhi }|{h_{o}}|{\varPhi }\rangle < \infty$, which we call finite-energy states. Since ${\tilde{b}}_{u, \alpha }^{\dag }{\tilde{b}}_{u, \alpha }$ and $n_{x, \alpha } n_{x, \beta }$ are positive semidefinite, the condition that $|\varPhi \rangle$ is a finite-energy state is equivalent to the following:

$$\begin{aligned} {\tilde{b}}_{u, \alpha } |\varPhi \rangle&= 0 \ \ \text {for all} \ u\in {\mathcal {I}}_{o} \end{aligned}$$

(58)

$$\begin{aligned} c_{x, \alpha } c_{x, \beta } |\varPhi \rangle&= 0 \ \ \text {for all} \ x \in \varLambda _{o} \ \text {and} \ \alpha \ne \beta . \end{aligned}$$

(59)

Repeating the same argument in the proof of Theorem 1, we find that the conditions ${\tilde{b}}_{u, \alpha } |\varPhi \rangle = 0$ for any $u \in {\mathcal {I}}_{o}$ and $c_{p, \alpha } c_{p, \beta } |\varPhi \rangle = 0$ for $p \in {\mathcal {E}}_{o}$ imply that any finite-energy state $|\varPhi \rangle$ is in the form of

$$\begin{aligned} |\varPhi \rangle = \sum _{\begin{array}{c} A_{1} ,\ldots , A_{n} \subset \{o\} \cup {\mathcal {E}}_{o} \\ A_{\alpha } \cap A_{\beta } = \emptyset \ \text {for all} \ \alpha \ne \beta \\ \sum _{\alpha =1}^{n} |A_{\alpha }| = F \end{array}} f' \left( \{A_{\alpha }\}\right) \left( \prod _{p \in A_{1}} {\tilde{a}}_{p, 1}^{\dag }\right) \ldots \left( \prod _{p \in A_{n}} {\tilde{a}}_{p, n}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(60)

where $A_{\alpha }$ is a subset of $\{o\} \cup {\mathcal {E}}_{o}$ such that $A_{\alpha } \cap A_{\beta } = \emptyset$ for $\alpha \ne \beta$ and $\sum _{\alpha =1}^{n}|A_{\alpha }| = F$, and $f'(\{A_{\alpha }\})$ is a certain coefficient. We divide the collections of subsets $\{A_{\alpha }\}$ into two groups: those that contain the origin o and those that do not, and then we express the finite-energy state as

$$\begin{aligned} |\varPhi \rangle = |\varPhi _{o}\rangle + |\tilde{\varPhi }\rangle , \end{aligned}$$

(61)

where

$$\begin{aligned} |\varPhi _{o}\rangle&= \sum _{\begin{array}{c} A_{1}, \ldots , A_{n} \in \{o\} \cup {\mathcal {E}}_{o} \\ A_{\alpha } \cap A_{\beta } = \emptyset \ \text {for all} \ \alpha \ne \beta \\ o \in \cup _{\alpha =1}^{n} A_{\alpha }, \sum _{\alpha =1}^{n}|A_{\alpha }| = F \end{array}} f' \left( \{A_{\alpha }\}\right) \left( \prod _{p \in A_{1}} {\tilde{a}}_{p, 1}^{\dag }\right) \ldots \left( \prod _{p \in A_{n}} {\tilde{a}}_{p, n}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(62)

$$\begin{aligned} |\tilde{\varPhi }\rangle&=\sum _{\begin{array}{c} A_{1}, \ldots , A_{n} \in {\mathcal {E}}_{o} \\ A_{\alpha } \cap A_{\beta } = \emptyset \ \text {for all} \ \alpha \ne \beta \\ \sum _{\alpha =1}^{n}|A_{\alpha }| = F \end{array}} f' \left( \{A_{\alpha }\}\right) \left( \prod _{p \in A_{1}}{\tilde{a}}_{p, 1}^{\dag }\right) \ldots \left( \prod _{p \in A_{n}}{\tilde{a}}_{p, n}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle . \end{aligned}$$

(63)

Using the conditions $c_{u, \alpha } c_{u, \beta } |\varPhi \rangle = 0$ for $u \in {\mathcal {I}}_{o}$, we see that $|\varPhi _{o}\rangle$ should be a fully polarized state and $h_{o} |\varPhi _{o}\rangle = 0$. Therefore, the expectation value $\langle {\varPhi }|{h_{o}}|{\varPhi }\rangle$ is the same as $\langle {{\tilde{\varPhi }}}|{h_{o}}|{{\tilde{\varPhi }}}\rangle$ and it is enough to examine normalized states in the form of Eq. (63). We can verify that for any states written as Eq. (63) it holds that

$$\begin{aligned} {\tilde{a}}_{o, \alpha }^{\dag } |\tilde{\varPhi }\rangle = \frac{\nu ^{2}}{\sqrt{\nu ^{2} + 1}} \sum _{p \in {\mathcal {E}}_{o}} {\tilde{a}}_{p, \alpha }^{\dag } |\tilde{\varPhi }\rangle , \end{aligned}$$

(64)

and hence

$$\begin{aligned} \langle {{\tilde{\varPhi }}}|{h_{o}}|{{\tilde{\varPhi }}}\rangle = s(2d\nu ^{2} + 1) - \frac{s\nu ^{4}}{\nu ^{2} + 1} \sum _{\alpha =1}^{n} \sum _{p, q \in {\mathcal {E}}_{o}} \langle {\tilde{\varPhi }}| {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{q, \alpha } |\tilde{\varPhi }\rangle . \end{aligned}$$

(65)

Since the state $|\tilde{\varPhi }\rangle$ does not have any double occupancy of ${\tilde{a}}_{p, \alpha }$ states, to evaluate Eq. (65), it is sufficient to consider the effective Hamiltonian

$$\begin{aligned} h_{\mathrm {eff}} = E_{0} P - s' \sum _{\alpha =1}^{n} \sum _{p, q \in {\mathcal {E}}_{o}} P {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{q, \alpha } P, \end{aligned}$$

(66)

where $E_{0} = s(2d\nu ^{2} + 1)$ and $s' = s\nu ^{4}/(\nu ^{2} + 1)$. The operator $P = \prod _{\alpha < \beta } \prod _{x \in {\mathcal {E}}_{o}} \Big (1 - {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{p, \alpha } {\tilde{a}}_{p, \beta }^{\dag } {\tilde{a}}_{p, \beta }\Big )$ denotes the projection operator onto the space with no double occupancy of ${\tilde{a}}_{p, \alpha }$ states, which is known as the Gutzwiller projection [48]. Here, we note that the operators ${\tilde{a}}_{p, \alpha }$ obey the usual anticommutation relations, $\{{\tilde{a}}_{p, \alpha }^{\dag }, {\tilde{a}}_{q, \beta }\} = \delta _{\alpha , \beta } \delta _{p, q}$. Noting that the operator P excludes the doubly occupied states, we do not have to consider the F-fermion Hilbert space where $F > 2d$ because doubly occupied states always appear in such a sector. Therefore, in the following, we restrict ourselves to the subspaces with $2 \le F \le 2d$ and evaluate $h_{\mathrm {eff}}$. Using the following relations

$$\begin{aligned} - P {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{q, \alpha } P = {\tilde{a}}_{q, \alpha } P {\tilde{a}}_{p, \alpha }^{\dag } \ \ \text {for} \ p ,q \in {\mathcal {E}}_{o}, p\ne q, \end{aligned}$$

(67)

$$\begin{aligned} P \left( 1 - \sum _{\alpha =1}^{n} {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{p, \alpha } \right) P = {\tilde{a}}_{p, \alpha } P {\tilde{a}}_{p, \alpha }^{\dag } \ \ \text {for} \ p \in {\mathcal {E}}_{o}, \end{aligned}$$

(68)

we have

$$\begin{aligned} h_{\mathrm {eff}} = \left( E_{0} - 2s' nd + s'F (n-1) \right) P + s' \sum _{\alpha =1}^{n} \varPsi _{\alpha } P \varPsi _{\alpha }^{\dag }, \end{aligned}$$

(69)

where $\varPsi _{\alpha }$ is defined as

$$\begin{aligned} \varPsi _{\alpha } = \sum _{p \in {\mathcal {E}}_{o}} {\tilde{a}}_{p, \alpha }. \end{aligned}$$

(70)

We treat the two cases, $n \le d$ and $n > d$, separately. We first discuss the former. Let us denote the local ground state energy of $h_{\mathrm {eff}}$ in ${\mathcal {H}}_{F}(\varLambda _{o})$ by $E_{F}^{\mathrm {GS}}$. Since the second term in Eq. (69) is positive semidefinite, we see that

$$\begin{aligned} E_{F}^{\mathrm {GS}} \ge E_{0} - 2s' nd + s' F(n-1). \end{aligned}$$

(71)

Therefore, we obtain

$$\begin{aligned} \min _{n \le F \le 2d} E_{F}^{\mathrm {GS}} \ge E_{0} - 2s' nd + s' n(n-1). \end{aligned}$$

(72)

In fact, we can write down the local ground state of Eq. (69) in ${\mathcal {H}}_{F}(\varLambda _{o})$ for $n \le F \le 2d$ using the method found in [49, 50] and find that $E_{F}^{\mathrm {GS}} = E_{0} - 2s' nd + s' F(n-1)$ for $n \le F \le 2d$. See Appendix B for details. For $F < n$, we note that in the direct sum decomposition (54), ${\mathcal {J}}_{M_{1}, \ldots , M_{n}}$ has at least $(M-F)$ zero subscripts because $\sum _{\alpha =1}^{n} M_{\alpha } =F < n$. Without loss of generality, we assume that $M_{F+1} = M_{F+2} = \cdots = M_{n} = 0$. For any state $|\psi \rangle \in {\mathcal {J}}_{M_{1}, \ldots , M_{F}, 0, \ldots , 0}$, we see that

$$\begin{aligned} \sum _{\alpha =1}^{n} \sum _{p, q \in {\mathcal {E}}_{o}} P {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{q, \alpha } P |\psi \rangle&= \sum _{\alpha =1}^{F} \sum _{p, q \in {\mathcal {E}}_{o}} P {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{q, \alpha } P |\psi \rangle , \end{aligned}$$

(73)

and

$$\begin{aligned} P\left( 1 - \sum _{\alpha =1}^{F} {\tilde{a}}_{p, \alpha }^{\dag } {\tilde{a}}_{p, \alpha }\right) P |\psi \rangle&= {\tilde{a}}_{p, \alpha } P {\tilde{a}}^{\dag }_{p, \alpha } |\psi \rangle \ \ \text {for} \ \alpha =1, \ldots , F \ \text {and} \ p \in {\mathcal {E}}_{o}. \end{aligned}$$

(74)

Using these relations, one finds that

$$\begin{aligned} h_{\mathrm {eff}} |\psi \rangle =\left\{ \left( E_{0} - 2s' F d + s' F (F-1) \right) P + s' \sum _{\alpha =1}^{F} \varPsi _{\alpha } P \varPsi _{\alpha }^{\dag } \right\} |\psi \rangle . \end{aligned}$$

(75)

Therefore, the local ground state energy $E_{M_{1}, \ldots , M_{F}, 0, \ldots , 0}^{\mathrm {GS}}$ in ${\mathcal {J}}_{M_{1}, \ldots , M_{F}, 0, \ldots , 0}$ satisfies

$$\begin{aligned} E_{M_{1}, \ldots , M_{F}, 0, \ldots , 0}^{\mathrm {GS}} \ge E_{0} - 2s' F d + s' F(F-1). \end{aligned}$$

(76)

Here, we note that the lower bound does not depend on the choice of ${\mathcal {J}}_{M_{1}, \ldots , M_{n}}$. From Eq. (76), we obtain

$$\begin{aligned} E_{F}^{\mathrm {GS}} \ge E_{0} - 2s' Fd + s' F(F-1) \ \ \text {for} \ 2 \le F < n. \end{aligned}$$

(77)

Furthermore, the lower bound of Eq. (76) is indeed saturated in ${\mathcal {H}}_{F}(\varLambda _{o})$ for $F < n$. See Appendix B. Noting that $F < n \le d$, the right-hand side of Eq. (77) takes the minimum value when $F = n-1$, and we find that

$$\begin{aligned} \min _{2 \le F < n} E_{F}^{\mathrm {GS}} \ge E_{0} - 2s' (n-1) d + s'(n-1)(n-2). \end{aligned}$$

(78)

Combining Eq. (72) with Eq. (78), we get

$$\begin{aligned} \min _{2 \le F \le 2d} E_{F}^{\mathrm {GS}} \ge E_{0} - 2s' n d + s' n (n-1), \end{aligned}$$

(79)

and its lower bound is non-negative when $0< \nu < \nu _{\mathrm {c}}(n, d)$, where $\nu _{\mathrm {c}}(n, d)$ is defined as

$$\begin{aligned} \nu _{\mathrm {c}}(n,d) = \sqrt{\frac{2d+1 + \sqrt{4(2d-n)(n-1) + (2d+1)^{2}}}{2(2d-n)(n-1)}} \ \ \text {for} \ \ n \le d. \end{aligned}$$

(80)

We next consider the case where $n > d$. In this case, we further discuss the following two cases, $d < n \le 2d$ and $2d < n$. When $d < n \le 2d$ and $n \le F \le 2d$, by repeating the previous argument, we obtain

$$\begin{aligned} \min _{n \le F \le 2d} E_{F}^{\mathrm {GS}} \ge E_{0} - 2s' n d + s' n(n-1). \end{aligned}$$

(81)

On the other hand, for $d < n \le 2d$ and $2 \le F < n$, we have

$$\begin{aligned} E_{F}^{\mathrm {GS}} \ge E_{0} - 2 s' F d + s' F(F-1) \ \ \text {for} \ 2 \le F < n . \end{aligned}$$

(82)

Since $n > d$, the lower bound of Eq. (82) has the minimum value at $F=d$. Therefore, we find that

$$\begin{aligned} \min _{2 \le F < n} E_{F}^{\mathrm {GS}} \ge E_{0} - s' d(d+1). \end{aligned}$$

(83)

Subtraction of the right-hand side of Eq. (83) from the right-hand side of Eq. (81) yields $s'(n-d)(n-d-1)$, which is non-negative because $n > d$. Thus, we obtain the following inequality

$$\begin{aligned} \min _{2 \le F \le 2d} E_{F}^{\mathrm {GS}} \ge E_{0} - s' d(d+1) \ \ \text {for} \ d < n \le 2d. \end{aligned}$$

(84)

When $n >2d$, for all F such that $2 \le F \le 2d$, it holds that $F < n$, and so we have

$$\begin{aligned} \min _{2 \le F \le 2d} E_{F}^{\mathrm {GS}} \ge E_{0} - s' d(d+1) \ \ \text {for} \ n> 2d. \end{aligned}$$

(85)

Summarizing the above inequalities, we obtain

$$\begin{aligned} \min _{2 \le F \le 2d} E_{F}^{\mathrm {GS}} \ge E_{0} - s' d(d+1) \ \ \text {for} \ n > d, \end{aligned}$$

(86)

and find that the lower bound of Eq. (86) is non-negative if $0< \nu < \nu _{\mathrm {c}}(n, d)$, where

$$\begin{aligned} \nu _{\mathrm {c}}(n, d) = \sqrt{\frac{2d+1 + \sqrt{8d^{2} + 1}}{2d(d-1)}}. \end{aligned}$$

(87)

Thus, we have proved Lemma 2. $\square$

Finally, Theorem 2 can be proved using Lemma 1 and 2. We note that $h_{p}$ can be regarded as a finite-dimensional matrix independent of the system size since the local Hamiltonian acts nontrivially only on $4d+1$ sites. This means that the energy levels of $h_{p}$ depend continuously on the parameters. Therefore, Lemma 2 ensures that $h_{p}$ is positive semidefinite when t/s and U/s are finite but sufficiently large. Lemma 1 implies that the ground states of $H_{2}$ are fully polarized states, which proves Theorem 2. $\square$

Appendices

Appendix A: The Linear Independence of Many-Body States

Here, we show the linear independence of the states of the form

$$\begin{aligned} |\varPhi _{\{A_{\alpha }\}, \{B_{\alpha }\}}\rangle = \left( \prod _{p \in A_{1}} a_{p,1}^{\dag }\right) \cdots \left( \prod _{p \in A_{n}} a_{p,n}^{\dag }\right) \left( \prod _{u \in B_{1}} b_{u, 1}^{\dag }\right) \cdots \left( \prod _{u \in B_{n}} b_{u, n}^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(88)

where $\{A_{\alpha }\}$ and $\{B_{\alpha }\}$ are arbitrary subsets of ${\mathcal {E}}$ and ${\mathcal {I}}$, respectively. Let us denote by ${\mathfrak {h}} \cong {\mathbb {C}}^{|\varLambda |}$ the single-particle Hilbert space. We define two types of localized states $\mathsf {a}_{p} = (\mathsf {a}_{p}(x))_{x \in \varLambda }$ and $\mathsf {b}_{u} = (\mathsf {b}_{u}(x))_{x \in \varLambda }$ in ${\mathfrak {h}}$ as

$$\begin{aligned} \mathsf {a}_{p}(x) = {\left\{ \begin{array}{ll} 1 &{} \ \ \text {if} \ x=p, \\ -\nu &{} \ \ \text {if} \ |x-p|=1/2, \\ 0 &{} \ \ \text {otherwise}, \end{array}\right. } \end{aligned}$$

(89)

for each $p \in {\mathcal {E}}$, and

$$\begin{aligned} \mathsf {b}_{u}(x) = {\left\{ \begin{array}{ll} 1 &{} \ \ \text {if} \ x=u, \\ \nu &{} \ \ \text {if} \ |x-u|=1/2, \\ 0 &{} \ \ \text {otherwise}, \end{array}\right. } \end{aligned}$$

(90)

for each $u \in {\mathcal {I}}$. In terms of $\mathsf {a}_{p}(x)$ and $\mathsf {b}_{u}(x)$, the operators $a_{p, \alpha }$ and $b_{u, \alpha }$ can be written as

$$\begin{aligned} a_{p, \alpha }&= \sum _{x \in \varLambda } \mathsf {a}_{p}(x) c_{x, \alpha }, \end{aligned}$$

(91)

$$\begin{aligned} b_{u, \alpha }&= \sum _{x \in \varLambda } \mathsf {b}_{u}(x) c_{x, \alpha }. \end{aligned}$$

(92)

The states $\{\mathsf {a}_{p}\}_{p \in {\mathcal {E}}}$ are linearly independent because for an arbitrary $p \in {\mathcal {E}}$, only $\mathsf {a}_{p}$ has a nonzero component on site p. Therefore, the Gram matrix ${\mathsf {A}}= \left( {\mathsf {A}}_{p, q}\right) _{p, q \in {\mathcal {E}}}$ defined by

$$\begin{aligned} {\mathsf {A}}_{p, q} = \sum _{x \in \varLambda } \mathsf {a}_{p}(x) \mathsf {a}_{q}(x) \end{aligned}$$

(93)

is a regular matrix. Now we introduce the dual operator of $a_{p, \alpha }$ as

$$\begin{aligned} {\bar{a}}_{p, \alpha } = \sum _{q \in {\mathcal {E}}} \left( {\mathsf {A}}^{-1}\right) _{p, q} a_{q, \alpha }, \end{aligned}$$

(94)

for which we can check that

$$\begin{aligned} \{{\bar{a}}_{p, \alpha }^{\dag }, a_{q, \beta }\} = \delta _{\alpha , \beta } \delta _{p, q}. \end{aligned}$$

(95)

Similarly, since the states $\{\mathsf {b}_{u}\}_{u \in {\mathcal {I}}}$ are also linearly independent, the Gram matrix ${\mathsf {B}} = \left( {\mathsf {B}}_{u, v}\right) _{u, v \in {\mathcal {I}}}$ defined by

$$\begin{aligned} {\mathsf {B}}_{u, v} = \sum _{x \in \varLambda } \mathsf {b}_{u}(x) \mathsf {b}_{v}(x) \end{aligned}$$

(96)

is a regular matrix. We define the dual operator of $b_{u, \alpha }$ as

$$\begin{aligned} {\bar{b}}_{u, \alpha } = \sum _{v \in {\mathcal {I}}} \left( {\mathsf {B}}^{-1}\right) _{u, v} b_{v, \alpha }, \end{aligned}$$

(97)

and see that

$$\begin{aligned} \{{\bar{b}}_{u, \alpha }^{\dag }, b_{v, \beta }\} = \delta _{\alpha , \beta } \delta _{u, v}. \end{aligned}$$

(98)

Having introduced the dual operators, we now prove the linear independence of the states Eq. (88). Suppose that

$$\begin{aligned} \sum _{\{A_{\alpha }\}, \{B_{\alpha }\}} f\left( \{A_{\alpha }\}, \{B_{\alpha }\}\right) |\varPhi _{\{A_{\alpha }\}, \{B_{\alpha }\}}\rangle = 0. \end{aligned}$$

(99)

For arbitrary $\{A_{\alpha }\}$ and $\{B_{\alpha }\}$, we operate Eq. (99) from the left with

$$\begin{aligned} \left( \prod _{p \in A_{1}} {\bar{a}}_{p, 1}\right) \cdots \left( \prod _{p \in A_{n}} {\bar{a}}_{p, n}\right) \left( \prod _{u \in B_{1}} {\bar{b}}_{u, 1}\right) \cdots \left( \prod _{u \in B_{n}} {\bar{b}}_{u, n}\right) , \end{aligned}$$

(100)

and then, using the anticommutation relations (95) and (98), we have $f\left( \{A_{\alpha }\}, \{B_{\alpha }\}\right) = 0$. Thus, the states $|\varPhi _{\{A_{\alpha }\}, \{B_{\alpha }\}}\rangle$ are linearly independent.

Appendix B: The Local Ground State

In this section, we construct the local ground state of Eq. (69) in ${\mathcal {H}}_{F}(\varLambda _{o})$ explicitly. To do this, we employ the method used in [49, 50]. For $n \le F \le 2d$, we define a state $|\varPsi _{F}\rangle \in {\mathcal {H}}_{F}(\varLambda _{o})$, as

$$\begin{aligned} |\varPsi _{F}\rangle = {\mathcal {N}} P \left( \prod _{\alpha =1}^{n} \varPsi _{\alpha }^{\dag } \right) \left( \prod _{q \in \{q_{1}, \ldots , q_{F-n}\} \subset {\mathcal {E}}_{o}} {\tilde{a}}_{q, \alpha _{q}}^{\dag } \right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(101)

where $q_{1}, \ldots , q_{F-n}$ are distinct lattice sites in ${\mathcal {E}}_{o}$ and $\alpha _{q}$ can take an arbitrary color. Here ${\mathcal {N}}$ denotes a normalization constant. From the operator identity

$$\begin{aligned} P \varPsi _{\alpha }^{\dag } P \varPsi _{\alpha }^{\dag } = 0, \end{aligned}$$

(102)

one finds that

$$\begin{aligned} s' \sum _{\alpha =1}^{n} \varPsi _{\alpha } P \varPsi _{\alpha }^{\dag } |\varPsi _{F}\rangle = 0, \end{aligned}$$

(103)

which means that the state $|\varPsi _{F}\rangle$ is the local ground state of Eq. (69) in ${\mathcal {H}}_{F}(\varLambda _{o})$ when $F \ge n$. We thus conclude that

$$\begin{aligned} E_{F}^{\mathrm {GS}} = E_{0} - 2 s' nd + s' F(n-1) \ \ \text {for} \ n \le F \le 2d. \end{aligned}$$

(104)

Next, we discuss the case when $2 \le F < n$. In this case, it suffices to consider the effective Hamiltonian $h_{\mathrm {eff}}$ within the subspace ${\mathcal {J}}_{M_{1}, \ldots , M_{F}, 0, \ldots , 0}$, in which $h_{\mathrm {eff}}$ acts as Eq. (75). As in the previous case, we can construct the local ground state in ${\mathcal {J}}_{M_{1}=1, \ldots , M_{F}=1, 0, \ldots , 0}$ as

$$\begin{aligned} |\varPsi \rangle = {\mathcal {N}}' P \left( \prod _{\alpha =1}^{F} \varPsi _{\alpha }^{\dag }\right) |\varPhi _{\mathrm {vac}}\rangle , \end{aligned}$$

(105)

where ${\mathcal {N}}'$ is a normalization constant. Again, from the operator identity (102), we have $s' \sum _{\alpha =1}^{F} \varPsi _{\alpha } P \varPsi _{\alpha }^{\dag } |\varPsi \rangle = 0$ and then obtain

$$\begin{aligned} E_{F}^{\mathrm {GS}} = E_{0} - 2s' F d + s' F(F-1). \end{aligned}$$

(106)