<figure><img src="image_1.png" /></figure>

HAL
open science

# Partial melt in mesoscale upper mantle upwellings beneath ocean basins

Isabelle Panet, Marianne Greff-Lefftz, Barbara Romanowicz

► To cite this version:

Isabelle Panet, Marianne Greff-Lefftz, Barbara Romanowicz. Partial melt in mesoscale upper mantle upwellings beneath ocean basins. Earth and Planetary Science Letters, 2024, 639, 10.1016/j.epsl.2024.118763. insu-05081568

HAL Id: insu-05081568

https://insu.hal.science/insu-05081568v1

Submitted on 23 May 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# Partial melt in mesoscale upper mantle upwellings beneath ocean basins

Isabelle Panet<sup>a, b, *</sup>, Marianne Greff-Lefftz<sup>a</sup>, Barbara Romanowicz<sup>a, c, d</sup>

<figure><img src="image_4.png" /></figure>

## ARTICLE INFO

**Keywords:**

Mantle mass structure

Gravity gradients

Seafloor topography

Seismic tomography

## ABSTRACT

How tectonic plate motions are coupled with mantle flows remains an open question. Quasi-periodic 2000 km wavelength undulations aligned with absolute plate motion in the Pacific and Indian Ocean basins, observed in gravity and seafloor topography and coinciding with seismically imaged low shear velocity fingers in the upper mantle, suggest the presence of meso-scale convection below the lithosphere. However, the correspondence of sub-lithospheric mantle mass excess, seafloor lows and slow upper mantle seismic velocities cannot be explained by temperature variations alone. Here we introduce a simplified system of bi-dimensional convective cells of width ~1000 km, extending from the base of the lithosphere through the extended mantle transition zone (down to 1000 km depth), at least partly driven from below. From mass balance considerations in a viscous Earth, we show that the density excess required in the hot upwelling limbs may reflect the formation of stable dense lenses of dehydration-induced partial melt atop the 410 km discontinuity, and upward entrainment of a small fraction of quasi-buoyant partially molten and recrystallized material across the upper mantle. Our model provides an explanation for the thin low shear velocity layer detected intermittently above the 410 km discontinuity in some parts of ocean basins away from subducted slabs, and supports the presence of water in the transition zone.

## 1. Introduction

Mapping the seismic velocity and density structure of the earth's mantle at high spatial resolution is a crucial step towards constraining the nature of mantle heterogeneity and the patterns of its convective flow at different scales. Global 3D seismic tomography has revealed striking differences between the shallow and deep mantle, as the patterns of seismic structure appear to be de-correlated above and below ~1000 km depth (Rudolph et al., 2015). Above 1000 km depth, enigmatic quasi-periodic, ~1800–2000 km wavelength, horizontally elongated slow shear velocity "fingers" have been found below all ocean basins (Lekic and Romanowicz, 2011; French et al., 2013), aligned with present-day absolute plate-motion (APM). The dynamic processes at the origin of these features remain unresolved, as are their potential interactions with other forms of mantle convection, such as deep mantle plumes (French and Romanowicz, 2015; Lei et al., 2020).

Owing to their sensitivity to the internal density structure, lateral variations in the gravity field and surface dynamic topography can be used to interpret the seismic velocity variations in terms of thermal and/or compositional heterogeneity, and obtain insights on the

corresponding dynamics. With a globally high spatial resolution, satellite gravity and altimetry data also provide unique information to better delineate the lateral shape and size of the convective structures (Panet et al., 2014). In previous work (Panet et al., 2022), we carried out an analysis of the second-order horizontal derivatives of the seafloor topography and of the gravity potential over the Pacific and Northern Indian ocean basins, and detected periodic undulations in these observables (Figure 1bc), at the same 2000-km wavelength and along the same orientation (i.e. present-day APM) as found for the slow seismic velocity fingers (Fig. 1a). Our results confirmed and extended those of Hayn et al. (2012), showing the presence of energy in the geoid at these wavelengths and orientations in the Central Pacific. Furthermore, we showed that, at least in some areas, these gravity and topography signals require alignments of mass excess in the sub-lithospheric convective mantle beneath seafloor lows, coinciding geographically with the slow seismic fingers (Fig. 1). This is opposite to what would be expected in the case of purely thermal convection, where slow seismic anomalies would correspond to hot and light upwellings beneath seafloor highs, and thus suggests the presence of dense compositional anomalies below the lithosphere in areas of slow seismic velocities, which is puzzling, but

<footnote><sup>a</sup> Université Paris Cité, Institut de physique du globe de Paris, CNRS, IGN, Paris F-75005, France</footnote>

<footnote><sup>b</sup> ENSG Géomatique, IGN, Marne-la-Vallée F-77455, France</footnote>

<footnote><sup>c</sup> University of California, Berkeley, CA 94720, USA</footnote>

<footnote><sup>d</sup> Collège de France, Paris, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail address: panet@ipgp.fr (I. Panet).</footnote>
confirms and generalizes an earlier inference made on a 2D profile between Fiji-Tonga and Hawaii (Katzman et al., 1998).

Here, we seek to quantify and interpret the corresponding thermal and compositional density anomalies in a dynamically plausible fashion, by building a regional-scale mantle mass model able to explain the three sets of periodic, 2000 km wavelength APM-oriented signals detected in the seismic, gravity and topography observables. In our approach, the geometry of the 3D sources at depth is initially constrained using the spatial properties of the seismic and the medium-scale geodetic signals, taking advantage of their complementary sensitivities. The densities of the sources are estimated in a second step from the geodetic data only, under an assumption of long-term vertical stability of this system of sources. We show that the geodetic and seismic observations can be reconciled by modifying a bi-dimensional thermal convection model to take into account the effect of dehydration melting in the upwellings, due to passage across the 410 km depth phase transition in the olivine system, in a moderately hydrous mantle.

## 2. The investigated signals and their relation to the mass sources

### 2.1. Data: periodic, APM-oriented anomalies in the Central Pacific

Fig. 1 shows the investigated signals over and below the Pacific ocean basin (see Panet et al., 2022, for a detailed description). To highlight the medium-scale APM-oriented elongated anomalies, we first filter the geoid and the seafloor topography at different spatial scales using wavelets. Then, for each spatial scale we calculate the second-order horizontal gradients of the wavelet-filtered geopotential and seafloor topography, expressed in local spherical frames. These gradients emphasize signals orthogonal to the East-West or North-South differentiation directions of the horizontal frame vectors; thus, we rotate the frame around the radial axis to investigate APM-oriented structures (Panet, 2018). We finally stack the gradient maps over a range of frame rotation angles in order to increase the signal-to-noise ratio.

These geodetic signals are based on analyses of GRACE and GOCE satellite gravity data, and of the Smith & Sandwell V16.1 bathymetric grid corrected for the isostatic contribution of the sediment load. This bathymetric grid is constrained from ship soundings, densely distributed

<figure><img src="image_5.png" /><figcaption>**Fig. 1.** APM-oriented, medium-scale periodic signals in the upper mantle shear velocities, seafloor slope gradients, and GRACE gravity gradients for the Pacific plate. Top: SEMum2 (French et al., 2013) slow velocity anomalies (dVs/Vs or dlnVs) at 260 km depth (a) with location of the cross-section shown in (d,e); GRACE gravity gradients (b) and seafloor slope gradients (c) along the APM direction (1100 km scale, 00 second-order gradients of the geopotential or the seafloor topography with a 10–40° clockwise rotation of the spherical frame around the radial axis to emphasize APM-oriented signals, see inset). Seafloor lows (resp. GRACE gravity lows) correspond to the positive 00 second-order gradients of the seafloor topography (resp. of the geopotential). For periodic signals, the spatial scale corresponds to the half-wavelength. Orange lines: plate boundaries (Bird, 2003); black arrows: present-day APM direction (Seton et al., 2012); green dots with red circles: hotspots (Steinberger, 2000). Bottom panel (d): Variations along the profile indicated in the top panels, for all 4 observables (shear velocity anomalies at 260 km depth). Signs of dVs/Vs and Bouguer gravity gradients have been reversed. Green rectangles: lower than average shear velocity. Pink rectangles: regions of seafloor lows. Bottom panel (e): depth cross-section in model SEMUCB-WM1 (French and Romanowicz, 2015), showing dVs/Vs with respect to the global mean, along the same profile. Redrawn from Panet et al. (2022).</figcaption></figure>
at the investigated spatial scales, validating its use as an independent dataset (SI Section A in Panet et al., 2022). Our analysis differs from those of Hoggard et al. (2017) and Holdt et al. (2022), who, in addition, remove a lithospheric plate cooling model and correct for crustal contributions by excising areas of anomalous crust and/or applying isostatic crustal corrections. Note that we have verified the robustness of our results to these corrections and the consistency of the discussed signals with the model of dynamic topography of Hoggard et al. (2017). The amplitudes of the geodetic signals over the younger and the older sides of the Pacific plate are summarized in Suppl. Table A1.

## 2.2. From the geodetic data to the present-day mass distribution in a dynamic Earth

A mass anomaly in the convective mantle induces a viscous flow, hence dynamic deflections of the internal and surface interfaces, and may also impact the depth of the mineral phase transitions by a predominantly thermal effect. These changes in the density structure generate an additional mass redistribution, in order to conserve mass in a compressible Earth. The detected gravity and surface dynamic topography signals result from the combined contributions of all these mass anomalies. This is why a widely used approach to relate observed geodetic signals to the corresponding present-day mass structure is to solve for the instantaneous flow driven by mantle mass anomalies in a self-gravitating Earth (Richards and Hager, 1984; Forte and Mitrovica, 2001), and we apply this approach (see SI Section B).

