BCS in the sky: signatures of inflationary fermion condensation

We consider a Bardeen-Cooper-Schrieffer (BCS)-like model in the inflationary background. We show that with an axial chemical potential, the attractive quartic fermion self-interaction can lead to a BCS-like condensation. In the rigid-de Sitter (dS) limit of inflation where backreaction from the inflaton and graviton is neglected, we perform the first computation of the non-perturbative effective potential that includes the full spacetime curvature effects in the presence of the chemical potential, subject to the mean-field approximation whose validity has been checked via the Ginzburg criterion. The corresponding BCS phase transition is always first-order, when the varying Hubble is interpreted as an effective Gibbons-Hawking temperature of dS spacetime. In the condensed phase, the theory can be understood from UV and IR sides as fermionic and bosonic, respectively. This leads to distinctive signatures in the primordial non-Gaussianity of curvature perturbations. Namely, the oscillatory cosmological collider signal is smoothly turned off at a finite momentum ratio, since different momentum ratios effectively probe different energy scales. In addition, such BCS phase transitions can also source stochastic gravitational waves, which are feasible for future experiments.


Introduction
The phenomenon of superconductivity is one of the most important discoveries in modern physics.In 1911, H. Onnes first observed that the mercury electrical resistance drops to zero below 4K, an effect he called superconductivity [1].Since then, many insightful attempts have been made to give a theoretical understanding of superconductivity, notably the phenomenological London equations in 1935 [2] and the Landau-Ginzburg model in 1950 [3].The successful microscopic theory of superconductors was finally put forward by J. Bardeen, L.N. Cooper and J.R. Schrieffer in 1957.In the Bardeen-Cooper-Schrieffer (BCS) theory [4], the electrons moving in the ionic lattice can develop an attractive interaction by exchanging phonons.At low temperatures, this attractive force binds electrons into Cooper pairs.These bosonic Cooper pairs then condense into the BCS vacuum, spontaneously breaking electromagnetic U (1) EM symmetry, introducing the Meissner effect and a non-zero supercurrent.
At the core of the BCS theory is the instability of electrons near the Fermi surface.Namely, as long as the interaction is attractive, and the fermion number is supported by a finite chemical potential, the Fermi sea is unstable against the formation of Cooper pairs, no matter how weak the attractive interaction is.The universality of the BCS model assumption makes it a flexible framework to generalize to other contexts.Over the past decades, the BCS theory and its variations played crucial roles in many aspects of modern physics from condensed matter to fundamental high energy theory.To name a few, BCS-like models are used for QCD chiral dynamics [5][6][7][8] and color superconductivity [9][10][11], early universe model building [12,13], neutrino condensation [14], dark energy/dark matter candidate [15][16][17][18][19][20][21][22], and fermion dynamics around black holes [23] .
Cosmology, on the other hand, is in a sense very much akin to condensed matter physics.For instance, Lorentz symmetry is always spontaneously broken (by the lattice in condensed matter systems and by the cosmic expansion in cosmology).There is, in both cases, a separation of a non-trivial homogeneous background and inhomogeneous perturbations.The observables both entail correlation functions of perturbations.In fact, both condensed matter systems and cosmology can be classified according to a common set of symmetry breaking patterns [24].Such a link between the two systems motivates us to look for the analogy of BCS condensation in a cosmological setup.
For instance, let us consider the fermion dynamics during inflation.It has been shown that the rolling inflaton can naturally provide an external chemical potential that supports a number of fermions [25][26][27][28][29][30].Then after turning on a weak attractive interaction, the fermions may pair up and condense, forming a BCS-like state.The crucial difference here is, of course, that the universe is expanding and the BCS mechanism must be modified by the spacetime curvature effects.We will see that the spacetime curvature can be effectively considered as a thermal effect, with a temperature set by the Hubble parameter during inflation.From the experience of condensed matter systems, we know that superconductivity typically requires low temperature.If the spacetime expansion is too fast, the effective temperature is too high, and the Cooper pairs may be teared apart, breaking the condensate.It is then possible to find a "critical temperature" for an inflationary BCS system.We also note that in the previous analysis of BCS in cosmology [12][13][14][15][31][32][33][34][35][36][37][38], the curvature effect has not yet been fully studied, and only leading-order results in Hubble constant expansion have been achieved so far [15,32,33].In the dS-invariant case, all-order results have been achieved but in order to generate a non-zero gap, either one needs to go to lower-dimensional spacetime [34,35] or introduce super-Planckian physics [12].It is thus doubtful whether any exponentially small gap can truly exist in an expanding spacetime.
To answer this theoretical question and explore the phenomenology of BCS in the sky, we set out in this paper to compute the effective action of a class of BCS-like models during inflation.We find an analytical result for the effective potential exact to all orders in the Hubble expansion.Inspecting this effective potential, we show that the BCS phase transition is always first-order under the naturalness assumption, suggesting the removal of any exponentially small gap by cosmic expansion.Since any expanding flat FRW universe can be locally approximated as dS when the Hubble parameter is slowly varying, our analysis can be thought of as an approximate IR completion of BCS in any FRW spacetime.In addition, we show that the squeezed limit of primordial bispectrum probes different phases of the BCS theory, leading to distinctive cosmological collider [39][40][41][42] signals that are feasible for future non-Gaussianity surveys [43][44][45][46][47][48][49].If the phase transition happens during inflation, the bubble collision process gives rise to stochastic gravitational waves that can be probed by a variety of detection methods, such as SKA [50], LISA [51] and DECIGO [52].
This paper is structured as follows.In Sect.2, we start with a brief review of BCS in the LAB and try to generalize it to a cosmological model, while giving some qualitative understandings.In Sect.3, we quantitatively analyze our BCS model first in flat spacetime and then in dS.After computing the non-perturbative effective potential, we examine the phase transition properties in Sect. 4. The impact on observables is discussed in Sect. 5.At last, we conclude and give future prospects in Sect.6.

