Phase separations induced by a trapping potential in one-dimensional fermionic systems as a source of core-shell structures

Ultracold fermionic gases in optical lattices give a great opportunity for creating different types of novel states. One of them is phase separation induced by a trapping potential between different types of superfluid phases. The core-shell structures, occurring in systems with a trapping potential, are a good example of such separations. The types and the sequences of phases which emerge in such structures can depend on spin-imbalance, shape of the trap and on-site interaction strength. In this work, we investigate the properties of such structures within an attractive Fermi gas loaded in the optical lattice, in the presence of the trapping potential and their relations to the phase diagram of the homogeneous system. Moreover, we show how external and internal parameters of the system and parameters of the trap influence their properties. In particular, we show a possible occurrence of the core-shell structure in a system with a harmonic trap, containing the BCS and FFLO states. Additionally, we find a spatial separation of two superfuild states in the system, one in the BCS limit as well as the other one in the tightly bound local pairs (BEC) regime.


Model and Method
In this paper, we study a one-dimensional system with the s-wave superconductivity. The attractive Hubbard Hamiltonian (i.e., with on-site pairing, < U 0) in presence of the magnetic field (h) has a following form: where σ † c i (c iσ ) denotes the creation (annihilation) operator of the particles at site i and spin { , } σ = ↑ ↓ . t is the hopping integral between nearest-neighbor sites, < U 0 is the on-site pairing interaction, h is the external magnetic Zeeman field, while μ i is an effective on-site chemical potential. The interaction term is treated within the mean-field broken-symmetry approximation: i n in n i n n where γ nσ and γ σ † n are the new quasi-particle fermionic operators, whereas u and v are the Bogoliubov-de Gennes (BdG) eigenvectors. One can obtain the self-consistent BdG equations in real space: is the single-particle Hamiltonian and δ Δ = Δ U ij i ij are the on-site SOPs. Using transformation (4), the SOPs can be found as: i i i n in in n i n in n Equations (5) can be solved self-consistently with respect to the distribution of Δ i . In this case, one can find the grand canonical potential for a given state as: www.nature.com/scientificreports www.nature.com/scientificreports/ From several solutions of the BdG equations, only those with a minimal value of grand canonical potential Ω (at fixed μ and h) indicate a thermodynamically stable state of the system. In the absence of a trap (i.e., μ μ = i at each site), the distribution of Δ i can be rewritten in momentum space, using the Fourier transform: where q, which are restricted to the first Brillouin zone, are total momenta of the Cooper pairs.