# 3. Step 1: constraints on the geometry of the sources

## 3.1. Principle

To investigate the thermal and compositional sources of the discussed signals, we develop a two-step approach. We first constrain the 3D geometry of the mass sources at depth from the spatial characteristics of the seismic and geodetic signals, interpreted in the light of a convection model able to fulfill the observation-based constraints on the sources (step 1). Then, the densities are estimated independently from the seismic velocities, under the dynamical assumption of a globally equilibrated system of sources (step 2). While the high spatial resolution of the seismic tomography model makes it possible to detect the slow velocity fingers and their correlation with the geodetic signals, we found that the shapes of the seismic anomalies are less well resolved laterally than those of the geodetic gradient signals. On the other hand, the seismic data have better depth sensitivity. Instead of directly converting the seismic velocity anomalies to density, we can combine the complementary sensitivities of the seismic and the geodetic data to depth and to lateral variations of the sources respectively, together with a simplified prior model, to define a common domain for sources of mass anomalies within the mantle.

## 3.2. Depth and lateral extent of the thermal sources

The depth range of the sources is primarily constrained by the seismic tomography model. The alignments of slow seismic velocities are best detected in the 200–350 km depth range (French et al., 2013), but they extend in depth from the asthenosphere through the transition zone (French and Romanowicz, 2015) and in some areas, down to ~1000 km depth (e.g. Fig. 1). We thus bound the radial extent of the sources from the top of the sub-lithospheric low velocity zone to the base of an extended transition zone (ETZ), around 1000 km depth. We note that the latter coincides with the maximum depth of sensitivity of the gravity gradient signals, which is commensurate with their ~1000 km scale (Panet, 2018).

Information on the lateral extent of the sources of excess density in the upper mantle can independently be derived from the characteristic scales of the observed gravity gradient and seafloor slope gradient

signals, which correspond to the horizontal length scale of the signal along the direction of the considered frame vector, e.g., $e_0$ for the $\theta\theta$ gradients in the local ($e_r, e_\theta, e_\varphi$) frame. To illustrate this, we perform synthetic sensitivity tests and investigate the characteristic scales of the geodetic signals generated by 100 km thick, APM-oriented elongated thermal or compositional mass anomalies of varying width and width-to-length aspect ratios between 1:3 and 1.10, located at different depths in a deformable earth, considering two reference viscosity profiles with or without a low viscosity channel (LVC) in the 60–220 km depth range (SI Section B, Suppl. Fig. B1; note that in all the tests of this manuscript, the LVC is chosen to globally coincide with the sub-lithospheric low velocity zone (LVZ)). We compare the cases of single sources (SI Section C) and of a periodic alignment of elongated sources (Fig. 2). We find that a 2000 km periodic spacing of sources, of widths close to the half-period of the alignment, generates a gravity gradient or a seafloor slope gradient signal of constant 1000 km lateral scale whatever the source depth (Fig. 2), consistent with the lateral pattern of the slow seismic fingers and the geodetic anomalies.

We conclude that internal mass loads able to explain the spatial properties of all observed signals require alignments of elongated sources in a depth range spanning from the top of the asthenosphere to the bottom of the ETZ, with a period of ~2000 km laterally. Above the transition zone, their width is constrained to ~800–1000 km (Fig. 2). These combined 3D geometrical constraints are reminiscent of Rayleigh-Benard convection with a ~1:1 cell aspect ratio, as would be expected in the case of convection in a fluid of constant viscosity heated from below (Bergé and Dubois, 1984). Even if the actual flow is probably three-dimensional, we use this bi-dimensional thermal model as a guidance, and thus consider a simplified model in which we schematically represent each upwelling or downwelling limb by a rectangular mass anomaly of width 1000-km extending from 170 to 1000 km depth, elongated horizontally in the APM direction, and carrying a constant thermal anomaly which takes opposite values for the upwellings and the downwellings: +/-$\Delta T$. This representation of the sources is very simplified, because the limited resolution obtained at depth from the data only allows us to estimate average densities over simple domains. Nevertheless, the above geometrical constraints bound the extent of the sources to first order, and the derived volume averages can be used in the future for comparisons with the predictions of more realistic convection models, averaged over the same domains. To account for uncertainties in the depth range of the sources, we vary the top and the bottom depths of the convection cells, from 60 to 170 km depth and from 800 to 1200 km depth, respectively.

However, this distribution of mass sources of thermal origin needs to be modified in order to account for the compositional variations required to jointly explain the seismic and geodetic observations (Section 1).

## 3.3. Asymmetry of the positive vs negative geodetic signals and compositional sources

A coherency analysis between the observed gravity gradient and seafloor slope gradient undulations provides further insights on the structure of the corresponding density anomalies. We estimated this coherency separately for the positive and the negative lobes of the signals (SI Section D). Over most of the Central Pacific, the coherency between the positive observed gravity and seafloor slope gradients ($C^+$) is higher than 0.5. In contrast, the coherency between the negative observed gravity and seafloor slope gradients ($C^-$) drops below 0.5 in a wide zone around Polynesia and on the older side of the Pacific plate (Fig. 3ab). Over the Indian ocean, we also find a higher coherency between the positive GRACE gravity gradient and seafloor slope gradients than between the negative ones over a large area (Suppl. Fig. D1). This asymmetry cannot be explained by a symmetric system of dynamic mantle sources, as would be the case for the purely thermal convection model introduced in Section 3.2.
<figure><img src="image_6.png" /><figcaption>**Fig. 2.** Characteristic scale of gravity gradient signals and seafloor slope gradient signals generated by a periodic alignment of APM-oriented elongated mass excess sources at increasing depths for different source widths (00 second-order gradients in a spherical frame aligned with the sources, see panel a). The characteristic scale is calculated above the mass excess sources (at the location of the red star in the map). Widths W of the sources: 300, 600, 800 and 1000 km, shown in different colors (see legend in panel d). The width-to-length aspect ratio of the sources is fixed to 1:6, their radial extent is between D-50 km and D + 50 km, where D is the depth. Panel a: location of the sources, for W = 800 km, in map view (top) and in a cross-section perpendicular to the APM direction (bottom). The period of the alignment is always 2000 km across the APM. Panel b (resp. c): case of gravity gradients computed in a dynamic Earth, for the viscosity profile 2 with a LVC (resp. profile 1 without a LVC), as described in Suppl. Fig. B1. Panel d (resp. e): case of seafloor slope gradients from surface dynamic topography, for viscosity profile 2 (resp. 1). Solid lines: scales for chemical mass sources; dots: scales for thermal mass sources crossing the 410 or the 660 km depth phase transitions. Note that these results do not change for aspect ratios of the sources 1:3 and above.</figcaption></figure>

This is informed by synthetic tests in which we consider a system of sources consisting of an alignment of horizontally elongated cells between the base of the ETZ and the base of the sub-lithospheric LVC (coinciding with the LVZ for simplicity). We focus on the mass anomalies representing thermal convection cells in the 170–1000 km depth range with a +/−75 K temperature anomaly (Section 3.2). We calculate the corresponding gravity and seafloor slope gradients signals in a dynamical modelling framework (Section 2.2), considering the viscosity profiles with or without a low viscosity channel (LVC) in the 60–180 km depth range (profiles 1 and 3 in Suppl. Fig. B1). As expected, we find that the coherency remains high and symmetrical, with C⁺ and C⁻ nearly equal to 1 for both viscosity profiles 1 and 3 (Figures not shown). Indeed, this system of sources is symmetrical, with a similar topographic response over the upwellings and the downwellings.

We then add 600 to 1000 km wide layers of density excess (respectively deficit) in the 60–160 km or 160–380 km depth range, within or above the upwelling limbs (respectively downwelling limbs) of this thermal convection model. These dense additional compositional sources are represented in a simplified manner by rectangular mass anomalies centered within the thermal mass anomalies of the upwelling limbs. A clear asymmetry is found in the presence of a sublithospheric LVC if we locate a mass excess in the LVC above the upwellings, with moderate

excess densities compared to those of the thermal rolls. The coherency between the positive lobes of the modelled gravity and seafloor slope gradients is then higher than that between the negative lobes when the width of the mass anomaly is smaller than the ~1000 km width of the convection cells (Fig. 3cd, Suppl. Figs. D2-D4). If, on the contrary, we consider a mass deficit above the downwellings, the respective coherency maps of the positive and the negative signals are reversed compared to Fig. 3cd, with a higher coherency between the negative modelled gravity and seafloor slope gradients than between the positive lobes. Thus, the observed asymmetry requires an additional source of mass excess on the upwelling side of the convective system, as opposed to a mass deficit on the downwelling side. The asymmetry increases for decreasing widths of the dense mass sources (Suppl. Figs. D2-D4), which, however, have to remain wider than ~800 km in order to explain the 1000 km scale of the signals (Section 3.2). We note that only a 100 km width difference compared to the 1000 km half-period is sufficient to break the symmetry (Suppl. Fig. D4, Figure 3 cd). Finally, mass excess sources located below the LVC within the thermal upwellings do not contribute to the observed asymmetry but do not prevent it either for a range of excess densities, as they simply reduce the thermal signal. However, as we describe below, they are necessary to jointly explain the topography and gravity gradient signals.
<figure><img src="image_7.png" /><figcaption>**Fig. 3.** Coherency analysis of the 1100-km scale, 0° modelled gravity gradients and seafloor slope gradients in the Pacific ocean (10° to 40° stacked frame rotation). Left (resp. right) panels: positive (resp. negative) coherency signals, see SI Section D. Coherency is computed over a 25° side window around each point. Panel a (resp. b): observed coherency. Panels c-d: predicted coherency for a system of 1000 km wide thermal convection cells with additional 900 km wide dense layers in the 60–160 km depth range, in the upwellings limbs of the convection system (see location in the inset map). Viscosity profile assumed as in Figure 4a with a sub-lithospheric LVC; temperature anomalies in the convection cells: ±75 K; density anomalies in the dense layer: 6 kg/m³. Panels e-f: same as panels c-d, for 800 km wide dense anomaly layers in the 60–160 km and 160–380 km depth ranges (density excesses, respectively: 6 kg/m³ and 12 kg/m³). For the considered viscosity profiles, high coherencies are obtained when the sources of the gravity and seafloor topography signals span the same depth range. This is not the case over the downwellings, in the presence of dense layers of width less than 1000 km in the LVC above the upwellings. In this case, side lobes of the large dynamic topography created by these dense layers strongly control the shape of the topographic signal above the downwellings. This degrades the coherency with the gravity signal, mostly produced by the deeper thermal convection cells.</figcaption></figure>
# 4. Step 2: estimation of the upper mantle mass excess

