{'fraction_non_alphanumeric': 0.039741426625739436, 'fraction_numerical': 0.013396214908625212, 'mean_word_length': 4.539241076455354, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 10, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Correlated Tunneling Dynamics of an Ultracold Fermi-Fermi Mixture Confined in a Double-Well

We unravel the correlated tunneling dynamics of a mass imbalanced few-body fermi-fermi mixture upon quenching the tilt of a double-well. The non-equilibrium dynamics of both species changes from Rabi-oscillations close to the non-interacting limit to a delayed tunneling dynamics for moderate interspecies repulsions. Considering strong interspecies interactions the lighter species experiences quantum self-trapping due to the heavier species which acts as an effective material barrier, while performing almost perfect Rabi-oscillations. The degree of entanglement, inherent in the system, is analyzed and found to be significant both at moderate and strong repulsions. To relate our findings with possible experimental realizations we simulate in-situ single-shot measurements and discuss how a sampling of such images dictates the observed dynamics. Finally, the dependence of the tunneling behavior on the mass ratio, the particle number in each species and the height of the barrier of the double-well is showcased.

We unravel the correlated tunneling dynamics of a mass imbalanced few-body fermi-fermi mixture upon quenching the tilt of a double-well. The non-equilibrium dynamics of both species changes from Rabi-oscillations close to the non-interacting limit to a delayed tunneling dynamics for moderate interspecies repulsions. Considering strong interspecies interactions the lighter species experiences quantum self-trapping due to the heavier species which acts as an effective material barrier, while performing almost perfect Rabi-oscillations. The degree of entanglement, inherent in the system, is analyzed and found to be significant both at moderate and strong repulsions. To relate our findings with possible experimental realizations we simulate in-situ single-shot measurements and discuss how a sampling of such images dictates the observed dynamics. Finally, the dependence of the tunneling behavior on the mass ratio, the particle number in each species and the height of the barrier of the double-well is showcased.
Regarding the non-equilibrium quantum dynamics, atoms trapped in a double-well potential constitute a prototype system to study the correlated tunneling dynamics in a controllable manner. For bosons such a system represents the analog of the well-known superconducting Josephson junction [44][45][46][47]. The bosonic Josephson junction exhibits various experimentally observed [44] intriguing phenomena, such as Josephson oscillations, π modes, macroscopic quantum self-trapping [48][49][50][51][52][53][54] and correlated pair tunneling [45,[55][56][57][58][59]]. An extension is provided by multicomponent bosonic setups trapped in a double-well evincing for instance coherent quantum tunneling and self-trapping [60], collapse and revival of population dynamics [61,62], symmetry breaking and restoring scenarios [63] and counterflow superfluidity [64]. Turning to the fermionic tunneling properties, the spinpolarized fermionic gas in a double-well has also been intensively studied mainly in three-dimensions unveiling coherent Josephson dynamics of two superfluid samples [65,66], temperature dependent deformations of the density profile [67] and incoherent single-particle oscillations [68]. However, the tunneling dynamics of multicomponent fermionic systems has received much less attention [69][70][71][72][73]. This includes the investigation of the dynamical properties of few atom FF mixtures confined in a double-well potential in relation to the few-body eigenspectrum [69], the necessity of the renormalization of the tunneling frequency in few Fermi polaron systems [70] and the observation of Josephson oscillations in threedimensional setups [72,73]. Besides the above-mentioned studies a systematic study of the FF mixture tunneling dynamics in a one-dimensional double-well covering the weak to strong interspecies correlation regimes has not been reported. In this setup, it would be particularly interesting to examine how the interspecies correlations modify the tunneling properties of the mixture i.e. the Josephson-like oscillations, quantum self-trapping or a more complex motion. Since the fermions of each species are spin-polarized the question of induced intraspecies correlations occurs. Another intriguing prospect is to analyze the interplay between the two fermionic clouds during the evolution and to examine whether one species can act as an effective material barrier for the other one, a result that is already known for bosons [74,75]. To track the non-equilibrium dynamics of the FF mixture we utilize the Multi-Layer Multi-Configurational Time-Dependent Hartree Method for Atomic Mixtures (ML-MCTDHX) [76], being a variational method that enables us to capture all the important particle correlations.
Motivated by the recent few fermion experiments [2][3][4] we investigate the correlated tunneling dynamics of a FF mixture upon quenching an initially tilted double-well to a symmetric one, thus favoring the tunneling of both components to the other well once the relative energy arXiv:1809.00908v1 [quant-ph] 4 Sep 2018 offset vanishes. Inspecting the single-and two-particle probabilities for the fermions of each species to reside in a certain well we unveil the occurrence of three distinct dynamical tunneling regimes with respect to the interspecies repulsion. In particular, close to the noninteracting limit both components perform almost perfect Rabi-oscillations [69] with the heavier species exhibiting either single-or two-particle tunneling and the lighter one three-particle transport. At intermediate interactions the tunneling behavior of both species is significantly altered. The dynamics of the heavier component is characterized by a higher-order quantum superposition while the lighter species undergoes single-particle tunneling. For strong interactions the lighter species experiences quantum self-trapping due to the heavier species which acts as a material barrier [74,75]. The heavier component exhibits almost perfect Rabi-oscillations performing either three-particle or pair tunneling. The degree of both inter-and intraspecies correlations is found to be overall significant especially for intermediate interactions. To provide possible experimental evidences of the observed dynamics we simulate in-situ single-shot measurements and showcase that utilizing a sampling of such images the fermionic tunneling behavior can be retrieved. Finally, the dependence of the tunneling dynamics on the mass ratio of the two components, the particle number and the height of the barrier of the double-well is discussed providing ways to further control it.
This work is structured as follows. In Section II we introduce our setup and many-body treatment as well as the relevant observables. Our results on the FF tunneling dynamics in a double-well are presented in Section III. Section IV provides the simulation of in-situ singleshot images. In Section V we discuss the robustness of the tunneling behavior for different system parameters. We provide a summary of our findings and an outlook in Sec. VI. In Appendix A, the influence of the quench strength on the tunneling dynamics is briefly discussed. Finally, in Appendix B we provide details on the numerical implementation of the single-shot procedure, while in Appendix C we address the convergence of our manybody simulations.

A. Setup
We consider a FF mixture consisting of N A and N B spin polarized fermions and mass ratio M B = 6M A for the species A and B respectively. This mass-imbalanced system can be experimentally realized by considering e.g. a mixture of isotopes of 40 19 K and 6 3 Li [19]. The mixture is confined in an one-dimensional tilted double-well external potential [77,78] which is comprised by a harmonic oscillator possessing a frequency ω and a centered Gaussian with height V 0 and width w as well as a linear tilt with tilt parameter d. Since the two species correspond to different atomic elements they possess distinct polarizations. This means that, experimentally using optical trapping [79,80], they experience different double-well potentials. For convenience here we choose the frequencies of the imposed oscillators to be the same for both species and therefore the two fermionic components experience different double-wells due to their different masses, see Fig. 1 (a) and the discussion below. The corresponding many-body Hamiltonian reads Operating in the ultracold regime s-wave scattering is the dominant interaction process and therefore the interspecies interactions can be adequately modeled by contact interactions that scale with the effective onedimensional coupling strength g AB for the distinct fermionic species. Note here that since s-wave scattering is forbidden for spin-polarized fermions due to the antisymmetry of the fermionic wavefunction [81,82], the fermions of the same species are considered to be non-interacting and thus only interspecies interactions are involved. The effective interspecies onedimensional coupling strength [83] is given by where ζ is the Riemann zeta function and µ = M A M B M A +M B denotes the corresponding reduced mass. The transversal length scale is a ⊥ = /µω ⊥ , where ω ⊥ is the frequency of the transversal confinement, and a s AB refers to the threedimensional s-wave scattering length between the two species. Experimentally g AB is tunable either by a s AB via Feshbach resonances [6,84] or by ω ⊥ and the resulting confinement-induced resonances [83,85].
In the following we shall rescale our Hamiltonian in units of ω ⊥ . Then, the corresponding length, time, and interaction strength scales are given in units of M A ω ⊥ , ω −1 ⊥ and 3 ω ⊥ M A respectively. Accordingly, the amplitude of the Gaussian barrier V 0 , its width w, the tilt parameter d and the frequency of the harmonic oscillator ω are expressed in terms of and ω ⊥ . Finally, in order to limit the spatial extension of our system we impose hard-wall boundary conditions at x ± = ±40.
Our system is initially prepared in the many-body ground state of the tilted double-well with a complete population imbalance, namely all fermions of both species reside in their corresponding left well (−40 < x < 0), see also Fig. 1 (a). Note that throughout this work, unless it is stated otherwise, we use V 0 = 1 and w = 0.1 thus having two doublets below the maximum of the barrier. The left well is, of course, significantly energetically favorable due to the presence of the linear external tilt dx for d sufficiently large. As an illustration Fig. 2 shows the ground state population of the left well of each species [see also Eq. (4)] for varying tilt magnitude d for the cases N σ = 5 (N A = N B ) when g AB = 0.1 and for N σ = 3 with g AB = 0.1 and g AB = 4.0. We observe that the population imbalance depends on both N σ and g AB . However, for d = 0.2 each species is fully localized in its left well independently of the interspecies interaction and the particle number. This independence of g AB is caused by the fact that for d = 0.2 the distinct fermionic clouds are non-overlaping, see also Fig. 1 (a). To ensure that both species reside initially in their left well we therefore use d = 0.2. To induce the dynamics in the symmetric double-well we quench at t = 0 the asymmetry to d = 0 and let the system evolve in time. Quenching the tilt to zero favors the tunneling of both components to the right well as the corresponding energy offset between the two distinct wells vanishes, see also Fig. 1 (b). It is worth mentioning that the case of different harmonic oscillator frequencies for each species does not yield fundamentally different tunneling phenomena but rather results in distinct tunneling frequencies of the two species and renders the tunneling regions to be presented below wider with respect to the corresponding interspecies interaction. We have checked this for the experimentally relevant [80] values ω B = 0.6ω A and M B = (40/7)M A (results not shown here for brevity). Note that the tunneling dynamics can also be induced for a final asymmetry, i.e. finite values of d, where the left well is energetically favorable (see Appendix A).

B. Many-Body Wavefunction Ansatz
To examine the quench induced tunneling dynamics of the FF mixture within the double-well we resort to ML-MCTDHX [76]. Within this approach, the many-body wavefunction is expanded with respect to a time-dependent and variationally optimized basis, allowing us to take into account both the inter and the intraspecies correlations. To incorporate inter and intraspecies correlations M distinct species functions, Nσ being the spatial σ = A, B species coordinates, for each component consisting of N σ fermions are firstly introduced. It holds that M ≤ min(dim(H A ), dim(H B )) with H σ being the Hilbert space of the σ-species. Accordingly, the manybody wavefunction Ψ M B is expressed as a truncated Schmidt decomposition [86] of rank M The Schmidt coefficients λ k (t) in decreasing order are referred to as the natural species populations of the k-th species function Ψ σ k of the σ-species and provide a measure of the system's entanglement (see also Sec. II C). Indeed, the system is termed entangled [87] or interspecies correlated when at least two distinct λ k (t) are nonzero, thus preventing the total many-body state [Eq. (2)] to be a direct product of two states.
Moreover in order to include interparticle correlations each of the species functions Ψ σ k ( x σ ; t) is expanded using the determinants of m σ distinct time-dependent fermionic single-particle functions (SPFs), ϕ 1 , . . . , ϕ mσ , namely In the latter expression, P denotes the permutation operator exchanging the particle configuration within the SPFs and sign(P i ) refers to the sign of the corresponding permutation. Also, C k,(n1,...,nm σ ) (t) are the timedependent expansion coefficients of a particular determinant and n i (t) is the occupation number of the SPF ϕ i ( x; t). In the following, we shall define that each fermionic species possesses intraspecies correlations if more than N σ SPFs are substantially occupied, otherwise the state of a species reduces to the Hartree-Fock ansatz [81,88,89]. Following e.g. the Dirac-Frenkel variational principle [90,91] for the many-body ansatz [see Eqs.
(2), (3)] yields the ML-MCTDHX equations of motion [76] for a binary fermionic mixture. These consist of M 2 linear differential equations of motion for the coefficients λ i (t), which are coupled to a set of M [ m A N A + m B N B ] non-linear integro-differential equations for the species functions and m A + m B integro-differential equations for the SPFs. Finally, let us mention that ML-MCTDHX is able to operate within different approximation orders, e.g. it reduces to the Hartree-Fock equation for M = 1 and m σ = N σ .

C. Observables of Interest
In this section we briefly introduce the main observables that will be subsequently employed for the interpretation of the tunneling dynamics. On the one-body level the tunneling dynamics can be examined by inspecting the population imbalance between the left and right wells during the time evolution. To this end, we measure the expectation value of the one-body density [56] e.g. in the left well Here, ρ σ (x; t) = Ψ M B (t)|Ψ σ † (x)Ψ σ (x)|Ψ M B (t) denotes the σ-species one-body density, which is normalized to the corresponding particle number N σ . Ψ σ † (x) [Ψ σ (x)] is the fermionic field operator that creates (annihilates) a σ-species fermion at position x. As it is evident from Eq. (4), ρ (1) σ (t) L essentially measures how many σ-species particles reside in the left well and it is normalized to the corresponding number of N σ particles such that ρ (1) To unveil the two-body correlation mechanisms that are responsible for the observed tunneling dynamics we employ the pair probability p (2) σσ [56,92] in the course of the evolution. It is based on the diagonal two-body reduced density matrix ρ which refers to the probability of measuring two fermions of the same (σ = σ ) or different species (σ = σ ) located at positions x 1 , x 2 at time t respectively [93]. In our case p (2) σσ provides the probability to find two fermions of the same or different species within the same well and it is given by The normalization used corresponds to σσ can be seen as a measure for the number of pairs that can be found in both wells and it is normalized to unity.
To expose the degree of both inter-and intraspecies correlations during the time-evolution for increasing interspecies repulsions we measure the entanglement and the fragmentation of the FF mixture [34,94,95].
In this way, also the corresponding departure from a Hartree-Fock state can be deduced. To quantify the presence of interspecies correlations or entanglement we calculate the eigenvalues λ k of the species reduced density matrix ρ Nσ ( , and σ = σ [see also Eq. (2)]. It is known that when multiple eigenvalues of ρ Nσ are macroscopically large the system is referred to as species entangled or interspecies correlated, otherwise it is said to be non-entangled. A well-known measure to quantify the degree of the system's entanglement is the Von-Neumann entropy [96,97], S, which is based on the natural populations λ i (t) of the species functions [see Eq. (2)] Indeed within the Hartree-Fock, non-entangled, limit S(t) = 0 since λ 1 (t) = 1 while for a beyond Hartree-Fock state where more than a single λ i is nonzero S(t) = 0.
To reveal the fragmented nature of each species we rely on the eigenvalues of the one-body reduced density matrix of the σ-species ρ Schematic representation of the quench protocol consisting of (a) an initially tilted double-well with V tilt = dx > 0 to (b) a symmetric one with V tilt = d = 0. The quench leads to the tunneling dynamics of the two fermionic species being initially localized in the left well of their corresponding double-well. In the sketch the double-wells V A DW , V B DW for each species are shown as well as the initial tilt V tilt (see legend). [98]. The eigenfunctions of ρ (1) σ (x, x ), are the so-called σ-species natural orbitals, φ σ i (x; t), which we consider to be normalized here to their corresponding eigenvalues n σ It can be shown that when = 0 and the corresponding Hartree-Fock wavefunction is retrieved. Therefore, if more than N σ natural orbitals are occupied, the system is said to be fragmented and the corresponding degree of fragmentation can be quantified via This quantity serves as a theoretical measure for the occupation of the m σ − N σ least occupied natural orbitals and thus for the deviation from a Hartree-Fock state when F σ > 0.

III. QUENCH INDUCED TUNNELING DYNAMICS
We consider a mass imbalanced repulsively interacting, g AB , FF mixture consisting of N A = N B = 3 fermions with M A = 1 and M B = 6. The system is initially prepared in the ground state of the tilted double-well described in Eq. (1) with tilt parameter d = 0.2, frequency ω = 0.1, barrier height V 0 = 1 and width w = 1.
For these values and since M B > M A the heavier Bspecies experiences a much more localized double-well around x = 0 when compared to the double-well of the lighter A-species. Moreover due to the tilt, each species is found to be fully localized in its respective double-well, see Fig. 1 (a) for a schematic respresentation. Then the distinct fermionic clouds are non-overlaping for t = 0, since d = 0.2 see Fig. 2, and therefore phase separated at x < 0. Concluding, this non-overlaping behavior between the two initial fermionic clouds is caused by their mass imbalance [13,39] and renders the ground state to be essentially independent of g AB , see Fig. 2. To examine the tunneling dynamics of the mass imbalanced FF mixture we perform a quench of the initially titled double-well with d = 0.2 to a symmetric one i.e. d = 0, see for instance Fig. 1 (b), keeping fixed all other parameters for a specific g AB and covering the regime from weak to strong interactions namely g AB ∈ [0, 5].

A. One-and Two-Body Tunneling Probabilities
To monitor the overall tunneling dynamics on the onebody level for different values of g AB we rely on the expectation value of the single-particle density within the left well ρ  σ (t) L essentially provides the percentage of the σspecies fermions within the left well or alternatively speaking the population imbalance between the two-sides of the double-well and offers therefore an adequate measure for the tunneling dynamics on the one-body level. Our primary aim is to investigate whether distinct tunneling regions for the fermions can be observed for a varying interspecies repulsion [69]. To conclude upon a certain tunneling regime we study ρ (1) σ (t) L as well as its spectrum, P    where almost all fermions of each species oscillate back and forth between their left-and right wells [69,94]. No-tice that ρ (1) σ (t) L takes values between 0 and 3, see also Table I. For increasing repulsion, 0.2 < g AB < 3, we enter region II where the tunneling process is modified with respect to region I, especially for the A-species. In particular the tunneling oscillations of the A-and B-species are characterized by one and two distinct frequencies respectively, see Figs. 3 (c) and (d). Regarding the tunneling behavior of the fermions this region can also be seen as a transition region between the weak and strong interaction regimes. For even stronger interactions 3 < g AB < 5, we realize region III where the A-species remains mainly in the left well, while the B-species still exhibits a strong amplitude two-frequency [see Figs. 3 (c), (d)] tunneling dynamics. As we shall argue below this latter interplay between the two species is caused by the fact that in this strongly interacting regime the heavier B-species acts as a material barrier for the lighter A-species, thus supressing the tunneling of the latter (see also Sec. III C). Note also that for the A-species this region III is reminiscent to the few-body analogue of quantum self-trapping exhibited for strongly interacting bosons trapped in a doublewell [45,[55][56][57][58]. We actually observe a small amplitude tunneling for times t > 800 (not shown here), while the observed tunneling amplitude of the B-species is even higher than in region I.
To gain deeper insights into the tunneling motion within and between the species and in order to reveal the underlying correlation mechanisms we invoke the pair probability p  σ (t) L enables us to identify the dominant particle configurations |N L , N R σ occuring in the course of the dynamics. Here, N L (N R ) refers to the number of σ-species fermions in the left (right) well. We remark that due to the Pauli exclusion principle [88,89] different fermions of the same species that reside in the same well populate energetically consecutive single-particle bands e.g. for |2, 1 A ≡ |1 0 ⊗ 1 1 , 1 0 where the upper index refers to the band number. The same holds for more than two fermions of different species that are in the same well. Referring to a certain time instant t = t 1 when e.g. p It is also worth stressing at this point that in order to systematically conclude upon certain number state configurations one needs to rely on a projection of the numerically obtained many-body wavefunction to a multiband Wannier number state basis as it has been demonstrated in [99][100][101][102][103] e.g. for the single component bosonic case. However, such a numerical analysis in our FF mixture case is computationally challenging, since we need to take into account a large sample of number states in order to form a complete basis. Below we discuss the behavior of p (2) σσ (t) within each of the above-mentioned tunneling regions I, II and III.
Regarding the A-species we observe that for both      weak (region I) and strong interspecies interactions (region III) all three fermions reside in the same well since p  Table I. The origin of the latter mechanism, in each region, can be better understood in combination with the corresponding behavior of ρ (1) Fig. 3 (a)]. In particular, within the weak interaction regime (region I) all three particles perform a simultaneous tunneling from the left to the right-hand well and vice versa, since ρ A (t) L oscillates between 0 and 3. Then, the system predominantly oscillates between the configurations |3, 0 A and |0, 3 A . However, for strong interactions (region III) the A-species fermions remain trapped in the left well as 2.9 < ρ   A (t) L < 2.7. As a consequence, at least one fermion tunnels from the left to the right well or better it is delocalized over both wells, see also Sec. III C. This indicates that the states |2, 1 A and |1, 2 A dominantly contribute.
Next, we turn our attention to the dynamics of the B-species and show p     Fig. 3 (b)]. Thus, the B-species cloud mainly tunnels between the |1, 2 B state (e.g. at t ≈ 30) and the superposition |2, 1 B + |3, 0 B (e.g. at t ≈ 60) resulting in two and single particle tunneling respectively. In region II, p  Region I Region II Region III σσ (t)], the Von-Neuman entropy S(t) and the degree of fragmentantion Fσ(t) for each species.
nates during the time evolution. Indeed, at the initial stages of the dynamics (10 < t < 50) all three B-species fermions tunnel to the right-hand well, see in particular that p B (t ≈ 100) L ≈ 0). In this way we can infer that the fermions in this regime of interactions tunnel as pairs [99,100,104,105].
To understand further the correlated dynamics within the different tunneling regions we additionally analyze the interspecies pair probability p  Table I. Indeed within the region I p   AB (t ≈ 60) ≈ 0.85. The latter indicate that here the most significant number state configurations of the system correspond to the |0, 3 A ⊗|1, 2 B and |3, 0 A ⊗ |2, 1 B + |3, 0 A ⊗ |3, 0 B respectively. Turning to region II, 0.2 < p (2) AB (t) < 0.7 and therefore a clear assignment in terms of number states can not be unavoided. However this suggests, as it has already been observed in p  BB (t), the participation of a multitude of number states and as a result the occurrence of a multimode dynamics. Is is also worth stressing again that this interaction regime II exhibits a strongly alternating behavior with respect to g AB , e.g. see ρ σσ (t) for 0.2 < g AB < 1.7 and 1.7 < g AB < 3.0 in Fig. 3. Moreover, as we shall argue below the dynamics in this region is characterized by an enhanced degree of both entanglement and fragmentation, see Sec. III B. Finally, for strong interactions (region III) p

