Mass-Imbalanced Ionic Hubbard Chain

TL;DR

This paper studies a one-dimensional mass-imbalanced ionic Hubbard model, revealing how charge and spin orderings evolve with interaction strength and identifying a weakly first-order phase transition.

Contribution

It introduces a mean-field analysis of a combined mass-imbalanced and ionic Hubbard model, detailing phase transitions and order parameters in the ground state.

Findings

Charge order dominates at low U

Magnetic order dominates at high U

Phase transition at critical U_c

Abstract

A repulsive Hubbard model with both spin-asymmetric hopping (${t_\uparrow\neq t_\downarrow}$) and a staggered potential (of strength $\Delta$) is studied in one dimension. The model is a compound of the mass-imbalanced (${t_\uparrow\neq t_\downarrow}$, ${\Delta=0}$) and ionic (${t_\uparrow = t_\downarrow}$, ${\Delta>0}$) Hubbard models, and may be realized by cold atoms in engineered optical lattices. We use mostly mean-field theory to determine the phases and phase transitions in the ground state for a half-filled band (one particle per site). We find that a period-two modulation of the particle (or charge) density and an alternating spin density coexist for arbitrary Hubbard interaction strength, ${U\geqslant 0}$. The amplitude of the charge modulation is largest at ${U=0}$, decreases with increasing $U$ and tends to zero for ${U\rightarrow\infty}$. The amplitude for spin alternation increases with $U$ and tends to saturation for ${U\rightarrow\infty}$. Charge order dominates below a critical value $U_c$, whereas magnetic order dominates above. The mean-field Hamiltonian has two gap parameters, $\Delta_\uparrow$ and $\Delta_\downarrow$, which have to be determined self-consistently. For ${U<U_c}$ both parameters are positive, for ${U>U_c}$ they have different signs, and for ${U=U_c}$ one gap parameter jumps from a positive to a negative value. The weakly first-order phase transition at $U_c$ can be interpreted in terms of an avoided criticality (or metallicity). The system is reluctant to restore a symmetry that has been broken explicitly.

Figures

Equations

$$H=$$

$$+\frac{\Delta }{2}\sum _{i\sigma }(-1{)}^i{n}_{i\sigma }+U\sum _i{n}_{i\uparrow }{n}_{i\downarrow },$$

$$\begin{array}{rl}\delta {\rho }_c& =-\frac{1}{L}\sum _{i\sigma }(-1{)}^i{n}_{i\sigma },\\ \delta {\rho }_s& =\frac{1}{L}\sum _{i\sigma }(-1{)}^i\sigma n_{i\sigma }.\end{array}$$

$$c_{i\uparrow }^{}\to c_{i\uparrow }^{},c_{i\downarrow }^{}\to (-1{)}^i{c}_{i\downarrow }^{†}$$

$$a_{m\sigma }=c_{2m-1\sigma },b_{m\sigma }=c_{2m\sigma },m=1,\dots ,\frac{L}{2}.$$

$$H=$$

$$-\frac{\Delta }{2}\sum _{m\sigma }(n_{m\sigma }^A-n_{m\sigma }^B)$$

$$+U\sum _m\left(n_{m\uparrow }^A{n}_{m\downarrow }^A+n_{m\uparrow }^B{n}_{m\downarrow }^B\right),$$

$$\begin{array}{rl}{a}_{m\sigma }^{}& =\sqrt{\frac{2}{L}}\sum _k{e}^{ikm}{a}_{k\sigma }^{},\\ b_{m\sigma }^{}& =\sqrt{\frac{2}{L}}\sum _k{e}^{ik(m+\frac{1}{2})}{b}_{k\sigma }^{},\end{array}$$

$$\displaystyle H_{0}=\sum_{k\sigma}\left(a_{k\sigma}^{\dagger},\,b_{k\sigma}^{\dagger}\right)\left(\begin{array}[]{cc}-\Delta/2&\varepsilon_{k\sigma}\\ \varepsilon_{k\sigma}&\Delta/2\end{array}\right)\left(\begin{array}[]{c}a_{k\sigma}^{\phantom{}}\\ b_{k\sigma}^{\phantom{}}\end{array}\right),$$

$${\epsilon }_{k\sigma }^{}=-2t_{\sigma }\cos \frac{k}{2}.$$

$$\begin{array}{rl}{a}_{k\sigma }^{}& =\cos {\phi }_{k\sigma }^{}{\alpha }_{k\sigma }^{}+\sin {\phi }_{k\sigma }^{}{\beta }_{k\sigma }^{},\\ b_{k\sigma }^{}& =-\sin {\phi }_{k\sigma }^{}{\alpha }_{k\sigma }^{}+\cos {\phi }_{k\sigma }^{}{\beta }_{k\sigma }^{},\end{array}$$

$$\tan 2{\phi }_{k\sigma }=\frac{2{\epsilon }_{k\sigma }}{\Delta },\cos 2{\phi }_{k\sigma }=\frac{\Delta }{2E_{k\sigma }},$$

$$E_{k\sigma }=\sqrt{{\epsilon }_{k\sigma }^2+(\Delta /2{)}^2}.$$

$$H_0=\sum _{k\sigma }{E}_{k\sigma }\left(-{\alpha }_{k\sigma }^{†}{\alpha }_{k\sigma }^{}+{\beta }_{k\sigma }^{†}{\beta }_{k\sigma }^{}\right)$$

$$\begin{array}{rl}\delta {\rho }_c& =\frac{1}{L}\sum _{k\sigma }\cos 2{\phi }_{k\sigma }\to \frac{\Delta }{4\pi }\sum _{\sigma }\frac{{\kappa }_{\sigma }K({\kappa }_{\sigma })}{t_{\sigma }},\\ \delta {\rho }_s& =-\frac{1}{L}\sum _{k\sigma }\sigma \cos 2{\phi }_{k\sigma }\to -\frac{\Delta }{4\pi }\sum _{\sigma }\sigma \frac{{\kappa }_{\sigma }K({\kappa }_{\sigma })}{t_{\sigma }},\end{array}$$

$${\kappa }_{\sigma }={\left[1+{\left(\frac{\Delta }{4t_{\sigma }}\right)}^2\right]}^{-\frac{1}{2}}.$$