Notations and conventions
Throughout the paper, we use the (−, +, +, +) metric sign convention.The Pauli matrices with Lorentz indices are defined as σ µ = (1, σ), σµ = (1, −σ), where σ are the usual three Pauli matrices.The Dirac matrices and its algebra are given by γ µ = 0 σ µ σµ 0 , γ 5 = −1 0 0 1 , {γ µ , γ ν } = −2η µν , {γ 5 , γ µ } = 0 . (1.1) We will also adopt the Van der Waerden notation for two-component spinors and the conventions follow from [53] and [54].In the dS context, we use dimensionless tilde variables for any parameters measured in Hubble units, for example, κ ≡ κ/H, ∆ ≡ ∆/H, μ ≡ µ/H, etc.A prime on quantum average values denotes the 2 From superconductors to the universe In the original BCS theory, the conducting-band electrons move in an ionic lattice and couple to the quanta of lattice vibration known as the phonon.In the zero temperature and decoupling limit, the electrons fill up the lowest energy states and form a sharp Fermi sphere.However, the coupling to the phonon introduces a non-local attractive interaction for pairs of electrons with opposite spin (s-wave pairing channel).Consequently, it becomes energetically favorable for the electrons to pair up and redistribute themselves, causing the spreading of the Fermi sphere.The width of the spread is given as the gap: where ω D is the Debye frequency, g is the coupling strength and N is the density of states near the Fermi surface.
The BCS ground state is then a two-mode squeezed state where ↑↓ denotes the spin-z eigenvalue, not the helicity.Any fermionic excitation on this ground state is gapped by ∆ in the energy spectrum.Since the electrons carry electromagnetic charge, this condensate of electrons necessarily breaks U (1) EM and leads to superconductivity.In contrast, in a cosmological setup, the BCS model structure must be slightly modified.This is because the universe must stay electrically neutral in the absence of any positively charged lattice.Therefore it is not possible to have a condensate of fermions charged under U (1) EM 1. Instead, we require fermions and anti-fermions to pair up.Consider a massless Dirac fermion model, where / ∇ = γ µ ∇ µ is the covariant derivative with the spacetime spin connection.This fermion may be understood as a new degree of freedom active during inflation, for example those from Grand Unification model building.Now the global U (1) V symmetry is exact, and the vector fermion current J µ V = J µ R + J µ L is conserved, where J µ R and J µ L are the right-and left-handed fermion current.The global U (1) A symmetry, however, may not be exact.This corresponds to the non-conservation of the axial current, Namely, there will be an asymmetry in the number of Cooper pairs with different constituent fermion helicities.To achieve this, we can introduce a parity-odd chemical potential via the operator where κ is the axial chemical potential.This non-zero chemical potential κ supports a net axial charge, while keeping a neutral vector charge, Here +, − denotes the helicity while N, N denotes the fermion and anti-fermion, respectively.In the absence of any helical chemical potential [55][56][57], the net helical charge is also zero, we have Combining (2.5), (2.6) and (2.7), we see that Thus our system contains an equal number of fermions and anti-fermions with the same helicity, but there is an imbalance of particle numbers for different helicities.Note that another technical reason for choosing an axial chemical potential as opposed to a vector chemical potential is that fermion chemical potential does not produce particles dynamically during inflation, and can be rotated away as shown in [58].In other words, a net vector fermion number density is not preserved by inflationary expansion, and if such a state is prepared initially, the inflationary expansion will dilute the vector fermion number away and drive the state toward the dS-invariant attractor solution of a Bunch-Davies vacuum.In contrast, an axial chemical potential is compatible with Bunch-Davies vacuum and cannot be rotated away when the condensation forms.It can assist dynamical particle production in a dS-invariant way, keeping a constant net axial fermion number density when the production rate and dilution rate achieve a balance.
We are now ready to introduce the attractive self-interaction that gives rise to BCS condensation.The simplest four-fermion interaction that binds fermions with anti-fermions is given by (2.9) The plus sign in the front guarantees the attractive nature of the interaction, since having a larger fermionic field bilinear decreases the potential energy.Thus fermion fields tend to be drawn closer.Such an interaction can originate from the exchange of a heavy scalar particle of mass scale M , or from the s-wave channel of the torsion-induced four-fermion interaction [59,60].Since the energy scale of inflation is not far below that of Grand Unified Theories (GUT) such as the Pati-Salam model [61][62][63][64][65], one can treat the four-fermion term as EFT operators of a broken GUT group, in a similar way that these operators appear in the SMEFT [66].For simplicity, we shall first focus on (2.9) and comment on other types of interaction in Sect.3.3.We also note that different from the non-local interaction of electrons in real BCS systems, the quartic self-interaction here appears to be local.Of course from an EFT perspective, local/non-local is a relative notion depending on the hierarchy between the energy scale of the probe and the underlying mechanism.The UV completion of L int can actually be a non-local exchange of massive scalars at the scale of M .This attractive interaction responsible for condensation is then an effective description at the gap energy scale ∆.As long as ∆ ≪ M , one can integrate out the heavy meson, and the non-local UV picture can be reduced to a local vertex in the IR.After introducing this self-interaction, one expects the fermions and anti-fermions of the same helicity will form Cooper pairs in the s-wave channel (see FIG. 2 for a cartoon depiction of the physical picture).Now although the helicity imbalance is still present, it becomes hidden in the internal structure of Cooper pairs, and is not visible in the IR picture at leading order.This fact will affect the parity property of the resulting cosmological collider signals, as we will demonstrate later.In summary, our BCS Lagrangian reads where g µν = a(τ ) 2 η µν is the background spacetime metric and D µ is the covariant derivative.For clarity, we will analyze our BCS model progressively, first in flat spacetime and then in dS and actual inflation.

Flat spacetime case
Here, we review the case of flat spacetime that has already been well studied, setting the footage to the dS case.
Reducing the covariant derivative to normal derivatives, the Lagrangian reads In the full theory partition function, the four-fermion interaction term can be resolved to a quadratic interaction by introducing an auxiliary field ∆, where we have introduced as the partition function for the fermion field with "effective mass" ∆.Then the full partition function can be written as a bosonic path integral weighed by an effective action of the auxiliary field ∆, Notice that this effective theory of the bosonic field ∆ should be understood as the IR face of the original UV theory of interacting fermions.
To determine the vacuum state of the effective theory in the IR, it suffices to take the mean-field description of the auxiliary field ∆(x) = ∆ ≈ const and integrate out the fermions, so that the partition function can be easily computed as an infinite determinant as or Finishing the loop momentum integral in (3.6) and inserting it back to (3.4), we obtain the effective potential defined as the leading order effective action in the uniform limit, with Here Λ is an explicit loop momentum cutoff, which should be replaced by renormalized parameters and scales after renormalization.However, we are going to consider embedding the BCS model in a time-dependent cosmological set-up later, and renormalization in a time-dependent system (such as cosmology) is by no means a simple task.For example, it is not clear to us what renormalization conditions should be posed to eliminate the cutoffs (Λ) and bare variables (M ), nor at what time should they be imposed.One possible choice is to fix the second-and fourth-order derivative at the gap: where µ, λ are some constants to be measured by experiments.Of course, the question in cosmology is that unlike in particle scatterings where pole masses and threshold couplings are known from experiments, there is no such measurement that tells us the value of µ, λ at ∆ 0 .In addition, some parameters entering the renormalization conditions are time-dependent, so it is not useful to stick to a fixed condition like (3.9).As a result, in favor of practicality, we do not pose any explicit renormalization conditions, but understand the cutoffs Λ, M as defining the renormalization conditions.Namely, we trade the functional dependence by and study the effective potential by varying M, Λ as an alternative to µ, λ.We also remark that the effective potential V eff is formulated as a tool to determine the ground state of the system, which is why we care about the functional dependence on the gap ∆.However, just like the Coleman-Weinberg effective potential, V eff is always defined up to a ∆-independent O(κ 4 ) constant.This O(κ 4 ) constant remainder does not affect the determination of the gap, but can be understood to represent the free energy of the free relativistic Fermi gas.
The lowest energy state of the theory is determined by the minimum of the effective potential.This gives the gap equation, The location of the minimum ∆ = ∆ 0 corresponds to the gap of BCS condensation.A way to see this is through computing the average value of the fermion bilinear in the UV theory of fermions, it is equal to the average value of the gap field in the IR theory of bosons: where in the third step we have integrated by parts.For this reason, we call the auxiliary field ∆ the gap field, and regard it as the field of Cooper pairs, which is the analogy of BCS ground state (2.2).The k ̸ = 0 modes are understood as bosonic excitations of the ground state.We will discuss these perturbations and their interesting phenomenology later in Sect. 5.
To solve for the gap, we use the approximation ∆ 0 ≪ M, Λ, κ, then the gap equation (3.11) yields a very similar result as the BCS gap (compare with (2.1)), Naturalness requires M ∼ M are of a similar magnitude.Perturbative unitarity of the four-fermion interaction also yields the condition Λ ≲ πM .In general, the chemical potential must also be below the cutoff, hence the hierarchy κ < Λ ≲ πM .As apparently seen from (3.13), the gap decreases with decreasing chemical potential and weaker coupling.Still, the gap equation always admits an exponentially small solution.In other words, condensation can appear no matter how small the chemical potential is or how weak the coupling is.This is, of course, an artifact of our zero temperature assumption.At a finite temperature, we expect thermal fluctuations to destabilize the condensation and unbound the fermions, leading to the disappearance of the gap [67].In cosmology, in addition to thermal effects present in the early universe, the expansion of spacetime is another important factor.It may thus be dangerous to rely on flat spacetime BCS computations in an expanding universe.Namely, the spacetime curvature may easily remove the exponentially small gap, unbinding the Cooper pairs.To fully understand the BCS condensation in an expanding universe such as inflation, we need to perform an "IR completion" of the BCS model regarding the effects of spacetime curvature.

