The 2PM Hamiltonian for binary Kerr to quartic in spin

From the S-matrix of spinning particles, we extract the 2 PM conservative potential for binary spinning black holes up to quartic order in spin operators. An important ingredient is the exponentiated gravitational Compton amplitude in the classical spin-limit for all graviton helicity sectors. The validity of the resulting Hamiltonian is verified by matching to known lower spin order results, as well as direct computation of the 2PM impulse and spin kicks from the eikonal phase and that from the test black hole scattering based on Mathisson-Papapetrou-Dixon equations.

B.4.2 Triangle integral coefficients 47 1 Introduction The detection of gravitational waves by LIGO and Virgo collaboration [1,2] has brought upon rich opportunity for deeper understanding of gravitational dynamics, galaxy and star evolutions, as well as opening up new alleys in the era of multi-messenger astronomy. The potential for great discoveries also calls for a vigorous pursuit of precision theoretical predictions. The dynamics of the in-spiraling phase are conventionally captured by the post-Newtonian (PN) expansion, which describes the dynamics as Newtonian gravity perturbed by special (v 2 ) and general ( GM r ) relativistic corrections [3][4][5][6][7][8][9][10], where the two dimensionless numbers controlling the expansion have similar order of magnitude (v 2 ∼ GM r ) due to the virial theorem. The approach has recently been augmented by effective field theory methods [11][12][13][14][15][16][17][18][19].
An alternative organization of the perturbative expansion is to consider special relativistic dynamics as the unperturbed theory, which is corrected by general relativity (GM/r). This is known as the post-Minkowskian (PM) expansion [20][21][22][23][24][25][26][27][28][29]. The revival of interest in PM dynamics is partly due to their application to effective one-body formalism [30,31], which is one of the major frameworks for constructing gravitational waveform templates used in LIGO [32]. Another reason the PM expansion has gathered interest is because this is a natural context for applying scattering amplitude methods to classical gravitational dynamics .
One major hurdle for analytic high-precision calculations in general relativity is the proliferation of high rank tensors; there are simply too many cranks and handles to coordinate in a computation. On the other hand, recent progress in scattering amplitudes has shown that conventional approaches contain redundant calculations which can be avoided by a judicious choice of computation methods. It was soon realised that the insights of modern amplitude techniques can also be applied to classical gravitational dynamics [70], and the idea was subsequently developed [36,37,39] and led to the breakthrough computation of 3 PM gravitational scattering [42,43]. The current state-of-the-art result is classical gravitational scattering at 4 PM order [60] based on various tools of modern amplitude computations. Our aim is to add another item to this arsenal of tools; efficient reorganization of spin degrees of freedom.
In the EFT approach, spins are introduced as extra degrees of freedom on the worldline [71], accompanied by an infinite number of spin-induced multipole operators whose strengths are parametrised by the corresponding Wilson coefficients [72,73]. These coefficients can be fixed to all orders in spin at 1 PM by matching to the 1 PM potential [74], or solving for the effective stress-energy tensor through the linearized field equations [75]. The on-shell S-matrix viewpoint alluded to the possibility of a simple organisation scheme for these seemingly complex infinitely many operators, through the observation that these infinitely many operators for a Kerr black hole corresponds to the minimally coupled gravitational three-point amplitude [40,41,76,77] which only consists of a single term. The relation between black holes and minimal coupling was confirmed through reproduction of the full 1 PM potential, impulse and spin-deflections [45,50,67,78].
The kinematics we are interested in has three length scales separated by large orders of magnitude; the de Broglie wavelength λ dB ∼ | p| −1 and the Compton wavelength λ C ∼ m −1 controlling the quantum properties of the bodies, 1 the spin length a ∼ | S|m −1 controlling the "size" of the bodies, 2 and the impact parameter | b| controlling the separation between the bodies. (1.1) Collectively denoting the wavelengths as λ, this separation of scales gives us two expansion parameters; λ/a (λ/| b|) and a/| b|. It is more convenient to recast these small numbers in the form where L = b × p is the orbital angular momentum, s is the spin quantum number, and l is the orbital angular momentum quantum number, with powers of restored to match the dimensions. The scaling of small numbers eq.(1.2) shows that the classical limit corresponds to leading order kinematics for the expansion parameter(s) λ/a ∼ 1/s (λ/| b| ∼ 1/l), while the spin expansion a/| b| is an expansion in the classical regime. An additional source of is provided by virtual gravitons inside loops, which is tracked by restoring to the coupling constant G → G . The inverse power of is cancelled by transferred momentum q µ , which should be counted as wavenumbersq µ ≡ q µ / in the classical limit [39]. Thus the classical limit corresponds to the following dependence for the variables: Of the one-loop master integrals relevant to O(G 2 ) order, only the triangle integral (∝ q −1 ) has the right scaling to have direct classical contributions: The bubble integral (∝ log q 2 ) has too many powers of q and corresponds to quantum corrections, while the box integral (∝ q −2 log q 2 ) has too less powers of q and the behaviour was termed superclassical in ref. [43]. For classical dynamics the bubble integral is irrelevant for obvious reasons, while the box integrals partially cancel between themselves and the remaining IR divergent part is given an interpretation as the iteration term from the lower perturbation order [37,43]. For the spin expansion we perform the subsitution | b| → | q| −1 , since the impact parameter b is interpreted as the dual Fourier variable of transferred momentum q. Therefore the classical spin expansion is an expansion by the small number which also could be understood as the ratio of spin to orbital angular momentum | S|/| L|. We perform the expansion in kinematic expansion parameter eq.(1.4) on top of the classical G expansion provided by the scaling eq.(1.3).
In computing the 2 PM potential from scattering amplitudes, one of the important ingredients for the on-shell method is the gravitational compton amplitude, which contributes to the one-loop unitarity cut. In ref. [41] the massive spinor helicity form of the amplitude was utilized accompanied by the holomorphic classical limit (HCL) [34], which misses certain spin coupling terms that contribute to the classical potential. This is not a huge limitation as missing terms become irrelevant in the phenomenologically interesting case of aligned-spin configuration [40], but it would be preferable to have a complete answer. Recent results at linear and quadratic order in spin [61,79] uses Feynman rules in the effective field theory framework [79]. For higher order in spin, comparable (one-loop) order results for effective action are available in the context of PN EFT [80,81] (for one-loop order results for lower orders in spin, see refs. [71,[82][83][84][85][86]).
To simplify the complexity of the computation beyond quadratic in spin order, we instead utilize the exponentiated form of the Compton amplitude. Indeed it was first observed that the three-point amplitude for gravitational minimally coupled spinning particle exponentiates in the classical spin limit and reproduces the spin multipoles of the Kerr solution [40,41,76]. Similarly, using soft factors such exponentiated form was presented for the opposite helicity graviton compton amplitude in [40,50,87]. Here we instead utilize BCFW recursion [88] on the three-point exponentiated form to derive the gravitational Compton amplitude for all helicity sectors, which also matches with previously known results. Note that the non-uniqueness of the Compton amplitude for minimal couplings with spins greater than 2 is reflected in the nonlocality of the exponentiated amplitude when expanded beyond quartic order in spins. This result reaches the limit constrained by our knowledge of the gravitational Compton that matches with Kerr. Equipped with the gravitational compton amplitude, we proceed to compute the oneloop amplitude of gravitationally coupled spinning particles and extract its classical contributions. The one-loop amplitude is presented in scalar integral basis with its coefficients computed via unitarity methods. Classical contributions only arrises from scalar box and triangles, with the former carrying infrared divergences. We expand the integral coefficients to quartic order in total spin operators and compute the effective Hamiltonian based on the effective field theory (EFT) matching technique [37] extended to spinning objects [79]. In the EFT approach both 1 PM and 2 PM potentials contribute to the O(G 2 ) amplitude, where the former appears as iteration terms containing IR divergences. The matching of IR divergence serves as an initial consistency test, and subtracting off the iteration terms gives us the 2 PM potential up to quartic order in spin. For the 1 PM potential, we use the all order in spin 1 PM potential presented in ref. [41] in a compact form.
Besides matching to earlier results at lower order in spin, we compute gauge invariant physical observables such as the impulse ∆ p and spin-kick ∆ S to test the validity of the Hamiltonian. Based on earlier works at 1 PM [39,45], Bern et. al. [79] proposed a general relation between these physical observables and the derivatives of the eikonal phase with respect to impact parameters and spin-vectors. These relations were verified up to quadratic order in spin at O(G 2 ) in ref. [61] (analogous relations in slightly different contexts appear in refs. [89,90]), which we extend to the quartic order by matching to our Hamiltonian calculations. Furthermore, we compute these observables in the probe limit using the Mathisson-Papapetrou-Dixon (MPD) equations, which describe motion of a spinning body on a curved background. We compute observables to quartic order in background spin and linear order in test body spin with orthogonal spin orientations to verify our results beyond the well-studied aligned-spin scattering.

The classical spinning Compton
In this section, we review the exponentiated gravitational Compton amplitude valid to quartic spin order O(S 4 ), which is equivalent to the known spinor-helicity expressions for minimally coupled particles with quantized spin s ≤ 2 [91,92]. While similar expressions have already appeared in the literature [40,50,77,87], our approach differs in that i) we use BCFW recursion in place of soft theorems, ii) we seek an expression independent of spin representations, bypassing the need to introduce generalized expectation value (GEV) used in refs. [40,50], and iii) we present the same helicity Compton amplitude in exponentiated form, which to the best of authors' knowledge has only appeared in the context of heavy particle effective theory [77]. 3 We also present extension of the exponentiated form to non-minimal couplings in appendix B.

