A&A 428, 255-260 (2004)
DOI: 10.1051/0004-6361:20035952

Effects of zonal deformations and the Earth's rotation rate variations on precession-nutation

S. Lambert1 - N. Capitaine2

1 - Department of the Navy, US Naval Observatory, Earth Orientation Department, 3450 Massachusetts Avenue NW, Washington DC 20392, USA
On leave from 2
2 - SYRTE - UMR 8630/CNRS - Observatoire de Paris, 61 avenue de l'Observatoire, 75014 Paris, France

Received 29 December 2003 / Accepted 9 August 2004

Because of their importance in the accurate modeling of the Earth's orientation in space, some non-negligible predictable effects on precession-nutation are investigated. This paper considers the coupling effects between the axial and the equatorial components of the Earth's rotation vector in the dynamical equations, and the effects of the second order lunisolar torque due to the Earth's zonal deformations. Firstly, the coupling effects are shown to contribute for less than 0.1 $\mu $as and are therefore negligible. Secondly, we demonstrate that the 0.7 mas contribution of the rotation rate variations due to zonal tides to the nutation in obliquity deduced by Bretagnon et al. (2000, Proc. IAU Coll., 180, 230; 2001, Celest. Mech. Dyn. Astr., 80, 177) is an artefact which comes from an incomplete way of taking into account the effect of the rotation rate variations. The net contribution is shown to be negligibly small. Thirdly, for an Earth model with an elastic mantle and decoupled liquid core, the contribution of the second-order lunisolar torque due to the Earth's zonal deformations is shown to be 207.9 $\mu $as and -9.7 $\mu $as on the 18.6-year nutation respectively in longitude and in obliquity, and a correction of -4925.9 $\mu $as/century on the precession in longitude.

Key words: reference systems - astrometry - time

1 Introduction

Very long baseline radio interferometry (VLBI) observations provide on a regular basis "celestial pole offsets'' that contain errors in the models for the position of the celestial pole in the celestial reference system. The most accurate available precession-nutation model is the IAU 2000A model recommended by the International Astronomical Union (IAU Resolution B1.8) adopted at the IAU 24th General Assembly (2000) and implemented by the IERS beginning on 1 January 2003.

Comparison of this new model against VLBI observations shows differences of the order of 200 $\mu $as. These differences are the consequence of various influences as well as geophysical processes that are summarized below. The main oscillation appearing in these residuals comes from the Free Core Nutation (FCN). It gives rise to a 430-day periodic term, the amplitude of which varies strongly and is therefore not predictable but can be described using an empirical model (Herring et al. 2002). As the nutation modeling contains only permanent predictable effects, the FCN is not included in IAU 2000A and therefore appears in the residuals. Other unpredictable effects come from the variable atmospheric and oceanic forcing which induces time varying nutation amplitudes of about 80 $\mu $as for the annual term (Bizouard et al. 1998). Some theoretical predictable effects are taken into account in the IAU 2000A precession-nutation model whereas they were not included in the previous models. These effects include:

Further discussions of these effects are still necessary to clarify how they can be related to previous studies. Coupling effects between the Earth's rotation rate and precession-nutation were neglected because their contributions are expected to be negligible compared to the accuracy of the observations. However, recent studies concluded that variations in Earth rotation due to zonal deformations were responsible for contributions of the order of 712 $\mu $as to the 18.6-year nutation in obliquity and a change of the precession in longitude of 4 mas per century (Bretagnon et al. 2000, 2001).

The deformations due to zonal tides have also direct effects on precession and nutation. The effect on nutation was first studied by Souchay & Folgueira (2000). It was concluded that, for a basic Earth model, this should modify the amplitude of the 18.6-year nutation in longitude by 168 $\mu $as, perturbations on other nutation terms all being less than 10 $\mu $as. Another work by Krasinsky (1999) showed larger corrections based on a computation of the gravitational torque exerted on the Earth. The largest correction provided is 5.11 mas on the 13.66-day term, whereas the semi-annual term and the 18.6-year term are changed respectively by 2.38 mas and 1.07 mas. Note that VLBI observations do not show such large discrepancies in the nutation model.

