# Third order corrections to the ground state energy of the gas of spin $s$
fermions with arbitrary densities of different spin projections

Piotr H. Chankowski, Jacek Wojtkiewicz and Szymon Augustynowicz111Emails:
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
Faculty of Physics, University of Warsaw,
Pasteura 5, 02-093 Warszawa, Poland

###### Abstract

Recently we have computed the third order corrections to the ground state
energy of the arbitrarily polarized diluted gas of spin 1/2 fermions
interacting through a spin-independent repulsive two-body potential. Here we
extend this result to the gas of spin $s$ fermions - a system the Hamiltonian
of which has an accidental $SU(2s+1)$ symmetry - with arbitrary densities of
fermions having different spin projections. The corrections are computed semi-
analytically using the effective field theory approach and are parametrized by
the $s$\- and $p$\- wave scattering lengths $a_{0}$ and $a_{1}$ and the
$s$-wave effective radius $r_{0}$, measurable in the low energy fermion-
fermion elastic scattering. The result is used to study the impact the higher
order corrections can have on the characteristics of the phase transition (at
zero temperature) to the ordered phase (on the emergence of the itinerant
ferromagnetism).

Keywords: Diluted gas of interacting fermions, effective field theory,
scattering lengths, itinerant ferromagnetism, phase transitions

##  1 Introduction

Although a clear experimental evidence is still lacking [1], there is a strong
conviction that in a finite density system of fermions a transition to the
ordered phase in which the densities of different spin components are not all
equal should be induced by a repulsive spin-independent interaction, if the
system is sufficiently dense and/or the interaction is strong enough. Such a
transition can be conveniently quantified by a nonzero value of an order
parameter $P$ which, in the usually considered case of spin $s=1/2$ fermions,
can be defined to have the meaning of the polarization (magnetization)
$P=(N_{\uparrow}-N_{\downarrow})/N$, $N=N_{\uparrow}+N_{\downarrow}$ of the
system. (An analogous order parameter can be defined in the general case of
spin $s>1/2$ if a specific pattern of the ordering is assumed.)

Theoretical investigation of this transition at zero temperature reduces to
computing the energy density of the system of $N$ fermions enclosed in the
volume $V$ as a function of the densities of different spin projections (or as
a function of the chosen order parameter) and taking the thermodynamic limit.
Such computations are most easily performed using the effective field theory
approach within which the underlying spatially non-local binary interaction of
fermions is replaced by a (in principle infinite) set of local operators of
decreasing length dimension and which yields the expansion of the computed
energy density in powers (in higher orders modified also by logarithms) of the
(overall) Fermi wave vector

$\displaystyle k_{\rm F}=\left({6\pi^{2}\over g_{s}}~{}\\!{N\over
V}\right)^{1/3},\phantom{aaa}g_{s}=2s+1~{}\\!,$ (1)

of the system. The expansion obtained in this way is naturally parametrized by
the scattering lengths $a_{\ell}$, effective radii $r_{\ell}$, etc.,
$\ell=0,1,\dots,$ which characterize the elastic scattering of low energy
fermions “in vacuum” and which are taken to specify the underlying interaction
potential. This approach, pioneered in [2] (see also [3, 4]), allowed to
complete recently [5] the order $k_{\rm F}^{4}$ (with respect to the energy of
the system of completely non-interacting fermions) corrections to the ground
state of the system of spin $s$ fermions with equal densities of different
spin projections. Using this approach it was also easy to reproduce (semi-
analytically) but in the universal setting the old (obtained by considering
the specific model of hard spheres interaction) result of Kanno [6] who
computed (fully analytically) the order $k_{\rm F}^{2}$ correction to the
ground state of the system of spin $1/2$ fermions [7] for an arbitrary value
of the system’s polarization $P$. The effective theory approach allows to
immediately extend [8] the result of Kanno to the case of spin $s$ fermions:
the complete up to the order $k_{\rm F}^{2}$ formula for the ground state
energy density $E_{\Omega}/V$ as a function of the densities of fermions
having different spin projections can be now written in the form222The order
$k_{\rm F}$ correction is the extension (obvious from the effective field
theory point of view) of the textbook [9, 10] mean-field result to the case of
spin $s$ fermions.

$\displaystyle{E_{\Omega}\over V}={k_{\rm F}^{3}\over 6\pi^{2}}~{}\\!{3\over
5}~{}\\!{\hbar^{2}k_{\rm F}^{2}\over
2m_{f}}\left\\{\sum_{\sigma=1}^{g_{s}}x_{\sigma}^{5}+{20\over
9\pi}~{}\\!k_{\rm
F}a_{0}\sum_{\sigma^{\prime}<\sigma}x_{\sigma^{\prime}}^{3}x_{\sigma}^{3}+{2\over\pi^{2}}~{}\\!(k_{\rm
F}a_{0})^{2}\sum_{\sigma^{\prime}<\sigma}J_{K}(x_{\sigma^{\prime}},x_{\sigma},)\right\\},$
(2)

in which

$\displaystyle x_{\sigma}={p_{{\rm F}\sigma}\over k_{\rm
F}}~{}\\!,\phantom{aaaa}p_{{\rm
F}\sigma}=\left({6\pi^{2}}~{}\\!{N_{\sigma}\over V}\right)^{1/3},$ (3)

and

$\displaystyle J_{K}(x,y)={1\over
21}\left\\{(22~{}\\!x^{3}y^{3}(x+y)-4~{}\\!x^{7}\ln\\!\left(1+{y\over
x}\right)-4~{}\\!y^{7}\ln\\!\left(1+{x\over y}\right)\right.$
$\displaystyle+{1\over
2}~{}\\!x~{}\\!y~{}\\!(x-y)^{2}(x+y)[15~{}\\!(x^{2}+y^{2})+11~{}\\!x~{}\\!y]\phantom{aaaaa}~{}$
$\displaystyle\left.-{7\over
4}(x-y)^{4}(x+y)((x+y)^{2}+x~{}\\!y)\ln\\!\left({x+y\over
x-y}\right)\right\\},\phantom{a}~{}$ (4)

is the function obtained by Kanno in his computation [6]. The extension to
$s>1/2$ is of obvious interest in view of the fact that the systems of
interacting fermions (and bosons) are nowadays investigated in experiments
with cold atoms in which the interaction strength can be tuned by exploiting
the physics of the so-called Feshbach resonance.

The system of spin of $1/2$ fermions with the (spin independent) repulsive
interaction constitutes the simple, textbook example [9, 10] of the continuous
phase transition because such a character of the paramagnetic to ferromagnetic
transition to the ordered state (signaled by the emergence of spontaneous
magnetization), both when it is induced by decreasing the temperature and by
increasing the density and/or the interaction strengths (that is when the
value of the parameter $k_{\rm F}a_{0}$ is increased) is found when the
computations are restricted to the mean field approximation (only the order
$k_{\rm F}a_{0}$ corrections to the thermodynamic functions are taken into
account). However when the order $(k_{\rm F}a_{0})^{2}$ corrections are
included this transition turns at low temperatures into the ordinary first
order one as has been found in [11] and at exactly $T=0$ can be easily
confirmed by using the formula (2) [12, 8]. Although in this approximation the
transition at zero temperature occurs at $k_{\rm F}a_{0}=1.054$ at which value
the perturbative expansion cannot most probably be trusted, such a character
of the transition seemed to be in line with the arguments [13] based on
general principles of the Ginzburg-Landau theory.

On the other hand in the series of papers [14] (overlooked in [7, 12]) the
ground state energy of the system of spin $1/2$ fermions as a function of the
polarization $P$ has been computed beyond a fixed order, by performing a
resummation of an infinite subclass of the effective field theory diagrams
contributing to it. As the authors of these papers claim, the transition to
the ordered state at $T=0$ becomes then again of the continuous type.
Moreover, the critical value of the parameter $k_{\rm F}a_{0}$ found by them
is in a surprisingly good agreement with the estimates based on the Quantum
Monte Carlo computations reported in [15]. While assessing the reliability of
the approximations made in [14] is in general difficult, it might be
instructive to see, using the complete order $k_{\rm F}^{3}$ corrections,
which we have recently computed in [16] in the case of spin $1/2$, how the
third order terms included in the resummation compare with the ones neglected
and how the character of the phase transition changes when it is analyzed
using the complete third order formula for the ground state energy. In view of
the mentioned circumstance that in experimental setups the role of fermions
are played by cold atomic gases which can have spins greater than $1/2$, it is
also of interest to investigate the transition in the case of higher spins.
Finally, according to the power counting rules which organize the expansion
within the effective field theory [17, 2], in the order $k_{\rm F}^{3}$ to the
energy contribute density for the first time interaction operators involving
in addition to four fermion field operators also two spatial derivatives which
bring in the dependence on the $p$-wave scattering length $a_{1}$ and the
$s$-wave effective radius $r_{0}$ which parametrize the departure of the
result from strict unversality. It is therefore natural to check their
potential impact on the phase transition.

Therefore the plan of the paper is as follows. In the next section, to prepare
the ground, we recall the effective theory approach allowing to easily recover
the order $k_{\rm F}a_{0}$ and $(k_{\rm F}a_{0})^{2}$ (this one semi-
analytically) corrections in (2). Then, in Section 3 and 4 we determine the
next order correction to the formula (2) adapting appropriately the
computation [16] of the third order correction to the ground state energy of
the arbitrarily polarized gas of spin $1/2$ fermions. The discussion of the
approach of [14] is included in Section 3. The character of the transition to
the ordered state is discussed in Section 5.

##  2 Second order corrections

As said, the effective theory naturally yields an expansion of the computed
quantities in powers of the product $kR$ of the wave vector $k$ corresponding
to the characteristic energy scale of the problem and of the characteristic
length scale $R$ of the “fundamental” interaction potential. In the case of
the corrections to the ground state energy of the gas of fermions the role of
the characteristic wave vector plays the Fermi wave vector $k_{\rm F}$ defined
in (3) and they are most conveniently computed (using the standard Feynman
rules formulated e.g. in [18]) as the sum of connected “vacuum” diagrams using
the Dyson expansion of the right hand side of the formula

$\displaystyle\lim_{T\rightarrow\infty}\exp(-iT(E_{\Omega}-E_{\Omega_{0}})/\hbar)=\lim_{T\rightarrow\infty}\langle\Omega_{0}|{\rm
T}\exp\\!\left(-{i\over\hbar}\\!\int_{-T/2}^{T/2}\\!dt~{}\\!V^{I}_{\rm
int}(t)\right)\\!|\Omega_{0}\rangle~{}\\!,$ (5)

in which $|\Omega_{0}\rangle$ is the ground state of the free Hamiltonian in
the $N$-fermion subspace of the Fock space, $V^{I}_{\rm int}(t)$ is the theory
interaction operator $V_{\rm int}$ taken in the interaction picture and the
limits $T\rightarrow\infty$, $V\rightarrow\infty$ are implicit (T stands for
the chronological product.

If, as assumed, the underlying ”fundamental” interaction of spin $s$ fermions
is spin-independent and (in the infinite volume limit) Galileo invariant,
$V_{\rm int}$ in (5) is the effective theory interaction operator333The
absence of the $\sigma^{\prime}=\sigma$ term in the sum is the immediate
consequence of the anticommutativity of the fermionic field operators; in [7],
where only the spin $s=1/2$ was treated, the sum in (6) consisted a single
term only. of the form [2]

$\displaystyle V_{\rm
int}=C_{0}\\!\int\\!d^{3}\mathbf{x}\sum_{\sigma^{\prime}<\sigma}\psi_{\sigma^{\prime}}^{\dagger}\psi_{\sigma^{\prime}}\psi_{\sigma}^{\dagger}\psi_{\sigma}+V_{\rm
int}^{(C_{2})}+V_{\rm int}^{(C_{2}^{\prime})}+\dots,$ (6)

The terms $V_{\rm int}^{(C_{2})}$ and $V_{\rm int}^{(C_{2}^{\prime})}$, which
will be specified in Section 4, involve two spatial derivatives and are of
lower length dimension than the first one; the ellipsis stands for other
operators of yet lower length dimension. Despite the infinite number of
interaction terms in (6), to a term proportional to a given power444Terms
proportional to higher powers of $(kR)$ are also modified by logarithms [2].
of $(kR)$ in the computed quantity contributes only a finite number of Feynman
diagrams (constructed using only a finite subset of the interaction terms of
(6)) - only those dictated by the so-called power counting rules [17, 2]. The
couplings $C_{0}$, $C_{2}$, $C_{2}^{\prime},\dots$ of the effective theory are
fixed by matching the amplitude of the elastic fermion-fermion scattering
computed using the effective theory onto the general form of such an amplitude
expanded in powers of the relative particle momenta and parametrized by the
scattering lengths $a_{\ell}$, the effective radii $r_{\ell}$, etc. which with
the sufficient accuracy characterize the “fundamental” binary interaction. Of
course, since the effective interaction is local, divergences appear and must
be regularized (in [7] a cutoff $\Lambda$ on the wave vectors circulating in
loops have been imposed); the relations between computed physical quantities
like $a_{\ell}$, $r_{\ell}$, $E_{\Omega},\dots$ and the effective theory
parameters $C_{0}$, $C_{2},\dots$ depend then on the regularization but this
dependence disappears when the computed energy is expressed in terms of the
scattering lengths and radii. To the order which is needed to compute the
third order corrections to the ground state energy of the gas of fermions the
relations between the couplings $C_{0}$, $C_{2}$ and $C_{2}^{\prime}$ and the
scattering lengths $a_{0}$, $a_{1}$ and the radius $r_{0}$ are known [5] and
(when the cutoff $\Lambda$ is used as the regulator) read

$\displaystyle C_{0}(\Lambda)={4\pi\hbar^{2}\over
m_{f}}~{}\\!a_{0}\left(1+{2\over\pi}~{}\\!a_{0}\Lambda+{4\over\pi^{2}}~{}\\!a^{2}_{0}\Lambda^{2}+\dots\right),\phantom{aaaaaaa}\\!$
(7) $\displaystyle C_{2}(\Lambda)={4\pi\hbar^{2}\over m_{f}}~{}\\!{1\over
2}~{}\\!a_{0}^{2}r_{0}+\dots~{}\\!,\phantom{aaaa}C_{2}^{\prime}(\Lambda)={4\pi\hbar^{2}\over
m_{f}}~{}\\!a_{1}^{3}+\dots~{}\\!.$ (8)

(30,20)(25,70,290) (30,20)(25,290,70)2 (50,20)(25,110,250)2
(50,20)(25,250,110) (40,-2.5)2 (40,42.5)2
(155,30)[]$A_{\sigma^{\prime},\sigma}(q)=$(190,30)2 (210,35)(20,195,345)
(210,25)(20,165,15)2 (230,30)2
(215,50)[]${}^{-k+q~{}\sigma^{\prime}}$(215,8)[]${}_{k~{}\sigma}$(295,30)[]$B_{\sigma^{\prime},\sigma}(q)=$(330,30)2
(350,35)(20,195,345) (350,25)(20,15,165)2 (370,30)2
(355,50)[]${}^{k~{}\sigma^{\prime}}$(355,8)[]${}_{k+q~{}\sigma}$ Figure 1: The
only type of diagrams contributing to the order $(k_{\rm F}R)^{2}$ correction
to the ground state energy of the gas of spin $s$ fermions and two
“elementary” one-loop diagrams out of which the second order and the third
order corrections with the $C_{0}$ couplings can be constructed. Solid and
dashed lines denote propagators of fermions with different spin projections.

The power counting rules tell that the order $k_{\rm F}R$ and $(k_{\rm
F}R)^{2}$ corrections to the ground state energy of the gas of fermions that
is, the order $k_{\rm F}a_{0}$ and $(k_{\rm F}a_{0})^{2}$ terms in the formula
(2), are generated solely by the term in (6) proportional to $C_{0}$. The
correction of order $k_{\rm F}a_{0}$ is therefore given by the sum over all
possible choices of the (different) pairs of spin projections
$\sigma^{\prime}$ and $\sigma$ circulating in the two loops of a simple two-
loop diagram that can be formed by sewing together the corresponding pairs of
lines representing the interaction vertex proportional to $C_{0}$ while that
of order $(k_{\rm F}a_{0})^{2}$ are given by the similar sum over spin
projections circulating in the three-loop diagram shown in Figure 1 on the
left. This diagram can be composed in three equivalent ways (corresponding to
different assignments of the momenta to its internal lines):

$\displaystyle{1\over i\hbar}~{}\\!{E^{(2)}_{\Omega}\over V}={1\over
2!}\left({C_{0}\over
i\hbar}\right)^{2}\sum_{\sigma^{\prime}<\sigma}\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\![A_{\sigma^{\prime},\sigma}(q)]^{2}={1\over
2!}\left({C_{0}\over
i\hbar}\right)^{2}\sum_{\sigma^{\prime}<\sigma}\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\![B_{\sigma^{\prime},\sigma}(q)]^{2}$
$\displaystyle={1\over 2!}\left({C_{0}\over
i\hbar}\right)^{2}\sum_{\sigma^{\prime}<\sigma}\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\!B_{\sigma,\sigma}(q)~{}\\!B_{\sigma^{\prime},\sigma^{\prime}}(q)~{}\\!,\phantom{aaaaaaaaaaaaaaaaa}$
(9)

out of the two “elementary blocks” (corresponding to the two “elementary”
loops shown in Figure 1)

$\displaystyle
A_{\sigma^{\prime},\sigma}(q)=\int\\!{d^{3}\mathbf{k}\over(2\pi)^{3}}~{}\\!i\\!\left[{\theta(|\mathbf{k}-\mathbf{q}|-p_{{\rm
F}\sigma^{\prime}})~{}\\!\theta(|\mathbf{k}|-p_{{\rm F}\sigma})\over
q^{0}-\omega_{\mathbf{k}}-\omega_{\mathbf{k}-\mathbf{q}}+i0}-{\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{k}-\mathbf{q}|)~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{k}|)\over
q^{0}-\omega_{\mathbf{k}}-\omega_{\mathbf{k}-\mathbf{q}}-i0}\right],\phantom{a}$
(10) $\displaystyle
B_{\sigma^{\prime},\sigma}(q)=\int\\!{d^{3}\mathbf{k}\over(2\pi)^{3}}~{}\\!i\\!\left[-{\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{k}|)~{}\\!\theta(|\mathbf{k}+\mathbf{q}|-p_{{\rm
F}\sigma})\over
q^{0}+\omega_{\mathbf{k}}-\omega_{\mathbf{k}+\mathbf{q}}+i0}+{\theta(|\mathbf{k}|-p_{{\rm
F}\sigma^{\prime}})~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{k}+\mathbf{q}|)\over
q^{0}+\omega_{\mathbf{k}}-\omega_{\mathbf{k}+\mathbf{q}}-i0}\right],$ (11)

