% IMPORTANT NOTES:
% - figures available on request from the authors

% =====================================================================
%Definitions:
\def \A{{\rm A}}
\def \B{{\rm B}}
\def \AB{{\rm AB}}
\def \E{{\rm E}}
\def \M{{\rm M}}
\def \EB{{\rm EB}}
\def \MB{{\rm MB}}
\def \n{{\bf n}}
\def \vv{{\bf v}}
\def \etal{{\it et~al. }}
\def \O{{\cal O}}
\def \etal{{et~al. }}
\input aa.cmm
%\input aa.mtm
%\refereelayout
  \MAINTITLE={ Relativistic spin effects in the Earth-Moon system}
  \SUBTITLE={ I. Formulation}
  \AUTHOR={ E. Bois$^1$ and D. Vokrouhlick\'y$^{\,1,\,2}$}
  \INSTITUTE={ @1  Observatoire de la C\^ote d'Azur, Dept. CERGA,
                 URA CNRS 1360, Av. N. Copernic, F-06130 Grasse, France
               @2 Institute of Astronomy, Charles University,
                 \v Sv\'edsk\'a 8, CZ-15000 Prague 5, Czech Republic }
  \ABSTRACT={ An accurate theory of the Moon's spin-orbital
              motion has been constructed by numerical integration.
              Several phenomena capable of producing effects of at least
              $10^{-4}$ arc seconds in the lunar physical librations have
              been included, analysed, and for the most part presented 
              in previous papers. The present work deals with
              the relativistic spin contributions acting on the
              Earth-Moon system: (i) de Sitter precession of the
              Earth reference frame, (ii) quasi-Newtonian torques
              acting on the lunar physical librations. The latter effect 
              does not seem to have been studied yet. The global behaviour 
              of the resulting librations is presented and described.
              The resulting magnitude of these librations is of the order
              of a few $10^{-4}$ arc seconds and may reach one milliarcsecond.
              Unlike the very faint relativistic contribution to the
              Earth's rotation, the corresponding relativistic terms in
              the Moon's rotational motion are not negligible with
              respect to the present observational accuracy of the
              lunar laser data. Consequently it would be appropriate to
              consider them in models adjusted to the observations.}
  \KEYWORDS={ Celestial mechanics -- Moon -- Gravitation -- Relativity }
  \THESAURUS={ 11(05.03.1; 07.13.3; 02.08.1; 02.18.8) }
  \OFFPRINTS={Bois "bois{\at}ocar01.obs-azur.fr"}
  \DATE={ ????? }
\maketitle
%\MAINTITLERUNNINGHEAD{ Relativistic spin effects in the Earth-Moon system}
%\AUTHORRUNNINGHEAD{ E. Bois and D. Vokrouhlick\'y}
\medskip
\titlea{Introduction}
The theory of lunar motion has represented one of the most important
keystones of celestial mechanics throughout its history. Having this
status, it has often served for testing the currently accepted 
theory of gravitation. Thus it is not surprising, that we already meet
important relativistic applications to lunar motion in the first years 
of Einstein's theory.
Namely, de Sitter (1916) discussed the appearance of the combined
geodetic and Thomas precessions of the Earth reference frame in the
context of lunar motion (nowadays known as `de Sitter precession'). 
Secondly, Nordtvedt recognized the possible importance of lunar motion for 
testing his strong equivalence breaking proposal (Nordtvedt 1968, 1973). 
Finally, first attempts at lunar laser ranging (LLR) in the late sixties 
were accompanied by a partial goal consisting in the measurement of
the gravitational constant variations. Such variations occur in several
tensor-scalar gravitational theories [see a review by Mulholland (1980)].

On the other hand, it should be understood that currently the importance 
of lunar motion for testing theories of gravitation has slightly changed.
The underlying reasons are as follows:
\smallskip

\item{i)} the Earth-Moon(-Sun) system is very `dirty' from
 the relativistic point of view because of many Newtonian physical
 phenomena related to the internal structure of the bodies that are
 only poorly modelled;
\item{ii)} relativistic (PN-)parameters are determined with
 considerable precision via other solar system experiments
 (see e.g. Will 1992; the Eddingtonian $\beta$ and $\gamma$ parameters
 currently coincide with the Einsteinian values with standard errors at
 the 0.2\% level);
\item{iii)} the `true frontiers' of the experimental testing of viable
 classes of theories of gravitation have moved to systems with a 
 strong field regime, represented mainly by binary pulsars 
 (e.g. Damour \& Taylor 1992; Damour \& Esposito-Far\`ese 1992).
\smallskip

\noindent Concerning the second point, we can mention that
the Einsteinian theory has been largely verified in the weak field
regime, at least with the solar system tools currently at our disposal 
(we are aware of planned Earth vicinity or solar system
missions designed to extend the upper limits of the parameters of the
post-Newtonian theories; it however appears that currently known
constraints on those parameters are sufficiently strong for the
possibilities of the LLR experiment). In this situation, one even
faces the question: {\it should the formulation of the lunar motion still
serve for testing the Einsteinian theory of gravitation or simply
incorporate it as a valid theory?}

It has been demonstrated that there is still a place for good routine work
in testing the weak-field limit of the theories of gravitation via lunar data
(e.g. Dickey \etal 1989; M\" uller \etal 1991).
Notably, investigations of the Nordtvedt effect still possess interesting
theoretical implications (e.g. Nordtvedt 1988; Damour \& Esposito-Far\`ese
1992). Secondly, de Sitter precession of the Earth reference
frame shows practical importance well beyond the pure theory
of lunar motion.

However, independently of the opinion concerning the answer to the previous
items, one has to ask: {\it which Einsteinian terms are to be retained in
the theory for a given precision?} We recall the important methodological
fact, that the Einsteinian theory of gravitation is embodied into the
general set of post-Newtonian theories (Damour 1993; Will 1993).
Testing relativity, one always distinguishes between
null and non-null tests. In the following, we shall have 
in mind non-null tests of general relativity.

Current state of the art LLR data reach the $3$ centimeter precision
level, averaged over data sets of the last few years (Dickey \etal 1994). 
It should be however emphasized that a careful choice of the best performing 
nights over the last years reaches even better properties 
-- about $1$ to $2$ centimeters
and improvement to a few millimeter level is already envisaged
(Ch. Veillet; private communication). In terms of lunar
librations this current best data quality means the level of a few
milliarcseconds (later referred to as {\it mas}) and better can be reached
in the future.

Generally, one can divide
relativistic dynamical contributions into those operating on the 
(i) translational motion of the bodies, (ii) rotational motion
of the bodies with respect to the locally transported frames
with the bodies, (iii) orientation of the local (reference)
frames transported by the bodies. Up to now, most of the attention
(both theoretical and observational) has been focused on the
relativistic contributions to the (lunar) translational motion and 
local frame definitions \footnote{$\dag$}{We leave here aside 
relativistic definition of the time scales in the solar system; 
Brumberg \& Kopejkin 1990, Seidelmann \& Fukushima 1992.}. 
Recent analytical insights can be found e.g. in Soffel (1989) and 
Nordtvedt (1994). It is interesting to note, that in the case
of the Moon, the dominant relativistic range terms (with the amplitude
of a few centimeters) are due to the relativistic tidal field of the
Sun in the Earth reference frame (e.g. Nordtvedt 1991, 1994l; 
Shahid-Saless 1992).
The Earth causes principally an advance of lunar perigee
(a so far unmeasurable value of about $0.6$ {\it mas} per year).
To our knowledge, none of the previous works touched on 
relativistic terms related to the lunar rotation. 
Even more generally: only a small fraction of all tests of relativity
theory relates to the spin motion of the bodies (for example
see the evidence for the spin-orbital coupling in the binary pulsar
systems; Weisberg \etal 1989). From this point of view, and the
fact of the existence of high-quality LLR data, it is an 
interesting task to investigate relativistic contributions
to the lunar rotational motion.

In a recent brief communication (Bois \& Vokrouhlick\'y 1994), we
reported on the possible importance of a relativistic phenomenon connected
with the lunar rotational motion not used so far in the LLR data
reduction programs. In this paper, we wish to elaborate more on the
topic. 

\titlea{ Relativistic spin terms in the Earth-Moon system}

The most suitable formulation of the post-Newtonian (PN) theory of motion
of a system of N weakly self-gravitating extended bodies for purposes
of celestial mechanics has been recently presented in a series
of papers by Damour, Soffel and Xu (Damour \etal 1991, 1992, 1993, 1994;
hereafter quoted as DSX). We shall essentially
follow their approach (the same applies for notation).
Success of DSX is based on precise distinguishing
between global (barycentric) translational motion of the bodies
and their local rotational motion resulting in careful usage of
corresponding global and local coordinate systems. Here, we shall
deal only with the local rotational motion. As the dynamical state of the
local frames is intimately connected with the rotational motion of the
corresponding body,
we shall also report on the results concerning de Sitter precession of the
Earth (dynamically non-rotating) reference frame.

Damour \etal (1993) have succeeded in defining the local spin vector 
$S^\A_a$ characterizing rotation of body A (embodied in the isolated
system of N extended bodies). Employing first principles they
showed that it satisfies Eulerian dynamical equations
$$ {dS^\A_a \over dT_\A} = \epsilon_{abc} \left[\sum_{k \geq 0} {1
 \over k!} M^\A_{bK} G^\A_{cK} + {1 \over 2c^2}
 S_b^\A H^\A_c \right] \; , \eqno(1) $$
\noindent where $(M^\A_L,S^\A_L)$ are mass and spin
(Blanchet-Damour; BD) multipoles characterizing the PN gravitational field 
of the extended bodies while $(G^\A_L,H^\A_L)$ are tidal gravitoelectric and
gravitomagnetic PN fields, $c$ is the light velocity. All these quantities
are carefully defined and discussed in Damour \etal (1991, 1992, 1993). Let
us just note that we adopt definition (4.24) in Damour \etal (1993) for spin
multipole moments possessing the simpler form (1) for the dynamical
equations. Superindex A denotes a corresponding body and $T_\A$ stands
for the local time measured along the worldline ${\cal L}_\A$ of the
origin of the local chart associated with body A (to be identified
later with its mass center). Local coordinate time $T_\A$ is related to
the global coordinate time $t$ by the formula
$$ {d T_\A \over dt} = 1 - {1 \over c^2} \left({1 \over 2} \vv_\A^2 +
 {\bar w}^\A\right) + \O(4) \; , \eqno(2) $$
\noindent where $\vv_\A$ is the velocity of the local
frame origin expressed in the global system and ${\bar w}^\A$
is the external gravitational potential influencing the motion of 
body A [see Eq.~(4)].

One learns from the DSX series that the first principles lead to dynamical 
laws for the lowest multipoles only (zero and first-order mass
multipoles and first-order spin multipoles). If one wanted to
develop a complete theory of the motion one would need to close the
laws by some physical assumptions. As the rela\-ti\-vi\-ty-attri\-bu\-ted
terms in the dynamical equations contain small parameters
(squares of the orbital velocities divided by the light speed,
Schwarz\-schild radii divided by the orbital  distances), 
we shall restrict ourselves to the 
pole-dipole-quad\-ru\-pole (PDQ) truncated model (clarified
in Damour \etal 1991, 1992). This means that we shall simply neglect all
higher multipole moments in the corresponding relativistic terms.
Moreover, because we do not dispose of dynamical equations
for the quadrupole moments $M^\A_{ab}$ we shall adopt
the `rigid' model of the extended bodies.
Although the notion of rigidity faces conceptual
problems in the theory of relativity, we shall
formally keep the rigid scheme as known from the Newtonian
approach. Practically this is acceptable since the
relativistic quadrupole contributions are very small
(see the Appendix for more details).

\titleb{ De Sitter precession of the local dynamically non-rotating
 reference frames}
We are free to fix several data concerning local commoving systems of the 
extended bodies. Namely, we shall follow the DSX recommendation to identify
the origin of the local system with the mass center of the corresponding
body: $M^\E_a = 0$ (in the following, we shall use index E when specifying
the Earth for the previously arbitrary body A; similarly, index M will
be reserved for the Moon).
This choice advantageously discards the first term in
Eq.~(1) for $k=0$. Secondly, in this study we shall use dynamically 
non-rotating local systems which are characterized by the condition:
$H^\E_a = 0$ (Coriolis effacing). Let us recall, that local dynamically
non-rotating frames show a slow (de Sitter) rotation with respect to the
kinematically non-rotating frames. The latter ones are practically
realized by the VLBI measurements, so that the
de Sitter rotation of the body A local chart is measured
with respect to the VLBI reference system (see e.g. Soffel \&
Brumberg 1991; Jacobs \etal 1993).

Damour \etal (1993) showed that the total first order
gravitomagnetic tidal moment reads
$$\eqalignno{ H_a\left(T_\E\right) =&\; -4\epsilon_{ijk} R^i_{{\rm E} a}
 \left[ \partial_j {\bar w}^\E_k + v_\E^j \partial_j {\bar w}^\E
 \right]_{\cal L_\E} \cr & +\epsilon_{abc}\left[V^\E_b A^\E_c
 + c^2 {dR^i_{{\rm E} b} \over dT_\E} R^i_{{\rm E} c}\right]
 _{\cal L_\E} \; , & (3)\cr} $$

\noindent where $R^i_{{\rm E} a}$ denotes the orthogonal rotational matrix
involved in the transformation between local and global systems.
Notice the clear theoretical interpretation of the first two terms 
in (3) as the (PN) gravitational influence of the external bodies, 
while the last two terms have coordinate origin in the transformation
between global and local systems.

In the numerical implementation of this theory (Sec. 3), we shall adopt
only the monopole part of the external metric, so that
$$ {\bar w}^\E = \sum_{\B \neq \E} {GM_\B \over r_\EB} \quad ,
 \qquad  {\bar w}^\E_i = \sum_{\B \neq \E} {GM_\B \over r_\EB}
 v_\B^i \; , \eqno(4)$$

\noindent neglecting `higher order terms'. It can be easily verified
that the contribution of these higher order terms is about 4 orders of 
magnitude smaller than the one of terms arising from the monopole part (4).

Introducing components of the rotational vector $\omega^\E_i$ characterizing
rotation of the dynamically non-rotating frame with respect to the
kinematically non-rotating (barycentric) one in the form
$$ \omega_i^\E = -{1 \over 2} \epsilon_{ijk} {dR^j_{{\rm E} a}
 \over dT_\E} R^k_{{\rm E} a} \; , \eqno(5) $$

\noindent ($\epsilon_{ijk}$ is the Levi-Civita fully antisymmetric tensor)
we arrive at
$$ {\vec \omega}^\E = -{1 \over 2c^2} {\bf \vec v}_\E\times{\bf \vec a}_\E
 + {2 \over c^2} {\bf \vec v}_\E\times{\nabla {\bar w}^\E} + {2 \over c^2}
 \nabla \times {\bf \vec {\bar w}}^\E\; , \eqno(6) $$

\noindent which is to be compared with results of Kopejkin (1988),
Will (1993). Note that some authors (e.g. Will 1993) prefer to include
the `gravitational part' of the barycentric acceleration ${\bf \vec a}_\E$
in the second term, arriving then at a factor $3/2$ rather than $2$.

It should be also noted that Damour \etal (1993) proved that local
frames induced by a dynamically non-rotating local chart are Fermi-Walker
transported along the worldline ${\cal L}_\A$ of the body A mass center
with respect to the external metric. This property shows that DSX
also succeeded in constructing Fermi (or generalized Fermi) frames 
around massive bodies (e.g. Ashby \& Bertotti 1986; Fukushima 1988).

\titleb{ Relativistic terms in the quasi-Newtonian torque}
Damour \etal (1993) succeeded in deriving an explicit `quasi-Newtonian' 
form for the first relativistic contribution in the rotational motion
arising from the first term in Eq.~(1) -- the only term remaining
in the PDQ model. Laborious evaluation of the tidal
gravitoelectric multipole $G^\M_{ab}$ leads to (Damour \etal 1993, 1994)
$$\eqalignno{G^\M_{ab} =&\; R^i_{{\rm M} a} R^j_{{\rm M} b} \sum_{\B \neq
 \M} {3GM_\B \over r_\MB^3} {\rm STF}_{ij} \biggl\{ n_\MB^{ij} +
 {1 \over c^2}\Bigl[n_\MB^{ij}\bigl(2\vv_\MB^2  \cr
 &\; \hbox{~~~~~~~}-2 {\bar w}^\M -{\bar w}^\B -{5 \over 2}(\n_\MB .\vv_\B)^2 -
 {1 \over 2} {\bf a}_\B.{\bf r}_\MB\bigr) \cr
 &\; \hbox{~~~~~~~}  + a_\MB^i r_\MB^j + v_\MB^i
  v_\MB^j -2(\n_\MB.\vv_\MB) n_\MB^i v_\MB^j \cr
 &\; \hbox{~~~~~~~} -  (\n_\MB.\vv_\M) n_\MB^i\bigl(v_\M^j-
 2v_\B^j \bigr)\Bigr]\biggr\}  +  \O(4) \; , & (7) \cr}$$
$$\eqalignno{ n_\AB^i =\; & (r_\M^i - r_\B^i)/r_\MB\; , \quad 
 r_\MB^i = r_\M^i - r_\B^i \; ,\cr
 v_\MB^i =\;& v_\M^i - v_\B^i\;  , \quad a_\MB^i = a_\M^i - a_\B^i\; ,
 \cr}$$

\noindent where the symbol ${\rm STF}_{ij}$ denotes the symmetric and 
trace-free part of the corresponding tensor quantity (see e.g. Damour \etal 
1991). In the first term one might recognize `the Newtonian torque'. 
However, one subtle point should be reminded here. Mass (BD) multipoles 
$M^\A_L$ are autonomous quantities defined in the PN level and have nothing
to do with analogous Newtonian quantities\footnote{$\ddag$}{ From this
point of view all harmonic coefficients corresponding to the
multipoles of the Earth-Moon system should be theoretically redefined
when considering higher relativistic terms. In practice, due to the 
smallness of the studied effects, we have accepted values of the
Earth-Moon gravity field coefficients as given by previous studies.}
(e.g. zeroth multipole $M_\A$
representing the mass of the body is trivially constant in the Newtonian
approach while variable on the PN level as soon as mass quadrupoles
are included; however, in what follows we shall neglect such variations,
as they can be shown to be extremely small for the solar system
bodies). It should be understood that, conceptually, relativity does not
provide {\it corrections} to the Newtonian solution. Rather one
can {\it identify} some terms, after developing a consistent PN approach, 
as {\it Newtonian-like}. This is exactly the case of the similarity
of the first term in Eq.~(7) with Newtonian torque. Clear appearance
of this point is also one important achievement of the DSX theory
(more discussion about this point might be found in Damour \etal 1993).

Formally, the remaining terms in the right hand sides of Eq.~(7) represent
the `relativistic' contribution. The aim of this paper is to investigate 
numerically their possible importance.

\titlea{ Numerical results}
We have constructed a model of complete lunar motion
including relativistic effects discussed previously (more details
about our model will be published elsewhere).
Our goal is to include all phenomena up to the precision level resulting 
from the LLR technology, and if possible even better for reasons of 
consistency
(e.g. at least 1 cm for the distance, 1 {\it mas} for the librations).
In particular, several phenomena capable of producing effects of at least
0.1 {\it mas} in the lunar physical librations have been included and 
analysed. The effects resulting from planetary actions (essentially Venus),
the Earth-Moon figure-figure interactions, and non-rigidities of the Moon 
and the Earth have been presented in previous papers (Bois \etal 1992;
Bois \& Journet 1993).
In the subsequent sections, we present the results of the relativistic
contributions related to the rotational motion of the Moon, as well as
the relativistic precession of the Earth reference frame.

The relativistic part of the model has been introduced as a 
pole-dipole-qua\-dru\-pole in the previous discussion. 
However, such truncation would not be sufficient
to meet the prescribed precision as mentioned earlier. One has to
include higher order mass multipoles (up to $l=5$ in the Moon case
while $l=4$ in the Earth case) of the bodies. Corresponding
`Newtonian-like' terms do appear in the systematic construction of
the PN theory of bodies endowed with higher order multipoles as can
be verified using a general scheme in the Appendix of Damour
\etal (1992). However, we neglected any terms of purely relativistic
origin associated with multipoles higher than second order.
Similarly, we took into account several terms arising from
the tidal deformations of the bodies (both elastic and anelastic).
Again, one can adopt a very broad model (including
such tidal effects on the relativistic base; see Damour \etal 1991),
but subsequently retain only necessary `Newtonian-like' terms
corresponding to the required precision.

Summarizing previous ideas, we can state that our model is not
Newtonian but rather `Newtonian-like', resulting from truncation
of the fully post-Newtonian (DSX) framework.
We shall focus on reporting the results related to the relativistic
effects mentioned in section 2.
\begtabfullwid
\tabcap{1}{ De Sitter precession -- results and comparison with previous
 works. The first row relates to the principal secular term, the other
 two to the periodic terms of year and half-year variations.}
\halign{#\hfil&&\quad#\hfil\cr
\noalign{\hrule\medskip}
 term & Fukushima & Bizouard \etal & our results & units \cr
\noalign{\medskip\hrule\medskip}
secular & 1.919 & 1.914 & 1.9193 & {\it as/cy} \cr
periodic {\it (1 year)} & 0.153 & 0.1525 & 0.1530 & {\it mas} \cr
periodic {\it ($1\over 2$ year)} & $2\times10^{-3}$ & $1.9\times10^{-3}$
                                     & $1.89\times10^{-3}$ & {\it mas} \cr
\noalign{\medskip\hrule}}
\endtab

\titleb{ De Sitter precession of the Earth reference frame}
We shall parametrize the vector $\omega^\E_i$ [see Eq.~(5)] corresponding to
the de Sitter rotation of the Earth reference frame by three Eulerian
angles $(\psi,\theta,\phi)$ referred to the ecliptic barycentric frame
($\psi$ being the precession angle as in usual astronomical practice).
As the phenomenon demonstrates precession, we
shall focus our attention on the $\psi$ angle. Investigating its
behaviour by a Fourier filter, we arrived at the results summarized
in table 1 [for comparison we include also results of
Fukushima (1991) and Bizouard \etal (1992)].

Comparison of our results with those obtained by Fukushima (1991) and 
Bizouard \etal (1992) shows a very good agreement, or at least a very close
one. 
It can be seen, that the half-year term amplitude would in reality be 
affected by higher than monopole contributions in the external potentials 
$({\bar w}, {\bar w}_i)$ in Eq.~(6). However, its practical importance is 
small because such a term is below the precision threshold of current 
astrometric methods (e.g. Kovalevsky 1990; Jacobs \etal 1993).

It should be noticed that the lunar-reference frame undergoes
a similar de Sitter precession to the Earth one. An alternative way of 
representing the two effects is to introduce the de Sitter precession 
of the common Earth-Moon barycentric reference frame (e.g. Dickey \etal 1989),
as it can be easily verified that the principal effects in table~1
originate in the solar action. However, due to the mutual
Earth-Moon action, the de Sitter precession of the two reference
frames differs slightly. Detailed inspection shows that the lunar
reference frame undergoes an additional precession of the order of
$30\, mas/cy$, so far not a measurable value. Moreover, bringing out 
this secular effect is a complicated task due to the other
(phenomenologically modeled) secular precessional effects related
to the selenophysical phenomena.

Finally, it is to be pointed that we discuss de Sitter (geodetic) precession
in `the narrow sense' following the terminology of Brumberg \etal
[1991]. We thus involve in our calculation only first two terms
in the fundamental formula (6). Discrepancy between our value and 
$1.9198\, as/cy$ given by Brumberg \etal [1991] is due to the fact
that the latter value includes also contribution from the Lense-Thirring
term due to coupling of the solar gravity field with the Earth-Moon
angular momentum. 

\titleb{ Relativistic terms in the lunar rotational motion}
In what follows, we shall report on the influence of `the
relativistic' terms in (7) on lunar rotational motion. It should
be pointed out, that as for as the authors know such terms have not yet
been applied in lunar theory. The main reason is that the theoretical description
of the relativistic contributions in the rotational motion of bodies
was unsatisfactory (see Damour \etal 1991).

Figure 1 shows the behaviour of resulting librations due to the
relativistic contributions in Eq.~(7). The three Eulerian angles 
$(\psi,\theta,\phi)$ of the classical 3-1-3 angular sequence represent 
the physical librations of the Moon's rotational motion with respect to a
reference system given by a terrestrial equatorial frame (J2000).
As the data drawn in figure 1 have been obtained
as a difference of the two numerical integrations of the lunar
rotation, first disregarding the relativistic terms, secondly
accounting for them, we call them $(\Delta\psi,\Delta\theta,
\Delta\phi)$.
The global frame time $t$ on the $x$ axis is the solar system barycentric 
time TCB. This choice is clearly governed by the fact that
the solar system barycentric time belongs to the set of the two most
important coordinate time scales defined in the solar system, contrary
to the lunar (local) coordinate time (e.g. Soffel \& Brumberg 1991;
Seidelmann \& Fukushima 1992). We have thus combined equations
(1) and (2) together with (7) in order to obtain the final form
of the Euler dynamical equations for lunar librations.
Let us also note, that the initial date of the integration
presented in Fig.~1 coincides with that of the JPL DE303 ephemeris, as
we have adjusted our model, with a few centimeter residuals, to it on 
the interval spanning the first $24$ years (interval of the DE303
validity).
\begfigwid 12.8 cm
\figure{1a-c)}{ Relativistic lunar librations as they formally appear
 in the series of the Euler angles characterizing lunar rotation --
 differences $(\Delta\psi,\Delta\theta,\Delta\phi)$ with respect to
 the solution free of the relativistic torque (see also Fig.~2).
 Milliarseconds are on the vertical axis and years on the horizontal axis.
 Dominant terms of $18.6$ and $80.1$ year periodicities are in all
 Euler angles. The angle of lunar proper rotation contains also a
 component related to the ordinary resonant frequency of $2.9$ years 
 for physical librations in longitude.}
\endfig

Let us recall the major changes in complete libration angles
during the considered time span: the lunar nutation angle
$\theta$ changes between $22$ and $25$ degrees with the major periodicity
of $18.6$ years being driven by the nodal period of the lunar orbit; the 
proper rotation of the Moon around its axis of figure is effected by the 
complete angle $\phi$ with the mean period of $27.3$ days; the precession 
angle $\psi$ is also locked in resonance with the lunar node period of
$18.6$ years.
Further forced behaviour due to the indirect action of the planets and 
higher harmonics of the Moon can be found in e.g. Eckhardt (1982) or 
Moons (1984).

We remark long-term components with periods of $18.6$ years (nodal
rate of the lunar orbit), $80.1$ years and $2.9$ years (in $\phi$) 
which are the two resonant frequencies of the spin-orbit motion of the 
Moon (Kopal 1969; Yoder 1981).
Appearance of these terms is not surprising in view of the form of the 
relativistic torque (7). The magnitude of the relativistic
librations is of the order of one milliarcsecond in all Euler
angles. The ratio of this magnitude with respect to the corresponding 
one coming from analysis of the Newtonian-like pole-quadrupole torque
[i.e. the first term in (7)]
is about $5\times 10^{-8}$. Such a value corresponds to the order of the
small (PN-) relativistic parameters $(v_A/c)^2$ and ${\bar w}^B/c^2$ both
of them being of the order of $10^{-8}$. The latter factor also clearly
determines the amplitude of the ratio of the relativistic and Newtonian
torques associated with pole-quadrupole interaction [see Eq.~(7)].

The non-linear features of the differential equations, the correlation
degree of the studied effect with respect to its neighbours (in Fourier
space) and the spin-orbit resonance, in the lunar case, make it hardly 
possible to speak about `pure' effects with their proper behaviour 
(even after fitting of the initial conditions). 
The differentiation method may give the right qualitative behaviour of
an effect and a good quantification of this effect relative to its
neighbours (see Table 2). The effects are not absolutely de-correlated
but relatively isolated. It is also possible to make some comparisons
with analogous problems (here for instance, with the analogous torques
acting on the Earth's rotation). When a rotational effect is simply periodic,
a fit of the initial conditions for a set of given parameters only refines
without changing completely the effect's behaviour. Nevertheless, the 
particular status of the relativistic effects may lead us to a precise 
quantification and will be discussed in a forthcoming paper
(see also discussion in the next section).

In order to get an idea about the importance of the relativistic
contribution, Table 2 presents a comparison of various effects
acting on the Moon's spin motion, obtained by Bois \etal (1992), Bois \&
Journet (1993), and in the present work. The amplitudes correspond to the
maxima reached by one of the three Eulerian angles $(\Delta\psi,\Delta
\theta,\Delta\phi)$ in each case.
Results are given respectively for two periods of librations when they exist,
namely 2.9 years and a period running around 18.6 to 24 years.
The acronym $ExMy$ signifies orders of multipoles respectively in
interaction, namely the order $x$ of the Earth acting on the order $y$
of the Moon. $V$ means Venus; $LTL$ means lunar tidal librations;
$RCL$ refers to the relativistic contributions of (7) represented by
the librations given in figure 1. We note that the amplitude of the
relativistic librations is comparable with the quasi-Newtonian
interaction of the Earth quadrupole and the Moon sextupole or the
influence of the Venus pole.
\begtabfull
\tabcap{2}{ Comparison of various effects on the Moon's librations. See
 the text for terminology.}
\halign{#\hfil&&\quad#\hfil\cr
\noalign{\hrule\medskip}
          effect & 2.9 years & 18.6/24 years & units \cr
\noalign{\medskip\hrule\medskip}
V0M2   &           & $\le 3$    & {\it mas} \cr
E0M5   & 10        & $\le 20$   & {\it mas} \cr
E2M2   &           & $\le 80$   & {\it mas} \cr
E2M3   & 0.5       & $\le 1$    & {\it mas} \cr
LTL    & 1         & $\le 4$    & {\it mas} \cr
RCL    & 0.1       & $\le 2$  & {\it mas} \cr
\noalign{\medskip\hrule}}
\endtab

\titlea{ Discussion and conclusion}
The main results of this communication can be summarized as
follows:
\smallskip

\item{i)} numerical simulation of the de Sitter precession of the 
 Earth reference frame yields results very close to those of Fukushima 
 (1991) or Bizouard \etal (1992);
\item{ii)} preliminary tests of importance of the relativistic torque 
 in the case of the Moon's rotational motion suggest that this effect ought 
 to be included in the lunar models.
\smallskip

\noindent Let us point out that the first result (i) is to be understood
more as an inner test of our model (program), because this phenomenon has 
already been studied both theoretically (e.g. Fukushima 1991; Bizouard \etal
1992) and even when interpreting the LLR data (e.g. Bertotti \etal 1987;
Shapiro \etal 1988; Dickey \etal 1989; M\" uller \etal 1991). On the contrary,
the second result (ii) does not seem to have been studied yet.
Bizouard \etal (1992) recently investigated the role of analogous terms in 
the case of the Earth's rotation. They concluded that the magnitude of this 
effect was negligible. Interestingly, in the lunar case the resulting 
contribution does not seem to be totally negligible.
Similar results were reported for the influence of the higher order 
multipole coupling on the `Newtonian level' (dominated
by the $J_2-J_2$ interaction). Souchay (1994; private communication) has
studied such effects for the Earth's rotation, while Bois \etal (1992) 
considered the lunar case. The resulting behaviour is mainly explained
by the ratio of the resulting torque intensity to the angular velocity
via the degree of asymmetry of the body. This ratio is larger for the Moon
than for the Earth. 

