Gravitational Spin-Orbit Hamiltonian at NNNLO in the post-Newtonian framework

We present the result of the spin-orbit interaction Hamiltonian for binary systems of rotating compact objects with generic spins, up to NNNLO corrections within the post-Newtonian expansion. The calculation is performed by employing the effective field theory diagrammatic approach, and it involves Feynman integrals up to three loops, evaluated within the dimensional regularization scheme. We apply canonical transformations to eliminate the non-physical divergences and spurious logarithmic behaviours of the Hamiltonian, and use the latter to derive the gauge-invariant binding energy and the scattering angle, in special kinematic regimes.

In particular, the PN analysis follows the system during the so-called inspiral phase, where the components of the binary move at non-relativistic velocities and their orbital separation is slowly decaying.At this stage, the non-relativistic regime of the motion allows for a perturbative treatment of the problem, that can be studied by a series expansion in powers of v/c, where v is the orbital velocity of the compact binary and c is the speed of light.For bound state systems like binaries of compact objects, in this phase, the kinetic and potential energy must share the same order of magnitude according to the virial theorem, namely v 2 = G N M/r, where G N is the Newton's constant, M is the typical mass of the system and r is the separation between the two bodies.From this perspective, the PN approach is an expansion in two parameters v and G N , whose n th order is referred to as nPN (or usually referred to as N n LO), has contributions from terms proportional to G l N (v 2 ) n+1−l where l = 1, 2, ..., n + 1.The analysis then is separated in two sectors, the conservative sector, where the outgoing radiation is ignored and orbits don't decay, and the radiation sector, where the emitted radiation is analysed.In the radiative sector, the emitted radiation carries energy and angular momentum to infinity.The radiations could also scatter off the background curvature to interact with the orbital dynamics giving rise to the so-called tail effects, which contributes to both the conservative and the radiative sector.See [113][114][115] for reviews.
The EFT approach is very well suited for this problem because it can be well separated into three different length scales, namely, the length scale associated with the compact object R s (Schwarzschild radius), the radius of the orbit r, and the wavelength of the emitted gravitational wave λ.Then an EFT at each scale exactly does the job of determining the essential degrees of freedom.In this method, the determination of any observable quantity at a given PN is equivalent to the computation of corresponding scattering amplitudes, which can be systematically computed through the evaluation of corresponding Feynman diagrams.The leading contribution to scattering amplitudes is given by tree diagrams, whereas higher order corrections are determined by diagrams with multiple loops.Tree diagrams represent rational functions of the kinematic variables (energy, momenta, and masses of the particles), therefore they are easy to compute.Loop diagrams, instead, represent very challenging integrals.This means that computing a scattering amplitude for each new level of precision requires drawing exponentially more Feynman diagrams and solving a vastly more complicated mathematical formula.
In this article we focus on the EFT techniques developed in [40,88] to carry the post-Newtonian analysis in the conservative sector of the two body problem.
The analysis of compact objects without spin up to the 1PN correction in the conservative orbital dynamics was reported in [40] for the first time using EFT techniques.Following this, the 2PN was computed in [41,42], the 3PN was computed in [43] and the 4PN was computed in [44][45][46][47][48].The current state of the art for the post-Newtonian calculation for the conservative potential is the 5PN correction computed in [50,51].Several partial results of 6PN [52,53] have been recently added to the catalogue.
Compact objects are expected to be rotating rapidly and hence the analysis of the spin becomes important for systems with such objects.The analysis of spin in the formalism of EFTs was initiated in [88].For systems with spinning compact objects, the post-Newtonian analysis is usually separated into different sectors with different powers of the spin variable appearing in the effective Lagrangian.The spin-orbit sector is analogous to the fine structure correction obtained to the hydrogen atom, which describe the interaction between the orbital angular momentum of the binary and the spin of one of its constituents.Whereas, the quadratic in spin are analogous to the hyperfine structure corrections to the hydrogen atom and describe the interaction between each other/self of the spins of the binary constituents.
For the spin-orbit sector, the leading order (LO) effective potential was first computed in [88].The next-to-leading order (NLO) potential in [89,90], and N 2 LO in [116].Partial results for the N 3 LO correction were recently presented in [91].Similarly, for the quadratic in spin sector, the LO effective potential was also computed in [88].The NLO in [92][93][94][95], N 2 LO in [96][97][98] and N 3 LO was computed in [99,100].The computation for cubic and higher orders in the spin variables can be found in [101][102][103], whereas the finite size effects are described in detail in [95,102,103].The spin-orbit N 3 LO was then completed in [117] using the techniques of first-order self-force.Here, for maximally rotating compact objects, the N n LO spin-orbit sector contributes at (3/2 + n)PN order and N n LO spin-squared sector contributes at (2 + n)PN order.
In this article, we present the complete conservative N 3 LO spin-orbit interaction, using the framework of EFT [40,88] and an extension of the diagrammatic approach presented in [46].In particular, we begin with deriving the required Feynman rules and then the Feynman diagrams: 1 (0PN) + 4 (1PN) + 21 (2PN) + 130 (3PN) = 156, for the non-spinning sector; 2 (LO) + 13 (NLO) + 100 (N 2 LO ) + 894 (N 3 LO ) = 1009, for the spinning sector, The corresponding tensor integrals are decomposed in terms of scalar loop integrals, up to three loops, which are decomposed in terms of a minimal basis of master integrals (MIs), by means of integration-by-parts (IBP) identities [118,119].The contribution of each diagram to the scattering amplitude is then obtained after substituting the analytic expression of the MIs.Finally, the Fourier transform of the amplitude generates the the effective Lagrangian, later converted into a Hamiltonian which is the main result of this article.The computational setup has been fully automated within an in-house Mathematica routine interfaced to QGRAF [120], used for the diagram generation, xTensor [121], for tensor algebra manipulation, LiteRed [122], for the IBP decomposition, elaborating on some of the ideas implemented in EFTofPNG [123].
The effective Lagrangian obtained in this way contain higher order derivatives of the position and the spin, which are eliminated by employing suitable coordinate transformations.The EFT Hamiltonian is then obtained by applying the Legendre transform, and it is found to contains poles in ǫ = d − 3 (d being the number of the continuous space dimensions), and logarithmic terms depending on the size of the binary system, which are eliminated by suitable canonical transformations.
Our novel result for the spin-orbit interaction Hamiltonian are used to derive the binding energy for circular orbits and the scattering angle, and we are glad to report that both expressions agree with the results given in [117].
The paper is organised as follows.In section 2, we review the description of the spinning binaries within the EFT formalism.In section 3, we present the computation for the N 3 LO spin-orbit potential employing the Feynman diagrammatic approach within the EFT framework.Then, in section 4, we describe the procedure of removing the residual divergences and logarithms to derive the EFT Hamil-tonian.We provide our main result of the spin-orbit Hamiltonian up to N 3 LO in section 5.In section 6, we compute two observable from the EFT Hamiltonian, namely, the binding energy of the binary system in circular orbits with aligned spin configuration and scattering angle for two spinning compact objects with aligned spins.We summarize our main results in section 7. The article contain three appendices: in appendix A, we describe the notations used in this article; in appendix B, we provide the required master integrals and their expression; in appendix C, we provide the Hamiltonians till 3PN in the non-spinning sector and till NNLO in the spin-orbit sector.
We provide the required EFT Feynman rules in the ancillary file Feynman Rules.m and the analytic results of the Hamiltonian till N 3 LO in the ancillary file Hamiltonian.m.

