IVS 2004 General Meeting Proceedings

Ring Downs of Free Core Nutations in the GSFC and USNO VLBI Nutation Series

D. E. Smylie(1), Andrew Palmer(2)

(1)Department of Earth and Space Science and Engineering, York University; (2)Department of Physics and Astronomy, York University

Abstract:

The VLBI nutation series from GSFC and the USNO, now both in excess of 23 years in length, allow the exploration of the temporal behaviour of the retrograde (RFCN) and prograde (PFCN) free core nutations.

Our analysis implements Singular Value Decomposition to solve the least squares problem of fitting the Discrete Fourier Transform to the non-equispaced VLBI nutation observations. A novel feature of our procedure is to use the Parseval relation to determine the number of singular values of the coefficient matrix to be eliminated.

We report the observation for the first time of the prograde mode predicted by Jiang (1993). We find periods of $-417\pm 5$ and $-413\pm 4$ days for the retrograde mode and $388\pm 2$ and $387\pm 5$ days for the prograde mode. Spectral analyses at 500 day intervals down the two records show the RFCN in both the GSFC and USNO data to be in free ring down. The PFCN in the GFSC series appears to be in a similar free ring down. The PFCN in the USNO series does not show similar behaviour, possibly because it is at the noise level.

Free decay of the nutation modes allow the direct measurement of the viscosity at the top of the core. We find a value of $795\;Pa\cdot s$ in contrast to the commonly cited value of $0.008\;Pa\cdot s$ found by Gans (1972) by the extrapolation of laboratory measurements.

   
1. Introduction

Earth's outer liquid core is capable of nearly rigid-body rotations relative to the shell (mantle and crust). Although the motions are mainly confined to the core (Jiang, 1993), conservation of angular momentum produces very small nutations of the rotation axis in space. The theoretical possibility of such motions has a long history (Rochester, Jensen and Smylie, 1974), going back to the classical literature. Only with the emergence of the VLBI technique (Rogers et al., 1983), in providing a close approximation to an inertial space frame, have these motions been brought into the realm of observation (Herring, Mathews and Buffet, 2003).

The VLBI nutation observations are inherently non-equispaced, since measurements can only be made when the radio sources are visible. Working in the frequency domain, we find the Discrete Fourier Transform (DFT) of the non-equispaced record by minimizing an objective function that weights the error power between the DFT and the measured values in inverse proportion to the squares of their standard errors. The Singular Value Decomposition (SVD) method (Golub and van Loan, 1989) is used in solving the resulting large least squares system of conditional equations. An innovation in our application of the SVD technique is the use of the Parseval relation in determining the number of singular values to be eliminated.

The long VLBI nutation series, now in excess of 23 years, allows us to explore the temporal behaviour in the frequency domain of both free core nutations. We discover that they are in free ring down. This opens the possibility of an accurate measure of the viscosity just below the core-mantle boundary. Since the motion there is largely a rotational shear, simple Ekman boundary layer theory permits the decay rates to be directly related to viscosity.

   
2. Fourier Analysis in the Non-Equispaced Case

When the time sampling is uneven and $g_j$ samples are at arbitrary times $t_j$, we make a least squares fit, $g_j^\prime$, of $\left(2L+1\right)$ points to the sum

$$g_j^\prime =\frac{1}{M}\sum_{k=-N}^{N} G_k e^{i2\pi k\Delta ft_j},$$ (1)

where $N\le L$. The $\left(2N+1\right)$ values of $G_k$, once calculated, give a frequency domain representation that is equispaced with sample interval $\Delta f=1/M$.

The VLBI nutations series have accompanying standard errors $\sigma_j$. If the error series in the representation is $\epsilon_j=g_j-g_j^\prime$, we fit to the $G_k$ that minimize the objective function

$$I=\sum_{j=-L}^L \frac{\epsilon_j \epsilon_j^*}{\sigma_j^2},$$ (2)

where the superscript asterisk indicates complex conjugation. The resulting conditional equations are

$$\sum_{k=-N}^N C_{m-k} G_k=d_m,\;\;m=-N,\cdots,N,$$ (3)

with

