Influence of the layered Moon and Earth's orientation on lunar rotation

Knowledge of the lunar internal structure gives us an opportunity to infer the origin as well as the evolution of our nearest celestial body (Matsumoto et al. 2015; Williams et al. 2001), in particular from the dynamical, rotational and internal structure. Given this situation, study of the lunar rotation has a long history, and many mathematicians and astronomers have spent a lot of time focused on the associated rotational theories. The Moon's rotational differential equations about its center of mass (CoM) were first published by Euler, and the Moon was treated as a perfectly rigid body. However, due to the complexity of the forcing terms, they were too difficult to solve exactly in an analytical way (Cappallo et al. 1981). The advent of digit computers in the 20th century changed this situation greatly. The solution of these complex equations of motion via numerical methods has became a feasible alternative to the analytical theories.

Since exploration of the Moon commenced during the Apollo era, seismological profiles have been obtained (Weber et al. 2011; Garcia et al. 2011), but left the structure near the core uncertain (Matsuyama et al. 2016; Matsumoto et al. 2015). Besides the seismic equipment, reflectors were also installed on the lunar surface. Over the years, Lunar Laser Ranging (hereafter, LLR) has benefited from various improvements both in data modeling and observing technology, which enabled the current accuracy of a few centimeters (Viswanathan et al. 2018; Pavlov et al. 2016; Yagudina 2009). Its wealth information and their analysis over the years provided a new effective method to study the dynamics and structure of the Moon.

On the other hand, during the second phase of the Chinese Lunar Exploration Program, the Chinese National Space Administration launched an unmanned lunar exploration mission called Chang'E-3 on 2013 December 14. It landed on the lunar surface and is located at about 44.12 ° N, 19.15° W with an elevation of −2640 (Tang et al. 2017), becoming the only spacecraft that soft-landed on theMoon after the Cold War.

During the Chang'E-3 extensional phase mission, this lander not only conducted several experiments for characterizing the Moon's surface, but also set up a radio transponder powerful enough to maintain contact with the ground. This Lunar Radio Data (hereafter LRD) in Xband with sub-millicycle accuracy was firstly implemented based on the planetary radio science receiver (Ping et al. 2014). Every full moon (when there is a lack of LLR dataset), the lander is observed continuously for about three to five hours per day during about seven days by two Chinese Deep Space Network (CDSN) ground stations, located at Jiamusi and Kashi. After applying the necessary corrections, such as Earth's atmospheric aberration, ionospheric disturbance and equipment thermal noise, high accuracy range and range ratemeasurements of the Chang'E-3 lander (root mean square, RMS) of about 80 cm and 0.2 mm sec⁻¹ were performed (Ping et al. 2017). The position of Chang'E-3 on the lunar surface is between the US reflector A15 and Russian L1 (see in Fig. 1(a) and 1(b)), which should be helpful to improve the spatial coverage of the data samples if we combine the LLR and LRD observation data together to study lunar rotation.

Zoom In Zoom Out Reset image size

The high precision of current ranging data imposes exacting demands upon theoreticians in this field (Williams 2018). It has become necessary to includemany small, previously ignorable effects in models of lunar rotation. For example, the gravity fields ofMoon must be described to a higher order, a more detailed structural model of the Moon should be refined, and the interaction of their figures between Earth and Moon should be considered.

Besides the Earth, the Moon is the only celestial body for which seismic data have been successfully employed to infer its internal structure. But unfortunately, the deepest regions which can be sounded by the travel time data from Apollo seismic data were limited to depths of ∼1378 km depth for P wave and ∼1138 km for S wave, leaving the structure near the core uncertain (Matsumoto et al. 2015). With the help of LLR, the Moon was modeled as a twolayer body with a molten core (Williams et al. 2001) and applied in the development of NASA's Planetary and Lunar Ephemeris 421 (DE421) (Folkner et al. 2009).

