Weak-coupling to unitarity crossover in Bose-Fermi mixtures: Mixing-demixing and spontaneous symmetry breaking in trapped systems

TL;DR

This paper develops a non-perturbative theoretical framework for Bose-Fermi mixtures that accurately describes the crossover from weak to strong interactions, revealing new phenomena like spontaneous symmetry breaking and demixing not seen in traditional models.

Contribution

It introduces analytic beyond-mean-field Lagrangians valid across the interaction crossover, improving upon the LHY correction approach and capturing unitarity effects in Bose-Fermi mixtures.

Findings

Identification of spontaneous symmetry breaking in the mixture.

Observation of demixing phenomena.

Qualitative differences from traditional models in density distributions.

Abstract

The usual treatment of a Bose-Fermi mixture relies on weak-coupling Gross-Pitaevskii (GP) and density-functional (DF) Lagrangians, often including the more realistic {perturbative} Lee-Huang-Yang (LHY) corrections. We suggest analytic {non-perturbative} beyond-mean-field Bose and Fermi Lagrangians valid along the crossover from weak- to strong-coupling limits of intra-species interactions consistent with the LHY corrections and the strong-coupling (unitarity) limit for small and large scattering lengths $|a|$, respectively, and use these to study the Bose-Fermi mixture. We study numerically mixing-demixing and spontaneous symmetry breaking in Bose-Fermi mixtures in spherically-symmetric and quasi-one-dimensional traps while the intra-species Bose and Fermi interactions are varied from weak-coupling to strong-coupling limits.The LHY correction is appropriate for medium to weak atomic interactions and diverges for stronger interactions (large scattering length $|a|$), whereas the present beyond-mean-field Lagrangian is finite in the unitarity limit ($|a|\to \infty$). We illustrate our results using the Bose-Fermi $^7$Li-$^6$Li mixture under a spherically-symmetric and a quasi-one-dimensional trap. The results obtained with the present model for density distribution of the Bose-Fermi mixture along the crossover could be qualitatively different from the usual GP-DF Lagrangian with or without LHY corrections. Specifically, we identified spontaneous symmetry breaking and demixing in the present model not found in the usual model with the same values of the parameters.

Figures

Equations

$$\frac{\partial n_i}{\partial t}$$

$$m_i\frac{\partial {\mathbf{v}}_i}{\partial t}$$

$$\displaystyle{\mbox{i}}\hbar\frac{\partial\phi_{i}}{\partial t}=\Big{(}-\frac{\hbar^{2}}{2m_{i}}\nabla^{2}+V_{i}+\mu_{i}\Big{)}\phi_{i},$$

$$\displaystyle{\mbox{i}}\hbar\frac{\partial\phi_{F}}{\partial t}=\Big{(}-\frac{\hbar^{2}}{8m_{F}}\nabla^{2}+V_{F}+\mu_{F}\Big{)}\phi_{F}.$$

$$\underset{|a_i|\to \infty }{\lim }{\mu }_i(n_i,a_i)=\frac{{\hslash }^2}{m_i}{\eta }_i{n}_i^{2/3},$$

$${\mu }_B(n_B,a_B)=\frac{{\hslash }^2}{m_B}(4\pi n_B{a}_B+2\pi \alpha n_B^{3/2}{a}_B^{5/2}+...),$$

$${\mu }_B(n_B,a_B)$$

$$f(x)$$

$${\mathcal{E}}_B(n_B,a_B)$$

$${\mathcal{E}}_B(n_B,a_B)=\frac{{\hslash }^2}{m_B}\frac{3}{5}{\eta }_B{n}_B^{2/3}.$$

$${\mathcal{E}}_B(n_B,a_B)$$

$${\mathcal{E}}_F(n_F,a_F)$$

$$c_1$$

$$\underset{|a_F|\to \infty }{\lim }{\mathcal{E}}_F(n_F,a_F)=\frac{3}{5}{E}_F{\eta }_F.$$

$${\mathcal{E}}_F(n_F,a_F)=\frac{3}{5}{E}_F\left[1+\frac{c_{1}y+(c_2-2c_1^2)y^2}{1-2c_{1}y+\frac{(c_2-2c_1^2)y^2}{{\eta }_F-1}}\right].$$

$${\mu }_F(n_F,a_F)$$

$$d_1$$

$${\mu }_F(n_F,a_F)\equiv$$

$$g(y)=$$

$$\underset{|a_F|\to \infty }{\lim }{\mu }_F(n_F,a_F)=E_F.$$

$$\mathcal{L}$$

$$+\frac{N_B{\hslash }^2}{2m_B}|\nabla {\varphi }_B{|}^2+\frac{N_F{\hslash }^2}{8m_F}|\nabla {\varphi }_F{|}^2$$

