Lens effect of remnant blocks on deep mantle upwelling causing anomalous subsidence

Domain: earth_science
Published: 2025-08-15
Editorial outcome: minor_revision

Abstract
Classical mantle convection models predict a broad surface uplift over a lower mantle upwelling. However, recent studies have identified anomalously localized surface subsidence above seismically imaged lower mantle upwellings, particularly in regions where upwellings are impeded by subducted/delaminated blocks not currently connected to a subducting/delaminating lithosphere (‘remnant blocks’ for simplicity), e.g., in the western USA and the South China Sea. Known geological processes cannot fully explain the observed localized subsidence, and its locality implies a strong association with the underneath upwelling. Here, we use numerical models to quantitatively explore the contribution of lower mantle upwelling to surface topography evolution. Our results demonstrate that the divergent mantle flow caused by lower-mantle upwelling can stretch the overlying lithosphere, inducing broad subsidence that can be reversed when the upwelling reaches the lithospheric bottom. Moreover, the interaction between the remnant block and upwelling can extend the duration of subsidence and focus the subsidence into a narrow region, which we call the lens effect of remnant blocks. Similar abnormal subsidence events in northeastern Asia further show the potential broad applicability of the proposed mechanism, although, in specific research areas, other shallow mechanisms could also have played important roles.

Introduction
Unlike the “ideal” scenario of plume-oceanic plate interaction where the surface develops a long-wavelength topographical swell1 (Fig. 1a), the rise of a plume toward continental lithosphere with a complex and layered rheological structure could cause “repetitive” topography variation across a broadly uplifted region2. Depending on the stress state of the lithosphere, ripple-shaped or Basin-and-Range-like topography variations (Fig. 1b) could potentially develop3. A previous numerical study further raised the possibility that the ponding of plume materials beneath the Mantle Transition Zone (MTZ) could induce broad surface subsidence (Fig. 1c)4. Therefore, the topographical response to deep-seated plume activities may potentially be diverse and might vary between different tectonic environments.

Recent seismic tomography models and geological studies potentially reveal an additional type of topographical response to mantle-plume-related activity, which features anomalous subsidence in elongated regions of thin lithosphere. As shown in Fig. 2, extended but narrow (tens of a few hundred kilometers wide) subsidence has developed in the vicinity of Yellowstone-Snake-River Plain (YSRP)5–7 and the South China Sea (SCS)8. Notably, these subsidence centers appear to lie above seismically imaged lower-mantle upwelling that may be impeded by cold anomalies9–11 (Fig. 2b, e), i.e., subducted/delaminated blocks not currently connected to a subducting/delaminating lithosphere10–12 (here, we will call these “remnant blocks”). In both scenarios, geological and geophysical observations suggest the presence of local ~1–2 km anomalous subsidence that cannot be easily explained by the regional tectonic evolution5–7,13.

The observed extra surface subsidence has been traditionally believed to reflect the dynamic effects of nearby oceanic subduction14. However, because the dynamic effect of subduction tends to affect a broad region (1000s of km)15, it cannot easily explain the intense subsidence in the ~100–200 km banded region of the western USA and the SCS7,8. In the YSRP, the depressed topography was suggested to be due to the existence of infiltrated mafic melts at the middle-upper crust depths, assuming that the YSRP’s crustal and lithospheric thicknesses are roughly equal to the surroundings16. However, more recent studies reveal that the YSRP crust is ~10–15 km thicker than the crust in the adjacent Northern Basin and Range (NBR) (Fig. S1)17. Even considering its basalt-enriched composition, the predicted elevation of the YSRP is ~1 km higher than observed18. In addition, anomalous subsidence could also be caused by the localized crustal thinning due to the lower crustal flow during the plume’s melt infiltration5. However, subsidence has been inferred to have happened nearly simultaneously throughout the YSRP and/or NBR region (Fig. 2a)7,19; In both YSRP and SCS, the Moho beneath the subsidence center is not shallower than beneath its surroundings20,21. Therefore, other deep-seated processes, like the previously discussed plume-upwelling-related activities (Fig. 1b, c) nowadays near both regions (Fig. 2b–e), could have played a role in causing these anomalous subsidence events.

As noted above, the topographic responses predicted by existing models (Fig. 1) do not appear to easily match the scenario of localized subsidence (Fig. 2a, c). In addition, as the temporal effect of a plume-remnant block interaction (Fig. 2b, e) was not considered, the specific mechanisms operating in this scenario remain to be analyzed.

Here, using numerical modeling, we study the surface topographic responses to interactions between lower mantle upwelling and remnant blocks. Our modeling results reveal that lower mantle upwelling can temporarily induce broad surface subsidence, whose magnitude and duration can be extended due to the existence of sinking remnant blocks above the upwelling’s path.

