Zero temperature momentum distribution of an impurity in a polaron state of one-dimensional Fermi and Tonks-Girardeau gases

We investigate the momentum distribution function of a single distinguishable impurity particle which formed a polaron state in a gas of either free fermions or Tonks-Girardeau bosons in one spatial dimension. We obtain a Fredholm determinant representation of the distribution function for the Bethe ansatz solvable model of an impurity-gas $\delta$-function interaction potential at zero temperature, in both repulsive and attractive regimes. We deduce from this representation the fourth power decay at a large momentum, and a weakly divergent (quasi-condensate) peak at a finite momentum. We also demonstrate that the momentum distribution function in the limiting case of infinitely strong interaction can be expressed through a correlation function of the one-dimensional impenetrable anyons.

1 Introduction 2 2 Model 4 2.1 Defining |min Q for impurity-gas repulsion 7 2.2 Defining |min Q for impurity-gas attraction: gas state 9 2.3 Defining |min Q for impurity-gas attraction: bound state 9 3 Fredholm determinant representation in the thermodynamic limit 10 3.1 Impurity-gas repulsion 10 3.2 Impurity-gas attraction: gas state 11 3.3 Impurity-gas attraction: bound state 11 1 Introduction Non-interacting Bose and Fermi systems have markedly different momentum distribution functions at low temperature. Bosons tend towards a macroscopic occupation of the zeromomentum state, and fermions spread over the volume of the Fermi sphere. When interparticle interactions are present, the distinction becomes not at all evident. We know from exactly solvable models in one spatial dimension that some observables evolve smoothly from bosonto fermion-like behavior, as a function of the inter-particle interaction strength. An example is provided by the Lieb-Liniger model, representing a gas of bosons interacting through a δ-function potential of an arbitrary strength g [1,2]. The excitation spectrum of the model in the g → ∞ limit, the Tonks-Girardeau gas, is the same as the one of a free Fermi gas [3,4]. Furthermore, any excitation in the Lieb-Liniger gas is parametrized by a set of distinct integers, same way as for a free Fermi gas, giving rise to the notion of the Pauli principle for one-dimensional interacting bosons [5]. This is consistent with the fact that the low-energy and momentum excitations of the interacting gapless one-dimensional Bose and Fermi systems can be interpreted as collective boson modes of a unique effective field theory, the procedure called the bosonization [6,7]. Despite of these similarities, the momentum distribution functions of the Tonks-Girardeau and free Fermi gases are radically different, which is seen from the exact [8,9], as well as asymptotic formulas [10,11]. How do interactions shape the momentum distribution function of a single distinguishable mobile particle, an impurity, interacting with a one-dimensional system? It has been demonstrated in Ref. [12] that the function n(k), defined as the probability to find the impurity in the state having the momentum k, does not have a single-particle delta-peak δ(k) in one spatial dimension, for any non-zero value of the impurity-gas coupling strength. Instead, n(k) ∼ k ν in the k → 0 limit. The value of ν was found only in the limit of the vanishing impurity-gas coupling strength [12]. Extending this result to an arbitrary coupling strength is a far-from-trivial problem. This is because the many-body spectrum of the whole system contains low-energy excitations with quadratic dispersion relation. The application of the bosonization technique is not straightforward for such a spectrum [13,14]. The recently developed paradigm of the non-linear Luttinger liquids [15] could perhaps be used to find ν for an arbitrary interaction strength. However, this has yet to be done. As for finding the exact shape of n(k) in the whole range of values of the momentum k, the Bethe ansatz solution remains the only non-perturbative analytical approach available thus far.
In the present paper we investigate the shape of the momentum distribution function n(k, Q) of an impurity interacting with a free Fermi gas in one spatial dimension. The system stays in the polaron state, defined as the minimum energy state at a given total momentum Q = P imp + N j=1 P j (Ref. [12] is dealing with Q = 0 only). The impurity has the same mass as the gas particle, and interacts with the gas through a δ-function potential of an arbitrary (positive or negative) strength g. The Hamiltonian reads Here, x j (P j ) is the coordinate (momentum) of a gas particle, j = 1, . . . , N, and x imp (P imp ) is the one of the impurity. Such a model is Bethe ansatz solvable; its eigenfunctions and spectrum have been found by McGuire [16,17]. McGuire's solution is a special case of the Bethe ansatz solution for the Gaudin-Yang model [18][19][20], having a peculiarity that any eigenfunction can be written as a single determinant resembling the Slater determinant for the free Fermi gas [21][22][23]. Such a representation, so far not available for any other interacting Bethe ansatz solvable model, enabled the derivation of an exact analytical expression for the time-dependent two-point impurity correlation function at zero [24] and arbitrary temperature [25]. Here, we present an exact analytical expression for n(k, Q) in the limit of infinite system size, L → ∞, valid for an arbitrary (positive or negative) coupling strength g and zero temperature. The answer is given in terms of the Fredholm determinant of a linear integral operator of integrable type (see, e.g, section XIV.1 of [5]). We use our exact analytical result (i) To obtain the largemomentum tails of n(k, Q), and the root mean-square uncertainty of the average momentum of the impurity. (ii) To extract a quasi-condensate-like divergence of n(k, Q) at k = Q. (iii) To establish the correspondence between n(k, Q) in the g → ∞ limit and a correlation function of the one-dimensional impenetrable anyons. The paper is organized as follows. In section 2 we define the model under consideration. In section 3 we summarize our exact analytical results expressed in terms of the Fredholm determinants. In sections 4 through 7 we analyze various limiting cases of the formulas from section 3. Section 8 explains principal steps of the calculation used to get the Fredholm determinant representation of section 3. We conclude in section 9. The appendices are selfexplanatory.

