## Abstract

Temperature gradients in a low-shear-velocity province in the lowermost mantle (*D″* region) beneath the central Pacific Ocean were inferred from the observation of a rapid *S*-wave velocity increase overlying a rapid decrease. These paired seismic discontinuities are attributed to a phase change from perovskite to post-perovskite and then back to perovskite as the temperature increases with depth. Iron enrichment could explain the occurrence of post-perovskite several hundred kilometers above the core-mantle boundary in this warm, chemically distinct province. The double phase-boundary crossing directly constrains the lowermost mantle temperature gradients. Assuming a standard but unconstrained choice of thermal conductivity, the regional core-mantle boundary heat flux (∼85 ± 25 milliwatts per square meter), comparable to the average at Earth's surface, was estimated, along with a lower bound on global core-mantle boundary heat flow in the range of 13 ± 4 terawatts. Mapped velocity-contrast variations indicate that the lens of post-perovskite minerals thins and vanishes over 1000 kilometers laterally toward the margin of the chemical distinct region as a result of a ∼500-kelvin temperature increase.

Heat transfer across Earth's core-mantle boundary (CMB) plays a central role in powering the core's magnetic-field–generating geodynamo and the configuration of mantle convection. Primary constraints on deep Earth temperatures are provided by the interpretation of seismologically detected velocity discontinuities as phase transitions, in which the associated pressures and temperatures can be experimentally and theoretically determined. Phase transitions that account for seismic discontinuities near 410- and 660-km depths and at the inner core–outer core boundary near 5149 km depth provide the few absolute temperature constraints available for the deep interior. There are large (±500 K) uncertainties in extrapolations along adiabats to the CMB at 2891 km depth (*1*). Mineral physics experiments have recently shown that a phase transition occurs in the primary mineral of the lower mantle, (Mg,Fe)SiO_{3} perovskite (Pv), yielding a post-perovskite (pPv) polymorph (*2*, *3*) at pressure and temperature (*P*-*T*) conditions close to those expected for a seismic velocity discontinuity observed several hundred kilometers above the CMB (*4*). Associating the phase transition with the discontinuity provides direct constraint on absolute temperature in the lowermost mantle, indicating ∼2500 K at 2700 km depth (125 GPa) (*2*), but there must be substantial lateral temperature variations.

The Pv-pPv phase transition is expected to produce a 2 to 4% *S*-wave velocity increase, little change in *P*-wave velocity, and an ∼1 to 2% increase in density, which are generally compatible with seismic observations (*2*, *3*, *5*–*7*). Experimental determinations of the Clapeyron *P*-*T* slope Γ for the Pv-pPv phase transition range from 5 to 13 MPa/K (*8*, *9*), and theoretical estimates fall within this range (*3*, *10*). The phase transition pressure may decrease with increasing bulk FeO content, although the magnitude of this effect is debated (*8*–*13*), and there is uncertainty in the partitioning coefficients for Fe, which might be influenced by a high- to low-spin transition in Fe^{2+} (*14*). The presence of Al_{2}O_{3} appears to have an opposite pressure effect to that of FeO and may broaden the two-phase region, affecting the seismic reflectivity (*15*–*18*).

The large positive Γ of the Pv-pPv phase change could potentially allow a rapid increase in temperature at greater depths to reverse the transformation. A strong temperature increase is expected in a thermal boundary layer at the base of the mantle (*1*), and the possibility of a second intersection of the geotherm with the phase boundary has been proposed (*19*). A shear velocity increase overlying a decrease may thus exist in the *D″* region; determining the discontinuity's depths would give two absolute temperature estimates and a temperature gradient based on the phase change *P*-*T* behavior. Initial investigations relevant to this question (*20*, *21*) are subject to other interpretations (*22*). A velocity reduction is intrinsically more difficult to detect than a velocity increase (*23*) and requires wave-form stacking of many signals. We conducted a detailed *S*-wave stacking-and-modeling analysis, resolving velocity structures that were consistent with local geotherms having a double intersection with the Pv-pPv phase boundary under the central Pacific. We obtained robust estimates of thermal gradients near the CMB.

