Fourth post-Newtonian Hamiltonian dynamics of two-body systems from an effective field theory approach

We calculate the motion of binary mass systems in gravity up to the fourth post--Newtonian order. We use momentum expansions within an effective field theory approach based on Feynman amplitudes in harmonic coordinates by applying dimensional regularization. We construct the canonical transformations to ADM coordinates and to effective one body theory (EOB) to compare with other approaches. We show that intermediate poles in the dimensional regularization parameter $\varepsilon$ vanish in the observables and the classical theory is not renormalized. The results are illustrated for a series of observables for which we agree with the literature.


Introduction
Precise predictions for observables describing the merging of two heavy astrophysical objects like black holes or neutron stars are very important [1]. In particular the spiral-in phase of these processes can be described by analytic calculations. Methods of non-relativistic effective field theory [2][3][4][5][6][7][8][9] allow to calculate the equations of motion of a binary mass system within the post-Newtonian (PN) approach. We will limit ourselves to the non-spinning case in the following. Until now the corrections have been computed to the 4th post-Newtonian order [2,4,5,8,9]. There are first corrections due to the static potential at 5PN [6,7], see also [10]. The corrections up to 4PN were derived by using other methods, cf. Refs. [11,12], before. Very recently first partial results up to the 5th post-Newtonian order have been obtained in [13]. In the post-Newtonian approach one retains all terms in the velocity being of the same order as the pure potential terms by the virial theorem [14]. Recently also important progress has been made in the post-Minkowskian (PM) approach reaching 3PM cf. [15,16], see also [17].
In this paper we present the corrections up to 4PN obtained using effective field theory methods. Here we apply dimensional regularization in calculating the Feynman integrals. At intermediary steps the dimensional parameter D = 4 − 2ε occurs. We perform multiple comparisons to results in the literature and find agreement. Starting at 3PN singular contributions of O(1/ε) occur when working in harmonic coordinates. From 4PN onward singularities of O(1/ε) also result from the tail terms, cf. also [5,8]. By applying canonical transformations one may map the Hamiltonian in harmonic coordinates into a class of pole-free Hamiltonians at 3PN, to which the ADM and EOB Hamiltonians belong. This is also the case at 4PN when accounting for the tail terms before. Therefore, the classical observables are free of the intermediary regularization parameter.
The paper is organized as follows. In Section 2 we describe the calculation method. We outline the general framework and present the results up to 4PN in harmonic coordinates by deriving the effective Hamiltonian using the Hamiltonian formalism [18] in the center of momentum frame. The cancellation of the pole contributions is discussed. We also present an associated pole-free Hamiltonian obtained after applying a canonical transformation. At 3PN one also obtains a ln(r) free Hamiltonian, with r = | x 1 − x 2 |. In Section 3 we construct the canonical transformation to ADM coordinates, to compare with these results. In Section 4 we study the canonical transformation to EOB coordinates, which are widely used in the literature and present the EOB Hamiltonian up to 4PN in explicit form. We present numerical results for the energy of the innermost stable orbit and the angular momentum for circular coordinates in Section 5. Section 6 contains the conclusions. In an appendix we present the general Ddimensional Hamiltonian in harmonic coordinates up to O(ε) in explicit form, which is also important for higher post-Newtonian calculations.

The general framework
The main steps of our calculation have already been described in Ref. [7] calculating the static potential to 5PN ab initio. Now we add the different velocity contributions up to 4PN. Starting from the Einstein-Hilbert Lagrangian, one parameterizes the metric g µν as proposed in [2] in terms of one scalar, a three-component vector and a six-component tensor field. Within this parameterization one can derive the contributing Feynman rules using the path integral. 1 The corresponding Feynman diagrams are generated using QGRAF [19], after providing the corresponding Feynman rules. We consider only the classical contributions and calculate the effective two-body potential. Here the necessary velocity contributions range up to O((G N M /r) k (v 2 ) n−k+1 ), with G N Newton's constant, M a mass scale, r the distance of the point masses and v ∈ {v 1 , v 2 } their velocities. We will set c = 1 in many places and keep it only when needed as an order parameter. The Lorentz algebra is carried out using Form [20] and the reduction to master integrals is performed using the code Crusher [21]. The contributing master integrals are known from the calculation of the static potential already, cf. [7] and references therein. One first obtains a Lagrange function of mth order containing also the accelerations a i and time derivatives thereof.
In Table 1 we give an overview on the complexity of the calculation up to 4PN and give the numbers of generated diagrams, the number of non factorizing diagrams, those with no world-line loops and no tadpoles by loop order. The initial number of 25324 diagrams reduces to 9266 contributing to the present result. Here we do not sort them in addition into equivalence classes and do thus determine the final combinatorial factors by the calculation itself.
The computation time for the complete project from Feynman diagram generation, IBP reduction to the final results amounts to a few hours using an Intel(R) Core(TM) i7-8650U CPU.
The present results and those from [4,8,22] are given in form of an mth order Lagrange-density, from which one may derive the associated equation of motion [23] with L the Lagrangian density and U Ω ∈ {x 1 , x 2 }, the coordinates of the two point masses and U Ω|t = ∂ t U Ω . The present results, although widely different in form and probably based on different Feynman rules, can be compared using the equation of motion (1) finding agreement with the corresponding equation of motion for the Lagrangians given in [4,8,22]. This applies to the potential contributions. In addition the tail terms have to be considered, cf. [8,24,33].
The next step is to eliminate the accelerations and their time derivatives. For this we use double zero insertions [25], up to linear terms in the acceleration. These are eliminated by a coordinate shift [25][26][27] implying the contribution of variational derivatives of the Lagrangian and a total time derivative. This procedure is applied to the D-dimensional Lagrangian for all post-Newtonian orders. The first order Lagrangian density is then obtained in somewhat different coordinates from harmonic coordinates [27]. Finally, we perform the Legendre transformation to obtain the Hamiltonian.

