Renormalization of pionless effective field theory in the A-body sector

Abstract

Current models of inter-nucleon interactions are built within the frame of Effective Field Theories (EFTs). Contrary to traditional nuclear potentials, EFT interactions require a renormalization of their parameters in order to derive meaningful estimations of the observables. In this paper, a renormalization procedure is designed in connection with many-body approximations applicable to large-$\mathrm {A}$ systems and formulated within the frame of many-body perturbation theory. The procedure is shown to generate counterterms that are independent of the targeted $\mathrm {A}$-body sector. As an example, the procedure is applied to the random phase approximation. This work constitutes one step towards the design of a practical EFT for many-body systems.

1 Introduction

The problem of describing accurately and systematically nuclear systems from their $\mathrm {A}$ constitutive protons and neutrons is now almost a century old. Starting with the discovery of the neutron in the 1930s, the nuclear many-body problem remains a great challenge in spite of the development of a large portfolio of theoretical methods (see Ref. [1] for a recent historical recapitulation).

Almost 30 years ago, the seminal papers of Weinberg [2, 3] and Rho [4] initiated a shift of paradigm. Rather than keeping on refining interaction potentials between nucleons by fine-tuning their short-range behavior, the focus has shifted towards the development and the study of various effective field theories (see Ref. [5] for a recent and pedagogical introduction to EFTs in the nuclear physics context, Ref. [6] for a global review, Refs. [7, 8] for reviews dedicated to chiral forces and their applications, and Ref [9] for a specific review on Halo EFT). Within this theoretical framework, the computation of nuclear observables has to face the demand of renormalizability.

In this article, the example of pionless effective field theory ($\pi \!\!{/}\,$EFT) [10] at leading order (LO) is used to demonstrate how the procedure of renormalization of the nuclear Hamiltonian affects the calculation of many-body nuclear observables. One key aspect in this respect relates to the fact that, while the Schrödinger equation can be solved exactly in few-body systems (as required by the power counting at LO), it is and will remain practically impossible to do so in large-$\mathrm {A}$ sectors. This feature must be explicitly considered when applying EFT to arbitrary nuclear systems and requiring an order-by-order renormalizability. In Sect. 2, the renormalizability issue and its consequences on $\mathrm {A}$-body calculations are introduced. A renormalization tailored to a given many-body approximation is advocated as a crucial step toward the construction of EFTs for nuclear many-body systems at low energy. The problem is then addressed in a generic way for many-body approximations formulated within the frame of many-body perturbation theory (MBPT). In Sect. 3, some essential theoretical tools, i.e. Weinberg’s asymptotic theorem [11] and the Bogoliubov–Parasiuk–Hepp–Zimmermann (BPHZ) theorem [12,13,14,15,16,17] are briefly recalled. In Sect. 4, a general procedure to renormalize the Hamiltonian is derived for a given many-body approximation and independently of the $\mathrm {A}$-body sector of interest. In Sect. 5, the procedure is applied to the Random Phase Approximation (RPA). Eventually, extensions of this approach and the potential impact on the shape of future many-body approximations are discussed in Sect. 6.

2 Generalities

In this section, the problem of renormalization of the nuclear Hamiltonian is introduced. First, emphasis is put on why it is becoming an important problem to be explicitly addressed in many-body calculations. Second, the formalism to study renormalization in a given many-body approximation is set up.

2.1 The nuclear many-body problem

Traditional many-body approaches rely on the given of a nuclear Hamiltonian H and aim at computing its exact eigenstates $\vert \Psi ^\mathrm {A}_m \rangle$ and their associated eigenvalues $E^\mathrm {A}_m$ in all $\mathrm {A}$-body sectors of interest, i.e. the goal is to solve the Schrödinger equation

$$\begin{aligned} H \vert \Psi ^\mathrm {A}_m \rangle = E^\mathrm {A}_m \vert \Psi ^\mathrm {A}_m \rangle \ \end{aligned}$$

(1)

where m labels the solutions, to the best accuracy possible.

In this context, the Hamiltonian can be modeled in various ways. The current paradigm consists of building H within the framework of chiral effective field theory ($\chi$EFT) [2,3,4, 18, 19] such that it takes the form of a series

$$\begin{aligned} H_{\chi } \equiv H^{\text {LO}}_{\chi } + H^{\text {SL}O}_{\chi } = H^{\text {LO}}_{\chi } + \sum ^{\infty }_{p=1} H^{{N^{p}LO}}_{\chi }\ \end{aligned}$$

(2)

where the leading-order (LO) and the sub-leading orders (SLOs) are organized according to a set of power-counting (PC) rules. First to be proposed historically, Weinberg’s power counting [2, 3, 19] happens to fit traditional many-body calculations, i.e. independently of the order at which SLOs are truncated, Eq. (1) is meant to be solved exactly to access observables such as $E^\mathrm {A}_m$.

However, Weinberg’s PC has since been shown to violate the demand that the EFT is (order-by-order) renormalizable [20]. A new PC addressing this problem and implying a perturbative treatment of pion exchanges has been developed [21]. Such a PC happens to fail in predicting accurately phase shifts at momenta near the pion mass in certain partial waves [22], which was interpreted as the consequence of the perturbative treatment of pion exchanges. Consequently, a modified PC relying on the non-perturbative treatment of pion exchanges was motivated by the study of an expansion around the chiral limit [23]. Based on such a PC, a careful numerical study [24] of the cut-off dependence of the phase shifts in different partial waves revealed the necessity to promote to LO one term per partial wave where the tensor force is attractive and singular. To avoid the promotion of an infinite number of terms, a refinement of the PC taking into account angular-momentum suppression was thus advocated in Ref. [24]. In this case, only a finite number of low partial waves together with their promoted counterterm are kept at LO. A careful analytical study of such promotion/demotion of terms due to the singular tensor force can be found in Ref. [25]. For a general review of these evolutions of the PC in $\chi$EFT, see for example Ref. [26]. In addition to modifying the order at which certain contributions enter the Hamiltonian, such new PCs generally stipulate that, while LO is to be solved exactly according to traditional many-body calculations, SLOs must be computed in perturbation with respect to the LO solution. It happens that the same computational scheme underlies the PC at play in $\pi \!\!{/}\,$EFT [6].

Generally speaking, the order-by-order renormalizability means that, at each order, observables can display a dependence on the regularization, e.g. in the form of a dependence on the ultraviolet cut-off that is not larger than the intrinsic uncertainty carried at the working order. To achieve that, observables are enforced to depend, at each order, on inverse powers of the cut-off and reach a finite limit as the cut-off is sent to infinity. In this context, $\pi \!\!{/}\,$EFT has been shown analytically to satisfy renormalizability at LO up to three-body systems [27] and numerically in the case of four-body systems [28].

2.2 Renormalization and many-body approximations

To this day, solving exactly Eq. (1) for $H^{\text {LO}}_{\pi \!\!{/}\,}$ remains in general numerically intractable for $\mathrm {A} \gg 10$. Consequently, one must design an additional expansion and truncation when attempting to solve

$$\begin{aligned} H^{\text {LO}}_{\pi \!\!{/}\,} \vert \Psi ^\mathrm {A}_m \rangle ^{\text {(LO)}} = E^{\mathrm {A}\text {(LO)}}_m \vert \Psi ^\mathrm {A}_m \rangle ^{\text {(LO)}} \ . \end{aligned}$$

(3)

Typical truncations applicable to $\mathrm {A}$-body systems with $\mathrm {A} \gg 10$ are nowadays implemented on the basis of a non-perturbative self-consistent Green’s function (SCGF) [29,30,31], coupled cluster (CC) [32, 33] and in-medium similarity renormalization group (IM-SRG) [34, 35] methods but also on the basis of MBPT [36,37,38].

Traditionally, and in agreement with power-counting rules of $\pi \!\!{/}\,$EFT, $H^{\text {LO}}_{\pi \!\!{/}\,}$ is renormalized on the basis of an all-order calculation in two- and three-body sectors. This happens to be technically feasible. Given, however, the impossibility to generate exact calculations in large $\mathrm {A}$ sectors, investigating how observables are impacted by the mandatory approximations appears to be a necessary task to validate the practical use of the current form of $\pi \!\!{/}\,$EFT across a large fraction of the nuclear chart. More specifically, one must question to which extent the approximate solving of Eq. (3), on the basis of a previously renormalized potential via an exact calculation in two- and three-body sectors, compromises renormalization invariance. In the next step, the rationale must be further extended to SLOs in agreement with the PC at play.

In this article, the problem is attacked via reverse engineering, i.e. instead of checking whether the renormalization invariance of observables obtained from $H^{\text {LO}}_{\pi \!\!{/}\,}$ via the exact solution of Eq. (3) in two- and three-body systems extends to large-$\mathrm {A}$ systems, one attempts to design renormalization prescriptions for a given many-body approximation, many-body observables thus being renormalization-invariant by construction. While this approach departs from the original scheme pursued in $\pi \!\!{/}\,$EFT, it is meant to serve as a first step towards rooting a successful (yet hypothetic) many-body approximation into an EFT that is suited to a large range of nuclear systems.

Finally, let us mention that a similar direction was explored for the dilute homogeneous fermionic matter [39, 40] and for the (1/N) approximation (N being here the spin–isospin degeneracy) [41]. Both works considered dimensional regularization, which is particularly useful for analytical calculations. In this context, the virtue of the present work is to lay down a general framework enabling the systematic study of the renormalization for (as much as possible) generic many-body approximations.

2.3 Set-up of $\mathrm {A}$-body calculations

2.3.1 Hamiltonian

Consider the Hamiltonian

$$\begin{aligned} H&\equiv \sum _{\mathbf {p}\sigma } \frac{p^{2}}{2m} a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma } \nonumber \\&\quad + \frac{1}{2!} \sum _{\sigma _1 \sigma _2} \sum _{\begin{array}{c} \mathbf {p}_1 \mathbf {p}_2 \\ \mathbf {p}_1^{\, \prime }\mathbf {p}_2^{\, \prime } \end{array}} \ (2\pi )^3 \delta (\mathbf {p}_1^{\, \prime }+ \mathbf {p}_2^{\, \prime }- \mathbf {p}_1 - \mathbf {p}_2) \ C_0 \nonumber \\&\qquad \qquad \qquad \qquad \times a^\dagger _{\mathbf {p}_1^{\, \prime }\sigma _1} a^\dagger _{\mathbf {p}_2^{\, \prime }\sigma _2} a_{\mathbf {p}_1 \sigma _1} a_{\mathbf {p}_2 \sigma _2} \end{aligned}$$

