Seismological observation of Earth's oscillating inner core.

We investigate the differential rotation of Earth's inner core relative to the mantle using pairs of precisely located nuclear explosions. We find that the inner core subrotated at least 0.1° from 1969 to 1971, in contrast to superrotation of ~0.29° from 1971 to 1974. These observations contradict models of steady inner core rotation and models that posit much faster rotation rates. The reversal of polarity, timing, and rotation rates is consistent with a model of oscillations about an equilibrium with gravitational locking of the mantle and inner core due to lateral density variations. The model, which has a 6-year period, can explain the variation in the length of day, which has oscillated fairly steadily for the past decades. Inner core oscillation would also allow interpretations of causal connections between inner core and mantle lateral variations, which are problematic if the inner core consistently superrotates.

INTRODUCTION

Earth’s inner core (IC), decoupled from the overlying mantle by the liquid iron outer core (OC), has been inferred from close observations to move and change over decades. However, the pattern, speed, and driving force of the change are still under debate.

The structure, motion, and evolution of the IC matter as the growth of the IC provides energy for the geodynamo, which is responsible for Earth’s critical magnetic field, with its well-known but poorly understood fluctuations and periodic reversals. It may or may not have formed, with its hemispherical division in structure, in alignment with the heterogeneity in the overlying mantle.

Initial models of fairly steady superrotation of Earth’s IC relative to mantle () have slowed from 1° per year to 0.05 to 0.15° per year (–) and have been challenged by models with only deformation of the IC boundary (ICB) (, ), seismological models with strong variations in rotation rate (, ), and geodynamic models with small equilibrium oscillations (–). It is fair to say that there is no current consensus.

RESULTS

We estimate the differential rotation of the IC by measuring changes in backscattering within the IC of P waves of a pair of nuclear explosions. We did this previously for another pair of explosions (–), and here, we follow the same procedure. Our data are two explosions recorded by the Large Aperture Seismographic Array (LASA) () (see fig. S1): Milrow in 1969 and Cannikin in 1971 (). LASA comprised more than 300 vertical-component stations across a 200-km aperture, emplaced 60 m deep in boreholes in the firm ground (). These data were heroically recovered from decaying magnetic tapes by Hedlin et al. (). The explosions Milrow and Cannikin (Fig. 1A) were an estimated 1 and 4 megatons and separated by 8 km under Amchitka Island ().

Fig. 1.

Map of LASA and nuclear tests.

(A) The locations of the two nuclear test pairs and LASA. The black triangle marks the location of LASA. The magenta stars show the locations of the pair of twin nuclear tests used in this study. (D) The raypaths of the PKiKP and IC scattered waves. (B and C) The distribution of the predicted time shifts for scattering mapped to its location at the IC boundary (ICB), for the case of 1° eastward rotation of the IC relative to the mantle. The R and T arrows in (B) and (C) illustrate the directions of the radial and transverse slownesses in Fig. 2, respectively. CMB, core-mantle boundary.

Fig. 2.

Data and prediction for ICS beams.

(A and D) Average beam amplitude from 1025 to 1250 s for the 1971–1974 pair and 1969–1971, respectively. The black dashed circles indicate the slownesses expected for the rims of the IC and OC observed at LASA. The magenta dashed line indicates the slownesses expected for the limit of the IC illuminated by direct waves from the events. The magenta dot is the predicted slowness of PKiKP, and the black dot is vertical incidence. (B and E) Average observed time shift. (C and F) The time shift predicted for 1° of superrotation. (G) Radial slice of the observed time shift averaged from transverse slowness of −0.01 to 0.01 s/km. Here, the slowness has an empirical correction ().

IC scattering (ICS) was originally identified by Vidale and Earle (), and models of the scatterers have been refined by Peng et al. () and Wu and Irving (). Recently, the three-dimensional heterogeneity in the IC has been further investigated ().

The ICS waves change over time if the IC rotates, with the waves scattered within receding regions arriving progressively later and waves scattered in approaching regions coming in earlier in the ICS wave train. Time shifts of the ICS waves observed from LASA for a pair of Novaya Zemlya (NZ) blasts in 1971 and 1974 conclusively showed 0.29° of IC superrotation in those 3 years ().

