Excitation of the Dissipationless Higgs Mode in a Fermionic Condensate
R. A. Barankov1 and L. S. Levitov2
Department of Physics, University of Illinois, 1110 West Green Street, Urbana, IL 61801
Department of Physics, Massachusetts Institute of Technology, 77 Massachusetts Ave, Cambridge, MA 02139
The amplitude mode of a fermionic superfluid, analogous to the Higgs Boson, becomes undamped
in the strong coupling regime when its frequency is pushed inside the BCS energy gap. We argue that
this is the case in cold gases due to the energy dispersion and nonlocality of the pairing interaction,
and propose to use the Feshbach resonance regime for parametric excitation of this mode. The results
presented for the BCS pairing dynamics indicate that even weak dispersion suppresses dephasing
and gives rise to persistent oscillations. The frequency of oscillations extracted from our simulation
of the BCS dynamics agrees with the prediction of the many-body theory.
The observation of resonance superfluidity in cold
atomic Fermi gases [1, 2] at magnetically tunable Fes-
hbach resonances [3] opened new avenue of exploring the
many-body phenomena. Similar to the earlier work on
cold Bose gases which triggered studies of fascinating col-
lective phenomena [5, 6], fermionic pairing at Feshbach
resonances [7, 8, 9, 10, 11, 12] presents new opportunities.
In particular, the high degree of coherence of trapped
atoms, and the possibility to control particle interaction
in situ on the shortest collective time-scale, the inverse
Fermi energy [4], can facilitate exploring new regimes
which are difficult to realize in solid state systems.
The theory of fermionic pairing predicts two princi-
pal collective modes intrinsic to the condensed state.
One is the massless Bogoliubov-Anderson mode related
to the order parameter phase dynamics. Being a Gold-
stone mode, it manifests itself in hydrodynamics in the
same way as in Bose systems, and was probed recently
in the experiments on gas expansion and oscillation in
traps [13, 14, 15, 16]. In addition, there exists a second
fundamental elementary excitation [17, 18, 19], related to
the dynamics of the order parameter modulus |∆|. No-
tably, this excitation is unique to fermionic pairing and
has no counterpart in Bose systems [19]. This massive
excitation, characterized by a finite frequency, is anal-
ogous to the Higgs Boson in particle physics. Like the
latter it remained elusive, for a long time evading direct
probes, although some indirect manifestations have been
discussed [20, 21]. The main obstacle to the detection
of the Higgs mode in superconductors is that it is es-
sentially decoupled from the phase mode responsible for
hydrodynamics and superfluidity.
In this work we propose to use the dynamical con-
trol of pairing interaction demonstrated in Refs.[1, 2] for
parametric excitation of the Higgs mode. We argue that
fermion superfluidity in the strong coupling regime real-
ized near Feshbach resonance represents a distinct ad-
vantage, since in this case the Higgs mode is pushed
inside the superconducting gap, h̄ω < 2∆, which elim-
inates damping due to coupling to quasiparticles. We
demonstrate that this mode can be excited by a time-
dependent pairing interaction, as illustrated in Fig. 1. In
FIG. 1: Non-decaying Higgs mode in a Fermi gas with
energy-dependent pairing interaction excited by the interac-
tion switching from gi at t < 0 to g at t > 0. Shown are
the time and energy dependence of the pairing amplitude (a),
the x-component of the pseudospin vector (b), and the pairing
amplitude at the Fermi energy (c) as obtained from the model
(2),(3) at g = 0.43, gi = 0.23, a1 = a2 = 0.5, γ/W = 0.01,
∆F /W = 0.016. Note the initial transient of few periods, ex-
hibiting some dephasing in pseudospin dynamics (b), followed
by synchronized collective oscillations of fermion states.
contrast, the BCS theory at weak coupling predicts the
Higgs mode frequency right at the edge of the quasiparti-
cle continuum, h̄ω = 2∆ [19], which leads to collisionless
damping of this mode [22, 28].
The departure from the behavior at weak coupling
arises from the change in the character of pairing inter-
action in the strong coupling regime, in particular due
http://arxiv.org/abs/0704.1292v1
to its finite spatial radius and frequency dispersion. Spa-
tial nonlocality of pairing interaction is known to lead
to discrete collective modes inside the BCS gap [17, 18].
Similarly, the energy dispersion of the pairing interac-
tion and pairing amplitude ∆p that becomes prominent
at strong coupling [23], leads to discrete collective mode
spectrum (see below). While the exact form of this dis-
persion is sensitive to the specifics of the strong cou-
pling problem, it is established in the literature that,
generally, both effects can occur near Feshbach reso-
nance [8, 10, 11, 12, 24]. Although our understanding
of the detailed microscopic picture may be hampered by
the nonpertubative nature of the strong coupling prob-
lem, we shall see that within a simplified model used be-
low the inequality h̄ω < 2∆ is fulfilled under very general
conditions.
The dissipationless BCS dynamics [22, 25] and the pos-
sibility to realize it in cold gases [26] attracted much at-
tention recently [24, 27, 28]. These investigations, with
the exception of Ref. [24], focused on the case of pair-
ing interaction which is constant in the entire fermion
energy band, concluding [29, 30] that several interesting
dynamical states, synchronized and desynchronized (or
dephased), can be realized by a sudden change in the
interaction strength (see the phase diagram in Ref. [29]).
In contrast, as we shall see below, the dephased be-
havior is suppressed in the strong coupling regime when
due to the energy dispersion of the pairing interaction
the Higgs mode falls inside the BCS energy gap. Un-
der these conditions an undamped Higgs mode can be
excited upon a sudden change in interaction. By analyz-
ing the limit when the interaction dispersion disappears
we show how the different regimes of Ref. [29] are recov-
ered. This correspondence suggests an interpretation of
the dephased oscillations discussed in Refs. [22, 29, 30]
as a manifestation of the Higgs mode, algebraically de-
phased at h̄ω = 2∆.
We shall analyze the pairing dynamics in a spatially
uniform system using the pseudospin representation [25]
of the BCS problem in which spin 1/2 operators s±
± isy
describe Cooper pairs (p,−p):
H = −
λpq(t)s
, (1)
where ǫp is the free particle spectrum. The interaction
λpq(t) that models the energy dispersion at strong cou-
pling is taken in the form of a sum of a dispersing and
nondispersing parts
λpq(t) =
(a1 + a2fpfq) , fp =
γ2 + ǫ2
, (2)
where the dimensionless parameter g(t) specifies the in-
teraction time-dependence, the constants a1,2 ≥ 0 satisfy
a1 + a2 = 1, and νF is the density of states at the Fermi
level. The second term in (2) features dispersion on the
energy scale γ. Our motivation for choosing the model
(2) was two-fold. Firstly, the form (2) is general enough
to provide insight into the role of different features, such
as the energy dispersion (which is controlled by the pa-
rameter γ) and separability (which is absent unless a1
or a2 vanishes). Secondly, our numerical method utilized
the rank two form of (2), allowing for substantial speedup
that could not be implemented for a more general inter-
action λpq. In addition, the model (2) is physically mo-
tivated by the theory of BCS pairing in the simultaneous
presence of a retarded and non-retarded interaction [31].
Within the mean-field approximation, the dynamical
equations derived from Eq.(1) assume a Bloch form:
= 2bp × rp, bp = −(∆
, ǫp), (3)
where rp = 2〈sp〉 are Bloch vectors, and the effective
magnetic field bp depends on the pairing amplitude ∆p.
The latter is defined self-consistently:
∆p = ∆
+ i∆y
λpq(t)
+ iry
. (4)
The interaction time dependence of interest is a step-like
change from the initial value gi to the final value g. With-
out loss of generality, the phase of the order parameter
can be chosen equal zero, allowing us to consider only the
x-component of the pairing amplitude, ∆p = ∆
. As an
initial state we take the paired ground state
(0) =
)2 + ǫ2
(0) =
)2 + ǫ2
. (5)
The equilibrium energy-dependent amplitude ∆p is de-
termined by the self-consistency equation
, (6)
in which λpq is given by (2) with the parameter values gi
and g for the initial and final state. The corresponding
equilibrium pairing gap values, ∆i
and ∆p, are found by
numerically solving the integral equation (6). Through-
out the paper we use the equilibrium value of the pairing
gap at the Fermi level, ∆F , at the final coupling g as a
natural energy scale to parameterize the dynamics.
We integrate Eqs.(3) using the Runge-Kutta method
of the 4-th order with a time step adjusted to achieve
sufficient precision of the calculation. In our simulation
we use N = 104, 105 equally spaced energy states within
bandwidthW , −W/2 < ǫp < W/2, with the level spacing
much smaller than all other energy scales in the problem.
We analyze the quantity ∆F (t) which at long times
oscillates between the maximum and minimum values
−2 −1 0 1
Inverse coupling 1/g−1/g
0 0.25 0.5
)−1/2
FIG. 2: Long-time behavior of ∆F (t), the pairing amplitude
at the Fermi level, oscillating between ∆+ and ∆−. Shown are
two examples of ∆± as a function of the initial state for a non-
separable (circles) and a separable (squares) interaction (2).
Parameters used: a1 = a2 = 0.5, g = 0.43, ∆F/W = 0.016,
γ/W = 0.01, and a1 = 0, a2 = 1, g = 0.61, ∆F /W = 0.005,
γ/W = 0.01, respectively. Inset: Linear fit of a sample trace
∆F (t) vs. t
−1/2 used to extract ∆±.
∆+ and ∆−. To find the asymptotic values ∆± we em-
ploy the numerical procedure sketched in Fig. 2 inset:
∆± are obtained from the linear fits to the maxima and
minima of ∆F vs. t
−1/2 intersection with the y-axis.
The t−1/2 time parameterization is motivated by the de-
phasing law δ∆(t) ∝ t−1/2 found in Refs.[22, 28] for the
energy-independent interaction. Should the dephasing
occur, the asymptotic values would coincide, ∆+ = ∆−.
In contrast to the above, for the interaction (2) the
dephased behavior is suppressed. Instead, as illustrated
in Fig. 2, we observe non-decaying periodic oscillations
for a wide range of initial states, both for the initial states
close to the normal state (gi ≪ g) as well as for the
initial states near equilibrium (gi ≈ g). At increasing gi
there is a critical point at which the asymptotic pairing
amplitude ∆± becomes zero.
To understand the origin of the oscillatory behavior for
the dispersive interaction, Eq.(2), we develop perturba-
tion theory near the point gi = g. Linearizing the Bloch
equations and taking a harmonic variation of the pairing
amplitude, δ∆x,y
(t) ∝ e−iωtδ∆x,y
pω , we find two collective
modes for the x and y components of ∆p corresponding
to the order parameter amplitude and phase variation
(see Ref.[19]). The amplitude (Higgs) mode with fre-
quency ω obeys the integral equation
λpqδ∆
− ω2/4
, (7)
where ∆p is the equilibrium gap obtained from Eq.(6).
The equation for δ∆y
(the phase mode) is similar to
Eq.(7) except for the denominator of the second frac-
−2 −1.5 −1 −0.5 0 0.5 1
Inverse coupling 1/g−1/g
0.02 0.04 0.06 0.08 0.1
Interaction dispersion γ/W
FIG. 3: (a): The Higgs mode frequency obtained from the
simulation with gi ≈ g (circles) and from Eq.(7) (solid line)
as a function of the dispersion parameter γ. The quasiparticle
energy minimum (dashed line) lies above the collective mode
frequency (parameters of the simulation: g = 0.61, a1 = a2 =
0.5). (b): Frequency of the Higgs mode as a function of the
initial state for non-separable and separable interactions with
the same parameters as in Fig.1. The frequency changes away
from gi = g as the amplitude of oscillations increases (Fig.1),
indicating unharmonicity of the Higgs mode.
tion which is ǫ2
− ω2/4. As expected from Goldstone
theorem, the equation for δ∆y
is solved by ω = 0.
To find the frequency ω of the Higgs mode, we note
that for the interaction λpq given by (2), which is an op-
erator of rank two, Eq.(7) turns into an algebraic equa-
tion involving a 2×2 determinant. Solving it we find that
for a2 > 0 the frequency ω lies within the BCS gap, as il-
lustrated in Fig.3a. To gain more insight, let us consider
a separable interaction, a1 = 0, a2 = 1, which yields
− ω2/4
, (8)
where ∆q ∝ fq. Balancing the factors under the sum in
order to obtain unity on the left hand side, and noting
that without the second factor Eq.(8) would be identical
to Eq.(6), it is easy to see that ω < 2∆F , i.e. the Higgs
mode is discrete.
Notably, as Fig.3a illustrates, the frequency obtained
from Eq.(7) coincides with the frequency of oscillations in
∆F (t) obtained by simulating BCS dynamics at g ≈ gi,
proving that the observed excitation is indeed the Higgs
mode. Furthermore, for g away from gi the frequency
extracted from ∆F (t) varies with g, decreasing below the
value at g ≈ gi and approaching zero at gi ≪ g and gi ≫
g (see Fig.3b). This indicates unharmonicity of the Higgs
mode that sets on at a large amplitude of oscillations.
To test these ideas further, we considered the regime
when the Higgs mode is strictly inside the quasiparticle
continuum, which can be realized in the model (2) with
−3 −2 −1 0 1 2 3
Inverse coupling 1/g−1/g
Onset of dephasing
FIG. 4: Quenching of dephasing for weakly dispersing in-
teraction. Asymptotic values of the pairing amplitude ∆±
for a dispersing (circles) and non-dispersing (red line) inter-
action for different initial states. The onset of dephasing is
marked by arrows. Parameters used: a1 = a2 = 0.5, g = 0.33,
∆F /W = 0.02, γ/W = 0.1, and a1 = 1, a2 = 0, g = 0.33,
∆F /W = 0.05, respectively.
the second term of a repulsive sign, a2 < 0. In this case
Eq.(7) has no real-valued solution in the region ω ≤ 2∆.
Simulating the BCS dynamics near gi ≈ g we find that
∆(t) exhibits exponentially decaying oscillations of the
form e−ηt cos(ω′t+φ) corresponding to a complex-valued
frequency ω. For a2 = 0 the collective mode frequency
ω = 2∆F lies at the edge of the quasiparticle continuum.
This property was linked to algebraic Landau damping
of this mode in Refs.[22, 28].
The discrete Higgs mode makes the BCS dynamics un-
damped for g near gi even for weakly dispersing interac-
tion λpq. It is interesting to connect this behavior to the
dephased BCS dynamics found in the case of constant in-
teraction. This is illustrated (Fig.4) by the dynamics at
weakly dispersing interaction γ ≫ ∆F , where we observe
that the region of dephased dynamics shrinks, with the
onset of dephasing shifting towards small g < gi. While
the oscillation amplitude 1
(∆+−∆−) is now finite, it re-
mains small due to dephasing in the transient region (see
Fig.1b). This behavior is consistent with the Higgs mode
approaching the quasiparticle continuum boundary.
In conclusion, we have shown that the energy disper-
sion of pairing interaction leads to quenching of dephas-
ing of the BCS dynamics, making the Higgs mode of the
pairing amplitude discrete. Parametric control of inter-
action in the strong coupling regime near a Feshbach res-
onance of cold atoms can be used to excite this mode.
This research was supported in part by the National
Science Foundation under Grant No. PHY05-51164.
[1] C. A. Regal, M. Greiner, and D. S. Jin, Phys. Rev. Lett.
92, 040403 (2004).
[2] M. W. Zwierlein, C. A. Stan, C. H. Schunck, S. M. F.
Raupach, A. J. Kerman, and W. Ketterle, Phys. Rev.
Lett. 92, 120403 (2004).
[3] E. Timmermans, P. Tommasini, M. Hussein, and A. Ker-
man, Physics Reports 315, 199 (1999).
[4] M. W. Zwierlein, C. H. Schunck, C. A. Stan, S. M. F.
Raupach, and W. Ketterle, Phys. Rev. Lett. 94, 180401
(2005).
[5] F. Dalfovo, S. Giorgini, L. P. Pitaesvkii, and S. Stringari,
Rev. Mod. Phys. 71, 463 (1999).
[6] C. J. Pethick and H. Smith, Bose-Einstein condensation
in dilute gases, Cambridge University Press, 2002.
[7] M. Holland, S. J. J. M. F. Kokkelmans, M. L. Chiofalo,
and R. Walser, Phys. Rev. Lett. 87, 120406 (2001).
[8] Y. Ohashi and A. Griffin, Phys. Rev. A 67, 063612
(2003).
[9] Tin-Lun Ho, Phys. Rev. Lett. 92, 090402 (2004).
[10] G. M. Bruun and C. J. Pethick, Phys. Rev. Lett. 92,
140404 (2004).
[11] Q. Chen, J. Stajic, S. Tan, and K. Levin, Physics Reports
412, 1 (2005).
[12] A. Schwenk and C. J. Pethick, Phys. Rev. Lett. 95,
160401 (2005).
[13] T. Bourdel, J. Cubizolles, L. Khaykovich, K. M. Maga-
lhes, S. J. Kokkelmans, G. V. Shlyapnikov, and C. Sa-
lomon, Phys. Rev. Lett. 91, 020402 (2003).
[14] F. Ferlaino, E. de Mirandes, G. Roati, G. Modugno, and
M. Inguscio, Phys. Rev. Lett. 92, 140405 (2004).
[15] J. Kinast, A. Turlapov, and J. E. Thomas, Phys. Rev.
Lett. 94, 170404 (2005).
[16] J. T. Stewart, J. P. Gaebler, C. A. Regal, and D. S. Jin,
Phys. Rev. Lett. 97, 220406 (2006).
[17] V. G. Vaks, V. M. Galitskii, A. I. Larkin, Zh. Eksp. Teor.
Fiz. 41, 1655 (1961) [Sov. Phys. JETP 14, 1177 (1962)].
[18] A. Bardasis and J. R. Schrieffer, Phys. Rev. 121, 1050
(1961).
[19] C. M. Varma, J. Low Temp. Phys. 126, 901 (2002).
[20] B. I. Ivlev, Pis’ma v Zh. Eksp. Teor. Fiz. 15, 441 (1972)
[JETP Lett. 15, 313 (1972)].
[21] P. B. Littlewood and C. M. Varma, Phys. Rev. B 26,
4883 (1982).
[22] A. F. Volkov and Sh. M. Kogan, Zh. Eksp. Teor. Fiz. 65,
2038, (1973) [Sov. Phys. JETP 38, 1018 (1974)].
[23] G. M. Eliashberg, Zh. Eksp. Teor. Fiz. 38, 966 (1960)
[Sov. Phys. JETP 11, 696 (1960)].
[24] M. H. Szymańska, B. D. Simons, and K. Burnett, Phys.
Rev. Lett. 94, 170402 (2005).
[25] R. A. Barankov, L. S. Levitov, and B. Z. Spivak, Phys.
Rev. Lett. 93, 160401 (2004).
[26] A. V. Andreev, V. Gurarie, and L. Radzihovsky, Phys.
Rev. Lett. 93, 130402 (2004); R. A. Barankov and L. S.
Levitov, Phys. Rev. Lett. 93, 130403 (2004).
[27] G. L. Warner and A. J. Leggett, Phys. Rev. B 71, 134514
(2005).
[28] E. A. Yuzbashyan, O. Tsyplyatyev, and B. L. Altshuler,
Phys. Rev. Lett. 96, 097005 (2006); Erratum: Phys. Rev.
Lett. 96, 179905 (2006).
[29] R. A. Barankov and L. S. Levitov, Phys. Rev. Lett. 96,
230403 (2006).
[30] E. A. Yuzbashyan and M. Dzero, Phys. Rev. Lett. 96,
230404 (2006).
[31] P. Morel and P. W. Anderson, Phys. Rev. 125, 1263
(1962).

  The amplitude mode of a fermionic superfluid, analogous to the Higgs Boson,
