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

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 observable. In this paper, a renormalization procedure is designed in connection with many-body approximations applicable to large-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 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.


Introduction
The problem of describing accurately and systematically nuclear systems from their A constitutive protons and neutrons is now almost a century old. Starting with the discovery of the neutron in the 1930s, the nuclear manybody 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 behaviour, 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 and Ref. [6] for a recent review). 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 (/ πEFT) [7] at leading order (LO) is used to demonstrate how the procedure of renormalization of the nuclear Hamiltonian affects the calculations 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-A sectors. This feature must be explicitly considered when applying EFT to arbitrary a Present address : Department of Physics, University of Surrey, Guildford GU2 7XH, United Kingdom. Mail : m.drissi@surrey.ac.uk nuclear systems and requiring an order-by-order renormalizability. In Sect. 2, the renormalizability issue and its consequences on 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 [8] and Bogoliubov-Parasiuk-Hepp-Zimmermann's (BPHZ) theorem [9,10,11,12,13,14] 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 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.

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.

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 Ψ A m and their associated eigenvalues E A m in all A-body sectors of interest, i.e. the goal is to solve the Schrödinger equation where m indexes the set of 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 (χEFT) [15,2,4,3,16] such that it takes the form of a series where the leading-order (LO) and the sub-leading orders (SLOs) are organised according to a set of power-counting (PC) rules. First to be proposed historically, Weinberg's power counting [2,3,16] happens to fit traditional manybody 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 A m . However, Weinberg's PC has since been shown to violate the demand that the EFT is (order-by-order) renormalizable [17] and alternative PCs have been proposed [18]. In addition to modifying the order at which certain contributions enter the Hamiltonian, new PCs 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 / π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, / πEFT has been shown analytically to satisfy renormalizability at LO up to three-body systems [19] and numerically in the case of four-body systems [20].

Renormalization and many-body approximations
To this day, solving exactly Eq. (1) for H LO / π remains in general numerically intractable for A 10. Consequently, one must design an additional expansion and truncation when attempting to solve Typical truncations applicable to A-body systems with A 10 are nowadays implemented on the basis of nonperturbative self-consistent Green's function (SCGF) [21,22,23], coupled cluster (CC) [24,25] and in-medium similarity renormalization group (IM-SRG) [26,27] methods but also on the basis of MBPT [28,29,30].
Traditionally, and in agreement with power-counting rules of / πEFT, H LO / π 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 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 / π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 LO / π via the exact solution of Eq. (3) in two-and three-body systems extends to large-A systems, one attempts to design renormalization prescriptions for a given many-body approximation, many-body observables thus being renormalizationinvariant by construction. While this approach departs from the original scheme pursued in / π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.

Consider the Hamiltonian
containing, for convenience, a sole spin-independent twobody interaction parametrised 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 / π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 consider at LO in / πEFT.

Perturbation theory
Approximations presently considered are formulated within MBPT on the basis of the particular partitioning of H 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 A-body sector, the unperturbed reference state, i.e. the ground state of H 0 , is given by the Slater determinant where |0 denotes the particle vacuum. Here, the index i is a shorthand notation to label hole states (σ i , p i ), i.e. one-body states that are occupied in Φ A 0 . 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 1 defined as where T denotes the time-ordering operator, creation and annihilation operators are in the Heisenberg picture, and Ψ A 0 is the exact ground-state of H connected to Φ A 0 via the adiabatic theorem of Gell-Mann and Low [31]. Exact Green's functions can be themselves expressed as a sum over the complete set of linked diagrams G (A,k) n carrying k incoming and k outgoing external lines i.e.
For example, the ground-state energy A-body energy is obtained from G (A,1) via the so-called Galitskii-Koltun sum rule. Based on an appropriate choice of the ordering of its time labels, the Lehmann's representation of the k-body Green's function (7) further accesses the eigenenergies of all systems whose number of particles is comprised between A−k and A+k. In practice, one has k A such that accessible observables are associated to states lying in the neighbourhood of the A-body sector used to set the reference state.
where S (A,k) Exact denotes the complete set of diagrams contributing to the k-body Green's function, n the number of interaction vertices and A 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 as where the limit η → 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 G Exact , the associated amplitude reads as 2πδ (±ω λ · · · − ω µ · · · + ω ν . . . ) where σ is an integer, λ, µ and ν denote generic singleparticle labels for internal, external outgoing and external incoming lines, respectively. Additionally, I + (I − ) denotes the set of internal particle (hole) lines, E + in (E − in ) the set of external incoming particle (hole) lines and E + out (E − out ) the set of external outgoing particle (hole) lines.
Many-body approximations considered here consist of selecting a subset S Exact so that the approximated k-body Green's function reads as