The NZ explosions that we previously studied were much closer, separated by only 1 km, and more closely matched in power. The greater separation, difference in yield and therefore burial depths, and spall timing () and the path through the Aleutian slab () make the seismograms of the Amchitka pair of explosions more difficult than those from the NZ pair to correlate to isolate waveform time shifts.

We take the same approach as Wang and Vidale () to estimate the IC differential rotation from 1969 to 1971 from the Amchitka explosions, simplifying it because of the reduced resolution arising from lower seismogram correlation. Predicted time shifts for an IC superrotating 1° for the Amchitka and NZ geometries (Fig. 1, B and C) are nearly orthogonal because of their pole-perpendicular and pole-parallel paths, respectively (Fig. 1A).

We filter the seismograms to retain the bandpass at 1 to 3 Hz. Most energy is less than 2 Hz, and the energy that is present above 2 Hz is too incoherent across the aperture of LASA to usefully stack. Arrays with much smaller apertures have examined signals up to 5 Hz (), but LASA’s large aperture, or perhaps its underlying scattering and attenuation structure, precludes coherence of such high-frequency waves.

We beamform the seismograms as described in the work of Wang and Vidale () using all available LASA stations. The 1969 and 1971 events had 323 and 331 stations, respectively, of which we used the 322 that were in common, while the 1971 and 1974 events had only 199 in common. Static shifts, on the order of tenths of seconds, correct for wavefront perturbations due to crust and mantle structures under LASA. The statics were estimated by aligning with cross-correlation according to iasp91 predictions of the observed PKIKP phases of an earthquake near the antipode ().

The beam amplitudes and the differential time shifts are shown as functions of slowness and time in Fig. 2. Both the Amchitka and NZ explosion pairs produced prominent ICS waves. The ICS arrivals start at the time and slowness of PKiKP and were bigger than the simultaneous late P-wave coda, which would appear off the plots at slownesses around 0.06 s/km. Both aspects are in contrast to the explosions in NZ, in which PKiKP was subtle and P coda dominated over ICS waves. The Amchitka’s ICS waves then expand to fill the range of slownesses possible for scattering across the IC. The scattering is strongest in the quadrant with positive radial and transverse slowness. This is the direction to a patch of IC strong scatterers mapped from an analysis of 73 earthquakes and explosions ().

The NZ explosion pair is one of the best-documented cases to prove the IC superrotation (–). The pattern of red at the positive transverse slowness and blue at the negative slowness (Fig. 2B) is consistent with that predicted for IC superrotation (Fig. 2C). The crucial observation in Fig. 2E is the observed predominance of red over blue areas at positive radial slowness. This pattern, although noisy, indicates that the ICS waves in those regions are arriving later in 1971 than in 1969, and the waves in the blue patches at lesser slownesses are, in contrast, arriving earlier, which is the opposite of the predicted pattern for superrotation (Fig. 2F).

The beamforming parameters used to generate Fig. 2 are critical as the signal is noisy. We use 10-s windows to capture about 10 to 15 cycles for the correlation; longer windows avoid locking into false minima and reduce the number of independent measurements but may incorrectly mingle scattering with different time shifts from different depths. We require a correlation greater than 0.4 for inclusion in the beam; too high a threshold only returns time shifts for a small fraction of the beam volume; too low a threshold returns more noise. We require an average beam amplitude greater than 10% of the maximum; very low amplitude parts of the beam did not show distinct signals and might be dominated by background coda. Last, we only searched for correlations to a maximum of ±0.25-s time shifts; most signals that appear plausible are considerably smaller than that value. Allowing noise filling a wider range of possible delay times degrades the stack of the faint signal that we report.

The PKiKP waveforms are prominent but complicated and different, partly because of the differences in source depth and slap-down phase timing () and possibly also because of complications from downdip propagation along the subducting Aleutian slab in conjunction with slightly different source locations. We tried numerous ways to deconvolve or cross-convolve to equalize the source time functions, without success. These explosions have simple waveforms for many global direct mantle and core-reflected P arrivals (, ). Our analysis assumes that the slab effects are localized close to the PKiKP raypath or that such effects average out over the wider swaths of ICS raypaths.