(4)

containing, for convenience, a sole spin-independent two-body interaction parametrized by $C_0$. All following derivations can be extended, with minimal modifications, to include spin–isospin dependence and three-body interactions (hence including the case of $\pi \!\!{/}\,$EFT at LO). Homogeneous neutron matter would typically be the first system of practical interest, in which case no three-body interaction and no isospin quantum number need to be considered at LO in $\pi \!\!{/}\,$EFT.

2.3.2 Perturbation theory

Approximations presently considered are formulated within MBPT on the basis of the particular partitioning of H

$$\begin{aligned} H&\equiv H_0 + H_1, \end{aligned}$$

(5a)

$$\begin{aligned} H_0&\equiv \sum _{\mathbf {p}\sigma } \frac{p^{2}}{2m} a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma }, \end{aligned}$$

(5b)

$$\begin{aligned} H_1&\equiv \frac{1}{2!} \sum _{\sigma _1 \sigma _2} \sum _{\begin{array}{c} \mathbf {p}_1 \mathbf {p}_2 \mathbf {p}_1^{\, \prime }\mathbf {p}_2^{\, \prime } \end{array}} \ (2\pi )^3 \delta (\mathbf {p}_1^{\, \prime }+ \mathbf {p}_2^{\, \prime }- \mathbf {p}_1 - \mathbf {p}_2) \ C_0\nonumber \\&\qquad \qquad \qquad \qquad \qquad \times a^\dagger _{\mathbf {p}_1^{\, \prime }\sigma _1} a^\dagger _{\mathbf {p}_2^{\, \prime }\sigma _2} a_{\mathbf {p}_1 \sigma _1} a_{\mathbf {p}_2 \sigma _2}, \end{aligned}$$

(5c)

i.e. the unperturbed Hamiltonian $H_0$ is taken to be the kinetic energy Hamiltonian. Because $H_0$ is invariant under spatial translations, this partitioning is convenient to study homogeneous nuclear matter. In the $\mathrm {A}$-body sector, the unperturbed reference state, i.e. the ground state of $H_0$, is given by the Slater determinant

$$\begin{aligned} \vert \Phi ^{\mathrm {A}}_0 \rangle \equiv \prod _{i=1}^{\mathrm {A}} a^\dagger _{i} \vert 0 \rangle , \end{aligned}$$

(6)

where $\vert 0 \rangle$ denotes the particle vacuum. Here, the index i is a shorthand notation to label hole states $(\sigma _i,\mathbf {p}_i)$, i.e. one-body states that are occupied in $\vert \Phi ^{\mathrm {A}}_0 \rangle$. Similarly a denotes particle states and Greek letters refer to generic states.

In this framework, observables of interest are obtained from the k-body Green’s functions¹ defined as

$$\begin{aligned}&i^k G^{(\mathrm {A},k)}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}} (t_{\mu _1},\dots ,t_{\mu _k}, t_{\nu _1}, \dots ,t_{\nu _k}) \nonumber \\&\quad \equiv \frac{ {\langle \Psi ^\mathrm {A}_0 | \mathrm {T}\left[ a_{\mu _k}(t_{\mu _k}) \dots a_{\mu _1}(t_{\mu _1}) a^\dagger _{\nu _1}(t_{\nu _1}) \dots a^\dagger _{\nu _k}(t_{\nu _k}) \right] | \Psi ^\mathrm {A}_0 \rangle } }{{\langle \Psi ^\mathrm {A}_0 | \Psi ^\mathrm {A}_0 \rangle }},\nonumber \\ \end{aligned}$$

(7)

where $\mathrm {T}$ denotes the time-ordering operator, creation and annihilation operators are in the Heisenberg picture, and $\vert \Psi ^\mathrm {A}_0 \rangle$ is the exact ground state of H connected to $\vert \Phi ^{\mathrm {A}}_0 \rangle$ via the adiabatic theorem of Gell-Mann and Low [42]. Exact Green’s functions can be themselves expressed as a sum over the complete set of linked diagrams $\mathcal {G}^{(\mathrm {A},k)}_n$ carrying k incoming and k outgoing external lines i.e.

$$\begin{aligned}&i^k G^{(\mathrm {A},k)}_{\begin{array}{c} \mu _1 \dots \mu _k\\ \nu _1 \dots \nu _k \end{array}} (t_{\mu _1} \dots t_{\mu _k} ,t_{\nu _1} \dots t_{\nu _k}) \nonumber \\&\quad =\, \sum ^{+\infty }_{n=0} \sum _{\mathcal {G}^{(\mathrm {A},k)}_n \in \mathcal {S}^{(\mathrm {A},k)}_{\text {Exact}}} \mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}_n}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}}(t_{\mu _1} \dots t_{\mu _k} ,t_{\nu _1} \dots t_{\nu _k})\nonumber \\ \end{aligned}$$

(8)

where $\mathcal {S}^{(\mathrm {A},k)}_{\text {Exact}}$ denotes the complete set of diagrams contributing to the k-body Green’s function, n the number of interaction vertices and $\mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}_n}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}} (t_{\mu _1} \dots t_{\mu _k} ,t_{\nu _1} \dots t_{\nu _k})$ the amplitude associated to $\mathcal {G}^{(\mathrm {A},k)}_n$.

Below, Green’s functions are expressed in the energy representation for convenience. A similar expression to Eq. (8) holds in that case. Further, the amplitudes can be expressed in terms of the particle and hole parts of the unperturbed one-body Green’s function associated to the partitioning (5) and the reference state (6), respectively reading

$$\begin{aligned} iG^{(\mathrm {A},1)0+}_{\mu \nu }&\equiv i \sum _{a>\mathrm {A}} \frac{\delta _{\mu a} \delta _{\mu \nu }}{\omega - \frac{\mathbf {p}_a^2}{2m} + i\eta }, \end{aligned}$$

(9a)

$$\begin{aligned} iG^{(\mathrm {A},1)0-}_{\mu \nu }&\equiv i \sum _{i=1}^{\mathrm {A}} \frac{\delta _{\mu i} \delta _{\mu \nu }}{\omega - \frac{\mathbf {p}_i^2}{2m} - i\eta }, \end{aligned}$$

(9b)

where the limit $\eta \rightarrow 0^+$ is implicit. Explicitly employing the two parts of the one-body Green’s function relates to using time-ordered diagrams. Focusing on one diagram $\mathcal {G}^{(\mathrm {A},k)}_n$ belonging to $\mathcal {S}^{(\mathrm {A},k)}_{\text {Exact}}$, the associated amplitude reads

$$\begin{aligned}&\mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}_n}_{\begin{array}{c} \mu _1 \dots \mu _k \nonumber \\ \nu _1 \dots \nu _k \end{array}} (\omega _{\mu _1} \dots \omega _{\mu _k} ,\omega _{\nu _1} \dots \omega _{\nu _k}) \nonumber \\&\quad =\, (-1)^\sigma \frac{(-i)^n}{n!} \sum _{\lambda } \frac{h^{22}_{\lambda \dots \lambda }}{(2!)^2} \dots \frac{h^{22}_{\lambda \dots \lambda }}{(2!)^2} \nonumber \\&\int \frac{\mathrm {d}\omega _\lambda }{2\pi } \dots \prod _{i=1}^{n} 2\pi \delta \left( \pm \omega _\lambda \dots - \omega _\mu \dots + \omega _\nu \dots \right) \nonumber \\&\qquad \times \prod _{e\in I^{+}} iG^{(\mathrm {A},1)0+}_{\lambda \lambda }(\omega _\lambda ) \prod _{e\in I^{-}} iG^{(\mathrm {A},1)0-}_{\lambda \lambda }(\omega _\lambda ) \nonumber \\&\qquad \times \prod _{e\in E_{\mathrm{in}}^{+}} iG^{(\mathrm {A},1)0+}_{\lambda \nu }(\omega _{\nu }) \prod _{e\in E_{\mathrm{in}}^{-}} iG^{(\mathrm {A},1)0-}_{\lambda \nu }(\omega _{\nu }) \nonumber \\&\qquad \times \prod _{e\in E_{\mathrm{out}}^{+}} iG^{(\mathrm {A},1)0+}_{\mu \lambda }(\omega _{\mu }) \prod _{e\in E_{\mathrm{out}}^{-}} iG^{(\mathrm {A},1)0-}_{\mu \lambda }(\omega _{\mu }), \end{aligned}$$

(10)

where $\sigma$ is an integer, $\lambda$, $\mu$ and $\nu$ denote generic single-particle labels for internal, external outgoing and external incoming lines, respectively. Additionally, $I^{+}$ ($I^{-}$) denotes the set of internal particle (hole) lines, $E_{\mathrm{in}}^{+}$ ($E_{\mathrm{in}}^{-}$) the set of external incoming particle (hole) lines and $E_{\mathrm{out}}^{+}$ ($E_{\mathrm{out}}^{-}$) the set of external outgoing particle (hole) lines.

Many-body approximations considered here consist of selecting a subset $\mathcal {S}^{(\mathrm {A},k)}_{\text {MB}} \subset \mathcal {S}^{(\mathrm {A},k)}_{\text {Exact}}$ so that the approximated k-body Green’s function reads

$$\begin{aligned}&i^k G^{(\mathrm {A},k)\text {MB}}_{\begin{array}{c} \mu _1 \dots \mu _k \nonumber \\ \nu _1 \dots \nu _k \end{array}} (\omega _{\mu _1} \dots \omega _{\mu _k} ,\omega _{\nu _1} \dots \omega _{\nu _k}) \nonumber \\&\quad \equiv \sum _{\mathcal {G}^{(\mathrm {A},k)} \in \mathcal {S}^{(\mathrm {A},k)}_{\text {MB}} } \mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}}_{\begin{array}{c} \mu _1 \dots \mu _k\\ \nu _1 \dots \nu _k \end{array}} (\omega _{\mu _1} \dots \omega _{\mu _k} ,\omega _{\nu _1} \dots \omega _{\nu _k}).\nonumber \\ \end{aligned}$$

(11)

2.3.3 Regularization

When considering local interactions, as those derived in any EFT, the amplitude $\mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}}$ contains, in general, ultraviolet (UV) divergences requiring the introduction of a regularization to suppress the high-momentum modes. In this work, the regularization is introduced via a momentum regulator $v_\Lambda (q)$ satisfying

