Guillaume Faye

Luc Blanchet

Alessandra Buonanno

Higher-order spin effects in the dynamics of compact binaries I. Equations of motion

Keywords: numbers: 04, 30, -w, 04, 25, -g

We derive the equations of motion of spinning compact binaries including the spin-orbit (SO) coupling terms one post-Newtonian (PN) order beyond the leading-order effect. For black holes maximally spinning this corresponds to 2.5PN order. Our result for the equations of motion essentially confirms the previous result by Tagoshi, Ohashi and Owen. We also compute the spinorbit effects up to 2.5PN order in the conserved (Noetherian) integrals of motion, namely the energy, the total angular momentum, the linear momentum and the center-of-mass integral. We obtain the spin precession equations at 1PN order beyond the leading term, as well. Those results will be used in a future paper to derive the time evolution of the binary orbital phase, providing more accurate templates for LIGO-Virgo-LISA type interferometric detectors.

I. INTRODUCTION

The laser interferometer gravitational-wave (GW) detectors LIGO, Virgo, GEO 600 and TAMA300 are currently searching for GWs emitted by inspiralling compact binaries composed of neutron stars and/or black holes. Analyzing the data using matched filtering technique requires a high-precision modelling of the inspiral waveform [1,2,3,4,5,6,7,8,9]. The post-Newtonian (PN) approximation to general relativity has been applied to build accurate theoretical templates up to the 3.5PN precision level1 for non-spinning compact bodies [10,11,12]. Post-Newtonian templates are currently used in analysing the data with ground-based detectors and in the future they will be used to detect GWs emitted by supermassive black-hole binaries with the space-based detector LISA.

Astrophysical observations suggest that black holes can have non-negligible spins, e.g., due to spin up driven by accretion from a companion during some earlier phase of the binary evolution. For a few black holes surrounded by matter, observations indicate a significant intrinsic angular momentum (see, e.g., Refs. [13,14,15] for stellar black holes and Refs. [START_REF] Fabian | Kerr Spacetime: Rotating Black Holes in General Relativity[END_REF][START_REF] Levin | [END_REF] for supermassive black holes). The spin may even be close to its maximal value [18]. Very little is known however about the black-hole spin magnitudes in binary systems [19].

To successfully detect GWs emitted by spinning, precessing binaries and to estimate the binary parameters, spin effects should be included in the templates. For maximally spinning compact bodies the spin-orbit coupling (linear in the spins) appears dominantly at the 1.5PN order, while the spin-spin one (which is quadratic) appears at 2PN order. The spin effect on the free motion of a test particle was first obtained in the form of a coupling to curvature by Papapetrou et al. [20,21,22]. Seminal works by Barker and O'Connell [23,24] yielded both the leading order spin-orbit and spin-spin contributions in the PN equations of motion. More recently, using an effective field theory approach [25], leading spin-orbit and spin-spin couplings in the two-body Hamiltonian were re-derived [26] and predictions for spin-spin couplings at 3PN order in the spin potential were obtained [27]. Based on the works [23,24], Kidder, Will and Wiseman [28,29] (see also Refs. [30,31]) computed the corresponding coupling terms in the radiation field, enabling thereby the derivation of the orbital phase evolution, the latter being the crucial quantity that determines the templates. Currently, only the leading order spin effects, i.e., the spin-orbit and spin-spin couplings have been implemented in the templates for spinning, precessing black-hole binaries [8,32,33,34,35].

More recently, Tagoshi, Ohashi and Owen [36,37] started the computation of the 1PN corrections to the leading spin-orbit coupling. Those corrections, linear in the spins, appear at 2.5PN order. However, their work has never been completed: the very important conserved integrals associated to the equations of motion at 2.5PN order and the mass quadrupole moment at the 2.5PN order were not computed.

The aim of the present paper together with its companion [38] is to complete the work of Refs. [36,37] and get the orbital phase evolution at 2.5PN order. In this paper we derive the equations of motion, confirming the main result of Ref. [37] (but correcting several important misprints) and compute the entire set of conserved Noetherian integrals of the motion associated with the Poincaré invariance, notably the energy and the total angular momentum. In Ref. [38] (henceforth paper II) we evaluate the multipole moments and the radiation field so as to deduce the orbital phase evolution.

The spin of a rotating body is of the order S true ∼ m a v spin , where m and a denote the mass and typical size of the body respectively, and where v spin represents the velocity of the body's surface. Here, by "true", we mean that the spin we are referring to is not rescaled [as in Eq. (1.1) below]. In this paper we shall consider bodies which are both compact, a ∼ G m c 2 , and maximally rotating, v spin ∼ c. For such objects the magnitude of the spin is roughly S true ∼ G m 2 c . The previous estimate shows that the spin goes as one power of 1/c, i.e., from the PN point of view, it is formally of order 0.5PN. Again, such a counting is appropriate for maximally rotating compact objects. It is then also customary to introduce a dimensionless spin parameter, generally denoted by χ, defined by S true = G m 2 c χ. In our computation the use of such parameter χ is not very convenient because it forces us to introduce some unwanted powers of the mass in front of the spins. On the other hand, it is useful to keep track of the correct PN order by counting all the powers of 1/c. Accordingly we shall "artificially" make explicit the factor 1/c in front of the spin by posing S true = S/c where S will be considered to be of "Newtonian" order. Hence, we shall denote the spin variable by S = c S true = G m 2 χ .

(1.1) Such a notation displays explicitly all powers of 1/c for maximally rotating compact objects. Notably, the spin-orbit (SO) effect always carries a factor 1/c 3 in front, so that it is regarded as being of order 1.5PN, while the spin-spin (SS) effect appears at order 2PN in our terminology, and the 1PN correction to the spin-orbit is 2.5PN order. This PN counting for spin effects corresponds to the standard practice when defining the templates of LIGO/Virgo and LISA detectors (see Refs. [7,8]).

For slowly rotating compact objects (v spin ≪ c) the spin is formally of higher order, namely

S true ∼ G m 2 v spin c 2
∼ 1/c 2 , hence the spin-orbit and spin-spin couplings are pushed at the 2PN and 3PN levels respectively. The 1PN correction to the spin-orbit manifests itself at the same level as the spin-spin coupling, namely 3PN. Of course all the computations in this paper and paper II [38] are still valid in the case of slow rotation, but in this case the spin terms are expected to be numerically smaller, and comparable to higher-order PN contributions.

This paper is organized as follows. In Sec. II we describe the stress-energy tensor of spinning point particles and review some relevant features of the spin formalism before defining our spin variables. In Sec. III we recall some general expressions of the PN metric and equations of motion, which are valid for arbitrary extended matter configurations. The PN metric is parametrized by certain elementary potentials computed in Sec. IV. Our final results for the spin-orbit terms in the equations of motion at the 2.5PN order are presented in Sec. V in a general frame. They are also specialized to the center-of-mass frame and reduced to circular orbits. The precessional equations for the spins including the 1PN relative correction are derived in Sec. VI. Finally, in Sec. VII, we obtain the spinorbit contributions to the conserved integrals associated with our 2.5PN dynamics. The two Appendices are devoted to some tests of our results.

II. STRESS-ENERGY TENSOR FOR SPINNING POINT-PARTICLES

Our calculations are based on the standard model of point-particles with spins [20,21,22,23,24,39,40,41,[START_REF] Dixon | Isolated systems in general relativity[END_REF][START_REF] Bailey | [END_REF]44,45,46]. In the Dixon formulation [START_REF] Dixon | Isolated systems in general relativity[END_REF], the stress-energy tensor,

T µν = T M µν + T S µν , (2.1)
is the sum of the "monopolar" (M) piece, which is a linear combination of monopole sources, i.e. made of Dirac delta-functions, plus the "dipolar" or spin (S) piece, made of gradients of Dirac delta-functions. The four-dimensional formulation of the monopolar part reads as

T M µν = c 2 A +∞ -∞ dτ A p (µ A u ν) A δ (4) (x -y A ) √ -g A , (2.2)
where δ (4) is the four-dimensional Dirac function. The world-line of particle A (A=1,2), denoted y µ A , is parametrized by the particle's proper time τ A . The four-velocity is given by c u µ A = dy µ A /dτ A , and normalized to g A µν u µ A u ν A = -1, where g A µν ≡ g µν (y A ) denotes the metric at the particle's location (the determinant of the metric at point A being denoted by g A ).

The four-vector p µ A is the particle's linear momentum satisfying Eqs. (2.4)-(2.5) below. The dipolar or spin part of the stress-energy tensor, which vanishes in the absence of spins, is2

T S µν = -c A ∇ ρ +∞ -∞ dτ A S ρ(µ A u ν) A δ (4) (x -y A ) √ -g A , (2.3)
where ∇ ρ is the covariant derivative associated with the metric g µν at the field point x, and the anti-symmetric tensor S µν A represents the spin angular momentum for particle A.

The momentum-like quantity p µ

A is a time-like solution of the equation

DS µν A dτ A ≡ cu ρ A ∇ ρ S µν A = c 2 (p µ A u ν A -p ν A u µ A ) . (2.4)
The equation of motion of the particle with spin, equivalent to the covariant conservation law of the total stress-energy tensor, namely ∇ ν T µν = 0, is given by the Papapetrou equation [20,21,22]

Dp µ A dτ A = - 1 2 S λρ A u ν A R µ A νλρ . (2.5)
The Riemann tensor is evaluated at the particle's position A, R µ A νλρ ≡ R µ νλρ (y A ). The equation of motion (2.5) can also be derived directly from the action principle of Bailey and Israel [START_REF] Bailey | [END_REF].

It is well-known that a choice must be made for a supplementary spin condition (SSC) in order to fix unphysical degrees of freedom associated with some arbitrariness in the definition of S µν . This arbitrariness can be interpreted, in the case of extended bodies, as a freedom in the choice for the location of the center-of-mass worldline of the body, with respect to which the angular momentum is defined (see e.g. [29] for discussion). In this paper we adopt the covariant supplementary spin condition

