Exploring phonon-like interactions in one-dimensional Bose-Fermi mixtures

With the objective of simulating the physical behavior of electrons in a dynamic background, we investigate a cold atomic Bose-Fermi mixture confined in an optical lattice potential solely affecting the bosons. The bosons, residing in the deep superfluid regime, inherit the periodicity of the optical lattice, subsequently serving as a dynamic potential for the polarized fermions. Owing to the atom-phonon interaction between the fermions and the condensate, the coupled system exhibits a Berezinskii-Kosterlitz-Thouless transition from a Luttinger liquid to a Peierls phase. However, under sufficiently strong Bose-Fermi interaction, the Peierls phase loses stability, leading to either a collapsed or a separated phase. We find that the primary function of the optical lattice is to stabilize the Peierls phase. Furthermore, the presence of a confining harmonic trap induces a diverse physical behavior, surpassing what is observed for either bosons or fermions individually trapped. Notably, under attractive Bose-Fermi interaction, the insulating phase may adopt a fermionic wedding-cake-like configuration, reflecting the dynamic nature of the underlying lattice potential. Conversely, for repulsive interaction, the trap destabilizes the Peierls phase, causing the two species to separate.


I. INTRODUCTION
Over the past decades, we have witnessed a rampant growth of experimental methods devised to cool and control dilute gases.The attainment of Bose-Einstein condensates (BECs) [1] was soon to be followed by in-depth explorations of BEC dynamics in light-induced periodic potentials [2,3], paving the way for the groundbreaking demonstration of the Mott-superfluid phase transition.This exemplary quantum phase transition was predicted by the Bose-Hubbard model [4,5] for an atomic gas loaded into an optical lattice [6].Since then, trapped and cooled dilute atomic gases have developed into a versatile laboratory where quantum matter can be studied in a controlled and detailed fashion, constituting one of the most promising platforms for the realization of analog quantum simulators [7][8][9].
Many of the open questions in the field of quantum phase transitions, especially those revolving around the formation and characterization of exotic quantum phases of matter, cannot, however, be directly addressed by the Bose-Hubbard model per se.Although the classical laser field forming the optical lattice is, in principle, dynamic, the back-action between the trapped particles (modeled as beam splitters) and the lattice is typically very weak [10,11].One may then, to a very good approximation, treat the lattice as the outcome of applying a periodic static potential.Certainly, this approximation would be sufficient to emulate many paradigmatic lattice models.However, such classical potentials do not take into account any back-action between the lattice and * e-mail address: Jonas.Larson@fysik.su.se the conducting matter, which we know is essential to explain several phenomena like the Peierls distortion [12] [see Figs. 1 (a) and (b)] and superconductivity [13].In solids, effects of that kind result from electron-phonon interactions emerging from the very nature of a dynamical lattice.
A central theme in our analysis is that assigning tractable degrees of freedom to the lattice renders its description dynamic and enables the simulation of some analog of phonon-like interactions.One possibility is to couple the atoms to the light field of an optical resonator [14,15], where substantial Stark shifts are observed even in the presence of a few photons.However, in current experiments, only a few modes of the resonator actively participate in the light-matter interaction, and it is, therefore, not possible to locally modify the dynamical lattice.Instead, one must turn to multi-mode cavities [16], which still pose difficulties in reaching configurations similar to those encountered in real solids [17].A viable alternative arises when considering atoms directly interacting with an ionic crystal [18].Here, harmonically trapped ions form a Wigner crystal, while additional neutral atoms move within this lattice.This ion-atom system bears obvious similarities to a real solid, albeit being experimentally challenging.In fact, the crystalline structure is not necessary for exploring phonon-like interactions.In mixtures of different atomic species -either different atoms or different internal states of the same atomic species -the interplay between subsystems can lead to intriguing effects.
Coming now to single atomic species, the use of Feshbach resonances allows experimenters to control the strength, and even the sign, of all involved interactions [19][20][21][22][23][24], aiming to probe the resulting phase dia-arXiv:2303.07083v2[cond-mat.quant-gas]15 Feb 2024 gram.In mixtures of bosons and spin-polarized fermions (hereinafter referred to as BF mixtures), it is well known that an attractive BF interaction leads to a so-called pairing phase in the strongly correlated regime.This phase has been studied for weak BF interactions, g bf , in one dimension (1D) [25][26][27], as well as in two [28,29] and three dimensions [30].The phase in question collapses if the interaction becomes too strong, resulting in clumping of the atoms and breaking of translational invariance.The effect of optical lattices on BF mixtures has also been investigated in Refs.[31][32][33][34][35][36][37].For deep lattices and/or very strong interaction, such systems can be described by a BF Hubbard model as they enter an insulating phase of composite fermions [38,39].The physical behavior is typically described within a Wannier-basis expansion for both species, where the bosons can be construed as agents of effective onsite energy shifts for the lighter fermions.The imposed approximations in such a scheme omit certain back-action between the two subsystems in comparison to the self-consistent analysis we employ here.Nevertheless, it is still possible to encounter Peierls phases, supersolids, and charge density waves [28,38].Furthermore, if the repulsive boson-boson interaction is weak, the system can enter a regime of phase separation where the bosons and fermions completely avoid each other [40].
Reporting from the experimental front, an early study focused on how the coherence of an atomic condensateheld in place by a cubic optical lattice -is affected by the presence of fermionic atoms [41].Quantum degeneracy for both bosons and fermions was attained in Ref. [23], where the two species were treated on equal footing, a trend followed in a series of subsequent papers [42][43][44][45].Further experimental investigations have also considered a particular regime emerging for a condensate in weak contact with much lighter fermionic atoms [24,44,46].Here, the fermions effectively induce a so-called Ruderman-Kittel-Kasuya-Yosida (RKKY) long-range bosonboson interaction [47].For attractive Bose-Fermi interaction, it was demonstrated that a self-sustained trap may emerge for those fermions located inside the bosonic condensate [44], an effect which may as well lead to the formation of soliton trains [24].The fermion-mediated spinspin interaction in a spinor condensate has also been recently observed by means of microwave spectroscopy [46].
Carrying on with the thread of recent investigations on collective phenomena brought about by a coherent interaction between Bose-Fermi atomic clouds, in this report, we study a 1D atomic BF mixture where only the bosons are subject to an optical lattice.The motivation for considering such a configuration is drawn from the resemblance to an actual solid -the optical lattice orders the bosons in a crystalline structure.Our spinpolarized fermions experience only the periodic potential arising from the boson atomic density but no externally imposed static potential, with the occasional exception of a harmonic trap.A similar idea, although experimentally more challenging, was proposed in Ref. [48], FIG. 1. (color online) Schematic description of the Peierls distortion.At a fermion filling corresponding to a wavenumber kF , opening up a gap of the dispersion at kF will lower the kinetic energy of the fermions, as depicted in (a).The gap opening results from a lattice distortion of the dynamical lattice, as shown in (b).Here, the period has been doubled through shifting every second atom by a distance δ, and consequently, the size of the Brillouin zone has been halved.Thus, for half filling, we encounter a realization of the Su-Schrieffer-Heeger (SSH) model, while for other fillings, the lattice distortion will generate another periodicity.In the normal phase, the condensate will share the same periodicity as the lattice and predominantly populate the lattice minima (c).In the Peierls phase, for half filling, the bosonic density will alternate between every second site (d).Upon comparison of (b) and (d), we note that the periodicity has been broken in two different ways: in (b), the densities are held fixed, but the locations of the sites have been shifted, while in (d), the locations are fixed and the densities have been altered.The former can be seen as a 'phase modulation,' and the latter as an 'amplitude modulation' of the densities.Finally, in (e), we depict the collapsed phase where all atoms, bosons, and fermions have been compressed to a small region of the lattice.
where the crystalline order for the bosons was established from a rapidly rotating condensate that created a triangular Abrikosov vortex lattice.We focus on the limit of a weak boson-boson (BB) interaction and an optical lattice of amplitude ≤ 20E R (E R is the recoil energy), where the bosonic gas is expected to form a superfluid, see Fig. 1(c).In this regime, we can work within a meanfield approximation for the condensed bosons, which act as a classical dynamical lattice felt by the fermions.The coupled system is solved self-consistently such that every back-action, at the mean-field level, is taken into account (without imposing, for instance, any single-band nor tight-binding approximations).At this hybrid mean-field level, the mixture displays a rich phase diagram.For strong BF interaction, with a coupling rate g bf , the system either experiences a collapse (g bf ≫ 0), in which the two species overlap and populate only a small fraction of the lattice, or a separation (g bf ≫ 0) where instead the two species avoid each other and populate different parts of the lattice.These are first-order transitions, even though translational invariance is spontaneously broken in both cases.For a non-zero interaction beyond a critical coupling g bf = g c bf , an instability occurs such that the period of the state is different from that of the underlying lattice; for example, for half-filling, a periodic doubling is found.The system then transitions from a Luttinger liquid (LL) phase into a Peierls phase via a Berezinskii-Kosterlitz-Thouless (BKT) transition.The Peierls phase manifests itself through a non-zero gap ∆ P at the Fermi wavenumber k F , i.e. ∆ P ̸ = 0 for g bf beyond the critical interaction strength g c bf .
We further present an account of the effects arising due to the presence of a confining harmonic trap.For a sufficiently strong repulsive boson-fermion interaction, it has already been demonstrated that, at the mean-field level, the trap induces a separation of the two species [49].This instance suggests ruling out a Peierls-like phase when a trap breaks the translational invariance.Inside a trap, the atomic densities vary in space, which, in a local density approximation, translates to a local Fermi wavenumber.A spatially varying wavenumber has a direct impact on the phases, such as by smoothing the discontinuous transition between the (possible) Peierls and collapsed/separated phases.Naturally, the presence of a trap entails a larger number of involved wavenumbers and thereby a broadening of the Bragg peaks in a timeof-flight detection.Additionally, we find that the interplay between the trap, the dynamical boson field, and the fermions results in an insulator bearing a clear fermionic "wedding cake structure", which differs from the density profile obtained for a static periodic potential.
Our narrative is structured as follows.Section II introduces the Bose-Fermi system under study and the model we employ for its description.There we also discuss some basic concepts and known results pertaining to the properties of the involved phases.The numerical methods are developed in Sec.III, while the associated results are presented in Sec.IV, first for the translationally invariant case IV A, and subsequently in the presence of a harmonic trap IV B, where we focus on the quantum coherence registered in the one-particle density matrix and on the decay of spatial correlations.Concluding remarks in Sec.V close our paper out, while the two appendixes present two possible paths to the derivation of effective model Hamiltonians taking quantum fluctuations into account for both species, alongside some details underlying the employed hybrid mean-field method.