B. Entanglement and Fragmentation
To expose the many-body nature of the dynamics of the FF mixture we next measure the degree of entanglement (or interspecies correlations) and the fragmentation (or intraspecies correlations) of the underlying manybody state. These features can be quantified by employing S(t) and F σ (t) respectively for varying g AB , see Fig. 4 and Table I. We remind here that when S(t) = 0 the system is called non-entangled, while S(t) = 0 signifies the presence of entanglement or interspecies correlations. Accordingly, for F σ (t) > 0 the σ-species is termed intraspecies correlated (see also Sec. II C). Most importantly, since the fermions of each species are spinpolarized, i.e. non-interacting, the occurrence of intraspecies correlations during the evolution is induced by the interspecies correlations. Recall also that in the initial (ground) state of the system the distinct fermionic clouds that reside in the left well, see Fig. 1 (a), are non-overlaping and therefore fully phase separated due to their mass imbalance [13,39]. Due to this non-overlaping behavior no inter-and intraspecies correlations occur for arbitrary g AB . As it can be easily seen in Figs. 4 (a), (b), (c) for weak interactions (region I) there is only a small amount of inter-and intraspecies correlations since both S(t) and F σ (t) are suppressed. Strikingly enough, for intermediate interactions (region II) S(t) as well as F σ (t) increase during the evolution and tend to saturate to a certain finite value, signifying that the many-body state is strongly entangled and fragmented. Inspecting region III we observe again that both S(t) and F σ (t) increase in the course of time but in a slower manner as compared to region II. Also they acquire moderate magnitudes with respect to region II. Interestingly enough within this strongly interacting regime III, F B (t) (heavier species) increases more rapidly during the evolution when compared to F A (t) (lighter species).

