Flat Subduction Versus Big Mantle Wedge: Contrasting Modes for Deep Hydration and Overriding Craton Modification

Subduction‐induced deep hydration and water cycling may strongly control the modification of overriding cratonic lithosphere. Two contrasting modes are generally proposed: (1) flat subduction (FS) regime with slab subducting subhorizontally beneath the overriding lithosphere and (2) big mantle wedge (BMW) regime with slab flattening at the bottom of mantle transition zone. Here, systematic numerical models are built to study the subduction‐induced deep hydration processes in the contrasting FS and BMW regimes as well as their effects on the modification of overriding lithosphere. The model results indicate that the dehydration process in the FS regime can significantly modify the overriding lithosphere for a region of about 600 km from the trench. During the progressive flat subduction, the partial melting and magmatism migrate toward the inner land of the overriding plate, which will be reversed and backward to the trench during the transition from flat to steep slab subduction. On the other hand, the deep hydration in the BMW regime is strongly dependent on the subcrustal serpentinite layer in the subducting slab, whereas the oceanic crust cannot carry water to the transition zone. The modification of the overriding lithosphere in the BMW regime occurs in a larger region of >1,000 km from trench, which is, however, generally slower and weaker. Thus, the flat subduction of Izanagi plate may play a more significant role in the modification of North China Craton in the late Jurassic to early Cretaceous, which could be accompanied by the effects of deep water cycling in the BMW regime.


Introduction
The water transportation from the surface to the deep mantle and its further circulation in the interior of the Earth play an important role in the evolution of the planet (e.g., Faccenda, 2014;Magni et al., 2014;Nakagawa & Nakakuki, 2019). There are several ways for the water to be transported upward, for example, through mid-ocean ridge (MOR), mantle plume, and island arc; however, the subducting slab is the only path to transport water downward to the mantle. The subduction-induced water cycling strongly affects various subduction-zone processes, for example, the widely studied arc magmatism (Schmidt & Poli, 1998) and intraslab earthquakes (Yamasaki & Seno, 2003), as well as the plausible overriding craton modification/destruction (Zhu et al., 2012).
Two subduction modes are generally proposed (Wu et al., 2019;Zhu & Xu, 2019), that is, flat subduction and big mantle wedge (Figure 1), which may result in contrasting deep hydration processes and thus play different roles in the overriding lithospheric modification. The "flat subduction" model indicates slab subducting subhorizontally beneath the overriding cratonic lithosphere, like the present-day flat slab beneath South America (Espurt et al., 2008). In this regime, the flat subducting slab carries and liberates water beneath the overriding lithosphere, which may result in the intense hydration and significant weakening of the cratonic root and further contribute to its modification ( Figure 1a). In addition, the mechanical bottom erosion of the overriding lithosphere by the flat slab may also play a certain role (Axen et al., 2018). On the other hand, the "big mantle wedge" model concerns the subducted and flatly stagnant slab in the mantle transition zone (MTZ). A certain amount of water may be carried by the sinking slab to the MTZ and could be released later on as the slab is being heated by the surrounding hot materials. The upward migration of the released water and its further interaction with the overriding lithosphere may lead to the craton modification. The key to compare and distinguish these two models is to constrain the subduction-induced deep hydration processes in the whole upper mantle, which is, however, far from better understanding (Liu & Li, 2018).
The oceanic crust can contain a certain amount of water in the hydrous minerals, which could be carried into the mantle by subduction and then liberate during the heating of sinking slab. However, the oceanic crust can only carry water to shallower depth (Maruyama & Okamoto, 2007), which thus prevent the deep hydration processes in the MTZ Nakao et al., 2016). Alternatively, the hydrated mantle rocks beneath the oceanic crust could be colder and then carry water to larger depths beyond the subarc (Hacker, 2008;Rüpke et al., 2004), which thus may contribute to the deep hydration and water cycling (Van Keken et al., 2011).
In this study, systematic thermomechanical models are conducted with the thermodynamic fluid-melt activity to investigate the subduction-induced deep water cycling as well as the contributions to the overriding craton modification. In particular, the deep hydration processes in the flat subduction and the big mantle wedge regimes are systematically compared, which shed new light on the modification/destruction of North China Craton (NCC) in the Mesozoic.