II. BACKGROUND: MODEL HAMILTONIAN AND PHASES
A. The Hamiltonian The physical system we study is a mixture of two atomic gases, one bosonic and one fermionic.Owing to tight transverse confinement (in the sense of ω ⊥ ≫ ω ∥ ), excited states of the transverse modes can be neglected and the low-energy physical description is quasi-one dimensional [26,50].To begin with, we will take the gas to be of infinite extent in the longitudinal direction, and later consider the effect of a trap.We work with dimensionless quantities, where energies are scaled by the boson recoil energy, E R = ℏ 2 k 2 /2m b with k the lattice wavenumber and m b the atomic boson mass, such that half the optical lattice wavelength sets the characteristic length scale of the problem.The behavior of such a system is modeled by the Hamiltonian where r = m f /m b is the mass ratio between the two species, g b is the effective 1D boson-boson (BB) interaction strength, and g bf is the effective 1D BF interaction strength.Later on, we will assume the bosons to be harmonically confined and subjected to an optical lattice, while the fermions will be harmonically confined but assumed to have negligible interaction with the optical lattice, The last term on the RHS of Eq. ( 2) is the optical-lattice potential, also denoted as , where the distance between potential minima is a = π in dimensionless units.At first, we remove the harmonic confinement and keep only the periodic optical lattice.To ensure that the lattice interacts only with the bosonic atoms, the light frequency should be such that it is not quasi-resonant with any of the dipole transitions -dictated by selection rules -of the fermionic atoms.The atomic field operators ψb,f (x) and ψ † b,f (x) obey the standard commutation relations for bosonic and fermionic creation and annihilation operators.Note that the Hamiltonian commutes with both particle number operators (α = b, f ), For further use, we define the average particle number of bosons/fermions as and the density ν α = an α = πn α as the filling of bosons/fermions per lattice site.
Hereinafter, we assume that the bosonic gas is dense enough to form a (quasi) condensed state, but not so dense that the interaction energy becomes dominant.In terms of the 1D and 3D scattering lengths a s,1D and a s,3D this implies where the healing length of the condensate is defined as ξ = (8πa s,3D ρ b ) −1/2 , with ρ b the condensate density.Furthermore, our approximations are consistent with limiting our attention to optical lattice strengths satisfying V 0 ≲ 20E R .We take the fermions to be lighter than the bosons, i.e. r = m f /m b < 1, which is typically the case for the relevant experiments, e.g. for a Li 6 -Li 7 mixture [23] one has r = 0.86 and for a Li 6 -Cs 133 mixture [44,51] one instead gets r ≈ 0.04.We will assume g b > 0 but allow any sign and magnitude of g bf .Both the BB and BF interactions should be possible to vary in the experiment by deploying Feshbach resonances [19][20][21][22][23].

