Correlations and synchronization in a Bose-Fermi mixture

We study a Bose-Fermi mixture within the framework of the mean-field theory, including three possible regimes for the fermionic species: fully polarized, BCS, and unitarity. Starting from the 3D description and using the variational approximation (VA), we derive 1D and 2D systems of equations, under the corresponding confining potentials. This method produces a pair of nonlinear Schr\"{o}dinger (NLS) equations coupled to algebraic equations for the transverse widths of the confined state. The equations incorporate interactions between atoms of the same species and between the species, assuming that the latter can be manipulated by means of the Feshbach resonance (FR). As an application, we explore spatial density correlations in the ground state (GS) between the species, concluding that they strongly depend on the sign and strength of the inter-species interaction. Also studied are the dynamics of the mixture in a vicinity of the GS and the corresponding spatiotemporal inter-species correlation. The correlations are strongly affected by the fermionic component, featuring the greatest variation in the unitary regime. Results produced by the VA are verified by comparison with full numerical solutions.


Introduction
Cold atom gases trapped in magneto-optical potentials have tremendously extended the range of experimental possibilities in condensed-matter physics. The possibility of controlling the intensity of the interactions, the depth of the trapping potential lattices, and their geometry enables emulating a broad variety of solid-state settings and the creation of novel states of the quantum matter [1]- [8].
In this context, a great variety of experimental and theoretical studies were stimulated by the possibility of mixing degenerate bosonic and fermionic gases. The first experimental works on Bose-Fermi mixtures (BFMs) were conducted with lithium isotopes [9,10], which was a precursor to the study of new combinations, such as 174 Yb-6 Li [11] and 87 Rb-40 K [12]. Particular attention has been paid to heavy-atom mixtures, e.g., 87 Sr-84 Sr [13] isotopes with a large nuclear spin, that have been proposed as prototypes for handling the quantum information. One of important points of the work in this direction is the use of Feshbach resonances (FRs) in the mixtures, as they make it possible to control the two-body interactions between the species. The FR has been observed for the 87 Rb-40 K mixture [14,15], and a giant FR has been reported for 85 Rb-6 Li [16]. Further, five FRs were observed in 6 Li-133 Cs [17], and over 30 resonances are known to occur in the 23 Na-40 K mixture [18]. In Ref. [19], a multiple heteronuclear FRs was reported in a triply quantum-degenerate mixture of bosonic 41 K and two fermionic species, 40 K and 6 Li. In the latter context, a wide s-wave resonance for combination 41 K-40 K transforming the system into a strongly interacting isotopic BFM immersed into a Fermi sea.
Parallel to the experiment, several theoretical works addressed the dynamics of BFMs [20]- [26]. An approach that has proven to be very useful for describing the ground state (GS) of the mixtures is provided by the mean field theory [27]- [31]. In this context, the FRs for mixtures of 23 Na- 6 Li, 87 Rb-40 K, 87 Rb- 6 Li, 3 He-4 He, 173 Yb-174 Yb, and 87 Sr-84 Sr have been studied [32,33].
The possibility of holding ultracold gases in magneto-optical traps makes it possible to create effectively one-or two-dimensional (1D or 2D) situations. A great number of works addressed such 1D settings [34]- [41]. In particular, the variational approximation makes it possible to reduce the 3D description to 1D or 2D [42,43]. In many cases, these approximations do not take into regard the fact that the transverse confinement width(s) may vary with the particle density in the unconfined direction. The variational approximation that does take this fact into account was developed for bosons in Ref. [44] and for fermions in Ref. [45], demonstrating the high accuracy of the 1D and 2D approximations in comparison with the full the 3D description.
In this paper we analyze the BFM at zero temperature confined to 1D or 2D by external potentials. First, we consider the cigar-shaped configuration produced by the strong confinement with the cylindrical symmetry acting in the transverse plane, while a generic weak potential may act in the axial direction. Second, we consider the disc-shaped BFM corresponding to the strong confinement applied in one direction, and a weak 2D potential acting in the perpendicular plane. We use the variational method developed in Refs. [44,45] to derive effective 1D and 2D nonlinear Schrödinger (NLS) equations for the cigar-shaped and disc-shaped configurations, respectively. In particular, our results apply to the 7 Li- 6 Li BFM. In both cases (1D and 2D), we analyze the influence of the density of both species on the width of the mixture in the confined dimensions, the main control parameter being the inter-species scattering length, a BF . The analysis is performed for three different regimes of the fermionic component: polarized, BCS, and unitarity. Using the dimensionally-reduced system, we consider spatial correlations between densities of particles of both species. Considering perturbations around the GS, we study the spatiotemporal correlation between the species too.
The paper is organized as follows: In Sec. 2, we derive 1D and 2D variational equations from the full 3D description. In Sec. 3, we apply our equations to the numerical analysis of the GS in the BFM, and consideration of the spatial correlation between the species. For the 1D setting, we calculate spatiotemporal correlation in the form of the so-called Pearson coefficient, considering perturbations around the GS. Conclusions are presented in Sec. 4.