Results
We design three types of numerical models to explore the effects of plume-like mantle upwelling: 1) an initial hot mantle portion in the lower mantle only to trigger plume-like upwelling22 (Fig. S3); 2) a remnant block only in the upper mantle; and 3) both initial structures (Fig. 10 in the methods and S4 in the supplements). For each type, we summarize 623, 120, and 120 test results, respectively, to show the model’s sensitivity to the critical variables considered here.

Topographic response to lower mantle upwelling
In general, most models with a deep mantle upwelling predict dome-like surface uplift before the upwelling reaches the lithosphere-asthenosphere boundary (LAB) (Figs. S5b–S5c). When the plume materials reach the LAB, they impose divergent horizontal drag and cause short-wave-length topography variations in the uplifted region (Figs. S5d–S5f vs. 1b)2.

However, in some instances, the models predict surface subsidence with horizontal stretching before the upwelling reaches the LAB. Run 4 represents a typical model exhibiting this behavior (Fig. 3). In this model, lower mantle thermal anomaly drives mantle upwelling and convection. The surface experiences ~1 km topography variation within a 1500 km wide horizontal range at ~2.5–3 Myr (Fig. 3b, c), consistent with previous studies (Fig. 1c). By further analyzing the deviatoric stress within the lithosphere, we find that a brief, but prominent extensional period happens just before the “plume” arrives at the LAB, which explains the crustal thinning and surface subsidence (Fig. 3g–i). This extensional episode is mainly due to the basal drag of divergent mantle flow caused by the lower mantle upwelling (Fig. 3c). However, when the plume meets the LAB, it uplift the LAB (Fig. 3g) and gradually reverses the topography (Fig. 3d–e, i). Towards the end of the model run, the topography changes back to a quasi-flat state (Fig. 3f).

More sensitivity tests for the models with lower mantle upwelling
Since the surface response to hot lower-mantle upwelling is affected by the thickness of the overriding plate and the phase transitions in the MTZ, in the models with lower mantle upwelling only, we mainly test seven variables to explore the full spectrum of possible dynamic regions, including the thickness of overriding plate (Fig. 5), the Clapeyron slope (γ410 and γ660, Fig. 4)23 and the density jump ratio (Δρ410 and Δρ660, Fig. 4) at the 410 km24,25 and 660 km26–28, and the belt width of the two mantle phase transitions (Fig. S7)29. We summarize model results in the form of checkbox diagrams to demonstrate the value of surface uplift or subsidence (Figs. 4–5, S7). Through the sensitivity test, we find that 1) surface subsidence mainly happens when Δρ660 is larger than 7%26–28, and γ660 is >-2.8 (Figs. 3–4); 2) the magnitude of topography uplift/subsidence increases with the decrease of the lithosphere thickness (Fig. 5a-d) while the wavelength of topographical variation does not change obviously, which is consistent with the previous 3D modeling results30 (Fig. 5e); and 3) the transition belt widths do not affect the modeling results, except when they are very narrow (i.e., 1 km) (Fig. S7).

Topographic response to a sinking remnant block
We also quantify the surface response to a remnant block in the deep upper mantle. In all models, the sinking cold block always causes surface subsidence (Fig. S8)15 with all tested values of initial length, thermal age (for calculating the thermal structure based on the half-space cooling model31), and the duration of the block’s stay within the mantle (e.g., Fig. S9). However, the subsidence magnitude is small (~200 m) (Fig. S8), not enough to explain the observed anomalous subsidence (Fig. 2).

Topographic response to upwelling-remnant block interactions
In the third model type, we further quantify the surface response to upwelling-remnant block interactions (Figs. 2 vs. 6; Figs. S6, S10). Because the number of free parameters is large, we only limitedly tested the effects of the scale, thermal structure, and compositional density of the remnant block and the temperature and scale of the lower mantle thermal anomaly (Figs. S6, S10).

As shown in Fig. 6, the combined effects of mantle-upwelling-induced stretching and remnant-block sinking can induce >2 km topographic variations, with the subsidence center at a depth deeper than -1 km (Fig. 6a, b). This enhanced subsidence compared to the Type-2 models (Fig. S8) is because the lower mantle upwelling in this model accelerates the upper mantle convection velocity. Thereafter, the reduced viscosity of the upper mantle (the strain rate weakening effects, Equation S4) allows the remnant block to sink rapidly to cause a more prominent surface subsidence.