EFT of spinning objects
In this section, we describe the action for the degrees of freedom of the gravitational field and the degrees of freedom of the compact objects, namely their center of mass and their spin.Then we describe the techniques of the Post-Newtonian formulation of GR in detail and briefly outline the procedure to compute the effective action.

Action
The dynamics of the gravitational field (g µν ) is given by the Einstein-Hilbert action along with a gauge fixing term, where Γ µ = Γ µ ρσ g ρσ (in the harmonic gauge Γ µ = 0), Γ µ ρσ is the Christoffel symbol, G N is the Newton's constant, R is the Ricci scalar, and g is the determinant of the g µν .
To model the dynamics of the compact object, we utilize worldlines x µ (a) (τ ) parametrized by an affine parameter τ and define tetrads Λ µ (a)A (τ ) along the worldlines which connects the body-fixed frame (denoted by upper case Latin indices) of the a th compact object and the general coordinate frame (denoted by Greek indices).Then the angular velocity tensor of the spinning object can be defined as and the corresponding conjugate momenta to the Λ µ (a)A is given by . (2.3) Then by demanding the reparameterization invariance of the point particle action, the dynamics for the compact objects is governed by the worldline point particle action given by [102], which, corresponds to Pryce, Newton, and Wigner gauge for spin supplementarity condition (SSC) given by S (a)µν (u ν (a) + u 2 (a) δ ν0 ) ≈ 0. In the above action, the center of mass of the compact object is modeled by the position of the point particle x µ (a) and the spin of the compact object is modeled by S (a)µν .The L (a)SI is given by equation (4.16) of [102] which corresponds to the spin-induced non minimal couplings that does not contribute up to the N 3 LO spin-obit sector at 4.5PN.Here u µ (a) = dx µ (a) /dτ is the four-velocity and u 2 (a) = g µν u µ (a) u ν (a) , whereas the proper time τ is related to the coordinate time t by dτ = c dt.
As we will be performing the computation using the techniques of multi-loop Feynman diagrams, it is necessary to write the gravitational coupling constant in d dimensions as where, R 0 is an arbitrary length scale.