The exponentiated representation of matter-graviton coupling
The 3-pt amplitude that describes a spin-s particle minimally coupled to helicity h boson is given as 4 M +|h| 3,min = g 0 x |h| 21 2s , M −|h| 3,min = g 0 x −|h| [21] 2s . (2.1) Here minimal coupling is defined kinematically by matching to the helicity conserving amplitude in the high energy limit [91]. Note that plus helicity couplings are simplest in the chiral basis, while negative helicity are simplest in the anti-chiral basis. The two expressions can be combined into a single expression 5 [50,87] where M 0 3,min is the minimal coupling of spinless matter (g 0 x |h| or g 0 x −|h| ), which we refer to as the scalar factor, and adopt the convention that massive leg 2 is interpreted as the final state particle with all momenta in the incoming convention. The subscript 3 The authors would like to thank Kays Haddad for bringing this paper to attention. 4 We will take the coupling constants to be real. Complex couplings correspond to dyons for electromagnetic couplings and Taub-NUT probes for gravitational coupling [93][94][95]. Spinor conventions follow appendix A of ref. [41].
is another valid expression, where ε * 3± = ε3∓ is the "complex conjugate" of the polarisation vector. This expression is singular in the m → 0 limit and 4pt amplitudes constructed from this expression are plagued by non-physical t-channel poles, so we do not consider this expression.
Rep(s) denotes the spin-s Lorentz representation for the polarisation spinors/tensors ε 1 and ε 2 , and the associated representation for the Lorentz generator J µν . We emphasise that the representation Rep(s) could be any spin-s Lorentz representation, 6 i.e. any of (s, 0), ( 2s−1 2 , 1 2 ), · · · , (0, s) can be used as Rep(s). The importance of representation independence in classical scattering of spinning particles has been emphasised in ref. [79], and this feature follows from the relations eq.(A.4) and eq.(A.5) reviewed in appendix A.
For finite spin representations, expansion of the exponential as a Taylor series terminates at finite orders [40,41] because the exponent is a nilpotent operator, reflecting the trace relations of generators for finite dimensional representations. In practice we are considering the classical spin limit, where the classical spin S is held fixed: In this limit the nilpotency of the exponent is pushed to infinite order and all terms in the expansion becomes relevant. We use the equality · = whenever the representation of polarisation spinors/tensors is irrelevant, and omit the unnecessary details. For example, eq.(2.2) is recast as An expression equivalent to eq.(2.2) was used in refs. [40,50] for studying binary black hole scattering, where generalized expectation value (GEV) was introduced to account for the chiral nature of the spin representation used in the computations. The main distinction in our approach is avoiding specific spin representations as in eq.(2.4) based on representation independence of eq.(2.2), which unnecessitates the introduction of GEV 7 and allows us to interpret the results using standard quantum field theory. The difference will become transparent for Compton amplitudes presented in the following sections.
The main advantage of an expression of the form eq.(2.4) is that the Lorentz generators appear in an exponentiated form. The little group of massive spinning particles in 4d is SU (2), where all finite dimensional representations are given as symmetric tensor products of the smallest nontrivial representation. Given a symmetric product spin-s state 8 |ψ s = |ψ 1/2 ⊗2s , the action of symmetry group element exp(−i θ · S) is given as; i.e. we can work with generators of the smallest representation instead of wrestling with generators of the actual representation our spinning matter has. This feature is especially 6 However, the spin representation of the spinning particle will matter when deviations from minimal coupling are considered as in appendix B.
7 For non-minimal coupling considered in appendix B the amplitudes will be dependent on the used spin representation, but GEV is avoided by use of non-chiral spin representations. 8 Spin coherent states form a subset of symmetric product states. Since coherent states are considered as "classical" states, we may restrict to symmetric product states without any loss of generality. Note that the bold notation of massive spinor helicity variables implicitly uses symmetric product states as wavefunctions for the particles. advantageous for large spin representations; working with 2 × 2 matrices is considerably simpler than working with (2s + 1) × (2s + 1) matrices with large s. For example, full evaluation using the Baker-Campbell-Hausdorff (BCH) formula often becomes complicated due to nested commutators, but 2 × 2 matrices are simple enough that matrix logarithms can be computed explicitly.
Moreover, exponentiation automatically symmetrises generator products when expanded into a polynomial of generators. This feature simplifies identification of spin products, which we adopt the conventions of refs. [61,79] where S µν (p 1 ) is the spin tensor for momentum p µ 1 in covariant spin supplementary condition (SSC) S µν (p)p ν = 0. The spin tensor can be rewritten using the spin vector S µ (p) for momentum p µ by where the spin vector is identified with the mass-rescaled Pauli-Lubanski pseudovector After extracting the spin contributions using eq.(2.6), the remaining polarisation spinor/tensor contraction ε 2 (p 2 ) · ε 1 (p 1 ) requires further processing due to the difference of little group frames between the initial and final states. Identification of this contribution is essential when constructing the Hamiltonian, and the procedure was named Hilbert space matching in ref. [78]. The frame difference causes the Thomas-Wigner rotation effect, which results in additional spin factors that were computed to full relativistic order in ref. [96]; see also ref. [79] for its connections to different choices of SSCs. The full form of the factor can be found in eq.(3.6).

