Testing general relativity with compact-body orbits: A modified Einstein-Infeld-Hoffmann framework
Clifford M. Will
gr-qc


Testing GR with compact-body orbits]Testing general relativity with compact-body orbits: A modified Einstein-Infeld-Hoffmann framework



cmw@phys.ufl.edu

^1 Department of Physics,
University of Florida, Gainesville FL 32611, USA
 

^2 GReCO, Institut d'Astrophysique de Paris, UMR 7095-CNRS,
Université Pierre et Marie Curie, 98^bis Bd. Arago, 75014 Paris, France




We describe a general framework for analyzing orbits of systems containing compact objects (neutron stars or black holes) in a class of Lagrangian-based alternative theories of gravity that also admit a global preferred reference frame.  The framework is based on a modified Einstein-Infeld-Hoffmann (EIH) formalism developed by Eardley and by Will, generalized to include the possibility of Lorentz-violating, preferred-frame effects.   It uses a post-Newtonian N-body Lagrangian with arbitrary parameters that depend on the theory of gravity and on “sensitivities” that encode the effects of the bodies' internal structure on their motion.    We determine the modified EIH parameters for the Einstein-Æther and Khronometric vector-tensor theories of gravity.
We find the effects of motion relative to a preferred universal frame on the orbital parameters of binary systems containing neutron stars, such as a class of ultra-circular pulsar-white dwarf binaries; the amplitudes of the effects depend upon “strong-field” preferred-frame parameters α̂_1 and α̂_2, which we relate to the fundamental modified EIH parameters.    We also determine the amplitude of the “Nordtvedt effect” in a triple system containing the pulsar J0337+1715 in terms of the modified EIH parameters.



Keywords: experimental gravity, general relativity, equations of motion, neutron stars, black holes, compact bodies, tests of general relativity