It should be emphasized that results presented in Sec.~3.2
(notably, direct differences of the Euler angles) place
only {\it the upper limit} for the magnitude of the
discussed relativistic lunar librations with respect to their
observability. Three reasons account for this: (i) 
a subclass of infinitesimal rotations is not suitably represented by
the 3-1-3 Eulerian parametrization as discussed in the next paragraph;
(ii) initial conditions of the integration have not been 
adjusted; (iii) interpretation of the relativistic effects
requires special care due to possible coordinate
dependence of the effects. We comment on these problems in
the following paragraphs.
\begfigwid 7 cm
\figure{2}{ The concepts of representation of the difference of
 the two rotated frames (1 and 2) with respect to the reference one
 (R): (a) simple difference of the Euler angles relating the two
 moving frames with respect to the reference one, (b) $(\delta,
 {\bf n})$ parametrization of the infinitesimal rotation of the
 frame 2 with respect to the frame 1. In our application, frame
 1 is represented by the lunar rotation not accounting for the
 effect of the relativistic lunar librations, while the
 rotation of the frame 2 contains contribution due to the
 relativistic librations. See the text for more details.}
\endfig

Differences of the Euler angles of the two configurations reported
in Sec.~3.2 (Fig.~1) show how the relativistic librations
formally appear in the lunar rotation. 
However, we must face the question of their possible observability through 
the LLR facility, the problem being qualitatively more difficult. 
The first step to reach the observability task would require 
reparametrization of the description
of the relativistic contribution to the lunar librations.
Instead of the difference of the two series of Euler
angles [accounting and not-accounting for the influence of the
relativistic torque (1)], we shall actually construct the rotation
of the reference frame defined by the lunar rotation solution,
without the relativistic term to the reference frame given
by the lunar rotation including the relativistic contribution.
It appears that a suitable representation of the rotation we are looking
for is obtained by specifying the instantaneous (virtual) rotation
axis and corresponding angle $\delta$ of rotation around the
axis (see figure~2 for this concept).
The position of the virtual axis is actually not important
for our order of magnitude argument, as we are basically seeking 
the amplitude of the additional lunar frame rotation imposed by
the appearance of the relativistic term in the right hand
sides of the Euler equations (1). Using elements of the
rotation group algebra (for information on the mathematical
techniques of rigid body rotation see e.g. Mignard 1989),
one arrives at the following formulas
$$\eqalignno{\cos\delta =&\, 2q_0^2 - 1 \;, & (8) \cr
 q_0 = &\, \cos{\theta_1 \over 2} \cos{\theta_2 \over 2}
  \cos\left(\psi_1-\psi_2+\phi_1-\phi_2\right) \cr
  &\, + \sin{\theta_1 \over 2} \sin{\theta_2 \over 2}
  \cos\left(\psi_1-\psi_2-\phi_1+\phi_2\right) \;,\cr}$$