$$\begin{aligned}&\lim _{\Lambda \rightarrow +\infty } v_\Lambda (q) = 1, \end{aligned}$$

(12a)

$$\begin{aligned}&\lim _{q \rightarrow +\infty } v_\Lambda (q) = 0, \end{aligned}$$

(12b)

$$\begin{aligned}&v_\Lambda (0) = 1, \end{aligned}$$

(12c)

$$\begin{aligned}&\int ^{+\infty }_{0} \mathrm {d}q \ v^2_\Lambda (q) < +\infty , \end{aligned}$$

(12d)

where $\Lambda$ represents the characteristic cut-off scale beyond which high-momentum modes are suppressed and $\mathbf {q}$ denotes the relative momentum between two nucleons. Eventually, one is interested in the limit where $\Lambda \gg Q$ with Q the scale characterizing the low-energy observables of interest. The limit of large $\Lambda$ ensures that the evaluation of observables is independent of the high-energy physics. The difference between the actual high-momentum interaction (which is unknown) and the regularized field theory one (which is arbitrary) is then captured effectively in its parameters referred to as low-energy constants (LECs). Correspondingly, the LECs explicitly depend on the regulator $v_\Lambda$ and in particular on the cut-off $\Lambda$. Eventually, the Hamiltonian to be actually considered reads

$$\begin{aligned} H&\equiv \sum _{\mathbf {p}\sigma } \frac{p^{2}}{2m} a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma } \nonumber \\&\quad + \frac{1}{2!} \sum _{\sigma _1 \sigma _2} \sum _{\begin{array}{c} \mathbf {p}_1 \mathbf {p}_2 \\ \mathbf {p}_1^{\, \prime }\mathbf {p}_2^{\, \prime } \end{array}} \ (2\pi )^3 \delta (\mathbf {p}_1^{\, \prime }+ \mathbf {p}_2^{\, \prime }- \mathbf {p}_1 - \mathbf {p}_2) \ C_0(\Lambda ) \ \nonumber \\&\qquad \qquad \qquad \qquad \quad \times v_\Lambda (\mathbf {q}) v_\Lambda (\mathbf {q}^{\, \prime }) \ a^\dagger _{\mathbf {p}_1^{\, \prime }\sigma _1} a^\dagger _{\mathbf {p}_2^{\, \prime }\sigma _2} a_{\mathbf {p}_1 \sigma _1} a_{\mathbf {p}_2 \sigma _2},\nonumber \\ \end{aligned}$$

(13)

where the dependence of the LEC $C_0(\Lambda )$ on the regularization is made explicit and where incoming and outgoing relative momenta are, respectively, defined as

$$\begin{aligned} \mathbf {q}&\equiv \frac{\mathbf {p}_1 - \mathbf {p}_2}{2}, \end{aligned}$$

(14a)

$$\begin{aligned} \mathbf {q}^{\, \prime }&\equiv \frac{\mathbf {p}_1^{\, \prime }- \mathbf {p}_2^{\, \prime }}{2}. \end{aligned}$$

(14b)

To make the cancellation of UV divergences explicit, LECs are usually decomposed, for a given renormalization scheme, into a renormalized component and a counterterm [43], i.e. in the case of (13) as

$$\begin{aligned} C_0(\Lambda ) = C_0^R + \delta C_0(\Lambda ), \end{aligned}$$

(15)

where $C^R_0$ is the finite renormalized LEC and $\delta C_0(\Lambda )$ encapsulates all counterterms that will be necessary to cancel UV divergences. For a given set of diagrams $\mathcal {S}^{(\mathrm {A},k)}_{\text {MB}}$ with renormalized vertices $C^R_0$, the aim is thus to derive an additional set of diagrams $\mathcal {S}^{(\mathrm {A},k)}_{\text {MB,ct}}$ with renormalized and counterterm vertices (to be derived as well) such that the total sum

$$\begin{aligned}&i^k G^{(\mathrm {A},k)\text {RMB}}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}} (\omega _{\mu _1} \dots \omega _{\mu _k} ,\omega _{\nu _1} \dots \omega _{\nu _k})\nonumber \\&\quad \equiv \sum _{\mathcal {G}^{(\mathrm {A},k)} \in \mathcal {S}^{(\mathrm {A},k)}_{\text {MB}} \cup \mathcal {S}^{(\mathrm {A},k)}_{\text {MB,ct}} } \mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}} (\omega _{\mu _1} \dots \omega _{\mu _k} ,\omega _{\nu _1} \dots \omega _{\nu _k}) \end{aligned}$$

(16)

converges to a finite value when $\Lambda \gg Q$.

In the perturbation theory employing the particle vacuum $\vert 0 \rangle$ as a reference state, the counterterms needed to achieve renormalization can be identified systematically via the application of the BPHZ theorem [12,13,14,15,16,17]. Applying a similar procedure to time-ordered diagrams stemming from a perturbative expansion around an in-medium reference state $\vert \Phi ^{\mathrm {A}}_0 \rangle$ poses some challenges. One must overcome these challenges given that perturbation theory, or any other expansion many-body method for that matter, can only be efficiently formulated around an in-medium reference state as soon as one targets systems containing more than a few particles. One example of apparent difficulty relates to the fact that the renormalization procedure, if possible, might have to be achieved for each $\mathrm {A}$ given that the particle and hole propagators do depend on $\mathrm {A}$. In numerical calculations, observables used to renormalize the LECs would have to be computed for each $\mathrm {A}$. This situation would penalize both the numerical efficiency and the predictive power of the theory.

Therefore, the goal of the present manuscript is to design a renormalization procedure for observables computed on the basis of a given MBPT approximation such that necessary k-body counterterms, hopefully limited to small k, are independent of $\mathrm {A} \ge k$.

3 Basic tools

Weinberg’s asymptotic theorem [11] and the BPHZ [12,13,14,15,16,17] theorem represent key tools to develop the renormalization procedure exposed in Sect. 4. Consequently, the two theorems are briefly recalled in the present section.

3.1 Weinberg’s asymptotic theorem

The UV convergence of the amplitude of a time-ordered diagram is analyzed as the convergence problem of a multidimensional integral of a multivariate function. The main ingredients to understand Weinberg’s asymptotic theorem are presently introduced for a generic multivariate function. Furthermore, the theorem is recast in terms of diagrams. For the complete proof and further discussion of Weinberg’s asymptotic theorem, see Ref. [11].

3.1.1 Asymptotic coefficients and multidimensional integrals

Let us first introduce the definition of asymptotic coefficients (provided they exist) of a function $f : \mathbb {R}^n \rightarrow \mathbb {C}$ as given in [11]. For any vector subspace²$S = \left\{ \mathbf {L}_1, \dots , \mathbf {L}_m \right\}$ of $\mathbb {R}^n$ with $m \le n$ and $\mathbf {L}_1, \dots , \mathbf {L}_m$ being m independent $\mathbb {R}^n$-vectors, and any compact region $W \subset \mathbb {R}^n$, the asymptotic coefficients are defined as the numbers $\alpha \left( \left\{ \mathbf {L}_1, \dots , \mathbf {L}_r \right\} \right)$ and $\beta \left( \left\{ \mathbf {L}_1, \dots , \mathbf {L}_r \right\} \right)$ (with $1 \le r \le m$) such that for any $\mathbf {C} \in W$

$$\begin{aligned}&f\left( \eta _1 \dots \eta _m \mathbf {L}_1 + \eta _2 \dots \eta _m \mathbf {L}_2 + \cdots + \eta _m \mathbf {L}_m + \mathbf {C} \right) \nonumber \\&\quad = O\left( \eta _1^{\alpha \left( \left\{ \mathbf {L}_1\right\} \right) } \left( \ln \eta _1 \right) ^{\beta \left( \left\{ \mathbf {L}_1\right\} \right) } \right. \nonumber \\&\qquad \qquad \left. \eta _2^{\alpha \left( \left\{ \mathbf {L}_1,\mathbf {L}_2\right\} \right) } \left( \ln \eta _2 \right) ^{\beta \left( \left\{ \mathbf {L}_1,\mathbf {L}_2\right\} \right) } \times \cdots \right. \nonumber \\&\qquad \qquad \left. \cdots \times \eta _m^{\alpha \left( \left\{ \mathbf {L}_1, \dots , \mathbf {L}_m\right\} \right) } \left( \ln \eta _m \right) ^{\beta \left( \left\{ \mathbf {L}_1, \ldots , \mathbf {L}_m\right\} \right) } \right) , \end{aligned}$$

(17)

if $\eta _1 \dots \eta _m$ go independently to infinity. The asymptotic coefficients $\alpha \left( S \right)$ and $\beta \left( S \right)$ can be interpreted as the asymptotic coefficients $\alpha \left( \left\{ \mathbf {L} \right\} \right)$ and $\beta \left( \left\{ \mathbf {L} \right\} \right)$ for $\mathbf {L}$ a ‘typical’ vector in S i.e. fixing $\eta _1 \dots \eta _{m-1}$ sufficiently large and $\mathbf {C} \in W$,

$$\begin{aligned}&f\left( \left[ \eta _1 \dots \eta _{m-1} \mathbf {L}_1 + \eta _2 \cdots \eta _{m-1} \mathbf {L}_2 + \cdots \right. \right. \nonumber \\&\quad \left. \left. \cdots + \eta _{m-1} \mathbf {L}_{m-1} + \mathbf {L}_m \right] \eta _m + \mathbf {C} \right) \nonumber \\&\quad = O\left( \eta _m^{\alpha \left( S\right) } \left( \ln \eta _m \right) ^{\beta \left( S\right) } \right) , \end{aligned}$$

(18)

when $\eta _m$ goes to infinity.

Considering now integrals of a function f, the integration along the directions $\mathbf {L}_1, \ldots , \mathbf {L}_r$ is defined as

$$\begin{aligned}&f_{\mathbf {L}_1, \cdots , \mathbf {L}_r}\left( \mathbf {X}\right) \nonumber \\&\quad \equiv \int _{-\infty }^{+\infty } \mathrm {d}y_1 \ \dots \int _{-\infty }^{+\infty } \mathrm {d}y_r\, f\left( \mathbf {X} + y_1 \mathbf {L}_1 + \cdots + y_r \mathbf {L}_r \right) , \nonumber \\ \end{aligned}$$

(19)