General Model Configuration and Boundary Conditions
Large-scale models are configured in a 2-D spatial domain of 5,000 × 1,400 km ( Figure 2). A 60-Ma-old oceanic lithosphere is set on the left part of the model domain, which is composed of an upper crustal layer (3 km), a lower crustal layer (5 km), and a lithospheric mantle layer. The thermal structure and whole thickness of the oceanic lithosphere are constrained by the infinite half-space cooling model, with the detailed calculations shown in Turcotte and Schubert (2002). A continental lithosphere is set on the right part of the model domain, which includes an upper felsic crust (20 km), a lower mafic crust (15 km), and an ultramafic lithospheric mantle layer (100 or 150 km). For the temperature configuration of the continental lithosphere, a linear gradient is applied with constant temperatures of 0°C and 1350°C at the top and bottom boundaries of the lithosphere, respectively. An initial weak zone is set between the oceanic and continental plates, The zoomed-in subduction zone showing clearly the transition between the oceanic and continental plates. The white lines demonstrate the isotherms with specific temperatures shown in the figure. The colors of rock types as specified in (c): 1 and 2, sticky air/water; 3 and 4, sediment during the surface processes; 5, sediment with partial melting; 6 and 7, solid and partially molten continental upper crust, respectively; 8 and 9, solid and partially molten continental lower crust, respectively; 10, oceanic upper crust; 11, oceanic lower crust; 12, oceanic crust with partial melting; 13, lithospheric mantle; 14, sublithospheric mantle; 15, hydrated mantle; 16, hydrated mantle with serpentinization; and 17, mantle with partial melting. (d and e) The density and effective viscosity structures of the model domain, which are validated in Li et al. (2019). used for localizing the initial subduction. The dynamics of subduction initiation at oceanic-continental plate boundary is not the focus of this study but has been systematically investigated in Li (2019, 2020). The initial temperature of the sublithospheric mantle is set with a constant adiabatic gradient of 0.5°C/km. A "sticky air" layer is set between the crustal surface and the top model boundary, which has a low density and viscosity (Schmeling et al., 2008). The algorithm of sticky air is systematically calibrated in Crameri et al. (2012), which can be used to mimic a free surface boundary condition. Detailed numerical parameters are shown in the supporting information (Tables S1 and S2).
For the velocity boundary conditions, free slip is satisfied for all boundaries. The convergence velocity is applied during the initial stages of model convergence ( Figure 2a) and removed after 20 Myr, leaving the subduction driven purely by the internal buoyancy. For the thermal boundary condition, fixed values of 0°C and 1975°C are applied for the top and bottom boundaries, respectively. The horizontal heat flux across the vertical boundaries is 0.

Initial Water Content
The water contents of the upper and lower oceanic crusts are constrained by the pressure-and temperature-dependent phase diagrams of basalt and gabbro, respectively ( Figure 3). In the initial model, the oceanic crust has low pressure and temperature conditions, which thus indicates high water contents of 2.78 wt% for the upper crust ( Figure 3b) and 1.47 wt% for the lower crust ( Figure 3c).
The water content of the mantle is zero in the initial model, except the topmost lithospheric mantle of oceanic plate that could be hydrated before subduction (Deschamps et al., 2013;Evans et al., 2013), the processes of which are not directly simulated here. Previous studies indicate that the subcrustal hydration of oceanic lithosphere may occur (i) in the MOR region, where the upwelling hot mantle materials are interacting with the ocean water (Sauter et al., 2013); (ii) around the scarps and transform fault regions, where the ocean water may penetrate downward to the subcrustal depth (Bideau et al., 1991;Morishita et al., 2009);and (iii) in the outer-rise region of the subducting slab close to the trench, where the ocean water could migrate along the bending-related faults into the topmost lithospheric mantle beneath the oceanic crust (Faccenda et al., 2009;Key et al., 2012;Ranero et al., 2003). The properties of the subcrustal hydrous layer of oceanic lithosphere are widely investigated with multiple geophysical observations ( Table 1). The results indicate that the thickness (H serp ) and water content (W serp ) of this hydrous mantle layer are quite variable among different subduction zones in nature, that is, H serp ∈ (0, 30) km and W serp ∈ (0, 4) wt% as compiled in Table 1. According to these observations, a subcrustal serpentinite layer with variable thickness (0-25 km) is applied for the oceanic lithosphere in this study, with the water content constrained by the phase diagram of hydrated peridotite ( Figure 3d). In the initial model with low pressure and temperature conditions, the subcrustal serpentinite layer has a high water content of 2.0 wt% ( Figure 3d).