in which $\omega_{\mathbf{k}}=\hbar\mathbf{k}^{2}/2m_{f}$, obtained by the
standard integration over $k^{0}$ using the residue method. Each of these
“blocks” is a sum of two parts - those in the “elementary block”
$A_{\sigma^{\prime},\sigma}(q)$ which will be needed in the discussion of the
approach of [14], can be given the interpretation as corresponding to the
propagation of two particles (the first part) and the propagation of two holes
(the other part). In either form of (9) after the integration over $dq^{0}$
only two (out of the possible four) terms are nonzero (and equal to one
another), because in the two others all poles are on the same side of the real
$q^{0}$ axis.

The final form of the order $(k_{\rm F}a_{0})^{2}$ contribution to the energy
density is obtained (see [7]) by making the appropriate change of the
remaining integration variables. In this way one arrives at

$\displaystyle{E^{(2)}_{\Omega}\over
V}={64m_{f}C_{0}^{2}\over\hbar^{2}(2\pi)^{6}}\sum_{\sigma^{\prime}<\sigma}J(p_{{\rm
F}\sigma^{\prime}},~{}\\!p_{{\rm F}\sigma})~{}\\!,$ (12)

with the function $J(p_{{\rm F}\sigma},~{}\\!p_{{\rm
F}\sigma^{\prime}})=J(p_{{\rm F}\sigma^{\prime}},~{}\\!p_{{\rm F}\sigma})$
given by the integral

$\displaystyle J(p_{{\rm F}\sigma^{\prime}},~{}\\!p_{{\rm
F}\sigma})=\int_{0}^{s_{\rm max}}\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{u}+\mathbf{s}|)~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{t}-\mathbf{s}|)~{}\\!g_{\sigma^{\prime},\sigma}(t,s)~{}\\!,$
(13)

in which $s_{\rm max}={1\over 2}(p_{{\rm F}\sigma^{\prime}}+p_{{\rm
F}\sigma})$. The function $g_{\sigma^{\prime},\sigma}(t,s)\equiv
g_{\sigma^{\prime},\sigma}(|\mathbf{t}|,~{}\\!|\mathbf{s}|)=g_{\sigma,\sigma^{\prime}}(t,s)$
given by the integral (symmetric in the labels $\sigma$ and $\sigma^{\prime}$)

$\displaystyle g_{\sigma^{\prime},\sigma}(t,s)\equiv{1\over
4\pi}\\!\int\\!d^{3}\mathbf{u}~{}\\!{\theta(|\mathbf{u}+\mathbf{s}|-p_{{\rm
F}\sigma^{\prime}})\theta(|\mathbf{u}-\mathbf{s}|-p_{{\rm
F}\sigma})\over\mathbf{t}^{2}-\mathbf{u}^{2}+i0}$
$\displaystyle=-\Lambda+g^{\rm
fin}_{\sigma^{\prime},\sigma}(t,s)+{t^{2}\over\Lambda}+\dots,\phantom{aaaaaaaaaaaa}$
(14)

has been obtained in [7] in the analytic form which assumes (without loosing
generality) that $p_{{\rm F}\sigma^{\prime}}\leq p_{{\rm F}\sigma}$:

$\displaystyle g^{\rm fin}_{\sigma^{\prime},\sigma}(t,s)={1\over
2}~{}\\!p_{{\rm F}\sigma}+{t\over 4}\ln{(p_{{\rm
F}\sigma}-t)^{2}-s^{2}\over(p_{{\rm F}\sigma}+t)^{2}-s^{2}}+{p^{2}_{{\rm
F}\sigma}-s^{2}-t^{2}\over 8s}\ln{(p_{{\rm F}\sigma}+s)^{2}-t^{2}\over(p_{{\rm
F}\sigma}-s)^{2}-t^{2}}~{}\\!,$ (15)

when $0<s\leq{1\over 2}(p_{{\rm F}\sigma}-p_{{\rm F}\sigma^{\prime}})$ and

$\displaystyle g^{\rm fin}_{\sigma^{\prime},\sigma}(t,s)={1\over 4}(p_{{\rm
F}\sigma^{\prime}}+p_{{\rm F}\sigma}+2s)+{t\over 4}\ln{p_{{\rm
F}\sigma^{\prime}}+s-t\over p_{{\rm F}\sigma^{\prime}}+s+t}+{t\over
4}\ln{p_{{\rm F}\sigma}+s-t\over p_{{\rm F}\sigma}+s+t}\phantom{aaaaaaaaa}$
$\displaystyle+~{}\\!{p^{2}_{{\rm F}\sigma^{\prime}}-t^{2}-s^{2}\over
8s}\ln{(p_{{\rm F}\sigma^{\prime}}+s)^{2}-t^{2}\over
u_{0}^{2}-t^{2}}+{p^{2}_{{\rm F}\sigma}-t^{2}-s^{2}\over 8s}\ln{(p_{{\rm
F}\sigma}+s)^{2}-t^{2}\over u_{0}^{2}-t^{2}}~{}\\!,$ (16)

where

$\displaystyle u_{0}^{2}={1\over 2}(p_{{\rm F}\sigma^{\prime}}^{2}+p_{{\rm
F}\sigma}^{2})-s^{2}~{}\\!,$ (17)

when ${1\over 2}(p_{{\rm F}\sigma}-p_{{\rm F}\sigma^{\prime}})<s\leq s_{\rm
max}$.

It is convenient to write the function $J(p_{{\rm
F}\sigma^{\prime}},~{}\\!p_{{\rm F}\sigma})$ as the sum

$\displaystyle J=J_{\rm div}+J_{\rm fin}+J_{1/\Lambda}+\dots,$

where the successive terms directly correspond to the terms in (14). $J_{\rm
fin}$ is therefore independent of the cutoff $\Lambda$ while $J_{1/\Lambda}$
and the ellipsis stand for terms vanishing in the limit
$\Lambda\rightarrow\infty$. As

$\displaystyle\int_{0}^{s_{\rm max}}\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{u}+\mathbf{s}|)~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{t}-\mathbf{s}|)={p^{3}_{{\rm F}\sigma^{\prime}}p^{3}_{{\rm
F}\sigma}\over 72}~{}\\!,$ (18)

[7], $J_{\rm div}=-(\Lambda/72)p^{3}_{{\rm F}\sigma^{\prime}}p^{3}_{{\rm
F}\sigma}$ and it is straightforward to check that the divergences arising
from $J_{\rm div}$ cancel against similar terms proportional to $\Lambda$
arising when $C_{0}$ in the correction

$\displaystyle{E^{(1)}_{\Omega}\over V}={C_{0}\over
36\pi^{4}}\sum_{\sigma^{\prime}<\sigma}p^{3}_{{\rm
F}\sigma^{\prime}}p^{3}_{{\rm F}\sigma}~{}\\!,$ (19)

arising from the mentioned two-loop diagrams is expressed in terms of the
$s$-wave scattering length $a_{0}$ as in (7). The cutoff independent terms of
(19) give then the term proportional to $k_{\rm F}a_{0}$ in (2) and the terms
proportional to inverse powers of $\Lambda$ (which would be absent, if the
dimensional regularization was used instead of the cutoff $\Lambda$) can be,
if the computation is restricted to the second order, simply discarded. In the
computation including higher order corrections they cancel against similar
spurious contributions arising in higher orders.