\noindent where $(\psi_1,\theta_1,\phi_1)$ and $(\psi_2,\theta_2,\phi_2)$
are two sets of Eulerian angles related to the two reference frames 
under consideration.
In the case of small differences in precession angles $\Delta\psi=
\psi_2-\psi_1$ and proper rotation angles $\Delta\phi = \phi_2-\phi_1$
(which is the case of curves in Fig.~1), one can simplify the previous
formula to the form
$$\delta^2 \approx \Delta\theta^2 + \left(\Delta\psi+\Delta\phi
  \right)^2 + 2\left(\cos{\bar \theta} -1\right)
  \Delta\psi\Delta\phi \; , \eqno(9)$$

\noindent where $\Delta\theta = \theta_2-\theta_1$ and ${\bar \theta}$ is
the averaged nutation angle of the two discussed configurations. 
Components of the unit vector ${\bf n} = (n_1,n_2,n_3)$ of the rotation 
axis direction then read
$$\eqalignno{ n_1 = &\, \left(\Delta\theta \cos{\bar \phi} +
 \Delta\psi \sin{\bar \theta} \sin{\bar \phi}\right)/\delta \; ,&(10a)\cr
 n_2 = &\, \left(\Delta\theta \sin{\bar \phi} -
  \Delta\psi \sin{\bar \theta} \cos{\bar \phi}\right)/\delta \; ,&(10b)\cr
 n_3 = &\, \left(\Delta\psi \cos{\bar \theta} + \Delta\phi\right)/\delta
  \; , & (10c) \cr}$$