Post-Newtonian formulation of General Relativity
In the bound state of two compact objects, we have three length scales, namely the length scale associated with the compact object R s (Schwarzschild radius), the radius of the orbit r, and the wavelength of the emitted gravitational wave λ.We assume the velocities of the particles to be small as compared to the velocity of light and the particles are far from each other, hence propagate on a flat background (g µν = η µν + h µν ), where the gravitational interaction between the two particles is governed by the gravitons h µν .Then we have a hierarchy of length scales As we are only interested in the long-distance physics at the scales of λ, we first decompose the graviton fields in short distance modes -potential gravitons H µν with scaling (k 0 , k) ∼ (v/r, 1/r) and long-distance modes -radiation gravitons hµν with scaling (k 0 , k) ∼ (v/r, v/r) [40].
Noting that v 2 ∼ 1/r for bound orbits due to the virial theorem (or the third Kepler law), the dimensionless expansion parameter can be taken as v 2 /c 2 ∼ G N M/rc 2 , which formally scales as 1/c 2 .Hence, following the majority of the PN literature, we equivalently adopt a formal expansion in 1/c with one PN order corresponding to 1/c 2 .For the spin variables, it holds S (a) = Gm 2 (a) χ (a) /c where the dimensionless spins χ (a) are at most O(1) for black holes and (realistic) neutron stars, so that S (a) ∼ 1/c.Henceforth we rescale the spins as S (a) → S (a) /c in order to make the PN counting in 1/c manifest.Now to compute the conservative binding potential of the two-body system, we ignore the radiation modes and decompose the potential modes in the Kaluza-Klein (KK) parameterization [124,125].In this, the different components of metric g µν (= η µν + H µν ) are encoded in three fields, a scalar φ, a 3-dimensional vector A i and a 3-dimensional symmetric rank two tensor σ ij .The decomposition is given by, where, The 2-body effective action is then given by integrating out the gravitational degrees of freedom from the above derived actions as where, L eff is the effective Lagrangian further decomposed as where, K eff is the kinetic term and V eff is the effective contribution due to gravitational interactions between the two objects.
The effective potential can be represented in the form of connected, classical, 1 particle irreducible (1PI) scattering amplitudes as where p is the momentum transfer between the two particles and the box diagram in the above equation refers to all possible Feynman diagrams with gravitons (φ, A i , and σ ij ) mediating the gravitational interaction between the two point particle represented by the two solid black lines.Our aim in this article is to compute the spin-orbit effective potential up to N 3 LO .For this, we further decompose the kinetic and potential terms as K eff = K pp + K spin and V eff = V pp + V spin where K pp and V pp represent the kinetic and potential terms for center of mass degrees of freedom for the point particle and K spin and V spin represent the kinetic and potential terms for the spin degrees of freedom.The expression for the kinetic terms are given by and the decomposition of the potential terms is defined as follows where, V N stands for the Newtonian potential and V j with j = {1PN, 2PN, 3PN} refers to the corresponding PN correction for the non-spinning part of the potential.The V SO j with j = {LO, NLO, N 2 LO, N 3 LO} refers to the corresponding correction to the spin-orbit coupling of the binary system.Then our aim in this article would be to compute the V j and V SO j using the techniques of multi-loop scattering amplitude, as further described in the next sections.