The variational approximation
We consider a diluted superfluid mixture formed by N B bosonic atoms of mass m B , and N F fermionic atoms of mass m F and spin s F .
The atoms interact through the pseudopotential, δ(r) [1]. We assume that bosons form a Bose-Einstein condensate (BEC), described by the Gross-Pitaevskii equation [1], while the local density approximation [1] applies to the description of the weakly interacting fermionic component. Accordingly, the dynamical equations for the BFM can be derived from action S, where L B and L F are the Lagrangian densities of the Bose and Fermi components, while L BF accounts for the interaction between the them: Here C F ≡ [6π 2 / (2s F + 1)] 2/3 is a constant that depends on spin s F ; g B ≡ 4π 2 a B /m B , , and g BF ≡ 4π 2 a BF /m BF are the three interaction  (3), along with s F , determine three different regimes [46]- [49] listed in the Table 1.
We apply the formalism developed below to the 7 Li-6 Li mixture, with the same scattering parameter for both species, a B/F = 5nm. The use of isotopes of the same alkali element is suggested by the similarity of their electric polarizability, thus implying similar external potentials induced by an optical trap. Unless specified otherwise, in what follows below we consider configurations with fully polarized fermions. Values of parameters for all fermionic regimes are collected in Table 1. The BCS and unitarity regimes involve more than one spin state of fermions, hence the magnetic trap will split the respective spin energy levels. For this reason, we assume the presence of the optical trap, which supports equal energy levels for all the spin states, making it possible to discriminate different regimes of the interaction in the Bose-Fermi mixture. In the BCS and unitarity regimes, we assume balanced populations of the two spin components.
Finally, varying action S with respect to Ψ * B and to Ψ * F , we derive the following system of nonlinear Schrödinger equations for bosons and fermions: As is known from previous works, the mean-field description of the fermionic component, based on the effective NLS equation (6) is valid in the hydrodynamic approximation, which implies that the Fermi component forms a correlated superfluid. In static situations, this approximation is relevant for sufficiently smooth configurations, whose characteristic spatial scale is much greater than the de Broglie wavelength at the Fermi surface [50,51]. In various dynamical settings, both pure fermionic and mixed Bose-Fermi ones, the hydrodynamic approach remains relevant for slow quasi-adiabatic evolution [42]- [49]; [52,53].
While numerical integration of this system in the 3D form system is very heavy computationally, the effective dimension may be reduced to 1D or 2D when the system is tightly confined by a trapping potential. To this end, the variational method is employed, making use of the factorization of the 3D wave function, which includes a Gaussian ansatz in the tightly confined transverse directions. The factorization has been widely used for Bose and Fermi system separately, as it shown in Refs. [44] and [45], respectively. In next two sections we apply the variational technique to a Bose-Fermi mixture, obtaining results that agree well with full 3D simulations. On the other hand, the applicability of the Gaussian ansatz for approximating of the transverse part of the 3D wave function is determined by a relation of the transverse-confinement strength, and the effective strength of the bosonic or fermionic nonlinearity. In case the latter factor is stronger, the transverse wave function may be used in the form of the Thomas-Fermi (TF) approximation [53].