**Seismic data analysis.** We analyzed 736 transverse horizontal-component *S*-wave seismograms from 46 intermediate- and deep-focus Tonga-Fiji earthquakes recorded at broadband seismic stations in California (Fig. 1). We added ∼300 seismograms to an earlier data set that sampled the lowermost mantle beneath the central Pacific Ocean (*24*, *25*) and used one-dimensional (1D) and 2.5D (axisymmetric spherical) modeling procedures to resolve the regional structure. Our *S* waves traversed a large low-shear-velocity province (LLSVP) imaged by global seismic tomography (*26*) (Fig. 1). Tomography models indicated that the CMB reflected phase (*ScS*) paths sample within the northern margin of the LLSVP, southeast of the Hawaiian hot spot (figs. S1 and S2).

The data have variable waveforms, with amplitude and timing fluctuations over small spatial apertures (figs. S3 to S5), which is consistent with prior work that proposed a laterally varying *D″* discontinuity at an average height near 230 km above the CMB in this region (*4*, *27*). Strong lateral gradients in seismic velocity structure are observed in nearby regions of the Pacific LLSVP (*28*), which suggest an abrupt northern margin similar to its southern margin (*29*) and the African LLSVP margins (*30*). In contrast to the southwest-to-northeast (SW-to-NE) regional trend from lower to higher *S*-wave velocities found in tomographic models, the *ScS* anomalies involve increasingly delayed arrivals from SW to NE over the small-scale region sampled by the data, along with a gradient in anisotropy (*31*). Our large data set agrees with previous analyses of various data attributes across the study area, which established that the primary variations are SW-to-NE, parallel to the Tonga-to-California ray paths.

Source complexity variations were accounted for by deconvolving the signals by averaged source wavelets for each event (fig. S3), which allowed data from different events to be combined. The deconvolved signals had bandwidth from 0.01 to 0.3 Hz. The data were subdivided into three parallel bins, separating *ScS* CMB reflection points at varying distances from California (Fig. 1). Double-array stacking of all station-receiver combinations sampling each bin characterized the reflectivity relative to height above the CMB. We normalized and aligned the seismograms so that *ScS* had unity reflectivity at the CMB (figs. S6 to S9). Each bin had several positive and negative peaks in reflectivity (Fig. 2 and fig. S7). Reflectivity more than 400 km above the CMB could not be confidently constrained as a result of contamination by the direct *S* phase, which turns at shallower depths.

The data stacks were first modeled with synthetics for localized 1D-layered structures to characterize the basic attributes of the quasi-2D structure. For each bin, a four-layer model, motivated by previous studies of the region (*24*, *25*, *27*), was developed to match salient attributes of the data stacks. The synthetics were processed identically to the data, accounting for bandwidth and deconvolution side-lobe effects (figs. S3 to S5).

The resulting models (Fig. 2) indicated a positive velocity increase (reflector U) at depths increasing from 2546 km (bin 1) to 2655 km (bin 2) to 2730 km (bin 3), with shear velocity increases of 1.1, 0.6, and 0.7%, respectively. This feature is similar to the 2 to 3% shear velocity increase detected below circum-Pacific regions commonly associated with the Pv-pPv transition. A small (–0.5 to –0.7%) velocity discontinuity (reflector A) was modeled ∼135 to 70 km above reflector U in bins 2 and 3. Although the data did not resolve reflector A in bin 1 because of *S*-wave interference, a similar structure was assumed to be present with a very small velocity decrease. This feature has no counterpart in models for circum-Pacific regions, and we associated it with the top of the LLSVP in each bin. Below the shear velocity increase, two abrupt velocity decreases were modeled: The shallower decrease (reflector L) was a weak –0.2% jump at 2805 km depth in bin 1, a –0.6% jump at 2800 km depth in bin 2, and a –0.9% jump at 2788 km depth in bin 3; and the deeper decrease (reflector B) was a –4.0% drop at 2874 km in bin 1, a –1.1% drop at 2859 km depth in bin 2, and a –0.6% drop at 2847 km in bin 3. The reflector L velocity decrease was the most unexpected feature in the structures. It is generally consistent with the hypothesis of a second intersection of the geotherm with the pPv phase boundary, involving the conversion of pPv back to Pv, effectively forming a lens of pPv material in *D″* above the CMB (Fig. 2). The discontinuity closest to the CMB may be associated with the top of a weak ultra-low-velocity zone (ULVZ), a thin low-velocity layer extensively observed just above the CMB (*4*, *25*, *31*). The predicted arrivals for these discontinuities are small (fig. S10); stacking is required to observe them clearly and to suppress scatter and the effects of localized lateral variations.