S µν A p A ν = 0 , (2.6)
which allows the natural definition of the spin four-vector S A µ in such a way that

S µν A = - 1 √ -g A ε µνρσ p A ρ m A c S A σ , (2.7)
where ε µνρσ is the four-dimensional antisymmetric Levi-Civita symbol such that ε 0123 = 1.

For the spin vector S A µ itself, we choose a four-vector which is purely spatial in the particle's instantaneous rest frame, where u µ A = (1, 0), hence the components of S A µ are (0, S A ) in that frame. Therefore, in any frame, 3 S A µ u µ A = 0 .

(2.8)

As a consequence of the supplementary spin condition (2.6), we easily verify that d(S µν A S A µν )/dτ A = 0 hence the spin scalar is conserved along the trajectories: S µν A S A µν = const. Furthermore, we can check, using (2.6) and also the Papapetrou law of motion (2.5), that the mass defined by m 2 A c 2 = -p µ A p A µ is indeed constant along the trajectories: m A = const. Finally, the relation linking the four-momentum p µ A and the four-velocity u µ A is readily deduced from the contraction of (2.4) with the four-momentum, which results in

p µ A (pu) A + m 2 A c 2 u µ A = 1 2c 2 S µν A S λρ A u σ A R A νσλρ , (2.9)
where (pu) A ≡ p A ν u ν A . Contracting further this relation with the four-velocity one deduces the expression of (pu) A and inserting it back into (2.9) yields the desired relation between p µ A and u µ A . Let us from now on focus our attention on spin-orbit interactions, which are linear in the spins, and therefore neglect all quadratic and higher corrections in the spins, say O(S 2 ). Drastic simplifications of the formalism occur in the linear case. Since the right-hand-side (RHS) of Eq. (2.9) is quadratic in the spins, we find that the four-momentum is linked to the four-velocity by the simple proportionality relation

p µ A = m A cu µ A + O(S 2 ) .
(2.10) Hence, Eq. (2.6) becomes

S µν A u A ν = O(S 3 ) . (2.11)
On the other hand, the equation of evolution for the spin, also sometimes referred to as the precessional equation, follows immediately from the relationship (2.4) together with the law (2.10) as DS µν A /dτ A = O(S 2 ), or equivalently

DS A µ dτ A = O(S 2 ) . (2.12)
This is simply the equation of parallel transport, which means that the spin vector S µ A remains constant in a freely falling frame, as could have been expected beforehand. Of course Eq. (2.12) preserves the norm of the spin vector, S A µ S µ A = const. When performing PN expansions it is necessary to use three-dimensional like expressions (instead of four-dimensional) for the stress-energy tensor. The field point is accordingly 3 The alternative choice S A µ p µ A = 0 is equivalent to S A µ u µ A = 0 modulo cubic terms in the spins O(S 3 ) (see below) which are neglected in the present paper. Such choices are also adopted in Refs. [28,29,36,37].

denoted by x = (c t, x), and similarly the source points are denoted y A = (c t, y A ). The particle trajectories are considered as functions of the coordinate time t = x 0 /c, say y A (t), and we introduce the ordinary (coordinate) velocity v µ A (t) = dy µ A /dt, also a function of coordinate time. Using Eq. (2.10) we can write the monopolar part (2.2) of the stressenergy tensor as

T M µν = T NS µν + O(S 2 ) , (2.13)
where NS T µν is just the standard piece appropriate to point masses without spins, which reads, in three-dimensional form,

T NS µν = A m A v µ A v ν A -g A ρσ v ρ A v σ A /c 2 δ(x -y A ) √ -g A . (2.14)
We have referred to this part of the stress-energy tensor as the "non-spin" contribution (NS) in spite of its implicit dependence on the spins through the metric tensor. Here δ ≡ δ (3) is the three-dimensional Dirac function. Similarly, the spin part of the stress-energy tensor, Eq. (2.3), can be re-written as

T S µν = - 1 c A ∇ ρ S ρ(µ A v ν) A δ(x -y A ) √ -g A , (2.15)
where the spin tensor S µν A (t) is now considered to be a function of coordinate time, like for the ordinary velocity v µ A (t). The covariant derivative ∇ ρ acts on x, which appears in the argument of the delta-function as shown in (2.15), and on time t through the timedependence of the positions y µ A (t), velocities v µ A (t) and spins S µν A (t). It is easy to further obtain the more explicit expression

√ -g T S µν = - 1 c A ∂ ρ S ρ(µ A v ν) A δ(x -y A ) + S ρ(µ A Γ ν)A ρσ v σ A δ(x -y A ) , (2.16)
where Γ νA ρσ ≡ Γ ν ρσ (y A ) denotes the Christoffel symbol evaluated at the source point A, and where one should notice that the square-root of the determinant √ -g in the left-hand-side (LHS) is to be evaluated at the field point (t, x), contrarily to the factor 1/ √ -g A in the RHS of Eq. (2.15) which is to be computed at the source point y A = (ct, y A ). The explicit form (2.16) of the spin stress-energy tensor is used in all our practical calculations.

In terms of three-dimensional variables the spin tensor reads [after taking into account the spin condition (2.8), namely

S A 0 = -S A i v i A /c] S 0i A = - 1 √ -g A ε ijk u A j S A k , (2.17a)
S ij A = - 1 √ -g A ε ijk u A 0 S A k + u A k v l A c S A l , (2.17b)
where ε ijk is the ordinary Levi-Civita symbol such that ε 123 = 1. Here, we have

u A 0 = u 0 A g A 00 + g A 0i v i A c , (2.18a)
u A j = u 0 A g A j0 + g A jk v k A c , (2.18b
)

with u 0 A = 1 -g A ρσ v ρ A v σ A /c 2 .
(2.18c)

In principle we could adopt as the basic spin variable the covariant vector (or covector) S A i . However, we shall instead use systematically the contravariant components of the vector S i A , which are obtained by raising the index on S A k by means of the spatial metric γ ik A , which denotes the inverse of the covariant spatial metric evaluated at point A, γ A kj ≡ g A kj (i.e. such that γ ik A γ A kj = δ i j ). Hence we define (and systematically use in all our computations)

S i A ≡ γ ik A S A k ⇐⇒ S A i ≡ γ A ij S j A . (2.19)
Beware of the fact that the latter definition of the contravariant spin variable S i A differs from the possible alternative choice g iν A S A ν . The spin vector S i A as defined by (2.19) agrees with the choice already made in Refs. [36,37].

III. POST-NEWTONIAN METRIC AND EQUATIONS OF MOTION

The starting point is the general formulation, i.e. valid for any matter stress-energy tensor T µν with spatially compact support, of the PN metric and equations of motion at 2.5PN order, as worked out in Ref. [47]. In harmonic (or De Donder) coordinates4 the 2.5PN metric is expressed in terms of certain "elementary" potentials as

g 00 = -1 + 2 c 2 V - 2 c 4 V 2 + 8 c 6 X + V i V i + V 3 6 + O 1 c 8 , (3.1a
)

g 0i = - 4 c 3 V i - 8 c 5 Ri + O 1 c 7 , (3.1b
)

g ij = δ ij 1 + 2 c 2 V + 2 c 4 V 2 + 4 c 4 Ŵij + O 1 c 6 . (3.1c)
These potentials, V , V i , • • • , are defined by some retarded integrals of appropriate PN iterated sources. To define them it is convenient to introduce the matter source densities

σ = T 00 + T kk c 2 , (3.2a
)

σ i = T 0i c , (3.2b
)

σ ij = T ij (3.2c) (with T kk ≡ δ ij T ij ).
Then, with -1 R denoting the usual flat space-time retarded operator, we have for the Newtonian like potential V ,

V = -1 R -4πGσ ≡ G d 3 x ′ |x -x ′ | σ (x ′ , t -|x -x ′ |/c) . (3.3a)
The higher-order PN potentials read

V i = -1 R -4πGσ i , (3.3b) Ŵij = -1 R -4πG(σ ij -δ ij σ kk ) -∂ i V ∂ j V , (3.3c) Ri = -1 R -4πG(V σ i -V i σ) -2∂ k V ∂ i V k - 3 2 ∂ t V ∂ i V , (3.3d) X = -1 R -4πGV σ ii + 2V i ∂ t ∂ i V + V ∂ 2 t V + 3 2 (∂ t V ) 2 -2∂ i V j ∂ j V i + Ŵij ∂ 2 ij V . (3.3e)
All these potentials are subject, up to the required PN order, to the differential identities

∂ t V + 1 c 2 1 2 Ŵii + 2V 2 + ∂ i V i + 2 c 2 Ri + V V i = O 1 c 4 , (3.4a)
∂ t V i + ∂ j Ŵij - 1 2 δ ij Ŵkk = O 1 c 2 , (3.4b)
which are consequences of the harmonic coordinate conditions; see Ref. [47].

In this paper we shall specialize the latter PN metric to systems of particles with spin. In this case, as we have reviewed in Sec. II, the stress tensor is the sum of the non-spin piece given by (2.14) and of the spin part (2.15), thus T µν = NS T µν + S T µν . Henceforth we often do not indicate the neglected O(S 2 ) terms. Hence, the source densities (3.2) will be of the form σ µν = NS σ µν + S σ µν , and all the potentials will thus admit similar decompositions, say

V = V NS + V S , • • • , (3.5a) Ŵij = Ŵ NS ij + Ŵ S ij , • • • . (3.5b)
The equations of motion of spinning particles are obtained from the covariant conservation of the total stress-energy tensor,

0 = ∇ ν T µν = ∇ ν T NS µν + ∇ ν T S µν + O(S 2 ) . (3.6)
To get the acceleration of the A-th particle, we insert into the conservation law (3.6) the expressions (2.14) and (2.15) of the stress tensor, integrate over a small volume surrounding the particle A (excluding the other particles B), and use the properties of the Dirac deltafunction. More precisely, in order to handle the delta-function, we systematically apply the rules appropriate to Hadamard's partie finie regularization and given by Eq. (4.6) below. As a result we obtain the equations of motion of the particle A and find useful to write them in the form