where $\mathbf {X}$ is a vector in $\mathbb {R}^n$. Thanks to Fubini’s theorem, if $f_{\mathbf {L}_1, \dots , \mathbf {L}_r}\left( \mathbf {X}\right)$ exists (in the sense that the integral is absolutely convergent), the integration depends only on the vector space $I = \left\{ \mathbf {L}_1, \dots , \mathbf {L}_r \right\}$ so that one defines

$$\begin{aligned} f_{I}\left( \mathbf {X}\right) \equiv \int _{\mathbf {Y} \in I} \mathrm {d}^r\mathbf {Y} \ f(\mathbf {X} + \mathbf {Y}) \equiv f_{\mathbf {L}_1, \dots , \mathbf {L}_r}\left( \mathbf {X} \right) . \end{aligned}$$

(20)

Furthermore, $f_{I}\left( \mathbf {X}\right)$ depends only on the projection of $\mathbf {X}$ along I³. Choosing a subspace E such that $\mathbb {R}^n = I \oplus E$, the domain of definition of the function $f_{I}\left( \mathbf {X}\right)$ can be restricted to E.

In the case of the amplitude of a time-ordered diagram, the general integrand depends on $(\omega _1, \mathbf {p}_1, \dots )$ and is integrated on the internal energies and momenta. Therefore in this case, I denotes the vector space of internal (one-body) energies and momenta whereas the vector space E denotes the space of external (one-body) energies and momenta⁴. The general vector space $\mathbb {R}^n = I \oplus E$ denotes the vector space of all (one-body) energies and momenta (internal and external).

As an example, the asymptotic coefficients $\alpha \left( S\right)$ of the in-vacuum propagator

$$\begin{aligned} iG^{(0,1)0}_{\mathbf {p}\sigma }(\omega ) = \frac{i}{\omega - \frac{p^2}{2m} + i \eta } \end{aligned}$$

(21)

are now extracted. The in-vacuum propagator $iG^{(0,1)0}_{\mathbf {p}\sigma }(\omega )$ is interpreted as a multivariate function on the vector space $\mathbb {R}^4 = \{ \mathbf {e}_\omega , \mathbf {e}_{p_x}, \mathbf {e}_{p_y}, \mathbf {e}_{p_z}\}$, so that⁵

$$\begin{aligned} f(\omega \mathbf {e}_\omega + p_x \mathbf {e}_{p_x} + p_y \mathbf {e}_{p_y} + p_z \mathbf {e}_{p_z} ) \equiv iG^{(0,1)0}_{\mathbf {p}\sigma }(\omega ). \end{aligned}$$

(22)

One can show that, in this case, the asymptotic coefficients of the in-vacuum propagator $\alpha ^0$ read

$$\begin{aligned} \alpha ^0\left( S\right) = {\left\{ \begin{array}{ll} -1 &{}\text { if } S = \{ \mathbf {e}_\omega \} ,\\ -2 &{}\text { if } S = \{ \mathbf {L} \} \text { with } \mathbf {L} \notin \{\mathbf {e}_\omega \}, \\ -2 &{}\text { if } \dim {S} \ge 2. \end{array}\right. } \end{aligned}$$

(23)

3.1.2 Asymptotic theorem

With all the notations introduced before, the general asymptotic theorem follows.

If a function $f\left( \mathbf {X}\right)$ possesses asymptotic coefficients $\alpha \left( S\right) ,\beta \left( S\right)$ for any non-null subspace $S \subset \mathbb {R}^n$, if $f(\mathbf {X})$ is integrable for any finite region of $\mathbb {R}^n$ and if $D_I < 0$ where

$$\begin{aligned} D_I \equiv \max \limits _{{S^{{\, \prime }}} \subseteq I} \left[ \alpha \left( {S^{{\, \prime }}}\right) + \dim {S^{\, \prime }} \right] \ , \end{aligned}$$

(24)

then $f_I\left( \mathbf {X}\right)$ exists i.e. is absolutely convergent.

To apply this theorem to a time-ordered diagram $\mathcal {G}$, it is reformulated in terms of sub-diagrams as detailed in Ref. [11]. This is done by associating to any sub-space of integration $S^{\, \prime }\subseteq I$ a sub-diagram $\gamma \subseteq \mathcal {G}$ made of internal lines. In particular $S^{\, \prime }= I$ corresponds to the sub-diagram $\gamma$ made of $\mathcal {G}$ itself without its external lines. The quantity $\alpha \left( S^{\, \prime }\right) + \dim {S^{\, \prime }}$ corresponds to the superficial degree of divergence of the associated sub-diagram $\gamma$ of $\mathcal {G}$ and is denoted $D(\gamma )$. For a diagram $\mathcal {G}^{(0,k)}_n$ where lines denote the unperturbed in-vacuum propagator (21), n denotes the number of vertices (5b) and k is the number of incoming (and of outgoing) external lines,

$$\begin{aligned} D\left( \mathcal {G}^{(0,k)}_n\right) = 5 - 3k + n. \end{aligned}$$

(25)

Having $D_I < 0$ is thus equivalent to having $D(\gamma ) < 0$ for all $\gamma \subseteq \mathcal {G}$ made of internal lines. Consequently, from Weinberg’s asymptotic theorem, the well-known power-counting theorem follows. The amplitude associated to $\mathcal {G}$ is finite if $D(\gamma ) < 0$ for any sub-diagram $\gamma \subseteq \mathcal {G}$ made of internal lines of $\mathcal {G}$.

Weinberg’s asymptotic theorem is very powerful, as it proves the convergence of the amplitude associated to a diagram with the sole knowledge of the asymptotic coefficients $\alpha \left( S\right)$ associated to the propagator. This particular property is fundamental for the development of the renormalization prescription exposed in Sect. 4.

3.2 Subtractions and BPHZ theorem

Whilst the asymptotic theorem allows one to prove the convergence of the calculations with adequate counterterms, the BPHZ theorem offers a systematic way to generate sufficient counterterms to lower the superficial degree of divergence of the internal sub-diagrams $\gamma \subseteq \mathcal {G}$. The procedure to do so is now briefly detailed.

3.2.1 Definitions

A sub-diagram $\gamma$ of $\mathcal {G}$ is defined as a subset of lines and vertices contained in $\mathcal {G}$ where end points of the lines of $\gamma$ belong to its vertices. A diagram $\mathcal {G}_1$ is said to be included in $\mathcal {G}_2$ and denoted as $\mathcal {G}_1 \subseteq \mathcal {G}_2$ if their set of lines verify the same inclusion relation. In particular, a sub-diagram $\gamma$ of $\mathcal {G}$ verifies $\gamma \subseteq \mathcal {G}$. The sub-diagram generated by the intersection of lines of two sub-diagrams $\gamma _1$ and $\gamma _2$ defines a sub-diagram $\gamma$ and is denoted as

$$\begin{aligned} \gamma \equiv \gamma _1 \cap \gamma _2. \end{aligned}$$

(26)

Two sub-diagrams $\gamma _1$ and $\gamma _2$ that have neither lines nor vertices in common are said to be disjoint and the result is denoted as

$$\begin{aligned} \gamma _1 \cap \gamma _2 = \emptyset . \end{aligned}$$

(27)

If neither $\gamma _1 \subseteq \gamma _2$ nor $\gamma _2 \subseteq \gamma _1$ and $\gamma _1 \cap \gamma _2 \ne \emptyset$ they are said to be overlapping. Otherwise they are said to be non-overlapping. For a set of non-overlapping sub-diagrams $\gamma _1,\gamma _2,\ldots ,\gamma _n$ of $\mathcal {G}$, the reduced diagram $\mathcal {G} {\setminus } \{\gamma _1,\gamma _2,\ldots ,\gamma _n\}$ is defined by the diagram resulting from $\mathcal {G}$ after contracting all lines of $\gamma _1,\gamma _2,\ldots ,\gamma _n$ to a point.

The amplitude associated to $\mathcal {G}$ is denoted as $\mathcal {A}^{\mathcal {G}}$. Given a set of mutually disjoint sub-diagrams $\gamma _1,\gamma _2,\ldots ,\gamma _n$ of $\mathcal {G}$, the corresponding amplitude is expressed as the product of the amplitudes $\mathcal {A}^{\gamma _j}$ of the sub-diagrams, while the remainder is denoted as $\mathcal {A}^{\mathcal {G} {\setminus } \{\gamma _1,\ldots ,\gamma _n \} }$ such that

$$\begin{aligned} \mathcal {A}^{\mathcal {G}} = \mathcal {A}^{\mathcal {G} {\setminus } \{\gamma _1,\ldots ,\gamma _n \} } \prod _{j=1}^{n} \mathcal {A}^{\gamma _j}. \end{aligned}$$

(28)

A sub-diagram $\gamma$ of $\mathcal {G}$ is referred to as a renormalization part if it is a one-particle irreducible (1PI) diagram with a superficial degree of divergence greater than or equal to 0, i.e. if

$$\begin{aligned} D(\gamma ) \ge 0. \end{aligned}$$

(29)

3.2.2 Recursive subtractions of UV divergences

The BPHZ procedure defines recursively a renormalized amplitude $R^{\mathcal {G}}$ associated to the diagram $\mathcal {G}$ [12, 13]. If the amplitude associated to $\mathcal {G}$ is convergent to begin with, then

$$\begin{aligned} R^{\mathcal {G}} \equiv \mathcal {A}^{\mathcal {G}}. \end{aligned}$$

(30)

If the diagram does not contain any renormalization part but is superficially divergent, it is called primitively divergent. In that case the renormalized amplitude is defined by

$$\begin{aligned} R^{\mathcal {G}} \equiv (1-t_{\mathcal {G}}) \mathcal {A}^{\mathcal {G}}, \end{aligned}$$

(31)

where $t_{\mathcal {G}}$ is the operator of the Taylor expansion with respect to the external momenta⁶ around 0 up to the order of the superficial degree of divergence $D(\gamma )$ of the diagram, i.e.

$$\begin{aligned}&t_{\gamma } \mathcal {A}^{\gamma }_{\begin{array}{c} p_{1}, \dots , p_{k}\\ p_{1}^{\, \prime }, \dots , p_{k}^{\, \prime } \end{array}} \equiv \sum _{j=0}^{D(\gamma )} \frac{1}{j!} \sum _{\begin{array}{c} s_{1} + \dots + s_{k}^{\, \prime }\ge 0 ,\nonumber \\ s_{1} +\dots + s_{k}^{\, \prime }=j \end{array}} ,\nonumber \\&\left. \frac{\partial ^{j}\mathcal {A}^{\gamma }}{\partial p_{1}^{s_{1}} \cdots \partial p_{k}^{s_{k}} \partial p_{1}^{{\, \prime }s_{1}^{\, \prime }} \cdots \partial p_{k}^{{\, \prime }s_{k}^{\, \prime }} } \right| _{ p_{1} = \dots = p_{k}^{\, \prime }= 0 } p_{1}^{s_{1}} \cdots p_{k}^{{\, \prime }s_{k}^{\, \prime }}.\nonumber \\ \end{aligned}$$

(32)

If $\mathcal {G}$ is superficially divergent and contains divergent sub-diagrams, the renormalized amplitude is defined recursively as

$$\begin{aligned} R^{\mathcal {G}} \equiv (1-t_{\mathcal {G}})\bar{R}^{\mathcal {G}}, \end{aligned}$$

(33)

where $\bar{R}^{\mathcal {G}}$ corresponds to the amplitude where all sub-divergences have already been subtracted. The subtraction by the Taylor operator $t_{\gamma }$ corresponds to adding the amplitude associated to a diagram where the divergent sub-diagram $\gamma$ has been replaced by a so-called counterterm vertex. Here $t_{\gamma }$ corresponds to a zero-momentum subtraction of UV divergences. Different subtractions can be used by modifying the definition of $t_{\gamma }$ [43].

3.2.3 Forest formula

The recursion relation (33) was solved explicitly by the forest formula [16], which is based on the concept of i-forest (for inclusion-forest). An i-forest is defined as any set of sub-diagrams (including the empty set) that are mutually non-overlapping. This way, the Hasse diagram⁷, for the order relation $\subseteq$ on the mutually non-overlapping sub-diagrams, represents a forest i.e. a set of disconnected trees (see the right panel of Fig. 1 for an example where the Hasse diagram of the i-forest is made of only one tree). An i-forest is said to be connected if its Hasse diagram is connected. A connected i-forest is also referred to as an i-tree. As for usual forests, an i-forest can be decomposed as the set of its connected components (i.e. as a set of i-trees). An i-forest is usually depicted by drawing boxes around the sub-diagrams as illustrated in the left panel of Fig. 1. The boxes are, thus, not allowed to overlap but can be nested.

Fig. 1

On the left panel, an example of an i-forest (depicted by boxes) for a three-loop ladder diagram contributing to $G^{(0,2)}$. The middle panel pictures the diagram with the counterterm associated to this i-forest (the vertex with an empty circle denotes the two-body counterterm appearing for this particular i-forest). The associated Hasse diagram is depicted on the right panel

An i-forest is restricted if each of its boxes contains only renormalization parts. To each restricted i-forest $\mathscr {F}$ one associates again an amplitude, namely

$$\begin{aligned} \Omega ^{\mathscr {F}} \equiv \tilde{\prod _{\gamma \in \mathscr {F}}} (-t_{\gamma }) \mathcal {A}^{\mathcal {G}}, \end{aligned}$$

(34)

where the tilde over the product sign stands for the fact that in case of nested diagrams within the i-forest one has to apply the Taylor operators from the innermost to the outermost diagram while for disjoint sub-diagrams the expressions are naturally independent of the order of the Taylor operators by virtue of Eq. (28). Each i-forest corresponds to a particular diagram with counterterms. The topology of the resulting diagram consists of the original diagram $\mathcal {G}$ where the sub-diagrams $\gamma$ of the i-forest have been contracted into vertices corresponding to counterterms. The nature of the counterterms depends on the i-forest. From now on, as a shorthand notation, a diagram where an i-forest $\mathscr {F}$ is pictured with boxes will represent directly the amplitude $\Omega ^{\mathscr {F}}$. See Fig. 2 for an example based on a one-loop diagram.

Fig. 2

Representation, in the case of a one-loop diagram, of the amplitude $\Omega ^{\mathscr {F}}$ both with the i-forest pictured on the original diagram and with an explicit counter-term vertex. The filled dot represents $C^R_0$ whereas the hatched vertex represents the counterterm associated to the i-forest represented by a box on the left-hand side diagram

Eventually, the forest formula states that the renormalized amplitude of the diagram $\mathcal {G}$ is given by the sum over all restricted i-forests, i.e.

$$\begin{aligned} R^{\mathcal {G}}=\sum _{\mathscr {F} \in \mathcal {F}_{R}(\mathcal {G})} \Omega ^{\mathscr {F}}, \end{aligned}$$

(35)

where $\mathcal {F}_{R}(\mathcal {G})$ denotes the set of restricted i-forests and where it is understood that the empty i-forest (i.e. without any box around a sub-diagram) stands for the diagram $\mathcal {G}$ itself. The term with the empty i-forest corresponds to the UV-divergent diagram, while all the other terms correspond to diagrams including counterterm vertices cancelling the original UV divergences.

4 Renormalization of in-medium diagrams

Theorems introduced in Sect. 3 are valid regardless of the peculiar perturbation theory at stake i.e. they can be used to analyze UV divergencies and derive sufficient counterterms whether the reference state is $\vert 0 \rangle$ or $\vert \Phi ^\mathrm {A}_0 \rangle$ [44]. However, as argued in Sect. 2.3, the direct use of BPHZ will generate $\mathrm {A}$-dependent counterterms that can and should be avoided.

A key component of the renormalization procedure exposed in this section is to relate UV divergences occurring in the calculation of the approximated in-medium Green’s functions (i.e. with $\vert \Phi ^\mathrm {A}_0 \rangle$ as a reference state) to UV divergences occurring in the calculation of in-vacuum Green’s functions (i.e. with $\vert 0 \rangle$ as a reference state) in a related approximation.

4.1 Cutting procedure

Let $\mathcal {G}^{(\mathrm {A},k)}_n \in \mathcal {S}^{(\mathrm {A},k)}_{\text {MB}}$ and $\mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}_n}$ be its associated amplitude given in Eq. (10). For any finite $\mathrm {A}$, the set of states $(\mathbf {p}_i,\sigma _i)$ labelling unperturbed hole propagators lie in a compact space.⁸ Consequently, to prove the UV convergence of Eq. (10), it is sufficient to prove the convergence of the sub-integral related to the sole particle propagators. The UV behavior of $\mathcal {G}^{(\mathrm {A},k)}_n$ is the same as the UV behavior of an associated diagram made only of particle propagators. Such a diagram can be defined as the one made out of the same n vertices but with incoming particle lines $E_{\mathrm{in}}^{+} \cup I^{-}$, outgoing particle lines $E_{\mathrm{out}}^{+} \cup I^{-}$ and internal particle lines $I^{+}$. As each line in the aforementioned diagram corresponds to an unperturbed particle propagator $iG^{(\mathrm {A},1)0+}$ and possesses $(k+p)$ outgoing and incoming external lines, where $p \equiv \# I^{-}$, the diagram is denoted in the following as $\mathcal {G}^{(\mathrm {A},k+p)}_{n}$. Diagrammatically, $\mathcal {G}^{(\mathrm {A},k+p)}_{n}$ is obtained from the original diagram by cutting all internal unperturbed hole lines in $\mathcal {G}^{(\mathrm {A},k)}_{n}$ and replacing each of them by an incoming and an outgoing external unperturbed particle propagator. Examples of the cutting procedure are displayed in Fig. 3.