Going to dS
Now let us move on to dS spacetime.The background geometry is given by where −∞ < t < ∞ is the cosmic time and −∞ < τ < 0 is the conformal time.Note that spacetime expansion is more transparent in the FRW coordinates (t, x), but the conformal coordinates (τ, x) are technically cleaner.Therefore, we will mostly work with the conformal coordinates in the computations below.An empty dS spacetime is maximally symmetric with 10 isometries.However, to implement BCS condensation, the 3 dS boost symmetries must be broken by the chemical potential.Hence, we are left with 7 symmetries, including 3 spatial translations, 3 rotations and 1 dilation: (τ, x) → (λτ, λx).These are the symmetries that constrain the form of our effective action, as we will see later.
We start by spelling out the covariant Lagrangian in the conformal coordinates: Here the covariant derivative is ∇ µ ≡ ∂ µ + i 4 ω µab σ ab , with the spin connection ω µab = −e ν b ∇ µ e c ν η ac .The dS vierbein can be chosen to be simply e a µ = aδ a µ , e µ a = a −1 δ µ a .The chemical potential is naturally introduced through the dimension-5 coupling to a rolling scalar background ϕ = ϕ 0 (t) (usually the inflaton during inflation).Λ c is an new physics scale arising from inflaton-fermion dynamics and is independent of the scale M of fermion self-interactions.It characterizes the strength of the Yukawa-like coupling between the inflaton and the fermion.Then the rolling speed of the scalar ∂µϕ Λc = a(τ )κδ 0 µ provides the chemical potential, where φ0 ≡ ∂ t ϕ 0 ≈ const.The size of the chemical potential is bounded by κ 2 < φ0 < Λ 2 c due to perturbative unitarity [26,58].Following the treatment in the preceding section, we introduce the gap field to resolve the quartic self-interaction, and define an effective action for the gap field by integrating out the fermions, (3.17) with and Unlike in flat spacetime, tracing over the fermion spectrum in dS is more non-trivial.We need to inspect its dynamics more clearly before performing the path integral over ψ.Due to the loss of time translational symmetry and the dS boosts, it is technically more convenient to work with the two-component-Weyl-spinor representation of the Dirac field, In the mean-field limit, ∆ ≈ const, the fermions can be quantized and expanded into orthogonal modes with different spatial momenta and helicities, with the spinor-helicity basis The orthonormal condition is The (anti)-commutation relations are while all other anticommutators vanish.The mode functions can be solved exactly in the form of the Whittaker W -function2 where we have denoted By either the IR expansion [25] or the Stokes-line method [29], one can find the occupation number for these fermions to be Now we are in a position to compute the effective potential.To reiterate, we wish to compute In the uniform limit, the effective action is simply a potential function integrated over spacetime volume: with Now take a derivative with respect to the gap field ∆, we obtain where ⟨• • • ⟩ f stands for the average value computed using the theory action S f [ψ] in (3.19).Switching to the 2-component spinor variables, the potential derivative becomes We can change the loop momentum variable into a dimensionless variable that corresponds to the physical momentum measured in Hubble units, The potential derivative now takes a more compact form, 4π 2 e −πκ I + + e πκ I − + c.c. , (3.34) 2Notice that two of our mode functions differ from those of [28] by a sign.This is because our "Dirac mass" is −∆ in (3.19).This difference is purely a matter of convention, and physical results do not depend on it.
where the key integrals are given by In the high-energy limit with z → ∞, we have , suggesting these integrals (and the potential derivative itself) are quadratically divergent, as one should expect from the flat spacetime loop integral.As in the flat spacetime case, we can keep the explicit cutoff but understand it as an alternative way of expressing the renormalized quantities.Thus the divergent part of the loop integral is trivial in this sense.As we shall see later, the non-trivial aspect of the dS loop integral is hidden in the finite part, encoding the information about the non-perturbative spacetime curvature effects that is not visible at any orders in the H n -expansion.Using a regularization trick based on the Whittaker product integral of Arthur Erdelyi [68], we are able to compute the finite remainder of the loop integral exactly for all parameter sizes.We spell out the detailed computation of the loop integrals I ± in Appendix A, and quote the final result of the effective potential here: The analytical expression for the finite remainder f ϵ ± and the coefficient functions A, B are given by (A.13) and (A.15) in Appendix A. Notice that the only divergent dS correction is the logarithmic term ∼ H 2 ∆ 2 ln Λ in the first line, all other corrections come from the finite remainder in the second line.
It is useful to take an H → 0 limit of V eff to examine the corrections of slow spacetime expansion.Note that the finite part of the effective potential contains terms non-perturbative in the small-H expansion, as seen from the trigonometric functions in f ϵ ± .To extract the leading non-perturbative correction, we expand the trigonometric functions into exponential functions and only keep terms with the leading exponent.For example, with μ = Expanding to the leading non-trivial order in this way, we obtain Here the first line reproduces the flat-spacetime result in (3.8) after taking H → 0. The second and third lines give perturbative corrections analytic in H, while the fourth and fifth lines correspond to the non-perturbative corrections that are non-analytic in H. Inspecting the form of the O(H 2 ) corrections, we find that there is a superficial dip dominating the small-gap region where ∆ 2 ≪ Hκ.Naively, this dip seems to suggest that the effective potential is unbounded from below in the gapless regime, and that no stable BCS condensation can occur for an arbitrarily small expansion rate.However, one must recall that the expansion (3.39) is organized in powers of H, and it would be inconsistent to take ∆ 2 ≪ Hκ without considering all the higher order contributions.
In fact, a full inspection using the exact effective potential (3.37) shows the effect of such terms is actually to lower the potential barrier at the origin.Therefore, raising the Hubble parameter leads to the formation of a true vacuum at ∆ = 0, and eventually to a symmetry-restoring phase transition and the disintegration of BCS condensation.We will come back to the details of BCS phase transition in Sect. 4. Note that with the benefit of hindsight, we have dressed the non-analytic Hubble corrections in (3.39) in the form of the dS background temperature Thus the non-analytic pieces are essentially suppressed by the corresponding Boltzmann factors.Take the first non-analytic piece for example, its Boltzmann factor reads which is the occupation number of the fermionic excitation upon the BCS ground state (more precisely, that of the enhanced left-handed fermions/anti-fermions as in (3.27)).Here µ = √ ∆ 2 + κ 2 is the effective mass of the fermions.From the analogy of BCS theory at a finite temperature, we recognize this non-analytic piece as the on-shell particle production effect of the dS background temperature.Of course, one should not take the thermal analogy too seriously, because there are power-law terms such as H n ∼ T n dS in (3.39), which should be absent in a gapped thermal system at low temperatures.These are more naturally viewed as generic spacetime curvature corrections that are not necessarily dS-like.