If the densities of fermions having different spin projections are all equal
($N_{\sigma}=N/g_{s}$, $p_{{\rm F}\sigma}=k_{\rm F}$ for all
$\sigma=1,\dots,g_{s}$), one numerically finds that

$\displaystyle J_{\rm fin}(k_{\rm F},k_{\rm F})={k_{\rm F}^{7}\over
840}~{}\\!(11-2\ln 2)~{}\\!.$

The sum over pairs of different spin projections gives in this case the factor
${1\over 2}g_{s}(g_{s}-1)$ and one recovers the well known old result re-
derived in [2] using the effective field theory. If the densities of fermions
with different spin projections are not equal, numerical evaluation of the
integrals defining the function $J(p_{{\rm F}\sigma^{\prime}},~{}\\!p_{{\rm
F}\sigma})$ shows it is proportional to the function $J_{K}$ (4) introduced by
Kanno the precise relation being555From the formulae (13) and (14) it is clear
that $J(p_{{\rm F}\sigma^{\prime}},~{}\\!p_{{\rm F}\sigma})\equiv
J(x_{\sigma^{\prime}}k_{\rm F},~{}\\!x_{\sigma}k_{\rm F})=k_{\rm
F}^{7}J(x_{\sigma^{\prime}},~{}\\!x_{\sigma})$.

$\displaystyle 160~{}\\!J_{\rm fin}(p_{{\rm F}\sigma^{\prime}},~{}\\!p_{{\rm
F}\sigma})=k_{\rm F}^{7}J_{K}(x_{\sigma^{\prime}},~{}\\!x_{\sigma})~{}\\!.$

This leads to the result (2) obtained recently in a somewhat different way in
[8].

##  3 Third order corrections proportional to $C_{0}^{3}$

According to the power counting rules the order $(k_{\rm F}R)^{3}$ corrections
to the ground state energy density arise from the four-loop diagrams with
three interaction vertices proportional to $C_{0}$ and from the two-loop
diagrams with a single $C_{2}$ or $C_{2}^{\prime}$ interaction vertex (the
terms $V_{\rm int}^{(C_{2})}$ and $V_{\rm int}^{(C_{2}^{\prime})}$ in (6)). In
this section we work out the contribution of the first class of diagrams; the
contribution of diagrams with the vertices generated by the interaction terms
proportional to $C_{2}$ and $C_{2}^{\prime}$ will be computed in Section 4.

(40,40)(40,330,90) (40,40)(40,90,210) (40,40)(40,180,360) (74,20)2 (6,20)2
(40,80)2 (40,-53)(80,115,65)2 (120,86)(80,235,185)2 (-40,86)(80,355,305)2
(200,40)(40,330,90) (200,40)(40,90,210) (200,40)(40,180,360) (234,20)2
(166,20)2 (200,80)2 (200,-53)(80,65,115)2 (280,86)(80,185,235)2
(120,86)(80,305,355)2 (320,40)(40,330,90)3 (320,40)(40,90,210)1
(320,40)(40,210,330) (354,20)2 (286,20)2 (320,80)2 (320,-53)(80,65,115)
(400,86)(80,185,235)3 (240,86)(80,305,355)1  Figure 2: The single particle-
particle and two different particle-hole vacuum diagrams contributing in the
order $C^{3}_{0}$ to $E_{\Omega}$. Different types of lines represent
propagators of fermions having different spin projections. The second
particle-hole diagram exists only if there are more than two different spin
projections ($s>{1\over 2}$).

If $s>{1\over 2}$ three kinds of nonvanishing four-loop diagrams with three
interaction vertices proportional to the $C_{0}$ coupling which can be formed.
They are shown in Figure 2. The first one is the so-called particle-particle
diagram while the remaining two are called particle-hole diagrams.

The contribution of the order $C_{0}^{3}$ particle-particle diagram shown in
Figure 2 on the left is (2 is the combinatoric factor)

$\displaystyle{1\over i\hbar}~{}\\!{E^{(3){\rm pp}}_{\Omega}\over V}={1\over
3!}\left({C_{0}\over
i\hbar}\right)^{3}2\sum_{\sigma^{\prime}<\sigma}\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\![A_{\sigma^{\prime},\sigma}(q)]^{3}~{}\\!,$
(20)

As it is clear from the spin structure of the interaction term proportional to
$C_{0}$ and the flow of the spin projections in the diagram, the sum runs only
over pairs of different spin projection labels $\sigma^{\prime}\neq\sigma$.
There are $g_{s}(g_{s}-1)/2$ such pairs and, if densities of fermions with
different spin projections are all equal, the summation reproduces the usual
spin factor associated with this contribution [2]. For the sake of further
discussion we note that the product of three $A$-”blocks” gives rise to
$2^{3}=8$ terms out of two are killed by the integration over $dq^{0}$ (they
have all poles on the same side of the real $q^{0}$ axis) while the remaining
six split into three identical terms which arise from the products of two
terms corresponding in (10) \- appealing to the interpretation of the two
terms of the $A$-”block” - to the propagation of “particles” and one
corresponding to the propagation of holes and another three identical ones in
which this composition is reversed. As a result this contribution to the
energy density can be written in terms of two functions, $G^{(1)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})$, and $G^{(2)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})$ (corresponding respectively to these
two kinds of compositions):

$\displaystyle{E^{(3){\rm pp}}_{\Omega}\over
V}={128m_{f}^{2}C_{0}^{3}\over(2\pi)^{8}\hbar^{4}}\sum_{\sigma^{\prime}<\sigma}\left[G^{(1)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})+G^{(2)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})\right].$ (21)

The functions $G^{(1)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$ and
$G^{(2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$ are given by the
integrals

$\displaystyle G^{(1)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})=\\!\int_{0}^{s_{\rm max}}\\!\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{t}+\mathbf{s}|)~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{t}-\mathbf{s}|)~{}\\!(g_{\sigma^{\prime},\sigma}(t,s))^{2}~{}\\!,\phantom{a}\\!$
(22) $\displaystyle G^{(2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})=\\!\int_{0}^{s_{\rm max}}\\!\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(|\mathbf{t}+\mathbf{s}|-p_{{\rm
F}\sigma^{\prime}})~{}\\!\theta(|\mathbf{t}-\mathbf{s}|-p_{{\rm
F}\sigma})~{}\\!(h_{\sigma^{\prime},\sigma}(t,s))^{2}~{}\\!,\phantom{a}\\!$
(23)

in which $g_{\sigma^{\prime},\sigma}(t,s)$ is the function defined by (14) and
$h_{\sigma^{\prime},\sigma}(t,s)\equiv
h_{\sigma^{\prime},\sigma}(|\mathbf{t}|,|\mathbf{s}|)=h_{\sigma,\sigma^{\prime}}(t,s)$
is given by the finite integral

$\displaystyle h_{\sigma^{\prime},\sigma}(t,s)={1\over
4\pi}\\!\int\\!d^{3}\mathbf{u}~{}\\!{\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{u}+\mathbf{s}|)\theta(p_{{\rm
F}\sigma}-|\mathbf{u}-\mathbf{s}|)\over\mathbf{t}^{2}-\mathbf{u}^{2}-i0}~{}\\!.$
(24)

The explicit form of the function $h_{\sigma^{\prime},\sigma}(t,s)$ has been
obtained666Analysis of the integration domains of (24) shows that the pole at
$\mathbf{t}^{2}=\mathbf{u}^{2}$ is never reached. using the same method (based
on the trick introduced in Appendix C of [20]) by which the function
$g_{\sigma^{\prime},\sigma}(t,s)$ has been computed in [7]. Assuming without
loss of generality that $p_{{\rm F}\sigma^{\prime}}\leq p_{{\rm F}\sigma}$, it
reads

$\displaystyle h_{\sigma^{\prime},\sigma}(t,s)=-{1\over 2}~{}\\!p_{{\rm
F}\sigma^{\prime}}-{t\over 4}\ln{t-(p_{{\rm F}\sigma^{\prime}}-s)\over
t+(p_{{\rm F}\sigma^{\prime}}-s)}-{t\over 4}\ln{t-(p_{{\rm
F}\sigma^{\prime}}+s)\over t+(p_{{\rm
F}\sigma^{\prime}}+s)}\phantom{aaaaaaaa}$ $\displaystyle+~{}\\!{t^{2}-(p_{{\rm
F}\sigma^{\prime}}^{2}-s^{2})\over 8s}\ln{t^{2}-(p_{{\rm
F}\sigma^{\prime}}+s)^{2}\over t^{2}-(p_{{\rm
F}\sigma^{\prime}}-s)^{2}}~{}\\!,\phantom{aaaaaaaaaaaaaaaaa}~{}$ (25)

if $0\leq s\leq{1\over 2}(p_{{\rm F}\sigma}-p_{{\rm F}\sigma^{\prime}})$, and
($u^{2}_{0}$ is again given by (17))

$\displaystyle h_{\sigma^{\prime},\sigma}(t,s)={1\over 2}\left(2s-p_{{\rm
F}\sigma^{\prime}}-p_{{\rm F}\sigma}\right)-{t\over 4}\ln{t-(p_{{\rm
F}\sigma^{\prime}}-s)\over t+(p_{{\rm F}\sigma^{\prime}}-s)}-{t\over
4}\ln{t-(p_{{\rm F}\sigma}-s)\over t+(p_{{\rm F}\sigma}-s)}\phantom{aaaaaa}$
$\displaystyle-{1\over 8s}\left[(p_{{\rm F}\sigma}-s)^{2}+(p_{{\rm
F}\sigma^{\prime}}-s)^{2}-2u_{0}^{2}\right]\phantom{aaaaaaaaaaaaaaaaaaaaaaaaaaa}$
(26) $\displaystyle-~{}\\!{t^{2}-p_{{\rm F}\sigma}^{2}+s^{2}\over
8s}\ln{t^{2}-(p_{{\rm F}\sigma}-s)^{2}\over t^{2}-u_{0}^{2}}-{t^{2}-p_{{\rm
F}\sigma^{\prime}}^{2}+s^{2}\over 8s}\ln{t^{2}-(p_{{\rm
F}\sigma^{\prime}}-s)^{2}\over t^{2}-u_{0}^{2}}~{}\\!,$

if ${1\over 2}(p_{{\rm F}\sigma}-p_{{\rm F}\sigma^{\prime}})\leq s\leq s_{\rm
max}={1\over 2}(p_{{\rm F}\sigma}+p_{{\rm F}\sigma^{\prime}})$. Expanding both
these expression one finds that the function $h_{\sigma^{\prime},\sigma}(t,s)$
behaves as $1/t^{2}$ when $t=|\mathbf{t}|\rightarrow\infty$. Therefore,
although the integral over $d^{3}\mathbf{t}$ in (23) covers the infinite
exterior of the two Fermi spheres, it is convergent - unlike the function
$G^{(1)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$ which depends on the
divergent function $g_{\sigma^{\prime},\sigma}(t,s)$, the function
$G^{(2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$ is finite.

Since the remaining order $C_{0}^{3}$ contributions to $E_{\Omega}/V$ will
turn out finite, already at this point one can check the cancellation of
ultraviolet divergences in the expression

$\displaystyle{\Delta E_{\Omega}\over
V}=\sum_{\sigma^{\prime}<\sigma}\left[{p_{{\rm F}\sigma^{\prime}}^{3}p_{{\rm
F}\sigma}^{3}\over
36\pi^{4}}~{}\\!C_{0}+{64m_{f}C_{0}^{2}\over(2\pi)^{6}\hbar^{2}}~{}\\!J(p_{{\rm
F}\sigma^{\prime}},p_{{\rm
F}\sigma})+{128m_{f}^{2}C_{0}^{3}\over(2\pi)^{8}\hbar^{4}}~{}\\!G^{(1)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})\right]+\dots$

Inserting here the expansions (7) of $C_{0}$ (truncated to the appropriate
order), writing

$\displaystyle(g_{\sigma^{\prime},\sigma}(t,s))^{2}=\Lambda^{2}-2\Lambda~{}\\!g^{\rm
fin}_{\sigma^{\prime},\sigma}+(g^{\rm
fin}_{\sigma^{\prime},\sigma})^{2}-2t^{2}+\dots$ (27)

and using result (18) it is a matter of a simple algebra to check that all
terms diverging with $\Lambda\rightarrow\infty$ as well as the spurious terms
specific for the adopted regularization (that is the finite term arising from
the product of the term proportional to $a_{0}\Lambda$ and the term
$t^{2}/\Lambda$ from $g_{\sigma^{\prime},\sigma}(t,s)$ in the contribution
proportional to $C_{0}^{2}$ and the term $-2t^{2}$ in (27)) completely cancel
out as they should. As a result, to obtain the relevant contribution to
$E_{\Omega}/V$ one can simply drop all $\Lambda$ dependent terms in the
function $g_{\sigma^{\prime},\sigma}(t,s)$ and set
$C_{0}=(4\pi\hbar^{2}/m_{f})a_{0}$ (the remaining part of the function
$G^{(1)}$ will be denoted $G^{(1)}_{\rm fin}$).

The first type of the particle-hole diagrams (the middle one in Figure 2) has
been computed in [16] in the case of spin $s={1\over 2}$ fermions. The
generalization of the formulae obtained there is straightforward and reduces
to taking again the sum over pairs of different spin projections:

$\displaystyle{1\over i\hbar}~{}\\!{E^{(3){\rm ph1}}_{\Omega}\over V}={1\over
3!}\left({C_{0}\over
i\hbar}\right)^{3}2\sum_{\sigma^{\prime}<\sigma}\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\![B_{\sigma^{\prime},\sigma}(q)]^{3}~{}\\!.$
(28)

(Again 2 is the combinatoric factor). This gives

$\displaystyle{E_{\Omega}^{(3){\rm ph1}}\over
V}=-{32m_{f}^{2}C_{0}^{3}\over(2\pi)^{8}\hbar^{4}}\sum_{\sigma^{\prime}<\sigma}\left[K^{(1)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})+K^{(2)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})\right],$ (29)

