Probing the bond order wave phase transitions of the ionic Hubbard model by superlattice modulation spectroscopy

TL;DR

This paper proposes using superlattice modulation spectroscopy to detect the elusive bond order wave phase in the ionic Hubbard model, leveraging advanced computational methods to identify phase transition signatures.

Contribution

It introduces a novel spectroscopic detection scheme for the bond order wave phase in the ionic Hubbard model, supported by theoretical analysis and simulations.

Findings

Spectroscopy reveals characteristics of phase transitions indicating the bond order wave.

The approach distinguishes between band and Mott insulators.

Insights into excitation spectra of different insulating phases.

Abstract

An exotic phase, the bond order wave, characterized by the spontaneous dimerization of the hopping, has been predicted to exist sandwiched between the band and Mott insulators in systems described by the ionic Hubbard model. Despite growing theoretical evidences, this phase still evades experimental detection. Given the recent realization of the ionic Hubbard model in ultracold atomic gases, we propose here to detect the bond order wave using superlattice modulation spectroscopy. We demonstrate, with the help of time-dependent density-matrix renormalization group and bosonization, that this spectroscopic approach reveals characteristics of both the Ising and Kosterlitz-Thouless transitions signaling the presence of the bond order wave phase. This scheme also provides insights into the excitation spectra of both the band and Mott insulators.

Figures

Equations

$${\hat{H}}_0$$

$$\hat{B}=\sum _{j=1,\sigma }^{L-1}(-1{)}^j\left(c_{j\sigma }^{†}{c}_{j+1\sigma }+\text{h.c.}\right).$$

$$\text{Im}{\chi }_{\hat{B}{\hat{B}}^{†}}(\omega )\sim L^{\frac{7}{4}}\sum _{n=0}^{\infty }{\left(\frac{\Gamma \left(n+1/8\right)}{\Gamma \left(n+1\right)}\right)}^2\delta \left(\hslash (\omega -\Omega (n))\right)$$

$$\text{Im}{\chi }_{\hat{B}{\hat{B}}^{†}}(\omega )\sim \frac{1}{J}\sqrt{{\left(\frac{\hslash \omega }{2\mu }\right)}^2-1}.$$

Full text

In solid state materials, the combination of strong interactions, quantum fluctuations and finely tuned energy scales gives rise to rich physics. For example, in a large class of materials including transition metal oxides Imada et al. (1998), organics Powell and McKenzie (2011) and iridates Witczak-Krempa et al. (2014), the presence of strong on-site repulsion between fermions leads to the suppression of charge motion and to the formation of Mott insulating states. Inducing charge fluctuations around these Mott insulators, e.g. by doping, reveals intricate phase diagrams highlighting the presence of multiple competing orders. Perhaps one of the best known examples is the emergence of $d$-wave superconductivity in high-temperature cuprates at the interface between antiferromagnetic and Fermi liquid phases BednorzMuelller1986 .

Complex states also arise near phase transitions when competing insulating effects are present. In the neighborhood of such transitions, where the strength of the insulating terms is comparable, the effect of smaller terms, such as the kinetic energy, leads to the emergence of metallic phases or exotic correlations. For example, at the interface between the Mott and band insulators, such a competition is believed to play an important role in the ionic to neutral transitions in organic charge-transfer solids Torrance et al. (1981); Nagaosa and ichi Takimoto (1986) and at ferroelectric transitions in perovskites Egami et al. (1993). The ionic Hubbard model, which gained prominence over the last decade, was first developed to explain the physics near these transitions. In this model, the on-site Hubbard repulsion and staggered potential terms induce insulating behavior when taken separately, but when taken together they can compete and give rise to a region of increased charge fluctuations. This region, occuring where these two terms are of comparable strength, is of great interest as this competition leads to the emergence of the bond order wave phase signaled by a spontaneous dimerization of the hopping.

This model is described by the standard Hubbard Hamiltonian to which a staggered potential, with an energy offset $\Delta$ between neighbouring sites, is added

[TABLE]

