# Free energy and specific heat near a quantum critical point of a metal

Shang-Shun Zhang  School of Physics and Astronomy and William I. Fine
Theoretical Physics Institute, University of Minnesota, Minneapolis, MN 55455,
USA  Erez Berg  Department of Condensed Matter Physics, Weizmann Institute of
Sciences, Rehovot, Israel  Andrey V. Chubukov  School of Physics and Astronomy
and William I. Fine Theoretical Physics Institute, University of Minnesota,
Minneapolis, MN 55455, USA

###### Abstract

We analyze free energy and specific heat for fermions interacting with gapless
bosons at a quantum-critical point (QCP) in a metal. We use the Luttinger-
Ward-Eliashberg formula for the free energy in the normal state, which
includes contributions from bosons, fermions, and their interaction, all
expressed via fully dressed fermionic and bosonic propagators. The sum of the
last two contributions is the free energy $F_{\gamma}$ of an effective low-
energy model of fermions with boson-mediated dynamical 4-fermion interaction
$V(\Omega_{m})\propto 1/|\Omega_{m}|^{\gamma}$ (the $\gamma-$model). This
purely electronic model has been used to analyze the interplay between non-
Fermi liquid (non-FL) behavior and pairing near a QCP, which are both
independent of the upper energy cutoff $\Lambda$. However, the specific heat
$C_{\gamma}(T)$, obtained from $F_{\gamma}$, does depend on $\Lambda$. We
argue that this dependence is spurious and cancels out, once we include the
contribution from bosons. We compare our $C(T)$ with the one obtained within
the $\gamma-$model using recently proposed regularization of $F_{\gamma}$. We
argue that for $\gamma<1$, the full $C(T)$ and the regularized $C_{\gamma}(T)$
differ by a $\gamma-$dependent prefactor, while for $\gamma>1$, the full
$C(T)$ is the sum of $C_{\gamma}(T)$ and the specific heat of free bosons with
fully dressed mass. For these $\gamma$, $C_{\gamma}(T)$ is negative. The
authors of Ref. [1] argued that a negative $C_{\gamma}(T)$ implies that the
normal state becomes unstable at some distance to a QCP. In our calculation,
both terms in $C(T)$ come from the same source, and $C_{\gamma}(T)$ is smaller
as long as vertex corrections can be safely neglected. We then argue that the
normal state remains stable even at a QCP.

##  I  Introduction.

In this work we analyze in detail the free energy and specific heat of a metal
near a critical point towards a spontaneous particle-hole order (Ising-
nematic, antiferromagnetic, etc), and of an electron-phonon system at
vanishing dressed Debye frequency of an optical phonon. In all these cases,
the low-energy physics is described by a model of fermions with Luttinger
Fermi surface, coupled by Yukawa-type interaction to a near-massless boson,
which represents either a critical fluctuation of a particle-hole order
parameter or a soft optical phonon [2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,
14, 15, 16, 17, 18, 19, *rech_2006_1, 21, 22, *Wang_H_17, *Wang_H_18,
*Fitzpatrick_15, 26, *sslee2, *lunts_2017]. The key motivation for our study
is current interest in a non-Fermi liquid (non-FL) behavior near a quantum-
critical point (QCP). Numerous previous studies have shown [15, 29, 16, 17,
18, 30, *metzner_new, 32, 33, 19, 34, 35, 21, 36, 22, *Wang_H_17, *Wang_H_18,
*Fitzpatrick_15, 37, 38, 26, *sslee2, *lunts_2017, 39, 40, 41, 42, *torroba_2,
44, 45, 46, *paper_2, *paper_3, *paper_4, *paper_5, *paper_6, *paper_odd, 53]
that at a QCP the self-energy at $T=0$ is singular in the frequency domain and
scales as $\Sigma(\omega)\propto\omega^{1-\gamma}$, where the exponent
$\gamma\ll 1$ in weakly anisotropic 3D systems, $\gamma=1/3$ at an Ising-
nematic and Ising-ferromagnetic QCP in 2D, $\gamma\approx 1/2$ at a 2D QCP
towards spin or charge density-wave order with a finite momentum, and
$\gamma=2$ for an electron-phonon problem. It is tempting to associate
$1+d\Sigma/d\omega$ with $m^{*}/m$ and associate $\omega$ with $T$. By this
reasoning, the leading term in the specific heat at small $T$,
$C(T)\propto(m^{*}/m)T$, should scale as $T^{1-\gamma}$, i.e. as $T^{2/3}$ at
an Ising-nematic QCP, as $T^{1/2}$ at a density-wave QCP, and as $1/T$ for
critical electron-phonon problem (although this last behavior obviously cannot
extend to $T=0$). Our goal is to check whether these formulas hold in
microscopic calculations.

A more specific motivation for our work is to clarify recent studies of the
free energy for critical fermion-boson systems [54, 55, 1, 56, 53, 57]. Some
of us and others recently analyzed [46, *paper_2, *paper_3, *paper_4,
*paper_5, *paper_6, *paper_odd] the interplay between non-FL in the normal
state and superconductivity within an effective low-energy model of fermions
with boson-mediated dynamical 4-fermion interaction $V(\Omega)\propto
1/|\Omega|^{\gamma}$ (the $\gamma-$model [58]). This model describes non-FL in
the normal state and superconductivity. Both are universal phenomena in the
sense that they come from fermions with energies well below the upper energy
cutoff of the model $\Lambda$. The condensation energy – the difference
between the free energy of a superconductor and of a would be normal state at
the same $T$ , is also independent on $\Lambda$ (Ref. [46, *paper_2, *paper_3,
*paper_4, *paper_5, *paper_6, *paper_odd]). However, the free energy of the
$\gamma-$model in the normal state is non-universal, even if we subtract its
value at $T=0$. Namely, its leading $T-$dependent term scales as $\Lambda
T^{1-\gamma}$ (Ref. [53]). The corresponding specific heat is then
$C(T)\propto\Lambda/T^{\gamma}$, in variation with the estimate based on the
self-energy. The authors of Refs. [1, 56] argued that the dependence of the
free energy and the specific heat on $\Lambda$ is spurious and has to be
regularized by adding the counter term to the free energy [1, 56], which
cancels out $\Lambda$ dependence. Once this is done, the regularized specific
heat becomes independent on $\Lambda$ and scales as $T^{1-\gamma}$, as
expected based on the self-energy. However, the regularization comes with the
cost: the prefactor in $C(T)\propto T^{1-\gamma}$ turns out to be negative for
$\gamma\geq 1$. 111The authors of [1, 56] found the prefactor to be negative
for $\gamma\geq 2$, which they only considered. The authors of [53] argued
that the prefactor is negative for $\gamma>1$.

Taken at a face value, a negative $C(T)$ would imply that the system becomes
unstable below a certain $T_{cr}$, when a negative $T^{1-\gamma}$ term, coming
from fermion-boson interaction, exceeds a positive $O(T)$ contribution to
$C(T)$ from free fermions. A potential resolution would be that this
instability is preempted by superconductivity, but it turns out that
$T_{cr}>T_{c}$ (Refs [1, 56, 53]).

In this work, we analyze free energy and specific heat within the full
fermion-boson model using the Luttinger-Ward-Eliashberg formula [60, 61] for
the variational free energy in the normal state. We assume that
superconductivity is suppressed, and extend the normal state analysis down to
small $T$. Luttinger and Ward argued [60] that the free energy of a system of
fermions with 4-fermion interaction can be expressed diagrammatically by
collecting skeleton diagrams with fully dressed fermionic propagators and
using conventional rules of the diagrammatic technique, but one has to add to
free energy the term $F_{el}$, which explicitly contains fermionic self-energy
$\Sigma$ (see below). This additional term is constructed such that the
stationary condition $\delta F/\delta\Sigma=0$ reproduces the diagrammatic
series for the self-energy. Eliashberg extended Luttinger-Ward approach to the
case of electron-phonon interaction. He argued that the free energy for such a
system is obtained by collecting skeleton diagrams with fully dressed
fermionic and bosonic propagators, and contains a second extra term $F_{bos}$,
which depends on the bosonic polarization operator $\Pi$ (the bosonic self-
energy) and is constructed such that the stationary condition $\delta
F/\delta\Pi=0$ reproduces the conventional diagrammatic series for $\Pi$.

We analyze the electron-phonon model and different electronic models in which
a certain collective bosonic mode becomes massless at a QCP. For these models,
the low-energy behavior of fermions and their soft collective excitations is
captured within an effective fermion-boson model, in which a collective mode
becomes an independent degree of freedom, coupled to fermions.

The full variational free energy of fermion-boson model is
$F=F_{bos}+F_{el}+F_{int}$, where $F_{int}$ is the sum of skeleton diagrams.
We assume, following earlier works, that both phonons and soft collective
modes are slow compared to dressed fermions, either because a velocity of a
boson is small compared to that of a dressed fermion, or because collective
bosons are Landau overdamped, and that the smallness of an (effective)
velocity of a boson is controlled by a dimensionless parameter $\lambda_{E}$,
often called Migdal-Eliashberg parameter (more on this below). In practical
terms, the fact that the bosons are slow compared to fermions means that
corrections to fermion-boson vertex are small as in the processes identified
with vertex corrections fermions are forced to vibrate at boson frequencies,
far away from their own resonance. This makes higher-loop terms in the
skeleton loop expansion of $F_{int}$ small compared to the one-loop term, and
we keep only this term in $F_{int}$.

The free energy of the $\gamma-$model is $F_{\gamma}=F_{el}+F_{int}$. Like we
said, the specific heat obtained from $F_{\gamma}$ depends on the cutoff. Our
goal is to understand the role of $F_{bos}$, specifically (i) whether it acts
as the counter term, which eliminated the cutoff dependence of $F_{\gamma}$,
and (ii) whether it also affects the universal part of $C(T)$.

We show below that for any $\gamma$, the full free energy $F$ near a QCP is

$F=-2\pi
TN_{F}\sum_{m}\rvert\omega_{m}\rvert+\frac{T}{2}\sum_{q}\log[{-D^{-1}_{q}}]$
(1)

where the first term is the contribution from free fermions, and in the second
$q\equiv({\bf q},\Omega_{m})$, $\Omega_{m}=2\pi Tm$, and $D_{q}$ is the fully
dressed bosonic propagator. This result holds even if we include thermal
fermionic self-energy, which near a QCP has to be computed self-consistently
beyond Eliashberg theory [37]. The fermion-boson interaction is present in
$D_{q}$ as it contains the bosonic polarization bubble. We compute the
specific heat from (1) and compare it with the one of the regularized
$\gamma-$model. We then address the issues (i) and (ii). Regarding (i), we
find that $F_{bos}$ cancels the cutoff-dependent terms in $F_{\gamma}$, i.e.
it provides the physical realization of the counter term. On (ii), the result
depends on whether $\gamma<1$ or $\gamma>1$. For $\gamma<1$ (Ising-nematic and
related models), the contribution from $F_{bos}$ to the universal part of the
specific heat is of the same order as the one from the regularized
$F_{\gamma}$. The two contributions differ by a $\gamma-$dependent factor,
e.g., by $3/2$ for $\gamma=1/3$. For $\gamma>1$, including the electron-phonon
case ($\gamma=2$), the contribution from $F_{bos}$ to the universal part of
$C(T)$ coincides with that from free bosons (a $T$-independent term for
$\gamma=2$). The full $F$ in this case (Eq. (1)) is the sum of contributions
from free bosons with the dressed mass and from the regularized $\gamma$
model. The last contribution is negative at small $T$, in agreement with Refs.
[1, 56]. However, the negative term appears in Eq. (1) as the subleading term
in the expansion $\log[{-D^{-1}_{q}}]$ to first order in the dynamical part of
the bosonic polarization, while the positive contribution from free bosons is
the leading term. The expansion holds in powers of the Migdal-Eliashberg
parameter $\lambda_{E}$ (defined below), and we argue that as long as
$\lambda_{E}\leq 1$, i.e., as long as the theory is under control, the
specific heat is positive. Based on this, we argue that the normal state
remains stable at a QCP and at any distance away from it. In this last respect
our conclusions are different from those in Refs. [1, 56].

The structure of the paper is as follows. In Sec. II we present the generic
Luttinger-Ward-Eliashberg expression for the free energy, briefly discuss the
Eliashberg theory, and use it to obtain the expressions for the full $F$, Eq.
(1), and for $F_{\gamma}$ in the purely fermionic $\gamma-$model. In Sec. III
we compare the two expressions for the Ising-nematic model in 2D. Here we also
show that the result for $F$ does not change if we include thermal self-
energy, which has to be calculated outside the Eliashberg theory, and estimate
the strength of vertex corrections once we include thermal self-energy. In
Sec. IV we consider antiferromagnetic QCP in 2D. In Sec. V we consider an
electron-phonon system near a QCP. Here we also discuss, in Sec. V.1, the
regularization of $F_{\gamma}$ from physical perspective. In Sec. VI we extend
the $\gamma=2$ model to arbitrary $\gamma$ between 1 and 2 and compute the
specific heat. We show that the full specific heat is positive, as long as
$\lambda_{E}\leq 1$, despite that the contribution from the regularized
$F_{\gamma}$ is negative. We present our conclusions in Sec. VII. Some
technical details of the calculations are presented in the Appendices.

##  II Free energy and specific heat

The variational free energy for interacting fermions has been derived by
Luttinger and Ward [60] and extended to fermion-boson systems by Eliashberg
[61] (see also [62, 63]). For more recent studies of variational free energy,
see Refs. [64, 65, *cmg_long, 67, 68, 54, 55, 53]. The free energy per unit
volume is the sum of the fermionic contribution, the bosonic contribution, and
the contribution due to fermion-boson interaction:

$F=F_{el}+F_{bos}+F_{int}.$ (2)

The fermionic part is

$F_{el}=-2T\sum_{k}\log{G^{-1}_{k}}+2iT\sum_{k}\Sigma_{k}G_{k}$ (3)

where $k\equiv({\bf k},\omega_{m})$, $\omega_{m}=\pi T(2m+1)$,
$T\sum_{k}=T\sum_{m}\int d{\bf k}/(2\pi)^{d}$ ($d$ is the spatial dimension),
${\tilde{\Sigma}}_{k}=\omega_{m}+\Sigma_{k}$, where $\Sigma_{k}$ is the self-
energy, and $G_{k}=(i{\tilde{\Sigma}}_{k}-\epsilon_{\bf k})^{-1}$ is the
Green’s function.

The bosonic part is

$F_{bos}=\frac{T}{2}\sum_{q}\left(\log[{-D^{-1}_{q}}]+\Pi_{q}D_{q}\right)$ (4)

where $q\equiv({\bf q},\Omega_{m})$, $\Omega_{m}=2\pi Tm$, $\Pi_{q}$ is the
bosonic self-energy, and $D_{q}=((D^{0}_{q})^{-1}-\Pi_{q})^{-1}$ is the
dressed bosonic propagator. For the bare bosonic propagator we set
$D^{0}_{q}=-D_{0}/(\Omega^{2}_{m}+\omega^{2}_{D})$ for the electron-phonon
case, where $\omega_{D}$ is a bare Debye frequency,
$D^{0}_{q}=-D_{0}/((\Omega_{m}/c)^{2}+{\bm{q}}^{2}+m^{2})$ for the Ising-
nematic case, and
$D^{0}_{q}=-D_{0}/((\Omega_{m}/c)^{2}+({\bm{q}}-{\bm{Q}})^{2}+m^{2})$ for the
antiferromagnetic case, where ${\bm{Q}}=(\pi,\pi)$, $c$ is of order of Fermi
velocity $v_{F}$, and $m$ is a bare boson mass. We set the lattice constant
$a=1$. For the last two cases, the $\Omega^{2}_{m}$ term in the bosonic
propagator can be neglected as for relevant $\Omega_{m}$ it is parametrically
smaller than the Landau damping term from $\Pi_{q}$. On the contrary, for the
electron-phonon case, $\Omega^{2}_{m}$ term is more relevant than the Landau
damping.

Finally, the interaction part is

$F_{int}=-T^{2}\sum_{k,k^{\prime}}g^{2}_{|{\bf k}-{\bf
k}^{\prime}|}G_{k}D(k-k^{\prime})G_{k^{\prime}}+....$ (5)

where $g_{\bm{q}}$ is the Yukawa coupling. The dots in (5) stand for higher-
order contributions, which account for vertex corrections (see Ref. [65,
*cmg_long, 69] for the discussion on higher-order terms in the loop expansion
of $F_{int}$) We assume, following [61], that vertex corrections can be
neglected (more on this below). For simplicity, we also approximate
$g_{\bm{q}}$ by $g$. We refer to the free energy described by Eqs. (2) - (5)
as the Eliashberg free energy.

Figure 1: Self-energy of (a) electron $i\Sigma(k)$ and (b) boson $-\Pi(q)$
(the polarization bubble). The solid (waggle) lines denote the dressed
electron (boson) Green’s functions. The polarization bubble contains factor of
$2$ from the spin degeneracy.

The stationary solutions for $\Sigma_{k}$ and $\Pi_{k}$ are obtained from
$\delta F/\delta\Sigma_{k}=0$ and $\delta F/\delta\Pi_{q}=0$. They give rise
to two Eliashberg equations for fermionic and bosonic self-energies [60, 61,
64, 65, *cmg_long, 69, 54, 55] (see Fig. 1)

$\displaystyle\Sigma_{k}$ $\displaystyle=$ $\displaystyle-
iT\sum_{q}g^{2}G_{k-q}D_{q},$ (6) $\displaystyle\Pi_{q}$ $\displaystyle=$
$\displaystyle 2g^{2}T\sum_{k}G_{k}G_{k-q},$ (7)

where the factor $2$ in Eq. (7) accounts for the spin degeneracy. These
equations are the same as one obtains diagrammatically, without invoking the
free energy. We emphasize in this regard that the diagrammatic loop expansion
with full $G$ and full $D$ holds only for $F_{int}$. The terms $F_{el}$ and
$F_{bos}$ are additional contributions to the free energy, constructed to
reproduce Eqs. (6) and (7) as stationary conditions for the full $F$.

Below we will analyze free energy in equilibrium, when $\Sigma_{k}$ and
$\Pi_{q}$ obey Eqs. (6) and (7). One can easily check that in this situation

$T/2\sum_{q}\Pi_{q}D_{q}=iT\sum_{k}\Sigma_{k}G_{k}$ (8)

because both expressions describe the same skeleton diagram, see Fig. 2. Along
the same lines,

$F_{int}=-iT\sum_{k}\Sigma_{k}G_{k}.$ (9)

Using these two expressions, we obtain

$F=-2T\sum_{k}\log{G^{-1}_{k}}+2iT\sum_{k}\Sigma_{k}G_{k}+\frac{T}{2}\sum_{q}\log[-{D^{-1}_{q}}]$
(10)

and separately

$F_{el}+F_{int}=-2T\sum_{k}\log{G^{-1}_{k}}+iT\sum_{k}\Sigma_{k}G_{k}$ (11)

Figure 2: Thermodynamic potential due to fermion-boson interaction. The three
diagrams are equivalent if we substitute the self-energy $\Sigma$ and the
polarization $\Pi$ in Eqs. (6) and (7).

###  II.1 Eliashberg theory

The Eliashberg formula for the free energy is valid when bosons are slow
compared to the fermions, either because $\omega_{D}\ll E_{F}$ in the
electron-phonon problem, or because the collective boson is Landau overdamped.
An extension to $N\gg 1$ fermionic flavors, which individually interact with a
boson, enhances the magnitude of the Landau damping term and increases the
applicability range of the Eliashberg theory [15, 4, 5, 22, *Wang_H_17,
*Wang_H_18, *Fitzpatrick_15, 19, *rech_2006_1, 70, 71, 72].

The condition that the bosons are slow compared to the fermions allows one to
factorize the momentum integration along and transverse to the Fermi surface
because in all three cases that we consider, the typical transverse momenta
are much smaller than typical longitudinal momenta (we illustrate this in Fig.
3). To obtain the leading contribution to the r.h.s. of (6) one can then
integrate over transverse $q_{\perp}$ in the fermionic propagator and over
$q_{\parallel}$ in the bosonic propagator with ${\bf q}$ connecting points on
the Fermi surface [64]. Integrating over momenta this way and extending both
integrations to infinity 222We assume, as in previous works on metallic QCP,
that fermionic bandwidth $W$ is the largest scale of the problem, and neglect
terms, which are small in $g/W$., one obtains a purely dynamical self-energy

