Resonantly interacting $p$-wave Fermi superfluid in two dimensions: Tan's contact and breathing mode

TL;DR

This paper theoretically investigates the properties of a two-dimensional resonantly interacting p-wave Fermi superfluid, focusing on Tan's contact parameters and breathing mode frequency, revealing how effective range influences these observables near resonance.

Contribution

It provides the first detailed analysis of Tan's contact and breathing mode behavior in 2D p-wave Fermi gases, including the effects of finite effective range.

Findings

The many-body component of Tan's contact peaks above resonance.

Breathing mode frequency dips near resonance due to contact behavior.

Frequency deviation from scale invariance depends logarithmically on effective range.

Abstract

Inspired by the renewed experimental activities on $p$-wave resonantly interacting atomic Fermi gases, we theoretically investigate some experimental observables of such systems at zero temperature in two dimensions, using both mean-field theory and Gaussian pair fluctuation theory. These observables include the two $p$-wave contact parameters and the breathing mode frequency, which can be readily measured in current cold-atom setups with $^{40}$K and $^{6}$Li atoms. We find that the many-body component of the two contact parameters exhibits a pronounced peak slightly above the resonance and consequently leads to a dip in the breathing mode frequency. In the resonance limit, we discuss the dependence of the equation of state and the breathing mode frequency on the dimensionless effective range of the interaction, $k_{F}R_{p}\ll1$, where $k_{F}$ is the Fermi wavevector and $R_{p}$ is the effective range. The breathing mode frequency $\omega_{B}$ deviates from the scale-invariant prediction of $\omega_{c}=2\omega_{0}$, where $\omega_{0}$ is the trapping frequency of the harmonic potential. This frequency shift is caused by the necessary existence of the effective range. In the small range limit, we predict that the mode frequency deviation at the leading order is given by, $\delta\omega_{B}\simeq-(\omega_{0}/4)\ln^{-1}(k_{F}R_{p})$.

Equations

$$\mathcal{H}=\sum _{\mathbf{k}}{\xi }_{\mathbf{k}}{\psi }_{\mathbf{k}}^{†}{\psi }_{\mathbf{k}}+\frac{1}{2}\sum _{\mathbf{k},{\mathbf{k}}^{\prime },\mathbf{q}}{V}_{\mathbf{k}{\mathbf{k}}^{\prime }}{b}_{\mathbf{kq}}^{†}{b}_{{\mathbf{k}}^{\prime }\mathbf{q}},$$

$$V_{\mathbf{k}{\mathbf{k}}^{\prime }}=\lambda \Gamma \left(\mathbf{k}\right){\Gamma }^{\ast }\left({\mathbf{k}}^{\prime }\right).$$

$$\Gamma \left(\mathbf{k}\right)=\frac{\left(k/k_F\right)}{{\left[1+{\left(k/k_0\right)}^{2n}\right]}^{3/2}}{e}^{i{\phi }_{\mathbf{k}}}$$

$$2{ϵ}_{\mathbf{k}}{\Psi }_{\mathbf{k}}+\sum _{{\mathbf{k}}^{\prime }}{V}_{{\mathbf{kk}}^{\prime }}{\Psi }_{{\mathbf{k}}^{\prime }}=E_b{\Psi }_{\mathbf{k}},$$

$$\frac{1}{\lambda }=-\mathcal{P}\sum _{\mathbf{p}}\frac{{\left|\Gamma \left(\mathbf{p}\right)\right|}^2}{2{ϵ}_{\mathbf{p}}-E_b},$$

$$k^2\cot {\delta }_p\left(k\right)=-\frac{1}{a_p}+\frac{2k^2}{\pi }\ln \left(R_{p}k\right)+\cdots .$$

$$T\left(\mathbf{k},\mathbf{k};E\right)={\left|\Gamma \left(\mathbf{k}\right)\right|}^2{\left[\frac{1}{\lambda }+\sum _{\mathbf{p}}\frac{{\left|\Gamma \left(\mathbf{p}\right)\right|}^2}{2{ϵ}_{\mathbf{p}}-E-i{0}^{+}}\right]}^{-1},$$

$$\frac{1}{\lambda }+\mathcal{P}\sum _{\mathbf{p}}\frac{{\left|\Gamma \left(\mathbf{p}\right)\right|}^2}{2{ϵ}_{\mathbf{p}}-E}=\frac{M{\left|\Gamma \left(\mathbf{k}\right)\right|}^2}{4{\hslash }^2{k}^2}\left[\frac{1}{a_p}-\frac{2k^2}{\pi }\ln \left(R_{p}k\right)\right].$$

$$\mathcal{P}\sum _{\mathbf{p}}\frac{{\left|\Gamma \left(\mathbf{p}\right)\right|}^2}{2{ϵ}_{\mathbf{p}}-E}$$

$$\frac{1}{a_p}$$

$$R_p$$

$$\frac{1}{a_p}=\frac{M{E}_b}{\pi {\hslash }^2}\left[\ln \frac{M\left|E_b\right|}{{\hslash }^2{k}_0^2}+\frac{3}{2n}\right].$$

$$-\frac{1}{a_p}-\frac{2{\kappa }^2}{\pi }\ln \left(R_p\kappa \right)=0.$$

$${\hat{\rho }}_{\mathbf{q}}=\Delta {\delta }_{\mathbf{q},0}+{\Delta }_{\mathbf{q}}.$$

$${\Omega }_{\text{MF}}=\frac{1}{2}\frac{{\Delta }^2}{\lambda }+\frac{1}{2}\sum _{\mathbf{k}}\left({\xi }_{\mathbf{k}}-E_{\mathbf{k}}\right),$$

$${\left|u_{\mathbf{k}}\right|}^2$$

$${\left|v_{\mathbf{k}}\right|}^2$$

$$u_{\mathbf{k}}{v}_{\mathbf{k}}^{\ast }$$

$$\frac{1}{\lambda }+\sum _{\mathbf{k}}\frac{{\left|\Gamma \left(\mathbf{k}\right)\right|}^2}{2E_{\mathbf{k}}}=0,$$