$$\displaystyle K(\kappa)\sim\left\{\begin{array}[]{ll}\frac{\pi}{2}\left[1+\frac{1}{4}\kappa^{2}+\frac{9}{64}\kappa^{4}+{\cal O}(\kappa^{6})\right],&\kappa\rightarrow 0\,,\\ \\ \frac{1}{2}\log\frac{16}{1-\kappa^{2}},&\kappa\rightarrow 1\,.\end{array}\right.$$

$$\delta {\rho }_c(U)$$

$$\times \left[2K({\kappa }_{\uparrow })K({\kappa }_{\downarrow })-K({\kappa }_{\uparrow })E({\kappa }_{\downarrow })-E({\kappa }_{\uparrow })K({\kappa }_{\downarrow })\right],$$

$$\delta {\rho }_s(U)$$

$$\times \left[K({\kappa }_{\uparrow })E({\kappa }_{\downarrow })-E({\kappa }_{\uparrow })K({\kappa }_{\downarrow })\right],$$

$$H$$

$$-h\sum _i(-1{)}^i{S}_i^{(z)},$$

$$J=\frac{4U{t}_{\uparrow }{t}_{\downarrow }}{U^2-{\Delta }^2},\gamma =\frac{t_{\uparrow }^2+t_{\downarrow }^2}{2t_{\uparrow }{t}_{\downarrow }},h=\frac{2\Delta (t_{\uparrow }^2-t_{\downarrow }^2)}{U^2-{\Delta }^2}.$$

$$H_{\text{mbox}mf}=-\sum _{i\sigma }\left[t_{\sigma }(c_{i\sigma }^{†}{c}_{i+1\sigma }^{}+H.c.)-\frac{{\Delta }_{\sigma }}{2}(-1{)}^i{n}_{i\sigma }\right],$$

$$H_{\text{mbox}mf}=$$

$$\displaystyle-\frac{L}{4}\big{[}(\Delta_{\uparrow}+\Delta_{\downarrow})\delta\rho_{c}+(\Delta_{\downarrow}-\Delta_{\uparrow})\delta\rho_{s}\big{]}.$$

$$\tan 2{\phi }_{k\sigma }=\frac{2{\epsilon }_{k\sigma }}{{\Delta }_{\sigma }},\cos 2{\phi }_{k\sigma }=\frac{{\Delta }_{\sigma }}{2E_{k\sigma }},$$

$$E_{k\sigma }=\sqrt{{\epsilon }_{k\sigma }^2+({\Delta }_{\sigma }/2{)}^2}.$$

$$|{\Psi }_0⟩$$

$$=\prod _{k\sigma }\left(\cos {\phi }_{k\sigma }{a}_{k\sigma }^{†}-\sin {\phi }_{k\sigma }{b}_{k\sigma }^{†}\right)|0⟩.$$

$$\begin{array}{rl}\delta {\rho }_c& =\frac{1}{4\pi }\sum _{\sigma }\frac{{\Delta }_{\sigma }{\kappa }_{\sigma }K({\kappa }_{\sigma })}{t_{\sigma }},\\ \delta {\rho }_s& =-\frac{1}{4\pi }\sum _{\sigma }\sigma \frac{{\Delta }_{\sigma }{\kappa }_{\sigma }K({\kappa }_{\sigma })}{t_{\sigma }},\end{array}$$

$${\kappa }_{\sigma }={\left[1+{\left(\frac{{\Delta }_{\sigma }}{4t_{\sigma }}\right)}^2\right]}^{-1/2}.$$

$$\frac{E({\Delta }_{\uparrow },{\Delta }_{\downarrow })}{L}=\sum _{\sigma }\left[\frac{{\Delta }_{\sigma }({\Delta }_{\sigma }-\Delta )}{8\pi t_{\sigma }}{\kappa }_{\sigma }K({\kappa }_{\sigma })-\frac{2t_{\sigma }}{\pi }\frac{E({\kappa }_{\sigma })}{{\kappa }_{\sigma }}\right]$$

Full text

I Introduction

During the last decade, the study of ultracold atoms in optical lattices has spawned new insight into the complex behavior of quantum many-body systems [Lewenstein_07, ; UCA-Review-MP1, ; UCA-Review-MP2, ]. Atomic gases stored in artificially engineered optical lattices constitute structures beyond those currently achievable in actual materials and, thanks to the easy manipulation of parameters, serve as a playground for the simulation of condensed-matter systems with unconventional properties [Wilczek04, ]. While in the solid-state context the use of simplified models cannot always be justified, because they may neglect relevant degrees of freedom, the clean and precisely controlled environment of ultracold atoms in optical lattices allows a direct mapping of a physical reality onto such models. Moreover, the possibility of manipulating the interaction strength using the Feshbach resonance [UCA_Feshbach_Rev, ] enables the observation of many-body phenomena from weak to strong coupling. Both bosonic [Jaksch_98, ; Greiner_02, ] and fermionic [Hofstetter_02, ; Esslinger_08, ; Esslinger_10, ] Hubbard models have been realized and investigated experimentally.

Optical lattices can be generated in various geometries, including bipartite lattices with different potential minima on the two sublattices [Hemmerich_91, ]. For bosonic atoms a checkerboard potential has been used to study metastable Bose-Einstein condensates with unconventional order parameters [Olschlager_11, ; Olschlager_13, ; Diliberto_14, ]. The relaxation of a bosonic gas, initially prepared in a state with alternating site occupancies, has been investigated for a one-dimensional optical lattice, by switching off and on a period-two potential [Trotzky_12, ]. Fermionic atomic gases in an optical honeycomb lattice – artificial graphene – have been tuned to the state of a Mott insulator for large enough on-site coupling $U$ [Tarruel_12, ]. With the addition of a disparity of well depths on the $A$ and $B$ sites of the honeycomb lattice a competition between the Mott insulator with homogeneous particle density and a band insulator with a modulated density appears [Esslinger_13, ; Messer_15, ]. In the limit of deep wells, such a system can be described by the ionic Hubbard model, where the single-particle levels on neighboring sites differ by an energy $\Delta$. In condensed-matter physics this model has appeared a long time ago in the context of the neutral-to-ionic transition [Nagaosa_86, ], observed in organic mixed-stack compounds [Torrance_81, ; Hubbard_81, ]. Later similar models have been used for describing metal-halogen chains [Gammel_90, ] and transition-metal oxides [Egami_93, ; Egami_94, ]. In these cases the $A$ and $B$ sites are associated with $d$- and $p$-orbitals, respectively.

The conventional Hubbard Hamiltonian is invariant with respect to translations by arbitrary lattice vectors and also has spin SU(2) symmetry. In the ionic Hubbard model the translational symmetry is reduced – the unit cell contains two sites – but the Hamiltonian preserves the SU(2) invariance. SU(2) symmetry can be broken in more than one way. One option for electronic systems is the coupling to a magnetic field, which induces an imbalance between up and down spins. In systems of cold atoms, where the spin degree of freedom may represent two hyperfine levels, an uneven mixture of the two components leads to a similar population imbalance [Zwierlein_06, ]. Another option for cold atoms is to use two atom species with different masses [Wille-et-al-08, ; Taglieber-et-al-08, ]. This mass imbalance leads to different hopping parameters $t_{\uparrow}$ and $t_{\downarrow}$ in the Hubbard model. State-dependent tunneling has also been realized for gases with one atomic species in two hyperfine levels [Mandel_04, ; Messer_15a, ]. In condensed-matter physics mixtures of fermions with different effective masses are realized in rare-earth-metal compounds where a localized $f$ level crosses a wide conduction band [Fazekas_99, ; Riseborough_16, ]. However, in the models advocated for describing these materials the two bands are usually strongly hybridized, and at least the light electrons have an additional spin degree of freedom.

In this paper we study the mass-imbalanced ionic Hubbard model, where translational and spin SU(2) symmetries are both explicitly broken. The Hamiltonian and the order parameters are presented in Sec. II. In Sec. III several limits are discussed, including the non-interacting case (${U=0}$), the small- and large-$U$ limits, the conventional ionic Hubbard model (${t_{\uparrow}=t_{\downarrow}}$, ${\Delta>0}$) and the mass-imbalanced Hubbard model (${t_{\uparrow}\neq t_{\downarrow}}$, ${\Delta=0}$). In Secs. IV and V we concentrate on mean-field predictions for the ground state. The method is explained in Sec. IV and used in Sec. V for describing the phase transition from (dominant) charge to spin order. It is found that the phase diagram for the mass-imbalanced case differs qualitatively from that of the conventional ionic Hubbard model. A brief summary is presented in Sec. VI, which also lists some questions to be addressed in the future.

II Hamiltonian and order parameters

The mass-imbalanced ionic Hubbard chain is defined by the Hamiltonian

[TABLE]

where $c_{i\sigma}^{\dagger}$ and $c_{i\sigma}^{\phantom{}}$ create and annihilate, respectively, fermionic particles at sites $i$ ($i=1,\dots,L$) with spin projections ${\sigma=\uparrow,\downarrow}$, and ${n_{i\sigma}=c_{i\sigma}^{\dagger}c_{i\sigma}^{\phantom{}}}$. The Hamiltonian commutes with the number of particles with spin $\sigma$, ${N_{\sigma}=\sum_{i}n_{i\sigma}}$, and hence with the total particle number, ${N=N_{\uparrow}+N_{\downarrow}}$. We restrict ourselves to the case of half filling, where $N$ is equal to the number of sites $L$. Moreover, we assume that the total magnetization vanishes, i.e., ${N_{\uparrow}=N_{\downarrow}}$. The parameters of the Hamiltonian will be chosen in the range ${t_{\uparrow}\geqslant t_{\downarrow}}$, ${\Delta\geqslant 0}$ and ${U\geqslant 0}$.

A finite ionic term (${\Delta>0}$) breaks translational invariance and induces a density imbalance between neighboring sites, whereas the Hubbard term (${U>0}$) suppresses density inhomogeneities and favors antiferromagnetic ordering. For different values of the hopping amplitudes (${t_{\uparrow}\neq t_{\downarrow}}$) the spin SU(2) symmetry is also explicitly broken. In this case charge and spin modulations are expected to be non-zero for any finite values of $U$ and $\Delta$. To characterize the two broken symmetries we introduce the order parameters $\delta\rho_{c}$ and $\delta\rho_{s}$, defined by the relations

[TABLE]

For ${|\delta\rho_{c}|>|\delta\rho_{s}|}$ charge ordering dominates (ionic phase), while for ${|\delta\rho_{s}|>|\delta\rho_{c}|}$ spin ordering prevails (antiferromagnetic phase).

The canonical transformation

[TABLE]

leaves the hopping term of the Hamiltonian (1) invariant, but exchanges the order parameters, ${\delta\rho_{c}\leftrightarrow-\delta\rho_{s}}$. The ionic potential is replaced by a Zeeman coupling to an alternating magnetic field and the interaction term changes sign. Therefore this transformation maps the repulsive ionic Hubbard model onto the attractive Hubbard model in an alternating magnetic field. The phase diagram for positive $U$ is then readily converted into the phase diagram for negative $U$ by interchanging ionic and antiferromagnetic phases and by associating the parameter $\Delta$ with the amplitude of a staggered magnetic field. It has been argued that the mapping (3) can be very useful for understanding the repulsive Hubbard model by using cold atoms in optical lattices as quantum simulators in the attractive regime [Ho_09, ].

The in-plane alternating spin $(-1)^{i}c_{i\uparrow}^{\dagger}c_{i\downarrow}^{\phantom{}}$ is transformed by Eq. (3) to a pair operator $c_{i\uparrow}^{\dagger}c_{i\downarrow}^{\dagger}$, which implies an intimate relationship between in-plane antiferromagnetic ordering for the repulsive ionic Hubbard model and superconductivity for the attractive Hubbard model in an alternating magnetic field. This field is detrimental for superconductivity and thus favors an alternating charge density for the attractive Hubbard model. Correspondingly, we expect in-plane antiferromagnetism to be suppressed by the ionic potential, and this is indeed found both within mean-field theory and in the large $U$ limit.

It is convenient to introduce unit cells with two sites and operators

[TABLE]

The Hamiltonian (1) then reads

[TABLE]

where $n_{m\sigma}^{A}=a_{m\sigma}^{\dagger}a_{m\sigma}^{\phantom{}},\,n_{m\sigma}^{B}=b_{m\sigma}^{\dagger}b_{m\sigma}^{\phantom{}}$.

III Limiting cases

Before embarking on a discussion of the mean-field ground state of the Hamiltonian (5), we explore certain limiting cases, namely $U=0$, small and large $U$, ${t_{\uparrow}=t_{\downarrow}}$, and ${\Delta=0}$.

III.1 Non-interacting particles: $U=0$

To diagonalize the Hamiltonian (5) for ${U=0}$, we first represent the Wannier operators $a_{m\sigma},b_{m\sigma}$ by Bloch operators $a_{k\sigma},b_{k\sigma}$,

[TABLE]

where $k=\frac{4\pi}{N}\nu$, ${-\frac{L}{4}<\nu\leqslant\frac{L}{4}}$. The Hamiltonian then reads

[TABLE]

where

[TABLE]

The Bogoliubov transformation

[TABLE]

diagonalizes $H_{0}$ if the angles $\varphi_{k\sigma}$ are chosen as

[TABLE]

where

[TABLE]

The transformed Hamiltonian

[TABLE]

has conduction and valence bands separated from each other by an energy gap $\Delta$. For an average number of one particle per site, the case considered here, the occupation numbers in the ground state are ${\langle\alpha_{k\sigma}^{\dagger}\alpha_{k\sigma}^{\phantom{}}\rangle=1}$ and ${\langle\beta_{k\sigma}^{\dagger}\beta_{k\sigma}^{\phantom{}}\rangle=0}$, as in a conventional semiconductor.

For this ground state the order parameters (2) are easily evaluated in the thermodynamic limit, ${L\rightarrow\infty}$, where $\frac{2}{L}\sum_{k}$ is replaced by $\frac{1}{2\pi}\int_{-\pi}^{\pi}dk$. We find

[TABLE]

where $K(\kappa_{\sigma})$ is the complete elliptic integral of the first kind with modulus

[TABLE]

Figure 1 shows the amplitudes of both charge and spin modulations as functions of $\Delta$ for various ratios of hopping parameters. The limiting cases ${\Delta\rightarrow 0}$ and ${\Delta\rightarrow\infty}$ are readily obtained from Eq. (17) using the asymptotic behavior of the elliptic integral,

[TABLE]

If both $t_{\uparrow}$ and $t_{\downarrow}$ are finite $\delta\rho_{c}$ vanishes for ${\Delta\rightarrow 0}$ and tends to 1 for ${\Delta\rightarrow\infty}$, while $\delta\rho_{s}$ vanishes in both limits. The result for ${\Delta\rightarrow 0}$ is different if one of the hopping amplitudes, say $t_{\downarrow}$, vanishes (as in the Falicov-Kimball model). In this case both $\delta\rho_{c}$ and $\delta\rho_{s}$ tend to $\frac{1}{2}$ for ${\Delta\rightarrow 0}$. However, this is a quite singular limit. In fact, for ${\Delta=0}$ (and ${U=0}$) the total energy is completely independent of the distribution of the infinitely heavy particles. This degeneracy is removed by a finite value of $U$.

III.2 Small and large $U$ limits

For ${U\ll\Delta}$ we can use perturbation theory for calculating interaction effects. This is justified, because for ${U=0}$ there is a gap $\Delta$ in the excitation spectrum. To leading order in $U$ we find (for ${t_{\uparrow}\geqslant t_{\downarrow}}$)

[TABLE]

where $E(\kappa_{\sigma})$ is the complete elliptic integral of the second kind. The density modulation $\delta\rho_{c}$ decreases linearly with $U$, since ${E(x)<K(x)}$ for ${x>0}$. For our choice of hopping parameters we have ${\kappa_{\uparrow}\geqslant\kappa_{\downarrow}}$. The function $K(x)$ is monotonically increasing, while $E(x)$ decreases with $x$. Therefore ${\delta\rho_{s}(U)-\delta\rho_{s}(0)\geqslant 0}$, i.e. the staggered magnetization increases as a function of $U$.

If the interaction strength is much larger than the other parameters, i.e. for ${U\gg t_{\uparrow},t_{\downarrow},\Delta}$, the ground-state configuration has essentially one particle at each site (${\delta\rho_{c}\approx 0}$), and the low-energy degrees of freedom are spin $\frac{1}{2}$ operators ${\bf S}_{i},i=1,\dots,L$. Degenerate perturbation theory yields the Hamiltonian [GMJ_16, ]

[TABLE]

with parameters

[TABLE]

For ${t_{\uparrow}=t_{\downarrow}}$ we recover the isotropic Heisenberg chain, where the only effect of $\Delta$ is a slight enhancement of the exchange constant $J$. For ${t_{\uparrow}>t_{\downarrow}}$ and ${\Delta=0}$ we obtain an $XXZ$ chain with anisotropy parameter ${\gamma>1}$, i.e., a uniaxial antiferromagnet [Kocharian_94, ; Fath95, ]. For ${t_{\uparrow}>t_{\downarrow}}$ and ${\Delta>0}$ the $z$-components of the spins are coupled to a staggered longitudinal field of strength $h$. It is interesting to note that this Hamiltonian has been proposed in the context of quasi-one-dimensional easy-axis spin $\frac{1}{2}$ antiferromagnets (such as CsCoCl3), where the staggered field arises in a mean-field treatment of interchain coupling [Shiba_80, ].

The Hamiltonian (24) has been studied with different techniques, but mostly for isotropic exchange (${\gamma=0}$) or for the easy-plane case (${\gamma<1}$), where the interaction and the magnetic field compete [Affleck_99, ; Lou_05, ; Mahdavifar_07, ; Paul_14, ]. In the present case (${\gamma>1}$) the staggered field and the interaction reinforce each other in producing long-range antiferromagnetic order, with up-spins predominantly on even sites and down-spins mostly on odd sites.

III.3 Ionic Hubbard model: $t_{\uparrow}=t_{\downarrow}$

For spin-independent hopping (${t_{\uparrow}=t_{\downarrow}}$) the Hamiltonian (1) embodies the one-dimensional ionic Hubbard model, which has been studied intensively during recent decades [Resta_Sorella_95, ; Ortiz_96, ; fab, ; tosatti_00, ; tor, ; Takada_Kido_01, ; Lou_03, ; brune, ; Manmana_04, ; Zhang_04, ; Solyom_06, ; Tincani_09, ; Hafez_10, ; Uhrig_2014, ]. Early studies found two phases (at half filling), a band insulator for ${U\ll\Delta}$ and a Mott insulator for ${U\gg\Delta}$, with a single quantum phase transition as a function of $U/\Delta$ [Resta_Sorella_95, ; Ortiz_96, ; tosatti_00, ]. The interest increased substantially when a field-theoretic treatment of the ionic Hubbard chain came up with two quantum critical points $U_{c1}$ and $U_{c2}$ [fab, ]. The system was found to be a band insulator (with finite and almost equal spin and charge gaps) for ${U<U_{c1}}$ and a Mott insulator (with gapped charge and gapless spin excitations) for ${U>U_{c2}}$, as expected, but a new phase was found to sneak in for ${U_{c1}<U<U_{c2}}$. At ${U=U_{c1}}$ the charge gap vanishes and one finds metallic behavior, the intermediate phase is a spontaneously dimerized insulator (with finite spin and charge gaps), and at ${U=U_{c2}}$ the spin degrees of freedom exhibit a Kosterlitz-Thouless transition [fab, ]. Subsequent numerical studies [tor, ; brune, ; Manmana_04, ; Tincani_09, ] clarified details of these phases and unambiguously confirmed the scenario of two quantum phase transitions [Manmana_04, ; Tincani_09, ].

The ionic Hubbard chain has been generalized by adding next-nearest-neighbor hopping [Japaridze_07, ] and by enlarging the unit cell to model for instance $MMX$ chains, where $M$ stands for metal atoms and $X$ for halogens [Torio_06, ]. In both cases rich phase diagrams have been found, where the band- to Mott-insulator transition goes through a sequence of unconventional insulating and/or metallic phases. In higher dimensions the phase diagram is more controversial. The ionic Hubbard model for ${d>1}$ at half filling has been studied by various methods, such as dynamical mean-field theory (DMFT) [Jabben_04, ; Randeria_06, ; Craco_07, ; Byczuk_09, ; Wang_14, ; Randeria_14, ; Krishnamurthy_15, ], quantum Monte Carlo [Scalettar_07a, ; Scalettar_07b, ], cluster DMFT [Kancharla_07, ], a variational cluster approach [Chen_10, ] and the coherent potential approximation [Hoang_10, ]. Intermediate phases are routinely found, but depending on applied constraints and computational details different types of order emerge, from metallic [Randeria_06, ; Craco_07, ; Scalettar_07a, ; Scalettar_07b, ; Hoang_10, ], to half-metallic and antiferromagnetic [Randeria_14, ; Krishnamurthy_15, ], as well as insulating and antiferromagnetic [Kancharla_07, ; Byczuk_09, ; Chen_10, ].

III.4 Mass-imbalanced Hubbard model: $\Delta=0$

In the absence of an alternating potential (${\Delta=0}$) the Hamiltonian (1) represents the mass-imbalanced Hubbard model, which was introduced in the early nineties [Brandt91, ] to interpolate between the standard Hubbard model (${t_{\uparrow}=t_{\downarrow}}$) and the Falicov-Kimball model (${t_{\downarrow}=0}$) [Falicov_69, ].

The one-dimensional mass-imbalanced Hubbard model is well understood at half filling in the large $U$ limit where it becomes equivalent to the easy-axis Heisenberg antiferromagnet, with anisotropy parameter $\gamma$ as in Eq. (25). It is generally accepted [Giamarchi_04, ; Mikeska_04, ], that in this case the ground state has long-range antiferromagnetic order. Moreover, the excitations are gapped. The situation is less clear for small and moderate values of $U$. Is there a (Kosterlitz-Thouless) transition to a gapped phase for an infinitesimal anisotropy [Fath95, ] or is the gapless phase extended over a finite region away from the line ${t_{\uparrow}=t_{\downarrow}}$ [Chan_08, ]? Away from half filling a lot of effort has been spent on the problem of phase separation, using rigorous techniques [Ueltschi04, ], bosonization [Wang-Chen-Gu-07, ], weak- and strong-coupling expansions [Barbiero_10, ], or numerical methods [Silva-Valencia07, ; Gu-Fan-Lin-07, ; Farkasovsky08, ]. Several studies were also made for higher dimensions, for instance on the Mott transition [Winograd_11, ; Dao_12, ] and on magnetic correlations [Sotnikov_12, ; Liu_15, ].

For ${t_{\downarrow}=0}$ we recover the Hamiltonian of the Falicov-Kimball model [Falicov_69, ; Gruber_96, ; Freericks_03, ]. Its ground state at half filling has been shown to consist of the most homogeneous configuration of immobile particles, which occupy the sites of one of the two sublattices [Kennedy_86, ; Brandt_86, ]. For other densities spin-up and spin-down particles are segregated for large enough repulsion [Brandt_91, ; Freericks_02, ]. In one dimension a more subtle phase separation can occur already for weak repulsion [Freericks_96, ].

IV Mean-field theory

In the rest of this paper we seek the ground state of the mass-imbalanced ionic Hubbard chain in mean-field approximation. We only consider the broken symmetries included in the Hamiltonian (1), which lead to a period-two density modulation and spin alternation. If we would admit arbitrary densities we would have to include more complicated spatial orderings as well as phase separation and even ferromagnetism [Farkasovsky_12, ], but we restrict ourselves to half filling.

IV.1 Mean-field Hamiltonian and its ground state

We define the mean-field state as the ground state of the single-particle Hamiltonian

[TABLE]

where $\Delta_{\uparrow}$ and $\Delta_{\downarrow}$ are variational parameters. This mean-field Hamiltonian can also be written as

[TABLE]

Therefore $\Delta_{\uparrow}+\Delta_{\downarrow}$ and $\Delta_{\uparrow}-\Delta_{\downarrow}$ can be interpreted as conjugate fields of the order parameters $\delta\rho_{c}$ and $\delta\rho_{s}$, respectively.

$H_{\mbox{\scriptsize mf}}$ is diagonalized by the Bogoliubov transformation (13), but now the ionic potential strength $\Delta$ has to be replaced by $\Delta_{\sigma}$ in Eqs. (14) and (15), i.e.,

[TABLE]

where

[TABLE]

The transformed Hamiltonian is again given by Eq. (16) and its ground state (for ${N=L}$ particles) is

[TABLE]

The order parameters are also calculated in exactly the same way as in Sec. III.1, giving

[TABLE]

where

[TABLE]

IV.2 Mean-field energy and gap equations

It is straightforward to calculate the expectation value of the Hamiltonian (1) with respect to the mean-field state (30). We find

[TABLE]

To calculate derivatives, we use the following relations

[TABLE]

and obtain

[TABLE]

with nonnegative coefficients

[TABLE]

At extrema of the energy the first derivatives vanish, and we get the gap equations

[TABLE]

They admit two different types of solutions, one where ${0<\Delta_{\sigma}<\Delta}$, and a second-one where the signs of $\Delta_{\uparrow}$ and $\Delta_{\downarrow}$ are different. These two possibilities correspond to different relative values of the order parameters (31). If both gap parameters are positive, charge modulations dominate, ${|\delta\rho_{c}|>|\delta\rho_{s}|}$, but if the gap parameters have different signs, the alternating spin is dominant, ${|\delta\rho_{s}|>|\delta\rho_{c}|}$. In general the gap equations have to be solved numerically, but for small on-site interaction we can use an expansion in $U$. To first order in $U$ we can obtain the mean-field corrections to the gap parameters by putting ${\Delta_{\sigma}=\Delta}$ on the r.h.s. of Eq. (37). Inserting these expressions into Eq. (31) we obtain the perturbative expressions (23), which implies that mean-field theory is exact to leading order in $U$.

V Phases and phase transitions

We have solved numerically the gap equations (37) to determine the various phases and phase transitions of the model. The solution is not always unique. In such a case one can compare the energies to find out which solution corresponds to the minimum. Another useful criterion is local stability, which requires the second derivatives of the energy to be positive definite. For values $\Delta_{\uparrow},\Delta_{\downarrow}$ satisfying the gap equations we obtain

[TABLE]

where the coefficients $A_{\sigma}$ are defined by Eq. (36). The eigenvalues of the “dynamical matrix” $\frac{1}{L}\left(\frac{\partial^{2}E}{\partial\Delta_{\sigma}\partial\Delta_{\sigma^{\prime}}}\right)$ are

[TABLE]

At a local minimum both eigenvalues have to be positive. Solutions that do not satisfy this condition have to be ruled out as unstable.

In some limiting cases analytical solutions of the gap equations can also be given. This will be particularly useful for discussing the behavior close to the transition. In the following we choose $t_{\uparrow}=1$.

V.1 From charge to spin order

Figure 2 shows the order parameters $\delta\rho_{c}$, $\delta\rho_{s}$ and the gap parameters $\Delta_{\uparrow}$, $\Delta_{\downarrow}$ for a potential strength ${\Delta=6}$ and two different anisotropies (${t_{\downarrow}=0.5}$ and $0.9$). At ${U=U_{c}}$ there is a weakly first-order transition from dominant charge order to prevailing antiferromagnetism. The steps of the order parameters at $U_{c}$ decrease with decreasing $t_{\downarrow}$ and finally disappear as ${t_{\downarrow}\rightarrow 0}$. The results for other values of $\Delta$ are qualitatively similar and lead to the phase diagram of Fig. 3. The transition line separating the ionic and antiferromagnetic phases depends only weakly on the mass imbalance $t_{\downarrow}/t_{\uparrow}$, especially for not too small values of $U$. For ${U\rightarrow\infty}$ the transition occurs at ${U=\Delta}$, as will be shown analytically below.

The barely visible discontinuities at $U_{c}$ (see Fig. 2) are linked to the fact that the gap parameter $\Delta_{\uparrow}$ changes sign and, instead of passing smoothly through zero, jumps from a small positive value to a small negative value. A state with ${\Delta_{\uparrow}=0}$ would be metallic for the particles with up-spins and semiconducting for the others. From the point of view of symmetry, the full translational invariance of the mean-field Hamiltonian would be restored for one spin component. A solution of Eq. (37) with ${\Delta_{\uparrow}=0}$ indeed exists, with ${\Delta_{\downarrow}=\Delta}$ (for a particular value $U^{*}$). However, when this point is approached, the quantity $A_{\uparrow}$, defined by Eq. (36), diverges, while $A_{\downarrow}$ remains finite. Therefore the eigenvalue $\lambda_{-}$ of Eq. (39) is negative and the solution is unstable. This result is a nice example for the resistance of many-body systems against the restoration of an explicitly broken symmetry.

To see how criticality is avoided, we show the solution of the gap equations close to $U_{c}$ in Fig. 4. The upper and lower lines of the $z$- and $s$-shaped curves correspond to locally stable solutions (${\lambda_{-}>0}$), while the middle section represents an unstable solution (${\lambda_{-}<0}$). In fact, we can easily show that the stability limit (${\lambda_{-}=0}$) coincides with the point where $d\Delta_{\sigma}/dU$ diverges. From Eq. (37) we obtain

[TABLE]

In view of Eq. (39) the prefactor diverges when ${\lambda_{-}\rightarrow 0}$. The value $U_{c}$ is calculated by comparing the energies of the two locally stable branches; the upper line has a smaller energy than the lower one for ${U<U_{c}}$ and the lower one is preferable for ${U>U_{c}}$.

We estimate $U_{c}$ for large values of $\Delta$ where the hysteresis is very small (it vanishes for ${\Delta\rightarrow\infty}$) and ${U_{c}\approx U^{*}}$. At the unstable critical point (with ${\Delta_{\uparrow}=0}$ and ${\Delta_{\downarrow}=\Delta}$) the gap equations boil down to

[TABLE]

where $\kappa_{\downarrow}^{0}=[1+(\Delta/4t_{\downarrow})^{2}]^{-\frac{1}{2}}$. For large values of $\Delta$ the parameter $\kappa_{\downarrow}^{0}$ is very small, $\kappa_{\downarrow}^{0}\approx 4t_{\downarrow}/\Delta\ll 1$, and we find

[TABLE]

and therefore $U_{c}$ indeed tends to $\Delta$ for ${\Delta\rightarrow\infty}$.

To estimate the discontinuity in the gap parameters at $U_{c}$ we determine the three solutions of the gap equations at $U^{*}$, again in the large $\Delta$ limit, where $\Delta_{\uparrow}$ is very small and $\Delta_{\downarrow}$ is close to $\Delta$. To leading orders of $\Delta_{\uparrow}$ and $({\Delta_{\downarrow}-\Delta})$ the gap equations read

[TABLE]

where Eq. (41) has been used. Besides the unstable solution there exist two other solutions, corresponding to the stable branches, namely

[TABLE]

The jump ${2|\Delta_{\uparrow}|}$ vanishes exponentially both for ${\Delta\rightarrow\infty}$ and for ${t_{\downarrow}\rightarrow 0}$. The same exponential behavior is found for $\Delta_{\downarrow}$.

V.2 Fidelity susceptibility

In addition to order parameters an interesting quantity characterizing the ground state is the fidelity susceptibility, which has been used successfully for quantum phase transitions [Zanardi_06, ; Gu_10, ], and also for the crossover from Bose-Einstein condensation of tightly bound pairs to BCS superconductivity [Khan_09, ]. For a given Hamiltonian, which depends on some parameter, typically a dimensionless coupling constant, the fidelity is defined as the overlap of the ground states for two different values of this parameter. In the present case we choose $U$ (for ${t_{\uparrow}=1}$) as characteristic parameter and use the (normalized) mean-field ground state ${|\Psi(U)\rangle}$ to calculate the fidelity

[TABLE]

The fidelity susceptibility is then defined as

[TABLE]

With our mean-field ground state (30) we obtain

[TABLE]

The Bogoliubov angles $\varphi_{k\sigma}$ are linked to the gap parameters $\Delta_{\sigma}$ by Eq. (28), i.e.

[TABLE]

and the derivatives of $\Delta_{\sigma}$ with respect to $U$ are given by Eq. (40). Working out the $k$-sum in Eq. (47) for ${L\rightarrow\infty}$ we finally get

[TABLE]

This expression is valid both for $U<U_{c}$ and for $U>U_{c}$. For a first-order transition at $U=U_{c}$, where the ground states above and below the transition point differ, the fidelity $F(U_{c}-\frac{1}{2}\delta U,U_{c}+\frac{1}{2}\delta U)$ is smaller than 1, even in the limit $\delta U\rightarrow 0$, and $\chi_{F}(U_{c})=\infty$. We expect therefore a sharp line emerging from a smooth background, in contrast to the case of a second-order transition, where $\chi_{F}$ diverges when a critical point is approached.

Figure 5 shows the fidelity susceptibility for ${\Delta=6}$ and various values of $t_{\downarrow}$. For ${t_{\downarrow}<0.99}$ a single peak is found, and it is located at $U_{c}$. We notice that this line is not simply an infinitely sharp peak as in an ordinary first-order transition, but $\chi_{F}$ is strongly enhanced upon approaching $U_{c}$, in agreement with the avoided criticality discussed above. Mathematically, this enhancement comes from the first factor of the r.h.s. of Eq. (40), which can also be written as $A_{\uparrow}A_{\downarrow}/(\lambda_{+}\lambda_{-})$. The eigenvalue $\lambda_{-}$ decreases when $U_{c}$ is approached and is very small at $U_{c}$. Surprisingly, a second peak appears if $t_{\downarrow}$ approaches 1. It has the form of a weak hump for ${t_{\downarrow}=0.999}$, which becomes sharper as $t_{\downarrow}$ approaches 1 and finally is singular at the isotropic point, ${t_{\downarrow}=1}$. We attribute the hump to a crossover due to a rapid change in gap parameters, whereas the singularity for ${t_{\downarrow}=1}$ marks a second-order phase transition, as will be discussed in more detail in the next section.

V.3 Mean-field theory for the ionic Hubbard model

Figure 7 shows the phase diagram for the ionic Hubbard chain. The main differences with respect to the mass-imbalanced case of Fig. 3 are the existence of a pure ionic phase, where ${\delta\rho_{s}=0}$, and the appearance of the intermediate phase. The bifurcation from a single transition line to two lines occurs at $\Delta\approx 3.3373088,\,U\approx 4.2398854$.

We now turn our attention to the ionic Hubbard model (${t_{\uparrow}=t_{\downarrow}}$), where the spin SU(2) symmetry is not explicitly broken. Figure 6 shows the results for the gaps $\Delta_{\sigma}$ and the order parameters $\delta\rho_{c}$, $\delta\rho_{s}$ as functions of $U$. For small $U$ the gap parameters are equal and, correspondingly, there is no antiferromagnetism. For large $U$ the gap parameters have different signs and therefore the antiferromagnetic order parameter $\delta\rho_{s}$ dominates. The nature of the transition between these two regimes depends on the strength of the ionic potential $\Delta$. For small $\Delta$ there is a single first-order transition, as in the mass-imbalanced case. For large $\Delta$ there are two transitions, first a continuous transition at $U_{c1}$ where the difference between gap parameters starts from zero and both of them remain first positive. In this intermediate phase the density modulation still dominates, but there is already antiferromagnetic order. A second transition occurs at a slightly larger value $U_{c2}$ and is of first order. For ${U>U_{c2}}$ antiferromagnetism dominates. It is worthwhile to mention that the two transitions at $U_{c1}$ and $U_{c2}$, respectively, coincide with the two singularities observed in the fidelity susceptibility for ${t_{\uparrow}=t_{\downarrow}}$ (see Fig. 5).

We focus first on the continuous transition at $U_{c1}$ and show that it is triggered by the softening of antiferromagnetic fluctuations. We start within the purely ionic phase, ${U<U_{c1}}$, where a single gap parameter $\Delta_{\uparrow}=\Delta_{\downarrow}=:\bar{\Delta}$ satisfies the equation

[TABLE]

We now enhance $U$ by an infinitesimal amount $\delta U$ and assume the gap parameters to change accordingly,

[TABLE]

Expanding the gap equations (37) up to first order in $\delta U$ and $\eta_{\sigma}$ and subtracting Eq. (50) we get the linear system

[TABLE]

Adding and subtracting these two equations we obtain

[TABLE]

The first line implies that the sum of the gap parameters varies linearly with $\delta U$, while the second line implies ${\eta_{\uparrow}=\eta_{\downarrow}}$, as long as the prefactor is finite. This prefactor is proportional to the eigenvalue $\lambda_{-}$ and therefore a new solution ${\eta_{\uparrow}\neq\eta_{\downarrow}}$ appears at the particular value of $U$ where ${\lambda_{-}=0}$, which we therefore identify with the critical point $U_{c1}$ for the transition to the intermediate phase.

We can also use the parametrization (51) for calculating the eigenvectors of the dynamical matrix (38). One readily finds ${\eta_{\uparrow}=\pm\eta_{\downarrow}}$ for $\lambda_{\pm}$. The eigenvector of $\lambda_{-}$ modifies primarily $\delta\rho_{s}$ and can therefore be associated with a magnetic excitation. This completes the soft-mode picture of the transition at the critical point $U_{c1}$, which is determined by the gap equation (50) (with $U=U_{c1}$) together with the soft-mode condition

[TABLE]

When $U_{c1}$ is approached from below, the order parameters are

[TABLE]

The treatment of the region above the critical point is slightly more tricky because we have to use a higher-order expansion. The calculation, detailed in the Appendix, gives the following result for the critical region

[TABLE]

with proportionality factors given by Eqs. (81) and (82), respectively. The behavior of Eqs. (55) and (56) agrees perfectly well with numerical calculations in the vicinity of $U_{c1}$, as shown in Fig. 8.

In Fig. 7 the intermediate phase appears as a tiny fjord, separating the dominating ionic and antiferromagnetic regions. To find out how this fjord widens further up we calculate the asymptotic behavior of the two transition lines. The transition at $U_{c2}$ is exactly of the same nature as the single transition at $U_{c}$ found in the case of unequal hopping parameters and therefore we can use the large $\Delta$ estimate of Eq. (42) for $U_{c2}$. The critical point $U_{c1}$ is determined by two equations, the gap equation (50) and the soft-mode condition (54). For large $\Delta$ $U_{c1}$ is also large and therefore the parameter $\bar{\kappa}$ has to be small. Expanding $\bar{\Delta}$ and the elliptic integrals in powers of $\bar{\kappa}$ we get

[TABLE]

together with

[TABLE]

and therefore

[TABLE]

We conclude that the size of the intermediate region steadily grows and diverges for ${\Delta\rightarrow\infty}$,

[TABLE]

Coming back to the fidelity susceptibility (see Fig. 5) we recall that for the ionic Hubbard model (${t_{\uparrow}=t_{\downarrow}}$) we found two peaks, one at $U_{c2}$ with a similar shape as that found for the mass-imbalanced case and one at $U_{c1}$ with a very different form. When approaching $U_{c1}$ from above we find a divergence characteristic of a second-order phase transition. This is not seen when approaching $U_{c1}$ from below. We can readily understand this interesting result from Eq. (40), because for ${U<U_{c1}}$ the singularity due to the first factor (triggered by the soft mode) is compensated by the second factor and we get

[TABLE]

Does the peculiar behavior of the fidelity susceptibility for ${U<U_{c1}}$ recur in other quantities, for instance in the staggered magnetic susceptibility? To calculate the response to a staggered magnetic field we add the perturbation

[TABLE]

to the Hamiltonian (1). Its expectation value with respect to the mean-field ground state is

[TABLE]

and the derivatives with respect to the gap parameters are

[TABLE]

Therefore the gap equations read

[TABLE]

We seek a solution as in Eq. (51) where $\bar{\Delta}$ is the gap parameter for ${h=0}$. For an infinitesimal value of $h$ we find

[TABLE]

The staggered magnetic response is measured by the susceptibility

[TABLE]

The eigenvalue $\lambda_{-}=\bar{A}(1-4U\bar{A})$ vanishes as ${(U_{c1}-U)}$ because $\bar{A}$ is finite and continuous close to $U_{c1}$. Therefore the staggered susceptibility diverges as ${(U_{c1}-U)^{-1}}$ and does not exhibit any anomalous behavior.

At this point we have to compare the mean-field phase diagram of Fig. 7 with that obtained using more elaborate methods (mentioned in Sec. III.3). There is qualitative agreement for the ionic phase, the location of the transition line and even the existence of an intermediate phase. However, the nature of both antiferromagnetic and intermediate phases differs. While we found long-range magnetic order in both phases, only quasi-long-range order is obtained by essentially exact treatments for the antiferromagnetic phase. This can be readily understood in terms of quantum fluctuations, which are neglected in mean-field theory. However, in the intermediate phase the ordering found both in numerical simulations and in analytical treatments for small $U$ is completely different from our mean-field result; this disparity cannot be simply attributed to order parameter fluctuations.

V.4 The ionic Falicov-Kimball model: $t_{\downarrow}=0$

For ${t_{\downarrow}\rightarrow 0}$ the gap equations are simplified considerably. With the limiting behavior

[TABLE]

the first gap equation yields ${\Delta_{\uparrow}=\Delta-U}$, and the order parameters are given by

[TABLE]

Therefore ${\delta\rho_{c}+\delta\rho_{s}=1}$, which implies ${n_{i\downarrow}=1}$ on odd sites and 0 on even sites. We have verified, using both analytical and numerical tools, that this mean-field solution reproduces the exact ground state of this ionic Falicov-Kimball model, namely the immobile particles occupy the potential minima and the mobile particles experience an effective ionic potential of strength ${\Delta-U}$. Both for ${U<\Delta}$, where ${\delta\rho_{c}>\delta\rho_{s}}$, and for ${U>\Delta}$, where ${\delta\rho_{s}>\delta\rho_{c}}$, the gap parameter $\Delta_{\uparrow}$ is finite, but for ${U=\Delta}$ it vanishes, and therefore the system is metallic just at the critical point ${U_{c}=\Delta}$. We can also easily show that the slopes of the order parameters diverge logarithmically by approaching the critical point. If $U$ tends to $U_{c}$ the parameter $\kappa_{\uparrow}$ tends to 1. Using the asymptotic behavior of the elliptic integrals together with Eq. (34) we find

[TABLE]

VI Summary and outlook

In our study of the one-dimensional mass-imbalanced ionic Hubbard model we have found a subtle interplay of explicit and spontaneous symmetry breakings. We have explored in detail the competition between two types of order, namely a modulation of the particle density induced by the breaking of translational symmetry and alternating magnetic order originating from the broken SU(2) symmetry. We have limited ourselves to the ground state at half filling (one particle per site), using mostly mean-field theory. For given parameters ${\Delta>0}$ (ionic potential strength) and ${t_{\downarrow}/t_{\uparrow}\neq 1}$ (hopping imbalance) we found a weakly first order transition as a function of the on-site interaction $U$ from a dominant density modulation to prevailing spin order. We have interpreted this transition in terms of an avoided criticality. Tiny steps appear in the two order parameters at a point $U_{c}$. The transition is very clearly seen as an infinitely sharp peak in the fidelity susceptibility.

Due to the broken symmetries, spin-dependent band gaps $\Delta_{\uparrow}$ and $\Delta_{\downarrow}$ appear which behave very differently as functions of $U$. At $U_{c}$ the gap for the light particles becomes very small, especially for large ionic potential strength $\Delta$, while the gap for the heavy particles remains large. This opens a very unusual avenue for spin-selected transport.

For symmetric hopping the existence of an intermediate phase is well established on the basis of advanced analytical and numerical methods. We have verified that such an intermediate phase also appears in mean-field theory (for large values of $U$), with a second-order transition at $U_{c1}$ and a (weakly) first-order transition at $U_{c2}$. However, for an arbitrarily small mass imbalance the second-order transition is replaced by a crossover and the intermediate phase disappears.

Mean-field theory is quite generally expected to be a bad approximation for one-dimensional systems, but in the present case fluctuations are strongly suppressed due to pre-existing gaps in the charge and spin excitation spectra and we expect our phase diagram to be essentially correct. Our mean-field result for the ground-state energy agrees to leading order in $U$ with perturbation theory, and it also reproduces the exact ground state in the Falicov-Kimball limit (${t_{\downarrow}\rightarrow 0}$).

The mean-field parameters $\Delta_{\uparrow}$ and $\Delta_{\downarrow}$, which determine the order parameters and also the gaps in the excitation spectrum behave very differently as $U_{c}$ is approached. The gap parameter for more mobile particles, $\Delta_{\uparrow}$, becomes very small and jumps to a negative value at $U_{c}$, while $\Delta_{\downarrow}$ remains large. Interestingly, the larger the strength of the ionic potential $\Delta$, the smaller the minimum gap reached at $U_{c}$. This value also strongly decreases with increasing mass imbalance and tends to zero in the Falicov-Kimball limit.

We have also discussed several limiting cases of the model, some of which have been thoroughly investigated in the recent past. In the limit of spin-independent hopping our mean-field results differ markedly from results obtained previously by large-scale numerical simulations. We attribute this discrepancy to the fact that the single-particle term of the Hamiltonian has ungapped spin excitations. In this case simple perturbation theory is not applicable and mean-field theory cannot be trusted.

For spin-dependent hopping, fluctuation effects are weak in most parts of the phase diagram. However, they may still be relevant close to the transition where one of the gap parameters is very small. Therefore it would be interesting to study the model in the region close to $U_{c}$ using more sophisticated methods, such as bosonization for small $U$ and the mapping onto a spin-1 model for large $U$ (as used in Ref. [Tincani_09, ] for the ionic Hubbard model). This would shed more light on the nature of this peculiar order-order transition and, in particular, on the avoided criticality. We could also learn how the intermediate phase found in the ionic Hubbard model is affected by a small mass imbalance.

Acknowledgements.

M.S. and G.I.J acknowledge support from the Georgian National Science Foundation through the Grant FR/265/6-100/14. D.B. and G.I.J. also acknowledge support from the Swiss National Science Foundation through the Project IZK0Z2-160773. L.J. acknowledges support from the World Federation of Scientist through the National Stipend Program.

Appendix A Critical behavior for the ionic Hubbard model

Here we present the solution of the gap equations for the ionic Hubbard model for $U>U_{c1}$. We use the notation $\bar{\Delta}$ for the gap parameter at $U_{c1}$. Expanding $\Delta_{\sigma}\kappa_{\sigma}K(\kappa_{\sigma})$ in powers of $\eta_{\sigma}=(\Delta_{\sigma}-\bar{\Delta})/\bar{\Delta}$ we find

[TABLE]

where the first three coefficients are given by

[TABLE]

It is convenient to introduce the dimensionless coupling constant $g=U_{c1}/(2\pi t)$. The soft-mode condition (54) then implies $c_{1}=1/g$. Using Eq. (50) we can rewrite the gap equations (37) for $U=U_{c1}(1+\varepsilon)$ as

[TABLE]

Subtracting the two equations we get

[TABLE]

We now use the factorization

[TABLE]

where $p_{n}(x,y)$ are polynomials in $x$ and $y$, in particular

[TABLE]

For $\eta_{\uparrow}\neq\eta_{\downarrow}$ Eq. (74) becomes

[TABLE]

The leading contributions in powers of $\eta_{\uparrow}+\eta_{\downarrow}$ and $\eta_{\uparrow}-\eta_{\downarrow}$ are

[TABLE]

A second relation is obtained by adding the two lines of Eq. (73). The leading terms are

[TABLE]

The last two equations yield the desired relations linking the gap parameters and $\varepsilon=(U-U_{c1})/U_{c1}$ slightly above the critical point $U_{c1}$,

[TABLE]

The order parameters are closely linked to the gap parameters. To leading order in $(U-U_{c1})$ we find

[TABLE]

for the magnetic order and

[TABLE]

for the charge modulation.