C. Single-Particle Density Evolution
To visualize the tunneling motion in the different interaction regimes identified above, we finally invoke the evolution of the single-particle density, ρ oscillates as a whole back and forth resulting in a threebody tunneling process. Only for large evolution times, t > 100, a small fraction of the A-species gets reflected from the barrier, see the dashed rectangle in Fig. 5 (a). For the B-species a substantial reflection already occurs during the first tunneling process and a corresponding fraction remains trapped in the left well. Consequently a delocalization of ρ One of these fragments travels back to the barrier and the other one continues towards the outer part of the right well. The effect of this collisional dynamics is, of course, much more pronounced for the A-species since it is the lighter one [76] and it is manifested by the delocalization of both species over the wells. This delocalization becomes even more apparent for large evolution times, t > 100, and it is a consequence of the fact that the distinct species suffer multiple collisions with one another, rendering the one-body density blurred. Recall also that the dynamics in this regime of interactions is strongly entangled and fragmented, see Fig. 4. Let us also comment at this point that within the Hartree-Fock approximation the above-mentioned delocalization of both species during the evolution becomes suppressed (not shown here for brevity) since it is essentially caused by the contribution of the higher-lying orbitals. Turning to the strongly interacting regime (region III), we observe that the heavier B-species acts as a hard-wall or material barrier for the lighter A-species [74,75], as shown in Figs. 5 (e), (f). In this way, the A-species remains fully trapped within the left well during the evolution, while the B-species still performs a tunneling motion. The only effect of the collision between the two species that is imprinted in the density of the B-species is a splitting of its density into two parts, one of which undergoes tunneling between the two wells and the other remains trapped in the right well where it performs tiny amplitude oscillations.