Figure S2 demonstrates that other highly scattered and prolonged phases do not show an organized pattern as clearly as ICS does, especially for the slownesses in which the energy is greatest. PKKP raypaths, which reflect from the far side of the core-mantle boundary, do at times touch the IC, but most energy goes around, as may be seen by the minimum in energy near vertical incidence. P′P′, which arrives noticeably off-azimuth, has the same properties. Both wave trains are shorter than the 225-s ICS waves, especially P′P′, which peaks for only about 50 s, and so benefit less from averaging than ICS.

The predicted pattern of Amchitka for 1° superrotation is shown in Figs. 1 and 2. It has the opposite polarity, blue near 0.02, which is also where the strongest energy is seen. The amplitude of the time shift seen is roughly ±0.2 s, while the prediction for 1° of rotation is ±1.8 s, so we estimate 0.1° of subrotation, which is an average of 0.05° per year. As a comparison, the same estimation for NZ is superrotation of about 0.1° per year (see Materials and Methods) (, ). This measurement is likely to be an underestimate, as we expect some damping of the signal by our processing and overprinting from unchanging waves traveling other paths in Earth, which would not be expected to show a time shift. A similar period of subrotation is also inferred by another seismological study of core-transmitted waves from earthquake doublets, which generally interpreted mostly superrotation, and at much greater rates ().

The Amchitka time shift measurements are not as clear as the shifts that we found for the NZ pair because, for the reasons given above, those explosions were much more similar repetitions. We made dozens of well-resolved time shift measurements in the NZ ICS waves, leading to a ±10% precision in rotation rate estimate. We could also evaluate the pole of rotation, as well as confirm that other phases interacting less with the IC showed very little change (). However, the polarity and rough amplitude of the movement are still apparent for the Amchitka explosions and indicate reversed motion at a comparable speed.

DISCUSSION

Our observations are consistent with a relatively simple model based on the mutual gravitational attraction between the lateral variations in the density in Earth’s mantle and IC and topography of the core-mantle and IC boundaries (, –). With reasonable parameters, this model can explain the clear 6-year oscillation in the length of day (LOD), although this resonance has alternative explanations, i.e., fluid OC motions (). The model predicts 0.5° peak-to-peak oscillation in the case that the OC does not differentially move with the IC, as is most likely (, ). Conservation of angular momentum would result in smaller IC oscillation if some of the OC participates in the oscillation.

This gravitational coupling model makes the specific prediction for IC rotation shown in Fig. 3 based on LOD variation observed back to the 1960s. As the LOD variation is proportional to the IC differential rotation rate, the integrated LOD variation is hypothesized to match the differential rotation angle (∆ΩIC; see Materials and Methods). During the interval from the first to the second Amchitka explosion, the integrated LOD shortened by 0.01 ms, corresponding to subrotation of 0.01° per year. From the first to the second NZ explosion, the integrated LOD lengthened by 0.12 ms, corresponding to superrotation at 0.08° per year. Our observations of IC motion accordingly show that it spun down by at least 0.05° per year in the first interval and spun up by 0.1° per year in the second. The IC motion’s amplitude and phase match the gravitational coupling model to within uncertainties.

Fig. 3.

Comparison with LOD data.

(A) Filtered LOD variation [gray line, with 6-month running average, following Holme and De Viron ()] and the predicted IC differential rotation angle (∆ΩIC, black line), computed from the integrated LOD data, assuming no coupling of IC differential rotation with the OC (, ). The blue and red lines delimit the interval between the twin explosive events pairs 1969–1971 and 1971–1974, respectively. The black dashed lines show the IC differential rotation angle trend between dates of the pairs of explosions, assuming that the LOD variance is only due to IC rotation. The blue and red dashed lines are the observed IC differential rotation angle trends. (B) Expanded time range between 1968 and 1977 with illustration of sub- and superrotation.

Alternating slow sub- and superrotation contradicts interpretations of fairly steady, progressive superrotation (, ), as well as rotation at much higher rates (, , ). Other studies acknowledge temporal changes in seismic waves but attribute them to nonrotational processes such as locally raising or lowering of the IC boundary topography (), material crossing the IC boundary either up or down (), or perhaps heterogeneity traveling in the flow of the OC ().