\noindent where similar quantities have been introduced for
the remaining Euler angles: $\Delta\psi = \psi_2-\psi_1$, 
$\Delta\phi = \phi_2-\phi_1$ and $({\bar \psi},{\bar \phi})$ for
the averaged quantities. One easily recognizes the vector character 
of composition of the infinitesimal rotations. We
have passed data shown in Fig.~1 through this test computing
$\delta$-angle and corresponding rotation axis motion.

Figure~3 shows the time dependence of the $\delta$-rotational
angle. The resulting amplitude decreased
to about $0.8$ {\it mas}, as can be expected from the series of the
$\Delta\psi$ and $\Delta\phi$ angles. Notice that there essentially
appears only the sum of the two quantities in Eq.~(9) corresponding to 
the fact that the infinitesimal rotation of the reference frame can be
expressed by arbitrary precession and proper rotation, provided that 
these two elementary rotations are nearly opposite 
(this is a well understood inconvenience of the 3-1-3 Euler sequence, 
Mignard 1989).
One remarks a pure $18.6$ year periodicity of this quantity
forced by the nodal motion of the Moon. It is to be reminded here
that the residual value of about one milliarcsecond is on the
level of the observability through LLR if the best quality nights
are selected (moreover a nearly half order improvement is
planned for the forthcoming years; Ch. Veillet, private communication).
The motion of the axis ${\bf n}$, showing precession in the lunar 
equatorial plane, can be decomposed into terms with periods of
$2.9$, $18.6$ and $80.1$ years. The third component of ${\bf n}$ shows 
oscillations with a period of $2.9$ years.