$$+\frac{1}{2}4\pi a_{BF}{N}_B{N}_F\frac{{\hslash }^2}{m_R}|{\varphi }_B{|}^2|{\varphi }_F{|}^2,$$

$$\displaystyle\,{\mbox{i}}\hbar\frac{\partial\phi_{B}({\bf r},t)}{\partial t}={\Big{[}}-\frac{\hbar^{2}\nabla^{2}}{2m_{B}}+m_{B}V({\bf r})+\mu_{B}(n_{B},a_{B})$$

$$\displaystyle\,+\frac{2\hbar^{2}\pi a_{BF}N_{F}}{m_{R}}|\phi_{F}|^{2}{\Big{]}}\phi_{B}({\bf r},t),$$

$$\displaystyle\,{\mbox{i}}\hbar\frac{\partial\phi_{F}({\bf r},t)}{\partial t}={\Big{[}}-\frac{\hbar^{2}\nabla^{2}}{8m_{F}}+m_{F}V({\bf r})+\mu_{F}(n_{F},a_{F})$$

$$\displaystyle\,+\frac{2\hbar^{2}\pi a_{BF}N_{B}}{m_{R}}|\phi_{B}|^{2}{\Big{]}}\phi_{F}({\bf r},t),$$

$$V(\mathbf{r})=\frac{1}{2}\left[{\omega }_x^2(x^2+y^2)+{\omega }_z^2{z}^2\right]$$

$$\text{mbox}i\frac{\partial {\varphi }_B(\mathbf{r},t)}{\partial t}=$$

$$\displaystyle+\frac{2\pi m_{B}a_{BF}n_{F}}{m_{R}}{\Big{]}}\phi_{B}({\bf r},t),$$

$$\text{mbox}i\frac{\partial {\varphi }_F(\mathbf{r},t)}{\partial t}=$$

$$\displaystyle+\frac{2\pi m_{B}a_{BF}n_{B}}{m_{R}}{\Big{]}}\phi_{F}({\bf r},t),$$

$$\mathcal{V}(\mathbf{r})$$

$$f(x)$$

Full text

I Introduction

Soon after the observation of Bose-Einstein condensates (BEC) in ultra-cold, ultra-dilute, harmonically trapped alkali atom vapors BEC , several groups were able to create and study trapped super-fluid Fermi gas FERMI and Bose-Fermi mixture bfmixture in a laboratory. Of these, a study of trapped super-fluid Bose-Fermi mixture is of interest because of a rich variety of phenomena it can exhibit rmpb ; rmpf . Such study can provide information about different intra- and inter-species interactions acting in this mixture. Phase-separation $-$ a typical feature of binary super-fluids $-$ in mixtures of quantum degenerate gases has been investigated in Bose-Fermi systems Molmer ; as . With the advent of experimental techniques, now it is possible to change the different inter- and intra-species interactions in a super-fluid Bose-Fermi mixture by manipulating an external electromagnetic interaction near a Feshbach resonance fesh . Hence, it is of natural interest to see how the different mixed and demixed phases of a super-fluid Bose-Fermi mixture change as the inter- and intra-species interactions are varied. It is also of interest to see if such phases could spontaneously break the symmetry of the underlying Lagrangian.

In the present paper, we study the mixing-demixing transition in Bose-Fermi super-fluid mixtures in three-dimensional isotropic and quasi-one-dimensional (quasi-1D) harmonic traps as the intra-species Bose and Fermi interactions are increased from weak to strong coupling. The hyperfine spin-1/2 Fermi super-fluid is considered to be in a fully paired (spin-up-down) state as ; Barbut rather than in a spin-polarized single hyperfine state Molmer ; bfmixture . We consider the intra-species interaction for Bose (Fermi) component as repulsive (attractive) which can be changed continuously from the weak-coupling limit to unitarity. The inter-species interaction between the Bose and the Fermi components is considered to be repulsive within the weakly coupling region. The strong-coupling unitary limit, where the gas parameter $x=|a|n^{1/3}\to+\infty$ with $a$ the $s$-wave scattering length, and $n$ the density, has recently drawn a great-deal of attention as it is characterized by universal laws arising from scale invariance. This limit is of great interest in different areas, such as, Bose and Fermi super-fluids univ ; castin ; salomon ; salomon2 ; hulet ; jinF , superconductivity string , string theory asoke , neutron star and Bose bstar stars, and quark-gluon plasma quark , and can be achieved in a laboratory fesh in a Bose-Fermi super-fluid mixture.