Consequently, the cutting procedure recasts the analysis of the UV behavior of a diagram contributing to the in-medium k-body Green’s function $i^k G^{(\mathrm {A},k)}$ into the analysis of a diagram solely made of unperturbed particle propagators contributing to the in-medium $(k+p)$-body Green’s function $i^{k+p}G^{(\mathrm {A},k+p)}$.

Fig. 3

Examples of the cutting procedure applied to diagrams $\mathcal {G}^{(\mathrm {A},k)}_2$ and resulting diagrams $\mathcal {G}^{(\mathrm {A},k+p)}_2$, with $p=2$ (top) and $p=1$ (bottom)

4.2 Relation to in-vacuum diagrams

At this stage if BPHZ is applied to the diagrams $\mathcal {G}^{(\mathrm {A},k+p)}_n$ resulting from the cutting procedure, $\mathrm {A}$-dependent counterterms will be generated. To avoid this, each diagram $\mathcal {G}^{(\mathrm {A},k+p)}_n$ is associated to a diagram $\mathcal {G}^{(0,k+p)}_n$ with the same UV behavior but that is made only of in-vacuum propagators.

This is done by noticing that the unperturbed particle one-body Green’s function $iG^{(\mathrm {A},1)0+}$ possesses the same asymptotic coefficients $\alpha \left( S\right)$ as the in-vacuum one-body Green’s function $iG^{(0,1)0}$, i.e. for any sub-space S of $\{ \mathbf {e}_\omega , \mathbf {e}_{p_x}, \mathbf {e}_{p_y}, \mathbf {e}_{p_z}\}$

$$\begin{aligned} \alpha ^+\left( S\right) = \alpha ^0\left( S\right) , \end{aligned}$$

(36)

where $\alpha ^+\left( S\right)$ correspond to the asymptotic coefficients of $iG^{(\mathrm {A},1)0+}_{\mathbf {p}\sigma }(\omega )$. Applying the asymptotic theorem stated in Sect. 3.1.2, one deduces that any diagram made solely of unperturbed particle propagators is UV-convergent if and only if the same diagram made of in-vacuum propagators is UV-convergent.

One can thus focus on the study of the diagram $\mathcal {G}^{(0,k+p)}_n$ associated to $\mathcal {G}^{(\mathrm {A},k+p)}_n$. As detailed in Sect. 3.2, the UV divergences are canceled out by additional diagrams, containing counterterms, generated by the BPHZ procedure. The renormalized amplitude of $\mathcal {G}^{(0,k+p)}_n$ reads

$$\begin{aligned} R^{\mathcal {G}^{(0,k+p)}_n} = \sum _{\mathscr {F} \in \mathcal {F}_{R}(\mathcal {G}^{(0,k+p)}_n)} \Omega ^{\mathscr {F}}. \end{aligned}$$

(37)

In order to transport the cancellation of UV divergences back to $\mathcal {G}^{(\mathrm {A},k+p)}_n$, one introduces the amplitude

$$\begin{aligned} R^{\mathcal {G}^{(\mathrm {A},k+p)}_n} \equiv \sum _{\mathscr {F} \in \mathcal {F}_{R}(\mathcal {G}^{(\mathrm {A},k+p)}_n)} \Omega ^{(\mathrm {A})\mathscr {F}}, \end{aligned}$$

