Quantum correlations and spatial localization in one-dimensional ultracold bosonic mixtures

We present the complete phase diagram for one-dimensional binary mixtures of bosonic ultracold atomic gases in a harmonic trap. We obtain exact results with direct numerical diagonalization for small number of atoms, which permits us to quantify quantum many-body correlations. The quantum Monte Carlo method is used to calculate energies and density profiles for larger system sizes. We study the system properties for a wide range of interaction parameters. For the extreme values of these parameters, different correlation limits can be identified, where the correlations are either weak or strong. We investigate in detail how the correlation evolve between the limits. For balanced mixtures in the number of atoms in each species, the transition between the different limits involves sophisticated changes in the one- and two-body correlations. Particularly, we quantify the entanglement between the two components by means of the von Neumann entropy. We show that the limits equally exist when the number of atoms is increased, for balanced mixtures. Also, the changes in the correlations along the transitions among these limits are qualitatively similar. We also show that, for imbalanced mixtures, the same limits with similar transitions exist. Finally, for strongly imbalanced systems, only two limits survive, i.e., a miscible limit and a phase-separated one, resembling those expected with a mean-field approach.


Introduction
The fascinating physics of interpenetrating superfluids has recently become a topic of large interest due to the experimental realisation of multi-component, atomic Bose-Einstein condensates [1,2,3,4,5]. In the weakly interacting regime, these mixtures are well described by coupled mean-field Gross-Pitaevskii equations (GPEs), and within this framework processes that lead to phase separation are well described [6,7,8,9,10,11,12,13,14] While mean-field theories allow to study weakly correlated systems, it is also important and interesting to examine quantum mixtures in strongly correlated regimes. In these regimes, analytic solutions can often only be obtained in limiting cases. Rather appealing results occur in strongly correlated regimes when the dimensionality is reduced. For quasi one-dimensional (1D) gas mixtures one finds that Luttinger liquid theory predicts many interesting effects, which include de-mixing for repulsive interactions or spin-charge separation analogous to that found in 1D electronic quantum systems [15,16,17,18]. Other relevant effects include the presence of polarized ground states, which allow to view the relative spatial oscillations as spin waves [19,20,21,22] and which have been experimentally observed [23,24,25].
Very strong correlations for single component bosons are realized in the Tonks-Girardeau (TG) gas [26,27,28], which was recently observed experimentally [33,34]. Bosonic mixtures in the strongly interacting limit have features common with the TG gas, and their ground-state wavefunction can be obtained analytically in certain interaction limits [35,36,37]. Experimental advances on Feshbach and confined induced resonances in recent years have made it possible to control both, the intra-species interactions and the inter-species interactions, over a wide range of parameters [38,39,40]. In the strongly interacting limit a number of relevant phenomena have been described including phase separation [15,16,17,41], composite fermionization [42,43,44], a sharp crossover between both limits [45], and quantum magnetism [46].
In this work we focus on mixtures where the number of atoms is small. The recent successful experimental trapping of ensembles of few atoms [47,48,49,50] has inspired an intense theoretical effort in few-atom systems [51,52,53,54,55,56,57,58,59,60,61,62,63,64]. For mixtures of few atoms, direct diagonalization methods [30,44,41], can be used together with other numerical methods efficient for larger numbers of atoms, like multiconfigurational Hartree-Fock methods (MCTDH) [66], density functional theory (DFT) [43], or quantum diffusion Monte Carlo (DMC) [65]. In the present work, we use direct numerical diagonalization to study the ground-state properties of a mixture of ultracold bosons confined in a 1D trap over a wide range of correlations regimes, determined by the scattering properties between the atoms. These are supplemented by DMC calculations to confirm trends for systems with larger particle numbers. While the extreme cases in which all correlations are either weak or strong are well known, here we calculate and discuss the full phase diagram and especially the transitions between the different regimes.
We study the ground-state wavefunction, and pay particular attention to the oneand two-body correlations in the extreme limits, and across the transitions between them. The quantum correlations between both components are characterized by means of the von Neumann entropy. This allows us to show that close to the crossover between the composite fermionization and phase separation, the ground state exhibits strong correlations between the two bosonic components.
Our manuscript is organized as follows. In Sec. 2 we introduce the model Hamiltonian and a general analytical ansatz for the ground-state wavefunction. Focusing first on balanced mixtures, we discuss in Sec. 3 the ground-state properties in terms of the densities, the coherence, the energies, the one-and two-body correlations, and the von Neumann entropy. In Sec. 4, we then present results on how the ground-state properties change when one component is larger than the other and finally summarize all our results in Sec. 4.