$\Sigma_{k}=\Sigma(\omega_{n})=\pi
T\sum_{m}{\text{sgn}}(\omega_{n}+\Omega_{m})D_{loc}(\Omega_{m}),$ (12)

At $T=0$,

$\Sigma(\omega_{n})=\int_{0}^{\omega_{n}}D_{loc}(\Omega_{m})d\Omega_{m}$ (13)

The form of $D_{loc}$ is model-specific, but in all cases we have at a QCP

$D_{loc}(\Omega_{m})=\left(\frac{\bar{g}}{|\Omega_{m}|}\right)^{\gamma}$ (14)

where $\gamma=1/3$ for the Ising-nematic case, $\gamma=1/2$ for the
antiferromagnetic case, and $\gamma=2$ for the electron-phonon case. The
coupling ${\bar{g}}$ is expressed via $g$ ( see Secs. III,IV,V below). Away
from the QCP, Eq. (14) is modified to

$D_{loc}(\Omega_{m})=\left(\frac{\bar{g}^{2}}{\Omega^{2}_{m}+M^{2}}\right)^{\gamma/2}$
(15)

where $M\sim m^{3}$ for Ising-nematic case, $M\sim m^{2}$ for
antiferromagnetic case, and $M=\bar{\omega}_{D}$ (renormalized Debye
frequency) for the electron-phonon case. We assume that $M$ and $\omega_{D}$
do not depend on temperature, or, more accurately, that their temperature
dependence yields smaller $C(T)$ compared to what we find below.

Because $\Sigma_{k}$ in (12) does not depend on $\epsilon_{\bm{k}}$, one can
explicitly integrate over momentum in Eq. (3) using $\int
d^{d}k/(2\pi)^{d}=N_{F}\int d\epsilon_{\bm{k}}$, where $N_{F}$ is the density
of states at the Fermi level per spin component. The integration yields

$\displaystyle-2T\sum_{k}\log{G^{-1}_{k}}=-2\pi
TN_{F}\sum_{m}\left(|\omega_{m}|+|\Sigma(\omega_{m})|\right)$ $\displaystyle
2iT\sum_{k}\Sigma_{k}G_{k}=2\pi TN_{F}\sum_{m}|\Sigma(\omega_{m})|$ (16)

Combing the two contributions, we find that the self-energy cancels out and
$F_{el}$ retains the same as for free fermions:

$F_{el}=-2\pi TN_{F}\sum_{m}|\omega_{m}|$ (17)

We emphasize that this holds only if $\Sigma_{k}$ does not depend on
$\epsilon_{\bm{k}}$. For a generic momentum and frequency dependent
$\Sigma_{k}$, $F_{el}$ does depend on the fermionic self-energy.

Applying the same procedure to Eqs. (10) and (11) we obtain [64, 65,
*cmg_long, 55, 53]

$F=-2\pi
T\sum_{m}|\omega_{m}|+\frac{T}{2}\sum_{q}\log[-{D^{-1}_{q}}]=F_{free}+\frac{T}{2}\sum_{q}\log[-{D^{-1}_{q}}]$
(18)

and

$F_{el}+F_{int}=-2\pi TN_{F}\sum_{m}|\omega_{m}|-\pi
TN_{F}\sum_{m}|\Sigma(\omega_{m})|$ (19)

Note that the self-energy $\Sigma(\omega_{m})$ cancels out in $F$, and that
the dependence on fermion-boson interaction comes about because $D_{q}$
depends on the polarization $\Pi(q)$.

At $T=0$,
$\Sigma(\omega_{m})=({\bar{g}}^{\gamma}/(1-\gamma))|\omega_{m}|^{1-\gamma}{\text{sgn}}\omega_{m}\equiv\omega_{0}^{\gamma}|\omega_{m}|^{1-\gamma}{\text{sgn}}\omega_{m}$,
where $\omega_{0}={\bar{g}}/(1-\gamma)^{1/\gamma}$. This holds for $\gamma<1$.
For $\gamma>1$, one has to add the contribution from the lower limit in (13).
This last contribution scales as $1/M^{\gamma-1}$ and diverges at $M\to 0$.
However, it does not contribute to the specific heat, as one can explicitly
verify. At a finite $T$, the self-energy becomes a function of a Matsubara
number, and there appears a separate singular contribution $O(1/M^{\gamma})$
from zero bosonic Matsubara frequency. This last contribution requires special
attention, and we discuss it in some detail in Sec. III.

Figure 3: Typical transverse and longitudinal momenta, $q_{\bot}$ and
$q_{\parallel}$, for the Ising-nematic case. At small ${\bf q}$,
$q_{\parallel}\gg q_{\bot}$.

###  II.2 A purely electronic $\gamma$-model

The $\gamma$-model is designated to reproduce some low-energy properties of
the fermion-boson system (more specifically, non-FL and superconductivity). It
is a fermion-only model in which $D_{loc}(\Omega_{m})$ from (15) plays the
role of an effective dynamical 4-fermion interaction [46, *paper_2, *paper_3,
*paper_4, *paper_5, *paper_6, *paper_odd]. The model allows one to analyze the
interplay between non-FL and pairing by solving coupled Eliashberg equations
for the dynamic fermionic self-energy and the dynamic pairing vertex
$\Phi(\omega_{m})$. In a more common and convenient formulation, these
equations are re-expressed in terms of the superconducting gap function
$\Delta(\omega_{m})$ and the inverse quasiparticle residue $Z(\omega_{m})$. By
construction, the model contains only the fermions, and its free energy in the
normal state is $F_{\gamma}=F_{el}+F_{int}$, given by (19):

$F_{\gamma}=-2\pi
TN_{F}\sum_{m}|\omega_{m}|-\pi^{2}T^{2}N_{F}{\bar{g}}^{\gamma}\sum_{m,m^{\prime}}\frac{{\text{sgn}}(\omega_{m}\omega_{m^{\prime}})}{((\omega_{m}-\omega_{m^{\prime}})^{2}+M^{2})^{\gamma/2}}$
(20)

The summation over $m$ is confined to frequencies below the upper energy
cutoff $\Lambda$ of the $\gamma-$model. In practice, this implies that the
summation holds over $M_{f}$ positive and $M_{f}$ negative fermionic Matsubara
frequencies $\omega_{n}=\pi T(2n+1)$ ($-M_{f}<n<M_{f}-1)$. The relation
between $M_{f}$ and $\Lambda$ can be obtained by comparing the exact sum of
$|\omega_{m}|$ with Euler-Maclauren formula, in which the integral is cut by
$\Lambda$. The comparison yields [53]
$4\pi^{2}T^{2}M^{2}_{f}=\Lambda^{2}+\pi^{2}T^{2}/3$, hence

$M_{f}={\tilde{\Lambda}}\left(1+\frac{1}{24{\tilde{\Lambda}}^{2}}+..\right)$
(21)

where ${\tilde{\Lambda}}=\Lambda/(2\pi T)$.

Applying this procedure to both terms in (20), we obtain [53] at $T\gg M$

$\displaystyle F_{\gamma}$ $\displaystyle=$ $\displaystyle-
N_{F}\left(\Lambda^{2}-{\bar{g}}^{\gamma}\Lambda^{2-\gamma}\frac{2(1-2^{-\gamma})}{(1-\gamma)(2-\gamma)}\right)$
(22) $\displaystyle-$ $\displaystyle N_{F}\pi
T\Lambda\left(\frac{\bar{g}}{M}\right)^{\gamma}-N_{F}\Lambda{\bar{g}}^{\gamma}(2\pi
T)^{1-\gamma}\zeta(\gamma)$ $\displaystyle+$ $\displaystyle
N_{F}\left(\frac{3}{2}{\bar{g}}^{\gamma}(2\pi
T)^{2-\gamma}\zeta(\gamma-1)-\frac{1}{3}\pi^{2}T^{2}\right),$

where $\zeta(s)$ is the Riemann zeta function. The first two terms in (22)
constitute the free energy at $T=0$. The next one, with $M$ in the
denominator, comes from the thermal piece in $\Sigma(\omega_{m})$ in (19), or,
equivalently, from the $m=m^{\prime}$ term in (20). The next term comes from
$\omega_{m},\omega_{m^{\prime}}\sim\Lambda$, but
$\omega_{m}-\omega_{m^{\prime}}=O(T)$. The last term is the combination of
cutoff independent contributions from both terms in (20).

The specific heat $C_{\gamma}(T)=-Td^{2}F_{\gamma}/d^{2}T$ is

$\displaystyle C_{\gamma}(T)$ $\displaystyle=$ $\displaystyle 2\pi\Lambda
N_{F}\left(\frac{\bar{g}}{2\pi
T}\right)^{\gamma}\gamma(\gamma-1)\zeta(\gamma)$ (23)
$\displaystyle+\frac{2}{3}\pi N_{F}\left(\pi
T-\frac{9}{4}{\bar{g}}^{\gamma}(2\pi
T)^{1-\gamma}(\gamma-2)(\gamma-1)\zeta(\gamma-1)\right)$

The first term in (23) is parametrically larger than the other two since it is
proportional to $\Lambda$. This term is positive, but depends linearly on the
upper energy cutoff. The second term is a universal contribution to $C(T)$.
This term is positive for $\gamma<1$, but becomes negative at small
$T<T_{0}=[{\bar{g}}/(2\pi)][9(\gamma-2)(\gamma-1)\zeta(\gamma-1)/2]^{1/\gamma}$
for $\gamma>1$, when $(\gamma-2)(\gamma-1)\zeta(\gamma-1)>0$. The temperature
$T_{0}$ increases with $\gamma$ up to $\gamma=3$.

The authors of [55, 1, 56] argued that the dependence of $C(T)$ on the cutoff
is spurious and must be eliminated by a proper regularization. They suggested
that this is achieved by adding to the r.h.s. of (20) the term

$\pi^{2}T^{2}N_{F}{\bar{g}}^{\gamma}\sum_{m,m^{\prime}}\frac{1}{((\omega_{m}-\omega_{m^{\prime}})^{2}+M^{2})^{\gamma/2}}$
(24)

This additional term cancels out all $\Lambda$-dependent terms in $F_{\gamma}$
in (22) and changes the prefactor for the universal $T^{2-\gamma}$ term. The
regularized free energy, which we label as ${\bar{F}}_{\gamma}$, is

$\displaystyle{\bar{F}}_{\gamma}=N_{F}\left({\bar{g}}^{\gamma}(2\pi
T)^{2-\gamma}\zeta(\gamma-1)-\frac{1}{3}\pi^{2}T^{2}\right)$ (25)

This yields a universal, cutoff-independent specific heat
${\bar{C}}_{\gamma}(T)$. At a QCP,

${\bar{C}}_{\gamma}(T)=2\pi N_{F}\left(\frac{\pi}{3}T-{\bar{g}}^{\gamma}(2\pi
T)^{1-\gamma}(\gamma-2)(\gamma-1)\zeta(\gamma-1)\right)$ (26)

For $\gamma<1$, all terms in (26) are positive as $\zeta(\gamma-1)$ is
negative. For $\gamma=1/3$,

${\bar{C}}_{1/3}(T)=2\pi N_{F}\left(\frac{\pi}{3}T+0.172{\bar{g}}^{1/3}(2\pi
T)^{2/3}\right)$ (27)

For $\gamma>1$, $C(T)$ given by (26) is still negative at small $T$ because
$(\gamma-2)(\gamma-1)\zeta(\gamma-1)>0$. For $\gamma\to 2$,

${\bar{C}}_{2}(T)=2\pi N_{F}\left(\frac{\pi}{3}T-\frac{{\bar{g}}^{2}}{2\pi
T}\right)$ (28)

###  II.3 Underlying fermion-boson model

We now return back to the underlying fermion-boson model, in which there is an
additional bosonic contribution to the free energy, and check whether the
effect of $F_{bos}$ is the same as of the extra term (24), which regularizes
$F_{\gamma}$.

The full free energy $F=F_{el}+F_{bos}+F_{int}=F_{\gamma}+F_{bos}$ is given by
Eq. (18) as the sum of the free-fermion contribution and the one expressed via
the full bosonic propagator. In contrast, $F_{\gamma}$, given by Eq. (19),
depends explicitly on the fermionic self-energy. We now study the relation
between these two expressions. We show that the outcome depends on the type of
a QCP. To see this, we consider separately Ising-nematic QCP,
antiferromagnetic QCP, and a QCP of an electron-phonon system.

##  III Ising-nematic QCP

We consider a 2D system. The bare bosonic propagator has the Ornstein- Zernike
form $D^{0}_{q}=-D_{0}/(q^{2}+m^{2})$. The static part of $\Pi(q)$
renormalizes $D_{0}$ and $m$. We assume that these renormalizations are
already incorporated into $D^{0}_{q}$. The dynamical part of $\Pi(q)$ accounts
for the Landau damping: $\Pi({\bf q},\Omega_{m})-\Pi({\bf
q},0)=(1/D_{0})\alpha|\Omega_{m}|/|{\bf q}|$. For a circular Fermi surface,
$\alpha=g^{*}k_{F}/(\pi v^{2}_{F})$, where $g^{*}=g^{2}D_{0}$ has the
dimension of energy and plays the role of an effective fermion-boson
interaction. The dressed bosonic propagator is

$D_{q}=-\frac{D_{0}}{\rvert{\bm{q}}\rvert^{2}+m^{2}+\alpha\frac{|\Omega_{m}|}{\rvert{\bm{q}}\rvert}}.$
(29)

Integrating over one momentum component and comparing with
$D_{{\text{loc}}}(\Omega_{m})$ from (14) for $\gamma=1/3$, we obtain
${\bar{g}}=(g^{*})^{2}/(162\sqrt{3}\pi^{2}E_{F})$ and
$\omega_{0}=(27/8){\bar{g}}=(g^{*})^{2}/(48\sqrt{3}\pi^{2}E_{F})$. The mass
$m$ is related to $M$ in the $\gamma-$model by
$M=32\pi/(81\sqrt{3})(mv_{F})^{3}/(g^{*}E_{F})$.

Substituting $D_{q}$ into (18), subtracting from $\log[-{D^{-1}_{q}}]$ its
static part, which does not contribute to the specific heat, and integrating
over momentum (the integral converges), we obtain at a QCP (i.e., at $m=0$)

$\displaystyle F=F_{I-N}$ $\displaystyle=$ $\displaystyle
F_{\text{free}}+\frac{\alpha^{2/3}}{4\pi\sqrt{3}}(2\pi
T)^{5/3}\sum_{1}^{M_{b}}n^{2/3}$ (30) $\displaystyle=$ $\displaystyle
F_{\text{free}}+\frac{\alpha^{2/3}}{4\pi\sqrt{3}}(2\pi
T)^{5/3}H_{-2/3}(M_{b})$

where $F_{\text{free}}=-N_{F}(\Lambda^{2}+\pi^{2}T^{2}/3)$ is free energy of a
gas of free fermions, and $H_{p}(M_{b})=\sum_{1}^{M_{b}}1/k^{p}$ is the
Harmonic number. The asymptotic expansion of $H_{p}(M_{b})$ at large $M_{b}$
is

$H_{p}(M_{b})=\frac{\left(M_{b}+\frac{1}{2}\right)^{1-p}}{1-p}+\zeta(p)+O(1/M_{b}^{p+1})$
(31)

The relation between $M_{b}$ and and the cutoff $\Lambda$ can be established
in a way similar to the procedure described above for fermions, by evaluating
$\sum_{m=1}^{M_{b}}m$ directly and using Eular-Maclaurin formula with
$\Lambda$ as the upper cutoff of frequency integration. This yields

$M_{b}+\frac{1}{2}={\tilde{\Lambda}}\left(1+\frac{1}{24{\tilde{\Lambda}}^{2}}+..\right).$
(32)

Substituting into (30), we obtain

$\displaystyle F_{I-N}$ $\displaystyle=$ $\displaystyle-
N_{F}\Lambda^{2}+\frac{\sqrt{3}(\alpha A)^{2/3}}{20\pi}\Lambda^{5/3}$ (33)
$\displaystyle-\frac{\pi^{2}}{3}N_{F}T^{2}+\frac{\alpha^{2/3}}{4\pi\sqrt{3}}(2\pi
T)^{5/3}\zeta(-2/3)$

where $\zeta(-2/3)\simeq-0.155$. Differentiating with respect to $T$, we find
that both the entropy $S_{I-N}(T)=-dF_{I-N}/dT$ and the specific heat
$C_{I-N}(T)=-Td^{2}F_{I-N}/dT^{2}=(2/3)S_{I-N}(T)$ are independent on
$\Lambda$. The specific heat is

$C_{I-N}(T)=\frac{2\pi^{2}}{3}N_{F}T-\frac{5\alpha^{2/3}}{9\sqrt{3}}(2\pi
T)^{2/3}\zeta(-2/3)$ (34)

Re-expressing the result in terms of ${\bar{g}}$ from Eq. (14), we obtain

$C_{I-N}(T)=\frac{2\pi^{2}}{3}N_{F}T-\frac{5}{3}{\bar{g}}^{1/3}(2\pi
T)^{2/3}\zeta(-2/3)$ (35)

Comparing this $C_{I-N}(T)$ with the ${\tilde{C}}_{1/3}(T)$ from (27) (a
regularized specific heat in the $\gamma-$model), we see that they agree up to
a numeric prefactor in the $T^{2/3}$ term (the one in $C_{I-N}(T)$ is larger
by $3/2$. The factor $3/2$ is the difference between the momentum integral of
$\log(-{D^{-1}_{q}})$ with static term subtracted and of $\Pi_{q}D_{q}$, i.e.,
between $\int_{0}^{\infty}dxx\log({1+1/x^{3}})=\pi/\sqrt{3}$ and
$\int_{0}^{\infty}dx/(x^{3}+1)=(2/3)(\pi/\sqrt{3})$).

The analysis at Ising-nematic QCP can be formally extended to other values of
$\gamma$ if we replace $q^{2}$ in $D_{q}$ by $q^{a}$ with some $a>1$. The
exponent $\gamma$ then changes from $1/3$ to $\gamma=(a-1)/(a+1)$, which
ranges between $0$ and $1$. One can easily verify (see Appendix D) that the
interaction contributions to $C_{I-N}(T)$ in the Ising-nematic model and in
the regularized $\gamma$-model have the same structure and just differ by
$1-\gamma$ (the prefactor is larger in $C_{I-N}(T)$). The conclusion here is
that for the Ising-nematic case the result of keeping the bosonic contribution
to the specific heat is almost entirely reproduced by either regularizing the
fermionic part of the free energy, as it was done in [55, 1, 56], or just
eliminating the cutoff-dependent term in the specific heat in (23).

###  III.1 Role of thermal fluctuations

Eq. (18) for the free energy is obtained under the assumption that the
momentum dependence of the self-energy can be neglected for
$\epsilon_{\bm{k}}\sim{\tilde{\Sigma}}(\omega)$. As mentioned above, this is
the case when the typical momenta transverse to the Fermi surface in Eq. (6)
are much smaller than typical momenta along the Fermi surface for the same
frequency. At $T=0$, this holds both at the QCP and away from it. At the QCP
we have $q^{\text{typ}}_{\perp}\sim{\tilde{\Sigma}}(\omega)/v_{F}$ and
$q^{\text{typ}}_{\parallel}\sim(\alpha\omega)^{1/3}$. We use
${\tilde{\Sigma}}(\omega)=\omega+\omega^{1/3}_{0}\omega^{2/3}$, where
$\omega_{0}\simeq(g^{*})^{2}/E_{F}\sim\bar{g}$. A simple analysis shows that
$q^{\text{typ}}_{\parallel}\gg q^{\text{typ}}_{\perp}$ up to
$\omega_{\text{max}}\sim(g^{*}E_{F})^{1/2}\sim{\bar{g}}^{1/4}E^{3/4}_{F}$.
This scale is much larger than the upper cutoff for non-FL behavior,
$\omega_{0}\sim{\bar{g}}$, which is also a typical scale for
superconductivity. Away from a QCP, typical
$q_{\parallel}^{\text{typ}}\sim\text{max}\\{m,(\alpha\omega)^{1/3}\\}$ are
even larger.