where the functions $K^{(1)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$
and $K^{(2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$ are given by

$\displaystyle K^{(1)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})=\int_{0}^{\infty}\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(|\mathbf{t}+\mathbf{s}|-p_{{\rm
F}\sigma^{\prime}})~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{t}-\mathbf{s}|)~{}\\!(f_{\sigma^{\prime},\sigma}^{(1)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s))^{2}~{}\\!,\phantom{a}$
(30) $\displaystyle K^{(2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})=\int_{0}^{\infty}\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{t}+\mathbf{s}|)~{}\\!\theta(|\mathbf{t}-\mathbf{s}|-p_{{\rm
F}\sigma})~{}\\!(f_{\sigma^{\prime},\sigma}^{(2)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s))^{2}~{}\\!,\phantom{a}$
(31)

and the functions
$f_{\sigma^{\prime},\sigma}^{(1)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)$
and
$f_{\sigma^{\prime},\sigma}^{(2)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)$
are defined by the integrals

$\displaystyle
f_{\sigma^{\prime},\sigma}^{(1)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)={1\over
4\pi}\\!\int\\!d^{3}\mathbf{u}~{}\\!{\theta(p_{{\rm
F}\sigma^{\prime}}-|\mathbf{u}+\mathbf{s}|)~{}\\!\theta(|\mathbf{u}-\mathbf{s}|-p_{{\rm
F}\sigma})\over(\mathbf{u}-\mathbf{t})\cdot\mathbf{s}+i0}~{}\\!,$ (32)

and

$\displaystyle
f_{\sigma^{\prime},\sigma}^{(2)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)={1\over
4\pi}\\!\int\\!d^{3}\mathbf{u}~{}\\!{\theta(|\mathbf{u}+\mathbf{s}|-p_{{\rm
F}\sigma^{\prime}})~{}\\!\theta(p_{{\rm
F}\sigma}-|\mathbf{u}-\mathbf{s}|)\over(\mathbf{u}-\mathbf{t})\cdot\mathbf{s}-i0}~{}\\!.$
(33)

Finally the contribution of the last type of diagrams (the rightmost one in
Figure 2) is given by

$\displaystyle{1\over i\hbar}~{}\\!{E^{(3){\rm ph2}}_{\Omega}\over V}={1\over
3!}\left({C_{0}\over
i\hbar}\right)^{3}(-3!)\\!\sum_{\sigma^{\prime\prime},\sigma^{\prime},\sigma}\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\!B_{\sigma^{\prime\prime},\sigma^{\prime\prime}}(q)~{}\\!B_{\sigma^{\prime},\sigma^{\prime}}(q)~{}\\!B_{\sigma,\sigma}(q)~{}\\!,$
(34)

where $-3!$ is the sign/combinatoric factor and the sum is over all triplets
of different spin projections
($\sigma\neq\sigma^{\prime}\neq\sigma^{\prime\prime}$,
$\sigma\neq\sigma^{\prime\prime}$). After the standard steps this leads to

$\displaystyle{E_{\Omega}^{(3){\rm ph2}}\over
V}={32m_{f}^{2}C_{0}^{3}\over(2\pi)^{8}\hbar^{4}}\sum_{\sigma^{\prime\prime},\sigma^{\prime},\sigma}\left\\{\tilde{K}^{(1)}(p_{{\rm
F}\sigma^{\prime\prime}};~{}\\!p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})+\tilde{K}^{(2)}(p_{{\rm F}\sigma^{\prime\prime}};~{}\\!p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})\right.\phantom{a}~{}$
$\displaystyle+\tilde{K}^{(1)}(p_{{\rm F}\sigma^{\prime}};~{}\\!p_{{\rm
F}\sigma},p_{{\rm F}\sigma^{\prime\prime}})+\tilde{K}^{(2)}(p_{{\rm
F}\sigma^{\prime}};~{}\\!p_{{\rm F}\sigma},p_{{\rm
F}\sigma^{\prime\prime}})\phantom{aa}$ (35)
$\displaystyle\left.+\tilde{K}^{(1)}(p_{{\rm F}\sigma};~{}\\!p_{{\rm
F}\sigma^{\prime\prime}},p_{{\rm F}\sigma^{\prime}})+\tilde{K}^{(2)}(p_{{\rm
F}\sigma};~{}\\!p_{{\rm F}\sigma^{\prime\prime}},p_{{\rm
F}\sigma^{\prime}})\right\\},$

where the functions $\tilde{K}^{(1)}(p_{{\rm
F}\sigma^{\prime\prime}};~{}\\!p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$,
and $\tilde{K}^{(2)}(p_{{\rm F}\sigma^{\prime\prime}};~{}\\!p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})$, symmetric in their second and third
Fermi momentum labels, are given by the following nested integrals

$\displaystyle\tilde{K}^{(1)}(p_{{\rm F}\sigma^{\prime\prime}};~{}\\!p_{{\rm
F}\sigma^{\prime}},p_{{\rm
F}\sigma})\phantom{aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa}$
$\displaystyle=\int_{0}^{\infty}\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(|\mathbf{t}+\mathbf{s}|-p_{{\rm
F}\sigma^{\prime\prime}})~{}\\!\theta(p_{{\rm
F}\sigma^{\prime\prime}}-|\mathbf{t}-\mathbf{s}|)~{}\\!f_{\sigma^{\prime},\sigma^{\prime}}^{(1)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)~{}\\!f_{\sigma,\sigma}^{(1)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)~{}\\!,\phantom{a}$
(36) $\displaystyle\tilde{K}^{(2)}(p_{{\rm
F}\sigma^{\prime\prime}};~{}\\!p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})\phantom{aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa}$
$\displaystyle=\int_{0}^{\infty}\\!ds~{}\\!s^{2}~{}\\!{1\over
4\pi}\\!\int\\!d^{3}\mathbf{t}~{}\\!\theta(p_{{\rm
F}\sigma^{\prime\prime}}-|\mathbf{t}+\mathbf{s}|)~{}\\!\theta(|\mathbf{t}-\mathbf{s}|-p_{{\rm
F}\sigma^{\prime\prime}})~{}\\!f_{\sigma^{\prime},\sigma^{\prime}}^{(2)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)~{}\\!f_{\sigma,\sigma}^{(2)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)~{}\\!,\phantom{a}$
(37)

However, the change $\mathbf{s}\rightarrow-\mathbf{s}$ in the integrals
defining the functions $\tilde{K}^{(1)}$ and $\tilde{K}^{(2)}$ makes it
obvious that

$\displaystyle\tilde{K}^{(2)}(p_{{\rm F}\sigma^{\prime\prime}};~{}\\!p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})=\tilde{K}^{(1)}(p_{{\rm
F}\sigma^{\prime\prime}};~{}\\!p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})~{}\\!.$ (38)

This allows to simplify the right hand side of (35) by retaining in it (and
doubling) only, say, the functions $\tilde{K}^{(2)}$.

If the densities $N_{\sigma}/V$ of fermions having different spin projection
are all equal ($p_{{\rm F}\sigma}=p_{{\rm F}\sigma^{\prime}}=p_{{\rm
F}\sigma^{\prime\prime}}=k_{\rm F}$), then $\tilde{K}^{(2)}(k_{\rm
F};~{}\\!k_{\rm F},k_{\rm F})=K^{(1)}(k_{\rm F},k_{\rm F})=K^{(2)}(k_{\rm
F},k_{\rm F})$. The spin factors of the sum of the contributions (29) and (35)
combine in this case into

$\displaystyle-{g_{s}(g_{s}-1)\over 2}+3~{}\\!{g_{s}(g_{s}-1)(g_{s}-2)\over
6}={1\over 2}~{}\\!g_{s}(g_{s}-1)(g_{s}-3)~{}\\!,$

which is the spin factor associated with the particle-hole diagram in this
case [2].

Analysis of the integration domains in the formulae defining the $K$-functions
one finds that the poles of the integrands at
$\mathbf{u}\cdot\mathbf{s}=\mathbf{t}\cdot\mathbf{s}$ are never reached, so
the prescriptions $\pm i0$ can be dropped. Still, obtaining the analytic forms
of the functions
$f_{\sigma^{\prime},\sigma}^{(1,2)}(\mathbf{t}\\!\cdot\\!\mathbf{s},~{}\\!s)$
is rather cumbersome - hence we do not give the technical details of the
computation. The final formulae (different in the regimes of small,
intermediate and large $s$) for these functions can be found in [16]; the
formulae for the function
$f^{(2)}_{\sigma,\sigma}(\mathbf{s}\cdot\mathbf{t},s)$ (not given in [16]) can
in principle be obtained as a $p_{{\rm F}\sigma^{\prime}}=p_{{\rm F}\sigma}$
limit of the function
$f^{(2)}_{\sigma^{\prime},\sigma}(\mathbf{s}\cdot\mathbf{t},s)$ but for this
particular case a much simpler formula can be obtained.777We have, of course,
checked numerically that
$f^{(2)}_{\sigma^{\prime},\sigma}(\mathbf{s}\cdot\mathbf{t},s)$ agrees with
$f^{(2)}_{\sigma,\sigma}(\mathbf{s}\cdot\mathbf{t},s)$ in the limit $p_{{\rm
F}\sigma^{\prime}}=p_{{\rm F}\sigma}$. It reads:

$\displaystyle
f^{(2)}_{\sigma\sigma}(\mathbf{t}\cdot\mathbf{s},s)=s~{}\\!(p_{{\rm
F}\sigma}-t)\phantom{aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa}$
$\displaystyle+{1\over 4}(p_{{\rm
F}\sigma}^{2}-s^{2}+2st)\ln\\!\left[{(t^{2}-p_{{\rm
F}\sigma}^{2}+s^{2}-2st)(tp_{{\rm F}\sigma}-p_{{\rm
F}\sigma}^{2}+s^{2}-st)\over t^{2}(tp_{{\rm F}\sigma}+k^{2}-s^{2}+st)}\right]$
$\displaystyle-{1\over 4}(p_{{\rm
F}\sigma}^{2}-s^{2}-2st)\ln\\!\left[{(t^{2}-p_{{\rm
F}\sigma}^{2}+s^{2}+2st)(tp_{{\rm F}\sigma}-p_{{\rm
F}\sigma}^{2}+s^{2}+st)\over t^{2}(tp_{{\rm F}\sigma}+p_{{\rm
F}\sigma}^{2}-s^{2}-st)}\right]$ $\displaystyle+{1\over 2}(p_{{\rm
F}\sigma}^{2}-s^{2})\ln{p_{{\rm F}\sigma}+s\over p_{{\rm
F}\sigma}-s}-t^{2}\ln{t-s-p_{{\rm F}\sigma}\over t+s-p_{{\rm
F}\sigma}}~{}\\!\phantom{aaaaaaaaaaaaaaaaaaaa}$

$f^{(2)}_{\sigma\sigma}(\mathbf{t}\cdot\mathbf{s},s)=-f^{(1)}_{\sigma\sigma}(\mathbf{t}\cdot\mathbf{s},s)$.

With the analytical formulae for these functions the remaining integrals in
(30), (31) and (37) can be evaluated numerically using the standard
Mathematica built-in routine for multidimensional integration over a
prescribed domain.

The resulting, order $(k_{\rm F}a_{0})^{3}$ correction, that adds up to the
terms in the curly bracket in (2) is888From the defining formulae it is clear
that all the functions $G^{(1,2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm
F}\sigma})$, $K^{(1,2)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$ and
$\tilde{K}^{2)}(p_{{\rm F}\sigma^{\prime\prime}};p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})$ scale like the eight power of the wave
vector $k_{\rm F}$ which effectively reduces them to functions of the single
variable $r=x_{\sigma^{\prime}}/x_{\sigma}$: $G^{(1)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})=p_{{\rm
F}\sigma}^{8}G^{(1)}(x_{\sigma^{\prime}}/x_{\sigma},1)$, etc.