Regularization
When considering local interactions, as those derived in any EFT, the amplitude A G (A,k) µ1...µ k ν1...ν k contains, in general, ultraviolet (UV) divergences requiring the introduction of a regularisation to suppress the high-momentum modes. In this work, the regularisation is introduced via a momentum regulator v Λ (q) satisfying where Λ represents the characteristic cut-off scale beyond which high-momentum modes are suppressed and q denotes the relative momentum between two nucleons. Eventually, one is interested in the limit where Λ Q with Q the scale characterising the low-energy observables of interest. The limit of large Λ 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 regularised 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 Λ and in particular on the cut-off Λ. Eventually, the Hamiltonian to be actually considered reads as where the dependence of the LEC C 0 (Λ) on the regularisation is made explicit and where incoming and outgoing relative momenta are respectively defined as To make the cancellation of UV divergences explicit, LECs are usually decomposed, for a given renormalization scheme, into a renormalized component and a counterterm [32], i.e. in the case of (13) as where C R 0 is the finite renormalized LEC and δC 0 (Λ) encapsulates all counterterms that will be necessary to cancel UV divergences. For a given set of diagrams S (A,k) MB with renormalized vertices C R 0 , the aim is thus to derive an additional set of diagrams S (A,k) MB,ct with renormalized and counterterm vertices (to be derived as well) such that the total sum converges to a finite value when Λ Q. In the perturbation theory employing the particle vacuum |0 as a reference state, the counterterms needed to achieve renormalization can be identified systematically via the application of the BPHZ theorem [9,10,11,12,13,14]. Applying a similar procedure to time-ordered diagrams stemming from a perturbative expansion around an inmedium reference state Φ A 0 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 A given that the particle and hole propagators do depend on A. In numerical calculations, obervables used to renormalize the LECs would have to be computed for each A. This situation would penalise both 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 A ≥ k.

Basic tools
Weinberg's asymptotic theorem [8] and the BPHZ [9,10,11,12,13,14] theorem represent key tools to develop the renormalization procedure exposed in Sec. 4 below. Consequently, the two theorems are briefly recalled in the present section.

Weinberg's Asymptotic theorem
The UV convergence of the amplitude of a time-ordered diagram is analysed 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 on Weinberg's asymptotic theorem, see Ref. [8].

Asymptotic coefficients and multidimensional integrals
Let us first introduce the definition of asymptotic coefficients (provided they exist) of a function f : R n → C as given in [8]. For any vector subspace 2 S = {L 1 , . . . , L m } of R n with m ≤ n and L 1 , . . . , L m being m independent R nvectors, and any compact region W ⊂ R n , the asymptotic coefficients are defined as the numbers α ({L 1 , . . . , L r }) and if η 1 . . . η m go independently to infinity. The asymptotic coefficients α (S) and β (S) can be interpreted as the asymptotic coefficients α ({L}) and β ({L}) for L a 'typical' vector in S i.e. fixing η 1 . . . η m−1 sufficiently large and when η m goes to infinity. Considering now integrals of a function f , the integration along the directions L 1 , . . . , L r is defined as where X is a vector in R n . Thanks to Fubini's theorem, if f L1,...,Lr (X) exists (in the sense that the integral is absolutely convergent), the integration depends only on the vector space I = {L 1 , . . . , L r } so that one defines Furthermore, f I (X) depends only on the projection of X along I 3 . Choosing a subspace E such that R n = I ⊕ E, the domain of definition of the function f I (X) can be restricted to E.
In the case of the amplitude of a time-ordered diagram, the general integrand depends on (ω 1 , p 1 , . . . ) 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 4 . The general vector space R n = I ⊕ E denotes the vector space of all (one-body) energies and momenta (internal and external).
As an example, the asymptotic coefficients α (S) of the in-vacuum propagator are now extracted. The in-vacuum propagator iG (ω) is interpreted as a multivariate function on the vector space One can show that, in this case, the asymptotic coefficients of the in-vacuum propagator α 0 read as

