The new lunar ephemeris INPOP17a and its application to fundamental physics

Abstract

We present here the new INPOP lunar ephemeris, INPOP17a. This ephemeris is obtained through the numerical integration of the equations of motion and of rotation of the Moon, fitted over 48 yr of lunar laser ranging (LLR) data. We also include the 2 yr of infrared LLR data acquired at the Grasse station between 2015 and 2017. Tests of the universality of free-fall are performed. We find no violation of the principle of equivalence at the (−3.8 ± 7.1) × 10⁻¹⁴ level. A new interpretation in the frame of dilaton theories is also proposed.

gravitation, ephemerides, Moon

1 INTRODUCTION

The Earth–Moon system is an ideal tool for carrying out tests of general relativity (GR) and more particularly the test of the universality of free-fall (UFF; Nordtvedt 1968a; Anderson et al. 1996). Since 1969, the lunar laser ranging (LLR) observations are obtained on a regular basis by a network of laser ranging stations Faller et al. (1969); Bender et al. (1973), and currently with a millimetre-level accuracy Samain et al. (1998); Murphy (2013). Thanks to this level of accuracy at the Solar system scale, the principle of the UFF can in theory be tested. However, at these accuracies (of 1 cm or below), the tidal interactions between the Earth and the Moon are complex to model, especially when considering that the inner structure of the Moon is poorly known Wieczorek (2007); Williams & Boggs (2015). This explains why the UFF test is only possible after an improvement of the dynamical modelling of the Earth–Moon interactions.

Recently, thanks to the GRAIL mission, an unprecedented description of the shape of the lunar gravity field and its variations were obtained for the 6 months of the duration of the mission Konopliv et al. (2014); Lemoine et al. (2014). This information is crucial for a better understanding of the dissipation mechanism over longer time span Matsumoto et al. (2015); Williams & Boggs (2015); Matsuyama et al. (2016). Furthermore, since 2015, the Grasse station, which produces more than 50  per cent of the LLR data, has installed a new detection path at 1064 nm (IR) ranging wavelength leading to a significant increase of the number of observations and of the signal-to-noise ratio Courde et al. (2017).

Together with these new instrumental and GRAIL developments, the Moon modelling of the INPOP planetary ephemeris was improved. Since 2006, INPOP has become a reference in the field of the dynamics of the Solar system objects and in fundamental physics Fienga et al. (2011, 2017).

The INPOP17a version presented here also benefits some of the planetary improvements brought by the use of updated Cassini-deduced positions of Saturn. The planetary and lunar Chebyshev polynomials built from INPOP17a have been made available on the INPOP website¹ together with a detailed technical documentation Viswanathan et al. (2017).

Since 2010, thanks to the millimetre-level accuracy of the LLR measurements and the developments in the dynamical modelling of the Earth–Moon tidal interactions, differences in acceleration of Earth and Moon in free-fall towards the direction of the Sun could reach an accuracy of the order of 10⁻¹⁴ Merkowitz (2010); Williams, Turyshev & Boggs (2012). With the improvement brought by GRAIL, addition of IR LLR observations, and the recent improvement of the dynamical modelling of INPOP17a, one can expect to confirm or improve this limit.

In this paper, we first present (see Section 2.1) the statistics related to the IR data set obtained at the Grasse station since 2015. In Section 2.2, we introduce the updated dynamical model of the Moon as implemented in the INPOP planetary ephemeris including contributions from the shape of the fluid core. In Section 2.4, we explain how we use the IR data to fit the lunar dynamical model parameters with the GRAIL gravity field coefficients as a supplementary constraint for the fluid core description.

Finally in Section 3, we describe how we test the UFF and give new constraints. In addition, we present a generalization of the interpretation in terms of gravitational-to-inertial mass ratios of UFF constraints, based on recent developments in dilaton theories Hees & Minazzoli (2015); Minazzoli & Hees (2016). Hinged on this generalization, we deduce that from a pure phenomenological point of view, one cannot interpret UFF violation tests in the Earth–Moon system as tests of the difference between gravitational and inertial masses only.

2 LUNAR EPHEMERIDES

The new INPOP planetary ephemeride INPOP17a (Viswanathan et al. 2017) is fitted to LLR observations from 1969 to 2017, including the new IR LLR data obtained at the Grasse station.

2.1 Lunar laser ranging

The principle of the LLR observations is well documented (Murphy et al. 2012; Murphy 2013). Besides the lunar applications, the laser ranging technique is still intensively used for tracking Earth orbiting satellites, especially for very accurate orbital (Peron 2013; Lucchesi et al. 2015) and geophysical studies (Jeon et al. 2011; Matsuo et al. 2013).

Non-uniform distributions in the data set are one contributor to correlations between solution parameters (Williams, Turyshev & Boggs 2009). Like one can see in Figs 1–3, about 70  per cent of the data are obtained after reflection on A15 reflector and on an average 40  per cent of the data are acquired within 30° of the quarter Moons.

Figure 1.

Histogram of annual frequency of LLR data with relative contribution from each LRR array including Grasse IR (1064 nm) observations. Points indicate the annual mean of post-fit residuals (in cm) obtained with INPOP17a. The dominance of range observations to A15 is evident. A change can be noticed after 2014 due to the contribution from IR at Grasse. The values in percentage indicate the LLR data contribution from each participating station.

Figure 2.

Histogram of synodic distribution of normal points obtained at Apache Point (c), at the Grasse station from 2012 from 2014 at 542 nm (b) and from 2014 to 2016 at 1064 nm (a). Q indicates the quarter Moon phase.

Figure 3.

Grasse LLR data retroreflector distribution at 532 and 1064 nm from 2015 to 2017. A and L indicate Apollo and Lunokhod retroreflectors, numbered by their respective lunar missions.

In this study, we show how the IR LLR observations acquired at the Grasse station between 2015 and 2017 (corresponding to 7  per cent of the total LLR observations obtained between 1969 and 2017 from all known ILRS ground stations) can help to reduce the presence of such heterogeneity.

2.1.1 Spatial distribution

Statistics drawn from the historical LLR data set (1969–2015) show an observer bias to range to the larger Apollo reflector arrays (mainly Apollo 15). This trend (see Figs 1 and 3) is also present on statistics taken during time periods after the re-discovery of Lunokhod 1 (L1) by Murphy et al. (2011). This is due to the higher return rate and thermal stability over a lunar day on the Apollo reflectors, thereby contributing to the higher likelihood of success.

With the installation of the 1064 nm detection path (see Fig. 3), as explained in Courde et al. (2017), the detection of photon reflected on all reflectors is facilitated, especially for Lunokhod 2 (L2): about 17  per cent of IR data are obtained with L2 when only 2  per cent were detected at 532 nm.

Owing to the spatial distribution of the reflectors on the Moon, Apollo 11 and 14 give sensitivity to longitude librations, Apollo 15 gives sensitivity to latitude librations, and the Lunokhod reflectors give sensitivity both in the latitude and longitude libration of the Moon. The heterogeneity in the reflector distribution of LLR data affects then the sensitivity of the lunar modelling adjustment Viswanathan et al. (2016). By acquiring a better uniformity in the reflector sampling, IR contributes to improve the adjustment of the Moon dynamical and rotational modelling (see Section 2.5).

2.1.2 Temporal distribution

The full and new Moon periods are the most favourable for testing gravity, as the gravitational and tidal effects are maximum. This was partially demonstrated by Nordtvedt (1998). In Fig. 2 are plotted the distributions of normal points relative to the synodic angle for APOLLO (in capitals, abbreviation for Apache Point Observatory Lunar Laser-ranging Operation, while Apollo refers to the US manned lunar missions) and Grasse station obtained at 532 and 1064 nm. About 25  per cent of the APOLLO data sample and almost 45  per cent of the Grasse 532 nm data sample are obtained within 30° of the quarter Moons. This can be explained by two factors.