BCFW construction of Compton amplitude (opposite helicity)
The 3pt amplitude eq.(2.4) can be used as the seed amplitude to construct 4pt amplitudes. We aim for the 4pt amplitude M 4 (p 1 , k + 2 , k − 3 , p 4 ) using the BCFW shift |2 → |2 +w|3 and |3] → |3]− w|2] (the 23] shift). Definitions for the Mandelstam variables are s = ( The new ingredient in this construction is the treatment of spin factors appearing in eq.(2.4). These factors are operators, and in general their order of composition matters because spin operators are non-abelian. We adopt the convention that leg 4 is the final state particle, so the spin factor associated with the leg 4 is always put to the left of the spin factor associated with the leg 1. On the s-channel cut (w s = − 2|p 1 |2] 3|p 1 |2] ), the spin factors combine into the following simple expression 9 where hatted variables denote BCFW-shifted variables andP = −p 1 −k 2 (w s ). The above result can be checked by explicitly comparing the matrix elements for chiral spin-1 2 representation and anti-chiral spin- 1 2 representation, which we provide the details in appendix A. Since matrix elements of other spin representations follow from their tensor products, checking consistency of the two cases is sufficient for generic spin representations. Interestingly, the u-channel cut (w u = − 2|p 4 |2] 3|p 4 |2] ) results in the same spin factor whereP ′ = −p 4 −k 2 (w u ), and the relevant shifted variables marked by hats are primed to denote that they correspond to the u-channel factorisation conditions. We do not present composition of scalar factors as it is not different from the usual BCFW construction [77,78,92]. Adding up all contributions, the 4pt amplitude M 4 (p 1 , k + 2 , k − 3 , p 4 ) is given as [87] M s , |h| = 1 , , |h| = 2 , with definitions for variables in eq.(2.9), which we present an extension to non-minimal coupling in appendix B.2. The expression eq.(2.11) appears in a slightly different form as eq.(56) in ref. [87] or as eq.(5.14) in ref. [77], and differs from the analogous expresison eq.(2.34) in ref. [40] by the spin factor eq.(2.9). For s ≤ 2, all the expressions are equivalent to the spinor-helicity expression for the same amplitude first presented as eq.(5.36) in ref. [91]. We make the following remarks regarding the expression eq.(2.11): 3. While ref. [87] used the BCH formula to arrive at the results, eq.(2.9) was obtained by solving matrix equations. The advantage of not resorting to the BCH formula will become clear for the same helicity Compton amplitude.
Once again, as we are only interested up to spin-quartic order, the spurious pole appearing at higher orders in J will be irrelevant to us.

BCFW construction of Compton amplitude (same helicity)
Following the previous section, we construct the 4pt amplitude M 4 (p 1 , k + 2 , k + 3 , p 4 ) using the 23] shift, with the same definitions for shifted variables used in the previous section. The spin factors turn out to be quite involved for this set-up, where the s-channel spin factor is given as This infinite sum have a radius of convergence of |x| < 4, reflecting the threshold condition t = 4m 2 . When restoring factors of , in the classical limit one simply takes f (x) = 1. The spin factor composition eq.(2.12) is obtained as the minimal coupling limit for the generic case presented in appendix B.3, and the details of the computation are given in appendix A. Similar to the opposite helicity case, the u-channel spin factor also turns out to be the same; which is a property specific to minimal coupling. Summing up, we get the following expression for M 4 (p 1 , k + 2 , k + 3 , p 4 ) using BCFW construction.
, |h| = 2 . (2.14) In contrast to opposite helicity case eq.(2.11), the expression eq.(2.14) has no spurious poles. This is consistent with spinor-helicity expression for the same amplitude [92], which is equivalent due to the relation which holds when m −1 is treated as a formal expansion parameter. We make following remarks regarding the expression eq.(2.14): 1. The non-local t-channel pole in a i coefficients actually cancel on the t-channel, so all physical channel residues are correctly captured by this expression to all orders in J.
2. Similar to eq.(2.11), the expression eq.(2.14) has a factorisation into the scalar factor and the spin factor. This is not true for generic couplings, as shown in appendix B.3.
3. The spin factor in eq.(2.12) was obtained by solving the matrix equations explicitly without resorting to the BCH formula, as explained in appendix A. The all-order resummation would have been impossibly complicated using the BCH formula.

4.
A similar exponentiated form exists in the context of heavy particle effective theory [77] where representation independence is attained by off-shell continuation of polarisation spinors, in contrast to purely on-shell expression eq.(2.14).
The same helicity Compton amplitudes do not contribute to the triple cut integrand at one-loop for minimal coupling [34], but they do contribute to the quadruple cut integrand. The quadruple cuts are needed for checking IR cancellations.