At a finite $T$ the self-energy $\Sigma_{k}$ can be split into two parts [5,
30, *metzner_new, 39, 42, *torroba_2, 37]. One is the thermal part,
$\Sigma^{th}_{k}$, which comes from zero bosonic Matsubara frequency, and the
other is the quantum part, $\Sigma^{q}_{k}$, which comes from all non-zero
Matsubara frequencies. For the quantum part, the condition
$q^{\text{typ}}_{\parallel}\gg q^{\text{typ}}_{\perp}$ holds up to
$\omega_{\text{max}}$, and one can evaluate $\Sigma^{q}$ using Eq. (12). For
$T\gg M\sim m^{3}/\alpha$,

$\Sigma^{q}_{k}=\text{sgn}(\omega_{m})\left[{3\over
2}{\bar{g}}^{1/3}\rvert\omega_{m}\rvert^{2/3}(1+O(T/\omega_{m}))+\zeta(1/3){\bar{g}}^{1/3}(2\pi
T)^{2/3}\right].$ (36)

(see Ref. [37] for the analysis of $\Sigma^{q}_{k}$ at all $T$).

For the thermal part, the situation is different: the condition
$q^{\text{typ}}_{\parallel}\gg q^{\text{typ}}_{\perp}$ holds only away from a
QCP, at a finite bosonic mass $m>{\tilde{\Sigma}}_{k}/v_{F}$, where
${\tilde{\Sigma}}_{k}=\omega+\Sigma^{th}_{k}+\Sigma^{q}_{k}$. Under this
condition, we obtain

$\Sigma^{th}_{k}=\text{sgn}(\omega_{m}){g^{*}T\over
4mv_{F}}=\text{sgn}(\omega_{m})\pi T\left({{\bar{g}}\over M}\right)^{1/3}.$
(37)

A straightforward analysis shows that Eq. (37) is valid for
$T<T^{*}\sim(m/k_{F})^{2}E^{2}_{F}/g^{*}$. At the QCP, $T^{*}$ vanishes, and
at any finite $T$, $q^{\text{typ}}_{\parallel}\ll q^{\text{typ}}_{\perp}$. The
thermal contribution to the self-energy then has to be computed differently,
by integrating over both components of momenta in the fermionic propagator.
For the one-loop self-energy this yields

$\Sigma^{th}_{k}=iBG_{k}=\frac{B}{{\tilde{\Sigma}}_{k}+i\epsilon_{k}}$ (38)

where $B=g^{*}T\log({k_{F}/m})$ and
${\tilde{\Sigma}}_{k}=\omega+\Sigma^{q}_{k}+\Sigma^{th}_{k}={\tilde{\Sigma}}^{q}_{k}+\Sigma^{th}_{k}$.
The thermal self-energy diverges, but only logarithmically. It has been argued
[74, 42, *torroba_2, 72] that the renormalization of $m$ at a finite $T$ by
high-energy fermions makes it $T$-dependent, in which case $m$ under the
logarithm is cut by $T/v_{F}$, i.e., $\log({k_{F}/m})$ can be approximated by
$\log({E_{F}/T})$. We follow these works and set $B=g^{*}T\log({E_{F}/T})$.

Figure 4: (a) The imaginary part of $\Sigma_{k}$ for different values of
$\omega_{m}$ at $B=0.1$. (b) $\Sigma_{k}\equiv\Sigma^{\prime}_{k}$ at
$\epsilon_{\bm{k}}=0$ at $B=0$ and $B=0.1$ (dashed and solid lines,
respectively).

The key new feature of $\Sigma^{th}_{k}$ in (38) is that it now depends on
both $\omega_{m}$ and $\epsilon_{\bm{k}}$. Then one has to redo the
integration over $\epsilon_{\bm{k}}$ in the fermionic part of the free energy
in Eq. (10). To do this, we solve Eq. (38) for $\Sigma^{th}$ in terms of
${\tilde{\Sigma}}^{q}=\omega_{m}+{\Sigma}^{q}(\omega_{m})$ and
$\epsilon_{\bm{k}}$. We obtain

$\Sigma^{th}_{k}=\sqrt{B+\left(\frac{{\tilde{\Sigma}}^{q}_{k}+i\epsilon_{\bm{k}}}{2}\right)^{2}}-\frac{{\tilde{\Sigma}}^{q}_{k}+i\epsilon_{\bm{k}}}{2}$
(39)

where we choose the branch cut of the square root along the negative real
axis. A similar expression, but at $\epsilon_{\bm{k}}=0$ and at
${\tilde{\Sigma}}^{q}_{k}\approx\omega_{m}$ has been obtained in [37]. In Fig.
4 (a), we plot the imaginary part of the total self-energy
$\Sigma^{th}_{k}+\Sigma^{q}_{k}$ from Eqs. (36) and (39) as a function of
$\epsilon_{\bm{k}}$ for different $\omega_{m}$. The dependence is linear in
$\epsilon_{\bm{k}}$ at the smallest $\omega_{m}$, with the universal slope
$-1/2$. This renormalizes the dispersion to $\epsilon_{\bm{k}}/2$. At larger
$\omega_{m}$, the renormalization of $\epsilon_{\bm{k}}$ becomes negligible.
The crossover between the two regimes is at $\omega_{m}\sim B^{3/4}$ at the
smallest $T$, and at $\omega_{m}\sim B$ at $T>(g^{*})^{3}/E^{2}_{F}$. In Fig.
4 (b) we plot $\Sigma_{k}$ at $\epsilon_{\bm{k}}=0$, where it is necessarily
real. At the smallest $\omega_{m}\sim T$, $\Sigma_{k}$ scales as
$\omega_{m}^{2/3}$ at $B=0$ and tends to a finite value
$\Sigma_{k}\approx\sqrt{B}$ at a finite $B$.

We now substitute $\Sigma^{th}(\omega_{m},\epsilon_{\bm{k}})$ from Eq. (39)
into ${\tilde{\Sigma}}_{k}={\tilde{\Sigma}}^{q}_{k}+\Sigma^{th}_{k}$ and then
into Eq. (10) and explicitly integrate over $\epsilon_{\bm{k}}$. After some
algebra (see Appendix B for details), we obtain that $\Sigma^{th}$ cancels out
from both terms in Eq. (10) which contain fermionic self-energy. Namely,

$\displaystyle-T\sum_{k}\log{[\epsilon^{2}_{\bf
k}+{\tilde{\Sigma}}^{2}_{k}]}=-2\pi
TN_{F}\sum_{n}|\omega_{m}+\Sigma^{q}(\omega_{m})|$ $\displaystyle
2iT\sum_{k}\Sigma_{k}G_{k}=2\pi TN_{F}\sum_{n}|\Sigma^{q}(\omega_{m})|$ (40)

As a result Eq. (18) holds despite that the thermal self-energy depends on
$\epsilon_{\bm{k}}$. For completeness, we verified explicitly that the
momentum dependence of $\Sigma^{th}$ does not generate a significant momentum
dependence of $\Sigma^{q}$ (which still contains the full self-energy in the
Green’s function in the r.h.s. of (6)).

###  III.2 Strength of vertex corrections

The free energy in (10) is obtained within a self-consistent one-loop
approximation, which neglects the vertex corrections. Several authors argued
[15, 19, *rech_2006_1, 26, 75, 76, 77, 78] that at $T=0$ lowest-order vertex
corrections are generally of order one (or of order $1/N$ in large $N$
theories), but higher-order corrections are $O(1)$ even at large $N$ (Ref.
[26]) and furthermore are logarithmically singular [75], except for special
cases [79]. The logarithms, however, likely modify the quasiparticle residue
but not the exponent $\gamma=1/3$, and hence do not affect the $T^{2/3}$
behavior.

In this section, we estimate the strength of vertex corrections at a finite
$T$. We recall that at a finite bosonic mass $m$, there is a range $T<T^{*}$,
where $q^{\text{typ}}_{\parallel}\gg q^{\text{typ}}_{\perp}$ and the thermal
self-energy $\Sigma^{th}\sim g^{*}T/(mv_{F})$ is obtained by factorizing the
momentum integration between fermionic and bosonic propagators, and a range
$T>T^{*}$, where the factorization does not hold. In this last regime,
$\Sigma^{th}$ generally depends on ${\tilde{\Sigma}}^{q}$ and
$\epsilon_{\bm{k}}$, and is of order $[g^{*}T\log({E_{F}/T)}]^{1/2}$ when it
is larger than ${\text{max}}\\{{\tilde{\Sigma}}^{q},|\epsilon_{\bm{k}}|\\}$.
At $m=0$, $T^{*}=0$, and the last regime holds for all $T$.

For $T<T^{*}$, a simple analysis shows that the leading vertex correction to
fermion-boson coupling $\delta g^{*}/g^{*}\sim T/T^{*}$, i.e., it remains
small in the same $T$ range where one can factorize the momentum integration.
For $T>T^{*}$, a similar analysis shows that $\delta g^{*}/g^{*}\sim
g^{*}T\log({E_{F}/T})/({\tilde{\Sigma}_{k}})^{2}$. This vertex correction is
at most of order one. It is then reasonable to expect that thermal vertex
corrections do not modify $C(T)\propto T^{2/3}$ at a QCP, and, moreover, the
prefactor differs from the one in Eq. (34) at most by a factor $O(1)$.

##  IV Antiferromagnetic QCP in 2D

Figure 5: Fermi surface (blue lines) and hot spots (black dots), connected by
the ordering wave vector ${\bf Q}$ (red arrows). The coordinate frame
$(k_{x},k_{y})$ is used in the main text.

The non-FL physics at a QCP towards spin order with momentum ${\bf
Q}=(\pi,\pi)$ is described by the $\gamma$-model with $\gamma=1/2$, as the
self-energy in hot regions on the Fermi surface (the ones in which both
$\epsilon_{\bm{k}}$ and $\epsilon_{{\bm{k}}+{\bm{Q}}}$ are small) scales as
$\Sigma(\omega_{m})\propto\omega^{1/2}_{m}$ [4, 5] 333The actual value of
$\gamma$ is somewhat different from $1/2$ as corrections to fermion-boson
vertex are logarithmical, and series of these corrections change the exponent
$\gamma$ to $1/2+\epsilon$, where $\epsilon$ is positive, but small
numerically [5, morr_97]. Besides, the dynamical exponent $z$ also flows
exponentially from $z=2$ to a smaller value [34]. It has been argued [28] that
in the absence of a superconducting instability, this flow eventually, at the
lowest energies, brings the system into the basin of attraction of a stable
fixed point with $z=1$. Our analysis is valid at energies where the dynamical
exponent is still $z\approx 2$.. The specific heat in the non-regularized
$\gamma$-model scales as $\Lambda/\sqrt{T}$, and the one in the regularized
$\gamma$-model scales as $\sqrt{T}$. Both are inconsistent with the specific
heat of the underlying fermion-boson model $C_{afm}\propto T\log T$, as we
show below. The reason for the inconsistency is, however, rather banal – the
non-FL behavior, described by the $\gamma=1/2$ model, holds only in hot
regions. Away from these regions the self-energy has a Fermi-liquid form at
the smallest frequencies. Hot fermions contribute most to superconductivity,
and the $\gamma=1/2$ model of hot fermions adequately describes the interplay
between non-FL and pairing. However, the free energy in the normal state is
the combined contribution from fermions over the whole Fermi surface, and the
one from hot fermions is proportional to the total width of the hot regions,
which is small compared to the circumference of the Fermi surface boundary.

To obtain the specific heat, we compute the free energy using Eq. (10). We
assume that $\Sigma_{k}$ depends on frequency and on the position on the Fermi
surface, but not on $\epsilon_{k}$. In this situation, one can still
explicitly integrate over the dispersion in the first two terms in (10). The
result is that the self-energy cancels out, even if it depends on the momentum
along the Fermi surface, and the free energy is given by Eq. (18). As before,
we assume that $D^{0}_{q}$ has Ornstein-Zernike form $D^{0}_{q}=-D_{0}/(({\bf
q}-{\bf Q})^{2}+m^{2})$ and incorporate the renormalizations from the static
polarization bubble $\Pi(0)$ into $m$ and $D_{0}$. We evaluate the dynamical
Landau damping term $\Pi(\Omega_{m})-\Pi(0)=\alpha|\Omega_{m}|$ right at
${\bm{q}}={\bm{Q}}$. The full propagator is

$D_{q}={D_{0}\over({\bf q}-{\bf Q})^{2}+m^{2}+\alpha\rvert\Omega_{m}\rvert}.$
(41)

Integrating over the component of ${\bf q}-{\bf Q}$ along the FS, we obtain

$D_{loc}(\Omega_{m})=\frac{{\bar{g}}^{1/2}}{\left(\Omega_{m}^{2}+M^{2}\right)^{1/4}},$
(42)

where ${\bar{g}}^{1/2}=g^{2}D_{0}/(4\pi v_{F}\sqrt{\alpha})$ and $M\sim
m^{2}$. This is approximately Eq. (15).

Substituting $D_{q}$ from (41) into (18) and integrating over ${\bm{q}}$, we
obtain

$F_{afm}=F_{free}-\frac{3\alpha}{2}T^{2}\sum_{n=1}^{M_{B}}n\log{\frac{nT}{T_{0}}}$
(43)

where the factor of $3$ is due to summation over spin components and
$T_{0}\sim k^{2}_{F}/\alpha$ is a non-universal scale related to the upper
cutoff of the integral over ${\bm{q}}$. Evaluating the frequency sum (see
Appendix C.2.2 for details) and using (32) to relate $M_{b}$ and the energy
cutoff $\Lambda$, we obtain

$\displaystyle F_{afm}$ $\displaystyle=$
$\displaystyle-\Lambda^{2}\left(N_{F}+\frac{3\alpha}{16\pi^{2}}\log{\frac{\Lambda}{2\pi
T_{0}\sqrt{e}}}\right)$ (44)
$\displaystyle+\frac{\pi^{2}}{3}T^{2}\left(-N_{F}+\frac{3\alpha}{8\pi^{2}}\log{\frac{T}{T^{*}_{0}}}\right)$

where $T^{*}_{0}$ differs from $T_{0}$ by a factor $O(1)$. Differentiating
over $T$, we obtain

$C_{afm}(T)=S_{afm}(T)=\frac{2\pi^{2}}{3}T\left(N_{F}+\frac{3}{8\pi^{2}}\alpha\log{\frac{T^{*}_{0}}{Te^{3/2}}}\right)$
(45)

We see that at small $T$ the interaction contribution to specific heat is
larger than the one from free fermions, hence $C_{afm}(T)$ scales as
$T\log(T^{*}_{0}/T)$. This behavior has been extensively discussed in the
context of non-FL behavior of cuprates and heavy fermion materials (see, e.g.,
41 and references therein).

The prefactor $\alpha$ in (44), (56) can be expressed in terms of the
effective electron-boson coupling $g^{*}$ and Fermi velocities at hot spots
${\bm{k}}_{hs}$ and ${\bm{k}}_{hs}+{\bm{Q}}$. We define
$\epsilon_{{\bm{k}}+{\bm{k}}_{hs}}=v_{x}k_{x}+v_{y}k_{y}$,
$\epsilon_{{\bm{k}}+{\bm{k}}_{hs}+{\bm{Q}}}=v_{x}k_{x}-v_{y}k_{y}$
($v^{2}_{F}=v^{2}_{x}+v^{2}_{y}$), see Fig. 5. In these notations [5],

$\alpha=\frac{4g^{*}}{\pi v^{2}_{F}\beta}$ (46)

where $\beta=2v_{x}v_{y}/v^{2}_{F}$. Substituting into (44), we obtain

$C_{afm}(T)=S_{afm}(T)=\frac{2\pi^{2}}{3}N_{F}T\left(1+\frac{3g^{*}}{2\pi^{2}E_{F}\beta}\log{\frac{T^{*}_{0}}{Te^{3/2}}}\right)$
(47)

It is instructive to compare this result with the specific heat in a purely
fermionic model, with and without regularization, but with the self-energy
averaged over the full Fermi surface. The self-energy at a Fermi point
$k=k_{F}$, located at $\delta k_{\parallel}=\delta k$ from a hot spot, is [34,
81]

$\Sigma(\delta
k,\omega_{m})=\frac{3g^{*}}{4v_{F}}T\sum_{\Omega_{m}}\frac{{\text{sign}}(\omega_{m}+\Omega_{m})}{\sqrt{\alpha|\Omega_{m}|+(\beta\delta
k)^{2}}}.$ (48)

At $T=0$,

$\Sigma(\delta k,\omega_{m})=\frac{3g^{*}}{2\pi
v_{F}\alpha}\left(\sqrt{\alpha|\omega_{m}|+(\beta\delta k)^{2}}-\beta|\delta
k|\right){\text{sign}}\omega_{m}$ (49)

At a hot spot, $\Sigma(0,\omega_{m})\propto|\omega_{m}|^{1/2}$, as in the
$\gamma$ model with $\gamma=1/2$. At the same time, the self-energy, averaged
over $\delta k$, scales as $\log({T_{0}/|\omega_{m}|})$. Such a self-energy
emerges in the $\gamma$-model with $\gamma=0+$ (Refs. [32, 33, 35, 38, 82,
83]) and hence proper comparison should be with this model. Indeed,
substituting $\Sigma(\delta k,\omega_{m})$ from (48) into (11), we obtain the
free energy of the $\gamma=0+$ model:

$F_{0+}=F_{free}-\frac{3g^{*}}{2\pi
v^{2}_{F}\beta}T^{2}\sum_{n,n^{\prime}}{\text{sign}}(\omega_{n}\omega_{n^{\prime}})\log{\frac{T^{**}_{0}}{|\omega_{n}-\omega_{n^{\prime}}|}}$
(50)

where $T^{**}_{0}$ is of the same order as $T_{0}$. The regularized free
energy is

${\bar{F}}_{0+}=F_{free}-\frac{3g^{*}}{2\pi
v^{2}_{F}\beta}T^{2}\sum_{n,n^{\prime}}\left({\text{sign}}(\omega_{n}\omega_{n^{\prime}})-1\right)\log{\frac{T^{**}_{0}}{|\omega_{n}-\omega_{n^{\prime}}|}}$
(51)

In Eqs. (50) and (51) the summation is over $-M_{f}<n,n^{\prime}<M_{f}-1$.
Evaluating the sum and relating $M_{F}$ to energy cutoff $\Lambda$ via (21),
we obtain

$\displaystyle F_{0+}$ $\displaystyle=$ $\displaystyle-
N_{F}\Lambda^{2}\left(1+\frac{3g^{*}}{2\pi^{2}E_{F}\beta}\log{2}\right)$ (52)
$\displaystyle+N_{F}\frac{3g^{*}}{2\pi E_{F}\beta}\Lambda
T\log{\frac{T^{**}_{0}}{T}}$
$\displaystyle-\frac{\pi^{2}}{3}N_{F}T^{2}\left(1+\frac{9g^{*}}{4\pi^{2}E_{F}\beta}\log\left(\frac{T_{1}}{T}\right)\right)$

and

$\displaystyle{\bar{F}}_{0+}$ $\displaystyle=$ $\displaystyle-
N_{F}\Lambda^{2}\left(1+\frac{3g^{*}}{2\pi^{2}E_{F}\beta}\log{\frac{0.89\Lambda}{T^{**}_{0}}}\right)$
(53) $\displaystyle-$
$\displaystyle\frac{\pi^{2}}{3}N_{F}T^{2}\left(1+\frac{3g^{*}}{2\pi^{2}E_{F}\beta}\log\left(\frac{T^{*}_{1}}{T}\right)\right),$

where $T_{1}\sim\Lambda$ and $T^{*}_{1}\sim\Lambda^{2}/T_{0}$. Differentiating
with respect to temperature, we obtain

$C_{0+}(T)=\frac{2\pi^{2}}{3}N_{F}\left(\frac{9g^{*}}{4\pi^{3}E_{F}\beta}\Lambda+T\left(1+\frac{9g^{*}}{4\pi^{2}E_{F}\beta}\log\left(\frac{T_{1}}{Te^{3/2}}\right)\right)\right)$
(54)

and

${\bar{C}}_{0+}(T)=S_{0+}(T)=\frac{2\pi^{2}}{3}N_{F}T\left(1+\frac{3g^{*}}{2\pi^{2}E_{F}\beta}\log{\frac{T^{*}_{1}}{Te^{3/2}}}\right)$
(55)