Our modeling includes the use of *ScS* and *S* separately as reference phases, with simultaneous modeling of both ensuring reliable velocity structures (fig. S8). The average *S*-wave velocity over the lowermost 350 km of the mantle decreases from bin 1 to bin 3 by anamount consistent with the lateral 3- to 4-s increase in *ScS* travel time delays (*31*), and this aspect of the data stacks is well modeled by our structures (fig. S8).

Uncertainties in the size of velocity contrasts and depths are ∼20% and ±5 km, respectively, for 1D models, but true uncertainties may be larger in the presence of lateral variations. Modeling with 2.5D finite-difference calculations confirmed that the depth estimates were realistic, but the velocity contrasts tended to be 50 to 100% higher in laterally varying models (figs. S11 and 12) as a result of the lateral termination of the reflecting interfaces. Unresolved small-scale heterogeneity can contribute to amplitude variability, affecting the inferred velocity contrasts (fig. S13).

The *S*-wave ray paths grazed the lowermost mantle (Fig. 1B). The three data bins spanned only ∼4° epicentral distance (240 km at the CMB) along the propagation direction, but the length scale over which the wave fields interact with the structure was several times larger, resulting from the spatial extent of the effective Fresnel zone (300 to 500 km) for the 3- to 4-s dominant periods of the data. 2.5D axisymmetric finite-difference modeling demonstrated that the spatial extent controlling signals in each bin was ∼5° epicentral distance; the entire data set sampled structure over about a 15°-wide region in epicentral distance (fig. S11). Energy in the data stacks between reflectors U and L in bins 1 and 2 not reproduced by 1D modeling was due to lateral sampling of adjacent bins. Although the precise structure was not definitively constrained, the basic layering and associated reflectivity are well characterized by the models in Fig. 2, particularly the discontinuity depths, which play a critical role in the following discussion because they provide the connection to pressure.

**Thermal models.** We focused on the depths of the paired velocity increase and decrease in each bin (reflectors U and L), which we interpreted as double crossings of the Pv-pPv phase boundary. We considered thermal models with mid-mantle adiabats overlying a conductive thermal boundary layer of variable thickness parameterized by an error-function decrease from a peak temperature at the CMB. This choice of model parameterization is discussed in the supporting online material (SOM) text S1. The phase-boundary position in *P*-*T* space is not precisely known because of uncertainties in the effects of mixed phases, Fe and Al partitioning in the deep mantle minerals, and experimental error. The position of the phase boundary was parameterized by Γ [ranging from 8 to 13 MPa/K, which spans the range of recent estimates, including those for pyrolitic composition (*8*, *9*, *32*)] and Δ*T*_{CMB,pPv} [the (positive) temperature difference between the CMB temperature and the pPv phase-boundary temperature at the CMB, ranging from 100 to 400 K].

The seismic models provide three samples of upper and lower intersections with the phase boundary in a laterally varying thermal regime (as required to produce regional topography on the phase boundary). Reflector U varies systematically relative to reflector L as the geotherm changes spatially. The seismic discontinuity depths yield reasonable geotherm fits for the range of Γ and Δ*T*_{CMB,pPv} considered, enhancing the plausibility of the double-intersection models. Example geotherms for two choices of Γ and two choices of Δ*T*_{CMB,pPv} are shown in Fig. 3 (see also fig. S14). As Γ increases, the radial thermal gradients decrease and, as Δ*T*_{CMB,pPv} increases, the gradients increase. The corresponding boundary-layer thickness estimates needed to match seismic observations range from 150 to 400 km (fig. S15).