Integral coefficients and the classical potential
We evaluate the 2PM conservative potential in this section. We start by computing the one-loop amplitude relevant to classical contribution, which is conveniently expanded in terms of scalar integral basis, where we only include the terms relevant to classical physics. The explicit forms of the box and triangle integrals are , (3.2) , .
We set the kinematics in the centre of momentum (COM) frame as: where p ′ = p− q and q is the momentum transfer between the two spinning particles. Expanding in | q| the basis integrals yields [98,99] The infrared divergent terms with ǫ −1 from the box and the crossed box integrals cancel up to a term proportional to | q| 2 |p|(E a +E b ) −1 , and the remaining IR divergent terms should cancel against iteration terms. This cancellation serves as a consistency check to our final result.
In summary, we use the formula: 5) where U a , U b are the Thomas-Wigner rotation factors [79,96] discussed in section 2.1. , U b is obtained by a simple exchange of labels a ↔ b. This factor ensures that the incoming states and the outgoing states are described in the same Hilbert space, where particles are labelled by their momenta and their little group indices in the COM frame; the iteration of the Hamiltonian makes sense only when the domain and the range are the same.V 1PM is the 1PM potential reported in ref. [96]. We first derive the relevant scalar integral coefficients using unitarity methods, and then compute the iteration terms. Using eq.(3.5) we derive the 2PM Hamiltonian up to quartic order in total spin vectors, where the covariant spin vector is defined as [79]: We test consistency of the results by matching physical observables in the next section.  Here are the unitraity cuts used to extract the corresponding integral coefficients. The dashed lines represents cuts from gluing two four-point tree Compton amplitudes (the double cut) and the dotted line is the additional cuts imposed to get the desired cuts. The solid and wavy lines represent spinning bodies and gravitons, respectively.

Integral Coefficients
We evaluate these integral coefficients in eq.(3.1) using the unitarity cut methods provided in ref. [100]. The unitarity cuts of loop amplitudes can be obtained by sewing tree amplitudes. For example, triple and quadruple cut are given by sewing of three and four corresponding tree amplitudes, respectively. In order to generate the triple and quadruple cuts, we start by sewing the t-channel double cut C 2 with classical Compton amplitudes eq.(2.11) and eq.(2.14) as our basic input, 10 The sewing imposes two cut conditions on the loop momentum, leaving two free parameters in four dimensions. Imposing another cut condition yields the triple cut and the number of free parameters reduces to one. For a suitable choice of loop momentum parametrisation the triangle coefficient is determined as the part of the triple cut independent of the free parameter, which can be extracted as the constant term of the Laurent expansion at infinity. Imposing a further cut condition yields the quadruple cut, fixing all parameters of the loop momentum. The (crossed) box coefficient is determined by summing over all possible quadruple cut solutions of the loop momentum.

The triangle coefficients
The triple cuts C 3,# can be obtained from the double cut C 2 , eq.(3.8), as We describe the evaluation of the triangle coefficient c △ from the cut C 3,△ . The other coefficient c ▽ and its associated cut C 3,▽ can be easily obtained by the substitutions The triple cut conditions relevant to c △ , are solved by the loop parametrisation with one free parameter z as [79], where the ± label the two solutions that satisfies the triple cut conditions and (3.13) Inserting eq.(3.9) into eq.(3.12) yields C 3,△ as a function of z. The triangle coefficient is determined as [100] c where Inf z C 3,△ (l) is the constant term of the Laurent expansion of C 3,△ at z = ∞. The other coefficient c ▽ is obtained via eq.(3.10), and the triangle contributions to the 1-loop amplitude is given as The classical limit of the triangle coefficients The triangle coefficients are functions of external kinematics. To extract classical contributions, factors of are restored through the map and only the terms with O( 0 ) scaling are kept. We expand the triangle coefficient by the following basis of spin operators that scales as O( 0 ): We write the amplitude as a sum over various terms with differing powers of the spin operator: △+▽ denotes the amplitude with S i a S j b .

The box and crossed box coefficients
Two more cut conditions are imposed on eq.(3.8) to get the quadruple cut, We only present the computation of c from C 4, , since the crossed box coefficients can be obtained by a simple substitution rule The quadruple cut conditions are solved by four solutions for the loop momentum [79], The null vectors η ± need to satisfy the properties Inserting the solutions into the quadruple cut C 4 yields the box coefficient as (3.36)

The classical limit of the box and crossed box coefficients
We adopt the normalization: • Scalar coupling The crossed box coefficient is obtained by simply replacing p 1 → −p 4 and p 4 → −p 1 , which is equivalent to σ → −σ + O(t), E i → −E i . The crossed box coefficient turns out to be the same as the box coefficient in the classical limit: The box the integrals in eq.(3.4) contain a classical piece that is infrared divergent, but combining the box and the crossed box contribution results in a purely imaginary IR divergent piece because they have the same coefficients: This needs to be canceled by the iteration terms from the 1 PM potential, which we compute next.

Iteration terms
We match the full theory 2 PM classical amplitude to that computed from the EFT to extract the classical potential following refs. [37,79] (for earlier works at 1 PM see refs. [96,101]). For computing iteration terms the on-shell potentials are continued off-shell with a hat to denote the continuationV ( p ′ , p), where p ( p ′ ) is the in-coming (out-going) momentum in the COM frame with | p| 2 = | p ′ | 2 . Spin operators S n satisfying the SU(2) commutation relations are included in the potential for spinning cases, where n = a, b is the particle label.
We begin with the potential expanded in the spin-vector basis [79]: where n a , n b denotes the degree of spin for particles a, b respectively, and α labels distinct tensor structures. The scalar function V na,n b denotes a tensorial combination of momentum p, Kronecker delta, and the 3D Levi-Civita tensor. This is a basis of operators with classical scaling O( −3 ). 11 The symmetrised product of spin operators S can be algebraically reduced as Each term in the potential can be represented as a contact term in the EFT Lagrangian, coupled to appropriate spin states. The EFT 2 PM amplitude is obtained as operators acting on the spin coherent states. Thus the bracket Ô denotes sandwiching the operatorÔ between the spin coherent state. This momentum space amplitude can then be matched to the classical piece of the QFT computation at 1-loop eq(3.1). By requiring M 1−loop classical = M EFT we find the 2PM potential V 2 = V 2 is given by: The second term which stems from contributions from the 1 PM potential is called the iteration terms, and needs to be subtracted from the one-loop amplitude to obtain the 2 PM potential.
As an explicit example, the 1PM potential vertex that goes into the iteration integral is:V where the coefficients V na,n b can be obtained by expanding the closed form 1PM potential [96]: and using eq.(3.17) to rewrite E in terms of spin vector S.