$\displaystyle{160\over\pi^{3}}~{}\\!(k_{\rm
F}a_{0})^{3}\left(\sum_{\sigma^{\prime}<\sigma}\\!\left[4G^{(1)}_{\rm
fin}(x_{\sigma^{\prime}},x_{\sigma})+4G^{(2)}(x_{\sigma^{\prime}},x_{\sigma})-K^{(1)}(x_{\sigma^{\prime}},x_{\sigma})-K^{(2)}(x_{\sigma^{\prime}},x_{\sigma})\right]\right.$
$\displaystyle\left.+\sum_{\sigma^{\prime\prime},\sigma^{\prime},\sigma}\\!\\!\left[2\tilde{K}^{(2)}(x_{\sigma^{\prime\prime}};x_{\sigma^{\prime}},x_{\sigma})+2\tilde{K}^{(2)}(x_{\sigma^{\prime}};x_{\sigma},x_{\sigma^{\prime\prime}})+2\tilde{K}^{(2)}(x_{\sigma};x_{\sigma^{\prime\prime}},x_{\sigma^{\prime}})\right]\right),\phantom{a}$
(39)

where the sum in the last line is over all triplets of different indices
$\sigma^{\prime\prime},\sigma^{\prime},\sigma$.

Figure 3: Plot of the functions $G^{(1)}_{\rm fin}(r,1)$ (the highest, blue,
curve), $G^{(2)}(r,1)$ (the lowest, red, curve) and ${1\over
4}[K^{(1)}(r,1)+K^{(2)}(r,1)]$ (the middle, green, curve).

With these results one can try to verify the reliability of the approach
developed in [14]. Its authors have considered only the case of spin $1/2$
fermions and only the contributions to the ground state energy density which,
like the order $C_{0}^{2}$ contribution (9) and the order $C_{0}^{3}$
contribution (20), can be composed out of the $A$-”blocks”. Of the
corresponding $n$-th order contribution ($(n-1)!$ is the combinatoric factor)

$\displaystyle{1\over i\hbar}~{}\\!{E^{(n)\rm pp}_{\Omega}\over V}={1\over
n!}\left({C_{0}\over
i\hbar}\right)^{n}(n-1)!\sum_{\sigma^{\prime}<\sigma}\\!\int\\!{d^{4}q\over(2\pi)^{4}}~{}\\![A_{\sigma^{\prime},\sigma}(q)]^{n}~{}\\!,$
(40)

which gives rise to $2^{n}-2$ nonvanishing terms they retain only $n$
identical ones obtained by taking $n-1$ parts of the $A$-”block” corresponding
to the propagation of particles and one part corresponding to the propagation
of holes (again appealing to the interpretation introduced earlier). Thus, out
of the complete order $C_{0}^{3}$ contribution obtained above they take only
the function $G^{(1)}_{\rm fin}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})$
discarding the contributions of the function $G^{(2)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})$ and of the sum ${1\over
4}[K^{(1)}(p_{{\rm F}\sigma^{\prime}},p_{{\rm F}\sigma})+K^{(2)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})]$. As the plot 3 shows, the first of
these two functions is indeed one order of magnitude smaller than
$G^{(1)}_{\rm fin}$ in the entire range of the ratio $r=p_{{\rm
F}\sigma^{\prime}}/p_{{\rm F}\sigma}$ (both these functions, similarly as the
functions $K^{(1,2)}$ and $\tilde{K}^{(2)}$ vanish at $r=0$ as a result of the
absence of $s$-wave interaction of fermions obeying the Pauli exclusion
principle) but the same is not true of the second function which is always
comparable with $G^{(1)}_{\rm fin}$ and enters (39) with the opposite sign
strongly reducing the contribution of the function $G^{(1)}_{\rm fin}$.
Moreover, in higher orders the number of discarded terms arising from (40)
rapidly grows and if each of the discarded $(2^{n}-2)-n$ terms of the
contribution (40) is only one order of magnitude smaller than each of the $n$
terms taken into account, the domination of the latter can be invalidated (as
the transition to the ordered state is expected to occur at $k_{\rm
F}a_{0}\stackrel{{\scriptstyle<}}{{{}_{\sim}}}1$, the higher order
contributions are not much suppressed by the powers of the expansion parameter
$k_{\rm F}a_{0}$). Moreover, in higher orders the number of other types of
diagrams (with topologies different than those shown in Figure 2) grows and
their contribution, like the one of $K^{(1)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})+K^{(2)}(p_{{\rm
F}\sigma^{\prime}},p_{{\rm F}\sigma})$, may dominate the one of the retained
terms. Thus the reliability of the approximations made in [14] may be
questioned even though the numerical results obtained there agree quite well
with those obtained by the Quantum Monte Carlo methods.

##  4 Contribution of the interaction terms proportional to $C_{2}$ and
$C_{2}^{\prime}$

According to the power counting rule in the order $k^{3}_{\rm
F}(\hbar^{2}k_{\rm F}^{2}/2m_{f})(k_{\rm F}R)^{3}$ one has to take into
account also the operators:

$\displaystyle V^{(C_{2})}_{\rm int}=-{C_{2}\over
8}\int\\!d^{3}\mathbf{x}\sum_{\sigma<\sigma^{\prime}}\left[\psi^{\dagger}_{\sigma}\psi^{\dagger}_{\sigma^{\prime}}(\psi_{\sigma^{\prime}}\mbox{\boldmath{$\nabla$}}^{2}\psi_{\sigma}-2\mbox{\boldmath{$\nabla$}}\psi_{\sigma^{\prime}}\\!\cdot\\!\mbox{\boldmath{$\nabla$}}\psi_{\sigma}+\mbox{\boldmath{$\nabla$}}^{2}\psi_{\sigma^{\prime}}~{}\\!\psi_{\sigma})\right.\phantom{a}~{}$
$\displaystyle\left.+(\mbox{\boldmath{$\nabla$}}^{2}\psi^{\dagger}_{\sigma}~{}\\!\psi^{\dagger}_{\sigma^{\prime}}-2\mbox{\boldmath{$\nabla$}}\psi^{\dagger}_{\sigma}\\!\cdot\\!\mbox{\boldmath{$\nabla$}}\psi^{\dagger}_{\sigma^{\prime}}+\psi^{\dagger}_{\sigma}\mbox{\boldmath{$\nabla$}}^{2}\psi^{\dagger}_{\sigma^{\prime}})\psi_{\sigma^{\prime}}\psi_{\sigma}\right],$
(41)

and

$\displaystyle V_{\rm int}^{(C_{2}^{\prime})}=-{C^{\prime}_{2}\over
8}\\!\int\\!d^{3}\mathbf{x}\sum_{\sigma,{\sigma^{\prime}}}(\mbox{\boldmath{$\nabla$}}\psi^{\dagger}_{\sigma}\psi^{\dagger}_{\sigma^{\prime}}-\psi^{\dagger}_{\sigma}\mbox{\boldmath{$\nabla$}}\psi^{\dagger}_{\sigma^{\prime}})\\!\cdot\\!(\mbox{\boldmath{$\nabla$}}\psi_{\sigma^{\prime}}\psi_{\sigma}-\psi_{\sigma^{\prime}}{\mbox{\boldmath{$\nabla$}}}\psi_{\sigma})~{}\\!.$
(42)

(Notice that in the second operator the sum involves also terms with
$\sigma={\sigma^{\prime}}$). Their contribution to the ground state energy
density in the case of nonequal densities of spin $1/2$ fermions having
different spin projections has been explicitly written down in [16].
Generalization of these formulae to the case of spin $s$ fermions with
different densities of their $g_{s}$ spin projections is straightforward. They
read

$\displaystyle{E^{(C_{2})}_{\Omega}\over V}={C_{2}\over
240\pi^{4}}\sum_{\sigma^{\prime}<\sigma}p_{{\rm F}\sigma^{\prime}}^{3}p_{{\rm
F}\sigma}^{3}(p_{{\rm F}\sigma^{\prime}}^{2}+p_{{\rm
F}\sigma}^{2})\phantom{aaaaaaaaaaaaaaaaaaaaaa}~{}$ $\displaystyle={k_{\rm
F}^{3}\over 6\pi^{2}}~{}\\!{3\over 5}~{}\\!{\hbar^{2}k_{\rm F}^{2}\over
2m_{f}}~{}\\!{1\over 6\pi}~{}\\!k_{\rm
F}^{3}a_{0}^{2}r_{0}\sum_{\sigma^{\prime}<\sigma}x^{3}_{\sigma^{\prime}}x^{3}_{\sigma}(x^{2}_{\sigma^{\prime}}+x^{2}_{\sigma})~{}\\!.\phantom{aaaaaaaaaaaa}$
(43) $\displaystyle{E^{(C^{\prime}_{2})}_{\Omega}\over V}={C^{\prime}_{2}\over
120\pi^{4}}\left\\{\sum_{\sigma=1}^{g_{s}}p_{{\rm F}\sigma}^{8}+{1\over
2}\sum_{\sigma^{\prime}<\sigma}p_{{\rm F}\sigma^{\prime}}^{3}p_{{\rm
F}\sigma}^{3}(p_{{\rm F}\sigma^{\prime}}^{2}+p_{{\rm
F}\sigma}^{2})\right\\}\phantom{aaaaaaaaaa}$ $\displaystyle={k_{\rm
F}^{3}\over 6\pi^{2}}~{}\\!{3\over 5}~{}\\!{\hbar^{2}k_{\rm F}^{2}\over
2m_{f}}~{}\\!{2\over 3\pi}~{}\\!(k_{\rm
F}a_{1})^{3}\left\\{\sum_{\sigma=1}^{g_{s}}x_{\sigma}^{8}+{1\over
2}\sum_{\sigma^{\prime}<\sigma}x^{3}_{\sigma^{\prime}}x^{3}_{\sigma}(x^{2}_{\sigma^{\prime}}+x^{2}_{\sigma})\right\\}.$
(44)

For equal all densities of all spin projections this reduces to the known
results given e.g. in [2]. The contribution $E^{(C^{\prime}_{2})}_{\Omega}$ is
special in that unlike the remaining ones (to this order) it does not vanish
if all fermions have the same spin projection. This has important consequences
for the transition to the ordered state.

##  5 Phase transition to the ordered state

A short range repulsive interaction, cooperating with the Pauli exclusion rule
can, if its strength is large enough and the gas density is not too low, cause
the transition to the ordered state. The qualitative explanation of this
phenomenon is that because of the Pauli exclusion, the dominant $s$-wave
interaction energy of two fermions with the same spin projection must vanish
(it is positive if the spin projections are different); hence assuming by the
majority of fermions the same spin projection and decreasing thereby the
interaction energy may become energetically more favourable than minimizing
the kinetic energy by having all spin projections equally populated. Since the
effect is due to the competition between the energy of the free system and the
interaction energy, it can occur only if the interaction is sufficiently
strong, which implies that it hardly can be treated perturbatively. In the
mean field approximation, equivalent to taking in (2) into account only the
order $k_{\rm F}a_{0}$ correction, and for spin $s=1/2$ fermions this happens
when $k_{\rm F}a_{0}>\pi/2$ (see e.g. [19, 9, 10]).

In this simplest case the role of the order parameter naturally plays the
polarization $P=(N_{+}-N_{-})/N$ and seeking the minimum of the ground state
energy $E_{\Omega}(P)$ as a function of $P$ is straightforward. In the case of
spin $s$ fermions there are in principle many possible configurations with
non-equal densities of different spin projections but from the heuristic
argument given above it follows that in this case the minimum of the ground
state energy should be in the configuration in which only one spin projection
has the density larger than the remaining ones which all have equal
densities.999The authors of [8] claim to have proved this mathematically but a
cursory look at their argument which relies on the obvious symmetry of the
equations determining the densities of the spin projections different than the
distinguished one shows it is not valid. Of course, in the ordered phase there
are then $g_{s}=2s+1$ degenerate ground states related to one another by the
$SU(g_{s})$ symmetry because in the absence of an external field any of the
$g_{s}$ spin projections can be that one which is populated by the majority of
fermions. Therefore the natural order parameter $P$ defined by setting
(without loss of generality we assume that it is the first spin component
which is populated differently than the remaining ones)

$\displaystyle x_{1}=(1+P)^{1/3}~{}\\!,\phantom{aaa}x_{\sigma}=\left(1-{P\over
g_{s}-1}\right)^{1/3}~{}\\!,\phantom{aaa}\sigma=2,\dots,g_{s}~{}\\!,$ (45)
$\displaystyle
x_{1}^{3}+x_{2}^{3}+\dots+x_{g_{s}}^{3}=g_{s}~{}\\!,\phantom{aaaa}-1\leq P\leq
g_{s}-1~{}\\!,\phantom{aaaaaaa}$ (46)

which corresponds to

$\displaystyle P={g_{s}N_{1}-N\over N}=(g_{s}-1)~{}\\!{N_{1}-N_{2}\over
N}~{}\\!,$ (47)

ceases to have the simple interpretation of the system’s
polarization.101010Obviously, if the spontaneous ordering is approached in the
usual thermodynamic limit by switching off an external magnetic field coupled
to the fermion spin, it will be the largest spin projection in the direction
of the field that will be singled out. In this case the magnetization will be
given by $N\mu P/(g_{s}-1)$, where $\mu$ is the magnetic moment associated
with a single spin. In terms of the factors $x_{1}(P)$ and $x_{2}(P,g_{s})$
given by (45) the curly bracket in the formula (2) takes the form

