RUSSIAN JOURNAL OF EARTH SCIENCES, English Translation, VOL 1, NO 1, DECEMBER 1998
Russian Edition: JULY 1998
On the effects of the inertia ellipsoid triaxiality in the theory of nutation
S. M. Molodensky
Joint Institute of the Physics of the Earth, 10 Bolshaja Gruzinskaja, Moscow 123810
Abstract. Modern nutation theories of the Earth and planets use the assumption, that the planetTs inertia ellipsoid is symmetrical about the axes of rotation. Below we estimate the effects of the inertia ellipsoid triaxiality accurate to ((B — A) / (C — A))A, where A, B, C are equatorial and polar principal moments of inertia, respectively.
1. Introduction
As known, the nutational motion of the Earth and planets in space is fully determined by Euler angles ‘3e,4>e,'4,e, which describe the position of the body-fixed (Tesserand’s) axes (x,y,z) with respect to an immobile reference frame. At the same time, dynamic equations of motion of an actual planet model define only the components (uix, ujy, uiz) of the angular velocity vector w with respect to a moving reference frame. In the case of an axially symmetrical planet (equatorial moments of inertia satisfy the condition A = B) the Euler angles are connected with these components through well-known Poinsot’s formulae (see, for example, [Landau and Lifshitz, 1964]), but in the case A ^ B Euler’s relations between ‘3e,4>e,'4,e and (u>x,u>y,u>z) present a rather complicated (nonlinear) system of ordinary differential equations which has not been investigated in all details as yet.
As was shown by Zharkov et al. [1996], in the case of small violation of the axial symmetry (if the parameter
r,= (A-B)/(C -A)
satisfies the condition \rj <C 1|) Euler’s equations can be solved by using the perturbation method in powers of this parameter; in such an approach, first-order (in rj) corrections describe only short-period (semidiurnal) perturbations of the values $e, 4>e, ipE- Nevertheless, this result cannot be considered as a final solution of the problem under consideration, because second-order corrections of the order of rf can contain also long-period terms, which perturb not only instant values of the nutational amplitudes, but also their mean values.
©1998 Russian Journal of Earth Sciences.
Paper No. TJE98003.
Online version of this paper was published on July 20 1998. URL: http://eos.wdcb.rssi.ru/tjes/TJE98003/TJE98003.htm
Modern VLBI measurements are most sensible just to these terms, and their analysis is most interesting for any practical purposes. Below we present a simple analytical expression that describes these corrections with the accuracy of the order of r/4.
2. Statement of the Problem
As mentioned above, modern astrometric measurements determine the motion of fixed points of a planetary surface R with respect to an immobile reference frame usually connected with extragalactic objects. At the same time, well-known Liouville’s equation
M + wxM = L, (1)
(where w is the angular velocity of the Tesserand’s reference frame (x,y,z), M is the angular momentum of the planet, L is the torque of external forces, and the dot above symbol denotes the time derivative) determines only the motion of the vector with respect to the moving system of coordinates (x,y,z). To express the components of the vector R in terms of the known values (u>x,u>y,u>z), it is convenient to use the well-known Poinsot’s kinematic relations (see, for example, [Zharkov et al., 1996]) which connect the trajectory of the vector u> with respect to the moving reference frame (x,y,z) (polhode) and the trajectory of the same vector in space (herpolhode).
In the case of an axially symmetrical rigid planet (A = B and the products of inertia 1^ vanish), the components M and w are connected by the well-known relation
Mx = Auix; My = Bujy] Mz = Cuiz; (2) and the torque L may be represented in the form
L = /(i sin(crt) + j cos(crt)), (3)
where i, j are unit vectors oriented along the axes x and y, respectively,
l=(C-A)vt/o
(4)
vt is the amplitude of the nearly diurnal tidal potential,
Vi = -^(xz cos(crt) + yz sin(crt)) =
—f-r2 sin 7? COS 7? cosicrt — <f>) , cr
(5)
a is the mean radius of the planet, (r, <j>) are spherical
coordinates, and cr is tidal frequency.
In the most general case, Liouville equation (1) represents a system of three nonlinear ordinary differential equations with respect to three unknown functions u>x(t),u>y(t),u>z(t). In the case rj = 0, it has a very simple exact solution. In fact, substituting (2)-(4) into (1), we obtain
= w(iecos(<7t) + jesin(crt) + k),
(6)
where k is a unit vector that coincides with the z axis and
£ =
vt{A - C){Aa - ui{A - C)) a2 A2a2ui + (A — C)2ui2
(7)
is the amplitude of nutational motion with respect to the moving reference frame (x, y, z).
The fact that the amplitudes the x and y components of vector w in (6) are equivalent is a direct consequence of the axial symmetry of the problem under consideration. In view of this, the surfaces of polhode and herpol-hode are exactly conical; in accordance with Poinsot’s kinematic relations, the ratio of apex angles of these cones is defined by the relation
= arcsin(—ea / (a + w));
(8)
2. Inverting Euler’s kinematic relations between the known components uix (t), u)y (t), ujz (t) of the vector w in the moving reference frame, the components
ui2 (t), uis(t) of the same vector in space, and Euler’s angles ‘3e,4>e,'4,e, these angles can be expressed in terms of uix (t), uiy (t), ujz (t). Taking into account, that the kinematic relations are nonlinear, we use the first-order method of perturbations with respect to the small parameter ?y.
3. To connect short-period first-order perturbations
of the functions 0Js(t) and long-period
second-order perturbations of the same functions, we use the well-known Poinsot’s theorem stating that the rotation of an arbitrary triaxial rigid body may be represented as the rolling of the polhode over the herpolhode without sliding. In accordance with this theorem, the ratio of trajectory lengths of the vector w in the moving and fixed reference frames L\ and L2 is equal to the ratio of the corresponding periods. The perturbation of the trajectory L2 in turn consists of (1) the known terms which describe semidiurnal oscillations and (2) unknown terms of higher orders of smallness which describe the mean radius of this trajectory. As is shown below, this permits the mean radius of the trajectory L2 to be found with sufficient accuracy (of the order of r/4).
3. Tidal Torque for a Triaxial Body
Using expression (5) for the components of tidal potential in Cartesian coordinates, it is easy to find the tidal torques which act on an axially nonsymmetrical planet
L = / pr x VVjrfr =
vt
(9)
the amplitude of nutational motion of a vector R in the body-fixed reference frame (x, y, z) is equal to a + e in the case of prograde nutational motion (if the directions of the oj motion and nutation coincide) and to in the opposite case of retrograde motion.
At rj 0 the solutions of dynamic and kinematic nonlinear equations can be rather complicated. Because long-period perturbations are most important in practical applications, below we use the method of its approximate solution based on the following procedure:
1. Dynamic equations (1) provide expressions for the components u>x(t), u>y(t), u>z(t) of the vector w with respect to the moving reference frame, in which the effects of the inertia ellipsoid triaxiality are taken into account;
= —« [i(C — B) sin(crt) + j(A — C) cos(crt)],
where r is the total volume of the Earth. Replacing the relation between uiy and My in equations (2) by the formula My = Bu)y and then substituting (2) and (9) in Liouville’s equation (1), we obtain
Auix + (C — B)uiyuiz = —|(C — B) sin(crt); (10a)
Buiy + (A - C)ujxujz = ^j(A — C) cos(crt); (10b)
Cdiz + (B - A)ujxujy = 0. (10c)
Using relations (6) and (7), it is easy to see that the
ratio of terms Cloz and (B — A)uixuiy in (10c) is of the order of rjei’t/(ga), where e is the flattening of the planet and g is the gravity on its surface. This ratio is negligibly small for all planets of the solar system, and the approximation uiz = u> = const is valid with a very high accuracy. In this approximation, the solution of equations (10a) and (10b) describing forced nutation and precession is
where
£l =
vt {B - C){B(T - lo{A - C)) a2 ABa2uj + (A — C)(B — C>3 ’
eui a + ui ’
(13d)
= <jj(l£l cos(crt) +J62 sin(crt) + k) (1 lcl)
(lib)
(if this amplitude satisfies the condition || 1)
4. First-order Approximation for a Slightly Triaxial Body
To estimate the effects of triaxiality of the planet in the first approximation, we first perform some simple transformations. Multiplying (12a) and (12b) by cos ipE and sin^B and then subtracting the results, we express {>E in terms of 6\, £2 as follows:
§E = w(ei cos ipE ~ £2 sin iPe)-
(14a
vt (A — C)(Aa — uj(B — C)) .
2 a2ABa2Lo + (A-C)(B-C)Lo3’ 1 ’
Substituting relations (11) in the well-known Euler relations between Euler angles (^e,<Pe,'Pe) and (u>x,u>y,u>z) [Landau and Lifshitz, 1964]
uix = (j)E sin $e sin ipE + cos ipE; (12a)
uiy = <j>E sin {>e cos ipE + $e sin ipE] (12b)
Wz = <Pe COS&E COS IpE + IpE, (12c)
we obtain a system of three nonlinear ordinary differential equations with respect to three unknown functions
$E , <pE, IpE •
In the case of an axially symmetrical body (A = B and Si = e2 = £, where e is defined by (7)) these equations have an elementary (Poinsot) solution in the form
Analogously, multiplying (12a) and (12b) by sin^_E and cos ipE and then summing the results, we have
<j>E sin §e = w(ei sin ipE + £2 cos iPe) (14b)
To linearize these equations, their solutions are represented by the superposition
<Pe = <Po + <Pi; IpE = Ipo + Ipl, $E = (I5)
where <po,ipo, tfo are the solutions of axially symmetrical problem (13a)-(13c), and <pi,ipi, are small corrections due to the triaxiality effects of the inertia ellipsoid.
Substituting (15) and (11a) into (14a), (14b), and (12c) and neglecting the terms of the order of r/2, we obtain
t>0 sin i90 + <p 1 sin ^0 + $i<Po cos i90 = U!£ 1 COs(<7t)(sin Ipo + Ipl COS ^0) + U}£2 sin(<7t)(c0s Ipo — Ipl sin ^0);
(16a
$0 = arctan
a + ui
£UJ
(13a) cos(<rt)(cos ipo ~ ipi sin ipo) —
UJ£2 sin(crt)(sin Ipo + Ipl COS ^0);
(16b)
and
£U1
sin i90
t + const,
(13b)
(13c
■ SCO
^0 = ----------7- « CT +CJ
Sin »n
i cos i90 + <Pi cos i90 — tfipo sin +ipo + 4>i = 0.
(16c
which show, that the angular frequency and amplitude of nutational motion in space are
Transforming, in accordance with (13b), the main terms in the right-hand sides of these equations,
U!£l COs((7t) sin Ipo + 0J£2 sin((Tt) COS Ipo =
^-H-£2 sin(V’C) + vt) + ^£l e'2 sin(V’C) - o-t) =
£l + £2 £1 — £2 x
ui-----------1- ui--------cos(2crt);
and
U!£ 1 COs((7t) COS Ipo ~ U!£2 sin((Tt) sin Ipo =
Si + £2
COs(lpo + <rt)
£1 — £2
^£]_£2 cos^Q _ ^ _ sin(2crt)
Since the system of equations (19) is of the third order, it has three linearly independent homogeneous solutions. It is easy to show, however, that these solutions are not interesting for the problem under consideration, and taking into account, that the zero-order values of because they describe the small variations in Euler an-<f>o,^o,ipo are connected with £,ui and a by relations gles caused by an arbitrary infinitesimal rotation of the (13), we obtain, after differentiation of these relations fixed reference frame.
with respect to the time,
The inhomogeneous solution of (19) has the form
60 sin $0 = Ul£ = U1
^o = 0;
Si + £2
(17a
(17b)
where
61 = —ip 1 = asin(2crt); ^ = bcos(2crt),
(20)
I cos^o + ^0 = 0.
(17c
Using these relations, equations (16) can be rewritten as follows:
' ' £ 1 & 2
<t> 1 sin 7?0 + ^1^0 COS 7?0 — IplW--------- sin(2crt) =
2 (18a)
£l _£2 ^
= u)---------cos(zcriJ;
and
{si - S2)(3(7 + uj)(a + lo)) £(ui2 + 2 (tui — 3 a2)
b = —a-
£U1 (T + u>
+ ipiui£l ~^£'2 =u£l 2 £2 sin(2(Tt); (18b)
4>1 COS $0 — (^0 sin $0 + i’l- (1®C)
Eliminating i90 and <p0 from these relations with the help of (13a) and (17c), we finally obtain the system of three inhomogeneous linear ordinary differential equations
Ci4>i + 02^1 = a cos(2<rt);
+ c3ipi = asin(2crt); (19)
ipi + cA(j) 1 - c3$i = 0.
with the coefficients
_ • q _______ __________su______________.
i i \2 1 i 1 \2\l/2 ’
((£UJ)2 + (a + uj)2)1/2 C2 = U!£ cot ^0=0"+^;
C3 = UJ£\
_ 0 _ ___________a + U!___________
C4 COS Vq \2\l/2 ’
((£uj)2 + (a + uj)2)1/2
and
a = cj(£ 1 - e2)/2.
These expressions fully determine first-order short-period perturbations of the nutational and precessional motion caused by the triaxiality of the planetTs figure.
To find short-period motion perturbations of the vector w, we introduce a fixed Cartesian reference frame (ei, &2, 63) in such a way that the axis e3 coincides with the axis of the unperturbed cone of herpolhode, and (pE is equal to the angle between the nodal line and vector ei. Then the relations between components (u>i,u>2,u>3) and Euler angles ($e, 4>e, iPe) are completely analogous to (12a)-(12c) after replacement of <p>E by ipE and ipE
by <Pe-
uii = ipE sin $e sin <Pe + $e cos <pE] (21a)
U>2 = IpE sin {>E cos (pE — $E sin <j>E] (21b)
UI3 = IpE COS $E sin (pE + (pE- (21c)
In the case of nutational motion of any actual planet, its spatial amplitude £u>/(cr + u>) does not exceed a few seconds of arc, and we can use, instead of (13a), a more simple expression
Substituting (15), (13b), (13c), (22), and (17) into (21), we have
eui
uji = (—a — 2acr cos(2crt))(----------------b
cr + ui
+b cos(2crt)) sin((cr + uj)t + a sin(2crt)) — —2crb sin(2crt) cos ((cr + uj)t + a sin(2crt));
SCO
ui2 = (—cr — 2acr cos(2<rt))(------------------b
(T + w
+& cos(2crt)) cos((cr + ui)t + a sin(2crt))-|-+2crb sin(2crt) cos ((cr + uj)t + a sin(2crt));
W3 = w.
(23a
(23b)
(23c
(24a
L2 = 40)(l + /2(5-A)),
(24b)
part of this function (B — A) which is connected with
short-period (semidiurnal) perturbations of herpolhode. Using this value, the total perturbation of herpolhode length can be represented as a sum of the known short-
(B — A) and unknown long-(B — A) that describes variation in the “mean” herpolhode radius:
period perturbation f2S
period perturbation f^1
L^f2(B-A) =
= Li°)(f^(B-A) + f21>(B-A))
(i)
(25a
Like the first-order corrections to Euler angles (20), these expressions contain only short-period (semidiurnal) perturbations.
5. Long-period terms in the second-order approximation
Now, we will use the well-known Poinsot’s theorem stating that the rolling of polhode over herpolhode takes place without sliding. On the strength of this theorem, a general method for the calculation of long-period herpolhode perturbations may be formulated as follows. Let the total length of the vector w trajectory in the moving reference frame (x,y,z) is defined by the expression
where Lis the length of the vector in the case of an axially symmetrical planet (A = B) and fi(B — A) is a small correction due to +effects of triaxiality. Because this trajectory coincides with circle (6) if A = B and with ellipse (11a) if A ^ B, the function /i actually describes the known difference between the length of ellipse (11a) with the semiaxes u>Si,u>S2 and the length 27TUJS of circle (6).
Likewise, using relations (23), the total length of the w trajectory in the fixed reference frame (ei, e2, 63) can be represented in the form
In accordance with Poinsot’s theorem, the condition of the polhode rolling over herpolhode without sliding is that the trajectories described by the angular velocity vector during equal time intervals have the same lengths in moving and fixed reference frames:
Li/Ti = L[°\l + f1(B-A))/T1 = = Li°\l + f2(B-A))/T2 = = \(T/(cr + ui)\l[ ^(1 + /2 \b — A) + +f(1\B-Aj)/T2,
where T\ and T2 are full revolution periods of the vector w in the moving and fixed reference frames, respectively. Since these values are connected through the relation
T2 = \(t/((t + ui)\T1,
the unknown value /^(B — A) can be expressed in terms of the known functions /1 (B — A) and (B — A) in the following, very simple form:
(B - A) = h (B — A) — f2s> (B - A) (25b)
Then, variation in the mean radius of herpolhode (which coincides with variation in the nutation amplitude SA multiplied by ui) may be represented in the form
SA = A0fi1)(B-A) =
= A0(f1(B-A)-f^(B-A)),
(26)
where L^ = |cr/(cr + w)|_Llj0'1 is the total length of the herpolhode in the biaxial case and f2(B — A) describes the difference between the lengths in the triaxial and biaxial cases.
As mentioned before, the absence of exact expressions for the Euler angles in terms of known components uix(t), uiy(t) precludes direct determination of this function. However, expressions (23) readily provide the
where Ao is the unperturbed nutation amplitude corresponding to the axially symmetrical approximation.
In accordance with (11), the complete trajectory length of the vector w in the moving reference frame
(X, y, Z) IS
/>27T
Li = ui (ef cos2 A + e\ sin2 A)1/,2c?A =
Jo
r2V?+e!
c2 c2 £1 - £2
cos(2A))1/,2c?A.
ui1 = ki cos(<f>o + (f)i) + k2 sin(c/>0 + 4>i);
Developing the integrand function as the TaylorTs series in small parameter k = (e\ — s2)/(s{ + £2) with the help of (lib) and (11c) and neglecting terms of the order of k4, we have
uj2 = ki sin(</>0 + <j)i) + k2 cos(<j>0 + </>i); (30b)
r2lT c2 _L c2
ii="i «£4£i»1'1-
-a-3L.»!(2A))<iA =
8V2(e2 + e2)3/2 1 ”
[Sl-62)2 (£?-£!)%
(27a
= 27tuj(s +
= 27xoj(e +
8e 64e3
(ei - S'j)2 x
16e
Comparing this expression and relation (24a) yields
where
u>3 = 0,
k\ = —ecruj — (bcr(cr + w)) +
+2acr
, eui
' (T + UJ
i(2crt) —
— (2aba(2a + w) + 4a2———) cos2(2<rt);
(30c
(31c
(e 1 ~ S'j)2 16e2
(27b)
ecruj
where the nutation amplitudes e\, e2 and e are defined by expressions (lib), (11c), and (7), respectively.
Strictly speaking, the definition of the complete trajectory length of the vector ui in the fixed reference frame (ei,e2,e3 needs some refinement, because this trajectory is generally nonclosed. In fact, it is obvious that, since the period of short-period perturbations is equal to 7r/<r and the period of nutational motion is equal to 27t/|<t + ui\, this trajectory is closed only if cr = k(cr + ui), where k is any integer number. Of course, nutation frequencies of actual planets do not satisfy this condition, and the herpolhode is not a closed curve.
In order to apply Poinsot’s theorem to the case of nonclosed trajectory of the angular velocity vector in space, we define the “mean length” of this trajectory as follows:
k2 = (4a------------------b 2crb(2cr + ui)) sin(2<rt) =
<7 + UJ
9
(7UJ
= —2 as-----------sin(2<rt).
CT+UJ [ ’
(31b)
2nN
N-
N
(28)
Substituting (23) into (28) and expanding, as before, the integrand function in Taylor’s series with the accuracy of the order of k4, we obtain
N
I £(TbJ I
L 2 = I———|x a + ui
v 1 [2ft n s s2 s3.,
bin — I (1 + - - — + rrr)dt
(29a
N
'o
8 16y
where
Substituting these expressions into (29), (30a), and (30b) and neglecting the terms of the fourth and higher orders, we have
(uji)2 + (uj2)2 + (W3)2 = kf + k2 =
= (ewer)2 (1 — 2 a-----------cos(2crt) +
y y uj + cr y ’
+a'
(ui — cr)2 — 4 a (a + ui)
(a + ui)2
+4 az
(cr + ui)2
sin2 (2 crt));
2ttN
ECU) 1 f <7 + u' (J
l-^-l lim — (1 + -)dt
'ct + uj'n^00 N J0 v 2’
, 5uj — 3cr — 6uia 4 (a + uj)2 >
<7 + UJ JV—>-oo N !<»>(„=
, ,3*»
1 ["+“ S2 , — dt
4(cr + ui)2 '
6 dt = 0.
(32a
(32b)
(32c
(33)
[ecruj
r[(wi)2 + (w2)2 + (w3)2] - 1.
Substituting (20a-c) into (26a, b), we have
(29b)
Summing these expressions, we finally obtain
L2 = \a/((T + cj)\L{0\l + fis\B-A)), (34a)
cos
where
(si — 82)^ (3(7 + 0j)2(0j2 — UtCT — <72) 4 £2(oj2 + 2crbj — 3<t2)2
(si — £2)2(u2 — UtCT — (72)
(34b)
4e2(cr — ui)2
Substituting (34b) and (26b) into (23), we obtain the final value of the correction to the mean nutation amplitudes
A = A0 + SA = A0(l + f(1}(B - A)) =
— j4q(1 +
(e 1 — S2)2 ,1 (w2 — uicr — <72)
— j4q(1 +
4e2 4 (<7 — ui):
(£1 ~ £2)2 (<r + gj)(5ct - 3u), 16e2 (<7 — Ul)2 '
)) =
(35)
respectively, where To = 2'k/ui is one sidereal day). Nu-tational motion with the period Ts = To/2 takes place only in the case of synchronised diurnal and orbital rotation (if at the points of equinoxes the orientation of the planet’s principal equatorial moments of inertia is invariable in space).
In this case, the nutation amplitudes are obviously determined not by the mean values of A, B averaged over the total period, but by their “effective” values in the vicinity of points at which the tidal torque is maximum. In this specific case, the corrections to the mean nutation amplitudes are proportional to (£1 £2)/£ rather than (£1 £2)2/£2 • This qualitative consideration explains, why second-order correction (35) tends to infinity in the limiting case a —?> u>.
Acknowledgement. This work was supported by the Russian Foundation for Basic Research, project No. 98-0564419.
Note that this correction vanishes if a = —u> or cr = 3w/5. The first root corresponds to the precession frequency. Vanishing of this expression for the precession frequency implies that the precessional constant depends only on the value C — (A + B)/2 but not on the difference of equatorial moments of inertia B — A.
The second root, as well as the denominator root cr = —ui, corresponds to the case of high nutation frequencies in space as compared with the period of diurnal rotation (these frequencies correspond to the spatial periods of nutational motion Ts = 5Tq/8 and Ts = To/2,
References
Landau, L. D. and E. M. Lifshitz, Mechanics, Moscow, Nauka, 1964 (in Russian).
Molodensky, S. M., Tides, Nutation, and Internal Structure of the Earth, Moscow, Nauka, 1980, 214 pp. (in Russian).
Zharkov V. N., S. M. Molodensky, E. Groten, A. Brzezinski, and P. Varga. The Earth and its rotation. Low-frequency geodynamics, Heidelberg, Wichman Verlag, 1996, 531 pp.
(Received May 15 1998.)