The usual mean-field treatment of super-fluid Bose-Fermi mixture as is confined to the weak-coupling limit of Bose and Fermi interactions described by the Gross-Pitaevskii (GP) gross ; rmpb Lagrangian for bosons and the density-functional (DF) rmpf Lagrangian for fermions. As the interaction strength is increased, we need a non-perturbative beyond-mean-field description. Lee, Huang and Yang (LHY) provided a perturbative Lagrangian for bosons LHYb and fermions LHYf . Although the LHY Lagrangian is valid for slightly stronger interactions, it is not appropriate for very strong interactions in the unitarity limit, where it diverges. For this investigation, we proposed minimal analytic forms of Bose and Fermi Lagrangians valid from weak-coupling to unitarity with proper LHY LHYb ; LHYf and unitarity limits without fitting parameter(s). Most of previous suggestions for such crossover functions for both bosons bose ; as and fermions fermi ; as were numerical with fitting parameter(s) and hence did not have the analytic LHY limits. With the present analytic weak-coupling to unitarity crossover functions, we write the dynamical beyond-mean-field equations for the Bose-Fermi system, which we use in this study of Bose-Fermi mixture in spherically-symmetric and quasi-1D traps. For a quasi-1D confinement, we use this beyond-mean-field 3D model rather than a strict 1D model obtained from a quantum mechanical many-body 1D Hamiltonian. Such a strict 1D model has novel properties like fermionization of bosons TG . The nonlinearities of the strict 1D model are also different from the present 3D model. For a large finite transverse trap in the present quasi-1D case, as in experiments, the present 3D model should be appropriate. How the present results will approximate the results of the strict 1D model under infinitely strong transverse trap is an open question beyond the scope of the present study.

In the present study, we find that the density distribution of the Bose-Fermi mixture along the weak- to strong-coupling crossover could be qualitatively different from that obtained employing the usual GP-DF Lagrangian with or without the LHY corrections. For example, we found demixing in the Bose-Fermi mixture obtained using the present model where the GP-DF Lagrangian predicted mixing of the components. We also found spontaneous symmetry breaking in the present model not found in the GP-DF model.

The plan of the paper is as follows. In Sec. II, we present the analytic expressions for energy and chemical potential of uniform Bose and Fermi super-fluids along the weak coupling to unitarity crossover and derive the nonlinear equations to study the mixing-demixing transition and spontaneous symmetry breaking in Bose-Fermi super-fluid mixtures. In Sec. III, we present the numerical results along the weak- to strong-coupling crossover for the Bose-Fermi mixture under spherically-symmetric and quasi-1D traps. We compared our results with those obtained from the usual weak-coupling GP-DF Lagrangian for the Bose-Fermi mixture and found that the present results could be qualitatively different. A summary of our findings is given in Sec. IV.

II Analytical model along the weak to strong coupling crossover

Bose and fully-paired Fermi super-fluids in the weak-coupling limit are well described by mean-field GP gross and DF DF rmpf equations, respectively. These equations for bosons rmpb and fermions as are equivalent to the super-fluid hydrodynamic equations. The Bose and Fermi super-fluids are described by a macroscopic order parameter. In case of bosons, the order parameter is also the single-particle wave function in the Hartree approximation of the many-body dynamics. In case of fermions, the macroscopic order parameter refers to a fully-paired bosonic entity known as Cooper pair lnc . Hence, the macroscopic hydrodynamic description of a Fermi super-fluid is formulated in terms of paired fermions in the form of Cooper pairs lnc and not in terms of single-particle Fermi wave function as . Such a description has led to excellent results for many collective xxx ; collective phenomena in many-fermion super-fluids, such as density density distribution or frequency of oscillation xxx , where a many-body description becomes unmanageable.

The macroscopic behavior of a Bose rmpb or Fermi rmpf super-fluid is governed by the classical Landau LANDAU equations of irrotational hydrodynamics with the velocity field ${\bf v}_{i}=\hbar\nabla S_{i}/2m_{i}$, where $S_{i}$ is the phase of the order parameter $\phi_{i}({{\bf r},t})=\sqrt{n_{i}({\bf r},t)}e^{{\mbox{i}}S_{i}({\bf r},t)}$ where ${\mbox{i}}=\sqrt{-1}$, $n_{i}$ is the density, $m_{i}$ is the mass of the fundamental entity responsible for super-fluidity: a bosonic atom or a Fermi pair. Here $i=B$ stands for bosons and $i=P$ for paired fermions. The continuity equation and the irrotational flow equation in this case, e. g., LANDAU

[TABLE]

are entirely equivalent to the following dynamical equation for the order parameter $\phi_{i}({{\bf r},t})$

[TABLE]

where $V$ is an external potential and $\mu_{i},$ is the bulk chemical potential for the uniform Bose or Fermi gas. The quantum pressure term $-{\hbar^{2}\nabla^{2}\sqrt{n_{i}}}/{2m_{i}\sqrt{n_{i}}}$ was not present in Eq. (2) in the original classical flow equations but were introduced later for an accurate description of the dynamics. This term leads to the proper kinetic energy term in the Galilean invariant as dynamical equation (3). For bosons Eq. (3) is the GP equation. For fermions we can relate the fermion and pair variables by $m_{P}=2m_{F},V_{P}=2V_{F},\mu_{P}=2\mu_{F}$ and Eq. (3) becomes the following equivalent density functional (DF) equation