Backreaction on inflation
At the end of this subsection, we point out that although non-perturbative information is encoded in (3.37), we are nevertheless working in the mean-field approximation (whose validity will be checked later in Sect.5.1) and have so far neglected the backreaction on the inflaton and graviton sector.Thus (3.37) is exact in a restricted sense.To estimate how much the BCS matter sector backreacts to the inflationary background, we can compare the typical energy deposit in the fermion condensation and the vacuum energy stored in the inflaton background.The former is estimated from the change of effective potential (3.39) during a typical gap field traverse δ ln ∆ = O(1), where we have chosen the second term in the expanded effective potential (3.39) as a benchmark of its dependence on the gap.In slow-roll inflation, the vacuum energy in the inflaton background is given in terms of the Hubble parameter and the Planck scale as where Here in the second step we have used ∆ 0 ≲ κ/2, which is found to be satisfied in most of the parameter region of the condensed phase.The last step comes from requiring κ 2 < φ0 < Λ 2 c , which is itself a consequence of perturbative unitarity [26,58].Henceforth we see that the phase transition in the fermion sector is negligible compared to the background inflaton dynamics, consistent with neglecting backreaction.

Other self-interactions
We have so far focused on the four-fermion interaction ( ψψ) 2 in the s-channel pairing.In principle, one can also consider other types of fermion self-interactions.At the level of dimension-6 operators, the most general set of Lorentz-invariant and CP T -invariant operators that couple fermion and anti-fermions is given by squaring the five basic bilinears3: ( ψψ) 2 , (i ψγ 5 ψ) 2 , ( ψγ µ ψ) 2 , ( ψγ µ γ 5 ψ) 2 , ( ψσ µν ψ) 2 . (3.45) However, due to Fierz identities, only three of the five combinations are linearly independent [69,70].We are free to choose the scalar, pseudoscalar and vector combinations to span the space of quartic interactions.We need to keep the spatial isotropy in any BCS condensations, so the vector condensate will be excluded4, leaving us with only the scalar and pseudoscalar combinations.Note that in the inflationary context, we do not necessarily impose so many symmetries.The rolling speed of the inflaton spontaneously breaks the boosts (potentially also P and T ) and introduces a distinction between the temporal and spatial components of fields.However, we do not wish to introduce such a distinction at the level of dimension-6 operators, because this would reintroduce the explicit inflaton derivative, i.e. ( ψγ 0 ψ) 2 ⊂ (∂ µ ϕ ψγ µ ψ) 2 , which dramatically increases the mass dimension of the operator.Henceforth, we will simply generalize our model to include the pseudoscalar self-coupling, In the equal-coupling case with M = M 5 , the theory enjoys a chiral U (1) A symmetry: ψ → e iαγ 5 ψ, (3.47) ψ → ψe iαγ 5 . (3.48) This U (1) A will be spontaneously broken whenever the fermions form a BCS condensate.Now, we need to introduce two auxiliary scalar fields to resolve the self-interactions with and In the uniform limit, Σ(x), Π(x) = const, we can always perform a U (1) A rotation ψ → e −iΘγ 5 /2 ψ to rotate away the mass phase, 3The Lorentz structure dictates the contraction of Lorentz indices, while the mixing between scalar (vector) and pseudoscalar (pseudovector) is forbidden by P .4In principle, one can also introduce a vector BCS condensate while preserving isotropy, by including an SU (2)I internal symmetry [24,71].Then the expectation value of the order parameter breaks the SU (2)I × SO(3) to the diagonal subgroup of effective rotations, which is still a symmetry.We leave this interesting scenario to future studies.
In doing so, one naturally performs a field redefinition Now the fermion partition function Z f takes exactly the same form of (3.19).The effective IR theory is thus turned into where We have also introduced the auxiliary ghost fields c, c to replace the Jacobian The ground state solution is obtained at the stationary points of the path integral (3.55).For the non-dynamical ghost fields, this is simply attained at c = c = 0.A subsequent variation with respect to the uniform gap fields ∆, Θ then selects the BCS vacuum, according to the minima of the effective potential where V 1−loop (∆) is identical to the case computed in the previous section.Inspecting the dependence of the effective potential on the coupling strengths, we find three limiting scenarios: 1. M = M 5 .In this equal-coupling case, the effective potential reduces to that of (3.30) in the single-coupling case, which is a function of the amplitude field ∆ only.The phase field Θ now enjoys a shift symmetry Θ → Θ + C, a legacy of the UV chiral U (1) A .Upon acquiring a non-zero average value with ∆ 0 ̸ = 0, the phase mode Θ becomes the massless Goldstone boson of the spontaneous broken shift symmetry.
2. |M − M 5 | ≪ M .As one deviates from the equal-coupling limit, the phase mode acquires a soft mass proportional to the coupling difference: 3. M ≪ M 5 (or M ≫ M 5 ).When there is a hierarchy between the couplings, the gap fields automatically align along the maximal symmetry breaking direction.In other words, the phase field becomes frozen at Θ = 0 (Θ = π/2) and decouples from the spectrum due to its heavy mass.Thus the theory reduces to the single-coupling model in the previous section.