Model
Our objective is to compute the momentum distribution function of an impurity, interacting with a free one-dimensional spinless Fermi gas at zero temperature. Here, |min Q is a polaron state, defined as the minimum energy state of the system having the total momentum Q and containing only one impurity. We discuss the properties of the polaron state later in this section. Note that our result for the function (2) is also valid for the impurity immersed into the Tonks-Girardeau gas. This can be explained using the arguments given in the end of section 2 in Ref. [25]. The Hamiltonian of the entire system is where is the Hamiltonian of the free Fermi gas, m is the particle mass, and The creation (annihilation) operators ψ † σ (ψ σ ) carry the subscript σ =↑ for the spinless Fermi gas, and σ =↓ for the impurity. We have The Hamiltonian (3) defines the fermionic Gaudin-Yang model [18][19][20], in which the number of the impurity particles, is arbitrary. However, the states with N imp > 1 do not contribute to the function (2). The first-quantized form of the Hamiltonian (3) with N imp = 1 and N particles from the Fermi gas is given by Eq. (1). The Planck constant, , is equal to one in our units. A commonly used dimensionless form of the impurity-gas coupling strength g is where is the gas density. To further simplify notations, we let m = 1 (10) and measure all momenta in the units of the Fermi momentum, We restore m and k F in the captions to the figures. Equation (2) can be written as where is the Q-dependent reduced density matrix of the impurity. The normalization condition k n(k, Q) = L 2π (14) is equivalent to (0) = 1.
For the system in a finite volume L, periodic boundary conditions are imposed. That n(k, Q) is real implies the involution (−y) = * (y), where the star stands for the complex conjugation. The symmetry n(−k, Q) = n(k, −Q).
The coordinate representation for |N is the Slater determinant All eigenstates of the Hamiltonian (1), |min Q being one of them, have been found in Refs. [16,17]. Let |Q be an eigenstate having total momentum Q. Such a state is parametrized by the quasi-momenta k 1 , . . . , k N +1 satisfying The energy of the state |Q reads Each k j should satisfy the equation where Here, where γ is given by Eq. (8). Thus, one has a system of N + 1 equations (26) for the variables k 1 , . . . , k N +1 and Λ. These equations, called the Bethe equations, are coupled through Eq. (24). Any solution to this system has the following properties [17]: (i) Λ is real. (ii) If α ≥ 0 all k j 's are real. (iii) If α < 0 either all k j 's are real, or k 1 , . . . , k N −1 are real, while k N and k N +1 have a non-zero imaginary part, and k N = k * N +1 . We will often use the following representation of the Bethe equations (26): where and Taking the derivative of Eq. (29) with respect to Λ we get , j = 1, . . . , N + 1.
The point of focus of our paper is n(k, Q) in the thermodynamic limit, defined as the limit of infinite system size, L → ∞, at a constant density In what follows, we use L → ∞ in place of L, N → ∞ for simplicity of the notations. The choice of the boundary conditions should play no role for n(k, Q) in the thermodynamic limit. The sum over momenta turns into the integral, and the normalization condition (14) becomes In sections 2.1 through 2.3 we proceed with solving the system of Eqs. (24) and (26) in the thermodynamic limit for the state |min Q entering Eq. (2).