Numerical Results
In this section, the obtained numerical results are presented. First, the homogeneous system (without the trap) is investigated. The magnetic field h versus chemical potential μ phase diagram is determined. Next, we show that phase separation occurs in the presence of a harmonic trap. Such a phase separation induced by a trapping potential we will call (artificially) enforced phase separation further in the text. For chosen parameters (i.e., μ and trapping potential), the phase separation between different phases (e.g., a state with the BCS core and the FFLO shell as well as the FFLO core and the BCS shell) can be realized in the system. phase diagram: homogeneous system. Numerical calculations presented in this section have been performed for a one-dimensional chain with the periodic boundary conditions and = N 200 sites. For the homogeneous system, i.e., μ μ = = . const i , the magnetic field h versus chemical potential μ phase diagram is shown in Fig. 1. The results for different values of U/t are presented (cf. also refs 24 and 49 ). In each case, they consist of three regions. For low values of magnetic field, the conventional BCS phase is stable (with Δ = const i and Δ ≠ 0 i ). With increasing h, one finds a discontinuous phase transition from the BCS to the FFLO phase (with Δ i changing from site to site). In the third region the normal (non-ordered, NO) phase is stable (with Δ = 0 i ). The transition from the BCS to the NO phase is continuous (for changing h and fixed μ) or discontinuous (at the vertical boundary in the phase diagram). Notice also that for small μ/t and large h/t, the so-called η phase can exist 50 (the FFLO phase with maximal q = π, at least for the FF ansatz 24,49 ). With increasing U/t, the region of the BCS phase extends, whereas the continuous FFLO-NO boundary weakly depends on U/t. In Fig. 1(a), the solid lines correspond to the case in which only the FF state is considered, i.e., all Cooper pairs have the same momentum q 25 . Including the fact that generally in the "full" FFLO state, Cooper pairs with different momenta can contribute (i.e., given by Eq. (6), one gets that the region of the FFLO phase occurrence is slightly extended in comparison to the one obtained only for the FF ansatz. This is an expected and well known relation between superconducting phases with different numbers of allowed q, i.e., a phase with larger number of q's is more stable than the FF phase with only one q 27,51 .
It should be also mentioned that for low filling (i.e., for |μ| ~ 2t), the BCS-BEC crossover can be realized. The detailed discussion of this issue has been presented in ref. 24 . In the context of the present work (particularly for the problem of the phase separation region in the system), the following property of the phase diagram is also very important. Namely, depending on values of the magnetic field h and on-site attraction U, with increasing |μ| the www.nature.com/scientificreports www.nature.com/scientificreports/ transitions from the BCS phase to the FFLO phase or from the FFLO phase to the BCS phase can occur [red and blue stars, respectively, in Fig. 1 However, it is important to keep in mind that the use of the mean-field approximation is generally restricted to the weak coupling limit and ground state properties. The limitations of the mean-field method affect the one-dimensional system the most because the pair fluctuations become very important in this case 52 . For such system geometry, the nature of phase transition between the BCS and FFLO phases can be faultily predicted. However, the mean-field approximation can give some useful description in the weak and intermediate couplings, which are comparable with the Bethe ansatz results 52 . While the mean-field FF-type calculations do not predict the correct type of the phase transition from the BCS phase to the FFLO phase, the self-consistent Bogoliubov-de Gennes results are in a good agreement with those obtained from the Bethe ansatz. In the continuum model, the first order phase transition is simply an artifact of the FF-type calculation 52 . This suggests a similar problem in the case of a lattice in the thermodynamic limit, when average particle concentration and lattice constant go to zero (i.e., → n 0 and → a 0 -the dilute limit).
phase separation and its types. Generally, the phase separation is a state of the system where two or more uniform phases (e.g., those which have been defined previously) occur in different parts (so-called domains) of the system. In this section, we distinguish two different types of phase separation, which can emerge in the system, i.e.: (a) a spontaneous phase separation, which can occur in a homogeneous system and (b) an artificially enforced phase separation, which can emerge in inhomogeneous system. In the present work these cases correspond to the system with μ = const i at each site and to the system with inhomogeneous spatial distribution of μ i , respectively. Below, we characterize briefly these two types of phase separations.
Spontaneous (macroscopic) phase separation. The discontinuous transitions between two (homogeneous) phases in the diagram, as a function of μ are usually related to the discontinuous change of particle concentration from value n + to value n − ( > + − n n ). Such a transition can be associated with the occurrence of the (macroscopic) phase separation in a defined range < < − + n n n of particle concentration 53 . In such a phase-separated state two domains with different particle concentrations n − and n + coexist (there can be also regions differing in the magnitude of the order parameter as well as thermodynamic phases). In this approach, because of neglecting the interface energy at the boundaries of the domains, such states can exist only in the thermodynamic limit (i.e., when → ∞ N ) [54][55][56][57][58] . In a finite system, the interface energy can lead to an occurrence of states with other textures 59-63 , besides the homogeneous states and the phase separated states discussed above. Due to the fact that the BCS-FFLO, BCS-NO, and FFLO-NO boundaries in Fig. 1(a) can be discontinuous in the diagram as a function 2 ) exist, respectively [i.e., for the model parameters denoted by the red stars in Fig. 1(b)]. (c) The inhomogeneous state (the enforced phase separation) for the system with two different chemical potentials: μ μ = i 1 for the left half of the system, and μ μ = i 2 for the right half. The red doted lines present the results from panels (a and b), for the comparison. (d) The local polarizations in the real space in the inhomogeneous state. The background colors (yellow and blue) correspond to the stable phases for fixed μ (lower value μ 1 and higher value μ 2 , respectively). [i.e., the model parameters denoted by the blue stars in Fig. 1(b)]. The background colors (yellow and blue) correspond to the stable phases for fixed μ (lower value μ 1 and higher value μ 2 , respectively). www.nature.com/scientificreports www.nature.com/scientificreports/ of n, one expects an occurrence of the following (macroscopic) phase separation regions: BCS/FFLO, BCS/NO, and FFLO/NO 24,49 . Notice that only the BCS-FFLO boundary is discontinuous for all model parameters. To distinguish the phase separation from other states discussed above, we will call it a spontaneous one, because it can occur spontaneously in the homogeneous system. [the red stars in Fig. 1(b)]. The SOP at each site for these two solutions is presented in Fig. 2(a,b). Analogously, we also choose the following parameters for U/ = − . t 3 0 and h/ = . t 0 55 [set  of the parameters, the blue stars in Fig. 1(b)]: μ = . t 0 75 1 (the BCS phase) and μ = . t 1 5 2 (the FFLO phase), and present the local dependence of the SOP in Fig. 3(a,b). One can notice that in the FFLO phase Δ i changes from site to site periodically, whereas for the BCS solution Δ i is homogeneous in the whole system. Now, we investigate the system with a particular distribution of the chemical potential: μ μ = i 1 in one half of system (favoring the FFLO phase or the BCS phase for set  or , respectively) and μ μ = , for  and , respectively. As can be expected, the spin polarization has a modulation, which is two times faster than that of the SOP. Additionally, maximal values of the spin polarization are located at the nodal points of the SOP 64,65 . This is a consequence of the interplay between the unpaired polarized particles and those in the superconducting state (i.e., the Cooper pairs). A small number of pairs (at sites with Δ ≈ 0 i ) supports the occurrence of polarized states (i.e., a large value of m i ). Note also that there is no coexistence of the superconducting and magnetic ordering in the BCS phase.
The nature of this inhomogeneous state is, however, different than the origin of the spontaneous (macroscopic) phase separation discussed in the previous point. In the present setup, the system parameters are inhomogeneous (different values of chemical potential in two parts of the system). We will call such separated state an (artificially) enforced phase separation in contrary to the spontaneous phase separation, which could occur in the homogeneous system.
Note, that even if μ μ μ = = c 1 2 (where μ c is a critical value of the chemical potential at the FFLO-BCS boundary for particular U/t and h/t), due to the fact that the interface has small but still finite energy, the finite system cannot exhibit spontaneous phase separation. In such a case, the whole system is in the FFLO phase or in the BCS phase (both states have equal energy and both solutions correspond to local minima of the grand canonical potential). system with a trap. After studying the model systems, let us investigate more realistic situations, in which the on-site potential μ i changes from site to site. An example of experimental realization of such systems are ultracold atomic gases in optical lattices. In the following, we consider two types of traps: (i) a linear trap with where r i is a location of i-th site, whereas r c is a location of the center of the trap (i.e., = i N / = 2 500). In this part, we performed the calculations for = N 1000 sites with open boundary conditions. V 0 was chosen in such a way to ensure a value much larger than the chemical potential corresponding to the BCS-NO transition in the homogenous system at = h 0 and fixed U. The self-consistent solutions of the BdG equations for the case of a linear trap are presented in Fig. 4. In this case, we show SOP UΔ i , magnetization m i , particle concentration n i , and the trapping potential ) as a function of the lattice site i (blue, green, red and dashed black lines, respectively). Due to the fact that μ i changes from site to site, n i changes at each site and one does not observe the solutions with, e.g., constant value of Δ i , corresponding to the BCS phase (but the region without oscillations of the SOP can be identified as corresponding to this phase, see below). However, one can indicate the regions along the chain which correspond to the phases indicated in Fig. 1(a). Namely, starting from the left side of Fig. 4 ). Note, that the location i of the boundaries between the regions along the chain corresponds, approximately, to the values of μ i for which the phase transition occurs in the homogeneous system [cf. Fig. 1(b)]. It is attributed to the fact that the interactions in the system are short-ranged and the interface energy is relatively small. Moreover, (2019) 9:6719 | https://doi.org/10.1038/s41598-019-42044-w www.nature.com/scientificreports www.nature.com/scientificreports/ for set , only the NO phase with = n 0 can exist, whereas for set  also the NO phase with ≠ n 0 is possible (cf. ref. 24 ).
The results for a harmonic trap, presented in Fig. 5, do not differ qualitatively from those obtained for a linear trap. We The results presented above show that in the system with a trap at each site, the state of a homogeneous system corresponding to μ μ = i is approximately realized. The trap shapes considered above are chosen in such a way to from a range that all of possible phases occurring in the system are present, at fixed h cf. Fig. 1 (from the NO phase at large |μ| at the boundaries of the system, i.e., at = i 0, − N 1, to the half-filled phases at μ = 0 in the center, i.e., at = i N /2). Thus, in other words, all phases for a homogeneous system at fixed μ can be realized in the setup with the trap. However, by modification of parameters of the trap for instance, one can realize a system which corresponds only to a part of the phase diagram (i.e., some limited part of μ range). A sample of the results is shown in Fig. 6. They are obtained for the system with a harmonic trap in the form of V t 1 25 . In this case, the minimum value of the trapping potential is shifted to some finite value, i.e., μ = V N /2 (black doted line). The chosen parameters allow to realize a part of the phase diagram in which two superconducting domains are separated by the NO state [cf. Fig. 1(b)]. Depending on a fixed value of magnetic field h, the states with BCS or BCS-FFLO core surrounded by the BCS shell, separated by the NO sub-shell-part can be realized in the system.
The results shown in Fig. 6 correspond to a part of the phase diagram [cf. Fig. 1(b)], where the BCS-BEC crossover is realized. One should notice that the system under consideration is in a trap, which can change the state at the site with μ i in comparison to the homogeneous solution with μ μ = i . Indeed, numerical results show that this state can be realized outside the BCS core state [ Fig. 6(a)], for h smaller than the one at which the FFLO  www.nature.com/scientificreports www.nature.com/scientificreports/ phase occurs in a homogeneous system. The realization of the two BCS-like domains, separated by the NO state, is possible by the decrease of h [ Fig. 6(b)]. At the boundary of the internal BCS core, one observes the FFLO-type gap oscillation with changing of its sign [66][67][68] . It is important to emphasize that, as a consequence of low filling at the external BCS shell state, we expect realization of the BEC condensate (the tightly bound local pairs region) of the Cooper pairs. According to the Leggett criterion 69 , the BEC (i.e., Bose-Einstein condensate) begins when the effective chemical potential is smaller than the lower band edge.