• New Moon phase. As the pointing of the telescope on to the reflectors is calibrated with respect to a nearby topographical feature on the surface of the Moon, the pointing itself becomes a challenge when the reference points lie in the unlit areas of the Moon. Also, as the new Moon phase occurs in the daylight sky, the noise floor increases and the detector electronics become vulnerable due to ranging at a very close angle to the Sun Williams et al. (2009); Courde et al. (2017).

• Full Moon phase. During this phase, thermal distortions remain as the primary challenge, arising due to the overhead Sun heating of the retroreflector arrays. This induces refractive index gradients within each corner cube causing a spread in the return beam, which makes detection more difficult. The proportion of this effect is partially linked to the thermal stability of the arrays. Since the A11, A14, and A15 arrays have a better thermal stability compared to the L1 and L2 arrays Murphy et al. (2014), observations to the latter become sparse during the full Moon phase (where A and L indicate Apollo and Lunokhod retroreflectors, respectively).

Despite these challenges, LLR observations during the above-mentioned phases of the Moon have been acquired with the IR detection.

After the first 2 yr of 1064 nm detection path at the Grasse station, the observations obtained within the 30° of the quarter Moons are reduced to 32  per cent, effectively increasing by around 10  per cent the portion of data sample close from the most favourable periods (new and full Moon) for tides and UFF studies.

This is primarily achieved due to the improved signal-to-noise ratio resulting from an improved transmission efficiency of the atmosphere at the IR wavelength of 1064 nm. In addition, high-precision data have also been acquired on the two Lunokhod reflector arrays during full Moon phase.

In Section 3, we will see how the IR LLR data help to improve the results related to the UFF tests.

2.1.3 Observational accuracy of the LLR observations

APOLLO observations are obtained with a 3.5 m telescope (under time sharing) at the Apache Point Observatory, while Grasse observations are obtained with a 1.5 m telescope dedicated for Satellite Laser Ranging (SLR) and LLR. A larger aperture is beneficial for statistically reducing the uncertainty of the observation Murphy (2013), which translates to millimetre-level accuracies for APOLLO. One can notice in Fig. 4 that the current lunar ephemerides have a post-fit residual scatter (rms) of about 1–2 cm for the recent observations while the LLR normal point accuracy is given to be at least two times smaller. This calls for an improvement of the Earth–Moon dynamical models within highly accurate numerically integrated ephemerides (see Section 2.5).

Figure 4.

APOLLO and Grasse LLR observations in terms of (i) observational accuracy as given by the annual mean of normal point uncertainty [converted from ps to one-way light time (LT) in cm] and (ii) annual weighted root mean square of post-fit residuals (one-way LT in cm) obtained with INPOP17a.

2.2 Lunar dynamical model

2.2.1 Lunar orbit interactions

In our model, we include the following accelerations perturbing the Moon's orbit.

• Point mass mutual relativistic interactions, in the parametrized post-Newtonian (pN) formalism, from the Sun, planets, and asteroids through Folkner et al. (2014, equation 27).

• Extended bodies’ mutual interactions, through Folkner et al. (2014, equation 28), which include

  • the interaction of the zonal harmonics of the Earth through degree 6;

  • the interaction between zonal, sectoral, and tesseral harmonics of the Moon through degree 6 and the point mass Earth, Sun, Jupiter, Saturn, Venus, and Mars;

  • the interaction of degree 2 zonal harmonic of the Sun.

• Interaction from the Earth tides, through Folkner et al. (2014, equation 32).

The tidal accelerations from the tides due to the Moon and the Sun are separated into three frequency bands (zonal, diurnal, and semi-diurnal). Each band is represented by a potential Love number k_{2m, E} with a matching pair of time delays τ_{Xm, E} (where subscript X is either associated with the daily Earth rotation τ_{Rm, E} or orbital motion τ_{Om, E}) to account for frequency-dependent phase shifts from an anelastic Earth with oceans. Here the time delay represents the phase lag induced by the tidal components. Although the time delay method inherently assumes that the imaginary component of k_{2m, E} varies linearly with frequency, it reduces the complexity of the dynamical model. The diurnal τ_{R1, E} and semi-diurnal τ_{R2, E} are included as solution parameters in the LLR analysis, while model values for potential Love numbers for a solid Earth are fixed to that from Petit & Luzum (2010, table 6.3) followed by corrections from the ocean model FES2004 Lyard et al. (2006). A detailed explanation about the most influential tides relevant to the Earth–Moon orbit integration can be found in Williams & Boggs (2016, table 6).

2.2.2 Lunar orientation and inertia tensor

• Lunar frame and orientation. The mantle coordinate system is defined by the principal axes of the undistorted mantle, whose moments of inertia matrix are diagonal. The time-varying mantle Euler angles (ϕ_m(t),θ_m(t),ψ_m(t)) define the orientation of the principal axis frame with respect to the inertial ICRF2 frame (see Folkner et al. 2014 for details). The time derivatives of the Euler angles are defined through Folkner et al. (2014, equation 14).

• Lunar moment of inertia tensor. The undistorted total moment of inertia of the Moon Ĩ_T is given by

  \begin{eqnarray} \tilde{I}_{\rm T} &=& \frac{\tilde{C}_{\rm T}}{m_{{\rm M}}R^{2}_{{\rm M}}}{\left[\begin{array}{ccc}1 &\quad 0 &\quad 0 \\ 0 &\quad 1 &\quad 0 \\ 0 &\quad 0 &\quad 1\end{array}\right]}\nonumber \\ &&+ {\left[\begin{array}{ccc}\tilde{C}_{2,0,M} - 2\tilde{C}_{2,2,M} &\quad 0 &\quad 0 \\ 0 &\quad \tilde{C}_{2,0,M} + 2\tilde{C}_{2,2,M} &\quad 0 \\ 0 &\quad 0 &\quad 0 \end{array}\right]}, \end{eqnarray}

  (1)

  where |$\tilde{C}$|_{n, m, M} is the unnormalized degree n, order m of the Stokes coefficient C_{n, m} for the spherical harmonic model of the undistorted Moon and |$\tilde{C}_{\rm T}$| is the undistorted polar moment of inertia of the Moon normalized by its mass m_M and radius squared |$R^{2}_{\rm M}$|⁠. Through equation (1), we are able to directly use the undistorted value of C₂₂ Manche (2011) from GRAIL-derived spherical harmonic model of Konopliv et al. (2013).
  The moment of inertia of the fluid core I_c is given by

  \begin{equation} I_{\rm c} = \alpha {_{\rm c}}\tilde{C_{\rm T}}{\left[\begin{array}{ccc}1-f_{\rm c} &\quad 0 &\quad 0 \\ 0 &\quad 1-f_{\rm c} &\quad 0 \\ 0 &\quad 0 &\quad 1 \end{array}\right]} ={\left[\begin{array}{ccc}A_{\rm c} &\quad 0 &\quad 0 \\ 0 &\quad B_{\rm c} &\quad 0 \\ 0 &\quad 0 &\quad C_{\rm c} \end{array}\right]}, \end{equation}

  (2a)

  where α_c is the ratio of the fluid core polar moment of inertia C_c to the undistorted polar moment of inertia of the Moon C_T, f_c is the fluid core polar flattening, and A_c and B_c are the equatorial moments of the fluid core. This study assumes an axisymmetric fluid core with A_c = B_c.
  The moment of inertia of the mantle I_m has a rigid-body contribution |$\tilde{I}_{\rm m}$| and two time-varying contributions due to the tidal distortion of the Earth and spin distortion as given in Folkner et al. (2014, equation 41). The single time delay model (characterized by τ_M) allows for dissipation when flexing the Moon Standish & Williams (1992); Williams et al. (2001); Folkner et al. (2014),

  \begin{equation} \tilde{I}_{\rm m} = \tilde{I}_{\rm T} - I_{\rm c}. \end{equation}

  (2b)