Asymptotic theorem
With all the notations introduced before, the general asymptotic theorem follows If a function f (X) possesses asymptotic coefficients α (S) , β (S) for any non-null subspace S ⊂ R n , if f (X) is integrable for any finite region of R n and if D I < 0 where then f I (X) exists i.e. is absolutely convergent.
To apply this theorem to a time-ordered diagram G, it is reformulated in terms of sub-diagrams as detailed in Ref. [8]. This is done by associating to any sub-space of integration S ⊆ I a sub-diagram γ ⊆ G made of internal lines. In particular S = I corresponds to the sub-diagram γ made of G itself without its external lines. The quantity α (S ) + dim S corresponds to the superficial degree of divergence of the associated sub-diagram γ of G and is denoted D(γ). For a diagram G (0,k) n where lines denote the unperturbed in-vacuum propagator (21), n denotes the number of vertices (5c) and k is the number of incoming (and of outgoing) external lines, 4 Finite sum on spin-isospin indices are omitted here as they introduce only a finite linear combination of integrals on energies and momenta so that the conclusion on the UV behaviour is not impacted. 5 The label σ is just considered as a fixed parameter so that it is dropped in the definition of f (ωeω + pxep x + pyep y + pzep z ).
Having D I < 0 is thus equivalent to having D(γ) < 0 for all γ ⊆ G made of internal lines. Consequently, from Weinberg's asymptotic theorem, the well-known powercounting theorem follows. The amplitude associated to G is finite if D(γ) < 0 for any sub-diagram γ ⊆ G made of internal lines of 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 α (S) associated to the propagator. This particular property is fundamental for the development of the renormalization prescription exposed in Sect. 4.

Subtractions and BPHZ theorem
Whilst the asymptotic theorem allows to prove the convergence of the calculations with adequate counterterms, BPHZ theorem offers a systematic way to generate sufficient counterterms to lower the superficial degree of divergence of the internal sub-diagrams γ ⊆ G. The procedure to do so is now briefly detailed.

Definitions
A sub-diagram γ of G is defined as a subset of lines and vertices contained in G where end points of the lines of γ belong to its vertices. A diagram G 1 is said to be included in G 2 and denoted as G 1 ⊆ G 2 if their set of lines verify the same inclusion relation. In particular, a sub-diagram γ of G verifies γ ⊆ G. The sub-diagram generated by the intersection of lines of two sub-diagrams γ 1 and γ 2 defines a sub-diagram γ and is denoted as Two sub-diagrams γ 1 and γ 2 that have neither lines nor vertices in common are said to be disjoint and the result is denoted as If neither γ 1 ⊆ γ 2 nor γ 2 ⊆ γ 1 and γ 1 ∩ γ 2 = ∅ they are said to be overlapping. Otherwise they are said to be non-overlapping. For a set of non-overlapping sub-diagrams γ 1 , γ 2 , . . . , γ n of G, the reduced diagram G \ {γ 1 , γ 2 , . . . , γ n } is defined by the diagram resulting from G after contracting all lines of γ 1 , γ 2 , . . . , γ n to a point. The amplitude associated to G is denoted as A G . Given a set of mutually disjoint sub-diagrams γ 1 , γ 2 , . . . , γ n of G, the corresponding amplitude is expressed as the product of the amplitudes A γj of the sub-diagrams while the remainder is denoted as A G\{γ1,...,γn} such that A sub-diagram γ of 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 D(γ) ≥ 0 .

Recursive subtractions of UV divergences
The BPHZ procedure defines recursively a renormalized amplitude R G associated to the diagram G [9,10]. If the amplitude associated to G is convergent to begin with, then 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 where t G is the operator of the Taylor expansion with respect to the external momenta 6 around 0 up to the order of the superficial degree of divergence D(γ) of the diagram, i.e.
If G is superficially divergent and contains divergent subdiagrams, the renormalized amplitude is defined recursively as whereR G corresponds to the amplitude where all subdivergences have already been subtracted. The subtraction by the Taylor operator t γ corresponds to adding the amplitude associated to a diagram where the divergent subdiagram γ has been replaced by a so-called counterterm vertex. Here t γ corresponds to a zero momentum subtraction of UV divergences. Different subtractions can be used by modifying the definition of t γ [32].