[TABLE]

In the strong-coupling regime, the scattering length $a_{i}$ is much larger than all length scales ($|a_{i}|\to\infty$) and consequently, the system shows universal behavior univ ; gior ; castin determined by the density $n_{i}$ independent of the parameter $a_{i}$. Here $i=B$ stands for bosons and $i=F$ for paired fermions. By dimensional arguments, the bulk chemical potential of a uniform Bose or Fermi gas at unitarity is given by as ; rmpf

[TABLE]

where $\eta_{i}$ is a universal parameter and $m_{i}$ the mass of an atom.

In the weak-coupling limit the bulk chemical potential of a uniform Bose gas is given by LHYb

[TABLE]

where $\alpha=64/(3\sqrt{\pi}$), and the first term $4\pi n_{B}a_{B}$ on the right-hand side is the weak-coupling mean-field GP gross result and the second term is the perturbative LHY contribution LHYb , which becomes important for moderate values of the scattering length $a_{B}$ $(>0)$. The LHY contribution to $\mu_{B}(n_{B},a_{B})$ has limited validity as it diverges as $a_{B}\to\infty$ at unitarity, whereas the correct $\mu_{B}(n_{B},a_{B})$ should remain finite at unitarity as given by Eq. (5). The minimal analytic form of the non-perturbative bulk chemical potential consistent with the LHY correction (6) and the unitarity limit (5) and also valid along the crossover from weak to strong coupling is scirep

[TABLE]

Equation (7) with (8) is a Padé approximant to the bulk chemical potential with the proper weak-coupling LHY and unitarity limits.

The LHY bosonic energy density consistent with Eq. (6) is

[TABLE]

The same at unitarity consistent with Eq. (5) is

[TABLE]

These two limiting values can be combined to give the following minimal energy density valid from weak coupling to unitarity

[TABLE]

Although there is no experimental estimate of the parameter $\eta_{B}$ for bosons despite some attempts univ , there are several microscopic many-body calculations of this parameter lying in the range from 3 to 9 otherB ; DG . Of these, Ding and Greene (DG) DG performed a microscopic Hartree calculation along the crossover and in addition to the value of $\eta_{B}=4.7$, they provided a reliable estimate of the universal function $f(x)$ along the crossover. In Fig. 1, we illustrate the present universal function $f(x)$ of Eq. (8) for $\eta_{B}=4.7$ and compare with the same from the microscopic calculation of DG DG and also with the GP functional $f(x)=4\pi x$ and the LHY functional $f(x)=4\pi(x+\alpha x^{5/2}/2)$. For very small $x$ or for small values $a_{B}$, both the GP and LHY functionals are in reasonable agreement with the present crossover functional as can be seen in Fig. 1. However, for larger $x$, near unitarity, the GP and the LHY contribution cannot describe the actual state of affairs.

For a fully-paired uniform super-fluid of spin-1/2 Fermi gas, the energy density is given by rmpf ; LHYf

[TABLE]

with Fermi momentum $k_{F}=(3\pi^{2}n_{F})^{1/3}$, Fermi energy $E_{F}=\hbar^{2}k_{F}^{2}/2m_{F}$, $a_{F}$ $(<0)$ the scattering length of spin up-down fermions. In Eq. (12), the first term $3E_{F}/5$ is the DF term as ; rmpf valid in the Bardeen-Cooper-Schrieffer (BCS) bcs weak-coupling limit. The next two terms represent the perturbative LHY contribution LHYf . The energy density at unitarity is written as xxx ; as

[TABLE]

The minimal analytic energy density along the weak-coupling to unitarity crossover consistent with the weak-coupling LHY limit (12) and the unitarity limit (14) is:

[TABLE]

The following expression for the bulk chemical potential of a uniform Fermi gas in the weak-coupling LHY limit can be obtained from Eq. (12):

[TABLE]

A non-perturbative expression for the chemical potential along the crossover is written in an analogous fashion

[TABLE]

This bulk chemical potential has the correct LHY limit (16) and the unitarity limit:

[TABLE]