Here $c_{j,\sigma}^{({\dagger})}$ are the fermionic annihilation (creation) operators and $n_{j,\sigma}$ is the particle number operator on site $j$ with spin $\sigma=\{\uparrow,\downarrow\}$. The amplitude $J$ is the hopping matrix element, $U$ the repulsive on-site interaction strength, and $L$ the number of sites. In one dimension, this model has been the subject of a large number of studies using a variety of techniques including bozonization, density renormalization group, exact diagonalization and quantum Monte Carlo methods Fabrizio et al. (1999, 2000); Gidopoulos, N. et al. (2000); Wilkens and Martin (2001); Torio et al. (2001); Kampf et al. (2003); Zhang et al. (2003); Batista and Aligia (2004); Manmana et al. (2004); Otsuka and Nakamura (2005); Legeza et al. (2006); Tincani et al. (2009). A smaller number of studies have also focused on the excitations of the ionic Hubbard model Aligia and Batista (2005); Hafez and Jafari (2010); Go and Jeon (2011); Hafez Torbati et al. (2014); Hafez-Torbati, Mohsen et al. (2015). Despite initial controversy, theoretical investigations point to the existence, at half filling, of a bond order wave phase occuring around $U\sim\Delta$ characterized by the spontaneous dimerization of the hopping, i.e. the order parameter $B=|\langle\hat{B}\rangle|$,

[TABLE]

The bond order wave spontaneously breaks site-inversion symmetry, and, in the limit of infinite system size, the state is two-fold degenerate with restored bond-inversion symmetry on either even or odd bonds. The bond order wave phase possesses finite charge and spin gaps and is separated, on the one side, from a band insulating state by an Ising quantum phase transition and, on the other side, from a Mott insulator by a Kosterlitz-Thouless (KT) transition. The ground state phase diagram as a function of the onsite interaction strength, $U$, is shown in Fig. 1(a). Despite strong theoretical evidence, the bond order wave has yet to be experimentally detected.

Ultracold fermionic gases provide an appealing novel avenue to detect this state. The ionic Hubbard model was recently realized using fermions loaded into an optical superlattice potential Tarruell et al. (2012). Furthermore, over the last decade, powerful detection techniques were developed to probe the intricate physics of quantum gases. For instance, characteristic excitations of correlated systems can be probed using various spectroscopic techniques such as radio frequency, Raman, Bragg and lattice modulation spectroscopy Bloch et al. (2008); ToermaeTarruell2015 . In particular, for fermionic systems, modulating the lattice potential was proposed as an approach to detect the pairing gap in a superfluid state and the spin ordering in Mott insulating states Kollath et al. (2006), and as a possible thermometer Loida et al. (2015). Numerous theoretical studies have also considered the response of double occupancy to the modulation of the lattice amplitude Kollath et al. (2006); Huber and Rüegg (2009); Sensarma et al. (2009); Massel et al. (2009); Korolyuk et al. (2010); Xu et al. (2011); Tokuno et al. (2012); Tokuno and Giamarchi (2012); Dirks et al. (2014). Experimentally, this response was used to investigate the fermionic Mott insulator Jördens et al. (2008), to probe nearest-neighbor correlations Greif et al. (2011), and to determine the lifetime of doublon excitations Strohmaier et al. (2010). In addition, mapping out higher Bloch bands in a quasimomentum-resolved fashion Heinze et al. (2011) and studying interband dynamics Heinze et al. (2013) were also carried out using lattice modulation. Directional lattice modulation spectroscopy was used to study charge density wave order in the two-dimensional ionic Hubbard model on an honeycomb lattice Messer et al. (2015).

Detecting the bond order wave and characterizing the nature of the neighboring phase transitions requires the development of a technique which couples directly to the order parameter. We demonstrate here, using time-dependent matrix renormalization group (t-DMRG) DaleyVidal2004 ; WhiteFeiguin2004 ; Schollwöck (2011) and bosonization techniques Delfino and Mussardo (1998); Bajnok et al. (2001); Takacs and Wagner (2006), that superlattice amplitude modulation spectroscopy reveals features of both the Ising and KT transitions in one-dimensional systems signaling the presence of the bond order wave phase. This modulation, also provides insights into the excitation spectra of both the band and Mott insulators, the two phases bordering the bond order wave. On the band insulating side, close to the Ising transition, and in the bond order wave, close to the KT transition, one can follow the closing of the excitation gaps.

The proposed detection method relies on inducing a small time-periodic modulation of the superlattice potential described by the perturbation $\hat{H}_{\text{pert}}=A~{}\sin(\omega t)~{}\hat{B}$, with $A$ the strength and $\omega$ the frequency of the modulation directly coupling to the order parameter. By contrast with standard lattice amplitude modulation, here the lattice is modulated in a dimerized fashion. Experimentally such a perturbation is achievable by time-periodically tuning the phase between the two laser waves generating the optical superlattice. The number of excitations created through the modulation is then assessed by monitoring the absorbed energy.