IV. SINGLE-SHOT SIMULATIONS
To provide possible experimental evidence for the many-body tunneling dynamics of the FF mixture we simulate in-situ single-shot absorption measurements [93,106,107]. This type of measurements probe the spatial configuration of the atoms which can be inferred by the  many-body probability distribution. Such experimental images refer to a convolution of the spatial particle configuration with a point spread function which essentially dictates the corresponding experimental resolution. For our calculations, to be presented below, we use a point spread function of Gaussian shape possessing a width w P SF = 1 l ≈ 3.2, where l = 1/ω denotes the corresponding harmonic oscillator length. It is also important to note here that in few-body experiments [3,4] besides the in-situ imaging of the cloud another technique to probe the state of the system and get rid of unavoidable noise sources that might destroy the experimental signal is fluorescence imaging [2,34]. However, the simulation of this technique which has been proven important for few-body experiments lies beyond our current scope.
Here we aim at demonstrating how in-situ imaging can be used to adequately monitor the above-described fewbody FF mixture dynamics.
Relying on the system's many-body wavefunction obtained within ML-MCTDHX we simulate in-situ single shot images for both species A, A A (x; t), and species B, A B (x |A A (x); t im ), at each instant t of the evolution when we consecutively image first the A and then the B species. For details regarding the simulation process of this experimental technique we refer to Appendix B. Below we focus on the FF dynamics, possessing g AB = 1.2 (region II), within the double-well upon quenching the tilt parameter from d = 0.2 to d = 0. Note that a similar procedure has been followed also for other values of g AB (not shown here). Figures 6 (a), (b) show the first simulated in-situ single-shot images for both species during the evolution, namely A A (x; t), and A B (x |A A (x); t) respectively. It can be deduced that the two species exhibit a tunneling behavior, resembling this way the overall tendency observed in the one-body density [see also Figs. 5 (c), (d)]. However, it is important to mention that a direct correspondence between the one-body density and one single-shot image is not possible due to the small particle number of the considered FF mixture, N A = N B = 3. Another source of the inability to explicitly observe the one-body density within a single-shot image is the presence of multiple orbitals in the system. Indeed, the many-body state builds upon a superposition of multiple orbitals [see Eqs. (2) and (3)] and therefore imaging an atom alters the many-body state of the other atoms and hence their one-body density. For a more elaborated discussion on this topic see [93,108,109]. To retrieve the one-body density of the system we next rely on an average of several single-shot images for each species, namelyĀ A (x; t) = 1/N shots