$$n_{2D}=-\frac{\partial {\Omega }_{\text{MF}}}{\partial \mu }=\frac{1}{2}\sum _{\mathbf{k}}\left(1-\frac{{\xi }_{\mathbf{k}}}{E_{\mathbf{k}}}\right)\equiv n_F.$$

$$\Gamma\left(\mathcal{Q}\right)=-\left[\begin{array}[]{ll}M_{11}\left(\mathcal{Q}\right)&M_{12}\left(\mathcal{Q}\right)\\ M_{21}\left(\mathcal{Q}\right)&M_{22}\left(\mathcal{Q}\right)\end{array}\right]^{-1},$$

$$M_{11}\left(\mathcal{Q}\right)$$

$$M_{12}\left(\mathcal{Q}\right)$$

$${\Omega }_{\text{GF}}\left[\mu ,\Delta \left(\mu \right)\right]=\frac{k_{B}T}{2}\sum _{i{\nu }_n}\sum _{\mathbf{q}}\ln \det \left[-{\Gamma }^{-1}\left(\mathcal{Q}\right)\right].$$

$$2n_B=-\frac{\partial {\Omega }_{\text{GF}}\left[\mu ,\Delta \left(\mu \right)\right]}{\partial \mu }.$$

$$n_{2D}=n_F+2n_B.$$

$$-\left[\frac{{\tilde{\Delta }}^2}{4}-\tilde{\mu }\Theta \left(\tilde{\mu }\right)\right]+\frac{{\tilde{\Delta }}^2}{2}\ln \frac{(k_F{R}_p{)}^{-1}}{\sqrt{\frac{{\tilde{\Delta }}^2}{4}-\tilde{\mu }\Theta \left(-\tilde{\mu }\right)}}=1$$

$$\tilde{\mu }\ln \frac{(k_F{R}_p{)}^{-1}}{\sqrt{\frac{{\tilde{\Delta }}^2}{4}-\tilde{\mu }\Theta \left(-\tilde{\mu }\right)}}=1-\frac{\pi }{2}\frac{1}{k_F^2{a}_p},$$

$$\frac{\Delta }{{\epsilon }_F}\simeq \frac{2}{e}\left(\frac{1}{k_F{R}_p}\right)\exp \left[\frac{\pi }{2}\frac{1}{k_F^2{a}_p}\right].$$

$$\frac{\mu }{{\epsilon }_F}\simeq 1+\frac{\pi }{e^2}{\left(\frac{1}{k_F{R}_p}\right)}^2\frac{1}{k_F^2{a}_p}\exp \left[\frac{\pi }{k_F^2{a}_p}\right].$$

$$\frac{\Delta }{{\epsilon }_F}\simeq {\left[-\frac{1}{2}\ln \left(R_p\kappa \right)-\frac{1}{4}\right]}^{-1/2}.$$

$$\frac{{\mu }_B}{{\epsilon }_F}\simeq {\left[-\frac{1}{2}\ln \left(R_p\kappa \right)-\frac{1}{4}\right]}^{-1}\simeq {\left(\frac{\Delta }{{\epsilon }_F}\right)}^2.$$

Full text

I Introduction

The realization of Feshbach resonances in ultracold atoms provides a unique opportunity to explore fascinating quantum many-body phenomena Chin2010 . By precisely tuning the $s$-wave scattering length using an external magnetic field, one can now routinely produce a stable cloud of strongly interacting fermions and observe novel Fermi superfluidity at the cusp of the crossover from a Bardeen-Cooper-Schrieffer (BCS) superfluid to a Bose-Einstein condensate (BEC) Eagles1969 ; Leggett1980 ; Nozieres1985 ; SadeMelo1993 ; Regal2004 ; Zwierlein2004 ; Hu2006 ; Diener2008 ; Bloch2008 ; Giorgini2008 . The manipulation of high-partial-wave interatomic interactions is also possible. In particular, about fifteen years ago the experimental demonstration of $p$-wave Feshbach resonances in 40K and 6Li atoms Regal2003 ; Zhang2004 ; Gunter2005 ; Schunck2005 ; Gaebler2007 ; Fuchs2008 ; Inada2008 opened the exciting perspective of creating a topological $p$-wave Fermi superfluid, which hosts non-trivial non-Abelian excitations at its edges or in its vortex cores - the so-called Majorana fermions

• that could enable topological quantum computation Read2000 ; Ivanov2001 ; Kitaev2003 ; Nayak2008 . Unfortunately, unlike a strongly interacting Fermi gas at BEC-BCS crossover, the $p$-wave resonantly interacting systems generally suffers from a heating problem due to serious loss in atom number and can hardly reach a low-temperature equilibrium state. Therefore, there is no significant experimental progress, in spite of a lot of interesting theoretical investigations at the early stage Gurarie2007 , exploring different aspects of a strongly interacting $p$-wave Fermi superfluid, such as the zero-temperature phase diagram Botelho2005 ; Ho2005 ; Gurarie2005 ; Cheng2005 ; Iskin2006 ; Cao2013 , the superfluid transition temperature in three dimensions (3D) Ohashi2005 ; Inotani2012 ; Inotani2015 and the Berezinskii-Kosterlitz-Thouless (BKT) phase transition in two dimensions (2D) Cao2017 .

This situation is much improved over the past few years Luciuk2016 ; Waseem2017 ; Yoshida2018 ; Wassem2018 . After a quench in the external magnetic field to the resonance limit, a quasi-equilibrium state of a 3D strongly interacting $p$-wave Fermi gas has been observed Luciuk2016 , and the contact parameters, which characterize the universal short-distance and large-momentum behavior of the system Tan2008a ; Tan2008b ; Tan2008c ; Yoshida2015 ; Yu2015 ; He2016 ; Peng2016 ; Zhang2017 ; Yao2018 ; Inotani2018 , have been measured using radio-frequency (rf) spectroscopy Luciuk2016 . Most recently, the atom loss close to the resonance has been found to reduce significantly in lower dimensions Waseem2017 , as theoretically predicted Levinsen2008 ; Fedorov2017 . These experimental advances suggest the possibility of realizing a 2D strongly interacting $p$-wave Fermi superfluid in future experiments.