## 4.1. System of sources

Based on these observational constraints, we can now fix the geometry of the considered system of sources to first order, by modifying a bi-dimensional thermal convection model. We thus consider a system of sources comprising a series of convective rolls with a nearly 1:1 aspect ratio, with thermal anomalies ranging from +/- 75 to +/- 200 K, and dense fingers within the 1000 km wide upwelling limbs or above them, just below the lithosphere. As in Section 3.3, these dense additional

compositional sources are represented by rectangular mass anomalies centered within the thermal mass anomalies of the upwelling limbs. The compositional density anomaly is assumed to remain constant within each source. Based on the geometrical constraints described above, we assume that the lateral extent of these dense fingers is 10% smaller than that of the thermal upwellings, corresponding to a 900 km width of the rectangular sources. In depth, we distribute the density excess in the upper mantle according to structural boundaries derived from seismology. Thus we consider 3 regions, or “layers”, within which we seek the distribution of excess density (Fig. 4a): a top layer spanning most of the upper mantle seismic low-velocity zone (60–160 km, layer 1), a

<figure><img src="image_8.png" /></figure>

**Fig. 4.** Best-fitting total density excess in layers 1 and 2, for convective cells extending in depth from 170 km to 1000 km and carrying ±75 to ±200 K thermal anomalies. (a): Two viscosity profiles and the depth ranges of the 3 layer model are considered. The LVC extends down to 180 km in the tests shown in Figures 3, 4 and SI Section E, and down to 220 km in those shown in Figure 2 and SI Sections C, D, H, but these differences do not change our conclusions. (b): cross-section of the thermal (blue and red rectangles) and compositional (black rectangles in the upwellings only) sources considered in the modelling. (c): Cumulative misfit on the Bouguer gravity and seafloor slope gradients, and global vertical equilibrium (Section 4, rms from the variance shown in Eq. (1)) for the case without an LVC. Contours represent 2.5% isolines up to 30% of residual standard deviation (external contour). A smaller number of contours indicates a degradation of the fit to the data. Colors indicate thermal anomaly in the convection cells: ±75 K (dark blue lines), ±100 K (violet lines), ±150 K (pink lines) or ±200 K (light pink lines). The compositional density anomaly in layer 3 is fixed to 1.5 kg/m³. (d-g): Same as panel c, in the case with an LVC, for different thicknesses and viscosities of the LVC. An LVC is needed to avoid too large anomalies in layer 1, leading to a well constrained mass anomaly in this layer. (h-k): Same as panel e, with layer 2 divided into two (see panel (a)), for different values of the total (thermal + compositional) density anomalies δρ₂b in layer 2b. δρ₂b = 5 kg/m³ (resp. 35 kg/m²) corresponds to 2.9 vol % (resp. 11.5 vol%) partial melting over 60 km, for a melt of 350 kg/m³ excess with respect to the PREM reference density model (Suppl. Fig. B1) and a 75 K thermal anomaly.
second layer extending down to the top of the transition zone (160–410 km, layer 2), and a third, bottom layer spanning the ETZ, from 410 down to 1000 km depth (layer 3). In a later stage of analysis, we refine the estimates in the layer 2 by dividing this layer into two sub-layers: a thicker layer 2a (160–350 km depth) and a thinner basal layer 2b (350–410 km depth), based on the average thickness (~60 km) of the seismic low velocity layer (LVL) reported above the transition zone in multiple locations around the globe (Revenaugh and Sipkin, 1994; Vinnik and Farra, 2007; Tauzin et al., 2010; Wei and Shearer, 2017).

We systematically vary the thermal anomalies within the rolls (taking values of +/-75, 100, 150, 200 K) in order to study the trade-off between the thermal and the compositional mass sources, and compute the corresponding negative density anomaly using a value of 2 × 10⁻⁵ K⁻¹ for the thermal expansivity (Stacey and Davis, 2008; SI Section E). The considered range of temperatures is globally consistent with available constraints on the sub-Pacific mantle thermal variations. As an upper bound, the excess temperature of plumes in the upper mantle is estimated to be ~200–300 K (Lin and Van Keken, 2006; Bao et al., 2022), and reaches 300 K in the case of the Hawaiian plume (Ballmer et al., 2011). At the lower end of the estimated variations, sub-lithospheric secondary convection patterns able to explain the 200 km wavelength gravity undulations reported from satellite altimetry involve a ~50–80 K temperature difference between the upwellings and the downwellings (Buck and Parmentier, 1986). On the other hand, low asthenospheric densities and seismic anomalies are associated with thermal anomalies of 150 K west of the East Pacific Rise (Harmon et al., 2011).

## 4.2. Estimation of the densities under an assumption of global equilibrium

To constrain the thermo-chemical nature of the mass sources, we investigate how the temperature (T) anomalies of the hypothesized convection cells trade-off with density excesses in the upwelling limbs. For a given vertical density distribution and a temperature excess ΔT in the upwelling limb, we calculate the Bouguer gravity gradients δg and the seafloor slope gradients δh from the dynamic topography, as predicted in the framework of a dynamic computation (Richards and Hager, 1984; Section 2.2). Thus we solve the system of differential equations that describe the Earth's response to the effects of a mantle internal load in spherical geometry, considering different radial viscosity profiles without (Figure 4ab, Suppl. Figs. E1a-d, E2a-c, E3, E5) and with a LVC in the 60–180 km (Figure 4acd,h-k, Suppl. Figs. E1e-h, E2d-f, E4, E6, E7, E8) or 60–160 km (Fig. 4e-g) depth range.

Furthermore, to evaluate the solutions, we hypothesize that the system of sources is fully equilibrated above the base of the ETZ. Indeed, as there is evidence for decoupling between flow above and below ~1000 km depth from the decorrelation of the long wavelength seismic structure (Rudolph et al., 2015), the deflection of plumes (French and Romanowicz, 2015) and stagnation of slabs (Fukao and Obayashi, 2013), and evidence for an intermittent seismic reflector at ~1000 km (Waszek et al., 2021; Zhang et al., 2023), we believe it is reasonable to assume that the horizontally structured upper mantle flows are stably decoupled from lower mantle circulation below the extended transition zone. Thus we set a compensation depth of ~1000 km for the system of masses, so that the dense mass sources are fully supported within the upper mantle and extended transition zone, independently from the lower mantle.

Next, focusing on the younger Pacific plate, we perform a parameter search for a set of models in which we systematically vary the average temperature within the convection cells and the excess density structure. A periodic signal is entirely characterized by two quantities, its period and its amplitude. Here, the geometry of the sources ensures a similar periodicity for the modelled and observed signals in the upwelling limbs, and only their amplitude needs to be matched. We thus assess the model quality by subtracting the modelled amplitude of this periodic pattern of signals from the observed one (Suppl. Table A1, based on

Table E2 in Panet et al., 2022).

In practice, we compare the modelled values at the center above the upwellings to the observed amplitude of the undulations using a cost function C defined as follows:

$$C^2(\delta\rho_1, \delta\rho_2, \delta\rho_3, \Delta T) = \left( \frac{\delta\rho_{\text{Bouguer}}^{\text{Obs}} - \delta\rho_{\text{Bouguer}}^{\text{Model}}}{\delta\rho_{\text{Bouguer}}^{\text{Obs}}} \right) + \left( \frac{\delta h_{\text{Seafloor Topo}}^{\text{Obs}} - \delta h_{\text{Seafloor Topo}}^{\text{Model}}}{\delta h_{\text{Seafloor Topo}}^{\text{Obs}}} \right) + \left( \frac{B^r - |B^f|}{B^r} \right)^2 \quad (1)$$

C² is a cumulated residual variance. It includes residuals in Bouguer gravity gradients, in seafloor slope gradients and an additional residual buoyancy term expressing the global vertical equilibrium of the system of sources. This last term is the relative difference of the dense compositional mass anomalies with respect to the light mass anomalies associated to the upwellings. It compares the negative buoyancy of the dense region Bf with the overall positive buoyancy Bf of the upwellings (including the deflections of the transition zone interfaces) and seafloor lows. If there is a large difference between the two, either the system is unstable and moving vertically, or it is dynamically supported by viscous stresses from a deeper flow that shares the same horizontal, plate-dependent geometry as the rolls considered. Because the signals are observed over broad oceanic domains, we consider that the system of sources should be stable over the long-term, and fully equilibrated above the base of the ETZ, as discussed above. We compute the total positive