• Lunar angular momentum and torques. The time derivative of the angular momentum vector is equal to the sum of torques (N) acting on the body. In the rotating mantle frame, the angular momentum differential equation for the mantle is given by

  \begin{equation} \frac{{\rm d}}{{\rm d}t}I_{\rm m}\omega _{\rm m} + \omega _{\rm m} \times I_{\rm m}\omega _{\rm m} = N, \end{equation}

  (2c)

  where N is the sum of torques on the lunar mantle from the point mass body A (N_{M, figM-pmA}), figure–figure interaction between the Moon and the Earth [N_{M, figM-figE}, using Folkner et al. (2014, equation 44)] and the viscous interaction between the fluid core and the mantle (N_CMB).
  The motion of the uniform fluid core is controlled by the mantle interior, with the fluid core moment of inertia (I_c) constant in the frame of the mantle. The angular momentum differential equation of the fluid core in the mantle frame is then given by

  \begin{equation} \frac{{\rm d}}{{\rm d}t}I_{\rm c}\omega _{\rm c} + \omega _{\rm m} \times I_{\rm c}\omega _{\rm c} = -N_{{\rm CMB}}, \end{equation}

  (2d)

  \begin{equation} N_{{\rm CMB}} = k_{\rm v} \big ( \omega _{\rm c} - \omega _{\rm m} \big ) + \big ( C_{\rm c} - A_{\rm c} \big ) \big (\hat{z}_{\rm m} \cdot \omega _{\rm c}\big )\big (\hat{z}_{\rm m} \times \omega _{\rm c} \big ), \end{equation}

  (2e)

  where k_v is the coefficient of viscous friction at the core-mantle boundary (CMB) and |$\hat{z}_{\rm m}$| is a unit vector aligned with the polar axis of the mantle frame. The second part on the right-hand side of equation (2e) is the inertial torque on the axisymmetric fluid core.

2.3 Reduction model

The reduction model for the LLR data analysis has been implemented within a precise orbit determination and geodetic software: gins Marty et al. (2011); Viswanathan et al. (2015) maintained by space geodesy teams at GRGS/OCA/CNES and written in fortran90. The subroutines for the LLR data reduction within gins are vetted through a step-wise comparison study conducted among the LLR analysis teams in OCA-Nice (this study), IMCCE-Paris, and IfE-Hannover, by using simulated LLR data and DE421 Folkner, Williams & Boggs (2009) as the planetary and lunar ephemeris. The modelling follows the recommendations of IERS 2010 Petit & Luzum (2010). To avoid any systematics in the reduction model, the upper limit on the discrepancy between the teams was fixed to 1 mm in one-way light time.

From each normal point, the emission time (in UTC) and the round-trip time (in seconds) are used to iteratively solve for the reflection time in the light-time equations. A detailed description is available in Moyer (2003, sections 8 and 11) for a precise round-trip light-time computation.

A detailed description of the reduction model used for this study is provided in Manche (2011).

2.4 Fitting procedure