B. Pairing phase and Peierls instability
For low bosonic densities, the occurrence of a socalled pairing phase for BF mixtures with attractive BF interaction is well known.For weak interaction strength and light fermions, the development of modulated densities can be understood from the fact that the fermion-mediated long-range interaction between two bosons is an ultracold atomic analog to the RKKY interaction [47,52,53].These spatial modulations of the boson density alter the potential landscape for the fermions, and subsequently instigate back-action between the two species.The spatial modulations will occur at the wavenumber 2k F if we assume an RKKY interaction from a "flat" background of fermions with a Fermi wavevector k F .In [25], a 1D BF mixture without optical lattice was investigated using bosonization and a transition was identified between a two-component gapless LL phase and a gapped pairing phase with periodic density modulations.The transition was found to "flow" towards the BKT fixed point with the gap opening as where Θ(x) is the Heaviside step function, and g c bf is the critical Bose-Fermi coupling strength.A study using discretization of the spatial coordinate and the numerical DMRG method also predicted a stable LL phase for small negative g bf [27].For repulsive interaction (g bf > 0), the aforementioned works found no pairing phase.However, in Ref. [26], a pairing phase was indeed reported upon employing the method of a random phase approximation.
A gap opening at 2k F and the corresponding periodic modulation in the density are reminiscent of the Peierls distortion.As first demonstrated by Peierls [12], a metal is unstable towards lattice distortions, i.e. a small displacement of the atoms.This phenomenon can be understood by treating the atomic displacements on the meanfield level in the well-known SSH model.Opening up a gap ∆ P in the fermion dispersion at the Fermi wavevector k F , as depicted in Fig. 1 (a), lowers the total energy and corresponds to a lattice modulation of wavenumber k P = 2k F , as sketched in Fig. 1 (b).There are however also some discrepancies.The gap ∆ P was derived from a mean-field treatment of the ions, while the pairing gap ∆ is the gap of excitations in the full quantum system.Still, it is interesting to observe that the gap in the SSH model also has an exponential form similar to (7).Moreover, in the original work of Peierls, the instability was demonstrated for any non-zero interaction [12], while the pairing effect arises as a consequence of a phase transition occurring at finite interaction strength.
An optical lattice is, on intuitive grounds, expected to stabilize the pairing phase against a collapse, since it renders the bosons less mobile by generating a larger effective mass.On the other hand, using the same argument, the lattice may also extend a possible metallic phase to larger regions in the phase diagram.We may note that for BF mixtures where both species are subject to an optical lattice, it has been reported that the pairing phase appears for both repulsive as well as attractive BF interactions [35].Furthermore, one could in principle expect a beating between the involved length scales, the Fermi wavenumber k F and the optical lattice wavenumber k, that could, in principle, give rise to socalled Devil's staircase structures [54].
A last remark on the periodically modulated Peierls phase is in order.We have envisioned the boson superfluid as the agent creating a dynamical lattice for the fermions.We may consider the opposite viewpoint of a condensate living in a partially dynamic background.Within this perspective, the Peierls phase is reminiscent of a supersolid [28,30], i.e., a superfluid state that has spontaneously broken the periodicity of the Hamiltonian.Somewhat similar scenarios have been studied in atomic condensates confined within optical resonators [55,56].