In this paper we will use Eq. (4) to describe Fermi dynamics with $\mu_{F}$ of Eq. (18) correct in the weak- and strong-coupling limits employing Fermi variables, and not pair variables. Actually, the chemical potential $\mu_{F}$ and the associated energy for fermions was used by von Weizsäcker von to describe properties of atomic nucleus without any super-fluid properties long before the work of Landau LANDAU . The expression for chemical potential $\mu_{F}$ (18) remains valid for any value of super-fluid fraction. Hence, Eq. (4) can be used to describe many macroscopic properties of equal mixture of spin-up and -down fermions, like density distribution and frequency of oscillation, independent of its super-fluid nature. This is pertinent as for attractive interaction between spin-up and -down fermions the super-fluid fraction could be small 2005 . For small values of super-fluid fraction, the use of Eq. (4) to describe super-fluid properties of attractive fermions, such as the generation of vortex lattice in a rotating Fermi super-fluid, may lead to qualitatively wrong result because the unpaired fermions remain in a normal state and do not contribute to vortex lattice formation. For describing the formation of vortex lattice, a more fundamental set of dynamical equations 2003 should be used. The simple DF equation (4) can, however, be used to describe non-super-fluid properties of the Fermi gas, such as density distribution and phase separation, even for small values of super-fluid fraction. There have been numerous successful applications of similar crossover models to study density distribution density ; xxx and collective dynamics collective of a Fermi gas along the crossover.

There are results of energy density ${\cal E}_{F}(n_{F},a_{F})$ from several theoretical microscopic calculations microF1 ; microF2 ; microF3 ; FNDMC ; LW ; NSR ; QMC ; ETMA and experimental estimates salomon ; prx ; ku . Most of the theoretical and experimental estimates for $\eta_{F}$ lie in the range $\eta_{F}=0.4\pm.05$ microF1 ; salomon ; jinF ; hulet ; salomon2 ; prx ; ku ; FNDMC ; LW ; NSR ; QMC ; ETMA and in this paper we will employ $\eta_{F}=0.4$. In Fig. 2(a), we plot present energy density ${\cal E}_{F}$ of Eq. (15) versus $(k_{F}a_{F})^{-1}$, its LHY limit (12) and the theoretical microF1 ; microF2 ; microF3 estimates of energy density. In 2(b) we compare the present result with several other experimental salomon ; prx ; ku and theoretical FNDMC ; LW ; NSR ; QMC ; ETMA estimates near unitarity. In Fig. 2(c), we plot the chemical potential (18) versus $(k_{F}a_{F})^{-1}$, its LHY limit (12), and the numerically calculated chemical potential from the energy density (15). From Figs. 2, we find that, for large values of $(|a_{F}|k_{F})^{-1}$ in the weak-coupling limit ($|a_{F}|\to 0$), the present crossover results for ${\cal E}_{F}$ and $\mu_{F}$ agree well with the LHY contribution. Along the whole crossover, the agreement of the present ${\cal E}_{F}$ (15) and $\mu_{F}$ (18) with other estimates is good, and we will use these in the present study with the crossover model.

There have been previous attempts to parameterize the chemical potential of a uniform Bose as ; bose and Fermi as ; fermi systems. These previous attempts heavily relied on fitting parameters and/or the agreement with known experimental and theoretical data was poor. The present analytic bulk chemical potentials for uniform bosons (7) and fermions (18) have no fitting parameters and have excellent agreement with known data as illustrated in Figs. 1 and 2, and we will use these in this study.

The Lagrangian density of the localized super-fluid Bose-Fermi mixture is written as as

[TABLE]

where $\phi_{i},N_{i}$ are the order parameter and number of atoms of the Bose or Fermi component, $m_{R}\equiv m_{B}m_{F}/(m_{B}+m_{F})$ is the reduced mass, $a_{BF}$ is the Bose-Fermi scattering length to characterize the spin-independent interaction between a boson and a fermion, and the energies ${{\cal E}_{i}}(n_{i},a_{i})$ are given by Eqs. (11) and (15) and the density $n_{i}\equiv N_{i}|\phi_{i}|^{2}.$

The Euler-Lagrange equations for a spherically symmetric trap corresponding to Lagrangian (II) are

[TABLE]

where $\mu_{B}$ and $\mu_{F}$ are given by Eqs. (7) and (18), respectively, valid along the crossover, and the confining trap is taken as

[TABLE]

where for a spherically-symmetric confinement $\omega_{x}=\omega_{z}$ and for a quasi-1D confinement $\omega_{x}\gg\omega_{z},$ with $\omega_{x}$ and $\omega_{z}$ the trapping frequencies along $x$ and $z$ axes, respectively. In the absence of the interaction between bosons and fermions ($a_{BF}=0$), Eqs. (22) and (23) reduce to Eqs. (3) and (4) for bosons and fermions, respectively.

In Eqs. (22) and (23) the Bose and Fermi chemical potentials $\mu_{B}$ and $\mu_{F}$ are valid along the crossover from weak to strong coupling, whereas for Bose-Fermi interaction we are using its value in the weak-coupling limit as in the GP equation, as there is no universally accepted form of the Bose-Fermi interaction in the weak and strong couplings. This is acceptable if the calculation is limited to only the weak coupling limit of Bose-Fermi interaction. From Fig. 1, we find that the GP functional agrees with the crossover formula for chemical potential for values of the gas parameter $x\lessapprox 0.2$ and the present study will be limited in this domain.