Our data do not preclude other, nonrotational processes. If the IC oscillates by only such a small angle, then such processes may well be necessary to explain most waveform change observations. The observation of different rates of changes for adjacent raypaths (, ) is more easily explicable if processes local to the IC raypaths are responsible than blaming IC rigid rotation.

While we cannot answer the question of what other processes are also driving changes in the IC seismograms, the separation of those changes from the changes due to oscillation will aid in their elucidation. While it is theoretically possible that the IC could precess as it oscillates if it is so inviscid that it could deform on the same time scale as it rotates, the range of possible models and their plausibility is sharply restricted.

Last, an IC oscillating about an equilibrium position would retain a consistent orientation to overlying mantle structures. This stability would permit models in which the lateral variation in IC structure, including the existence of an inclined boundary of the IC’s hemispherical structure (), might evolve guided by the influence on OC flow of strong thermal variations at the base of the mantle (–).

MATERIALS AND METHODS

Estimation of IC motion

The differential time shifts from the beamformed seismograms of the pair of nuclear tests recorded at LASA can be used to infer the IC’s rotation (–). Wang and Vidale () applied an improved backprojection approach to locating the strong scatterers within the IC based on the slownesses and lapse times of the backscattered energy in ICS waves. The time shifts of the ICS waves from the same scatterers for the twin nuclear tests at different times are compared with the predicted time shifts assuming 1° of IC eastward rotation. The linearity with a positive slope between the observed time shifts and predicted ones proves the IC superrotation, and the slope of the linear regression illustrates the amount of the rotation. On the basis of the above method, we measure the IC superrotation rate as 0.1° per year from 1971 to 1974 ().

In this study, the beams are noisier for the Amchitka explosions in 1969 and 1971. Instead of locating each scatterer from the beamformed ICS waves, we only use the average time shifts within the 225-s ICS waves, following the method in the work of Vidale () and Vidale et al. (). We next assume that all the scatterers are located close to the IC boundary and compute the predicted time shifts for 1° superrotation, using the source location in Amchitka (see Fig. 2 and fig. S3). The maximum and minimum values of the predicted and observed time shifts are used to compute the rotation amount. The predicted peak time shifts are ±1.8 s, compared to the observed ±0.2 s, with a reversed sign (Fig. 2, D to F). Thus, we estimate 0.1° subrotation. Wang and Vidale () found that if the scatterers are deeper than modeled, then the inferred rotation rate will be greater. On the basis of the observed distribution of IC scatterers (, ), we also consider an IC scatterer depth of 300 km as a lower boundary. The predicted time shifts are then ±1.25 s (see fig. S3), so the subrotation of the IC would be about 0.16°. The deep scatterers are not estimated to contribute at the larger slownesses close to grazing the IC, which correspond to some of the larger predicted time shifts (see fig. S3). This result indicates that the depths of IC scatterers lie below the IC boundary, and the IC’s subrotation rate is more consistent with the prediction from the LOD variation.

Connection between IC oscillation and LOD variation

The angular momentum of Earth is conservative and may be expressed aswhere IIC,OC,M and ωIC,OC,M are the moment inertia and angular velocity of the IC, OC, and mantle (including the crust) of Earth, respectively, and C is a constant. The mantle and IC are assumed to be rigid, so no viscous relaxation of the IC will be considered. The angular momentum of the system includes the OC, as well as the IC and mantle. However, the OC is sufficiently inviscid to have negligible viscous drag from small 6-year oscillations of the mantle and IC. Some electromagnetic coupling of the IC and OC is possible (), but following Ding and Chao (), we present the case of no coupling and note that the amplitude of the IC differential rotation would be reduced if much of the OC was entrained in the oscillation; therefore, the angular velocity change of the mantle can be rewritten as

Considering that ωM is almost constant, the angular velocity change ΔωM is proportional to the LOD variation (ΔLOD) ()

Combining Eqs. 2 and 3, we see that the IC differential rotation rate (ΔωIC) is proportional to the LOD variation (ΔLOD)

To obtain the IC differential rotation angle (∆ΩIC), we integrate the differential rotation rate (ΔωIC)

In summary, we obtain the IC differential rotation angle that is proportional to the integrated LOD variation.