The one-dimensional variational system
The common form of the confinement is provided by the tight harmonic-oscillator potential acting in the trapping dimensions. Thus, the confinement to 1D (the cigarshaped configuration along the z-axis) corresponds to the following potentials, which include weak axial components, V 1D,B/F : Assuming that the transverse trapping potential is strong enough, the dimensional reduction is carried out by means of the usual factorized ansatz for the wave functions, Ψ B/F (r, t) = v B/F x, y; σ B/F (z, t) f B/F (z, t), where the transverse groundstate Gaussians, v B/F ,with widths σ B/F , and the axial functions, f B/F , are normalized to 1 and N B/F respectively: For both species, we define the axial density as n 1D,B/F ≡ f B/F 2 . The relevance of the use of the Gaussian form in the confined directions for bosons and fermions was demonstrated in in Refs. [44] and [45], respectively, by showing that the respective factorization produces GSs close to their exact 3D counterparts. As demonstrated below, the Pauli principle makes the transverse width of the Fermi gas larger, in comparison with the Bose gas in the presence of a similar external potential. Replacing Ψ B/F by the factorized density in the Lagrangian density and integrating in the transverse plane of (x, y), the expression for action (1) is cast into takes the following form: where the effective 1D densities are and e 1D,B and e 1D,F are the 1D energy densities of the boson and fermion species, respectively: with C 1D,F ≡ (3/5)(6π(2s F + 1)) 2/3 in Eq. (14). Next, varying action S given by Eq. (9) with respect to f * B and f * F yields the respective Euler-Lagrange equations, i.e., the motion equations for the BFM in the 1D approximation: Relationships between σ B/F and f B/F are obtained by varying the 1D action (9) with respect to σ B/F : . Thus, Eqs. (15)-( 18) constitute a system of four 1D coupled equations produced by the reduction of the underlying 3D system (5), (6). A numerical solution of Eqs. (17) and (18) by means of the Newton's method produces the dependences of the transverse widths σ B and σ F on the 1D densities, n 1D,B and n 1D,F , which are shown, respectively, in the left and right columns of Fig. 1, for the BFM with the attractive and repulsive interactions (the top and bottom rows, respectively). In all the cases, the Fermi species has, generally, a greater transverse width than its Bose counterpart. We consider for the bosonic density a range greater than for its fermionic counterpart because our calculations correspond to BFMs in which the number of bosons exceeds the number of fermions. For the attractive mixture (a BF = −25 nm), it is observed that, for the bosons (see Fig. 1(a)) the width, σ B , slightly increases with n 1D,B and n 1D,F . For the fermions (see Fig. 1(b)), the width, σ F , strongly increases with n 1D,F , due to the Pauli principle, and decreases with n 1D,B , because of the attractive boson-fermion interaction. In the repulsive mixture (a BF = 25 nm), the situation for bosons is almost the same as in the case of the attraction, with the difference that σ B slightly decreases with the fermionic density, n 1D,F . For the fermions (see Fig. 1(d)) the increase in the boson density (n 1D,B ) generates a significant increase in σ F when the fermionic density (n 1D,F ) is low. These four diagrams demonstrate the significant dependence of the widths of both species on the densities, thus showing the importance of treating σ B and σ F as the variational variables.
The above results can be explained by the dependence of the interaction between the species on their densities. To do this, we define a single scattering parameter, and its counterpart for the fermions, obtained by replacing subscript B with F. Making use of the so defined interaction parameter 19, in Eqs. 15 and 16, these equations assume the form given by and where the Eqs. 20 and 21 represented the motion equations for bosons and fermions separately, with effective scattering parameter g e,B and g e,F respectively. Figure 2 shows the variation of the effective scattering parameters, g e,B and g e,F , with the inter-species scattering length, a BF . Figure 2(a) shows that the interaction between the bosons is suppressed when the interaction is attractive, and increases when it is repulsive. When the fermionic density is much lower that the bosonic density (circles), the influence of the variation of the interaction between the species on the effective interaction strength is minimal. Figure 2(b) shows the variation of the effective scattering parameter of the fermions (g e,F ) under the same conditions as in Fig. 2(a), demonstrating the strong influence of the inter-species interaction on the fermions, due to the higher bosonic density. It is worthy to note that, for the range of a BF presented in Fig. 2(b), g e,F is about three orders of magnitude higher than g F . The inset in Fig. 2(b) additionally shows the dependence of g e,F on the fermionic density n 1D,F (keeping n 1D,B fixed). A much larger number of bosons than that of fermions (the situation that that we address in this work), in the case of an attractive mixture, implies that the fermionic profile will be close to that of the bosonic gas. Accordingly, the Gaussian ansatz is appropriate for in both species, which is confirmed below by the comparison between 3D numerical and variational solutions.

