Dark-Bright Solitons in a Superfluid Bose-Fermi Mixture

The recent experimental realization of Bose-Fermi superfluid mixtures of dilute ultracold atomic gases has opened new perspectives in the study of quantum many-body systems. Depending on the values of the scattering lengths and the amount of bosons and fermions, a uniform Bose-Fermi mixture is predicted to exhibit a fully mixed phase, a fully separated phase or, in addition, a purely fermionic phase coexisting with a mixed phase. The occurrence of this intermediate configuration has interesting consequences when the system is nonuniform. In this work we theoretically investigate the case of solitonic solutions of coupled Bogoliubov-de Gennes and Gross-Pitaevskii equations for the fermionic and bosonic components, respectively. We show that, in the partially separated phase, a dark soliton in Fermi superfluid is accompanied by a broad bosonic component in the soliton, forming a dark-bright soliton which keeps full spatial coherence.


I. INTRODUCTION
A long standing problem in the context of quantum fluids is the description of mixtures composed of two kinds of interacting superfluids belonging to different statistics. The first theoretical analysis of superfluid 4 He mixed with superfluid 3 He dates back to the 70's (see, for example, [1,2] and references therein). In experiments, however, the simultaneous superfluidity of the two components of liquid 3 He-4 He mixtures has never been realized, since the miscibility of 3 He in 4 He is very small (a few percent) and the temperature needed to reach superfluidity of fermions in the mixture is too low to be reached with available cryogenic techniques. Dilute ultracold atomic gases are instead excellent candidates for studying superfluid properties of mixtures. Superfluidity has been recently obtained experimentally in a mixture of a Bose condensed gas and a superfluid Fermi gas of two lithium isotopes, 6 Li and 7 Li [3,4], where a new mechanism for superfluid instability was observed, related to the dynamical instability of the supercurrent counterflow rather than to the more standard Landau criterion [5,6]. In addition, in ultracold atomic gases the strength of the interspecies and intraspecies interaction can be varied by means of an external magnetic field, thanks to the occurrence of Feshbach resonances. One can thus foresee the exploration of the whole phase diagram of the mixture, which is expected to be very rich [7][8][9][10][11][12][13][14][15].
In the present work, we focus on the case where the bosonic superfluid is the minority component, while the fermionic superfluid exhibits a dark soliton. We choose this case as a paradigmatic configuration in which the interplay between miscibility and immiscibility, together with superfluidity, gives rise to a peculiar behaviour * Electronic address: marek.tylutki@ino.it which reveals the crucial effects of nonlinearity caused by interactions. The opposite limit of bright solitons in a mixture with fermions as the minority component was discussed in [16].
In pure Fermi superfluids of dilute atomic gases, theoretical predictions of the structure and dynamics of dark solitons [17][18][19][20][21][22][23][24] have recently stimulated experimental investigations [25]. The experiments confirm the theoretical expectation that dark solitons in a three-dimensional fermionic superfluid quickly decay into vortical excitations due to snaking instability [26], as it was earlier observed with bosons [27,28]. With bosons very long-lived dark solitons have been generated by filling the soliton with atoms in another hyperfine state [29][30][31], thus creating a dark-bright solitonic structure of a two-component Bose-Bose superfluid [32]. Here we theoretically investigate the analogue dark-bright soliton in a Bose-Fermi superfluid mixture, the main difference between the two cases being that the Bose-Fermi phase diagram is known to admit, in addition to a fully mixed phase and a fully separated phase, also a third phase consisting of pure fermions in equilibrium with a mixture of fermions and bosons. The stability conditions of such an intermediate phase in a uniform system were studied in [7] by using the equation of state of an ideal Fermi gas weakly interacting with a dilute Bose gas. Such a phase is predicted to occur also in the strongly interacting regime [11]. A more refined equation of state, including the interaction among fermions, was later applied also to nonuniform configurations by treating the interaction energy in local density approximation [8]. However, since a dark soliton is localized on the length scale of the healing length of the superfluid, which is of the order of the inverse Fermi wave vector for fermions at unitarity, its characterization requires a theory which properly includes nonlocal effects, beyond the local density approximation. For this purpose, we use coupled Bogoliubov-de Gennes and Gross-Pitaevskii equations for the fermionic and bosonic components respectively.
Our paper is organized as follows: in § II we discuss the stability condition of the uniform phase of a Bose-Fermi mixture when the Fermi gas is at unitarity (infinite scattering length); in § III we write the mean-field equations, which are subsequently used to find the stationary solitonic configurations of the system; in § IV we analyse the behaviour of solitons of a Bose-Fermi mixture when fermions are at unitarity. We pay special attention to the transition from the miscible state to the so called partially separated phase and to the fully separated state. We find that, in the partially separated phase, the density depletion of fermions possesses a solitonic character and the phase coherence between the left and right hand sides is maintained; conversely, in the fully separated phase, the depletion in the Fermi density is completely filled by bosons (the density of the fermions vanishes) and the phase coherence between the two sides is lost. Finally in § V we investigate the behaviour of solitons along the crossover from the Bardeen-Cooper-Schrieffer phase to Bose-Einstein condensation (BCS-BEC crossover), and we conclude that, while on the BCS side the system exhibits a behaviour similar to that at unitarity, on the BEC side the partially separated configuration disappears as expected for a Bose-Bose mixture.