The Earth-figure torque was originally referred to by Breedlove (1977) as figure-figure interaction. This torque results from the interaction of the 2^nd zonal harmonic (J₂) of the gravity field of the Earth with the complete 2^nd degree gravity field of the Moon. Since the calculation of this torque needs the position of Earth's pole in the inertial system, we should acquire this position as accurately as possible. The transformation between Geocentric Celestial Reference System (GCRS) and International Terrestrial Reference System (ITRS) can be done strictly in accordance with International Astronomical Union (IAU) resolutions. IERS Conventions 2010 (Petit & Luzum 2010) provide a detailed description of the algorithms for the transformation. Besides the well-modeled transformation procedure, there are also unmodeled minor offsets due to the complexity of Earth's rotation. They can only be corrected by observation. Unmodeled offsets in the Earth's north pole, dX(t) and dY(t), were adopted from the published IERS C04 solution (Pavlov et al. 2016; Bizouard & Gambis 2009), which are combined from Global Positioning System (GPS), Satellite Laser Ranging (SLR), and Very Long Baseline Interferometry (VLBI) observational data, and observations from the Quasar network are also included (Pavlov et al. 2016; Finkelstein et al. 2008).

In the ephemeris integration, only long-term change in the orientation of Earth is modeled. The model of Earth's orientation implemented in the ephemeris integration is a combination of the IAU 1976 precession model with an estimated linear correction and a modified IAU 1980 nutation model including only an 18.6 yr period (Folkner et al. 2014). In previous literature (Pavlov et al. 2016), Earth Orientation Parameters (EOP) must be a central element that is considered in the LLR data reduction. However, they were not taken into account during the ephemeris integration, and this discrepancy motivates us to explore this factor.

In this paper, the rotational equations of the Euler angles are applied to describe the Moon's orientation with respect to (w.r.t.) the inertial coordinate system (Yang et al. 2018; Viswanathan et al. 2017; Pavlov et al. 2016; Folkner et al. 2014; Cappallo et al. 1981). Based on the integration approach implemented independently, the difference between a one-layer and two-layer lunar model and precise position of Earth's pole in calculating the figure interaction between Earth and Moon is analyzed. This work will help us to understand lunar rotationmore clearly and lay a foundation for our continuedwork of combining LRD and LLR data together in the near future, with the prospect making a distinct contribution to study of the Earth-Moon system and also to evaluate the impact of adding new types of observations on the accuracy of lunar ephemerides.

The outline of this article is as follows: In Section 2, the numerical model of lunar rotation is described. The Euler angle results based on different lunar inner models and Earth orientation are calculated in Section 3. Some brief summaries are then provided to conclude this paper in the last section.

In this work, our numerical model is described by Euler's rotational equations which were also often used in investigating the rotation of other terrestrial planets (van Hoolst 2007). The origin of the coordinates is at the lunar CoM. As the Moon is a not perfectly rigid body, we consider three coordinate axes to be along the axes of principal moments of inertia of the mantle and the instantaneous angular velocity ω utilized to describe it in an inertial reference system. We refer to Yang et al. (2018, 2017) for full descriptions. As displayed in Figure 2, OX_i (i = 1, 2, 3) is defined, where OX₁ is along the minimum moment of inertia axis and OX₃ is the maximum moment of inertia axis. The OX₁ axis is along the positive direction of Mean Earth; OX₃ points to the positive direction of the lunar spin pole. The International Celestial Reference System (ICRS) is adopted as our inertial reference system Oξ_i (i = 1, 2, 3). The relations between these two systems are defined by three angles ϕ, θ and ψ, i.e. the Euler angles. Here ϕ is the angle along the Mean Earth Equatorial plane, θ is the angle between the Mean Earth plane and the lunar equator plane, i.e. the precession and nutation angle respectively, and ψ is the rotation angle along the lunar equator plane (see Fig. 2). The components of the angular velocity vector ω in the principal axes (PA) system can be described by Euler angles (Goldstein 2011; Danby 1992)

$$\begin{eqnarray}&&\left[\begin{array}{l}{\omega }_{1}\\ {\omega }_{2}\\ {\omega }_{3}\end{array}\right]=\left[\begin{array}{ccc}\sin \theta \sin \psi & \cos \psi & 0\\ \sin \theta \cos \psi & -\sin \psi & 0\\ \cos \theta & 0 & 1\end{array}\right]\left[\begin{array}{l}\dot{\phi }\\ \dot{\theta }\\ \dot{\psi }\end{array}\right].\end{eqnarray} \tag{ 1 }$$

Zoom In Zoom Out Reset image size