The iteration integral term in the 3D integral basis
Evaluation of the full integral is not needed for determining classical contributions, and considering the integrand in the limit | k| 2 ∼ | p| 2 is enough [37]. Expanding the integrand of iteration term to O(| k| 2 − | p| 2 ) 0 order yields In the first equality V na,n b has been expanded at | k| = | p| so that V na,n b becomes a function of | p|, and the following identity was used.
Iter △ corresponds to the first term in the penultimate line of eq.(3.48) having three propagators, and Iter corresponds to the remaining terms with two propagators. Note that the operators O i ( k) are tensors of k contracted to other vectors, therefore the last line of eq.(3.48) is a tensor integral of k. We expand the tensor integrals by scalar intgrals are the tensors constructed from 3D Kronecker deltas and a scalar function depending on | p| 2 , | q| 2 . Ellipsis denotes terms irrelevant for classical physics. The 3D scalar integrals I, J are given as: (3.51) The full theory 1-loop amplitude of section 3.1 should also be recast using three-dimensional variables. The integral coefficients are converted using eq.(3.3) and eq.(3.7). The 4D scalar integrals are reduced to the 3D scalar integrals as (3.52) The 4D triangle integral is proportional to the 3D scalar bubble integral I, and the imaginary part of the 4D box integral is proportional to the 3D scalar triangle integral J.

The 3D bubble coefficient
We present the 3D bubble part of the iteration integrals up to quartic order in spin based on the discussions of section 3.2. (3.53) The variables are defined as where V i,j are the 1PM coefficients defined in eq.(3.46).
• Scalar couplingĨ This result is consistent with the one given in eq.(6.24) of [79] • S aĨ This result is consistent with the one given in eq.(6.24) of [79] • S a S bĨ 1,1 =Ĩ The corresponding coefficients are: 1,0;0,1 + 1,0;0,1 + This result is consistent with the one given in eq.(6.24) of [79] • S 2 The corresponding coefficients are: This result is consistent with the one given in eq.(4.11) of [61] • S 2 a S bĨ Note that this basis is not linearly independent, since different tensor structures may be related by the Schouten identity: where we used the on-shell condition p · q = | q| 2 2 on the last line. When restoring the factors, the last line of eq.(3.61) carries one more power of , so it can be neglected. The coorsponding coefficients we derived are: 2,1;0,0 + A 2,0;0,1 ) + (3.63) The corresponding coefficientsĨ (α) 3,0 are: 3,0;0,0 + A The expressions are too lengthy to be shown here. The explicit results can be found in the ancillary Mathematica file "3D Integral Coefficients.wl".

The 3D triangle coefficient: Cancellation of IR divergence
IR divergences should not appear in the classical potential. This requires cancellation between the 3D triangle iteration integral and the sum of 4D box and crossed box integrals.
Decomposing Iter △ into a sum over different power of spin: whereJ na,n b are the 3D triangle coefficients for various powers of spin, we obtain the following consistency condition:J whereĉ B denotes the 4D box coefficient with rotation factors U a U b incorporated. This is the primary consistency check of our computations.
• S a S bJ The coorsponding coefficientsJ This result is consistent with the one presented in [79].
• Quartic Spin Order: These coefficients also cancel with their corresponding 4D box coefficients. The expressions are too lengthy to be shown here and they can be found in the ancillary Mathematica file "3D Integral Coefficients.wl".

Conservative spin Hamiltonian
The full spinning Hamiltonian is: where H 1PM is the Fourier transform of (3.47). The components of the 2 PM Hamiltonian can be constructed from the components presented in sections 3.1 and 3.2. Only the 4D triangle integrals and the 3D bubble integrals contribute to eq.(3.5) and the momentum space Hamiltonian isH where M △+▽ was given in section 3.1.1 and Iter was given in section 3.2.1.
We report the full momentum dependence in the ancillary file "2PM Hamiltonian.wl" and only present the Hamiltonian to the first few orders in the | p| expansion. The cancellation between 1/| p| n poles in the iteration term and the triangle part of the 1-loop amplitude provides another consistency check of our results.
Putting eq.(3.23b) and eq.(3.62) into eq.(3.75), one would get: The superficial 1/| p| n poles are removed by applying eq.(3.61) so that they pushed to higher order in , leaving only: 12 (3.84) Note that these numbers are not unique due to the Schouten identity eq.(3.61), but physical observables are not affected by the identity.
The corresponding coefficients H α 3,0 are: (3.86) 12 In the ancillary file "2PM Hamiltonian.wl", we introduce a parameter δ to turn on and off the Schouten identity. When δ = 0, it gives theH α 2,1 in eq.(3.83). When δ = 1, it gives theH Note that this basis is also not linearly independent since different contractions can be related by the Schouten identity: where we use the on-shell condition p· q = | q| 2 2 on the second line. Application of eq.(3.75) results in the expressioñ  13 In the ancillary file "2PM Hamiltonian.wl", we introduce a parameter δ to turn on and off the Schouten identity. When δ = 0, it gives theH α 2,2 in eq.(3.88). When δ = 1, it gives theH (3.90a) The corresponding coefficients H (α) 3,1 are: 3,1 = −
In general, the canonical transformation is given by eq.(3.101). However at NLO, ∆H NLO na,n b takes the form: 102) where H N is the Newtonian Hamiltonian: and H EFT LO,(na,n b ) are the 1PM spinning potential given by eq.(3.47) with | p| = 0.
Starting from the scalar sector, one needs the generator: to show the equivalence of our leading order 2PM scalar Hamiltonian with that from ref. [102]. For the spinning sectors, the NLO O(G) generators g NLO na,n b are already fixed by the NLO O(G) Hamiltonian matching in [96], and g NLO 0,0 is fixed by eq.(3.104), leaving no more room for additional generators at NLO. Indeed, combining the results in [96] and eq.(3.104), we match the results from EFT formalism up to NLO quadratic in spin order.

Observables
We compute gauge invariant observables from Hamiltonian equations of motion (EOM) in the COM frame, focussing on the impulse ∆ p and the spin kick ∆ S a/b . Alternative method for computing these observables is to use the eikonal phase [39,45,79], or in the test body limit to use the Mathisson-Papapetrou-Dixon (MPD) equations. Matching with these results provide further consistency checks of our Hamiltonian.