Regarding the well defined frequencies and the stable character of the
signal due to relativistic librations, we hope that our results
represent well the true behaviour of the phenomena.
However, we plan to investigate in some detail the effect of the
initial condition adjustment on the results referred to here.
In any case we can conclude, that the
relativistic lunar librations introduced in this study are
among a class of phenomena influencing lunar rotation, which
`are standing just behind the door of observability in the LLR
data' and thus would require greater care in the near future.
\begfig 7 cm
\figure{3}{ Rotational angle $\delta$ (in {\it mas}) vs. time in
 years. The signal essentially decomposes into the $18.6$ years
 long-term component and the short-term $27.3$ day component.}
\endfig

Concerning to the third warning mentioned above [item (iii)], we
recall that the final proof of the observability of the relativistic 
librations is based upon their influence on the two-way laser beam timing,
expressed in the proper time of the observer (this being a coordinate 
independent quantity). Moreover, suitable de-correlation with the other 
stronger disturbing effects has to be met. Full analysis of such
problems lies out of the scope of this paper. 

We hope to devote a forthcoming study to the two topics mentioned
in the two preceding paragraphs [also items (ii) and (iii) mentioned 
previously].

Let us finally note that another interesting topic to be tackled in the
future concerns the generalization of the DSX formula (7) to the 
parametrized framework of the PN theory of gravitation
standardly used in the LLR data interpretation programs
(Dickey \etal 1989; M\" uller \etal 1991)\footnote{$\P$}{Methodological
importance of embodying of the Einsteinian general
relativity into a wider field of parametrized PN theories for purposes
of testing its validity is stressed e.g. by Damour (1992).}.

\acknow{We are grateful to J. Bryant for careful reviewing of the first
 version of this paper and the unknown referee for valuable 
 comments. D.V. has finished this paper while staying at the
 Observatoire de la C\^ote d'Azur (Dept. CERGA, Grasse) thanks to the
 H. Poincar\'e fellowship.}

\appendix { } In this appendix, we shall comment on our rigid-body 
Ansatz assumed for the relativistic part of the lunar rotation.
As noted in Sec.~2, conceptually we do not dispose of the equations 
for the quadrupoles $M^\A_{ab}$ unless a physical concept of the bodies
is settled. The most restrictive is rigidity, leaving
no free parameter in the description of the body mass multipoles.
Necessarily a set of three parameters of the rotation group
should be introduced (it is thus rather a geometric concept, 
being extreme in its physical content). However, in the sense of
a formal decomposition, the rigidity notion expresses successfully 
the behaviour of the mass multipoles of the solar system bodies.
It is then conventional to describe more complicated physical body-phenomena
(e.g. tides, rotational deformation etc.) as additional terms to
those due to the rigid concept (Munk \& MacDonald 1975).

As the relativistic contribution is small due to the presence of the 
first-order relativistic parameters (squares of the orbital velocities 
divided by the light speed, Schwarz\-schild radii divided by the orbital 
distances), it is advisable to couple it with the rigid concept
quadrupoles [see the first term in Eq.~(1)].
In this case, the mass quadrupoles $M^\A_{ab}$ satisfy equations
(see e.g. Dixon 1979)
$$ {d M^\A_{ab} \over d T_\A} = - {\rm STF}_{ab}\left\{ \Omega^\star_{{\rm A} ca}
 M^\A_{bc}\right\} \; , \eqno(A1)$$

\noindent where $\Omega^\star_{{\rm A} ab} \equiv \epsilon_{abc}
\Omega^\A_c$ is a dual antisymmetric tensor to the vector of the
angular velocity of the body A. Though theoretically elegant, this
approach is not suitable for astronomical calculations. Luckily
the solution of equation (A1) is known explicitly and allows us
to express quadrupole components directly in terms
of the parameters of the proper rotation of the bodies. Moreover,
because it is conventional in astronomical and geodynamical research
to use spherical harmonics analysis of the gravitational fields
with the corresponding notion of harmonic coefficients $(C^\A_{lm},
S^\A_{lm})$, we seek expressions of the $M^\A_{ab}$ in those terms.
If the system co-rotating with a particular body is 
defined by the principal axis of the moment of inertia, then all
quadrupole terms vanish with the exception of two: 
$J^\A_2 \equiv -C^\A_{20}$ and $C^\A_{22}$. 
Next, assuming that the rotation of body A is
parametrized by the set of Euler angles $(\psi_\A,\theta_\A,
\phi_\A)$, one gets after expansion the following
expressions of the mass quadrupoles [to be inserted in Eq.~(1)]
$$ M^\A_{ab}(T_\A) = J_2^\A M_\A a_A^2 {\bf A}_{ab}(\theta_\A,\psi_\A)
 \; , \eqno(A2)$$

\noindent where
$$\eqalignno{
 {\bf A}_{11} = & {1 \over \sqrt{2}}\sin^2\theta_\A\cos 2\psi_\A +
  {1 \over 3}  \left(1 - {3 \over 2}\sin^2 \theta_\A\right) \; , \cr
 {\bf A}_{22} = & -{1 \over \sqrt{2}}\sin^2\theta_\A\cos 2\psi_\A +
  {1 \over 3} \left(1 - {3 \over 2}\sin^2 \theta_\A\right) \; , \cr
 {\bf A}_{33} = & - {2 \over 3} \left(1 - {3 \over 2}\sin^2 \theta_\A
  \right) \; , \cr
 {\bf A}_{12} = & {1 \over \sqrt{2}} \sin^2\theta_\A \sin 2\psi_\A \; ,
   \cr
 {\bf A}_{13} = & -{1 \over 2} \sin 2\theta_\A \sin\psi_\A \; , \cr
 {\bf A}_{23} = & {1 \over 2} \sin 2\theta_\A \sin\psi_\A \; , \cr}$$

\noindent for those related to the principal quadrupole terms
$J_2^\A$, and
$$ M^\A_{ab}(T_A) = C_{22}^\A M_\A a_\A^2 {\bf B}_{ab}(\psi_\A,\theta_\A,
 \phi_\A) \; , \eqno(A3)$$

\noindent where
$$\eqalignno{
 {\bf B}_{11} = & \left(1+\cos^2\theta_\A\right)\cos\psi_\A\cos 2\phi_\A\cr
  &\; \hbox{~~~}+ 2\cos\theta_\A\sin\psi_\A\sin 2\phi_\A + 6\sin^2\theta_\A
   \cos 2\phi_\A \; ,  \cr
 {\bf B}_{22} = & -\left(1+\cos^2\theta_\A\right)\cos\psi_\A\cos 2\phi_\A\cr
  &\; \hbox{~~~}- 2\cos\theta_\A\sin\psi_\A\sin 2\phi_\A + 6\sin^2\theta_\A`
   \cos 2\phi_\A \; ,  \cr
 {\bf B}_{33} = & -12 \sin^2\theta_\A \cos 2\phi_\A \; , \cr
 {\bf B}_{12} = & \left(1+\cos^2\theta_\A\right)\sin\psi_\A\cos 2\phi_\A\cr
  & \; \hbox{~~~} + 2\cos\theta_\A\cos\psi_\A\sin 2\phi_\A \; , \cr
 {\bf B}_{13} = &\, 2 \sin\theta_\A \left(\sin\psi_\A\sin 2\phi_\A +
  \sin\psi_\A\cos\theta_\A\cos 2\phi_\A\right) \; , \cr
 {\bf B}_{23} = &\, 2 \sin\theta_\A \left(\sin\psi_\A\sin 2\phi_\A -
  \sin\psi_\A\cos\theta_\A\cos 2\phi_\A\right) \; , \cr}$$