Defining |min Q for impurity-gas repulsion
In the case of the repulsive interaction, γ ≥ 0, the L → ∞ limit of Eq. (27) reads We adopt the convention that the distinct integers n j are enumerated in the increasing order, n 1 < · · · < n N +1 . Equation (24) turns into the algebraic relation between Λ and Q: where Z = arctan(α − Λ) + arctan(α + Λ) απ (38) and ϕ = 1 The function Q D encompasses all n j 's: The energy (25) turns into where Let Such a choice leads to Q D = 0, and corresponds to the minimum energy state |min Q for −1 ≤ Q ≤ 1. Equation (37) turns into The parameter Λ runs from −∞ to ∞ when Q runs from −1 to 1. Equations (42) and (44) determine E min as a function of Q for −1 ≤ Q ≤ 1. The minimum energy state for Q outside of that interval is parametrized by consecutive sets of n j 's other than given by Eq. (43). The result is a smooth periodic function of Q, plotted in the left panel of Fig. 1. Note that and Therefore, decreases from 1/2 to zero when γ increases from zero to infinity.

Defining |min Q for impurity-gas attraction: gas state
The gas state is defined for the attractive interaction, γ < 0, as the minimum energy state for all k j 's being real. Such a state has been realized experimentally for the Lieb-Liniger gas in the experiment with ultracold atoms [26]. The analysis following the steps from section 2.1 leads to Eqs. (42) and (44) in which γ is now negative. This results in E min (Q) being an odd function of γ. Therefore, the function [E min (Q) − E min (0)]/[E min (1) − E min (0)] coincide with the one for the repulsive case, plotted in the left panel of Fig. 1. The function decreases from zero to −1/2 when γ increases from minus infinity to zero. This means that the minimum energy state for a weak repulsion, γ 1, does not go continuosly to the gas state for a weak attraction, −γ 1. Rather, it turns into the weakly attractive bound state, discussed in section 2.3.

Defining |min Q for impurity-gas attraction: bound state
The bound state is the true minimum energy state for the attractive interaction, γ < 0. That is, k j 's are not required to be real, as it was for the gas state, section 2.2. As a result, the phase shifts take the form (36) for the real k 1 , . . . , k N −1 , and [17] where k ± is defined by Eq. (31). Therefore, Eq. (24) takes the form where Q D is given by Eq. (40) with j running from 1 to N − 1. Like in the case γ > 0, we have Q D = 0 for the minimum energy states in the interval −1 ≤ Q ≤ 1: where The function entering Eq. (41) is plotted in the right panel of Fig. 1, and is a periodic function of Q. Unlike for the repulsive and the attractive gas state, (i) Λ runs through the finite interval, −Λ F ≤ Λ ≤ Λ F , when Q runs from 1 to −1 in Eq. (51); (ii) E min (Q) has cusps at Q = ±1, ±2, . . ..

One has
The function E min (1) is obtained by substituting Λ F into Eq. (53), and increases from 1/4 to 1/2 when γ goes from minus infinity to zero.

Fredholm determinant representation in the thermodynamic limit
In this section we show the main results of our paper: exact analytic formulas for the impurity momentum distribution function n(k, Q) at zero temperature and an arbitrary positive and negative impurity-gas interaction strength g. These formulas contain Fredholm determinants of linear integral operators. Let V be an M × M matrix with the entries V jl = V (k j , k l ), I be the identity matrix, and Then the Fredholm determinant is The right hand side of Eq. (57) taken for a large but finite M can be used to evaluate the Fredholm determinant numerically [27]. An equivalent definition, appears in the mathematical literature on the linear integral operators theory (see, e.g., [28], vol IV, p.24). Naturally,V can be recognized as a linear integral operator with the kernel V (q, q ) on the domain [−1, 1] × [−1, 1]. The necessary existence and convergence conditions are fulfilled for the operators encountered in our paper. The energy of the state |min Q is a periodic function of Q, and n(k, Q), defined by Eq. (2), inherits this periodicity. We rewrite Eq. (12) as In what follows, we write (y) explicitly for the positive values of y, and use the involution (16) to get it for the negative values.

Impurity-gas repulsion
The Fredholm determinant representation in the case of the repulsive impurity-gas interaction, The identity operator is denoted byÎ. The kernels of the linear integral operatorsK andŴ , on the domain [−1, 1] × [−1, 1], are defined by and respectively. The kernel (61) belongs to a class of integrable kernels [5,29]. The functions e ± are defined as and Here, the phase shift δ(q) is defined as and the value of Λ can be found as a function of Q by using Eq. (44). The behavior of the momentum distribution function is illustrated in Fig. 2(a).

Impurity-gas attraction: gas state
All formulas from the section 3.1 are valid for the gas state after letting γ be negative. The behavior of the momentum distribution function is illustrated in Fig. 2

Impurity-gas attraction: bound state
The presence of the bound state qualitatively affects the Fredholm determinant representation for the function (y), as compared with Eq. (60): Here, the kernels of the linear integral operatorsK b andŴ b are defined by and respectively. The function f is defined as k + and k − are defined by Eq. (31), and The other functions entering Eqs. (66)-(71) are defined in section 3.1. The typical behavior of the momentum distribution function is shown in Fig. 2(c).