Observales in the COM frame
We compute the impulse and spin kick in the COM frame using two methods; by integrating the Hamiltonian EOM iteratively to O(G 2 ), and by applying differential operators to 2 PM eikonal phase [79]. The results were shown to be consistent for all computed sectors.

Observables from the Hamiltonian EOM
The Hamiltonian EOM of an observable O i iṡ The EOM for the observables of our interest are We write O i as a series expansion in G, i.e.
The scattering variables are integrals of the EOM from the past to future infinity, which often suffers from IR divergences that can be avoided by a suitable change of integration variable; 14 we divideȮ i byż and use dz as integration variable instead The suitable boundary conditions are: We compute the 2PM impulse ∆ p and spin kicks ∆ S a/b for general initial spin directions

Observables from the eikonal phase
The 1PM and 2PM eikonal phases are defined as [79]: where b is the impact parameter. It was shown that the impact parameter and the spin-kick can be obtained as suitable derivatives of eikonal phase at O(G) order [39,45], which was generalized to O(G 2 ) in ref. [79]. The latter is the version we utilise: where χ = χ 1PM + χ 2PM + · · · is the full eikonal phase, and the Poisson bracket and D SL are defined by: (4.9) 14 We acknowledge Chia-Hsien Shen for discussion of the following method.
p ⊥ denotes the momentum components perpendicular to the incoming momentum p (0) = (0, 0, p ∞ ), and corresponds to x and y components in our set-up. The z component of the impulse ∆p z is determined by the energy conservation condition We present our impulse ∆ p and the spin kicks ∆ S a , ∆ S b from eq.(4.8) up to quartic order in spin in the ancillary file "Observables.wl". We have verified that the results match with that from Hamiltonian EOM for all computed sectors, to all orders in momentum | p|.
As an additional check we compared and matched with the results of ref. [40] where the 2 PM aligned-spin scattering angle was computed up to quartic order in spin. The aligned-spin kinematics imposes the condition b · S = p · S = 0 and the scattering angle simplifies to [61]: whereχ is the eikonal phase without the rotation factors U a,b . The rotation factors are dropped because ref. [40] uses the proper impact parameter, which is different from the canonical impact parameter [75,89] corresponding to the canonical (or Pryce-Newton-Wigner) SSC [104][105][106] suited for Hamiltonian dynamics; 15 see ref. [79] for the rotation factor's connection to SSC choice. The difference between covariant and canonical impact parameter will be important when comparing with test black hole scattering results in the next section.

The background-probe limit
In the extreme mass ratio limit the dynamics of the binary system reduces to that of a test spinning black hole moving on a Kerr background. Scattering observables in this set-up can be calculated using the Mathisson-Papapetrou-Dixon (MPD) equations [108][109][110], which describes a test body moving on a curved background. We use the pole-dipole approximation which describes a test particle endowed with a monopole moment (momentum) p µ and a dipole moment (spin tensor) S µν . The equations can be solved iteratively to compute scattering variables to linear order in test-body spin, which we provide the details in appendix C. This method provides an independent computation that can be used as consistency checks of our results.
We consider scattering on a Kerr background with initial data

12)
15 A review on consequences for different choices of SSCs can be found in ref. [107].
where b µ is the impact parameter, u µ 0 is the initial 4-velocity, and s µ 0 is the mass-rescaled initial Pauli-Lubanski spin vector s µ 0 = S µ 0 /m. The black hole spin is directed along the z-axis and covariant SSC is implied. We choose this configuration to test genericity of our results beyond well-studied aligned-spin cases while keeping some simplicity. The full results are presented in appendix C.4.
To make contact with the observables computed in section 4.1.1 and 4.1.2 we take the test body limit, by keeping only the leading contributions in the large mass ratio M m → ∞ expansion, where we identify m a = M as the black hole's mass and m b = m as the spinning test body's mass. The difference between canonical and covariant variables need to be accounted for; Hamiltonian EOM uses canonical variables, while MPD EOM uses covariant variables. This would require us to do a shift to the impact parameter [75] .
For the spin kick, we compare the spin 4-vector since it is free of "frame dependence" as demonstrated in ref. [67]. In practice, the spin kick ∆S b computed from the Hamiltonian EOM is boosted as p , (4.14) where Up to the explained technicalities, the observables computed using three different methods were confirmed to be consistent with each other.
For the simpler case of aligned-spins, the scattering angle for test black hole scattering was computed in ref. [111] to G 2 order with quadrupole moment J µνλσ included in the analysis. The computations amount to quadratic-in-spin order computation for the scattering angle, which was shown to match with amplitude results eq.(1.12) of ref. [40], which was also compared with quartic-in-spin order results for the test black hole in ref. [112]. Therefore matching with these results would be redundant, since the scattering angle of aligned spin case presented in ref. [40] was already compared in section 4.1.2.

Conclusions and outlook
In this paper, we've computed the 2 PM conservative Hamiltonian up to quadratic order in total spin for binary spinning black holes without resorting to Feynman diagrams. As NLO corrections to the effective action that are cubic and quartic in spin PN expansion are given in [80,81], it will be interesting to compare at the observable level. It will certainly be interesting to extend this to general spinning two bodies, which has been explored to quadratic order in spin [61,89,90] based on Feynman diagrams. As discussed in appendix B, direct implementation of the BCFW recursion for non-minimal coupling results in an expression for the Compton amplitudes with spurious poles, obstructing application of on-shell methods to the gravitational two-body problem. It would be interesting to develop a systematic way of cancelling the spurious poles in the Compton amplitude, and to compare with that based on EFT Feynman rules.
Certainly the most interesting theoretical question is the Kerr Compton amplitude beyond spin-2. Recently a new spin- 5 2 amplitude was proposed based on double copy structure [113]. It would be interesting to see if there are any further properties that singles out this amplitude. For example: When compared to the BCFW result, is the difference minimal in the basis of local contact terms constructed along the lines of ref. [114]? Does it have the best Regge behaviour at high energies?
On the other hand, it may very well be that these contact term ambiguities, which by definition is of order G, will not contribute classically. Indeed tidal operators for scalars, which appear as contact terms in the Compton amplitude, leads to modification to the triangle coefficients that has distinct q-scalings [54,55]. Assuming the same behaviour carries over to spinning cases, unless the contact terms are allowed to have extra inverse powers of they will lead to quantum effects and become irrelevant for classical observables. If this is the case, then one only needs a systematic, and preferably economic, way of removing the spurious poles.