buoyancy Bf as follows: $B^r = \int_{r_0}^{r_1} \delta\rho(r)dr$ where $r_0$ and $r_1$ are the radii of the base of the convective system and its top, respectively. Likewise:

$$B^f = \int_{r_2}^{r_3} \delta\rho(r)dr, \text{ where } r_2 \text{ and } r_3 \text{ are the radii of the bottom and top of the dense layer, respectively. Thus the cost function allows for a trade-off between fitting the observations and fulfilling the prior of vertical equilibrium.}$$

Finally, we require that each term of this cumulative residual variance remains “small”: more than ~90% of the gravity and the seafloor slope gradients variance (~30% of their rms) has to be explained by the model, and the vertical equilibrium should be attained within 30% of error:

$$\left\{ \begin{array}{l} \sqrt{\left( \frac{\delta g_{\text{Bouguer}}^{\text{Obs}} - \delta g_{\text{Bouguer}}^{\text{Model}}}{\delta g_{\text{Bouguer}}^{\text{Obs}}} \right)^2} < 0.3 \\ \\ \sqrt{\left( \frac{\delta h_{\text{Seafloor Topo}}^{\text{Obs}} - \delta h_{\text{Seafloor Topo}}^{\text{Model}}}{\delta h_{\text{Seafloor Topo}}^{\text{Obs}}} \right)^2} < 0.3 \\ \\ \sqrt{\left( \frac{B^r - |B^f|}{B^r} \right)^2} < 0.3 \end{array} \right. \quad (2)$$

We now present the results obtained for the 3- and the 4-layers models.

# 5. Results: trade-offs between thermal and chemical density anomalies

## 5.1. Sensitivity to geometrical and physical parameters (3-layers model)

We first solve for density anomalies in each of the three layers 1, 2, 3, and test the sensitivity of the estimated densities to different geometrical and physical parameters of the sources. The results are presented in detail in SI Section E, where we assess the sensitivity of the estimated densities to the viscosity profile and the presence of a sub-lithospheric low viscosity channel (Suppl. Fig. E1), to an increase in the 660 km
depth Clapeyron slope (Suppl. Fig. E2), and to variations in the top and bottom depths of the convective range (Suppl. Figs. E3-E6). To account for uncertainties in the depth range of the convective system, we varied the depth to the top of the convective system within the sub-lithospheric LVC (Suppl. Figs. E3-E4), and we varied its bottom depth by ± 200 km (Suppl. Figs. E5-E6). A summary of these results is provided below in Section 5.3.

## 5.2. Impact of mass accumulation in a thin layer (4-layers model)

We then refine the estimate of the density anomalies in layer 2 by assessing the trade-off between the total (thermal + compositional) excess densities in layers 2a and 2b (Fig. 4h-k). In particular, we investigate the smallest values of density excess in layer 2a when mass accumulates in the layer 2b. Tables 1 and 2 summarize the compositional density anomalies obtained in layers 2a and 2b for the two cases shown in Figs. 4h and 4j respectively, which correspond to a small 1.5 kg/m³ compositional anomaly in the deep layer 3, and 5 or 35 kg/m³ total anomaly in layer 2b. In the presence of a small density excess in the deeper layer 3, the layer 2a total density decreases to ~ 2 kg/m³, with +/- 75 K thermal anomalies in the convection cells and a 35 kg/m³ total anomaly in layer 2b (Fig. 4j, Table 2): in that case, the compositional contribution of layer 2a (about 7 kg/m³) still corresponds to a mass excess. Note that the deep mass excess in layer 3 must be limited in order for the ascending currents to remain buoyant. Also, in our tests, we find that it cannot exceed ~ 4 kg/m³, while the corresponding total density anomalies in layer 2a still remain positive (Suppl. Fig. E7) or only slightly negative (Suppl. Fig. E7d, but with a degraded fit compared to the positive solutions), such that the associated compositional anomaly is always positive. These conclusions remain valid when decreasing the thickness of layer 2b by a factor of 2 (Suppl. Fig. E8, Tables G1-G2). We obtain similar density anomalies in layer 2a, for approximately twice larger density anomalies in layer 2b compared to Main Fig. 4h-k.

## 5.3. Summary: constraints on the mantle mass anomalies

In Fig. 4c-g we summarize the results for the 170–1000 km depth convective range. Our tests have shown that the presence of a sub-lithospheric LVC (Fig. 4d-g) is required to avoid too large (above 20 kg/m³) densities in layer 1 (Fig. 4, Suppl. Figs. E1–E6), and that the density excess in the LVC (layer 1) is well constrained from the surface dynamic topography, to ~10 kg/m³ with respect to the PREM reference density profile. The layer 2 density anomaly increases as a function of the upwellings’ thermal anomaly within a range of strictly positive values (Fig. 4e), even when increasing the layer 3 density excess in the limited 0–6 kg/m³ range of its possible values (Suppl. Fig. E1e-h). Actually, as mentioned above, the surface dynamic topography constrains the densities in the LVC (layer 1), but the gravity signals require a larger amount of mass excess. This is why some of the excess mass must be distributed deeper within the upper mantle. This system of sources is able to produce seafloor lows above the hot upwellings because of the downward flow induced by the excess mass in layer 1, which overcomes,

in the sub-lithospheric LVC, the upward contributions from the deeper upwellings (Suppl. Fig. F1). Our calculations suggest that the thermal anomaly should not exceed +/- 100 K within the convection cells, in order to avoid too large upper mantle density excesses (Ballmer et al., 2013b), i.e. above 20 kg/m³ in layer 2 (Fig. 4, Suppl. Fig. E1). Further, we find that the bottom of the convective range should not be much shallower than 1000 km, as the light buoyancy of the upwellings would become insufficient to support the mass excess needed in the upper mantle to fit the geodetic data (Suppl. Figs. E5-E6). The tests show that the top of the convection cells can reach the base of the lithosphere in the presence of a sub-lithospheric LVC and for small density anomalies in layer 3 (Suppl. Figs. E3-E4), with a relatively limited impact on the estimated densities in layers 1 and 2, as is also the case when increasing the bottom depth (Suppl. Fig. E6d-f). Whatever the physical source of this mass excess, explaining the geodetic signals thus requires both an LVC beneath the lithosphere and excess mass distributed across the upper mantle in and above the upwellings, with a small potential contribution of the deeper layer 3. This conclusion remains valid also in the case of a mass accumulation in layer 2b (Fig. 4h-k, Suppl. Fig. E7), still calling for near-neutral or preferably positive densities in layer 2a, while the moderate temperatures prevent too large mass accumulation in layer 2b (Fig. 4, Suppl. Fig. E7a-d). Let us recall that the estimated densities correspond to average values over the fixed depth extents of layers 1 to 3 (Section 4.1), thereby showing a limited radial resolution of the sources. Whether these averages result from a radial density gradient, local concentrations of mass or homogeneous sources, cannot be determined from our analysis. We finally note that in these tests, we did not consider the 520 km discontinuity as it would only produce a mass deficit in the upwellings, further stressing the need for additional mass excess therein (SI Section B).

Our modelling is based on a purely radial viscosity structure, but the obtained density estimates should remain valid, to first order, also in the presence of lateral viscosity variations. The lateral viscosity decrease caused by the presence of a 300 K hot, 1000 km wide plume in the upper mantle and transition zone, leads to a 12% decrease of the corresponding dynamic topography (Steinberger et al., 2023). This is less than the impact of lateral viscosity variations associated to the whole mantle thermal structure, stiff slabs, cratonic keels, weak plate boundaries and low viscosity post-perovskite, estimated to add up to 20% of the dynamic topography (Ghosh et al., 2010). Thus, neglecting these variations in our modelling should lead to an over-estimation of the dynamic topography above the hot upwellings, hence a corresponding under-estimation of the mass excess in layer 1 in the presence of a sub-lithospheric LVC, probably by up to ~10–15%. At the same time, the densities below the LVC are not expected to be significantly impacted because of their limited contribution to the surface dynamic topography, as discussed above. The opposite holds over the cold downwellings, where a purely radial viscosity profile produces an over-estimation of the mass sources; however the error should remain small for convection cells located below the LVC, and not larger than ~10–15% in the depth range of the LVC for a shallower top of the convection cells. In this latter case, the volume of the convection cells within the LVC represents only a

**Table 1**

Melt fractions and iron enrichment, corresponding to Figure 4 h. Estimation of melt fraction in layers 2a and 2b in the case of a 350 kg/m³ dense melt, for the best-fitting total density anomalies shown in Figure 4 h. The best-fitting total density anomalies (left column) are split into their thermal (column 2) and compositional component (column 3), purely attributed to melt (column 4) or to iron enrichment of the solid olivine matrix with respect to magnesium (column 5).

<table><thead><tr><th>Case of Figure 4h</th><th>Total density δρ<sub>thermal</sub> + δρ<sub>chemical</sub></th><th>Thermal density δρ<sub>thermal</sub></th><th>Compositional density δρ<sub>chemical</sub></th><th>Melt fraction for a 350 kg/m<sup>3</sup> dense melt (without solid enrichment)</th><th>Percentage of iron enrichment of olivine (without melt)</th></tr></thead><tbody><tr><td rowspan="2">Layer 2a<br/>±75 K</td><td>11 kg/m<sup>3</sup></td><td>-5.25 kg/m<sup>3</sup></td><td>16,25 kg/m<sup>3</sup></td><td>4,6%</td><td>1,1%</td></tr><tr><td>17 kg/m<sup>3</sup></td><td>-7 kg/m<sup>3</sup></td><td>24 kg/m<sup>3</sup></td><td>6,9%</td><td>1,6%</td></tr><tr><td rowspan="2">Layer 2b<br/>±75 K</td><td>5 kg/m<sup>3</sup></td><td>-5.25 kg/m<sup>3</sup></td><td>10.25 kg/m<sup>3</sup></td><td>2,9%</td><td>0,7%</td></tr><tr><td>5 kg/m<sup>3</sup></td><td>-7 kg/m<sup>3</sup></td><td>12 kg/m<sup>3</sup></td><td>3,4%</td><td>0,8%</td></tr></tbody></table>
<caption>**Table 2**</caption>

Melt fractions and iron enrichment, corresponding to Figure 4j. Same as Table 1, for the best-fitting total density anomalies shown in Figure 4j.

<table><thead><tr><th>Case of<br>Figure 4j</th><th>Total density<br>δρ<sub>thermal</sub> +<br>δρ<sub>chemical</sub></th><th>Thermal<br>density<br>δρ<sub>thermal</sub></th><th>Compositional<br>density<br>δρ<sub>chemical</sub></th><th>Melt fraction for a 350 kg/m<sup>3</sup> dense melt (without<br>solid enrichment)</th><th>Percentage of iron enrichment of olivine<br>(without melt)</th></tr></thead><tbody><tr><td>Layer 2a</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>±75 K</td><td>2 kg/m<sup>3</sup></td><td>-5.25 kg/m<sup>3</sup></td><td>7.25 kg/m<sup>3</sup></td><td>2,1%</td><td>0,5%</td></tr><tr><td>±100 K</td><td>9 kg/m<sup>3</sup></td><td>-7 kg/m<sup>3</sup></td><td>16 kg/m<sup>3</sup></td><td>4,6%</td><td>1,1%</td></tr><tr><td>Layer 2b</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>±75 K</td><td>35 kg/m<sup>3</sup></td><td>-5.25 kg/m<sup>3</sup></td><td>40.25 kg/m<sup>3</sup></td><td>11,5%</td><td>2,7%</td></tr><tr><td>±100 K</td><td>35 kg/m<sup>3</sup></td><td>-7 kg/m<sup>3</sup></td><td>42 kg/m<sup>3</sup></td><td>12%</td><td>2,8%</td></tr></tbody></table>

small fraction of their total volume, so that the total error is probably much smaller than 10–15%.

**6. Discussion: convection across phase transitions in a hydrated mantle**

**6.1. Melting and recrystallization in upwellings crossing the transition zone**

What could be the nature of this dense material? A model of eclogite entrainment from subducted oceanic crust, as proposed for the Hawaiian plume (Ballmer et al., 2013b), is ruled out here. Indeed, in the region of Hawaii it involves an overall mass excess only below 300 km depth, due to the phase transition of stishovite into lighter coesite which yields an overall mass default above 300 km depth. Furthermore, eclogite is expected to increase the seismic velocities (Hu et al., 2019). Considering the spatial distribution of the mass excess, mostly above the transition zone and up to the sub-lithospheric LVC, a possible explanation involves

the effect of a moderately hydrous transition zone on the upwelling limbs of the meso-scale convection system (Fig. 5). Indeed, hydrated mantle in the wadsleyite phase, upwelling across the 410 km discontinuity, is expected to partially melt when transforming into olivine (Bercovici and Karato, 2003), and this could be enhanced by the presence of carbonatites (Karato et al., 2007; Dasgupta and Hirschmann, 2006). Depending on the relative content in iron and water, a wide range of densities can be obtained for these deep melts, a topic of current debate among mineral physicists. Some experimental and theoretical studies indicate that such melts can be neutrally buoyant above the transition zone (e.g. Matsukage et al., 2005; Sakamaki et al., 2006; Freitas et al., 2017; Karki et al., 2021; Andrault and Bolfan-Casanova, 2022), while other investigations yield positively buoyant melts except for very iron-rich environments (Drewitt et al., 2022). Independent geodetic and seismological observations as presented here could provide further constraints on this issue. Note that a stable melt layer has been invoked to explain the presence of a ~20 to 100 km thick seismic low velocity layer above the 410 km discontinuity, as detected in some

<figure><img src="image_9.png" /></figure>

**Fig. 5.** Sketch of seafloor topography and geoid signals and their proposed sources. (a) : The seafloor topography and geoid undulations, both referred to a "background" (horizontal lines) above zones of upwellings and downwellings, and elongated along the absolute plate motion (APM) direction (see Fig. 3 in Panet et al., 2022, for the corresponding effects on the seafloor slope gradients and the GRACE gravity gradients). (b): Proposed dynamic model, where hot upwelling limbs coincide with seafloor and geoid lows, and slow seismic velocities in model SEMUCB-WM1. In the presence of a moderate amount of water, stable lenses of Fe-enriched partial melt form within the upwellings as they cross the wadsleyite → olivine phase transition. The upward transport of some of this partially molten and recrystallized material may accumulate in the asthenosphere, at the top of the upwellings, and then be recycled by lateral transport towards ridges or subduction zones. This secondary thermo-chemical convection system could be fed by wet, iron-enriched mantle plumes originating from the lower mantle (dark pink arrows). Vertical scale is not the same in panels a and b. The fluctuations represented here are relative, contrasting the behavior in upwelling limbs (red arrows) and downwelling limbs (blue arrows), with limited sensitivity to the precise depth of the upper mantle discontinuity in a reference 1D model, at the scale of our analysis.
areas and attributed to partial melt (e.g. Revenaugh and Sipkin, 1994; Vinnik and Farra, 2007; Tauzin et al., 2010; Wei and Shearer, 2017). Enriched in FeO and H₂O, the obtained melt fraction can be both denser than the surrounding upper mantle and lighter than transition zone minerals. We propose that a small fraction of the melt produced at the 410 km discontinuity could be entrained upward by the motion of the upwellings, potentially crystallizing and/or promoting chemical exchange by re-crystallization of the surrounding mineral matrix (Das-gupta and Hirschmann, 2006). Partially molten and recrystallized material would finally accumulate in the asthenosphere at the top of the convective upwellings, as detected below the lithosphere in some areas sampled away from ridges (e.g. Kawakatsu et al., 2009; Mehouachi and Singh, 2018; Hua et al., 2023). This could lead to a non-uniform distribution of melt in the asthenosphere. The accumulated dense material would be finally recycled by lateral asthenospheric flow towards ridges or subduction zones (e.g. Toomey et al., 2002; Hoink and Lenardic, 2010).

In what follows, we test whether the estimated densities in the 4-layer model and the corresponding instantaneous flow field support this scenario, and what constraints they bring on the melt fraction and iron-enrichment of the upper mantle. Because the densities remain constant over the chosen domains of the sources, all the buoyancies, compositions and melt fractions discussed should be considered as average values over the same domains.

## 6.2. Gravitational stability and upward motion of the dense material

We first verify the gravitational stability of the dense material in the top and bottom upper mantle layers 1 and 2b. With moderate +/- 75 K thermal anomalies in the meso-scale convection system and a limited compositional excess density in the deep layer 3, the density excess in layer 2b (Fig. 4h-k) is always much smaller than $\Delta\rho_{410} + \delta\rho_3$, where $\Delta\rho_{410}$ is the 180 kg/m³ density jump at the 410 km discontinuity in the PREM reference density profile and $\delta\rho_3$ the density excess in layer 3. This means that the absolute density in layer 2b is smaller than that in layer 3. Consequently, in steady state, the dense lense in layer 2b would largely remain stagnant (Bercovici and Karato, 2003; Karato et al., 2007) within the convective upwellings. Likewise, the top layer 1 can be stable, given its small density excess compared to the ~75 kg/m³ density jump at the base of the sub-lithospheric LVZ in the PREM model, which results in an absolute density smaller than in the convective upwelling limb below it, taking into account its excess density $\delta\rho_{2a}$.

While the deflections of the 410 km interface and the presence of a light layer 3 (in total density) compensate the upper mantle mass excess and ensure the global stability of the whole secondary convection system in the mantle, the local upward entrainment by upwellings of the dense material formed in layer 2b across layer 2a depends on the buoyancy ratio B between the stabilizing compositional buoyancy and the destabilizing thermal buoyancy. Upward entrainment requires B smaller or close to 1; larger values are also possible for a low viscosity contrast of the dense basal layer (e.g. Lin and van Keken, 2006; Davaille, 1999). Here, the density excess in layer 2a can decrease to near-neutral buoyancy for a large enough density excess in layer 2b (Fig. 4j). Thus, we find

$B = 1.4$ in a convection system associated with ± 75 K thermal anomalies, for a 35 kg/m³ total density contrast $\delta\rho_{2b}$ in layer 2b and a small compositional density excess $\delta\rho_3 = 1.5$ kg/m³ in layer 3. This corresponds to a total density anomaly in layer 2a $\delta\rho_{2a} \sim 2$ kg/m³ (Fig. 4j, Table 2). Note that with these values for $\delta\rho_{2a}$ and $\delta\rho_{2b}$, the average density over both layers 2a + 2b is ~9.5 kg/m³, which is consistent with the data (e.g. Suppl. Fig. E1). The obtained compositional density excess in layer 2a is compatible with upward entrainment up to the base of the LVC. Note that low buoyancy ratios in layer 2a require mass accumulation in layer 2b, which supports the considered hypothesis of partial melting across the 410 km discontinuity. The total density excess in layer 2a could be further decreased to zero, corresponding to a case of neutral buoyancy, for a 42 kg/m³ total density anomaly in layer 2b (Table 3). However, the quality of the fit degrades as we progressively reduce the mass excess in layer 2a to zero, as shown by the decreasing number of cost function isolines (Fig. 4h-k), so that there is no solution with a compositional mass default in layer 2a, as also illustrated in Suppl. Fig. E7. Thus, in the context of the studied scenario, the geodetic data favor solutions with a buoyancy ratio slightly above 1 in layer 2a. To confirm that upward entrainment of dense material across layers 2a and 2b is possible, we have computed the instantaneous flow field associated with the system of sources shown in Fig. 4j and Table 2, for ± 75 K thermal anomalies. The results indeed feature positive upward velocities in the upwellings, and downward flow in the downwellings (Suppl. Fig. F1).

## 6.3. Density constraints on upper mantle melt fraction and iron-enrichment

According to our modeling, the compositional density excess of layer 2b should be larger than ~40 kg/m³ in order to not violate the possibility of upward entrainment of the dense material across layer 2a. At the same time, it cannot increase much above ~65 kg/m³ to avoid a degradation of the fit to the geodetic data (Fig. 4k, Suppl. Fig. E7). Theoretical and experimental investigations (Freitas et al., 2017; Wei and Shearer, 2017; Hammond and Humphreys, 2000) suggest at most ~2% melt in the LVL above the transition zone, based on the observed seismic velocity reduction in that layer. Attributing the entire density excess in layer 2b to melt would require too much melt and is thus incompatible with current estimates of melt densities and melt fractions (SI Section G). The melt fraction can be reduced, and the seismic and geodetic data reconciled, if the detected source of mass excess also results from crystallized iron-enriched minerals in both layers 2a and 2b. Because of the large difference in their respective atomic masses, and assuming iron dominates the computed density variations (e.g. Drewitt et al., 2022), only 1% iron enrichment of olivine with respect to magnesium is indeed able to produce a 15 kg/m³ density anomaly. Thus, the 40 kg/m³ compositional density excess in layer 2b (Table 2) could be explained by 2.7% iron enrichment of olivine with respect to magnesium (Fig. 6a, for a 0% melt fraction abscissa), or a combination of not more than 5% of melt and ~1.5 to 2.7% iron enrichment of the solid matrix (Fig. 6a). In this layer, 5% of partial melting with a gravitationally stable melt of excess density < ~200 kg/m³, could be consistent with 15–18 wt

<caption>**Table 3**</caption>

Melt fractions and iron enrichment, in the case of neutral buoyancy of the layer 2a. Same as Table 1, for the best-fitting total density anomalies corresponding to a neutral buoyancy of the layer 2a.

<table><thead><tr><th>Neutral buoyancy in layer 2a</th><th>Total density<br>δρ<sub>thermal</sub> + δρ<sub>chemical</sub></th><th>Thermal density<br>δρ<sub>thermal</sub></th><th>Compositional density<br>δρ<sub>chemical</sub></th><th>Melt fraction for a 350 kg/m<sup>3</sup> dense melt<br>(without solid enrichment)</th><th>Percentage of iron enrichment of olivine<br>(without melt)</th></tr></thead><tbody><tr><td>Layer 2a</td><td>scheme</td><td></td><td></td><td></td><td></td></tr><tr><td>±75 K</td><td>0 kg/m<sup>3</sup></td><td>-5.25 kg/m<sup>3</sup></td><td>5.25 kg/m<sup>3</sup></td><td>1.5%</td><td>0.35%</td></tr><tr><td>±100 K</td><td>0 kg/m<sup>3</sup></td><td>-7 kg/m<sup>3</sup></td><td>7 kg/m<sup>3</sup></td><td>2%</td><td>0.5%</td></tr><tr><td>Layer 2b</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>±75 K</td><td>42 kg/m<sup>3</sup></td><td>-5.25 kg/m<sup>3</sup></td><td>47.25 kg/m<sup>3</sup></td><td>13.5%</td><td>3.1%</td></tr><tr><td>±100 K</td><td>66 kg/m<sup>3</sup></td><td>-7 kg/m<sup>3</sup></td><td>73 kg/m<sup>3</sup></td><td>20.9%</td><td>4.9%</td></tr></tbody></table>
<figure><img src="image_10.png" /><figcaption>**Fig. 6.** Iron enrichment and melting in layers 2a and 2b. Panel a: percentage of iron enrichment of solid olivine in layer 2b (right colorbar), corresponding to the solution of Figure 4j for a ± 75 K thermal anomaly (40.25 kg/m³ compositional density anomaly in this layer), as a function of melt fraction and melt average density excess. Panel b: possible compositions of melt in layer 2b (with FeO content in color), as a function of melt average density excess (right vertical axis), irrespective of melt fraction. Compositions have been calculated based on Andrault & Bolfan-Casanova (2022, their Figure 5a). Note that a higher iron content corresponds to lower degrees of melting and lower melt fractions (case of batch melting). Two example values of melt fractions as a function of its composition are indicated, based on Freitas et al. (2017, Suppl. Info); the discussed melt with 350 kg/m³ excess density is obtained for a 0.7% melt fraction in their experiments. Panel c: same as panel a, for layer 2a (7.25 kg/m³ compositional density anomaly). The 0% iso-contour corresponds to the case where the entire compositional density anomaly in layer 2a or 2b is explained by melt. The case where it is entirely explained by enriched solid olivine corresponds to the 0% melt fraction (abscissa). Dashed gray line in panels b and c: 180 kg/m³ limit of gravitational stability for the melt above the 410 km discontinuity, without any dynamic support by upwellings (see SI Section G). Blue ellipses in panels a-c indicate the space of acceptable parameter values.</figcaption></figure>

% FeO and a few wt% H₂O in the melt (Fig. 6b). At the same time, the 7.25 kg/m³ compositional anomaly of layer 2a could reflect up to 0.5% iron enrichment of the solid olivine with respect to magnesium (Table 2, Fig. 6c), and only a very small amount of melt. Such small variations in the Mg/(Mg+Fe) ratio of solid olivine in both layers 2a and 2b are within the natural range of variability for the forsteritic composition, which comprises an 8–12% Fe component (Pearson et al., 2003). Our estimates of iron enrichment of the solid matrix are also well below the variations reported when silicate melts crystallize or react with the surrounding rocks, which may lead to an >15% iron component of the solid olivine (Pearson et al., 2003). These results would also be consistent with electromagnetic and magnetotelluric data, which imply an almost melt-free upper mantle above the LVL (Freitas and Manthilake, 2019), while allowing for the presence of melt across the whole upper mantle (Hirschmann et al., 2005).

Note that entrainment of dense, crystallizing melts by mantle upwellings has been proposed to explain how dense lunar magmas can rise for hundreds of kilometers up to the surface of the Moon (Hess, 2000). If the upward velocity of the upwelling instability exceeds the downward velocity of the melt, the melt can rise with respect to the background mantle while sinking with respect to the upwelling: stirring overcomes the buoyancy effect (Hess, 1991). As mentioned earlier, such a model is also reinforced by the failure of a scenario where, instead of dense sources in the upwellings, light sources would be present on top of the downwellings. In addition to failing to produce the asymmetry in the observed coherency between gravity and topography data, the light density anomalies would be too large and unable to jointly explain the gravity and the topography signals (SI Section H and Suppl. Table H1).

## 6.4. An active role of the deep upwelling flow

The thermo-chemical model we present is able to explain our observations on the younger side of the Pacific plate. Actually, meso-scale gravity and seafloor topography undulations are also detected on the older side of this plate (Adam et al., 2021; Panet et al., 2022), and the Pacific ocean basin shows a clear dichotomy. The younger side of the plate exhibits larger seismic anomalies, lower amplitude gravity and topography signals, and a concentration of deeper plumes. The opposite holds on the older side, where the amplitude of the gravity and

topography undulations increases, the slow seismic fingers disappear (Fig. 1abc), and fewer, cooler hotspots are found (Bao et al., 2022). A model of lithospheric drips has been proposed for the older side of the Pacific plate (Adam et al., 2021), where 1000-km spaced, cold instabilities dripping from the base of a thickened lithosphere sink and are sheared away by the large-scale mantle flow (see also Coltice et al., 2018). However, this model is not able to explain our observations on the younger side of the Pacific (SI Section I). There, cold lithospheric drips beneath seafloor lows would be located above the hot and slow seismic fingers, which is unlikely and underlines an important difference of our observations compared to those of Adam et al. (2021) on older seafloor, which is the puzzling correspondence of the slow seismic fingers with seafloor lows we find on younger seafloor. These conclusions suggest a change in the form of the meso-scale convection between the younger/hotter, and older/colder sides of the Pacific ocean basin, stressing the difference of our model for the younger side, compared to purely plate-driven secondary convection as proposed by Adam et al. (2021) for the older side.

Indeed, the model of deep meso-scale convective flow we investigate is probably not everywhere imposed by plate motions, in contrast to Richter rolls, but rather originates at depth, because its geodetic and seismic fingerprints are detected under both fast-moving and slow-moving plates and across spreading ridges (see discussion in Panet et al., 2022). Thus, it should be driven, at least partially, by thermal and maybe also compositional anomalies at the base of the extended transition zone. In particular, thermal anomalies can be expected as this region of the Pacific mantle contains some of the hottest mantle plumes (Bao et al., 2022). Our model then implies that the bottom boundary layer of this convective flow can be rooted at the base of the extended transition zone near 1000-km depth. This is consistent with our assumption of decoupling of the flow above and below this depth, and supports inferences of a viscosity increase around this depth, thought to explain the stagnation of slabs, deflections of plumes and changes in the large-scale mantle structure (Rudolph et al., 2015). Finally, we hypothesize that wet, iron-rich plumes (e.g. Simons et al., 2022) in the lower mantle (Tenner et al., 2012) as well as return flow associated with subduction might contribute to feeding this meso-scale circulation and the compositional anomalies in its active upwellings.

Deciphering how this meso-scale flow, at least partly driven from
below, interacts with plate motions probably requires a better understanding of its interplay with convective patterns at smaller spatial scales, such as are reflected in APM-oriented geoid and seafloor topography undulations, detected over a broad range of wavelengths, from 200 to 1400 km, some commensurate with the depths of seismic discontinuities in the upper to mid-mantle (e.g. Cazenave et al., 1992; Wessel et al., 1994). At 200 km wavelength, they have been interpreted as evidence for secondary-scale thermal convection in the upper mantle below the cooling lithosphere of fast-moving plates, organized by the vertical shear due to plate motions and promoted by the existence of a sub-lithospheric low viscosity layer (e.g. Richter and Parsons, 1975; Buck and Parmentier, 1986). Alternatively, models of viscous fingering, where the alignment results from the injection of a low viscosity fluid into a higher viscosity one between two close parallel plates, have also been proposed (e.g. Weeraratne et al., 2007, Harmon et al., 2011). Three-dimensional dynamic models of the interacting deeper and shallower mantle circulations at different scales, constrained by the three sets of seismic, gravity and topography data, are needed to further describe how the circulation in the upper mantle is coupled to the Earth's deeper interior dynamics and to the tectonic plates (Coltice et al., 2017), and how it evolves in time.

## 7. Conclusions and outlook

This study shows how the combination of seismic and geodetic observations with mass balance models in a viscous earth can bring insights on the present-day meso-scale structure of the convective mantle. Even if the actual flow is most likely more complex and three-dimensional, the simplified two-dimensional models explored here provide first order evidence for the thermo-chemical nature of the corresponding dynamics, constraining the depth extent of the meso-scale convection from the lithosphere down to ~1000 km depth, with cells of aspect ratio ~1:1. We have shown that the observations can be reconciled by taking into account the effect of the 410 km depth olivine phase transition on bi-dimensional convection cells crossing a moderately hydrous transition zone, leading to the formation of dense, iron-enriched melts and recrystallized minerals in the upwellings. This suggests that the thin LVL atop the 410 km discontinuity observed in the study of seismic receiver functions and attributed to dehydration melting may only be observed intermittently in ocean basins away from subduction zones (Wei and Shearer, 2017). While seismic data provide constraints on the velocity reduction, our results provide new and complementary constraints on its density.

These results are based on a present-day mantle mass model, able to explain the three sets of periodic, 2000 km wavelength APM-oriented signals detected in the seismic, gravity and topography observables, thereby addressing the joint analysis of these observables at a rarely considered intermediate scale of mantle dynamics. At this stage, a coupling of the seismic tomography model and the geodetic signals using density-to-velocity scaling factors (e.g. Forte and Mitrovica, 2001; Tiberi et al., 2003, Simmons et al., 2010), or by converting a thermal and compositional structure into density and velocities (Trampert et al., 2004), is outside of the scope of this study but would be important to address in future work. For that, it would be advantageous to enhance the lateral resolution of the seismic tomography model using the gravity gradient signals, in a joint inversion of these signals with the seismic waveform and wavespeed data. Nevertheless, the average densities over the considered source domains obtained here can be used for comparison with similar averages built from 3-dimensional mantle convection models.

## Data statement

The seismic and geodetic signals studied in this work have been published in French et al. (2013), Panet et al. (2022). The software developed by the authors to perform the multiscale gravity gradient and

seafloor slope gradient analysis was introduced in Panet (2018), Panet et al. (2022); the corresponding codes are available at https://doi.org/10.6084/m9.fishshare.20255679. The software used by the authors to perform the dynamical modelling of the signals is available upon request.

## CRediT authorship contribution statement

Isabelle Panet: Conceptualization, Methodology, Software, Writing – original draft, Writing – review & editing. Marianne Greff-Lefftz: Methodology, Software, Writing – review & editing. Barbara Romanowicz: Conceptualization, Methodology, Writing – review & editing.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Data availability

Data will be made available on request.

## Acknowledgments

This work was supported by CNES. It is based on observations from the GRACE satellite mission as well as observations from Cryosat-2, Jason-1 and Envisat. We thank Denis Andrault and Jie Deng for helpful discussions about the properties of deep melts, and Mark Richards for helpful suggestions on the dynamics. We thank an anonymous reviewer and Manuele Faccenda for their constructive reviews which contributed to improve the manuscript. BR acknowledges support from NSF grant EAR-1800324. The authors used the GMT software (Wessel and Smith, 1995) for plotting. This study contributes to the IdEx Université de Paris ANR-18-IDEX-0001.

## Supplementary materials

Supplementary material associated with this article can be found, in the online version, at doi:10.1016/j.epsl.2024.118763.

## References

* Adam, C., Vidal, V., Pandit, B., Davaille, A., Kempton, P.D., 2021. Lithosphere destabilization and small-scale convection constrained from geophysical data and analogical models. *Geochem., Geophys., Geosyst.* 22, e2020GC009462.

* Andrault, D., Bolfan-Casanova, N., 2022. Mantle rain toward the Earth's surface: a model for the internal cycle of water. *Phys. Earth Planet. Interiors* 322, 106815.

* Ballmer, M.D., Ito, G., van Hunen, J., Tackley, P.J., 2011. Spatial and temporal variability in Hawaiian hotspot volcanism induced by small-scale convection. *Nat. Geosci.* 4 (7), 457–460.

* Ballmer, M., Ito, G., Wolfe, C., Solomon, S.C., 2013b. Double layering of a thermochemical plume in the upper mantle beneath Hawaii. *Earth Planet. Sci. Lett.* 376, 155–164.

* Bao, X., Lithgow-Bertelloni, C.R., Jackson, M.G., Romanowicz, B., 2022. On the relative temperatures of Earth's volcanic hotspots and mid-ocean ridges. *Science* (1979) 375 (6576), 57–61.

* Bercovici, D., Karato, S.-I., 2003. Whole-mantle convection and the transition-zone water filter. *Nature* 425, 39–44.

* Bergé, P., Dubois, M., 1984. Rayleigh-Bénard convection. *Contemp. Phys. 25* (6), 535–582.

* Bird, P., 2003. An updated digital model of plate boundaries. *Geochem. Geophys. Geosyst.* 4, 1027–1078.

* Buck, W.R., Parmentier, E.M., 1986. Convection beneath young oceanic lithosphere: implications for thermal structure and gravity. *J. geophys. Res.* 91, 1961–1974.

* Cazenave, A., Houry, S., Lago, B., Dominh, K., 1992. Geosat-derived geoid anomalies at medium wavelength. *J. Geophys. Res.* 97, 7081–7096.

* Coltice, N., Graault, M., Ulvrova, M., 2017. A mantle convection perspective on global tectonics. " *Earth. Sci. Rev.* 165, 120–150. https://www.sciencedirect.com/science/journal/001282520/".
Coltice, N., Larrouturou, G., Debayle, E., Garnero, E.J., 2018. Interactions of scales of convection in the Earth's mantle. Tectonophysics. 746, 669-677.

Davaille, A., 1999. Simultaneous generation of hotspots and superswells by convection in a heterogeneous planetary mantle. Nature 402, 756-760.

Dasgupta, R., Hirschmann, M., 2006. Melting in the Earth's deep upper mantle caused by carbon dioxide. Nature 440, 659-662.

Drewitt, J.W.E., Walter, M.J., Brodholt, J.P., Muir, J.M.R., Lord, O.T., 2022. Hydroous silicate melts and the deep mantle H₂O cycle. Earth Planet. Sci. Lett. 581, 117408.

Forte, A., Mitrovica, J.X., 2001. Deep-mantle high-viscosity flow and thermochemical structure inferred from seismic and geodynamic data. Nature 410, 10491053.

Freitas, D., Manthilake, G., Schiavi, F., Chantel, J., Bolfan-Casanova, N., Bouhifd, M.A., Andrault, D., 2017. Experimental evidence supporting a global melt layer at the base of the Earth's upper mantle. Nat. Commun. 8, 2186.

Freitas, D., Manthilake, G., 2019. Electrical conductivity of hydrous silicate melts: implications for the bottom-up hydration of Earth's upper mantle. Earth Planet. Sci. Lett. 523, 115712.

French, S.W., Lekic, V., Romanowicz, B., 2013. Waveform Tomography Reveals Channeled Flow at the Base of the Oceanic Asthenosphere. Science (1979) 342, 227-230.

French, S.W., Romanowicz, B., 2015. Broad plumes rooted at the base of the Earth's mantle beneath major hotspots. Nature 525, 95-99.

Fukao, Y., Obayashi, M., 2013. Subducted slabs stagnant above, penetrating through, and trapped below the 660 km discontinuity. J. Geophys. Res. 118, 5920-5938.

Ghosh, A., Becker, T.W., Zhong, S.J., 2010. Effects of lateral viscosity variations on the geoid. Geophys. Res. Lett. 37, L01301.

Hammond, W.C., Humphreys, E.D., 2000. Upper mantle seismic wave velocity: effects of realistic partial melt geometries. J. Geophys. Res. 105 (B5), 10975-10986.

Harmon, N., Forsyth, D.W., Weeraratne, D.S., Yang, Y., Webb, S.C., 2011. Mantle heterogeneity and off axis volcanism on young Pacific lithosphere. Earth Planet. Sci. Lett. 311, 306-315.

Hayn, M., Panet, I., Diament, M., Holschneider, M., Mandea, M., Davaille, A., 2012. Wavelet-based directional analysis of the gravity field: evidence for large-scale geoid undulations. Geophys. J. Int. 189 (3), 1430-1456.

Hess, P.C., 1991. Diapirism and the origin of high TiO₂ mare glasses. Geophys. Res. Lett. 18 (11), 2069-2072.

Hess, P.C., 2000. On the source regions for mare picrite glasses. J. Geophys. Res. 105 (E2), 4347-4360.

Hirschmann, M., Aubaud, C., Withers, A., 2005. Storage capacity of H₂O in nominally anhydrous minerals in the upper mantle. Earth Planet. Sci. Lett. 236, 167-181.

Hoggard, M.J., Winterbourne, J., Czarnota, K., White, N., 2017. Oceanic residual depth measurements, the plate cooling model, and global dynamic topography. J. Geophys. Res. Solid Earth 122, 2328-2372.

Holdt, M.C., White, N.J., Stephenson, S.N., Conway-Jones, B.W., 2022. Densely sampled global dynamic topographic observations and their significance. J. Geophys. Res.: Solid Earth 127, e2022JB024391.

Hoink, T., Lenardic, A., 2010. Long wavelength convection, Poiseuille-Couette flow in the low-viscosity asthenosphere and the strength of plate margins. Geophys. J. Int. 180, 23-33.

Hu, J., Liu, L., Faccenda, M., Zhou, Q., Fischer, K., Marshak, S., Lundstrom, C., 2019. Modification of the Western Gondwana craton by plume-lithosphere interaction. Nat. Geosci. 11, 203-210.

Hua, J., Fischer, K., Becker, T., Gazel, E., Hirth, G., 2023. Asthenospheric low-velocity zone consistent with globally prevalent partial melting. Nat. Geosci. 16, 175-181.

Karato, S.-I., Bercovici, D., Leahy, G., Richard, G., Jing, Z., 2007. The transition-zone water filter model for global material circulation: where do we stand? Earth's deep water cycle. Geophys. Monogr. Ser. 168, 289-313.

Karki, B.B., Ghosh, D.B., Karato, S.-I., 2021. Behavior and properties of water in silicate melts under deep mantle conditions. Sci. Rep. 11, 10588.

Katzman, R., Zhao, L., Jordan, T.H., 1998. High-resolution, two-dimensional vertical tomography of the central Pacific mantle using ScS reverberations and frequency-dependent travel times. J. Geophys. Res. 103 (B8), 17933-17971.

Kawakatsu, H., Kumar, P., Takei, Y., Shinohara, M., Kanazawa, T., Araki, E., Suyehiro, K., 2009. Seismic evidence for sharp lithosphere-asthenosphere boundaries of oceanic plates. Science (1979) 324, 499-502.

Lei, W., Ruan, Y., Bozdag, E., Peter, D., Lefebvre, M., Komatitsch, D., Tromp, J., Hill, J., Podhorszki, N., Pugmire, D., 2020. Global adjoint tomography-model GLAD-M25. Geophys. J. Int. 223, 1-21.

Lekic, V., Romanowicz, B., 2011. Inferring uppermantle structure by full waveform tomography with the spectral element method. Geophys. J. Int. 185 (2), 799-831.

Lin, S.-C., van Keken, P.E., 2006. Dynamics of thermochemical plumes: 1. Plume formation and entrainment of a dense layer. Geochem., Geophys., Geosyst. 7 (2), Q02006.

Matsukage, K.N., Jing, Z., Karato, S.-I., 2005. Density of hydrous silicate melt at the conditions of Earth's deep upper mantle. Nature 438, 488-491.

Mehouachi, F., Singh, S.C., 2018. Water-rich sublithospheric melt channel in the equatorial Atlantic ocean. Nat. Geosci. 11, 65-69.

Panet, I., Pajot-Métivier, G., Greff-Leffitz, M., Métivier, L., Diament, M., Mandea, M., 2014. Mapping the mass distribution of Earth's mantle using satellite-derived gravity gradients. Nat. Geosci. 7, 131-135.

Panet, I., 2018. An analysis of gravitational gradients in rotated frames and their relation to oriented mass sources. J. Geophys. Res. 123, 11062-11090.

Panet, I., Greff-Leffitz, M., Romanowicz, B., 2022. Dense mantle flows periodically spaced below ocean basins. Earth Planet. Sci. Lett. 594, 117745.

Pearson, D.G., Canil, D., & Shirey, S.B. (2003). Mantle samples included in volcanic rocks: xenoliths and diamonds, In Treatise on Geochemistry, Volume 2: The Mantle and the Core, Holland H.D. & Turekian, K.K., Elsevier, pp 171-275 ISBN 978-0-08-043751-4.

Revenaugh, J., Sipkin, S.A., 1994. Seismic evidence for silicate melt atop the 410-km mantle discontinuity. Nature 369, 474-476.

Richards, M.A., Hager, B.H., 1984. Geoid anomalies in a dynamic Earth. J. Geophys. Res. 89, 5987-6002.

Richter, F., Parsons, B., 1975. On the interaction of two scales of convection in the mantle. J. Geophys. Res. 80, 2529-2541.

Rudolph, M.L., Lekic, V., Lithgow-Bertelloni, C., 2015. Viscosity jump in Earths mid-mantle. Science (1979) 6266, 1349-1352.

Sakamaki, T., Suzuki, A., Ohtani, E., 2006. Stability of hydrous melt at the base of the Earth's upper mantle. Nature 439, 192-194.

Seton, M., et al., 2012. Global continental and ocean basin reconstructions since 200 Ma. Earth Sci. Rev. 113, 212-270.

Simmons, N.A., Forte, A.M., Boschi, L., Grand, S.P., 2010. GyPSuM: a joint tomographic model of mantle density and seismic wave speeds. J. Geophys. Res. 115, B12310.

Simons, K., Dixon, J., Schilling, J.G., Kingsley, R., Poreda, R., 2022. Volatiles in basaltic glasses from the Easter-Salas y Gomez Seamount Chain and Easter Microplate: implications for geochemical cycling of volatile elements. Geochem, Geophys., Geosys 3, 1039.

Stacey, F.D., Davis, P.M., 2008. Physics of the Earth. Cambridge University Press, Cambridge, UK. Fourth Edition.

Steinberger, B., 2000. Plumes in a convecting mantle: models and observations for individual hotspots. J. Geophys. Res. 105, 11127-11152.

Steinberger, B., Grasnick, M.-L., Ludwig, R., 2023. Exploring the origin of geoid low and topography high in West Antarctica: insights from density anomalies and mantle convection models. Tektonika 1 (2), 118-139.

Tauzin, B., Debayle, E., Vittlinger, G., 2010. Seismic evidence for a global low-velocity layer within the Earths upper mantle. Nat. Geosci. 3, 718-721.

Tenner, T.J., Hirschmann, M.H., Withers, A.C., Ardia, P., 2012. H₂O storage capacity of olivine and 406 low-Ca pyroxene from 10 to 13 GPa: consequences for dehydration melting above the transition zone. Contrib. Mineral. Petrol. 163, 297-316.

Trampert, J., Deschamps, F., Resovsky, J., Yuen, D., 2004. Probabilistic tomography maps chemical heterogeneities throughout the lower mantle. Science (1979) 306, 853-586.

Tiberi, C., Diament, M., Déverchère, J., Petit-Mariani, C., Mikhailov, V., Tikhotsky, S., Achauer, U., 2003. Deep structure of the Baikal rift zone revealed by joint inversionof gravity and seismology. J. Geophys. Res. 108 (B3), 2133.

Toomey, D.R., Wilcock, W.S.D., Conder, J.A., Forsyth, D.W., Blundy, J.D., Parmentier, E. M., Hammond, W.C., 2002. Asymmetric mantle dynamics in the MELT region of the East Pacific Rise. Earth Planet. Sci. Lett. 200, 287-295.

Vinnik, L., Farra, F., 2007. Low S velocity atop the 410-km discontinuity and mantle plumes. Earth Planet. Sci. Lett. 262, 398-412.

Waszek, L., Tauzin, B., Schmerr, N.C., Ballmer, M.D., Afonso, J.C., 2021. A poorly mixed mantle transition zone and its thermal state inferred from seismic waves. Nat. Geosci. 14, 949-955.

Weeraratne, D.S., Forsyth, D.W., Yang, Y., Webb, S.C., 2007. Rayleigh wave tomography beneath intraplate volcanic ridges in the South Pacific. J. Geophys. Res. 112, B06303.

Wei, S.S., Shearer, P.M., 2017. A sporadic low-velocity layer atop the 410 km discontinuity beneath the Pacific Ocean. J. Geophys. Res. Solid Earth 122, 5144-5159.

Wessel, P., Bercovici, D., Kroenke, L., 1994. The possible reflection of mantle discontinuities in Pacific geoid and bathymetry. Geophys. Res. Lett. 21 (18), 1943-1946.

Wessel, P., Smith, W.H.F., 1995. New version of the generic mapping tool. EOS Trans. Am. Geophys. Un. 76 (33), 329.

Zhang, Z., Irving, J.C.E., Simons, F.J., Alkhalifah, T., 2023. Seismic evidence for a 1000 km mantle discontinuity under the Pacific. Nat. Commun. 14, 1714.