Estimating heat flux for a given thermal structure requires knowledge of lowermost mantle thermal conductivity *K*. We adopted the most commonly used estimate, *K* =10W/(m·K) (*33*), noting that this parameter has not been directly constrained and both higher and lower values have been advocated (uncertainty in *K* is discussed in SOM text S2). The suite of models that fit the data yields average heat flux estimates shown in Fig. 4 and depends directly on the assumed *K*. The lateral variation in temperature at the depth of reflector U in bin 2 is shown for each set of geotherms. The lateral variation in temperature between the bins is controlled by Δ*T*_{CMB,pPv}. Fitting the data is difficult for Δ*T*_{CMB,pPv} values less than 100 K, whereas unrealistically high heat fluxes are predicted for values more than 300 K. The heat flux estimates, ranging from 60 to 110 mW/m^{2}, are on the order of Earth'saverage surface heat flux (∼86 mW/m^{2}) (*34*). Reducing the uncertainty of Γ and *K* will directly narrow the range of heat flux estimates.

The core's surface area is ∼30% of Earth's surface area, so our result suggests that heat flow through the CMB may be a corresponding fraction of the 44 TW at the surface (*34*). Assuming that the study region is relatively warm, as indicated by the low *S*-wave velocities, we estimated a lower bound of 13 ± 4 TW for CMB heat flow for our choice of *K*. Substantially higher heat flux would be estimated if radiative *K* is high, whereas lower heat flux could result from acutely heterogeneous grain-boundary structures and interstitial impurities that impede phonon conductivity (SOM text S2). Modeling the discontinuity depths for individual bins gives very consistent results (fig. S16). Exploring the effects of the ±5-km uncertainties in the discontinuity depths indicates that the results are most sensitive to the reflector L depths.

The sensitivity to reflector L depths raises the issue of distortion of the phase boundary by dynamical instability. The layer of pPv minerals is expected to be denser (∼1 to 2%) than the underlying Pv, and the latter is likely to have strong viscosity reduction resulting from the rapid increase in temperature in the thermal boundary layer. The lower edge of the pPv lens will thus tend to sink, heat up, and convert to Pv, establishing a dynamical equilibrium. This process is difficult to simulate, lacking constraints on the overall size and geometry of the pPv lens and on boundary-layer viscosity structure.

Another implication of these thermal models is that there must be horizontal conduction of substantial heat (on the order of 10 mW/m^{2}), given the temperature variations between bin models and the horizontal scale of ∼1000 km sampled by the seismic waves. This lateral thermal gradient suggests dynamical flow with relatively hot upwelling material toward the northeast (near bin 3).

**Chemical piles.** We consider the foregoing analysis of thermal gradients under the central Pacific in a global context. The study area is near the northern edge of the LLSVP under the Pacific Ocean, far from recent subduction zones (*35*), and is generally considered to be warm (relative to cooler circum-Pacific locations) because of the low shear velocities. Lacking a chemical effect, the Pv-pPv phase boundary should be deeper or absent in warm mantle (*36*). However, this region appears to be anomalously dense, possibly because of the presence of excess FeO (*37*, *38*), favoring the notion of chemically distinct material swept into a large pile by larger-scale subduction-driven flow in the surrounding mantle. There may be alternate explanations for the density excess, with the accumulation of mid-ocean ridge basalt–enriched mantle being one possibility (*39*), but this suggestion still invokes a dense, chemically distinct pile of material. Numerical modeling indicates that a dense pile is likely to have an increase in temperature toward its edges as a result of boundary layer separation (*40*).