Model Hamiltonian
Let us consider a mixture of two bosonic components, A and B, with a small, fixed number of atoms in each component, N A and N B . We assume that the two components are two different hyperfine states of the same atomic species of mass m, and that they are trapped in the same, one-dimensional parabolic potential V (x) = 1 2 mω 2 x 2 . At low temperatures, all scattering processes between the atoms are assumed to be described by contact interactions , where the positions of atoms of species A(B) are given by the coordinates x j (y j ). The 1D intra-and inter-species coupling constants g A(B) and g AB are assumed to be tunable independently by means of confinement induced resonances [38]. We will restrict our study to repulsive interactions. The many-body Hamiltonian isĤ =Ĥ A +Ĥ B +Ĥ int , withĤ There are three coupling constants g A , g B , g AB each of them ranging from g = 0 for ideal Bose gas interaction to g → ∞ for strong Tonks-Girardeau interaction. This defines eight limits schematically shown in Fig. 1. The composite fermionization limit is reached when g AB → ∞ with the other coupling constants vanishing [42,43,44]. We termed TG-BEC gas a system with one of the intra-species coupling constants large, while other coupling constants vanish [41]. If one of the intra-species coupling constants together with the inter-species coupling constant are large, the phase separation limit is reached [15,16,17,41]. Finally, if all coupling constants tend to infinity, the wavefunction is known exactly and can be mapped to the one of an ideal Fermi gas [35]. We call this limit full fermionization. In the following we will calculate and discuss the complete phase diagram, which includes the transitions between these limits. To restrict the large number of free parameters, we note the transition between TG and phase separation limit is symmetric when switching the values of g A and g B and we can therefore circumscribe the discussion to the situation where g B is small and change g A . In the following, we will use harmonic oscillator units and scale all lengths in units of oscillator length a 0 = /(mω) and all energies in units of level spacing ω. To solve the Hamiltonian (1) we use two different numerical approaches: direct diagonalization [41] and DMC [65]. The former allows us to calculate the full density matrix of the system and therefore gives us access to all single and multi-particle correlations. However, since it is limited to small particle numbers, the latter will be used to check for trends when the number of particles becomes larger. While DMC is well described in the literature, let us briefly explain our approach to direct diagonalization. For this we expand the second quantised field operators into eigenfunctions, φ n (x), of the single-particle (SP) Hamiltonian for the harmonic oscillator where the creation and annihilation operators,â † k andâ k , satisfy the bosonic = 0, and similarly forb † k andb k , while all commutators between operators belonging to different species vanish. Here, n A(B) is the number of modes used in the expansion. The Hamiltonian can then be written as [45] where The ground state can be expressed in terms of Fock vectors where D n A(B) ,i !) − 1 2 and Φ 0 is the vacuum. The occupation numbers of the n A (n B ) modes for each component are given by Note that Ω increases exponentially with the number of particles and modes, which is the reason why the numerical solution using this approach is restricted to small numbers of atoms.
A good ansatz for the unnormalized ground-state wavefunction of the mixture when g B = 0 and outside of the phase-separated regime can be constructed using the solution for non-interacting atoms in the harmonic trap, Φ(X) = exp Here the 1D s-wave scattering length a A for the intra-species interactions and a AB for the inter-species interactions are related to the 1D coupling constants as g A = −2 2 /(ma A ) and g AB = −2 2 /(ma AB ) and we assume that both coupling constants are non-negative corresponding to repulsive interactions. For practical purposes, we find that the coupling constant g = 20 is close enough to the infinite limit, and therefore we use this value in the direct diagonalization method in describing the large coupling constant limits.

Balanced Mixtures
In the following we will first concentrate on systems in which both components have the same particle number. Unless otherwise stated, we will use N A = N B = 2.