For APOLLO station observations, scaling the uncertainties of the normal points depending on the change of equipments, or a change in the normal point computation algorithm, is advised (see http://physics.ucsd.edu/~tmurphy/apollo/151201_notes.txt). Unrealistic uncertainties present in observations from Grasse, McDonald MLRS2, and Matera between time periods 1998–1999, 1996, and 2010–2012, respectively, are rescaled.

During the fitting procedure, bounds are used Stark & Parker (1995) for limiting the variability of the estimated parameters, while considering the parameter correlation and variance within the normal matrix. For the gravity field coefficients (including C_{2, 0, M} and C_{2, 2, M}), the bounds are placed using the uncertainties provided by GRAIL [after scaling the formal uncertainties by a factor of 40, following the recommendation by Konopliv et al. (2013)] with their values centred on the GRAIL gravity field estimates.

Additional details of the weighting scheme and the fitting procedure used for the construction of INPOP17a solution can be found in Viswanathan et al. (2017). A filtering scheme is enforced during the iterative fit of the parameters. At each iteration, the residuals are passed through a 3σ filter (where σ is recomputed at each iteration).

2.4.1 Biases

Changes in the ground station introduce biases in the residuals. These biases correspond either with a known technical development at the station (new equipment, change of optical fibre cables) or systematics. Any estimated bias can be correlated with a corresponding change in the ground station, provided the incidents have been logged. A list of known and detected biases is given in Viswanathan et al. (2017).

2.5 Results

Table 1 gives the comparison of post-fit residuals of LLR observations from different ground stations, obtained with the previous solution INPOP13c and the new solution INPOP17a. Table 2 and 3 show the improvement brought by the IR LLR observations on the post-fit residuals of Grasse and APOLLO stations, respectively. Table 4 provides a list of the fixed parameters while Table 5 gives the list of the adjusted parameters related to the lunar interior. The fitted coordinates of the Moon reflectors and of the LLR stations can be found in Viswanathan et al. (2017). As the LLR observations are not included in the construction of the ITRF Altamimi et al. (2016), small corrections to the LLR station coordinates help for the improvement of LLR residuals during the construction of the lunar ephemerides. The Earth orientation parameters and the modelling of the Earth rotation are however kept fixed to the IERS convention (see Section 2.3).

The solution INPOP_G with an axisymmetric core fitted to LLR observations serves as a validation of our lunar model and analysis procedure, against the DE430 Jet Propulsion Laboratory planetary and lunar ephemeris analysis described in Folkner et al. (2014) and Ephemeris of Planets and the Moon (EPM) Institute of Applied Astronomy Russian Academy of Sciences ephemeris in Pavlov, Williams & Suvorkin (2016). Only 532 nm wavelength LLR data are used for matching with the DE430 and EPM ephemeris. In Folkner et al. (2014), Pavlov et al. (2016), and INPOP_G, gravity field coefficients up to degree and order 6 are used for the Moon (GL0660b from Konopliv et al. 2013) and the Earth [GGM05C from Ries et al. (2016) for INPOP17a ephemeris and EGM2008 from Pavlis et al. (2012, 2013) for DE/EPM ephemerides]. Coefficients C₃₂, S₃₂, and C₃₃ are then included in the fit parameters as they improve the overall post-fit residuals. For INPOP_G, the improvement of the formal uncertainty compared to Pavlov et al. (2016), especially in the estimation of parameter k_v/C_T, indicates a strong dissipation mechanism within the Moon, through viscous torques at the fluid core–mantle boundary. Overall, INPOP uncertainties are consistent with EPM Pavlov et al. (2016) published values. DE Williams, Boggs & Folkner (2013); Folkner et al. (2014) uncertainties are greater than INPOP and EPM, and should therefore be considered as more realistic.

Differences between GL0660b values and fitted C₃₂, S₃₂, and C₃₃ from Folkner et al. (2014), Pavlov et al. (2016), or in INPOP_G are several orders of magnitude greater than the mean GRAIL uncertainties (see Konopliv et al. 2013). These results suggest that some significant effects impacting the LLR observations are absorbed by the adjustment of the degree 3 of the full Moon gravity field.

The solution INPOP_G + IR refers to the addition of 2 yr of IR LLR observations Courde et al. (2017) described in Section 2.1 and built in following the same specification as of INPOP_G.

This data set is weighted at the same level as the APOLLO station normal points within the estimation procedure (see Section 2.4).

The first outcome from the introduction of the IR data sets is the improvement of the post-fit residuals obtained for L1 reflector as one can see in Tables 2 and 3 and in Figs 5 –8. This is due to the increase of normal points obtained for this reflector as discussed in Section 2.1.1.

Figure 5.

Post-fit residuals in (cm) versus time (year) obtained with INPOP_G + IR specification (Section 2.5) for McDonald, MLRS1, MLRS2, Haleakala, and Matera stations.

Figure 6.

Post-fit residuals in (cm) versus time (year) obtained with INPOP_G + IR specification (Section 2.5) for Grasse station with the green wavelength.

Figure 7.

Post-fit residuals in (cm) versus time (year) obtained with INPOP_G + IR specification (Section 2.5) for Grasse station with the IR wavelength.

Figure 8.

Post-fit residuals in (cm) versus time (year) obtained with INPOP_G + IR specification (Section 2.5) for APOLLO station.

The second conclusion is that because of only 2 yr on data, the improvement brought by the addition of IR data on the estimated parameters characterizing the Moon and its inner structure is significant, especially for those quantifying the dissipation mechanism such as Q_27.212 and τ_M with a decreasing uncertainty or |$\frac{k_{{\rm v}}}{C_{{\rm T}}}$| and f_c with a significant change in the fitted value (see Table 5).

A significant global improvement is noticeable when one compares post-fit residuals obtained with INPOP_G and with INPOP_G + IR with those obtained with INPOP13c as presented in Fienga et al. (2014) or in Tables 2 and 3. Finally, one should notice in Table 1 the 1.15 cm obtained for the post-fit weighted rms obtained for the 3 yr of the last period of the APOLLO data (group D) as well as that for the IR Grasse station.

Table 1.

Comparison of post-fit residuals of LLR observations from ground stations with corresponding time span, number of normal points available, and number of normal points used in each solution after a 3σ rejection filter. The WRMS (in cm) is obtained with solutions INPOP13c (1969–2013) and INPOP17a (1969–2017). INPOP13c statistics are drawn from Fienga et al. (2014).

				INPOP13c		INPOP17a
Code	Station	Time span	Available	Used	WRMS	Used	WRMS
					(cm)		(cm)
70610	APOLLO, NM, USA (group A)	2006–2010	941	940	4.92	929	1.27
70610	APOLLO, NM, USA (group B)	2010–2012	506	414	6.61	486	1.95
70610	APOLLO, NM, USA (group C)	2012–2013	361	359	7.62	345	1.52
70610	APOLLO, NM, USA (group D)	2013–2016	832	–	–	800	1.15
01910	Grasse, FR	1984–1986	1187	1161	16.02	1161	14.01
01910	Grasse, FR	1987–1995	3443	3411	6.58	3407	4.11
01910	Grasse, FR	1995–2006	4881	4845	3.97	4754	2.86
01910	Grasse, FR	2009–2013	999	990	6.08	982	1.41
01910	Grasse, FR	2013–2017	3351	–	–	3320	1.51
56610	Haleakala, HI, USA	1984–1990	770	739	8.63	728	4.80
07941	Matera, IT	2003–2013	83	70	7.62	37	2.37
07941	Matera, IT	2013–2015	30	–	–	28	2.93
71110	McDonald, TX, USA	1969–1983	3410	3302	31.86	3246	18.87
71110	McDonald, TX, USA	1983–1986	194	182	20.60	148	16.77
71111	MLRS1, TX, USA	1983–1984	44	44	29.43	44	32.73
71111	MLRS1, TX, USA	1984–1985	368	358	77.25	356	62.58
71111	MLRS1, TX, USA	1985–1988	219	207	7.79	202	11.07
71112	MLRS2, TX, USA	1988–1996	1199	1166	5.36	1162	3.81
71112	MLRS2, TX, USA	1996–2012	2454	1972	5.81	1939	3.72
71112	MLRS2, TX, USA	2012–2015	17	–	–	15	2.59
TOTAL		1969–2017	25 289	20 160		24 089

Table 2.

Grasse LLR data retroreflector statistics computed using post-fit residuals obtained with INPOP_G and INPOP_G + IR, within the fit intervals 01/01/2015 to 01/01/2017 (with a 3σ filter), with the WRMS in m (rms weighted by the number of normal points from each reflector).

Grasse
LRRR	INPOPG	INPOPG + IR	per cent change	NPTs
A15	0.0183	0.0181	1.1	1018
A14	0.0203	0.0177	12.8	172
A11	0.0267	0.0239	10.5	215
L1	0.0215	0.0166	22.8	265
L2	0.0246	0.0215	12.6	256
WRMS	0.0207	0.0189	9.5	1926

Table 3.

APOLLO LLR data retroreflector statistics computed using post-fit residuals obtained with INPOP_G and INPOP_G + IR, within the fit intervals 01/01/2015 to 01/01/2017 (with a 3σ filter), with the WRMS in m (rms weighted by the number of normal points from each reflector).

APOLLO
LRRR	INPOPG	INPOPG + IR	per cent change	NPTs
A15	0.0127	0.0127	0.2	344
A14	0.0192	0.0177	7.8	176
A11	0.0185	0.0169	8.7	164
L1	0.0186	0.0157	15.6	89
L2	0.0136	0.0137	−0.7	64
WRMS	0.0159	0.0149	6.7	837

Table 4.

Fixed parameters for the Earth–Moon system.

Parameter	Units	INPOP	DE430	EPM
(EMRATa − 81.300 570) × 106		1.87	−0.92	−0.92c
(RE − 6378.1366) × 104	km	0.0	−3	0.0
|$(\dot{J}_{2E} -2.6\times 10^{-11})$|	yr−1	0.0	0.0	0.0
(k20, E − 0.335)		0.0	0.0	0.0
(k21, E − 0.32)		0.0	0.0	0.0
(k22, E − 0.301 02)		−0.019 02	0.018 98	−0.019 02
(τO0, E − 7.8 × 10−2) × 102	d	0.0	−1.4	0.0
(τO1, E + 4.4 × 10−2)	d	0.0	0.0b	0.0
τO2, E + 1.13 × 10−1) × 101	d	0.0	0.13	0.0
(RM − 1738.0)	km	0.0	0.0	0.0
(αC − 7.0 × 10−4)		0.0	0.0	0.0
(k2, M − 0.024 059)		0.0	0.0	0.0
(l2 − 0.0107)		0.0	0.0	0.0

Notes.^a EMRAT is fitted during the joint analysis between the lunar and planetary part.

^bτ_{O1, E} in Folkner et al. (2014) given as −0.0044 is a typographical error.

^cEMRAT in the EPM solution Pavlov et al., 2016, is fixed to a value obtained from DE430

3 TEST OF THE EQUIVALENCE PRINCIPLE

3.1 Context

Among all possibilities to test GR, the tests of the motion of massive bodies as well as the propagation of light in the Solar system were historically the first ones, and still provide the highest accuracies for several aspects of gravity tests (see Berti et al. 2015; Joyce et al. 2015; Yunes, Yagi & Pretorius 2016 for recent overviews of constraints on alternative theories from many different types of observations). This is in part due to the fact that the dynamics of the Solar system is well understood and supported by a long history of observational data.

In GR, not only do test particles with different compositions fall equally in a given gravitational field, but also extended bodies with different gravitational self-energies. While a deviation from the former case would indicate a violation of the weak equivalence principle (WEP), a deviation from the latter case would be a sign of a violation of the strong equivalence principle (SEP; Will 2014). Violations of the equivalence principles are predicted by a number of modifications of GR, often intending to suggest a solution for the problems of dark energy and dark matter (Capozziello & de Laurentis 2011; Berti et al. 2015; Joyce et al. 2015) and/or to put gravity in the context of quantum field theory (Kostelecký 2004; Woodard 2009; Donoghue 2017). The UFF, an important part of the equivalence principle, is currently tested at a level of about 10⁻¹³ with torsion balances Adelberger et al. (2003) and LLR analyses Williams et al. (2012).