Multiply the two sides of Equation (1) by the inverse matrix of the first term on the right side of Equation (1), then differentiate Equation (1) w.r.t. time and we can derive the following second order differential equations for the Euler angles

$$\begin{eqnarray}&&\left[\begin{array}{l}\ddot{\phi }\\ \ddot{\theta }\\ \ddot{\psi }\end{array}\right]=\left[\begin{array}{ccc}\displaystyle \frac{\sin \psi }{\sin \theta } & \displaystyle \frac{\cos \psi }{\sin \theta } & 0\\ \cos \psi & -\sin \psi & 0\\ 0 & 0 & 1\end{array}\right]\left[\begin{array}{l}{\dot{\omega }}_{1}\\ {\dot{\omega }}_{2}\\ {\dot{\omega }}_{3}\end{array}\right]+\left[\begin{array}{c}\displaystyle \frac{\dot{\theta }(\dot{\psi }-\dot{\phi }\cos \theta )}{\sin \theta }\\ -\dot{\phi }\dot{\psi }\sin \theta \\ -\ddot{\phi }\cos \theta +\dot{\phi }\dot{\theta }\sin \theta \end{array}\right].\end{eqnarray} \tag{ 2 }$$

$\mathop{\omega }\limits^{.}$, through the theory of Euler's equations, depends on the total torques Γ from external bodies. The relation governing the rate of change of angular momentum of a body in an inertial coordinate system can be written succinctly as

$$\begin{eqnarray}&&\displaystyle \frac{d}{dt}({\rm{I}}\omega )=\Gamma,\end{eqnarray} \tag{ 3 }$$

where I represents the moment of inertia tensor, Γ is the torques from external bodies and t is time. But in a frame rotating at angular velocity ω, the equation should be written in the following form

$$\begin{eqnarray}&&\displaystyle \frac{d}{dt}({\rm{I}}\omega )+\omega \times {\rm{I}}\omega =\Gamma .\end{eqnarray} \tag{ 4 }$$

Combining the mathematical relationship (Eqs. (1) and (2)) between ω and Euler angles together with the associated physical law (Eqs. (3) and (4): principal axis theorem), we can construct the lunar rotational numericalmodel. In these formulas, the key question is how to model the inertia tensor I and total torques Γ from external bodies.

The Moon is treated as an anelastic body with a liquid core. The core is constrained by the shape of the coremantle boundary (CMB) at the interior of the mantle and assumed to be rotating like a solid (Pavlov et al. 2016). The effect of the existing liquid core will be investigated in the next section. The contributions to total torques Γ contained in our integration of these equations include all planets in the solar system and Sun, with each of these bodies in the lunar gravitational field being modeled as a point mass. The additional figure interaction between Earth and Moon is a torque corresponding to the Earth's oblateness.

The rotation of the Moon is strongly coupled to the orbital motion of the Moon, so called spin-orbit synchronism. Any deviation from uniform orbital motion that is not accompanied by a similar excursion from a uniform rotation motion state will lead to a torque which forces lunar physical libration, i.e. influence the Euler angles. For this reason, the orbital model applied in generating a description of the lunar rotation is of central importance. We present the following accelerations perturbing the Moon's orbit in our numerical model (Yang et al. 2018):

• Interactions among the Moon, planets, asteroids and the Sun, all of which are treated as point masses (Folkner et al. 2014) (eq. (27)). i.e. the Parametrized Post-Newtonian (PPN) effects.
• Considering the Earth's figure, and the Moon is considered as a point mass;
• Considering the theMoon's figure, the influences from Earth, Sun, Venus, Mars and Jupiter are included;
• Considering the effects from the distorted part of Earth (raised by the Sun and theMoon) and the Moon (raised by the Sun and the Earth) from Folkner et al. (2014) (eq. 32).

In this section, based on our independently realized numerical model and integrated approach to lunar motion (Yang et al. 2018), we analyze the influence of the layered models and precise Earth's pole position on lunar rotation, respectively.

3.1. One-layer and Two-layer Models of the Moon

To understand the contributions of the liquid core to the rotation of the Moon, i.e. the CMB effect, we compare the Euler angles calculated from the one-layer model and twolayer model. Firstly, based on previous study (Newhall & Williams 1996; Cappallo et al. 1981), the moment of inertia and equations of rotational motion when the Moon is considered as a one-layer model can be written as