Comparing (47) and (55) we see that the prefactor for the $T\log T$ term is
the same, i.e., for an antiferromagnetic QCP regularization of the free energy
of a purely electronic $\gamma=0+$ model yields the same $C(T)\sim T\log{T}$
as from the bosonic term. This agrees with the analysis of Sec. III of the
$\gamma=1/3$ model, extended to arbitrary $0<\gamma<1$, where we found that
the leading $T^{1-\gamma}$ terms in ${\bar{C}}_{\gamma}(T)$ and in the full
$C(T)$ differ by a factor $1-\gamma$, which tends to one at $\gamma\to 0$. The
specific heat in the non-regularized $\gamma=0+$ model has a parasitic
temperature-independent piece that scales with $\Lambda$. The prefactor for
the universal $T\log{T}$ term in (54) is larger than the one in (55) by the
factor $3/2$ – the same number as we found in Sec.III.

Away from the critical point, $m$ is finite, and at the smallest $T$ the
$\log({1/T})$ dependence in (44) is replaced by $\log({1/m^{2}})$. The total
$C_{afm}(T)$ can then be cast in the form

$C_{afm}(T)=\frac{2\pi^{2}}{3}N_{F}T\left(1+\frac{3g^{*}}{\pi^{2}E_{F}\beta}\log{\frac{k_{F}}{m}}\right)$
(56)

In this Fermi liquid regime, the self-energy, averaged along the Fermi
surface, is $\Sigma_{av}=\lambda_{av}\omega$, where

$\lambda_{av}=\frac{3g^{*}}{\pi^{2}E_{F}\beta}\log{\frac{k_{F}}{m}}$ (57)

Comparing (56) and (57), we see that in a Fermi liquid regime at a finite m,
$C_{afm}=C_{free}(1+\lambda_{av})$, as is expected.

##  V QCP in electron-phonon system

We now analyze the free energy for the case of electrons interacting with an
Einstein boson. We use Eq. (18) as an input and compute the bosonic
contribution to the specific heat. The propagator of an Einstein boson is
$D_{q}=-D_{0}/(\Omega^{2}_{m}+\omega^{2}_{D}+\Pi_{q})$, where $\omega_{D}$ is
the bare Debye frequency and $\Pi_{q}$ (which incorporates the overall factor
$D_{0}$) comes from the interaction with electrons. We set $\omega_{D}$ to be
finite, but much smaller than the Fermi energy $E_{F}=v_{F}k_{F}/2$. We define
the dimensionless coupling $\lambda$ via

$\lambda=\frac{{\bar{g}}^{2}}{\omega^{2}_{D}},~{}~{}{\bar{g}}^{2}=g^{2}N_{F}D_{0}$
(58)

where $g$ is the same as in (5). We consider temperatures smaller than
$\omega_{D}$. At such $T$, the contribution to the specific heat from free
bosons, $C_{bos}\propto e^{-\omega_{D}/T}$, is exponentially small.

For definiteness we consider the 2D case. The form of the 2D polarization
operator of free fermions at small momentum and frequency is well-known:

$\Pi_{q}=2{\bar{g}}^{2}\left(1-\frac{\Omega_{m}}{\sqrt{\Omega^{2}_{m}+(v_{F}q)^{2}}}\right)$
(59)

We assume and then verify that typical $v_{F}q$ are of order $E_{F}$, while
typical $\Omega_{m}$ for the specific heat are of order $T$. For such $v_{F}q$
and $\Omega_{m}$, we can compute $\Pi_{q}$ to linear order in $\Omega_{m}$,
but need a more accurate dependence on $q$. In 2D, the static part of
$\Pi_{q}$ remains equal to $2g^{2}$ for all momenta up to $2k_{F}$ and drops
at larger momentum. The dynamical part changes between small $q$ and $q\sim
k_{F}$, and for arbitrary $q<2k_{F}$ is

$-2{\bar{g}}^{2}\frac{|\Omega_{m}|}{v_{F}q}\frac{2k_{F}}{\sqrt{4k^{2}_{F}-q^{2}}}.$
(60)

Substituting $\Pi_{q}$ at small $\Omega_{m}$ and arbitrary $q<2k_{F}$ into
$D_{q}$, we obtain

$D^{-1}_{q}=\Omega^{2}_{m}+{\bar{\omega}}^{2}_{D}+2{\bar{g}}^{2}\frac{|\Omega_{m}|}{v_{F}q}\frac{2k_{F}}{\sqrt{4k^{2}_{F}-q^{2}}}$
(61)

where ${\bar{\omega}}_{D}=\omega_{D}(1-2\lambda)^{1/2}$ is the dressed Debye
frequency. The dressed ${\bar{\omega}}_{D}$ vanishes at $\lambda=1/2$ for all
momenta $q<2k_{F}$. At larger $q$, the static $\Pi_{q}$ decreases and
${\bar{\omega}}_{D}$ remains finite even at $\lambda=1/2$.

Strictly speaking, the polarization operator has to be computed using full
fermionic propagators, which include the self-energy. This does affect the
static $\Pi$, which is generally different from $2g^{2}$ and has contributions
from fermions with energies of order $E_{F}$, of the order of the upper cutoff
$\Lambda$ in the $\gamma-$model (Ref. [84]). To simplify the discussion, below
we keep the free-fermion result with the understanding that the actual
renormalization of $\omega_{D}$ likely differs somewhat from
$(1-2\lambda)^{1/2}$. The corrections to the Landau damping term
$|\Omega_{m}|/v_{F}q$ is of order of $\lambda_{E}$ and hence are small. We
assume without proof that this holds even when we extend the Landau damping
formula to $q\sim k_{F}$.

We show below that within the regime of validity of the Eliashberg theory,
$\lambda_{E}\ll 1$, the last term in (61) is small compared to the first two.
The corresponding $\gamma$-model then has $\gamma=2$.

The vanishing of the dressed Debye frequency at some finite $\lambda$
($\lambda=1/2$ if we use free-fermion expression for static $\Pi_{q}$) has
been noticed before (see e.g. [13, 1, 85] and references therein), both in 2D
and 3D systems. However, in 3D, ${\bar{\omega}}_{D}$ is not flat for
$q<2k_{F}$, and the dressed ${\bar{\omega}}_{D}$ vanishes at a critical
$\lambda$ only at $q=0$ and scales as $q^{2}$ at small $q$. In this situation
the full bosonic propagator has the same form as in the 3D Ising-nematic
model, and the corresponding $\gamma$-model has $\gamma=0+$, with the
effective interaction

$D_{loc}(\Omega_{m})=\log{\frac{\bar{g}}{|\Omega_{m}|}}$ (62)

In 2D, corrections to free-fermion form of $\Pi_{q}$ also introduce quadratic
momentum dependence of ${\bar{\omega}}_{D}$ around $q=0$, even for an
isotropic fermionic dispersion [86], such that very near QCP critical theory
becomes the same as in a 2D Ising-nematic case. Alternatively, a non-parabolic
fermionic dispersion can also introduce a quadratic term in the bosonic
dispersion. However, the momentum dependence may be weak, resulting in a wide
range around a QCP, where ${\bar{\omega}}_{D}$ can be approximated by
momentum-independent constant, $\omega_{D}(1-2\lambda)^{1/2}$. For a system on
a square lattice, quantum Monte Carlo data show that the minimum of
${\bar{\omega}}_{D}$ is at $M=(\pi,\pi)$ (Ref. [87]). The dispersion is flat
around the minimum, and the overall variation of ${\bar{\omega}}_{D}$ with
momentum is quite small. At the minimum, ${\bar{\omega}}_{D}$ displays
$(1-2\lambda)^{1/2}$ dependence up to $\lambda\sim 0.4$ (Ref. [13]).

In our analysis of the free energy we focus on the regime where the momentum
dependence of ${\bar{\omega}}_{D}$ can be neglected. In this regime
$F=F_{free}+(T/2)\sum_{q}\log(-{D^{-1}_{q}})$, where
$F_{free}=-\pi^{2}T^{2}N_{F}/3$ is the free energy of a free Fermi gas. We
assume and then verify that the largest contribution to the specific heat
comes from the ${\bm{q}}-$independent term in $D_{q}$ and approximate
$\log(-{D^{-1}_{q}})$ by expanding to leading order in the Landau damping term
(which we shall later show is a small correction for an $\gamma>1$)

$\log(-{D^{-1}_{q}})=\log{\left(\Omega^{2}_{m}+{\bar{\omega}}^{2}_{D}\right)}+2\frac{{\bar{g}}^{2}}{v_{F}q}\frac{|\Omega_{m}|}{\Omega^{2}_{m}+{\bar{\omega}}^{2}_{D}}\frac{2k_{F}}{\sqrt{4k^{2}_{F}-q^{2}}}$
(63)

Substituting into the free energy and integrating over $|{\bf q}|<2k_{F}$, we
obtain

$\displaystyle\frac{T}{2}\sum_{q}\log(-{D^{-1}_{q}})=\frac{k^{2}_{F}}{\pi}T\sum_{n=1}^{M_{b}}\log{(4\pi^{2}T^{2}n^{2}+{\bar{\omega}}^{2}_{D})}$
$\displaystyle+\frac{{\bar{g}}^{2}k_{F}}{2\pi
v_{F}}\sum_{n=1}^{M_{b}}\frac{n}{n^{2}+({\bar{\omega}}_{D}/(2\pi T))^{2}}$
(64)

The first term is the free energy of a free Einstein boson with the
renormalized Debye frequency ${\bar{\omega}}_{D}$, the second one is the
contribution from fermion-boson interaction. Evaluating the frequency sums
(see the Appendix C.2.3 for details) and using the relation between $M_{b}$
and the upper energy cutoff $\Lambda$, Eq. (32), we obtain

$\displaystyle
F=N_{F}\left(-\Lambda^{2}+\frac{4E_{F}\Lambda}{\pi}\log{\frac{\Lambda}{e}}+{\bar{g}}^{2}\log{\frac{\Lambda}{{\bar{\omega}}_{D}}}\right)$
$\displaystyle+N_{F}\left[-\frac{\pi^{2}T^{2}}{3}+4TE_{F}\log{\left(1-e^{-\frac{{\bar{\omega}}_{D}}{T}}\right)}+{\bar{g}}^{2}f\left(\frac{{\bar{\omega}}_{D}}{2\pi
T}\right)\right]$ (65)

where

$f(x)=\log{x}-\frac{1}{2}\left(\psi(1+ix)+\psi(1-ix)\right)$ (66)

and $\psi(y)$ is di-Gamma function. We see that the $\Lambda$-dependent terms
in (65) are independent of $T$, and hence do not contribute to the entropy and
the specific heat. Differentiating twice with respect to temperature, we
obtain the total specific heat for the isotropic electron-phonon system

$C_{ep}(T)=\frac{2\pi^{2}}{3}N_{F}\left[T+\frac{6E_{F}}{\pi^{2}}Q\left(\frac{{\bar{\omega}}_{D}}{2\pi
T}\right)\right]$ (67)

where

$\displaystyle Q(x)=\left(\frac{\pi x}{\sinh{\pi
x}}\right)^{2}-\frac{\pi}{2}\lambda_{E}x^{2}\left(x\frac{d^{2}f}{dx^{2}}+2\frac{df}{dx}\right)$
(68)

and $\lambda_{E}={\bar{g}}^{2}/({\bar{\omega}}_{D}E_{F})$ is the Eliashberg
parameter. The Elishberg theory, which neglects vertex corrections, is valid
when $\lambda_{E}$ is small. In Fig. 6, we plot $Q(x)$ for different
$\lambda_{E}$. We see that this function is positive for all $x$. Accordingly,
$C_{ep}(T)$ given by (66) is also positive for all temperatures. We plot
$C(T)$ in Fig. 7.

Figure 6: Function $Q(x)$, given by Eq. 68, for different $\lambda_{E}<1$.

The limiting forms of $C_{ep}(T)$ are

$C_{ep}(T)=\frac{2\pi^{2}}{3}N_{F}T\left(1+\frac{\lambda}{1-2\lambda}\right)$
(69)

at $2\pi T\ll{\bar{\omega}}_{D}$, and

$\displaystyle C_{ep}(T)$ $\displaystyle=$
$\displaystyle\frac{2\pi^{2}}{3}N_{F}\left(T+\frac{6E_{F}}{\pi^{2}}\left(1-\lambda_{E}\frac{{\bar{\omega}}_{D}}{4T}\right)\right)$
(70) $\displaystyle=$
$\displaystyle\frac{2\pi^{2}}{3}N_{F}\left(T+\frac{6E_{F}}{\pi^{2}}-\frac{3{\bar{g}}^{2}}{2\pi^{2}T}\right)$

at $2\pi T\gg{\bar{\omega}}_{D}$, which includes the case
${\bar{\omega}}_{D}\to 0$ at finite $T$. In the two limits, the entropy
$S_{ep}(T)=C_{ep}(T)$ at $2\pi T\ll{\bar{\omega}}_{D}$, and

$S_{ep}(T)=\frac{2\pi^{2}}{3}N_{F}\left(T+\frac{6E_{F}}{\pi^{2}}\log{\frac{T}{{\bar{\omega}}_{D}}}+\frac{3{\bar{g}}^{2}}{2\pi^{2}T}\right)$
(71)

at $2\pi T\gg{\bar{\omega}}_{D}$. Note that in this last limit the entropy is
always positive. It diverges logarithmically at ${\bar{\omega}}_{D}\to 0$ at a
finite $T$.

We now take a more careful look at the expression for the specific heat at
$2\pi T\gg{\bar{\omega}}_{D}$. The first term in the second line in (70) is
the contribution from free fermions, the second is the contribution from free
bosons, but with an effective Debye frequency, renormalized by the interaction
with fermions, and the third term is the direct contribution from the
electron-phonon interaction. This last term is negative and is the same as the
interaction contribution to the specific heat in the regularized
$\gamma$-model, Eq. (28). Without the middle term, the specific heat would
become negative below a certain temperature,
$T=(3/(2\pi^{2}))^{1/2}{\bar{g}}\simeq 0.39{\bar{g}}$, which exceeds the onset
temperature for superconductivity $T_{c}\simeq 0.18{\bar{g}}$ (Ref. [88, 89,
44]). Because of the middle term, however, the full $C(T)$ remains positive.
This holds even at $T\to 0$, as one can see from the first line in (70). The
key here is the condition
$\lambda_{E}={\bar{g}}^{2}/({\bar{\omega}}_{D}E_{F})<1$, which requires one to
treat the case of vanishing dressed Debye frequency as a double limit, in
which $E_{F}$ tends to infinity simultaneously with ${\bar{\omega}}_{D}\to 0$
(Refs. [50, 53]).

Figure 7:  The ratio $C_{ep}(T)/T$, where $C_{ep}(T)$ is the specific heat of
the electron-phonon system, given by Eq. (67). We plot $C(T)/T$ a function
$T/{\bar{\omega}}_{D}$, where ${\bar{\omega}}_{D}$ is the dressed Debye
frequency, for different values of the Eliashberg parameter $\lambda_{E}$. We
set $E_{F}=10\bar{g}$, in which case
$\lambda_{E}=0.1({\bar{g}}/{\bar{\omega}}_{D})^{2}$. The specific heat is
positive at all $T$. At $T\ll\bar{\omega}_{D}$, $C_{ep}(T)/T$ saturates at
$2\pi^{2}/3N_{F}(1+\lambda/(1-2\lambda))$ at $T\gg\bar{\omega}_{D}$,
$C_{ep}(T)/T$ asymptotically approaches its value for free bosons.

The authors of ([55, 1]) argued that the negative prefactor for the $1/T$ term
in Eq. (28) indicates that the normal state becomes unstable below a certain
$T$ despite that the total $C_{ep}(T)$ is positive. Their argument is that the
$T-$independent term in (70), which renders the total $C_{ep}(T)$ positive, is
the contribution from free bosons and as such does not affect the electrons.
Our counter-argument is that both positive and negative parts of $C_{el}(T)$
come from the term in the free energy $(T/2)\sum_{q}\log{[-D^{-1}_{q}]}$, once
one expands it in the dynamical part of $\Pi_{q}$: the positive contribution
is the the zeroth order term and the negative $1/T$ contribution comes from
the first order in the expansion. In our view, this shows that both terms
should be treated on equal footings. Besides, despite the fact that the
leading $T-$ independent part of $C_{el}(T)$ has the same form as the specific
heat of a free massless boson, this term does depend on fermion-boson
interaction as the latter renormalizes the bare $\omega_{D}$ into
${\bar{\omega}}_{D}=\omega_{D}\sqrt{1-2\lambda}$. For $\lambda\approx 1/2$,
$\omega_{D}={\bar{g}}/\sqrt{\lambda}\approx\sqrt{2}{\bar{g}}$ is comparable to
${\bar{g}}$, and without interaction-driven renormalization of $\omega_{D}$
into ${\bar{\omega}}_{D}$ the specific heat of free bosons would be
exponentially small at $T\leq{\bar{g}}$. In this respect, the fermion-boson
coupling gives rise to two effects: it generates a negative $T$-dependent
contribution to $C_{ep}$, and simultaneously gives rise to a much larger,
positive $T-$independent contribution.

For completeness, we also compute the specific heat within our model for
larger $\lambda_{E}$, using the full formula for $\log{D^{-1}_{q}}$ rather
than expanding in the Landau damping in Eq. (63). We find that the specific
heat is positive for all $\lambda_{E}$. We plot $C_{ep}/T$ in Fig. 8. This
result is of limited validity, however, as in the free energy we didn’t
include higher-order terms in the skeleton expansion in $F_{int}$. These terms
are of higher order in $\lambda_{E}$, when $\lambda_{E}$ is small, but are not
small when $\lambda_{E}>1$. Still, we emphasize that within the model we used
here, $C_{ep}(T)$ is positive for all $T$.

Figure 8:  Normalized specific heat of the electron-phonon system,
$C_{ep}(T)/T$, obtained by using the full expression for the free energy
$F=F_{free}+(T/2)\sum_{q}\log(-D_{q}^{-1})$, without expanding
$\log(-D_{q}^{-1})$ in the Landau damping. We used the same parameters as in
Fig. 7. The specific heat remains positive for all values of $\lambda_{E}$.

###  V.1 Physical origin of the regularization of $F_{\gamma}$

We now argue that the interaction-driven renormalization of $\omega_{D}$ is
related to the issue of the regularization of $F_{\gamma}$ in the
$\gamma$-model. To relate the two, we recall that $iT\sum_{k}\Sigma_{k}G_{k}$,
which is the interaction part of $F_{\gamma}$, can be re-expressed as
$(T/2)\sum_{q}\Pi_{q}D_{q}$ (see Eq. (8), where $q=({\bf q},\Omega_{m})$ and
$\Pi_{q}=2g^{2}T\sum_{k}G_{k}G_{k+q}$. In the analysis above, we computed this
last term neglecting in $D_{q}$ the dynamical part of $\Pi_{q}$, which also
depends on momentum ${\bf q}$. Without this term, $D_{q}$ depends only on
frequency, and the momentum integration involves only $\Pi_{q}$. The double
integral over ${\bf q}$ and ${\bf k}$ can be transformed into the integration
over the two fermionic momenta ${\bf k}$ and ${\bf k}+{\bf q}$ and then into
the integration over the two dispersions $\epsilon_{k}$ and $\epsilon_{k+q}$.
Each integral is proportional to ${\text{sign}}\omega_{m}$, where $\omega_{m}$
is the Matsubara frequency in the corresponding Green’s function, hence the
momentum integration gives rise to the factor
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})$, where
$\omega_{m}-\omega_{m^{\prime}}=\Omega_{m}$. This is the same factor as in the
second term in Eq. (20) for the free energy $F_{\gamma}$ of the non-
regularized $\gamma$-model. The same holds for the interaction term in the
free energy in the $\gamma$ model:
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})$ in the interaction term in (20)
has been obtained by integrating independently over two fermionic dispersions:
one of $G_{k-q}$ in Eq. (6) and the other of $G_{k}$ in
$iT\sum_{k}\Sigma_{k}G_{k}$. Using now
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})=1+({\text{sign}}(\omega_{m}\omega_{m^{\prime}})-1)$,
we immediately see that the first and second terms correspond to contributions
from the static and dynamical parts of $\Pi_{q}$, respectively.