A Spin factor composition
The main objectives of this appendix are to prove representation independence of the expression eq.(2.4) and to evaluate composition of the spin factors in amplitude products of the form M 3 × M 3 , which should be evaluated as operator products.
The two relations eq.(A.4) and eq.(A.5) establish eq.(2.4), because all finite spin representations are constructed from tensor products of spin-1 2 representations. Consider the gluing configuration as depicted in figure 3, where all lines are considered on-shell. We attach auxiliary variables z i to the spin factor of eq.(2.4), bearing in mind the 3pt amplitude representation eq.(B.1). The composition of spin factors is an operator product which we work out explicitly in this section. The composition can be used in the BCFW construction of 4pt amplitudes and in triple/quadruple cut computation of the spin factor for amplitude products of the form M 3 × M 3 .
As emphasised in the beginning of this appendix, the LHS is a product of operators and the RHS are 2 × 2 matrix expressions for the operator in chiral and aniti-chiral spin-1 2 representations. The task is to find the b µν that solves the operator relation . (A.8) in both spin- 1 2 representations. The solution to the two matrix equations is (A.9) Setting z 1 = z 4 = 1, the above formula interpolates the two expressions of (2.35) in ref. [40] to arbitrary spin representations.
The rest is straightforward; we only exchange the roles of angles and squares. (A.11)

A.2 Equal helicity case
A.2.1 Helicity configuration 2 + 3 + Without loss of generality we may set the reference spinors of ε 2+ (ε 3+ ) as |3 (|2 ). Working out the factor eq.(A.6) in chiral and aniti-chiral spin-1 2 representations yields First we find a basis for b µν that solves the relation for the two matrix equations. This is just a pair of 2 × 2 matrix equations, which is simple enough to be solved without the Baker-Campbell-Hausdorff formula. After some algebra, one finds the ansatz where the coefficients a i are not simple (rational) functions of kinematics. A formal series expansion in z −1 i yields the solution where the function f (x) is defined after eq.(2.12), reproduced below.
Restoring to massless momenta, the function f (x) can be simplified to f (x) = 1 + O( 2 ) for computing the classical limit.

B Classical spinning Compton amplitude for generic Wilson coefficients
We construct a generalisation of matter-graviton coupling eq.(2.4) which maintains the exponential structure while allowing for non-minimal Wilson coefficients. The generalisation is used to construct Compton amplitudes and to compute scalar integral coefficients for non-minimal Wilson coefficients.

B.1 The residue integral representation and the spin representation of spinning matter
Motivated by the residue integral representation introduced in ref. [78], we deform the 3pt amplitude eq.(2.2) to incorporate the Wilson coefficients C S n of ref. [73] parametrising how different the coupling is from black holes. We use the notations of ref. [41] for the Wilson coefficients where C S n = 1 is the Kerr black hole limit.
The parameter η denotes the helicity of the massless boson, and the scalar factor M 0,η|h| 3 denotes the coupling to spinless matter; g 0 x |h| or g 0 x −|h| . The integration variable z is a bookkeeping device where the integration contour encircles the origin CCW, and for minimal coupling the Wilson polynomial C S n z n resums to z n = 1 1−z , localising the variable to z = 1 and yielding eq.(2.4). To get to the second line, we substitute the Lorentz generator J µν by the spin generator S µν = − 1 m ǫ µνλσ p 1λ S σ and use 3pt kinematics to write The second line of eq.(B.1) was used in ref. [41] to argue that Kerr black holes couple minimally to gravitons. We remark that while ref. [90] considers the O(S 2 ) order operator corresponding to deviations from the black hole limit to be distinct from the O(S 2 ) order operator corresponding to the black hole limit, they seem to be on an equal footing for eq.(B.1) and ref. [61].
We are interested in the classical spin limit, s → ∞ with S = s fixed, of the expression eq.(B.1). Unlike eq.(2.4), the amplitude eq.(B.1) depends on the spin representation of the polarisation tensor for the spinning matter. Therefore the classical spin limit depends on how the limit is approached, and we have restored polarisation tensors and spin representation in eq.(B.1). For example, consider the following two extreme cases for Rep(s) when coupling to a positive helicity massless particle. where the first line was inserted into eq.(B.1) to obtain the second line. For the polarisation tensors eq.(B.5) the same computation leads to (B.8) While the coefficients g i,CA and g i,PT are distinct, their ratios converge to unity in the limit s → ∞.
Therefore the two sequences of reprentations can be considered equivalent in the classical limit. We use the equality · = for the amplitudes in the following sections and omit unnecessary details.

B.2 BCFW construction of Compton amplitude (opposite helicity)
We use the set-up described in section 2.2 to construct the Compton amplitude. For the s-channel, the product of spin factors is given by inserting shifted variables into eq.(A.9) (B.10) The bookkeeping parameter z 1 is associated with leg 1 and z 4 is associated with leg 4. The s-channel propagator pole is absorbed by the scalar factor as in section 2.2. Similar to the minimal coupling case, the u-channel spin factor turns out to be the same.
Note that the bookkeeping parameters z 1 and z 4 have been exchanged; z 1 is the parameter associated with leg 4 and z 4 is the parameter associated with leg 1. Adding up all contributions, the 4pt amplitude M 4 (p 1 , k + 2 , k − 3 , p 4 ) is given as with the relevant definitions in eq.(B.10). This expression should be understood as a sum over symmetric polynomials of J with coefficients determined by residue integral of bookkeeping parameters z 1 and z 4 . In other words, the exponential factor containing J should be expanded before performing the auxiliary residue integrals. Note that for minimal coupling the residue integrals localise to z 1 = z 4 = 1, yielding eq.(2.11) upon resummation.
An important feature of the amplitude eq.(B.12) is that all physical channel residues including the t-channel for graviton coupling are correctly captured to all orders in J, just like the minimally coupled case eq.(2.11). Also, the order of J at which spurious poles start to appear is the same for eq.(B.12) and eq.(2.11).