Model Limitations
The current numerical model includes the most critical processes of hydration, partial melting, and multiple phase transitions . There are still some uncertainties and limitations. One uncertainty of the deep hydration model is the water capacity of nominally anhydrous minerals (NAMs) of the mantle, which is not included in the thermodynamic database ( Figure 3). It is also not well constrained by the laboratory experiments, some of which have shown the water capacity of the mantle NAMs is generally <0.1-0.2 wt%, whereas the minerals (primarily Wadsleyite and Ringwoodite) in the MTZ may contain more water, for example, on the order of 1.0 wt%, or even as much as 3.0 wt% (Bercovici & Karato, 2003;Bolfan-Casanova et al., 2000). It is worth noting that the Wadsleyite and Ringwoodite still belong to NAMs, although their water capacity could be even higher than hydrous minerals. Another large uncertainty lies in the preexisting water condition of the MTZ. Nakagawa (2017) suggested that the MTZ could be nearly water-saturated during the long timescale mantle convection, which means that the MTZ minerals cannot absorb more water. These previous studies indicate that either the actual water capacity or the preexisting water content of the NAMs in the MTZ is not well constrained for the natural Earth. For simplicity, a reference water capacity of 0.1 wt% and no preexisting water content are generally applied for the mantle NAMs in the current model. Their effects are further tested with additional models.
Another limitation lies in the partial melting of mantle rocks, which is only applied in the regions above 300-km depth, according to the parameterization of Katz et al. (2003). Thereby, the partial melting of mantle rocks is neglected below 300-km depth. It thus prevents the test of partial melting or the water filter model of the MTZ (Bercovici & Karato, 2003), which requires further studies.

General Model Without Initial Serpentinite Layer in the Oceanic Plate
In the general model setup (Figure 2), the initial water is only present in the oceanic crust which is constrained by the thermodynamic database (Figure 3), whereas the lithospheric mantle is absent of water ( Figure 4a′).
The converging plates are initially pushed with a fixed velocity of 6 cm/year (Figure 4e), which includes a subducting plate velocity of 4 cm/year and an overriding plate velocity of −2 cm/year. During this stage with constant pushing velocity, the oceanic slab subducts and stagnates in the MTZ (Figure 4b). At 20 Myr from the beginning of the model, the prescribed convergent velocities are removed, after which the purely buoyancy-driven free subduction is applied. Accompanying with trench retreat, the slab stagnates and flattens at the bottom of the MTZ for >1,000 km ( Figure 4c). Finally, the steep subduction with negligible trench retreat ( Figure 4e) results in the foundering of stagnant slab into the topmost lower mantle (Figure 4d). During the whole subduction processes, most water of the sinking slab is lost in the shallower subarc depth, with some water carried to a deeper depth of 200-300 km (Figures 4b′-4d′). However, the water taken to the MTZ is very limited (≤0.1 wt% in the present study), because all the hydrous phases broke down and the residual water is only present in the mantle NAMs, which have a low water capacity of 0.1 wt%. Figure 3. Computed water capacity for four rock types in the subduction model, that is, sediment (a), basalt for oceanic upper crust (b), gabbro for oceanic lower crust (c), and the hydrated peridotite for mantle (d). The representative compositions of these four rock types are, respectively, according to the global subducting sediment (GLOSS) (Plank & Langmuir, 1998), the mean composition of global ocean ridge basalts (Gale et al., 2013), the synthetic composition of oceanic lower crustal gabbro (Behn & Kelemen, 2003), and the LOSIMAG-C1 composition of peridotite (Hart & Zindler, 1986). The phase diagram sections are computed by Perple_X using Gibbs free energy minimization (Connolly, 2005). The maximum water content for each rock type is limited to account for initial incomplete or heterogeneous hydration (Gerya et al., 2006): (a) sediment ¼ 7.60 wt%; (b) oceanic upper crust ¼ 2.78 wt%; (c) oceanic lower crust ¼ 1.47 wt%; and (d) hydrated mantle ¼ 2.00 wt%.