Computational Algorithm
To obtain the effective potential from the diagrammatic approach as shown in equation (2.10), we begin by generating all the relevant generic topologies contributing at different orders of G N .It could be easily seen from the virial theorem that N n LO has contributions from terms proportional to G l N where l = 1, 2, ..., n + 1, and we consider all the topologies at l − 1 loops contributing to the specific order l.So, for the computation of the N the order G 4 N (3-loop) using QGRAF [120].There is 1 topology at order G N (tree-level), 2 topologies at order G2 N (one-loop), 9 topologies at G 3 N (two-loop), and 32 topologies at order G 4 N (three-loop).Then we dress these topologies with the KK field and Feynman rules derived from the action of PN1 expansion of GR given in (2.1) and (2.4), to obtain all the Feynman diagrams that contribute to the given order of G N and v depending on the specific perturbation order.The number of diagrams that contribute at particular order in 1/c and of particular loop topology are given in table 1a and 1b 2 .

Multi-loop Diagrams
Within the EFT framework, the sources remain static and as a result the generated Feynman diagrams are mapped to two-point multi-loop Feynman diagrams with mass-less internal lines and an external momentum (the momentum transferred between two sources) as shown in figure 1.We translate these Feynman diagrams to their corresponding Feynman amplitudes after performing the tensor algebra using xTensor [121].The generic form of the effective potential corresponding to any l-loop Feynman graph G can be expressed as ) where, (i) N C is a tensor polynomial depending on the world-line coordinates (x µ (a) ), the spin tensor (S (a)µν ), and their higher-order time derivatives, (ii) p is the momentum transfer between the sources (Fourier momentum), (iii) N F is the tensor polynomial built out of momenta p, (iv) k i are the loop momentums, (v) D σ denotes the set of denominators corresponding to the internal lines of G, (vi) N M stands for a tensor polynomial built out of external momenta p and loop momenta k i .
We perform the reduction of the multi-loop tensor integrals to scalar integrals by applying a set of projectors exploiting the Lorentz invariance.We build the projectors depending on the Lorentz invariant external momentum(p) and the background metric.After applying the projector the numerator (N M ) of the multi-loop integral is transformed as where, ÑM is a polynomial depending on the scalar products between the external momentum (p) and the loop momentums (k i ).We use projectors up to rank 6 for the complete evaluation of the N 3 LO Spin-Orbit effective potential.The generated scalar integrals are not all independent and there exist linear relations between these integrals, which are envisaged via the Integration-By-Parts (IBP) relations.We employ LiteRED [122] to generate these IBP relations to obtain a smaller set of independent scalar integrals, known as Master integrals (MIs).For the specific computation of N 3 LO Spin-Orbit effective potential, we obtain 1 MI at one-loop, 2 MIs at two-loop, and 3 MIs at three-loop.These MIs are well known and admit closed analytic expressions in d dimension.We provide the explicit expressions in B.1 for convenience.
After the evaluation of the multi-loop integral, we perform the Fourier transform of the tensor polynomials, which have the form where, N F and ÑF are the tensor polynomial built out of Fourier momenta p and f (p) is a function of the external momentum arising from the multi-loop integral.The Fourier integrals up to rank 8 are required, which are obtained by iterative differentiation of the expression given in section B.2.Then, we obtain the effective potential that only depends on the orbital variables x (a) , S (a) , and their higher-order time derivatives by expanding the complete scattering amplitude as Laurent series in ǫ around d = 3.

Integrand generation
Multi-loop methods The computation of the effective potential starting from the generation of the required Feynman diagrams, expressing them in multi-loop integrands, performing IBP reduction, and then applying the Fourier transformations have been automated through an in-house code, elaborating on some of the ideas implemented in EFTofPNG [123], and using xTensor [121] for tensor algebra manipulations as well as successful interface to LiteRed [122], Reduze [126], KIRA [127] for the IBP reduction.A flow chart for the complete computational algorithm for the effective potential as implemented in our in-house code is shown in Fig. 2.

