Tests of Lorentz symmetry in the gravitational sector

Lorentz symmetry is one of the pillars of both General Relativity and the Standard Model of particle physics. Motivated by ideas about quantum gravity, unification theories and violations of CPT symmetry, a significant effort has been put the last decades into testing Lorentz symmetry. This review focuses on Lorentz symmetry tests performed in the gravitational sector. We briefly review the basics of the pure gravitational sector of the Standard-Model Extension (SME) framework, a formalism developed in order to systematically parametrize hypothetical violations of the Lorentz invariance. Furthermore, we discuss the latest constraints obtained within this formalism including analyses of the following measurements: atomic gravimetry, Lunar Laser Ranging, Very Long Baseline Interferometry, planetary ephemerides, Gravity Probe B, binary pulsars, high energy cosmic rays,... In addition, we propose a combined analysis of all these results. We also discuss possible improvements on current analyses and present some sensitivity analyses for future observations.


Introduction
The year 2015 was the centenary of the theory of General Relativity (GR), the current paradigm for describing the gravitational interaction (see e.g. the Editorial of this special issue [1]). Since its creation, this theory has passed all experimental tests with flying colors [2,3] ; the last recent success was the discovery of gravitational waves [4], summarized in [5]. On the other hand, the three other fundamental interactions of Nature are described within the Standard Model of particle physics, a framework based on relativistic quantum field theory. Although very successful so far, it is commonly admitted that these two theories are not the ultimate description of Nature but rather some effective theories. This assumption is motivated by the construction of a quantum theory of gravitation that has not been successfully developed so far and by the development of a theory that would unify all the fundamental interactions. Moreover, observations requiring the introduction of Dark Matter and Dark Energy also challenge GR and the Standard Model of particle physics since they cannot be explained by these two paradigms altogether [6]. It is therefore extremely important to test our current description of the four fundamental interactions [7].
Lorentz invariance is one of the fundamental symmetry of relativity, one of the corner stones of both GR and the Standard Model of particle physics. It states that the outcome of any local experiment is independent of the velocity and of the orientation of the laboratory in which the experiment is performed. If one considers non-gravitational experiments, Lorentz symmetry is part of the Einstein Equivalence Principle (EEP). A breaking of Lorentz symmetry implies that the equations of motion, the particle thresholds, etc... may be different when the experiment is boosted or rotated with respect to a background field [8]. More precisely, it is related to a violation of the invariance under "particle Lorentz transformations" [8] which are the boosts and rotations that relate the properties of two systems within a specific oriented inertial frame (or in other words they are boosts and rotations on localized fields but not on background fields). On the other hand, the invariance under coordinates transformations known as "observer Lorentz transformations" [8] which relate observations made in two inertial frames with different orientations and velocities is always preserved. Considering the broad field of applicability of this symmetry, searches for Lorentz symmetry breaking provide a powerful test of fundamental physics. Moreover, it has been suggested that Lorentz symmetry may not be a fundamental symmetry of Nature and may be broken at some level. While some early motivations came from string theories [9][10][11], breaking of Lorentz symmetry also appears in loop quantum gravity [12][13][14][15], non commutative geometry [16,17], multiverses [18], brane-world scenarios [19][20][21] and others (see for example [22,23]). Tests of Lorentz symmetry have been performed since the time of Einstein but the last decades have seen the number of tests increased significantly [24] in all fields of physics. In particular, a dedicated effective field theory has been developed in order to systematically consider all hypothetical violations of the Lorentz invariance. This framework is known as the Standard-Model Extension (SME) [8,25] and covers all fields of physics. It contains the Standard Model of particle physics, GR and all possible Lorentz-violating terms that can be constructed at the level of the Lagrangian, introducing a large numbers of new coefficients that can be constrained experimentally.
In this review, we focus on the gravitational sector of the SME which parametrizes deviations from GR. GR is built upon two principles [2,26,27]: (i) the EEP and (ii) the Einstein field equations that derive from the Einstein-Hilbert action. The EEP gives a geometric nature to gravitation allowing this interaction to be described by spacetime curvature. From a theoretical point of view, the EEP implies the existence of a spacetime metric to which all matter minimally couples [28]. A modification of the matter part of the action will lead to a breaking of the EEP. In SME, such a breaking of the EEP is parametrized (amongst others) by the matter-gravity coupling coefficientsā µ andc µν [29,30]. From a phenomenological point of view, the EEP states that [2,27]: (i) the universality of free fall (also known as the weak equivalence principle) is valid, (ii) the outcome of any local non-gravitational experiment is independent of the velocity of the free-falling reference frame in which it is performed and (iii) the outcome of any local non-gravitational experiment is independent of where and when in the universe it is performed. The second part of Einstein theory concerns the purely gravitational part of the action (the Einstein-Hilbert action) which is modified in SME to introduce hypothetical Lorentz violations in the gravitational sector. This review focuses exclusively on this kind of Lorentz violations and not on breaking of the EEP.
A lot of tests of GR have been performed in the last decades (see [2] for a review). These tests rely mainly on two formalisms: the parametrized post-Newtonian (PPN) framework and the fifth force formalism. In the former one, the weak gravitational field spacetime metric is parametrized by 10 dimensionless coefficients [27] that encode deviations from GR. This formalism therefore provides a nice interface between theory and experiments. The PPN parameters have been constrained by a lot of different observations [2] confirming the validity of GR. In particular, three PPN parameters encode violations of the Lorentz symmetry: the α 1,2,3 PPN coefficients. In the fifth force formalism, one is looking for a deviation from Newtonian gravity where the gravitational potential takes the form of a Yukawa potential characterized by a length λ and a strength α of interaction [31][32][33][34]. These two parameters are very well constrained as well except at very small and large distances (see [35]).
The gravitational sector of SME offers a new framework to test GR by parametrizing deviations from GR at the level of the action, introducing new terms that are breaking Lorentz symmetry. The idea is to extend the standard Einstein-Hilbert action by including Lorentz-violating terms constructed by contracting new fields with some operators built from curvature tensors and covariant derivatives with increasing mass dimension [36]. The lower mass dimension (dimension 4) term is known as the minimal SME and its related new fields can be split into a scalar part u, a symmetric trace free part s µν and a traceless piece t κλµν . In order to avoid conflicts with the underlying Riemann geometry, the Lorentz violating coefficients can be assumed to be dynamical fields and the Lorentz violation to arise from a spontaneous symmetry breaking [37][38][39][40][41][42]. The Lorentz violating fields therefore acquire a non-vanishing vacuum expectation value (denoted by a bar). It has been shown that in the linearized gravity limit the fluctuations around the vacuum values can be integrated out so that only the vacuum expectation values of the SME coefficients influence observations [39]. In the minimal SME, the coefficientū corresponds to a rescaling of the gravitational constant and is therefore unobservable and the coefficientst κλµν do not play any role at the post-Newtonian level, a surprising phenomenon known as the t-puzzle [43,44]. Thes µν coefficients lead to modifications from GR that have thoroughly been investigated in [39]. In particular, the SME framework extends standard frameworks such as the PPN or fifth force formalisms meaning that "standard" tests of GR cannot directly be translated into this formalism.
In the last decade, several measurements have been analyzed within the gravitational sector of the minimal SME framework: Lunar Laser Ranging (LLR) analysis [45,46], atom interferometry [47,48], planetary ephemerides analysis [49,50], short-range gravity [51], Gravity Probe B (GPB) analysis [52], binary pulsars timing [53,54], Very Long Baseline Interferometry (VLBI) analysis [55] andČerenkov radiation [56]. In addition to the minimal SME, there exist some higher order Lorentz-violating curvature couplings in the gravity sector [43] that are constrained by short-range experiments [57][58][59],Čerenkov radiation [30,56] and gravitational waves analysis [60,61]. Finally, some SME experiments have been used to derive bounds on spacetime torsion [62,63]. A review for these measurements can be found in [30]. The classic idea to search for or to constrain Lorentz violations in the gravitational sector is to search for orientation or boost dependence of an observation. Typically, one will take advantage of modulations that will occur through an orientation dependence of the observations due to the Earth's rotation, the motion of satellites around Earth (the Moon or artificial satellites), the motion of the Earth (or other planets) around the Sun, the motion of binary pulsars, . . . The main goal of this communication is to review all the current analyses performed in order to constrain Lorentz violation in the pure gravitational sector.
Two distinct procedures have been used to analyze data within the SME framework. The first procedure consists in deriving analytically the signatures produced by the SME coefficients on some observations. Then, the idea is to fit these signatures within residuals obtained by a data analysis performed in pure GR. This approach has the advantage to be relatively easy and fast to perform. Nevertheless, when using this postfit approach, correlations with other parameters fitted in the data reduction are completely neglected and may lead to overoptimistic results. A second way to analyze data consists in introducing the Lorentz violating terms directly in the modeling of observables and in the global data reduction. In this review, we highlight the differences between the two approaches.
In this communication, a brief theoretical review of the SME framework in the gravitational sector is presented in Section 2. The two different approaches to analyze data within the SME framework (postfit analysis versus full modeling of observables within the SME framework) are discussed and compared in Section 3. Section 4 is devoted to a discussion of the current measurements analyzed within the SME framework. This discussion includes a general presentation of the measurements, a brief review of the effects of Lorentz violation on each of them, the current analyses performed with real data and a critical discussion. A "grand fit" combining all existing analyses is also presented. In Section 5, some future measurements that are expected to improve the current analyses are developed. Finally, our conclusion is presented in Section 6.