As the Earth and the Moon both fall in the gravitational field of the Sun – and because they neither have the same compositions nor the same gravitational self-energies – the Earth–Moon system is an ideal probe of both the WEP and the SEP, while torsion balance Adelberger et al. (2003) or MICROSCOPE Liorzou et al. (2014) is only sensitive to violations of the WEP.

In this paper, we implemented the equations given in Williams et al. (2012) and introduce in the INPOP fit the differences between the accelerations of the Moon and the Earth.

The aim of this work is first to give the most general constraint in terms of acceleration differences without assuming metric theories or other types of alternative theories (Section 3.3). In a second step (Section 3.4), we propose two interpretations: one following the usual formalism proposed by Nordtvedt (see e.g. Nordtvedt 2014 and references therein), and the other following the dilaton theory Damour & Polyakov (1994); Hees & Minazzoli (2015); Minazzoli & Hees (2016).

3.2 Method

In order to test possible violations of GR in terms of UFF, a supplementary acceleration is introduced in the geocentric equation of motion of the Moon, such that the UFF violation-related difference between the Moon and the Earth accelerations reads Nordtvedt (1968b)

\begin{eqnarray} \Delta \boldsymbol a^{\overline{{\rm UFF}}}\equiv (\boldsymbol a_{\rm M} - \boldsymbol a_{\rm E})^{\overline{{\rm UFF}}} = \boldsymbol a_{{\rm E}} \Delta _{{\rm ESM}}. \end{eqnarray}

(3)

Δ_ESM is estimated in the LLR adjustment together with the other parameters of the lunar ephemerides given in Table 5. In what follows, we shall name Δ_ESM ‘UFF violation parameter’. ESM stands for the three bodies involved, namely the Earth, the Sun, and the Moon, respectively. As we shall see in Section 3.4.2, some theoretical models induce a dependence of the UFF violation parameter on the composition of the Sun, in addition to the ‘more usual’ dependence on the compositions and on the gravitational binding energies of the Moon and the Earth.

Table 5.

Extended body parameters for the Earth and the Moon. Uncertainties for INPOP_G and INPOP_G + IR (1σ) are obtained from a 5  per cent jackknife, while other solutions (DE430 and EPM) are assumed as (1σ) formal uncertainties.

Parameter	Units	INPOPG	INPOPG + IR	DE430	EPM
(GMEMB − 8.997 011 400 × 10−10) × 1019	au3/d2	4 ± 2	4 ± 2	−10	10 ± 5
(τR1, E − 7.3 × 10−3) × 105	d	0 ± 4	6 ± 3	6 ± 30	57 ± 5
(τR2, E − 2.8 × 10−3) × 105	d	9.2 ± 0.4	8.7 ± 0.3	−27 ± 2	5.5 ± 0.4
(CT/(mMR2) − 0.393 140) × 106		6.9 ± 0.2	8.2 ± 0.2	2a	2a
(C32 − 4.840 4981 × 10−6b) × 109		4.1 ± 0.3	3.9 ± 0.3	4.4	4.4 ± 0.1
(S32 − 1.666 1414 × 10−6b) × 108		1.707 ± 0.006	1.666 ± 0.006	1.84	1.84 ± 0.02
(C33 − 1.711 6596 × 10−6b) × 108		−1.19 ± 0.04	−2.40 ± 0.04	−3.6	−4.2 ± 0.2
(τM − 9 × 10−2) × 104	d	−14 ± 5	−35 ± 3	58.0 ± 100	60 ± 10
|$(\frac{k_{{\rm v}}}{C_{{\rm T}}}-1.6\times 10^{-8}) \times 10^{10}$|	d−1	12.7 ± 0.4	15.3 ± 0.5	4.0 ± 10.0	3.0 ± 2.0
(fc − 2.1 × 10−4) × 106		37 ± 3	42 ± 3	36 ± 28	37 ± 4
(h2 − 3.71 × 10−2c) × 103		6.3 ± 0.2	6.8± 0.2	11.0 ± 6	6 ± 1
Q27.212 − 45 (derived)		3.9 ± 0.5	5.0 ± 0.2	0 ± 5	0 ± 1

Notes.^aDerived quantity.

^bC₃₂, S₃₂, and C₃₃ are reference values from the GRAIL analysis by Konopliv et al. (2013).

^ch₂ reference value from LRO-LOLA analysis by Mazarico et al. (2014).

In order to estimate Δ_ESM with the appropriate accuracy, one should correct for supplementary effects such as the solar radiation pressure and thermal expansion of the retroreflectors Vokrouhlický (1997); Williams et al. (2012). An empirical correction on the radial perturbation (Δr_EM) induced by the UFF test has to be applied. For instance, with some simplifying approximations (Nordtvedt 2014), one can show that the UFF additional acceleration would indeed lead to an additional radial perturbation (Δr_EM) of the Moon's orbit towards the direction of the Sun given by

\begin{eqnarray} \Delta r_{{\rm EM}}= S \Delta _{{\rm ESM}} \cos D, \end{eqnarray}

(4)

where S is a scaling factor of about −3 × 10¹⁰ m Williams et al. (2012) and D is the synodic angle. A correction Δr = 3.0 ± 0.5 mm Vokrouhlický (1997); Williams et al. (2012) is then applied in order to correct for solar radiation pressure and thermal radiation of the retroreflectors, and a new corrected value of Δ_ESM is then deduced (see Table 6).

3.3 Results

Fits were performed including in addition to the previous fitted parameters presented in Table 5 the UFF violation parameter Δ_ESM given in equation (3). Two different fits were considered including 532 and 1064 nm data sets (solution labelled INPOP_G + IR), or just the 532 nm data sets (solution labelled INPOP_G). A supplementary adjustment was also performed for a better comparison to the previous determination from other LLR analysis groups, which were limited to a data sample up to 2011 (labelled as limited data). Results are given in Table 6.

Table 6.

Comparison of results for the value of Δ_ESM (column 4) estimated with the solution INPOP17A fitted to LLR data set between (1) 1969–2011 (for comparison with Müller et al. 2012; Williams et al. 2012); (2) 1969–2017 with data obtained only in green wavelength, (3) 1969–2017 with data obtained with both green and IR wavelength. Column 5 empirically corrects the radial perturbation from effects related to solar radiation pressure and thermal expansion of retroreflectors using equation (4), with a value Δr = 3.0 ± 0.5 mm Williams et al. (2012). Column 6 contains the value of Δ_ESM after applying the corrections of column 5. Column 7 contains the parameter η obtained using equation (13). See discussion in Section 4.

Reference	Data	Uncertainty	Estimated	Corrected	Corrected	Parameter
	time span		ΔESM	cos D	ΔESM	ηc
	(Year)		(× 10−14)	(mm)	(× 10−14)	(× 10−4)
Williams et al. (2009)a	1969–2004	N/A	3.0 ± 14.2	2.8 ± 4.1	−9.6 ± 14.2	2.24 ± 3.14
Williams et al. (2012)	1969–2011	N/A	0.3 ± 12.8	2.9 ± 3.8	−9.9 ± 12.9	2.25 ± 2.90
Müller et al. (2012)a, b	1969–2011	3σ	−14 ± 16	–	–	–
INPOP17A (limited data)	1969–2011	3σ	−3.3 ± 17.7	4.0 ± 5.2	−13.5 ± 17.8	3.03 ± 4.00
Hofmann & Müller (2016)a	1969–2016	3σ	–	–	−3.0 ± 6.6	0.67 ± 1.48
INPOP17A (green only)	1969–2017	3σ	5.2 ± 8.7	1.5 ± 2.6	−5.0 ± 8.9	1.12 ± 2.00
INPOP17A (green and IR)	1969–2017	3σ	6.4 ± 6.9	1.1 ± 2.1	−3.8 ± 7.1	0.85 ± 1.59