Processing the effective Lagrangian
The effective potential obtained in this way usually contains higher-order time derivatives of the position (a (a) , ȧ(a) , ä(a) ,• • • ) and the spin ( Ṡ(a) , S(a) ,• • • ).In our computation of the N3 LO spin-orbit potential, we have 6th order time derivative of position and 3rd order of time derivative in spin.We need to eliminate these higher-order time derivatives from the potential to facilitate the computation of the Hamiltonian, to obtain gauge-invariant quantities as well as to pave the way for the implementation within the effective one-body formalism.Additionally, the potential in the non-spinning sector at 3PN and the N 3 LO spin-orbit potential contains poles in the dimensional regularization parameter ǫ and logarithmic terms of the form log( r R0 ).These terms may originate either from the choice of coordinates or because of tail effects, and they can be removed, after combining the potential and the tail contribution to conservative effect, by a suitable coordinate transformation [47].However, the tail effects do not arise either upto 3PN order in the non-spinning case, or up to N 3 LO , alias 4.5PN, in spin-orbit sector [128], therefore all the divergent terms and logarithms in our computation are removed by finding a proper coordinate transformation.First, we discuss the procedure to eliminate the higher-order time derivatives and then we show how we remove the divergent pieces and the logarithmic terms to obtain the EFT Hamiltonian free of any divergences and logarithmic terms.