(38)

where $\Omega ^{(\mathrm {A})\mathscr {F}}$ denotes the amplitude associated to the same diagram as the one originating from $\Omega ^{\mathscr {F}}$ except that lines denoting in-vacuum propagators are replaced by unperturbed particle propagators. Therefore, the amplitudes $\Omega ^{(\mathrm {A})\mathscr {F}}$ contains k-body counterterms as in $\Omega ^{\mathscr {F}}$ which are, by construction, the same for all $\mathrm {A}\ge k$.

4.3 General procedure

The procedure to derive UV-finite k-body Green’s functions computed with respect to $\vert \Phi ^\mathrm {A}_0 \rangle$ on the basis of a many-body approximation defined by a truncated set of particle–hole diagrams $\mathcal {S}^{(\mathrm {A},k)}_{\text {MB}}$ is now recapitulated. For any diagram $\mathcal {G}^{(\mathrm {A},k)} \in \mathcal {S}^{(\mathrm {A},k)}_{\text {MB}}$ with p internal hole lines:

• 1. Apply the cutting procedure to $\mathcal {G}^{(\mathrm {A},k)}$, i.e. cut internal hole lines into external particle propagators and replace external hole lines by particle propagators. From this first step the associated diagram $\mathcal {G}^{(\mathrm {A},k+p)}$ is obtained.

• 2. Replace in $\mathcal {G}^{(\mathrm {A},k+p)}$ all particle propagators by in-vacuum propagators, thus generating the diagram $\mathcal {G}^{(0,k+p)}$ contributing to the $(k+p)$-body Green’s function with respect to $\vert 0 \rangle$.

The set of diagrams thus obtained is denoted as $\mathcal {S}^{(0,k+p)}_{\text {MB}}$.

• 3. The renormalization of LECs can be carried out as usual on $\mathcal {G}^{(0,k+p)}$, leading to the introduction of an additional set of diagrams with counterterms denoted as $\mathcal {S}^{(0,k+p)}_{\text {MB,ct}}$. The LECs are typically fixed by matching (a subset of) the $(k+p)$-body Green’s functions to observables in their approximations defined by $\mathcal {S}^{(0,k+p)}_{\text {MB}} \cup \mathcal {S}^{(0,k+p)}_{\text {MB,ct}}$.

Then, for any diagram $\mathcal {G}^{(0,k+p)}_{\text {ct}} \in \mathcal {S}^{(0,k+p)}_{\text {MB,ct}}$:

• 4. Replace each in-vacuum propagator with a particle propagator. This generates the diagram $\mathcal {G}^{(\mathrm {A},k+p)}_{\text {ct}}$.

• 5. External lines obtained via the cut (replacement) of internal hole lines are closed (replaced) by a hole line. This leads to the diagram $\mathcal {G}^{(\mathrm {A},k)}_{\text {ct}}$.

The set of all diagrams $\mathcal {G}^{(\mathrm {A},k)}_{\text {ct}}$ thus obtained is denoted as $\mathcal {S}^{(\mathrm {A},k)}_{\text {MB,ct}}$.

• 6. The approximated UV-finite k-body Green’s function eventually can be considered as Eq. (16), i.e.

  $$\begin{aligned}&i^k G^{(\mathrm {A},k)\text {RMB}}_{\begin{array}{c} \mu _1 \dots \mu _k \\ \nu _1 \dots \nu _k \end{array}} (\omega _{\mu _1} \cdots \omega _{\mu _k} ,\omega _{\nu _1} \cdots \omega _{\nu _k}) \nonumber \\&\equiv \sum _{\mathcal {G}^{(\mathrm {A},k)} \in \mathcal {S}^{(\mathrm {A},k)}_{\text {MB}} \cup \mathcal {S}^{(\mathrm {A},k)}_{\text {MB,ct}} } \mathcal {A}^{\mathcal {G}^{(\mathrm {A},k)}}_{\begin{array}{c} \mu _1 \cdots \mu _k \nonumber \\ \nu _1 \cdots \nu _k \end{array}} (\omega _{\mu _1} \cdots \omega _{\mu _k} ,\omega _{\nu _1} \cdots \omega _{\nu _k}) \ . \end{aligned}$$

4.4 Discussion

Section 4 introduced a general procedure to transpose a renormalization scheme formulated for a perturbative expansion around the particle vacuum $\vert 0 \rangle$ to a perturbative expansion using an $\mathrm {A}$-body reference state $\vert \Phi ^\mathrm {A}_0 \rangle$ whilst keeping the same partitioning of H i.e. taking $H_0$ as the kinetic energy. Given an approximation to the in-medium k-body Green’s function, $\mathrm {A}$-independent counterterms can be fixed by matching in-vacuum $(k+p)$-body Green’s functions $i^{k+p}G^{(0,k+p)}$ to appropriate observables instead of in-medium k-body ones $i^{k}G^{(\mathrm {A},k)}$. Depending on the approximation, the index p may span a large range of values, e.g. one may consider diagrams containing an arbitrary large number of hole lines p. A key practical point relates thus to which subset of those $(k+p)$-body Green’s functions $i^{k+p}G^{(0,k+p)}$ actually need to be used to fix the LECs. Indeed, if the range of $(k+p)$-body sectors effectively needed to renormalize a given set of many-body diagrams is not small, the procedure will be laborious. In practice, however, the BPHZ procedure stipulates that only 1PI sub-diagrams that are superficially divergent, the so-called renormalization parts, must be considered. This hopefully limits the possible topology of the renormalization parts and of the associated counterterms such that the set of in-vacuum Green’s functions that needs to be eventually considered to achieve renormalization remains very limited.

The hope is thus that the set of identified counterterms is gentle. For pure neutron systems, where the three-body contact term entering $H^{\text {LO}}_{\pi \!\!{/}\,}$ is inactive, renormalization was shown to be achievable with only $\delta C_0(\Lambda )$ counterterms on the basis of exact (no approximation) vacuum two- and three-body Green’s functions [27]. If both neutrons and protons are present, the $D_0(\Lambda )$ three-body contact term is necessary to achieve renormalization. Numerical calculations seem to indicate that this conclusion extends to the vacuum four-body Green’s function [28].

As a matter of fact, the set of counterterms needed to handle a given many-body approximation is mainly driven by two important features, i.e. (i) the topology of the diagrams controlling their ultraviolet character and (ii) the degeneracy factor g of the interacting fermions, e.g. $g=2$ for spin one-half neutrons and $g=4$ if protons are added. While both considerations are of different nature, the second one must not be overlooked in practice given that Pauli blocking can forbid certain topologies and their otherwise divergent character, thus drastically limiting the number of potentially needed counterterms compared to the case of an arbitrary g. As discussed in detail for the dilute Fermi gas in [39, 45], a sole⁹$\delta C_0(\Lambda )$ counterterm is needed to achieve renormalization up to 4th order in perturbation theory for $g=2$. In the case where $g=4$, an additional $D_0(\Lambda )$ three-body contact term becomes necessary at 4th order, as first pointed out in [39]¹⁰. Going to higher orders, diagrams with a larger number of particle lines propagating simultaneously are allowed, at least if the particle number $\mathrm {A}$ is large enough¹¹ potentially leading to renormalization parts with more external legs and thus requiring counterterms of higher ranks. Given that k-body contact terms are Pauli blocked beyond $k=2$ ($k=4$) for $g=2$ ($g=4$), this situation would require derivative terms. While being speculative, this situation cannot be excluded at this point as there exists no general proof, to the best of our knowledge, that it is not the case.

Starting from a bare $C_0$ and conjecturing that $\delta C_0(\Lambda )$ counterterms are sufficient to obtain UV-finite $\mathcal {G}^{(0,k)}_{n}$ in pure neutron systems for any value of k and n, this property would remain true for any $\mathcal {G}^{(\mathrm {A},k)}_n$ thanks to the procedure discussed in Sect. 4.3. In terms of the subtracting operator $t_\gamma$, the conjecture consists in assuming that there exists an operator with $t_\gamma \ne 0$ only for two-body diagrams $\gamma$. This would mean that the set of renormalization parts could be restricted to two-body sub-diagrams such that the renormalization in the two-body sector would be sufficient to ensure it for any perturbative approximation and any particle number. As mentioned above, and to the best of our knowledge, this conjecture remains, however, unproven. Proving it would be an important step forward in the understanding of $\pi \!\!{/}\,$EFT and its application to all $\mathrm {A}$-body sectors.

5 Application and extension

In this section, the procedure described in Sect. 4.3 is applied to the calculation of the in-medium one-body propagator within the random phase approximation (RPA). The celebrated RPA truncation scheme acts as an example of practical interest for many-body calculations. Then, an extension of the general procedure is discussed in connection with the use of a more general partitioning of H, namely the Hartree–Fock (HF) partitioning, to define the perturbation theory and the associated one-body unperturbed propagators.

5.1 Random phase approximation

Historically, the RPA was first introduced in Ref. [47] as a way to deal with collective phenomena such as charge screening effect in the electronic gas. Later, it was reformulated in Refs. [48, 49] as a particular resummation of perturbation theory diagrams, namely forward and backward particle–hole excitations. Regarding nuclear systems, RPA and its extensions play an important role to tackle collective excitations [50].

Here the Hamiltonian is partitioned as in Eq. (5) and the reference state is chosen to be $| \Phi ^\mathrm {A}_0 \rangle$ as defined in Eq. (6). The RPA is formulated as an approximation to the self-energy in the Dyson equation, namely

$$\begin{aligned}&G^{(\mathrm {A},1)\text {RPA}}_{\mu \nu }(\omega ) = G^{(\mathrm {A},1)0}_{\mu \nu }(\omega ) \nonumber \\&\quad + \sum _{\lambda _1 \lambda _2} G^{(\mathrm {A},1)\text {RPA}}_{\mu \lambda _1}(\omega ) \ \Sigma ^{(\mathrm {A})\text {RPA}}_{\lambda _1\lambda _2}(\omega ) \, G^{(\mathrm {A},1)0}_{\lambda _2\nu }(\omega ), \end{aligned}$$

(39)

where $G^{(\mathrm {A},1)\text {RPA}}_{\mu \nu }(\omega )$ and $\Sigma ^{(\mathrm {A})\text {RPA}}_{\lambda _1\lambda _2}(\omega )$ denote, respectively the one-body propagator and the self-energy in the RPA, and the unperturbed propagator reads