\noindent for the part associated with the $C^\A_{22}$ coefficient.
Parameter $a_\A$ denotes a characteristic dimension of body
A related to the definition of the quadrupole harmonic coefficients
(it finally drops out from the definition of
quadrupole tensor $M^\A_{ab}$). Notice that the $J_2$-associated
part does not depend on the angle $\phi_\A$ of the proper rotation
due to axial symmetry.

Finally, let us recall the definition of the symmetric and trace-free
part of a second-rank tensor $t_{ij}$ [Eq.~(1)]
$$ {\rm STF}_{ij}\left[t_{ij}\right] = {1 \over 2}\left(t_{ij}+t_{ji}
 \right) - {1 \over 3} \delta_{ij} {\rm Tr}\left({\bf t}\right) \; ,
 \eqno(A4)$$

\noindent where $\delta_{ij}$ is the Kronecker symbol and ${\rm Tr}
\left({\bf t}\right)$ stands for the tensor trace. In particular, we have
$${\rm STF}_{ij}\left[n_\AB^{ij}\right] = {\rm STF}_{ij} \left[
 n_\AB^i n_\AB^j\right] = n_\AB^{ij} - {1 \over 3}\delta_{ij}\; .$$

\begref{References}

\ref Ashby N., Bertotti B., 1986, Phys. Rev. D34, 2246