Elimination of higher-order time derivatives
We briefly describe the procedure to remove all the higher-order time derivatives, using coordinate transformations [97,[129][130][131][132].The Lagrangian under a coordinate transformation x (a) → x (a) + δx (a) changes by So, when the equation of motion (EOM) is linear in a (a) at LO, we can construct a perturbatively small δx (a) such that terms depending on a (a) drop out in L + δL.Similarly, for terms involving higher-order time derivatives of a (a) , one can take δx (a) to be a total time derivative such that the higher-order time derivatives in L + δL cancel upon partial integration.In general, the O δx 2 (a) contributions have to be kept, but will turn out to be negligible for the explicit steps outlined below, making the procedure equivalent to insertion of EOM [130].Following the same approach, when we apply a small transformation simultaneously to the rotation matrix . The variation of the rotation matrix, including quadratic terms in ω (a) , reads as, where the ω kj (a) is the antisymmetric generator of the rotation matrix, 3 and similarly the spin transforms as Due to such transformation, the Lagrangian changes by So, when the above equation is linear in Ṡ(a) at LO, we can construct the ω ij (a) such that all the terms depending on Ṡ(a) and higher-order derivatives drops out in L + δL.Here, it is important to keep terms quadratic in ω (a) , as first pointed out in [133], merely inserting EOM would be missing these contributions.We apply this procedure iteratively to eliminate the higher-order time derivatives from the N 3 LO spin-orbit potential.Specifically, we perform 4 iterations, where 1. we remove the terms with a (a) and its higher-order time derivatives from the LO and NLO spin-orbit potentials, 2. we remove Ṡ(a) and its higher-order time derivatives from the NLO spin-orbit potential, 3. we remove the a (a) and its higher-order time derivatives from non-spinning 2PN and 3PN potentials, 4. we remove the a (a) , Ṡ(a) and their higher-order time derivatives from the N 2 LO and N 3 LO spin-orbit potentials.
After each iteration, we obtain a new Lagrangian (with new EOM) to be used in the next iteration.
One can check that at each step contributions quadratic in δx (a) and cubic in ω (a) are negligible (higher order in spin or the PN approximation).Following these steps, we obtain the effective Lagrangian, which depends on the position, velocity, and spin only.

Computation of the EFT Hamiltonians
We apply the Legendre transformation on the effective Lagrangian obtained in the previous step to derive the EFT Hamiltonian where the canonical momenta p i is defined as .
We invert this relation (order by order in 1/c) to express ẋi (a) in terms of p i (a) .We use this relation between ẋi (a) and p i (a) in the equation (4.5) to obtain the required Hamiltonian H(x, p, S).

Removal of the poles and logarithms
The Hamiltonian derived in the previous step, and the Lagrangian both contain divergent pieces and logarithmic terms, which can be removed by finding suitable coordinate transformations.While using the Lagrangian description, we can add total time derivative terms with arbitrary coefficients, and fix these coefficients to a set of values such that the divergent piece drops out while removing the higher-order time derivatives as described in section 4.1.
On the other hand, in the Hamiltonian description, we can define a canonical transformation given by 5 where G is the infinitesimal generator of the canonical transformation.Usually, the ansatz for the terms in the generic total time derivative to be added to the Lagrangian contains a large number of terms. 5The Poisson bracket is defined as where the last term is ignored in the context of the current analysis since it contributes to higher PN orders.
Here, we remove the poles and the logarithms in the Hamiltonian description by defining a suitable canonical transformation.We put up an ansatz for the infinitesimal generator G with arbitrary coefficients and build a system of linear equations in these coefficients by demanding the cancellation of the divergent pieces in H ′ .We use the solution to fix the arbitrary coefficients to a set of values, thereby obtaining the final Hamiltonian H ′ free of poles and logarithms.In the following sections, we describe the procedures for the removal of the poles and logarithms from the 3PN non-spinning sector and the Spin-orbit N 3 LO sector.

3PN non-spinning sector
The 3PN non-spinning Lagrangian has been studied extensively [43] and we follow the steps suggested to remove the divergent pieces.Specifically, following [43] we add a total derivative term with the complete 3PN Lagrangian where, and n ≡ r/r is the separation unit vector for the binary.We follow the usual procedure of removing a i and Ṡij and their higher-order time derivatives as described in Sec.4.1 including the total derivative term with the 3PN Lagrangian.It turns out that the divergent pieces drop out and we obtain a finite Lagrangian till 3PN.Then we remove the terms involving logarithms by finding a canonical transformation to the Hamiltonian obtained from the above Lagrangian.For this purpose, we use the following ansatz for the infinitesimal generator, where, We use this generator to build a canonical transformation following Eq.4.8.By requiring the removal of the logarithm terms in the new Hamiltonian, we build a system of equations and solving them we obtain

Spin-orbit NNNLO sector
For the Spin-orbit Hamiltonian at N 3 LO , we remove both the divergent terms as well as terms involving logarithms by finding a suitable canonical transformation.Inspired by [97], we build the following ansatz for the infinitesimal generator where the arbitrary coefficients g n are defined as Following Eqn.4.8, we define a canonical transformation using the above ansatz for the infinitesimal generator and build a system of linear equations by requiring the removal of the divergent as well as the logarithmic terms.We solve these set of equations and find the arbitrary coefficients as 10 , 10 , Following these steps, we ultimately obtain the full effective Hamiltonian free of divergences and logarithms for the non-spinning sector till 3PN and for the spin-orbit sector till N 3 LO .

Results
To begin, we first define a set of dimensionless variables to express the Hamiltonian in a compact form.We introduce the total mass M = m (1) + m (2) , the reduced mass of the two body system µ = m (1) m (2) /M , the mass ratio q = m (1) /m (2) , and the symmetric mass ratio ν = µ/M .These are related as follows (5.1) Additionally, we express the results in the center of mass (COM) frame of reference and define the momentum in COM frame as p ≡ p (1) = −p (2) .In the COM, the orbital angular momentum is defined as L = (r × p).So, we can write p 2 = p 2 r + L 2 /r 2 , where p r = p • n, p ≡ |p| and L ≡ |L|.All the variables are rescaled to write the expressions in the form of dimensionless parameters The total EFT Hamiltonian can be written as where, In the non-spinning part, the Hamiltonian till 3PN is known in the literature [43] and in the spinning sector the Hamiltonian is known till NNLO [116].The novel result of the HN 3 LO in the EFT approach in the COM frame is presented here + (1 ↔ 2) . (5.6) Here, (1 ↔ 2) implies changing the label 1 and 2 on the spin variables along with q ↔ 1/q.The computation of the previously known Hamiltonians, together with the novel computation of the HN 3 LO , allows us to obtain the complete expression for the spin-orbit N 3 LO Hamiltonian within the EFT framework.We provide the known Hamiltonians in appendix C. The Hamiltonian in the generic frame is also provided in the ancillary file Hamiltonian.m with this article.

Computation of observables with spin
The obtained Hamiltonian is still gauge dependent as it depends on the radial coordinate, we compute gauge invariant observable to compare with the literature [117].In this section, we describe the procedure for the computation of two gauge invariant observable: the binding energy for circular orbit with aligned spin configuration and the scattering angle with aligned spin configuration.
For this, we work in the COM defined by p (1) + p (2) = 0 as given in section 5.With this, we assume that the spins are aligned to the direction of the orbital angular momentum of the binary system.Such aligned spin configuration is given by where, L = |L| and S (a) = |S (a) |.

Binding energy for circular orbits with aligned spins
We compute the gauge invariant relation between the binding energy and the orbital frequency for circular orbits by removing the dependence on the radial coordinate.For circular orbits, we know We invert the above relation again to express L as a function of ω.Additionally, we define a gauge invariant PN parameter x = ω 2/3 .Using the above two equations, we express L as a function of x, and substituting this in the Hamiltonian, we obtain the gauge invariant binding energy of the circular orbits in the aligned spin configuration Ẽ(x, S(a) ).Following the above procedure with the Hamiltonian given in section 5 we obtain where, , and S ⋆ = S (1) /q + S (2) q.Our result given in equation (6.7) agrees with the result given in equation ( 10) of [117].

Scattering angle with aligned spins
In this section, we study the scattering angle considering an aligned spin binary system following [105].First, we re-scale the spin variables as a (a) = S (a) /m (a) .In the COM, the Hamiltonian H is expressed as a function of p r , L, r, and S(a) and inverting that we obtain p r = p r (H, L, r, S (a) ) .(6.8) Then the scattering angle χ is given by χ(H, L, S (a) ) = − dr ∂p r (H, L, r, S (a) ) ∂L − π .(6.9) Now, the total center of mass energy (E) and the total energy per total rest mass (Γ) is given by where γ is the Lorentz factor.We invert the above relation to express the Lorentz factor γ in terms of Γ, and we obtain where, v ≡ |ṙ| is the relative velocity of the compact objects.Moreover, the total angular momentum L can be expressed in terms of the impact parameter b and the aligned spin configurations a (+) and a (−) as where, δ = (m (1) − m (2) )/M , a (+) = a (1) + a (2) and a (−) = a (1) − a (2) .We use equation (6.10) and (6.12) to trade H for v and using equation ( 6 Now, applying the above procedure with the Hamiltonian given in section 5, we obtain the scattering angle computed in the COM for aligned spins, which can be expressed as where, and Our result given in equation (6.17) agrees with the result given in equation ( 7) of [117].We note that the N 3 LO spin-orbit scattering angle contains all gauge-invariant information of the corresponding Lagrangian or Hamiltonian [128], albeit restricted to the center-of-mass frame.

Conclusion
In this work, we presented the complete evaluation of N 3 LO Post-Newtonian (PN) correction to the the spin-orbit Hamiltonian for rapidly rotating compact objects, within the effective field theory diagrammatic approach of General Relativity.The required Feynman diagrams in momentum space were automatically generated using EFT Feynman rules for the interaction of spinning compact objects.They were algebraically decomposed in terms of master integrals, by means of integration by parts identities for dimensionally regulated integrals, whose expressions were available in the literature.The contribution of each diagram to the Lagrangian, namely to the effective potential, was found after taking the Fourier transform to position space, and series expanding around d = 3 + ǫ space dimension.The Hamiltonian is finally found by means of Legendre transform, followed by suitable canonical transformations, required for eliminating residual non-physical divergences and spurious logarithmic behaviours.
We computed the gauge invariant relation between binding energy and the orbital frequency in the circular orbit with aligned spins as well as the the scattering angle with aligned spins and they were found in agreement with the results available in the literature, previously obtained using the self-force formalism.
Our results complete the description of coalescing rapidly rotating binary systems up to 4.5 PN order within the spin-orbit sector.
The techniques and the in-house code developed for the current project are very flexible, and can be applied to other interesting problems, at higher order in the PN expansion, for both non-spinning and spinning compact objects.
Note added During the completion of this work, a similar study has appeared in [133]

B.2 Fourier Integrals
The expression for the scalar Fourier integral is given by which is for generic dimensions but only in Euclidean signature.The required tensor Fourier integrals are then obtained by taking derivatives on both sides with respect to x i .

C Lower-order Hamiltonians
In this appendix we give the results for all the Hamiltonians given in eq.( 5.4) and eq.(5.5).
C.1 Non-spinning sector up to 3PN

d = 3 Figure 2 :
Figure 2: Flowchart of the computational algorithm

3 )
We invert the above relation to express r as a function of L and we define the orbital frequency asω = ∂ H( L, S (a) ) .13) we trade L for b.This allows us to express the scattering angle as χ(v, b, S (a) ) = − γ µγv dr ∂p r (v, b, r, S (a) ) ∂b − π .(6.14)

Table 1 :
3LO spin-orbit potential, we generate all the topologies till Number of Feynman diagrams contributing different sectors.