Densities
The main feature of the density evolution in this system is the occurrence of phase separation for increasing inter-species interactions. However this process takes two, fundamentally different forms: in the composite fermionization limit atoms of different species avoid each other even though the species' densities still occupy the same space, whereas in the phase separation limit the overlap of the respective densities goes to zero.
The density along the transition from the BEC-BEC limit (all couplings small) to the composite fermionization limit (g AB large) is shown in Figs. 2(a-b). There are crucial differences in the evolution of the density along the transition from the TG-BEC to the phase separation limit (Figs. 2(c-d)). One immediately notices that the transition into the composite fermionization state happens at a finite value of g AB ∼ 2, whereas the transition to the phase-separated regime happens already for very small values of g AB . Also the final state reached in the composite fermionization or the phase separation limit are very different.
This difference in the final states can be understood by looking at the one-body density matrix (OBDM) given by with a similar expression for ρ B 1 (x, x ). The decomposition in terms of natural orbitals f k (x) of the OBDM and their corresponding occupations λ k A is given in Eq. (11). The densities shown in Fig. 2 are the diagonals of these matrices, calculated with direct diagonalization. As discussed in Ref. [45], the OBDM of both components in the composite fermionization limit are identical and show two peaks. Contrary, in the phase separation limit the OBDM of B shows a single peak located at the center of the trap, while the OBDM for A shows two peaks at the edges. The largest used value of the coupling constant g = 8 is big enough, so that the density profiles shown in Figs. 2 are practically the same as for the infinite coupling constant.
Finally, the transition from the composite fermionization to the phase-separated regime is shown in Figs. 2(e) and (f). One can see that the spatial separation of the clouds happens for a finite value of g A . At the transition between both limits, the OBDM of both species show a complicated structure, which we discuss in detail in subsec. 3.4.

Coherence and Entanglement
Since increasing the coupling constant will drive the system from the weakly to the strongly correlated regime, the coherence is a good quantity for identifying different regions in the phase diagram. It can be characterised by the largest eigenvalue of the OBDMs (11), λ , which provides the largest occupation of a natural orbital. In our numerical calculations with direct diagonalization we normalize the OBDM to 1 instead of the number of atoms. In Figs. 3(a) and (b) we show the largest occupation numbers for the A and the B species, respectively, over the whole range of interactions. Note that the sum of all eigenvalues of each component sum up to 1, in accordance with the chosen normalization.
One can see from Fig. 3(a) that the coherence in the A species decreases monotonically along the transition from the BEC-BEC (λ A 0 = 1) to the TG-BEC (λ A 0 ∼ 0.7) limit, as well as to the composite fermionization limit (λ A 0 ∼ 0.55). However, the transition for increasing g A at a finite g AB shows that a maximum of coherence is reached for finite values of g A ∼ 5, which corresponds roughly to the value where the cloud de-mixing happens (see Figs. 2(e) and (f)). This maximum in coherence within species A is very surprising, as usually the presence of interactions is thought of as detrimental to coherence. Here, however the presence of interactions within the A component to a certain degree "counterbalance" the interactions between the species and therefore allows to re-establish a higher degree of coherence again. Note that after the de-mixing transition the coherence within species A goes down again, which is a clear indication that the enhancement is somehow mediated using the overlap with species B. As expected, species B shows a large degree of coherence in all limits, except the composite fermionization one (see Fig. 3(b)) . However, the re-establishment of coherence along the transition from composite fermionization to the phase-separated limit happens over a definite and narrow region, which corresponds to the area in which the coherence in species A shows a maximum.
One might, at this point wonder how the transition to phase separation manifests itself during the transition from the TG-BEC to the phase-separated limit, as no obvious signature is visible in the coherence phase diagram. The answer is that phase separation happens already for small values of g AB , which can be seen in Figs. 2(c).
It is important to observe that there are no phase transitions in the whole phase diagram. The ground-state energy is always a continuous and smooth function of the parameters, so that the transition between the different regimes is of crossover type.
Closely related to the coherence in the sample is the entanglement between the two components. This can be quantified by calculating the von Neumann entropy, S A = −Tr (ρ A ln ρ A ), which is a function of the reduced density matrix for a single component Here ρ = |Ψ 0 Ψ 0 | is the density matrix, Ψ 0 is the system ground state, and is the Fock vector for species B only. This matrix is obtained by means of direct diagonalization. In Fig. 3(c) we show the von Neumann entropy S A along the transition between BEC-BEC and composite fermionization. S A can be seen to approach a constant value as g AB is increased, corresponding to the large inter-species correlations present in the composite fermionization. The same plot also shows S A along the transition between TG-BEC and phase-separated limit. The two species are less correlated throughout this transition, but still S A saturates to a constant value in the phase-separated limit. In Fig. 3(d) we plot S A for different values of g AB when g A is tuned from zero to a large value. When g AB = 20, this corresponds to the transition between composite fermionization and a phase-separated gas. We observe a peak which coincides with the crossover between both limits. This peak disappears as g AB is reduced, as observed in the curves for g AB = 4, 2 in Fig. 3(d) (for g AB = 0, S A is zero for every value of g A ).