While the remnant block is cold and strong, it can divert the path of upwelling (Figs. 6c–d and 7 vs. Fig. 2b)11. Depending on the scale of the remnant block, the latter can either be pushed to be nearly vertical (Fig. 6e) or only mildly tilted by the lower mantle upwelling (Fig. 7c). Then, the upwelling reaches the LAB and uplifts the left shoulder of the subsidence region (Fig. 6d–g). Subsequently, the left shoulder keeps rising, while to the right, because this side is less affected by the upwelling, the topography remains low (Fig. 6e, f). As a result, the subsidence center shrinks to a ~200 km wide V-shape at the bottom (Fig. 6f).

Unlike the scenario shown in Fig. 3, the dynamic topography caused by the remnant block’s sinking, which induces downwelling in the upper mantle, enhances the initial subsidence magnitude (Fig. 6b). Later, the upper crust of the subsidence center experiences quick thinning after ~3 Myr (Fig. 6h). At the same time, the total crustal thickness gradually increases (Fig. 6d, i) as weak lower crustal materials flow from the compressive left shoulder to the extensive subsidence center (Fig. 6g)32. Because the lithosphere at the subsidence center is (comparatively) thinned more significantly than the crust (Fig. 6h) and is no longer affected by the mantle upwelling, its topography remains low during the following evolution of the model (Fig. 6f, i).

Instead, in models where the remnant block is large and weak (Fig. 7f), the plume can penetrate through it7,33. In this kind of model, the topography of the subsidence center can be partially reversed.

Methods
We use a finite element code written in MATLAB to conduct these 2D numerical experiments38,73–76. This code uses well-vectorized matrix-assembly techniques popularized by MILAMIN77. It includes a free surface on its top surface (a surface with no shear or normal stress)78 and can simulate mantle melting behavior that depends on the temperature, pressure, and melting degree (Fig. S2, Equations S10–S12)49,50. Because the largest density change in current models is <10%, the mantle should be well approximated by the Boussinesq approximation75. Therefore, we assume an incompressible viscoplastic rheology (Fig. 10c, e) that is governed by the conservation of mass, moment, and energy (Equations S1–S7). The phase transitions of mantle materials (olivine-to-wadsleyite and ringwoodite-to-bridgmanite mainly) are calculated from Equation S8–S9 (Fig. 10d)79,80. More details on our numerical method can be found in the supplements.