The two-dimensional variational system
To derive 2D equations for the disc-shaped gas, the shape of the confinement potential is taken as (cf. Eq. (7) for the cigar-shaped configuration): where the second term corresponds to the strong harmonic trap acting along the z direction.
The corresponding factorized ansatz is adopted as where u B/F and φ B/F are normalized to 1 and N B/F , respectively, and u B/F is represented by the Gaussian wave function, ξ B/F (x, y, t) being the widths of the gas in the confined direction. Substituting the factorized ansatz (23) into action (1) and integrating over z, we arrive at the following expression for the effective 2D action: where n 2D,B/F ≡ φ B/F (x, y) 2 are the 2D particle densities of the boson and fermion species, and e 2D,B and e 2D,F are their energy densities: with C 2D,F ≡ (3/5) 1/2 (6/(2s F + 1)) 2/3 π in Eq. (29).
The motion equations of the 2D system are obtained by the variation of the action S given by Eq. (24) with respect to variables φ B and φ F : Relations between ξ B/F and φ B/F are produced by the Euler-Lagrange equations associated to ξ B/F : (32) and (33) for ξ B/F were solved numerically. The four graphs in Fig. 3 show the dependence of ξ B and ξ F (in the left and right columns, respectively) on n 2D,B/F , for the attractive and repulsive mixtures (top and bottom panels, respectively). The analysis for the attractive case can be divided into two parts: when the boson density is high, the fermions gas tends to compress, so that at low densities its width is found to be even less than that of the boson gas, a situation that tends to reverse for high fermionic densities; for low bosonic densities, the width of the boson gas rapidly decreases with the increase of the fermionic density, until becoming nearly constant, while, on the contrary, the width of the fermion gas increases progressively when its density is higher. Conversely, when the mixture is repulsive, the situation is similar to that outlined above for the 1D setting: low densities of the mixtures cause the gas to compress, while high densities cause the mixture to expand, this trend being much stronger for the fermion component.
To determine the influence that one species exerts on the other in 2D, we define the effective scattering parameters, cf. a similar definition (19) adopted in the 1D case: g e,F being produced by replacing subscript B with F. As shown in detail below, our calculations were performed assuming the same order of magnitude for the number of particles for both 1D and 2D cases (for this reason that the effective 1D density is about 10 times the magnitude of its 2D counterpart). With regard to this, conclusions concerning the effective interaction strength, which follow in the 1D setting from Fig.  2, apply to 2D case as well.