Interaction Energies
An interesting question is how the interaction energy changes across the transitions between the different limits. The average interaction energy in species A is We display this energy in Fig. 4(a). For zero g A there are no interactions between A atoms and U A is equal to zero. By increasing g A the energy U A first grows as correlations are being introduced. For larger repulsion, particles avoid each other which leads to very strong correlation and the interaction energy drops down to zero. Starting from the BEC-BEC region, this is a long drawn process, however for a finite value of g AB this happens over a very well defined domain of the parameter g A , located at small values of g A . Note that for g A = 0 and in the presence of interaction with species B the particles in species A are much more localised than for g AB = 0. Therefore, small increases in the interaction strength g A leads to strong increases in the interaction energy U A . This is also consistent with the maximum found in the correlation strength within component A. The interaction energy goes to zero in the TG-BEC limit, which is the behaviour expected for a single component gas [29,30,31,32], as the increased energy is now stored in the single particle harmonic oscillator energies. During the whole process the total energy is increased from to Quantum correlations and spatial localization in 1D ultracold bosonic mixtures . . . The energy for N A = N B = 2 is shown in Fig. 4(c). The energy obtained by the direct diagonalization and DMC methods coincides. For no interactions between different species, g AB = 0, the energy can be expressed as E = ω + E 2 (g A ), where E 2 (g A ) is the energy of two trapped particles interacting with the coupling constant g A [51]. In order to prove that the described limits exist in larger systems, we calculate the energy for N A = N B = 10 particles with DMC method. The energy per particle in the BEC-BEC limit (15) does not depend on the number of particles, E BECBEC /N = 1/2. We show it in Figs. 4(c,e) with green circles for N A = N B = 2 and 10.
In Figs. 4(d,f) we depict the energy per particle as a function of g AB , starting from the BEC-BEC (solid line) and the TG-BEC (dahsed line) limits. Here the green circles (red squares) indicate the energy per atom in the BEC-BEC (TG-BEC) limit. The  Figure 5. (Color online) The first and second columns show the OBDMs, the third and fourth columns the TBDFs, each time for species A and B, respectively, and the last column shows the CTBDF. g AB is large in all panels. The first and last row displays the numerical result obtained for a value of g A just before and after the crossover, respectively, and the middle row shows the results obtained from calculating the OBDMs and TBDF directly using the ansatz given in Eq. (9) with a A = a AB = 0. Good agreement is clearly visible with the numerical results before the crossover. energy in the TG-BEC limit given by Eq. (16) is E TGBEC /N = 3/4 for N A = N B = 2 and E TGBEC /N = 11/4 for N A = N B = 10. In the transition from the BEC-BEC limit to the composite fermionization one, the energy saturates to certain value, for which we do not have an analytical prediction. As well, a monotonic behavior is observed in the transition from the TG-BEC to the phase separation limit (Figs. 4 (d) and (f)). The average interaction energy between both species, given by is important to quantify the transition to the composite fermionization or the phaseseparated regime. The interaction energy rapidly increases from zero to a maximum at g AB ≈ 2 (see Fig. 4(b)) and decreases again towards zero for g AB → ∞. For g A = 0 this corresponds to building up strong correlations between the particles of different species in the composite fermionization limit, whereas in the limit of large g A this reflects the transition to macroscopic phase separation of the two components.