$\displaystyle f_{s}(P)=x_{1}^{5}+(g_{s}-1)x_{2}^{5}+{20\over
9\pi}~{}\\!k_{\rm F}a_{0}~{}\\!(g_{s}-1)~{}\\!x_{2}^{3}\left[x_{1}^{3}+{1\over
2}~{}\\!(g_{s}-2)~{}\\!x_{2}^{3}\right]$
$\displaystyle+~{}\\!{2\over\pi^{2}}~{}\\!(k_{\rm
F}a_{0})^{2}(g_{s}-1)\left[J_{K}(x_{1},x_{2})+{1\over
2}~{}\\!(g_{s}-2)~{}\\!{44-8\ln 2\over 21}\right].\phantom{a}~{}\\!$ (48)

In the case of spin $1/2$ fermions the phase transition to the ordered state
at $k_{\rm F}a_{0}=\pi/2$ in the mean field approximation, equivalent to
retaining in (2) only the first two terms in the curly bracket, is known to be
continuous. This results however from a numerical coincidence which makes the
second derivative of $E^{\rm MF}_{\Omega}(P)$ with respect to $P$ to vanish at
$P=0$ for $k_{\rm F}a_{0}=\pi/2$, precisely when the nontrivial minimum first
appears. This does not happen if $s>1/2$ and in the mean field approximation
the phase transitions to the ordered states are in all these cases of first
order.111111Taking the effects of the interaction in the first nontrivial
order and restricting oneself to the configurations specified by (45) one can
repeat the steps taken in [9] and obtain the full low temperature profile of
the transition to the ordered phase for any spin $s$. At zero magnetic field
$H$, the equation  $\displaystyle g_{s}\mu H+{4g_{s}\over 3\pi}~{}\\!k_{\rm
F}a_{0}~{}\\!P=(g_{s}-1)\left\\{(1+P)^{2/3}-\left(1-{P\over
g_{s}-1}\right)^{2/3}\right.\phantom{aaaaaaaaaaaaaaaaaaaaaaa}$
$\displaystyle\left.-{\pi^{2}\over 12}\left({k_{\rm
B}T\over\varepsilon^{(s)}_{\rm
F}}\right)^{2}\left[(1+P)^{-2/3}-\left(1-{P\over
g_{s}-1}\right)^{-2/3}\right]+\dots\right\\},$ which determines the parameter
$P$ starts, for $s>1/2$ to have a nontrivial solution before the slopes at
$P=0$ of its both sides equalize - this shows that the transition is of first
order. It can be also easily checked that the critical value of $k_{\rm
F}a_{0}$ at which it occurs decreases with the spin $s$ but rather slowly -
$(k_{\rm F}a_{0})_{\rm cr}^{\rm MF}\approx 1.43$ for $s=3/2$, $(k_{\rm
F}a_{0})_{\rm cr}^{\rm MF}\approx 1.21$ for $s=7/2$ and $(k_{\rm F}a_{0})_{\rm
cr}^{\rm MF}\approx 1.14$ for $s=9/2$. In view of the qualitative explanation
of the reasons for the transition such a behaviour is easy to understand: the
gas kinetic energy is linear in the number $g_{s}$ of the possible spin
projections while the interaction depends on it (in the mean field
approximation) quadratically. One can also observe that only for $s=3/2$ (and,
of course, $s=1/2$) is the spontaneous polarization $P$ not maximal at
$(k_{\rm F}a_{0})_{\rm cr}^{\rm MF}$; for all higher spins it becomes maximal
($P=g_{s}-1$) already at the transition point.

If the corrections of order $(k_{\rm F}a_{0})^{2}$ are included, i.e. the full
formula (48) is used, the transition to the ordered state (at zero
temperature) becomes of first order also for $s={1\over 2}$ (and retains, of
course, this character for $s>{1\over 2}$). The relevant plots can be found in
[12] and in [8] for $s={1\over 2}$ and in [8] for $s\geq{1\over 2}$. The
transition occurs now at $(k_{\rm F}a_{0})_{\rm cr}^{\rm 2nd}\approx 1.054$
for $s=1/2$, $(k_{\rm F}a_{0})_{\rm cr}^{\rm 2nd}\approx 0.954$ for $s=3/2$,
$(k_{\rm F}a_{0})_{\rm cr}^{\rm 2nd}\approx 0.840$ for $s=7/2$ and $(k_{\rm
F}a_{0})_{\rm cr}^{\rm 2nd}\approx 0.804$ for $s=9/2$. Moreover, in this
approximation only for $s=1/2$ the polarization gradually approaches the
maximal value $P=1$ ($P\approx 0.58$ exactly at the critical coupling $k_{\rm
F}a_{0}=(k_{\rm F}a_{0})_{\rm cr}^{\rm 2nd}$); in the remaining cases it jumps
immediately to the maximal possible value $P_{\rm max}=g_{s}-1$.

Figure 4: Plots of the ground state energy of the system of spin $1/2$
fermions as a function of the order parameter $P$ in units of $(k_{\rm
F}^{3}/6\pi^{2})(3/5)(\hbar^{2}k_{\rm F}^{2}/2m_{f})$, i.e. of the function
$f_{s}(P)$. Left: without the order $(k_{\rm F}a_{0})^{3}$ corrections - the
lines from below correspond to $k_{\rm F}a_{0}=1.05404$ (blue), $1.05409$
(red) and 1.105414 (green). Right: with the order $(k_{\rm F}a_{0})^{3}$
corrections \- the lines from below correspond to $k_{\rm F}a_{0}=0.99091$
(blue), $0.99096$ (red) and 0.99101 (green). In both panels the middle line
corresponds to the respective critical value of $k_{\rm F}a_{0}$.

Inclusion of the order $(k_{\rm F}a_{0})^{3}$ corrections to the energy
density treated as a function of the order parameter $P$ amounts to adding to
the function $f_{s}(P)$ given by (48) the term

$\displaystyle{160\over\pi^{3}}~{}\\!(k_{\rm
F}a_{0})^{3}(g_{s}-1)\left[4G^{(1)}_{\rm
fin}(x_{2},x_{1})+4G^{(2)}(x_{2},x_{1})-K^{(1)}(x_{2},x_{1})-K^{(2)}(x_{2},x_{1})\right.\phantom{aaaa}$
$\displaystyle+{1\over 2}(g_{s}-2)\left(4G^{(1)}_{\rm
fin}(x_{2},x_{2})+4G^{(2)}(x_{2},x_{2})-K^{(1)}(x_{2},x_{2})-K^{(2)}(x_{2},x_{2})\right)\phantom{aa}$
(49)
$\displaystyle\left.+(g_{s}-2)\left(\tilde{K}^{(2)}(x_{1};x_{2},x_{2})+2\tilde{K}^{(2)}(x_{2};x_{2},x_{1})+(g_{s}-3)\tilde{K}^{(2)}(x_{2};x_{2},x_{2})\right)\right].$

This has the following consequences. Firstly, the critical value of the
expansion parameter $k_{\rm F}a_{0}$ is further reduced: the transition occurs
now at $(k_{\rm F}a_{0})_{\rm cr}^{\rm 3rd}=0.99096$ for $s=1/2$ (at $(k_{\rm
F}a_{0})_{\rm cr}^{\rm 3rd}=0.74190$ and $0.49445$ for $s=3/2$ and $s=9/2$,
respectively). Secondly, the first order character of the transition, quite
clear without this correction, becomes now much less pronounced - while for
$s=1/2$ the height of the hill separating the minimum at $P=0$ from the one at
$P\neq 0$ remains almost unchanged, greatly reduced is the nonzero value of
the order parameter right at the transition. This is clear from the comparison
of the two panels of Figure 4 showing the shape of the function $f_{s}(P)$ for
three values of $k_{\rm F}a_{0}$ close to the critical one (corresponding to
the approximation used): in the left one the order $(k_{\rm F}a_{0})^{3}$
corrections are not included, in the right one they are. This aspect - the
reduction of the the nonzero value of the order parameter right at the
transition - is even more dramatic in the case of higher spins: without the
third order corrections this value was maximal $P=g_{s}-1$; with these
corrections it is reduced to $\sim 0.12$ (practically independent of the value
of $s$)! Moreover, for $s>1/2$ also the height of the hill separating the
symmetry breaking minimum from the symmetry preserving one at $P=0$ is reduced
by roughly three for $s=3/2$ up to four for $s=9/2$ orders of magnitude. In
the case of spin $3/2$ this is shown in Figure 5. The results for different
values of the spin $s$ are presented in Table 1.

These purely perturbative results strongly point towards the possibility that
the transition to the ordered state is, when all effects associated with the
interactions other than that proportional to the coupling $C_{0}$ in (6) are
ignored (which is equivalent to setting to zero all the higher partial waves
scattering lengths $a_{\ell}$ and all the effective radii $r_{\ell}$), indeed
of the continuous type. Most probably with the inclusion of higher and higher
corrections its character becomes less and less first order and becomes truly
continuous when a resummation of the sort performed in [14] is made but
becomes practically indistinguishable from such already at a finite (not very
high, as our computation shows) order of the ordinary perturbative expansion.

Figure 5: Plots of the ground state energy of the system of spin $3/2$
fermions as a function of the order parameter $P$ in units of $(k_{\rm
F}^{3}/6\pi^{2})(3/5)(\hbar^{2}k_{\rm F}^{2}/2m_{f})$, i.e. of the function
$f_{s}(P)$. Left: without the order $(k_{\rm F}a_{0})^{3}$ corrections \- the
lines from below correspond to $k_{\rm F}a_{0}=0.9532$ (blue), $0.9542$ (red)
and $0.9552$ (green). Right: with the order $(k_{\rm F}a_{0})^{3}$ corrections
\- the lines from below correspond to $k_{\rm F}a_{0}=0.74188$ (blue),
$0.74190$ (red) and $0.74192$ (green). In both panels the middle line
corresponds to the respective critical value of $k_{\rm F}a_{0}$.

Finally we can briefly discuss effects induced by the potential presence of
nonnegligibly small $s$-wave effective radius $r_{0}$ and/or $p$-wave
scattering length $a_{1}$. They are introduced through the two lower length
dimension operators in (6) which are given explicitly by the formulae (41) and
(42). In the configuration of densities assumed in this Section these
interactions add to the function $f_{s}(P)$ the terms

$\displaystyle{1\over 6\pi}~{}\\!k_{\rm
F}^{3}a_{0}^{2}r_{0}~{}\\!(g_{s}-1)\left[x_{1}^{3}x_{2}^{3}(x_{1}^{2}+x_{2}^{2})+(g_{s}-2)~{}\\!x_{2}^{8}\right]\phantom{aaaaaaaaa}~{}$
$\displaystyle+{2\over 3\pi}~{}\\!(k_{\rm F}a_{1})^{3}\left[x_{1}^{8}+{1\over
2}g_{s}(g_{s}-1)~{}\\!x_{2}^{8}+{1\over
2}(g_{s}-1)~{}\\!x_{1}^{3}x_{2}^{3}(x_{1}^{2}+x_{2}^{2})\right].$ (50)

If the effects of these operators are subleading compared to the effects
caused by the order $k_{F}^{3}$ corrections of the operator proportional to
$C_{0}$ (that is, when that $|r_{0}|,|a_{1}|\ll a_{0}$) they modify the
general picture described above only slightly. For instance in the case of
$s=1/2$ the presence of $r_{0}\approx a_{0}/10$ decreases by less than one
percent the critical value of the parameter $k_{\rm F}a_{0}$, while
$r_{0}\approx-a_{0}/10$ increases it by the same amount in agreements with the
observation made in [14]. In both cases the character of the phase transition
is not appreciably changed. Somewhat surprising is, however, the observation
that for larger spin values this patern is reversed: negative $r_{0}$
decreases slightly the critical value of $k_{\rm F}a_{0}$, while positive
$r_{0}$ increases it. The reason for this is that while for $s=1/2$ the first
line of (50) is a function of $P$ monotonically decreasing to zero (the
contribution to the energy density of the interaction proportional to $C_{2}$
is by the Pauli exclusion principle also bound to vanish at $P=1$), for higher
spins it a is slightly increasing function, precisely in the range in which
the new minimum of the energy density forms.

In general, the correction induced by the operator (42) is larger than the one
induced by (41). Moreover in the perturbation series for $E_{\Omega}/V$ it is
the first correction which is not bound to vanish at $P=1$ by the Pauli
exclusion principle. Nevertheless, if one assumes that $|a_{1}|$ is of the
same order as $|r_{0}|$ and much smaller than $a_{0}$, its contribution is not
larger than that of (41) simply because it is proportional to $a_{1}^{3}$. In
the configuration of the densities considered here the expression in the
second line of (50) is for all values of the spin $s$ a monotonically
increasing function. For this reason the contribution increases the value of
$(k_{\rm F}a_{0})_{\rm cr}$ for $a_{1}>0$ and lowers if for $a_{1}<0$ by
roughly one percent so long as
$|a_{1}|\stackrel{{\scriptstyle<}}{{{}_{\sim}}}a_{0}/10$ without affecting
significantly the character of the transition. Thus, it seems likely, that if
all the scattering lengths $a_{\ell}$, $\ell=1,2,\dots$, all the radii
$r_{\ell}$, etc. are subleading with respect to $a_{0}$, the continuos
character of the transition will emerge after all corrections are taken into
account (being, from the practical point of view indistinguishable from
continuos already starting from some fixed order of the perturbative
expansion).