$$\begin{eqnarray}&&\begin{array}{l}{\rm {I}}={\rm {I}}_{{\rm{rigid}}}+{\rm {I}}_{{\rm{tidal}}}+{\rm {I}}_{{\rm{spin}}}\\ \mathop{\omega }\limits^{.}={\rm {I}}^{-1}\{\displaystyle \sum _{i}{\Gamma }_{{\rm{fig-pm}}}+{\Gamma }_{{\rm{fig-fig}}}-\mathop{\rm {I}}\limits^{.}\omega -\omega \times {\rm {I}}\omega \}\end{array}\end{eqnarray} \tag{ 5 }$$

where I_rigid, I_tidal and I_spin represent the rigid-body component, the effect of tidal deformation from the Earth and lunar deformation due to lunar rotation (Newhall & Williams 1996), respectively Γ_fig–pm is the torque acting from the point-mass in the lunar gravitational field and Γ_{fig – fig} is the torque due to the Earth's oblateness.

As to the two-layer model (i.e. mantle + fluid core) (Folkner et al. 2014; Viswanathan et al. 2017; Pavlov et al. 2016), they are

$$\begin{eqnarray}&&\begin{array}{l}{\rm {I}}={\rm {I}}_{{\rm{rigid}}}+{\rm {I}}_{{\rm{tidal}}}+{\rm {I}}_{{\rm{spin}}}-{\rm {I}}_{c}\\ \mathop{\omega }\limits^{.}={\rm {I}}^{-1}\{\displaystyle \sum _{i}{\Gamma }_{{\rm{fig-pm}}}+{\Gamma }_{{\rm{fig-fig}}}+{\Gamma }_{{\rm{cmb}}}-\mathop{\rm {I}}\limits^{.}\omega -\omega \times {\rm {I}}\omega \}\end{array}\end{eqnarray} \tag{ 6 }$$

where Γ_cmb represents the torque on the mantle due to the interaction between mantle and core. It should be noticed that I and ω represent the moment of inertia and angular velocity for the lunar mantle, respectively. But different with the I and ω in Equation (5) where they mean the moment of inertia and angular velocity for the whole Moon respectively, I_c is the fluid core's moment of inertia. When the fluid core is treated in the mantle frame, it can be expressed as a constant matrix as (Folkner et al. 2014)

$$\begin{eqnarray}{\rm {I}}_{c}={\alpha }_{c}C\left[\begin{array}{ccc}1-{f}_{c} & 0 & 0\\ 0 & 1-{f}_{c} & 0\\ 0 & 0 & 1\end{array}\right]\end{eqnarray} \tag{ 7 }$$

where α_c is a dimensionless coefficient for the ratio of core to total polar moments of inertia, C is the rigid polar moment of inertia and f_c is the oblateness of the fluid core.

The Γ_cmb and the core's angular velocity ${\mathop{\omega }\limits^{.}}_{c}$ (referring to lunar mantle) can be written as

$$\begin{eqnarray}&&\begin{array}{l}{\Gamma }_{{\rm{cmb}}}=C\left[\displaystyle \frac{{k}_{\nu }}{C}({\omega }_{c}-\omega )+{\alpha }_{c}{f}_{c}(\hat{z}\cdot {\omega }_{c})(\hat{z}\times {\omega }_{c})\right]\\ {\mathop{\omega }\limits^{.}}_{c}={({\rm {I}}_{c})}^{-1}[-\omega \times {\rm {I}}_{c}{\omega }_{c}-{\Gamma }_{{\rm{cmb}}}]\end{array}\end{eqnarray} \tag{ 8 }$$

where $\displaystyle \frac{{k}_{\nu }}{C}$ is a friction parameter and $\hat{\rm {z}}=(0,0,1)$. In this paper, $\displaystyle \frac{{k}_{\nu }}{C}$ and f_c refer to DE430 (Folkner et al. 2014), fixed as 1.6365616533709530 × 10⁻⁸ and 2.4623904789198150 × 10⁻⁴, respectively.

