Collective oscillations of a two-component Fermi gas on the repulsive branch

We calculate frequencies of collective oscillations of two-component Fermi gas that is kept on the repulsive branch of its energy spectrum. Not only is a paramagnetic phase explored, but also a ferromagnetically separated one. Both in-, and out-of-phase perturbations are investigated, showing contributions from various gas excitations. Additionally, we compare results coming from both timedependent Hartree-Fock and density-functional approaches.


I. INTRODUCTION
Theoretical and experimental studies of multicomponent quantum mixtures have always played a crucial role in our attempts to understand the quantum theory [1,2]. For several decades, the central point of such considerations was reserved for mixtures of different states of helium [3]. However, the first realizations of quantum degenerate gases in 1990s [4][5][6] have created a playground for scientists that not only allows one to combine different fermionic and bosonic ingredients of the mixture, but also to freely tune the interatomic interactions by the means of Feshbach resonances [7].
One of the main techniques to study the excitations of a trapped, interacting quantum mixture is to examine its collective oscillations [8,9]. Depending of the particular scenario of an experimental excitation procedure, such oscillations can be classified into different categories [2]. Firstly, if the sample does not undergo a change in the volume, but its geometric shape is altered, the mode is called a surface one. Such modes were employed to study a transition from collisionless to hydrodynamic regimes in both bosonic [10,11] and fermionic [12,13] gases. On the other hand, if a volume change happens, one deals with a so called breathing or compression mode [14,15]. Those types of excitations have proved useful in the investigation of equations of state for strongly interacting fermionic gases [16][17][18].
The richness of effects increases even more when decays through different channels and collapses [41,42] [47,48] and Bose-Fermi [49] mixtures in the presence of phase-separating domain wall were studied experimentally. We focus on the former case, in which phase separation can be interpreted as a manifestation of the long-standing problem of ferromagnetic (Stoner) instability [47,48,[50][51][52][53][54][55][56][57][58][59]. In such a system, two equally populated spin species (or in experimental terms -two hyperfine states) of Fermi gas interact only via repulsive short-range potential. However, the ground state of such a mixture is a superfluid of paired atoms of opposite spins (so called lower or attractive branch of the energy spectrum), rather than a ferromagnet (upper or repulsive branch) [57]. It stems from the fact that true zero-range interactions tuned by Feshbach resonances necessarily need an underlying attractive potential with a weakly bound molecular state [7].
Due to phase-separated state being intrinsically unstable towards decay into such molecules, early experiments dealt with high rates of losses [51]. To overcome these problems, in recent experiments artificial domain structure was initially prepared, making the phase-separated state stable for a finite time [47,48].
As for now, theoretical studies have skipped analysis of small oscillations in this system while excited onto repulsive branch, instead focusing mainly on the attractive one [33,[60][61][62][63]. In this work, we fill this void by considering in-and out-of-phase radial and breathing modes of purely repulsive two-component fermionic mixture for both overlapping and phase-separated states. We compare our results with the previously available weakly-interacting ones and show how to refine them by means of the renormalization of the interaction. Moreover, we compare the results obtained by the means of time-dependent Hartree-Fock calculations with the hydrodynamical approach [64], pointing out the regimes of applicability of the latter method.
The paper is structured as follows. In Sec. II, we analyze excitation frequencies by the means of timedependent Hartree-Fock methods and compare them to past theoretical techniques. Sec. III is dedicated to analysis of applicability of hydrodynamical methods in similar by comparing such results to atomic orbital calculations. In the final Sec. IV, we briefly review obtained results and provide suggestions for further work.

II. TIME-DEPENDENT HARTREE-FOCK CALCULATIONS
We start our considerations by describing the method used to analyze the statics and dynamics of the mixture. We employ time-dependent Hartree-Fock (or atomicorbital) approach [64,65] in which we assume that manybody wave function Ψ (x 1 , ..., x N ) of N indistinguishable fermionic atoms is given by a single Slater determinant: The coordinates x i of an atom comprise both spatial and spin variables and ϕ i (x), i = 1, ..., N denote different, orthonormal spin-orbitals. As in scenario that we describe both spin states are equally populated, we further assume that the spin-dependent part of spin-orbitals is twofold and that exactly half of the atoms occupy each spin state.
We restrict ourselves to only contact interaction between different species, as such a description is very well accurate for systems with realistic short-range potentials (e.g. coming from broad Feshbach resonances). Such an interaction can be described by a single parameter, s-wave scattering length a, which can be connected to the coupling constant g ≥ 0 by g = 4π a 2 /m. Therefore, the time-dependent Hartree-Fock equations for the spatial parts of the spin-orbitals, i.e. spatial orbitals representing the first component, ϕ i,+ (r, t), and the second component, ϕ i,− (r, t), are given by for i = 1, ..., N/2 and with n ± (r, t) = Here, m is the mass of fermionic atom, V tr (r) is the trapping potential taken to be a radially symmetric harmonic trap, V tr (r) = 1 2 mω 2 r 2 , and n ± are single-particle densities of both species.
To further refine the interspecies interaction and include many-body quantum corrections, we locally renormalize the coupling strength by replacing the bare scattering length, a, by the effective (and symmetrized) one: where k + (r) and k − (r) are the local Fermi momenta for the first and the second component, respectively, and ζ(k F a) is the renormalization function [64,66]. The renormalization function ζ(k F a) is expanded perturbatively into powers ofk F a: For the physical interpretation and numerical values of consecutive perturbative terms, see e.g. [64]. Taking only first-order terms in the above expansion results in mean-field Eqs. (2). However, taking into account the second-and third-order terms changes the timedependent Hartree-Fock Eqs. (2) in the following way: where The HF equations can describe both statics and dynamics of the mixture -the usual method is employment of the split-step operator to propagate the equations in time. To obtain a ground state, this time is imaginary, and for the dynamics -real.
In Fig. 1 we present the results of propagation of the HF equations in imaginary time, yielding ground states of the mixture for given interactions. We can see that for small interactions the clouds fully overlap, but become larger as the interaction increases. At some critical interaction strength, k F a c , phase separation occurs. k F is the Fermi wave number in the center of a trap in Thomas-Fermi approximation, k F = (24N ) 1/6 /a HO , and a HO is the harmonic oscillator length. In order to get value of k F a c consistent with both more sophisticated theoretical approaches (Quantum Monte Carlo, LOCV, large N expansion) and experimental results, we renormalize perturbatively the scattering length in the way described above. In the trap, this critical value is k F a c ≈ 1.0. For weak interaction, clouds stay in one place, but for the strong one, they separate, entering ferromagnetic phase.

A. Radial compression mode
We start with the case of radial compression (monopole) mode. It is usually performed by symmet-rically perturbing both clouds, by either slightly expanding or deepening the trap, exciting oscillations with radial symmetry preserved. The gas changes its volume, therefore this mode is classified as a compression one. The frequency of oscillations is taken from analysis of the cloud's width in a radial or axial direction.
Firstly, we compare previous theoretical approaches in a weakly interacting regime without a renormalization of the interaction strength. From here onwards, we will use atomic units, = ω = m = 1. The starting point is a sum rule approach, firstly used for the oscillations in spin half Fermi gas by Vichi and Stringari [33]. The expression for a radial mode in a harmonic trap reads We are interested in weakly interacting case, for which phase separation does not occur, and Thomas-Fermi ground-state densities provide a reliable approximation. Under the assumption that spherical symmetry is not broken and both clouds fully overlap, n(r) = n(r) = n + (r) = n − (r), consecutive mean field energies read In the initial works, noninteracting Thomas-Fermi profiles are used, where A = 6 5/3 π 4/3 12 . The result is then quite easily calculated analytically, yielding simple expression This result is shown in Fig. 2 as a blue line. However, one can refine this result by considering ground-states profiles that come from Thomas-Fermi equations with mean-field interactions included (see Eq. 12 in Ref. [59]): where µ is a global chemical potential, that is implicitly given by the number of atoms through normalization procedure. This equation can be readily solved (see Appendix in Ref. [59]) by means of e.g. Cardano formulas, yielding ground state profiles that have analytical, however quite complicated form. Using these density profiles, one gets an orange line in Fig. 2 Apart from curves from sum rule approach, there are two more results presented in Fig. 2: oscillation frequencies extracted from analysis of time evolution of TDHF equations, and frequencies from linearization of hydrodynamical formulation of the problem (further discussed in Sec. III). One can see that all of them give very similar results for this particular mode.
For the further analysis, we decide to use TDHF approach over other two for the following reasons. Firstly, sum rule approach is not well suited to account for a system that undergoes phase separation and moreover has to be provided with ground-state profiles from external theory. As for the hydrodynamical approach, it can be readily utilized for both statics and dynamics of considered system, however it is not obvious for which situation the system really behaves hydrodynamically. On the other hand, atomic orbital approach does not differentiate between collisionless and collisionally hydrodynamical regimes, as it should work well in both of them. The connection between hydrodynamical and atomic orbital methods will be further elaborated in Sec. III.
a. In-phase oscillations Firstly, we analyze radial oscillations with fully renormalized interaction for both weakly and strongly interacting gas. Fig. 3 provides the results. 1 Inclusion of the higher order terms in the interaction shifts the phase transition towards smaller values of k F a and increases the maximum value of the oscillation frequency, which is situated just before the transition occurs. It is now ca. ω/ω 0 ≈ 2.2 in comparison to ω/ω 0 ≈ 2.1 for the bare mean-field model.
After the transition, again there is only one dominating frequency that decreases with growing interacting and decreasing overlap of the clouds, achieving noninteracting value ω/ω 0 ≈ 2.0 for a very strong repulsion. It suggests that in the fully separated regime, the domain wall has no effect on the frequency of the radial oscillations as compared to the noninteracting, polarized gas of N fermionic atoms. This is the only case in which the oscillations can be fitted into a sum of sinuses for the interaction above the critical value. For all the other cases, the analysis is based only on the Fourier transform for each interaction strength.  b. Out-of-phase oscillations Additionally, we consider out-of-phase radial oscillations. It is achieved by perturbing both clouds differently -when one of them undergoes deepening of the trap, the other feels the trap to be widening. Specifically, the small deviations of the densities ought to behave like δn + = −δn − . Potential experimental realization of such a scenario could involve imposing magnetic trapping on top of the optical one. The magnetic field would repel one of the species from the center of the trap, and attract the other one. Moreover, it is not unthinkable to find such a out-of-phase contribution in some other excitation schemes (e.g. quenching from one value interaction to another).
The results are presented in Fig. 4. Below the critical value of the interaction, two branches are clearly visible. The upper one constitutes a contribution from in-phase oscillations that are excited by this particular scheme. The lower branch consists of out-of-phase contribution. Contrary to the upper one, the frequency decreases with the interaction, achieving ca. ω/ω 0 ≈ 1.2 just before the phase transition. The Fourier transform shown as a background confirms that the lower branch is much stronger, with only small part of oscillations excited on the upper one.
After the transition, oscillations tend not to have a clear decomposition into small finite number of frequencies. However, the Fourier transforms of each of the time evolutions show strong contributions from the frequencies that can be considered extensions of two branches from the paramagnetic phase. The upper branch follows the trend from the in-phase perturbation scheme, decreasing and achieving nonincteracting value. The same thing happens for the loweer branch -it increases and goes to ω/ω 0 ≈ 2.0 for very strong interaction. Contrary to the pre-transition regime, upper branch appears to be excited stronger. Additionally, contributions from other, lower-and higher-lying frequencies are clearly visible.

B. Quadrupole modes
The next excitation that we analyze is a compression quadrupole (breathing) mode. The perturbing scheme involves slightly compressing the gas in two out of three directions, leaving one of the axes unperturbed. Again, we perform in- (Fig. 5) and out-of-phase (Fig. 6) calculations.
As for the in-phase case, we again can distinguish between two branches before the transition, however this time they are equally strong. The upper branch again achieves ω/ω 0 ≈ 2.2 near the transition, but the lower one falls down to only ω/ω 0 ≈ 1.8, giving rise to an-other distinctive excitation. After the transition, upper branch appears to be clearly more dominating, falling to noninteracting value for strongly interacting gas.
The out-of-phase case exhibits three distinctive branches -the upper (ω/ω 0 ≈ 2.2 near the transition), the middle (ω/ω 0 ≈ 1.8) and the lower (ω/ω 0 ≈ 1.2) one. The latter clearly dominates. Again, after the transition, each of them is visible, going back to noninteracting regime, with the upper one staying the most pronounced. The last excitation we consider is a type of a surface mode, called a radial quadrupole mode, in which gas is again excited only in two axes, however alternately, giving rise only to change of the shape of the clouds, but not their volume. The results for in-phase perturbing scheme are presented in Fig. 7 and for the out-of-phase one in Fig. 8.
The in-phase excitation shows two distinctive branches -the usual upper one, and the lower one that near the transition achieves ω/ω 0 ≈ 1.8. Contrary to previous cases, the lower branch stays stronger than the upper one after the transition. The out-of-phase case is characterized by the same previous three branches, with the lowest being the most dominating. This time, the upper branch stays the most pronounced over the transition. Frequency of radial quadrupole mode for out-ofphase perturbation of both fermionic species for varying interaction strength.

III. COMPARISON TO HYDRODYNAMICAL APPROACH
We now proceed to briefly compare the atomic orbital approach and the hydrodynamical (in this case sometimes called time-dependent Thomas-Fermi [67]) one. The most important assumption that underlies a hydrodynamical behavior is a fast relaxation of any dynamical distortions in momentum distribution towards a local Fermi sphere centered around local hydrodynamic momentum of the collective flow. Such a thermalization is expected when both clouds overlap and interact by two-body collisions, but is not present in collisionless, noninteracting regime that governs the nonoverlapping parts of the clouds 2 . However, previous studies showed that such a hydrodynamical description works well for a spin-dipole mode for which experimental results are retrieved even in a weakly interacting regime. It is partially due to the fact, that the noninteracting frequency of this particular oscillation coincides with the hydrodynamical one. Then, corrections due to the interaction (even the weak one) seem to stem from the overlapping parts of the clouds, and can be described by the hydrodynamics. It is not obvious for which type of oscillations considered in this work hydrodynamical description can be used. We try to answer this question by comparing TDHF and hydrodynamical results.
The starting point for hydrodynamical methods is the density functional approach introduced for such a repulsive mixture in Ref. [59]. In this theory, the full Hamiltonian is given by H = T tot + E int + E pot , and each of the terms can be written as a functional of one-particle density, making use of local density approximation. The first term, the total kinetic energy T tot = T + T c consists of the intrinsic kinetic energy T , which is approximated by the Thomas-Fermi functional, T = j=± dr 3 5 An 5/3 j , and the kinetic energy of the collective motion, T c = j=± dr m 2 n j v 2 j . The potential energy is takes the form E pot = j=± dr V j n j . The contact interaction term is given as an overlap between the density profiles of the components, E int = g dr n + n − .
To connect this DFT description to hydrodynamics, classical pseudo-wavefunction in the spirit of early Madelung's works [69] is introduced: where n + + n − = ψ † ψ retrieves one-particle density, and ∇χ ± = v ± are the irrotational velocity fields of the collective motion. The Euler-Lagrange equation for these 2 However, some works suggested that Pauli principle can act as a substitute to interactions, and the hydrodynamic description can be working even for very weakly interacting Fermi gas [68].
fields under the Hamiltonian H then take hydrodynamical form: To obtain frequency of small oscillations one can either linearize these equations or perform real-time evolution of the pseudo-Schrödinger equation that appears as a consequence of the inverse Madelung transformation: Both methods generally yield the same results, and we use either of them, depending on the difficulty of performing the calculations in a given case. We confirm that hydrodynamics retrieves the TDHF results for radial monopole in both, in-and out-of-phase, cases. This statement is true for both types of interaction, renormalized and not renormalized. In the case of other modes, three branches were characterized -the upper (growing from ω/ω 0 = 2.0 to ω/ω 0 ≈ 2.2), the middle (decreasing from ω/ω 0 = 2.0 to ω/ω 0 ≈ 1.8), and the lower one (decreasing from ω/ω 0 = 2.0 to ω/ω 0 ≈ 1.2). The first and the last one are retrieved quantitatively with TDHD approach in all cases. However, the middle branch is always characterized by ω/ω 0 ≈ 1.41, independent of the interaction. It is not surprising, as such a solution of linearized hydrodynamical equations (in basic form, without renormalization) is explicitly independent from the coupling constant, yielding ω/ω 0 = √ 2. These results strongly suggest that overlapping regions of fermionic clouds can be described hydrodynamically. It shows that if the hydrodynamical excitation frequency of the noninteracting gas is as it should be (coincidentally, as noninteracting gas is not hydrodynamic), the corrections due to interaction can be evaluated by the means of hydrodynamical description. It proves to be useful, when the TDHF description becomes intractable -namely when there are too many atoms in the system. Therefore, TDHF and TDHD seem to be complementary methods for such a use, TDHF being suitable for small systems, and TDHD being utilized to consider settings closer to experimental setups.

IV. CONCLUSIONS
Recapitulating, we have filled the gap in the literature by analyzing small oscillations of repulsively interacting Fermi-Fermi mixture. Going beyond previous considerations, we investigated the gas on the repulsive energy branch, in contrast to attractive one of superfluid dimers. We did not limit ourselves to the paramagnetic phase in which the density profiles of each species are equal, but also covered the phase-separated state. We calculated radial monopole and two quadrupole modes by the means of time-dependent Hartree-Fock methods. Moreover, we compared these results to the hydrodynamical approach, validating the applicability of the latter to investigation of the dynamics of a weakly interacting Fermi gas. As a future line of work, in the presence of ongoing experiments on quantum mixtures, further investigation of va-lidity of hydrodynamical formalism in such settings can be proposed.