BCS phase transition during inflation
Equipped with the exact formula of the effective potential, let us now analyze the parameter dependence of this effective potential and the phase diagram of the inflationary BCS condensation.Without loss of generality, we will focus on the single-coupling model (3.15), and the result can be directly translated to the more general two-coupling model (3.46).The effective potential V eff (∆) is dependent on four parameters: the chemical potential κ, the Hubble scale of inflation H, the coupling strength M , and the cutoff scale Λ which is understood as playing the role of mass renormalization parameter.Without referring to any external energy scales such as the Planck scale M p ∼ 10 18 GeV, the property of the BCS system only depends on the relative size of these four parameters.We need to determine which of them should be fixed and which of them are to be considered tunable.
In a realistic inflationary spacetime, the background has a weak time dependence characterized by the slow roll parameters ε, η, introduce a weak time dependence in the chemical potential and the Hubble parameter, where t * is a reference time during inflation.In the slow-roll phase of inflation, ε, η > 0, thus the chemical potential tends to increase, while the Hubble parameter tends to decrease during inflation.With this observation in mind, we will mainly study our BCS model dependence on the chemical potential and Hubble parameter, while holding the coupling and cutoff fixed.In all cases, the Hubble scale turns out to be the smallest scale of the problem, thus we will typically quantify other scales in units of Hubble.When we need to vary the Hubble scale, we will use a reference Hubble scale H * to quantify other dimensional parameters.Before a quantitative discussion, recall that in the flat spacetime case, the gap equation always admits a finite but exponentially small solution.Therefore, the system is always in the condensed phase.However, with the spacetime expansion effects of inflation, it is not true that the gap equation always admits a non-zero solution.We have seen in the previous section that the dS correction looks thermal at the leading order, with a Hubble temperature T dS ∼ H.For BCS superconductors, when the temperature rises above a critical temperature T c , superconductivity will be lost and the BCS condensation will disintegrate Borrowing this analogy, during inflation, this should happen at and the explanation here is that spacetime expansion tears the Cooper pairs apart, destroying the condensate.This estimation gives a loose lower bound on the gap size, ∆ 0 ≳ H/(2π).We will see later the actual lower bound is tighter.We plot the effective potential and the gap for different sizes of the chemical potential and Hubble parameter in FIG. 2. Clearly, at a constant chemical potential, the gap only exists when Hubble is below a critical value.Meanwhile, for a fixed Hubble parameter, the chemical potential has to be large enough to support the existence of the gap.More interestingly, we find that the gap changes in a discontinuous way, and that the phase transition is first order, in agreement with our expectation above.We are also able to draw the phase diagram in FIG. 3, where the normal phase and condensed phase are discontinuously separated by a curve that extends all the way from the origin to the cutoff.
In fact, if we assume several reasonable conditions Naturalness: we can give a numerical proof that the phase transition is always first order.Here naturalness means the loop correction to the mass operator M 2 ∆ 2 is not much larger than the tree-level value5, typically not by a factor of order two for M .Also, the fermion self-interactions must be established as an EFT operator above the Hubble 5As mentioned in Sect.3.1, Λ is understood as a physical parameter determined through implicit renormalization conditions.Therefore it is not meaningful to directly compare Λ to M and conclude about EFT validity.Instead, one should look at the loop correction of the potential in the large-∆0 regime and require it to be less than the tree-level contribution.Therefore, EFT validity is still preserved as along as Λ < πM .This leaves room for Λ being greater than M without breaking EFT validity.On the other hand, if one goes beyond Λ = πM , the effective potential is overturned and a tachyonic instability appears signalling the breakdown of EFT.To avoid fine tuning, we have also imposed a naturalness bound here to avoid being dangerously close to instability.scale.The proof is simple: second-order phase transitions happen when the second-order derivative of the effective potential vanishes, However, it is a straightforward although technical task to show that the above equation has no solution under the set of reasonable conditions (4.4).More explicitly, we can numerically solve the minimum of the second-order derivative of the effective potential, min κ,M,Λ,H∈(4.4) Therefore, V ′′ eff (0) = 0 has no solution in the domain (4.4) defined by naturalness and EFT validity, and second-order phase transition is impossible.Note that the above proof is a sufficient condition for first-order In the exact flat spacetime case with H = 0, the system always stays in the condensed phase, while an exponentially small Hubble expansion rate can destabilize the condensate.phase transition, but is sensitive to the naturalness bound tightness.It will be interesting to investigate the phase transition type more systematically under even more relaxed constraints in future works.
At last, we briefly comment on the two-coupling model case.As mentioned in the previous section, the order parameter Θ automatically picks the smaller of {M 2 , M 2  5 } to align along.All the discussions above are then directly translated via the replacement M 2 → min{M 2 , M 2 5 }.

Observational signatures
Our BCS model in inflation also brings characteristic observational signatures in the sky (hence the title).First, in the condensed phase, the excitation modes of the BCS vacuum couple to the inflaton through the chemical potential operator, and lead to cosmological collider signals that probe the particle spectrum of the system.Such oscillatory signals can be observed in the squeezed limit of the primordial curvature bispectrum.Second, if the universe starts out in the normal phase, due to the slow time dependence of the inflaton rolling speed, the increasing chemical potential tends to drive the fermions to the condensed phase via a first-order phase transition.Such a phase transition is accompanied by bubble nucleations and collisions, producing stochastic gravitational waves.In this section, we will directly adopt the two-coupling model and investigate these two potential observables.