[
    Clifford M. Will^1,2
    December 30, 2023
========================





§ INTRODUCTION



The N-body equations of motion in the post-Newtonian limit of general relativity have been derived by numerous authors <cit.>.
A key assumption that
went into those analyses was that the weak-field, slow-motion limit of gravitational theory applied everywhere, in the interiors of the bodies as well
as between them.  This assumption restricted the applicability of the equations of motion to systems such as the solar system.
However, when dealing with systems
containing a neutron star or a black hole with highly relativistic spacetimes near or inside them, one can no longer
apply the assumptions of the post-Newtonian limit everywhere, except
possibly in the interbody region between the relativistic bodies.  Instead,
one must employ a method for deriving equations of motion for compact
objects that, within a chosen theory of gravity, involves solving the full,
relativistic equations for the regions inside and near each body, solving
the post-Newtonian equations in the interbody region, and matching
these solutions in an appropriate way in an “overlap region” surrounding
each body.  This matching  leads to constraints on the motions
of the bodies
which constitute the sought-after equations of motion.  In general relativity, a related procedure that avoided explicit inner solutions by exploiting surface integrals around each body  was carried out in the classic work of Einstein, Infeld and Hoffmann
(EIH) <cit.>.

In the fully weak-field post-Newtonian (PN) limit, it is also known that the
motion of self-gravitating, post-Newtonian fluid or dust bodies is independent of their internal structure, i.e., there is no Nordtvedt effect <cit.>.  Each body moves on a geodesic of
the post-Newtonian interbody metric generated by the other bodies, with
proper allowance for post-Newtonian terms contributed by its own interbody field.  The equations of motion obtained by this method are identical to the EIH equations.   This structure independence has also been verified explicitly to second post-Newtonian (2PN) order for self-gravitating fluid bodies <cit.>. It is therefore reasonable to assume that the EIH equations are valid for systems of compact bodies (neutron stars or
black holes).  The only restriction is that they be quasistatic, nearly spherical, and sufficiently small compared to their separations that tidal interactions may be neglected (throughout this discussion, we will ignore the effects of spin).  
This was verified in a restricted context for nonrotating black holes in a seminal paper 
by d'Eath <cit.>, and was subsequently extended to a variety of contexts  <cit.>.   

Key to this result is the validity of the Strong Equivalence
Principle within general relativity, which
guarantees that the structure of each body is independent of the surrounding gravitational environment.  By contrast, most alternative theories of gravity possess additional gravitational fields, whose values in the
matching region can influence the structure of each body, and thereby can
affect its motion.    Using a cyclic gedanken experiment within a freely falling frame that encompasses a given body, and that assumes only conservation of energy (see Sec. 2.5 of <cit.> for discussion), it can be shown that, if the
binding energy E_ B of a body varies as a result of a variation in an external parameter or field ψ, the body experiences an additional acceleration given by
δa∼m^-1∇ E_ B (x, v) ∼ (∂ln m/∂ψ) ∇ψ.
Thus, the bodies need not follow geodesics of any metric, but instead their
motion may depend on their internal structure.

In practice, the EIH-inspired matching procedure  is very cumbersome 
<cit.>.  Within general relativity, a simpler method for
obtaining the EIH equations of motion is to treat each body as a “point” mass of inertial mass m_a and to solve Einstein's equations using a point-mass matter action or energy-momentum tensor, with proper care to neglect or regularize infinite “self” fields (at high PN orders, this regularization can become very complicated <cit.>).   In the action for general relativity, we thus write

    I = 1/16 π G∫ R √(-g)   d^4x - ∑_a m_a ∫ dτ_a  ,

where τ_a is proper time along the world line of the ath body.  By solving the field equations to 
1PN order, it is then possible to derive straightforwardly from the matter action an N-body EIH action in the form

    I_ EIH = ∫ L(_1, …_N, _1, …_N ) dt  ,

with a Lagrangian L
written purely in terms of the variables (_a, _a) of the bodies.  The result is the 1PN N-body Lagrangian of general relativity.  The N-body EIH equations of motion are then given by 

    d/dt∂ L/∂ v_a^j -  ∂ L/∂ x_a^j  = 0 ,    a = 1, … N  .


In alternative theories of gravity, we assume that the only difference is the possible
dependence of the mass on the boundary values of the auxiliary fields.  Thus, following the suggestion of Eardley <cit.>, we merely replace the
constant inertial mass m_a in the matter action with the variable inertial
mass m_a (ψ_A), where ψ_A represents the values of the external auxiliary
fields, evaluated at the  body (we neglect their variation across the interior of the matching region), with infinite self-field contributions
excluded.  
The functional dependence of m_a upon the variable ψ_A will
depend on the nature and structure of the body.  Thus, we write the action of the alternative theory in the form

    I = I_G - ∑_a ∫ m_a  (ψ_A  ) dτ_a  ,

where I_G is the action for the metric and auxiliary fields ψ_A.  
In varying the action with respect to the fields g_μν and ψ_A, the variation of
m_a must now be taken into account.   In the post-Newtonian limit, where the fields ψ_A are expanded about asymptotic values ψ_A^(0) according to
ψ_A = ψ_A^(0) + δψ_A, it is generally sufficient to expand m_a (ψ_A) in the form

    m_a (ψ_A) = m_a (ψ_A^(0) ) + ∑_A ∂ m_a/∂ψ_A^(0)δψ_A + 1/2∑_A,B∂^2 m_a/∂ψ_A^(0)∂ψ_B^(0)δψ_Aδψ_B + … .

Thus, the final form of the metric and of the N-body Lagrangian will
depend on m_a ≡ m_a (ψ_A^(0) ) and on the parameters ∂ m_a /∂ψ_A^(0), and so on. 
It is conventional to use
the term “sensitivities” to describe these parameters, since they measure the
sensitivity of the inertial mass to changes in the fields ψ_A. 
Thus, we 
define

    s_a^(A)   ≡∂ln m_a/∂lnψ_A^(0) ,   s'_a^(AB)≡∂^2 ln m_a/∂lnψ_A^(0)∂lnψ_B^(0) ,

and so on, where the derivatives are typically taken holding the total baryon number of the body fixed, for bodies made of matter.   For black holes, some other method must be used to identify the sensitivities.     Similar sensitivities can be defined for the radius and moment of inertia of the compact body. 
Gralla <cit.> developed a more general theory of the motion of “small” bodies characterized by such parameters as mass, spin and charge in an environment of external fields, and
 argued that Eardley's ansatz is a special case of that general framework.

These ideas were incorporated into a “modified EIH formalism” in Sec. 11.3 of Theory and Experiment in Gravitational Physics <cit.> (based in part on unpublished notes by Eardley).  However, that formalism was restricted by requiring that the N-body Lagrangian be invariant under suitable low-velocity Lorentz transformations.   This precluded the possibility of preferred-frame effects in compact-body dynamics.   A number of developments since that time (1981) have made it desirable to relax that assumption and to create a more general modified EIH formalism.  The first was the formulation of a variety of alternative theories of gravity with auxiliary vector and tensor fields, whose presence establishes a preferred frame of reference where the components of those fields take on privileged values.  Examples include vector-tensor theories such as Einstein-Æther theory and Khronometric theory <cit.>, and scalar-vector-tensor theories, such as TeVeS and STV <cit.>.   A second development was a number of very precise tests of preferred-frame effects in the orbits of ultralow-eccentricity binary pulsar systems <cit.>.   These analyses quoted bounds on preferred-frame parameters α̂_1 and α̂_2, where the hats are meant to denote some unspecified generalization of the usual parametrized post-Newtonian (PPN) preferred-frame parameters α_1 and α_2 to strong-field situations.   One goal of this paper is to relate these “hatted” parameters explicitly to parameters of the modified EIH formalism.  In the special cases of Einstein-Æther and Khronometric theories, we will obtain the modified EIH parameters explicitly, making it possible to express α̂_1 and α̂_2 in terms of the fundamental parameters of these theories and of the sensitivities of the compact bodies in the system.  We will also express the “Nordtvedt” parameter describing a failure of the universality of free fall in a compact-body system such as the pulsar J0337+1715 in a triple system with two white dwarf companions <cit.> in terms of parameters of this modified EIH formalism.

In Sec. <ref>, we describe this generalization of the modified EIH formalism, and in Sec. <ref>, we obtain the parameters of the formalism in a selection of alternative theories.  Section <ref> applies the formalism to two-body and three-body dynamics, obtaining a number of observable effects in terms of the modified EIH parameters.   In Sec. <ref> we make concluding remarks.   We use units in which the locally-measured gravitational constant G and the speed of light c are unity.  



§ MODIFIED EINSTEIN-INFELD-HOFFMANN FORMALISM



We construct a general EIH formalism  using arbitrary parameters whose values depend both on the theory under study and on the nature of the bodies in the system.   In this case, however, the parameters appear in the N-body 1PN Lagrangian rather than in the metric.
This implies that we are restricting attention to Lagrangian-based metric theories of
gravity, known as semiconservative theories (see Sec. 4.4 of <cit.> or <cit.>  for discussion).
We will also restrict attention to theories that have no
Whitehead term in the post-Newtonian limit; all of the currently popular theories as described in <cit.> satisfy this constraint.   Our goal is to generalize the formalism presented in <cit.>, which made a restriction to fully conservative theories.  The distinction between the two classes of theories is that semiconservative theories admit preferred-frame effects at 1PN order, i.e. one or more of their PPN parameters α_1 and  α_2 is nonzero, while fully conservative theories have α_1=α_2=0 and do not have preferred-frame effects at 1PN order.    Nordtvedt <cit.> developed a similar formalism for compact-body dynamics, but did not address the consequences of preferred-frame effects.

Each body is characterized by an inertial mass m_a, defined to be the
quantity that appears in the conservation laws for energy and momentum
that emerge from the EIH Lagrangian.  We then write for the metric,
valid in the interbody region and far from the system,

    g_00   = -1 + 2 ∑_a α^*_a m_a/| - _a| + O(ϵ^2)  ,
    
    
    g_0j    = O(ϵ^3/2)  ,
    
    
    g_jk   =  ( 1 + 2 ∑_a γ^*_a m_a/| - _a| ) δ_jk + O(ϵ^2)  ,

where α^*_a  and γ^*_a are functions of the parameters of the theory and of the
structure of the ath body, and ϵ∼ m/r ∼ v^2.   For test-body geodesics in this metric, the quantities α^*_a m_a, and ∑_a α^*_a m_a are the Kepler-measured active gravitational
masses of the individual bodies and of the system as a whole. The metric (<ref>) is used mainly to discuss the propagation of photons in systems with compact bodies, and plays a role, for example,  in obtaining timing formulae for binary pulsars.  In general
relativity, α^*_a ≡γ^*_a ≡ 1.


To obtain the modified EIH Lagrangian, we first generalize the post-Newtonian
semiconservative N-body Lagrangian, given by Eq. (6.80) of <cit.> or Eq. (6.81) of <cit.>:

    L     = - ∑_a m_a  ( 1 - 1/2 v_a^2 - 1/8 v_a^4  )
    
          + 1/2∑_a ∑_b  am_a m_b/r_ab [ 1 + (2γ + 1) v_a^2
    -(2β-1) ∑_c  am_c/r_ac
            
    - 1/2 (4γ + 3 +α_1 - α_2) _a ·_b 
    - 1/2 (1 + α_2 ) (_a ·n_ab)(_b ·n_ab)
    
           
    - ξ_ab/r_ab^2·∑_c  ab m_c  ( _bc/r_ac
    -  _ac/r_bc )
      ]  ,

where r_ab = |_a - _b| and n_ab = _ab/r_ab.   
The Lagrangian is expressed in a coordinate system at rest with respect to the preferred-frame singled out by the theory of gravity in question; this is generally assumed to coincide with the frame in which the cosmic background radiation is isotropic.  

We set the Whitehead parameter ξ = 0, and replace PPN parameters with parameters dependent upon each body, according to

    L_ EIH   = - ∑_a m_a  [ 1 - 1/2  v_a^2 - 1/8 (1+ A_a  ) v_a^4  ]
    
           +1/2∑_a  ∑_b  am_a m_b/r_ab [  G_ab + 3  B_ab v_a^2 
    - 1/2 ( G_ab + 6  B_(ab) +  C_ab ) _a ·_b 
     .
    
        .
         
    - 1/2 ( G_ab +  E_ab ) (_a ·n_ab)(_b ·n_ab)  ] 
    
          
    - 1/2∑_a ∑_b  a∑_c  a D_abcm_a m_b/r_abm_c/r_ac .

The quantities
A_a, G_ab, B_ab, C_ab, E_ab and D_abc
are functions of the parameters of the theory and of the structure of each
body, and satisfy

    G_ab =  G_(ab)  ,    C_ab =  C_(ab) ,    E_ab =  E_(ab) ,    D_abc =  D_a(bc) .

Note that B_ab has no special symmetry, in general.

Notice that we did not introduce a parameter in front of the kinetic v_a^2 term in Eq. (<ref>).  Any such parameter can always be absorbed into a new definition of the inertial mass m'_a of body a.  We are then free to change the constant term -∑_a m_a to be the sum of the new inertial masses -∑_a m'_a.  This has no effect on the equations of motion, but does allow the Hamiltonian derived from L_ EIH to be the sum of the new inertial masses at lowest order.   We also did not include a term of the form
(m_a m_b/r_ab)(_b ·n_ab)^2; such a term can be associated (via a total time derivative in the Lagrangian) with the Whitehead term in Eq. (<ref>), which we have chosen to reject.   

In general relativity, G_ab=  B_ab =  D_abc = 1, while A_a =  C_ab = E_ab = 0.  
In the  post-Newtonian limit of semiconservative theories (with ξ=0), for structureless
masses (no self-gravity), the parameters have the values [compare Eq. (<ref>)]

    G_ab   =1  , 
       B_ab = 1/3 (2γ +1)  ,    D_abc = 2β -1  ,
    
     A_a    = 0  ,    C_ab =  α_1 - α_2  ,    E_ab =  α_2 
      .

In the fully conservative case, including contributions of the 1PN-order self-gravitational binding energies of the bodies, the parameters can be shown to have the values

    G_ab   = 1 + (4β-γ-3)  ( Ω_a/m_a +  Ω_b/m_b )  ,
    
     B_ab   = 1/3 (2γ +1)  ,    D_abc = 2β -1  ,
    
     A_a      =  C_ab =   E_ab = 0 
      ,

where Ω_a is the self-gravitational energy of the ath body.  The connection between these parameters and those introduced by Nordtvedt <cit.> is detailed in an Appendix.

To obtain the Lagrangian in a moving frame, we make a Lorentz transformation from the original preferred
frame to a new frame which moves at velocity w relative to the old frame.
In order to preserve the post-Newtonian character of the Lagrangian, we
assume that w ≡ |w| is small, i.e. of O(ϵ^1/2). This transformation from 
rest coordinates x^α=(t, x) to moving coordinates ξ^μ = (τ, ξ)  can be expanded in powers of
w to the required order.  This approximate form of the Lorentz transformation is sometimes called a post-Galilean transformation <cit.>, and has the form

    = ξ +  ( 1+ 1/2 w^2  ) wτ + 1/2 (ξ·w ) w + ξ× O(ϵ^2)  ,
    
    
    t    = τ ( 1 + 1/2 w^2 + 3/8 w^4  ) +  ( 1+ 1/2 w^2  ) ξ·w + τ× O(ϵ^3)  ,

where wτ is assumed to be O(ϵ^0).  
From the transformation (<ref>), we have

    _a    = ν_a + w - dν_a/dτ (ξ_a ·w) - 1/2 w^2 ν_a - (ν_a ·w)  (ν_a + 1/2w )  ,
    
    1/r_ab   = 1/ξ_ab [ 1 + 1/2 (w·n'_ab )^2 
    + 1/ξ_ab  (w·ξ_a )(ν_a ·n'_ab )
     .
    
            .
    -  1/ξ_ab(w·ξ_b )(ν_b ·n'_ab )  ]  ,

where ξ_a, ξ_b and  ν_a ≡ dξ_a/dτ are to be evaluated at the same time τ, given by a clock at the spatial origin (ξ=0) of the moving coordinate system, and n'_ab≡ξ_ab/ξ_ab.  
The new Lagrangian is given by

    L( ξ, τ) = L(, t) dt/dτ - df/dτ ,

where dt/dτ is evaluated at ξ=0, and where we are free to subtract a total time derivative of a function f to simplify the new Lagrangian.
Substituting these results into Eqs. (<ref>) and (<ref>), 
dropping constants and total time derivatives, and replacing ξ_ab and ν_a with r_ab and _a, we obtain the Lagrangian in the moving frame

    L   =  L_ EIH  + 1/4∑_a m_a  A_a  [ v_a^2 w^2 + 2v_a^2 (v_a ·w) + 2(v_a ·w)^2  ]
    
          - 1/4∑_a  bm_a m_b/r_ab{ C_ab w^2 - 2 ( 6  B_[ab] -  C_ab ) (v_a ·w)
      .
     
           . 
     +  E_ab [ (w·n_ab )^2 + 2 (w·n_ab )(v_a ·n_ab )  ]
     } .

Notice that the Lagrangian is post-Galilean invariant if and only if

    A_a  ≡ B_[ab]≡ C_ab≡ E_ab≡ 0  .

These quantities are then the preferred-frame parameters of our modified EIH formalism.



§ MODIFIED EIH PARAMETERS IN ALTERNATIVE THEORIES



As an illustration of this modified EIH framework for compact bodies, we will focus on specific theories where calculations have been carried out.   As we have already discussed, a variety of approaches have shown that  the EIH equations of motion for compact objects within general relativity are identical to those of the post-Newtonian limit with weak fields everywhere.  In other words, in general relativity, G_ab≡ B_ab≡ D_abc≡ 1, and the remaining coefficients vanish, independently of
the nature of the bodies.



 §.§ Scalar-tensor theories



The modified EIH formalism was first developed by Eardley <cit.> for application to the Jordan-Fierz-Brans-Dicke theory.  It makes use of the fact that
only the scalar fleld ϕ produces an external influence on the structure of
each compact body via its boundary values in the matching region. This boundary value of ϕ is related to the local value of the gravitational
constant as felt by the compact body by

    G_ local = G/ϕ (4 + 2ω/3+2ω )  ,

where G is the fundamental gravitational coupling constant.   Thus we will treat the inertial mass m_a of each body as a being a function of ϕ.   Then, by defining the deviation of ϕ from its asymptotic value ϕ_0 by ϕ≡ϕ_0 (1+ Ψ), we can write down the expansion 

    m_a(ϕ) =
     m_a [1  +  s_a Ψ + 1/2  (s_a^2 +s'_a   -s_a) Ψ^2 + O(Ψ^3)   ] ,
  
where m_a ≡ m_a (ϕ_0), and  we define the dimensionless sensitivities

    s_a    ≡ ( d ln m_a(ϕ)/d lnϕ )_0  ,
    
    
    s'_a    ≡ ( d^2 ln m_a(ϕ)/d (lnϕ)^2 )_0  .

The action for massless scalar-tensor theory is then written

    I = 1/16π G∫ [ ϕ R - ω(ϕ)/ϕ g^μνϕ_,μϕ_,ν ] √(-g) d^4x  - ∑_a ∫ m_a (ϕ) dτ_a  ,

where the integrals over proper time  τ_a are to be taken along the world line of each body a.  It is straightforward to vary the action with respect to g_μν and ϕ to obtain the field equations,

    G_μν   = 8π G/ϕ T_μν + ω(ϕ)/ϕ^2 ( ϕ_,μϕ_,ν - 1/2 g_μνϕ_,λϕ^,λ ) 
    
         + 1/ϕ ( ϕ_;μν - g_μν_g ϕ )  , 
    
    _g ϕ   = 1/3 + 2ω(ϕ) ( 8π G T  - 16π G ϕ∂ T/∂ϕ - dω/dϕϕ_,λϕ^,λ )  ,



where 

    T^μν = (-g)^-1/2∑_a m_a (ϕ) u_a^μ u_a^ν (u_a^0)^-1δ^3 ( x -  x_a)   ,

where u_a^μ is the four-velocity of body a.
The equations of motion take the form

    T^μν_;ν - ∂ T/∂ϕϕ^,ν = 0  .


Carrying out a post-Newtonian calculation of the metric as described, for example in <cit.> or in Sec. 5.3 of <cit.>, we obtain, to lowest order 

    Ψ    =2 ζ∑_b m_b/r_b (1 - 2s_b  ) + O(ϵ^2 )  ,
    
    
    g_00    = -1 + 2 ∑_b m_b/r_b ( 1 - 2 ζ s_b  ) + O(ϵ^2 )  ,
    
    
    g_0j   = - 4 (1-ζ) ∑_b m_b v_b^j/r_b + O(ϵ^5/2)  ,
    
    
    g_jk   = δ_jk [ 1 + 2 ∑_b m_b/r_b ( 1 - 2 ζ +2 ζ s_b  )  ]  ,

where ζ = 1/(4 + 2ω_0), r_b = |-_b|, and we have chosen units in which G_ local =1.  For the explicit O(ϵ^2) terms in g_00 and Ψ, see <cit.>.  Notice that the active gravitational mass as measured by test-body Keplerian orbits far from each body is given by 

    (m_ A)_a = m_a (1- 2 ζ s_a )  .


From the complete post-Newtonian solution for g_μν and Ψ, we can obtain the matter  action for the ath body, given by 

    I_a = - ∫ m_a(ϕ)  ( - g_00 - 2 g_0j v_a^j - g_jk v_a^j v_a^k  )^1/2 dt  .

To obtain an N-body action in the form of Eq. (<ref>), we first make
the gravitational terms in I_a, manifestly symmetric under interchange of
all pairs of particles, then take one of each such term generated in I_a, and
sum the result over a.  The resulting N-body Lagrangian then has the form of Eq. (<ref>) with

    G_ab   = 1 - 2ζ (s_a +s_b  - 2s_a s_b  )  ,
    
     B_ab   = 1/3 [  G_ab + 2(1-ζ)  ]  ,
    
     D_abc   =  G_ab G_ac + 2ζ (1-2s_b)(1-2s_c)  [ λ (1-2s_a) + 2 ζ s'_a  ]  ,
    
     A_a    =  B_[ab] =   C_ab =  E_ab = 0  .




 §.§ Vector-tensor theories



Because the norm of the vector field K in Einstein-Æther and Khronometric theories is constrained to be -1,
the structure of a spherically symmetric compact body at rest with rest to the preferred rest frame does not depend on it.  However it could depend on the time component K^0, or more properly on the invariant quantity K· u, where u is the body's four-velocity.  (The structure of a  rotating body could also depend on the projection of K along the body's spin axis, but here we will focus on non-rotating bodies.)   We define for a body with four-velocity u, 

    γ≡ - K^μ  u_μ≡ 1 + Ψ ,

where we assume that far from the system, for a test body at rest, K· u = -1.    We define the sensitivities 

    s_a    ≡ ( d ln m_a (γ)/d lnγ )_γ =1 ,
    
    
    s'_a    ≡ ( d^2 ln m_a (γ)/d (lnγ)^2  )_γ =1 ,
  
where the derivatives are to be taken holding baryon number fixed.  Then the expansion of m_a(γ) in powers of Ψ is again given by Eq. (<ref>).
With this assumption, Foster <cit.> and Yagi et al. <cit.>  derived the metric and equations of motion to post-Newtonian order for systems of compact bodies in Einstein-Æther theory.  In the preferred rest frame, the metric is given by

    g_00   = -1 + 2U - 2U^2 - 2Φ_2 + 3Φ_1s + O(ϵ^3)  ,
    
    
    g_0j   = g^j + O(ϵ^5/2)  ,
    
    
    g_jk   = (1+2U) δ_jk + O(ϵ^2)  ,

and the vector field is given by

    K^0    = (-g_00)^-1/2 + O(ϵ^3)  ,
    
    
    K^j    = k^j + O(ϵ^5/2)  ,

where

    U     = ∑_b G_N m_b/r_b ,   Φ_2 = ∑_b,cG_N^2 m_b m_c/r_b r_bc ,   Φ_1s = ∑_b G_N m_b/r_b v_b^2 (1-s_b)  ,
    
    
    g^j    = ∑_b  G_N m_b/r_b [ B_b^- v_b^j + B_b^+ n_b^j (n_b ·_b)  ]  ,
    
    
    k^j    = ∑_b  G_N m_b/r_b [ C_b^- v_b^j + C_b^+ n_b^j (n_b ·_b)  ]  ,

where G_N = 2G/(2-c_14), G is the gravitational coupling constant of the theory,  n_b = (x - x_b)/r_b, and c_14 = c_1 +c_4.   
The quantities B_b^± and C_b^± are complicated expressions involving the constants c_1, c_2, c_3 and c_4 of Einstein-Æther theory and the sensitivities s_b (see Eqs. (23) and (29) of <cit.>).   

From Eqs. (<ref>) and (<ref>) we obtain

    Ψ (_a)  = 1/2 v_a^2 + 3/8 v_a^4 + 2 v_a^2 U(_a) - v_a^j k^j(_a)  + O(ϵ^3)  .

Note that, when v_a^j =0, Ψ = 0, resulting in no dependence of the inertial mass on the vector field.  
Writing the action for the ath body as 

    I_a = - ∫ m_a(Ψ)  ( - g_00 - 2 g_0j v_a^j - g_jk v_a^j v_a^k  )^1/2 dt  ,

we expand to post-Newtonian order, make the action manifestly symmetric under interchange of all pairs of particles, select one of each term, and sum over a.   After rescaling each mass by m_a → m_a/(1-s_a) and replacing the constant term in the Lagrangian by the sum of the rescaled masses, we obtain the modified EIH Lagrangian in the form  of Eq. (<ref>), with

    G_ab   = G_N/(1-s_a)(1-s_b) , 
     
     B_ab   =  G_ab (1-s_a)  ,
    
     D_abc   = G_ab G_ac (1-s_a)  ,
    
     A_a    =  s_a - s'_a/1-s_a ,
    
     C_ab   =   G_ab [  α_1 - α_2 + 3  (s_a + s_b  )  -  Q_ab - R_ab ]  ,
    
     E_ab   =  G_ab [  α_2 +  Q_ab -  R_ab ]  ,

where α_1 and α_2 are the PPN preferred-frame parameters of Einstein-Æther theory, given by 


    α_1     =  -8 (c_3^2 + c_1 c_4)/2 c_1 - c_1^2 + c_3^2 ,
    
    α_2     =  1/2α_1 - (2c_+-c_14)(c_++c_14+3c_2)/c_123 (2-c_14) ,

and 

    Q_ab   = 1/2 (2-c_14/2c_+ - c_14 ) (α_1 - 2α_2) (s_a +s_b) +  2 - c_14/c_123 s_a s_b   ,
    
     R_ab   = 8+α_1/4c_1 [ c_- (s_a + s_b) + (1-c_-) s_as_b  ] .

Here 
 c_± = c_1 ± c_3, and c_123 = c_1 + c_2 +c_3.   The two-body equations of motion that follow from the Lagrangian with these coefficients agree with Eq. (33) of <cit.> (after correcting a sign and a parenthesis in Eqs. (34) and (35) of that paper).

Note that the Lagrangian is in general not Lorentz invariant, and therefore will exhibit preferred-frame effects.  Even when the parameters c_i are constrained so as to enforce α_1 = α_2 = 0, making the dynamics Lorentz invariant in the fully weak field post-Newtonian limit, the dynamics of compact bodies can still be dependent on the overall motion of the system via the motion-induced sensitivities of the bodies.

The corresponding equations in Khronometric theory can be obtained from these by setting c_1 = - ϵ, c_2 = λ_K, c_3 = β_K + ϵ, and c_4 = α_K + ϵ, and taking the limit ϵ→∞ <cit.>.  The parameters of the modified EIH Lagrangian are given by Eqs. (<ref>), but now with G_N = 2G/(2-α_K), and 

    Q_ab   =  1/β_K + λ_K [(α_K + β_K + 3 λ_K) (s_a +s_b) + (2-α_K)  s_a s_b  ]  ,
    
     R_ab   = 1/2 (8+α_1)  [ s_a + s_b - s_as_b  ] ,

where α_1 is the PPN parameter of Khronometric theory, given by 

    α_1 = 4(α_K - 2 β_K)/β_K -1 .


Yagi et al. <cit.>  calculated neutron star sensitivities in both Einstein-Æther  and Khronometric theories.  In order to do so, it was necessary to construct models for neutron stars moving uniformly relative to the preferred frame.  From
Eqs. (<ref>) and (<ref>), assuming uniform motion with no external bodies (U = k^j =0), the sensitivity is given by s =  v^-1 d ln m/dv.   They chose the coefficients c_2 and c_4 so that the PPN parameters α_1 and α_2 saturate the bounds from solar-system measurements, and obtained fitting formulae for the sensitivities as a function of c_+, c_- and the compactness M/R of the neutron star.



§ APPLICATION TO TWO- AND THREE-BODY DYNAMICS



Since one of our goals is to apply this formalism to binary systems
containing compact objects, let us now restrict
attention to two-body systems.  
We obtain from L the two-body equations of motion

    a_1    =  -m_2 n/r^2{ G_12 
    -  (3 G_12 B_12 +  D_122 ) m_2/r .
    
        .
        
    - 1/2 [ 2 G_12^2 + 6 G_12 B_(12) + 2 D_211 +  G_12( C_12 +  E_12)  ] m_1/r .
    
        .
        
    + 1/2 [3 B_12 -  G_12 (1+  A_1)  ] v_1^2
    + 1/2(3 B_21 +  G_12+  E_12) v_2^2 
     .
    
        .
        
    - 1/2 ( 6 B_(12) + 2  G_12 +  C_12 +  E_12 ) v_1 ·v_2  
    - 3/2 ( G_12+  E_12 ) (n·_2)^2 
     .
    
        .
        
    + 1/2 (  C_12 +  G_12 A_1  ) w^2
    + 1/2 (  C_12 - 6 B_[12]+  E_12 + 2  G_12 A_1  ) _1 ·w .
    
        .
        
    + 1/2 (  C_12 + 6 B_[12]-  E_12 ) _2 ·w .
    
        .
        
    + 3/2 E_12 [ (w·n )^2 + 2(w·n )(v_2 ·n )  ]
    }
         
    +  m_2 v_1/r^2n·{ [3 B_12 +  G_12 (1+  A_1)  ] _1 - 3 B_12_2  +  G_12 A_1 w}
         
    - 1/2m_2 v_2/r^2n·{ ( 6 B_(12) + 2  G_12 +  C_12 +  E_12 ) _1
     .
    
        .
        
    - ( 6 B_(12)  +  C_12 -  E_12 ) _2
    +2  E_12w}
         
    - 1/2m_2 w/r^2n·{ (   C_12 -6 B_[12] +  E_12 -2 G_12 A_1  ) _1
     .
    
        .
        
    -  ( C_12 - 6 B_[12] -  E_12 ) _2
    -2  (  G_12 A_1 -  E_12 ) w} ,
    
    a_2    = { 1⇌ 2;  n→ -n} ,

where a_a ≡ d_a/dt, ≡_1 - _2, r ≡ ||, and n≡/r.  If the conditions of Eq. (<ref>) hold, then the preferred-frame terms vanish.

To Newtonian order, the center of mass X = (m_1 _1 + m_2 _2)/m of the system is unaccelerated, thus, to sufficient accuracy in the post-Newtonian terms, we can set X = Ẋ = 0 and 
write

    _1    =  [ m_2/m + O(ϵ)  ]  ,
    
    _2    = - [ m_1/m + O(ϵ)  ]  .

We define

    ≡_1 - _2  ,   a≡a_1 - a_2  ,
    
    
    m    ≡ m_1 +m_2  ,   η≡m_1m_2/m^2 ,   Δ = m_2 - m_1/m .

We also define the two-body coefficients

    G   ≡ G_12 ,    B_+≡ B_(12) ,    B_-≡ B_[12] ,    D≡m_2/m D_122 + m_1/m D_211 , 
    
     C   ≡ C_12 ,    E≡ E_12 ,    A^(n)≡ ( m_2/m )^n  A_1 -  (- m_1/m )^n  A_2   .

Then the equation of motion for the relative orbit takes the form 

    a = a_ L + a_ PF ,

where the purely two-body, or “local” contributions have the form (we use “hats” to denote parameters associated with compact bodies)

    a_ L =  - G m n/r^2 
    + m/r^2 [ n ( Â_1 v^2 + Â_2 ^2 + Â_3 m/r ) 
    + B̂_1  ]  ,

where

    Â_1    = 1/2{ G ( 1- 6η ) - 3 B_+ - 3 Δ B_- - η ( C + 2 E) +   G A^(3)} ,
    
    Â_2    = 3/2η ( G +  E)  ,
    
    Â_3    =  G [ 2 η G + 3 B_+ + η (  C +  E ) + 3 Δ B_- ] +  D ,
    
    B̂_1    =  G  (1-2η)  + 3 B_+ + 3 Δ B_- + η C +   G A^(3) ,

and the preferred-frame contributions have the form

    a_ PF   =  
     m/r^2{n [  ( 1/2α̂_1 + 2 G A^(2) ) (w· ) + 3/2 (α̂_2+  G A^(1) ) (w·n)^2  ]
     .
    
          . 
    - w [ 1/2α̂_1 (n· )
    + α̂_2 (n·w)  ]
    +  G A^(2) (n·w)
    } ,

where 

    α̂_1     = Δ (  C +  E)
      - 6  B_-  - 2  G A^(2) ,
      
    α̂_2     =  E -  G A^(1) .

These two parameters play the role of compact-body analogues of the PPN parameters α_1 and α_2 in the context of binary systems.

In the Newtonian limit of the orbital motion, we have a =  Gm/r^3, with Keplerian orbit solutions given by x= r  n, r=p/(1+e cos f), where

    n   ≡ ( cosΩcosϕ- cosιsinΩsinϕ )  e_X 
    
         
     +  ( sinΩcosϕ + cosιcosΩsinϕ ) e_Y
    + sinιsinϕ e_Z  ,

where the orbit elements are: semilatus rectum p, eccentricity e, inclination ι, longitude of ascending node Ω and pericenter angle ω; 
f ≡ϕ - ω is the  true anomaly, ϕ is the orbital phase measured from the ascending node and 
e_A are chosen reference basis vectors.   We also have that  v = ṙ n + (h/r) λ and ṙ = (he/p) sin f, where λ = ∂n/∂ϕ.  The orbital angular momentum per unit reduced mass is given by h = x×v = ( Gmp)^1/2ĥ, where ĥ = n×λ.   The semilatus rectum p is related to the semimajor axis a by p = a(1-e^2).  

In the presence of perturbations, these equations define the “osculating” orbit, with variable orbital elements.   
Folowing the standard methods for the perturbed Kepler problem as described, for example, in Sec. 3.3 of <cit.>, we obtain the secular changes in the orbit elements. 
From the post-Newtonian terms in Eq. (<ref>), we find that the pericenter advance per orbit is given by

    Δω = 6 π m/p P G^-1 ,

where 

    P   =  G B_+ + 1/6 ( G^2 -    D )
    
    
    + 1/6 G [  6  Δ B_-  +η G (2 C +  E )  +  G A^(3) ] .

 This is the only secular perturbation produced by the post-Newtonian terms in Eq. (<ref>).  
 
We now calculate the secular changes in the orbit elements resulting from the preferred-frame perturbations in Eq. (<ref>).  We define

    e_ P   ≡n|_ϕ = ω = e_Ωcosω + e_⊥sinω ,
    
    e_ Q   ≡λ|_ϕ = ω = - e_Ωsinω + e_⊥cosω ,
    
    ĥ   ≡e_ P×e_ Q =  e_Ω×e_⊥ ,

where e_ P is a unit vector pointing toward the pericenter and e_ Q=ĥ×e_ P; e_Ω is a unit vector pointing along the ascending node, and e_⊥ = ĥ×e_Ω.   For any vector A, we then define components A_ P, A_ Q, A_ h, A_Ω and A_⊥ accordingly.  

The secular changes in the orbit elements are then given by

    Δ a    = 0  ,
    
    Δ e    = -πα̂_1 Δ (m/p )^1/2 w_ P (1-e^2) F(e) 
    +2πα̂_2  w_ Pw_ Q e√(1-e^2) F(e)^2 
      ,
    
    Δϖ   = - πα̂_1 Δ (m/p )^1/2 w_ Q√(1-e^2) F(e)/e 
    -πα̂_2   (w_ P^2 - w_ Q^2  )  F(e)^2 
      ,
    
    Δι   =  πα̂_1 Δ (m/p )^1/2 w_ hsin(ω) e F(e)
    -2πα̂_2 w_ h w_ RF(e)/√(1-e^2) ,
    
    ΔΩ   = - πα̂_1 Δ (m/p )^1/2w_ h/sinιcos(ω) e F(e)
    -2πα̂_2 w_ h w_ S/sinιF(e)/√(1-e^2) ,

where Δϖ = Δω + cosι ΔΩ, 
F(e) ≡ (1+ √(1-e^2))^-1,
and for any vector A,

    A_ R   ≡ A_ Pcos(ω) - A_ Qsin(ω) √(1-e^2) ,
    
    
    A_ S   ≡ A_ Psin(ω) + A_ Qcos(ω) √(1-e^2) .

Note that the secular perturbations do not depend on the coefficients G A^(n).

It is tempting to interpret these secular changes as implying linearly growing values of the orbital elements.  However, the expressions in Eqs. (<ref>) depend on ω, both from explicit ω dependence, and via the P, Q, R and S components of w.  The pericenter angle is actually advancing at an average rate dω/dϕ≡ω' ≃ 3(m/p)  P G^-1 [see Eq. (<ref>)], which we anticipate is much larger than the preferred-frame effects shown in Eq. (<ref>) – the goal is to set strong upper bounds on such effects.   Thus the variations in the orbit elements will be modulated on a pericenter precession timescale and could even change sign.    So in order to find the proper long-term evolution of the elements, we define, for a given element X_a, dX_a/dϕ≡Δ X_a /2π, insert ω = ω_0 + ω' ϕ in the right-hand-sides of Eqs. (<ref>), including in the P, Q, R and S components of various vectors, and integrate with respect to ϕ.   As we will be interested in low-eccentricity binary pulsars, we will assume that e ≪ 1.   Inserting the resulting integrals for e and  ϖ, into the expression

    r = a [1-ecos(ϕ - ω) + O(e^2)]  ,

and expanding to first order in the preferred-frame perturbations, we obtain

    r/a   = 1 - e_0 cos (ϕ - ω_0 - ω' ϕ ) - 1/4α̂_1 Δ (m/a )^1/2w/ω' ( ŵ_⊥cosϕ - ŵ_Ωsinϕ )
    
         
    + e_0/4α̂_2 w^2 sinω' ϕ/ω' [ 2 ŵ_Ωŵ_⊥cos( ϕ + ω_0) 
    +  ( ŵ_⊥^2 - ŵ_Ω^2  ) sin(ϕ + ω_0)  ] .
 
The first term in Eq. (<ref>) is the normal contribution to r/a resulting from the small eccentricity e_0, with the pericenter advancing at a rate ω'.  The second term is a forced eccentricity of the orbit, with an amplitude proportional to ( ŵ_⊥^2 + ŵ_Ω^2)^1/2 = (1 - ŵ_h^2)^1/2≡sinψ, where ψ is the angle between the orbital angular momentum ĥ and the velocity w relative to the preferred frame, and a phase given by tan^-1 (-ŵ_Ω / ŵ_⊥ ) <cit.>.   This effect is present even in the limit e_0 → 0.  The final term is also a polarization of the orbit, proportional to e_0, with an amplitude modulated by the factor sinω' ϕ/ω'.   However, it vanishes in the limit e_0 → 0.  

The other important effect of the preferred-frame perturbations is to cause the orbital angular momentum to precess.  Since ĥ = sinι (sinΩe_x - cosΩe_y ) + cosιe_z, variations in ĥ are given by

    Δĥ  = sinιΔΩe_Ω - Διe_⊥ .

Inserting the expressions for Δι and ΔΩ from Eqs. (<ref>), taking the small e limit, and noting that e_Ωcosω + e_⊥sinω = e_ P and that A_⊥e_Ω - A_Ωe_⊥ = A×ĥ we obtain

    Δĥ   = - π/2α̂_1 Δ (m/a )^1/2 w_ h e e_ P
    -πα̂_2 w_h (w×ĥ )   ,

leading to a precession of the angular momentum vector h.  

Searches for eccentricities induced by the α̂_1 term in Eq. (<ref>) resulted in bounds on α̂_1  as small as a few parts in 10^5 <cit.>.   The tightest bound used a specific binary pulsar J1738+0333, whose orbit around its white-dwarf companion has an eccentricity 3.4 × 10^-7.  The analysis was helped by the fact that the white dwarf is bright enough to be observed spectroscopically, leading to accurate determinations of the key orbital parameters.  Furthermore, because the pericenter advances at a rate of about 1.6   deg   yr^-1, the decade-long data span made it possible to partially separate any induced eccentricity, whose direction is fixed by the direction of w, from the natural eccentricity, which rotates with the pericenter.  For this system the result was |α̂_1| < 3.4 × 10^-5 <cit.>.

Limits on α̂_2  were obtained by looking for the precession of the orbital plane of a binary  system [see Eq. (<ref>)].   Such a precession would lead to a variation in the “projected semimajor axis” of the pulsar, a_p sinι, a quantity that is measured very accurately in binary pulsar timing.   Combining data from the two wide-binary millisecond pulsar systems J1738+0333 and 
J1012+5307, Shao and Wex <cit.> obtained the bound |α̂_2 | < 1.8 × 10^-4.    

A given theory of gravity can be constrained or ruled out by combining these bounds on α̂_1 and α̂_2 with Eqs. (<ref>), together with estimates of the sensitivities of any compact bodies in the system.  Those sensitivities will depend on both the theory and the equation of state of nuclear matter.

Returning to the full N-body Lagrangian (<ref>), and working at quasi-Newtonian order, we can derive the leading signal of the failure of the universality of free fall (Nordtvedt effect) in a hierarchical three-body system, such as the pulsar J0337+1715 in a triple system with two white dwarf companions.   For a  two-body system in the
presence of a third body, the equations of motion become

    a_1    = -
     G_12 m_2 _12/r_12^3 
    -  G_13 m_3 _13/r_13^3 ,
    
     a_2    =  G_12 m_1 _12/r_12^3 
    -  G_23 m_3 _23/r_23^3 .
  
Following the method described for example in Sec. 13.3.3 of <cit.>, it is straightforward to show that, for nearly circular coplanar orbits, the perturbation of the inner orbit induced by the Nordtvedt effect is given by 

    δ r = - η̂_N R/a_0ω_b^2(1+2ω_b/Λ)/ω_b^2 - Λ^2cos(Λ t + Φ)   ,

where a_0 and R are the semimajor axes of the inner and outer orbits, respectively, ω_b is the angular frequency of the inner orbit, and Λ≡ω_b - ω_3 is the difference between the inner and outer orbit angular frequencies.   The “strong-field” Nordtvedt parameter η̂_N is given by

    η̂_N ≡ G_12 -  G_13 .

For the pulsar J0337+1715 we can ignore the sensitivities of the two white-dwarf companions, and write, for the specific cases of scalar-tensor and Einstein-Æther theories,

    η̂_N =  {[             - ζ s_1    :  scalar-tensor;         s_1/(1-s_1) :  Einstein-Æther , ]

where we have set G_N=1 in Einstein-Æther theory.   Bounds on the Nordtvedt effect signal obtained from the data  can then be used to constrain these specific theories. 
  
 


§ CONCLUDING REMARKS



We have extended the modified EIH framework to incorporate the possibility of preferred-frame effects, giving a direct link between bounds on such effects derived from observations of binary pulsar systems and the fundamental parameters of alternative theories of gravity.  In contrast to the simplicity of the PPN formalism, the link between observation and theory here depends on the internal structure (sensitivities) of the compact bodies in the system.     

This work was supported in part by the National Science Foundation,
Grant No. PHY 16–00188.   We are grateful to Enrico Barausse, Guillaume Faye and Viraj Sanghai for useful discussions.






§ THE COMPACT-BODY LAGRANGIAN FORMALISM OF NORDTVEDT


Nordtvedt's <cit.> compact body formalism introduced two kinetic parameters δ M_a^(2), and δ M_a^(4), six symmetric parameters, Γ_ab, Θ_ab, τ_ab,  σ_ab, and ξ_ab, two antisymmetric parameters Φ_ab and Ψ_ab, and one three-body parameter, Γ_abc, symmetric on the second two indices.   Our mass rescaling argument implies that δ M_a^(2)=0, and our rejection of Whitehead terms implies that 
ξ_ab = Ψ_ab = 0.   The remaining parameters are related to our modified EIH parameters by

    G_ab    = (m_a m_b)^-1Γ_ab ,
    
     B_ab   =1/3  (m_a m_b)^-1 ( Γ_ab + 2γΘ_ab - 4 Φ_ab )  ,
    
     D_abc   = (m_a m_b m_c)^-1 (2β-1)Γ_abc ,
    
     A_a     = m_a^-1δ M_a^(4) ,
    
     C_ab   = (m_a m_b)^-1 (α_1 τ_ab - 2α_2 σ_ab ) ,
    
     E_ab   = 2 (m_a m_b)^-1α_2 σ_ab ,

where γ, β, α_1 and α_2 are the standard PPN parameters.



§ REFERENCES


iopart-num