V. FURTHER CHARACTERISTICS OF THE TUNNELING DYNAMICS
Having discussed in detail the properties of the tunneling dynamics of the FF mixture in the double-well, we finally showcase how the overall dynamics depends on certain system parameters. To this end, we first study    the effect of the barrier height V 0 on the quench induced tunneling dynamics. As a reference system we consider the double-well with V 0 = 1 that has already been discussed above, while as an indicator for the tunneling dynamics we employ ρ   B (t) L respectively for the case V 0 = 4 for varying interspecies repulsion g AB . We observe that the tunneling dynamics of both species, and especially the B-species, is strongly influenced by the barrier height and in particular it is overall suppressed [69,95,101], compare Figs. 3 (a), (b) and Figs. 7 (a), (b). Regarding the A-species the three above-identified tunneling regions are shifted to weaker interactions and in particular region I becomes negligible in size and exists only very close to the non-interacting limit [hardly visible in Fig. 7 (a)]. Therefore, an increasing barrier height shows a similar effect on the A-species as a stronger interaction strength where the B-species acts as a material barrier [74,75], see for instance Fig. 5 (e). For the heavier B-species the barrier height has the most prominent effect as the tunneling dynamics is significantly suppressed for every value of g AB when compared to the V 0 = 1 case where tunneling oscillations occur independently of the interaction strength. Here the tunneling takes place for g AB > 1 and it exhibits a small amplitude.
As a next step, we examine the influence of the mass ratio on the tunneling behavior of the FF mixture by inspecting again ρ        [74,75] for the lighter A-species and its collision with the B-species pushes the latter to the right well, enforcing its tunneling motion.
To shed light on the particle number dependence of the tunneling motion, we finally utilize the normalized expectation value ρ (1) σ (t) L /N σ for a varying particle number N σ ∈ {2, 3, 5} within the different interaction regimes, see Figs. 7 (i)-(n). Overall, we observe that for all three interaction regimes (I, II, III) an increasing particle number results in a stronger damping of ρ  σ (t) L /N σ oscillation, remain essentially insensitive. An inspection of the two-body tunneling characteristics unveils more pronounced differences since for larger particle numbers more modes are triggered which are imprinted e.g. in the p (2) σσ evolution (results not shown here).