We obtain a set of coupled dimensionless equations from Eqs. (22) and (23) by expressing length in units of $l_{0}\equiv\sqrt{\hbar/m_{B}\omega_{x}}$, time in units of $t_{0}=m_{B}l_{0}^{2}/\hbar$, $|\phi_{i}|^{2}$ in units of $l_{0}^{-3}$, and energy in units of $\hbar^{2}/m_{B}l_{0}^{2}$, etc:

[TABLE]

where $k_{F}=(3\pi^{2}N_{F}|\phi_{F}|^{2})^{1/3}$, and functions $f$ and $g$ are given by Eqs. (8) and (19).

III Numerical Results

Equations (25) and (26) do not have analytic solution and different numerical methods, such as split time-step Crank-Nicolson CN ; CN2 and Fourier pseudo-spectral FS methods, are usually used to obtain their solution. The ground state of the Bose-Fermi mixture is obtained by solving Eqs. (25)-(26) in imaginary time CN .

We consider a Bose-Fermi mixture of 7Li and 6Li atoms in a three-dimensional isotropic trap ($\omega_{x}=\omega_{z}$), where we use the radial symmetry of the system, trapping potential and emergent solutions, to cast the Eqs. (25)-(26) in terms of a single spatial coordinate $r$, i.e. the radial coordinate of the spherical polar coordinate system. The radial spatial step size $\Delta r$ and time step $\Delta t$ used to solve the Eqs. (26)-(25) numerically are 0.05 and 0.0001, respectively, in dimensionless units.

The solutions obtained with the present crossover model which smoothly connects the weak-coupling regime with the unitarity regime are compared with the two models applicable in the weak-coupling regime, i.e. the GP-DF and LHY models. Before we proceed, let us precisely state what we mean by these three models. Solutions of Eqs. (25) and (26) with $f(x=a_{B}N_{B}^{1/3}|\phi_{B}|^{2/3})$ and $g(y=a_{F}k_{F})$ given, respectively, by Eqs. (8) and (19) are termed as solutions obtained by the present model (denoted by symbol $P$ for present). In the LHY and GP-DF models, $f(x)$ and $g(y)$ in Eqs. (25) and (26) are given, respectively, by

[TABLE]

We find that the ground state solution of the present model can be different quantitatively as well as qualitatively from the LHY and GP-DF models.

We find that without the inter-species interaction ($a_{BF}=0$) the density of both components is maximum at the center resulting in a mixed phase. With the increase of repulsive inter-species interaction the density of one of the components may reduce at the center. With further increase of inter-species repulsion the density of one of the components could be zero at the center and when that happens we will call the resultant Bose-Fermi state a demixed state.

We first consider a Bose-Fermi 7Li-6Li mixture with $N_{\rm B}=50000,N_{\rm F}=1000$, $a_{\rm F}=0,a_{\rm BF}=1500a_{0}$ for different $a_{\rm B}$. In Figs. 3 we plot component densities $|\phi_{i}(r)|^{2}$ for the three models normalized as $4\pi\int r^{2}dr|\phi_{i}(r)|^{2}=1$. For $a_{\rm B}=500a_{0}$, there is a good agreement between the component densities obtained from the three models as is shown in Fig. 3(a), which indicates that the system is in the weak-coupling regime. As $a_{\rm B}$ is increased progressively to (b) $1400a_{0}$ and (c) $2000a_{0}$, in both the LHY and GP-DF models, the system slowly changes from demixed state to mixed state with the transition first occurring for the LHY model, viz. Figs. 3 (a)-(b), and then for the GP-DF model, viz. Figs. 3(b)-(c); whereas in the present model the system remains always demixed. On further increase in $a_{\rm B}$, the system remains demixed as per the present model as shown in Fig. 4(a) for $a_{\rm B}=\infty$. However, a large $a_{BF}$ ($=1500a_{0}$) is necessary for demixing in the present model and for a slightly smaller $a_{BF}$ ($=1000a_{0}$), we have mixing in the present model, viz. Fig. 4(b). Hence the qualitative difference among the results of the present model on one hand and the LHY and the GP-DF models on the other hand as found in Figs. 3, as $a_{B}$ is changed from weak to strong coupling, is caused by a relatively large value of $a_{BF}$ ($=1500a_{0}$) used in numerical simulation. If we used the value $a_{BF}=1000a_{0}$ instead, keeping all other parameters unchanged, we verified that there will not be any qualitative difference in the results of the three models: as $a_{B}$ is increased from weak to strong coupling in all three models there will be transition from demixed to mixed configuration (figure not presented). The same will be true for the physical value $a_{BF}=40a_{0}$ for the Bose-Fermi scattering length in the 7Li-6Li system expbf .