becomes undamped in the strong coupling regime when its frequency is pushed
inside the BCS energy gap. We argue that this is the case in cold gases due to
the energy dispersion and nonlocality of the pairing interaction, and propose
to use the Feshbach resonance regime for parametric excitation of this mode.
The results presented for the BCS pairing dynamics indicate that even weak
dispersion suppresses dephasing and gives rise to persistent oscillations. The
frequency of oscillations extracted from our simulation of the BCS dynamics
agrees with the prediction of the many-body theory.

Excitation of the Dissipationless Higgs Mode in a Fermionic Condensate
R. A. Barankov1 and L. S. Levitov2
Department of Physics, University of Illinois, 1110 West Green Street, Urbana, IL 61801
Department of Physics, Massachusetts Institute of Technology, 77 Massachusetts Ave, Cambridge, MA 02139
The amplitude mode of a fermionic superfluid, analogous to the Higgs Boson, becomes undamped
in the strong coupling regime when its frequency is pushed inside the BCS energy gap. We argue that
this is the case in cold gases due to the energy dispersion and nonlocality of the pairing interaction,
and propose to use the Feshbach resonance regime for parametric excitation of this mode. The results
presented for the BCS pairing dynamics indicate that even weak dispersion suppresses dephasing
and gives rise to persistent oscillations. The frequency of oscillations extracted from our simulation
of the BCS dynamics agrees with the prediction of the many-body theory.
The observation of resonance superfluidity in cold
atomic Fermi gases [1, 2] at magnetically tunable Fes-
hbach resonances [3] opened new avenue of exploring the
many-body phenomena. Similar to the earlier work on
cold Bose gases which triggered studies of fascinating col-
lective phenomena [5, 6], fermionic pairing at Feshbach
resonances [7, 8, 9, 10, 11, 12] presents new opportunities.
In particular, the high degree of coherence of trapped
atoms, and the possibility to control particle interaction
in situ on the shortest collective time-scale, the inverse
Fermi energy [4], can facilitate exploring new regimes
which are difficult to realize in solid state systems.
The theory of fermionic pairing predicts two princi-
pal collective modes intrinsic to the condensed state.
One is the massless Bogoliubov-Anderson mode related
to the order parameter phase dynamics. Being a Gold-
stone mode, it manifests itself in hydrodynamics in the
same way as in Bose systems, and was probed recently
in the experiments on gas expansion and oscillation in
traps [13, 14, 15, 16]. In addition, there exists a second
fundamental elementary excitation [17, 18, 19], related to
the dynamics of the order parameter modulus |∆|. No-
tably, this excitation is unique to fermionic pairing and
has no counterpart in Bose systems [19]. This massive
excitation, characterized by a finite frequency, is anal-
ogous to the Higgs Boson in particle physics. Like the
latter it remained elusive, for a long time evading direct
probes, although some indirect manifestations have been
discussed [20, 21]. The main obstacle to the detection
of the Higgs mode in superconductors is that it is es-
sentially decoupled from the phase mode responsible for
hydrodynamics and superfluidity.
In this work we propose to use the dynamical con-
trol of pairing interaction demonstrated in Refs.[1, 2] for
parametric excitation of the Higgs mode. We argue that
fermion superfluidity in the strong coupling regime real-
ized near Feshbach resonance represents a distinct ad-
vantage, since in this case the Higgs mode is pushed
inside the superconducting gap, h̄ω < 2∆, which elim-
inates damping due to coupling to quasiparticles. We
demonstrate that this mode can be excited by a time-
dependent pairing interaction, as illustrated in Fig. 1. In
FIG. 1: Non-decaying Higgs mode in a Fermi gas with
energy-dependent pairing interaction excited by the interac-
tion switching from gi at t < 0 to g at t > 0. Shown are
the time and energy dependence of the pairing amplitude (a),
the x-component of the pseudospin vector (b), and the pairing
amplitude at the Fermi energy (c) as obtained from the model
(2),(3) at g = 0.43, gi = 0.23, a1 = a2 = 0.5, γ/W = 0.01,
∆F /W = 0.016. Note the initial transient of few periods, ex-
hibiting some dephasing in pseudospin dynamics (b), followed
by synchronized collective oscillations of fermion states.
contrast, the BCS theory at weak coupling predicts the
Higgs mode frequency right at the edge of the quasiparti-
cle continuum, h̄ω = 2∆ [19], which leads to collisionless
damping of this mode [22, 28].
The departure from the behavior at weak coupling
arises from the change in the character of pairing inter-
action in the strong coupling regime, in particular due
http://arxiv.org/abs/0704.1292v1
to its finite spatial radius and frequency dispersion. Spa-
tial nonlocality of pairing interaction is known to lead
to discrete collective modes inside the BCS gap [17, 18].
Similarly, the energy dispersion of the pairing interac-
tion and pairing amplitude ∆p that becomes prominent
at strong coupling [23], leads to discrete collective mode
spectrum (see below). While the exact form of this dis-
persion is sensitive to the specifics of the strong cou-
pling problem, it is established in the literature that,
generally, both effects can occur near Feshbach reso-
nance [8, 10, 11, 12, 24]. Although our understanding
of the detailed microscopic picture may be hampered by
the nonpertubative nature of the strong coupling prob-
lem, we shall see that within a simplified model used be-
low the inequality h̄ω < 2∆ is fulfilled under very general
conditions.
The dissipationless BCS dynamics [22, 25] and the pos-
sibility to realize it in cold gases [26] attracted much at-
tention recently [24, 27, 28]. These investigations, with
the exception of Ref. [24], focused on the case of pair-
ing interaction which is constant in the entire fermion
energy band, concluding [29, 30] that several interesting
dynamical states, synchronized and desynchronized (or
dephased), can be realized by a sudden change in the
interaction strength (see the phase diagram in Ref. [29]).
In contrast, as we shall see below, the dephased be-
havior is suppressed in the strong coupling regime when
due to the energy dispersion of the pairing interaction
the Higgs mode falls inside the BCS energy gap. Un-
der these conditions an undamped Higgs mode can be
excited upon a sudden change in interaction. By analyz-
ing the limit when the interaction dispersion disappears
we show how the different regimes of Ref. [29] are recov-
ered. This correspondence suggests an interpretation of
the dephased oscillations discussed in Refs. [22, 29, 30]
as a manifestation of the Higgs mode, algebraically de-
phased at h̄ω = 2∆.
We shall analyze the pairing dynamics in a spatially
uniform system using the pseudospin representation [25]
of the BCS problem in which spin 1/2 operators s±
± isy
describe Cooper pairs (p,−p):
H = −
λpq(t)s
, (1)
where ǫp is the free particle spectrum. The interaction
λpq(t) that models the energy dispersion at strong cou-
pling is taken in the form of a sum of a dispersing and
nondispersing parts
λpq(t) =
(a1 + a2fpfq) , fp =
γ2 + ǫ2
, (2)
where the dimensionless parameter g(t) specifies the in-
teraction time-dependence, the constants a1,2 ≥ 0 satisfy
a1 + a2 = 1, and νF is the density of states at the Fermi
level. The second term in (2) features dispersion on the
energy scale γ. Our motivation for choosing the model
(2) was two-fold. Firstly, the form (2) is general enough
to provide insight into the role of different features, such
as the energy dispersion (which is controlled by the pa-
rameter γ) and separability (which is absent unless a1
or a2 vanishes). Secondly, our numerical method utilized
the rank two form of (2), allowing for substantial speedup
that could not be implemented for a more general inter-
action λpq. In addition, the model (2) is physically mo-
tivated by the theory of BCS pairing in the simultaneous
presence of a retarded and non-retarded interaction [31].
Within the mean-field approximation, the dynamical
equations derived from Eq.(1) assume a Bloch form:
= 2bp × rp, bp = −(∆
, ǫp), (3)
where rp = 2〈sp〉 are Bloch vectors, and the effective
magnetic field bp depends on the pairing amplitude ∆p.
The latter is defined self-consistently:
∆p = ∆
+ i∆y
λpq(t)
+ iry
. (4)
The interaction time dependence of interest is a step-like
change from the initial value gi to the final value g. With-
out loss of generality, the phase of the order parameter
can be chosen equal zero, allowing us to consider only the
x-component of the pairing amplitude, ∆p = ∆
. As an
initial state we take the paired ground state
(0) =
)2 + ǫ2
(0) =
)2 + ǫ2
. (5)
The equilibrium energy-dependent amplitude ∆p is de-
termined by the self-consistency equation
, (6)
in which λpq is given by (2) with the parameter values gi
and g for the initial and final state. The corresponding
equilibrium pairing gap values, ∆i
and ∆p, are found by
numerically solving the integral equation (6). Through-
out the paper we use the equilibrium value of the pairing
gap at the Fermi level, ∆F , at the final coupling g as a
natural energy scale to parameterize the dynamics.
We integrate Eqs.(3) using the Runge-Kutta method
of the 4-th order with a time step adjusted to achieve
sufficient precision of the calculation. In our simulation
we use N = 104, 105 equally spaced energy states within
bandwidthW , −W/2 < ǫp < W/2, with the level spacing
much smaller than all other energy scales in the problem.
We analyze the quantity ∆F (t) which at long times
oscillates between the maximum and minimum values
−2 −1 0 1
Inverse coupling 1/g−1/g
0 0.25 0.5
)−1/2
FIG. 2: Long-time behavior of ∆F (t), the pairing amplitude
at the Fermi level, oscillating between ∆+ and ∆−. Shown are
two examples of ∆± as a function of the initial state for a non-
separable (circles) and a separable (squares) interaction (2).
Parameters used: a1 = a2 = 0.5, g = 0.43, ∆F/W = 0.016,
γ/W = 0.01, and a1 = 0, a2 = 1, g = 0.61, ∆F /W = 0.005,
γ/W = 0.01, respectively. Inset: Linear fit of a sample trace
∆F (t) vs. t
−1/2 used to extract ∆±.
∆+ and ∆−. To find the asymptotic values ∆± we em-
ploy the numerical procedure sketched in Fig. 2 inset:
∆± are obtained from the linear fits to the maxima and
minima of ∆F vs. t
−1/2 intersection with the y-axis.
The t−1/2 time parameterization is motivated by the de-
phasing law δ∆(t) ∝ t−1/2 found in Refs.[22, 28] for the
energy-independent interaction. Should the dephasing
occur, the asymptotic values would coincide, ∆+ = ∆−.
In contrast to the above, for the interaction (2) the
dephased behavior is suppressed. Instead, as illustrated
in Fig. 2, we observe non-decaying periodic oscillations
for a wide range of initial states, both for the initial states
close to the normal state (gi ≪ g) as well as for the
initial states near equilibrium (gi ≈ g). At increasing gi
there is a critical point at which the asymptotic pairing
amplitude ∆± becomes zero.
To understand the origin of the oscillatory behavior for
the dispersive interaction, Eq.(2), we develop perturba-
tion theory near the point gi = g. Linearizing the Bloch
equations and taking a harmonic variation of the pairing
amplitude, δ∆x,y
(t) ∝ e−iωtδ∆x,y
pω , we find two collective
modes for the x and y components of ∆p corresponding
to the order parameter amplitude and phase variation
(see Ref.[19]). The amplitude (Higgs) mode with fre-
quency ω obeys the integral equation
λpqδ∆
− ω2/4
, (7)
where ∆p is the equilibrium gap obtained from Eq.(6).
The equation for δ∆y
(the phase mode) is similar to
Eq.(7) except for the denominator of the second frac-
−2 −1.5 −1 −0.5 0 0.5 1
Inverse coupling 1/g−1/g
0.02 0.04 0.06 0.08 0.1
Interaction dispersion γ/W
FIG. 3: (a): The Higgs mode frequency obtained from the
simulation with gi ≈ g (circles) and from Eq.(7) (solid line)
as a function of the dispersion parameter γ. The quasiparticle
energy minimum (dashed line) lies above the collective mode
frequency (parameters of the simulation: g = 0.61, a1 = a2 =
0.5). (b): Frequency of the Higgs mode as a function of the
initial state for non-separable and separable interactions with
the same parameters as in Fig.1. The frequency changes away
from gi = g as the amplitude of oscillations increases (Fig.1),
indicating unharmonicity of the Higgs mode.
tion which is ǫ2
− ω2/4. As expected from Goldstone
theorem, the equation for δ∆y
is solved by ω = 0.
To find the frequency ω of the Higgs mode, we note
that for the interaction λpq given by (2), which is an op-
erator of rank two, Eq.(7) turns into an algebraic equa-
tion involving a 2×2 determinant. Solving it we find that
for a2 > 0 the frequency ω lies within the BCS gap, as il-
lustrated in Fig.3a. To gain more insight, let us consider
a separable interaction, a1 = 0, a2 = 1, which yields
− ω2/4
, (8)
where ∆q ∝ fq. Balancing the factors under the sum in
order to obtain unity on the left hand side, and noting
that without the second factor Eq.(8) would be identical
to Eq.(6), it is easy to see that ω < 2∆F , i.e. the Higgs
mode is discrete.
Notably, as Fig.3a illustrates, the frequency obtained
from Eq.(7) coincides with the frequency of oscillations in
∆F (t) obtained by simulating BCS dynamics at g ≈ gi,
proving that the observed excitation is indeed the Higgs
mode. Furthermore, for g away from gi the frequency
extracted from ∆F (t) varies with g, decreasing below the
value at g ≈ gi and approaching zero at gi ≪ g and gi ≫
g (see Fig.3b). This indicates unharmonicity of the Higgs
mode that sets on at a large amplitude of oscillations.
To test these ideas further, we considered the regime
when the Higgs mode is strictly inside the quasiparticle
continuum, which can be realized in the model (2) with
−3 −2 −1 0 1 2 3
Inverse coupling 1/g−1/g
Onset of dephasing
FIG. 4: Quenching of dephasing for weakly dispersing in-
teraction. Asymptotic values of the pairing amplitude ∆±
for a dispersing (circles) and non-dispersing (red line) inter-
action for different initial states. The onset of dephasing is
marked by arrows. Parameters used: a1 = a2 = 0.5, g = 0.33,
∆F /W = 0.02, γ/W = 0.1, and a1 = 1, a2 = 0, g = 0.33,
∆F /W = 0.05, respectively.
the second term of a repulsive sign, a2 < 0. In this case
Eq.(7) has no real-valued solution in the region ω ≤ 2∆.
Simulating the BCS dynamics near gi ≈ g we find that
∆(t) exhibits exponentially decaying oscillations of the
form e−ηt cos(ω′t+φ) corresponding to a complex-valued
frequency ω. For a2 = 0 the collective mode frequency
ω = 2∆F lies at the edge of the quasiparticle continuum.
This property was linked to algebraic Landau damping
of this mode in Refs.[22, 28].
The discrete Higgs mode makes the BCS dynamics un-
damped for g near gi even for weakly dispersing interac-
tion λpq. It is interesting to connect this behavior to the
dephased BCS dynamics found in the case of constant in-
teraction. This is illustrated (Fig.4) by the dynamics at
weakly dispersing interaction γ ≫ ∆F , where we observe
that the region of dephased dynamics shrinks, with the
onset of dephasing shifting towards small g < gi. While
the oscillation amplitude 1
(∆+−∆−) is now finite, it re-
mains small due to dephasing in the transient region (see
Fig.1b). This behavior is consistent with the Higgs mode
approaching the quasiparticle continuum boundary.
In conclusion, we have shown that the energy disper-
sion of pairing interaction leads to quenching of dephas-
ing of the BCS dynamics, making the Higgs mode of the
pairing amplitude discrete. Parametric control of inter-
action in the strong coupling regime near a Feshbach res-
onance of cold atoms can be used to excite this mode.
This research was supported in part by the National
Science Foundation under Grant No. PHY05-51164.
[1] C. A. Regal, M. Greiner, and D. S. Jin, Phys. Rev. Lett.
92, 040403 (2004).
[2] M. W. Zwierlein, C. A. Stan, C. H. Schunck, S. M. F.
Raupach, A. J. Kerman, and W. Ketterle, Phys. Rev.
Lett. 92, 120403 (2004).
[3] E. Timmermans, P. Tommasini, M. Hussein, and A. Ker-
man, Physics Reports 315, 199 (1999).
[4] M. W. Zwierlein, C. H. Schunck, C. A. Stan, S. M. F.
Raupach, and W. Ketterle, Phys. Rev. Lett. 94, 180401
(2005).
[5] F. Dalfovo, S. Giorgini, L. P. Pitaesvkii, and S. Stringari,
Rev. Mod. Phys. 71, 463 (1999).
[6] C. J. Pethick and H. Smith, Bose-Einstein condensation
in dilute gases, Cambridge University Press, 2002.
[7] M. Holland, S. J. J. M. F. Kokkelmans, M. L. Chiofalo,
and R. Walser, Phys. Rev. Lett. 87, 120406 (2001).
[8] Y. Ohashi and A. Griffin, Phys. Rev. A 67, 063612
(2003).
[9] Tin-Lun Ho, Phys. Rev. Lett. 92, 090402 (2004).
[10] G. M. Bruun and C. J. Pethick, Phys. Rev. Lett. 92,
140404 (2004).
[11] Q. Chen, J. Stajic, S. Tan, and K. Levin, Physics Reports
412, 1 (2005).
[12] A. Schwenk and C. J. Pethick, Phys. Rev. Lett. 95,
160401 (2005).
[13] T. Bourdel, J. Cubizolles, L. Khaykovich, K. M. Maga-
lhes, S. J. Kokkelmans, G. V. Shlyapnikov, and C. Sa-
lomon, Phys. Rev. Lett. 91, 020402 (2003).
[14] F. Ferlaino, E. de Mirandes, G. Roati, G. Modugno, and
M. Inguscio, Phys. Rev. Lett. 92, 140405 (2004).
[15] J. Kinast, A. Turlapov, and J. E. Thomas, Phys. Rev.
Lett. 94, 170404 (2005).
[16] J. T. Stewart, J. P. Gaebler, C. A. Regal, and D. S. Jin,
Phys. Rev. Lett. 97, 220406 (2006).
[17] V. G. Vaks, V. M. Galitskii, A. I. Larkin, Zh. Eksp. Teor.
Fiz. 41, 1655 (1961) [Sov. Phys. JETP 14, 1177 (1962)].
[18] A. Bardasis and J. R. Schrieffer, Phys. Rev. 121, 1050
(1961).
[19] C. M. Varma, J. Low Temp. Phys. 126, 901 (2002).
[20] B. I. Ivlev, Pis’ma v Zh. Eksp. Teor. Fiz. 15, 441 (1972)
[JETP Lett. 15, 313 (1972)].
[21] P. B. Littlewood and C. M. Varma, Phys. Rev. B 26,
4883 (1982).
[22] A. F. Volkov and Sh. M. Kogan, Zh. Eksp. Teor. Fiz. 65,
2038, (1973) [Sov. Phys. JETP 38, 1018 (1974)].
[23] G. M. Eliashberg, Zh. Eksp. Teor. Fiz. 38, 966 (1960)
[Sov. Phys. JETP 11, 696 (1960)].
[24] M. H. Szymańska, B. D. Simons, and K. Burnett, Phys.
Rev. Lett. 94, 170402 (2005).
[25] R. A. Barankov, L. S. Levitov, and B. Z. Spivak, Phys.
Rev. Lett. 93, 160401 (2004).
[26] A. V. Andreev, V. Gurarie, and L. Radzihovsky, Phys.
Rev. Lett. 93, 130402 (2004); R. A. Barankov and L. S.
Levitov, Phys. Rev. Lett. 93, 130403 (2004).
[27] G. L. Warner and A. J. Leggett, Phys. Rev. B 71, 134514
(2005).
[28] E. A. Yuzbashyan, O. Tsyplyatyev, and B. L. Altshuler,
Phys. Rev. Lett. 96, 097005 (2006); Erratum: Phys. Rev.
Lett. 96, 179905 (2006).
[29] R. A. Barankov and L. S. Levitov, Phys. Rev. Lett. 96,
230403 (2006).
[30] E. A. Yuzbashyan and M. Dzero, Phys. Rev. Lett. 96,
230404 (2006).
[31] P. Morel and P. W. Anderson, Phys. Rev. 125, 1263
(1962).