C. Stability of the Bose-Fermi mixture
For sufficiently strong BF interaction, the system is unstable towards long-wavelength fluctuations.The stability of BF mixtures without an optical lattice has already been investigated in, e.g., Refs.[26,50].Applying a hydrodynamic (mean-field) description of both fermions and bosons, a linear stability analysis yields the condition where the constant g f = π 2 /r appears as an effective fermion interaction strength.Note that the different definition compared to Ref. [50] is due to our use of rescaled dimensionless units.The fact that a high fermion density stabilizes the mixture is particular to the 1D case.
Attractive interaction leads to a collapse of both species, while a repulsive interaction induces a phase separation where the bosons and fermions either avoid each other or lump together to form BF soliton mixtures [34,57,58].More precisely, for repulsive boson-fermion interaction, the bosons can form a dark soliton (density dip), filled by a bright fermionic soliton.For attractive interaction, both species form a bright soliton.This is reminiscent of the Townes solitons predicted for electromagnetic waves [59], and recently demonstrated in ultracold atomic Bose-Bose mixtures [60,61].Soliton solutions are known to be unstable beyond the mean-field approach [62,63] but BF mixtures can form stable selfbound systems (e.g., in Ref. [64],) which in 1D are stable within the mean-field approximation as well as to higher order [65].In the course of our analysis, we will find out that an optical lattice increases the range of g bf values for which the mixture is stable.

A. A hybrid mean-field approximation
To investigate the pairing phase, we employ a "hybrid" approach that intertwines quantum and classical dynamics, a method which has already been followed in Refs.[66,67].We note that this path bears similarities to the mean-field approximation employed to derive the SSH model.It is also akin to the DFT method [66].Due to its affinity with the SSH model, we will refer to the gapped phase in the hybrid approximation as the Peierls phase.For the ground state, we adopt the ansatz where |Ψ f ⟩ is a general state of N f fermions, and |ψ b ⟩ is a (generalized) coherent state of the bosons satisfying [68] ψb with the complex variable ψ b (x) called the condensate wavefunction or order parameter.We impose the normalization 1 for a finite system and thereby factor out the total boson number N b .If m f /m b = r ≪ 1, it is justified to assume an adiabatic evolution, similar to the Born-Oppenheimer approximation in molecular physics [69], where the lighter fermions adjust approximately instantaneously according to the heavier bosons.Under such conditions, any gauge potential [70], characterizing non-adiabatic corrections and arising due to back-action between the condensate and the fermions, can be neglected.The fermionic part of the Hamiltonian can then be diagonalized while keeping the bosonic degrees of freedom fixed, thus acting as an "adiabatic potential".The approximation is equivalent to replacing the bosonic field operators with the condensate wavefunction in the Hamiltonian, yielding the hybrid operator where the mean-field energy functional of the bosonic part is and the second part is just the fermionic part of the Hamiltonian written in diagonal form, ϵ f,n being the eigenvalues of the Hartree equation for the fermion orbitals, Note that the solution of this Hartree equation exhibits a functional dependence on the boson density, hence from now on we will be writing 11), a nonlinear Schrödinger equation for ψ b can be derived, where the expectation value of the fermion density can be calculated as To determine the ground state, we resort to Eq. ( 14) and numerically propagate an initial condensate wavefunction in imaginary time, employing the split-operator method [71].This procedure amounts to solving the eigenvalue problem (13) at each time step.
To assess the Peierls phase for a system of infinite extent, we will use the following prescription: Assume a filling ν f = 1/2, so that we expect density modulations with a period of two sites.Since the modulations are commensurate with the optical lattice, the entire problem is periodic and the solutions of the Hartree equation are Bloch waves, with a reduced unit cell in reciprocal space due to the doubled periodicity in the Peierls phase.Again, for the system evolution, we propagate in imaginary time, solving the Hartree equation at each time instant.The case of incommensurate fermion filling is interesting in its own right, as it may lead to coexisting spatially separated regions of commensurable and incommensurate phases [72].However, our method does not allow for an incommensurate filling at present, whence we leave the problem aside for later investigation.
For the Peierls phase, we appeal to the same order parameter as in the SSH model, namely the energy gap of fermionic excitations above the Fermi surface, which in turn is identical to the gap of the hybrid Hamiltonian (11).As we discuss in Sec.II B, this quantity should approximate the excitation gap of the entire system.

B. Thomas-Fermi model for the fermions
To ascertain the transition to the collapsed or separated phases, we resort to a Thomas-Fermi model for the fermions (not to be confused with the Thomas-Fermi approximation); for an introduction, see Refs.[73,74] and for applications, see Refs.[64,75,76].We define a classical field ψ f (x) = n f (x)/N f , called the Thomas-Fermi wavefunction, where n f (x) has been defined in Eq. ( 15) under the normalization dx |ψ f | 2 = 1.The kinetic energy is approximated as a functional This form can be motivated by means of a dimensional analysis, and is derived in detail in [73,74].We obtain the Thomas-Fermi energy functional as where the first term on the right is given by ( 12) and the second by Taking the functional derivative, we find an equation of a two-component Gross-Pitaevskii type, (20) Employing once more the split-operator method [71], we then evolve an initial state in imaginary time to determine the ground state.Due to the mass separation between the two species, one may argue that an adiabatic elimination of the fermion field should be justified.This would result in an effective bosonic model where the fermion-mediated RKKY-like boson-boson interactions would appear.Such an approximation, however, has a negligible benefit to the numerical solution, whence we solve the full mean-field model formulated by Eqs.(20).