Correlation Matrices
Since in the presence of strong interactions the system has non-trivial many-body correlations, it is interesting to look not only at single-particle densities, but also at pair-wise correlation functions. The single-particle densities are quantified by the OBDM, Eq. (10). For the particles of the same species, the two-particle correlations are quantified by the two-body distribution function (TBDF) with an analogous expression for B. If the two atoms stem from different species, their pair-wise correlations are captured in the cross two-body distribution function (CTBDF) given by Both functions are proportional to the joint probability for finding two atoms at two given positions. It was shown in Ref. [42,45] that the correlation functions are very useful for a description of the composite fermionization and the phase separation limits. In the following we will carefully look at the transition between these two limits. The phase separation occurs for g A and g AB large and implies a density distribution with atoms of species B are localized at the center of the trap, while the atoms of species A gather at the edges of the density of B. As discussed above, the de-mixing point can also be identified in the coherence, the interaction energies and the entanglement.
In Fig. 5 we show the OBDMs, TBDFs and CTBDFs just before (g A = 5) and just after (g A = 7) the crossover. The upper row and lower row show numerical results while the middle row represents the analytical results obtained from ansatz (9) with a A = a AB = 0. One can see that just before the crossover the densities of both species, i.e. the diagonals of the OBDMs, significantly overlap (panel (a) and (b)), whereas the overlapping is greatly reduced after the crossover (panels (k) and (l)). The TBDFs and CTBDFs before and after the crossover (panels (c) to (e) and (m) to (o), respectively) demonstrate that the atoms of species A are anticorrelated with themselves and with the atoms of species B, as both functions vanish along the diagonal. Note that at the same time atoms of species B are not strongly correlated. This is also captured by ansatz (9), where strong correlations are induced by zeros whenever A-A or A-B atoms overlap (see panels (f) to (j)). All densities and pair correlations computed with this ansatz qualitatively resemble the exact correlation functions just before demixing. However, the ansatz fails to describe the ground state of the system once the system has phase separated.
Let us note that the TBDF for the A species shown in Fig. 5(c) corresponding to the crossover for N A = N B = 2 look similar to those obtained for N A = 4 and a very heavy atom in component B (discussed in [54,55]) or a large number of atoms in B (discussed in [41]). Those cases belong to the phase-separated limit, in which B formed Figure 6. (Color online) Densities with N A = N B = 2, 4, 6, 8, 10 atoms. a) densities plotted in the composite fermionization limit, showing that the two peaks appear farther appart as N is increased. b) densities for B in the phase-separated limit. The atoms tend to localize more and more in the center as N is increased. c) densities for A in the phase-separated limit. The atoms of A are in the edges of B, forming two TG gases with N/2 atoms. a material barrier. Therefore, the two atoms of A stay at each side of B. Very differently in this case, there are only two atoms of A, and they can be localized in either side of B.
For N A = N B > 2 the results discussed above remain qualitatively valid. We show in Fig. 6(a) the densities for the composite fermionization limit when N A = N B = 2, 4, 6, 8, 10 calculated with DMC. In this situation, the OBDMs are equal for both species. The two peaks present in the density tend to spatially separate as N is increased, as a consequence of the large repulsion between both species, which increases with the number of atoms. In Fig. 6(b) and (c) we show the densities for B and A, respectively, in the phase-separated limit. As N is increased, the atoms of B have a greater tendency to localize in the center of the trap. The numerically calculated density for A shows that this component is localized at each side of B, forming two TG gases with N/2 atoms in each side.
The difference in the energy between BEC-TG and TG-TG regimes is further increased in balanced systems of a larger size, N A = N B 1. Indeed, according to Eq. (15), the energy in the BEC-TG scales linearly with the number of particles N , which is a typical behavior of weakly interacting bosons. Instead, in TG-TG limit according to Eq. (16) the dependence on N is quadratic. The resembles the behavior of the energy of fermionic particles and is a manifestation of Girardeau mapping. Comparing the results for N A = N B = 2 with N A = N B = 10 we already observe how the difference in the energy between limits increases, see Fig. 4.