$$\begin{aligned} iG^{(\mathrm {A},1)0}_{\mu \nu }(\omega ) = iG^{(\mathrm {A},1)0+}_{\mu \nu }(\omega ) + iG^{(\mathrm {A},1)0-}_{\mu \nu }(\omega ). \end{aligned}$$

(40)

The 1PI time-unordered Feynman diagrams contributing to the self-energy in this approximation consist of the so-called ring diagrams. Examples of first, second, third and fourth order (in terms of number of vertices) contributions to the self-energy are displayed in Fig. 4. Once the counterterms are correctly taken into account for the 1PI part of the one-body Green’s function, no additional UV divergences appear in the full (1PR) one-body Green’s function. Thus, in the following, we focus only on 1PI diagrams.

Fig. 4

Examples of 1PI time-unordered diagrams contributing to the one-body Green’s function in the RPA. Oriented lines denote here the unperturbed propagator (40)

Fig. 5

Example of decomposition of a time-unordered diagram (left-hand side) into a sum of time-ordered diagrams (right-hand side). In a time-unordered diagram, any line refers to the complete propagator $iG^{(\mathrm {A},1)0}$. In a time-ordered diagram, any ascending (descending) line refers to the unperturbed particle (hole) propagator $iG^{(\mathrm {A},1)0+}$ ($iG^{(\mathrm {A},1)0-}$). Time-ordered diagrams in the first two rows contain one loop made of unperturbed particle propagators whereas time-ordered diagrams in the last two rows contain zero loop made of unperturbed particle propagators. Due to conservation of momentum, diagrams with one particle and one hole external leg vanish so that they are not represented. For a general partitioning of the Hamiltonian, however, this would not be necessarily the case

Each 1PI time-unordered diagram contributing to the one-body Green’s function is decomposed in a sum of time-ordered diagrams. In a time-unordered diagram, each line refers to the complete unperturbed propagator $iG^{(\mathrm {A},1)0}$. In a time-ordered diagram, an ascending (descending) line refers to the unperturbed particle (hole) propagator $iG^{(\mathrm {A},1)0+}$ ($iG^{(\mathrm {A},1)0-}$). Example of this decomposition is represented in Fig. 5. Consequently, the RPA can be recast as a truncation on the sum of particle–hole diagrams.

Table 1 Examples of diagrams $\mathcal {G}^{(\mathrm {A},1)}_n$, with p hole lines and one-particle loop, contributing to $G^{(\mathrm {A},1)}$ in the RPA. Associated cut diagrams $\mathcal {G}^{(0,1+p)}_n$ and their superficial degree of divergence $D({\mathcal {G}}^{{(0,1+p))}}_n)$ are given

Table 2 Examples of diagrams $\mathcal {G}^{(\mathrm {A},1)}_n$, with p hole lines and no particle loop, contributing to $G^{(\mathrm {A},1)}$ in the RPA. Associated cut diagrams $\mathcal {G}^{(0,1+p)}_n$ are given

Applying the procedure designed in Sect. 4.3 to 1PI time-ordered diagrams, the set of diagrams with p hole lines belonging to $\mathcal {S}^{(\mathrm {A},1+p)}_{\text {RPA}}$ is explicitly pictured at second, third and fourth orders in Tables 1 and 2. The resulting diagrams displayed in Table 2 contain no loops and, thus, are free of any UV divergence. Those displayed in Table 1 can be considered as one-loop diagrams (made of n internal lines) contributing to the n-body Green’s function defined with respect to $\vert 0 \rangle$. With the notations of Sect. 4.3, one has

$$\begin{aligned} k&= 1, \end{aligned}$$

(41a)

$$\begin{aligned} p&= n-1. \end{aligned}$$

(41b)

Following the BPHZ procedure for a diagram $\mathcal {G}^{(0,n)}_{n}~\in ~\mathcal {S}^{(0,n)}_{\text {RPA}}$ is straightforward. As $\mathcal {G}^{(0,n)}_{n}$ contains only one loop, any potential renormalization part $\gamma$ must contain at least all n internal lines building the loop so that

$$\begin{aligned} L^\gamma&= 1, \end{aligned}$$

(42a)

$$\begin{aligned} I^\gamma&\ge n, \end{aligned}$$

(42b)

$$\begin{aligned} n^\gamma&\le n, \end{aligned}$$

(42c)

where $L^{\gamma }$ is the number of loops, $I^{\gamma }$ the number of internal lines and $n^{\gamma }$ the number of vertices of $\gamma$. Using the topological identity $L^{\gamma } = I^{\gamma } - n^{\gamma } + 1$ implies that $I^{\gamma } = n^{\gamma }$ so that

$$\begin{aligned} L^\gamma&= 1, \end{aligned}$$

(43a)

$$\begin{aligned} I^\gamma&= n, \end{aligned}$$

(43b)

$$\begin{aligned} n^\gamma&= n. \end{aligned}$$

(43c)

Eventually, the only potential renormalization part of $\mathcal {G}^{(0,n)}_{n}$ is $\mathcal {G}^{(0,n)}_{n}$ itself. From Eq. (25), the superficial degree of divergence of $\mathcal {G}^{(0,n)}_{n}$ for $n\ge 2$ reads

$$\begin{aligned} D(\mathcal {G}^{(0,n)}_{n}) = 5 - 2 n. \end{aligned}$$

(44)

The only solution to $D(\mathcal {G}^{(0,n)}_{n}) \ge 0$ is obtained for $n=2$. Consequently, there is only one UV-divergent diagram which is pictured in Fig. 6 along with the additional counterterm generated via the application of BPHZ.

Fig. 6

The only UV-divergent diagram appearing in the RPA. Its associated counterterm is added. The filled vertex represents the renormalized finite LEC $C^R_0$ while the hatched vertex represents the counterterm $\delta C_0^{\text {RPA}}(\Lambda )$

The counterterm depends on the chosen renormalization scheme. Starting from the regularization of Eq. (12), the subtraction of the UV divergence can be achieved by a zero-derivative contact counterterm $\delta C_0^{\text {RPA}}(\Lambda )$. Matching the one-loop in-vacuum two-body Green’s function to the experimental S-wave scattering length $a_0$ leads to

$$\begin{aligned} C^R_0&= \frac{4\pi }{m} a_0, \end{aligned}$$

(45a)

$$\begin{aligned} \delta C_0^{\text {RPA}}(\Lambda )&= \frac{4\pi }{m} \ \frac{2}{\pi } \ \left( \int ^{+\infty }_{0} \mathrm {d}q \ v^2_\Lambda (2q) \right) \ a_0^2. \end{aligned}$$

(45b)

One thus concludes that an additional diagram is required with the pure contact counterterm (45b) to achieve renormalization. The final set of RPA 1PI time-unordered diagrams contributing to the in-medium propagator are pictured in Fig. 7. The scheme is applicable to any A.

The conclusions regarding the counterterm needed to renormalize the RPA one-body propagator are in agreement with explicit calculations of the RPA energy [41]. In [41], calculations were made with dimensional regularization and minimal subtractions for infinite homogeneous matter. Only the second order contribution to the energy was noted to be UV divergent and has been taken care of in their specific renormalization scheme. Such a UV-divergent contribution is the counterpart of the UV-divergent diagram identified in the computation of the RPA self-energy (up to an extra hole-line connecting its two external amputated legs).

The advantage of the approach used here is that no explicit calculations are required to determine the structure of counterterms. This is crucial for its application to many-body approximations that are intrinsically numerical. For example, RPA calculations with a generic cut-off regularization are in general not explicitly computable analytically because of the momentum dependence of the interaction introduced by the regulator.¹² Even in the absence of explicit calculations, counterterms can sometimes be explicitly computed, as illustrated by the RPA with the counterterm given in (45b).

Fig. 7

Time-unordered 1PI diagrams contributing to the one-body Green’s function in the RPA with the additional counterterm derived following the procedure in Sect. 4.3

5.2 Hartree–Fock partitioning

So far, MBPT has been formulated on the basis of choosing the kinetic energy as the unperturbed Hamiltonian. In the present section, a different partitioning of the Hamiltonian is considered to illustrate the flexibility of the procedure derived in Sect. 4.3.

One can typically consider a partitioning leading to dressed one-body propagators associated with a non-trivial on-shell self-energy. The first idea is to use the Hartree–Fock (HF) self-energy, which leads to the partitioning

$$\begin{aligned} H&= H_0 + H_1, \end{aligned}$$

(46a)

$$\begin{aligned} H_0&\equiv \sum _{\mathbf {p}\sigma } e^{\text {HF}}_{\mathbf {p}\sigma }(\Lambda ) \ a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma }, \end{aligned}$$

(46b)

$$\begin{aligned} H_1&\equiv \sum _{\mathbf {p}\sigma } \left( \frac{p^2}{2m} - e^{\text {HF}}_{\mathbf {p}\sigma }(\Lambda )\right) a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma } \nonumber \\&\quad + \frac{1}{2!} \sum _{\sigma _1 \sigma _2} \sum _{\begin{array}{c} \mathbf {p}_1 \mathbf {p}_2 \mathbf {p}_1^{\, \prime }\mathbf {p}_2^{\, \prime } \end{array}} \ (2\pi )^3 \delta (\mathbf {p}_1^{\, \prime }+ \mathbf {p}_2^{\, \prime }- \mathbf {p}_1 - \mathbf {p}_2) \ C_0(\Lambda ) \nonumber \\&\qquad \qquad \qquad \qquad \qquad \times v_\Lambda (\mathbf {p}_1 - \mathbf {p}_1^{\, \prime }) \ a^\dagger _{\mathbf {p}_1^{\, \prime }\sigma _1} a^\dagger _{\mathbf {p}_2^{\, \prime }\sigma _2} a_{\mathbf {p}_1 \sigma _1} a_{\mathbf {p}_2 \sigma _2}, \nonumber \\ \end{aligned}$$

(46c)

where

$$\begin{aligned} e^{\text {HF}}_{\mathbf {p}\sigma }(\Lambda ) \equiv \frac{p^2}{2m} + \Sigma ^{(\rho )\text {HF}}_{\mathbf {p}\sigma }(\Lambda ) \end{aligned}$$

(47)

denotes HF single-particle energies. Note that in Eqs. (46) the chosen regularization is different from Eq. (13). As the regularization must not impact many-body calculations, it should be chosen based on convenience. In this case, the regularized potential only depends on momentum transfer $\mathbf {p}_1^{\, \prime }- \mathbf {p}_1$.