A. Translation-invariant BF mixture
We first focus on the periodic case in the absence of a trap, where we set V trap = 0.This idealization targets the thermodynamic limit in the system response, which is meaningful when assessing the universal properties of the phase transitions reported herein.A Bose-Fermi mixture with no lattice was studied in Ref. [25], where it was found that, for negative g bf , a BKT transition separates a LL from a Peierls phase; that gap had the form dictated by Eq. ( 7).
In the upper two frames of Fig. 2, we depict the numerically extracted phase diagrams; frame (a) shows a region in the (g bf , g b ) plane when keeping V 0 = 0 and r = 0.04 constant, while frame (b) displays a region in the (g bf , V 0 ) plane for g b = 0.4 and r = 0.04 constant.Both (a) and (b) show the region of stability of the Peierls phase towards collapse or separation (black solid lines), located via a calculation using the Thomas-Fermi approximation and compared to the mean-field prediction (8) (black dotted lines).The transition is identified by detecting the discontinuity in the overlap between boson and fermion densities.The LL-to-Peierls transition (red solid lines) is calculated via the hybrid method, and the critical coupling g c bf is obtained from a least square fit of the expression (7) to the numerical data.From frame (a) we conclude that the region of stability has an extent which abides by the inequality ( 8), but more importantly that our method predicts a Peierls (i.e.pairing) phase for g bf > 0, in disagreement with Ref. [25].The reason for such a disparity could be simply attributed to the fact that the hybrid method entails the unphysical assumption that the solution is periodic, which hinders the detection of long-wavelength fluctuations.From the literature on Bose-Fermi-Hubbard systems, where both gases are subject to an optical lattice, we expect however an optical lattice to stabilize the Peierls phase for repulsive interaction [77].From frame (b) we infer that this is indeed the case: the optical lattice stabilizes the Peierls phase beyond the boundary set by condition (8) for a homogeneous configuration.More precisely, the extent of this phase is growing with the optical lattice potential while the LL, collapsed, and separated phases are shrinking for the shown parameter ranges.Furthermore, the Peierls gap also grows with increasing lattice amplitude, meaning that it should be easier to observe in an experiment with a stronger optical lattice, provided the full system does not enter an insulating phase due to the inclusion of an optical lattice.
To illustrate the type of transition occurring between the LL and Peierls phases, Fig. 2(c) shows three examples of the Peierls gaps as a function of g bf , for V 0 = 0, 4, and 8, keeping g b = 0.4 constant.As demonstrated by the dotted curves, close to the critical value, the numerically extracted gaps fit very well the analytical expression (7) for the gap of a BKT transition.Away from the  8) is plotted with dotted black lines.The critical coupling g c bf for the LL-to-Peierls transition is calculated from a least square fit of the Peierls gap to the expression (7), using the hybrid mean-field method.(b) Same as for (a) but in the (g bf , V0) plane.Here it becomes clear that the presence of the optical lattice alters the phase diagram quantitatively, in particular through stabilizing the Peierls phase.(c) The Peierls gap as a function of g bf , once more calculated using the hybrid code, for V0 = 0 (yellow), V0 = 4 (green), and V0 = 16 (blue), fitted to the BKT formula (7) (dotted lines).(d) Boson and fermion densities in the Peierls phase for the parameters V0 = 1, g b = 0.4, and g bf = −1.In all plots, we consider half filling, ν f = 1/2.critical point, there is a notable "knee" feature, which occurs when the bosons only populate every other site, and solely for V 0 ̸ = 0.
Figure 2(d) displays the density of bosons and fermions within the unit cell, as a typical example for a particular point picked in the phase diagram: V 0 = 1, g b = 0.4, and g bf = −1.The doubling of the periodicity is in evidence.Contrary to most earlier studies, our method is capable of capturing both amplitude and phase modulations, however, we only find modulations in amplitude.For the repulsive case, g bf > 0, the situation is similar to that shown in (d), but with the two densities being "out-of-phase" instead.
In order to experimentally probe the phase diagram, one must detect the order parameter or any other quantity capable of telling the different phases apart.The gap could, in principle, be extracted via pump-probe experiments [78,79], while one could also imagine transport experiments to differentiate insulating and conduction phases, i.e.Peierls vs. LL.As will be evident in the next section, the trap induces a varying Fermi wave number k F which makes the Peierls gap opening not that pronounced.In a time-of-flight detection, for instance, the Bragg peaks are smeared out, and the particular one corresponding to the gap is almost invisible.Hence, in a realistic setup, including the trap, time-of-flight measurements might not turn out to be the most optimal scheme.Instead, the onsite densities, which display al- ternating variations, could be accessed via single-site resolution detection [80].This technique has been successfully implemented for cold fermionic gases, including even detection of higher order correlators [81][82][83].