Using t-DMRG, we determine the evolution of the full system $\hat{H}_{0}+\hat{H}_{\text{pert}}$, with open boundary conditions, and compute the energy of the system as a function of time. Typically, the energy shows a quadratic rise at initial times before entering a linear regime and then saturates at later times. We extract the slope of the linear energy increase, a quantity we identify as the energy absorption rate. For the time evolution using t-DMRG we typically keep $120$ states (except at the Ising critical point where we keep $160$). The error analysis are performed by increasing the matrix dimension to $160$ states (and $240$ at the Ising critical point). The time-step of the Trotter-Suzuki time evolution is set to $J\Delta t=0.001\hbar$ and we use $J\Delta t=0.0005\hbar$ to conduct the error analysis. The error bars provided in the figures represent the maximal uncertainties due to the matrix dimension, the time-step and the fit error (as the fit range has been varied).

In the limit of sufficiently weak modulation strength, where linear response theory applies, the energy absorption rate is related to the imaginary part of the dynamic susceptibility associated with the order parameter of the bond order wave, i.e. $\overline{dE}/dt\sim\omega A^{2}\text{Im}\chi_{\hat{B}\hat{B}^{\dagger}}(\omega)$. Here the imaginary part of the dynamic susceptibility $\chi_{\hat{B}\hat{B}^{\dagger}}(\omega)$ is determined by the Fourier transform of the retarded correlation function $\chi_{\hat{B}\hat{B}^{\dagger}}(t)=-i\theta(t)\langle[\hat{B}(t),\hat{B}^{\dagger}(0)]\rangle_{0}$ where $\langle\cdot\rangle_{0}$ denotes the expectation value with respect to the ground state.

We first consider the structure of the energy absorption rate under the effect of superlattice modulation at the Ising critical point $U_{c1}$. Within bosonization, $\hat{B}$ is proportional to the Ising order parameter in the vicinity of the Ising transition Fabrizio et al. (1999, 2000). At the Ising quantum critical point, the scaling dimension of the order parameter is $1/8$. Estimating its dynamic susceptibility at criticality by a scaling argument, we find the absorbed power to diverge as $L~{}\omega^{-3/4}$ in the thermodynamic limit as the modulation frequency decreases to zero signaling the Ising transition. At frequencies lower than the spin gap, the charge fluctuations dominate causing the divergence. One should note that normal lattice modulation fails to detect the existence of the Ising phase transition. In this case, the energy absorption rate does not present a divergence at the Ising transition as this latter modulation scheme does not couple to the bond order wave.

The form of the divergence is affected by the system size. In a finite system, the imaginary part of the dynamic susceptibility is

[TABLE]

with $\hbar\Omega(n)=4\pi\frac{\hbar u}{aL}(n+\frac{1}{16})$ where $a$ is the lattice spacing and $u$ the sound velocity of the low energy excitations. Thus, for a finite system, the divergence will be signaled by the presence of a series of peaks occuring at $\hbar\Omega(n)$ with spectral weight scaling as $\omega^{-3/4}$.

In order to test these low energy predictions, we time-evolve systems of different sizes at the Ising critical point for a range of modulation frequencies and extract the energy absorption rates. At low energies ($\hbar\omega/J<0.15$), our numerical results are in good agreement with bosonization and in particular the peak height follows well the predicted $(\hbar\omega/J)^{-3/4}$ divergence. To identify values of $U/J$ near the Ising critical point $(U/J)_{c1}$, we find for each system size the location of the maximum of the derivative of the order parameter, $\partial|\langle\hat{B}\rangle|/\partial U$ (Fig. 1(b)). Analyzing the energy absorption rate at these values, we then identify the $n=1$ peak footnote2 . We estimate the peak position by fitting a Gaussian to the t-DMRG data for the system of size $L=128$ and modulation amplitude $A=0.0005J$ finding $\hbar\omega_{\text{peak}}/J\approx 0.098\pm 0.005$. From this we extract the sound velocity of excitations $\frac{\hbar u}{aJ}\approx 0.94\pm 0.05$, a typical value for lattice systems. Using the value of the sound velocity $u$ extracted for $L=128$, we then determine the $n=1$ peak positions for the other system sizes which we indicate in Fig. 1(d) by shaded regions. The agreement between these peaks predicted from bosonization and the t-DMRG ones is good. We attribute the disagreement in the peak position and height for $L=64$ to the breakdown of bosonization in the energy range where the corresponding $n=1$ peak appears. While the asymmetry of the $n=1$ peak at $L=192$ and $A=0.001J$ is due to saturation effects in the numerics.