Motivated by this possibility, here we present a detailed theoretical study of two important experimental observables of a $p$-wave Fermi superfluid at zero temperature: the two $p$-wave contact parameters and the breathing mode frequency. The investigation is based on our recent results of the zero-temperature equations of state Hu2018 , which are reliably calculated using the Gaussian pair fluctuation (GPF) theory beyond mean-field Hu2006 ; Diener2008 ; Hu2007 ; He2015 . A finite-temperature investigation is also possible, by applying the Noziéres and Schmitt-Rink theory above the superfluid phase transition Nozieres1985 ; SadeMelo1993 . We note that, at sufficient high temperatures close to the Fermi degenerate temperature, the calculations of $p$-wave contact parameters and breathing mode frequency were recently performed by Yi-Cai Zhang and Shizhong Zhang Zhang2017 , using the virial expansion theory Yu2009 ; Liu2009 ; Hu2011 ; Liu2013 .

In this work, we are particularly interested in the breathing mode frequency right at the resonance. In three dimensions, an $s$-wave resonantly interacting Fermi gas acquires scale-invariant zero-energy wave-functions, which are eigenstates of the dilation operator Werner2006 . In the presence of an isotropic harmonic trap with frequency $\omega_{0}$, there is a hidden symmetry $SO(2,1)$, yielding a scale-invariant breathing mode frequency $\omega_{c}=2\omega_{0}$ Werner2006 . In two dimensions, this hidden symmetry was nicely explained by Pitaevskii and Rosch using the same contact $s$-wave interatomic interaction, which is scale-invariant classically Pitaevskii1997 . However, the quantum renormalization of the $s$-wave contact interaction necessarily introduces a new length scale of the 2D scattering length and explicitly breaks the scale-invariance of the interaction. Therefore, the breathing mode frequency deviates from the classically invariant value of $\omega_{c}$, i.e., $\delta\omega_{B}=\omega_{B}-\omega_{c}\neq 0$. This frequency shift is now referred to as quantum anomaly Hofmann2012 ; Taylor2012 ; Cao2012 ; Vogt2012 ; Mulkerin2018 ; Holten2018 ; Peppler2018 ; Hu2019 . In our case of a resonantly interacting $p$-wave interaction, where the 2D scattering area disappears, we anticipate that the system may also have scale-invariant zero-energy wave-functions if there is no length scale set by interactions, and in the presence of isotropic harmonic trap it has the scale-invariant breathing mode frequency $\omega_{c}$. This is unfortunately not true. The renormalization of $p$-wave interaction necessarily gives a length scale of the effective range of interactions. The breathing mode frequency then deviates from $\omega_{c}$. We find that this frequency shift in the $p$-wave channel is much larger than its $s$-wave counterpart of quantum anomaly and could be more easily measured in experiments.

The rest of the paper is set as follows. In the next section (Sec. II), we provide the Hamiltonian of a 2D spinless Fermi gas near a $p$-wave Feshbach resonance described by a separable interaction potential. We show how to calculate the scattering area $a_{p}$ and the effective range $R_{p}$ for the separable potential. As discussed in Sec. III, this enables us to obtain the equations of state of the system as functions of $a_{p}$ and $R_{p}$. We derive the analytic mean-field equations and present the numerical GPF results beyond mean-field for the equations of state. In Sec. IV, we discuss the two $p$-wave contacts and the related breathing mode frequency. In Sec. V, we focus on the resonance limit and discuss the significant frequency shift from the scale-invariant frequency $\omega_{c}$, due to the existence of the effective range. Finally, Sec. VI is devoted to conclusions.

II Model Hamiltonian and two-body scattering

A spinless 2D $p$-wave interacting Fermi gas of $N$ atoms can be described by the Hamiltonian (with the area $A=1$) Botelho2005 ; Hu2018 ,

[TABLE]

where $\psi_{{\bf k}}^{\dagger}$ ($\psi_{{\bf k}}$) is the creation (annihilation) field operator for atoms, $\xi_{{\bf k}}\equiv\epsilon_{{\bf k}}-\mu=\hbar^{2}\mathbf{k}^{2}/(2M)-\mu$ is the single-particle dispersion with mass $M$ and chemical potential $\mu$, and $b_{\mathbf{k}\mathbf{q}}^{\dagger}\equiv\psi_{{\bf k}{\bf+q}/2}^{\dagger}\psi_{-{\bf k}{\bf+q}/2}^{\dagger}$ is the composite operator that creates a pair of atoms with center-of-mass momentum $\mathbf{q}$. The inter-particle interaction takes a separable form with the chiral $p_{x}+ip_{y}$ symmetry Nozieres1985 ; Botelho2005 ; Ho2005 ; Hu2018 ,

[TABLE]

Here, $\lambda$ is the bare interaction strength and

[TABLE]

is a dimensionless regularization function with the cut-off momentum $k_{0}$, polar angle $\varphi_{{\bf k}}$ and exponent $n$ that is introduced for the convenience of numerical calculations. The Fermi wavevector $k_{F}$ is related to the number density of atoms $n_{2D}=N/A$ by the relation, $k_{F}=\sqrt{4\pi n_{2D}}$. Our choice of the chiral $p_{x}+ip_{y}$ channel is motivated by the phase diagram established by Gurarie et al. Gurarie2005 . Although experimentally the Feshbach resonances for $m=0$ and $m=\pm 1$ are nearly degenerate, at low temperature the system will spontaneously break the spin-rotational symmetry and condense into the $p_{x}+ip_{y}$ superfluid state.

We may replace $\lambda$ with a characteristic energy $E_{b}$ by solving the following two-body problem at zero center-of-mass momentum Botelho2005 ; Hu2018 ,

[TABLE]

where $\Psi_{\mathbf{k}}$ is the two-body wave-function in momentum space and $\mathbf{k}$ is the relative momentum of two particles. By using the separability of the interaction potential, after some algebra, it is easy to find that,

[TABLE]