Limit of strong interaction, |γ| → ∞
Correlation functions of the model (1) in the γ → ∞ limit has been represented as Fredholm determinants in the works [30,31]. Using the Fredholm determinant representation we demonstrate that the one-body density matrix (y) in the γ → ∞ limit can be written as a correlation function of the one-dimensional impenetrable anyons. Such a correspondence remains valid for the gas state in the γ → −∞ limit.

Impurity-gas repulsion
We begin with discussing the γ → ∞ limit of the impurity-gas repulsion. The kernels (61) and (62) simplify significantly when compared to arbitrary γ. Using that we have in the leading order in α This gives us for the kernel (61), and for the kernel (62). The γ → ∞ limit of Eq. (44) reads Substituting this formula into Eq. (75) we get Let us now show how (y) emerges in the model of one-dimensional impenetrable anyons [32]. Recall that the anyon field operators satisfy the commutation relations and Here, sgn(x) = |x|/x, sgn(0) = 0, and κ is the statistics parameter. The correlation function ψ † A (y)ψ A (0) has the Fredholm determinant representation, given by Eq. (4) from Ref. [32]. The transformation explained in Ref. [5] (see the discussion of the equivalence between Eqs. (3.12) and (3.13) in Ch. XIII therein) leads us to the equality where λ entering the kernel (74) is related to the statistical parameter κ as follows: Comparing Eqs. (78) and (82) we get for κ and Q in the interval between minus one and one. The left hand side of Eq. (81) has also been extensively evaluated numerically [33,34]. However, no connection between the mobile impurity and anyon correlation functions, as suggested by Eqs. (81) and (83), has been given in the literature. Furthermore, the Jordan-Wigner transformation connects the anyon field operators and the fermion operators. Therefore, the right hand side of Eq. (81) is a correlation function of a free spinless Fermi gas: where |FS stands for the Fermi sea. Since Eq. (2.19) from Ch. XIII in Ref. [5] gives it is ψ † F and ψ F that lead to the emergence of the rank-one operatorŴ in Eq. (60). Note that the evaluation of the right hand side in Eqs. (85) and (86) can be done by using the Wick's theorem (for Eq. (86) see, e.g., Ref. [35]), without any use of the coordinate representation of the wave functions of the model. Interestingly, in a recent work [36] a two-dimensional impurity model has been linked to anyons, albeit in a different manner: There the statistical parameter is related to the impurity-phonon coupling.

Impurity-gas attraction: gas state
We now turn to the case of the gas state for the impurity-gas attraction, introduced in section 2.2. The γ → −∞ limit of the Fredholm determinant representation introduced in section (3.2), leads to the same formulas as the γ → ∞ limit, discussed in section (4.1).

Impurity-gas attraction: bound state
Finally, we consider the bound state for the impurity-gas attraction. We take the γ → −∞ limit in the formulas of section (3.3) and get in the leading order Furthermore, it follows from Eq. (51) Therefore, we write the following asymptotic expression: Substituting this into Eq. (59) we get The γ → ∞ expansion (90) is not a uniform estimate of the exact result for n(k, Q), since it misses the divergence at k = Q, discussed in detail in section 7. Still, it conveys an important message: the impurity momentum distribution becomes completely flat, and infinitely broad, in the γ → −∞ limit.

Total momentum
is particular (recall that k F = 1 everywhere but in the captions to the figures). One finds that n(k, 1) for the repulsive ground state, section 3.1, and the attractive gas state, section 3.2, coincide with the momentum distribution of a free Fermi gas. It follows from Eq. (44) that Λ = ∞ at Q = 1. We have in the leading order in Λ and Eq. (60) takes the form Plugging this function into Eq. (59) we indeed get the momentum distribution of a free Fermi gas. This result can also be obtained without using the Fredholm determinants. For the Gaudin-Yang model, Eq. (3), all three of the Hamiltonian, the spin-ladder operator and the total momentum P , commute with each other. Therefore, any state |Ψ FF of a free Fermi gas with N + 1 particles can be turned into an eigenstate of the Hamiltonian (1), containing N host particles and one impurity, and having the same energy and momentum as |Ψ FF . Furthermore, The state (96) is the minimum energy state |min Q for Q given by Eq. (91) and |Ψ FF being the minimum energy state of a free Fermi gas for the same Q. This can be shown very straightforwardly by examining the exact eigenfunctions and spectrum of the model (1), see, for example, section 5 of Supplementary Information in Ref. [37]. Equation (97) gives the momentum distribution of a free Fermi gas immediately. The case of the bound state for the attractive interaction, sections 2.3 and 3.3, is different. The shape of n(k, Q) is qualitatively the same at Q given by Eq. (91) in comparison with any other value of Q. This is because the state (96) is not the minimum energy state of the Hamiltonian at any value of Q. We plot n(k, 1) in Fig. 2(c).
6 n(k, Q) in the k → ∞ limit The large k limit of n(k, Q), following Eq. (59), is determined by an expansion of (y) in the vicinity of y = 0. It turns out that , ∂ y , and ∂ 2 y are continuous at y = 0. Therefore and The third derivative of (y) has a discontinuity at y = 0. This implies for the leading term of the large k expansion n(k, Q) = 1 2π Taking into account the involution (16) we arrive at where Each of Eqs. (98), (99), and (101) has a lot of physics behind. We discuss them one-by-one.