II. PHASE SEPARATION
Let us first discuss the conditions for the stability of the homogeneous phase of the mixture, which is a crucial point in order to understand the numerical results of the following sections. The phase diagram of a weakly interacting Bose-Fermi mixture, at zero temperature, can be derived starting from the following expression for the energy density: Here n b and n f are the densities of the Bose and Fermi gas respectively, and E F =h 2 k 2 F /2m f is the Fermi energy, with k F = (3π 2 n f ) 1/3 . In the above expression we assume that the Fermi gas is at unitarity (infinite scattering length) and η is the dimensionless Bertsch parameter [33][34][35], which simply rescales the energy density of the Fermi gas with respect to the ideal gas expression. The quantities g bb and g bf are the bosonic intraspecies and the Bose-Fermi interspecies coupling constants, which are related to the corresponding scattering lengths according to g bb = 4πh 2 a bb /m b and g bf = 4πh 2 a bf (m b + m f )/(2m b m f ), where m f and m b are the masses of fermions and bosons. We assume that the Bose-Fermi scattering length does not depend on the internal state of the Fermi atoms, as in the case of the recent experiments with lithium atoms [3,4]. Equation (1) was first used in [7] and [36] to describe the phase diagram of a dilute Bose gas interacting with an ideal Fermi gas (η = 1). In [7] the phase diagram was explored as a function of the densities of the two components and the existence of three phases was predicted for a positive Bose-Fermi scattering length: i) a uniform mixture, where both components occupy the entire space at constant densities; ii) a partially separated phase, where part of the space is occupied by pure fermions and part by a Bose-Fermi mixture; iii) a fully separated phase, where bosons and fermions are completely separated. The same happens at unitarity, with the only difference that the Fermi energy is renormalized by the universal Bertsch parameter η. It is worth noticing that the existence of three phases is peculiar of the Bose-Fermi mixture. In fact a Bose-Bose mixture only admits the mixed uniform phase and the fully separated phase, because of the different power-law dependence on the densities in the equation of state.
The stability condition predicted by the energy density (1) for the uniform mixture is ( For n f larger than this critical value the uniform mixture is unstable and the system exhibits either partial or full phase separation.

III. MEAN-FIELD EQUATIONS
In order to describe both uniform and nonuniform configurations of the mixture, including solitons, we use the Bogoliubov-de Gennes (BdG) equations [35] for the interacting superfluid fermions and the Gross-Pitaevskii (GP) equation [37] for the Bose-condensed bosons, which are coupled via the interspecies interaction term, fixed by g bf . We also use a simple geometry consisting of a threedimensional rectangular box with hard walls in the longitudinal direction z (box size L) and periodic boundary conditions in the transverse directions x and y (box size L ⊥ ). The system is assumed to be uniform in the transverse directions and all spatial variations are along z.
The Fermi gas is described in terms of a set of quasiparticle amplitudes {u j } and {v j }, solutions of the BdG equations, which in our case take the form: and j is the generic quantum number associated to the longitudinal degrees of freedom. With this choice, the BdG equations become are the order parameter (gap function) and the density of the superfluid fermionic component. The quantity µ f is the chemical potential of fermions, and the eigenvalues ε j,k ⊥ are the quasiparticle energies. The Bose-condensed component is instead described by the Gross-Pitaevskii equation for the bosonic order parameter (macroscopic wave function) ψ b : where µ b is the chemical potential of bosons and the den- The four equations (3a)-(3d) must be solved selfconsistently. An energy cutoff E c = 50E F is used to solve the BdG equations; this implies a proper renormalization of the Fermi-Fermi coupling constant, for which we use the relation 1/(k F a f f ) = 8πE F /(g f f k 3 F ) + 2 E c /E F /π [18]. The key parameter characterizing the interaction among fermions is 1/k F a f f , and we perform calculations in the range −1 ≤ 1/k F a f f ≤ 1, in the crossover from the BCS regime (negative values) to the BEC regime (positive values), passing through unitarity (1/k F a f f = 0).

IV. SOLITONIC SOLUTIONS IN THE UNITARITY REGIME
In the following we will consider a mixture where the number of bosons is about 10% of the number of fermions. For simplicity we also impose m b = m f , which is a reasonable approximation for mixtures of two isotopes of the same atomic species, but this assumption does not affect the main results of the work. We typically solve Eqs (3a)-(3d) with N f ≈ 500 fermions and N b ≈ 50 bosons in a box with L = 70k −1 F and L ⊥ = 15k −1 F . We start from trial functions ∆(z) and ψ b (z) and iterate till convergence to a stationary solution which does not depend on the initial choice.
Examples are given in Fig. 1  values of the parameter g 2 bf n 0 /(g bb E F ). The three values are chosen as representative of the different phases: a fully mixed phase for the smallest value of g 2 bf n 0 /(g bb E F ) (top panel), a fully separated phase for the largest one (bottom panel), and an intermediate partially separated phase (middle panel). In the latter case fermions occupy the whole volume, partly as a pure phase and partly in a mixed phase, and the interface between the two regions is significantly wider than the domain wall found for the fully separated phase (bottom panel), where the width of the interface is of the order of the healing lengths of the two superfluids. The healing length of bosons is much smaller than the size of the box in our case, because we have chosen the value of g bb such that the solution of the GP equation is well approximated by the Thomas-Fermi approximation n b (z) = (µ b − g bf n f (z)) /g bb , except near the box boundaries.
To gain further insight, we perform calculations for several values of g 2 bf n 0 /(g bb E F ) and, in each case, we calculate the overlap integral dz n f (z)n b (z). The results are shown as markers in Fig. 2. The overlap integral is 1 in the fully mixed phase and vanishes in the fully separated phase. The instability of the uniform mixture is clearly visible as a sharp deviation from 1, which occurs at the critical value of g 2 bf n 0 /(g bb E F ) 0.4. In an infinite system (i.e., in the thermodynamic limit, where surface and interface effects are ignored) this value is expected to be well approximated by the linear stability condition (2); by using the mean-field value of the Bertsch parameter, η = 0.59 [35], this threshold is g 2 bf n 0 /(g bb E F ) 0.6 (vertical green dashed line). In the same figure we also show the percentage of space occupied by the pure Fermi gas (thick dashed line) and the pure Bose gas (solid line). One can see that, by increasing g 2 bf n 0 /(g bb E F ) above the critical value 0.4, pure fermions start occupying a significant part of the box while bosons remain still mixed with fermions, which is another indication of the occurrence of the intermediate phase.
Having discussed the ground state as a test case, let us now concentrate on the second column of Fig. 1, where we show the density profiles in the presence of a dark soliton in the Fermi component. The parameters are the same as in the first column. The dark soliton is imprinted in the Fermi superfluid by imposing a node (∆ = 0) of the order parameter at the centre of the box (z = 0) and a π phase difference between the two sides. In a purely fermionic superfluid, the stationary solution of the BdG equation with such constraints exhibits a deep density depletion of width of the order of k −1 F [18]. If bosons are present, depending on the values of the interaction strengths g bf and g bb , they can be attracted into the soliton, thus changing its structure. The figure shows that, for small values of g bf (top panel), the system favours a uniform mixed configuration in agreement with the condition (2) and only a small fraction of bosons is attracted by the soliton, which remains almost unaltered. The resulting structure is analog to the dark-anti-dark soliton pairs that are predicted for Bose-Bose mixtures [38,39]. In the opposite limit of large g bf (bottom panel) all bosons form a pure phase at the centre, pushing fermions apart, and the soliton is replaced by two domain walls separating pure bosons from pure fermions. The central case is the most interesting: bosons like to stay mixed with fermions, but fermions like to form a pure phase. The net effect is that all bosons are pushed into the soliton but with a broad overlap between the two components. The overall structure keeps its solitonic character and becomes a bright-dark soliton. The shaded area in Fig. 2 is the interval of g 2 bf n 0 /(g bb E F ) where we find stationary solutions having such a brightdark soliton structure, with both the node of the order parameter and the minimum of the fermionic density at a single point, z = 0. The resulting scenario shares interesting analogies with the dark-bright solitonic structure exhibited by two-component Bose superfluids [32]. We find similar solutions also for different values of the bosonic fraction: N b = 20%, 30% and 50% of N f .
For the three dark soliton solutions in the second column of Fig. 1 the phase of the order parameter is π for negative z and 0 for positive z, which implies that  Fig. 1, with the same interaction parameters, but with fermions in the BCS regime (1/kF a f f = −1). The structure is qualitatively the same as at unitarity, but the soliton is shallower and the fermionic density exhibits Friedel oscillations. The order parameter, not shown, is also similar to the one in the right column of Fig. 1, but with ∆0 = 0.22EF . ∆(z) is a real function. The third column of the same figure shows the corresponding profiles of ∆(z). However, in the case of complete separation (bottom panel) the phase difference between the two sides is irrelevant, because the two regions of pure fermions, separated by the central Bose gas, behave as independent superfluids with no phase coherence; in fact, we have numerically checked that, by arbitrarily changing the phase difference of the order parameter, both n f and |∆| remain unchanged in the solution of the coupled BdG and GP equations. Moreover, adding more bosons would simply result in further separating the fermionic superfluids, keeping the same structure of the domain walls. This is not the case of the partially mixed phase shown in the central panel where, if we suddenly change the sign of ∆ on the left side, for instance, the solitonic structure is quickly lost and the solution converges to one without soliton, as in the first column. This proves that the dark-bright soliton in the intermediate partially mixed phase has indeed solitonic character, being an effect of nonlinearity and phase coherence.

V. FERMIONS IN THE BCS AND BEC REGIMES
Now we discuss what happens when the fermionic component of the mixture is in the BCS-BEC crossover, away from unitarity.
In the BCS regime (1/k F a f f < 0) fermions can also phase-separate from bosons either partially or completely. This is not surprising, as the BCS and the unitary regimes share many similar features as suggested, for example, by the fact that the chemical potential is positive in both cases. The qualitative picture of phase separation is the same: a uniform mixture is stable for small values of the the parameter g 2 bf n 0 /(g bb E F ), a partial phase separation occurs for larger values, and a complete phase separation for an even larger coupling. As a consequence, also the structure of the solitonic solutions is expected to be similar. In Fig. 3 we show the results for 1/k F a f f = −1, for the same parameters used at unitarity in the central column of Fig. 1. Apart from more pronounced Friedel oscillations and a shallower solitonic depletion in the density distribution of fermions [18], we find that the results look indeed very similar. If one moves further into the BCS regime, however, the solitonic structure becomes broader and broader, and less robust against instability mechanisms associated to the fermionic degrees of freedom [23].
In order to discuss the BEC regime (1/k F a f f > 0) we have to recall that the fermionic superfluid of ultracold atoms in the BCS-BEC crossover is actually made of an equally populated two-spin-component Fermi gas and the s-wave scattering length a f f accounts for the interaction between atoms with different spins. Due to pairing, this two-component gas behaves as superfluid described by an order parameter ∆ and the total density n f , as we have done so far. However, in the BEC regime fermions with opposite spins form tightly bound bosonic dimers, which in turn form a Bose-Einstein condensate. In this situation neither Eq. (1) nor Eq. (2) hold. The Bose-Fermi mixture thus behaves as a Bose-Bose mixture, where one of the two bosonic components is the condensate of dimers, and the system can be described by two coupled GP equations where the mass of the dimer is m d = 2m f and the coupling constants of the dimer-dimer and dimer-boson interactions, in the first Born approximation, are given by g dd = 2g f f andg = 2g bf (note that, though the exact many-body values for these parameters are different [35], we use the mean-field expressions for consistency with the BdG equations). Replacing the BdG equation for the order parameter ∆ of paired fermions with a GP equation for the order parameter ψ d of bosonic dimers is expected to be a good approximation for 1/(k F a f f ) of the order of, or larger than one. An important consequence is that the region of the phase diagram of the mixture where one finds the intermediate, partially mixed phase becomes narrower when moving from unitarity towards the BEC regime [8] and eventually disappears in the deep BEC limit. An example is given in Fig. 4, where we show the density profiles of fermions (blue lines) and bosons (red lines) of the Fermi-Bose mixture, obtained by solving the coupled BdG and GP equations (3a)-(3d) for 1/k F a f f = 1, the other parameters remaining the same as in the central row of Fig. 1. The solution in the top panel corresponds to the case where the phase of the order parameter ∆ is constant, while the one in the bottom panel is obtained by imposing a π phase difference and a node at the centre. Two main comments are in order: i) the comparison with Fig. 1 shows that the partially mixed phase, which was present at unitarity for the same value of g 2 bf n 0 /(g bb E F ), is lost in the BEC regime in favour of a fully separated phase; ii) if we calculate the density profile by solving the coupled GP equations (4) (dashed lines) instead of the the coupled BdG and GP equations (3a)-(3d) (solid lines), we find very similar results.
These considerations suggest that the unitary regime is the most suitable and interesting for the investigation of dark-bright solitons in Fermi-Bose mixtures.