Numerical Results
Our analysis was developed for the GS and dynamics of perturbations around it in the presence of the 3D harmonic-oscillator trap of the following form: In particular, for the GS of the 1D and 2D settings we focused on determining the spatial correlation C s between the spatial particle densities in both species, defined as wheren B/F = n B/F − n B/F , standing for the spatial average. For dynamical perturbations around the GS, a spatiotemporal correlation, which is defined by replacing the spatial average with the spatiotemporal average, is known as the Pearson coefficient C s−t [54]. The GSs in the 1D and 2D cases were found by means of the imaginary-time integrations of Eqs. (15), (16) and (30) Fig. 4 3.1.1. Accuracy of the variational method as a function of scattering parameter a BF for the ground states. As a starting point in the 1D case with free coordinate z, we analyze the effect of the magnitude and sign of the interaction parameter on the spatial profile of both species, and the accuracy of the variational method (see Eqs. 15 -18) compared to the 3D solution, by varying the scattering length a D,B . The plots in Fig. 4 shows this situation, with panels 4(a) and 4(b) corresponding to the profiles of n 1D,B and n 1D,F , respectively, where the inset in panel 4(b) specifies whether the profile was produced by the VA or by means of the 3D solution, while 4(c) and 4(d) correspond to the profiles of σ B and σ F . The mixture with many more bosons than fermions is considered here: N B = 5 × 10 4 , N F = 2.5 × 10 3 . Because of this condition, the bosonic profile is mainly determined by its self-interaction and the external potential. First, we deal with the influence of the interspecies interaction (g BF ) on the density profile, and then address an error resulting from the variational approach. Figure 4(a) shows that variations of the bosonic density profile are very small in comparison to the significant changes of the inter-species scattering length. The situation is opposite for the fermionic species. As the repulsive scattering length increases, the fermions tend to be pushed to the periphery of the bosonic-gas density profile. This phenomenon is known as demixing [27,32,35,52]. On the other hand, for the attractive case, fermions are, naturally, concentrated in the same region where the bosons are located. Figures Fig. 4(c) and 4(d) show that the width of the bosonic profile (Fig. 4(c)) slightly increases as going from the interspecies attraction to repulsion. A similar trend is observed for fermions, see Fig. 4(d), but amplified in the spatial zone of the interaction with the bosons, where the gas is compressed in the case of the attraction and expands in the case of the repulsion. It is noted that the fermionic component expands in the confined direction much more than its bosonic counterpart, and that the fermionic width fluctuates markedly with changes in density. Now, to analyze the accuracy of the method we solved Eqs. 5 and 6 and obtained the 1D density associated through of ρ 1D,B/F = Ψ B/F 2 dxdy. The inset in panel 4(a) shows that the difference between the bosonic profiles obtained by both methods (VA or 3D) is ∼ 2% of the maximum density for all cases (the fact that the error changes very little with variations in a BF is a consequence of the greater number of bosons). Figure 4(b) shows that, for the case of the attractive mixture, the variational profile is very close to the 3D result, in particular for the case of a BF = −6nm. For the repulsive mixture, it is observed that the error increases, which is a consequence of the lower fermionic density at the center of the 3D harmonic potential, dominated for the bosons, hence a monotonously decreasing function in the transverse direction, such as the Gaussian, is not a sufficiently good approximation. We define the global error for the variational method as E %,1D = |ρ 1D − n 1D | dz (for both species). Figure 5 shows the global error for both species as a function of the scattering parameter, a BF . This picture demonstrates that the error for the bosonic species is around 2%, and it does not change much (as already noted in the comment to Fig. 4(a)). For the fermions, the error attains a minimum value ∼ 0.5% at a BF ≈ 4nm, thus corroborating that the Gaussian is a very good approximation. This is a consequence of the fact that, for this value of a BF , the interspecies interaction term practically compensates the Pauli repulsion term, making the dynamics of the fully polarized Fermi gas close to that giverned by the Schödinger equation (recall that the Gaussian is the solution for the ground state). When the mixture becomes more attractive, the fermionic dynamics is dominated by the bosons, producing a similar error for both species, while for the repulsive mixture the Gaussian approximation is not appropriate. Interestingly, for the non-interacting mixture, the error for the fermions is smaller than for the bosons, because the fermionic density is very low, making the self-interaction terms weak in comparison to the external potential, therefore it is appropriate to use the Gaussian ansatz for describing the 1D dynamics. Using the TF profile may produce, in principle, better results, but the respective analysis would be cumbersome.
3.1.2. Spatial correlations in the ground state. As the first application, we calculate the spatial correlation between densities of the two species (C s ), as defined by Eq. (36) for a wide range of parameters, ranging from very attractive to to very repulsive mixtures. We keep in mind that the variational calculation produces a minor error for attractive mixtures, therefore the values obtained for the repulsive case will be considered as help to asses the validity of our analysis. The eventual objective is to produce a parameter which determines the mutual influence of both species. In the Fig. 6, the spatial correlation C s is shown versus N F and a BF , which shows that, for an attractive interaction, C s reaches values greater than 0.8 (in the mixed state), whereas for the repulsive interaction C s decreases until reaching values around C s = 0.3 (in the demixed state). It is also noted that the number of fermions strongly affects the dependence of the correlation on the inter-species scattering length. When N F increases, the difference between the maximum and the minimum value of C s decreases, as a BF varies from −25 nm to +25 The initial conditions are the same in all cases and are described in the text. Other parameters are the same as in Fig. 4 nm, tending to reach values close to C s = 0.8, which is the value attained in the absence of the interaction, depending only on the presence of the same trap acting on both components.