Reference Model With an Initial Serpentinite Layer of 20 km Thick
In this model, an initial serpentinite layer is prescribed beneath the oceanic crust, which is rheologically weak and has a thickness of 20 km ( Figure 5). All the other parameters remain the same as the general model ( Figure 4).
At the initial stage of the model driven by the constant pushing velocity, the oceanic slab is subhorizontally subducted beneath the overriding continental lithosphere (Figure 5a). The water is carried by the slab to the inner continent within about 600 km from the trench, during which a certain quantity of water is liberated and migrates upward to hydrate the continental lithosphere ( Figure 5a′). The flat subduction region is characterized by cold thermal condition (Figure 5a), which thus prevents the partial melting of the hydrous mantle rocks above the flatly underthrusting slab. The partially molten rocks (red color) are only predicted in the leading end of the flat subduction region where relatively high temperature is present due to the contact with hot asthenosphere (Figure 5a). The increasing of slab pull leads to gradual transition from flat to steep subduction (Figures 5b and 5c). After 20 Myr, the prescribed convergent velocity is removed, leading to further slab rollback and long slab stagnation at the bottom of MTZ (Figure 5d). The previously hydrated mantle rocks above the flat slab will be heated by the incoming asthenospheric flow, which leads to the lateral migration of partial melting toward the trench (Figures 5a-5d). On the other hand, significant water is carried mainly by the subcrustal hydrated mantle rocks to the MTZ (Figures 5b′-5d′).
In the big mantle wedge regime, the stagnant and flatten slab at the bottom of MTZ is conductively heated by the surrounding hot materials, as shown by the evolution of isothermal contours (Figures 5d and 5e). The increasing temperature of the slab leads to the decomposition of high-pressure hydrosilicates. The released water migrates upward, hydrates the upper mantle as well as the overriding lithosphere (Figures 5d′ and 5e′). The hydration process leads to vigorous partial melting around the bottom of the continental lithosphere (Figures 5e and 5e′).
The complex processes of melt extraction and further emplacement of the magmatic rocks are not well constrained and are, thus, not directly simulated in the current study. In order to give implications for the potential magmatic activity, the spatiotemporal distributions of partial melting are plotted (Figures 6a and 6b). It shows that the partial melting migrates far away from the trench during the progressive flat slab subduction, which goes instead toward the trench during the transition from flat to steep subduction, that is, flat slab rollback. At the late stages from about 50 Myr (Figures 6a and 6b), the partial melting occurs at around the bottom of the overriding lithosphere, locating far from the trench (i.e., >1,000 km), due to the hydration in the big mantle wedge regime.  (2015) Mariana 150 6.5 ± 1.5 24 ± 5 2 Cai et al. (2018) a H serp defines the thickness of hydrated lithospheric mantle layer in the subducting slab. b W serp is the estimated water content of the hydrated mantle layer.
The overriding continental lithosphere is significantly modified (i.e., hydrated and/or partially molten) by the fluid and melt activities in both the flat subduction and the big mantle wedge regimes (Figures 6c and  6d). The flat subduction can affect the regions of about 600 km from the trench, where the cratonic lithosphere is significantly modified. The thickness of unmodified lithospheric mantle is changing from the original 100 km to ≤35 km after 60 Myr (Figure 6c). Consequently, the unmodified continental lithosphere (including crust) is changing from 135 to 35-70 km (Figure 6d). The flattened slab at the bottom of the MTZ, that is, the big mantle wedge regime, can affect larger regions of >1,000 km from the trench; however, the degree of hydration and related cratonic modification is lower (Figures 6c and 6d).