B. BF mixture in a harmonic trap
Current experiments rely on the application of harmonic potentials for the attainment of sufficiently long trapping times.However, the trap may have a relatively small frequency yet still confine the atoms, such that the system locally experiences a periodic potential.To study these more realistic situations, we include a harmonic trap in our analysis but still consider a very tight trans-verse confinement and consequently a quasi-1D configuration.
We have used the hybrid method to investigate the different phases in the aforementioned setup.Figure 3 presents the numerical results obtained for the densities, along the main diagonal of the one-particle density matrix.In the upper row, we present the real-space densities, and in the lower row, the fermion density per site is depicted, which could, in principle, be accessed in an experiment, as mentioned in the last section.Illustrated in the left column is a novel insulating state not present in the translationally invariant system.We refer to this profile as a "wedding cake", drawing an analogy to the wedding cake-like structure found in the harmonically confined Bose-Hubbard model [5] and also in FIG. 4. From a Mott to a Peierls insulator with extended eigenstates in a trap.One-particle density matrix, |⟨ ψf (x) ψ † f (x ′ )⟩| 2 , where only the top-right quadrant is plotted owing to the symmetry about x = 0. Two cases are shown: (a) a trapped mixture in the wedding cake phase, for the same state as the one presented in the left column of Fig. 3; (b) Peierls phase for the same state as in the central column of Fig. 3.The non-diagonal part of the one-particle density matrix is significantly suppressed in the wedding cake phase, where the decay is Gaussian and on a length scale similar to the characteristic length of the optical lattice.In the Peierls phase, there are instead periodic oscillations signifying coherence.
fermionic systems [84,85].In this state, the fermions form plateaus (or "layers of the cake") with an integer number of fermions per site.In the figure, a plateau of three fermions per site is observed in the center of the trap, transitioning in a step-like fashion to two and finally one fermion per site as we move further out from the center.The wedding cake pattern results from the interplay between kinetic, potential, and interaction energies and only appears when an optical lattice is present.It cannot be accurately modeled using a local density approximation due to the abrupt changes in density.
In the central column of Fig. 3, we elaborate on the correspondence to the Peierls phase occurring in the absence of a trap.Generally, the boson and fermion clouds form a region of a size determined by the trapping strength alongside the BB and BF interactions.The optical lattice induces a spatial variation in the bosonic profile, imparted to the fermions via the BF interaction.The Friedel oscillations of the fermions in this potential are then imprinted on the bosons via the RKKY interaction.Because the Peierls phase is a compressible supersolid, the wavelength varies continuously over the width of the trap.Therefore, in the trapped system, the Peierls modulations are generally not commensurate with the lattice.These modulations may be detected using single-site addressing, leading to a spatial variation in density that cannot be attributed to the profile of the trap or that of the optical lattice.We note here that the beating between the fermion and the optical lattice wavenumbers may give rise to local commensurate-incommensurate transitions in which the fermionic wavenumber of an atomic density wave spatially adjusts to the optical lattice period [72].
For any positive value of g bf and any strength of the optical lattice, we have observed that the system phases separate.The rightmost column of Fig. 3 displays the separated state, where the heavier atoms occupy the center of the trap, surrounded by regions of the lighter species.Perhaps the most surprising feature of Fig. 3 is the emergence of the wedding cake phase where we would expect a collapsed phase.The integer filling per site clearly indicates a fermionic Mott insulator state [86][87][88].One way to ascertain whether this instance is something more than a coincidence is to examine the off-diagonal elements of the one-particle density matrix, given by To derive the above expression, we have written where ψf,n are fermionic operators and the orbitals ϕ f,n (x) are solutions of the Hartree equation (see Appendix B) An insulating state is expected to have an exponentially decaying one-particle density matrix [89,90].However, it turns out that this is not the case for the BF mixture in a harmonic trap.
In Fig. 4, we focus into a more detailed comparison between the wedding cake [frame (a)] and Peierls [frame (b)] phases, showcasing the corresponding one-body density matrix.As seen in Fig. 3(a), the wedding cake is distinguished by having an integer number of fermions per site.The one-particle density matrix provides an additional way of distinguishing the phases.In the wedding-cake phase, it displays a Gaussian decay on a length scale similar to the characteristic length of the optical lattice.This feature, in turn, signals the opening of a more pronounced gap in the spectrum and the presence of an insulating phase.On the other hand, the spatial decay of the matrix elements in the Peierls phase reveals a sustained coherence originating from the individual trapped-particle wavefunctions, which we will attempt to elucidate.
A characteristic example of the spatially decaying one-point density matrix n(x, x ′ ) is plotted in Fig. 5 for It is evident that in the normal phase, the "wavelength" of the fitted sinc function scales with √ k, while the relation is not so direct in the Peierls phase.To motivate the mechanism underlying the periodically decaying oscillations, we start by setting g bf = 0, in which case the fermions are free from their interaction with the bosons in the same trap.In this case, we find (see Appendix C for further details) Equation ( 24) explicitly reveals an algebraically decaying envelope and a N f -dependence of the spatial frequency, with N f being the number of fermions confined in the trap, both attributes evidenced in the left plots of Fig. 5.While gapped phases might suggest an exponential decay of the single-particle density [91], the spatial confinement imposed by the trap implies that the extent of excitations is typically of the same order as the characteristic trap length (unless very weak trap frequencies are considered).Thus, we observe an algebraic decay instead of an exponential one.Turning on an attractive Bose-Fermi interaction g bf to access the regime of Peierls instability in the right frames of Fig. 5 marks the departure from the sinc profile and the N f scaling, yet the periodicity is still evident, and the frequency carries on increasing with the extent of the trap.The counterdiagonal now defines a particular cut through ripples visible in Fig. 4(d) -remnants of quantum interference developing symmetrically to the main diagonal and signifying a Peierls phase with a larger spatial extent than the normal.The overall amplitude of the density matrix along the counterdiagonal, depicted in the right frames of the figure, also increases with N f , in line with the prediction of Eq. (C2).Taking N f → ∞, we obtain the delta function as the familiar limit of the sinc function, writing n(x, −x) ∝ δ(x).