Dynamics near the ground state
We begin the study of conservative dynamics, considering a mixture with arbitrary initial conditions for the 1D fields: we assume a Gaussian shape along the z-axes with widths (standard deviation) of 100µm and 200µm for bosons and fermions, respectively. Figure 7 shows three cases of the temporal evolution with these initial conditions for a BF = −6nm, a BF = −12nm, and a BF = −14nm. In the first case (panels 7(a) and 7(b)), it is observed that the densities converge towards the center of the potential in an expected pattern of oscillations around the potential minimum; in addition, the fermions are affected by bosons, as can be seen for the mark left by the bosons in the fermionic density. For the second case (panels 7(c) and 7(d)), it is observed that the increase in the magnitude of the attractive interaction generates dark solitons in the fermionic density, some of which shows oscillatory dynamics very similar to that observed, experimentally and theoretically, in Refs. [55]- [60]. Finally, the last case (panels 7(e) and 7(f)) shows that further increase in the magnitude of the interspecies interaction generates a larger number of dark solitons. In other words, we show that the attractive interaction of fermions with bosons in a state different from the ground state eventually generates a gas of dark solitons. Now we aim to address the system's dynamics in the vicinity of the GS. We start with the GS found in the absence of the inter-species interaction (a BF = 0). Then, at t = 0, we switch the interaction on, which may imply the application of the magnetic field, that gives rise to a BF = 0 via the FR.  observed in Fig. 6, Fig. 9 shows that the spatial-correlation parameter C s , defined as per given Eq. (36) rapidly decreases following the transition from negative to positive values of a BF (the mixing-demixing transition). When the inter-species interaction is more attractive, it increases the correlations, imposing spatiotemporal synchronization of density fluctuations. In the case of strong repulsion, a BF = 50 nm, correlation fluctuations are suppressed because the system falls into a demixed state. Figure 10 presents the analysis of the spatiotemporal synchronization of the mixture, where three possible regimes are considered for the fermions: fully polarized, BCS, and unitarity. Parameters of the the Lagrangian density for each fermionic regime are given in Table 1. The synchronization was calculated through the Person coefficient C s−t defined above. For each of the three cases, the spatial correlation, Cs, was calculated too. When the interaction is attractive, there is no discrepancy between the synchronization and correlation curves. On the other hand, in all the three regimes it is observed that, as the interaction becomes more repulsive, the synchronization of the mixture decreases in comparison to the spatial correlation. This trend enhances as the system approaches the demixing transition. Another generic feature is that stronger correlations between the species near the GS imply a stronger synchronization on the spatiotemporal dynamics too.
The correlation properties are different for the distinct fermionic regimes. C s−t is always positive for the system with fully polarized fermions and in the BCS regime, decreasing as the interaction gets more repulsive. On the other hand, the unitarity regime features significant differences: at first, C s−t grows to a maximum value at a BF = −10 nm; then it decreases, reaching negative values of C s−t /C s in the demixed state. In all the cases, no significant difference are observed between C st and C s in the case of the attractive Bose-Fermi interaction, i.e., a highly correlated GS supports the dynamical synchronization too.   Fig. 11 the spatio-temporal dynamics of the 1D density is displayed as produced by the solution of the 3D equations 5 and 6, using the same procedure as in Section 3, for initial conditions similar to those used in Fig. 7 but with a BF = −10nm. The initial conditions for the 3D dynamics are given by the ansatz in Eq. 8, with the Gaussian profile along the z-axes, similar to those used in the Fig. 7. Panels 11(a) and 11(b) show spatio-temporal diagrams of the bosonic and fermionic density, respectively, making the emergence of dark solitons obvious, see also Fig. 7. This result corroborates that the dark solitons will emerge too in the 3D dynamics, which is approximated by the present 1D model. The other panels show a comparison of the 1D spatial profiles, as obtained from 3D simulations and the 1D VA, for three instants of time: t = 0ms (panels 11(c) and 11(d)), t = 25ms (panels 11(e) and 11(f)), and t = 50ms (panels 11(g) and 11(h)). The results demonstrate that the VA profiles are very similar to their counterparts produced by the 3D simulations, hence the present approximation provides good accuracy and allows one to study dynamical features of the Bose-Fermi mixture in a sufficiently simple form.