Effect of the Thickness of Initial Serpentinite Layer
The thickness of the serpentinized mantle in the subducting slab is changing among different subduction zones in nature (Table 1), the effects of which are studied in this section (Figure 7). In the presented models, the serpentinite layer has no effect on the buoyancy, but it does have an effect (weakening) on the viscosity. The model results indicate that flat subduction is not predicted if the serpentinite layer is thin, for example, H serp ¼ 5 km (Figures 7a and 7a′) or H serp ¼ 10 km (Figures 7b and 7b′). The steeply subducting slab flattens in the MTZ, with forming a big mantle wedge. The quantity of water carried by the subducting slab to the MTZ is limited, which only affects the neighboring regions of the stagnant slab but does not contribute to the modification of overriding continental lithosphere (Figures 7a and 7b).
In contrast, if the initial serpentinite layer is thicker, for example, H serp ¼ 15 or 25 km (Figures 7c and 7d), the evolution of the model is similar to the reference model with H serp ¼ 20 km ( Figure 5). The flat subduction is formed at the initial stage, which gradually changes to the steep subduction due to the increasing of slab pull. Then the slab stagnates and flattens in the MTZ with forming a big mantle wedge. A large

Journal of Geophysical Research: Solid Earth
LI quantity of water is carried by the sinking slab to the MTZ, which liberates later and significantly contributes to the hydration of the upper mantle and the modification of the overriding continental lithosphere (Figures 7c and 7d).
The subduction-induced hydration and modification of the overriding lithosphere is further compared in Figure 8. In the models with thinner initial serpentinite layer, the overriding lithospheric modification is restricted to a narrower region of 200-300 km with H serp ¼ 5 km (Figure 8a) and 300-400 km with H serp ¼ 10 km (Figure 8b). In contrast, a wider region of about 600 km in the overriding lithosphere is modified by the flat subduction-induced hydration in the models with thicker initial serpentinite layers (Figures 8c and 8d). In addition, the hydration in the big mantle wedge regime may lead to an additional overriding lithospheric modification of about 600 km (e.g., Figure 8d).

Effect of a Thick Overriding Continental Lithosphere
In the previous models, the overriding continental lithosphere has a normal thickness of 135 km, which may be thinner than the stable craton (Peslier et al., 2010;Sleep, 2005). In this section, an additional model with a thick overriding plate of 185 km is further conducted (Figure 9). All the other parameters keep the same as the reference model ( Figure 5). The general model evolution with a thick overriding plate is similar to the reference model (cf. Figures 9 and  5). The flat slab subduction is resulting at the initial stages, which changes to steep subduction style and forms a big mantle wedge. The dehydration occurs during either the flat subduction or the slab stagnation in the MTZ. The length scale of the flat subduction is similarly around 600 km from the trench, with intense hydration and modification of the overriding lithosphere (Figures 9d and 9d′). In addition, the dehydration from stagnant slab in the MTZ may influence a further region from the trench, although the effects are much weaker (Figures 9d and 9d′).
The different phenomenon in this model with a thick overriding lithosphere is the limited partial melting in either the flat subduction or the big mantle wedge regimes (cf. Figures 9 and 5). This is due to the increased solidus temperature according to the larger pressure at the bottom of the thick overriding lithosphere. Thus, the conditions of partial melting are more difficult to be achieved with the similar degree of hydration. The partial melting in this model could be promoted by increasing the water capacity of the mantle NAMs ( Figure S1).

End-Member Models With High or Low Water Capacity of MTZ Rocks
In the above models, the potential high water capacity and high preexisting water content of the rocks in the MTZ are neglected, with applying a similar low capacity of 0.1 wt% and no initial water content. It thus can be considered as a 0.1 wt% vacancy for the MTZ rocks to be saturated. In order to better understand its effects, two additional end-member models are further conducted: (1) The MTZ NAMs have a high water

Journal of Geophysical Research: Solid Earth
LI capacity of 1.0 wt% but still no preexisting water content, which thus need to absorb a large amount of water to be saturated; and (2) the entire mantle NAMs, including the MTZ ones, have zero water capacity and again no preexisting water content, which is like a saturated condition for the mantle NAMs.
In the first model with a high water capacity of 1.0 wt% for the mantle NAMs in the MTZ (Figures 10a and  10a′), the result shows that the water liberated from the slab in the MTZ is totally absorbed by the neighboring mantle rocks which have a high water capacity. Thus, it differs significantly from the reference model with upward migration of the hydration front (e.g., Figure 5). This model represents an end-member regime with an initially dry MTZ, which instead has a high water capacity. Not surprisingly, the limited water carried by a single subducting slab can only be used to feed the "thirsty" rocks in the MTZ.
In another model with zero water capacity of NAMs in the entire mantle (Figures 10b and 10b′), the water liberated from the stagnant slab in the MTZ migrates upward and is finally absorbed by the cold lithospheric mantle, which has a large water capacity. Instead, the sublithospheric mantle rocks have negligible ability to store water because of the zero water capacity of NAMs and the high temperature-induced breaking down of all the hydrous minerals. Thus, the water carried by the subducting slab to the MTZ depth will be mostly, if not all, transported to the overriding continental lithosphere and contribute to the craton modification. However, the degree of lithospheric modification is still lower in the big mantle wedge regime than the flat subduction regime (Figures 10b and 10b′). It is worth noting that the mantle rocks along the pathway of the upward water migration are still marked as hydrated mantle (with light blue color in Figure 10b), although they cannot absorb any water.