where $\mathcal{\mathscr{P}}$ stands for taking Cauchy principal value. As we shall see in the next subsection, $E_{b}$ is related to the 2D scattering area $a_{p}$ (see Eq. (12) below). It can be either negative or positive Botelho2005 ; Hu2018 . In the former case, the so-called BEC side, it is simply the ground-state energy of a two-body bound state and the associated binding energy $\varepsilon_{B}=-E_{b}>0$. In the latter, it may be viewed as a scattering energy $E_{b}=\hbar^{2}\mathbf{k}_{b}^{2}/M>0$ for two particles colliding with a characteristic relative momentum $\mathbf{k}_{b}$ within the two-particle continuum. In this case, the sum at the right-hand-side of Eq. (5) is not well-defined and we have taken the Cauchy principal value of the sum to remove possible ambiguity. From now on, for convenience we name $E_{b}$ as the scattering energy, in spite of the fact that it can take negative values on the BEC side.

In the previous work Hu2018 , we determined the equations of state of the 2D $p$-wave Fermi superfluid as functions of the parameters ($E_{b},k_{0},n=1$). Here, for the purpose of calculating the $p$-wave contacts and breathing mode frequency, it is more useful to parameterize the inter-particle interaction by using the 2D scattering area $a_{p}$ and the effective range of interactions $R_{p}$, which are formally defined through the $p$-wave phase shift $\delta_{p}\left(k\right)$ Zhang2017 ; Levinsen2008 ; Hu2018 ,

[TABLE]

As shown in Appendix of our previous work Hu2018 , we find that $R_{p}=k_{0}^{-1}$ in the limit of $n\rightarrow\infty$. In the following, we derive the general expressions of $a_{p}$ and $R_{p}$ for an arbitrary exponent $n$. This is necessary, since we have to take a finite value of $n$ in actual numerical calculations. At low energy, all the physical results of interest should be functions of $a_{p}$ and $R_{p}$, independent of the different choice for $n$.

II.1 The expressions of $a_{p}$ and $R_{p}$

To relate the bare interaction strength $\lambda$ to the scattering parameters, we calculate the two-body $T$-matrix in vacuum Hu2018 ,

[TABLE]

where $E\equiv\hbar^{2}k^{2}/M$. Using the relation $T^{-1}({\bf k},{\bf k};E)=-M[\cot\delta_{p}(k)-i]/(4\hbar^{2})$, we find that in the limit $k\rightarrow 0$,

[TABLE]

As shown in Appendix A, for arbitrary exponent $n$ we have

[TABLE]

Therefore, we obtain

[TABLE]

In the limit $n\rightarrow\infty$, we recover the known relation $R_{p}=k_{0}^{-1}$ Hu2018 . For the expression of the scattering area $a_{p}$, we may replace $\lambda$ in favor of the scattering energy $E_{b}$. In the low-energy limit, i.e., $\left|E_{b}\right|\ll\hbar^{2}k_{0}^{2}/M$, we find that,

[TABLE]

It is easy to see that the scattering energy $E_{b}$ changes sign in the unitary limit $a_{p}=\pm\infty$. On the BEC side with $a_{p}>0$, we may write $E_{b}=-\varepsilon_{B}$, where $\varepsilon_{B}\equiv\hbar^{2}\kappa^{2}/M$ is the binding energy, and obtain

[TABLE]

This equation agrees with the low-energy expansion of the phase shift in Eq. (6), where $k=i\kappa$ is simply the pole of the $p$-wave scattering amplitude $f_{p}(k)=\sqrt{2/(\pi k)}[\cot\delta_{p}(k)-i]^{-1}$.

III Zero temperature theory

The zero-temperature mean-field and GPF theories of a 2D chiral $p$-wave Fermi superfluid were laid out in our previous work Hu2018 . Here, for self-containedness, we give a brief summary. In the superfluid phase, two fermions can pair up via the separable attraction $V_{\mathbf{k}\mathbf{k}^{\prime}}$ to form a Copper pair, described by a generalized density operator $\hat{\rho}_{\mathbf{q}}\equiv\lambda\sum_{{\bf k}}\Gamma^{*}(\mathbf{k})b_{\mathbf{kq}}$. The pairs then condense into the zero center-of-mass momentum state, as described by a nonzero pairing order parameter $\Delta$, i.e.,

[TABLE]

On top of this condensate are strong pair fluctuations, represented by the field operator $\Delta_{\mathbf{q}}$ for the non-condensed Cooper pairs.

Neglecting $\Delta_{\mathbf{q}}$ leads to the mean-field description. At a given chemical potential, the zero-temperature thermodynamic potential takes the form Hu2018 ,

[TABLE]

where $E_{{\bf k}}=[\xi_{{\bf k}}^{2}+\Delta^{2}\left|\Gamma({\bf k})\right|^{2}]^{1/2}$ is the energy of fermionic Bogoliubov quasi-particles. The associated quasi-particle wave functions are given by,

[TABLE]

By minimizing the mean-field thermodynamic potential with respect to $\Delta$ and $\mu$, we obtain the mean-field gap equation,

[TABLE]

and the mean-field number equation,

[TABLE]

The contribution of strong pair fluctuations to the thermodynamic potential can be accounted for, by taking an approximate Green function $\Gamma(\mathcal{Q}\equiv\left\{\mathbf{q},i\nu_{n}\right\})$ for non-condensed Copper pairs at the Gaussian level Hu2006 ,

[TABLE]

where the matrix elements are given by,

[TABLE]

$M_{21}(\mathcal{Q})=M_{12}^{*}(\mathcal{Q})$, and $M_{22}(\mathcal{Q})=M_{11}^{*}(\mathcal{Q})$. Here, $\nu_{n}\equiv 2n\pi k_{B}T$ with integer $n=0,\pm 1,\pm 2,\cdots$ are bosonic Matsubara frequencies, and the abbreviations $u_{\pm}\equiv u_{{\bf q}/2\pm{\bf k}}$, $v_{\pm}\equiv v_{{\bf q}/2\pm{\bf k}}$, and $E_{\pm}\equiv E_{{\bf q}/2\pm{\bf k}}$ are used. At the Gaussian level, the effective interaction between non-condensed Cooper pairs is treated within the Bogoliubov approximation, so there is no residual interaction between bosonic quasi-particles. Therefore, it is straightforward to write down the fluctuation part of the thermodynamic potential for non-interacting quasi-particles Abrikosov1963 ,