Model setup
In all runs, we utilized a free surface boundary condition78 as most of the models would be applied to stretched continental margins; in submarine regions of the marginal sea (e.g., SCS), the topographic relief should be increased to reflect the extra load associated with water infill. Other boundary sections are set to be free slip. The overriding lithosphere has varied thicknesses (except for tests in Fig. 5 and S10, the lithosphere is ~65 km thick) and contains a 30 km thick crust (20 km lower crust + 10 km upper crust) (Fig. 10a), and its initial temperature structure is controlled by 1-D steady-state conductive thermal profile and the mantle solidus (Figs. 10b and S2)31,81. In most models with lower mantle upwelling, the geometry of the hot anomaly in the lower mantle comprises a 400 km × 200 km square and a 400 km × 200 km equilateral triangle (Fig. 10b), which could cause a ~1000 km wide and ~1 km surface uplift after the plume contacts with the lithosphere (Fig. S5), mimicking the scale of hotspot swells82,83; the effect of the hot anomaly scale is also tested (Fig. S6). The hot center (orange color in Fig. 10a) is initially 1650 °C84 (except for the tests in Fig. S6), and the temperature gradually reduces to 1350 °C towards the edge (within the gray region in Fig. 10a). The temperature and material structures are reinforced every time step within the 400 km × 200 km square region to approximate the continuous heat supply from the core-mantle boundary (e.g., Shi and Morgan, 2024). In models with the remnant block, the block initially lies above the lower mantle (600 km deep), and its initial temperature follows the half-space cooling model (Fig. 10b)31. We further assume that the remnant block has stayed within the mantle for a varied period and, thus, is warmer than the prediction of the half-space cooling model (Fig. 10a, b). Besides, the surface and bottom temperatures remain at 0 and 1350 °C, respectively. Key model parameters are shown in Table S1.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Review of the Manuscript NCOMMS-25-20191, entitled Deciphering anomalous subsidence: lens effect of subducted blocks\non deep mantle upwelling, by Liu et al.\nThis manuscript investigates the dynamics and topographic effects of the interaction of a mantle upwelling with a sinking\nslab. The conceptual model is based on two case studies and is tested using 2D numerical models. The results show that\nsuch an interaction can lead to complex topographic patterns, including the presence of a topographic low (subsidence)\nabove a mantle upwelling, in the zone where it interacts with the sinking slab. This result is somewhat counterintuitive, but it\nis a natural outcome of the models.\nThe manuscript is very well written, well-illustrated, and it is a pleasure to read. The implications are well-extracted, and the\ncomparison with natural examples is informative. The limitations of the model are also well discussed and leave the door\nopen for further studies with more complex 3D geometries. The parametric sensitivity tests are also a major plus.\nI just would like to disclose that I am not a numerical modelling developer, and I would suggest having the manuscript seen\nby a numerical modelling expert. From my side, the models are intuitive and are doing what they are supposed to. I\nrecommend this manuscript for publication. The presented results will interest a broad range of geoscientists.\nJoão Duarte\nUniversity of Lisbon", "reviewer_2": "Overall, I think the topic is interesting, the approach undertaken by the authors is sound and\nthe structure of the paper is logical. The discussion is also of interest. However, I think there\nare a few weaknesses in the way they present their numerical results, and some of their choice\nof natural application cases are questionable. That being said, I believe that after major\nrevisions, the authors can provide a manuscript that is worth it publication in Nature\nCommunications.\nBelow I provide three major comments followed by other minor comments by line.\n1) Description of the numerical results\nI find that the results are not always correctly described and some physical explanation/interpretations are\nmissing.\nFor instance, L. 144, the authors describe a period of “prominent extensional period” before the plume arrival.\nBut there is no explanation of why this happens. Please provide seom more detailed description of this\nobservation.\nLikewise, L171: “the subducted block’s sinking enhances the initial subsidence magnitude”-> Why?\nFigures 3 and 4 are a bit hard to read given the scale, specially flow field. Would it be possible to have only the initial panel encompassing the whole model domain but later in time to focus only in regions of interest?\n(whole-domain figures could appear in the Supp Material).\nIn doing, so, I would strongly recommend to plot the deviatoric stress field within the\nlithosphere to help understanding the meaning of the average value (sub-panels g, see below)\nSome sub-panel descriptions need improvement. For instance, Fig. 3g is not clear to me what is the width over\nwhich the stress is averaged (lithosphere thickness L 131 or crust thickness L 143)? What does an average\nmeans if the sign change with depth?\nAlso, Fig 3i show an “elevation difference between the subsidence center and its left flank”\nbut I’m not sure what is the subsidence center. On sub-panels Fig3 d o f, it seems that there is\nno subsidence at all. Please explain what you mean.\n2) Model sensitivity to parameters\nThe authors ran approximately 800 simulations as said in L 118 (that is indeed a large number!) to span a\nrange of possible model parameter space. I think this is an important thing to report more clearly in the main\ntext.\nL148, the authors should find a way of summarizing their main findings over the sensitivity tests (same for L\n152) even if they do not show the results in the main text.\nMost importantly, I think some of the sensitivity tests should be moved into the main text, especially those\nthat could more directly affect the model results that they want to compare to natural observations. From my\npoint of view there are at least two modeling choices that have to be tested, shown in main text and discussed.\nThese are :\n• The lithospheric strength/rheology: This modeling assumption is never clearly explained.\nNotably, given that the strength/thickness of the lithosphere may impact the wavelength of\nthe topographic signal, this should be tested. Could a weak/thin lithosphere reduce the\nwavelength of topographic variations?\n• The length of the sinking block: It sounds like the narrowness of the subsidence observed in models\ncombining an upwelling and a block is controlled by the apparent width of the block\n(its length if it lies horizontally or its width if it lies vertically). What if the initial length is much\nlarger ? (as it is the case of the fast anomalies discussed for the natural cases)\n3) Relevance of models and/or the application cases\nTwo of the natural cases that are discussed by the authors : Western US and Northeast China, are cases where\nthe stagnant/sinking block is still attached to a subducting plate (Gorda plate and Pacific plate, respectively).\nThis is not the case in the models. For me, this is a strong weakness of the study and I think that the authors\nshould better justify either their model setup in regard to the regions that they want to apply the results to.\nOther comments per line number\nL27 (and in many other part): related to my main comment above, the use of “subducted block” is sometimes\nconfusing as the authors use it for both stagnant slab still attached to a surface plate, or to some isolated block\nlying in the mantle. May the use of “remnant block” could be more appropriate for the case where it is not\nlonger attached to an active subduction zone?\nL105: Does the MILAMIN code really has a true free-surface? If, so then the description in L339 of an “air\nlayer above the top surface” is confusing. If there is a sticky layer on top of the surface domain, is it instead a\nfree-slip boundary condition on top of the air layer? Please clarify.\nL111: Please be more explicit on the phase transition that you take into account. In general, I would said,\nthat the solved equations could/should go into the Method section.\nFigure3: Displaying the deviatoric stresses within the lithosphere should be considered. In sub-panel g) it is\nhard to understand what has been averaged. Also when you say deviatoric stress, I assume that it is the\nhorizontal deviatoric stresses? It is neverd said, I think.\nL141: When does the surface experience ~1 km topography variation? Please provide a little bit more details\nwhen you describe your results.\nL148: A brief summary of the parametric tests should be provided in the main text and not only referred to the\nSupp Material. In particular, emphasis should be put in this summary on the sensitivity of the wavelength and\nmagnitude of vertical motions.\nL154 and 163: Why does the topographic variation is so different in between the two. Please provide some\nphysical explanation. Is the sinking velocity of the block the same in the two type of models?\nL171: Why does the initial magnitude increases?\nL197: “direction peperndicular to the page” it’s rather a surprising wording. Please correct.", "reviewer_3": "The authors have explored the dynamics of slab-plume interactions using high-resolution 2D models of mantle convection.\nTheir results demonstrate the influence of proximal slabs in guiding plume-induced evolution of surface topography --\nspecifically, the observed subsidence patterns in two well-studied regions that cannot be adequately explained within the\ntraditional conceptual framework.\nThe results appear well-supported by detailed numerical models, within the specific set of model parameterisations\nemployed. I imagine this study to be well-received within the community and spur further studies exploring other modalities\nof complex slab-plume interactions and their influence on the evolution of surface topography.\nHowever, I have some concerns regarding the model setup (Fig E1). The base of the plume is represented as a virtual\nsource of heat/materials, hanging 'mid-air', fixed in space, some 650 km above the bottom of the computational domain. With\nthis approach, a key element of dynamic interaction with descending slabs is compromised, because the downward flux of\nmaterial is unable to push plume sources away from its path. The authors should add a new section to the manuscript\nexploring the implications of such alternate parameterisations, where the base of the plume is able to move more freely, as is\nthe case within the Earth's mantle. If producing additional models to explore such possibilities is prohibitive, authors\nshould include speculative, plausible scenarios for such and further soften the claims with regard to this conceptual framework\nbeing widely applicable in explaining anomalous subsidence patterns in other regions."}, "author_response": "We really appreciate your recognition of our efforts and the potential interest among\ngeoscientists raised in our study! With your encouragement, we will try to explore more\npossibilities of this scenario in the future, including 3D modeling approaches applied to more\ncomplex geometries.\nWe really appreciate your recognition of our efforts in this study! We acknowledge your\nconstructive suggestions for improving our manuscript, and we have carefully considered all\nof them and modified the text accordingly. Following your suggestions, we have tried to\nrestructure the manuscript by merging the online method into the main text and moving the\nsensitivity test summaries regarding the mantle phase transition and lithospheric thickness\ninto the main text, too. By making these modifications, our manuscript should better match\nthe style preferred by Nature Communications.\n1) Description of the numerical results\nFor instance, L. 144, the authors describe a period of “prominent extensional period” before\nthe plume arrival. But there is no explanation of why this happens. Please provide some more\ndetailed description of this observation.\nThanks. We have added the following explanation: This extensional episode is mainly due to\nthe basal drag of divergent mantle flow caused by the lower mantle upwelling (Figs. 4c).\nLikewise, L171: “the subducted block’s sinking enhances the initial subsidence magnitude”->\nWhy?\nThanks. Because there is additional subsidence caused by the block sinking (or dynamic\ntopography for simplicity), the superposed effects of block sinking and lower mantle upwelling\ncaused more subsidence.\nFigures 3 and 4 are a bit hard to read given the scale, specially flow field. Would it be possible\nto have only the initial panel encompassing the whole model domain but later in time to focus\nonly in regions of interest? (whole-domain figures could appear in the Supp Material).\nIn doing, so, I would strongly recommend to plot the deviatoric stress field within the\nlithosphere to help understanding the meaning of the average value (sub-panels g, see below)\nWe apologize for the confusing figures. In the revision, we have changed the figure format\npartially based on the suggestions of the reviewer, hoping this helps to make the figure more\nstraightforward to interpret. Please see below:\nFig. 3 New. (now it is Fig. 4)\nFig. 4 New. (now it is Fig. 7)\nSome sub-panel descriptions need improvement. For instance, Fig. 3g is not clear to me what\nis the width over which the stress is averaged (lithosphere thickness L 131 or crust thickness L\n143)? What does an average means if the sign change with depth?\nWe apologize for this too-crude description. The horizontal deviatoric stress with the largest\nmagnitude within the lithosphere is sampled every 5 km to reveal the dominant lithospheric\ndeformation state.\nAlso, Fig 3i show an “elevation difference between the subsidence center and its left flank”\nbut I’m not sure what is the subsidence center. On sub-panels Fig3 d o f, it seems that there is\nno subsidence at all. Please explain what you mean.\nWe apologize for this too-crude description. We have changed the description to be:\nComparisons of the elevation between the model center (horizontal axis: 1750 km-2250 km)\nand its left flank (horizontal axis: 1250 km-1750 km).\n2) Model sensitivity to parameters\nThe authors ran approximately 800 simulations as said in L 118 (that is indeed a large number!)\nto span a range of possible model parameter space. I think this is an important thing to report\nmore clearly in the main text.\nThank you for recognizing our efforts!\nL148, the authors should find a way of summarizing their main findings over the sensitivity\ntests (same for L152) even if they do not show the results in the main text.\nMost importantly, I think some of the sensitivity tests should be moved into the main text,\nespecially those that could more directly affect the model results that they want to compare\nto natural observations. From my point of view there are at least two modeling choices that\nhave to be tested, shown in main text and discussed.\nBased on this suggestion, we have moved the sensitivity test summary that was in the method\nsection into the main text.\nThese are:\n• The lithospheric strength/rheology: This modeling assumption is never clearly explained.\nNotably, given that the strength/thickness of the lithosphere may impact the wavelength of\nthe topographic signal, this should be tested. Could a weak/thin lithosphere reduce the\nwavelength of topographic variations?\nBased on this comment, we have added a comparison of the topographic variations based on\nthe results summarized in Figure S6 of the initial supplements (now it is Fig. 6 in the main text).\nAccording to Fig. 6e, we find the lithosphere thickness does not clearly affect the subsidence\nwavelength, except when the subsidence is weak in the scenarios for >=100 km thick\nlithosphere.\n• The length of the sinking block: It sounds like the narrowness of the subsidence observed in\nmodels combining an upwelling and a block is controlled by the apparent width of the block\n(its length if it lies horizontally or its width if it lies vertically). What if the initial length is much\nlarger ? (as it is the case of the fast anomalies discussed for the natural cases)\nThanks for raising this point. If the slab is larger than the models in the first-version manuscript,\nthe path of the mantle upwelling would be clearly shifted (see the figure below). Figure S11\nwas in the supplements, and in the revision, we have moved it into the main text to form parts\nof Fig. 8. We also added the following discussion to explain our point of view towards a\nscenario with an even larger slab in the mid-mantle: Depending on the scale of the remnant\nblock, the latter can either be pushed to be nearly vertical (Fig. 7e) or only mildly tilted by the\nlower mantle upwelling (Fig. 8c). Please see more of our replies to this question below, where\nwe think we have better clarified the scenarios we wished to study here.\nFigure 8. Result of typical plume-remnant block interaction models (a-c: Run 821; d-f: Run 813).\n(a)-(f) Snapshots of model evolution. The remnant block is 1000 km wide; the seafloor age is\n60 and 20 Myr in Run 821 and 813, respectively, and the heating time is 20 Myr (Fig. S10).\nOther plotting habits are the same as in Fig. 4.3) Relevance of models and/or the application\ncases\nTwo of the natural cases that are discussed by the authors : Western US and Northeast China,\nare cases where the stagnant/sinking block is still attached to a subducting plate (Gorda plate\nand Pacific plate, respectively). This is not the case in the models. For me, this is a strong\nweakness of the study and I think that the authors should better justify either their model\nsetup in regard to the regions that they want to apply the results to.\nWe accept the criticism, and we think this might be a good place to better clarify the term\n‘subducted block’ that we had been using, so helping us to improve the presentation.\nGenerally, the stagnant/sinking block appears still attached to a subducting plate in both cases.\nHowever, based on more localized studies of these mantle structures, the connection between\nthe stagnant and subducting portions is seen to be locally broken. There are some low-wave-\nspeed ‘holes’/channels (or seismic gaps) within the subducting plate (see the figures below\nfrom recently published papers). These holes offer an opportunity for the deep mantle\nupwelling to rise, affecting the upper mantle convection and sometimes causing localized\nsubsidence (as our sensitivity tests demonstrate. This has been proposed to happen either\nwhen the MTZ is anomalously wet, which is proposed to need a long-lived nearby subduction\nzone in any case) and intraplate magmatism (e.g., (Camp et al., 2015; Camp et al., 2017; Q.\nZhou et al., 2018; Zhou et al., 2025). In addition, we still reserve our opinion that some of the\nhigh-wave-speed anomalies within both areas’ MTZ could contain the delaminated cratonic\nroot, given that significant cratonic lithosphere thinning has happened there (Liu et al., 2016;\n2018a, b; 2019). In this possibility, the remnant blocks do not necessarily need to be connected\nto a subducting plate at all (we admit that our opinion is highly debatable; we are open to\nremoving the related texts if they appear to be offensive).\nAll in all, as presented by these studies (see below), the stagnant fragments are generally\nsmaller than 1000 km in at least one axis, covered by our tested scenarios. However, following\nthe above point raised by the reviewer, if the slab is too broad and complete and so blocks\nfocused lower mantle upwelling (e.g., in regions away from the subsidence locus in both\nregions), the scenarios we discussed here would not be appropriate to allow any surface\nsubsidence. We have added more discussion regarding this potential weakness to the main\ntext:\nIn regions like the western United States and northeastern China, previous studies have\nshown that high-wave-speed anomalies within the mantle transition zone are typically\nassociated with actively subducting plates81,82. If this observation is consistently accurate\nacross these regions, the presence of a subducting plate would overshadow the lower mantle\nupwelling, leading to effects on upper mantle convection that differ from our model\npredictions (Fig. 7)83. However, recent localized seismological studies indicate that some, if not\nall, high-wave-speed anomalies in these areas are separated from the subducting plate by a\nseismic gap84 or a low-wave-speed ‘hole’ or channel11,68. Given that lithospheric delamination\nevents could have happened in these regions, some high-wave-speed anomalies can be\ndelaminated cratonic roots instead32,33,66. Indeed, the subducting plate itself adjacent to the\ntarget anomalies seems to have locally penetrated into the lower mantle independently12,85,\nparticularly beneath regions experiencing anomalous subsidence, where some lower mantle\nupwelling may occur (Figs. 9 vs. 1 & 10). From this perspective, the simultaneous occurrence\nof lower mantle upwelling, remnant blocks, and localized anomalous subsidence could be\ndynamically connected by our proposed scenario (Fig. 9b). Nonetheless, we emphasize that\nour model cannot be generalized to other areas within these regions that lie outside the\nseismic gap. Further numerical modeling incorporating more complex initial geometries than\nthose used here would be necessary.\nOther comments per line number\nL27 (and in many other part): related to my main comment above, the use of “subducted block”\nis sometimes confusing as the authors use it for both stagnant slab still attached to a surface\nplate, or to some isolated block lying in the mantle. May the use of “remnant block” could be\nmore appropriate for the case where it is not longer attached to an active subduction zone?\nThanks for pointing out this potential confusion. As we explained above, we did not intend to\nstudy a more complex scenario where the block is still connected to a subducting portion. We\nagree that remnant block is a better word (thank you for this suggestion!), and we have\nreplaced ‘subducted block’ with ‘remnant block’ throughout the text and the title, with its\ndefinition given in the abstract and the second paragraph: subducted/delaminated blocks not\ncurrently connected to a subducting/delaminating lithosphere.\nL105: Does the MILAMIN code really has a true free-surface? If, so then the description in L339\nof an “air layer above the top surface” is confusing. If there is a sticky layer on top of the surface\ndomain, is it instead a free-slip boundary condition on top of the air layer? Please clarify.\nOur version of the code does have a free surface (Andrés-Martínez et al., 2015). We have\nclarified the text as follows: In all runs, we utilized a free surface boundary condition(Andrés-\nMartínez et al., 2015) as most of the models would be applied to stretched continental margins\nL111: Please be more explicit on the phase transition that you take into account. In general, I\nwould said, that the solved equations could/should go into the Method section.\nThanks for the reminder. We have clarified the text: The phase transitions of mantle materials\n(olivine-to-wadsleyite and ringwoodite-to-bridgmanite mainly) are calculated from Equation\nS8-S9.\nFigure3: Displaying the deviatoric stresses within the lithosphere should be considered. In sub-\npanel g) it is hard to understand what has been averaged. Also when you say deviatoric stress,\nI assume that it is the horizontal deviatoric stresses? It is neverd said, I think.\nWe apologize for this crude description. We have added the following description: The\nhorizontal deviatoric stress with the largest magnitude within the lithosphere is sampled every\n5 km horizontally to reveal the dominant lithospheric deformation state.\nL141: When does the surface experience ~1 km topography variation? Please provide a little\nbit more details when you describe your results.\nWe apologize for this too-crude description. The largest topography variation happens at ~2.5-\nL148: A brief summary of the parametric tests should be provided in the main text and not\nonly referred to the Supp Material. In particular, emphasis should be put in this summary on\nthe sensitivity of the wavelength and magnitude of vertical motions.\nBased on this comment, we have added the following words to the main text: the magnitude\nof topography uplift/subsidence increases with the decrease of the lithosphere thickness (Figs.\n6a-6d) while the wavelength of topographical variation does not change obviously, which is\nconsistent with previous 3D modeling results (Wang and Li, 2021) (Fig. 6e).\nL154 and 163: Why does the topographic variation is so different in between the two. Please\nprovide some physical explanation. Is the sinking velocity of the block the same in the two\ntype of models?\nThe sinking is not completely the same in both models. Because we used a non-linear rheology,\nthe lower mantle upwelling would drive a fast mantle convection in the upper mantle,\ndecreasing its rheology and leading to a fast sinking of the remnant block. Based on your\ncomments, we have added the following explanation to the main text: This enhanced\nsubsidence compared to the Type-2 models (Fig. S8) is mainly because the lower mantle\nupwelling in this model accelerates the upper mantle convection velocity. Thereafter, the\nreduced viscosity of the upper mantle (the strain rate weakening effects, Equation S4) allows\nthe remnant block to sink rapidly to cause a more prominent surface subsidence.\nL171: Why does the initial magnitude increases?\nBecause there is additional subsidence caused by the block sinking (or dynamic topography\nfor simplicity), the superposed effects of block sinking and lower mantle upwelling caused\nmore subsidence.\nL197: “direction peperndicular to the page” it’s rather a surprising wording. Please correct.\nWe have changed the text to a direction perpendicular to the model slice.\nThank you very much for recognizing our start-up efforts towards some possibilities of slab-\nmantle-upwelling!\nThis is an important potential shortcoming of our presented results. Thank you for pointing\nthis out!\nAs mentioned by the reviewers, adding a free-moving mantle upwelling is not technically\npossible in our current method — it would be somewhat easier to add an ad-hoc moving\nplume source within a 3D model. To address this weakness, we have added the following\ndiscussion to the main text:\nModel Limitations Awaiting Future Investigation\nIn some of the models we present here, we assume that a lower mantle upwelling\nencounters a previously subducted/delaminated block that is not currently connected to an\nactively subducting/delaminating plate (Fig. 7). Our idealized model setups are mainly\ndesigned to explore the potential topographical response to interactions between deep\nmantle upwelling and remnant blocks, without accounting for additional interfering factors\nsuch as nearby active subduction or the stability of the lower mantle upwelling.\nIn regions like the western United States and northeastern China, previous studies have\nshown that high-wave-speed anomalies within the mantle transition zone are typically\nassociated with actively subducting plates81,82. If this observation is consistently accurate\nacross these regions, the presence of a subducting plate would overshadow the lower mantle\nupwelling, leading to effects on upper mantle convection that differ from our model\npredictions (Fig. 7)83. However, recent localized seismological studies indicate that some, if not\nall, high-wave-speed anomalies in these areas are separated from the subducting plate by a\nseismic gap84 or a low-wave-speed ‘hole’ or channels11,68. Given that lithospheric delamination\nevents could have happened in these regions, some high-wave-speed anomalies can be\ndelaminated cratonic roots instead32,33,66. Indeed, the subducting plate itself adjacent to the\ntarget anomalies seems to have locally penetrated into the lower mantle independently12,85,\nparticularly beneath regions experiencing anomalous subsidence, where some lower mantle\nupwelling may occur (Figs. 9 vs. 1 & 10). From this perspective, the simultaneous occurrence\nof lower mantle upwelling, remnant blocks, and localized anomalous subsidence could be\ndynamically connected by our proposed scenario (Fig. 9b). Nonetheless, we emphasize that\nour model cannot be generalized to other areas within these regions that lie outside the\nseismic gap. Further numerical modeling incorporating more complex initial geometries than\nthose used here would be necessary.\nFurthermore, our idealized model presupposes there exists an active lower upper mantle\nupwelling with a path in the lower mantle that remains relatively unchanged, similar to\nprevious studies41. In nature, such vertical upwellings may be more transient. Nearby\nsubducted slabs may induce changes in the lower mantle convection field and/or the position\nof thermal or chemical piles at the core-mantle boundary, potentially distorting the plume\npath over time53,86. This represents another potential complexity that our current modeling\napproach cannot readily address. However, given that our modeling timeframe, particularly\nfor the occurrence of abnormal subsidence, spans only several million years (Figs. 7-8), the\nsimplified model setup here may still be valid for a preliminary exploration of surface\nresponses to deep mantle processes (Fig. 9b)."}, {"round": 2, "reviewer_comments": {"reviewer_2": "Dear editor, dear authors,\nI think that all my main concerns have been addressed by the authors. I just find the numbers (topography change) on the\nnew Figure 5 difficult to read. Maybe the authors could consider making it bigger (maybe in 2 cols 3 rows?) for publication.\nApart from that, I have no more comments and I think the manuscript is worth publication as is.\nBest wishes"}, "author_response": "We really appreciate your constructive suggestions for improving our manuscript, and we have\nremade the figure based on your suggestions."}]

Decision letter