Hence, the static part of $\Pi_{q}$ accounts for the renormalization of the
bare $\omega_{D}$ into ${\bar{\omega}}_{D}$, which vanishes at the QCP. In the
underlying fermion-boson model, $\Pi_{q}$ is the full polarization operator,
with static and dynamics parts, and the renormalization
$\omega_{D}\to{\bar{\omega}}_{D}$ must be taken into consideration. This
implies that $(T/2)\sum_{q}\Pi_{q}D_{q}$ and $iT\sum_{k}\Sigma_{k}G_{k}$ have
to be computed without adding counter terms, and both depend on the upper
cutoff. Like we demonstrated, the two terms cancel out in the full free energy
$F$. The latter is expressed in terms of $D^{-1}_{q}$, which contains the
dressed ${\bar{\omega}}_{D}$ and the dynamical part of $\Pi_{q}$.

Then even when the renormalization of the bosonic mass does depend on the
cutoff (e.g., in the case of a lattice dispersion), the free energy is
expressed via the fully dressed mass, which vanishes at a QCP.

The $\gamma-$model is constructed differently. In this model, the
renormalization of $\omega_{D}$ into ${\bar{\omega}}_{D}$ is already absorbed
into $D_{loc}(q)$, which, by construction, depends on the dressed
${\bar{\omega}}_{D}$. Hence, the terms which renormalize $\omega_{D}$ must be
excluded to avoid double counting. The way to do this is to eliminate the
contribution from the static part of $\Pi_{q}$ by replacing
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})$ by
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})-1$. This is precisely the
counter term, which the authors of [55, 1, 56] suggested to add to regularize
the free energy of the $\gamma$-model.

The same reasoning holds for other values of $\gamma$. In each $\gamma-$model,
one has to subtract the renormalization of the bosonic mass to avoid double
counting. This is achieved by the same substitution
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})$ by
${\text{sign}}(\omega_{m}\omega_{m^{\prime}})-1$ in Eq. (20).

##  VI Extension to $\gamma<2$

It is instructive to verify how the $T-$independent and the $1/T$ term in Eq.
(70) evolve if we add a momentum-dependent term to the bosonic propagator
$D_{q}$ in (61) and gradually change the exponent $\gamma$ in the
corresponding $\gamma$-model to $\gamma<2$. A way to do this
phenomenologically is to consider a fermion-boson model with the bosonic
propagator

$D^{-1}_{q}=\Omega^{2}_{m}+(cq)^{2a}+{\bar{\omega}_{D}}^{2}+2{\bar{g}}^{2}\frac{|\Omega_{m}|}{v_{F}q}$
(72)

with $a>1$. We assume that the $q^{2a}$ term comes from fermions with energies
of order $E_{F}$, and set the prefactor $c$ to be of order
$E^{1/a}_{F}/k_{F}$. As before, we consider the double limit in which
${\bar{\omega}}_{D}$ tends to zero and simultaneously $E_{F}$ tends to
infinity.

We verified that the leading contribution to the fermionic self-energy
$\Sigma(\omega_{m})$ comes from the first two terms in (72), while the Landau
damping term accounts for a negative correction. Specifically,

$\Sigma(\omega_{m})\propto|\omega_{m}|^{1/a-1}\left(1-\left(\frac{T_{a}}{|\omega_{m}|}\right)^{\frac{a+1}{a}}\right)$
(73)

where

$T_{a}\sim{\bar{g}}\left(\frac{{\bar{g}}}{E_{F}}\right)^{\frac{a-1}{a+1}}$
(74)

For $E_{F}\to\infty$, $T_{a}$ tends to zero, hence $T_{a}/|\omega_{m}|$ is
vanishingly small for all $\omega_{m}$. Comparing (73) with
$\Sigma(\omega_{m})\propto|\omega_{m}|^{1-\gamma}$ in the $\gamma-$model, we
find $\gamma=2-1/a$. This exponent ranges between $1$ and $2$, when $a$ ranges
between $1$ and infinity. At $a=1+0$, a more accurate analysis shows that
$\Sigma(\omega_{m})\propto\log{|\omega_{m}|}$.

The free energy and the specific heat can be obtained in the same way as
above. For brevity, we skip the details of the calculations and just list the
results. We also neglect the free-fermion part of the specific heat and label
the specific heat due to fermion-boson interaction as $C_{int}(T)$. Up to a
positive overall factor,

$C_{int}(T)\propto
T^{\frac{2}{a}}\left(1-\left(\frac{T_{a}}{T}\right)^{\frac{a+1}{a}}+...\right)$
(75)

where dots stand for higher-order terms in the expansion in $T_{a}/T$. The
positive term in (75) comes from the $\Omega^{2}_{m}$ and $(cq)^{2a}$ terms in
the bosonic ptopagator, and the negative term comes from the Landau damping
term in (72). This negative term is vanishingly small as $T_{a}$ tends to zero
when $a>1$. The exponent $1/a$ equals to $2-\gamma$, hence $C_{int}(T)\propto
T^{2(2-\gamma)}$. For $\gamma=2$, $C_{int}(T)$ becomes temperature
independent. This is consistent with the result that we obtained in the
previous Section.

###  VI.1 Extension to $1/2<a<1$

For completeness, we also present the results for smaller values of the
exponent $a$: $1/2<a<1$. The condition $a>1/2$ is required for ultraviolet
convergence.

Evaluating the fermionic self-energy, we now obtain

$\Sigma(\omega_{m})\propto|\omega_{m}|^{\frac{2}{2a+1}}\left(1-\left(\frac{|\omega_{m}|}{T_{a}}\right)^{\frac{a+1}{a}}\right)$
(76)

where $T_{a}$ is the same as in (74). The dominant contribution to the self-
energy now comes from the Landau damping term and from the $(cq)^{2a}$ term in
$D_{q}$ in (72), while the $\Omega^{2}_{m}$ term accounts for a negative
correction. Because $T_{a}$ now tends to infinity at $E_{F}\to 0$, the second
term in (56) is vanishingly small for all $\omega_{m}$. Associating the
exponent $2/(2a+1)$ with $1-\gamma$, we find that for $a<1$,
$\gamma=(2a-1)/(2a+1)$.

For the specific heat we find

$C_{int}(T)\propto
T^{\frac{2}{2a+1}}\left(1-\left(\frac{T}{T_{a}}\right)^{\frac{a+1}{a}}+...\right)$
(77)

The positive contribution to $C(T)$ now comes from the Landau damping term and
the $(cq)^{2a}$ term in (72), while the negative contribution comes from the
$\Omega^{2}_{m}$ term. The dots stand for terms with higher powers of
$T/T_{a}$. Because for $a<1$, $T_{a}$ tends to infinity at $E_{F}\to\infty$,
the negative term is vanishingly small at any $T$. As a result $C_{int}(T)$ is
again positive. Using the relation $\gamma=(2a-1)/(2a+1)$, valid for $a<1$, we
find that $C_{int}(T)\propto T^{1-\gamma}$. This agrees with the results in
Sec. (III). At $a\to 1/2$, a more accurate analysis yields $C_{int}(T)\propto
T\log T$, as in Sec. (IV).

There is a discontinuity in $\gamma$ at $a=1$, i.e., the model with $a=1+0$
corresponds to $\gamma=1$, and the one with $a=1-0$ corresponds to
$\gamma=1/3$. This is the consequence of discontinuity of $T_{a}$ at $a=1$ and
$E_{F}\to\infty$: $T_{a}$ tends to zero at $a>1$, is of order ${\bar{g}}$ at
$a=1$, and tends to infinity at $a<1$. Right at $a=1$, the frequency
dependence of the self-energy and the temperature dependence of the specific
heat undergo a crossover from $\Sigma(\omega_{m})\propto|\omega_{m}|^{2/3}$
and $C_{int}(T)\sim T^{2/3}$ at $\omega_{m},T\ll{\bar{g}}$ to
$\Sigma(\omega_{m})\propto\log{|\omega_{m}|}$ and $C_{int}(T)\sim T^{2}$ at
$\omega_{m},T\gg{\bar{g}}$ (modulo logarithms). In both cases the specific
heat is positive. The low-temperature behavior of the model with $a=1$ is the
same as in the $\gamma=1/3$ model.

##  VII Conclusions

In this paper, we analyzed the free energy and specific heat for a system of
fermions interacting with nearly gapless bosons near a QCP in a metal. The
effective low-energy model for quantum-critical fermions is the one in which
bosons are integrated out, and the fermions are interacting via an effective,
purely dynamical interaction $V(\Omega_{m})\propto 1/|\Omega_{m}|^{\gamma}$.
This $\gamma-$model is adequate for the description of non-FL behavior and
pairing near a QCP, and the competition between tendencies towards non-FL and
pairing. This physics is fully determined by low-energy fermions and is
independent on the upper energy cutoff in the theory, $\Lambda$. The
condensation energy, associated with pairing, is also independent of
$\Lambda$. At the same time within the $\gamma$-model, the free energy in the
normal state does depend on $\Lambda$. Furthermore, the dependence on
$\Lambda$ extends to temperature-dependent terms in the free energy. As a
result, the specific heat in the $\gamma-$model also depends on the cutoff. In
recent papers [1, 56], the authors argued that the dependence on $\Lambda$ is
a spurious one and has to be eliminated by proper regularization. They added a
term to the free energy, which cancels out cutoff dependence of the free
energy. However, the regularized specific heat turns out to be negative for
$\gamma\geq 2$, considered in [1, 56]. Some of us and others [53] argued that
the specific heat becomes negative at small enough $T$ already at $\gamma>1$.

We analyzed the specific heat near the QCP by returning back to the underlying
fermion-boson model and collecting contributions to the free energy from
fermions, bosons, and their interaction. This allowed us to obtain the full
expression for the specific heat and compare it with the regularized specific
heat in fermions-only $\gamma-$model.

Our key result is that the specific heat in the full fermion-boson model is
independent on the cutoff and is positive all the way up to a QCP. This holds
within the Eliashberg theory, which we used in the calculations, in the
parameter range where the theory is rigorously justified, i.e., the
dimensionless Eliashberg parameter $\lambda_{E}$, which measures the strength
of vertex corrections, is small (number-wise, $C(T)$ remains positive even
when $\lambda_{E}=O(1)$).

We considered three cases, all in 2D: Ising-nematic QCP, antiferromagnetic
QCP, and QCP for electrons interacting with Einstein phononons. For the first
case, the exponent in the purely electronic model is $\gamma=1/3$. For the
second it is $\gamma=1/2$ for fermions near the hot spots, but is reduced to
$\gamma=0+$ in the effective model with the interaction averaged over the
Fermi surface. For electron-phonon case, the effective fermion-only model has
$\gamma=2$.

For the two cases with $\gamma<1$, where the specific heat in the regularized
$\gamma-$model is positive, we found that the regularization and the effect of
keeping the bosonic piece in the free energy is largely the same thing.
Specifically, the regularized specific heat has correct temperature dependence
($T^{2/3}$ for the Isng-nematic case, and $T\log T$ for the AFM case), and the
prefactor differs from the correct one only by a numerical factor, which,
moreover, is equal to one in the AFM case.

In the electron-phonon case ($\gamma=2$) the modified electronic specific heat
reproduces the temperature dependence of the actual $C(T)$. However, $C(T)$
has an additional temperature-independent piece, which also comes from the
electron-phonon interaction. Both terms originate from $\log{D_{q}}$, and the
temperature-independent term is the leading one. Because of this, the actual
$C(T)$ is positive for all values of the dressed Debye frequency, i.e., at any
distance from the QCP. The same holds for other models, whose fermionic part
is described by the $\gamma$ model with $\gamma>1$. We believe that this
result implies that the normal state of a critical fermion-boson model remains
stable at all $T$, as long as one neglects the pairing instability. In this,
our conclusions differ from the ones in Refs. [1, 56].

##  VIII Acknowledgement

##### Acknowledgment

We thank Ar. Abanov, B. Altshuler, A. Klein, A. Levchenko, D. Maslov, J.
Schmalian, G. Torroba, Y. Wang, Y. Wu, and E. Yuzbashyan for fruitful
discussions. This project was supported by the US-Israel Binational Science
Foundation (BSF). The work by A.V.C. was supported by U.S. Department of
Energy, Office of Science, Basic Energy Sciences, under Award No. DE-
SC0014402. E.B. was supported by the European Research Council (ERC) under
grant HQMAT (Grant Agreement No. 817799).

##  Appendix A Details about Ising-nematic case

###  A.1 Self-energy of an electron at a finite $T$

The one-loop self-energy of an electron is given by

$i\Sigma_{k}=-g^{*}T\sum_{\Omega_{m}}\int\frac{d^{2}\bm{q}}{(2\pi)^{2}}\frac{1}{i\tilde{\Sigma}_{k+q}-\epsilon_{\bm{k}+\bm{q}}}\frac{D_{0}}{q^{2}+m^{2}+\alpha\frac{\rvert\Omega_{m}\rvert}{\rvert\bm{q}\rvert}},$
(78)

where $\Sigma_{k}=\Sigma({\bf k},\omega_{m})$ and the notations are the same
as in the main text: $\tilde{\Sigma}_{k}\equiv\omega_{m}+\Sigma_{k}$ and
$\alpha=g^{*}k_{F}/(\pi v_{F}^{2})$, where $g^{*}$ is the effective fermion-
boson coupling.

At $T=0$, the sum is replaced by $T\sum_{\Omega_{n}}=(1/2\pi)\int
d\Omega_{n}$. The leading term in $\Sigma_{k}$ is obtained by factorizing the
momentum integration along and transverse to the Fermi surface (see Fig. 3 of
the main text). This leading term depends only on frequency, i.e., the self-
energy is local. At a QCP,

$\Sigma_{k}=\frac{3}{2}\bar{g}^{1/3}\rvert\omega_{m}\rvert^{2/3}\text{sgn}(\omega_{m}),$
(79)

where $\bar{g}$ \- the coupling constant of the corresponding fermionic
$\gamma=1/3$ model is

$\bar{g}=\frac{1}{3^{9/2}}\left(\frac{v_{F}}{k_{F}}\right)^{3}\alpha^{2}=\frac{1}{3^{9/2}\pi^{2}}\frac{(g^{*})^{2}}{E_{F}}.$
(80)

The factorization of momentum integration is valid as long as typical
fermionic momenta
$q^{\text{typ}}_{f}\sim\text{max}(\tilde{\Sigma}(\omega_{m}),\epsilon_{\bm{k}})/v_{F}$
(same as typical momenta transverse to the Fermi surface
$q^{\text{typ}}_{\perp}$) is much smaller than typical bosonic momentum
$q^{\text{typ}}_{b}\sim(\alpha|\Omega_{m}|)^{1/3}$ (same as typical momenta
along the Fermi surface $q^{\text{typ}}_{\parallel}$). The comparison of the
two scales shows that the factorization is valid in the whole range where
$\Sigma_{k}>\omega_{m}$ and at larger frequencies holds up to
$\omega_{max}\sim(g^{*}E_{F})^{1/2}\sim{\bar{g}}^{1/4}E_{F}^{3/4}$.

At a finite temperature, there are two types of bosonic fluctuations: the
thermal one with $\Omega_{m}=0$ and the quantum one with $\Omega_{m}\neq 0$.
This splits the self-energy into two parts [5, 30, *metzner_new, 39, 42,
*torroba_2, 37]

$\Sigma_{k}(k)=\Sigma^{th}_{k}+\Sigma^{q}_{k},$ (81)

where, we remind, $\Sigma_{k}=\Sigma({\bf k},\omega_{m})$. We have

$i\Sigma^{th}({\bf
k},\omega_{m})=-g^{*}T\int\frac{d^{2}\bm{q}}{(2\pi)^{2}}\frac{1}{i\tilde{\Sigma}(\bm{k}+\bm{q},\omega_{m})-\epsilon_{\bm{k}+\bm{q}}}\frac{1}{|{\bm{q}}|^{2}+m^{2}},$
(82)

and

$i\Sigma^{q}({\bf k},\omega_{m})=-g^{*}T\sum_{\Omega_{n}\neq
0}\int\frac{d^{2}\bm{q}}{(2\pi)^{2}}\frac{1}{i\tilde{\Sigma}({\bm{k}}+{\bm{q}},\omega_{m}+\Omega_{m})-\epsilon_{\bm{k}+\bm{q}}}\frac{1}{|{\bm{q}}|^{2}+m^{2}+\alpha\frac{\rvert\Omega_{n}\rvert}{\rvert\bm{q}\rvert}}.$
(83)