Summary of the Numerical Models
The influences of variable model parameters on the subduction-induced deep hydration and overriding craton modification are summarized in Table 2. It demonstrates that the flat subduction is resulting in the models with a relatively thick (≥15 km) subcrustal serpentinite layer in the oceanic plate. In contrast, if the

Journal of Geophysical Research: Solid Earth
LI serpentinite layer is absent or relatively thin (≤10 km), the flat subduction is not predicted. The mechanism of flat subduction is an important issue, although not the focus of this study, which has been systematically investigated previously (e.g., Huangfu et al., 2016;Li et al., 2011;Manea et al., 2017;van Hunen et al., 2004). It could be attributed to many factors, for example, the young subducting slab, the oceanic plateau subduction, and the seaward movement of the overriding plate. In this study, the complex density variation during the mantle hydration is not applied, which indicates the same reference density for the dry and hydrated mantle rocks, although the density will be decreased during partial melting (Table S2). The contribution of the thick serpentinite layer to the formation of flat subduction may be due to the hydration-induced rheological weakening of the slab. Thus, the slab unbending occurs more easily when subducting to the sublithospheric depth, which finally results in the flat slab subduction beneath the overriding lithosphere.
The hydration in the big mantle wedge and its contribution to the overriding lithospheric modification depend on not only the serpentinite layer that carries water to the deep mantle but also the preexisting water condition of the NAMs in the upper mantle and particularly the MTZ (Table 2). In order for the big mantle wedge regime to work on the modification of the overriding craton, a greatly hydrated or even saturated MTZ should have already been present; otherwise, the water carried by a subducting slab will be absorbed by the potential high water capacity of the MTZ rocks (Figures 10a and 10a′). Since the hydration of the MTZ is generally considered to be correlated with subduction, it thus indicates that the big mantle wedge regime may not take effect with a single subducting slab but with a complex history of many subduction events.

Subduction-Induced Deep Water Cycling
The numerical models indicate that the water carried by the oceanic crust to the MTZ is negligible (Figure 4), which is mainly attributed to the low temperature condition of the "choke point" in the water capacity diagrams of basalt and gabbro (Figures 3b and 3c). On the other hand, the oceanic crust of the subducting slab is generally thin, that is, 6-8 km, which can be easily heated by the thermal conduction during slab subduction. Even for an oceanic plateau with a thicker crust of up to 32 km, it is still hard to carry water to the MTZ ( Figure S2). It thus indicates that most of the water contained in the oceanic crust cannot pass through the "choke point" of the hydrous phase transitions at about 300-km depth (Figures 3b and 3c). The subcrustal serpentinite layer is an efficient way to carry water to the deeper mantle (e.g., Figure 5). The "choke point" in the water capacity diagram of the mantle rock locates at about 600°C with the pressure of about 6 GPa (Figure 3d). The temperature condition (~600°C) is higher than those of crustal rocks (~400-500°C), whereas the pressure condition (~6 GPa) is lower than the choke point of the crust (~10 GPa) (Figure 3). Both conditions help to keep the temperature lower than the threshold value of the "choke point" of the hydrated mantle rocks. Consequently, significant quantity of water can still be contained in the subcrustal hydrous layer of the subducting slab and carried to the deeper mantle.