B.3 BCFW construction of Compton amplitude (same helicity)
We construct the 4pt amplitude M 4 (p 1 , k + 2 , k + 3 , p 4 ) using BCFW recursion in two ways; i) the 23] shift used in section 2.2, and ii) the other viable shift, |2] → |2] + w|3] and |3 → |3 − w|2 (the [23 shift), which we only report the final outcome. The reason we include both shifts is because they lead to different expressions, in contrast to minimal coupling case.
For the 23] shift, the s-channel spin factor is given as where the function f (x) is given in eq.(A. 16). For minimal coupling the bookkeeping parameters localise to z 1 = z 4 = 1 and reduces to eq.(2.12). The spin factor composition eq.(B.13) follows from eq.(A.15) with shifted momentak 2 andk 3 , where z −1 i are used as formal expansion parameters 16 . The u-channel spin factor is given as (B.14) 16 Although the integration contours for zi are small circles encircling the origin, the formal expansion is justified because residue integrals in zi are bookkeeping devices and have no physical content. 1), unless we impose the minimal coupling condition z 1 = z 4 = 1 which removes all spurious poles as in eq.(2.14).
4. Similar to the Kerr case, the full form of f (x) factor is only relevant for the classical spin limit eq.(2.3) and the factor can be simplified to f (x) = 1 + O( 2 ) when computing the classical limit of the integral coefficients. This limit was used at the integrand level when computing the coefficients given in appedix B.4.2.
5. When computing triangle integral coefficients, the classical contributions seem to vanish when (J) n≥2 -order expansion term of eq.(B.15) is involved, subject to the condition x 1 + x 2 = 1. This behaviour was observed in computation of c Based on the last observation, we expect spurious pole cancellation terms to play an important role for classical scattering beyond binary black holes. Systematic removal of spurious poles will be treated in a future work.

B.4 Partial results for integral coefficients in the classical limit
The results presented in this section do not contain the Hilbert space matching factor. The integral coefficients were evaluated using the cut parametrisations of ref. [79] summarised in section 3. The notation for scalar integrals follow that of section 3.

B.4.1 Box integral coefficients
The integral coefficients are normalised by iM = (16πG) 2 c I + · · · where elipsis denotes triangle and quantum contributions. The coefficients are evaluated to the leading order, which is sufficient for checking IR divergence cancellation and eikonal exponentiation. The orthogonal vector n µ is defined as n µ = ǫ αβγµ p 1α p 2β q γ . (B.25)

B.4.2 Triangle integral coefficients
The integral coefficients are normalised by iM = (16πG) 2 c △ I △ + · · · where elipsis denotes other contributions. The part of the Compton amplitude eqs.(B.15-B.16) needed for the S n 1 S 0 2 and S n 1 S 1 2 sectors do not have spurious poles, and we compute full results for those sectors. The orthogonal vector n µ is defined as n µ = ǫ αβγµ p 1α p 2β q γ . Only the sectors containing nontrivial Wilson coefficients C S n≥2 are shown. The Although the results eqs.(B.31-B.36) are not full results, the coefficients provide full information for black holes C S n = 1 since the same helicity graviton exchange contributions c △,±± do not contribute at this order in [34]. The vanishing of same helicity exchange contributions was confirmed explicitly for the sectors S n 1 S 0 2 and S n 1 S 1 2 provided in eqs.(B.26-B.30). The following remarks could be of interest. 2. The Kerr limit C S n = 1 kills all terms with spin contractions of the form (q · S) in the triangle coefficients, where all spin contractions were expressed using the basis {(q · S) , (p · S) , (n · S)}. This partly explains the pattern observed in refs. [61,79], where terms of the form q 2 S 2 and (q · S) 2 appeared in the Hamiltonian only through the combination in the Kerr limit. This is the combination obtained by applying Levi-Civita identities to (n · S i )(n · S j ). The observation extends this pattern to higher spin orders.
C Solving MPD equations iteratively to linear order in test-body spin

C.1 Spin vector description of MPD equations
The Mathisson-Papapetrou-Dixon (MPD) equations describe motion of a spinning body on a curved background, and leading order equations are given as [108][109][110] where S µν is the rank-2 spin tensor, p µ is the 4-momentum, and u µ = Dx µ /Dτ is the 4-velocity of the particle moving on the worldline x µ (τ ) where τ is the proper time. D Dτ denotes parallel transport on the worldline which will be abbreviated by a dot above the tensor, i.e.ȧ µ = Da µ /Dτ = u ν ∇ ν a µ .
Using the covariant SSC (also known as the Tulczyjew-Dixon SSC [110,116]) S µν p ν = 0 and only keeping leading and sub-leading orders in S, we find which means p µ = mu µ + O(S 2 ) where m 2 = p 2 is the rest mass. For the covariant SSC the spin tensor can be described using the spin vector instead, which is given by the relations and the equations of motion for the spin vector becomes 17 which is just the Fermi-Walker transport 18 along the worldline x µ (τ ) when neglecting O(S 4 ) terms. Using the rescaled spin s µ = S µ /m, the equations of motion are given aṡ u µ = 1 2 R µ αβγ ǫ βγδσ u α u δ s σ + O(s 2 ) , s µ = −u µu ν s ν + O(s 4 ) . (C.5)
R µ ναβ = R µ νλσ ǫ λσ αβ . (C.9) The raising and lowering of indices are done with respect to the full metric g µν and g µν up to this point. The exact differential proper time dτ is given as The boundary conditions for the spin deviations σ µ (i) are simply σ µ (i) (τ 0 → −∞) = 0 . (C.11) However, we cannot use the same boundary conditions for the position deviations ξ µ (i) : There is an IR divergence due to long-ranged nature of Coulomb-like potentials coming from integrating velocity deviations at large distances ( v − v 0 ) ∝ 1 r . The resolution used in refs. [90,117] and related works was to fix the boundary conditions at τ 0 = 0 instead, which is defined as the time of closest encounter. This approach is rather inconvenient for our purposes, since we need to invert the relations between variables at τ 0 = 0 and variables at τ 0 → −∞ to specify the initial scattering conditions. The approach we take is to allow arbitrary "subtractions" (0) ξ µ (i) of logarithmic dependence on τ 0 which do not change the trajectory of the test particle at τ 0 → −∞; (C.12) The conditions ensure that initial scattering conditions are unchanged; dξ µ (i) dτ 0 (τ 0 → −∞) = 0. An analogous boundary condition was used in ref. [118].
For computing O(G 2 ) order scattering variables, we only need to perform the subtraction to ξ µ (1) . The logarithmically divergent terms coming from indefinite integration of dξ µ (1) dτ 0 were observed to be one of the two combinations, both of which we separate out as (0) ξ µ (i) .
For the subleading order accelerations, we use the condition u 0 ·ξ (1) = 0 to simplify the equations. The position acceleration is given as; while the spin evolution equation is given as; We only integrate eq.(C.20) once, since we are only interested in scattering variables at O(G 2 ) order.
20 For brevity, we drop the subscripts of initial velocity and impact parameter in this section.