Analysis of P imp
For the repulsive ground state, and the attractive gas state we have where Λ and Q are connected by Eq. (51). Using the Hellmann-Feynman theorem as explained in Ref. [38] gives the average momentum of the impurity in terms of the group velocity, Though the curves in Fig. 3(b) are not straight lines, the difference cannot be seen by a naked eye. It follows from Eq. (105) that the slope of P imp at Q = 0, is set by the value of the effective mass m * defined by the expansion of E(Q) at Q = 0: The explicit form of E(Q) is discussed in section 2. The analytic formula for m * corresponding to Eq. (103) is γ > 0 ground state, and γ < 0 gas state (108) (π + arctan α) 2 π + arctan α − α(1 + α 2 ) −1 , γ < 0 bound state.
The analytic expressions (108) and (109) for the effective mass were obtained for the first time in the works [16] and [17], respectively. The γ → ∞ limit of Eq. (108) is m/m * = 0: the impurity becomes infinitely heavy. This is contrasted with the γ → −∞ limit of Eq. (109), which is m/m * = 1/2: the mass of the impurity bound to the gas particles remains finite. A quantitative comparison between m * for γ > 0 from Eq. (108), and m * for γ < 0 from Eq. (109) is made in Fig. 4 .
6.2 Analysis of the coefficient C in the large k expansion n(k, Q) = C/k 4 In this section we give the explicit analytic formula for the coefficient C in Eq. (102). For the repulsive ground state, and the attractive gas state we have γ > 0 ground state, and γ < 0 gas state, where Z is defined in Eq. (64) and ϕ by Eq. (39). Recall that Λ and Q are connected by Eq. (44), and note that C in Eq. (110) is an even function of γ. For the attractive bound state Z is replaced with Z b , Eq. (71). Hence, The former line tends to zero, and the latter tends to 1/2, in the |γ| → ∞ limit.
To what extent C could be extracted numerically from the large momentum behavior of n(k, Q) is illustrated in Fig. (6). We evaluated n(k, Q) from the Fredholm determinant representation presented in section 3.

Analysis of P 2 imp
The average of P 2 imp , Eq. (99), is expressed through P imp and C: where, by definition, is a root-mean-square deviation. Equation (116) is valid for the repulsive ground state and attractive gas state. The result for the attractive bound state is obtained by replacing Z with Z b . Exemplary plots of σ are shown in Fig. 7. 7 n(k, Q) in the k → Q limit In this section we present the y → ∞ expansion of (y). We use it to prove the existence of the power-law singularity seen in Fig. 2, as well as to calculate the exponent ν, and the numerical prefactor. So far, ν has only been found at Q = 0 and γ → +0 in Ref. [12]; this result follows from our formulas as a particular case. for Q = 0. The dotted red, and dashed blue lines are for Q = 0.8k F : the former is for k > 0, and the latter is for k < 0. Note that the Fredholm determinants are numerics-friendly, but n(k, Q) decays very fast with increasing |k|, and this makes the numerical evaluation of C a challenge.