Collective modes
Up till now, we have been focusing on the background solution of the BCS model, namely, that of the vacuum selection.The state in reality is, however, not the exact BCS vacuum.One can view this as due to the on-shell particle production in an expanding spacetime.We have seen that the non-perturbative correction to the effective potential does include Boltzmann-suppressed real-particle effects.The collective modes describing fluctuations of the Cooper-pair condensate in magnitude and in phase are known as the amplitude mode (or Higgs mode) and phase mode (or Anderson-Bogoliubov mode), respectively [72][73][74][75].
To analyze these collective modes, we need to perturb the state around the BCS vacuum, and look at the effective action at next-to-leading order in the gradient expansion, The free dynamics of these modes will then be dictated by the quadratic Lagrangian expanded from the above effective action.Since here we are interested in probing the perturbation dynamics, we will assume a hierarchy between the gap and the Hubble scale, ∆ 0 ≫ H such that the flat spacetime approximation is reliable.This will significantly simplify our calculation.Due to the boost-breaking chemical potential, we do not expect Z µν ∆ , Z µν Θ ∝ g µν , and these scalar modes will a non-unit sound speed in general.In the flat spacetime limit O(H 0 ), the field strength renormalization coefficients are easily solved by computing 1-loop diagrams, Here the subscript stands for a = ∆, Θ. Regulating the logarithmically divergent integrals using a hard cutoff, we find (see more details in Appendix B) Inserting this back to the effective action, we perturb the vacuum by the canonically normalized modes, where σ and θ indicate the amplitude and phase modes, respectively.Θ 0 = 0 (Θ 0 = π/2) is the vacuum angle when M ⩽ M 5 (M > M 5 ).Expanding to quadratic order and restoring the scale factors adiabatically, we obtain This quadratic action is useful for reading out the sound speed and the mass information.Interestingly, we find that in the small gap limit with κ ≫ ∆ 0 > H, the sound speeds universally approach to that of a relativistic fluid, This is a reminiscence of the sound speed of non-relativistic BCS excitations in superconductors [73,76].The mass of the collective modes is given by (5.10) As mentioned in Sect.3.3, in the equal-coupling case, m 2 θ = 0, and the phase mode θ is the massless Goldstone of U (1) A symmetry breaking.Otherwise, the θ mode acquires a mass proportional to the coupling difference.Notice also that the amplitude mode mass √ 6∆ 0 is slightly larger than that of the non-relativistic case (2∆ 0 ), although the fermionic excitations still have a gap of ∆ 0 .Due to the cosmic expansion, the gap solution cannot be arbitrarily small.We find in most parameter regions, it is of order ∆ 0 ≳ O(10)H.Thus the amplitude mode is very heavy in Hubble units.This suggests a negligible cosmological collider signal of the amplitude mode.The phase mode can be light in the case of soft U (1) A breaking, but it turns out that its coupling to inflaton is in general rather weak, hence also a small signal strength, as we will see later.
A check on the mean-field approximation The calculation of gap equation and the phase analysis above are performed in the mean-field approximation where the gap field is assumed to take a uniform background with small fluctuations.Now that the linear dynamics of the fluctuations σ have been understood, we can estimate how well the mean-field approximation works.More specifically, we need to check the Ginzburg criterion for reaching a stable condensed vacuum, (5.11) Namely, we require the standard deviation of the quantum fluctuations of the order parameter to be smaller than its mean value.To the leading order in gradient expansion, the two-point function (power spectrum) of σ is given by that of a canonical massive scalar field in de Sitter spacetime6 [77,78], (5.12) Thus the left-hand side of (5.11) yields where we have applied (5.3) and the hierarchies H ≪ κ, H ≲ ∆ 0 .Therefore we conclude that the Ginzburg criterion is indeed met, and that the mean field approximation is valid whenever the gap is allowed to form.

Cosmological collider signals
Now let us turn on the inflaton fluctuations and examine the observational consequences of the collective modes.
Our goal is to estimate the resulting curvature bispectrum