The approach to the Ising phase transition is also detected from the signal obtained on the band insulating side as one can follow the linear closing of the charge gap as the system approaches the critical point. In Fig. 1(c), we present the response deep within the band insulator. The absorption peak is located approximately at $\hbar\omega\sim(\Delta-U)$, it is very sharp and its width scales approximately as $J^{2}/(\Delta-U)$. A strong rise occurs at its left boundary corresponding to a divergence at the lower excitation band edge as also seen in the non-interacting model. The location, width and shape of the peak deep in the band insulator (see solid line in Fig. 1(c)) is understood within an effective model using a Schrieffer-Wolff transformation for the excited states inpreparation . When approaching the Ising transition by increasing $U$, the peak broadens and shifts towards smaller energies. The position of the peak maximum as a function of $U$ is shown in Fig. 2, thus confirming that the charge gap closes when approaching the Ising phase transition.

By comparison to the Ising transition, the location of the KT quantum critical point is harder to pinpoint in finite size systems as the bond order wave order parameter falls off very slowly, see Fig. 1(b). Nevertheless, to the left of the KT transition on the bond order wave side, bosonization predicts a gapped response followed by a sharp increase of the energy absorption rate at twice the minimum energy (mass) of a soliton, $\mu$, corresponding to the creation of a soliton-antisoliton pair

[TABLE]

In the spontaneous dimerized phase, singlets form either on even or odd bonds giving rise to a doubly degenerate ground state and solitons are interpreted as domain walls between the two ground states. The operator $\hat{B}$ is SU($2$) invariant, and can only induce transitions from the singlet ground state to excited states within the same spin sector. The lowest available state is a pair of domain walls (solitons) of opposite spins giving a threshold at twice the soliton mass Lecheminant (2005).

To test this prediction, we time-evolve the system near this second phase boundary at $U=52J$ for a range of modulation frequencies and extract the energy absorption rate. As shown in Fig. 1(e), we observe sharp rises of the absorption rate at multiples of the spin gap value. The spin gap is obtained from static DMRG calculations in different $S_{z}=N_{\uparrow}-N_{\downarrow}$ sectors, $\Delta_{s}=E_{0}(N=N_{\uparrow}+N_{\downarrow}=L,S_{z}=2)-E_{0}(N=N_{\uparrow}+N_{\downarrow}=L,S_{z}=0)$. For $U=52J$ we find $\Delta_{s}\approx 0.049J$ for $L=64$ converged in the matrix dimension. Note, that the normal lattice modulation does not couple to the soliton-antisoliton excitation. Approaching the KT transition from the bond order wave side, the soliton mass becomes smaller and smaller until at the transition the gap closes and a low energy feature, associated with spin excitations, arises on the Mott-insulating side. Hence, superlattice modulation succeeds in signaling the proximity of the KT quantum critical point.

In the Mott insulating phase, the averaged energy absorption rate is predicted by bosonization to be constant at low modulation frequencies in the infinite size limit. In fact, in a finite system we expect equal weight peaks which blend into a constant spectrum when $L\rightarrow\infty$. This spectrum, whose amplitude decreases with increasing $U$, is bounded by the low energy cutoff. Our numerical results for the case of the Mott insulator show at low energies a broad excitation spectrum (see Fig. 3 (b)). We also find that the width and height of the spectrum decrease with increasing $U$. These different features corroborate the predictions from bosonization. Additionally, due to the dominant role played by the spin degrees of freedom in this region of the phase diagram, for $J\ll U-\Delta$, the ionic Hubbard model can be mapped at low energies onto an isotropic Heisenberg chain with exchange interaction $J_{XY}=J_{Z}=\frac{4J^{2}}{U}\frac{1}{1-(\Delta/U)^{2}}$ Nagaosa and ichi Takimoto (1986). As seen in Fig. 3, the width of the spectrum increases linearly with the strength of the Heisenberg exchange interaction confirming the spin nature of this excitation spectrum. Substructures in the spectrum might arise at longer timescales as revealed in Massel et al. (2009) for the homogeneous Hubbard model.

In summary, we demonstrated that superlattice modulation spectroscopy can be used to detect features of both the Ising and KT transitions signaling the presence of the bond order wave phase. This approach would provide a first experimental glimpse into a phase that has evaded detection in the solid state context and highlights the versatility of spectroscopic methods.

Acknowledgment

We thank A. Sheikhan for insightful discussions. We acknowledge financial support of the DFG (SFB 1238 project C05), ERC Phonton (648166).