dP A µ dt = F A µ , (3.7)
where both the "linear momentum density" P A µ and "force density" F A µ (per unit mass) involve a non-spin piece (NS) and the spin part (S),

P A µ = P NS A µ + P S A µ (3.8a) and F A µ = F NS A µ + F S A µ . (3.8b)
The non-spin parts correspond to the geodesic equations and read

P NS A µ = v ν A g A µν -g A ρσ v ρ A v σ A /c 2 , ( 3.9)
F NS A µ = 1 2 v ν A v λ A (∂ µ g νλ ) A -g A ρσ v ρ A v σ A /c 2 . (3.10)
Their complete expressions in terms of the elementary potentials (3.3) were given in Ref. [47].

We shall need them for a spatial index (µ = i) and for completeness we report here the result (see Eqs. (8.3) in [47])

P NS A i = v i A + 1 c 2 -4V i + 3V v i + 1 2 v 2 v i A + 1 c 4 -8 Ri + 9 2 V 2 v i + 4 Ŵij v j -4V V i + 7 2 V v 2 v i -2v 2 V i -4v i v j V j + 3 8 v i v 4 A + O 1 c 6 , (3.11a)
F NS A i = (∂ i V ) A + 1 c 2 -V ∂ i V + 3 2 v 2 ∂ i V -4v j ∂ i V j A + 1 c 4 4∂ i X + 8V j ∂ i V j -8v j ∂ i Rj + 9 2 v 2 V ∂ i V + 2v j v k ∂ i Ŵjk -2v 2 v j ∂ i V j + 7 8 v 4 ∂ i V + 1 2 V 2 ∂ i V -4v j V j ∂ i V -4v j V ∂ i V j A + O 1 c 6 . (3.11b)
These expressions are still valid in the present situation, but we have to remember that the elementary potentials therein do involve contributions from the spins, e.g. V = NS V + S V . Therefore it is crucial to compute the spin parts of the potentials and to insert them into the non-spin (geodesic-like) contributions to the equations of motion, Eqs. (3.11). Now the purely spin parts, S P A µ and S F A µ , will produce a deviation from the geodesic motion which is induced by the effect of spins. We have found that they admit the following expressions,

m A c P S A µ = - 1 2c d dt g A µν S 0ν A + 1 2 (∂ ρ g µν ) A S ρν A - 1 2 g A ρν Γ νA µσ S ρ0 A v σ A c - 1 2 g A µν Γ 0A ρσ S ρν A v σ A c , (3.12a) m A c F S A µ = 1 2 (∂ µρ g νσ ) A S ρν A v σ A - 1 2 (∂ µ g νλ ) A Γ νA ρσ S ρλ A v σ A . (3.12b)
To compute them is relatively straightforward because all the metric coefficients and Christoffel symbols therein take their standard non-spin expressions (since we are looking for an effect linear in the spins), and these have already been computed in Ref. [47].

As a check of our calculations we have also used an alternative formulation of the equations of motion, which is directly obtained from the Papapetrou equations of motion (2.5) and reads, at linearized order in the spins,

m A c Du µ A dτ A = - 1 2 S λρ A u ν A R µ A νλρ + O(S 2 ) . (3.13)
We lower the free index µ so as to use the convenient relation

Du A µ /dτ A = du A µ /dτ A - 1 2 u ν A u λ A (∂ µ g νλ ) A .
The resulting equation takes the same form as Eq. (3.7),

dP A µ dt = F A µ , (3.14)
but with some distincts linear momentum and force densities P A µ and F A µ . It is clear that the non-spin parts, corresponding to geodesic motion, can be taken to be exactly the same as in our previous formulation, namely Eqs. (3.11). However the spin parts are different; they are given in terms of the Riemann tensor R A µλστ ≡ R µλστ (y A ) as follows,

m A c P S A µ = 0 , (3.15a) m A c F S A µ = R A µλστ ε νρστ 2 g A g A πǫ v π A v ǫ A v λ A g A νω v ω A S A τ , (3.15b)
where S A τ is the covariant spin covector appearing in (2.19). The difference with Eqs. (5.1-3) in Tagoshi et al. [37] is due to the fact that these authors work on the contravariant version of the Papapetrou equation. The advantage of the formulation (3.15) over the previous one (3.12) is of course that it is manifestly covariant. This advantage is however relatively minor in practical PN calculations, since the manifest covariance of the equations is anyway broken from the start. It remains that the two formulations are very useful, and their joint use provides a very good check of the calculations.

The quantity (3.15b) can be computed from the 2.5PN metric, by inserting it into the curvature tensor R A µλστ , but we may also express them directly by means of the elementary potentials (3.3). Let us give here the complete result at the required PN order,

m A c F S A i = 1 c 3 ε ijk ∂ j ∂ t V + v l ∂ jl V S k + 2ε jkl ∂ il V v j -V j S k A + 1 c 5 ε ijk ∂ j ∂ t V + v l ∂ jl V S k V + 1 2 v 2 S k -(Sv)v k + ∂ 2 t V + v l ∂ l ∂ t V v j S k + 2∂ l V ∂ l V k S j + S k ∂ j V l + v j S k ∂ l V + ∂ j V ∂ t V -v l ∂ l V S k + ε jkl 2 2∂ j V ∂ l V i -2V j ∂ il V + 2∂ i V ∂ l V j + V ∂ il V j -2v j ∂ l V ∂ i V + v j V ∂ il V + v j v m ∂ lm V i -v j v m ∂ il V m + v j ∂ l ∂ t V i -v l ∂ i ∂ t V j + v m ∂ il Ŵjm -v m ∂ lm Ŵij -∂ l ∂ t Ŵij -2∂ il Rj S k + v 2 -∂ il V j + v j ∂ il V S k + 2(Sv)v k ∂ il V j A .
(3.16)

IV. COMPUTATION OF THE SPIN PARTS OF ELEMENTARY POTENTIALS

We shall compute all the spin parts of the elementary potentials listed in Eqs. (3.3), which are needed for insertion into the "non-spin" parts of the momentum and force densities as defined by Eq. (3.11). Here we do not compute the non-spin parts of the potentials since they are known from Ref. [47].

Let us start by deriving a few lowest-order results. First, it is immediate to see that the non-spin parts of the matter source densities σ, σ i and σ ij , Eqs. (3.2), start at Newtonian order, and that their spin parts start at 0.5PN order ∼ 1/c in the cases of the vectorial σ S i and tensorial densities σ S ij , and only at 1.5PN order ∼ 1/c 3 in the case of the scalar density σ S . Here we are using our counting for the PN order of spins [see Eq. (1.1)], which is physically appropriate to maximally rotating compact objects. With lowest-order precision the expressions of the source densities for two spinning particles read

σ S = - 2 c 3 ε ijk v i 1 S j 1 ∂ k δ 1 + 1 ↔ 2 + O 1 c 5 , (4.1a)
σ S i = - 1 2c ε ijk S j 1 ∂ k δ 1 + 1 ↔ 2 + O 1 c 3 , (4.1b)
σ S ij = - 1 c ε kl(i v j) 1 S k 1 ∂ l δ 1 + 1 ↔ 2 + O 1 c 3 . (4.1c)
The symbol 1 ↔ 2 means adding the same terms but corresponding to the other particle. The Dirac delta-function is denoted by δ 1 ≡ δ(xy 1 ), and ∂ k δ 1 means the spatial gradient of δ 1 with respect to the field point x. The lowest-order potentials are then straightforward to obtain from the fact that ∆(1/r 1 ) = -4π δ 1 (where r 1 ≡ |xy 1 |), and we get

V S = - 2G c 3 ε ijk v i 1 S j 1 ∂ k 1 r 1 + 1 ↔ 2 + O 1 c 5 , (4.2a)
V S i = - G 2c ε ijk S j 1 ∂ k 1 r 1 + 1 ↔ 2 + O 1 c 3 , (4.2b)
Ŵ S ij = - G c ε kl(i v j) 1 S k 1 ∂ l 1 r 1 + G c δ ij ε klm v k 1 S l 1 ∂ m 1 r 1 + 1 ↔ 2 + O 1 c 3 , (4.2c) Ŵ S kk = 2G c ε klm v k 1 S l 1 ∂ m 1 r 1 + 1 ↔ 2 + O 1 c 3 . (4.2d)
At the dominant level only contribute to the potentials some compact-support terms (proportional to the source densities σ S µν ) -notably the non-compact support term ∼ ∂V ∂V in the spin part of the potential Ŵij , Eq. (3.3c), turns out to be negligible.

To find all the spin terms in the equations of motion up to 2.5PN order, we see from Eq. (3.11) that we need V to 2.5PN order and V i at 1.5PN order [i.e. 1PN beyond what is given by (4.2b)], together with Ŵij , Ri and X at order 0.5PN. As we see, the potential Ŵij is already given by Eq. (4.2c) with the right precision. Our first problem is to obtain the compact-support "Newtonian" potential V to the 2.5PN order. Definition (3.3a) shows that the mass density σ, source of V , admits at an arbitrary high PN order the structure

σ = μ1 + μ1 S δ 1 + 1 √ -g ∂ t ν 1 S δ 1 + 1 √ -g ∂ i ν i 1 S δ 1 + 1 ↔ 2 . (4.
3)

The factors μA , S μA , S ν A and S ν i A are functions of the spins and the velocities v i A , and functionals of the metric components or, equivalently at 2.5PN, of the elementary potentials (3.3). Note that though g µν (x, t), by contrast to S µν A (t), depends on the field point, this is not the case of the moment-like quantities entering the square brackets of Eq. (2.15). Each of them, being multiplied by the Dirac distributions δ A , is indeed evaluated at point x = y A , after the Hadamard procedure described below. Thus, it depends on time only (via the point-mass positions y A and velocities v A ). The index S indicates an additional linear dependence in the spin components, but of course, the full spin dependence is more complicated due to the implicit occurrence of S µν A in the potentials themselves. Notably, the effective mass μ1 whose expression in terms of V , V i , Ŵii and v 2 1 can be found in Ref. [47] contains a net contribution due to the spin at the 2.5PN order and given by