Based on Equations (5)–(8), we can compare the difference in the Euler angles from the two models of the Moon, the one-layer model and two-layer model. Figure 3 shows the difference in the Euler angles from the two models for 250 yr. The precession angle ϕ and nutation angle θ exhibit a trend of −1.99 × 10⁻⁶ arcsec d⁻¹ and −2.74 × 10⁻⁷ arcsec d⁻¹, respectively. The difference in rotational angle ψ is depicted by a mirror image of ϕ, with a trend of about 1.82 × 10⁻⁶ arcsec d⁻¹. To analyze the difference between the two sets of Euler angles, frequency analysis was applied to Δ ϕ, Δ θ and Δ ψ. The main frequencies are listed in Table 1. As seen in Figure 3, the most obvious period of the three angles is 18.6 yr, exactly the signal of Γ_cmb, i.e. the core (Williams et al. 2001). On the other hand, in general, we need to re-fit the observation data such as LLR data when the physics of the model are changed to get the new initial conditions that are consistent with the lunar free librations. However, the fitting procedure is not included in this work, which means that we need to identify the free librations. Based on our frequency analysis, we identified 1056.12 d free libration in longitude, and 81 yr and 27.3 d free librations in pole position mode (in inertial space). In particular, the 1056.12 d free libration in longitude is even visible in Δ ψ of Figure 3. The 0.2'' in Δθ (latitude) will bring about a 1.2 m range change for Chang'E-3. It implies a potential contribution by using radio data from the Chinese Chang'E project to study the influence of the lunar fluid core (e.g. the 18.6 yr oscillation in Fig. 3).

Zoom In Zoom Out Reset image size

Table 1. Frequency Analysis of Δ ϕ, Δ θ and Δ ψ

	Δϕ	Δθ	Δψ
1	18.61 yr	18.59 yr	18.61 yr
2	27.29 d	27.29 d	27.29 d
3	-	-	1056.12 d
4	-	81.26 yr	-

3.2. EOP in the Earth-Moon Figure Torque

The distribution of mass within the Earth deviates from spherical symmetry enough so that the lowest-order effects of the Earth's figuremust be accounted for in calculation of the torque acting on the Moon. Yoder (1979) and Eckhardt (1982) demonstrated that the three most significant terms in the torque can be calculated from

$$\begin{eqnarray}&&\begin{array}{l}{\Gamma }_{{\rm{fig-fig}}}=\displaystyle \frac{15{\mu }_{e}{ {\mathcal R} }_{e}^{2}{J}_{2e}}{2{r}_{e}^{5}}\{(1-7{\sin }^{2}\eta )[{\hat{\rm {r}}}_{e}\times {\rm {I}}{\hat{\rm {r}}}_{e}]\\ +2\sin \eta [{\hat{\rm {r}}}_{e}\times {\rm {I}}{\hat{\rm {P}}}_{e}+{\hat{\rm {P}}}_{e}\times {\rm {I}}{\hat{\rm {r}}}_{e}]-\displaystyle \frac{2}{5}[{\hat{\rm {P}}}_{e}\times {\rm {I}}{\hat{\rm {P}}}_{e}]\}\end{array}\end{eqnarray} \tag{ 9 }$$

where J_2e is the Earth's zonal harmonics, ${ {\mathcal R} }_{e}$ is Earth's equatorial radius, ${\hat{\rm {r}}}_{e}$ is the normalized direction vector of the Earth from the Moon, ${\hat{\rm {P}}}_{e}$ is the normalized direction vector of the Earth's north pole and $\sin \eta ={\hat{\rm {r}}}_{e}\cdot {\hat{\rm {P}}}_{e}$. All the vectors in the equation are described in the lunar frame

From Equation (9) we can see that the direction of the Earth's pole ${\hat{\rm {P}}}_{e}$ in the inertial reference system plays an important role in calculating the figure-figure torque on lunar rotation. To analyze its contributions to lunar rotation, we first need to know the precise orientation of the Earth. The orientation of the Earth is defined as the rotation between a rotating geocentric system connected to the ITRF (the terrestrial system realized by observatories on Earth) and a non-rotating geocentric system connected to the ICRF (the celestial system realized by coordinates of quasars, stars and also objects in the solar system). The transformation matrix between the two systems is usually utilized to describe the Earth's rotation (Petit & Luzum 2010).