Considering such a large disparity in the amplitudes of the contributions to nutation due to zonal tides that are reported in various studies, a rigorous computation of these effects is necessary. In this paper, we compute the effects due to the torque produced by the action of the tesseral potential on the deformations produced by the zonal potential together with the effects on the rotation rate produced by these deformations. In Sect. 2 we recall the expressions of the Euler-Liouville equations developed up to the second order, emphasing the couplings between variations of the axial component and the equatorial component of the instantaneous rotation vector. In Sect. 3 we introduce the second order part of the lunisolar torque coming from Earth's zonal deformations, and in Sect. 4, we compute all the predictable effects due to zonal tides on the Earth's precession-nutation.

2 Dynamical equations of the Earth's rotation up to the second order

We consider a deformable Earth with fluid core. In a rotating frame with its axes oriented towards the Earth's mean axes of inertia, the angular momentum equation is:

\frac{{\rm d}\vec H}{{\rm d}t}+\vec\omega\times\vec H=\vec\Gamma
\end{displaymath} (1)

where $\vec H$ is the whole Earth's angular momentum vector, corresponding to the product of the inertia tensor $\underline{\overline I}$:

\underline{\overline I}=\left(\begin{array}{ccc}
A & 0 & 0 ...
..._{22} & c_{23} \\
c_{31} & c_{32} & c_{33}
\end{displaymath} (2)

wherein second-order quantities cij characterize the effects of the deformations, with the Earth's instantaneous vector of rotation $\vec\omega$:

m_1 \\
m_2 \\
\end{displaymath} (3)

$\Omega$ being the mean Earth's rotation rate.

$\vec\Gamma$ is the tidal torque expressed in the terrestrial frame. The equatorial component $\Gamma=\Gamma_1+{\rm i}\Gamma_2$ of the torque is computed from the expression $\Gamma'$ in the celestial frame by the relation:

\Gamma=\Gamma'{\rm e}^{-{\rm i}\Phi}
\end{displaymath} (4)

and the Euler's kinematical relations:
$\displaystyle \dot\theta+{\rm i}\dot\Psi\sin\theta$ = $\displaystyle -\Omega m{\rm e}^{{\rm i}\Phi}$ (5)
$\displaystyle \dot\Phi+\dot\Psi\cos\theta$ = $\displaystyle \Omega(1+m_3)$  

relate the terrestrial motion of the Earth's axis of rotation $m=m_1+{\rm i}m_2$ and the relative variations of the Earth's rotation rate m3 to the Euler's angles $\theta$, $\Psi$ and ${\Phi }$ between the terrestrial frame and the mean ecliptic frame of the epoch of reference.

Using the equations of Sasao et al. (1980), one obtains the following expressions:

$\displaystyle \dot m-{\rm i}\sigma_{\rm r}m+\frac{\dot c+{\rm i}\Omega c}{A}+\frac{A_{\rm f}}{A}(\dot m_{\rm f}+{\rm i}\Omega m_{\rm f})$     (6)
$\displaystyle -{\rm i}\sigma_{\rm r}mm_3-{\rm i}\frac{\Omega}{A}(A_{\rm f}m_{\rm f}m_3-C_{\rm f}m_{{\rm f},3}m)$      
$\displaystyle -{\rm i}\frac{\Omega}{A}cm_3+{\rm i}\frac{\Omega}{A}c_{33}m$ = $\displaystyle \frac{\Gamma}{A\Omega}$  
$\displaystyle \dot m+\dot m_{\rm f}-{\rm i}\sigma_{\rm f}m_{\rm f}+{\rm i}\Omega m_{\rm f}+\frac{\dot c_{\rm f}}{A_{\rm f}}$ = 0                                                     

in which $c=c_{13}+{\rm i}c_{23}$, $c_{\rm f}=c_{{\rm f},13}+{\rm i}c_{{\rm f},23}$, $\sigma_{\rm r}=\Omega(C-A)/A$ and $\sigma_{\rm f}=\Omega(C_{\rm f}-A_{\rm f})/A_{\rm f}$, quantities subscripted by f being relevant to the core.

The previous equations reflect two different couplings between the axial and the equatorial component of the rotation vector of the Earth. One coupling comes from the products $m\times m_3$ or $m_{\rm f}\times m_3$ in the left hand side of (6). Another coupling appears through the sidereal rotation angle ${\Phi }$ depending upon m3, in Eqs. (4) and (5) in the opposite sense to transform coordinates between the terrestrial to the celestial frame.

Variations of Euler's angles are related to the instantaneous vector of rotation by Euler's kinematical relations (5). One gets the nutation angles referred to the ecliptic of epoch:

$\displaystyle %
\Delta\psi_1\sin\epsilon_1$ = $\displaystyle \Delta\Psi\sin\theta$ (7)
$\displaystyle \Delta\epsilon_1$ = $\displaystyle -\Delta\theta.$  

Note that the periodic part of the nutation angles referred to the ecliptic of date ( $\Delta\psi$ and $\Delta\epsilon$) is the same as the one of $\Delta\psi_1$ and $\Delta\epsilon_1$, and that the precession is expressed along the fixed ecliptic. The contribution to the precession $\psi_A$ is equivalent to the secular term of $\Delta\psi_1$.

3 Torque including Earth's zonal deformations

The tidal gravitational torque exerted on the Earth is derived from the tesseral part of the lunisolar potential $\phi$ in the celestial reference frame (Sasao et al. 1980):

                                 $\displaystyle \Gamma'$ = $\displaystyle \Gamma^{\prime(1)}+\Gamma^{\prime(2)}$ (8)
  = $\displaystyle -{\rm i}\Omega^2(C-A)\phi-{\rm i}\Omega^2(c_{33}-c_{11})\phi$  

where the superscripts (1) and (2) stand for the first and second order parts of the torque respectively. The second order part includes the deformations c33 and c11 of the Earth along the directions of the mean axes of inertia (cf. relation (2)).

The complex-valued quantity $\phi$ depends on the right ascension $\alpha$ and declination $\delta$ of the perturbing body:

\phi=\frac{3K^2\mu}{\Omega^2r^3}\sin\delta\cos\delta(\cos\alpha+{\rm i}\sin\alpha)
\end{displaymath} (9)

in which K2 is the geocentric gravitational constant, $\mu $ is the ratio between the mass of the disturbing body and the mass of the Earth and r is its geocentric distance.

Table 1: Excess of LOD ($\mu $s) due to zonal tides from IERS Conventions 2003, relative variations of the Earth's rotation rate m3z (10-12 rad) and c33/C ($\mu $as) from Eqs. (11) and (12).

Assuming as generally (see e.g. Melchior 1978; Yoder et al. 1981) that the trace of the inertia tensor is constant, one has in the case of zonal deformations:

\end{displaymath} (10)

The variations of the dynamical ellipticity are related to corresponding changes in the excess of length-of-day (LOD) and to the axial component m3 of the instantaneous vector of rotation of the whole Earth:

m_3^z=-\frac{\delta LOD^z}{\overline{LOD}}
\end{displaymath} (11)

where the superscript z indicates that these rotation rate variations come from zonal tides and $\overline{LOD}=86~400$ s. According to Mathews et al. (2002) and Mathews (2004), the increment of inertia of the whole Earth produced by the deformations due to zonal tides is:
               c33 = $\displaystyle -C_{\rm eff}m_3^z$ (12)
$\displaystyle C_{\rm eff}$ = $\displaystyle \frac{C_{\rm m}}{1-\gamma C_{\rm f}/\kappa C}$  

$\gamma$ and $\kappa$ being undimensional parameters representing the deformabilities of the mantle and the core-mantle boundary under degree 2 tidal forcing (see e.g. Sasao et al. 1980; Mathews et al. 1991) and $C_{\rm f}$ the axial moment of inertia of the core.

The second order torque due to the zonal deformations is finally expressed using quantity c33 or equivalently m3z:

                                  $\displaystyle \Gamma^{\prime(2)}$ = $\displaystyle -{\rm i}\frac{9K^2\mu}{2r^3}c_{33}\sin\delta\cos\delta(\cos\alpha+{\rm i}\sin\alpha)$ (13)
  = $\displaystyle {\rm i}\frac{9K^2\mu C_{\rm eff}}{2r^3}m_3^z\sin\delta\cos\delta(\cos\alpha+{\rm i}\sin\alpha).$  

Direct computation of the changes in the dynamical ellipticity can be based on expression for c33 corresponding to a given Earth model. For example Souchay & Folgueira (2000) used a model provided by Melchior (1978) on the basis of a simplified elastic Earth model. Computing the effect on nutation with a microarcsecond accuracy requires a more sophisticated Earth model. The effect of the Earth model appears through the zonal response coefficient that can be estimated from observed Earth Orientation Parameters time series (Hefty & Capitaine 1990). The model of zonal variations in Earth's rotation of Yoder et al. (1981) considered an equilibrium ocean tide and was adopted in the IERS Conventions 1996 (McCarthy 1996). A more refined model provided by Defraigne & Smits (1999) was adopted in the IERS Conventions 2003. It includes a decoupled core and an inelastic mantle, and ocean effects included using a transfer function worked out by Mathews et al. (2002).

Our computation of c33 uses this model for LOD variations and the mean moments of inertia of the whole Earth provided by the IERS Conventions 2003:

$\displaystyle C=8.037\times 10^{37}\;{\rm kg}\;{\rm m^2}$     (14)
$\displaystyle A=8.011\times 10^{37}\;{\rm kg}\;{\rm m^2}$      
$\displaystyle C_{\rm m}=7.040\times 10^{37}\;{\rm kg}\;{\rm m^2}.$      

Using the compliances from Mathews et al. (2002) estimates, one gets:

C_{\rm eff}=9.213\times 10^{37}\;{\rm kg}\;{\rm m^2}.
\end{displaymath} (15)

Table 1 gives the excess of the length-of-day due to tidal deformations provided by the IERS Conventions 2003 together with the variations m3z and the relative variations of C from Eqs. (11) and (12).

4 Solution of the dynamical equations

4.1 Couplings with the rotation rate m3

All the terms showing products $m\times m_3$ or $m_{\rm f}\times m_3$ in the left-hand side of Eq. (6) express coupling between the motion of the axis of rotation and the Earth's rotation rate (note that the use of a decoupled core involves the relation $m_{{\rm f},3}=-m_3$). The largest contribution of these terms to nutation has been evaluated as being lower than 0.1 $\mu $as for the 18.6-year term. This confirms that this coupling effect is negligible and it will not be considered further in the following.

4.2 The coupling with the sidereal rotation angle ${\Phi }$

The effect of taking into account the Earth's rotation irregularities in the sidereal rotation angle ${\Phi }$ in computing the solution for precession-nutation has to be investigated.

Bretagnon et al. (2000, 2001) constructed non-rigid nutation series by applying the transfer function MHB 2000 (Mathews et al. 2002) to the SMART97 rigid nutation series (Bretagnon et al. 1998). Given that the transfer function of MHB 2000 does not take into account the variations in the Earth's rotation rate, these effects were added to the third component of the Earth's rotation vector and the Euler's angles were recomputed from Euler's kinematical relationships. This method gave rise in particular to a 18.6-year term in obliquity with an amplitude reaching 712 $\mu $as.

We investigate this effect considering the equations for a rigid ellipsoidal Earth (the quantity c is null) together with variations in Earth's rotation rate. Equation (6) is reduced to:

\dot m-{\rm i}\sigma_{\rm r}m=\frac{\Gamma}{A\Omega}
\end{displaymath} (16)

and therefore:

m(t)=K_{\rm r}{\rm e}^{{\rm i}\sigma_{\rm r}t}+{\rm e}^{{\r...
...\rm e}^{-{\rm i}\Phi}{\rm e}^{-{\rm i}\sigma_{\rm r}t}{\rm d}t
\end{displaymath} (17)

where $K_{\rm r}$ is a constant depending on initial values of the problem.

The first term of solution (17) is the free motion of the rotation axis of the Earth, known as the Euler motion. The second term of Eq. (17) corresponds to the forced motion of the rotation axis. Variations of the Euler's angles are therefore:

                      $\displaystyle \dot\theta+{\rm i}\dot\Psi\sin\theta$ = $\displaystyle -\Omega {\rm e}^{{\rm i}\Phi}{\rm e}^{{\rm i}\sigma_{\rm r}t}\int...
...Gamma'}{A\Omega}{\rm e}^{-{\rm i}\Phi}{\rm e}^{-{\rm i}\sigma_{\rm r}t}{\rm d}t$ (18)
$\displaystyle \dot\Phi+\dot\Psi\cos\theta$ = $\displaystyle \Omega(1+m_3).$  

Using an iterative process, the integration of the system (18) converges rapidly after a few iterations. Note that the sidereal rotation angle ${\Phi }$ is affected by the irregularities m3 in the Earth's rotation rate, according to the second equation of the system. In the first equation of Eq. (18), the quantity ${\Phi }$ appears twice inside and outside the time integral in opposite sense. First, it is used to transform the torque $\Gamma'$ from the celestial frame to the terrestrial frame (rotation of $-\Phi$). Second, it is used in Euler's kinematical relations to transform the rotation vector expressed in the terrestrial frame into the Euler's angles with respect to the space-fixed reference frame (rotation of $+\Phi$).

By a numerical evaluation, one can see that omitting the variations in the Earth's rotation rate in the two rotations of angle plus or minus ${\Phi }$, that is to say writing $\dot\Phi+\dot\Psi\cos\theta=\Omega$ instead of $\dot\Phi+\dot\Psi\cos\theta=\Omega(1+m_3)$, leads to an error below 0.1 $\mu $as on $\theta$ and $\Psi$. The effect of the former rotation is actually mostly cancelled out by the latter so that the net effect is negligibly small.

However, omitting to take into account the Earth's rotation rate irregularities in one of the rotations of plus or minus ${\Phi }$ leads to large terms which could be interpreted as being coupling effects. For example, if the computation of the torque in the terrestrial frame neglects these irregularities, one gets a substantial contribution in the in-phase 18.6-year term of the obliquity (-729 $\mu $as) and on the 9.3-year in-phase term of the longitude (71 $\mu $as). This 18.6-year oscillation of the angle $\theta$ is due to a coupling between the precession in longitude $\Psi=\Psi_0\times t$ and the 18.6-year term in the variations of the Earth's rotation rate of amplitude a18.6:

\theta_{18.6}=\Psi_0\sin\theta_0\frac{\Omega a_{18.6}}{\sigma_{18.6}^{\prime 2}}\cos\sigma_{18.6}'t.
\end{displaymath} (19)

The precession in longitude is also affected by a spurious effect of 4 mas/c. This explains the erroneous values provided in the above-mentioned studies that followed such an approach. These "artefacts" are displayed in Table 2 together with the terms provided by Bretagnon et al. (2000, 2001) and the net effect.

Table 2: Effect of zonal variations in Earth's rotation: artefacts induced by omitting the variations when transforming the torque from the celestial frame to the terrestrial frame compared to the values of Bretagnon et al. (2001), and net contribution.

4.3 Effect of the second order zonal tides-induced torque

The computation of the effects of the second order torque induced by zonal tides on nutation is done using Eq. (6) in the frequency domain once the torque is computed according to Eq. (13). The calculation procedure of the lunisolar torque made in the present paper is based on the lunar theory ELP2000 (Chapront-Touzé & Chapront 1983) and the solar system semi-analytical solution VSOP87 (Bretagnon & Francou 1988). Poisson series manipulations have been processed using the GREGOIRE software package developed by Chapront (2003). We use standard values (IERS Conventions 2003) for other orbital parameters involved in the computation of the torque.

The computed variations of the nutations angles are the following:

(i) Secular term.
The constant term of the torque leads to a secular variation in longitude of -4926 $\mu $as/c. Note that such a contribution would be automatically included in the observed value for the precession. The obliquity is not affected.
(ii) Periodical terms.
The main second-order contribution to the periodic terms comes from the 18.6-year term of the torque which is a coupling between the 18.6-year variation of the dynamical ellipticity and the constant term of the lunisolar potential. Other contributions are below 5 $\mu $as. It is interesting to note that, according to Table 1, the 18.6-year term is in the opposite sense compared to the others. Namely, it is negative and decreases the excess of LOD (or equivalently, decreases the difference between axial and equatorial moment of inertia so that the Earth becomes less elliptic). This means that the effect of the 18.6-year variation is to concentrate the masses of the Earth around the rotation axis, so that the Earth rotation speed increases to satisfy the angular momentum conservation law, contrarily to the other contributions which are positive and make the Earth oblateness increase. A consequence of the oblateness decrease is that the torque exerted on the equatorial bulge is lowered and the 18.6-year nutational response to this excitation will consequently be lowered. The 18.6-year nutation coefficient should be decreased in absolute value. In fact, this contribution reaches 208 $\mu $as in the nutation in longitude whereas it produces only a small effect of almost -10 $\mu $as on the nutation in obliquity. Since the 18.6-year nutation amplitude in longitude is negative (-17 $^{\prime\prime}$), such a result is in the right direction.
Table 3 summarizes the results and provides a comparison between our study and other works. The difference of the order of 40 $\mu $as for the 18.6-year nutation in longitude with respect to Souchay & Folgueira (2000) is due to the more realistic Earth model considered in this study. The differences with respect to Mathews et al. (2002) need some additional comments.

Mathews et al. (2002) investigated a refined Earth model (elastic mantle and decoupled liquid core) to evaluate the influence of the so-called "non-linear'' terms in the dynamical equations in order to remove them from the observations before fitting the geophysical parameters used in the transfer function. This gave a contribution to the 18.6-year term of 94 $\mu $as and -29 $\mu $as respectively in longitude and in obliquity, including the second order torque due to both zonal and sectorial tides. Note that the differences between the values obtained in this study and those of Mathews et al. (2002) appearing in Table 3 are in fact considerably reduced with respect to revised values obtained by Mathews (2004). Excluding the sectorial part, the revised contributions in the 18.6-year nutation are 194 $\mu $as and -10 $\mu $as respectively in longitude and in obliquity (Mathews 2004). The remaining difference with respect to our values is due to anelasticity and ocean tide effects taken into account in the MHB work. It should moreover be noted that since this work has been carried out, Mathews (2003) pointed out that the above computed effects could be cancelled out by the sectorial and zonal parts of the potential acting on the deformations due to the tesseral potential. This point needs further investigation.

Table 3: Effects of the second order contribution of the tidal potential (zonal tides contribution) on nutation angles ($\mu $as) and precession ($\mu $as/c) from different studies. MHB 2000Z: values of Mathews et al. (2002) including only the zonal tides contribution; SF99: Souchay & Folgueira (2000).

5 Discussion and conclusion

Previous studies (Bretagnon et al. 2000, 2001) pointed out that the Earth's rotation rate variations due to zonal tides could have a noticeable effect on precession-nutation, with an amplitude of the order of 700 $\mu $as in the 18.6-year nutation in obliquity, although such an effect was not detected in the observations.

The present paper has shown that this effect is in fact cancelled out by taking into account the variations of the rotation rate in the expression of the torque as seen from the terrestrial frame used in the dynamical equations. This point resolves a serious dilemma of several years concerning the magnitude of these effects. Moreover, this term does not contain the effects induced by the variations of the dynamical ellipticity, although it may sometimes have been understood as being so.

We also recomputed the effects due to the coupling between axial and equatorial components of the rotation vector in the dynamical equations of the Earth rotation and we took into account the second-order contribution of the torque due to zonal deformations of the Earth's shape. These results are obtained by solving the dynamical Eqs. (6) using the model for variations in the rotation rate due to zonal tides of the IERS Conventions 2003. For an elastic Earth with a decoupled liquid core, we have concluded that:

the coupling effects between axial and equatorial components of the rotation vector do not bring any contribution larger than 0.1 $\mu $as and are therefore negligible;

the second-order torque induced by Earth's zonal deformations produces larger effects. The effect on the 18.6-year nutation is of 207.9 $\mu $as in longitude and -9.7 $\mu $as in obliquity. Other contributions are about a few microarcseconds. The correction to the precession in longitude is -4925.9 $\mu $as/c.
Note that for comparison with VLBI observations, the other second order contributions of the torque on precession-nutation as mentioned in Mathews (2003) should be considered. Although the effect of zonal tides on nutations is already included in MHB 2000 and should be completed by other second order effects (not yet published), this computation gives a more realistic model than that of Souchay & Folgueira (2000) as well as a result independent of MHB 2000 and confirms an inaccuracy in the MHB 2000 tables.

The authors are grateful to J. Souchay and P. M. Mathews for valuable discussions and to J. Chapront for having kindly provided the software package GREGOIRE. They thank the anonymous referee for very useful suggestions in improving the manuscript.



Copyright ESO 2004