[TABLE]

For a given $\mu$, once $\Omega_{\textrm{GF}}$ is numerically calculated, we determine the number of Cooper pairs $n_{B}$ by using numerical differentiation,

[TABLE]

The number equation Eq. (20) is then updated to,

[TABLE]

This leads to an updated chemical potential in the GPF theory.

It is worth noting that within GPF the pairing gap $\Delta\left(\mu\right)$ is always calculated at the mean-field level, by solving the gap equation Eq. (19). This is necessary to ensure a gapless Goldstone phonon mode Hu2006 ; Diener2008 , i.e., $\det\Gamma^{-1}(\mathcal{Q}=0)=0$. In principle, it is possible to have a generalized approximation to improve the gap equation beyond mean-field. Accordingly, we could improve the vertex function $\Gamma(\mathcal{Q})$ beyond GPF. This possibility will be explored in future studies.

III.1 Analytic solutions from mean-field theory

In two dimensions, the integrals involved in mean-field equations can often be integrated out explicitly, leading to some nice analytic solutions. In the following, we take Fermi wave-vector $k_{F}$ and Fermi energy $\varepsilon_{F}=\hbar^{2}k_{F}^{2}/(2M)$ as the units of wave-vector and energy, respectively. In particular, we define the dimensionless pairing gap $\tilde{\Delta}=\Delta/\varepsilon_{F}$ and the dimensionless chemical potential $\tilde{\mu}=\mu/\varepsilon_{F}$. By setting $n\rightarrow\infty$ (i.e., taking a step-like function for the regularization function $\Gamma({\bf k})$) and performing the integrals in the gap and number equations, we arrive at two coupled equations:

[TABLE]

and

[TABLE]

where $\Theta(x)$ is the step function.

In the BCS limit, where $\tilde{\Delta}\rightarrow 0$ and $\tilde{\mu}\rightarrow 1$, we find from Eq. (28) that

[TABLE]

By substituting it into Eq. (27), we obtain

[TABLE]

In the BEC limit, the chemical potential becomes negative and approaches the half of the bound state energy, i.e., $\mu\rightarrow-\hbar^{2}\kappa^{2}/(2M)$, where $\kappa$ is the solution of Eq. (13). From Eq. (27), it is readily seen that,

[TABLE]

For the chemical potential, we rewrite it into the form, $\mu=-\hbar^{2}\kappa^{2}/(2M)+\mu_{B}/2$, where the molecular chemical potential $\mu_{B}$ is approximately equal to $g_{B}n_{2D}/2$ with $g_{B}$ being the strength of the interaction between two pairs. After some algebra, we find that

[TABLE]

The pair-pair interaction strength $g_{B}=\mu_{B}/(n_{2D}/2)$ is then given by,

[TABLE]

in agreement with the previous result (see Eq. (48) in Ref. Hu2018 ).

III.2 Numerical results on equation of state

At the level beyond mean-field, the GPF theory can only be solved numerically. In Fig. 2, we report the GPF chemical potential (a), total energy (b) and the pairing gap (i.e., the inset) as a function of the inverse scattering area $-1/(k_{F}^{2}a_{p})$ at a given effective range $k_{F}R_{p}=0.05$, using the black solid lines with circles. For comparison, we show also the corresponding mean-field results by the red dashed lines. These thermodynamic variables have been shown in the previous work as a function of the scattering energy $E_{b}$ Hu2018 .

Both chemical potential and total energy suppress significantly from their non-interacting values $\varepsilon_{F}$ and $E_{0}=N\varepsilon_{F}/2$, respectively. In particular, on the BEC side, we observe a flat molecular chemical potential and total energy, which are nearly independent on the scattering area $a_{p}$. As discussed in the previous work Hu2018 , this is an indication of the formation of an interacting Bose condensate of composite Copper pairs in two dimensions, with a constant pair-pair interaction strength $g_{B}\sim\hbar^{2}/M$.

III.3 Super-Efimov trimers

It is worth mentioning that, for 2D fermions with $p$-wave interaction, Nishida and co-workers recently discovered a series of three-particle bound states, namely super-Efimov states Nishida2013 . How would the many-body properties of the system (i.e., contact and breathing mode as addressed in this work) be affected by these super-Efimov trimers is an interesting research topic to explore Zhang2017SuperEfimov . Naïvely, due to the double exponential scaling of the super-Efimov trimers Nishida2013 , we anticipate that only one trimer with an emergent energy scale will appear under current experimental conditions. The neighboring trimer with smaller energy cannot exist due to its large spatial extent, while the one with larger energy is simply too deep to be experimentally observed. In this respect, the impact of super-Efimov states to the many-body physics could be less significant than that of conventional Efimov states.

IV Results and discussions

We are now ready to discuss the $p$-wave contacts and the related breathing mode frequency. There are two contact parameters, characterizing the short-distance and large-momentum behaviors of different correlation functions, such as momentum distribution and pair-pair correlation function Yoshida2015 ; Yu2015 . As shown by Yi-Cai Zhang and Shizhong Zhang Zhang2017 , these two contacts $C_{a}$ and $C_{R}$ satisfy the adiabatic relations,

[TABLE]

and therefore can be determined once the total energy is known at a given entropy $S$. At zero temperature, where the entropy is always zero, we simply take the two first-order derivatives.

IV.1 Tan’s $p$-wave contacts

For this purpose, we may write the zero-temperature total energy in a dimensionless form $\xi(x,y)$,

[TABLE]

and the two $p$-wave contacts can similarly be rewritten in the dimensionless way,

[TABLE]

where $\xi_{x}\equiv\partial\xi/\partial x$ and $\xi_{y}\equiv\partial\xi/\partial y$. Following the dimensionless form of the total energy, it is easy to find that the chemical potential $\mu=\partial E/\partial N$ and the pressure $P=\mu n_{2D}-E/A$,