VI. CONCLUSIONS
In this work we theoretically study a mixture of Bose and Fermi superfluids by using a mean-field theory based on the solutions of coupled BdG and GP equations. We tune the interaction between fermions in the BCS-BEC crossover, and we also vary the interspecies and intraspecies interactions in order to explore the three regimes of fully mixed, partially mixed and fully separated phases. The focus of the work is on solitonic solutions, and we find that, in the regime of partial mixing, a dark soliton in the Fermi component becomes wider and deeper, taking all bosons in, but maintaining phase coherence. This dark-bright solitonic structure can serve to stabilize dark solitons in Fermi superfuids against snaking instability, but it is also interesting in itself as an example of many-body state where nonlinearity and coherence play a relevant role. Here we have found that bosons favour and amplify an inhomogeneous (solitonlike) order parameter of fermions. A similar mechanism might also favour the realization of the Fulde-Ferrell-Larkin-Ovchinnikov (FFLO) phase in the unitary regime of the Fermi gas, which involves density modulations and could be enhanced by the presence of bosons with spindependent interaction, thus making the observation of such elusive phase easier.
We finally notice that the range of parameters in which we find partial mixing may be experimentally accessible using, for instance, mixtures of 6 Li-7 Li, 6 Li-87 Rb or 23 Na-87 Rb. When the mass of the two components is different, and particularly when m f /m b 1, the situation is expected to be even more favourable. In fact, by expressing Eq. (2) in terms of masses and scattering lengths, one finds that the condition for the stability of the miscible phase becomes (n