VI. CONCLUSIONS
We have investigated the tunneling dynamics of a FF mixture with spin polarized components confined in a double-well. The fermionic mixture is initially prepared within the left well of a tilted double-well and the tunneling dynamics is induced by removing the tilt and let the system evolve in a symmetric double-well. We particularly examined the impact of the interspecies interactions on the tunneling behavior of each species and unveiled the significant role of intra-and interspecies correlations. The emergent dynamics has been characterized on both the one-and two-body level by invoking the single-and two-particle probabilities for the fermions of each species to reside e.g. in the left-hand well in the course of the dynamics.
For very weak interactions close to the noninteracting limit both species perform almost perfect Rabi-oscillations. The heavier species exhibit small deviations from the perfect Rabi scenario since they either show single or two-particle tunneling i.e. during the evolution not all three particles tunnel together. Increasing the interspecies repulsion the tunneling dynamics is suppressed for both species. The lighter species undergoes single-particle tunneling, while for the heavier one a more complex dynamics takes place which is characterized by a higher-order quantum superposition further indicating the presence of strong correlations in the system. Turning to strong interactions we observe that the lighter species experiences a quantum self-trapping due to the heavier species which acts as a material barrier. In particular, the collision between the two species causes a splitting of the density of the heavier component into two parts. The first part performs tunneling between the two wells and the other one remains trapped in the right well ex-hibiting small amplitude oscillations. As a consequence, this heavier component undergoes almost perfect Rabioscillations being characterized either by three-particle or pair tunneling at different instants of the evolution. The degree of both inter-and intraspecies correlations for all interaction regimes has been exposed and found to be overall significant especially for intermediate interactions. To relate our findings with possible experimental realizations we simulate in-situ single-shot measurements and showcase how an increasing sampling of such images can be used to adequately reproduce the observed fermionic tunneling dynamics.
The dependence of the observed tunneling behavior on the mass ratio of the two components, the particle number and the height of the barrier of the double-well has been discussed. It is shown that the mass imbalance between the components possesses a strong influence on the tunneling dynamics depending on the interspecies interaction strength. Namely, for weaker interactions the tunneling amplitude of the heavier species becomes smaller for an increasing mass ratio, while the tunneling features of the lighter species remain essentially unaffected. Increasing the repulsion a decrease of the tunneling amplitude for both species occurs during evolution. For strong interactions the heavier species acts as a material barrier for the lighter one thus suppressing the tunneling motion of the latter. On the other hand, for fixed interspecies repulsion and larger particle numbers a damping of the tunneling oscillations takes place especially for the lighter species. Finally, we show that for a higher barrier the tunneling dynamics of both species, and in particular of the heavier one, is overall suppressed since an increasing barrier disfavors the tunneling process.
There is a multitude of interesting directions worth pursuing in future studies. A straightforward one would be to examine the tunneling dynamics of a dipolar FF mixture, under the quench protocol considered herein, and investigate how the long-range character of the interactions alters the emergent tunneling behavior including the inherent entanglement generation. Another interesting prospect is to unravel the collisional dynamics of two fermi components which are initially well separated, e.g. due to the presence of a high barrier and then are left to collide by removing this intermediate barrier in a sense of the counterflow experiment. Inspecting the many-body character of the spontaneously generated non-linear excitations such as dark solitons [110] or domain-wall structures [108,109] between the two species is a challenging future task. Finally, the generalization of our current findings utilizing time-dependent quench protocols would be desirable in order to design schemes for selective transport between the wells of each component.    σ (t) L occur in the sense that the oscillation frequency and amplitude decrease for a larger offset. Finally, for d > 0.09 both species remain in the left well without tunneling throughout the evolution as the energy offset is adequately large and suppresses any possible tunneling process. We remark that we have performed the same investigation for both intermediate and strong interspecies interactions g AB , observing a similar behavior of the emerging dynamics within the same postquench intervals of the tilt magnitude d as above (not shown here for brevity).