Here we are using a large value of $a_{BF}$ in the range $1000a_{0}\sim 1500a_{0}$, and are using a mean-field GP-type Bose-Fermi repulsion valid for small values of the gas parameter $x_{BF}\equiv a_{BF}(n_{BF})^{1/3}\lessapprox 0.2$, where $n_{BF}$ is the geometric mean of Bose and Fermi densities. From Figs. 3 and 4, we find that the typical average values of $|\phi_{i}|^{2}$ are less than 0.0005 and the densities are $n_{i}=N_{i}|\phi_{i}|^{2}$ with $N_{B}=50000$, $N_{F}=1000$ and $a_{BF}=1500a_{0}$; consequently, $x_{BF}\approx 0.12<0.2$, where the GP approximation for Bose-Fermi interaction is valid. In these figures the bosonic gas parameter with $N_{B}=50000$ and $a_{B}=1500a_{0}$ is larger than 0.2, thus requiring the present crossover formula for a proper description of dynamics.

Keeping $N_{\rm F}$ and $a_{\rm F}$ fixed, the $N_{\rm B}$-$a_{\rm B}$ phase plots showing the parameter domains of mixed and demixed states are illustrated in Figs. 5(a) and (b) for $a_{BF}=1500a_{0}$ and $1000a_{0}$, respectively. The demixed (mixed) states for the present, GP-DF, and LHY models appear on the left (right) side of the respective lines.

Let us consider another case with $N_{\rm B}=1000,N_{\rm F}=50000$, $a_{\rm B}=500a_{0},a_{\rm BF}=1500a_{0}$, while $a_{\rm F}$ is progressively decreased. In this case, for $-2000a_{0}\leq a_{F}<0$, the three models lead to qualitatively similar demixing with the Bose component forming a core with the Fermi component forming a shell surrounding it as shown in Fig. 6(a) for $a_{F}=-2000a_{0}$. The results of the GP-DF model are independent of $a_{\rm F}$, and hence in Fig. 6 we do not show the results of this model. As $a_{\rm F}$ is decreased slightly to $a_{F}=-2100a_{0}$, we find a Bose component forming a shell outside a Fermi core at the center as per the present model as is shown in Fig. 6(b) for $a_{F}=-2100a_{0}$. Upon further reduction of $a_{\rm F}$ to $-8000a_{0}$, the Bose shell moves further away from the center of the trap according to the present model, whereas LHY model ends up in the mixed phase as is shown in Fig. 6(c).

If we perform the analysis illustrated in Fig. 6 for $a_{BF}=1000a_{0}$, the inter-species repulsion is weaker and in the weak-coupling limit ($a_{F}=0$), we have overlapping states for all models as shown in Fig. 7(a). Again we do not show the results of the GP-DF model as the same do not change with $a_{F}$. In the strong-coupling limit ($a_{F}=-10000a_{0}$), the present model results in the demixed state, whereas the LHY model remains in the mixed state, viz. Fig. 7(b).

The mixing-demixing phenomena illustrated in Figs. 6 and 7 leads to the $N_{F}$-$|a_{F}|$ phase plots of Figs. 8. In Fig. 8(a) the fixed parameters are $N_{B}=1000,a_{B}=500a_{0}$ and $a_{BF}=1000a_{0}$. In this case mixing-demixing transition takes place only for the present model for $N_{F}\geq 270.$ In Fig. 8(b), with fixed parameters $N_{B}=1000,a_{B}=500a_{0},a_{BF}=1500a_{0}$, the present model, which results in a demixed ground state, predicts a position-swapping transition with the Bose component forming the core which is surrounded by the fermionic shell below a critical $N_{F}$ for a given value of $a_{F}$ (as is shown in Fig. 6(a)). Above this critical $N_{F}$, the Fermi component moves to the core with the Bose component forming the shell around it as in Fig. 6(c). These two qualitatively different states are marked as B(in) and B(out), respectively, in Fig. 8(b). For the same fixed parameters $N_{B}=1000,a_{B}=500a_{0},a_{BF}=1500a_{0}$, the LHY model leads to the demixing-mixing transition as one increases $N_{F}$ at the fixed value of $a_{F}$ as is shown in Fig. 8(c).