Considering homogeneous neutron matter at density $\rho$ such that the reference Slater determinant reads

$$\begin{aligned} \vert \Phi ^\rho _0 \rangle \equiv \prod _{\sigma , p < k_F} a^\dagger _{\mathbf {p}\sigma } \vert 0 \rangle , \end{aligned}$$

(48)

with $k_F \equiv \left( 3\pi ^2 \rho \right) ^{1/3}$ the associated Fermi momentum, the HF self-energy is given by

$$\begin{aligned} \Sigma ^{(\rho )\text {HF}}_{\mathbf {p}\sigma }(\Lambda ) = C_0(\Lambda )\left( \rho - \int _{\vert \mathbf {p}^{\, \prime }\vert < k_F} \frac{\mathrm {d}^3\mathbf {p}^{\, \prime }}{(2\pi )^3} \ v_\Lambda (\mathbf {p} - \mathbf {p}^{\, \prime }) \right) ,\nonumber \\ \end{aligned}$$

(49)

leading to the HF in-medium propagator

$$\begin{aligned} G^{(\rho ,1)\text {HF}}_{\mathbf {p}\sigma }(\omega ;\Lambda ) = \frac{\theta \left( p - k_F \right) }{\omega - e^{\text {HF}}_{\mathbf {p}\sigma }(\Lambda ) + i\eta } + \frac{\theta \left( k_F - p \right) }{\omega - e^{\text {HF}}_{\mathbf {p}\sigma }(\Lambda ) - i\eta }. \end{aligned}$$

(50)

The difficulty encountered in the above scheme is that the renormalization is inherently a self-consistent problem, i.e. the Hamiltonian obtained after regularization and renormalization characterized by $C_0(\Lambda )$ depends on the set of diagrams eventually included in the post-HF step of the calculation, i.e. on the chosen many-body approximation, whose renormalization itself depends on the UV behavior of the one-body propagator defined through the HF partitioning. This problem is beyond the task treated in the present work. It is actually part of a greater problem in SCGF theories [51,52,53]. To appreciate the extent of the problem, one can naively insert into Eq. (49) the $C_0(\Lambda )$ coupling needed to renormalize typical many-body approximations on the basis of the original partitioning, e.g.:

• Employing $\pi \!\!{/}\,$EFT at LO [6], the coupling constant is renormalized exactly so that typically $C_0(\Lambda ) \propto \Lambda ^{-1}$. From Eq. (49) clearly $\Sigma ^{(\rho )\text {HF}}_{\mathbf {p}\sigma }(\Lambda )$ converges to 0 and the unperturbed propagator to the one at play in the free kinetic energy partitioning, hence losing the potential benefits of using a different partitioning.

• Renormalizing for the RPA, one has a divergent coupling constant. Typically $\delta C^{\text {RPA}}_0(\Lambda ) \propto \Lambda$ such that $\Sigma ^{(\rho )\text {HF}}_{\mathbf {p}\sigma }(\Lambda )$ linearly diverges when $\Lambda \rightarrow +\infty$.

These two examples demonstrate that achieving a self-consistent renormalization is not straightforward and requires some attention in the future.

One way to circumvent such a difficulty relies on a slightly different partitioning that consists of setting aside the counterterms in the definition of the (pseudo) HF one-body self-energy, thus leading to

$$\begin{aligned} H&= H_0 + H_1, \end{aligned}$$

(51a)

$$\begin{aligned} H_0&\equiv \sum _{\mathbf {p}\sigma } e_{\mathbf {p}\sigma }(\Lambda ) \ a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma }, \end{aligned}$$

(51b)

$$\begin{aligned} H_1&\equiv \sum _{\mathbf {p}\sigma } \left( \frac{p^2}{2m} - e_{\mathbf {p}\sigma }(\Lambda )\right) a^\dagger _{\mathbf {p} \sigma } a_{\mathbf {p} \sigma } \nonumber \\&\quad + \frac{1}{2!} \sum _{\sigma _1 \sigma _2} \sum _{\begin{array}{c} \mathbf {p}_1 \mathbf {p}_2 \\ \mathbf {p}_1^{\, \prime }\mathbf {p}_2^{\, \prime } \end{array}} \ (2\pi )^3 \delta (\mathbf {p}_1^{\, \prime }+ \mathbf {p}_2^{\, \prime }- \mathbf {p}_1 - \mathbf {p}_2) \ C_0(\Lambda ) \nonumber \\&\qquad \qquad \qquad \qquad \, \times v_\Lambda (\mathbf {p}_1 - \mathbf {p}_1^{\, \prime }) \ a^\dagger _{\mathbf {p}_1^{\, \prime }\sigma _1} a^\dagger _{\mathbf {p}_2^{\, \prime }\sigma _2} a_{\mathbf {p}_1 \sigma _1} a_{\mathbf {p}_2 \sigma _2}, \nonumber \\ \end{aligned}$$

(51c)

where

$$\begin{aligned} e_{\mathbf {p}\sigma }(\Lambda ) \equiv \frac{p^2}{2m} + \Sigma ^{(\rho )}_{\mathbf {p}\sigma }(\Lambda ) \end{aligned}$$

(52)

with

$$\begin{aligned} \Sigma ^{(\rho )}_{\mathbf {p}\sigma }(\Lambda ) = C_0^R\left( \rho - \int _{\vert \mathbf {p}^{\, \prime }\vert < k_F} \frac{\mathrm {d}^3\mathbf {p}^{\, \prime }}{(2\pi )^3} \ v_\Lambda (\mathbf {p} - \mathbf {p}^{\, \prime }) \right) . \end{aligned}$$

(53)

Because $\Sigma ^{(\rho )}_{\mathbf {p}\sigma }(\Lambda )$ is now defined in terms of $C_0^R$ rather than $C_0(\Lambda )$, it converges for $\Lambda \rightarrow +\infty$ toward

$$\begin{aligned} \Sigma ^{(\rho )}_{\mathbf {p}\sigma } = \rho \left( 1-\frac{1}{2}\right) C_0^R. \end{aligned}$$

(54)

Thus, the unperturbed particle propagator converges toward

$$\begin{aligned} G^{(\rho ,1)+}_{\mathbf {p}\sigma }(\omega ) = \frac{\theta \left( p - k_F \right) }{\omega - (\frac{p^2}{2m} + \rho \left( 1-\frac{1}{2}\right) C_0^R) + i\eta } \end{aligned}$$

(55)

and its asymptotic coefficients satisfy¹³

$$\begin{aligned} \alpha \left( S\right) = {\left\{ \begin{array}{ll} -1 &{}\text { if } S = \{ \mathbf {e}_\omega \} ,\\ -2 &{}\text { if } S = \{ \mathbf {L} \} \text { with } \mathbf {L} \notin \{\mathbf {e}_\omega \} ,\\ -2 &{}\text { if } \dim {S} \ge 2, \end{array}\right. } \end{aligned}$$

(56)

such that the procedure described in Sect. 4.3 can be applied to identify the $C_0(\Lambda )$ needed to achieve renormalization on the basis of the post-HF step of the calculation, i.e. on the chosen many-body approximation.

6 Conclusions

Starting from the derivation of a Hamiltonian H describing the interaction between nucleons within an EFT approach (namely $\pi \!\!{/}\,$EFT), the power-counting rules proposed to compute observables at LO require one to exactly solve the $\mathrm {A}$-body Schrödinger equation for the truncated Hamiltonian $H^{\text {LO}}$ in such a way that renormalization is ensured. However, exact calculations remain intractable in large $\mathrm {A}$-body sectors ($\mathrm {A} \gg 10$), which may compromise the renormalization invariance of computed observables. In order to overcome this tension in the context of many-body methods that can be formulated in terms of time-ordered diagrams, the idea pursued in this article is to design a method to renormalize workable truncations and check whether or not this renormalization procedure is consistent with the one traditionally employed in $\pi \!\!{/}\,$EFT via the exact calculation of few-body systems. The answer to this question depends on the truncation scheme employed.

In this article, a method is thus designed to identify divergences for any given set of many-body diagrams generated through a perturbative expansion of the k-body Green’s function around an $\mathrm {A}$-body Slater determinant reference state. This set can be strictly perturbative or eventually correspond to the resummation of an infinite (subset) of diagrams, e.g. summing particle-particle or particle–hole ladders. The method involves a so-called ’cutting’ procedure allowing one to relate the ultraviolet divergences of the in-medium k-body Green’s function to those displayed by a set of in-vacuum$(k+p)$-body Green’s functions. Applying the BPHZ procedure to the diagrams making up the in-vacuum Green’s functions, counterterms necessary to renormalize the original in-medium k-body Green’s function are identified in a systematic fashion. This procedure delivers the desired property that k-body counterterms are independent of the $\mathrm {A}$-body sector ($\mathrm {A}\ge k$) one starts from. Eventually, the present development is similar to what has been done in QFT at finite temperature [43, 54].

This work only constitutes a first step forward and critical extensions remain to be carried out. First, one must go from perturbative to intrinsically non-perturbative methods based on in-medium diagrams such as CC, SCGF or IM-SRG. Regarding SCGF, let us mention the important work showing possible additional counterterms compared to the naive application of BPHZ to diagrams formulated in terms of fully-dressed propagators [51,52,53]. Second, additional steps are needed to extend the present developments to more general partitioning. Of importance are partitionings based on unperturbed Hamiltonians breaking exact symmetries of the Hamiltonian. The most trivial case consists of breaking translational invariance such that unperturbed propagators are no longer diagonal in momentum space, e.g. using a harmonic oscillator Hamiltonian. The use of Weinberg’s asymptotic theorem must be extended in such a case. Less trivial are partitionings breaking U(1) symmetry that are employed to tackle the superfluid character of nuclear matter and open-shell nuclei, e.g. MBPT [38, 55, 56], CC [57] or SCGF [31, 58] using a Bogoliubov vacuum as reference state. The associated diagrammatic relies on the use of anomalous propagators in addition to normal propagators. As a result, the analysis of ultraviolet divergencies is fundamentally different so that diagrams contributing to the mean-field (i.e. Hartree, Fock and Bogoliubov) already require counterterms [59, 60] contrary to those (i.e. Hartree and Fock) at play in the perturbation theory considered here. Let us mention some interesting work along this direction such as [53] in the case of SCGF applied to relativistic scalar field theories as well as [40] in the case of a dilute non-relativistic Fermi gas in a conventional BCS superfluid phase.