Notes.^aThermal expansion correction not applied.

^bSRP correction not applied.

^cDerived using |$\frac{|\Omega _{\rm E}|}{m_{\rm E} c^2} - \frac{|\Omega _{\rm M}|}{m_{\rm M} c^2}$| = −4.45× 10⁻¹⁰ (Williams et al. 2012, equation 6).

The additional acceleration of the Moon orbit in the direction of the Sun correlates with a coefficient of 0.95 and 0.90 with GM_EMB and the Earth–Moon mass ratio (EMRAT), respectively. In all the solutions w.r.t. LLR EP estimation, the gravitational mass of the Earth–Moon barycentre (GM_EMB) remains as a fit parameter due its high correlation with the EP parameter (Δ_ESM). EMRAT was estimated from a joint planetary solution and kept fixed during LLR EP tests (for all INPOP solutions in Table 6) due to its weak determination from LLR.

A test solution that fitted EMRAT, with GM_EMB as a fixed parameter, gives an estimate of Δ_ESM = (8 ± 7.0) × 10⁻¹⁴. However, the value of EMRAT estimated from an LLR-only solution has an uncertainty of one order of magnitude greater than that obtained from the joint planetary fit. This is also consistent with a similar result by Williams et al. (2009). As a result, EMRAT was not included as a fit parameter for the estimates provided in Table 6, as it resulted in a degraded fit of the overall solution.

Williams et al. (2012) show that including annual nutation components of the Earth pole direction in space, to the list of fitted parameters during the estimation of LLR EP solution, increases the uncertainty of the estimated UFF violation parameter (Δ_ESM) by 2.5 times. Moreover, it is to be noted that within Table 6, the solutions by Williams, Turyshev & Boggs (2009, 2012) and Müller, Hofmann & Biskupek (2012) use the IERS 2003 McCarthy & Petit (2004) recommendations within the reduction model, while all INPOP17 solutions use IERS 2010 Petit & Luzum (2010) recommendations. The notable difference between the two IERS models impacting the LLR EP estimation is expected to be from the precession nutation of the celestial intermediate pole within the ITRS–GCRS transformation (Petit & Luzum 2010, p. 8).

Equation (4) shows the dependence of Δ_ESM w.r.t. the cosine of the lunar orbit synodic angle, synonymous with the illumination cycle of the lunar phases. Due to the difficulties involved with ranging to the Moon during the lunar phases with the extreme values of cos D (new and full Moon) as described in Section 2.1.2, the LLR observations during these phases remain scarce. The availability of IR LLR observations from Grasse contributes to the improvement of this situation, as shown in Fig. 2. This is reflected in the improvement of the uncertainty of the estimated value of Δ_ESM by 14  per cent, with solutions including the IR LLR data.

Using both IR and green wavelength data, and empirically correcting for the radial perturbation for effects related to solar radiation pressure and thermal expansion, our final result on the UFF violation parameter is given by (see also, Table 6)

\begin{eqnarray} \Delta _{{\rm ESM}} = (-3.8 \pm 7.1) \times 10^{-14}. \end{eqnarray}

(5)

The continuation of the IR observational sessions at Grasse will help to continue the improvement in the Δ_ESM estimations.

An observable bias in the differential radial perturbation of the lunar orbit w.r.t. the Earth, towards the direction of the Sun, if significant and not accounted for within the dynamical model, would result in a false indication of the violation of the principle of equivalence estimated with the LLR observations. Oberst et al. (2012) show the distribution of meteoroid impacts with the lunar phase. Peaks within the histogram in Oberst et al. (2012, p. 186) indicate a non-uniform temporal distribution with a non-negligible increase in both small and large impacts during the new and full Moon phase. Future improvements to the LLR EP estimation must consider the impact of such a bias that could potentially be absorbed during the fit by the LLR UFF violation parameter Δ_ESM.

3.4 Theoretical interpretations

3.4.1 Nordtvedt's interpretation: gravitational versus inertial masses

Although equations of motion are developed at the pN level in INPOP Moyer (2003), violations of the UFF can be cast entirely in the Newtonian equation of motion with sufficient accuracy. As described by Nordtvedt (1968b), a difference of the inertial (m^I) and gravitational (m^G) masses would lead to an alteration of body trajectories in celestial mechanics according to the following equation:

\begin{eqnarray} \boldsymbol {a}_{\rm T}=-\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm T} \sum _{A\ne T} \frac{G m^{\rm G}_A}{r_{AT}^3}\boldsymbol r_{AT} , \end{eqnarray}

(6)

where |$\boldsymbol r_{AT}= \boldsymbol x_{\rm T} - \boldsymbol x_A$| and G is the constant of Newton.

Following Williams et al. (2012), the relative acceleration at the Newtonian level between the Earth and the Moon due to the attraction of the Sun reads

\begin{eqnarray} &&\!\!\!\!&&{\boldsymbol a_{\rm M} - \boldsymbol a_{\rm E}= - \frac{ G \mu }{r_{{\rm EM}}^3}\boldsymbol r_{{\rm EM}}+ G m^{\rm G}_{\rm S}\left[ \frac{\boldsymbol r_{{\rm SE}}}{r_{{\rm SE}}^3}-\frac{\boldsymbol r_{{\rm SM}}}{r_{{\rm SM}}^3} \right]} \nonumber \\ &&+\, G m^{\rm G}_{\rm S}\left[ \frac{\boldsymbol r_{{\rm SE}}}{r_{{\rm SE}}^3} \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm E} -1\right)- \frac{\boldsymbol r_{{\rm SM}}}{r_{{\rm SM}}^3} \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm M} -1\right) \right], \end{eqnarray}

(7)

with |$\mu \equiv m^{\rm G}_{\rm M} + m^{\rm G}_{\rm E}+ \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm E} -1\right)m^{\rm G}_{\rm M}+\left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm M} -1\right) m^{\rm G}_{\rm E}$|⁠. |$\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm E}$| and |$\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm M}$| are the ratios between the gravitational and the inertial masses of the Earth and Moon, respectively.

With ephemerides, the first term of equation (7) does not lead to a sensitive test of the UFF, because it is absorbed in the fit of the parameter |$m^{\rm G}_{\rm M}+m^{\rm G}_{\rm E}$| (e.g. Williams et al. 2012), while the last term does. At leading order, one can approximate both distances appearing in this last term as being approximately equal. One gets

\begin{eqnarray} \Delta \boldsymbol a^{\overline{{\rm UFF}}} &\equiv & (\boldsymbol a_{\rm M} - \boldsymbol a_{\rm E})^{\overline{{\rm UFF}}} \nonumber \\ &\approx & G m^{\rm G}_{\rm S}\left[ \frac{\boldsymbol r_{{\rm SE}}}{r_{{\rm SE}}^3} \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm E} -1\right)- \frac{\boldsymbol r_{{\rm SM}}}{r_{{\rm SM}}^3} \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm M} -1\right) \right]\nonumber \\ &\approx & \boldsymbol a_{{\rm E}}\left[ \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm E} -1\right)- \left(\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm M} -1\right) \right] \nonumber \\ &\equiv & \boldsymbol a_{{\rm E}} \Delta _{{\rm ESM}} \end{eqnarray}

(8)

with

\begin{eqnarray} \Delta _{{\rm ESM}} =\left[\left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm E}- \left(\frac{m^{\rm G}}{m^{\rm I}} \right)_{\rm M} \right]. \end{eqnarray}

(9)

One recovers equation (3). Therefore, in this context, constraints on Δ_ESM can be interpreted as constraints on the difference of the gravitational-to-inertial mass ratios between the Earth and the Moon.