\ref Bertotti B., Ciufolini I., Bender P. L., 1987, Phys. Rev. Lett. 58,
 1062

\ref Bizouard C., Schastok J., Soffel M., Souchay J., 1992, In: N. Capitaine
 (ed.), Journees 1992: Systemes de reference spatio-temporels,
 Observatoire de Paris

\ref Bois E., Journet A., 1993, Celest. Mech. 57, 295

\ref Bois E., Vokrouhlick\'y, D., 1994, In: A. Yu (ed.), Relativistic
 Gravitational Experiments in Space \& Related Theoretical Topics, 
 $2^{\rm nd}$ W. Fairbanks conference, Hong-Kong, in press

\ref Bois E., Wytrzyszczak I., Journet, A., 1992, Celest.
 Mech. 53, 185

\ref Brumberg V. A., Kopejkin S. M., 1990, Celest. Mech. 48, 23

\ref Brumberg V. A., Bretagnon P., Francou G., 1991, In: N. Capitaine,
 Journees 1991 Systemes de reference spatio-temporels, Observatoire
 de Paris, Paris

\ref Damour T., 1992, In: B. Julia, J. Zinn-Justin, Gravitation and
 Quatizations, North-Holland, Amsterdam

\ref Damour T., Taylor J. H., 1992, Phys. Rev. D45, 1840