Results in two dimensions.
3.2.1. Accuracy of the variational method as a function of scattering parameter a BF for the ground state. In the 2D setting, we studied the GS of the BFM, fixing the number of bosons to N B = 5 × 10 4 (similar to the value used in the 1D case), which was and much greater than the number of fermions. Fig. 12 shows the profile of the 2D bosonic and fermionic densities, n 2D,B/F (top), as obtained from the VA and 3D simulations, and Here a BF = −10nm, the initial conditions and other fixed parameters being the same as in Fig. 7.
the widths of the mixture in the transverse direction, ξ B/F (bottom), with respect to the radial coordinate. To obtain the 2D profile from the 3D simulations (results of solving Eqs. 5 and 6), we integrate the 3D density along the z so that ρ 2D,B/F = Ψ 2D,B/F dz.
The panels for the bosonic and fermionic components are displayed on the left and right, respectively. The results are very similar to what was observed in the 1D case in Fig. 4, for the same value of the physical parameters. Therefore, the following conclusions are also similar to what inferred in the 1D setting: the repulsive mixture concentrates the bosons at the center, while the attractive mixture concentrates both species at the center. Panels 12(c) and 12(d) show that only the width of the fermionic density profile varies significantly with the change of the scattering length of the inter-species interaction, which is consequence of a greater number of bosons than fermions. It is clearly seen that fermions are stronger confined when the interaction is attractive, and their spatial distribution significantly expands when the interaction is repulsive. Similar results have been observed in [27,28,32]. Now, to compare the results obtained from the VA with those produced by the 3D simulations, we note that both profiles are practically identical, except for the repulsive case in which a discrepancy is observed (smaller than in the 1D case). The inset in panel 12(a) shows that the difference between the two results has a magnitude of nearly three orders of magnitude lower than the density itself. The error is lower than in the 1D case because the 2D reduction is closer to the full 3D model. Now, similar to the 1D case, we define the overall percentage error of the VA as E %,2D = |ρ 2D − n 2D | dxdy (for both species). Figure 13 shows the error for both species as a function of interspecies scattering parameter (a BF ). For bosons it is around of 0.2% (one order of magnitude lower than in the 1D case), and does not change much (as shown in the inset to Fig. 12(a)). For fermions the error is greater than for bosons throughout the observed range, but it is quite small for the attractive mixture. Like in the 1D case, the error increases for the repulsive mixture, but remains lower than that calculated for bosons in the 1D case. In summary, the 2D approximation is very accurate, implying, as shown below, that the description of the dynamics is also more accurate in this case.