Furthermore, the LLR test of UFF captures a combined effect of the SEP, from the differences in the gravitational self-energies, and the WEP due to compositional differences, of the Earth–Moon system. In general, one has

\begin{eqnarray} \Delta _{{\rm ESM}} =\Delta _{{\rm ESM}}^{{\rm WEP}} + \Delta _{{\rm ESM}}^{{\rm SEP}}. \end{eqnarray}

(10)

In order to separate the effects of WEP, we rely on results from laboratory experiments that simulate the composition of the core and the mantle materials of the Earth–Moon system. One such estimate is provided by Adelberger (2001), which translates to the following mass ratio difference:

\begin{eqnarray} \Delta _{{\rm ESM}}^{{\rm WEP}}&=& \Big [\Big (\frac{m^{\rm G}}{m^{\rm I}}\Big )_{{\rm E}} - \Big (\frac{m^{\rm G}}{m^{\rm I}}\Big )_{{\rm M}}\Big ]_{{\rm WEP}} \end{eqnarray}

(11)

\begin{eqnarray} \hphantom{\Delta _{{\rm ESM}}^{{\rm WEP}}}&=& (1.0 \pm 1.4) \times 10^{-13}. \end{eqnarray}

(12)

It is also possible to deduce the Nordtvedt parameter (η) defined as

\begin{eqnarray} \Delta _{{\rm ESM}}^{{\rm SEP}} &=& \eta _{{\rm SEP}} \left[\left(\frac{|\Omega |}{m c^2}\right)_{\rm E} - \left(\frac{|\Omega |}{m c^2}\right)_{\rm M}\right] \end{eqnarray}

(13)

\begin{eqnarray} \hphantom{\Delta _{{\rm ESM}}^{{\rm SEP}}}&\approx & \eta _{{\rm SEP}} \times (-4.45 \times 10^{-10}), \end{eqnarray}

(14)

where Ω and mc² are the gravitational binding and rest mass energies, respectively, for the Earth and the Moon (subscripts E and M, respectively). The value of −4.45 × 10⁻¹⁰ is obtained from Williams et al. (2009, equation 7).

However, all metric theories lead to a violation of the SEP only. Therefore, for metric theories, it is irrelevant to try to separate violation effects of the WEP and SEP, as the WEP is intrinsically respected.

3.4.2 Dilaton theory and a generalization of the Nordtvedt interpretation

Starting from a general dilaton theory, a more general equation governing celestial mechanics than equation (6) has been found to be Hees & Minazzoli (2015); Minazzoli & Hees (2016)

\begin{eqnarray} \boldsymbol {a}_{\rm T}=-\sum _{A\ne T} \frac{G m^{\rm G}_A}{r_{AT}^3}\boldsymbol r_{AT}\left(1+\delta _{\rm T}+\delta _{AT}\right) . \end{eqnarray}

(15)

The coefficients δ_T and δ_AT parametrize the violation of the UFF. In this expression, the inertial mass |$m^{\rm I}_A$| writes in terms of the gravitational mass |$m^{\rm G}_A$| as |$m^{\rm G}_A=(1+\delta _A) m^{\rm I}_A$| Hees & Minazzoli (2015); Minazzoli & Hees (2016). Of course, since |$m_A^{\rm G}/m_A^{\rm I}=1+\delta _A$|⁠, one recovers equation (6) when δ_AB = 0 for all A and B. From equation (15), one can check that the gravitational force in this context still satisfies Newton's third law of motion:

\begin{equation} m^{\rm I}_A\boldsymbol {a}_A =\frac{G m^{\rm I}_A m^{\rm I}_B}{r_{AB}^3}\boldsymbol r_{AB}\left(1+\delta _A+\delta _B+\delta _{AB}\right)=- m^{\rm I}_B\boldsymbol {a}_B. \end{equation}

(16)

In the dilaton theory, the δ coefficients are functions of ‘dilatonic charges’ and of the fundamental parameters of the theory Damour & Donoghue (2010); Hees & Minazzoli (2015); Minazzoli & Hees (2016). However, in what follows, we will consider the phenomenology based on the δ parameters independently of its theoretical origin, as a similar phenomenology may occur in a different theoretical framework.

In general, δ_T can be decomposed into two contributions: one from a violation of the WEP and one from a violation of the SEP:

\begin{equation} \delta _{\rm T}=\delta _{\rm T}^{{\rm WEP}}+\delta _{\rm T}^{{\rm SEP}},\qquad \textrm{with}\qquad \delta _{\rm T}^{{\rm SEP}} = \eta \frac{|\Omega _{\rm T}|}{m_{\rm T} c^2}. \end{equation}

(17)

The quantity |$\delta _{\rm T}^{{\rm SEP}}$| depends only on the gravitational energy content of the body T. On the other hand, |$\delta _{\rm T}^{{\rm WEP}}$| depends on the composition of the falling body T (Damour & Donoghue 2010; Hees & Minazzoli 2015; Minazzoli & Hees 2016). In some theoretical situations (see e.g. Damour & Donoghue 2010), if |$\delta _{\rm T}^{{\rm WEP}} \ne 0$|⁠, then |$\delta _{\rm T}^{{\rm WEP}} \gg \delta _{\rm T}^{{\rm SEP}}$|⁠, such that one can have either a clean WEP violation or a clean SEP violation.

Like the parameter |$\delta _{\rm T}^{{\rm WEP}}$|⁠, δ_AT depends on the composition of the falling bodies. However, unlike |$\delta _{\rm T}^{{\rm WEP}}$|⁠, it also depends on the composition of the body A that is the source of the gravitational field in which the body T is falling (Hees & Minazzoli 2015; Minazzoli & Hees 2016). As a consequence, the relative acceleration of two test particles with different compositions cannot be related to the ratios between their gravitational-to-inertial masses in general (i.e. |$m_A^{\rm G}/m_A^{\rm I}=1+\delta _A$|⁠). This contrasts with the usual interpretation (see for instance Williams et al. 2012). However, with some theoretical models, |$\delta _{\rm T}^{{\rm WEP}}$| is much greater than δ_AT (Damour & Donoghue 2010; Hees & Minazzoli 2015; Minazzoli & Hees 2016).

At the Newtonian level, the relative acceleration between the Earth and the Moon reads

\begin{eqnarray} \boldsymbol a_{\rm M} - \boldsymbol a_{\rm E} &=& - \frac{ G \mu }{r_{{\rm EM}}^3}\boldsymbol r_{{\rm EM}}+ G m^{\rm G}_{\rm S}\left[ \frac{\boldsymbol r_{{\rm SE}}}{r_{{\rm SE}}^3}-\frac{\boldsymbol r_{{\rm SM}}}{r_{{\rm SM}}^3} \right] \nonumber \\ &&+\, G m^{\rm G}_{\rm S}\left[ \frac{\boldsymbol r_{{\rm SE}}}{r_{{\rm SE}}^3} (\delta _{\rm E}+\delta _{{\rm SE}})- \frac{\boldsymbol r_{{\rm SM}}}{r_{{\rm SM}}^3} (\delta _{\rm M}+\delta _{{\rm SM}}) \right], \end{eqnarray}

(18)

with |$\mu \equiv m^{\rm G}_{\rm M} + m^{\rm G}_{\rm E}+ (\delta _{\rm E}+\delta _{{\rm EM}})m^{\rm G}_{\rm M}+(\delta _{\rm M}+\delta _{{\rm EM}}) m^{\rm G}_{\rm E}$|⁠. As discussed already in the previous subsection, the first term of equation (18) does not lead to a sensitive test of the UFF, because it can be absorbed in the fit of the parameter |$m^{\rm G}_{\rm M}+m^{\rm G}_{\rm E}$| (e.g. Williams et al. 2012), while the last term does. At leading order, one can approximate both distances appearing in this last term as being approximately equal. One therefore has