summary
In conclusion, we have shown that properties of the system in the presence of a trap are strongly associated to the phase diagram of the homogeneous system (in the absence of a trap). It is clearly seen from the spatial dependence of different quantities in the system with linear and harmonic trap (Figs 4 and 5, respectively) and provides a confirmation about the validity of so-called local density approximation, that is very often employed in theoretical calculations for trapped system starting from previous results for homogeneous systems 22,48,70 . We have shown that in such a system the core-shell structures can be created due to the trapping potential. For chosen parameters, the order of states realized in such structures corresponds to the sequence of phases occurring in the phase diagram for the homogeneous system (Fig. 1). The states occurring in a particular sequence depend mainly on values of model parameters (mainly h and U). Generally, the shape of the trapping potential does not change the structure of phases occurrence. Moreover, even in the same type of trap, for different values of the model parameters (i.e., interaction U or magnetic field h), the structure with the BCS core and the FFLO shell as well as with the FFLO core and the BCS shell can be realized (depending on U and h).
An increase of the on-site interaction at low filling can lead to the BCS-BEC crossover in the system. By tuning of the trap parameters, we can realize the part of the phase diagram in which this crossover occurs. Thus, one can obtain the more complex core-shell structures. In particular, one can obtain two BCS states separated by the NO state (e.g., Fig. 6). Due to the low filling, at the outer BCS shell, the BEC should emerge in the system.
In this work, we have shown that different core-shell structures can occur in a trapped system. These structures are examples of so-called (artificially) enforced phase separation, occurring in spatially inhomogeneous systems, whose origin is different than that of the macroscopic phase separation in homogeneous systems. It is important to emphasize that such experimental setup with a trap allows to investigate phase diagrams of homogeneous systems with short-range interactions. The theoretical prediction presented in this work should be realizable experimentally in a relatively simple way. h t 0 375 . The dependence of (local) chemical potential μ i is marked with the dotted line. Insets present schematically the region of the phase diagram for homogeneous systems corresponding to the range of μ i changes [cf. with Fig. 1(b)]. The background colors (blue, yellow and pink) indicate the boundaries between the phases (NO, BCS and FFLO, respectively).