V. CONCLUDING REMARKS
In this work, we have studied the response of a dilute mixture of two atomic clouds of disparate masses, one fermionic and one bosonic, to an external lattice potential alongside an atomic trap.The gases are confined in a cylindrical trap so that the low-energy physical behavior of the system is effectively one-dimensional, while the optical lattice imprints a period structure onto the bosonic gas.We have assumed a weak optical lattice and a high density of bosons, such that the bosons form a condensate, which amounts to an effective interaction of the fermions with a coherent field.
We confirmed that a Peierls instability, manifested via a BKT-type phase transition to a Peierls phase, persists when the optical lattice potential is applied solely to the bosons.Moreover, the presence of this phase was found to stabilize the system against collapse and separation, as well as to enhance the significance of the Peierls effect, a property that is highly desirable for experimental explorations of the associated phonon-like physics.The presence of the Peierls phase was confirmed upon developing and applying a hybrid mean-field method, which retained an amount of coherence able to sustain sinctype oscillations along the counterdiagonal of the oneparticle density matrix.It is interesting to note that such a method can reveal the BKT nature of the phase transition, as well as the similarity to the SSH model of the Peierls instability.The crucial scaling is captured by the method we have developed owing to the fact that we take the deformation of the Fermi surface into account, while we also retain the quantum degrees of freedom for the fermions albeit within an adiabatic approximation.An algebraic decay of correlations is found for all phases occurring in the trapped BF mixture, in accordance with a BKT phase transition.
The natural step forwards would be to take quantum fluctuations of the boson degrees of freedom into account.This entails the development of a formulation and related numerical methods manifestly beyond any mean-field reduction, and to that aim we outline possible directions for deriving the corresponding many-body Hamiltonian in Appendix A below.Another promising direction for future research concerns the investigation of fermion fillings incommensurate with the optical lattice as well as studying the quenched time-evolution problem of a single fermion in a dynamical bosonic potential.Furthermore, an interesting extension replaces the optical lattice with a BB mixture, using either trapped ions or dipole bosons.Finally, systems of the kind can be analyzed in greater detail using matrix product states or bosonization techniques.
Appendix A: Routes towards known quantum many-body models In Sec.III A, we reported on the presence of a Peierls phase at a hybrid mean-field level.We argue that this description accounts for the pairing effect qualitatively.However, there are clear limitations to such an approach.Effectively, the model is a free fermion theory, and as long as there is a finite Peierls gap, ∆ ̸ = 0, however we find a nontrivial insulator with algebraically decaying order for the trapped BF mixture.This is in line with what one expects from a BKT transition, where one phase should exhibit algebraic decay.Furthermore, the classical field stemming from the boson condensation cannot build up quantum correlations with the fermions, while a quantum field can.To quantitatively analyze such situations one must go beyond mean-field.In this appendix we outline, without going into details, two possible approaches in order to derive effective Hamiltonians for which both species are treated quantum mechanically.
In Ref. [36], the condensate wavefunction was calculated for zero coupling between the two species, and subsequently, Bogoliubov excitations around the mean-field response were taken into account, arising from a weak interaction.Due to the presence of a deep external lattice potential in their model, a Hubbard-Holstein model is obtained which can be understood in terms of polarons.In our case, we instead find the fermion-phonon Hamiltonian where are given in terms of the u µ (x), v µ (x) of the Bogoliubov transformation and ω µ are the energies of the Bogoliubov modes.The third term will confine the fermions to their lowest Bloch band if the interaction is strong enough, and in this limit, the model approaches the Hubbard-Holstein model of Ref. [36].However, the approach above would have to be modified since it rests on the assumption that g bf is weak.For a deep optical lattice one could come up with another method: expand the boson field in the corresponding Wannier functions, w j (x) localized at site j, as [7,8] ψb (x) = j âj w j (x). (A3) If we expand the fermion field as ψf (x) = k ĉk e ikx and impose the single-band and tight-binding approximations, we arrive at the many-body Hamiltonian where is the Bose-Hubbard Hamiltonian, with J and U the tunneling rate and onsite interaction strengths, respectively, nj = â † j âj , and D(k − k ′ ) is the overlap integral (which is Gaussian in the harmonic approximation).Since Nb is preserved, the quadratic fermionic Hamiltonian can be readily diagonalized, and we find no backaction on the bosons due to the fermions.To incorporate such effects one would need to go beyond the single-band or tight-binding approximations.Alternatively, the approximations may be kept but other additional degrees of freedom should be introduced corresponding to the phonon modes.For heavy bosons, one can follow the idea of Ref. [92] to allow a variation in the position of the Wannier-function centers jπ, but keep the shape of the functions intact.Thus, we associate a quantized shift δj with each Wannier function, w j (x) → w j (x − δj ).Assuming small shifts δx j ≪ 1 one may then expand to linear order in them and derive an effective Fröhlich-like Hamiltonian [13] ĤF r = Ĥmb + ĤBP + ĤF BP , (A6) where the boson-phonon interaction is given by If we introduce the local phonon annihilation/creation operators dj / d † j , the phonon displacement is expressed as δj = dj + d † j /2, (A8) while the fermion-boson-phonon interaction term takes the form where the first term is the bare phonon energy (µ is the characteristic frequency and Nd = j d † h dj ).In principle, within the harmonic approximation, the coefficients J 1 , U 1 and D 1 can be analytically determined.Both ĤBP and ĤF BP describe phonon-assisted tunneling, either between neighboring lattice sites (bosons) or between different momentum modes (fermions).It may be noted that if Appendix B: The hybrid mean-field-quantum method In this appendix, we provide some further details about the employed hybrid mean-field ansatz.As already mentioned in Sec.III A, a gap of the form (7) opens up in the pairing phase.A superfluid of bosons is gapless and, to the lowest order, is described by the condensate wavefunction ψ b (x) -the mean-field order parameter.We define the generalized coherent state where |0⟩ is the bosonic vacuum.Taking the expectation value of the time-dependent Schrödinger equation, we obtain an effective hybrid Hamiltonian where the first term is given by ( 12), the second term is and the third term is the Berry connection [70].We will assume that the state of the condensate is slowly modified, so that this term can be neglected.The above is equivalent to a product-state ansatz for the ground state where |Ψ f ⟩ is a general fermion state of N f fermions.Such an ansatz neglects any entanglement built between the bosons and fermion subsystems; correlations may actually arise within the non-condensed fraction of bosons and, more importantly, different configurations of the condensate may get entangled with the fermions.In general, such mixed quantum-classical dynamics is interesting in a much wider context in condensed matter physics, high-energy physics, and quantum gravity, and has been investigated in detail in Refs.[93][94][95].In Ref. [93], a multi-configurational mean-field approximation based on hybrid quantum-classical theory was developed, with the central object being the quantum-classical distribution function -a map from the classical phase space to the set of quantum density operators, defined as ρ(q, p) = ij ϱ ij |Ψ i ⟩⟨Ψ j |δ(q − q ij )δ(p − p ij ), (B5) where q, p are the classical generalized coordinates and momenta, and |Ψ i ⟩ is a basis state of the quantum subsystem.In this formalism, "non-diagonal" contributions correspond to the quantum subsystem generating a coupling between different trajectories of the classical subsystems.When considering only a single trajectory, the quantum-classical distribution function can be written as a single delta-function term, where the brackets indicate a functional dependence.In other words, the system is represented by a single point in classical phase space, evolving with a quantum density matrix of a pure state.We find the ground state of Ĥeff self-consistently first observing that it can be readily diagonalized by selecting ψf (x) = n ϕ f,n (x) ψf,n , (B7) where the orbitals are solutions of the Hartree (singleparticle) equation Second, we define a non-linear Schrödinger equation for the condensate, where depends explicitly on the bosons through (B8).The ground state of (B9) is found from imaginary time propagation utilizing the split-operator method [71].
As mentioned in Sec.III A, the above method can also be easily adapted to study a translationally invariant system of infinite extent for fermion filling ν f = 1/2.In this case, we assume a periodic solution ψ b (x + 4π) = ψ b (x).The solutions of the Hartree equations are then Bloch waves e iqx/2 u f,n (q, x).Due to the assumed double periodicity of 4π, the first Brillouin zone is halved and the fermionic filling is in direct correspondence to the edges of the Brillouin zone.It is straightforward to find the fermion density as the integral n f (x) = 1BZ dq |u f,0 (q, x)| 2 . (B11) Another important difference here is that we have to diagonalize the Bloch/Hartree equation for each quasimomentum -in practice for a large number of sampled values.Apart from these discrepancies between the two methods, the solution to Eq. (B9) is also found using imaginary time propagation.
Appendix C: Decay of correlations in the absence of Bose-Fermi interaction For g bf = 0 (normal phase), the one-point density matrix can be calculated analytically via the Cristoffel-Darboux formula, involving summed products of individual orthogonal fermionic wavefunctions in terms of the Hermite polynomials H n (x) (with n = 0, 1, . . .N f ) derived from the quantum harmonic oscillator eigenstates.Letting N f = N for convenience, the result reads H N (x ′ )H N +1 (x) − H N (x)H N +1 (x ′ ) x − x ′ .