\begin{eqnarray} \Delta \boldsymbol a^{\overline{{\rm UFF}}} &\equiv & (\boldsymbol a_{\rm M} - \boldsymbol a_{\rm E})^{\overline{{\rm UFF}}} \nonumber \\ &\approx & G m^{\rm G}_{\rm S}\left[ \frac{\boldsymbol r_{{\rm SE}}}{r_{{\rm SE}}^3} (\delta _{\rm E}+\delta _{{\rm SE}})- \frac{\boldsymbol r_{{\rm SM}}}{r_{{\rm SM}}^3} (\delta _{\rm M}+\delta _{{\rm SM}}) \right] \nonumber \\ &\approx & \boldsymbol a_{{\rm E}} \left[(\delta _{\rm E}+\delta _{{\rm SE}})- (\delta _{\rm M}+\delta _{{\rm SM}})\right] \nonumber \\ &\equiv & \boldsymbol a_{{\rm E}} \Delta _{{\rm ESM}}, \end{eqnarray}

(19)

where |$\Delta \boldsymbol a^{\overline{{\rm UFF}}}$| is the part of the relative acceleration between the Earth and the Moon that violates the UFF. Once again, one recovers equation (3) – although its theoretical interpretation is different compared to the previous subsection.

When δ_SM = δ_SE, and especially when δ_SM = δ_SE = 0, one recovers the usual equation (9). But it is not the case in general because the composition of the Sun may affect the dynamics in some cases as well. Therefore, in a more general context than in Section 3.4.1, constraints on Δ_ESM cannot be uniquely interpreted as constraints on the difference of the gravitational-to-inertial mass ratios between the Earth and the Moon.

As a consequence, from a pure phenomenological point of view – or, equivalently, from an agnostic point of view – one should not interpret Δ_ESM in terms of gravitational-to-inertial mass ratios only. Indeed, a more general expression of the UFF violating parameter is given by

\begin{eqnarray} \Delta _{{\rm ESM}} =\left[(\delta _{\rm E}+\delta _{{\rm SE}})- (\delta _{\rm M}+\delta _{{\rm SM}}) \right], \end{eqnarray}

(20)

where one can see that the Sun's composition may affect the dynamics as well, through the coefficients δ_SE and δ_SM.

Otherwise, see a discussion on how to decorrelate the dilaton parameters from planetary ephemeris in Minazzoli et al. (2017).

4 DISCUSSION

As emphasized in Section 3.4.1, metric theories lead to a violation of the SEP only. Hence, it is tempting to use equation (13) in order to convert the result on Δ_ESM in equation (5) into a constraint on the Nordtvedt parameter η_SEP – when considering a metric theory prior.

However, such a conversion would not give a clean constraint on the actual Nordtvedt parameter η_SEP. The reason is that, since η_SEP depends on the pN parameters, one should also fit the extra pN parameters in the Einstein–Infeld–Hoffmann (EIH) equations of motion, at the same time in both the lunar and the planetary ephemeris – because the latter is used in the derivation of the former. Hence, unless a global fit of the various pN parameters and Δ_ESM is done at the same time for the whole Solar system solution, the conversion of Δ_ESM into η_SEP through equation (13) does not give a constraint on the actual Nordtvedt parameter η_SEP, but on another parameter that we shall call η instead – and that is simply defined by equation (13).

Despite this fact, the result on Δ_ESM that is given in equation (5) can nevertheless be interpreted in terms of fundamental physics, because a whole subset of theories predict a large domination of the WEP over the SEP in Δ_ESM Damour & Donoghue (2010); Minazzoli & Hees (2016) – meaning that one would have a violation of the UFF while the pN parameters would be either equal to their value in GR or their difference with respect to their value in GR would be negligible at the present level of experimental accuracy.

However, in order to separate the SEP and WEP contributions to Δ_ESM in a general case – or to determine the Nordtvedt parameter η_SEP when considering a metric theory prior – one would need to consider the whole Solar system simultaneously in a consistent parametrized pN framework. This interesting study is left for a future work.

Nevertheless, an internal test on the impact of the extra pN parameters γ and β in the EIH equations under their known limits [taken from Bertotti, Iess & Tortora (2003) and Fienga et al. (2015), respectively] shows no significant impact on our results, due to the little sensitivity of these parameters to the LLR data. Hence, η represents a good quantitative approximation of the Nordtvedt parameter η_SEP, as deduced from testing the UFF with LLR data only. Moreover, since UFF constraints are often reported in terms of η, this quantity can still be used in order to compare the sensitivity of the various lunar ephemeris solutions with respect to testing the UFF. The estimates of η are reported in Table 6.

5 CONCLUSIONS AND FUTURE WORK

In this paper, we present an improvement in the lunar dynamical model of INPOP ephemeris (version 17a) compared to the previous release (version 13c). The model is fitted to the LLR observations between 1969 and 2017, following the model recommendations from IERS 2010 Petit & Luzum (2010). The lunar parameter estimates obtained with the new solution are provided in Table 5 with comparisons to that obtained by other LLR analysis groups. The improvement brought by the new IR LLR data from Grasse station on the parameter estimates is characterized. The post-fit LLR residuals obtained with INPOP17a are between 1.15 and 1.95 cm over 10 yr of APOLLO data and 1.47 cm over 2 yr of the new IR LLR data from Grasse Viswanathan et al. (2017). Our solution benefits also the better spatial and temporal distribution of the IR Grasse data with an improvement of 14 per cent of the UFF tests and better estimations of the Moon dissipation parameters.

We take advantage of the lunar ephemeris improvements to perform new tests of the UFF. A general constraint is obtained using INPOP, in terms of the differences in the acceleration of the Earth and the Moon towards the Sun. In addition to the Nordtvedt interpretation Nordtvedt (1968b, provided in Section 3.4.1), we propose an alternative interpretation and a generalization of the usual interpretation from the point of view of the dilaton theory (Damour & Polyakov 1994; Hees & Minazzoli 2015; Minazzoli & Hees 2016), provided in Section 3.4.2. We obtain an estimate of the UFF violating parameter Δ_ESM = (− 3.8 ± 7.1) × 10⁻¹⁴, showing no violation of the principle of equivalence at this level. Future work may further allow us to separate between the SEP and the WEP contributions to Δ_ESM by studying the whole Solar system simultaneously in a consistent parametrized pN framework – see discussion in Section 4.

Thermal expansion of the retroreflectors and solar radiation pressure are currently employed as empirical corrections following Vokrouhlický (1997) and Williams et al. (2009). Future LLR analysis will consider an implementation of these effects within the reduction procedure, so as to improve the uncertainty of the EP test. Oberst et al. (2012) show the distribution of meteoroid impacts with the lunar phase, indicating a non-uniform temporal distribution during the new and full Moon phase that could impact the test of EP. The impact of this effect needs to be characterized during the EP test, to be considered as negligible at the present LLR accuracy.

The use of a strictly GRAIL-derived gravity field model Konopliv et al. (2013) highlights longitude libration signatures well above the LLR noise floor, arising from unmodelled effects in lunar ephemeris Viswanathan (2017). Other LLR analysis groups Folkner et al. (2009, 2014); Pavlov et al. (2016) prefer to fit the degree-3 components away from GRAIL-derived gravity field coefficients. Extra periodic terms on the longitude libration present in the DE430 lunar model are not considered within this paper. Instead, a work is in progress to further improve the lunar dynamical model and to identify the cause of the low-degree spacecraft-derived lunar gravity field inconsistency with that from the analysis of LLR data.

Acknowledgements

The authors extend their sincere gratitude to all the observers and engineers at Grasse, APOLLO, McDonald, Matera, and Haleakala LLR stations for providing timely and accurate observations over the past 48 years.

Footnotes

1

Available at: http://www.imcce.fr/inpop