The Standard-Model Extension in the gravitational sector
Many of the tests of Lorentz and CPT symmetry have been analyzed within an effective field theory framework which generically describes possible deviations from exact Lorentz and CPT invariance [8,25] and contains some traditional test frameworks as limiting cases [64,65]. This framework is called, for historical reasons, the Standard-Model Extension (SME). One part of the activity has been a resurgence of interest in tests of relativity in the Minkowski spacetime context, where global Lorentz symmetry is the key ingredient. Numerous experimental and observational constraints have been obtained on many different types of hypothetical Lorentz and CPT symmetry violations involving matter [24]. Another part, which has been developed more recently, has seen the SME framework extended to include the curved spacetime regime [37]. Recent work shows that there are many ways in which the spacetime symmetry foundations of GR can be tested [29,39].
In the context of effective field theory in curved spacetime, violations of these types can be described by an action that contains the usual Einstein-Hilbert term of GR, a matter action, plus a series of terms describing Lorentz violation for gravity and matter in a generic way. While the fully general coordinate invariant version of this action has been studied in the literature, we focus on a limiting case that is valid for weak-field gravity and can be compactly displayed. Using an expansion of the spacetime metric around flat spacetime, g µν = η µν + h µν , the effective Lagrange density to quadratic order in h µν can be written in a compact form as where L EH is the standard Einstein-Hilbert term, G αµνβ is the double dual of the Einstein tensor linearized in h µν , G the bare Newton constant and c the speed of light in a vacuum. The Lorentz-violating effects in this expression are controlled by the 9 independent coefficients in the traceless and dimensionlesss µν [39]. These coefficients are treated as constants in asymptotically flat cartesian coordinates. The ellipses represent additional terms in a series including terms that break CPT symmetry for gravity; such terms are detailed elsewhere [43,56,60] and are part of the so-called nonminimal SME expansion. Note that the process by which one arrives at the effective quadratic Lagrangian (1) is consistent with the assumption of the spontaneous breaking of local Lorentz symmetry, which is discussed below. Also of interest are the matter-gravity couplings. This form of Lorentz violation can be realized in the classical point-mass limit of the matter sector. In the minimal SME the point-particle action can be written as where the particle's worldline tangent is u µ = dx µ /dλ [29]. The coefficients controlling local Lorentz violation for matter are c µν and a µ . In contrast tos µν , these coefficients depend on the type of point mass (particle species) and so they can also violate the EEP. When the coefficientss µν , c µν , and a µ vanish perfect local Lorentz symmetry for gravity and matter is restored. It is also interesting to mention that this action with fixed (but not necessarily constant) a µ and c µν represents motion in a Finsler geometry [66,67]. It has been shown that explicit local Lorentz violation is generically incompatible with Riemann geometry [37]. One natural way around this is assumption of spontaneous Lorentz-symmetry breaking. In this scenario, the tensor fields in the underlying theory acquire vacuum expectation values through a dynamical process. Much of the literature has been devoted to studying this possibility in the last decades [9,38,[68][69][70][71][72][73][74][75][76][77][78], including some original work on spontaneous Lorentz-symmetry breaking in string field theory [10,11]. For the matter-gravity couplings in Eq. (2), the coefficient fields c µν , and a µ are then expanded around their background (or vacuum) valuesc µν , andā µ . Both a modified spacetime metric g µν and modified point-particle equations of motion result from the spontaneous breaking of Lorentz symmetry. In the linearized gravity limit these results rely only on the vacuum valuesc µν , andā µ . The dominant signals for Lorentz violation controlled by these coefficients are revealed in the calculation of observables in the post-Newtonian limit.
Several novel features of the post-Newtonian limit arise in the SME framework. It was shown in Ref. [39] that a subset of thes µν coefficients can be matched to the PPN formalism [2,27], but others lie outside it. For example, a dynamical model of spontaneous Lorentz symmetry breaking can be constructed from an antisymmetric tensor field B µν that producess µν coefficients that cannot be reduced to an isotropic diagonal form in any coordinate system, thus lying outside the PPN assumptions [78]. We can therefore see that the SME framework has a partial overlap with the PPN framework, revealing new directions to explore in analysis via thes µν ,c µν , andā µ coefficients. The equations of motion for matter are modified by the matter-gravity coefficients for Lorentz violation c µν andā µ , which can depend on particle species, thus implying that these coefficients also control EEP violations. One potentially important class of experiments from the action (2) concerns the Universality of Free Fall of antimatter whose predictions are discussed in [29,79]. In addition, the post-Newtonian metric itself receives contributions from the matter coefficientsc µν andā µ . So for example, two (chargeless) sources with the same total mass but differing composition will yield gravitational fields of different strength.
For solar-system gravity tests, the primary effects due to the nine coefficientss µν can be obtained from the post-Newtonian metric and the geodesic equation for test bodies. A variety of ground-based and space-based tests can measure these coefficients [80][81][82]. Such tests include Earth-laboratory tests with gravimeters, lunar and satellite laser ranging, studies of the secular precession of orbital elements in the solar system, and orbiting gyroscope experiments, and also classic effects such as the time delay and bending of light around the Sun and Jupiter. Furthermore, some effects described by the Lagrangian (1) can be probed by analyzing data from binary pulsars and measurements of cosmic rays [56].
For the matter-gravity coefficientsc µν andā µ , which break Lorentz symmetry and EEP, several experiments can be used for analysis in addition to the ones already mentioned above including ground-based gravimeter and WEP experiments. Dedicated satellite EEP tests are among the most sensitive where the relative acceleration of two test bodies of different composition is the observable of interest. Upon relating the satellite frame coefficients to the standard Sun-centered frame used for the SME, oscillations in the acceleration of the two masses occur at a number of different harmonics of the satellite orbital and rotational frequencies, as well as the Earth's orbital frequency. Future tests of particular interest include the currently flying MicroSCOPE experiment [83,84].
While the focus of the discussion to follow are the results for the minimal SME coefficientss µν , recent work has also involved the nonminimal SME coefficients in the pure-gravity sector associated with mass dimension 5 and 6 operators. One promising testing ground for these coefficients is sensitive short-range gravity experiments. The Newtonian force between two test masses becomes modified in the presence of local Lorentz violation by an anisotropic quartic force that is controlled by a subset of coefficients from the Lagrangian organized as the totally symmetric (k eff ) jklm , which has dimensions of length squared [43]. This contains 14 measurable quantities and any one short-range experiment is sensitive to 8 of them. Two key experiments, from Indiana University and Huazhong University of Science and Technology, have both reported analysis in the literature [57,58] . A recent work combines the two analyses to place new limits on all 14, a priori independent, (k eff ) jklm coefficients [59]. Other higher mass dimension coefficients play a role in gravitational wave propagation [60] and gravitationalČerenkov radiation [56].
To conclude this section, we ask: what can be said about the possible sizes of the coefficients for Lorentz violation? A broad class of hypothetical effects is described by the SME effective field theory framework, but it is a test framework and as such does not make specific predictions concerning the sizes of these coefficients. One intriguing suggestion is that there is room in nature for violations of spacetime symmetry that are large compared to other sectors due to the intrinsic weakness of gravity. Considering the current status of the coefficientss µν , the best laboratory limits are at the 10 −10 -10 −11 level, with improvements of four orders of magnitude in astrophysical tests on these coefficients [56]. However, the limits are at the 10 −8 m 2 level for the mass dimension 6 coefficients (k eff ) jklm mentioned above. Comparing this to the Planck length 10 −35 m, we see that symmetry breaking effects could still have escaped detection that are not Planck suppressed. This kind of "countershading" was first pointed out for theā µ coefficients [85], which, having dimensions of mass, can still be as large as a fraction of the electron mass and still lie within current limits.
In addition, any action-based model that breaks local Lorentz symmetry either explicitly or spontaneously can be matched to a subset of the SME coefficients. Therefore, constraints on SME coefficients can directly constrain these models. Matches between various toy models and coefficients in the SME have been achieved for models that produce effectives µν ,c µν ,ā µ , and other coefficients. This includes vector and tensor field models of spontaneous Lorentz-symmetry breaking [29,39,[75][76][77][78], models of quantum gravity [12,65] and noncommutative quantum field theory [17]. Furthermore, Lorentz violations may also arise in the context of string field theory models [86].

Postfit analysis versus full modeling
Since the last decade, several studies aimed to find upper limits on SME coefficients in the gravitational sector. A lot of these studies are based on the search of possible signals in post-fit residuals of experiments. This was done with LLR [45], GPB [52], binary pulsars [53,54] or Solar System planetary motions [49,50]. However, two new works focused on a direct fit to data with LLR [46] and VLBI [55], which are more satisfactory.
Indeed, in the case of a post-fit analysis, a simple modeling of extra terms containing SME coefficients are least square fitted in the residuals, attempting to constrain the SME coefficients of a testing function in residual noise obtained from a pure GR analysis, where of course Lorentz symmetry is assumed. It comes out correlations between SME coefficients and other global parameters previously fitted (masses, position and velocity. . . ) cannot be assessed in a proper way. In others words, searching hypothetical SME signals in residuals, i.e. in noise, can lead to an overestimated formal error on SME coefficients, as illustrated in the case of VLBI [55], and without any chance to learn something about correlations with other parameters, as for example demonstrated in the case of LLR [46]. Let us consider the VLBI example to illustrate this fact. The VLBI analysis is described in Section 4.2. Including the SME contribution within the full VLBI modeling and estimating the SME coefficients TT altogether with the other parameters fitted in standard VLBI data reduction leads to the estimates TT = (−5 ± 8) × 10 −5 . A postfit analysis performed by fitting the SME contribution within the VLBI residuals obtained after a pure GR analysis leads tō s TT = (−0.6 ± 2.1) × 10 −8 [55]. This example shows that a postfit analysis can lead to results with overoptimistic uncertainties and one needs to be extremely careful when using such results.

Data analysis
In this section, we will review the different measurements that have already been used in order to constrain the SME coefficients. The different analyses are based on quite different types of observations. In order to compare all the corresponding results, we need to report them in a canonical inertial frame. The standard canonical frame used in the SME framework is a Sun-centered celestial equatorial frame [64], which is approximately inertial over the time scales of most observations. This frame is asymptotically flat and comoving with the rest frame of the Solar System. The cartesian coordinates related to this frame are denoted by capital letters The Z axis is aligned with the rotation axis of the Earth, while the X axis points along the direction from the Earth to the Sun at vernal equinox. The origin of the coordinate time T is given by the time when the Earth crosses the Sun-centered X axis at the vernal equinox. These conventions are depicted in Figure 2 from [39]. In the following subsections, we will present the different measurements used to constrain the SME coefficients. Each subsection contains a brief description of the principle of the experiment, how it can be used to search for Lorentz symmetry violations, what are the current best constraints obtained with such measurements and eventually how it can be improved in the future.

Atomic gravimetry
The most sensitive experiments on Earth searching for Lorentz Invariance Violation (LIV) in the minimal SME gravity sector are gravimeter tests. As Earth rotates, the signal recorded in a gravimeter, i.e. the apparent local gravitational acceleration g of a laboratory test body, would be modulated in the presence of LIV in gravity. This was first noted by Nordtvedt and Will in 1972 [87] and used soon after with gravimeter data to constrain preferred-frame effects in the PPN formalism [88,89] at the level of 10 −3 .
This test used a superconducting gravimeter, based on a force comparison (the gravitational force is counter-balanced by an electromagnetic force maintaining the test mass at rest). While superconducting gravimeters nowadays reach the best sensitivity on Earth, force comparison gravimeters intrinsically suffer from drifts of their calibration factor (with e.g. aging of the system). Development of other types of gravimeters has evaded this drawback: free fall gravimeters. Monitoring the motion of a freely falling test mass, they provide an absolute measurement of g. State-of-the art free fall gravimeters use light to monitor the mass free fall. Beyond classical gravimeters that drop a corner cube, the development of atom cooling and trapping techniques and atom interferometry has led to a new generation of free fall gravimeters, based on a quantum measurement: atomic gravimeters.
Atomic gravimeters use atoms in gaseous phase as a test mass. The atoms are initially trapped with magneto-optical fields in vacuum, and laser cooled (down to 100 nK) in order to control their initial velocity (down to a few mm/s). The resulting cold atom gas, containing typically a million atoms, is then launched or dropped for a free fall measurement. Manipulating the electronic and motional state of the atoms with two counterpropagating lasers, it is possible to measure, using atom interferometry, their free fall acceleration with respect to the local frame defined by the two lasers [90]. This sensitive direction is aligned to be along the local gravitational acceleration notedẑ; the atom interferometer then measures the phase ϕ = kaẑT 2 , where T is half the interrogation time, k 2(2π/λ) with λ the laser wavelength, and aẑ is the free fall acceleration along the laser direction. The free fall time is typically on the order of 500 ms, corresponding to a free fall distance of about a meter. A new "atom preparation -free fall -detection" cycle is repeated every few seconds. Each measurement is affected by white noise, but averaging leads to a typical sensitivity on the order of or below 10 −9 g [91][92][93].
Such an interferometer has been used by H. Müller et al. in [47] and K. Y. Chung et al. in [48] for testing Lorentz invariance in the gravitational sector with Caesium atoms, leading to the best terrestrial constraints on thes µν coefficients. The analysis uses three data sets of respectively 2.5 days for the first two and 10 days for the third, stretched over 4 years, which allows one to observe sidereal and annual LIV signatures. The gravitational SME model used for this analysis can be found in [39,47,48]; its derivation will be summarized hereunder. Since the atoms in free fall are sensitive to the local phase of the lasers, LIV in the interferometer observable could also come from the pure electromagnetic sector. This contribution has been included in the experimental analysis in [48]. Focusing here on the gravitational part of SME, we ignore it in the following.
The gravitational LIV model adjusted in this test restricts to modifications of the Earth-atom two-body gravitational interaction. The Lagrangian describing the dynamics of a test particle at a point on the Earth's surface can be approximated by a post-Newtonian series as developed in [39]. At the Newtonian approximation, the two bodies Lagrangian is given by where R and V are the position and velocities expressed in the standard SME Sun-centered frame and R = R/R with R = |R|. In addition, we have introduced the G N is the observed Newton constant measured by considering the orbital motion of bodies and defined by (see also [39,50] or Section IV of [52]) and the 3-dimensional traceless tensors From this Lagrangian one can derive the equations of motion of the free fall mass in a laboratory frame (see the procedure in Section V.C.1. from [39]). It leads to the modified local acceleration in the presence of LIV [39] given by where is the Earth's boost, R ⊕ is the Earth radius, M ⊕ is the Earth mass and χ the colatitude of the lab whose reference frame'ŝ z direction is the sensitive axis of the instrument as previously defined here. This model includes the shape of the Earth through its spherical moment of inertia I ⊕ which appears in In [48], Earth has been approximated as spherical and homogeneous leading to i ⊕ = 1 2 , i 1 = 7 6 and i 4 = − 1 2 . The sensing direction of the experiment precesses around the Earth rotation axis with sidereal period, and the lab velocity varies with sidereal period and annual period. At first order in V ⊕ and ω ⊕ and as a function of the SME coefficients, the LIV signal takes the form of a harmonic series with sidereal and annual base frequencies (denoted resp. ω ⊕ and Ω) together with first harmonics. The time dependence of the measured acceleration aẑ from Eq. (7) arises from the terms involving theẑ indices. It can be decomposed in frequency according to [39] δaẑ aẑ The model contains seven frequencies l ∈ {Ω, ω ⊕ , 2ω ⊕ , ω ⊕ ± Ω, 2ω ⊕ ± Ω}. The 14 amplitudes C l and D l are linear combinations of 7s µν components:s JK ,s T J ands XX −s YY which can be found in Table 1 of [48] or Table IV from [39]. In order to look for tiny departures from the constant Earth-atom gravitational interaction, a tidal model for aẑ variations due to celestial bodies is removed from the data before fitting to Eq. (8). This tidal model consists of two parts. One part is based on a numerical calculation of the Newtonian tide-generating potential from the Moon and the Sun at Earth's surface based on ephemerides. It uses here the Tamura tidal catalog [94] which gives the frequency, amplitude and phase of 1200 harmonics of the tidal potential. These arguments are used by a software (ETGTAB) that calculates the time variation of the local acceleration in the lab and includes the elastic response of Earth's shape to the tides, called "solid Earth tides", also described analytically e.g. by the DDW model [95]. A previous SME analysis of the atom gravimeter data using only this analytical tidal correction had been done, but it led to a degraded sensitivity of the SME test [47]. Indeed, a non-negligible contribution to aẑ is not covered by this non-empirical tidal model: oceanic tide effects such as ocean loading, for which good global analytical models do not exist. They consequently need to be adjusted from measurements. For the second analysis, reported here, additional local tidal corrections fitted on altimetric data have been removed [96] allowing to improve the statistical uncertainty of the SME test by one order of magnitude.
After tidal subtraction, signal components are extracted from the data using a numerical Fourier transform (NFT). Due to the finite data length, Fourier components overlap, but the linear combinations of spectral lines that the NFT estimates can be expressed analytically. Since the annual component ω l = Ω has not been included in this analysis, the fit provides 12 measurements. From there, individual constraints on the 7 SME coefficients and their associated correlation coefficients can be estimated by a least square adjustment. The results obtained are presented in Table 1. Table 1. Atom-interferometry limits on Lorentz violation in gravity from [48]. The correlation coefficients can be derived from Table III of [48].
Correlation coefficients All results obtained are compatible with null Lorentz violation. As expected from boost suppressions in Eq. (7) and from the measurement uncertainty, on the order of a few 10 −9 g [97], typical limits obtained are in the 10 −9 range for purely spatials µν components and 4 orders of magnitude weaker for the spatio-temporal componentss T J . It can be seen e.g. with the purely spatial components that these constraints do not reach the intrinsic limit of acceleration resolution of the instrument (which has a short term stability of 11 × 10 −9 g/ √ Hz) because the coefficients are still correlated. Their marginalized uncertainty is broadened by their correlation.
Consequently, improving the uncertainty could be reached through a better decorrelation, by analyzing longer data series. In parallel, the resolution of these instruments keeps increasing and has nowadays improved by about a factor 10 since this experiment. However, increasing the instrument's resolution brings back to the question of possible accidental cancelling in treating "postfit" data. Indeed, it should be recalled here that local tidal corrections subtracted prior to analysis are based on adjusting a model of ocean surface from altimetry data. In principle, this observable would as well be affected by gravity LIV; fitting to these observations thus might remove part of SME signatures from the atom gravimeter data. This was mentioned in the first atom gravimeter SME analysis [47]. The adjustment process used to assess local corrections in gravimeters is not made directly on the instrument itself, but it always involves a form of tidal measurement (here altimetry data, or gravimetry data from another instrument in [98]). All LIV frequencies match to the main tidal frequencies. Further progress on SME analysis with atom gravimeters would thus benefit from addressing in more details the question of possible signal canceling.

Very Long Baseline Interferometry
VLBI is a geometric technique measuring the time difference in the arrival of a radio wavefront, emitted by a distant quasar, between at least two Earth-based radio-telescopes. VLBI observations are done daily since 1979 and the database contains nowadays almost 6000 24hours sessions, corresponding to 10 millions group-delay observations, with a present precision of a few picoseconds. One of the principal goals of VLBI observations is the kinematical monitoring of Earth rotation with respect to a global inertial frame realized by a set of defining quasars, the International Celestial Reference Frame [99], as defined by the International Astronomical Union [100]. The International VLBI Service for Geodesy and Astrometry (IVS) organizes sessions of observation, storage of data and distribution of products, in particular the Earth Orientation parameters. Because of this precision, VLBI is also a very interesting tool to test gravitation in the Solar System. Indeed, the gravitational fields of the Sun and the planets are responsible of relativistic effects on the quasar light beam through the propagation of the signal to the observing station and VLBI is able to detect these effects very accurately. By using the complete VLBI observations database, it was possible to obtain a constraint on the γ PPN parameter at the level of 1.2 × 10 −4 [101,102]. In its minimal gravitational sector, SME can also be investigated with VLBI and obtaining a constrain on thes TT coefficient is possible.
Indeed, the propagation time of a photon emitted at the event (cT e , X e ) and received at the position X r can be computed in the SME formalism using the time transfer function formalism [103][104][105][106][107] and is given by [39,80] T (X e , T e , where the terms a 1 and a 2 from [80] are taken as unity (which corresponds to using the harmonic gauge, which is the one used for VLBI data reduction), R e = |X e |, R r = |X r |, R er = |X r − X e | with the central body located at the origin and where we introduce the following vectors and where G N is the observed Newton constant measured by considering the orbital motion of bodies and is defined in Eq. (5). This equation is the generalization of the well-known Shapiro time delay including Lorentz violation. The VLBI is actually measuring the difference of the time of arrival of a signal received by two different stations. This observable is therefore sensitive to a differential time delay (see [108] for a calculation in GR). Assuming a radio-signal emitted by a quasar at event (T e , X e ) and received by two different VLBI stations at events (T 1 , X 1 ) and (T 2 , X 2 ) (all quantities being expressed in a barycentric reference frame), respectively, the VLBI group-delay ∆τ (SME) in SME formalism can be written [55] ∆τ (SME) = 2 where we only kept thes TT contribution (see Eq. (7) from [55] for the full expression) and we use the same notations as in [108] by introducing three unit vectors Ten million VLBI delay observations between August 1979 and mid-2015 have been used to estimate thes TT coefficient. First, VLBI observations are corrected from delay due to the radio wave crossing of dispersive media by using 2 GHz and 8 GHz recordings. Then, we used only the 8 GHz delays and the Calc/Solve geodetic VLBI analysis software, developed at NASA Goddard Space Flight Center and coherent with the latest standards of the International Earth Rotation and Reference Systems Service [109]. We added the partial derivative of the VLBI delay with respect tos TT from Eq. (11) to the software package using the USERPART module of Calc/Solve. We turned to a global solution in which we estimateds TT as a global parameter together with radio source coordinates. We obtaineds with a postfit root mean square of 28 picoseconds and a χ 2 per degree of freedom of 1.15. Correlations between radio source coordinates ands TT are lower than 0.02, the global estimate being consistent with the mean value obtained with the session-wise solution with a slightly lower error.
In conclusion, VLBI is an incredible tool to test Lorentz symmetry, especially thes TT coefficient. This coefficient has an isotropic impact on the propagation speed of gravitational waves as can be noticed from Eq. (27) below (or see Eq. (9) from [56] or Eq. (11) from [60]). The analysis performed in [55] includes the SME contribution in the modeling of VLBI observations and includes thes TT parameter in the global fit with other parameters. It is therefore a robust analysis that produces the current best estimate on thes TT parameter. In the future, the accumulation of VLBI data in the framework of the permanent geodetic monitoring program leads us expect improvement of this constraint.

Lunar Laser Ranging
On August, 20th 1969, after ranging to the lunar retro-reflector placed during the Apollo 11 mission, the first LLR echo was detected at the McDonald Observatory in Texas. Currently, there are five stations spread over the world which have realized laser shots on five lunar retro-reflectors. LLR is used to conduct high precision measurements of the light travel time of short laser pulses emitted at time t 1 by a LLR station, reflected at time t 2 by a lunar retro-reflector and finally received at time t 3 at a station receiver. The data are presented as normal points which combine time series of measured light travel time of photons, averaged over several minutes to achieve a higher signal-to-noise ratio measurement of the lunar range at some characteristic epoch. Each normal-point is characterized by one emission time (t 1 in universal time coordinate -UTC), one time delay (∆t c in international atomic time -TAI) and some additional observational parameters as laser wavelength, atmospheric temperature and pressure etc. According to [110], the theoretical pendent of the observed time delay (∆t c = t 3 − t 1 in TAI) is defined as where T 1 is the emission time expressed in barycentric dynamical time (TDB) and ∆τ t is a relativistic correction between the TDB and the terrestrial time (TT) at the level of the station. The reception time T 3 expressed in TDB is defined by the following two relations with T 2 the time in TDB at the reflection point X o and X o are respectively the barycentric position vector at the emitter and the reception point, X r is the barycentric position vector at the reflection point, ∆T (grav) is the one way gravitational time delay correction and ∆τ a is the one way tropospheric correction. LLR measurements are used to produce the Lunar ephemeris but also provide a unique opportunity to study the Moon's rotation, the Moon's tidal acceleration, the lunar rotational dissipation, etc [111]. In addition, LLR measurements have turn the Earth-Moon system into a laboratory to study fundamental physics and to conduct tests of the gravitation theory. Nordtvedt was the first to suggest that LLR can be used to test GR by testing one of its pillar: the Strong Equivalence Principle [112][113][114]. He showed that precise laser ranging to the Moon would be capable of measuring precisely the ratio of gravitational mass to inertial mass of the Earth to an accuracy sufficient to constrain a hypothetical dependence of this ratio on the gravitational self-energy. He concluded that such a measurement could be used to test Einstein's theory of gravity and others alternative theories as scalar tensor theories. The best current test of the Strong Equivalence Principle is provided by a combination of torsion balance measurements with LLR analysis and is given by [115][116][117] where η is the Nordtvedt parameter that is defined as m G /m I = 1 + ηU/mc 2 with m G the gravitational mass, m I the inertial mass and U the gravitational self-energy of the body. Using the Cassini constraint on the γ PPN parameter [118] and the relation η = 4β − γ − 3 leads to a constraint on β PPN parameter at the level β − 1 = (1.2 ± 1.1) × 10 −4 [116].
In addition to tests of the Strong Equivalence Principle, many other tests of fundamental physics were performed with LLR analysis. For instance, LLR data can be used to search for a temporal evolution of the gravitational constantĠ/G [115] and to constrain the fifth force parameters [119]. In addition, LLR has been used to constrain violation of the Lorentz symmetry in the PPN framework. Müller et al. [119] deduced from LLR data analysis constraints on the preferred frame parameters α 1 and α 2 at the level α 1 = (−7 ± 9) × 10 −5 and Considering all the successful GR tests performed with LLR observations, it is quite natural to use them to search for Lorentz violations in the gravitation sector. In the SME framework, Battat et al. [45] used the lunar orbit to provide estimates on the SME coefficients. Using a perturbative approach, the main signatures produced by SME on the lunar orbit have analytically been computed in [39]. These computations give a first idea of the amplitude of the signatures produced by a breaking of Lorentz symmetry. Nevertheless, these analytical signatures have been computed assuming the lunar orbit to be circular and fixed (i.e. neglecting the precession of the nodes for example). These analytical signatures have been fitted to LLR residuals obtained from a data reduction performed in pure GR [45]. They determined a "realistic" error on their estimates from a similar postfit analysis performed in the PPN framework. The results obtained by this analysis are presented in Table 2. It is important to note that this analysis uses projections of the SME coefficients into the lunar orbital planes 11 ,s 22 ,s 0i (see Section V.B.2 of [39]) while the standard SME analyses uses coefficients defined in a Sun-centered equatorial frame (and denoted by capital letters I J ). Table 2. Estimation of SME coefficients from LLR postfit data analysis from [45]. No correlations coefficients have been derived in this analysis. The coefficientss ij are projections of thes I J into the lunar orbital plane (see Eq. (107) from [39]) while the linear combinationss Ω ⊕ c ands Ω ⊕ s are given by Eq. (108) from [39].

Coefficient
However, as discussed in Section 3 and in [46,55], a postfit search for SME signatures into residuals of a data reduction previously performed in pure GR is not fully satisfactory. First of all, the uncertainties obtained by a postfit analysis based on a GR data reduction can be underestimated by up to two orders of magnitude. This is mainly due to correlations between SME coefficients and others global parameters (masses, positions and velocities, . . .) that are neglected in this kind of approach. In addition, in the case of LLR data analysis, the oscillating signatures derived in [39] and used in [45] to determine pseudo-constraints are computed only accounting for short periodic oscillations, typically at the order of magnitude of the mean motion of the Moon around the Earth. Therefore, this analytic solution remains only valid for few years while LLR data spans over 45 years (see also the discussions in footnote 2 from [50] and page 22 from [39]).
Regarding LLR data analysis, a more robust strategy consists in including the SME modeling in the complete data analysis and to estimate the SME coefficients in a global fit along with others parameters by taking into account short and long period terms and also correlations (see [46]). In order to perform such an analysis, a new numerical lunar ephemeris named "Éphéméride Lunaire Parisienne Numérique" (ELPN) has been developed within the SME framework. The dynamical model of ELPN is similar to the DE430 one [120] but includes the Lorentz symmetry breaking effects arising on the orbital motion of the Moon. The SME contribution to the lunar equation of motion has been derived in [39] and is given by where G N is the observed Newtonian constant defined by Eq. (5) In addition to the orbital motion, effects of a violation of Lorentz symmetry on the light travel time of photons is also considered. More precisely, the gravitational time delay ∆T (grav) appearing in Eq. (14) is given by the gravitational part of Eq. (9) [80].
Estimates on the SME coefficients are obtained by a standard chi-squared minimization: the LLR residuals are minimized by an iterative weighted least squares fit using partial derivatives previously computed from variational equations in ELPN. After an adjustment of 82 parameters including the SME coefficients a careful analysis of the covariance matrix shows that LLR data does not allow to estimate independently all the SME coefficients but that they are sensitive to the following three linear combinations:s XX −s YY ,s TY + 0.43s TZ ,s XX +s YY − 2s ZZ − 4.5s YZ .
The estimations on the 6 SME coefficients derived in [46] is summarized in Table 3. In particular, it is worth emphasizing that the quoted uncertainties are the sum of the statistical uncertainties obtained from the least-square fit with estimations of systematics uncertainties obtained with a Jackknife resampling method [121,122]. In summary, LLR is a powerful experiment to constrain gravitation theory and in particular hypothetical violation of the Lorentz symmetry. A first analysis based on a postfit estimations of the SME coefficients have been performed [45] which is not satisfactory regarding the neglected correlations with other global parameters as explained in Section 3. A full analysis including the integration of the SME equations of motion and the SME contribution to the gravitational time delay has been done in [46]. The resulting estimates on some SME coefficients are presented in Table 3. In addition, some SME coefficients are still correlated with parameters appearing in the rotational motion of the Moon as the principal moment of inertia, the quadrupole moment, the potential Stockes coefficient C 22 and the polar component of the velocity vector of the fluid core [46]. A very interesting improvement regarding this analysis would be to produce a joint GRAIL (Gravity Recovery And Interior Laboratory) [123][124][125] and LLR data analysis that would help in decorrelating the SME parameters from the lunar potential Stockes coefficients of degree 2 and therefore improve marginalized estimations of the SME coefficients. Finally, in [45] and [46], the effects of SME on the translational lunar equations of motion are considered and used to derive constraints on the SME coefficients. It would be also interesting to extend these analyses by considering the modifications due to SME on the rotation of the Moon. A first attempt has been proposed in Section V. A. 2. of [39] but needs to be extended.

Planetary ephemerides
The analysis of the motion of the planet Mercury around the Sun was historically the first evidence in favor of GR with the explanation of the famous advance of the perihelion in 1915. From there, planetary ephemerides have always been a very powerful tool to constrain GR and alternative theories of gravitation. Currently, three groups in the world are producing planetary ephemerides: the NASA Jet Propulsion Laboratory with the DE ephemerides [120,[126][127][128][129][130][131], the French INPOP (Intégrateur Numérique Planétaire de l'Observatoire de Paris) ephemerides [132][133][134][135][136][137] and the Russian EPM ephemerides [138][139][140][141][142]. These analyses use an impressive number of different observations to produce high accurate planetary and asteroid trajectories. The observations used to produce ephemerides comprise radioscience observations of spacecraft that orbited around Mercury, Venus, Mars and Saturn, flyby tracking of spacecraft close to Mercury, Jupiter, Uranus and Neptune and optical observations of all planets. This huge set of observations have been used to constrain the γ and β post-Newtonian parameter at the level of 10 −5 [136,137,[141][142][143], the fifth force interaction (see [32] and Figure 31 from [143]), the quantity of Dark Matter in our Solar System [144], the Modified Newtonian Dynamics [131,[145][146][147], . . .
A violation of Lorentz symmetry within the gravity sector of SME induces different types of effects that can have implications on planetary ephemerides analysis: effects on the orbital dynamics and effects on the light propagation. Simulations using the Time Transfer Formalism [104,106,107] based on the software presented in [148] have shown that only thes TT coefficients produce a non-negligible effect on the light propagation (while it has impact only at the next post-Newtonian level on the orbital dynamics [29,39]). On the other hand, the other coefficients produce non-negligible effects on the orbital dynamics [39] and can therefore be constrained using planetary ephemerides data. In the linearized gravity limit, the contribution from SME to the 2-body equations of motion within the gravitational sector of SME are given by the first line of Eq. (17) (i.e. for a vanishing V k ). The coefficients TT is completely unobservable in this context since absorbed in a rescaling of the gravitational constant (see the discussion in [39,52]).
Ideally, in order to perform a solid estimation of the SME coefficients using planetary ephemerides, one should include the full SME equations in the integration of the planets motion and fit them simultaneously with the other estimated parameters (positions and velocities of planets, J 2 of the Sun, . . . ). This solid analysis within the SME formalism has not been performed so far.
As a first step, a postfit analysis has been performed [49,50]. The idea of this analysis is to derive the analytical expression for the secular evolution of the orbital elements produced by the SME contribution to the equations of motion. Using the Gauss equations, secular perturbations induced by SME on the orbital elements have been computed in [39] (see also [49] for a similar calculations done for thes T J coefficients only). In particular, the secular evolution of the longitude of the ascending node Ω and the argument of the perihelion ω is given by where a is the semimajor axis, e the eccentricity, i the orbit inclination (with respect to the ecliptic), n = (G N m /a 3 ) 1/2 is the mean motion, ε = 1 − (1 − e 2 ) 1/2 , δm the difference between the two masses and M their sum (in the cases of planets orbiting the Sun, one has M ≈ δm). In all these expressions, the coefficients for Lorentz violation with subscripts P, Q, and k are understood to be appropriate projections ofs µν along the unit vectors P, Q, and k, respectively. For example,s k = k isTi , s PP = P i P jsij . The unit vectors P, Q and k define the orbital plane (see [39] or Eq. (8) from [50]). Instead of including the SME equations of motion in planetary ephemerides, the postfit analysis uses estimations of supplementary advances of perihelia and nodes derived from ephemerides analysis [135,140,144] to fit the SME coefficients through Eq. (19). In [50], estimations of supplementary advances of perihelia and longitude of nodes from INPOP (see Table 5 from [135]) are used to fit a posteriori the SME coefficients. This analysis suffers from large correlations due to the fact that the planetary orbits are very similar to each other: nearly eccentric orbit and very low inclination orbital planes. In order to deal properly with these correlations a Bayesian Monte Carlo inference has been used [50]. The posterior probability distribution function can be found on Figure 1 from [50]. The intervals corresponding to the 68% Bayesian confidence levels are given in Table 4 as well as the correlation matrix. It is interesting to mention that a decomposition of the normal matrix in eigenvectors allows one to find linear combinations of SME coefficients that are uncorrelated with the planetary ephemerides analysis (see Eq. (15) and Table IV from [50]). Table 4. Estimations of the SME coefficients from a postfit data analysis based on planetary ephemerides from [50]. The uncertainties correspond to the 68% Bayesian confidence levels of the marginal posterior probability distribution function. The associated correlation coefficients can be found in Table III from [50]. In summary, planetary ephemerides offer a great opportunity to constrain hypothetical violations of Lorentz symmetry. So far, only postfit estimations of the SME coefficients have been performed [49,50]. In this analysis, estimations of secular advances of perihelia and longitude of nodes obtained with the INPOP planetary ephemerides [135] are used to fit a posteriori the SME coefficients using the Eqs. (19). The 68% marginalized confidence intervals are given in Table 4. This analysis suffers highly from correlations due to the fact that the planetary orbits are very similar. A very interesting improvement regarding this analysis would be to perform a full analysis by integrating the planetary equations of motion directly within the SME framework and by fitting the SME coefficients simultaneously with the other parameters fitted during the ephemerides data reduction.

Gravity Probe B
In GR, a gyroscope in orbit around a central body undergoes two relativistic precessions with respect to a distant inertial frame: (i) a geodetic drift in the orbital plane due to the motion of the gyroscope in the curved spacetime [149] and (ii) a frame-dragging due to the spin of the central body [150]. In GR, the spin of a gyroscope is parallel transported, which at the post-Newtonian approximation gives the relativistic drift whereŜ is the unit vector pointing in the direction of the spin S of the gyroscope, r and v are the position and velocity of the gyroscope,r = r/r and J is the angular momentum of the central body.
In 1960, it has been suggested to use these two effects to perform a new test of GR [151,152]. In April 2004, GPB, a satellite carrying 4 cryogenic gyroscopes was launched in order to measure these two precessions. GPB was orbiting Earth on a polar orbit such that the two relativistic drifts are orthogonal to each other [153]: the geodetic effect is directed along the NS direction (North-South, i.e. parallel to the satellite motion) while the frame-dragging effect is directed on the WE direction (West-East, see [52,153] for further details about the axes conventions in the GPB data reduction Within the SME framework, if one considers only thes µν coefficients, the equation of parallel transport in term of the spacetime metric is not modified (see Eq. (143) from [39]). Nevertheless, the expression of the spacetime metric is modified leading to a modification of the relativistic drift given by Eq. (150) from [39]. In order to focus only on the dominant secular part of the evolution of the spin orientation, the relativistic drift equation has been averaged over a period. The SME contribution to the precession can be written as [39] where G N is the effective gravitational constant defined by Eq. (5), the coefficientsĩ are defined bỹ i (β) = 1 + βI ⊕ /(M ⊕ r 2 ),σ J is a unit vector normal to the gyroscope orbital plane, r and v are the norm of the position and velocity of the gyroscope ands JK t is the traceless part ofs JK as defined by Eq. (6). Using the geometry of GPB into the last equation and using Eq. (20a), one finds that the gyroscope anomalous drift is given by where the units are mas/yr. These are the SME modifications to the relativistic drift arising from the modification of the equations of evolution of the gyroscope axis (i.e. modification of the parallel transport equation due to the modification of the underlying spacetime metric). In addition to modifying the evolution of the spin axis, a breaking of Lorentz symmetry will impact the orbital motion of the gyroscope. As a result, the position and velocity of the gyroscope will depend on the SME coefficients and therefore also impact the evolution of the spin axis through the GR contribution given by Eq. (20b). The best way to deal with this effect is to use the GPB tracking measurements (GPS) in order to constrain the gyroscope orbital motion and eventually constrain the SME coefficients through the equations of motion. In [52], these tracking observations are not used and only the gyroscope drift is used in order to constrain the SME contributions coming from both the modification of the parallel transport and from the modification of GPB orbital motion. In order to do this, the contribution of SME on the evolution of the orbital elements given by Eq. (19) and (26) are used, averaged over a period and in the low eccentricity approximation. This secular evolution for the osculating elements is introduced in the relativistic drift equation for the gyroscope from Eq. (20b) and averaged over the measurement time using Eq. (20a). Using the GPB geometry, this contribution to the relativistic drift is given by ∆R NS = 5.7 × 10 6 (s XX −s YY ) + 1.7 × 10 7sXY − 1.9 × 10 7sXZ − 6.6 × 10 7sYZ , (24a) with units of mas/yr. The sum of the two SME contributions to the gyroscope relativistic drift given by Eq. (23) and (24) can be compared to the GPB estimations given by Eq. (21). The result is given in Table 5. The main advantage of GPB comes from the fact that it is sensitive to thes TT coefficient. The constraint on this coefficient is at the level of 10 −3 , a little bit less good than the one obtained with VLBI or with binary pulsars but relying on a totally different type of observations. The constraints on the spatial part of the SME coefficients (s I J ) are at the level of 10 −7 and are superseded by the other measurements. The constraints on these coefficients come mainly from the contribution arising from the orbital dynamics of GPB and not from a direct modification of the spin evolution. Constraining the orbital motion from GPB by using the gyroscope observations only is not optimal and tracking observations may help to improve the corresponding constraints (in this case, a dedicated satellite may be more appropriate as discussed in Section 5.3). Table 5. Estimations of the SME coefficients from a postfit data analysis based on GPB [52].
In summary, the GPB measurement of a gyroscope relativistic drifts due to geodetic precession or frame-dragging can be used to search for a breaking of Lorentz symmetry. The main advantage of this technique comes from its sensitivity tos TT . As already mentioned, this coefficient has an isotropic impact on the propagation velocity of gravitational waves as can be noticed from Eq. (27) below (see also Eq. (9) from [56] or Eq. (11) from [60]). A preliminary result based on a post-fit analysis performed after a GR data reduction of GPB measurements gives a constraint ons TT at the level of 10 −3 [52]. This should be investigated further since the Earth's quadrupole moment has been neglected and Lorentz-violating effects on the aberration terms can also change slightly the results. In addition, impacts from Lorentz violations on frame-dragging arising in other contexts such as satellite laser ranging (see Section 5.3) or signals from accretion disks around collapsed stars [154] would also be interesting to consider.

Binary pulsars
The discovery of the first binary pulsars PSR 1913+16 by Hulse and Taylor in 1975 [155] has opened a new window to test the theory of gravitation. Observations of this pulsar have allowed one to measure the relativistic advance of the periastron [156] and more importantly to measure the rate of orbital decay due to gravitational radiation [157]. Pulsars are rotating neutron stars that are emitting very strong radiation. The periods of pulsars are very stable which allows us to consider them as "clocks" that are moving in an external gravitational field (typically in the gravitational field generated by a companion). The measurements of the pulse time of arrivals can be used to infer several parameters by fitting an appropriate timing model (see for example Section 6.1 from [2]): (i) non-orbital parameter such as the pulsar period and its rate of change; (ii) five Keplerian parameters and (iii) some post-Keplerian parameters [158]. In GR, the expressions of these post-Keplerian parameters are related to the masses of the two bodies and to the Keplerian parameters. If more than 2 of these post-Keplerian parameters can be determined, they can be used to test GR [159]. Nowadays, more than 70 binary pulsars have been observed [160]. A description of the most interesting binary pulsars in order to test the gravitation theory can be found in Section 6.2 from [2] or in the supplemental material from [53].
The model fitted to the observations is based on a post-Newtonian analytical solution to the 2 body equations of motion [161] (see also [162]) and includes contribution from the Einstein time delay (i.e. the transformation between proper and coordinate time), the Shapiro time delay, the Roemer time delay [158]. The model also corrects for several systematics like atmospheric delay, Solar system dispersion, interstellar dispersion, motion of the Earth and the Solar System, . . . (see for example [163]).
Pulsars observations provide some of the best current constraints on alternative theories of gravitation (for a review, see [164,165]). In addition to the Hulse and Taylor pulsar, the double pulsar [166] now provides the best measurement of the pulsar orbital rate of change [165]. In addition, the post-Keplerian modeling has been fully derived in tensor-scalar theories [167][168][169] such that pulsars observations have provided some of the best constraints on this class of theory [165,170,171]. It is important to mention that non perturbative strong field effects may arise in binary pulsars system and needs to be taken into account [169,172].
In addition, binary pulsars have also been successfully used to test Lorentz symmetry. For example, analyses of the pulses time of arrivals provide a constraint on the α 1,2,3 PPN parameters. Since non perturbative strong field effects may arise in binary pulsars system (see for example [173] for strong field effects in Einstein-Aether theory), the obtained constraints are interpreted as strong field version of the PPN parameters denoted byα i . Estimates of these parameters should be compared carefully to the standard weak field constraints since they may depend on the gravitational binding energy of the neutron star. The best current constraint onα 1 = −0.4 +3.7 −3.1 × 10 −5 is obtained by considering the orbital dynamics of the binary pulsars PSR J1738+0333 [174,175]. The best current constraint onα 2 takes advantage from the fact that this parameter produces a precession of the spin axis of massive bodies [176]. The combination of observations of two solitary pulsars lead to the best current constraints on |α 2 | < 1.6 × 10 −9 [177]. Finally, the parameterα 3 produces a violation of the momentum conservation in addition to a violation of the Lorentz symmetry. This parameter will induce a self-acceleration for rotating body that can be constrained using binary pulsars [178]. The best current constraint uses a set of 5 pulsars (4 binary pulsars and one solitary pulsar) and is given byα 3 < 5.5 × 10 −20 [179].
Furthermore, specific Lorentz violating theories have also been constrained with binary pulsars. In [72,73], binary pulsars observations are used to constrain Einstein-Aether and khronometric theory. In these theories, the low-energy limit Lorentz violations can be parametrized by four parameters: the α 1 and α 2 PPN parameters and two other parameters. It has been shown [72,73,173] that the orbital period decay depends on these four parameters. Assuming the solar system constraints on α 1 and α 2 [2], measurements of the rate of change of the orbital period of binary pulsars have been used to constrain the two other parameters (see for example Fig. 2 from [72]). In this work, strong field effects have been taken into account by solving numerically the field equations in order to determine the neutron stars sensitivity [73].
Finally, binary pulsars have been used in order to derive constraints on the SME coefficients. As in the PPN formalism, constraints obtained from binary pulsars need to be considered as constraints on strong-field version of the SME coefficients that may include non perturbative effects. Two different types of effects have been used to determine estimates on the SME coefficients: (i) tests using the spin precession of solitary pulsars and (ii) tests using effects on the orbital dynamics of binary pulsars [53]. The SME contribution to the precession rate of an isolated spinning body has been derived in [39] and is given by where P is the spin period andŜ j is the unit vector pointing along the spin direction. The effects from the pulsar spin precession on the pulse width can be found in [177,180]. Two solitary pulsars have been used to constrain the SME coefficients with this effect. The second type of tests come from the orbital dynamics of binary pulsars. As mentioned in Sections 4.3 and 4.4, the SME will modify the two-body equations of motion by including the term from Eq. (17). At first order in the SME coefficients, this will produce several secular effects that have been computed in [39]. In particular, an additional advance in the argument of periastron and of the longitude of the nodes has been mentioned in Eq. (19) and used to constrain the SME with planetary ephemerides. For binary pulsars, it is possible to constrain a secular evolution of two other orbital elements: the eccentricity and the projected semi-major axis x. The secular SME contributions to these quantities have been computed in [39,53,54] and are given by where m P is the mass of the pulsar, m C is the mass of the companion and all others quantities have been introduced after Eqs. (19). For each binary pulsar, in principle 3 tests can be constructed by usingω,ė,ẋ. In [53], 13 pulsars have been used to derive estimates on the SME coefficients. The combination of the observations from the solitary pulsars and from the 13 binary pulsars are reported in Table 6. Both orbital dynamics and spin precession are completely independent ofs TT whose constraint will be discussed later. Several comments can be made about this analysis. First of all, it can be considered as a postfit analysis done after an initial fit performed in GR (or within the post-Keplerian formalism). In particular, correlations between the SME coefficients and other parameters (e.g. orbital parameters) are neglected. Secondly, for most of the pulsars,ẋω andė are not directly measured from the pulse time of arrivals but rather estimated from the uncertainties on x, ω and e divided by the time span of the observations. Further, it is important to mention that effects of Lorentz violations have been considered only for the orbital dynamics but never on the Einstein delay or on the Shapiro time delay in this analysis. The full timing model within SME can be found in Section V.E.3 from [39] (see also [181] for a similar derivation with the matter-gravity couplings). In addition, some parameters are not measured like for example the longitude of the ascending node Ω or the azimuthal angle of the spin. These parameters have been marginalized by using Monte Carlo simulations. It is unclear what type of prior probability distribution function has been used in this analysis and what is the impact of this choice. Nevertheless, the results obtained by this analysis (which does not include thes TT parameter) are amongst the best ones currently available demonstrating the power of pulsars observations. The main advantages of using binary pulsars come from the fact that their orbital orientation vary which allows one to disentangle the different SME coefficients and to end up with low correlations. Furthermore, they are so far the only constraints on the strong field version of the SME coefficients. Table 6. Estimation of SME coefficients from binary pulsars data analysis from [53,54]. No correlations coefficients have been derived in this analysis. These estimates should be considered as estimates on the strong field version of the SME coefficients that may include non perturbative strong field effects due to the gravitational binding energy.

Coefficient
In addition, a different analysis has been performed to constrain the parameters TT alone [54]. While the orbital dynamics and the spin precession is completely independent ofs 00 (i.e. the time component ofs µν in a local frame), the boost between the Solar System and the binary pulsar frame makes appear explicitly thes TT coefficient. In [54], the assumption that there exists a preferred frame where thes µν tensor is isotropic is made, which makes the results specific to that case (although the analysis can be done without this assumption). The analysis requires the knowledge of the pulsar velocity with respect to the preferred frame as well as the velocity of the Solar System with respect to the same frame. Three pulsars have their radial velocity measured, which combined with proper motion in the sky can be used to determine their velocity. The velocity of the Solar System is taken as its velocity with respect to the Cosmic Microwave Background (CMB) frame w (with |w | = 369 km/s). The analysis is completely similar to the ones performed for the other SME coefficients (see the discussion in the previous paragraph). It is known thats TT has a strong effect on the propagation of the light neglected in [54], which may impact the result. In addition, all correlations betweens TT and the other SME coefficients are neglected. Finally, two different scenarios have been considered regarding the preferred frame: (i) a scenario where the preferred frame is assumed to be the CMB frame and (ii) a scenario where the orientation of the preferred frame is left free and is marginalized over but the magnitude of the velocity of the Solar System with respect to that frame is still assumed to be the 369km/s. The general case corresponding to a completely free preferred frame has not been considered. If the CMB frame is assumed to be the preferred frame, the constraint ons TT is given by s TT < 1.6 × 10 −5 which is a bit better than the one obtained with VLBI (see Eq. (13)) although the VLBI analysis does not assume any preferred frame. The scenario where the orientation of the preferred frame is left as a free parameter leads to an upper bound on s TT < 2.8 × 10 −4 .
In summary, observations of binary pulsars are an incredible tool to test the gravitation theory. These tests are of the same order of magnitude (and sometimes better) than the ones performed in the Solar System. Moreover, observations of binary pulsars are sensitive to strong field effects. Observations of the pulse arrival times have been used to search for a breaking of Lorentz violation within the PPN framework by constraining the strong field version of the α i parameters. The parameterα 1 is constrained at the level of 10 −5 ,α 2 at the level of 10 −9 andα 3 at the level of 10 −20 [164]. In addition, constraints on Einstein-Aether and khronometric theory have also been done by combining Solar System constraints with binary pulsars observations [72,73]. Finally, within the SME framework, a postfit analysis has been done by considering the spin precession of solitary pulsars and the orbital dynamics of binary pulsars. The obtained results are given in Table 6 and constrain the strong field version of the SME coefficients. The main advantage of using binary pulsars comes from the fact that they proved an estimate of all the SME coefficients with reasonable correlations. It has to be noted that the modification of the orbital period due to gravitational waves emission has not been computed so far in the SME formalism. In addition, the constraint ons TT suffers from the assumption of the existence of a preferred frame. Moreover, the corresponding analysis has neglected all effects on the timing delay that may also impact the results and has neglected the other SME coefficients that may also impact this constraint.

4.7.Čerenkov radiation
GravitationalČerenkov radiation is an effect that occurs when the velocity of a particle exceeds the phase velocity of gravity. In this case, the particle will emit gravitational radiation until the particle loses enough energy to drop below the gravity speed [56]. In modified theory of gravity, the speed of gravity in a vacuum may be different from the speed of light andČerenkov radiation may occur and produces energy losses for particles traveling over long distances. Observations of high energy cosmic rays that have not lost all their energy throughČerenkov radiation can be used to put constraints on models of gravitation that predicts gravitational waves that are propagating slower than light. This effect has been used to constrain some alternative gravitation theories [182,183]: a class of tensor-vector theories [184], a class of tensor-scalar theories [185], extended theories of gravitation [186] and some ghost-free bigravity [187].
The propagation of gravitational waves within the SME framework has been derived in [56,60] (including nonminimal SME contributions). In particular, in the minimal SME, the dispersion relation for the gravitational waves is given by [56] where l α is the 4-momentum of the gravitational wave. A similar expression including nonminimal higher order SME terms can be found in [56,60]. If the minimal SME produces dispersion-free propagation, the higher order terms lead to dispersion and birefringence [60]. As can be directly inferred from the last equation, gravitationalČerenkov radiation can arise when the effective refractive index n is wherel µ = l µ / |l|. The expression for the energy loss rate due to Lorentz-violating gravitationaľ Cerenkov emission has been calculated from tree-level graviton emission for photons, fermions and scalar particles and is given by [56] where d is the dimension of the Lorentz violating operator (d = 4 for the minimal SME), F w (d) is a dimensionless factor depending on the flavor w of the particle emitting the radiation, p is the particle incoming momentum (withp = p/ |p|) ands (d) is a direction-dependent combination of SME coefficients. In the minimal SME,s (4) (p) is decomposed on spherical harmonics as where we explicitly indicated the (SH) to specify that these coefficients are spherical harmonic decomposition of the SME coefficients. The calculation of the dimensionless factor F w (d) for scalar particles, fermions and photons has been done in [56]. The integration of Eq. (29) shows that if a cosmic ray of specie w is observed on Earth with an energy E f after traveling a distance L along the directionp, this implies the following constraint on the SME coefficients where F w (d) = (2d − 5)/F w (d) is another dimensionless factor dependent on the matrix element of the tree-level process for graviton emission.
Using data for the energies and angular positions of 299 observed cosmic rays from different collaborations [188][189][190][191][192][193][194][195], Kostelecký and Tasson [56] derived lower and upper constraints on 80 SME coefficients, including the nine coefficients from the minimal SME whose constraints are given by the Table 7. In their analysis, they consider the coefficients from the different dimensions separately and did not fit all of them simultaneously. In addition, in the minimal SME, they did a fit for thē s TT parameter alone and another fit for the other 8 coefficients. The number of sources and their directional dependence across the sky allow one to disentangle the SME coefficients and to derive two-sided bounds from the Eq. (31). The only coefficient that is one sided iss TT because it produces isotropic effects. The bounds are severe for these coefficients, on the order of 10 −13 . However, this analysis assumes that the matter sector coefficients vanish. Furthermore, several assumptions have been made in order to derive the bounds from Table 7. It is assumed that the cosmic ray primaries are nuclei of atomic weight N = 56 (iron), that theČerenkov radiation is emitted by one of the fermionic partons in the nucleus that carries 10 % of the cosmic ray energy and that the travel distance of the cosmic ray is 10 Megaparsec (Mpc) [56]. Although only conservative assumptions are used for the astrophysical processes involved in the production of high-energy cosmic rays, the observations rely on the sources on the order of 10 Mpc distant, and thus the analysis is of a different nature than a controlled laboratory or even Solar-System test. Table 7. Lower and upper limits on the SME coefficients decomposed in spherical harmonics derived fromČerenkov radiation [56]. For the sake of completeness and to allow an easy comparison with the estimations of the other standard cartesians µν coefficients, the following relations give the links between the spherical harmonic decomposition and the standard cartesian decomposition of the SME coefficients:

Coefficient Lower bound Upper bound
In summary, observations of cosmic rays allow one to derive some stringent boundaries on the SME coefficients. The idea is that if Lorentz symmetry is broken, these high energy cosmic rays would have lost energy by emittingČerenkov radiation that has not been observed. The boundaries on the spherical harmonic decomposition of the SME coefficients are given in the Table 7 (in order to compare these boundaries to other constraints, they have been transformed into boundaries on standard cartesian SME coefficients in Table 8). For the minimal SME, one can limit the isotropics TT (one sided bound) or the other eight other coefficients ins µν , but not all the nine simultaneously. These boundaries are currently the best available in the literature at the exception ofs TT whose constraint is only one sided. Nevertheless, several assumptions have been made in this analysis and the observations rely on sources located at very high distances. This analysis is therefore of a different nature than the other ones where more control on the measurements is possible.

Summary and combined analysis
To summarize, several measurements have already successfully been used to constrain the minimal SME in the gravitational sector (i.e. thes µν coefficients): • Atom interferometry [47,48].
A detailed description of all these analyses is provided in the previous subsections and the Table 8 summarizes the current estimates. It is also interesting to combine all these estimations together to provide the best estimates on the SME coefficients. In order to do this, we perform a large least-square fit including all the results from the Table 8 including the covariance matrices quoted in the previous subsections. The results from theČerenkov radiation are not included since they rely on a very different type of observations. Two combined fits are presented: one without including the pulsars results and one including the pulsars results. This is due to the fact that pulsars are sensitive to a strong version of the SME coefficients that may include non perturbative strong field effects as described in Section 4.6. If this is the case, then the pulsars results cannot be directly combined with the weak gravitational field estimates on the SME coefficients. If no non perturbative strong field effect arises, then the right column from Table 9 presents a combined fit that includes these observations as well. The results from Table 9 include all the information currently available in the literature on thes µν (estimations and correlation matrices). It can also be noted that the pulsars results improve significantly the marginalized estimations ons TY ands TZ by reducing strongly the correlation between these two coefficients.
A review of these measurements can be found in [30]. GPB [52] s TT Table 9. Estimation of SME coefficients resulting from a fit combining results from: atomic gravimetry (see Table 1), VLBI (see Eq. 13), LLR (see Table 3), planetary ephemerides (see Table 4), Gravity Probe B (see Table 5). The correlation matrices from all these analyses have been used in the combined fit. The right column includes the pulsars results from Table 6 as well. The three estimates ons J J are obtained by using the traceless conditions TT =s XX +s YY +s ZZ .

The future
In addition to all the improvements related to existing analysis suggested in the previous sections, there are a couple of sensitivity analyses that have been done within the SME framework. First of all, a thorough and detailed analysis of a lot of observables related to gravitation can be found in [39]. In addition, we will present in the next subsections a couple of analyses and ideas that may improve the SME coefficients estimates in the future.

The Gaia mission
Launched in December 2013, the ESA Gaia mission [196] is scanning regularly the whole celestial sphere once every 6 months providing high precision astrometric data for a huge number (≈ 1 billion) of celestial bodies. In addition to stars, it is also observing Solar System Objects (SSO), in particular asteroids. The high precision astrometry (at sub-mas level) will allow us to perform competitive tests of gravitation and to provide new constraints on alternative theories of gravitation.
First of all, the Gaia mission is expected to provide an estimate of the γ PPN parameter at the level of 10 −6 [197] by measuring the deflection of the light on a 5 years timescale. Furthermore, in addition to this global determination of a global PPN parameter from observations of light deflection, it has been proposed to use Gaia observations to map the deflection angle in the sky and to look for a dependence of the γ PPN parameter with respect to the Sun impact parameter [198][199][200][201][202]. Such a dependence of the gravitational deflection with respect to the observation geometry is also a feature predicted by SME as shown in [82]. Therefore, the global mapping of the light deflection with Gaia can also be efficiently used to constrain some SME coefficients. A first sensitivity analysis can be found in [82] and is reported on Table 10. Note that proposals observations and missions like AGP [203] or LATOR [204] can in the long term improve these estimates further by improving the light deflection measurement. Table 10. Sensitivity of the SME coefficients to the measurement of the light deflection by several space missions or proposals (these estimates are based on Table I from [82]).

Missions TTsT JsI J
Gaia [196] 10 −6 10 −6 10 −5 AGP [203] 10 −7 10 −7 10 −6 LATOR [204] 10 −8 10 −8 10 −7 In addition to gravitation tests performed by measuring the light deflection, Gaia also provides a unique opportunity to test gravitation by considering the orbital dynamics of SSO. One can estimate that about 360 000 asteroids will be regularly observed by Gaia at the sub-mas level, which will allow us to perform various valuable tests of gravitation [205,206]. In particular, realistic simulations of more than 250 000 asteroids have shown that Gaia will be able to constrain the β PPN parameter at the level of 10 −3 [205]. The main advantage from Gaia is related to the huge number of bodies that will be observed with very different orbital parameters as illustrated on Figure 1. As a consequence, the huge correlations appearing in the planetary ephemerides analysis (see Section 4.4) will not appear when considering asteroids observations and the marginalized confidence intervals will be highly improved compared to planetary ephemerides analysis. A realistic sensitivity analysis of Gaia SSOs observations within the SME framework has been performed (see also [206] for preliminary results). In this analysis, 360 000 asteroids have been considered over the nominal mission duration (i.e. five years) and a match between the SSO trajectories with the Gaia scanning law is performed to find the observation times for each SSO. Simultaneously with the equations of motion, we integrate the variational equations, the simulated SSO trajectories being transformed into astrometric observables as well as their partial derivatives with respect to the parameters considered in the covariance analysis. The covariance analysis leads to the estimated uncertainties presented in Table 11. These uncertainties are incredibly good, which is due to the variety of the asteroids orbital parameters as discussed above. Using our set of asteroids, the correlation matrix for the SME parameters is very reasonable: the most important correlation coefficients are 0.71, -0.68 and 0.46. All the other correlations are below 0.3. Therefore, Gaia offers a unique opportunity to constrain Lorentz violation through the SME formalism. Finally, the Gaia mission is likely to be extended to 10 years, therefore doubling the measurements baseline which will also impact significantly the expected uncertainties. Finally, it is worth mentioning that the Gaia dataset can be combined with radar observations [207] that are complementary in the time frame and orthogonal to astrometric telescopic observations. In summary, the Gaia space mission offers two opportunities to test Lorentz symmetry in the Solar System by looking at the deflection of light and by considering the orbital dynamics of SSO. The second type of observations is extremely interesting in the sense that the high number and the variety of orbital parameters of the observed SSO leads to decorrelate the SME coefficients.

Analysis of Cassini conjunction data
The space mission Cassini is exploring the Saturnian system since July 2004. During its cruising phase while the spacecraft was on its interplanetary journey between Jupiter and Saturn, a measurement of the gravitational time delay was performed [118]. This measurement occurred during a Solar conjunction in June 2002 and was made possible thanks to a multi-frequency radioscience link (at X and Ka-band) which allows a cancellation of the solar plasma noise [118]. The related data spans over 30 days and has been analyzed in the PPN framework leading to the best estimation of the γ PPN parameter so far given by (2.1 ± 2.3) × 10 −5 [118].
The exact same set of data can be reduced within the SME framework and is expected to improve our currents TT estimation. The time delay within the SME framework has been derived in [80] and is given by Eq. (9).
A simulation of the Cassini link during the 2002 conjunction within the full SME framework has been realized using the software presented in [148] (see also [208,209]). The signature produced by thes TT coefficients on the 2-way Doppler link during the Solar conjunction is illustrated on Figure 2. In [80], a crude estimate of attainable sensitivities in estimate of the SME coefficients using the Cassini conjunction data is given (see Table I from [80]). It is shown that some combinations of thes I J coefficients can only be constrained at the level of 10 −4 , which is 7 to 8 orders of magnitude worse than the current best constraints on these coefficients. It is therefore safe to neglect these and to concentrate only on thes TT coefficient. A realistic covariance analysis performed over the 30 days of the Solar conjunction and assuming an uncertainty of the Cassini Doppler of 3 µm/s [210,211] shows that thes TT parameter can be constrained at the level of 2 × 10 −5 using the Cassini data allowing an improvement of a factor 4 with respect to the current best estimate coming from VLBI analysis (see Eq. (13)). Therefore, a reanalysis of the 2002 Cassini data within the SME framework would be highly valuable.

Satellite Laser Ranging (LAGEOS/LARES)
Searching for violations of Lorentz symmetry by using the orbital motion of planets (see Section 4.4), binary pulsars (see Section 4.6), the Moon (see Section 4.3) and asteroids (see Section 5.1) has turned out to be highly powerful. It is therefore logical to consider the motion of artificial satellite orbiting around Earth to search for Lorentz violations. In particular, laser ranging to the two LAGEOS and to the LARES satellites has successfully been used to test GR by measuring the impact of the Schwarzschild precession on the motion of the satellites [212][213][214]. It has also been claimed that the impact of the frame-dragging (or Lense-Thirring effect) due to the Earth's spin on the orbital motion of the satellites has been measured [215][216][217][218][219][220] although this claim remains controversial [221][222][223][224][225][226]. Similarly, the LAGEOS/LARES satellites can also be used to search for Lorentz violations. A sensitivity analysis has been done in [49] and it has been shown that the LAGEOS satellites are sensitive at the level of 10 −4 to thes T J coefficients. Using LARES should improve significantly this value. Further numerical simulations are required in order to determine exactly the SME linear combinations to which the ranging to these satellites is sensitive to. A data analysis within the full SME framework (i.e. including the integration of the SME equation of motion and including the SME coefficients with the other global parameters in the fit) would also be highly interesting. In addition, similar tests of Lorentz symmetry can also be included within the scientific goals of the LAser RAnged Satellites Experiment (LARASE) project [227] or within the OPTIS project [228].

Gravity-Matter coefficients and breaking of the Einstein equivalence principle
All the measurements mentioned in Section 4 can be analyzed by considering the gravity-matter coupling coefficientsā µ andc µν [29] that are breaking the EEP. Some atomic clocks measurements have already provided some constraints on theā µ coefficients [229][230][231]. In addition, in [50] the planetary ephemerides analysis is interpreted by considering theā µ coefficients and the atomic interferometry results from [48] and the LLR results from [45] are also reinterpreted by considering the gravity-coupling coefficients. Clearly this is a preliminary analysis that needs to be refined by more solid data reductions. Considering the increasing number of fitted parameters, it is of prime importance to increase the number of measurements used in the analysis and to produce combined analysis with as many types of observations as possible. The measurements developed in Section 4 are a first step in order to reach this goal. The gravity-coupling coefficients can also be constrained by more specific tests related to the EEP like for example tests of the Universality of Free Fall with MicroSCOPE [83,84], tests of the gravitational redshift with GNSS satellites [232], with the Atomic Clocks Ensemble in Space (ACES) project [233], or with the OPTIS project [228], . . .

Conclusions
Lorentz symmetry is at the heart of both GR and the Standard Model of particle physics. This symmetry is broken in various scenarios of unification, of quantum gravity and even in some models of Dark Matter and Dark Energy. Searching for violations of Lorentz symmetry is therefore a powerful tool to test fundamental physics. The last decades have seen the number of tests of Lorentz invariance arise dramatically in all sectors of physics [24]. In this review, we focused on searches for Lorentz symmetry breaking in the pure gravitational sector. Mainly two frameworks exist to parametrize violations of Lorentz invariance in the gravitation sector. First of all, the three α 1,2,3 PPN parameters phenomenologically encode a violation of Lorentz symmetry at the level of the spacetime metric [2]. These parameters are constrained by LLR (see Section 4.3) and by pulsars timing measurements (see Section 4.6). In addition, it is interesting to notice that the corresponding PPN metric parametrizes also Einstein-Aether and Khronometric theories in the weak gravitational field limit [72] while these theories have a more complex strong field limit (and can show non perturbative effects) that have been constrained by pulsars observations (see Section 4.6 and [72,73]).
In addition to the PPN formalism, the SME formalism has been developed by including systematically all possible Lorentz violations terms that can be constructed at the level of the action. In the pure gravitational sector, the gravitational action within the SME formalism contains the usual Einstein-Hilbert action but also new Lorentz violating terms constructed by contracting new fields with some operators built from curvature tensors and covariant derivatives with increasing mass dimension [36]. The lower mass dimension term is known as the minimal SME. In the limit of linearized gravity, the observations within the minimal SME formalism depend on 9 coefficients, thes µν symmetric traceless tensor. This formalism offers a new opportunity to search for deviations from GR in a framework different from the standard PPN formalism. We reviewed the different observations that have been used so far to constrain the SME coefficients. The main idea is to search for a signature (usually periodic) that arises from a dependence on the orientation of the system measured (the dependence on the orientation is typically due to the Earth's rotation, the orbital motion of the planets around the Sun, etc. . . ) or from a dependence on the boost of the system observed (so far, only the binary pulsarss TT constraint comes from this type of dependence [54]). Most of SME analyses are postfit analyses in the sense that analytical signatures due to SME are fitted in residuals noise obtained in a previous data reduction performed in pure GR. In Section 3, we showed that this approach can sometimes lead to overoptimistic constraint on the SME coefficients and that one should be careful in interpreting results obtained using such an approach.
In Section 4, we discussed in details the different measurements used so far to constrain thē s µν coefficients: atomic gravimetry (Section 4.1), VLBI (Section 4.2), LLR (Section 4.3), planetary ephemerides (Section 4.4), Gravity Probe B (Section 4.5), pulsars timing (Section 4.6) andČerenkov radiation (Section 4.7). In each of these subsections, we describe the current analyses performed in order to constrain the SME coefficients and provide a critical discussion from each of them. We also provide a summary of these constraints on Table 8. In addition, we used all these results to produce a combined analysis of the SME coefficients. This fit is done by taking into account the correlation matrices for each individual analysis. The results of this combined fit are presented in Table 9 and are the current best estimates of the SME coefficients that are possible to derive with all available analyses. In addition to the minimal SME, there exists higher order Lorentz violating terms that have been considered and constrained by short-range gravity experiments [57][58][59], gravitational waves analysis [60] andČerenkov radiation [30,56].
In Section 5, we discussed some opportunities to improve the current constraints on the SME coefficients. In particular, the European space mission Gaia offers an excellent opportunity to probe Lorentz symmetry through the measurement of light deflection and through the orbital motion of asteroids. The Cassini conjunction data also offers a way to constrain thes TT coefficient that impacts severely the propagation of light. Finally, existing satellite laser ranging data can also be analyzed within the SME framework.
In addition, as mentioned in Section 5.4, all the analyses presented in this review can include gravity-matter coefficients [29]. While considering these, the number of coefficients fitted increase significantly and it becomes crucial to produce a fit combining several kinds of experiments. A preliminary analysis considering these coefficients for planetary ephemerides, LLR and atomic gravimetry has been performed in [50] but needs to be refined. In addition, some atomic clocks experiments have already been used to constrain matter-gravity coefficients [229][230][231].
In conclusion, though no violation of Lorentz symmetry has been observed so far, an incredible number of opportunities still exists for additional investigations. There remains a large area of unexplored coefficients space that can be explored by improved measurements or by new projects aiming at searching for breaking of Lorentz symmetry. In addition, the increasing number of parameters fitted (by including the gravity-matter coupling coefficients simultaneously with the pure gravity coefficients in the analyses) will deter the marginalized estimates of each coefficient. This verdict emphasizes the need to increase the types of measurements that can be combined together to explore the vast parameters space as efficiently as possible. The current theoretical questions related to the quest for a unifying theory or for a quantum theory of gravitation suggests that Lorentz symmetry will play an important role in the search for new physics. Hopefully, future searches for Lorentz symmetry breaking will help theoreticians to unveil some of the mysteries about Planck-scale physics [22].