Next we consider a quasi-1D 7Li-6Li Bose-Fermi mixture Q1D along the $z$ axis with strong traps in the $x-y$ plane ($\omega_{x}=10\omega_{z}$) and consider a few illustrative examples to show the qualitatively different ground state solutions obtained from the present model as compared to the GP-DF and LHY models. The spatial step size $\Delta x=\Delta y=\Delta z$ and time step $\Delta t$ used to solve Eqs. (25)-(26) numerically are $0.15$ and $0.005625$, respectively, in dimensionless units. In such a quasi-1D Bose-Fermi mixture, the essential collective dynamics and mixing-demixing transition take place in the $z$ direction. Hence for quasi-1D 7Li-6Li Bose-Fermi mixtures, we will illustrate only the reduced 1D density along the $z$ direction $d_{i}(z)=\int dxdy|\phi_{i}(x,y,z)|^{2}$. Apart from mixing-demixing transition in the $z$ direction, we also find spontaneous symmetry-broken states in the quasi-1D 7Li-6Li Bose-Fermi mixtures. Two examples of spontaneous symmetry breaking in the present model are shown in Fig. 9(a)-(b) for the reduced 1D density along $z$ direction where one of the components moved away from the center breaking the parity symmetry. In this case the results of the GP-DF and LHY models lie very close to each other and the result of only the former model is shown. In Figs. 9 there is a demixing in the present crossover model with the density of one of the components being zero at the center while the other component having a density maximum at the center. The densities of the LHY model remain overlapping and parity symmetric.

In addition to the spontaneous symmetry broken states of Fig. 9, it is also possible to have partially to fully demixed states in quasi-1D Bose-Fermi mixtures. For this purpose, we consider (a) $\omega_{x}=10\omega_{z}$ and (b) $\omega_{x}=100\omega_{z}$ , where the spatial step size $\Delta x=\Delta y=\Delta z$ and time step $\Delta t$ used to solve the equations (25)-(26) numerically are taken to be $0.1$ and $0.0025$, respectively. For both (a) and (b), we consider $N_{\rm B}=1000$, $N_{\rm F}=100,a_{B}=500a_{0},a_{\rm F}=-20000a_{0}$, and $a_{\rm BF}=1500a_{0}$. The reduced 1D density $d_{i}(z)$ obtained in these two cases with the present and GP-DF models are shown in Figs. 10(a) and (b); it is evident that the bosons stay in the central region and the fermions are expelled symmetrically in two directions. With the present model, the ground state is partially demixed with $\omega_{x}=10\omega_{z}$, viz. Fig. 10(a) and fully demixed with $\omega_{x}=100\omega_{z}$, viz. Fig. 10(b). In the GP-DF models, there is partial demixing which becomes more pronounced as the trap becomes more confined along the radial direction; the same is the case with LHY model (not shown in the figure).

IV SUMMARY

Here we have demonstrated spontaneous symmetry breaking and mixing-demixing transition in trapped super-fluid Bose-Fermi mixture along the crossover from weak-coupling to unitarity for both intra-species Bose and Fermi interactions. For Bose-Fermi inter-species interaction, we have used the weak-coupling GP interaction. The usual description of the super-fluid Bose-Fermi mixture employs the weak-coupling GP Lagrangian for the bosons and DF Lagrangian for the fermions. The interaction term in the GP Lagrangian is essentially the same as the many-body Hartree interaction term; that in the DF Lagrangian is the total kinetic energy of the fermions in the Fermi sea. Usual treatment of the Bose-Fermi mixture employing GP Lagrangian for bosons and DF Lagrangian for fermions is termed GP-DF formulation. To study the Bose-Fermi mixture along the weak-to-strong coupling crossover, we suggested analytic non-perturbative Lagrangians for intra-species Bose and Fermi interactions with correct LHY limit(s) in the weak-coupling domain and with correct unitarity limit(s) in the strong-coupling domain. These analytic Lagrangians have a single parameter: the universal parameter(s) $\eta_{i}$ at unitarity for bosons and fermions. In this study we have used the most precise value(s) of this universal parameter, which should be updated in future applications, if possible. The Euler-Lagrange equations of the Bose-Fermi mixture describe the dynamics and have been used in this study.

In this study, we considered spherically-symmetric and quasi-1D traps. In both cases, we compared the results of the present model valid along the crossover with those of the usual weak-coupling GP-DF and LHY models and found that the two types of treatments may lead to qualitatively different results. For example, in spherically symmetric traps we identified cases of demixing in Bose-Fermi mixture using the present crossover model not found in the weak-coupling GP-DF and LHY models. In quasi-1D traps, we found spontaneous symmetry breaking in the present crossover model not found in GP-DF and LHY models. Hence we conclude that the perturbative LHY model is unable to properly describe the Bose-Fermi mixture away from the weak-coupling limit.

Acknowledgements.

SG thanks the Science & Engineering Research Board, Department of Science and Technology, Government of India (Project: ECR/2017/001436) and ISIRD of Indian Institute of Technology, Ropar (Project: 9-256/2016/IITRPR/823) for support. SKA thanks the Fundação de Amparo à Pesquisa do Estado de São Paulo (Brazil) (Project: 2016/01343-7) and the Conselho Nacional de Desenvolvimento Científico e Tecnológico (Brazil) (Project: 303280/2014-0) for support.