Below we consider the two components of the self-energy separately. We assume
for simplicity that the bosonic mass $m$ acquires some weak temperature
dependence via mode-mode coupling and cut $\log{m}$ singularity in the
formulas below by $\log{T}$ (for the analysis of $\Sigma_{k}$ for
$T-$independent mass see [37]. In this approximation, the quantum self-energy
$\Sigma^{q}$ can still be computed by factorizing the momentum integration and
remain local [37, 72]. The result is

$\displaystyle\Sigma^{q}(\omega_{m})$ $\displaystyle=\pi T\sum_{\Omega_{n}\neq
0}\left(\frac{\bar{g}}{\rvert\Omega_{n}\rvert}\right)^{1/3}\text{sgn}(\omega_{m}+\Omega_{n})$
$\displaystyle=\bar{g}^{1/3}(2\pi T)^{2/3}H_{1/3}(m),$ (84)

where $H_{\gamma}(m)=\sum_{n=1}^{m}1/n^{\gamma}$ is the Harmonic number. At
frequencies $\omega_{m}\gg T$, one can use the expansion of a Harmonic number
at large $m$: $H_{1/3}(m)\simeq 3/2(m+1/2)^{2/3}+\zeta(1/3)+...$ and obtain

$\Sigma^{q}(\omega_{m})\simeq\frac{3}{2}\bar{g}^{1/3}\text{sgn}(\omega_{m})\left(\rvert\omega_{m}\rvert^{2/3}+{2\over
3}\zeta(1/3)(2\pi T)^{2/3}+...\right),$ (85)

This formula is valid up to the same $\omega_{\text{max}}$ as at $T=0$.

For thermal self-energy, momentum integration can be factorized only in a
particular parameter range, which we identify below. Outside this range, the
leading contribution to $\Sigma^{th}_{k}$ in (82) is obtained by integrating
over both momentum components in the bosonic propagator.

Below we consider separately parameter ranges where $\Sigma^{th}_{k}$ is local
and where it is not.

###  A.2 Local self-energy: $\Sigma_{k}\equiv\Sigma(\omega_{m})$

In this section, we consider the situation when the momentum integration in
Eq. (82) can be factorized. The factorization implies that for the same
frequency, typical fermionic momentum (the one transverse to the Fermi
surface) is much smaller than typical bosonic momentum connecting points on
the Fermi surface. Typical fermionic momentum is
$q^{\text{typ}}_{f}\sim\text{max}(\tilde{\Sigma}(\omega_{m}),\epsilon_{\bm{k}})/v_{F}$,
while typical bosonic momentum is $q^{\text{typ}}_{b}\sim m$. Factorization is
justified when $q^{\text{typ}}_{f}\ll q^{\text{typ}}_{b}$. Under this
condition

$i\Sigma^{th}(\omega_{m})=-\frac{g^{*}T}{4\pi^{2}}\int_{-\Lambda_{q}}^{\Lambda_{q}}\frac{dq_{\bot}}{i\tilde{\Sigma}(\omega_{n})-\epsilon_{\bm{k}}-v_{F}q_{\bot}}\int_{-\Lambda_{q}}^{\Lambda_{q}}\frac{dq_{\parallel}}{q_{\parallel}^{2}+m^{2}}.$
(86)

where $\Lambda_{q}\sim k_{F}$ is the upper cutoff of momentum integration.
Assuming both $q_{f}$ and $q_{b}$ are far smaller than $\Lambda_{q}$, one can
set $\Lambda_{q}\to\infty$. Momentum integration then can be done explicitly,
and the result is

$\Sigma^{th}(\omega_{m})=\frac{g^{*}T}{4mv_{F}}\text{sgn}(\omega_{m})\equiv\pi
T\left(\frac{\bar{g}}{M}\right)^{1/3}\text{sgn}(\omega_{m}).$ (87)

where

$\displaystyle M$ $\displaystyle=\frac{64}{3^{9/2}}\frac{m^{3}}{\alpha}.$ (88)

The total self-energy $\Sigma(k)=\Sigma^{th}(k)+\Sigma^{q}(k)$ is

$\Sigma(\omega_{m})\simeq\left[\pi
T\left(\frac{\bar{g}}{M}\right)^{1/3}+\frac{3}{2}\bar{g}^{1/3}\rvert\omega_{m}\rvert^{2/3}\right]\text{sgn}(\omega_{m}),$
(89)

The two terms become comparable at

$\omega_{{\text{c}ross}}(T)\sim\frac{T^{3/2}}{M^{1/2}}.$ (90)

Thermal self-energy is larger at $\omega_{m}<\omega_{{\text{c}ross}}(T)$.

Eq. (89) is valid when $q^{\text{typ}}_{f}\ll q^{\text{typ}}_{b}$. i.e., when

$\tilde{\Sigma}(\omega_{m})/v_{F}\ll m.$ (91)

At $\omega_{m}<\omega_{{\text{c}ross}}$, $\Sigma^{th}>\Sigma^{q}$, and Eq.
(91) sets the condition on temperature

$M<T<T^{*}\sim\bar{g}\left(\frac{M}{\bar{g}}\right)^{2/3}\left(\frac{E_{F}}{\bar{g}}\right)^{1/2}.$
(92)

At $\omega_{m}>\omega_{{\text{c}ross}}$, $\Sigma^{th}<\Sigma^{q}$, and Eq.
(91) sets the condition on frequency

$\omega_{{\text{c}ross}}<\omega_{m}<\omega^{*}\sim\bar{g}\left(\frac{M}{\bar{g}}\right)^{1/2}\left(\frac{E_{F}}{\bar{g}}\right)^{3/4}.$
(93)

One can check that self-consistency condition
$\omega^{*}>\omega_{{\text{c}ross}}$ leads to the same condition on $T$ as Eq.
(92). Then, when Eq. (92) is satisfied, factorization of momentum integration
is valid for all frequencies up to $\omega^{*}$. We illustrate this in Fig. 9.

We emphasize that the $T$ range in Eq. (92) does exists at small but finite
$M$ simply because $M^{2/3}>M$, but collapses at a QCP, where $M=0$. In other
words, factorization of momentum integration in the integral for $\Sigma^{th}$
holds only away from a QCP.

There is one more condition. We assumed above that $\Sigma_{k}\gg\omega_{m}$.
A simple analysis shows that this condition is satisfied at arbitrary ratio of
$\Sigma^{th}$ and $\Sigma^{q}$ when $M<\bar{g}^{5/2}/E_{F}^{3/2}$. This
relation obviously holds for small $M$.

Figure 9: Parameter range where the self-energy given by Eq. (89) is local,
i.e., momentum-independent (marked by “Local” in the plot).

###  A.3 Non-local self-energy

At temperatures above $T^{*}$, the condition $q_{\parallel}^{typ}\gg
q_{\bot}^{typ}$ in the integral for $\Sigma^{th}$ is not satisfied. The
integration over $q_{\parallel}$ in (82) can be done explicitly:

$\int_{-\Lambda_{q}}^{\Lambda_{q}}\frac{dq_{\parallel}}{q_{\parallel}^{2}+q_{\bot}^{2}+m^{2}}\approx\frac{\pi}{\sqrt{q_{\bot}^{2}+m^{2}}}.$
(94)

Then

$i\Sigma^{th}(k)=-\frac{g^{2}AT}{4\pi}\int_{-\Lambda_{q}}^{\Lambda_{q}}\frac{dq_{\bot}}{i\tilde{\Sigma}(\omega_{m},\bm{k}+\bm{q})-\epsilon_{\bm{k}}-v_{F}q_{\bot}}\frac{1}{\sqrt{q_{\bot}^{2}+m^{2}}}.$
(95)

One can verify (see below) that at $T\gg T^{*}$, the leading term in this
integral is obtained by ignoring the $\bm{q}$ dependence in the ferminonic
propagator and pulling it out of the integral, i.e., by approximating

$\displaystyle\int_{-\Lambda}^{\Lambda}\frac{dq_{\bot}}{i\tilde{\Sigma}(k+q)-\epsilon_{\bm{k}}-v_{F}q_{\bot}}\frac{1}{\sqrt{q_{\bot}^{2}+m^{2}}}$
$\displaystyle\simeq\frac{2\log\left(\frac{1}{m}\right)}{i\tilde{\Sigma}(k)-\epsilon_{\bm{k}}}.$
(96)

This leads to an algebraic relation

$\Sigma^{th}(k)=\frac{B}{\left(\Sigma^{th}(k)+\tilde{\Sigma}^{q}(k)\right)+i\epsilon_{\bm{k}}},$
(97)

where $\tilde{\Sigma}^{q}(k)=\Sigma^{q}(k)+\omega_{m}$, and

$B=\frac{g^{*}T}{2\pi}\log\left(\frac{1}{m}\right).$ (98)

Eq. (97), viewed as quadratic equation on $\Sigma^{th}(k)$, has two solutions.
The physical one must satify the boundary condition $\Sigma^{th}=0$ at $B=0$.
This selects out the solution

$\Sigma^{th}(k)=-\frac{\tilde{\Sigma}^{q}(\omega_{n})+i\epsilon_{\bm{k}}}{2}+\text{sgn}(\omega_{n})\sqrt{\frac{\left(\tilde{\Sigma}^{q}(\omega_{n})+i\epsilon_{\bm{k}}\right)^{2}}{4}+B}.$
(99)

We remind that we define $\sqrt{z}$ with a branch cut along the negative real
axis of the complex variable $z$. One can verify that upon
$\omega_{m}\leftrightarrow-\omega_{m}$ and
$\epsilon_{\bm{k}}\leftrightarrow-\epsilon_{\bm{k}}$, $\Sigma^{th}(k)$
transforms as

$\displaystyle\text{Re}\Sigma^{th}(\omega_{m},\epsilon_{\bm{k}})$
$\displaystyle=-\text{Re}\Sigma^{th}(-\omega_{m},\epsilon_{\bm{k}})=+\text{Re}\Sigma^{th}(\omega_{m},-\epsilon_{\bm{k}}),$
(100) $\displaystyle\text{Im}\Sigma^{th}(\omega_{m},\epsilon_{\bm{k}})$
$\displaystyle=+\text{Im}\Sigma^{th}(-\omega_{m},\epsilon_{\bm{k}})=-\text{Im}\Sigma^{th}(\omega_{m},-\epsilon_{\bm{k}}).$
(101)

When $\Sigma^{th}(k)>\tilde{\Sigma}^{q}$,
$\Sigma^{th}(k)\approx\sqrt{B}\text{sgn}(\omega_{n})$.

Eq. (99) has been obtained in [37] for $\epsilon_{k}=0$. We will be chiefly
interested in the consequences of the dependence of $\Sigma^{th}(k)$ on
$\epsilon_{k}$.

Figure 10: Parameter range where the self-energy $\Sigma^{th}(k)$ is non-local
(marked by “Non-local” in the plot).

The total self-energy is given by $\Sigma(k)=\Sigma^{th}(k)+\Sigma^{q}(k)$,
with the quantum part given by Eq. (85). Expanding the self-energy to linear
order in $\epsilon_{\bm{k}}$ we find

$\Sigma(\omega_{m},\epsilon_{\bm{k}})\simeq\Sigma(\omega_{m},0)-\frac{i}{2}\left(1-\frac{\rvert\Sigma^{q}(\omega_{m})\rvert}{\sqrt{[\Sigma^{q}(\omega_{m})]^{2}+4B}}\right)\epsilon_{\bm{k}},$
(102)

where

$\Sigma(\omega_{m},0)=\frac{\tilde{\Sigma}^{q}(\omega_{n})}{2}+\text{sgn}(\tilde{\Sigma}^{q}(\omega_{n}))\sqrt{\frac{1}{4}[\tilde{\Sigma}^{q}(\omega_{n})]^{2}+B}.$
(103)

The first term renormalizes the frequency dependence of the Green’s function,
while the second term renormalizes the Fermi velocity into

$v_{F}^{*}=\frac{1}{2}\left(1+\frac{\rvert\Sigma^{q}(\omega_{m})\rvert}{\sqrt{[\Sigma^{q}(\omega_{m})]^{2}+4B}}\right)v_{F}.$
(104)

The renormalized velocity becomes $v_{F}/2$ when
$\Sigma^{th}>{\tilde{\Sigma}}^{q}$, i.e., when
$2\sqrt{B}\gg\Sigma^{q}(\omega_{m})$, and differs only slightly from $v_{F}$
when $\Sigma^{th}<{\tilde{\Sigma}}^{q}$. The crossover between the two regimes
is at frequency

${\tilde{\omega}}_{\text{cross}}\sim\frac{B^{3/4}}{\bar{g}^{1/2}}=\bar{g}\left(\frac{T}{\bar{g}}\right)^{3/4}\left(\frac{E_{F}}{\bar{g}}\right)^{3/8}\log^{3/4}\left(\frac{1}{m}\right).$
(105)

We next consider the applicability range for Eq. (101). Let’s set
$\epsilon_{k}=0$ to avoid unnecessary complications. In obtaining (101) we
assumed that

$\tilde{\Sigma}_{k}/v_{F}>m$ (106)

At $\omega_{m}\ll{\tilde{\omega}}_{\text{cross}}$,
$|\Sigma_{k}|\approx|\Sigma^{th}|=\sqrt{B}$, and the inequality in (106) sets
the condition on $T$:

$T>\bar{g}\left(\frac{M}{\bar{g}}\right)^{2/3}\left(\frac{E_{F}}{\bar{g}}\right)^{1/2}\frac{1}{\log(1/m)}\sim\frac{T^{*}}{\log{1/m}}$
(107)

Up to a logarithm, this is $T>T^{*}$, i.e., $T=T^{*}$ is a sharp boundary
between local and non-local forms of $\Sigma^{th}$. Keeping the logarithm one
obtains [37] an extended crossover regime. It formally becomes wide at $m\to
0$, but like we said, we assume that mode-mode coupling cuts $\log{m}$ at
$\log{T}$. Then the crossover regime is rather narrow. The upper limit on $T$,
at which

$T<T_{max}\sim\bar{g}\left(\frac{E_{F}}{\bar{g}}\right)^{1/2}$ (108)

is set by the boundary condition on the momentum-independence of the thermal
self-energy.

At $\omega_{m}\ll{\tilde{\omega}}_{\text{cross}}$, Eq. (101) is valid in the
same range of $T$, up to a frequency $\omega\sim{\bar{g}}$. We illustrate this
in Fig. 10.

##  Appendix B Cancellation of $\Sigma^{th}(k)$ in the free energy

In this Appendix, we show explicitly that the thermal self-energy
$\Sigma^{th}$ cancels out in the free energy $F_{el}$, Eq. (3). This holds
when $\Sigma^{th}$ is local, and when it is non-local and given by (101).

###  B.1 Case of local self-energy: $\Sigma(k)=\Sigma(\omega_{n})$

When the self-energy is independent to $\epsilon_{\bm{k}}$, the momentum
integration is straightforward,

$\sum_{\bm{k}}\log\left(\frac{i\tilde{\Sigma}(\omega)-\epsilon_{\bm{k}}}{\epsilon_{\bm{k}}}\right)=N_{F}\left(\pi\rvert\tilde{\Sigma}(\omega_{n})\rvert+i\pi\Lambda_{q}\text{sgn}(\omega_{n})\right),$
(109)
$\sum_{\bm{k}}\frac{\Sigma(\omega)}{\tilde{\Sigma}(\omega)+i\epsilon_{\bm{k}}}=N_{F}\pi\rvert\Sigma(\omega_{n})\rvert.$
(110)

Upon summation over $\omega_{n}$ we obtain

$F_{el}=-2T\sum_{\omega_{n}}\pi
N_{F}\left(\rvert\tilde{\Sigma}(\omega_{n})\rvert-\rvert\Sigma(\omega_{n})\rvert\right)\equiv-2\pi
TN_{F}\sum_{\omega_{n}}\rvert\omega_{n}\rvert,$ (111)

which is equal to the free energy of the non-interacting Fermi gas. The self-
energy cancels out from this expression.

###  B.2 Case of non-local $\Sigma^{th}(k)$

We now show that the cancellation holds even when $\Sigma^{th}$ depends on the
dispersion $\epsilon_{\bm{k}}$.

The electronic part of the free energy per volume is

$F_{el}=-T\sum_{k}\ln\left(\frac{(i\tilde{\Sigma}(k)-\epsilon_{\bm{k}})(i\tilde{\Sigma}(-k)-\epsilon_{-\bm{k}})}{\epsilon_{\bm{k}}^{2}}\right)+2T\sum_{k}\frac{\tilde{\Sigma}(k)-\omega_{m}}{\tilde{\Sigma}(k)+i\epsilon_{\bm{k}}}.$
(112)

We show below that the non-local $\Sigma^{th}(k)$ actually cancels out in each
of two contributions to $F_{el}$.

Substituting $\Sigma^{th}$ from (99) into the first term, we obtain after
simple algebra

$\displaystyle(i\tilde{\Sigma}(k)-\epsilon_{\bm{k}})(i\tilde{\Sigma}(-k)-\epsilon_{-\bm{k}})$
$\displaystyle=\left(\frac{\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert\pm
i\epsilon_{\bm{k}}}{2}+\sqrt{\frac{\left(\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert\pm
i\epsilon_{\bm{k}}\right)^{2}}{4}+B}\right)$
$\displaystyle\left(\frac{\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert\pm
i\epsilon_{\bm{k}}}{2}+\sqrt{\frac{\left(\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert\pm
i\epsilon_{\bm{k}}\right)^{2}}{4}+B}\right),$ (113)

where $\pm$ refers to $\text{sgn}(\omega_{m})$. Introducing
$\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert=2\sqrt{B}y$ and
$\epsilon_{\bm{k}}=2\sqrt{B}z$, we re-express (113) as

$\displaystyle\frac{(i\tilde{\Sigma}(k)-\epsilon_{\bm{k}})(i\tilde{\Sigma}(-k)-\epsilon_{-\bm{k}})}{\epsilon_{\bm{k}}^{2}}$
$\displaystyle=\frac{1}{4z^{2}}\left(y\pm iz+\sqrt{\frac{\left(y\pm
iz\right)^{2}}{4}+B}\right)^{2},$ (114)

To obtain the first term in $F_{el}$, we need to integrate this expression
over $\epsilon_{k}$ (i.e., over $z$) and sum up over Matsubara frequencies.
Combining contribution from positive and negative $z$, we obtain

$\displaystyle F_{el}^{(1)}=$
$\displaystyle-T\sum_{k}\ln\left(\frac{(i\tilde{\Sigma}(k)-\epsilon_{\bm{k}})(i\tilde{\Sigma}(-k)-\epsilon_{-\bm{k}})}{\epsilon_{\bm{k}}^{2}}\right)$
$\displaystyle=$
$\displaystyle-4\sqrt{B}N_{F}T\sum_{\omega_{m}}\int_{0}^{\Lambda}dz\ln\left[\frac{1}{4z^{2}}\left(y+iz+\sqrt{\frac{\left(y+iz\right)^{2}}{4}+B}\right)\left(y-iz+\sqrt{\frac{\left(y-iz\right)^{2}}{4}+B}\right)\right]$
$\displaystyle=$ $\displaystyle-4\pi\sqrt{B}N_{F}T\sum_{\omega_{m}}y=-2\pi
N_{F}T\sum_{\omega_{m}}\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert.$ (115)

We see that the result is the same as if $\Sigma^{th}$ was absent.

For the second term in $F_{el}$, we again use Eq. (99) and express
$\Sigma^{th}$ in terms of ${\tilde{\Sigma}}^{q}$ and $\epsilon_{k}$. This
yields

$\frac{\tilde{\Sigma}(k)-\omega_{m}}{\tilde{\Sigma}(k)+i\epsilon_{\bm{k}}}=\frac{\frac{\tilde{\Sigma}^{q}(k)-i\epsilon_{\bm{k}}}{2}+\sqrt{\frac{\left(\tilde{\Sigma}^{q}(k)+i\epsilon_{\bm{k}}\right)^{2}}{4}+B}\text{sgn}(\omega_{m})-\omega_{m}}{\frac{\tilde{\Sigma}^{q}(k)-i\epsilon_{\bm{k}}}{2}+\sqrt{\frac{\left(\tilde{\Sigma}^{q}(k)+i\epsilon_{\bm{k}}\right)^{2}}{4}+B}\text{sgn}(\omega_{m})+i\epsilon_{\bm{k}}}.$
(116)

Re-expressing in terms of $y$ and $z$, as did before, we obtain

$\displaystyle\sum_{\bm{k}}\frac{\tilde{\Sigma}(k)-\omega_{m}}{\tilde{\Sigma}(k)+i\epsilon_{\bm{k}}}$
$\displaystyle=$ $\displaystyle
2\sqrt{B}N_{F}\int_{-\Lambda_{q}/(2\sqrt{B})}^{\Lambda_{q}/(2\sqrt{B})}dz\frac{y-iz+\sqrt{\left(y+iz\right)^{2}+1}}{y+iz+\sqrt{\left(y+iz\right)^{2}+1}}$
$\displaystyle-2\rvert\omega_{m}\rvert
N_{F}\int_{-\Lambda_{q}/(2\sqrt{B})}^{\Lambda_{q}/(2\sqrt{B})}dz\frac{1}{y+iz+\sqrt{\left(y+iz\right)^{2}+1}}.$
(117)

This integral is convergent with typical $z=O(1)$. Given that
$\Lambda\gg\sqrt{B}$, the $z-$integration can be extended to infinite limits.
Integrating in infinite limits, we obtain

$\displaystyle\int_{-\infty}^{\infty}dz\frac{y-iz+\sqrt{\left(y+iz\right)^{2}+1}}{y+iz+\sqrt{\left(y+iz\right)^{2}+1}}$
$\displaystyle=y\int_{-\infty}^{\infty}dt\frac{1-it+\sqrt{\left(1+it\right)^{2}+y^{-2}}}{1+it+\sqrt{\left(1+it\right)^{2}+y^{-2}}}$
$\displaystyle=\pi y.$ (118)
$\displaystyle\int_{-\infty}^{\infty}dz\frac{1}{y+iz+\sqrt{\left(y+iz\right)^{2}+1}}$
$\displaystyle=\int_{-\infty}^{\infty}dt\frac{1}{1+it+\sqrt{\left(1+it\right)^{2}+1}}$
(119) $\displaystyle=\pi/2.$ (120)

Collecting contributions, we find

$\displaystyle
F_{el}^{(2)}=2T\sum_{k}\frac{\tilde{\Sigma}(k)-\omega_{m}}{\tilde{\Sigma}(k)+i\epsilon_{\bm{k}}}$
$\displaystyle=2\pi
TN_{F}\left(\rvert\tilde{\Sigma}^{q}(\omega_{m})\rvert-\rvert\omega_{m}\rvert\right)$
(121)

as if $\Sigma^{th}$ was absent. Combining $F_{el}^{(1)}$ and $F_{el}^{(2)}$,
we obtain

$F_{el}=-2\pi TN_{F}\sum_{\omega_{m}}\rvert\omega_{m}\rvert,$ (122)

which is the free energy of a non-interacting Fermi gas. We see that the self-
energy cancels out in $F_{el}$ even when $\Sigma^{th}$ depends on
$\epsilon_{\bm{k}}$.

##  Appendix C Evaluation of free energy

###  C.1 $\gamma$-model at $\gamma=0^{+}$

In the purely electronic $\gamma$-model, the free energy is
$F_{\gamma}=F_{el}+F_{int}$. For a generic non-zero $\gamma$, the free energy
has been evaluated in Ref. [53]. Here, we compute the free energy for the
special case $\gamma\to 0^{+}$, relevant to the analysis of the
antiferromagnetic QCP (see the main text). The case $\gamma\to 0^{+}$ requires
special care as the interaction
$V(\Omega_{m})\propto\log{{\bar{g}}/|\Omega_{m}|}$ For the free energy, we
have in this case $F_{0+}=F_{free}+F_{0+,int}$, where in the notations from
the main text

$F_{0+,int}=\frac{3g^{*}}{8\pi^{3}v_{F}^{2}\beta}S_{0+},$ (123)

$\beta=2v_{x}v_{y}/v_{F}^{2}$, and

$S_{0+}=(2\pi
T)^{2}\sum_{n,n^{\prime}=-M_{f}}^{M_{f}-1}\text{sgn}(2n+1)\text{sgn}(2n^{\prime}+1)\log\frac{\rvert
n-n^{\prime}\rvert 2\pi T}{T_{0}^{**}}.$ (124)

The thermal contribution, from $n=n^{\prime}$, has to be evaluated at a non-
zero bosonic mass. This contribution to $F_{0+}$ is linear in $T$ and does not
affect the specific heat. Summing over $n^{\prime}\neq n$, we obtain

$\displaystyle S_{0+}$ $\displaystyle=4(2\pi
T)^{2}\left(2\sum_{n=1}^{M_{f}-1}\log(n!)-\frac{1}{2}\sum_{n=1}^{2M_{f}-1}\log(n!)\right)-4\pi
T\Lambda\log{\frac{2\pi T}{T_{0}^{**}}}$ (125)

Contributions from $n\sim O(1)$ are of order $\sim T^{2}$ We show that the
summation over $n\gg 1$ yields a larger $\sim T^{2}\log(T)$ term. To evaluate
this contribution, we use the asymptotic formula

$\log(n!)=\left(n+\frac{1}{2}\right)\log(n)-n+\frac{1}{2}\log(2\pi)+\frac{1}{12n}+O(\frac{1}{n^{2}}),$
(126)

Substituting into (125) and using

$\displaystyle\sum_{n=1}^{M_{f}-1}\left(n+\frac{1}{2}\right)\log(n)={1\over
2}M_{f}^{2}\log M_{f}-{1\over 4}M_{f}^{2}-{1\over 2}M_{f}+O(1)$
$\displaystyle\sum_{n=1}^{M_{f}-1}\frac{1}{n}=\log\left(M_{f}\right)+O(1),$
$\displaystyle\sum_{n=1}^{M_{f}-1}n=M_{f}^{2}/2-M_{f}/2,~{}~{}~{}\sum_{n=1}^{M_{f}-1}1=M_{f}-1$
(127)

and the relation between $M_{f}$ and the upper theory cutoff $\Lambda$, we
obtain

$\displaystyle S_{0+}$ $\displaystyle=$ $\displaystyle 4(2\pi
T)^{2}\left[-M_{f}^{2}\log 2+{1\over 2}\log(2\pi)M_{f}-{1\over 8}\log
M_{f}+O(1)\right]-2(2\pi T)^{2}M_{f}\log\left({T\over T_{0}^{**}}\right)$
(128) $\displaystyle=$ $\displaystyle-\Lambda^{2}\log(16)+4\pi\Lambda
T\log\left(\frac{T_{0}^{**}}{2\pi
T}\right)-2\pi^{2}T^{2}\log\left(\frac{\Lambda}{T}\right)+O(T^{2}).$

Hence

$\displaystyle F_{0+}^{int}$ $\displaystyle=N_{F}\frac{3g^{*}}{2\pi^{2}\beta
E_{F}}\left[-\Lambda^{2}\log(2)\right.$ $\displaystyle\left.+\pi\Lambda
T\log\left(\frac{T_{0}^{**}}{T}\right)\right.$ $\displaystyle\left.-{1\over
2}\pi^{2}T^{2}\log\left(\frac{\Lambda}{T}\right)+O(T^{2})\right].$ (129)

Differentiating twice with respect to temperature, one obtains the specific
heat

$C_{0+}^{int}(T)=N_{F}\frac{3g^{*}}{2\pi\beta E_{F}}\left[\Lambda+\pi
T\log\left(\frac{\Lambda}{T}\right)+O(T)\right].$ (130)

It contains a constant $\propto\Lambda$ and a universal $T\log(1/T)$ term.

For comparison, we evaluate the free energy of the regularized $\gamma$-model,
$\bar{F}_{0+}=F_{frre}+\bar{F}_{0+}^{int}$, where

$\bar{F}_{0+}^{int}=\frac{3g^{*}}{8\pi^{3}v_{F}^{2}\beta}\bar{S}_{0+},$ (131)

and

$\bar{S}_{0+}=(2\pi
T)^{2}\sum_{n,n^{\prime}=-M_{f}}^{M_{f}-1}\left(\text{sgn}(2n+1)\text{sgn}(2n^{\prime}+1)-1\right)\log\frac{\rvert
n-n^{\prime}\rvert 2\pi T}{T_{0}^{**}}.$ (132)

Since the summand is non-zero only when $2n+1$ and $2n^{\prime}+1$ has
opposite signs, the thermal part with $n=n^{\prime}$ is avoided. The sum is
evaluated in the same way as for the original $\gamma$-model, and the result
is gives rise to

$\displaystyle\bar{S}_{0+}$ $\displaystyle=4(2\pi
T)^{2}\left(2\sum_{n=0}^{M_{f}-1}\log(n!)-\sum_{n=0}^{2M_{f}-1}\log(n!)\right)-4(2\pi
T)^{2}M_{f}^{2}\log\frac{2\pi T}{T_{0}^{**}}$ $\displaystyle=4(2\pi
T)^{2}\left[-M_{f}^{2}\log M_{f}+\log\left({e^{3/2}\over
4}\right)M_{f}^{2}-{1\over 12}\log M_{f}+O(1)\right]$ $\displaystyle-4(2\pi
T)^{2}M_{f}^{2}\log\frac{2\pi T}{T_{0}^{**}}$
$\displaystyle=4\Lambda^{2}\log\frac{e^{3/2}T_{0}^{**}}{4\Lambda}-\frac{4}{3}\pi^{2}T^{2}\log\left(\frac{\Lambda^{2}}{2\pi
TT_{0}^{**}}\right)+O\left(T^{2}\right).$ (133)

As expected, the cutoff-dependent $\Lambda T\log(1/T)$ term is removed. The
coefficient of the universal $T^{2}\log(1/T)$ term is $2/3$ of that in the
original $\gamma$-model. This is the same ratio as for a non-zero $\gamma$
(see the main text). The interaction part of the free energy is

$\displaystyle\bar{F}_{0+}^{int}$
$\displaystyle=-N_{F}\frac{3g^{*}}{2\pi^{2}\beta
E_{F}}\left[\Lambda^{2}\log\left(\frac{4\Lambda}{e^{3/2}T_{0}^{**}}\right)\right.$
$\displaystyle\left.+{1\over 3}\pi^{2}T^{2}\log\left(\frac{\Lambda^{2}}{2\pi
TT_{0}^{**}}\right)+O\left(T^{2}\right)\right].$ (134)

Differentiating twice with respect to temperature, we obtain the specific heat

$\bar{C}_{0+}^{int}(T)=N_{F}\frac{g^{*}}{\beta
E_{F}}T\log\left(\frac{\Lambda^{2}}{2\pi TT_{0}^{**}}\right)+O\left(T\right).$
(135)

###  C.2 Boson-fermion model

The free energy of the underlying boson-fermion model is given by
$F=F_{free}+F_{bos}$, where

$F_{bos}=\frac{k}{2}T\sum_{q}\log{\left(-D_{q}^{-1}\right)}$ (136)

and $k$ is the number of components of the bosonic fields: $k=1$ for Ising-
nematic and electron-phonon cases, and $k=3$ for an antiferromagnetic QCP. We
presented the results for $F_{bos}$ for the three cases in the main text. Here
we show the details of the evaluation of $F_{bos}^{*}$.

####  C.2.1 Ising-nematic QCP

Subtracting frequency-independent term from $\log\left(-D_{q}^{-1}\right)$ and
integrating over the momentum in Eq. (136) we obtain

$\displaystyle F_{bos}$
$\displaystyle=\frac{T}{2}\sum_{\Omega_{n}}\int\frac{d^{2}\bm{q}}{4\pi^{2}}\log\left(1+\frac{\alpha\rvert\Omega_{n}\rvert}{q^{3}}\right)=\frac{\alpha^{2/3}}{4\sqrt{3}}T\sum_{\Omega_{n}}\rvert\Omega_{n}\rvert^{2/3}.$
(137)

The frequency sum over $2M_{b}+1$ Matsubara frequencies is expressed via the
Harmonic number $\sum_{n=1}^{M_{b}}n^{2/3}=H_{-2/3}(M_{b})$. Then
$F_{bos}=\alpha^{2/3}(2\pi T)^{5/3}H_{-2/3}(M_{b})/4\sqrt{3}\pi$.

Using the expansion of Harmonic number at large argument,
$H_{-2/3}(M_{b})=(3/5)(M_{b}+1/2)^{5/3}+\zeta(-2/3)+O(1/(M_{b}+1/2)^{1/3})$,
and using the relation between $M_{b}$ and $\Lambda$, Eq. (32), we obtain

$F_{bos}={\alpha^{2/3}\over 4\sqrt{3}\pi}\left({3\over
5}\Lambda^{5/3}+\zeta(-{2\over 3})(2\pi T)^{5/3}\right).$ (138)

Differentiating twice over temperature and combining with free-fermion
contribution, we obtain $C_{I-N}(T)$, given by Eq. (35).

####  C.2.2 Antiferromagnetic QCP

For this case, the momentum integral in Eq. (136) is logarithmically singular
and depends on the upper momentum cutoff $\Lambda_{q}\sim k_{F}$. Integrating
over $q$, we obtain

$\displaystyle F_{bos}$ $\displaystyle=\frac{3\alpha
T}{8\pi}\sum_{\Omega_{n}}\rvert\Omega_{n}\rvert\log\frac{\Lambda^{2}_{q}}{\alpha\rvert\Omega_{n}\rvert}\equiv-\frac{3\alpha}{2}T^{2}\sum_{n=1}^{M_{b}}n\log\frac{nT}{T_{0}},$
(139)

where $T_{0}\sim\Lambda^{2}_{q}/\alpha$. The frequency sum over $2M_{b}+1$
Matsubara frequencies is expressed in terms of the hyperfactorial function
$H(x)$ as

$\sum_{n=1}^{M_{b}}n\log\frac{nT}{T_{0}}=\log\left[H(M_{b})\right]+\frac{M_{b}\left(M_{b}+1\right)}{2}\log\frac{T}{T_{0}}.$
(140)

At large $M_{b}\gg 1$, $\log\left(H(M_{b})\right)$ is expanded as

$\displaystyle\log\left(H(M_{b})\right)$
$\displaystyle=-\frac{1}{4}M_{b}^{2}+\left(\frac{1}{12}+\frac{1}{2}M_{b}(M_{b}+1)\right)\log(M_{b})$
$\displaystyle+{\cal O}(1).$ (141)

Using the relation between $M_{b}$ and $\Lambda$, Eq. (32), we obtain after
simple algebra

$\displaystyle(2\pi T)^{2}\sum_{n=1}^{M_{b}}n\log\frac{nT}{T_{0}}$
$\displaystyle=-\frac{1}{4}\Lambda^{2}+\frac{1}{2}\Lambda^{2}\log\frac{\Lambda}{2\pi
T_{0}}+\frac{1}{3}\pi^{2}T^{2}\log\frac{T_{0}}{T}+{\cal O}(T^{2}).$ (142)

Hence

$\displaystyle F_{bos}$
$\displaystyle=-\frac{3\alpha}{16\pi^{2}}\Lambda^{2}\log\frac{\Lambda}{2\pi
T_{0}\sqrt{e}}+\frac{\alpha}{8}T^{2}\log\frac{T}{T_{0}}+{\cal O}(T^{2}).$
(143)

####  C.2.3 QCP of an Einstein phonon

Near a QCP at which the dressed Debye frequency vanishes for $q<2k-F$, the
dressed phonon propagator takes the form
$D_{q}^{-1}=\Omega_{n}^{2}+\bar{\omega}_{D}^{2}+2\bar{g}^{2}\rvert\Omega_{n}\rvert/(v_{F}q)(2k_{F}/\sqrt{4k^{2}_{F}-q^{2}})$,
where $\omega_{D}$ and $\bar{\omega}_{D}=\omega_{D}(1-2\lambda)^{1/2}$ are
bare and dressed Debye frequencies, and $\lambda=\bar{g}^{2}/\omega_{D}^{2}$.
Substituting into (136) and treating the Landau damping term as perturbation,
we obtain

$\displaystyle F_{bos}$
$\displaystyle\simeq\frac{T}{2}\sum_{\Omega_{n}}\int\frac{d^{2}\bm{q}}{4\pi^{2}}\log\left(\Omega_{n}^{2}+\bar{\omega}_{D}^{2}\right)+\frac{T}{2}\sum_{\Omega_{n}}\int\frac{d^{2}\bm{q}}{4\pi^{2}}\frac{2\bar{g}^{2}}{v_{F}q}\frac{\rvert\Omega_{n}\rvert}{\Omega_{n}^{2}+\bar{\omega}_{D}^{2}}{2k_{F}\over\sqrt{4k_{F}^{2}-q^{2}}}.$
(144)

where the integration over $q$ is up to $2k_{F}$. The first term is the free
energy of a free Einstein phonon with the dressed Debye frequency
$\bar{\omega}_{D}$:

$F_{bos}^{(1)}=4N_{F}E_{F}T\left(\log{\bar{\omega}_{D}}+2\sum_{n=1}^{M_{B}}\log({2\pi
Tn})+\sum_{n=1}^{M_{B}}\log\left(1+\frac{\bar{\omega}_{D}^{2}}{4\pi^{2}T^{2}n^{2}}\right)\right)$
(145)

Using

$\displaystyle\sum_{1}^{M_{B}}\log{n}=(M_{b}+1/2)\log{(M_{b}+1/2)/e}+\frac{1}{2}\log{2\pi}$
$\displaystyle\sum_{1}^{M_{B}}\log{2\pi T}=(M_{b}+1/2)\log{2\pi
T}-\frac{1}{2}\log{2\pi T}$ (146)

and the relation between $M_{B}$ and $\Lambda$, we obtain

$F_{bos}^{(1)}=4N_{F}E_{F}\left[\frac{\Lambda}{\pi}\log{\Lambda}{e}+\sum_{1}^{M_{B}}\log\left(1+\frac{\bar{\omega}_{D}^{2}}{4\pi^{2}T^{2}n^{2}}\right)-log{T}\right]$
(147)

The first term is $T$-independent and does not contribute to entropy and
specific heat. In the second term, the sum over $m$ converges and the
summation can be extended to $M_{b}=\infty$. Evaluating the sum using Euler-
Maclauren formula and combining with the last term, we obtain

$F_{bos}^{(1)}=4N_{F}E_{F}\left[{\Lambda\over\pi}\log\left({\Lambda\over
e}\right)+T\log\left(1-e^{-\bar{\omega}_{D}/T}\right)\right].$ (148)

We note in passing that the exponential temperature dependence of
$F_{bos}^{(1)}$ at the smallest $T$ implies that all terms in Euler-Maclauren
series expansion in $T/{\bar{\omega}}_{D}$ vanish, as we explicitly verified.

Carrying out the momentum integration in the second term in (144), we obtain

$\displaystyle F_{bos}^{(2)}$
$\displaystyle=\pi\bar{g}^{2}N_{F}T\sum_{\Omega_{n}}\frac{\rvert\Omega_{n}\rvert}{\Omega_{n}^{2}+\bar{\omega}_{D}^{2}}=\bar{g}^{2}N_{F}\sum_{n=1}^{M_{b}}\frac{n}{n^{2}+\left(\frac{{\bar{\omega}}_{D}}{2\pi
T}\right)^{2}},$ (149)

The sum over Matsubara frequencies is expressed via di-Gamma functions as

$\sum_{n=1}^{M_{b}}\frac{n}{n^{2}+\left(\frac{{\bar{\omega}}_{D}}{2\pi
T}\right)^{2}}=\text{Re}\left[\psi\left(1+i\frac{{\bar{\omega}}_{D}}{2\pi
T}+M_{b}\right)-\psi\left(1+i\frac{{\bar{\omega}}_{D}}{2\pi T}\right)\right].$
(150)

Using the asymptotic expression $\psi(z)\simeq\log(z)$ at $\rvert z\rvert\gg
1$ and re-expressing $\log{\Lambda/(2\pi T)}$ as
$\log{\Lambda/{\bar{\omega}}_{D}}+\log{{\bar{\omega}}_{D}/(2\pi T)}$ we obtain

$\displaystyle F_{bos}$ $\displaystyle=$ $\displaystyle
N_{F}\left[4E_{F}{\Lambda\over\pi}\log\left({\Lambda\over
e}\right)+\bar{g}^{2}\log\left({\Lambda\over\bar{\omega}_{D}}\right)\right]$
(151) $\displaystyle+$ $\displaystyle
4N_{F}E_{F}T\left[\log{\left(1-e^{-\bar{\omega}_{D}/T}\right)}+\lambda_{E}{\bar{\omega}_{D}\over
4T}f\left(\frac{\bar{\omega}_{D}}{2\pi T}\right)\right].$

where the dimensionless function $f(x)$ is

$f(x)=\log x-\frac{1}{2}\psi(1+ix)-\frac{1}{2}\psi(1-ix).$ (152)

This is Eq. (65) in the main text.

##  Appendix D Phenomenological models that map to the $\gamma$-model with
$0<\gamma<1$

In this Appendix, we consider a phenomenological extension of the Ising-
nematic model, which maps to the $\gamma$ model with $\gamma=1/3$, to a family
of boson-fermion models that map to the $\gamma$-model with $0<\gamma<1$. The
boson propagator takes the form

$D_{q}^{-1}=-\left(q^{2-a}+\frac{\alpha\rvert\Omega\rvert}{q}\right)/D_{0},$
(153)

where the parameter $a$ is tunable. We assume that the Fermi surface is
circular, like in the Ising-nematic case.

To establish the relation with the $\gamma$-model, we compute the free energy,
$F=F_{el}+F_{int}$. As in the Ising-nematic case, it can be re-expressed as
$F=F_{free}+F_{int}$, where $F_{free}$ is the contribution of free Fermi gas,
and $F_{int}$ comes from fermion-boson interaction

$F_{int}=-g^{2}T^{2}\sum_{m,m^{\prime}}\int\frac{d^{2}\bm{k}d^{2}\bm{k}^{\prime}}{(4\pi^{2})^{2}}\frac{1}{i\tilde{\Sigma}(\omega_{m})-\epsilon_{\bm{k}}}\frac{1}{i\tilde{\Sigma}(\omega_{m}^{\prime})-\epsilon_{\bm{k}^{\prime}}}D_{q},$
(154)

where $\bm{q}=\bm{k}-\bm{k}^{\prime}$ by momentum conservation. We assume and
then verify that typical momentum scale in the boson propagator,
$\sim\omega^{1/(3-a)}$, is much larger than the one in the fermion propagator,
$\sim\tilde{\Sigma}(\omega)/v_{F}$. In this situation, the momentum
integration can be factorized as

$\displaystyle F_{int}$
$\displaystyle=g^{*}T^{2}\sum_{m,m^{\prime}}\int\frac{dk_{\bot}}{2\pi}\frac{1}{i\tilde{\Sigma}(\omega_{m})-v_{F}k_{\bot}}\int\frac{dk_{\bot}^{\prime}}{2\pi}\frac{1}{i\tilde{\Sigma}(\omega_{m}^{\prime})-v_{F}k_{\bot}^{\prime}}$
$\displaystyle\int\frac{dq_{\parallel}}{2\pi}\frac{1}{|q_{\parallel}|^{2-a}+\frac{\alpha\rvert\omega_{m}-\omega_{m}^{\prime}\rvert}{|q_{\parallel}|}}.$
(155)

Carrying out the momentum integration, we obtain

$\displaystyle F_{int}$
$\displaystyle=-\pi^{2}T^{2}N_{F}\bar{g}^{\frac{1-a}{3-a}}\sum_{m,m^{\prime}}\sum_{\bm{k}\bm{k}^{\prime}}\frac{\text{sgn}(\omega_{m}\omega_{m^{\prime}})}{\rvert\omega_{m}-\omega_{m^{\prime}}\rvert^{\frac{1-a}{3-a}}}.$
(156)

This is equivalent to the free energy of the $\gamma$-model with
$\gamma=(1-a)/(3-a)$ and the effective coupling constant

$\bar{g}=\left(\frac{1}{(3-a)\sin\frac{2\pi}{3-a}}\frac{g^{*}}{2\pi
v_{F}\alpha^{\frac{1-a}{3-a}}}\right)^{\frac{3-a}{1-a}}.$ (157)

The effective $\gamma$ changes continuously from $0$ to $1$ when $a$ is
changes between $1$ to $-\infty$. For all these $a$, the coupling constant
$\bar{g}^{\gamma}$ remains positive-defined. The sum in Eq. (156) has been
evaluated in the main text. It contains $\Lambda-$dependent terms and the
universal term of order $T^{(5-a)/(3-a)}$. In the regularized $\gamma$-model,
$\Lambda-$dependent terms cancel out. The free energy is

$\displaystyle\bar{F}_{\gamma}$
$\displaystyle=F_{free}+\frac{2}{4(3-a)\sin\frac{2\pi}{3-a}}\zeta\left(-\frac{2}{3-a}\right)\left(2\pi\alpha\right)^{\frac{2}{3-a}}T^{\frac{5-a}{3-a}}.$
(158)

The full free energy of the model includes the contribution from bosons.

$F_{full}=F_{full}(T=0)-{\pi^{2}\over
3}N_{F}T^{2}+\frac{1}{4\sin\frac{2\pi}{3-a}}\zeta\left(-\frac{2}{3-a}\right)\left(2\pi\alpha\right)^{\frac{2}{3-a}}T^{\frac{5-a}{3-a}},$
(159)

where $F_{full}(T=0)$ comes from the zero-temperature quantum fluctuations and
depends on cutoff $\Lambda$. Comparing the $T$-dependent terms in $F_{full}$
and $\bar{F}_{\gamma}$, we see that they have the same form, but the
prefactors for the $T^{(5-a)/(3-a)}$ term differ by $2/(3-a)$. The prefactors
agree at $a=1+0$, when $\gamma=0+$, as we also found in the explicit analysis
of the $\gamma=0+$ model in the main text.

## References

  * Yuzbashyan and Altshuler [2022a] E. A. Yuzbashyan and B. L. Altshuler, Breakdown of the migdal-eliashberg theory and a theory of lattice-fermionic superfluidity, Phys. Rev. B 106, 054518 (2022a).
  * Nayak and Wilczek [1994] C. Nayak and F. Wilczek, Non-fermi liquid fixed point in 2 + 1 dimensions, Nuclear Physics B 417, 359 (1994).
  * Sachdev  _et al._ [1995] S. Sachdev, A. V. Chubukov, and A. Sokol, Crossover and scaling in a nearly antiferromagnetic fermi liquid in two dimensions, Phys. Rev. B 51, 14874 (1995).
  * Millis [1992] A. J. Millis, Nearly antiferromagnetic fermi liquids: An analytic eliashberg approach, Phys. Rev. B 45, 13047 (1992).
  * Abanov  _et al._ [2003] A. Abanov, A. V. Chubukov, and J. Schmalian, Quantum-critical theory of the spin-fermion model and its application to cuprates: Normal state analysis, Advances in Physics 52, 119 (2003).
  * Abanov  _et al._ [2001a] A. Abanov, A. V. Chubukov, and J. Schmalian, Fingerprints of spin mediated pairing in cuprates, Journal of Electron spectroscopy and related phenomena 117, 129 (2001a).
  * Scalapino [2012a] D. J. Scalapino, A common thread: The pairing interaction for unconventional superconductors, Rev. Mod. Phys. 84, 1383 (2012a).
  * Bergeron  _et al._ [2012] D. Bergeron, D. Chowdhury, M. Punk, S. Sachdev, and A.-M. S. Tremblay, Breakdown of fermi liquid behavior at the $(\pi,\pi)=2{k}_{F}$ spin-density wave quantum-critical point: The case of electron-doped cuprates, Phys. Rev. B 86, 155123 (2012).
  * Efetov [2015] K. B. Efetov, Quantum criticality in two dimensions and marginal fermi liquid, Phys. Rev. B 91, 045110 (2015).
  * Tsvelik [2017] A. M. Tsvelik, Ladder physics in the spin fermion model, Phys. Rev. B 95, 201112 (2017).
  * Marsiglio [2020] F. Marsiglio, Eliashberg theory: A short review, Annals of Physics 417, 168102 (2020).
  * Chubukov  _et al._ [2020a] A. V. Chubukov, A. Abanov, Y. Wang, and Y.-M. Wu, The interplay between superconductivity and non-fermi liquid at a quantum-critical point in a metal, Annals of Physics 417, 168142 (2020a).
  * Chubukov  _et al._ [2020b] A. V. Chubukov, A. Abanov, I. Esterlis, and S. A. Kivelson, Eliashberg theory of phonon-mediated superconductivity – when it is valid and how it breaks down, Annals of Physics 417, 168190 (2020b).
  * Scalapino [2012b] D. J. Scalapino, A common thread: The pairing interaction for unconventional superconductors, Rev. Mod. Phys. 84, 1383 (2012b).
  * Altshuler  _et al._ [1994] B. L. Altshuler, L. B. Ioffe, and A. J. Millis, Low-energy properties of fermions with singular interactions, Phys. Rev. B 50, 14048 (1994).
  * Abanov  _et al._ [2001b] A. Abanov, A. V. Chubukov, and A. M. Finkel’stein, Coherent vs . incoherent pairing in 2d systems near magnetic instability, EPL (Europhysics Letters) 54, 488 (2001b).
  * Oganesyan  _et al._ [2001] V. Oganesyan, S. A. Kivelson, and E. Fradkin, Quantum theory of a nematic fermi fluid, Phys. Rev. B 64, 195109 (2001).
  * Metzner  _et al._ [2003] W. Metzner, D. Rohe, and S. Andergassen, Soft fermi surfaces and breakdown of fermi-liquid behavior, Phys. Rev. Lett. 91, 066402 (2003).
  * Rech  _et al._ [2006] J. Rech, C. Pépin, and A. V. Chubukov, Quantum critical behavior in itinerant electron systems: Eliashberg theory and instability of a ferromagnetic quantum critical point, Phys. Rev. B 74, 195126 (2006).
  * Chubukov  _et al._ [2004] A. V. Chubukov, C. Pépin, and J. Rech, Instability of the quantum-critical point of itinerant ferromagnets, Phys. Rev. Lett. 92, 147003 (2004).
  * Efetov  _et al._ [2013] K. B. Efetov, H. Meier, and C. Pepin, Pseudogap state near a quantum critical point, Nature Physics 9, 442 (2013).
  * Raghu  _et al._ [2015] S. Raghu, G. Torroba, and H. Wang, Metallic quantum critical points with finite bcs couplings, Phys. Rev. B 92, 205104 (2015).
  * Wang  _et al._ [2017] H. Wang, S. Raghu, and G. Torroba, Non-fermi-liquid superconductivity: Eliashberg approach versus the renormalization group, Phys. Rev. B 95, 165137 (2017).
  * Wang  _et al._ [2018] H. Wang, Y. Wang, and G. Torroba, Superconductivity versus quantum criticality: Effects of thermal fluctuations, Phys. Rev. B 97, 054502 (2018).
  * Fitzpatrick  _et al._ [2015] A. L. Fitzpatrick, S. Kachru, J. Kaplan, S. Raghu, G. Torroba, and H. Wang, Enhanced pairing of quantum critical metals near $d=3+1$, Phys. Rev. B 92, 045118 (2015).
  * Lee [2009] S.-S. Lee, Low-energy effective theory of fermi surface coupled with u(1) gauge field in $2+1$ dimensions, Phys. Rev. B 80, 165102 (2009).
  * Dalidovich and Lee [2013] D. Dalidovich and S.-S. Lee, Perturbative non-fermi liquids from dimensional regularization, Phys. Rev. B 88, 245106 (2013).
  * Schlief  _et al._ [2017] A. Schlief, P. Lunts, and S.-S. Lee, Exact critical exponents for the antiferromagnetic quantum critical metal in two dimensions, Phys. Rev. X 7, 021010 (2017).
  * Bonesteel  _et al._ [1996] N. E. Bonesteel, I. A. McDonald, and C. Nayak, Gauge fields and pairing in double-layer composite fermion metals, Phys. Rev. Lett. 77, 3009 (1996).
  * Dell’Anna and Metzner [2006] L. Dell’Anna and W. Metzner, Fermi surface fluctuations and single electron excitations near pomeranchuk instability in two dimensions, Phys. Rev. B 73, 045127 (2006).
  * Yamase and Metzner [2012] H. Yamase and W. Metzner, Fermi-surface truncation from thermal nematic fluctuations, Phys. Rev. Lett. 108, 186405 (2012).
  * Son [1999] D. T. Son, Superconductivity by long-range color magnetic interaction in high-density quark matter, Phys. Rev. D 59, 094019 (1999).
  * Chubukov and Schmalian [2005] A. V. Chubukov and J. Schmalian, Superconductivity due to massless boson exchange in the strong-coupling limit, Phys. Rev. B 72, 174520 (2005).
  * Metlitski and Sachdev [2010a] M. A. Metlitski and S. Sachdev, Quantum phase transitions of metals in two spatial dimensions. ii. spin density wave order, Physical Review B 82, 075128 (2010a).
  * Mross  _et al._ [2010] D. F. Mross, J. McGreevy, H. Liu, and T. Senthil, Controlled expansion for certain non-fermi-liquid metals, Phys. Rev. B 82, 045121 (2010).
  * Chubukov and Wölfle [2014] A. V. Chubukov and P. Wölfle, Quasiparticle interaction function in a two-dimensional fermi liquid near an antiferromagnetic critical point, Phys. Rev. B 89, 045108 (2014).
  * Klein  _et al._ [2020] A. Klein, A. V. Chubukov, Y. Schattner, and E. Berg, Normal state properties of quantum critical metals at finite temperature, Phys. Rev. X 10, 031053 (2020).
  * Metlitski  _et al._ [2015] M. A. Metlitski, D. F. Mross, S. Sachdev, and T. Senthil, Cooper pairing in non-fermi liquids, Phys. Rev. B 91, 115111 (2015).
  * Punk [2016] M. Punk, Finite-temperature scaling close to ising-nematic quantum critical points in two-dimensional metals, Phys. Rev. B 94, 195113 (2016).
  * Maslov and Chubukov [2010] D. L. Maslov and A. V. Chubukov, Fermi liquid near pomeranchuk quantum criticality, Phys. Rev. B 81, 045110 (2010).
  * Varma [2020] C. M. Varma, Colloquium: Linear in temperature resistivity and associated mysteries including high temperature superconductivity, Rev. Mod. Phys. 92, 031001 (2020).
  * Damia  _et al._ [2020] J. A. Damia, M. Solís, and G. Torroba, How non-fermi liquids cure their infrared divergences, Phys. Rev. B 102, 045147 (2020).
  * Wang and Torroba [2017] H. Wang and G. Torroba, Non-fermi liquids at finite temperature: Normal-state and infrared singularities, Phys. Rev. B 96, 144508 (2017).
  * Wang  _et al._ [2016] Y. Wang, A. Abanov, B. L. Altshuler, E. A. Yuzbashyan, and A. V. Chubukov, Superconductivity near a quantum-critical point: The special role of the first matsubara frequency, Phys. Rev. Lett. 117, 157001 (2016).
  * Lederer  _et al._ [2017] S. Lederer, Y. Schattner, E. Berg, and S. A. Kivelson, Superconductivity and non-fermi liquid behavior near a nematic quantum critical point, Proceedings of the National Academy of Sciences 114, 4905 (2017).
  * Abanov and Chubukov [2020] A. Abanov and A. V. Chubukov, Interplay between superconductivity and non-fermi liquid at a quantum critical point in a metal. i. the $\gamma$ model and its phase diagram at $T=0$: The case $0<\gamma<1$, Phys. Rev. B 102, 024524 (2020).
  * Wu  _et al._ [2020a] Y.-M. Wu, A. Abanov, Y. Wang, and A. V. Chubukov, Interplay between superconductivity and non-fermi liquid at a quantum critical point in a metal. ii. the $\gamma$ model at a finite $T$ for $0<\gamma<1$, Phys. Rev. B 102, 024525 (2020a).
  * Wu  _et al._ [2020b] Y.-M. Wu, A. Abanov, and A. V. Chubukov, Interplay between superconductivity and non-fermi liquid behavior at a quantum critical point in a metal. iii. the $\gamma$ model and its phase diagram across $\gamma=1$, Phys. Rev. B 102, 094516 (2020b).
  * Wu  _et al._ [2021a] Y.-M. Wu, S.-S. Zhang, A. Abanov, and A. V. Chubukov, Interplay between superconductivity and non-fermi liquid at a quantum critical point in a metal. iv. the $\gamma$ model and its phase diagram at $1<\gamma<2$, Phys. Rev. B 103, 024522 (2021a).
  * Wu  _et al._ [2021b] Y.-M. Wu, S.-S. Zhang, A. Abanov, and A. V. Chubukov, Interplay between superconductivity and non-fermi liquid behavior at a quantum-critical point in a metal. v. the $\gamma$ model and its phase diagram: The case $\gamma=2$, Phys. Rev. B 103, 184508 (2021b).
  * Zhang  _et al._ [2021] S.-S. Zhang, Y.-M. Wu, A. Abanov, and A. V. Chubukov, Interplay between superconductivity and non-fermi liquid at a quantum critical point in a metal. vi. the $\gamma$ model and its phase diagram at $2\gamma<3$, Phys. Rev. B 104, 144509 (2021).
  * Wu  _et al._ [2022] Y.-M. Wu, S.-S. Zhang, A. Abanov, and A. V. Chubukov, Odd frequency pairing in a quantum critical metal, Phys. Rev. B 106, 094506 (2022).
  * Zhang  _et al._ [2022a] S.-S. Zhang, Y.-M. Wu, A. Abanov, and A. V. Chubukov, Superconductivity out of a non-fermi liquid: Free energy analysis, Physical Review B 106, 144513 (2022a).
  * Protter  _et al._ [2021] M. Protter, R. Boyack, and F. Marsiglio, Functional-integral approach to gaussian fluctuations in eliashberg theory, Phys. Rev. B 104, 014513 (2021).
  * Yuzbashyan and Altshuler [2022b] E. A. Yuzbashyan and B. L. Altshuler, Migdal-eliashberg theory as a classical spin chain, Phys. Rev. B 106, 014512 (2022b).
  * Yuzbashyan  _et al._ [2022] E. A. Yuzbashyan, M. K.-H. Kiessling, and B. L. Altshuler, Superconductivity near a quantum critical point in the extreme retardation regime, Phys. Rev. B 106, 064502 (2022).
  * Grossman  _et al._ [2021] O. Grossman, J. S. Hofmann, T. Holder, and E. Berg, Specific heat of a quantum critical metal, Phys. Rev. Lett. 127, 017601 (2021).
  * Moon and Chubukov [2010] E.-G. Moon and A. Chubukov, Quantum-critical pairing with varying exponents, Journal of Low Temperature Physics 161, 263 (2010).
  * Note [1] The authors of [1, 56] found the prefactor to be negative for $\gamma\geq 2$, which they only considered. The authors of [53] argued that the prefactor is negative for $\gamma>1$.
  * Luttinger and Ward [1960] J. M. Luttinger and J. C. Ward, Ground-state energy of a many-fermion system. ii, Phys. Rev. 118, 1417 (1960).
  * Eliashberg [1960] G. M. Eliashberg, Interactions between electrons and lattice vibrations in a superconductor, JETP 11, 696 (1960), [ZhETF, 38, 966, (1960)].
  * Bardeen and Stephen [1964] J. Bardeen and M. Stephen, Free-energy difference between normal and superconducting states, Phys. Rev. 136, A1485 (1964).
  * Prange and Kadanoff [1964] R. E. Prange and L. P. Kadanoff, Transport theory for electron-phonon interactions in metals, Phys. Rev. 134, A566 (1964).
  * Haslinger and Chubukov [2003] R. Haslinger and A. V. Chubukov, Condensation energy in strongly coupled superconductors, Phys. Rev. B 68, 214508 (2003).
  * Chubukov  _et al._ [2005a] A. V. Chubukov, D. L. Maslov, S. Gangadharaiah, and L. I. Glazman, Thermodynamics of a fermi liquid beyond the low-energy limit, Phys. Rev. Lett. 95, 026402 (2005a).
  * Chubukov  _et al._ [2005b] A. V. Chubukov, D. L. Maslov, S. Gangadharaiah, and L. I. Glazman, Singular perturbation theory for interacting fermions in two dimensions, Phys. Rev. B 71, 205112 (2005b).
  * Secchi  _et al._ [2020] A. Secchi, M. Polini, and M. I. Katsnelson, Phonon-mediated superconductivity in strongly correlated electron systems: A luttinger–ward functional approach, Annals of Physics 417, 168100 (2020).
  * Benlagra  _et al._ [2011] A. Benlagra, K. Kim, and C. Pépin, The luttinger–ward functional approach in the eliashberg framework: a systematic derivation of scaling for thermodynamics near the quantum critical point, Journal of Physics: Condensed Matter 23, 145601 (2011).
  * Maslov and Chubukov [2009] D. L. Maslov and A. V. Chubukov, Nonanalytic paramagnetic response of itinerant fermions away and near a ferromagnetic quantum phase transition, Phys. Rev. B 79, 075112 (2009).
  * Chowdhury and Berg [2020] D. Chowdhury and E. Berg, Intrinsic superconducting instabilities of a solvable model for an incoherent metal, Phys. Rev. Research 2, 013301 (2020).
  * Esterlis  _et al._ [2021] I. Esterlis, H. Guo, A. A. Patel, and S. Sachdev, Large-$n$ theory of critical fermi surfaces, Phys. Rev. B 103, 235129 (2021).
  * Guo  _et al._ [2022] H. Guo, A. A. Patel, I. Esterlis, and S. Sachdev, Large-$n$ theory of critical fermi surfaces. ii. conductivity, Phys. Rev. B 106, 115151 (2022).
  * Note [2] We assume, as in previous works on metallic QCP, that fermionic bandwidth $W$ is the largest scale of the problem, and neglect terms, which are small in $g/W$.
  * Blaizot and Iancu [1997] J.-P. Blaizot and E. Iancu, Lifetimes of quasiparticles and collective excitations in hot qed plasmas, Phys. Rev. D 55, 973 (1997).
  * Metlitski and Sachdev [2010b] M. A. Metlitski and S. Sachdev, Quantum phase transitions of metals in two spatial dimensions. i. ising-nematic order, Physical Review B 82, 075127 (2010b).
  * Holder and Metzner [2015] T. Holder and W. Metzner, Fermion loops and improved power-counting in two-dimensional critical metals with singular forward scattering, Phys. Rev. B 92, 245128 (2015).
  * Eberlein  _et al._ [2016] A. Eberlein, I. Mandal, and S. Sachdev, Hyperscaling violation at the ising-nematic quantum critical point in two-dimensional metals, Phys. Rev. B 94, 045133 (2016).
  * Pimenov  _et al._ [2021] D. Pimenov, A. Kamenev, and A. V. Chubukov, One-dimensional scattering of two-dimensional fermions near quantum criticality, Phys. Rev. B 103, 214519 (2021).
  * Damia  _et al._ [2019] J. A. Damia, S. Kachru, S. Raghu, and G. Torroba, Two-dimensional non-fermi-liquid metals: A solvable large-$n$ limit, Phys. Rev. Lett. 123, 096402 (2019).
  * Note [3] The actual value of $\gamma$ is somewhat different from $1/2$ as corrections to fermion-boson vertex are logarithmical, and series of these corrections change the exponent $\gamma$ to $1/2+\epsilon$, where $\epsilon$ is positive, but small numerically [5, morr_97]. Besides, the dynamical exponent $z$ also flows exponentially from $z=2$ to a smaller value [34]. It has been argued [28] that in the absence of a superconducting instability, this flow eventually, at the lowest energies, brings the system into the basin of attraction of a stable fixed point with $z=1$. Our analysis is valid at energies where the dynamical exponent is still $z\approx 2$.
  * Wang and Chubukov [2013] Y. Wang and A. V. Chubukov, Superconductivity at the onset of spin-density-wave order in a metal, Phys. Rev. Lett. 110, 127001 (2013).
  * Fitzpatrick  _et al._ [2013] A. L. Fitzpatrick, S. Kachru, J. Kaplan, and S. Raghu, Non-fermi-liquid fixed point in a wilsonian theory of quantum critical metals, Phys. Rev. B 88, 125116 (2013).