Comparison of Slab Dehydration in Flat Subduction Versus Big Mantle Wedge Regimes
The numerical models indicate that the hydration processes in both the flat subduction and the big mantle wedge regimes can contribute to the overriding craton modification (e.g., Figure 5). In order to further compare their efficiencies and timescales, two simplified models are conducted (Figure 11), in which the subduction process is not simulated. In the flat subduction model, an oceanic slab of 60 Ma old and 500 km long, with a 15-km-thick subcrustal serpentinite layer, is prescribed directly beneath the overriding continental lithosphere, whereas a slab with the same configuration is put at the bottom of the MTZ in the big mantle wedge model.
In both models, the dehydration occurs first from the edges of the oceanic slab, which then migrates to the center ( Figure 11). Finally, most of water in the prescribed slab is transferred to the overriding lithosphere and the upper mantle. However, the total dehydration in the flat subduction model is much faster (~30 Myr) than that in the big mantle wedge model (~62 Myr), because the threshold temperature condition for dehydration at a shallower depth is much lower than that with higher pressure in the MTZ (Figure 3d). Consequently, the slab dehydration in the MTZ requires higher temperature condition and thus longer time for heating.
Finally, the resulting hydration of the overriding lithosphere is more intense in the flat subduction model than the big mantle wedge model, due to the much longer pathway of water transportation in the latter, which leads to the loss of water during the upper mantle hydration. Alternatively, most of water in the flat slab will be transferred into the cold core of the overriding lithosphere and absorbed in the hydrous minerals. Thus, the hydration-induced craton modification is more efficient in the flat subduction regime than the big mantle wedge regime, which agrees with the results of the complex subduction models (e.g., Figure 5).