Circum-Pacific regions have *D″* shear velocities that are 3 to 5% faster than under the central Pacific and stronger shear velocity discontinuities (*4*). It is widely inferred that the contrast in velocities is partly due to temperatures being 700 to 1200 K lower in the circum-Pacific because of the presence of cool downwellings and partly due to the compositional differences in the dense central Pacific chemical pile. A 10% molar increase of FeSiO_{3} in pPv can reduce shear velocities by 2.6% (*13*), which is comparable to the effect of a 1000-K temperature increase at CMB pressure.

Comparing the structure of the double crossing of the Pv-pPv phase boundary in the central Pacific with structure in circum-Pacific regions permits a specific assessment of the effect of chemical variations on the phase boundary. The only circum-Pacific region with evidence for a double crossing is under the eastern Cocos Plate (*41*): A 1.7% shear velocity increase at 2704 km at the upper Cocos discontinuity (UC) is accompanied by a 1.7% decrease at 2821 km at the lower Cocos discontinuity (LC). The average shear velocity between these depths is 7.55 km/s, ∼5% faster than in the pPv lens in the central Pacific. If we assume that the Cocos Plate and central Pacific environments have the same Γ of 11.5 MPa/K (*8*) and a 700-K contrast in temperature, we can find a geotherm compatible with the Cocos region velocity discontinuity depths, UC and LC (Fig. 5).

This scenario requires a 6-GPa pressure shift of the phase boundary resulting from an effect such as higher FeO content of the LLSVP, which is consistent with some predictions (*9*, *11*). A 10–mole percent increase in Fe is expected to give a 5- to 10-GPa reduction of the phase-boundary pressure (*9*). This additional Fe should broaden the phase transition loop, giving a velocity transition rather than a discontinuity and weakening the reflectivity (fig. S17). The velocity discontinuity also diminishes because of the greater reduction of shear velocity with added Fe for pPv (–0.26%/mol % Fe) than for Pv (–0.2%/mol % Fe) (*13*). For a 10% Fe–enriched region, the shear velocity contrast should be ∼0.6% weaker than in a normal region, depending on Fe partitioning coefficients with other phases. The combined effects of added Fe and lateral temperature variations relative to circum-Pacific regions can thus account for the 5% shift of absolute shear velocity and the reduced velocity discontinuity. If the UC phase transition occurs in relatively Fe-depleted MgSiO_{3} [resulting from the low-spin transition of Fe favoring partitioning into (Mg,Fe)O] at 2500 to 2600 K (*2*, *9*), the CMB temperature for this scenario is ∼4100 K.

**Conclusions**. Seismic reflectivity under the central Pacific is consistent with a laterally vanishing lens of pPv material within and near the margin of a chemically distinct pile (Fig. 6). The pile is probably undergoing internal flow, with shear coupling resulting in hot upwelling on the pile margins, giving a lateral increase in temperature that causes the pPv lens to laterally revert to Pv. The internal transition from horizontal to vertical flow near the edge of the pile can account for an observed gradient in shear wave anisotropy (*31*). The velocity decrease observed above the U horizon may correspond to the top of an Fe-rich pile. The deeper ULVZ feature varies on the edge of the pile and is thicker and not as extreme in velocity reduction as observed in some regions (*42*). The *S*-velocity reductions are comparable to the *P*-velocity reductions (*27*) in this ULVZ. This finding indicates that little or no partial melt accumulated under this edge of the pile, whereas strong ULVZs with partial melting are found on other margins (*43*).

Thermal modeling gives direct determinations of thermal gradients for the CMB region. Regional heat flux into the base of the LLSVP is 85 ± 25 mW/m^{2}, close to average surface heat flux, for *K* = 10 W/(m·K). Global extrapolation suggests a lower bound on CMB heat flow of 13 ± 4 TW, subject to large uncertainty in *K*. These relatively high values favor the sequestration of heat-producing radiogenic elements in the core and a relatively young age for the inner core.

**Supporting Online Material**

www.sciencemag.org/cgi/content/full/314/5803/1272/DC1

SOM Text S1 and S2

Figs. S1 to S17

References