.14)
6The two-point function of quantum fields is usually singular in the coincident limit and needs renormalization.The leading contributions to ⟨σ 2 (x)⟩ are of order O(m 2 σ ) and O(H 2 ), but they are renormalization-scheme dependent since they come from UV divergences.The leading contribution that does not depend on renormalization scheme turns out to be of order O(H 4 /m 2 σ ).It is straightforward to see that even if we evoke the renormalization-scheme dependent leading order contribution ⟨σ 2 (x)⟩ ∼ m 2 σ /(16π 2 ), the Ginzburg criterion is still met due to the hierarchy ∆ 2 0 < κ 2 .
where P ζ ≃ 2 × 10 −9 is the curvature power spectrum, f N L is the signal amplitude and S(k 1 , k 2 , k 3 ) is the normalized dimensionless shape function.In the cosmological collider paradigm [39][40][41][42], the squeezed limit of the shape function shows an oscillatory behavior where ω directly gives the mass of the particle exchanged by inflatons, hence cosmological "collider physics".
As mentioned in Sect.3.1, our model in the BCS condensed phase can be viewed from two perspectives at different energy scales.
• At energy scales above the gap E > ∆ 0 , it is more natural to take the perspective of a fermionic theory.
In this UV picture, the degrees of freedom are carried by interacting fermions with an effective mass generated by the condensate, as described by S ′ f [ψ] in (3.51) with ⟨Σ + iγ 5 Π⟩ = ∆ 0 e iγ 5 Θ 0 .
• At energy scales below the gap E < ∆ 0 , it is more natural to take the perspective of a bosonic theory.In this IR picture, the degrees of freedom are carried by the collective modes, as described by The inflaton fluctuation in the cosmological collider plays the role of a detector with variable energy scales.This is because in an expanding spacetime, the physical wavelength The effective energy scale it probes therefore scans through the theory it couples to.
Since physical resonances between the inflaton and different theory particles happen at different energy scales, varying the comoving momentum ratio7 is therefore tantamount to varying the center-of-mass energy of a particle collider.Hence the cosmological collider can in principle scan through the particle spectrum during inflation, i.e. a tomography in energy scales.In particular, for the inflaton 3-point function, the equilateral limit probes higher-energy physics, while the squeezed limit probes lower-energy physics8.This idea is further illustrated by the cartoon diagram in FIG. 4. We shall see this more clearly when we discuss the shape function below.More explicitly, we need to find the inflaton φ coupling to different fields in both the UV picture and the IR picture.The curvature correlator is then translated from inflaton correlators via a gauge transformation (5.16) Naturally, we expect the axial chemical potential to introduce the inflaton coupling via the Stückelberg trick Therefore, the inflaton perturbation φ(x) couples to fermions in the UV picture via Fermions must form closed loops to contribute to inflaton correlators.In the IR picture, φ(x) couples to both the amplitude mode and the phase mode.To find the explicit couplings, we need to expand S eff in powers of φ after doing the Stückelberg trick.We find the leading interactions to be 7It is crucial that only the momentum ratio encodes the physics at different energy scales, since the overall size of the momenta is irrelevant due to the scale invariance.8One can also view this cosmological collider tomography from the perspective of the cosmological flow [79].The evolution of the equal-time correlator shows that the formation of cosmological collider signals starts out in the equilateral limit, and extends to the squeezed limit later.This suggests that signals at more squeezed momentum ratios form at lower energy scales due to the redshift.Here the vertical axis is the conformal time τ , and the horizontal axis is the physical distance x ph ≡ a(τ The black curve traces the inflaton mode with a specific comoving momentum k 3 .The stars correspond to different resonant decay events producing the φ(k 1 ), φ(k 2 ) modes.As the universe expands, the energy scale of φ(k 3 ) decreases, crossing from the UV regime to the IR regime, where it effectively drives different massive degrees of freedom.Now changing the comoving momentum of the other two inflaton modes is equivalent to changing the resonant decay time of the massive degrees of freedom, which is also equivalent to changing the energy scale of the massive degrees of freedom.In this way, the more squeezed the momentum ratio is (k 1 → ∞), the later resonance happens, and the lower energy scale is probed by the resonance. with where we have assumed M = M 5 = 40H, Λ = Λ σ = Λ θ = 100H, κ = 45H for a typical numerical estimate.Notice that by construction (5.17), the inflaton φ is always derivatively coupled.This prevents a tadpole correction for the inflaton background.The phase mode θ is protected by a Z 2 ⊂ U (1) A symmetry whether U (1) A is explicit or not, and always couples to other fields in pairs.This fact brings two implications.First, there can be no mixing between the phase mode with the inflaton.Second, the phase mode contribution to inflaton correlators must starts at loop-level at leading order.The amplitude mode σ, however, is not protected by any symmetry.It indeed mixes with the inflaton and contributes to the inflaton correlators at tree-level.The leading order bispectrum comes from a fermion loop in the UV picture, a amplitude mode exchange and phase mode loop in the IR picture (see FIG. 5).To estimate the amplitude of signals, we follow the same strategy as in [58], Propagators are usually suppressed by 1/(mass).Additionally, the cosmological collider signals always contain the Boltzmann suppression factor e −π(mass)/H , attributed to either the resonance saddle point or from the particle vacuum production rates [80].The fermion loop contributes in the UV a signal strength [58] 1) , (5.22) where the Boltzmann factor is from the particle production rate (3.27).In the IR, the amplitude mode gives a negligible tree-level exchange signal e −πmσ/H ∼ O(10 −30 ). (5.23) Here the mass power enhancement is from the three-point vertex, and can be checked through previous works [80,81].The phase mode loop also yields a negligible signal, ( In all the above estimates, we have assumed Therefore, as we expected, the fermions yield observably large cosmological collider signals, while the collective modes are hardly excited and give essentially no cosmological collider signal.Since different momentum ratios in the squeezed limit probe physics at different energy scales, we expect the oscillatory signal of the UV fermions to damp out beyond a scale k 1 /k 3 ≳ (• • • ) × ∆−1 0 .To better understand this behavior and find out the waveform of the signal, we will adopt the saddle point method to be developed in details in a forthcoming paper [82].Our only assumption is that the condensation scale ∆ 0 , serves as a natural boundary between two distinct "phases"(more precisely, two distinct pictures of the same phase), and that unbounded fermions can only observed in the UV phase where E ∼ k a > ∆ 0 .First, let us start be noticing that cosmological collider signals measure the total phase cumulated between different events [80].For fermions with a comoving momentum close to that of φ(k 3 ) in the UV, the first relevant event is its particle production, at a time τ * given by [29] Production: The UV fermions then propagate and interact, before condensing at a time τ c with Condensation: The last possible event is the resonant decay of fermions into inflaton pairs with momenta k 1 , k 2 .To determine more precisely the resonance time τ • , we note that the dynamical phase of the fermion mode function takes the form e −i w(τ )dτ with dispersion relation w 2 = (k 3 − κa) 2 + ∆2 0 a 2 , while the inflaton oscillates as e ik 12 τ .Then the resonance event happens when the total dynamical phase strikes a saddle point, Resonance: The factor of two in the fermion dynamical phase arises from the loop where two internal fermion propagators are soft.We have also approximated the hard fermion in the loop as an EFT contact term.The cosmological collider signal records the cumulated phase from the production and the resonance, while the condensation time sets a cutoff for the latest resonance.Therefore, the necessary condition for the existence of the cosmological collider signal is the time hierarchy (5.28) This gives a momentum ratio cutoff for the cosmological collider signal, (5.29) For a momentum ratio beyond this cutoff scale, the resonance event cannot happen because the fermions would have paired and condensed due to the attractive self-interaction.The exact expression for the correlator is hard to calculate even at the tree-level [83][84][85] or for a bubble loop diagram [86], not to mention the triangle loop diagram here.Thus here we simply estimate the cosmological collider signal by the cumulated phase We plot this approximate waveform in FIG.6 with a schematic window function to indicate the damping effect beyond r max .We stress that this is can only be understood as a very rough description of what the true signal might be, due to the various approximations made above.A more refined computation is left for future works.
1 r max 10 50 Figure 6.The schematic cosmological collider signal in our BCS model.Here we have chosen the parameters M = M 5 = 40H, Λ = 100H, κ = 45H.We have applied an artificial window function centered around r max to smoothly turn off the fermion signals in the squeezed limit.We caution the readers that this figure is only a rough indication of the features in the true signal, by applying a window function by hand to connect the UV and the IR.We leave a detailed calculation to future works.

Gravitational waves from phase transition
The time dependence of the chemical potential κ = φ0 (t)/Λ c directly enters the shape of the effective potential.Therefore, if the fermion theory starts out in the normal phase with ⟨ ψψ⟩ = 0, the acceleration of the rolling inflaton may trigger a first-order phase transition that is associated with BCS condensation ⟨ ψψ⟩ = 2M 2 ∆ 0 .The bubble collisions during the phase transition then become the source of stochastic gravitational waves.
To estimate the gravitational wave energy density today, we follow the treatment of [87][88][89] and examine the changing rate of the bubble Euclidean action.In the Euclidean IR effective theory, the leading order action reads (5.31) Before solving for the bounce solution, we canonicalize the gap field by introducing and rewrite the action as with c 2 σ ≡ −Z 00 ∆ /Z 33 θ ≈ 1/3 and ∆ r being an arbitrary reference scale.To solve for an O(4)-invariant bounce, we rescale the spatial coordinates by defining x ≡ c s x and absorb the non-unit sound speed, (5.34) Now S4 gives the O(4)-invariant bounce action in the thin-wall limit as Here is the one-dimensional instanton action and δV eff is the vacuum energy difference.During the phase transition, the typical change of vacuum energy can be estimated as (5.37) The ρ field traverse is typically δρ = O(1)κ/ √ 6π.Therefore, we approximate .38)This leads to an estimated bounce action and its time derivative where η is the second slow-roll parameter.The spectrum of gravitational waves today can be written as [87][88][89]] Here Ω R is the present-day radiation abundance, while ρ inf is the inflaton potential energy density.The definition of S is provided in [87] and the flat space spectrum is given by [90] (5.43) Despite the constraints imposed by experimental data on the magnitude of the slow-roll parameter at CMB scales, it is worth noting that the generation of gravitational waves during inflation may occur in later stages of the inflation, where the slow-roll parameter ε may take a relatively larger value (e.g.ε ∼ 0.1), then the density ratio δρ r can be as large as O(10 −3 ) .where N e is the e-folding number before the end of inflation.The experiment curves are from IPTA [91], SKA [50], as well as space-based detectors such as LISA [51], DECIGO [52], TianQin [92], Taĳi [93].The gray line is dashed to indicate a risky parameter choice ε = 0.5.
Unfortunately, for ε ≲ 0.1 in general, we do not expect the signal strength (brown solid line) to reach the sensitivity of LISA or DECIGO.
We plot the resulting gravitational waves signals in the Fig. 7.The e-folding number is chosen as N e = 37 for the frequency region relevant to Pulsar Timing Array (PTA) experiments, whereas for space-based detectors such as LISA, it is set to N e = 23.We see that the gravitational waves associated with the BCS phase transition is feasible for the SKA experiment, while the detection probability for space-based detectors is less optimistic, unless the parameters are tuned to risky values such as ε ∼ 0.5.