Forest formula
The recursion relation (33) was solved explicitly by the forest formula [13] that 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 7 , for the order relation ⊆ on the mutually non-overlapping subdiagrams, 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 6 The renormalization point is chosen here to be 0 but could be chosen arbitrarily. 7 A Hasse diagram associated to an order relation is a diagrammatic representation of the ordering between the objects considered. 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. An i-forest is restricted if each of its boxes contains only renormalization parts. To each restricted i-forest F one associates again an amplitude, namely 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 G where the sub-diagrams γ 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 F is pictured with boxes will represent directly the amplitude Ω F . See Fig. 2 for an example based on a one-loop diagram. Eventually, the forest formula states that the renormalized amplitude of the diagram G is given by the sum over all restricted i-forests, i.e.
where F R (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 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.

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 analyse UV divergencies and derive sufficient counterterms whether the reference state is |0 or Φ A 0 [33]. However, as argued in Sect. 2.3, the direct use of BPHZ will generate 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 Φ A 0 as a reference state) to UV divergences occurring in the calculation of in-vacuum Green's functions (i.e. with |0 as a reference state) in a related approximation.

Cutting procedure
Let G (A,k) n ∈ S (A,k) MB and A G (A,k) n be its associated amplitude given in Eq. (10). For any finite A, the set of states (p i , σ i ) labelling unperturbed hole propagators lie in a compact space 8 . 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 behaviour of G (A,k) n is the same as the UV behaviour of an associated diagram made only of particle propagators. Such diagram can be defined as the one made out of the same n vertices but with incoming particle lines E + in ∪ I − , outgoing particle lines E + out ∪ I − and internal particle lines I + . As each line in the aforementioned diagram corresponds to an unperturbed particle propagator iG (A,1)0+ and possesses (k + p) outgoing and incoming external lines, where p ≡ #I − , the diagram is denoted in the following as G original diagram by cutting all internal unperturbed hole lines in G (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 behaviour of a diagram contributing to the in-medium k-body Green's function i k G (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 (A,k+p) .