Large y expansion of (y) in case of impurity-gas repulsion
The density matrix and the momentum distribution are related by Eq. (59). Both are 2k Fperiodic in Q (recall that k F = 1 everywhere but in the captions to the figures). This property together with Eq. (17) makes it sufficient to examine for 0 ≤ Q ≤ 1 only. The large y expansion of the determinant representation (60) can be obtained by a finite-size analysis of the form-factors followed by a resummation of the soft modes, along the lines of the works [15,[44][45][46][47][48]. We leave the details for a separate publication. The result is The numerical prefactor depends on γ and Q through the phase shift (65): Here, the coefficient Z is given by Eq. (64): and G stands for the Barnes G-function, defined by the functional equation with the normalization G(1) = 1, where Γ(z) is the Euler Gamma function. The functionF entering the second term on the right hand side of Eq. (119) is andÃ follows from A by replacing F withF in Eqs. (120), (122) and (123). The second term on the right hand side of Eq. (119) is, generally, subleading -it decays faster than the first one: However, the inequality turns into an equality at Q = 1, that is, the subleading term becomes of the same order as the leading one, and their sum in Eq. (119) reproduces the exact formula (y) = sin y y , Q = 1.
We show (y) evaluated from the exact expression (60), and the convergence of the asymptotic formula (119) to this exact expression in the panels (a) and (d) of Fig. (8), respectively. We would like to emphasize that the decay rates of the leading and the first subleading terms in Eq. (119) are close to each other when Q is close to one. Top panels: absolute value of (y) from the exact formula. Bottom panels: the absolute value of the ratio of (y) from the exact and large-y-asymptotic formulas.
7.2 Large y expansion of (y) in case of impurity-gas attraction: gas state All formulas from the section 7.1 are valid for the gas state after letting γ be negative. We show (y) evaluated from the exact expression (60), and the convergence of the asymptotic formula (119) to this exact expression in the panels (b) and (e) of Fig. (8), respectively.

Large y expansion of (y) in case of impurity-gas attraction: bound state
In case of the attractive bound state, the explicit expression for (y) is given by Eq. (66), and the leading term in the y → ∞ expansion reads where and Z b given by Eq. (71). The prefactors A andÃ, Eq. (120), depend on γ and Q through the phase shift only. By contrast, the prefactor A b , Eq. (129), depends on γ and Q explicitly. We show (y) evaluated from the exact expression (66), and the convergence of the asymptotic formula (128) to this exact expression in the panels (c) and (f) of Fig. (8), respectively.

The exponent ν and the prefactor in Eq. (118) for n(k, Q)
The singular part of the momentum distribution, Eq. (118), is fully characterized by the asymptotic expressions for (y). Equation (119) leads to the exponent γ > 0 ground state, and γ < 0 gas state (131) and Eq. (128) leads to Both Eqs. (131) and (132) tend to the same value in the |γ| → ∞ limit, which coincides with the result from Ref. [23]. This limiting value is indicated with the thin dotted line in Fig. 9. One can also see that ν = 0 when Q reaches the Fermi momentum for the γ > 0 ground state, and γ < 0 gas state. Recall that n(k, Q) turns into the Fermi function at Q = 1, as illustrated in the panels (a) and (b) of Fig. 2 and discussed in section 5. The case γ < 0 bound state is different, there ν is a non-trivial function of γ at Q = 1. Letting Q = 0 and γ → +0 in Eq. (131) we get This gives the same dependence on γ as in Ref. [12].

Determinant representation for finite N
In this section we present the impurity momentum distribution function n(k, Q) for a finite particle number N through determinants of finite-dimensional matrices. This result is crucial for deriving the Fredholm determinant representation of section 3. Recall that we stick to the notations of the paper [25], whenever possible.  Our starting point is Eq. (19). We write the form-factor as given by Eq. (5.23) from Ref. [25]: Here, ∂k j /∂Λ is defined by Eq. (32), and for the determinant of the (N +1)×(N +1) matrix. The momentum Q of the state |min Q is the sum of the quasi-momenta k 1 , . . . , k N +1 , Eq. (24). How these quasi-momenta are specified is discussed in sections 2.1 through 2.3. The momentum of the state |N is the sum of p 1 , . . . , p N . Combining Eqs. (21) and (24) implies the constraint for the sum over p 1 , . . . , p N in Eq. (19). We transform Eq. (19) by replacing the constraint (137) with the Kronecker delta: The summations over p 1 , . . . , p N on the right hand side of Eq. (138) run independently from each other. One can see from Eqs. (135) and (136) that N |ψ k↓ |min Q = 0 if p j = p l at j = l. The factor 1/N ! is to compensate counting the form-factor multiple times upon the permutations of p 1 , . . . , p N . Equations (12) and (16), and the representation The terms on the right hand side of Eq. (141) are determined by Eq. (135), and p 1 , . . . , p N are quantized as given by Eq. (22). We now take the sum over p 1 , . . . , p N in Eq. (141). Let us consider the function where det D is defined by Eq. (136), f is an arbitrary function, and p j s are quantized as given by Eq. (22). After some elementary transformations (used, for example, to get the identities in appendix B.3 from Ref. [25]) we come at the following representation for Eq. (142): Here, and p = 2πn/L, n = 0, ±1, ±2, . . .. For γ > 0 repulsive ground state and γ < 0 attractive gas state the quasi-momenta k 1 , . . . , k N +1 are real. This implies Furthermore, one can show that for any real-valued k j (see, for example, section 5.2 from Ref. [25]). We, therefore, can use the identity (143) for the function (141), and get where and The matrix B has rank one, and we can write Eq. (147) as We now turn to the γ < 0 bound state. Here, k 1 , . . . , k N −1 are real, and k N = k * N +1 are complex. This implies It follows from Eq. (24) that Since Q and Λ are connected by Eq. (51), we get Using the identity (143) for the function (141) we come at Eqs. (147)-(149). Later, we will use the following representation for the entries of the matrix (148): where The uncertainty in Eq. (154) at j = l can be resolved by L'Hôpital's rule, which amounts to making use of the expansion That is, where c(k) is given by Eq. (155) and ∂k j /∂Λ by Eq. (32). Let us represent the function c from Eq. (155) as where Γ is a union of counter-clockwise-oriented contours around the points z = 2πn/L. Assuming that k is real, we deform Γ into a contour encircling the point z = k, and two straight lines infinitesimally above and below the real axis: We assume 0 < y < L; the result for y = 0 and y = L follows from the continuity of (y). The first integral is equal to zero, which is seen by using Cauchy's residue theorem (the integration contour is extended to the closed one by adding a half-circle in the upper half-plane). The second integral is equal to zero for the same reason (the integration contour is extended to the lower half-plane). Therefore, we get for Eq. (155): We now introduce the function Substituting the Bethe equations (29) into Eq. (160) we find c(k j ) = e(k j ), j = 1, . . . , N + 1. Furthermore, and Using Eqs. (161)-(163) we get for Eq. (157) This expression can be represented as follows Thus, we can write the matrix (154) as where Here, where ∂k j /∂Λ is defined by the exact formula (32). The uncertainty in Eq. (169) at j = l can be resolved by L'Hôpital's rule. The matrix (149) can be written as where Using Eqs. (168)-(172) we get for Eq. (150) Recall that we are working at a finite constant density, Eq. (9). The expression (173) is valid in the interval 0 ≤ y ≤ L.
That the exact function (y) is L-periodic and satisfies the involution (16) implies the exact identity (L − y) = * (y).
We have verified numerically that Eq. (173) with the kernels (169)-(172) satisfies Eq. (174) for any N , and y in the interval 0 ≤ y ≤ L. We have also verified it by performing symbolic computations using Mathematica package for N = 2.
Let us now discuss the case of the complex quasi-momenta: Im(k N ) < 0 and Im(k N +1 ) > 0, Eq. (49). The representation (158) leads to The first (second) integral gives non-zero contribution for Im(k) > 0 (Im(k) < 0). In both cases one arrives at Eq. (160). Further analysis is the same as for the real quasi-momenta, it leads to Eqs. (169)-(173). Note that and the involution (174) holds true.
We plot (y) in Fig. 10. The top panels show that it oscillates if Q = 0. The bottom panels (d) and (e) demonstrate that the oscillations are largely, but not fully, suppressed for the function e iQy (y). Since the number of the gas particles, N = 40, used in the plot, is large, the residual oscillations seen in the bottom panels (d) and (e) can be attributed to the subleading term written explicitly on the right hand side of Eq. (119), valid in the thermodynamic limit. There are no visible oscillations in the bottom panel (f), consistent with the small contribution of the subleading terms to the asymptotic formula (128). Note that the oscillations of the function e iQy (y) can be seen in Fig. 4 from Ref. [23], though the thermodynamic limit have not been taken in the analytic formulas used therein, and the period of the oscillations has not been identified.

Conclusion
The main result of the present paper is the Fredholm determinant representation, Eqs. (60) and (66), for the momentum distribution function, n(k, Q), of an impurity which formed a polaron state with a free Fermi gas (or the Tonks-Girardeau gas [3,4]). Using this representation we examined how the properties of the impurity depend on the strength g of the impurity-gas δ-function interaction potential, and on the value of the total momentum Q of the system (which is the same as the momentum of the polaron). We have found that the formation of the bound state strongly affects the behavior of n(k, Q). In the absence of the bound state n(k, Q) turns into the Fermi function at Q = 1 (recall that the momenta are given in the units of the Fermi momentum k F everywhere except for the captions to the figures). This can be seen in Fig. 2(a) and 2(b). In the presence of the bound state n(k, Q) has a weak singularity at k = Q for any Q, including Q = 1, and the Fermi function does not emerge. This can be seen in Fig. 2(c). The distinct role the bound state plays in the behavior of the impurity's momentum distribution function is reflected at the level of the dispersion relation of the polaron. Indeed, the group velocity of the polaron vanishes at Q = 1 in the absence of the bound state, see Fig. (1)(a). Such a vanishing velocity is consistent with the impurity spreading over the Fermi sea and mimicking the distribution of the gas particles in the momentum space. In the presence of the bound state the group velocity of the polaron does not vanish at Q = 1, therefore the momentum distribution of the impurity cannot have the shape of the Fermi distribution function. Another distinct feature of the polaron in the presence of the bound state is almost linear dependence of the group velocity on Q for all values of the coupling strength, Fig. 3(b). That is, the impurity can be viewed as a free particle having the effective mass m * , and its momentum is P imp Q/m * . This is also seen from perturbative calculations, Ref. [49], not limited to the exactly solvable case considered in our paper.
We have used the exact wave functions and spectrum of the model. In a number of papers the mobile impurity problem is investigated by using approximate wave functions. Being constructed from a few particle-hole excitations, Refs. [50,51], these functions predict rather accurately some static properties, Ref. [52], and time dynamics, Fig S4 in Ref. [37], of the mobile impurity in one dimension. The momentum distribution function has not been treated using the aforementioned basis of the variation functions, to the best of our knowledge. Other natural ways to construct variation functions, by taking solely a product of coherent states [53][54][55], including Gaussian state correlations between different momentum modes [56], or correlations to an arbitrarily high order [57], are also promising. How to perform a resummation of the excitations containing arbitrary number of the particle-hole pairs for a weak impurity-gas coupling is discussed in Refs. [58][59][60].
An exciting development of ultracold atomic physics made it possible to setup experiments on diffusion and drag of quantum impurities embedded in a degenerate ultracold gas. Special to one dimension is the observation of the Bloch oscillations of a mobile impurity moving through a quantum fluid in the absence of a periodic lattice [61]. The momentum distribution function of the impurity has been measured in that experiment. However, there the impurity neither started out in the equilibrium ground state, nor reached such a state in the course of the temporal evolution. A The L → ∞ limit of Eq. (173): transformation to Eqs. (60) and (66) In this appendix we explain how we arrive at the Fredholm determinant representation (60) and (66), valid for N → ∞, starting from Eq. (173), valid for any finite N . Recall that we are working at a finite gas density, therefore N → ∞ implies L → ∞.
Combining the definitions (30) and (65) we write The L → ∞ limit of Eq. (32) reads for the real k 1 , . . . , k N +1 . This way, we get the kernel (61) Recall that Z b < 0. The L → ∞ limit of k N and k N +1 is given by Eq. (49). The leading term in the large L expansion of Eq. (160) in the interval 0 ≤ y ≤ L is and Substituting equation (49) into (32) we obtain in place of Eq. (178) for j = N, N + 1. Further, we limit y to the interval 0 ≤ y ≤ L/2, which implies e(k + ) = 0 for Eq. (181). This gives and for the L → ∞ limit of the functions e ± defined by Eq. (170). Evidently, Therefore, we get for the L → ∞ limit of the function (169) and and For the diagonal terms we use Eq. (165) and combine it with Eq. (168). This gives Using the identity where ξ is a number, and R is a rank one matrix, we write Eq. (173) as The two last rows and columns of this matrix are special because k N and k N +1 are complex: Here, e − (k j )e − (k l ) π , j, l = 1, . . . , N − 1, and where We calculate the determinant and the inverse of D omitting the terms which are higher than the first order in ξ: We use this identity for the determinant (194), where D is given by Eq. (198). We have This gives [BD −1 C] jl = 2π L e − (k j )e − (k l ) π e −ik − y −f 1 (k j )e ik l y − e ik j y f 1 (k l ) + 2ye ik j y e ik l y where f 2 (q) = αe ik − y − (α − 2y)e iqy − f 1 (q).
This is the desired representation (66).

B Small distance expansion of the reduced density matrix
In this appendix we derive formulas presented in section 6. We start from the finite-size expression for the reduced density matrix, Eq. (173). This way the repulsive ground state, the attractive gas state, and the attractive bound state are treated all at once. The expansion of the kernels (169) and (172) up to order three in y reads 2π L ∂k j ∂Λ ∂k l ∂Λ −1/2 K(k j , k l ) = −α + i(i + Λ)y − i 2 yα(k j + k l ) respectively. After substituting the expansions (218) and (219) into the determinants on the right hand side of Eq. (173) we use the following identity det N (I + U V T ) = det s (I + V T U ). (220) Here, U and V are N × s matrices with the columns formed by N + 1-component vectors u 1 , . . . , u s and v 1 , . . . v s , respectively. As a result (Mathematica package has been used to evaluate the determinants) we expanded Eq. (173) up to order three in y: We now take the thermodynamic limit in Eq. (222). For the repulsive ground state and the attractive gas state we have S 0 = Z, Here, ϕ = 1 2πα 2 ln 1 + (α − Λ) 2 1 + (α + Λ) 2 .
(227) and Z is given by Eq. (64). Notably, the result for the attractive bound state follows by just replacing Z with Z b , Eq. (71). Finally, the expansion (221) gives for Eqs. (98), (99), and (102) and respectively. This leads us to the results discussed in section 6.