[TABLE]

where $P_{0}\equiv n_{2D}\varepsilon_{F}/2$. By substituting the expressions of the dimensionless contact into the last equation for pressure, we obtain the pressure relation Zhang2017 ,

[TABLE]

It is useful to distinguish the two- and many-body contributions to the contact parameters. For the two-body contribution, we assume that the system can be viewed as an ideal, non-interacting gas of $N/2$ pairs, each of which has the energy,

[TABLE]

In other words, on the BEC side the pair takes the ground-state energy of the two-body bound state; while on the BCS side, the minimum energy of the pair should be zero (i.e., the lower threshold of the two-particle continuum). The two-body contribution to the total energy of the system can then be written as,

[TABLE]

On the BEC side, by using Eq. (13), we find that the two-body contribution to the contact from the energy $E_{2B}$, denoted by $C_{a,2B}$ and $C_{R,2B}$, is given by Zhang2017 ,

[TABLE]

On the BCS side, $C_{a,2B}=0$ and $C_{R,2B}=0$, as a result of $E_{2B}=0$. On both sides, either BEC or BCS, we obtain that,

[TABLE]

This equation is easy to understand from the pressure relation, since the two-body bound state does not contribute to the many-body observables such as pressure.

In Fig. 2, we plot the two dimensionless contact parameters as a function of $-1/(k_{F}^{2}a_{p})$ at a given effective range $k_{F}R_{p}=0.05$, calculated by using the GPF theory (black lines with circles) and the mean-field theory (red dashed lines). The two-body contribution from the bound state to the contacts is also shown by green dot-dashed lines. As we see in Fig. 2(a), the contact related to the scattering area $C_{a}$ is always positive. It increases with increasing interaction strength $(k_{F}^{2}a_{p})^{-1}$. On the BEC side with a positive scattering area, we find that the GPF result of $C_{a}$ is exhausted by the two-body contribution $C_{a,2B}$. The mean-field theory seems to under-estimate $C_{a}$, with the largest under-estimation occurs at $(k_{F}^{2}a_{p})^{-1}\sim 0.5$. On the other hand, the contact related to the effective range, $C_{R}$, has a non-monotonic dependence on the inverse scattering area, as shown in Fig. 2(b). As $(k_{F}^{2}a_{p})^{-1}$ increases, $C_{R}$ initially increases, reaches a maximum at $(k_{F}^{2}a_{p})^{-1}\sim-1$, and then decreases to zero at about the resonance limit. Towards the BEC limit, it decreases very rapidly. We find similarly that the GPF result of $C_{R}$ is almost exhausted by the two-body contribution $C_{R,2B}$. The mean-field theory generally over-estimates $C_{R}$ and the over-estimation becomes increasingly larger when we increase interaction strength. This is related to the unreliable prediction of the mean-field theory on the pair-pair interaction strength (see Eq. (33)).

We have separated out the many-body parts of the two contact parameters, $C_{a,MB}=C_{a}-C_{a,2B}$ and $C_{R,MB}=C_{R}-C_{R,2B}$, and show the GPF predictions in Fig. 3, for two effective ranges of interactions, $k_{F}R_{p}=0.05$ (black line with circles) and $k_{F}R_{p}=0.10$ (blue line with squares). On the BEC side, the many-body parts of both contact parameters are small, consistent with the observation in Fig. 2 that the contacts are exhausted by the two-body contribution. Across the resonance limit, they exhibit a pronounced peak. The peak in $C_{a,MB}$ is slightly above the resonance limit. The peak in $C_{R,MB}$ locates at $(k_{F}^{2}a_{p})^{-1}\sim-1$ and shifts towards the BCS limit with decreasing effective range. We note that the many-body parts of the two $p$-wave contacts are always positive.

IV.2 Breathing mode frequency

The interesting dependence of the many-body part of the contacts on the interaction strength may lead to a non-trivial breathing-type oscillation mode, when the 2D $p$-wave Fermi superfluid is confined in a harmonic trap with trapping frequency $\omega_{0}$. This is a mode excited by the perturbation $\lambda(t)\mathcal{O}\equiv\lambda(t)\sum_{i=1}^{N}r_{i}^{2}$, i.e., by slightly modulating the harmonic trapping frequency for a certain period. For non-interacting bosons or fermions, the breathing mode frequency is simply $\omega_{c}=2\omega_{0}$. The inter-particle interaction generally leads to a frequency shift. As shown by Yi-Cai Zhang and Shizhong Zhang Zhang2017 , the frequency shift at the leading order is proportional to the contact parameters. By using virial expansion, the frequency shift of the breathing mode at high temperatures was then theoretically studied Zhang2017 .

In our zero-temperature case, we calculate the breathing mode frequency using the well-known scaling approach Menotti2002 ; Hu2004 ; Hu2014 . This amounts to assuming a polytropic form for the pressure equation of state, $P\propto n_{2D}^{\gamma+1}$, where the polytropic index $\gamma$ may be calculated using

[TABLE]

at the center of the harmonic trap. The scaling approach then leads to a breathing mode frequency Hu2014 ,

[TABLE]

where at zero temperature we rewrite $\partial P/\partial n_{2D}$ in terms of the compressibility $\kappa_{T}=[n_{2D}^{2}(\partial\mu/\partial n_{2D})]^{-1}$ and $\kappa_{T}^{(0)}\equiv(n_{2D}\varepsilon_{F})^{-1}$ is its non-interacting value. This expression emphasizes the sound-wave nature of the breathing mode frequency. Qualitatively, the breathing mode frequency can be estimated as $c_{s}k_{\textrm{min}}$, where $c_{s}$ is the sound velocity and $k_{\textrm{min}}\propto A^{-1/2}$ is the minimum wavevector of the Fermi cloud with an area $A$. By recalling the relation $\kappa_{T}^{-1}\propto c_{s}^{2}$ and assuming the pressure $P\propto A^{-1}$ under the soft-wall confinement of the harmonic trap, we find $\omega_{B}^{2}\propto(P\kappa_{T})^{-1}$.