Implications for the Modification of NCC in the Mesozoic
The NCC used to have an old (~2.5 Ga), thick (180-200 km), and cold (~40 mW/m 2 ) lithosphere (e.g., Menzies et al., 1993;Wu et al., 2019;Xu, 2001;Zhu & Xu, 2019); however, a present thin and hot lithosphere is observed beneath the eastern NCC (Figure 12), which suggests that the ancient cratonic root has been thinned for~100 km. It is generally believed that the peak stages of the NCC modification occur in the early Cretaceous, or more precisely at~125-120 Ma (Zheng et al., 2018;Zhu et al., 2012).
It is generally proposed that the modification/destruction of the NCC is related to the fluid/melt activity during the Izanagi plate subduction in the Mesozoic (Liu & Li, 2018, and references therein). The previous numerical models are generally focusing on the big mantle wedge regime, which can be further divided into two types, that is, thermal convection (He, 2014) and hydrous MTZ upwelling . Both processes may contribute to the bottom erosion of the overriding lithosphere. In this study, the subduction-induced hydration and overriding lithospheric modification in the flat subduction and the big mantle wedge regimes are systematically compared, which are further compared to the geological records of the NCC.
One of the most important geological responses of the lithospheric modification is the resulting magmatism. The most prominent magmatism in the NCC occurs in early Cretaceous, which appears to have an eastward younging trend (Wu et al., 2019). In contrast, the magmatism in Jurassic, that is, before the NCC modification, is characterized by a reverse, westward younging trend. The numerical model shows that the progressive flat subduction leads to the partial melting migrating further away from the trench (Figures 6 and 9), which is thus consistent with the westward younging trend of Jurassic magmatism in the NCC. However, we need to keep in mind that it is just the potential agreement because the complex processes of magma migration as well as the further crustal-level partial melting and magmatic emplacement are not directly simulated. Afterward, the magmatism migrates toward the trench during the transition from flat to steep slab subduction (Figures 6 and 9), which agrees with the eastward younging trend of early Cretaceous magmatism in the NCC. On the other hand, the magmatism in the big mantle wedge regime does not show clear spatiotemporal trends, which may occur simultaneously in multiple cratonic regions above the stagnant slab ( Figure 6). If constrained by the distribution of magmatism, it indicates that the flat slab subduction beneath the overriding continental lithosphere may play a more important role in the NCC modification in the Mesozoic.
A problem is the temporal and spatial scales of progressive flat subduction, which are generally short in the current numerical models, for example,~15 Myr and~600 km from the trench (e.g., Figure 5a). However, if constrained by the westward younging distribution of magmatic rocks in the NCC, the progressive flat subduction may exist for the entire Jurassic of~50 Myr and~2,000 km from the trench (Wu et al., 2019). Such long flat subduction is not obtained in the current numerical models, which also does not exist in the present Earth. It is worth noting that Wu et al. (2019) finally estimated the flat subduction lasting for 20 Myr, that is, 160-140 Ma, which is comparable to the current models ( Figures 5 and 9). In addition, the spatial scale of 2,000 km includes the later back-arc extension and the resulting marginal sea stretching (Figure 12). If going backward in time to the Cretaceous, the width of east NCC is much narrower and may be comparable with not only the current numerical models but also the maximal length scale of the present-day flat slab beneath South America (Espurt et al., 2008).
In the big mantle wedge regime, the water liberates from the stagnant slab in the MTZ, which may further contribute to the hydration of overriding cratonic lithosphere ( Figure 5). In the region not far from the trench, that is, within the spatial scale of flat subduction, the effects of both regimes are overlapped. In contrast, the region beyond the flat subduction may be purely controlled by the dehydration of the flattened slab Colors represent topography as shown in the color bar at the bottom, which is produced by GMT (Wessel et al., 2013) with the data from ETOPO1 Global Relief Model (https://www.ngdc.noaa.gov/mgg/global/). The thin black lines are shown for the boundary faults and/or suture zones of the NCC and surrounding regions (Zhu et al., 2012). The thin colored lines are the contours of lithospheric thickness of NCC with the open-sourced data from "http://www.craton.cn/data" (Chen, 2010;Chen et al., 2014). The thick yellow lines denote the trenches of present Pacific and Philippine plate subduction zones, respectively. The thick blue line is the estimated trench position of the Izanagi subduction zone in the Jurassic before the major modification of eastern NCC (Wu et al., 2019). at the bottom of the MTZ. The spatial scale of the cratonic hydration in the big mantle wedge regime is >1,000 km in the current numerical models, which is however hard to quantify and is strongly dependent on the length and properties of the stagnant slab in the MTZ. The length scale of stagnant slab in the MTZ is much larger than the flatly subducted slab beneath the overriding lithosphere in both the numerical models (e.g., Figure 5) and the present-day natural Earth (Espurt et al., 2008;Fukao & Obayashi, 2013).

Conclusions
A series of numerical models are constructed to study the dynamics of subduction-induced deep hydration processes and the effects on the overriding craton modification. The main conclusions are shown below: 1. Subducting oceanic crust cannot carry water to the deeper mantle, for example the MTZ; however, the subcrustal serpentinite layer in the sinking slab is an efficient way for the deep water cycling. 2. The flat slab subduction can significantly hydrate the overriding cratonic lithosphere for a region within about 600 km from the trench. During the progressive flat subduction, the partial melting and magmatism migrate far away from the trench, which will be reversed and go backward to the trench during the transition from flat to steep subduction. 3. Subduction-induced deep hydration in the MTZ and the big mantle wedge is strongly dependent on the subcrustal serpentinite layer in the sinking slab as well as a nearly saturated preexisting MTZ. It can contribute to the overriding craton modification for a larger region of >1,000 km from the trench, which is, however, generally slower and weaker than the flat subduction regime. 4. The water capacity of MTZ rocks and their preexisting water content before subduction can strongly modify the upward water migration and the overriding craton modification in the big mantle wedge regime. 5. The modification of NCC is more likely to be controlled by the flat subduction of Izanagi plate in the late Jurassic to early Cretaceous. The flattened slab in the MTZ, that is, the big mantle wedge regime, may also play a certain role. This work was supported by the NSFC project (41688103) and the Strategic Priority Research Program (B) of Chinese Academy of Sciences (XDB18000000). Numerical simulations were run with the clusters of National Supercomputer Center in Guangzhou (Tianhe-II). T. Gerya, J. Connolly, S. Karato, and Ling Chen are greatly acknowledged for the helpful discussion. The constructive reviews by V. Magni, T. Nakagawa, and the Associate Editor F. Capitanio improved the paper. The topographic map of Figure 12 is produced with GMT (Wessel et al., 2013) with the data from ETOPO1 Global Relief Model (https:// www.ngdc.noaa.gov/mgg/global/). The contours of lithospheric thickness in Figure 12 are plotted with the open-sourced data from "http://www. craton.cn/data," which is also provided in the public repository of Zenodo (https://doi.org/10.5281/ zenodo.3946588) along with the model data. The figures of numerical models are produced by Matlab and further compiled by Adobe Illustrator.