μ1 S = m 1 - 1 c 2 V S + 1 c 4 -4 V i S v i 1 -2 Ŵii S 1 + O 1 c 6 , (4.4)
where the value at the particle's location is meant in the sense of Eq. (4.6a). The expressions of the other moments will not be provided here. It is in fact sufficient for our purpose to observe that, as shown by Eq. (4.4), we have (μ 1 ) S + S μ1 = O(1/c 5 ), and that S ν 1 is at least of order O(1/c 7 ) whereas S ν i 1 is of order O(1/c 3 ). As the spin contribution in σ, say S σ, is already of order 1.5PN ∼ 1/c 3 , see Eq. (4.1a), we need to expand the retardations in V only at relative 1PN order, hence

V S = G d 3 x ′ |x -x ′ | σ(x ′ , t) S - G c d 3 x ′ ∂ ∂t σ(x ′ , t) S + G 2c 2 d 3 x ′ |x -x ′ | ∂ 2 ∂t 2 σ(x ′ , t) S + O 1 c 5 . (4.5)
We then substitute the value of σ following from Eq. (4.3). The integrals are evaluated with the help of the formulas

d 3 x ′ F (x ′ ) δ(x ′ -y 1 ) = (F ) 1 , (4.6a)
d 3 x ′ F (x ′ ) ∂ ′ i δ(x ′ -y 1 ) = -(∂ i F ) 1 , (4.6b)
where the values at point y 1 are denoted by parenthesis like for (F ) 1 . These formulas extend the usual formulas of distribution theory, which are valid for a smooth function F with compact support, to singular functions with a finite number of singular points and deprived of essential singularities (see Ref. [48] for full explanations about this generalization). The formulas (4.6) are part of Hadamard's self-field regularization which is systematically employed in the present approach and the one of [49,50]. 5 In the end we are led to

V S (x, t) = G r 1 μ1 S (t) + μ1 S (t) -G ν i 1 S (t) ∂ ′ i 1 -g(x ′ , t) |x -x ′ | 1 + G 2c 2 ∂ 2 t ν i 1 S ∂ i r 1 - G 2c 2 m 1 a i 1 S ∂ i r 1 + 1 ↔ 2 + O 1 c 6 . (4.7)
The final result for V is obtained by replacing the moments and the determinant of the metric at the 2.5PN level by their explicit values derived from the lowest order approximation of the potentials. The computation of S V i is similar to that of S V , though slightly simpler since the counterpart of μ for σ i does not depend implicitly on the spin at the 1.5PN order.

Next we explain how to compute the non-compact (NC) support terms, and we take the example of the particular NC term in the potential Ri given by R

S (NC) i = ∆ -1 -2∂ k V ∂ i V S k + O 1 c 3 . (4.8)
In the source of this term we have to insert the Newtonian approximation of the potential V , which is simply

V = G m 1 r 1 + G m 2 r 2 + O(c -2
), together with the leading-order spin term S V k given previously in Eq. (4.2b). The source being known we are then able to integrate (using the same techniques as in Ref. [47]) and we get R

S (NC) i = G 2 m 1 8c ε ikl S k 1 ∂ l 1 r 2 1 - G 2 m 2 c ε klm S k 1 ∂ 1 il ∂ 2 m g + 1 ↔ 2 + O 1 c 3 , (4.9)
in which g = ln (r 1 + r 2 + r 12 ) (4.10a)

satisfies ∆g = 1 r 1 r 2 . ( 4

.10b)

The crucial fact which enables the latter integration in closed analytic form is the existence of the function g (first introduced by Fock [START_REF] Fock | Theory of space, time and gravitation[END_REF]). This function and its generalizations are extremely useful in the computation of the spinless equations of motion at 2PN and 3PN orders [47,49]. Finally all the necessary spin parts of the potentials are computed by PN iteration, ready for insertion into the non-spin contribution of the equations of motion as given by Eqs. (3.11). For all the potentials we are in agreement with the results reported by Tagoshi et al. [37] in their Appendix. 6 5 Hadamard's regularization is known to yield some ambiguous coefficients in the equations of motion and the radiation field of non-spinning point particles at 3PN order. When using dimensional regularization these ambiguities are seen to be associated with the appearance of poles ∝ 1/ε (or "cancelled" poles) in the dimension of space d = 3 + ε [12]. The PN order considered in the present paper is merely 1PN, since we are computing the 1PN correction to the leading spin-orbit effect. At this order there are no poles; therefore dimensional and Hadamard's regularizations are equivalent. 6 We have however noticed the following misprints in Ref. [37]: in Eq. (A1h) for S Ri , the third term in the first parenthesis of the first line should be +m 2 /(r 12 s 2 ); in Eq. (A1i) for S X, the first term in the parenthesis following (n 12 v 2 ) in the third line must be read -m 2 /(r 12 s 2 ).

V. THE 2.5PN EQUATIONS OF MOTION WITH SPIN-ORBIT EFFECTS A. Equations in a general frame

In addition to the spin parts of the potentials computed in Sec. IV and inserted into Eqs. (3.11), we add the required spin corrections to the geodesic motion as given by either the formulation of Eqs. (3.12) or that of (3.15)- (3.16). The latter corrections are computed by inserting into them the non-spin parts of the potentials taken from [47]. We find that the two formulations [respectively given by (3.12) and (3.15)-(3.16)] are equivalent and agree on the result. Finally the 2.5PN equations of motion with spin-orbit effects are obtained in the form

dv 1 dt = A N + 1 c 2 A 1PN + 1 c 3 A S 1.5PN + 1 c 4 A 2PN + A SS 2PN + 1 c 5 A 2.5PN + A S 2.5PN + O 1 c 6 .
(5.1) Here the Newtonian acceleration is

A i N = -G m 2 r 2 12 n i 12 ,

and we denote by

A i N , A i 1PN , A i 2PN and A i 2.
5PN the standard non-spin contributions (in harmonic coordinates) which are well-known, see Eqs. (8.4) in [47] and earlier works reviewed in [START_REF] Damour | Gravitational Radiation[END_REF]. In particular A i 2.5PN represents the standard radiation reaction damping term. (For simplicity we henceforth suppress the subscript NS on non-spin-type contributions.)

The leading-order spin effect is the 1.5PN spin-orbit term. For this term we recover the standard expression, known from Refs. [23,24] and given in [28,29] in the center-of-mass frame, and in [37] in a general frame. In the following we shall sometimes use some formulas relating the "mixed products" of three vectors in three dimensions,

(U 1 , U 2 , U 3 ) U = (UU 1 ) U 2 × U 3 + (UU 2 ) U 3 × U 1 + (UU 3 ) U 1 × U 2 (5.2a) = (U, U 2 , U 3 ) U 1 + (U 1 , U, U 3 ) U 2 + (U 1 , U 2 , U) U 3 , (5.2b)
valid for any vectors U, U 1 , U 2 , U 3 (in 3 dimensions). Here the vectorial product of ordinary Euclidean vectors is indicated with the × symbol, for instance (U 1 × U 2 ) i = ε ijk U j 1 U k 2 ; parenthesis denote the usual Euclidean scalar product, (UU 1 ) = U i U i 1 = U • U 1 ; and the mixed product, or determinant between three vectors, is denoted (U

1 , U 2 , U 3 ) ≡ U 1 • (U 2 × U 3 ) = ε ijk U i 1 U j 2 U k 3 . This yields A S 1.5PN = Gm 2 r 3 12 6 (S 1 , n 12 , v 12 ) m 1 + 6 (S 2 , n 12 , v 12 ) m 2 n 12 + 3(n 12 v 12 ) n 12 × S 1 m 1 + 6(n 12 v 12 ) n 12 × S 2 m 2 -3 v 12 × S 1 m 1 -4 v 12 × S 2 m 2 .
(5.3a)

We use, whenever convenient, the notation v 12 = v 1v 2 for the relative velocity. The next-order spin correction is the spin-spin (SS) at 2PN order. We do not give this term since we are concerned here with spin-orbit effects which are linear in the spins. The SS term is quadratic in the spins, O(S 2 ), and can be found in Refs. [23,24] and e.g. in Eq. (5.9) of [37]. Now the 1PN correction to the spin-orbit effect, which is the aim of this paper and the work [37], reads (5.3b) Surprisingly, we find that our expression has substantial differences with the result given in Eq. (5.10) of [37]. However, since we recovered in the last section exactly the same potentials as given in the Appendix of [37], and since as we shall see below we find perfect agreement with the equations of motion computed in [37] in the case of the center-of-mass frame, we believe that the latter differences can only be due to some trivial misprints (and most probably to some mixup of Mathematica files) in the last stage of the work [37]. 7 In Appendix A we shall prove that the equations of motion stay invariant under global Poincaré transformations. Such a verification is quite important to test the correctness of the equations (it played an important role during the computation of the 3PN non-spin terms in [49,50]). Furthermore, we show in Appendix B that the test mass limit of the 7 For completeness we indicate here all the misprints in Eq. (5.10) of [37]: in order to recover the correct acceleration, the last term before the closing curly brackets on the fifth line of Eq. (5.10),