Effect of a larger population in the weakly interacting species
In the imbalanced case, N B > N A , the wavefunction (9) can be equally used as an ansatz for the exact ground state of the systems. The four limits discussed above equally exist. Nevertheless, the weakly interacting species has now a greater tendency to localize in the center of the trap and condense, which modifies the boundaries between the different regimes associated to these limits. In Fig. 7(a)-(b) and (e)-(f) we report the largest eigenvalue of the OBDM for species A and B to quantify the coherence, covering the whole range of coupling constants, when N A = 2 and N B = 3, 4, respectively. As N B is increased we observe that the region in which B is not condensed is reduced (the light blue area in Figs. 7 (b) and (f)). Moreover, the minimum value of λ B 0 , which occurs in this non-condensed area, grows with N B for fixed N A . Notice also that the area in which λ A 0 approaches the largest possible value λ 0 = 1, i.e. close to the g A axis, is reduced as N B is increased. In Figs. 8 (a) and (b) we show the density profiles for A and B along the transition between the BEC-BEC limit and composite fermionization, for N A = 2 and N B = 4. The atoms of species B are now more concentrated in the center than when both populations were equal, even though species B is still not fully condensed. The two peaks in species A appear at a smaller value of g AB , and are more spatially separated than in the case N B = N A . We note that in the composite fermionization limit, the density of species A in the center for the balanced case is finite, while in the imbalanced case it vanishes (compare Figs. 2 (a) and Figs. 8 (a)). The density profiles along the transition between the TG-BEC and the phase-separated gas are presented in Figs. 8 (c) and (d). Comparing with the balanced case plotted in Figs. 2 (c) and (d) we notice that, in the phase separation limit, the two peaks in the density profile of A are now more separated and the squeezing in the density of B is smaller. The average interaction energy U AB (Fig. 7 (c) and (g)) tends to zero when phase separation occurs. Figs. 8 (e)  Fig. 2. The density for B in the composite fermionization limit is more similar to the one obtained in the phase separation limit. and (f) report the density along the transition between composite fermionization and phase separation. We observe that the position of the two peaks in the density profile of A in the phase-separated and the composite fermionization limit is closer than in the balanced case (compare with Figs. 2 (e) and (f)). Also, the crossover occurs now at a smaller value of g A . The average interaction energy U A (Figs. 7 (d) and (h)) decreases abruptly to zero after the crossover. We conclude that for larger imbalances, N B N A , the composite fermionization region is highly suppressed, and therefore the surviving limits are those associated to BEC-BEC, TG-BEC and the phase-separated mixtures.
If the macroscopic limit is reached in such a way that the number of atoms in one of the species is fixed, the minority species plays role of an impurity which perturbs the majority species. The relative contribution of the minority species to the energy becomes smaller and polaronic description might be applicable.
Current experimental advances in ultracold atomic physics allow one to scrutinize the onset and evolution of correlations in few-atom bosonic fluids. Small samples can be trapped, and their interactions can be largely tuned, thus providing a fantastic ground to understand how quantum many-body correlations build in small samples. Binary mixtures are specially appealing as they provide the first step towards understanding the effect of environments on quantum systems in a controlled way. To advance in that direction, we study the effect of embedding a quantum fluid (component A) within a second quantum fluid (B) with tunable intra-and inter-species interactions at zero temperature. We fix the coupling constant of B-B interactions to that of ideal bosons, g B = 0, and vary A-A and A-B interactions in a wide range, 0 ≤ g A < ∞, 0 ≤ g AB < ∞. This permits us to explore the phase diagram for a variety of regimes. The energy, one-and two-body correlation functions, density profiles and von Neumann entropy are calculated exactly using diagonalization method. For larger system sizes, the results are complemented with the energy and density profiles obtained by diffusion Monte Carlo method.
We have described the transition between the following four limits: a) BEC-BEC limit, where both components interact weakly and thus remain condensed, b) BEC-TG limit, where the two components interact weakly among each other and A has strong intra-species interaction, c) composite fermionization limit, where the interaction between both species is large, inducing strong correlations within both species, and d) a phase separation limit, where both the intra-species interaction in A and the inter-species interactions are large. We show that the transition between the different limits involves sophisticated changes in the one-and two-body correlations. The energetic properties change in a smooth way, with the energy and its derivatives remaining continuous, which implies a transition of a crossover type rather than a true phase transition. At the same time, the entanglement between the two components has a much sharper dependence on the interactions. This is demonstrated by reporting the von Neumann entropy, which manifests a sharp peak along the transition between composite fermionization and phase separation. The evolution of the density profiles of A and B components is studied in detail both for the balanced and the imbalanced case. The effect of a large number of particles on the energy and the density profiles is discussed. We analyze the coherence properties by expanding the one-body density matrix in natural orbitals and obtaining the occupation numbers. We demonstrate that full condensation (largest occupation number equal to one) for A species is reached only in the BEC-BEC regime, while the weakly interacting B species also remains fully condensed in the TG-BEC regime, and the condensation is almost complete in the phase separation regime. We argue that the described picture of the transition between four mentioned regimes remain valid also in a macroscopically large balanced mixtures, N A = N B → ∞. Contrarily, when the macroscopic limit is reached by increasing the number of atoms of the weakly-interacting species, N B → ∞, the composite fermionization limit is suppressed. Therefore the phase diagram in this highly imbalanced case resembles the one expected within a mean-field approach. The studied effects are relevant to ongoing and future experiments with small two-component systems.

Acknowledgments
This project was supported by Science Foundation Ireland under Project No. 10/IN.1/I2979. We acknowledge also partial financial support from the DGI (Spain)