Relation to in-vacuum diagrams
At this stage if BPHZ is applied to the diagrams G (A,k+p) n resulting from the cutting procedure, A-dependent counterterms will be generated. To avoid this, each diagram G (A,k+p) n is associated to a diagram G (0,k+p) n with the same UV behaviour but that is made only of in-vacuum propagators. This is done by noticing that the unperturbed particle one-body Green's function iG (A,1)0+ possesses the same asymptotic coefficients α (S) as the in-vacuum onebody Green's function iG (0,1)0 , i.e. for any sub-space S of {e ω , e px , e py , e pz } where α + (S) correspond to the asymptotic coefficients of iG (A,1)0+ pσ (ω). 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 G (0,k+p) n associated to G (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 G (0,k+p) n reads as In order to transport the cancellation of UV divergences back to G (A,k+p) n , one introduces the amplitude where Ω (A)F denotes the amplitude associated to the same diagram as the one originating from Ω F except that lines denoting in-vacuum propagators are replaced by unperturbed particle propagators. Therefore, the amplitudes Ω (A)F contains k-body counterterms as in Ω F which are, by construction, the same for all A ≥ k.

General procedure
The procedure to derive UV-finite k-body Green's functions computed with respect to Φ A 0 on the basis of a many-body approximation defined by a truncated set of particle-hole diagrams S

Discussion
Section 4 introduced a general procedure to transpose a renormalization scheme formulated for a perturbative expansion around the particle vacuum |0 to a perturbative expansion using an A-body reference state Φ A 0 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, 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 (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 LO / π is inactive, renormalization was shown to be achievable with only δC 0 (Λ) counterterms on the basis of exact (no approximation) vacuum twoand three-body Green's functions [19]. If both neutrons and protons are present, the D 0 (Λ) 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 [20].
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 details for the dilute fermi gas in Ref. [34], a sole 9 δC 0 (Λ) counterterm is needed to achieve renormalization up to 4 th order in perturbation theory for g = 2 whereas an additional D 0 (Λ) three-body contact term is necessary for g = 4. Going to higher orders, diagrams with a larger number of particle lines propagating simultaneously are allowed, at least if the particle number A is large enough 10 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 our knowledge, that it is not the case.
Starting from a bare C 0 and conjecturing that δC 0 (Λ) counterterms are sufficient to obtain UV-finite G (0,k) n in pure neutron systems for any value of k and n, this property would remain true for any G (A,k) n thanks to the procedure discussed in Sect. 4.3. In terms of the subtracting operator t γ , the conjecture consists in assuming that there exists an operator with t γ = 0 only for two-body diagrams γ. 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 / πEFT and its application to all A-body sectors.

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 well-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.

Random phase approximation
Historically, the RPA was first introduced in Ref. [35] as a way to deal with collective phenomena such as charge screening effect in the electronic gas. Later, it was reformulated in Refs. [36,37] 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 [38].
Here the Hamiltonian is partitioned as in Eqs. (5) and the reference state is chosen to be |Φ A 0 as defined in Eq. (6). The RPA is formulated as an approximation to the selfenergy in the Dyson equation, namely where G (ω) denote respectively the one-body propagator and the self-energy in the RPA approximation, and the unperturbed propagator reads as 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. Each 1PI time-unordered diagram contributing to the one-body Green's function is decomposed in a sum of timeordered diagrams. In a time-unordered diagram, each line refers to the complete unperturbed propagator iG (A,1)0 . In a time-ordered diagram, an ascending (descending) line refers to the unperturbed particle (hole) propagator iG (A,1)0+ (iG (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.
Applying the procedure designed in Sect. 4.3 to 1PI time-ordered diagrams, the set of diagrams with p hole lines belonging to S In a time-ordered diagram, any ascending (descending) line refers to the unperturbed particle (hole) propagator iG (A,1)0+ (iG (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.
displayed in Table 2 contain no loop and, thus, are free of any UV-divergence. Those displayed in Table 1 read as one-loop diagrams (made of n internal lines) contributing to the n-body Green's function defined with respect to |0 . With the notations of Sect. 4.3, one has Following the BPHZ procedure for a diagram G  , with p hole lines and one particle loop, contributing to G (A,1) in the RPA approximation. Associated cut diagrams G n internal lines building the loop so that where L γ is the number of loops, I γ the number of internal lines and n γ the number of vertices of γ. Using the topological identity L γ = I γ − n γ + 1 implies that I γ = n γ n G , with p hole lines and no particle loop, contributing to G (A,1) in the RPA approximation. Associated cut diagrams G (0,1+p) n are given. so that Eventually, the only potential renormalization part of G The only solution to D(G (0,n) n ) ≥ 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.
The counterterm depends on the chosen renormalization scheme. Starting from the regularisation of Eq. (12), the subtraction of the UV-divergence can be achieved by a zero-derivative contact counterterm δC RPA 0 (Λ). Matching the one-loop in-vacuum two-body Green's function to the experimental S-wave scattering length a 0 leads to When computing the one-body Green's function within / πEFT at LO in neutron matter at the RPA level, one 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.

Hartree-Fock partitioning
the renormalization procedure of Sect. 4.3 can be applied. In the case of Eq. (54), using Eq. (12b) is sufficient for Eqs. (55) to hold. For derivative contact interactions it becomes necessary to take v Λ decreasing sufficiently fast, e.g. exponentially.

Conclusions
Starting from the derivation of a Hamiltonian H describing the interaction between nucleons within an EFT approach (namely / πEFT), the power-counting rules proposed to compute observables at LO require to exactly solve the A-body Schrödinger equation for the truncated Hamiltonian H LO in such a way that renormalization is ensured. However, exact calculations remain intractable in large A-body sectors (A 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 / π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 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 A-body sector (A ≥ k) one starts from. Eventually, the present development is similar to what has been done in QFT at finite temperature [32,39].
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 [40,41,42]. Furthermore, additional steps are needed to extend the present developments to more general partitioning. Second, 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 [43,30,44], CC [45] or SCGF [46,23] 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 [47,48] contrary to those (i.e. Hartree and Fock) at play in the perturbation theory considered here.