$\begin{array}[]{|c||c|c||c|c||c|c|}\hline\cr s&(k_{F}a_{0})_{\rm cr}^{\rm
MF}&P^{\rm MF}_{\rm cr}&(k_{F}a_{0})_{\rm cr}^{\rm 2nd}&P^{\rm 2nd}_{\rm
cr}&(k_{F}a_{0})_{\rm cr}^{\rm 3rd}&P^{\rm 3rd}_{\rm cr}\\\ \hline\cr
1/2&\pi/2&0&1.0540&0.58&0.99097&0.11\\\ 3/2&1.4298&2.78&0.9543&P_{\rm
max}&0.74190&0.14\\\ 5/2&1.3016&P_{\rm max}&0.8877&P_{\rm
max}&0.62560&0.125\\\ 7/2&1.2118&P_{\rm max}&0.8405&P_{\rm
max}&0.54976&0.115\\\ 9/2&1.1439&P_{\rm max}&0.8040&P_{\rm
max}&0.49445&0.105\\\ 11/2&1.0901&P_{\rm max}&0.7753&P_{\rm
max}&0.45164&0.10\\\ \hline\cr\end{array}$

Table 1: Characteristics of the transitions to the ordered state at $T=0$ for
different values of the spin $s$ of fermions.

##  6 Conclusions

Computing the ground-state energy of a finite density system of fermions
interacting through a binary spin independent repulsive interaction is a
classic problem of many-body quantum mechanics. The modern effective field
theory approach greatly simplifies this task. In particular it allowed to
complete the computation of the fourth order (i.e. proportional to $k_{\rm
F}^{4}$) corrections in the case of spin $s$ fermions and equal densities of
different spin projections [5] and to extend to the third order the
computation of energy of spin $1/2$ fermions for an arbitrary value of the
system’s polarization [16]. Here we have extended the latter result to the
case of spin $s$ fermions and arbitrary densities of different spin
projections. The derived formulae allow to compute the ground state energy
semi-analytically exploiting only the built-in routines of the Mathematica
package.

We have used this result to discuss two issues. Firstly we have checked
numerically how the magnitude of the third order term included in the
resummation of an infinite subset of corrections done in [14] compares with
the magnitudes of the rejected terms and found that it is not obviously
dominant. Secondly, we have discussed the impact the third order corrections
have on the characteristics of the system’s transition to the ordered phase at
zero temperature. We have found that already the third order corrections tend
to erase the first order character of this transition independently of the
value $s$ of the spin of fermions. Although our observations are made on the
basis of the perturbative expansion used in the regime which is probably
beyond the domain of its applicability (the comparison of the results of the
Quantum Monte Carlo computations of the ground state energy of the unpolarized
gas of spin $1/2$ fermions with the perturbative coputation reveals [15] that
the expansion is reliable up to $k_{\rm
F}a_{0}\stackrel{{\scriptstyle<}}{{{}_{\sim}}}0.6$; in the case of $s>1/2$ the
value of the expansion parameter $k_{\rm F}a_{0}$ at which the transition to
the ordered state occurs is smaller but since at the same time the magnitude
of the successive terms of the perturbation series are increased by the
growing powers of the factor $g_{s}$, the limit of the reliability of the
expansion probably also decreases with $s$), they nevertheless seem to lend
some support to the claim made in [14] that the transition is continuous
rather than first order. It seems that the picture which emerges is quite
sensible: since continuous transition are associated with fluctuations at all
length scales, at any finite order of the perturbative expansion such a
transition should look like a first order one (the continuous transition in
the system of spin $1/2$ fermions obtained in the mean field approximation is
from this perspective a mere numerical coincidence), though its first order
character may rapidly disappear with increased order of the approximation (so
that it quickly may become indistinguishable from a continuous one from the
practical point of view); the truly continuous character of a transition can
really be revealed only by a (partial) resummation of of all orders
contributions, of the sort proposed in [14].

Acknowledgments. We would like to thank Dr Krzysztof Jachymski for a
discussion.

## References

  * [1] G.-B. Jo et al., Science 325 (2009) 1521; D. Pekker et al., Phys. Rev. Lett. 106 (2011), 050402; Y.-R. Lee et al., Phys. Rev. A85 (2012), 063615; C. Sanner et al., Phys. Rev. Lett. 108 (2012), 240404.
  * [2] H.-W. Hammer, R. J. Furnstahl, Nucl. Phys. A678, 277 (2000); arXiv:0004043 [nucl-th].
  * [3] R.J. Furnstahl and H.-W. Hammer, Phys. Lett. B531, 203 (2002); H.-W. Hammer, S. König and U. van Kolck, Rev. of Mod. Phys. 92 (2020), 025004; R.J. Furnstahl, V. Steele and N. Tirfessa, Nucl. Phys. A671 (2000), 396.
  * [4] see e.g. Proceedings of the Joint Caltech/INT Workshop Nuclear Physics with Effective Field Theory, ed. R. Seki, U. van Kolck and M. Savage (World Scientific, 1998); Proceedings of the INT Workshop Nuclear Physics with Effective Field Theory II, ed. P.F. Bedaque, M. Savage, R. Seki and U. van Kolck (World Scientific, 2000).
  * [5] C. Wellenhofer, C. Drischler and A. Schwenk, Phys. Lett. B802 135247 (2020) arXiv:1812.08444 [nucl-th]; Phys. Rev. C104, 014003 (2021) arXiv:2102.05966 [cond-mat.quant-gas].
  * [6] S. Kanno, Prog. Theor. Phys. 44, 813 (1970).
  * [7] P.H. Chankowski and J. Wojtkiewicz, Phys. Rev. B104 144425 (2021) arXiv:2108.00793 [cond-mat.quant-gas].
  * [8] J. Pera, J. Casulleras and J. Boronat, arXiv:2205.13837 [cond-mat.quant-gas].
  * [9] K. Huang, Statistical Mechanics, John Willey and Sons, Inc., New York 1963.
  * [10] R.K. Pathria, Statistical Mechanics, Pergamon Press, Oxford, 1972.
  * [11] R.A. Duine and A.H. MacDonald, Phys. Rev. Lett. 95, 230403 (2005).
  * [12] P.H. Chankowski and J. Wojtkiewicz, Acta Phys. Pol. B53 (2022), 1-A3, arXiv:2111.14927 [cond-mat.quant-gas].
  * [13] D. Belitz, T.R. Kirkpatrick and T. Vojta, Phys. Rev. Lett. 82, 4707 (1999).
  * [14] L. He and X.-G. Huang, Phys. Rev. A85, 043624 (2012), arXiv:1106.1345; L. He, Ann. of Phys. 351, 477 (2014), arXiv:1405.3338, Phys. Rev. A90, 053633, arXiv:1405.5242; L. He, X.-J. Liu, X.-G. Huang and H. Hui, Phys. Rev. A93, 063629 (2016), arXiv:1412.2412.
  * [15] S. Pilati, G. Bertaina, S. Giorgini and M. Troyer, Phys. Rev. Lett. 105, 030405 (2010), arXiv:1004/1169 [cond-mat.quant-gas].
  * [16] P.H. Chankowski, J. Wojtkiewicz and R. Bakhshizada, Acta. Phys. Pol. B53 (2022), 9-A4, arXiv:2206.05076 [cond-mat.quant-gas].
  * [17] S. Weinberg, The Quantum Theory of Fields, Vol. II, The Press Syndicate of the University of Cambridge, 1996.
  * [18] A.L. Fetter and J.D. Walecka, Quantum Theory of Many Particle Systems, McGraw Hill, 1971.
  * [19] E. Stoner, The London, Edinburgh and Dublin Phil. Mag. and Journal of Science 15 (1933), 1018.
  * [20] R.J. Furnstahl, H.-W. Hammer and N. Tirfessa, Nucl. Phys. A689, 846 (2001).