(C1)
For x ′ = −x, along the counterdiagonal of the density matrix that we select for Fig. 5, the above expression can be recast into the following form: C2) which follows from the familiar asymptotic N ≫ 1 expansion of quantum harmonic oscillator eigenstates [96] e −x 2 /2 H N (x) ∼ (C3) and the Legendre duplication formula for the gamma function.In fact, the number of fermions confined in the trap need not be large.Even for N f ≳ 10 we note a reasonable agreement between the exact formula of Eq. (C1), and the asymptotic sinc profile of (C2).

FIG. 2 .
FIG. 2. (color online)Phase properties in the absence of a trap.(a) Phase diagram in the (g bf , g b ) plane, marking the different phases: LL, Peierls, collapsed, and separated.The solid black lines mark the boundaries of first-order transitions, while the solid red lines represent the BKT transitions.The mean-field prediction (8) is plotted with dotted black lines.The critical coupling g c bf for the LL-to-Peierls transition is calculated from a least square fit of the Peierls gap to the expression (7), using the hybrid mean-field method.(b) Same as for (a) but in the (g bf , V0) plane.Here it becomes clear that the presence of the optical lattice alters the phase diagram quantitatively, in particular through stabilizing the Peierls phase.(c) The Peierls gap as a function of g bf , once more calculated using the hybrid code, for V0 = 0 (yellow), V0 = 4 (green), and V0 = 16 (blue), fitted to the BKT formula (7) (dotted lines).(d) Boson and fermion densities in the Peierls phase for the parameters V0 = 1, g b = 0.4, and g bf = −1.In all plots, we consider half filling, ν f = 1/2.

FIG. 3 .
FIG. 3. (Color online) Density profiles in the presence of a harmonic trap.The top row shows the boson density (magenta) and the fermion density (green), while the bottom row shows the fermion density per site.Only x > 0 is shown since the densities are symmetric around x = 0.Both plots pertain to the ground state in a BF mixture in a trap, obtained using the hybrid method.(a): "Wedding-cake" structured fermionic insulator, with parameters V0 = 4, g bf = −48, ω b = 0.03 and ω f = 0.3.(b): Peierls phase, with parameters V0 = 1, g bf = −4, ω b = 0.072 and ω f = 0.0072.The inset in the upper plot zooms on the density profile in the region 0 ≤ x/a ≤ 8, evincing density modulations in the bulk, with their origin in the Peierls instability.(c): Separated phase with parameters V0 = 4, g bf = 4, ω b = 0.08 and ω f = 0.8.In all simulations, we have set the remaining parameters to g b = 0.4, r = 0.04, N f = 32, and N b = 400.

FIG. 5 .
FIG. 5. (Color online) Scaling of the two-point density matrix with the width of the harmonic trap.The one-point density matrix n(x, x ′ ) is calculated for x ′ = −x and plotted with x scaled with a √ k.The oscillations of the one-point density matrix are fitted to the expression n(x, −x) = Csinc(bx) in line with Eq. (C2).Left column: normal phase (g bf = 0, ω 0,b = 0.005, ω 0,f = 0.07).Right column: Peierls phase (g bf = −1.5, ω 0,b = ω 0,f = 0.0036).The harmonic confinements have been set as ω b/f = kω 0,b/f , with k = 1, 2, 4 for the top, middle and bottom rows.In all simulations, we have fixed the remaining parameters to g b = 0.1 and V0 = 0.25.