Appendix B: Single-Shots Implementation in Fermi-Fermi Mixtures
To simulate the experimental single-shot procedure we perform a sampling of the many-body probability distribution [93,106,108,109] being accessible within the ML-MCTDHX framework. We remark that such a theoretical implementation of the experimental process has already been performed for single-species bosons [106][107][108][109] and fermions [34] as well as for binary bosonic ensembles [93] but not yet for FF mixtures. As in the twospecies bosonic case, the corresponding single-shot procedure for FF mixtures depends strongly on the interand intraspecies correlations inherent in the system. Indeed within a many-body treatment the presence of entanglement between the different species is of significant importance regarding the image ordering. This can be understood by inspecting the Schmidt decomposition [see Eq. (2)] and especially the involved Schmidt coefficients λ k 's. Below let us elaborate on the corresponding process where we image first the A and then the B-species, obtaining in this way the images A A (x) and A B (x |A A (x)). However, we note that the same overall procedure has to be followed in order to image first the B and then the A species, resulting in single-shots images A B (x) and A A (x |A B (x)) respectively.
To image first the A and then the B species we consecutively annihilate all the N A fermions. Referring to a specific time instant of the imaging process, e.g. t im , a random position is drawn that satisfies the constraint ρ (1) N A (x 1 ) > q 1 , where q 1 denotes a random number lying in the interval [0, max{ρ (1) N A (x; t im )}]. Next, we project the (N A + N B )-body wavefunction onto the (N A − 1 + N B )-body one. The latter is accomplished by the use of the projection operator 1 N (Ψ A (x 1 )⊗Î B ), wherê Ψ A (x 1 ) refers to the fermionic field operator that annihilates an A species fermion located at x 1 and N is the normalization constant. An important observation, here, is that the above process affects the Schmidt weights, λ k , and therefore besides that the B-species has not already been imaged, both ρ (1) To comprehend the latter we rely on the Schmidt decomposition according to which the many-body wavefunction after the first measurement reads In this expression, i refer to the Schmidt coefficients of the (N A − 1 + N B )-body wavefunction. Repeating the above steps N A − 1 times we obtain the following distribution of positions (x 1 , x 2 ,...,x N A −1 ) which is then convoluted with a point spread function. This for the A-species, wherex are the spatial coordinates within the image and w P SF is the width of the employed point spread function. It is also important to mention at this point that after annihilating all A-species fermions the many-body wavefunction reads where x|Φ A j,1 refers to the single-particle orbital of the j-th mode and the wavefunction of the B-species, |Ψ N B M B (t im ) , corresponds to the second term of the cross product on the right-hand side. Evidently, |Ψ N B M B (t im ) obtained after the annihilation of all N A fermions corresponds to a non-entangled N B -particle many-body wavefunction and its corresponding single-shot procedure reduces to that of the single-species case [106,108,109]. The latter has been extensively used in a variety of settings (see e.g. [106,108,109]) and therefore it is only briefly outlined below. For a time instant of the imaging t = t im we calculate ρ (1) N B (x; t im ) from the manybody wavefunction |Ψ N B ≡ |Ψ(t im ) . Next, a random position x 1 is drawn according to the constraint ρ Let us finally briefly discuss the basic aspects of our numerical method and then showcase the convergence of our results. ML-MCTDHX [76] consists of a variational method for solving the time-dependent manybody Schrödinger equation of atomic mixtures with constituents being either of Bose [93,111] or Fermi [34,94] type. To construct the many-body wavefunction a timedependent variationally optimized many-body basis is employed, enabling us to take into account the system's important correlation effects using a computationally feasible basis size. Therefore the system's relevant subspace of the Hilbert space is spanned at each time instant of the evolution in a more efficient manner when compared to expansions relying on a time-independent basis as the number of basis states can be significantly reduced. Finally, as a result of the multi-layer ansatz for the total wavefunction we are able to capture both the intra-and the interspecies correlations emerging during the nonequilibrium dynamics of a bipartite system.
Within our approach the used Hilbert space truncation, namely the order of the considered approxima- tion, corresponds to the employed numerical configuration space being denoted by C = M − (m A , m B ). In this notation M = M A = M B and m A , m B refer to the number of species and single-particle functions respectively for each of the species. We also note that for our simulations, a primitive basis consisting of a sine discrete variable representation is employed which involves 400 grid points. To infer about the convergence of our many-body simulations we systematically check that upon varying the numerical configuration space C = M − (m A , m B ) the observables of interest remain insensitive. We remark that all many-body calculations discussed in the main text rely on the numerical configuration C = 9−(9, 9) for N A = N B = 3 and C = 10−(10; 10) when N A = N B = 5. In this way, the available Hilbert space for the simulation includes 8037 (9140) coefficients when N σ = 3 (N σ = 5). This is in sharp contrast to an exact diagonalization procedure which should take into account 10586800 (83219 10 10 ) coefficients for the N σ = 3 (N σ = 5) case, rendering these simulations infeasible. To conclude, let us briefly showcase the convergence of our results upon varying the number of species functions and single-particle functions. To perform this investigation we employ the expectation value of the one-body density in the left well, ρ (1) σ;C (t) L , and calculate its absolute deviation during the time evolution, for each of the species, between the C = (9; 9; 9) and other numerical configurations where N σ denotes the particle number of the σ-species.
As it is evident from Eq. (C1) ∆ ρ (1) σ;C,C (t) L is normalized to unity. Fig. 9 (a) [(b)] shows ∆ ρ (1) A;C,C (t) L [∆ ρ (1) B;C,C (t) L ] following a quench of the tilt parameter from d = 0.2 to d = 0 for g AB = 1.2. For reasons of completeness we remark that this value of g AB lies within the region that the tunneling dynamics is characterized by strong inter-and intraspecies correlations, see also Fig. 4. We indeed observe an adequate convergence of both ∆ ρ (1) A;C,C (t) L and ∆ ρ (1) B;C,C (t) L . In particular, comparing the C = 6 − (6, 6) and C = 9 − (9, 9) approximations, the corresponding relative difference, for both species, is less than 5% throughout the evolution. The same observations can also be made for the case of strong interspecies interactions as illustrated in Figs. 9 (c) and (d) for species A and B respectively. Finally, we remark that a similar analysis has been performed for all other observables used in the main text and found to be ade-quately converged, i.e. their relative deviations between the C = 9 − (9, 9) and the C = 6 − (6, 6) configurations lie below 6% (not shown here for brevity).