$$C_m=\sum_{j=-L}^L \frac{1}{\sigma_j^2}e^{-i2\pi m\Delta ft_j... ...M\sum_{j=-L}^L \frac{g_j}{\sigma_j^2} e^{-i2\pi l\Delta ft_j}.$$ (4)

Although the coefficient matrix in the system of equations (3) is equidiagonal and Hermitean, and therefore of Toeplitz form, the Levinson algorithm is found inadequate for their solution even when implemented in double quad precision (Hida, Lee and Bailey, 2000). Instead, we employ the method of Singular Value Decomposition (SVD) considered in the next Section.

   
3. Singular Value Decomposition and Parseval's Theorem

The reason for the failure of the Levinson algorithm, when applied to the solution of the system (3), is that the coefficient matrix for large $N$ is nearly singular (Lawson and Hanson, 1974). The DFT values are typically very large, and the representation (1) of the original time sequence always has isolated values greatly in error. For such large, least squares systems, there is a wide range of literature which suggests the application of Singular Value Decomposition (Press et al., 1992, Chs. 2 and 15).

In matrix notation, the system (3) is ${\bf C}\cdot{\bf G}={\bf d}$. In SVD, the matrix ${\bf C}$ is factored into the triple product ${\bf C}={\bf U}\cdot{\bf W}\cdot{\bf V}^{*T}$. ${\bf U}$ and ${\bf V}$ are unitary matrices and ${\bf W}$ is a diagonal matrix with the singular values of ${\bf C}$ in descending order down the diagonal. In our analysis, we represent the position of the celestial pole as the complex valued sequence $g_j$. The system (3) is then complex as are the unitary matrices ${\bf U}$ and ${\bf V}$. The singular values and ${\bf W}$ are real as a consequence of the Hermitean symmetry of ${\bf C}$. We have implemented algorithm 358 of Businger and Golub (Businger and Golub, 1969) for the SVD of a complex matrix in both double and quadruple precision.

Since the unitary matrices ${\bf U}$ and ${\bf V}$ have inverses which are the complex conjugates of their transposes, the factored form of ${\bf C}$ allows the solution to be written ${\bf G}={\bf V}\cdot{\bf W}^{-1}\cdot{\bf U}^{*T}\cdot{\bf d}$.

Implementation of the SVD procedure involves the inverses of the smallest singular values in ${\bf W}^{-1}$ being set to zero, giving a modified diagonal matrix ${\bf W}^{\prime -1}$. This leads to an approximate DFT, ${\bf G}^\prime$, such that the residual error vector in the solution is minimized in length (Press et al., 1992, pp. 53-56). ${\bf G}^\prime$ is found as the product of a matrix and a vector as ${\bf G}^\prime={\bf V}{\bf W}^{\prime -1}\cdot{\bf U}^{*T}{\bf d}$.

The main difficulty in implementing the SVD procedure, as described in the existing literature, is choosing the number of singular values to discard. Our solution to this problem is to use the Parseval relation between the mean square in the time domain and the mean square in the frequency domain. Before any singular values are discarded, the DFT is many orders of magnitude too large. As the smallest singular values are discarded, its magnitude decreases monotonically until the Parseval relation is closely satisfied.

The sinusoids in the representation (1) are orthogonal in addition over equispaced time sample points, yielding the Parseval relation in the discrete case,

$$\frac{1}{2L+1}\sum_{j=-L}^L g_j^\prime g_j^{\prime *}=\left(\Delta f\right)^2 \sum_{k=-N}^N G_k G_k^*.$$ (5)

The left hand side of this expression is just the power in the time domain signal. For stationary sequences it is $E\{g_j^\prime g_j^{\prime *}\}$, the expected value of the mean square in the time domain. In the case of stationary sequences, the expected value of the mean square will be independent of whether or not the sampling is equispaced or non-equispaced. For stationary sequences, the Parseval relation holds for both cases.

Our procedure then is to discard the smallest singular values in $W^{-1}$ until the Parseval ratio $R$ given by

$$R=\frac{\left(\Delta f\right)^2\sum_{k=-N}^N G_k G_k^*}{E\{g_j^\prime g_j^ {\prime *}\}}$$ (6)

is closest to unity.

   
4. Spectral Analysis of the GSFC and USNO Nutation Series

The Goddard Space Flight Center (GSFC) nutation series were extracted from the file 2003b.eops on their website (http://gemini.gsfc.nasa.gov/solutions/2003b/2003b.eops). The file contains the residuals in longitude, $\Delta \psi$, and obliquity, $\Delta \epsilon$, compared to the IAU 1980 nutation series, and their standard errors, in milliarcseconds. The series runs from Julian day $2,444,089.993750$ (August 3, 1979) to Julian day $2,452,705.270139$ (March 6, 2003) and consists of 3369 points. Tests reveal that four of the points have null values and 23 have identical time tags. The former were removed and the values of the latter were averaged with weights in inverse proportion to the squares of their standard errors. A series of 3,342 points were left for analysis.

Nutation series from the United States Naval Observatory (USNO) were extracted from the file usn2002c.eops obtained by anonymous ftp (ftp://cddisa.gsfc.nasa.gov/vlbi /ivsproducts/eops). Again residuals in longitude and obliquity compared to the IAU 1980 nutation series are provided, along with standard errors in milliarcseconds. The series runs from Julian day $2,444,089.993750$ (August 3, 1979) to Julian day $2,452,719.269444$ (March 20, 2003) and consists of 2,944 points. Fifteen of the points were found to have identical time tags and these were again averaged with weights in inverse proportion to the squares of their standard errors. A series of 2,929 points were left for analysis.

Before performing spectral analysis on these two series, they were first converted to series in the Cartesian coordinates of the celestial pole, $x=\Delta\psi\sin\epsilon$, where $\epsilon$ is the obliquity with sine $0.39777716$, and $y=\Delta\epsilon$, with corresponding conversion of the standard errors. Linear trends and 18.6 yr (6,798.58 d), 9.3 yr (6,798.58/2 d), annual (365.2597 d) and semi-annual (365.2597/2 d) periodic terms were then removed by least squares fits. Spectral analyses were performed on the residuals of the complex valued series of positions of the celestial pole, $g_j=x_j+iy_j$.

The GSFC series covers a time span of $8,615.276 d$. The total record length was taken as $8,617 d$ and it was divided into four segments of length $M=4,924 d$ with 75% overlap. The USNO series covers a time span of $8,629.276 d$ and was divided into four segments of length $M=4,932 d$, again with 75% overlap. Before the DFTs of the segments were calculated, they were windowed for $-M/2<t<M/2$ with a Parzen window to suppress finite record effects.

Power spectral estimates were calculated as the average over the four segments, for each series, of the squared magnitudes of the DFTs, normalized by division with the integral of the square of the Parzen window (Smylie, Hinderer, Richter and Ducarme, 1993). They are shown plotted in Figure 1.

 

Figure 1: Average spectral estimates of the free core nutation modes. Those for Goddard Space Flight Center data are shown on the left, while those for data from the United States Naval Observatory are shown on the right.

Parameters for the Retrograde Free Core Nutation (RFCN) and the Prograde Free Core Nutation (PFCN) were recovered by fitting resonances of the form

$$\frac{A^2}{1+4Q^2\left(\frac{f-f_0}{f_0}\right)^2}.$$ (7)

Recovered parameters are listed in Table 1.

 

Table 1: Fitted parameters of the RFCN and PFCN resonances in the average spectral estimates of the GSFC and USNO series shown in Figure 1.

		RFCN
	$A^2\left(mas^2/cpd\right)$	$f_0\left(10^{-3}cpd\right)$	Period(days)	Q
GSFC	122.44	-2.39773	-417.061	6.4263
USNO	135.46	-2.42368	-412.596	5.4113
		PFCN
GSFC	22.932	2.57627	388.158	5.7932
USNO	22.547	2.58356	387.063	5.8856

The power spectral estimates for each segment were based on DFTs calculated by Singular Value Decomposition, eliminating the smallest singular values until the Parseval ratio $R$, defined in expression (6), is as close as possible to unity.

   
5. Poinsot Representation of Nutations and Wobbles

In the Earth frame, free core nutation modes have an associated nearly diurnal retrograde wobble. The relation between nutation and wobble is best understood by the Poinsot construction. In the case of the RFCN, a small body cone rolls, without slipping, once per sidereal day on the inside of a large space-fixed cone, as illustrated on the left of Figure 2. The PFCN can be described by a small body cone rolling, without slipping, once per sidereal day on the outside of a large space-fixed cone, as illustrated on the right of Figure 2.

The kinematics of the wobble-nutation modes, illustrated by the Poinsot constructions, allow the amplitudes and Q's to be related by

 

Figure 2: Poinsot constructions for the RFCN (left) and the PFCN (right).

$$a_W=\pm\frac{f_{N_0}}{f_{N_0}-1/T_s}a_N,$$ (8)

and

$$Q_W=\pm\frac{f_{N_0}-1/T_s}{f_{N_0}}Q_N,$$ (9)

with the positive sign referring to the RFCN and the negative sign referring to the PFCN.

The wobble and nutation frequencies are related by $f_W=f_N-1/T_s$. The wobble period $T_W$ and the nutation period $T_N$ are related by

$$T_W=\frac{T_N T_s}{T_s-T_N}.$$ (10)

Estimates of the squares of the nutation and wobble amplitudes can be made from the peak spectral densities $A_N^2$ and $A_W^2$, returned from fits to the resonances of the form (7), by integrating between half-power points. The estimates of the squares of the amplitudes are then

$$a_N^2=\mp\pi A_N^2 f_{N_0}/\left(2Q_N\right)\;and\;a_W^2=-\pi A_W^2 f_{W_0}/ \left(2Q_W\right).$$ (11)

Amplitudes, periods and $Q$'s for the average spectral estimates shown in Figure 1, as calculated by the foregoing relations, are listed in Table 2.

 

Table 2: Parameters recovered for the average spectral estimates shown in Figure 1.

			RFCN
	$a_N\left(\mu arcseconds\right)$	$T_N\left(days\right)$	$a_W\left( \mu arcseconds\right)$	$T_W\left(days\right)$	$Q_N$	$Q_W$
GSFC	267.881	-417.061	0.63902	-0.994891	6.4263	2,693.9
USNO	308.711	-412.596	0.74437	-0.994865	5.4113	2,244.2
			PFCN
GSFC	126.566	388.158	0.32602	-0.999839	5.7932	2,249.0
USNO	124.686	387.063	0.32208	-0.999846	5.8856	2,278.5

   
6. Decay of the Free Core Nutation Modes and Viscosity at the Top of the Core

Spectral densities of the four individual segments entering the average spectral density shown in Figure 1 suggest the modes are in free ring down.

In order to examine the temporal behaviour more fully, we have divided the GSFC and USNO nutation series into 2,000 day segments advancing down the time axis in 500 day steps. Spectral analyses were performed on each of the 14 segments for each series by Singular Value Decomposition, using the Parseval criterion to determine the number of singular values to be eliminated. Average spectral estimates were made, based on each of four successive segments with 75% overlap. Before the analysis of each segment, the series were windowed with a Parzen window. The result was 11 spectral estimates for each series, centered at 1,750 days into the record and at 500 day increments down the time axis thereafter. Fits to the RFCN and PFCN resonances of the form (7) were used to recover amplitudes through expression (11).

In free decay, the nutation amplitude follows the exponential decay scheme

$$a_N=a_{N_0} e^{-\pi t/\left(Q_N T_N\right)},$$ (12)

with $a_{N_0}$ denoting the amplitude at time $t=0$. On taking logarithms of both sides, if the mode is in free decay, we obtain the linear function of time, $\log a_N=ct+d$, with $c=-\pi\log e/\left(Q_N T_N\right)=-\log e/\tau$, $d=\log a_{N_0}$, $t_{1/2}=\tau\ln 2$, where $\tau$ is the e-folding time and $t_{1/2}$ is the half life of the decay. Plots of $\log a_N$ as functions of time for both nutation modes of the GSFC and USNO series are shown in Figure 3.

 

Figure 3: Logarithm of the amplitudes of the free core nutations as functions of time. The upper left plot shows the RFCN from the GSFC series while the upper right plot is that for the RFCN from the USNO series. The lower plots show the results for the prograde mode. Lower left is that for the PFCN from the GSFC series and the lower right plot is for the PFCN from the USNO series. Linear fits to the time dependencies are shown directly on the plots.

The free decay of the nutations allows direct measure of the $Q$s of the modes. In Table 3, we list the fitted parameters $c$ and $d$ and the recovered $Q$s, initial amplitudes $a_{N_0}$, half lives $t_{1/2}$, and Ekman numbers $E_k$, for the plots shown in Figure 3, using the periods of the average spectra given in Table 2.

 

Table 3: Free decay parameters of the retrograde and prograde nutations in the spectra of the GSFC and USNO series of VLBI nutation measurements.

			RFCN
	$c$	$d$	$a_{N_0}$	$Q_N$	$Q_W$	$t_{1/2}$	$E_k$
	$\left(10^{-4}\right)$		$\left(\mu arcseconds\right)$			$\left(days\right)$	$\left(10^{-11}\right)$
GSFC	-1.36446	2.9204	832.53	23.976	10,051	2206.24	10.8135
USNO	-1.16525	2.8525	712.03	28.379	11,769	2583.44	7.8869
			PFCN
GSFC	-1.18929	2.7142	517.85	29.556	11,474	2531.22	8.2976
USNO	0.126810	2.2224	168.88	-278.0	-107,608

The free ring downs of the nutations provide an opportunity to measure the viscosity just below the Core-Mantle boundary (CMB).

Even in realistic Earth models, the free core nutations are closely rigid-body rotations of the outer fluid core with respect to the nearly stationary mantle and crust (shell) (Jiang, 1993). They are both nearly diurnal retrograde wobbles with angular frequencies close to $-\Omega$. For a nearly diurnal retrograde wobble of amplitude $A$, the velocity field at radius $r$ is $-Ar\left(\hat{\bf {\theta}}\sin\Omega t + \hat{\bf {\phi}}\cos\theta\cos \Omega t\right)$. By comparison, the shell is in nearly uniform rotation and the adjustment of the velocity field to that of the shell takes place across an Ekman layer of thickness $\delta$ of $O\left(\sqrt{E_k}\right)$ with the Ekman number given by

$$E_k=\frac{\nu}{b^2 \Omega},$$ (13)

where $b$ is the radius of the core-mantle boundary and $\nu$ is the kinematic viscosity. The leading order boundary layer equations (Moore, 1978; Smylie and McMillan, 1998) yield a boundary layer thickness

$$\delta=\sqrt{\frac{\nu}{\Omega \left(1/2+\left\vert\cos\theta\right\vert\right)}}.$$ (14)

The leading order stresses are

$$\sigma_{r\theta}=\frac{\eta}{\delta}v_\theta,\;\;\sigma_{r\phi}=\frac{\eta} {\delta}v_\phi,$$ (15)

where $\eta$ is the dynamic viscosity. The rate of dissipation of energy per unit area in the motion against these stresses is

$$\frac{de}{dt}=v_\theta\sigma_{r\theta}+v_\phi\sigma_{r\phi}=... ...delta}\left[\sin^2 \Omega t+\cos^2\theta\cos^2\Omega t\right].$$ (16)

With $\rho_0$ representing the density at the top of the core, and integrating over the CMB and over one cycle, the energy dissipated per cycle is found to be

$$\frac{2}{105}\pi^2\rho_0 A^2 b^4 \sqrt{\frac{2\nu}{\Omega}}\left(138\sqrt{3}-37\right).$$ (17)

The energy of the motion is closely

$$\frac{1}{2}I_{oc}A^2,$$ (18)

where $I_{oc}=912\times 10^{34}\;kg\cdot m^2$ is the moment of inertia of the outer core.

The quality factor $Q_W$ is defined as $2\pi$ times the ratio of the energy $E$ to the energy dissipated per cylce, and is then

$$Q_W=\frac{105I_{oc}}{2\pi\rho_0 b^4\sqrt{2\nu/\Omega}\left(138\sqrt{3}-37 \right)}.$$ (19)

The kinematic viscosity just below the CMB can then be retrieved from the $Q_W$ of the ring down as

$$\nu=\frac{11025I^2_{oc}\Omega}{8\pi^2\rho_0^2 b^8\left(138\sqrt{3}-37\right)^2 Q_W^2}.$$ (20)

With $b=3,480$km as the radius of the core-mantle boundary, using the values of $Q_W$ in Table 3, the ring down of the RFCN from the GSFC series gives a kinematic viscosity of $0.0955\;m^2\;s^{-1}$, while the ring down of the RFCN in the USNO series gives $0.0696\;m^2\;s^{-1}$. The ring down of the PFCN in the GSFC series gives $0.0733\;m^2\;s^{-1}$. Since the density at the top of the core is close to $10^4\;kg\cdot m^{-3}$, the respective values of the dynamic viscosities are $955\;Pa\cdot s$, $696\;Pa\cdot s$ and $733\;Pa\cdot s$. These compare with the value of Gans (1972) of $8$ centipoise or $8\times 10^{-3}\;Pa\cdot s$, based on the extrapolation of laboratory data, and commonly used by dynamo theorists! Our values are closer to the value of $7\;m^2\;s^{-1}$, measured by Davis and Whaler (1997), from the spin-up of the core following the 1969 geomagnetic jerk.

   
7. Discussion

Two unusual circumstances have allowed us to make an important direct measurement of the viscosity at the top of the core. First, we now have available quite long VLBI nutation measurement series which allow us to explore the temporal behaviour of the free core nutations. Second, during the more than twenty-three year lengths of the records from Goddard Space Flight Center and the United States Naval Observatory, the free core nutations appear to have been in free decay.

Since, in the Earth frame, the free core nutations are closely rigid-body rotations of the outer fluid core with respect to a nearly stationary mantle at close to a retrograde diurnal period, their ring downs, in combination with simple Ekman boundary layer theory, provide direct measures of viscosity at the top of the core.

The average recovered viscosity of $795\;Pa\cdot s$ is very much larger than the value of $8\times 10^{-3}\;Pa\cdot s$ estimated by Gans (1972), on the basis of the extrapolation of laboratory measurements, and commonly used by dynamo theorists. We are fortunate in the fact that, at these short periods, the Lorentz force and electromagnetic damping effects are completely negligible (Crossley and Smylie, 1975).

The detailed spectral analyses involved have been made possible by the addition of the Parseval criterion for the determination of the number of smallest singular values to be eliminated, in the Singular Value Decomposition technique, used to solve the conditional equations for the calculation of the Discrete Fourier Transforms of the non-equispaced series of nutation observations.

While the free core nutations have now rung down to nearly unobservable levels, at the outset of the observations in 1979, the newly observed Prograde Free Core Nutation, predicted by Jiang (1993), had an amplitude in the GSFC data as high as 62% of the dominant RFCN. Should the modes be excited again, they promise to allow further refinement of Earth property measurements at the top of the core.

8. Acknowledgments

D.E.S. is grateful for financial support from the Natural Sciences and Engineering Research Council of Canada.

Bibliography

15

Businger, P. A. and G. H. Golub 1969. Singular Value Decomposition of a Complex Matrix, Comm. ACM,, 12, 564-565.

15

Crossley, D. J. and D. E. Smylie 1975. Electromagnetic and Viscous Damping of Core Oscillations, Geophys.J.R. astr.Soc., 42, 1011-1033.

15

Davis, R. G. and K. A. Whaler 1997. The 1969 geomagnetic impulse and spin-up of the Earth's liquid core, Phys. Earth Planet. Inter., 103, 181-194.

15

Gans, R. 1972. Viscosity of the Earth's core, J. Geophys. Res., 77, 360-366.

15

Golub, G. H. and C. F. van Loan 1989. Matrix Computations,, 3rd ed., Johns Hopkins University Press, Baltimore, Md.

15

Herring, T. A., P. M. Mathews and B. A. Buffett, 2002. Modeling of Nutation-Precession: Very long baseline interferometry results, J. Geophys. Res., 107, ETG4-1 to ETG4-13.

15

Hida, Y., X. Li and D. Bailey, 2000. Quad Double Arithmetic: Algorithms, Implementation, and Application, Technical Report LBNL-46996, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, October 2000. Available at http://www.nersc.gov/dhbailey/ mpist/mpdist.html.

15

Jiang, Xianhua and D. E. Smylie, 1996. Variational Calculation of the Free Core Nutation Mode, Phys. Earth Planet. Inter., 94, 159-182.

15

Jiang, Xianhua, 1993. Wobble-Nutation Modes of the Earth, Ph.D. thesis, York University, Toronto, Canada.

15

Lawson, C. L. and R. Hanson, 1974. Solving Least Squares Problems, Prentice-Hall Inc., Englewood Cliffs, New Jersey.

15

Moore, D. W., 1978. Homogeneous Fluids in Rotation. In: Roberts, P. H., Soward, A. M. (Eds.), Rotating Fluids in Geophysics, Academic Press, New York.

15

Press, W. H., S. A. Teukolsky, W. T. Vetterling and B. P. Flannery 1992. Numerical Recipes,, 2nd ed., Cambridge University Press, Cambridge UK.

15

Rochester, M. G., O. G. Jensen and D. E. Smylie, 1974. A Search for the Earth's `Nearly Diurnal Free Wobble', Geophys. J. R. astr. Soc., 38, 349-363.

15

Rogers, A., R. Cappallo, H. Hinteregger, J. Levine, E. Nesman, J. Webber and A. Whitney 1983. Very-Long-Baseline Radio Interferometry: the Mark III System for Geodesy, Astrometry, and Aperture Synthesis, Science, 219, 51-53.

15

Smylie, D. E. and D. G. McMillan, 1998. Viscous and Rotational Splitting of the Translational Oscillations of Earth's Solid Inner Core, Phys. Earth Planet Inter., 106, 1-18.

IVS 2004 General Meeting Proceedings