The polytropic index $\gamma$ can be directly calculated once the energy or pressure equation of state is known. By taking derivative with respect to density in Eq. (40) and neglecting all small second-order derivatives, we find that,

[TABLE]

where the bar over $\xi_{x}$ and $\xi_{y}$ indicates that we do not include the irrelevant two-body contribution. By substituting it into Eq. (47), we obtain the frequency shift $\delta\omega_{B}=\omega_{B}-\omega_{c}$,

[TABLE]

By replacing the two derivatives with the help of Eqs. (37) and (38), we finally arrive at

[TABLE]

where the virial theorem $PA=M\omega_{0}^{2}\left\langle\mathcal{O}\right\rangle/2$ in the presence of harmonic traps is used NoteVirialTheorem . We therefore recover Eq. (84) in Ref. Zhang2017 and explicitly show the relation between the many-body parts of the two $p$-wave contacts and the frequency shift in the breathing mode.

In Fig. 4, we present the frequency shift of the breathing mode as a function of the inverse scattering area, at two effective ranges $k_{F}R_{p}=0.05$ (a) and $k_{F}R_{p}=0.10$ (b). These results are calculated using Eq. (48) within the mean-field theory (dashed lines) and the GPF theory (lines with symbols). We find that the frequency shift is negative on the BCS side and exhibits a broad dip at $-(k_{F}^{2}a_{p})^{-1}\sim 1-2$. This dip structure is apparently related to the peak structure in the many-body part of the two $p$-wave contacts, according to Eq. (51). The two contacts contribute differently in opposite signs and the contribution from $C_{a,MB}$ seems to dominate. We note that, in a 1D harmonically trapped $p$-wave Fermi superfluid, the breathing mode frequency shows qualitatively similar dependence on the interacting strength in the weak-coupling regime Imambekov2010 ; Chen2016 .

On the BEC side, we see that the frequency shift predicted by the GPF theory becomes flat and small. This is associated with the formation of tight-binding molecules who interact via a nearly constant molecular scattering length, as we discussed earlier. The GPF frequency shift is positive and is about 5% at $k_{F}R_{p}=0.10$. In contrast, the mean-field frequency shift is negative and shows a non-trivial cusp at $(k_{F}^{2}a_{p})^{-1}\sim 0.5$. This mean-field behavior is unphysical, arising from the unreliable equations of state predicted by the mean-field theory. It is interesting to note that, the breathing mode frequency shift of a weakly-interacting 2D Bose gas was investigated both theoretically and experimentally Olshanii2010 ; Merloti2013 . In that case, the shift is too small to be experimentally observed. The moderately interacting molecular condensate formed in the strongly-interacting 2D $p$-wave Fermi superfluid could be a possible candidate to observe the breathing mode frequency shift due to beyond-mean-field effects.

V Frequency shift in the resonance limit

Here we focus on the breathing mode frequency in the resonance limit $a_{p}=\pm\infty$. If we neglect the dependence of the equations of state on the effective range, the dimensionless energy function $\xi$ is simply a constant. From the pressure $P=\xi P_{0}\propto n_{2D}^{2}$, we find a polytropic index $\gamma=2$ and hence $\omega_{B}=2\omega_{0}$. This could be an exact result ensured by the scale invariance of the system Werner2006 . However, the necessary existence of the effective range breaks the scale invariance and leads to a derivation of the breathing mode frequency away from the scale-invariant result of $\omega_{c}=2\omega_{0}$. A similar situation happens in an $s$-wave 2D Fermi superfluid Hofmann2012 ; Vogt2012 . While the superfluid with $s$-wave contact interaction is scale invariant in the classical treatment, i.e., the model Hamiltonian simply scales upon stretching the length of the system Pitaevskii1997 , the renormalization of the contact interaction necessarily introduces a 2D $s$-wave scattering length $a_{2D}$ that violates the scale-invariance. This leads to an up-shift in the breathing mode frequency, the so-called quantum anomaly, which is about 10% in the strongly-interacting regime Hofmann2012 ; Taylor2012 . It is not a surprise to see the similarity between the effective range $R_{p}$ in a $p$-wave Fermi superfluid and the 2D scattering length $a_{2D}$ in an $s$-wave Fermi superfluid. This is discussed in more detail in Appendix B.

V.1 Chemical potential in the resonance limit

Before we discuss the frequency shift in a resonantly interacting $p$-wave Fermi superfluid, it is useful to first understand the chemical potential in this limit. Using the mean-field equations Eqs. (19) and (20), we find that,

[TABLE]

By treating $-[\ln(k_{F}R_{p})]^{-1}=-y^{-1}$ as the small parameter, we obtain,

[TABLE]

Thus, towards the zero-range limit, the mean-field chemical potential at resonance vanishes linearly.

More accurate predictions from the GPF theory should be determined numerically. Empirically, we find that the GPF chemical potential at resonance can be nicely fitted by the formalism,

[TABLE]

where $A\simeq 0.078\ll 1$. While the GPF chemical potential at resonance still vanishes linearly in the zero-range limit, the slope (i.e., the value of $A$) is much slower than that of the mean-field chemical potential.

In Figs. (5) and (6), we report the mean-field and GPF predictions of the chemical potential at resonance as a function of $-\ln^{-1}(k_{F}R_{p})$, respectively. The analytic expressions and the empirical formalism discussed in the above are also shown. At small effective range, they agree well with the numerical results.

V.2 Frequency shift

The mean-field analytic expression and the GPF empirical formalism for the chemical potential at resonance are very useful to understand the shift of the breathing mode frequency. To see this, we may calculate the polytropic index related to the chemical potential, i.e., $\mu\propto n_{2D}^{\gamma}$, by using

[TABLE]

where we have rewritten $\mu=\varepsilon_{F}\xi^{(\mu)}(y)$ and have assumed $\xi_{y}^{(\mu)}\ll\xi^{(\mu)}$. By taking the chemical potential at the leading order, i.e., $\xi^{(\mu)}\propto-1/y$, we obtain immediately $\gamma\simeq 1-1/[2\ln(k_{F}R_{p})]$, and consequently,

[TABLE]