To analyze the influence of the unmodeled pole offset during the figure torque calculation on lunar rotation, firstly, we compare the position of Earth's pole from these two models. Figure 4 shows the Earth's pole position in the XY-plane of the inertial frame from the IERS suggested model (labeled as 'IauPole') and the ephemeris model used (labeled as 'EphPole'). Moreover, the difference between the angle of ${\hat{\rm {r}}}_{e}$ and ${\hat{\rm {P}}}_{e}$ is plotted in Figure 5, which demonstrates that the difference between these two models is less than 1''.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We compare the Euler angle result, with the IAU solution when calculating the figure-figure torque, to the Euler angle result from the ephemeris DE430 model. Figure 6 exhibits the difference between these two results.

Zoom In Zoom Out Reset image size

From Figure 6 we can find that the model of the transformation algorithms performs well in Earth orientation transformation from TRS to CRS. The maximum discrepancy between cases with C04 and those without C04 is only 2 × 10⁻⁷ arcsec, far smaller than not only the current precision of numerical integration (50 mas difference between DE430, EPM2017 and INPOP17a), but also the current RMS of observation minus calculation (O – C, 1 cm, about 3 mas in Δ θ) of LLR. In this situation, we can neglect its contribution to the calculation of the torque for figure-figure interaction. In addition, frequency analysis of the data found a period of 0.5 yr, which is exactly the largest period of Earth's nutation except for 18.6 yr. The difference in the free libration modes between the two models is too small to be detected.

Lunar rotation is linked to its interior structure. Deviation from uniform rotation motion will result in torque, forcing lunar physical librations because of the spin-orbit synchronism (Petrova et al. 2017; Yang et al. 2017; Rambaux & Williams 2011). For that reason, understanding the rotation of the Moon and elements influencing these processes is crucial for studying the lunar internal structure and its composition (Wieczorek 2006; Williams et al. 2001). The internal structure and composition of the Moon that accounts for the formation and evolution of such bodies makes the study of lunar rotation a hot research topic.

There is still considerable uncertainty about the size and properties of the lunar core (Williams et al. 2006). In this work, we studied the differences between a twolayer and one-layer lunar model in lunar rotational motion. The discrepancies between these two models confirm that the core structure of the Moon can be inferred from rotation data. A recently published paper developed a three-layer lunar rotation model including a solid inner core (Dumberry & Wieczorek 2016), and will contribute to lunar rotation studies; however, the feasibility of implementing this model into numerical calculations must be investigated.

The development of observation techniques and data reduction methods enables possibly observing a highly accurate unmodeled Earth orientation offset, stimulating our interest in the figure-figure torque from the Earth's shape on lunar rotation. We calculated the difference in our numerical integration framework, and found that the influence of this factor is negligible due to present observing ability and precision. Euler angle differences with and without the IERS C04 data are much smaller than the observational noise. Future lunar in-situ geodetic experiments, such as photographic zenith tubes, telescopes, active lasers or laser transponders on the lunar surface, would experience the same shortcoming.

The factors influencing the Moon's rotation are numerous and complex. In this paper, we compared two factors involved in lunar rotation to gain a clearer understanding of how these factors influence lunar rotation. It would be helpful to construct a more realistic Earth-Moon model in the future when employing LRD observation data together with LLR to infer lunar properties. This work is part of our ongoing work to build an independent processing platform for re-processing all radio tracking data from the Chinese Lunar Exploration Project, combining the radio data and existing laser ranging data together. This combination will improve the understanding of constraints on the lunar core. This represents a new way to study the Moon, and will provide reference for further radio science investigations during Chinese Deep Space Exploration missions.

The EOP 08 C04 solution is available at ftp://hpiers.obspm.fr/iers/eop/eopc04_08/. This research is supported by LIESMARS Special Research Funding, by the National Natural Science Foundation of China (U1831132, 41590851, 11373060, 10973030 and 10778635), the State Key Project for Science and Technology (2015CB857101), National Astronomical Observatories, Chinese Academy of Sciences, a grant from the Hubei Province Natural Science (2018CFA087), Open Project of Lunar and Planetary Science Laboratory, Macau University of Science and Technology (FDCT 119/2017/A3), and Open Funding of Guizhou Provincial Key Laboratory of Radio Astronomy and Data Processing (KF201813). We thank the anonymous referees for their valuable comments that improved the manuscript.