Conclusion
The BCS theory of superconductivity has long been a thriving paradigm that inspired many excellent works and exciting discoveries across the vast territories of modern physics.It is particularly interesting to consider its generalization to cosmology due to the similarity between the expanding universe and condensed matter systems.In this work, we considered embedding a class of BCS-like models with chemical potential into inflation, and studied in detail the properties of the fermion condensate and its phenomenology.Starting from the singlecoupling case, we computed analytically the effective potential for BCS condensation in dS non-perturbatively in the spacetime curvature, subject to the decoupling limit M p → ∞ and neglecting backreaction from the inflaton and graviton sector.We found that the spacetime curvature effect strongly impacts the fermion condensation.Namely, the spacetime expansion may rip Cooper pairs apart and forbid the formation of an exponentially small gap as would otherwise be present in the flat spacetime approximation.During the slow-roll phase of inflation, BCS condensate can form with a gap at least a few times the Hubble scale, with the assistance of chemical potential.Assuming some reasonable constraints, we were able to prove that the corresponding BCS phase transition is always first-order, and obtained the phase diagram of our system.Most features remain the same after we generalize to the two-coupling model, where an additional phase mode appears in the IR EFT.The analysis of the collective modes shows that they share a universal sound speed of relativistic fluids, and masses proportional to the gap.We checked the Ginzburg criterion and showed the validity of the mean-field approximation.Our BCS model in the sky also yields interesting (although challenging at the same time) observational signatures via the cosmological collider tomography or stochastic gravitational waves induced by the first-order BCS phase transition.There are certainly many aspects for further investigations in the future: • Throughout this work, we have limited ourselves to the adiabatic case with a slowly-varying yet non-zero Hubble parameter.Such is the case for inflation.However, in later stages of the universe, the Hubble parameter is itself time-dependent.The adiabatic result with an instantaneous Hubble parameter can only serve as a crude leading order approximation.It would be interesting to give a better analysis of BCS condensations in a general FRW spacetime with significant deviations from dS.
• As mentioned in Sect.3.3, we have excluded the simplest vector condensate model due to the spatial isotropy.
However, in more complicated cases where the internal symmetries mix with spacetime symmetries, it is possible to form a vector condensate while maintaining the (effective) spatial isotropy.We leave this interesting scenario for future studies.
• For simplicity, we only focused on the essential part of the BCS model in this work, namely its pairing of fermions into Cooper pairs and their condensation.The gauge fields have not yet been introduced, and there is no Higgs phase of superconductivity.Hence our model in this work resembles more that of a fermionic superfluid.It is thus interesting to introduce charged fermion condensate and gauge fields along with their spontaneous symmetry breaking, thereby achieving a cosmological superconducting phase.This is left for future work.
• Our proof for the BCS phase transition being first-order is made under a naturalness assumption.This is reasonable from a practical point of view, yet it may also be possible to relax such constraints and check the transition type in more general cases.In other words, if we relax certain constraints, are there any critical points where the phase boundary terminates, and a second-order phase transition takes over?
• Despite the arguments built upon strong physical intuitions, we did not give any precise computation of the cosmological collider signal and its attenuation behavior when crossing into the IR regime.Thus an important question for us to work on in the future is how to acquire a more quantitative grasp on this tomography of cosmological colliders.
where G 22  33 is the Meĳer G-function [94].To proceed, we deform I − by (A.9) The expansion of generalized hypergeometric function at unit argument is [95] Γ(a 1 )Γ(a 2 ) which is valid if s = p i=1 b i − p+1 i=1 a i is equal to an integer t such that −t ≤ 0. For a detailed expression of all types of coefficients, we refer readers to [95] Notice that now the finite remainder f ϵ − (∆, κ, H) is a known function.One can also compare (A.11) with (A.4) and see that the divergence structure is similar.In particular, the coefficients for the logarithmic divergence are identical.Repeating the procedure for f ϵ + and sum over helicities, we obtain .12)where the last term enjoys an exact analytical expression, Here γ is the Euler constant and ψ(z) = Γ ′ (z)/Γ(z) is the digamma function.As we have mentioned before, different regularization methods should yield the same physical result after renormalization.Any constantmismatch before renormalization should be regarded as the mismatch of divergent terms.To match with the Λ-reg result, we require (A.6) and set As a small caveat, this only recovers the leading order of A, B in H.But practically it does not matter, because one always needs to perform the renormalization.In addition, it does not influence the shape of the effective potential by much.One could have started with ϵ-reg directly.We prefer the Λ-reg only because it is conceptually simple and convenient for the numerical check.In summary, our exact result for the effective potential is To confirm its validity, we numerically compute the loop integrals, and find the result being consistent with our analytical result, as shown in FIG. 8.

B Field strength renormalization
In order to calculate the 1-loop diagrams in (5.2), we first determine the fermion propagator with chemical potential, which can be obtained by inverting its equation-of-motion, Integrate out p 0 first and set a cutoff for the spatial momentum integral, we obtain

Figure 1 .
Figure 1.A cartoon illustration of the Cooper pair constituents in our BCS-like model.Black and gray circles represent fermions and anti-fermions, with their momenta indicated by the arrow and spin labeled with respect to the momenta.Each Cooper pair consists of a fermion-anti-fermion pair with opposite momenta and same helicity.Due to the axial chemical potential, there are more Cooper pairs built from left-handed fermions (red dashed circle) than that from right-handed fermions (blue dashed circle).

Figure 2 .
Figure 2. Left column: The effective potential for different choices of chemical potential (first row) and Hubble parameter (second row), with the other parameters fixed to be M = 40H, Λ = 100H.Right column: The BCS gap for different choices of chemical potential (first row) and Hubble parameter (second row), with the other parameters fixed to be M = 40H * , Λ = 100H * , κ = 50H * .In the condensed phase, the BCS gap increases with chemical potential and decreases with the Hubble parameter.The discontinuous change of the gap signals a first-order phase transition, where the two vacua possess equal free energy but are separated by a potential barrier.Beyond the critical chemical potential/temperature (solid circles), the asymmetric vacuum is metastable and may tunnel back to the symmetric vacuum via bubble nucleation that dismantles the condensate.The empty circles correspond to the inflection point where the false vacuum disappears.

Figure 3 .
Figure 3.The phase diagram of our inflationary BCS model.The solid blue line separates the normal phase and the condensed phase, all the way from the origin to the cutoff scale Λ = 100H * .Here the coupling is chosen to be M = 40H * .In the exact flat spacetime case with H = 0, the system always stays in the condensed phase, while an exponentially small Hubble expansion rate can destabilize the condensate.

Figure 4 .
Figure 4.The cosmological collider tomography, i.e. different momentum ratios probe physics at different energy scales.

Figure 5 .
Figure 5.The leading diagrams that contribute to the bispectrum.

Figure 7 .
Figure 7.The gravitational waves energy density Ω GW from BCS phase transition as a function of frequency f .The legend insertion of the figure explicitly displays the various parameter selections, where N e is the e-folding number before the end of inflation.The experiment curves are from IPTA[91], SKA[50], as well as space-based detectors such as LISA[51], DECIGO[52], TianQin[92], Taĳi[93].The gray line is dashed to indicate a risky parameter choice ε = 0.5.Unfortunately, for ε ≲ 0.1 in general, we do not expect the signal strength (brown solid line) to reach the sensitivity of LISA or DECIGO.