\ref Damour T., Esposito-Far\`ese G., 1992, Class. Quantum Grav. 9, 2093

\ref Damour T., Soffel M., Xu Ch., 1991, Phys. Rev. D43, 3273

\ref Damour T., Soffel M., Xu Ch., 1992, Phys. Rev. D45, 1017

\ref Damour T., Soffel M., Xu Ch., 1993, Phys. Rev. D47, 3124

\ref Damour T., Soffel M., Xu Ch., 1994, Phys. Rev. D49, 618

\ref Dickey J. O., Newhall X. X., Williams J. G., 1989, Adv. Space Res.
 9, 975

\ref Dickey J. O., Bender P., Faller J., \etal, 1994, Science 265, 482

\ref Dixon W. G., 1979, In: J. Ehlers (ed), Isolating graviting systems
 in general relativity, North-Holland, Amsterdam

\ref Eckhardt, D.H., 1982, In: O. Calame (ed.), High-Precision Earth
 Rotation and Earth-Moon Dynamics, 193

\ref Fukushima T., 1988, Celest. Mech. 44, 61

\ref Fukushima T., 1991, A\&A 244, L11

\ref Jacobs C. S., Sovers O. J., Williams J. G., Standish E. M.,
 1993, Adv. Sp. Res. 13, (11)161

\ref Kopal Z. 1969, The Moon, D. Reidel Publ. Co., Dordrecht

\ref Kopejkin S. M., 1988, Cel. Mech. 44, 87

\ref Kovalevsky J., 1990, Astrom\'etrie Moderne, Springer-Verlag,
 Berlin

\ref Mignard F., 1989, In: D. Benest, Cl. Froeschl\'e (eds.), Modern
 methods in celestial mechanics, Editions Fronti\`eres, Gif-sur-Yvette

\ref Moons, M., 1984, Cel. Mech. 34, 263

\ref Mulholland J. D., 1980, Rev. Geoph. and Sp. Sci. 18, 549 

\ref Munk W. H., MacDonald G. J. F., 1975, The rotation of the
 Earth, Cambridge Univ Press, Cambridge

\ref M\" uller J., Schneider M., Soffel M., Ruder H., 1991, ApJLett
 382, L101

\ref Nordtvedt K., 1968, Phys. Rev. 170, 1186

\ref Nordtvedt K., 1973, Phys. Rev. D7, 2347

\ref Nordtvedt K., 1988, Phys. Rev. D37, 1070

\ref Nordtvedt K., 1991, Phys. Rev. D43, 3131

\ref Nordtvedt K., 1994, Icarus, submitted

\ref Seidelmann P. K.,Fukushima T., 1992, A\&A 265, 833

\ref Shahid-Saless B., 1992, Phys. Rev. D46, 5404

\ref Shapiro I. I., Reasenberg R. D., Chandler J. F., Babcock
 R. W., 1988, Phys. Rev. Lett. 61, 2643

\ref de Sitter W., 1916, MNRAS 77 (Secs. 24-26), 155

\ref Soffel M. H., 1989, Relativity in astrometry, celestial mechanics
 and geodesy, Springer-Verlag, Berlin

\ref Soffel M. H., Brumberg V. A., 1991, Celest. Mech. 52, 355

\ref Weisberg J. M., Romani R. W., Taylor J. H., 1989, ApJ 347, 1030

\ref Will C. M., 1992, Int. J. Mod. Phys. 1, 13

\ref Will C. M., 1993, Theory and experiment in gravitational physics,
 Cambridge Univ. Press, Cambridge

\ref Yoder C. F., 1981, Phil. Trans. R. Soc. London A303, 327

\bye