More careful treatments of the chemical potential to the next order in Eqs. (54) and Eq. (56) lead to the results,

[TABLE]

for the mean-field theory and

[TABLE]

for the GPF theory, respectively.

In Fig. 7, we show the up-shifts of the breathing mode frequency predicted by the mean-field theory and the GPF theory in the resonance limit, together with the asymptotic behaviors in the zero-range limit, as discussed in the above. According to the GPF theory, the shift of the breathing mode frequency can easily reach 10% at a relatively small effective range, i.e., $[\ln(k_{F}R_{p})]^{-1}\simeq-0.6$ or $k_{F}R_{p}\sim 0.2$.

At this point, it is interesting to compare the frequency shift exhibited by a resonantly interacting $p$-wave Fermi superfluid and by a strongly interacting $s$-wave Fermi superfluid, both in two dimensions. In the latter case, the theoretically predicted maximum quantum anomaly of 10% is yet to be experimentally confirmed Vogt2012 ; Holten2018 ; Peppler2018 . The main obstacle comes from the confinement-induced effective range $R_{s}$, which is significant under the current experimental condition. Indeed, our recent analysis indicates that the effective range $R_{s}$ in a 2D $s$-wave superfluid can strongly suppress the quantum anomaly down to 1-2% Hu2019 . In sharp contrast, for a resonantly interacting $p$-wave Fermi superfluid, the effective range $R_{p}$ enhances the frequency shift. Owing to the great feasibility in tuning $R_{p}$ in cold-atom experiment, therefore, we anticipate that a low-temperature $p$-wave Fermi gas at Feshbach resonances would be an ideal candidate to conclusively confirm the predicted frequency shift.

VI Conclusions and outlooks

In conclusions, we have theoretically determined two important experimental observables - Tan’s contact parameter and the breathing mode frequency

• of a resonantly interacting $p$-wave Fermi superfluid in two dimensions at the BEC-BCS evolution. Both observables can be easily accessed in current cold-atom experiment, as soon as a stable $p$-wave superfluid is realized in reduced dimensions. The two Tan’s contact parameters can be directly extracted from the tail of momentum distribution probed by radio-frequency spectroscopy Luciuk2016 , and the breathing mode measurement is now a routine tool in cold-atom laboratories Vogt2012 ; Holten2018 ; Peppler2018 .

We have proposed that, similar to an $s$-wave Fermi superfluid at the BEC-BCS crossover, the $p$-wave Fermi superfluid in the resonance limit experiences a frequency shift, due to the non-vanishing effective range of interactions that explicitly breaks the scale invariance of the system. The up-shift in the breathing mode frequency, away from the scale-invariant value $\omega_{c}=2\omega_{0}$, turns out to be significant. At the leading order, it is inversely proportional to the logarithm of the effective range. As a result of this slow-decay logarithmic dependence, the frequency shift can reach 5-10% over a wide range of the effective range.

Acknowledgements.

We thank Shizhong Zhang and Yi-Cai Zhang for stimulating discussions. This research was supported by Australian Research Council’s (ARC) Discovery Programs Grant No. DP170104008 (HH), Grant No. FT140100003 and Grant No. DP180102018 (XJL).

Appendix A An integral in the two-body $T$-matrix

Here we consider the integral,

[TABLE]

By introducing the variable $z\equiv(p/k_{0})^{2}$, we find that,

[TABLE]

where $z_{0}\equiv(k/k_{0})^{2}\ll 1$. To handle the operator $\mathcal{\mathscr{P}}$ for Cauchy principle value, we divide the whole integral into three parts $[0,z_{0})\cup[z_{0},2z_{0})\cup[2z_{0},\infty)$. Upon changing the dummy variable, the integral $\mathcal{I}$ can be rewritten in terms of $I_{1}$ and $I_{2}$,

[TABLE]

where

[TABLE]

and

[TABLE]

It is clear that $I_{1}=2z_{0}+o(z_{0})$. For $I_{2}$, by neglecting the higher contribution $o(z_{0})$, it can be separated into two parts,

[TABLE]

These two parts can be integrated out explicitly:

[TABLE]

Putting $I_{1}$ and $I_{2}$ together, up to the order $o(z_{0})$ we obtain the expression,

[TABLE]

which is Eq. (9) in the main text.

Appendix B Quantum anomaly in a strongly interacting $s$-wave Fermi superfluid

In an $s$-wave Fermi superfluid, Tan’s adiabatic relation is given by Werner2012 ,

[TABLE]

which takes exactly the same form as the adiabatic relation for the effective range of interactions, as given in Eq. (35), up to an unimportant pre-factor. This same form emphasizes the similar role played by the effective range $R_{p}$ in a $p$-wave Fermi superfluid and by the scattering length $a_{2D}$ in an $s$-wave Fermi superfluid.

Let us now write the total energy of the $s$-wave Fermi superfluid in a dimensionless form He2015 ,

[TABLE]

where for the two-component Fermi gas the Fermi wavevector $k_{F}=\sqrt{2\pi n_{2D}}.$ By using the adiabatic relation Eq. (70), we then find,

[TABLE]

where $\xi_{z}=\partial\xi/\partial z$. The dimensionless chemical potential and pressure are also easy to obtain,

[TABLE]

By calculating the polytropic index related to the pressure Hofmann2012 , we find that,

[TABLE]

where again the bar denotes the exclusion of the two-body bound-state contribution. The calculation of the polytropic index related to the chemical potential leads to the same expression at the same level of approximation. Thus, we obtain the quantum anomaly,

[TABLE]

According to the GPF calculation or quantum Monte Carlo simulations, at around the strongly interacting regime $\ln\left(k_{F}a_{2D}\right)\sim 0$, the many-body part of the contact shows a peak with $C_{MB}/k_{F}^{4}\sim 0.05$ He2015 . This is correlated with a total energy $(E+N\varepsilon_{B}/2)/E_{0}\sim 0.25$ He2015 . By using these two numbers, we find a quantum anomaly $\delta\omega_{B}/(2\omega_{0})\sim 0.1$ for a strongly interacting $s$-wave Fermi superfluid.