Spatial correlations of the two-dimensional density
Here we again apply the formalism to the calculation of the spatial correlation C s . The Fig. 14 shows the dependence of the spatial correlation C s , defined as per Eq. (36), on the scattering length a BF for three values of N F and fixed N B = 5 × 10 4 , in the regime of the fully polarized fermionic component. The case with a smaller number of fermions (N = 2500, the same as in Fig. 12) gives rise to the greatest contrast in the values of the correlation: from a value close to 1 for g BF = −25 nm to a value near zero for g BF = 25 nm. This transition is more drastic than that observed for the same parameters in the 1D case. This is mainly due to different factors multiplying the interaction terms in the 2D equations (31), ( 32), in comparison with their 1D counterparts (15), (16). The contrast between the mixing and demixing decreases as the number of fermions increases, due to the strengthening of the self-interaction of the fermions, hence the correlation gets less sensitive to the inter-species scattering length, a BF .
The Fig. 15 shows the spatial correlation, defined as per Eq. (36), for the three fermionic regimes: polarized, BCS and unitarity. The figure is a counterpart of Fig. 10 drawn for the 1D case, but only for the spatial correlation in the GS (C s ). One may conclude that these fermionic regimes show the behavior similar ti that observed in the 1D case, a difference being that the three curves demonstrate stronger demixing when a BF changes from positive to negative values. In the unitarity regime it is again observed that the correlation reaches a maximum close to 1 at a BF ≈ −10 nm, dropping to negative values when the mixture is strongly repulsive. Figure 16 shows, for the attractive mixture (g BF = −12nm), the comparison between the results obtained by the t3D integration (Eqs. 5 and 6) with the results obtained by the 2D integration of the variational equations (Eqs. 30 and 31). The approach is similar to that developed in Sec. 3.1.4 for the 1D case. It starts from arbitrary initial conditions in 2D, viz., 2S fields with Gaussian shapes in the radial direction, a standard deviation of 50µm for bosons and 100µm for fermions being assumed. The initial conditions for the 3D dynamics are obtained with the ansatz given by Eq. 23 and assuming a Gaussian in the xy-plane. Panels 16(a) and 16(b) display two spatiotemporal diagrams, obtained from the 3D integration, showing a cut along the x axis of the 2D density; since the dynamics, for the time range observed, obeys the cylindrical symmetry, this cut is representative. It is seen that the boson profile clearly affects the fermion profile. As the initial conditions do not correspond to the ground state, the fermions are attracted to the center by the bosonic pattern and are subject to the confinement imposed by the harmonic potential, generating strong interaction with the boson core. A gray soliton ring is generated (in panel 16(b) it is seen as a gray soliton). In panels 16(c) and 16(d), the radial profile obtained by the 3D integration and by VA for three instants of time are compared, showing that it is not possible to distinguish them. The inset in panel 16(c) shows that, in the observed range, the 2D overall error is less than 1% (at initial time, it is zero because the 3D initial conditions are obtained from the 2D configuration using the ansatz, Eq. 23). The approach delivers very good results. In addition, in panel 16(d) it was observed that the rings are gray (the minimum value of the fermion density is not zero). Panels 16(e) and 16(f) show diagrams of the 2D density in the xy-plane for the last moment of time, highlighting the gray ring soliton observed in the fermionic component, which expands wherein the boson density is zero. We note that, in this case, the approximation is very good and allows one to considerably reduce the computing time needed to observe the 2D dynamics. The dynamics of the mixture thus produces a vast variety of localized nonlinear structures.

Conclusions
The incorporation of the width of the confined species as a variational variable shows that the term accounting for the interaction between the bosonic and fermionic species in the lower-dimensional equations is inversely proportional to the width in the confined direction(s). The derived systems of equations, which include the algebraic equations for the widths as functions of the 1D or 2D bosonic and fermionic densities. These approximations yield solutions which agree well with the full 3D theory. It is worthy to note that the ensuing relation between the widths and densities significantly varies with the strength of the inter-species interaction, which can be tuned by means of the Feshbach resonance. For very low densities, our equations are similar to the simplest models which postulate constant confinement width in the transverse direction(s). As an application of the variational method, we have studied the spatial correlation between the species, in which a strong dependence on the interaction parameter was observed. This sensitivity becomes stronger for the repulsive interaction, which causes demixing of the species. It was observed too that, in 1D and 2D alike, the spatial correlation between the species is also strongly affected by the inter-species interaction strength. In addition, it was concluded that the greatest contrasts in the correlation, changing from the large correlation in the highly attractive mixtures to a weak correlation in the strongly repulsive ones, occur in the unitarity fermionic regime. Further, in the 1D setting the dynamics near the ground state demonstrates that the spatiotemporal correlation follows the same trends as the spatial correlation, the largest (but still small) difference being observed for strongly repulsive mixtures Finally, it is relevant to mention that the mixture can also give rise to stable dark solitons in the Fermi field, in both the 1D and in 2D cases. However, as the hydrodynamic conditions was assumed, the solitons may occur in some form but are unconfirmed and should not be believed in too strongly. A detailed characterization of these states and the validity of the hydrodynamic approximation for them will be presented elsewhere.