Remarks on the 4PN Tail Term
The quadrupole structure of the tail term at 4PN results from the leading term in the multipole expansion of the retarded radiation field from linearized Einstein gravity, cf. [28], with η µν the Minkowski metric and |g| the modulus of the determinant of the metric tensor. Let L be the multi-index L = i 1 , i 2 , ..., i l , then the components of h µν are given by whereĨ and The dots are Newton's notation for the time-derivative. I L (t) are the mass-type multipole moments in Cartesian coordinates in symmetric-trace-free (STF) tensor representation. We do not consider the current-type multipoles, since they do not contribute at the level of 4PN. The leading order tail contribution results from the scattering of the I ij (t)-generated quadrupole radiation field off the static field h 00 = − 4G c 2Ĩ0 from the monopole source l = 0, where M = m 1 + m 2 denotes the total mass of the radiating system. The source representation of I ij (t) in the case of a binary system reads in its center-of-mass system with µ = m 1 m 2 /M . More details are given in [28,34]. A derivation of the tail-term in ADMcoordinates was given in [29,30].
Since the derivation of the 4PN tail term is not thoroughly unique in the literature due to the use of different regularization schemes or combinations thereof, cf. [8,11,[31][32][33][34][35], we add a few clarifying remarks in the present treatment based on Feynman amplitudes. We follow Ref. [8], Eq. (23), performing the calculation in D dimensions. The gravitational constant is taken in D dimensions, G N → G N µ −2ε 1 , with a mass scale µ 1 which reappears in and γ E denotes the Euler-Mascheroni constant. 2 Adjusting the expression from [8] for the different choice of regulator −2ε = ε FS the representation for the tail contribution to the action S is given by where Q(k 0 ) denotes the Fourier transform of the quadrupole moment I ij To convert this expression to the time domain we have to make use of the well-known property of the Fourier transform [38] We finally obtain the following expression 3 and τ 0 is given by The square of the third time derivative of I ij evaluates to Note, that we need to keep terms O(ε) that are multiplied with the pole term in δS tail .
2 In other calculations [35] the regularization has been performed using the Hadamard symbol [36,37]. 3 Here and in the following summation over equal indexes is understood.
On circular orbits, setting p.n = 0, the τ integration in Eq. (14) can be performed analytically and one obtains with the angular frequency Ω given by with E circ the energy and J the angular momentum at the orbit. For later use we introduce the variable x In observables on circular orbits only the combination ln r r 0 + ln(2Ωτ 0 ) = ln(4Ωr) = 1 2 ln(16x) (20) appears, which unambiguously gives rise to the well-known result in the literature without the need to manually choose a particular scale. Aside from this logarithm, logarithms of the type ln(r/r 0 ) remain in the Hamiltonian in harmonic coordinates. They will be treated in Section 2.4.

The results in harmonic coordinates
In the following we list the different contributions to the Hamiltonian resulting from the first order Lagrange-function derived in Section 2.1 up to 4PN, which we write for center of mass system coordinates, normalizing the momentum p by µ. The coordinates are close to harmonic coordinates, modified due to shifts to eliminate the accelerations. We normalize the Hamiltonian H byĤ = H/M − 1, and refer to the parameter ν = µ/M . Furthermore, we use u = G N M/r. We account for the O(ε) terms up to 1PN since they play a role in the transformation to other coordinates discussed later. The full O(ε) up to 4PN is displayed in appendix A.
One obtainŝ and with p.n = p.r/r, r = | r| and ( ... I ij (t)) 2 is given in Eq. (16). The terms of first order in ν can be obtained in the Schwarzschild test-particle limit in harmonic coordinates [39,40] Both the 3PN and 4PN (effective) Hamiltonian contain pole-terms in the dimensional parameter ε and logarithmic contributions. This is not unusual since these quantities are not observables. 4

The cancellation of the pole contributions
Using the effective field theory approach of Einstein-Hilbert gravity in the post-Newtonian orders, the harmonic coordinates belong to a class in which pole terms already occur at the level of 3PN.
During the whole calculation one works in 3 − 2ε Euclidean space dimensions, starting at the Newtonian level. In various places one has to expand up to O(ε) because certain terms appear multiplicatively together with pole terms. Because of the possible diagrammatic structures, the 6PN terms can potentially contain double poles, the 9PN terms triple poles etc., which is also indicated by the structure of the canonical transformations given e.g. in (38).
On the other hand, ADM coordinates [11] and EOB coordinates [13] are pole-free at 3PN. 5 Since at 3PN a canonical transformation exists, cf. Section 2.5, transforming to pole-free coordinates, no observable will exhibit poles, and, at 3PN, associated to that, no logarithmic terms of the kind ln(r/r 0 ).
At 4PN for the complete Hamiltonian, including the tail terms, the poles arrange such, that they can be transformed away by a canonical transformation ending up with a pole-free Hamiltonian. However, logarithmic contributions remain. The latter ones are free of the dimensional mass scale µ 1 , cf. Section 2.2. Therefore, the potential renormalization group equation [41,42] of the bare Hamiltonian is trivial and no renormalization takes place; both the interacting point masses remain constant w.r.t. their rest mass and Newton's constant remains a constant, cf. also [43].

Hamiltonians in Pole-free Coordinates
Although first order bare Lagrangians and Hamiltonians of the effective field approach to classical gravity are allowed to have poles in the dimensional parameter ε, there exist classes of coordinates which are free of poles. This is well-known at 3PN and examples are the ADM and EOB coordinates. For convenience, we provide a canonical transformation from harmonic coordinates to this class in the following. These Hamiltonians have no real preference against Hamiltonians with poles, as e.g. obtained in the harmonic gauge, since all observables are pole free. In conjunction with the poles also associated logarithmic contributions emerge. These related logarithms are eliminated in the same manner, while others, stemming from the tail terms, may remain.
Canonical transformations form a group, cf. e.g. [44], and one may therefore also consider their products for some purpose. Let us consider the canonical transformation of a Hamiltonian H to a Hamiltonian H . Given the differential operator D g By using y k = q k for 1 ≤ k ≤ l and y k = p k for l + 1 ≤ k ≤ 2l, one can write the canonical transformation induced by the generator g as Lie-series [44,45] 6 in the following form, or with the linear differential operator By the aid of the fundamental relations one obtains The canonical transformations are associated to symplectic groups [47][48][49]. The generating functions for the canonical transformation ofĤ harm to the associated pole-free HamiltonianĤ polefree harm are given by The transformations read with the Poisson brackets One obtainŝ with Note that here the contributions to the Hamiltonian are used in D dimensions. Unlike the case in (51) a double commutator does not contribute to (41). It emerges for the first time at 6PN. Since canonical transformations leave the dynamics of a system invariant [50], one obtains that all observables calculated using (21)(22)(23)(24)(25) lead to the same results. Here the canonical transformations may contain also the regularization parameter ε = 0, which would formally imply singularities in the limit ε → 0 in the Hamiltonian. However, this is not relevant since all observables do not depend on this parameter. It only remains a necessary asset by applying dimensional regularization for the calculation.

ADM coordinates
We notice thatĤ N andĤ 1P N in harmonic coordinates are the same as in ADM coordinates [11]. Now we construct the canonical transformation from our results to those of Ref. [11] using the same variables at both sides. Both Hamiltonians are then equivalent and lead to the same observables, which are also frame-independent [51]. The generating functions for the canonical transformation fromĤ harm toĤ ADM are given by The Hamiltonians starting from 2PN are related bŷ

Effective one body coordinates
We will now consider the link to the Hamiltonian of effective one body (EOB) theory [52,53], which is widely used. Here we consider the effective Hamiltonian [52], Eq. (6.5), and transform from ADM to EOB coordinates. We definê Note that in [55], Eq. (6.1), α 3 = α 4 = 0 and with the functions A,D and Q given in [13], cf. [52], Eq. (6.10). Note that one also has to rescale in [13] in the root-term in (53), to be able to expand the contributions to (53) to compare the respective contributions. Here the coordinates in the r.h.s. are those of ADM and in the l.h.s. of EOB. In the following we will first consider the local contributions. Furthermore, it is convenient to consider the canonical transformation form ADM to the EOB coordinates for the quantity withĤ The different post-Newtonian contributions in both sides of (55) are now labeled by powers of 1/c 2 . We determine the canonical transformation mapping the expansion coefficients of T (k) ADM to T (k) EOB . 7 We write the canonical transformations using the respective terms of the associated Lie-series [44,45].
At the Newtonian order T T    EOB up to the 2nd post-Newtonian order, which is equivalent to the present approach at this level.
We rescale r → G N M r and obtain the generating functions g i | 4 i=1 of the following canonical transformations In this way we can directly compare to results obtained by using the effective one-body approach, which have been worked out to the second post-Newtonian order in [52], at third [54] and fourth post-Newtonian order [55], with first results at the fifth post-Newtonian order in [13].
For convenience we display the Hamiltonian of EOB in closed form. Its principal structure has been given in [55], Eq. (2.12), The local contribution to Eq. (65) has been derived from relations presented in Ref. [13]. The tail term J tail at 4PN is the same as in Eq. (45).

Observables at the fourth post-Newtonian level
In the following we calculate a series of observables up to 4PN starting from the Lagrangian in harmonic coordinates. 8 We present the energy of the innermost stable orbit E(Ω) as a function of the orbital frequency Ω and the angular momentum J(Ω) along circular orbits, since we would like to use the tail term in closed analytic form directly from harmonic coordinates. Otherwise one would have to perform expansions in the excentricity [30,55,56]. 9 We agree with the previous results given in [11], derived by using different calculation methods. The following kinematic decomposition holds For circular orbits, p.n = 0, the energy E is obtained by The relation between the angular momentum J and r is found by cf. [18]. One may express E = E(J, ν) using the variable x given in Eq. (19) and obtains The innermost stable circular orbit is defined by the requirement dE/dx = 0. The test particle solution reads [18] It is convenient to normalize the angular momentum J by j = J/(G N m 1 m 2 ), for which we obtain +O(x 5 ) .
(71) 8 By removing the acceleration contributions already deviations to the original harmonic coordinates have been implied, cf. [27]. 9 For related results see also [57]. The test particle solution reads [18] 1 Both relations (69,71) are well-known from the results of other calculations, cf. e.g. [18]. The position of the innermost stable circular orbit varies with ν [58]. In the limit ν → 0 it assumes x = 1/6 and in the equal mass case at 4PN it reaches x ≈ 0.245. We have obtained agreement with the above relations by explicit calculation, which are illustrated in Figures 1-2. Here r s = 2G N M/c 2 denotes the Schwarzschild radius. It is evident that even in the region of ISCO at 4PN convergence has not been reached and even higher order corrections are important. In these Figures we limited ourselves to the post-Newtonian approach. When extending the x-range, further matching to the results of other representations becomes necessary, cf. e.g. [59].
The final goal of the post-Newtonian calculations is the prediction of the oscillation signal due to gravitational waves in the detectors. This requires in addition to the Hamiltonian also the luminosity function L, cf. e.g. [60], which currently only allows to use the 3PN approximation on the Hamiltonian side, which we agree, cf. e.g. [60,61].

Conclusions
We have calculated the Hamiltonian of two-body systems in an effective field theory approach to gravity to the fourth post-Newtonian order in harmonic coordinates using dimensional regularization. Here, starting with the third post-Newtonian order divergences in the dimensional regularization parameter ε occur. However, these terms do not affect any observable. These are gauge-and coordinate system artifacts related to the regularization process. No renormalization of the observables related to the process under consideration is needed. We have also presented polefree Hamiltonians, obtained after a canonical transformation. Already before similar calculations in harmonic coordinates were performed in Refs. [8,12,62]. The present calculation, however, cannot be compared literally to these, because of probable differences in Feynman rules and differences in total time derivatives, which are partly large contributions. We have checked for [8] that we obtain the same equation of motion for the mth order Lagrange density. Other differences lay in the elimination of the acceleration terms and their time derivatives by forming the first order effective Lagrangian. The approach in [12] (and references therein) in harmonic coordinates is not an effective field theory approach, but is based on the Fokker action [63]. However, we have shown by an explicit calculation, that the Lagrangian given [12] yield the same result for E(Ω) which we obtain.
We compared to other results obtained working in harmonic coordinates [8,12], ADM coordinates [11] and EOB coordinates [13] and demonstrated equivalence on the Hamiltonian level by the construction of explicit canonical transformations. We have illustrated our results in the case of circular motion for E(Ω) and j(Ω). We also present the Hamiltonian resulting from harmonic coordinates up to O(ε) after eliminating the acceleration by a shift in D dimensions, which is useful for further higher order calculations.