A S 2.5PN = Gm 2 r 3 12 n 12 -6(n 12 , v 1 , v 2 ) (v 1 S 1 ) m 1 + (v 2 S 2 ) m 2 - (S
+ v 1 -3 (S 1 , n 12 , v 1 ) m 1 (n 12 v 1 ) + (n 12 v 2 ) + 6(n 12 v 1 ) (S 1 , n 12 , v 2 ) m 1 -3 (S 1 , v 1 , v 2 ) m 1 -6(n 12 v 1 ) (S 2 , n 12 , v 1 ) m 2 + (S 2 , n 12 , v 2 ) m 2 12(n 12 v 1 ) -6(n 12 v 2 ) -4 (S 2 , v 1 , v 2 ) m 2 + v 2 6(n 12 v 1 ) (S 1 , n 12 , v 12 ) m 1 + 6(n 12 v 1 ) (S 2 , n 12 , v 12 ) m 2 -n 12 × v 1 3(n 12 v 12 ) (v 1 S 1 ) m 1 + 4 Gm 1 r 12 (n 12 S 2 ) m 2 -n 12 × v 2 6(n 12 v 12 ) (v 2 S 2 ) m 2 -4 Gm 1 r 12 (n 12 S 2 ) m 2 + v 1 × v 2 3 (v 1 S 1 ) m 1 + 4 (v 2 S 2 ) m 2 + n 12 × S 1 m 1 - 15 2 (n 12 v 12 )(n 12 v 2 ) 2 + 3(n 12 v 2 )(v 12 v 2 ) -14 Gm 1 r 12 (n 12 v 12 ) -9 Gm 2 r 12 (n 12 v 12 ) + n 12 × S 2 m 2 -15(n 12 v 12 )(n 12 v 2 ) 2 -6(n 12 v 1 )(v 12 v 2 ) + 12(n 12 v 2 )(v 12 v 2 ) + Gm 1 r 12 - 35
-12ε jkl n j v k 1 v l 2 (v 1 S 1 )/m 1 , must be replaced by -6ε jkl n j v k 1 v l 2 [(v 1 S 1 )/m 1 +(v 2 S 2 )/m 2 ]
; the term before the last one in the seventh line has to be read -6n j v k 2 (nv 12 )(v 2 S 2 )/m 2 instead of -6n j v k 2 (nv 12 )(v 1 S 1 )/m 1 ; and the very last term +7v j 1 v k 2 (v 1 S 1 )/m 1 must be modified as +v

j 1 v k 2 [3(v 1 S 1 )/m 1 + 4(v 2 S 2 )/m 2 ].
of motion is identical with the geodesic equations around a Kerr black hole (for simplicity we restrict ourself to circular orbits). Both verifications have already been made in Ref. [37] but we present some alternative ways to do the checks.

B. Equations in the center-of-mass frame

Let us now present the result in the center-of-mass (CM) frame, defined by the nullity of the center-of-mass vector, equal to the conserved integral associated with the boost invariance of the equations of motion, which will be checked in Appendix A. We shall derive the center-of-mass integral at the 2.5PN order in the next section, however for the present computation we need it only at the 1.5PN order. When working in the CM frame we find it convenient to introduce the same spin variables as chosen by Kidder [29] (except that we denote by Σ what he calls ∆), namely

S ≡ S 1 + S 2 ,
(5.4a)

Σ ≡ m S 2 m 2 - S 1 m 1 . ( 5
v = dx/dt = v 1 -v 2 = v 12
(see e.g. Ref. [37])

y 1 = m 2 m + ν 2c 2 δm m v 2 - G m r x + ν m c 3 v × Σ , (5.5a
)

y 2 = - m 1 m + ν 2c 2 δm m v 2 - G m r x + ν m c 3 v × Σ .
(5.5b)

In addition to the spin-orbit effect at order 1.5PN ∼ 1/c 3 (last term in these relations), we have included the well-known 1PN ∼ 1/c 2 non-spin term. This term is obviously needed here because during the reduction of the equations of motion to the CM frame at order 2.5PN in the spins, we shall need to take into account the 1PN non-spin term coupled to the lowest-order 1.5PN spin term; such coupling evidently produces some 2.5PN spin terms.

In the CM frame the equation of the relative motion reads

dv dt = B N + 1 c 2 B 1PN + 1 c 3 B S 1.5PN + 1 c 4 B 2PN + B SS 2PN + 1 c 5 B 2.5PN + B S 2.5PN + O 1 c 6 ,
(5.6) where we recognize all the various terms similarly to Eq. (5.1). We find that the spin-orbit term and the 1PN correction to the spin-orbit are given, in terms of the spin variables (5.4), by

B S 1.5PN = G r 3 n 12(S, n, v) + 6 δm m (Σ, n, v) + 9(nv)n × S + 3 δm m (nv)n × Σ -7v × S -3 δm m v × Σ , (5.7a) B S 2.5PN = G r 3 n (S, n, v) -30ν(nv) 2 + 24νv 2 - Gm r (38 + 25ν) + δm m (Σ, n, v) 15ν(nv) 2 + 12νv 2 - Gm r (18 + 29 2 ν) + (nv)v (S, n, v)(-9 + 9ν) + δm m (Σ, n, v)(-3 + 6ν) + n × v (nv)(vS)(-3 + 3ν) -8 Gm r ν(nS) - δm m 4 Gm r ν(nΣ) + 3(nv)(vΣ) + (nv)n × S - 45 2 ν(nv) 2 + 21νv 2 - Gm r (25 + 15ν) + δm m (nv)n × Σ -15ν(nv) 2 + 12νv 2 - Gm r (9 + 17 2 ν) + v × S 33 2 ν(nv) 2 + Gm r (21 + 9ν) -14νv 2 + δm m v × Σ 9ν(nv) 2 -7νv 2 + Gm r (9 + 9 2 ν) . (5.7b)
We find perfect agreement with Eqs. (5.18) and (5.20) of Tagoshi et al. [37]. 8

C. Reduction to quasi circular orbits

Finally we present the case where the orbit is nearly circular, i.e. whose radius is constant apart from small perturbations induced by the spins (as usual we neglect the gravitational radiation damping at 2.5PN order). Following Ref. [29] we introduce an orthonormal triad {n, λ, ℓ} defined by n = x/r as before, ℓ = L N /|L N | where L N ≡ µ x × v denotes the Newtonian angular momentum, and λ = ℓ × n. The orbital frequency ω is defined for general, not necessarily circular orbits, v = ṙn + rωλ where ṙ = (nv). The components of the acceleration a = dv/dt along the basis {n, λ, ℓ} are then given by

n • a = r -rω 2 , (5.8a) λ • a = r ω + 2 ṙω , (5.8b)
ℓ • a = -rω λ • dℓ dt . (5.8c)
We project out the spins on this orthonormal basis, defining S = S n n + S λ λ + S ℓ ℓ and similarly for Σ. Next we impose the restriction to circular orbits which means r = 0 = ṙ and v 2 = r 2 ω 2 (neglecting radiation reaction damping terms). In this way we find that the equations of motion (5.6) with (5.7) are of the type

dv dt = -ω 2 r n + a ℓ ℓ + O 1 c 6 .
(5.9)

8 Note that the spin variables adopted in [37] are defined by

χ s ≡ 1 2 S1 m 2 1 + S2 m 2 2 and χ a ≡ 1 2 S1 m 2 1 -S2 m 2 2
and differ from our own. We have

S = m 2 (1 -2ν)χ s + δm m χ a and Σ = m 2 - δm m χ s -χ a .
is no component of the acceleration along λ. Comparing with Eqs. (5.8) in the case of circular orbits, we see that ω is indeed the orbital frequency, while a ℓ = -rω(λ • dℓ/dt) is proportional to the variation of ℓ in the direction of the velocity v = rωλ. We find that ω 2 is of the form

ω 2 = G m r 3 1 + 1 c 2 ζ 1PN + 1 c 3 ζ S 1.5PN + 1 c 4 ζ 2PN + ζ SS 2PN + 1 c 5 ζ S 2.5PN + O 1 c 6 , (5.10)
where ζ 1PN and ζ 2PN denote the standard non-spin contributions, 9 and where

ζ S 1.5PN = G m r 3/2 1 G m 2 -5S ℓ -3 δm m Σ ℓ , (5.11a)
ζ S 2.5PN = G m r 5/2 1 G m 2 39 2 - 23 2 ν S ℓ + 21 2 - 11 2 ν δm m Σ ℓ , (5.11b)
with, e.g. S ℓ ≡ (Sℓ) = S • ℓ. On the other hand, we get

a ℓ = 1 c 3 α S 1.5PN + 1 c 4 α SS 2PN + 1 c 5 α S 2.5PN + O 1 c 6 , (5.12)
with spin-orbit coefficients

α S 1.5PN = G m r 3/2 1 m r 2 7S n + 3 δm m Σ n , (5.13a) α
S 2.5PN = G m r 5/2 1 m r 2 - 63 2 + ν 2 S n - 27 2 δm m Σ n . (5.13b)
We see that the resulting motion cannot be exactly circular for general orientations of the spins. Let us show however that the time-averaged acceleration coincides with the acceleration of a particle that rotates uniformly about the origin. In a first step, we must make explicit the time dependence of the dynamical variables x, S and Σ. As the motion is uniformly circular in the absence of spin, the position x decomposed along a fixed orthonormal basis {e 1 , e 2 , ℓ} reads x(t) = e 1 r cos(ω NS t) + e 2 r sin(ω NS t) ,

with ω NS being the orbital frequency when the spins are turned off.

The spin variables are computed by means of the precession equations, which decouple in the case of a pure spin-orbit interaction. The spin 1, for instance, obeys an equation whose right-hand-side is polynomial in Gm/r = v 2 , (nS 1 ) and (vS 1 ). For dimensional reasons, it must then have the form (for circular orbits, up to say the 2PN order)

dS 1 dt = k=1,2 Gm rc 2 k a (k,n) S 1 (vS 1 )n + a (k,v) S 1 (nS 1 )v + O 1 c 5 , (5.15)
and similarly for dS 2 /dt. The functions of m 1 /m, m 2 /m denoted by a

(k,n) S 1
and a (k,v) S 1 may be obtained from the results of the next section (see also paper II). They allow us to define dimensionless coefficients like a

(n) S = k=1,2 (Gm/r) k a (k,n) S 1 .
The key point is that the latter coefficients are constant, which suggests to solve the above differential equations in the moving basis {n, λ, ℓ}. Indeed, the time derivative of a spin component in this basis, say S 1 n = (nS 1 ), is given by a relation of the type

dS 1 n dt = (n Ṡ1 ) + ( ṅS 1 ) (5.16)
with ṅ = ẋ/r = ω NS λ. This results, after eliminating Ṡ1 by means of Eq. (5.15), in a linear differential equation with constant coefficients for S 1 n . Proceeding in the same way for the other components of the first spin, we arrive at the following system:

dX S 1 dt = M S 1 .X S 1 , (5.17)
where M S 1 is a 3 × 3 constant matrix and X S 1 = (S 1 n , S 1 λ , S 1 ℓ ). The relations ℓ.dS 1 /dt = ℓ.dS 2 /dt = 0 (since ℓ is constant because we neglect the SS terms) imply that (0, 0, 1) is an eigenvector associated with the eigenvalue λ 0 = 0. There remain two eigenvalues, say λ + 1 and λ - 1 ; but since the trace of M S 1 vanishes because (nv) = 0, we have λ - 1 = -λ + 1 . In the end, we notice that the spins are almost constant at Newtonian order in the basis {e 1 , e 2 , ℓ}, which means that they precess about ℓ with angular velocity -ω NS in the moving frame. Therefore, λ ± 1 is purely imaginary and reduces to ±iω NS at Newtonian order. At higher order we shall have λ

± 1 = ±i(ω NS -Ω 1 )
where Ω 1 = O(1/c 2 ) represents the precession frequency. The components S 1 n and S 1 λ solving Eq. (5.17) are then linear combinations of cos[(-ω NS +Ω 1 )t] and sin[(-ω NS +Ω 1 )t]. 10 As for the component S 1 ℓ , it is constant neglecting terms quadratic in the spins. We complete our proof by time averaging the term a ℓ in the acceleration (5.9). We first observe that the conservative part of the dynamics involves three different angular frequencies (ω NS , Ω 1 and Ω 2 ), so that it cannot be periodic in general. Therefore, it is not appropriate to average the particle motion on the orbital period. Instead, the time-average will be achieved on infinite time. Defining

S 1 n = lim T →+∞ 1 T t+T t dt ′ S 1 n (t ′ ) , (5.18)
we find S 1 n = 0. We next notice that the orbital frequency ω is actually constant (neglecting SS terms), for it depends on the spin through S 1 ℓ and S 2 ℓ only, which are constant. The average of a ℓ is a linear combination of S n = (nS) = 0 and Σ n = (nΣ) = 0; hence it does not contribute: a ℓ = 0.

VI. THE 2PN SPIN-ORBIT EQUATIONS OF PRECESSION

In this Section we give the equations of evolution of the spins, or precession equations, at relative 2PN order, i.e. one PN order beyond the dominant term. The precession equations quite simple to derive from the equation of parallel transport (2.12), which we recall is valid at the linear order in the spins [neglecting O (S 2 )], but at any PN order in that term which is linear in the spins. The PN corrections are easily computed from the non-spin part of the metric and Christoffel symbols computed in Ref. [47]. The precession equations in a general frame take the form

dS 1 dt = 1 c 2 T S 1PN + 1 c 3 T SS 1.5PN + 1 c 4 T S 2PN + O 1 c 5 , (6.1)
together with the equation with 1 ↔ 2. At the lowest order we find

T S 1PN = Gm 2 r 2 12 S 1 (n 12 v 12 ) -2n 12 (v 12 S 1 ) + (v 1 -2v 2 )(n 12 S 1 ) . (6.2)
The above equation is already known [29,37]. See e.g. Eq. ( 4.3) in [37] and the paragraph afterward commenting about the difference with formulations based on an alternative definition for the spin, like that of Ref. [29]. The spin-spin (SS) term is also known but is out of the scope of the present paper (and the parallel transport equation we employ); it can be found elsewhere, see Eqs. ( 2)-( 3) of [8]. Then we find that the next-order spin-orbit term is For completeness and for the benefit of users of these formulas in the data analysis of detectors, we present also the precession equations in the CM frame, using our specific spin variables defined by (5.4). These are

T S 2PN = Gm 2 r 2 12 S 1 (n 12 v 2 )(v 12 v 2 ) - 3
dS dt = 1 c 2 U S 1PN + 1 c 3 U SS 1.5PN + 1 c 4 U S 2PN + O 1 c 5 , (6.4a
)

dΣ dt = 1 c 2 V S 1PN + 1 c 3 V SS 1.5PN + 1 c 4 V S 2PN + O 1 c 5 , (6.4b)
where all the spin-orbit coefficients are given by

U S 1PN = G m ν r 2 n -4(vS) -2 δm m (vΣ) + v 3(nS) + δm m (nΣ) + (nv) 2S + δm m Σ , (6.5a)
U S 2PN = G m ν r 2 n (vS) -2v 2 + 3(nv) 2 -6ν(nv) 2 + 7 Gm r -8ν Gm r -14 Gm r (nS)(nv) + δm m (vΣ) ν -3(nv) -4 Gm r + δm m Gm r (nΣ) (nv) 2 - ν 2 + v (nS) 2v 2 -4νv 2 -3(nv) 2 + 15 2 ν(nv) 2 + 4 Gm r -6ν Gm r + (vS)(nv) (2 -6ν) + δm m (nΣ) - 3 2 νv 2 + 3ν(nv) 2 - Gm r - 7 2 ν Gm r -3 δm m (vΣ) (nv) ν + S (nv) v 2 -2νv 2 - 3 2 (nv) 2 + 3ν(nv) 2 - Gm r + 2ν Gm r + δm m Σ (nv) -νv 2 + 3 2 ν (nv) 2 - Gm r + ν Gm r , (6.5b)
and

V S 1PN = G m r 2 n (vΣ) (-2 + 4ν) -2 δm m (vS) + v (nΣ) (1 -ν) + δm m (nS) + (nv) Σ (1 -2ν) + δm m S , (6.6a)
V S 2PN = G m r 2 n (vΣ) ν -2v 2 + 6ν(nv) 2 + 3 Gm r + 8ν Gm r + Gm r (nΣ) (nv) 2 - 45 2 ν + 2ν 2 + δm m (vS) ν -3(nv) 2 -4 Gm r + δm m Gm r (nS) (nv) 2 - ν 2 + v (nΣ) ν 2 v 2 + 2ν 2 v 2 - 9 2 ν 2 (nv) 2 - Gm r + 9 2 ν Gm r + 8ν 2 Gm r +(vΣ)(nv)ν (-1 + 6ν) -3 δm m (vS) (nv) ν + δm m (nS) - 3 2 νv 2 + 3ν(nv) 2 - Gm r - 7 2 ν Gm r + Σ (nv) 2ν 2 v 2 -3ν 2 (nv) 2 - Gm r + 4ν Gm r -2ν 2 Gm r + δm m S (nv) -νv 2 + 3 2 ν (nv) 2 - Gm r + ν Gm r . (6.6b)
To these expressions one may add the SS terms in the standard way (see Eqs. ( 2)-( 3) of [8]).

VII. SPIN EFFECTS IN THE CONSERVED INTEGRALS OF THE MOTION

Having obtained in Sec. V the equations of motion, the important task is now to deduce from them the complete set of conserved integrals of the motion associated with the global Poincaré invariance of these equations (which has been checked in Ref. [37] and Appendix A below). In principle, the conserved integrals of the motion, which generalize the usual notions of energy, angular and linear momenta, and center of mass position, should be best derived from a Lagrangian. In the present paper, however, we did not attempt to derive complete Lagrangian for the particles with spins (see [START_REF] Bailey | [END_REF] for a discussion on how to formulate Lagrangians with spins); rather, we have obtained the integrals of the motion by "guess-work", starting from their most general admissible form, and then imposing the conservation laws when the equations of motion are satisfied. 11 Here we simply state the results.

The PN expansion of the conserved integral of the energy, namely E such that dE/dt = 0, reads as

E = E N + 1 c 2 E 1PN + 1 c 3 E S 1.5PN + 1 c 4 E 2PN + E SS 2PN + 1 c 5 E S 2.5PN + O 1 c 6 , (7.1)
where the non-spin pieces, E N , E 1PN and E 2PN , are known and can be found e.g. in Ref. [54].

For instance we have

E N = 1 2 m 1 v 2 1 + 1 2 m 2 v 2 2 -G m 1 m 2 r 12
. For the lowest-order spin-orbit effect we find, in agreement with the standard result,

E S 1.5PN = Gm 2 r 2 12 (S 1 , n 12 , v 1 ) + 1 ↔ 2 , (7.2a)
where we employ a special notation for the totally anti-symmetric "mixed product" between three vectors, as given in (5.2). For the spin-orbit contribution at 2.5PN order we find

E S 2.5PN = Gm 2 r 2 12 (S 1 , n 12 , v 1 ) 1 2 v 2 1 -3v 2 2 + 3(n 12 v 1 )(n 12 v 2 ) + 3 2 (n 12 v 2 ) 2 -2 G m 1 r 12 + G m 2 r 12 + (S 1 , n 12 , v 2 ) 2v 2 1 -(v 1 v 2 ) + 2v 2 2 -3(n 12 v 1 ) 2 -3(n 12 v 1 )(n 12 v 2 ) -3 Gm 2 r 12 + (n 12 , v 1 , v 2 ) (v 1 S 1 ) + 2(v 2 S 1 ) + 1 ↔ 2 . (7.2b)
Notice that several equivalent forms can be given to this result. For instance if wished one could introduce the mixed product (S 1 , v 1 , v 2 ) in place of a (n 12 , v 1 , v 2 ) in the last term of (7.2b), making use of linear combinations such as (n

12 v 1 )(S 1 , v 1 , v 2 ) = (n 12 , v 1 , v 2 )(v 1 S 1 )+ (S 1 , n 12 , v 2 )v 2 1 -(S 1 , n 12 , v 1 )(v 1 v 2 )
[a consequence of Eq. (5.2)]. As before we do not give the SS contribution at 2PN order (see [29] for instance).

We give here the corresponding result for the conserved center-of-mass energy in the CM frame:

E = m ν c 2 e N + 1 c 2 e 1PN + 1 c 3 e S 1.5PN + 1 c 4 e 2PN + e SS 2PN + 1 c 5 e S 2.5PN + O 1 c 6 , (7.3)
where e N = 1 2 v 2 -Gm r (see [55] for the other non-spin contributions). The SO coupling terms in the CM frame are found to be e

S 1.5PN = G r 2 (S, n, v) + (Σ, n, v) δm m , (7.4a
) 11 As usual we neglect the radiation reaction effect at 2.5PN order. Indeed we know that such an effect does not depend on the spins. The contribution of the spins in the radiation reaction force comes in at 1.5PN order beyond the dominant effect, which means at the 4PN level, and has been computed in Ref. [START_REF] Will | [END_REF].

Radiation reaction effects will be included into the present formalism when we obtain the contributions of the spins in the GW flux [38].

e S 2.5PN = G r (S, n, v) - 3 2 (1 + ν)v 2 - 3 2 ν(nv) 2 + G m r +(Σ, n, v) 1 2 (1 -5ν)v 2 + 1 2 (2 + ν) G m r δm m . (7.4b)
Let us next deal with the conserved total angular momentum J, i.e. dJ/dt = 0, sum of orbital and spin contributions, which we write as

J = L + 1 c S 1 + 1 c S 2 , (7.5)
where L is the orbital angular momentum, and where S 1 and S 2 are the contravariant spin vectors defined following the specific choice made in Eq. (2.19) [recall also their peculiar dimension which follows from (1.1)]. The angular momentum L admits the PN expansion

L = L N + 1 c 2 L 1PN + 1 c 3 L S 1.5PN + 1 c 4 L 2PN + L SS 2PN + 1 c 5 L S 2.5PN + O 1 c 6 , (7.6)
where all the non-spin pieces are given by Eq. (4.4) of [54]. For instance,

L N = m 1 y 1 × v 1 + m 2 y 2 × v 2
. Now, in order to express in the best way the spin-orbit contributions in L, we find that they must be written in the following way,

L S 1.5PN = y 1 × p S 1 + y 2 × p S 2 + K S 1.5PN , (7.7a)
L S 2.5PN = y 1 × q S 1 + y 2 × q S 2 + K S 2.5PN , (7.7b)
in which we have introduced some convenient notions of the "individual linear momenta" of the particles, say S p 1 and S p 2 at 1.5PN order, and S q 1 and S q 2 at 2.5PN order. The extra terms in the RHS, S K 1.5PN and S K 2.5PN , incorporate all what remains, the point being that they depend on the positions of the particles only through their relative separation, i.e. r 12 = |y 1y 2 | and n 12 = (y 1y 2 )/r 12 . The only dependence of the conserved angular momentum on the individual positions y 1 and y 2 is the one which is given explicitly by the first terms of Eqs. (7.7). 12 The results we find for these momenta are Let us add a comment on the meaning of the conservation of the total angular momentum J at 2.5PN order [Eq. (7.5) with (7.6)]. When differentiating J with respect to time, we generate several spin contributions at 2.5PN order: (i) The "main" one is coming from the differentiation of the Newtonian term L N , and is due to the replacement of the acceleration by the equations of motion (5.1) with (5.3b); (ii) There is the one coming from the differentiation of the 1PN part L 1PN , since the replacement of the accelerations at order 1.5PN [Eq. (5.3a)] therein does also produce some terms at 2.5PN order; (iii) When differentiating the lowestorder spin-orbit term S L 1PN , the derivative of the spins gives other 2.5PN terms via the precessional equations; (iv) When differentiating the spin vectors themselves, S 1 and S 2 , one must make use of the precessional equations with their full 2PN accuracy 14 which are given by Eqs. (6.1)-(6.3). Only when account is taken of all these replacements (i)-(iv) of accelerations and spin precession, does one find that J is conserved, dJ/dt = 0, up to 2.5PN order (neglecting the 2.5PN non-spin radiation reaction damping).

p S 1 = - G m 2 r 2
The orbital angular momentum in the CM frame reads:

L = ν ℓ N + 1 c 2 ℓ 1PN + 1 c 3 ℓ S 1.5PN + 1 c 4 ℓ 2PN + ℓ SS 2PN + 1 c 5 ℓ S 2.5PN + O 1 c 6 , (7.10)
where ℓ N = m x × v; the non-spin contributions can be found in Refs. [54,55,56]. We have

+ (-2 -3ν) G m r 3 (x v) (v S) + v -4 G m r ν (v Σ) + - 1 2 + 2ν v 2 (v Σ) -4ν G m r (v Σ) + 5ν 2 G m r 3 (x Σ) (x v) δm m + v - 1 2 + 3ν 2 v 2 (v S) -7 ν G m r (v S) + (6 ν -3) G m r 3 (x S) (x v) . (7.11b)
Finally let us give the conserved integrals of the linear momentum P and center of mass position G, which are related to each other by dG/dt = P. Recall that the existence of the center-of-mass integral G is a consequence of the boost-invariance of the equations of motion (cf. Appendix A). Both P and G admit a PN expansion exactly like those of E and L. Quite naturally, we find that the spin-orbit contributions in P are simply given by the sum of the "individual" linear momenta for each particles that we found convenient to introduce in order to express the angular momentum in Eqs. (7.7) The derivation of the complete set of integrals of the motion gives us further confidence in the physical soundness of the equations of motion derived in this paper. Those results, together with the analyses performed in Appendices A and B, complete the resolution of the problem of linear spin-orbit effects in the binary's equations of motion at 2.5PN order.

  .4b)Mass parameters are denoted by m ≡ m 1 + m 2 , δm ≡ m 1 -m 2 and ν ≡ m 1 m 2 /m 2 (such that 0 < ν ≤ 1/4). At the leading order in the spins we have the following relation between the positions y 1 and y 2 in the CM frame and the relative position x = y 1y 2 and velocity

9

  They are given by ζ 1PN = Gm r (-3 + ν) and ζ 2PN = Gm r 2 6 + 41 4 ν + ν 2 in harmonic coordinates.

13c 2 r 12 S 1 - 1 2c 2 (v 1 S 1 ) v 1 . 4

 121114 Ref.[29] uses different definitions for the spin variables, which are related to ours by(S 1 ) Kidder = 1 + Gm 2 (x v) 2 (x S)

  1 , n 12 v 12 ) m 1 15(n 12 v 2 ) 2 + 6(v 12 v 2 ) + 26

					Gm 1 r 12	+ 18	Gm 2 r 12
	-	(S 2 , n 12 , v 12 ) m 2	15(n 12 v 2 ) 2 + 6(v 12 v 2 ) +	49 2	Gm 1 r 12	+ 20	Gm 2 r 12

  2 (n 12 v 2 ) 2 (n 12 v 12 ) + 12 (v 12 S 1 ) 3(n 12 v 2 ) 2 + 2(v 12 v 2 ) + Gm 1 r 12 -16(n 12 S 1 )(n 12 v 12 ) + 3(v 1 S 1 ) -7(v 2 S 1 ) 12 S 1 )(n 12 v 2 ) 2 + (v 12 S 1 )(n 12 v 2 ) -m 2 ) + v 2 (n 12 S 1 ) 2(v 12 v 2 ) + 3(n 12 v 2 ) 2+ 2(n 12 v 12 ) (v 1 S 1 ) + (v 2 S 1 ) -5(n 12 S 1 )

		Gm 1 r 12	(n 12 v 1 ) -	Gm 2 r 12	(n 12 v 12 )
	+ n + 2(n 12 S 1 ) (n -(n 12 S 1 ) Gm 2 r 12 (n 12 v 12 ) -v 1 3 2 G r 12 (6m 1 G r 12	(m 1 -m 2 ) .		(6.3)

  + m 2 ) + 3v 12 × S 1 (n 12 v 1 ) + n 12 × v 1 (v 1 S 1 ) + n 12 × v 12 3(n 12 S 1 ) (n 12 v 1 ) + (n 12 v 2 ) + (v 12 S 1 ) ,(7.8b)with the equations with 1 ↔ 2. The last terms in the RHS of Eqs. (7.7) are explicitly given by 12 S 1 )n 12 -S 1 + 1 2 v 2 1 S 1 -(v 1 S 1 )v 1 + 1 ↔ 2 , (7.9a) + 3(v 12 v 2 ) + 3(n 12 v 1 ) 2 -4(n 12 v 1 )(n 12 v 2 ) 12 S 1 ) -3(n 12 v 1 ) + (n 12 v 2 ) 12 S 1 )(n 12 v 2 ) -2(v 1 S 1 ) + Gm 2 v 2 r 12 -3(n 12 S 1 )(n 12 v 2 ) + 3(v 1 S 1 ) + 4(v 2 S 1 ) + 1 ↔ 2 . (7.9b)The 1.5PN term in the conserved angular momentum, Eq. (7.7a), agrees with the result ofKidder [29].13

	12 Gm 2 r 2 12 1.5PN = q S 1 = G m 2 n 12 × S 1 -n 12 × S 1 , r 12 2(n K S S 2.5PN = S 1 3 8 v 4 1 + Gm 2 r 12 1 + 5 2 v 2 1 + 4(v 1 v 2 ) -2v 2 2 + -1 2 v 2 3 2 (n 12 v 2 ) 2 + 1 2 Gm 1 r 12 + 3 2 Gm 2 r 12 + Gm 2 n 12 3 2 (n 12 v 2 ) 2 + r 12 (n 12 S 1 ) 3v 2 2G r 12 12 + (v 1 v 2 ) (7.8a) (m 1 K -3(n 12 v 1 ) 2 -4 Gm 1 r 12 -Gm 2 r 12 + (v + v 1 -1 2 (v 1 S 1 )v 2 1 + Gm 2 r 12 3(n

  . Thus, , n 12 , v 1 ) -2v 1 × S 1 + 3v 2 × S 1 + n 12 × v 1 + n 12 × v 2 (n 12 S 1 ) + 1 ↔ 2 .

					P S	1.5PN = p S	1 + p S	2 ,	(7.12a)
					P S	2.5PN = q S	1 + q S	2 ,	(7.12b)
	where the explicit expressions (7.8) hold. For G, we obtain rather simple expressions,
	G S	1.5PN = v 1 × S 1 + v 2 × S 2 ,			(7.13a)
	G S	2.5PN =	1 2	v 1 × S 1 v 2 1 -	Gm 2 r 12	-	y 1 r 12	(S 1 (7.13b)

As usual nPN refers to terms of order (v/c)

2n where v is the internal velocity and c the speed of light. In this paper we explicitely display all powers of c and of Newton's constant G.

Recall that with our convention the spin variable has the dimension of a true spin times c; the stress-energy tensor has the dimension of an energy density.

Thus ∂ ν ( √ -g g µν ) = 0, where g µν is the inverse of the usual covariant metric g µν , and g = det(g ρσ ).

This can also be deduced immediately from introducing a different spin variable S c 1 with constant magnitude (described in Sec. VII of paper II) and obeying dS c 1 /dt = Ω 1 × S c 1 ; noticing that the components of S c 1 in the basis {n, λ, ℓ} are linear combinations of those of S 1 , with constant coefficients.

However, let us stress that the definition of some individual momenta for the particles is merely introduced here as a convenient notation. In order to define in a meaningful way the notions of individual linear momenta of the particles with spins, we would need a Lagrangian, which as said before we did not compute, and the linear momenta would simply be the conjugate momenta of the ordinary positions.

This is the only place where one needs the precessional equations with 2PN accuracy.

A: LORENTZ INVARIANCE OF THE EQUATIONS OF MOTION

Because of the global Poincaré invariance of the Einstein equations (with bounded sources), and the manifest covariance of the De Donder harmonicity condition, it is not possible to physically distinguish between two harmonic-coordinate grids differing by a mere Lorentz transformation. As a result, the equations of motion must be of the same form in two such grids. In other words, up to an arbitrary PN order n, the link between the boosted acceleration a ′ 1 (y C , v C , a C ) and the boosted positions y ′ B (y C , v C ), velocities, v ′ B (y C , v C ) and spins S ′ 1 (y C , v C , S C ), must be given by the original equations of motion [i.e. Eq. ( 5.1) at the 2.5PN level] with the original variables being replaced by their primed counterparts. Note that the Euclidean metric and the totally antisymmetric tensors remain unchanged under Lorentz transformations. Schematically, we may write

The resulting relation between un-boosted quantities,

Equivalence with the equations of motion in the un-boosted frame:

up to negligible PN corrections. This property constitutes the so-called explicit Lorentz boost invariance of the equations of motion. It happens to be a very powerful check for the coefficients entering the functions A of Eq. (5.1), and in particular its contribution A 2.5PN [see Eq. (5.3b)].

In order to verify the validity of Eq. (A2), we need to determine the function A ′ explicitly, which requires to know how y B , v B , a 1 and S B transform under a Lorentz boost. Let us start with considering an arbitrary space-time event P with coordinates x µ in the current working frame (F ). Its coordinates in a boosted frame (F ′ ) of relative velocity V are related to the original ones by x ′µ = Λ µ ν (V)x ν , where the Lorentz matrix Λ µ ν (V) is given by

with γ being the Lorentz factor 1/ 1 -V 2 /c 2 . An event Q with coordinates y ′µ in (F ′ ) is simultaneous to P in the new frame if and only if y ′0 = x ′0 . There exist two such events located on the two world-lines of the binary companions. Their coordinates in (F ′ ) are denoted by y ′ µ 1 = (ct ′ , y ′ 1 ) and y ′ µ 2 = (ct ′ , y ′ 2 ) respectively. The mapping t ′ → y ′ 1 defines a function y ′ 1 (t ′ ), and similarly for the second body. The events having coordinates (ct ′ , y ′ 1 (t ′ )) and (ct ′ , y ′ 2 (t ′ )) in (F ′ ) do not generally appear as simultaneous in (F ). They may be referred to in components as (ct 1 , y 1 (t 1 )) and (ct 2 , y 2 (t 2 )) in that frame, the functions y 1 (t) and y 2 (t) being the original trajectories. By construction, we have

Let us express in the end RHS in terms of the coordinate time t. A derivation of the general formula linking y ′ 1 (t ′ ) to y 1 (t) in the PN scheme can be found in [57]. This relation reads, see Eqs. (3.20) in [57],

The velocity and acceleration follow from the partial derivation with respect to t ′ together with the formula

The spin components in the new frame cannot be obtained directly from the linear Lorentz transformation law. This is because the definition of S 1 and S 2 involves the inverse of the 3-metric γ ij induced by g µν on a slice t = const. Now, γ ij implicitly depends on the choice of the coordinate time and is generally singular because of the particle's self-gravitation. To avoid complications rising from this second issue, we shall first focus on the case of test particles on a fixed background.

In the frame (F ′ ), the spin components of the first test body read

Whereas the transformation law of γ ij 1 is more difficult, that of (γ ij ) 1 = (g ij ) 1 results straightforwardly from the transformation of the space-time metric:

with Λ µ λ (V) = Λ µ λ (-V) denoting the inverse transformation. Therefore, computing γ ′ij 1 (t ′ ) amounts to expressing the latter quantity as a function of (g ′ ij ) 1 . This is achieved by means of the relation det(γ ′ kl ) 1 (γ ′ij ) 1 = (Comγ ′ ) ji 1 , valid for any matrix (γ ′ ij ) 1 between its determinant det(γ ′ kl ) 1 , its comatrix (Comγ ′ ) ij 1 and its inverse. For 3-dimensional matrices, the determinant may be written in an Euclidean covariant form as

Similarly, we have for the comatrix

The inverse spatial metrics γ ′ij 1 is then given by the ratio of the RHS of Eqs. (A10) and (A9), where the primed metric relates to (g µν ) 1 after Eq. (A8).

We finally look at the determination of the covariant spin components S ′1 i . As S 1 µ is a Lorentzian vector, they are at once seen to be equal to

by virtue of the supplementary condition (2.8), S 1 0 = -S 1 i v i 1 /c. At this stage, we have expressed S ′ i 1 in terms of quantities evaluated at time t 1 , which has led us to a relation of the form S 1 = S 1 (t 1 ). It remains to rewrite S 1 (t 1 ) as a function of t. For the present purpose, we restrict ourselves to a perturbative approach, and resort to the convenient formula

generalizing in a straightforward way Eq. (3.16) of Ref. [57] to any smooth function f (see also the Appendix A of [57]). In the end, this yields the following identity for the spin "vector" S ′ 1 = (S ′i 1 ) defined in the frame (F ′ ):

where

These expressions are valid at any order in the boost velocity V. After specializing the above equation truncated at the PN level to the metric (3.1), we arrive at

Note that all powers of V consistent with the 1PN approximation beyond the leading spin-orbit term have been included. In principle, Eq. (A14a) holds only for test particles. Nonetheless, it turns out not to depend on any regularized field. It is thus legitimate to extend it to the conditions of the present problem.

With the previous transformation laws in hand, we are in position to check the Lorentz invariance as explained before. After a lengthy calculation, we arrive at the expected identity

APPENDIX B: TEST-MASS LIMIT OF THE EQUATIONS OF MOTION

In the limit where one of the objects, say the number 1, is nearly at rest while its companion has a very small mass for a finite ratio S 2 /m 2 , we must recover the dynamics of a spinning test particle in the background of a Kerr black hole of mass m 1 and spin S 1 = m 1 a (in this Appendix we pose G = c = 1). To allow direct comparison with the PN equations of motion for m 2 → 0 at S 2 /m 2 = const, we shall work with the Kerr metric in harmonic coordinates. The link between the Boyer-Lindquist grid (indicated by the label BL henceforth) and some spatial harmonic coordinates can be obtained from Eqs. ( 41) and (43) of Ref. [58]:

with r ± = m 1 ± m 2 1 -a 2 1 and i 2 = -1. Since ∇ µ ∇ µ t BL = 0, we may also choose t = t BL . The exact expression of the metric in the new grid is rather complicated, but we shall not need it beyond the linear order in the spin. Neglecting the quadratic terms O(S 2 1 ), the line element reduces to

which coincides with the one deriving from the metric (3.1) at the dominant order, hence the harmonic coordinates defined by Eqs. (B1) and t = t BL are the same as those of the PN formalism.

At this level, we may derive the equations of motion of a test particle with spin per unit mass S 2 /m 2 orbiting in the gravitational field (B2). For simplicity, we assume the trajectory to be circular and lie in the equatorial plane θ = π/2; the vector ∂ z points to the direction of the spin black hole, so that S 1 = S z 1 = m 1 a 1 ; the spherical coordinate basis is denoted by (∂ r , ∂ θ , ∂ φ ). The circularity conditions state in particular that r remains constant in time. The spatial components of the four-velocity are then

After taking these relations into account, the explicit form of the evolution equations (3.14) becomes

The harmonic gravitational field only depends on r and θ, both of which do not change with time. It is itself independent of t. Thus, Eqs. (B4a) and (B4d) imply that u 0 and ω = dφ/dt are constant, whereas (B4c) yields S 2 r = 0; (B4b) shows that S 2 θ = const and fixes the value of ω. We draw the time variation of the spin S 2 from the precession equation (2.12) specialized to the Kerr background (B2): Noticing that dS 2 r /dτ = 0, it is immediate to see from (B5a) together with the condition S 2 0 u 0 = -S 2 φ u φ that S 2 φ = 0. The remaining equations are identically satisfied. As a result, the spin of the small object is aligned (or anti-aligned) with the spin of the black hole, meaning that

up to possible quadratic contributions. In the test particle limit, the spin vectors are related to S and Σ as S