* (22) J. Gonzalez, F. Guinea and M. A. H. Vozmediano, Marginal-Fermi-liquid behavior from two-dimensional Coulomb interaction, Phys. Rev. B 59 (1999), 2474-2477 http://doi.org/10.1103/PhysRevB.59.R2474
  * (23) J. Gonzalez, F. Guinea and M. A. H. Vozmediano, Electron-electron interactions in graphene sheets Phys. Rev. B 63 (2001) 134421 http://doi.org/10.1103/PhysRevB.63.134421
  * (24) M. A. H. Vozmediano, Renormalization group aspects of graphene, Phil. Trans. Roy. Soc. Lond. A 369 (2011), 2625-2642 http://doi.org/10.1098/rsta.2010.0383
  * (25) L. Fernández, V. Alves, L. O. Nascimento, F. Peña, M. Gomes and E. C. Marino, Phys. Rev. D 102 (2020) no.1, 016020 doi:10.1103/PhysRevD.102.016020 [arXiv:2002.10027 [hep-th]].
  * (26) E. V. Gorbar, V. P. Gusynin and V. A. Miransky, Dynamical chiral symmetry breaking on a brane in reduced QED, Phys. Rev. D 64, (2001) 105028 https://doi.org/10.1103/PhysRevD.64.105028
  * (27) S. Teber, Electromagnetic current correlations in reduced quantum electrodynamics, Phys. Rev. D 86 (2012) 025005, https://doi.org/10.1103/PhysRevD.86.025005
  * (28) A. V. Kotikov and S. Teber, Two-loop fermion self-energy in reduced quantum electrodynamics and application to the ultrarelativistic limit of graphene, Phys. Rev. D 89, no. 6, (2014) 065038 https://doi.org/10.1103/PhysRevD.89.065038
  * (29) S. Teber, Two-loop fermion self-energy and propagator in reduced QED3,2, Phys. Rev. D 89, no. 6, (2014) 067702, https://doi.org/10.1103/PhysRevD.89.067702
  * (30) S. Teber and A. Kotikov, Field theoretic renormalization study of reduced quantum electrodynamics and applications to the ultrarelativistic limit of Dirac liquids, Phys. Rev. D 97, no.7, (2018) 074004. https://doi.org/10.1103/PhysRevD.97.074004
  * (31) V. Alves, W. S. Elias, L. O. Nascimento, V. Juric̆ić and F. Peña, “Chiral symmetry breaking in the pseudo-quantum electrodynamics,” Phys. Rev. D 87, no.12, (2013) 125002. https://doi.org/10.1103/PhysRevD.87.125002
  * (32) A. Kotikov and S. Teber, Critical behaviour of reduced QED4,3 and dynamical fermion gap generation in graphene, Phys. Rev. D 94, no.11, (2016) 114010. https://doi.org/10.1103/PhysRevD.94.114010
  * (33) L. O. Nascimento, V. Alves, F. Peña, C. M. Smith and E. Marino, Chiral-symmetry breaking in pseudoquantum electrodynamics at finite temperature, Phys. Rev. D 92, (2015) 025018. https://doi.org/10.1103/PhysRevD.92.025018
  * (34) J. Báez-Cuevas, A. Raya and J. C. Rojas, Chiral symmetry restoration in reduced QED at finite temperature in the supercritical coupling regime, Phys. Rev. D 102 (2020), 056020. https://doi.org/10.1103/PhysRevD.102.056020
  * (35) M. E. Carrington, Effect of a Chern-Simons term on dynamical gap generation in graphene, Phys. Rev. B 99 (2019) 115432. https://doi.org/10.1103/PhysRevB.99.115432
  * (36) J. A. C. Olivares, L. Albino, A. J. Mizher and A. Raya, Influence of a Chern-Simons term in the dynamical fermion masses in reduced or pseudo QED, Phys. Rev. D 102 (2020) no.9, 096023 https://doi.org/10.1103/PhysRevD.102.096023 [arXiv:2009.08484 [hep-ph]].
  * (37) G. C. Magalhães, V. Alves, E. C. Marino and L. O. Nascimento, Pseudo Quantum Electrodynamics and Chern-Simons theory Coupled to Two-dimensional Electrons, Phys. Rev. D 101 (2020) no.11, 116005 https://doi.org/10.1103/PhysRevD.101.116005.
  * (38) M. E. Carrington, A. R. Frey and B. A. Meggison, Effect of anisotropy on phase transitions in graphene, Phys. Rev. B 102 (2020) no.12, 125427. https://doi.org/10.1103/PhysRevB.102.125427
  * (39) P. I. C. Caneda and G. Menezes, Reduced quantum electrodynamics in curved space , Phys. Rev. D103, (2021) 065010. https://doi.org/10.1103/PhysRevD.103.065010
  * (40) L. D. Landau and I. M. Khalatnikov, The gauge transformation of the Green function for charged particles, Sov. Phys. JETP 2 (1956) 69 [Zh. Eksp. Teor. Fiz. 29 (1955) 89].
  * (41) E. S. Fradkin, Concerning some general relations of quantum electrodynamics, Zh. Eksp. Teor. Fiz. 29, 258 (1955) [Sov. Phys. JETP 2, 361 (1956)].
  * (42) A. Bashir, Non-perturbative Fermion Propagator for the Massless Quenched QED3 Phys. Lett. B491, (2000) 280. https://doi.org/10.1016/S0370-2693(00)01043-1
  * (43) A. Bashir and A. Raya, Landau-Khalatnikov-Fradkin transformations and the fermion propagator in quantum electrodynamics, Phys. Rev. D 66, (2002) 105005. https://doi.org/10.1103/PhysRevD.66.105005
  * (44) A. Bashir and A. Raya, Fermion propagator in quenched QED3 in the light of the Landau-Khalatnikov-Fradkin transformation, Nucl. Phys. Proc. Suppl. 141, (2005) 259 https://doi.org/10.1016/j.nuclphysbps.2004.12.039
  * (45) A. Ahmad, J. J. Cobos-Martínez, Y. Concha-Sánchez, and A. Raya, Landau-Khalatnikov-Fradkin transformations in reduced quantum electrodynamics, Phys. Rev. D 93, (2016) 094035. https://doi.org/10.1103/PhysRevD.93.094035
  * (46) D. Dudal, A. J. Mizher and P. Pais, Exact quantum scale invariance of three-dimensional reduced QED theories, Phys. Rev. D 99 (2019) no.4, 045017- https://doi.org/PhysRevD.99.045017
  * (47) C. G. Bollini and J. J. Giambiagi, Arbitrary powers of D’Alembertians and the Huygens’ principle, J. Math. Phys. 34 (1993), 610-621 https://doi.org/10.1063/1.530263
  * (48) A. Coste and M. Luscher, Parity Anomaly and Fermion Boson Transmutation in Three-dimensional Lattice QED, Nucl. Phys. B 323 (1989), 631-659 https://doi.org/10.1016/0550-3213(89)90127-2
  * (49) C. D. Roberts and A. G. Williams Dyson-Schwinger Equations and the Application to Hadronic Physics Prog. Part. Nucl. Phys. 33 (1994) 477. https://doi.org/10.1016/0146-6410(94)90049-3
  * (50) M. R. Pennington Swimming with Quarks, J. Phys. Conf. Ser. 18 (2005) 1. https://doi.org/10.1088/1742-6596/18/1/001
  * (51) J. C. Ward, An Identity in Quantum Electrodynamics, Phys. Rev. 78 (1950) 182 https://doi.org/10.1103/PhysRev.78.182
  * (52) H. S. Green, A Pre-renormalized quantum electrodynamics, Proc. Phys. Soc. (London) A66, (1953) 873. https://doi.org/10.1088/0370-1298/66/10/303
  * (53) Y. Takahashi, On the generalized ward identity, Nuovo Cimento 6, (1957) 371. https://doi.org/10.1007/BF02832514
  * (54) Y. Takahashi, Canonical quantization and generalized Ward relations: Foundation of nonperturbative approach, Print-85-0421 (Alberta), (1985).
  * (55) V. Miransky and K. Yamawaki, Conformal phase transition in gauge theories Phys. Rev. D 55, (1997) 5051. https://doi.org/10.1103/PhysRevD.55.5051.
  * (56) A. Bashir, A. Huet, and A. Raya Gauge dependence of mass and condensate in chirally asymmetric phase of quenched three-dimensional QED  Phys. Rev. D 66, (2002) 025029. https://doi.org/10.1103/PhysRevD.66.025029.
  * (57) K. Raya, A. Bashir, S. Hernández-Ortiz, A. Raya, and C. D. Roberts, Multiple solutions for the fermion mass function in QED3, Phys. Rev. D88, (2013) 096003. https://doi.org/10.1103/PhysRevD.88.096003
  * (58) T. Appelquist, M. J. Bowick, E. Cohler, and L. C. R. Wijewardhana Chiral-symmetry breaking in 2+1 dimensions Phys. Rev. Lett. 55 (1985) 1715. https://doi.org/10.1103/PhysRevLett.55.1715.
  * (59) T. W. Appelquist, M. Bowick, D. Karabali, and L. C. R. Wijewardhana Spontaneous chiral-symmetry breaking in three-dimensional QED Phys. Rev. D 33 (1986) 3704. https://doi.org/10.1103/PhysRevD.33.3704.
  * (60) T. Appelquist, D. Nash, and L. C. R. Wijewardhana, Critical Behavior in (2+1)-Dimensional QED  Phys. Rev. Lett. 60 (1988) 2575. https://doi.org/10.1103/PhysRevLett.60.2575
  * (61) A. Bashir, A. Raya, S. Sánchez-Madrigal and C.D. Roberts, Gauge Invariance of a Critical Number of Flavours in QED3. Few-Body Syst 46 (2009) 229. https://doi.org/10.1007/s00601-009-0069-9.
  * (62) A. V. Kotikov and S. Teber, Critical Behavior of (2+1)-Dimensional QED: $1/N$ Expansion, Particles 2020, 3, 345-354. https://doi.org/10.3390/particles3020026
  * (63) A. V. Kotikov and Sofian Teber, Critical behaviour of (2+1)-dimensional QED: 1/N-corrections, EPJ Web of Conferences 138 (2017), 06005. https://doi.org/10.1051/epjconf/201713806005
  * (64) A. V. Kotikov and S. Teber, Critical behavior of (2+1)-dimensional QED: $1/N_{f}$ corrections in an arbitrary nonlocal gauge, Phys. Rev. D94, (2016) 114011. https://doi.org/10.1103/PhysRevD.94.114011
  * (65) A. V. Kotikov, V.I. Shilin and S. Teber, Critical behavior of (2+1)-dimensional QED: $1/N_{f}$ corrections in the Landau gauge, Phys.Rev.D 94 (2016) 5, 056009. https://doi.org/10.1103/PhysRevD.94.056009
  * (66) J. A. Casimiro-Olivares Generación dinámica de Masas en R(P)QED M. Sc. Thesis (In Spanish) UMSNH (2019).
  * (67) J. S. Ball and T.-W. Chiu Analytic properties of the vertex function in gauge theories. II Phys. Rev. D 22 (1981) 2550. https://doi.org/10.1103/PhysRevD.22.2550
  * (68) S. Metayer and S. Teber Two-loop mass anomalous dimension in reduced quantum electrodynamics and application to dynamical fermion mass generation e-Print: 2107.07807 [hep-th] https://arxiv.org/abs/2107.07807
  * (69) A. Ahmad, A. Ayala, A. Bashir, E. Gutiérrez and A. Raya, QCD Phase Diagram and the Constant Mass Approximation, J. Phys. Conf. Ser. 651, no.1, (2015) 012018 https://doi.org/10.1088/1742-6596/651/1/012018
  * (70) A. Khare, Fractional Statistics and Quantum Theory (World Scientific, Singapore, 2005), 2nd ed.
  * (71) K. I. Kondo and P. Maris, First-order phase transition in three-dimensional QED with Chern-Simons term Phys. Rev. Lett. 74(1995) 18. https://doi.org/10.1103/PhysRevLett.74.18
  * (72) K. I. Kondo and P. Maris, Spontaneous chiral symmetry breaking in three-dimensional QED with a Chern-Simons term, Phys. Rev. D 52, (1995) 1212. https://doi.org/10.1103/PhysRevD.52.1212
  * (73) C. P. Hofmann, A. Raya and Saúl Sánchez-Madrigal, Confinement in Maxwell-Chern-Simons planar quantum electrodynamics and the $1/N$ approximation Phys. Rev. D 82, (2010) 096011. https://doi.org/10.1103/PhysRevD.82.096011
  * (74) S. Sánchez Madrigal, C. P. Hofmann and A Raya, Dynamical Mass Generation and Confinement in Maxwell-Chern-Simons Planar Quantum Electrodynamics J. Phys.Conf. Ser. 287 (2011) 012028. https://doi.org/10.1088/1742-6596/287/1/012028
  * (75) R. F. Ozela, V. Alves, G. C. Magalhães and L. O. Nascimento, Effects of the pseudo-Chern-Simons action for strongly correlated electrons in a plane Phys. Rev. D 105 (2022) no.5, 056004 https://doi.org/10.1103/PhysRevD.105.056004
  * (76) E. C. Marino and Leandro O. Nascimento and Van SÃ©rgio Alves and N. Menezes and C. Morais Smith, Quantum-electrodynamical approach to the exciton spectrum in transition-metal dichalcogenides, 2D Materials, 5, 4 (2018). https://doi.org/10.1088%2F2053-1583%2Faacc3f
  * (77) D. Dudal, A. J. Mizher and P. Pais, Remarks on the Chern-Simons photon term in the QED description of graphene Phys. Rev. D 98 (2018) 065008. https://doi.org/10.1103/PhysRevD.98.065008
  * (78) G. García Naumis, Electronic properties of 2D materials and its heterostructures: a minimal review , Rev. Mex. Fis. 67 (5), (2021) 1. https://doi.org/10.31349/RevMexFis.67.050102
  * (79) K.-I. Kondo, Transverse Ward-Takahashi Identity, Anomaly and Schwinger-Dyson Equation, Int. J. Mod. Phys. A12, (1997) 5651. https://doi.org/10.1142/S0217751X97002978
  * (80) H.-X. He, F. Khanna and Y. Takahashi, Transverse Ward-Takahashi identity for the fermion-boson vertex in gauge theories, Phys. Lett. B480, (2000) 222. https://doi.org/10.1016/S0370-2693(00)00353-1
  * (81) H.-X. He, Transverse vector vertex function and transverse Ward-Takahashi relations in QED, Commun. Theor. Phys. 46, (2006) 109. https://doi.org/10.1088/0253-6102/46/1/025
  * (82) H.-X. He, Transverse ward-takahashi relation for the fermion-boson vertex function in four-dimensional abelian gauge theory, Int. J. Mod. Phys. A22, (2007) 2119. https://doi.org/10.1142/S0217751X07036257
  * (83) K. Johnson and B. Zumino, Gauge Dependence of the Wave-Function Renormalization Constant in Quantum Electrodynamics, Phys. Rev. Lett.3, (1959) 351. https://doi.org/10.1103/PhysRevLett.3.351
  * (84) B. Zumino, Gauge properties of propagators in quantum electrodynamics, J. Math. Phys. 1, (1960) 1. https://doi.org/10.1063/1.1703632
  * (85) S. Okubo, The gauge properties of Green’s functions in quantum electrodynamics, Nuovo Cim. 15, (1960) 949. https://doi.org/10.1007/BF02860201
  * (86) I. Bialynicki-Birula, On the Gauge Properties of Green’s Functions, Nuovo Cim. 17, (1960) 951. https://doi.org/10.1007/BF02732140
  * (87) H. Sonoda, On the gauge parameter dependence of QED, Phys. Lett. B 499 (2001) 253. https://doi.org/10.1016/S0370-2693(01)00030-2
  * (88) L. A. Fernandez-Rangel, A. Bashir, L. X. Gutierrez-Guerrero and Y. Concha-Sanchez, Constructing Scalar-Photon Three Point Vertex in Massless Quenched Scalar QED, Phys. Rev. D 93, no. 6, (2016) 065022. https://doi.org/10.1103/PhysRevD.93.065022
  * (89) N. Ahmadiniaz, A. Bashir and C. Schubert, Multiphoton amplitudes and generalized LKF transformation in Scalar QED, Phys. Rev. D 93 (2016) 4, 045023 https://doi.org/10.1103/PhysRevD.93.045023
  * (90) A. Kizilersu and M. R. Pennington, Building the Full Fermion-Photon Vertex of QED by Imposing Multiplicative Renormalizability of the Schwinger-Dyson Equations for the Fermion and Photon Propagators, Phys. Rev. D 79, (2009) 125020. https://doi.org/10.1103/PhysRevD.79.125020
  * (91) A. Bashir and A. Raya, Constructing the fermion boson vertex in QED(3), Phys. Rev. D 64, (2001) 105001. https://doi.org/10.1103/PhysRevD.64.105001
  * (92) A.V. Kotikov and S. Teber Landau-Khalatnikov-Fradkin transformation and the mystery of even $\zeta$-values in Euclidean massless correlators Phys. Rev. D 100, (2019) 105017. https://doi.org/10.1103/PhysRevD.100.105017
  * (93) A. James, A. V. Kotikov, and S. Teber Landau-Khalatnikov-Fradkin transformation of the fermion propagator in massless reduced QED  Phys. Rev. D 101 (2020) 045011. https://doi.org/10.1103/PhysRevD.101.045011
  * (94) A. V. Kotikov and S. Teber, Landau-Khalatnikov-Fradkin Transformation and Hatted $\zeta$-Values Phys. Part.Nucl. 51, (2020) 562. https://doi.org/10.1134/S1063779620040425
  * (95) V. P. Gusynin, A. V. Kotikov, and S. Teber, Landau-Khalatnikov-Fradkin transformation in three-dimensional quenched QED Phys. Rev. D 102, (2020) 025013 https://doi.org/10.1103/PhysRevD.102.025013
  * (96) D. Dudal, F. Matusalem, A. J. Mizher, A. R. Rocha and C. Villavicencio, Half-integer anomalous currents in 2D materials from a QFT viewpoint Sci. Rep. 12 (2022) no.1, 5439 https://doi.org/10.1038/s41598-022-09483-4