Dynamics of the Subducted Izanagi-Pacific Plates Since the Mesozoic and Its Implications for the Formation of Big Mantle Wedge Beneath Eastern Asia

Introduction

The tectonic units in East Asia mainly consist of four terrain blocks: the Siberia Craton to the north, the North China Craton (NCC) in the middle, the South China Block (SCB) on the south and the Central Asian Orogenic Belt (CAOB) in between (Domeier and Torsvik, 2014; Figure 1). These four blocks were separately drifting on their own at least before the Triassic (e.g., Domeier and Torsvik, 2014; Matthews et al., 2016). Since the Late Paleozoic, the NCC amalgamated with the CAOB during the Late Permian to Early Triassic and collided with the SCB in the Triassic, while the northeast China combined with the Siberia Craton during the Jurassic to Early Cretaceous, and these four terrains converged together to become part of the Eurasian plate (e.g., Sorokin et al., 2020). To the east, the western Pacific plate is actively subducting along the Kuril, Japan and Mariana trenches, whereas the Philippine Sea Plate is descending along the Nankai trough and the Ryukyu trench (e.g., Matthews et al., 2016; Figure 1).

FIGURE 1

FIGURE 1. Major tectonic settings of Eastern Asia. The blue, red and cyan lines denote the shapes of South China Block (SCB), North China Craton (NCC) and Amuria Block, respectively. The red lines with triangles denote the subduction zones. The white dashed lines denote locations of the vertical cross-sections P1-P4 shown in Figure 2. The start, center (red star) and end coordinates with a great circle distance of 10^o for each profile are marked (white circle) on each cross-section.

Since the early Mesozoic, extensive tectonic deformations have occurred in Eastern Asia (Zhu et al., 2012; Wang et al., 2018; Li et al., 2019). The northeast China experienced intensive extensional deformations in the Early Cretaceous which are characterized by widespread volcanism and extensional basins (Meng, 2003; Wu et al., 2005; Wang et al., 2006). The interior of the NCC had undergone extensive modification and reactivation from Mesozoic to Cenozoic as indicated by large-scale deformation, magmatic activities, and basin formation (Zhu et al., 2012; Wu et al., 2019). The mantle xenoliths studies show that the ancient lithospheric mantle of the SCB was replaced by the newly accreted mantle through lithospheric extension and asthenosphere upwelling during the Meso-Cenozoic time (Li et al., 2014a; Li et al., 2014b). The tectonic evolution of the Eastern Asia area has been suggested to be closely related to the subduction history of the Paleo-Pacific plate (i.e., the Izanagi plate) and the west Pacific plate (e.g., Zhu and Xu, 2019 and references therein). The large-scale extension and volcanic activities in northeast China are supposed to be associated with the slab rollback of the subducted Izanagi plate (e.g., Xu and Zheng, 2017; Li et al., 2019; Ma and Xu, 2021), and many researchers regard Izanagi-Pacific subduction as one of the principal triggers for the reactivation and destruction of the NCC (e.g., Zhu and Xu, 2019; Liu et al., 2021). It has also been proposed that the Pacific subduction influenced the SCB mainly during the middle to late Mesozoic and produced large-scale granitoid and volcanic rocks (Sun et al., 2007; Liu et al., 2012).

However, the subduction history of the Izanagi-Pacific plate is poorly constrained and remains controversial (e.g., Sun et al., 2007; Yang, 2013; Domeier and Torsvik, 2014; Matthews et al., 2016; Müller et al., 2019). The link between the Izanagi-Pacific plate subduction and the tectonic evolution of Eastern Asia has not been fully explored, especially lack of geodynamical verifications for the Mesozoic time. Based on seismic tomography studies many scholars have argued that the Pacific slab, and perhaps the older Paleo-Pacific slab are stagnant in the mantle transition zone beneath Eastern Asia (Fukao et al., 2009; Zhao et al., 2011) to form a “big mantle wedge” (BMW) that the overriding lithosphere evolved in response to subduction dynamics of the BMW (Zhao et al., 2011; Xu et al., 2018; Zhu and Xu, 2019). Some recent studies on magmatic activities in Eastern Asia have also indicated the development of a big mantle wedge between 145–120 Ma (e.g., Ma and Xu, 2021). On the other hand, recent studies have proposed that the present-day stagnant slab in the mantle transition zone beneath Eastern Asia is the subducted Pacific slab rather than the Izanagi slab (Ma et al., 2019). In addition, global mantle convection models show that the stagnant slab under Eastern Asia largely results from subduction in the past 20–30 Myr (Mao and Zhong, 2018). The stagnation of the Izanagi-Pacific slab beneath Eastern Asia has been investigated in many numerical experiments (e.g., Seton et al., 2015; Peng et al., 2021a; Liu et al., 2021). However, most of these studies have focused on the stagnation of the Pacific slab in the Cenozoic time, while whether the Izanagi slab in the Mesozoic time could be stagnant or not remains unclear. Therefore, some remaining questions that 1) how does the Izanagi-Pacific plate evolve beneath Eastern Asia since the Mesozoic time, and 2) how its subduction dynamics would potentially affect the formation of the BMW, especially in the Mesozoic time, deserve to be well explored from a geodynamical perspective.

In this study, we perform 3-D global mantle convection models to study the slab dynamics of the subducted Izanagi-Pacific plate and its influence on the tectonic evolution of Eastern Asia. In numerical models, the plate motion history is incorporated to 1) constrain the boundary conditions and 2) track the Izanagi and the Pacific plate separately to distinguish their evolutionary processes.

Methods

The dynamical evolution of the mantle convection system is governed by partial differential equations for conservation of mass, momentum and energy and advection of composition (McNamara and Zhong, 2004). We assume the mantle is incompressible and solve the governing equations under the Boussinesq approximation with modified numerical code CitcomS (Zhong et al., 2008) and the non-dimensional governing equations are:

$\begin{array}{c}\nabla \cdot \overrightarrow{u}=0\end{array}(1)$

$\begin{array}{c}-\nabla P+\nabla \cdot \left(\eta \dot{\epsilon }\right)+\left(RaT-{\displaystyle \sum _{k}R{b}^k{\mathrm{\Gamma }}_k}\right)\widehat{e_r}=0\end{array}(2)$

$\begin{array}{c}\frac{\partial T}{\partial t}+\left(\overrightarrow{u}\cdot \nabla \right)T={\nabla }^{2}T+Q\end{array}(3)$

where $\overrightarrow{u}$ is the velocity, P is the dynamic pressure, $\eta$ is the viscosity, $\dot{\epsilon }$ is the strain rate, Ra is the Rayleigh number, T is the temperature, $\widehat{e_r}$ is the unit vector in the radial direction, t is the time, Q is the internal heating. The thermal Rayleigh number Ra and the phase-change Rayleigh number Rb^k are defined as:

$\begin{array}{c}Ra=\frac{\alpha \rho g_0\mathrm{\Delta }T{R}^3}{\kappa {\eta }_0}\end{array}(4)$

$\begin{array}{c}R{b}^k=\frac{\mathrm{\Delta }{\rho }_k{g}_0{R}^3}{\kappa {\eta }_0}\end{array}(5)$

where $\alpha ,\rho ,g_0,\mathrm{\Delta }T,\kappa ,{\eta }_0$ are dimensional parameters for the reference thermal expansivity, density, gravitational acceleration, temperature difference between the bottom and the surface, thermal diffusivity and viscosity, respectively. $\mathrm{\Delta }{\rho }_k$ is the density jump for the kth phase change and R is the radius of the Earth. Notice that the Rayleigh number is defined using the radius of the Earth and is ∼10 times larger than that defined using the mantle thickness.

A phase-change function formulation is used here to represent phase changes as in earlier studies (Christensen and Yuen, 1985; Zhong and Gurnis, 1994). Γ_k is defined in dimensionless form as:

$\begin{array}{c}{\mathrm{\Gamma }}_k=\frac{1}{2}\left[1+\tanh \left(\frac{{\pi }_k}{\delta }\right)\right]\end{array}(6)$

where δ is the phase change width that measures the depth segment of phase change, and π_k is the dimensionless “excess pressure” as

$\begin{array}{c}{\mathrm{\pi }}_k=d-d_k-{\gamma }_k\left(T-T_k\right)\end{array}(7)$

where d is the depth, d_k and T_k are the reference depth and temperature of phase change k, and γ_k is the Clapeyron slope. The Clapeyron slope is normalized by ρgR/ΔT. The phase change at 670 km depth from spinel to post-spinel changes is included in our models. All relevant physical parameters are listed in Table 1.

TABLE 1

TABLE 1. Physical parameters of numerical models.

The models have a dimensionless radius of 0.55 and 1.0 at the CMB and the surface, respectively. The computational domain is divided into 12 caps with each cap containing 128 × 128 × 80 elements and employs grid refinements near the surface (∼25 km), the phase transition (∼25 km) and the CMB (∼25 km) in the radial direction. The surface temperature is constant with nondimensional T = 0 (i.e., 273 K), whereas the CMB is thermally insulating (i.e., zero heat flux boundary condition) to prohibit the formation of upwelling plumes. The heat flux that is expected to come out of the core in our models is smaller due to the relatively small temperature drop (i.e., 2500 K, Table 1) than estimated constraints (Lay et al., 2008). However, since we focus on upper mantle slab dynamics and do not consider the bottom thermal boundary layer above the CMB, we, therefore, use a zero heat flux instead of isothermal condition for the bottom boundary and a relatively smaller temperature drop. Models are internally heated with a nondimensional internal heating rate of Q = 100 (i.e., ∼6 × 10^–12 W/kg). The velocity boundary condition for the top boundary is time-dependent by imposing the plate velocity history from the reconstruction model (e.g., Matthews et al., 2016), while the CMB is free-slip. The initial temperature for the lithosphere (e.g., above 150 km) is obtained by following the half-space cooling model with a plate age of 100 Ma and below the lithosphere is 0.52 (or 1573 K) everywhere, respectively (Supplementary Figure S1).

The viscosity structure in our models is both depth and temperature-dependent, expressed as $\eta \left(T,r\right)={\eta }_{r}exp\left[A\left(0.5-T\right)\right]$, where η_r is the depth-dependent pre-factor and A is the activation coefficient which is 9.21 in all cases (e.g., equivalent to ∼190 kJ/mol of activation energy dimensionalized by R_gasΔT), leading to a maximum of 4 orders of viscosity change due to variation of temperature. The viscosity pre-factor η_r is 1.0, 0.02, 2.0 for the lithosphere, asthenosphere and the lower mantle, respectively. All cases in this study include a thin (∼60 km) weak layer (η_r is set to 0.002) below the phase change boundary to simulate the rheological effects of the phase change (e.g., Mitrovica and Forte, 2004), which is suggested to play an important role in producing horizontally deflected slabs in the transition zone (Mao and Zhong, 2018) (Supplementary Figure S1). Our models also include a linear increase of thermal diffusivity and a decrease of thermal expansivity from the surface to the CMB, as indicated from mineral physics experiments (Gibert et al., 2003; Katsura et al., 2009), by a factor of 2.18 and 2.5, respectively (e.g., Mao and Zhong, 2018).

In this study, we use ∼315 million tracers (20 per element) in numerical models to track the Izanagi plate and the Pacific plate separately. We introduce a tracking technique by which we can explicitly distinguish the evolution of the Izanagi slab and the Pacific slab in the mantle. The general tracking procedure is as follows: firstly, a time-dependent plateID that represents different plates at the surface is obtained by extracting the plateID information from the reconstruction model with a time interval of 1 Myr using the GPlates software (Müller et al., 2018); secondly, at each timestep, 1) we choose the closest time point to assign the plateID for each surface node, 2) for each tracer in an element above a given depth (e.g., above 150 km), and determine the corresponding surface element that is in the same vertical column with this tracer’s element; if the nodal plateID of this surface element matches the plate that we want to track, then we remark the tracer as the target plate with a given flavor, 3) for tracers below the given depth, we do not remark the tracers and leave the flavors of the tracers as they are and move to next timestep.

Results

We present five models that differ in plate motion history model, initial setup and the presence of weak layer in the transition zone and use case 1 as reference (Table 2). Case 1 uses the plate model of Matthews et al. (2016) and starts to run from 410 Ma. Case 2 uses the plate reconstruction model of Young et al. (2019) and starts to run from 410 Ma. Case 3 is similar to case 1, except that the plate motion history after 90 Ma is replaced by that from Yang (2013), while the plate motion history before 90 Ma is the same as case 1 from Matthews et al. (2016). In case 3, the moving direction of the Izanagi plate was rotated by ∼20° clockwise after 90 Ma so that the plate moved north-northwestward at around 84 Ma. The average moving speed of the Izanagi plate during 90–77 Ma is ∼20 cm/yr. After 77 Ma the Izanagi plate subducted beneath northeast Asia and totally fell into the mantle around 55 Ma. Case 4 uses the plate model of Matthews et al. (2016) but starts to run from 300 Ma. Case 5 is the same as case 1 but without the weak layer at 670-km transition zone.

TABLE 2

TABLE 2. Numerical model parameters.

Slab Structure of the Izanagi-Pacific Plate

We first show the results of case 1 which uses the plate model of Matthews et al. (2016) started from 410 Ma. Figure 2 shows several cross-sections of the temperature and compositional fields at different locations for this case at t = 0.0 Ma, i.e., present-day time. The P-wave velocity perturbations from two seismic tomography models TX2019slab (Lu et al., 2019) and GAP_P4 (Fukao and Obayashi, 2013) are also shown in Figure 2. The slabs show a variety of morphologies at different locations beneath Eastern Asia. The slab subducts nearly vertically into the transition zone depth along the Mariana Trench in southeast Asia (Figures 2A,B, columns 1 and 2). However, beneath the north and northeast Asia, the slabs are horizontally deflected and tend to stagnate near the transition zone depth (Figures 2C,D, columns 1 and 2). The general slab structures in the upper mantle beneath Eastern Asia are quite similar for different seismic models (Figure 2, columns 1 and 2). For case 1, the structures of the subducted slabs as represented by the cold thermal anomalies and the compositional fields show good first-order consistency with the structures of fast velocity anomalies in seismic models, especially the existence of stagnant slabs beneath Eastern Asia in the mantle transition zone (Figure 2, columns 3 and 4).

FIGURE 2

FIGURE 2. Cross-sections of P-wave velocity perturbations from seismic tomography models TX2019slab (first column) (Lu et al., 2019) and GAP_P4 (second column) (Fukao and Obayashi, 2013), and the present-day non-dimensional residual temperature (third column) and composition (last column) from numerical model case 1. The black line in each cross-section corresponds to the 670-km depth. The locations of (A–D) correspond to profiles P1-P4 as shown in Figure 1, respectively. In figures of tomography models (i.e., first and second columns) the green lines indicate contours of non-dimensional residual temperature δT of −0.05 from the third columns, i.e., material 5% colder than the ambient mantle that supposed to be associated with subducted slabs. In figures of compositional field, the blue and red colors represent the Izanagi Plate and the Pacific Plate, respectively.

However, although the seismic observations show fast velocity anomalies beneath Eastern Asia, they do not tell whether the velocity anomalies are associated with the Izanagi plate or the Pacific plate. By tracking the Izanagi plate and the Pacific plate separately in our models (Methods), we find that at the present-day time, most of the Izanagi slabs beneath the Eastern Asia region have entered the lower mantle, except a few regional slab fragments left near the 670-km depth (Figure 2, column 4). The fast velocity anomalies at the 670-km transition zone are associated with the new Pacific plate rather than the old Izanagi plate (Figure 2, column 4), which is consistent with other previous studies (Ma et al., 2019; Liu et al., 2021).

Figure 3 shows some snapshots of the temperature and composition fields at 670 km depth for case 1 together with a cross-section showing slab structures beneath Eastern Asia, and a full-time evolution is provided in Supplementary Figure S2. In the early Mesozoic (e.g., 250 Ma), the mantle structure beneath Eastern Asia is rather complex and abundant Izanagi slabs are trapped and accumulate at 670 km depth before falling into the lower mantle (Figure 3A and Supplementary Figure S2). The subducted slabs distribute not only at the leading edge of the trench but also at the trailing edge due to slab rollback (Figure 3A and Supplementary Figure S2). From early to middle Mesozoic (∼150 Ma), the Izanagi plate continues to subduct westward beneath the Eurasian plate and the slabs mostly distribute at the leading edge of the trench, while many of the previously accumulated Izanagi slabs have fallen into the lower mantle (Figure 3B and Supplementay Figure S2). To the late Mesozoic, the subducted slabs generally first accumulate at the transition zone depth before entering the lower mantle, which results in a buckling morphology when they sink to greater depths (Figures 3C,D and Supplementary Figure S2). Meanwhile, the middle ocean ridge of the Izanagi-Pacific plate also moves westward closely to the Eurasian plate while no subduction initiates for the Pacific plate (Figure 3D and Supplementary Figure S2).

FIGURE 3

FIGURE 3. Slab structures of the Izanagi-Pacific plate from case 1 at (A) 250 Ma, (B) 150 Ma, (C) 80 Ma, (D) 60 Ma, (E) 40 Ma and (F) 10 Ma. The first column shows the 3D slab structures of the Izanagi plate (blue) and the Pacific plate (red) at the Eastern Aisa region. The second and third columns show the residual temperature and composition fields at 670 km depth, respectively. Convergent and ridge and transform boundaries are represented by magenta and black lines, respectively. The fourth and fifth columns show the cross-section of residual temperature and composition fields at section P3 as shown in Figure 1. In figures of compositional fields, the blue and red colors denote the Izanagi plate and the Pacific plate, respectively. A full-time evolution of the slab structures is provided in Supplementary Figure S2.

To the early Cenozoic (∼40 Ma), the entire Izanagi plate has subducted to the mantle including the ridge and the subduction of the Pacific plate has initiated (Figure 3E and Supplementary Figure S2). Although most Izanagi slabs have fallen to the lower mantle, some remnants are still preserved in the upper mantle at transition zone depth, while the slabs of the Pacific plate remain at a shallower depth (Figure 3E and Supplementary Figure S2). Since then, the residuals of the Izanagi slabs continue to sink from the transition zone until they have completely entered the lower mantle at the present-day time (Figure 3F and Supplementary Figure S2). Meanwhile, the subduction of the Pacific plate has fully developed, and the slabs are mostly trapped or flattened above the transition zone instead of penetrating to the lower mantle (Figure 3F and Supplementary Figure S2).

Time-Dependent Slab Structure Beneath the NCC and the SCB

Here, we present the slab structures of the Izanagi-Pacific plate beneath the NCC and the SCB. Figure 4 shows some snapshots of the temperature and composition fields beneath the NCC and the SCB from case 1, in which the locations of the cross-sections are not fixed but move along with the NCC and SCB over time. The fully-time evolution of slab structures is provided in Supplementary Figures S3, S4, respectively. From early to middle Mesozoic, the slab structures beneath the NCC and the SCB are quite similar, that is the subducted slabs of the Izanagi plate mostly directly pass through the 670-km boundary into the lower mantle and tend to be stirred into small pieces by mantle convection (e.g., 250–180 Ma, Figure 4A, Supplementary Figures S3, S4). From middle to late Mesozoic, instead of directly penetrating to the lower mantle, some of the subducted slabs of the Izanagi plate are trapped or stagnant above the mantle transition zone depth beneath the NCC and spread westward to the central part of the NCC (e.g., ∼130 Ma, Figure 4B, left columns). Similar slab structures are also observed beneath the SCB at this stage but with a different time of slab stagnation (e.g., ∼140 Ma, Figure 4B, right columns). The stagnation of the subducted slabs which extends with a horizontal length of 800–1,000 km facilitates the formation of a big mantle wedge beneath the lithosphere of the NCC and the SCB (Figure 4B).

FIGURE 4

FIGURE 4. Slab structures of the Izanagi-Pacific plate beneath the North China Craton (left columns) and the South China Block (right columns) at some time snapshots from case 1. In (A–F), the left and right columns are cross-sections of the residual temperature and composition fields across the NCC and the SCB, respectively. The fully-time evolution of the slab structures for the NCC and the SCB are provided in Supplementary Figures S3, S4, respectively.

With time, more and more slab materials accumulate at the transition zone depth and eventually they pass through the 670-km boundary to sink into the lower mantle (e.g., Figure 4C). The mantle wedge is thus destroyed due to this subsequent accumulation and sinking of Izanagi slabs in the middle Mesozoic (∼110 Ma, Figure 4C). From the middle to the late Mesozoic (∼70 Ma), another slab stagnation event has developed beneath both the NCC and the SCB, although the duration times are somewhat different beneath these two blocks (Figure 4D, Supplementary Figures S3, S4). The timescales of these multiple slab stagnation events, however, are quite similar with a magnitude on the order of tens of million years (Supplementary Figures S3, S4). In the early Cenozoic (∼40 Ma), most of the Izanagi slabs beneath the NCC and the SCB have completely fallen into the lower mantle (Figure 4E). From the middle Cenozoic to the present, the subducted slabs of the Pacific plate are mostly trapped and lie horizontally at the 670-km mantle transition zone depth (Figure 4F, Supplementary Figures S3, S4). However, abundant slabs of the Pacific plate have subducted and spread to the center part of the NCC which promotes the formation of a big mantle wedge beneath the lithosphere, while there are fewer Pacific slabs beneath the SCB due to the block of the Philippine plate (Supplementary Figures S3, S4).

Effect of Plate Reconstruction Models and Model Setup

The motion histories of many plates at the early stage have been destroyed, and the reconstructions of the plate motion history remain under debate. Here we test the influences of different plate motion history models on our results. Case 2 uses the same model setup as case 1 but with surface plate motion history from a different plate reconstruction model of Young et al. (2019). The present-day upper mantle slab structures beneath Eastern Asia in case 2 are quite similar to case 1, although some differences exist for slab morphologies in the lower mantle depths (Figures 2A–D, 5A–D). However, the slab dynamics of the Izanagi plate for case 2 reveal some significant differences in the Mesozoic time compared to case 1. For example, there is no slab stagnation in the middle Mesozoic (∼130 Ma) beneath the NCC for case 2 to facilitate the formation of a big mantle wedge (Figure 5E), while the horizontal extent of the stagnant slabs beneath the SCB (∼140 Ma) in case 2 seems to be also smaller than that in case 1 (Figures 4B, 5G). On the other hand, the slab dynamics of the Izanagi plate in the late Mesozoic time in case 2 is similar to case 1, in which the slabs tend to be temporarily stagnant above the transition zone beneath both the NCC and the SCB (Figures 5F,H).

FIGURE 5

FIGURE 5. Slab structures of the Izanagi-Pacific Plate beneath the Eastern Aisa for case 2 in which the plate motion history is from Young et al. (2019). (A–D) are present-day temperature and composition of slab structures for cross-sections P1-P4 in Figure 1, respectively. (E,F) are temperature and compositional fields from cross-sections through the NCC at 130 Ma and 70 Ma, respectively. (G,H) are temperature and compositional fields from cross-sections through the SCB at 140 Ma and 70 Ma, respectively. In panels of temperature fields, the green lines indicate contours of non-dimensional residual temperature δT of -0.05 from case 1 at the same time and cross-section. In panels of compositional fields, the magenta and green lines indicate compositional contours of -0.5 and 0.5 from case 1 which associated with the Izanagi slab and the Pacific slab, respectively.

There are also some debates on the subduction history on the Izanagi-Pacific plate: one view suggests that the ridge system of the Izanagi plate is parallel to the coastline and subducts to the west together with the Izanagi plate (e.g., Matthews et al., 2016), whereas another view suggests that the ridge system is perpendicular to the coastline, and the Izanagi plate subducts to the north and finally extinct (e.g., Yang, 2013). Case 3 uses the same model setup as case 1, except that the plate motion history after 90 Ma is replaced by that from Yang (2013), while the plate motion history before 90 Ma is the same as case 1 from Matthews et al. (2016). The general slab structures beneath Eastern Asia at present-day for case 3 are consistent with case 1, especially in the upper mantle depths (Figures 6A–D). However, from the late Mesozoic to early Cenozoic, the slabs of the Izanagi plate in case 3 tend to only subduct to the northern part of Eastern Asia, while in the southern part (e.g., beneath the SCB) there are no remnants of the Izanagi slabs left at the transition zone depths (Figures 6E,F). To the middle Cenozoic, there are no Izanagi slabs left in the upper mantle beneath the NCC and SCB in case 3 and the subducted slabs of the Pacific plate have reached the margin of the Eastern Asia region, which is earlier than that in case 1 (Figures 3E, 6G). The slabs of the Pacific plate in the late Cenozoic for case 3 are similar to case 1 while all Izanagi slabs have sunk to the lower mantle (Figures 3F, 6H).

FIGURE 6

FIGURE 6. Slab structures of the Izanagi-Pacific Plate beneath the Eastern Aisa for case 3 in which the plate motion history after 90 Ma is from Yang (2013). (A–D) are present-day temperature and composition of slab structures for cross-sections P1-P4 in Figure 1, respectively. (E–H) are some time snapshots of the residual temperature and composition fields at 670 km depth, respectively. Line descriptions are the same as that in Figure 5.

In this study, we use a simple initial temperature condition and let the model run very earlier to obtain well-developed mantle convection and subduction pattern before the Mesozoic. To justify this model setup we also run another case, case 4, which uses the same initial condition as case 1 except that the model runs starting at 300 Ma with plate motion history from Matthews et al. (2016). In the Mesozoic time, the slab structures of the subducted Izanagi plate beneath the NCC and SCB are quite similar between case 4 and case 1, both in the upper mantle and in the lower mantle, except that there are some remnants of the Izanagi slabs near the core-mantle boundary in case 1 rather than in case 4 (Figures 4, 7). The slab stagnation events of the Izanagi plate beneath the NCC (∼130 Ma) and the SCB (∼140 Ma) are also observed in case 4 which are consistent with that in case 1 (Figure 7). To the Cenozoic, the slab structures beneath the NCC and the SCB in case 4 are nearly identical to case 1 (Figures 4, 7).

FIGURE 7

FIGURE 7. Slab structures of the Izanagi-Pacific Plate beneath the North China Craton (left columns) and the South China Block (right columns) at some time snapshots for case 4. In (A–F), the left and right columns are cross-sections of the residual temperature and composition fields across the NCC and the SCB, respectively. Line descriptions are the same as those in Figure 5.

In most of our models, we include a weak viscosity layer that is associated with the phase transition at 670 km depth, which has been suggested to be crucial for slab stagnation in the mantle transition zone (e.g., Mao and Zhong, 2018). Case 5 uses the same model setup as case 1 except that the weak layer due to phase transition at 670 depth is removed. Although in case 5 the slab structures in the upper mantle beneath the NCC and SCB are similar with case 1, the slab stagnation of the subducted slabs has not well developed as prominent as that in case 1, especially for the Izanagi slabs in the Mesozoic time (Figure 8). However, the Pacific slabs are still mostly stagnant in the mantle transition zone even without the presence of the weak layer which is consistent with that in case 1 (Figures 4, 8).

FIGURE 8

FIGURE 8. Slab structures of the Izanagi-Pacific Plate beneath the North China Craton (left columns) and the South China Block (right columns) at some time snapshots for case 5. In (A–F), the left and right columns are cross-sections of the residual temperature and composition fields across the NCC and the SCB, respectively. Line descriptions are the same as those in Figure 5.

Discussions

In this study, we perform global mantle convection models to examine the behaviors of the Izanagi plate and the Pacific plate since the Mesozoic time. We introduce a tracking procedure in numerical models to explicitly distinguish the slab materials of the two plates, which provides a more robust demonstration of slab dynamics in their evolutionary processes than identifying slabs solely from thermal structures as in many previous studies (e.g., Seton et al., 2015; Peng et al., 2021a; Liu et al., 2021). Our results show a general match of present-day slab structures beneath Eastern Asia with those inferred from seismic tomography observations, especially in the upper-to-middle mantle depths. Although a quantitative comparison between the predicted and seismic imaged slabs (e.g., Seton et al., 2015; Flament, 2019) could be helpful, our models still show a first-order consistency which is valid for the current study. There may be uncertainties in plate reconstruction models which could affect the structures of the subducting slabs, such as the poor-constrained motion history of the Izanagi plate, but the one we used here (i.e., Matthews et al., 2016) has reproduced the consistent mantle structures. In addition, in our models we use a uniform upper thermal boundary layer as the initial condition and let the models run fully self-consistent, rather than constraining the thermal structure of the lithosphere and the slabs in the upper mantle through data assimilation (e.g., Bower et al., 2015). However, we have performed models with different initial conditions which show little difference in the time evolution of slab structures in the upper mantle beneath Eastern Asia. Moreover, although double-sided subduction could potentially develop in our models which can be restrained by imposing slab dip in shallower depths (e.g., Bower et al., 2015) or a lubricating oceanic crust (e.g., Hu et al., 2018) as in data assimilation method, our models have revealed a more fully self-consistent mantle convection system without introducing other artificial constraints.

We show that the subducted slabs of the Izanagi plate and the Pacific plate beneath Eastern Asia exhibit different slab behaviors since the Mesozoic time. In the Mesozoic, the subducted slabs of the Izanagi plate do not always directly pass through the 670-km boundary into the lower mantle, and sometimes they tend to first accumulate and spread laterally above the boundary before falling to the lower mantle. In some cases, the stagnation of the slabs promotes the formation of a big mantle wedge beneath the overriding lithosphere. In particular, we find that multiple stagnation events of the Izanagi slabs have occurred in the Mesozoic time beneath Eastern Asia, e.g., beneath the NCC at 150–120 Ma and 90–70 Ma, and beneath the SCB at 160–140 Ma and 90–70 Ma, respectively. These two stages of slab stagnation events correspond to the episodes of trench advance (170–150 Ma and 120–100 Ma) in plate reconstruction model (Matthews et al., 2016). Since the Mesozoic time, multistaged magmatism activities have been identified in Eastern Asia. For example, three main periods of magmatism formation are identified in NCC, e.g., 160–140 Ma, 130–110 Ma and 100–80 Ma, with the peak period at ∼125 Ma (Zhai et al., 2016). Regional igneous rocks in the SCB have been divided into four major emplacement episodes: 190–175 Ma, 165–155 Ma, 145–125 Ma and 105–95 Ma (Cao et al., 2021). The magmatism from Jurassic to Early Cretaceous in NCC has been attributed to flat-slab subduction followed by a rollback of the Izanagi plate during the late Mesozoic (Wu et al., 2019). However, the time periods of slab stagnation events beneath Eastern Asia in our numerical models also generally coincide with the episodes of magmatic activities in this region. The formation of such mantle wedge, which may lead to slab-triggered water release and vigorous mantle convective flow that could potentially weaken the lithosphere (e.g., Wang et al., 2016; Yang et al., 2017), is therefore supposed to be responsible for the Mesozoic evolution of the overriding plate that are characterized by magmatic activities and/or lithospheric instabilities in Eastern Asia (Ma and Xu, 2021). The intraplate magmatism along East Asia has disappeared in the Late Mesozoic (∼80 Ma) and a flat slab model has been proposed to explain this magmatic termination (Peng et al., 2021b). In contrast, our numerical results still show a stagnation event at 90–70 Ma, but the resultant mantle wedge seems to be smaller than that in the middle Mesozoic (∼130 Ma) which might be insufficient to cause large-scale magmatism. Moreover, the mantle wedge is later destroyed due to the subsequent sinking of the stagnant slabs within a duration timescale of ∼10–20 Myrs, suggesting a transient nature of slab stagnation which has also been revealed by previous global mantle convection studies (Mao and Zhong, 2018). To the middle Cenozoic, most slabs of the Izanagi plate have fallen into the lower mantle and currently, there are few remnants preserved in the upper mantle beneath the Eastern Asia region, which is consistent with recent studies (Liu et al., 2017). In comparison, the subducted slabs of the Pacific plate mostly lie horizontally in the mantle transition zone depth beneath Eastern Asia, which is supposed to be responsible for the Cenozoic magmatic activities and lithosphere evolution in this area (e.g., Zhu et al., 2012; Liu et al., 2017).

In this study, we apply a weak viscosity layer associated with the 670-km phase transition in our models, in which the stagnant slabs remain coherent as they spread horizontally along the 670-km boundary rather than dissolve in the weak viscosity layer (Supplementary Figure S5). Although some regional models have suggested that the weak layer has a small effect on slab stagnation (e.g., Li et al., 2019), global models have demonstrated that the weak layer is crucial to explain not only the stagnant slabs in the transition zone (e.g., East Asia) but also other slab structures (e.g., North America) that are observed in seismic tomography studies (Mao and Zhong, 2018). Recent numerical studies, on the other hand, have argued that a pressure-driven Cenozoic mantle wind is the dominant mechanism for the formation of the stagnant slabs beneath East Asia (Peng et al., 2021a). The linkage between the Izanagi plate subduction and the tectonic evolution of East Asia is still controversial because the slab dynamics of the Izanagi plate have not been well resolved. Previous studies have proposed flat-slab subduction of the Izanagi plate to be responsible for the East Asian lithospheric evolution in Mesozoic time (Liu et al., 2021; Peng et al., 2021b). Our numerical results, which reveal multiple stagnation events of the Izanagi slab in Mesozoic time as an alternative mechanism, provide new insight into the Izanagi slab dynamics and its relation with the East Asian tectonics. On the other hand, there are some potential limitations in our models. We assume an incompressible mantle under the Boussinesq approximation which may overpredict slab volumes particularly in the lower mantle (Flament, 2019). In this study, we focus on slab dynamics in the upper mantle and we suspect mantle compressibility may have a small effect on upper mantle slab dynamics.

Conclusion

Our 3-D global geodynamic model reproduces the typical slab structures as observed by seismic tomographic imaging beneath Eastern Asia. Our model also shows that multiple slab stagnation events, which facilitate the formation of big mantle wedge beneath Eastern Asia, have developed during the subduction of the Izanagi plate in Mesozoic time with a timescale on the order of tens of million years. Importantly, the duration times of those slab stagnation events generally coincide with the episodes of magmatic activities in Eastern Asia. Our model provides an attractive relationship between the subduction of the Izanagi-Pacific plate and the tectonic evolution of Eastern Asia.

Data Availability Statement

The seismic tomography models GAP_P4 and TX2019slab are available at the websites http://d-earth.jamstec.go.jp/ and IRIS website http://ds.iris.edu/ds/products/emc-tx2019slab/, respectively. The numerical models are performed using modified CitcomS (www.geodynamics.org) that is available at https://github.com/bchwu/plate_tracer.

Author Contributions

BW, JH, and YW contributed to the conception and design of the study. BW, JH, and YW designed the model. BW ran the model and analyzed the results. All authors contributed to writing the manuscript and approved the submitted version.

Funding

The work is supported by National Natural Science Foundation of China (42074105; 41674096; 41764004) and the National Key Scientific and Technological Infrastructure project “Earth System Science Numerical Simulator Facility” (EarthLab).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank Editor Lijun Liu, Jiashun Hu and Nicolas Flament for their insightful and constructive reviews and comments that significantly improve the paper. We also thank Profs. Mingming Li and Yang Li for helping improve the paper a lot. Figures are drawn using the Generic Mapping Tools (GMT, www.soest.hawaii.edu/gmt/).

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/feart.2022.829163/full#supplementary-material

References

Bower, D. J., Gurnis, M., and Flament, N. (2015). Assimilating Lithosphere and Slab History in 4-D Earth Models. Phys. Earth Planet. Interiors 238, 8–22. doi:10.1016/j.pepi.2014.10.013

CrossRef Full Text | Google Scholar

Cao, X., Flament, N., Li, S., and Müller, R. D. (2021). Spatio-temporal Evolution and Dynamic Origin of Jurassic-Cretaceous Magmatism in the South China Block. Earth-Science Rev. 217, 103605. doi:10.1016/j.earscirev.2021.103605

CrossRef Full Text | Google Scholar

Christensen, U. R., and Yuen, D. A. (1985). Layered Convection Induced by Phase Transitions. J. Geophys. Res. 90 (Nb12), 10291–10300. doi:10.1029/JB090iB12p10291

CrossRef Full Text | Google Scholar

Domeier, M., and Torsvik, T. H. (2014). Plate Tectonics in the Late Paleozoic. Geosci. Front. 5 (3), 303–350. doi:10.1016/j.gsf.2014.01.002

CrossRef Full Text | Google Scholar

Flament, N. (2019). Present-day Dynamic Topography and Lower-Mantle Structure from Palaeogeographically Constrained Mantle Flow Models. Geophys. J. Int. 216 (3), 2158–2182. doi:10.1093/gji/ggy526

CrossRef Full Text | Google Scholar

Fukao, Y., Obayashi, M., Nakakuki, T., and Grp, D. S. P. (2009). Stagnant Slab: A Review. Annu. Rev. Earth Planet. Sci. 37, 19–46. doi:10.1146/annurev.earth.36.031207.124224

CrossRef Full Text | Google Scholar

Fukao, Y., and Obayashi, M. (2013). Subducted Slabs Stagnant above, Penetrating through, and Trapped below the 660 Km Discontinuity. J. Geophys. Res. Solid Earth 118 (11), 5920–5938. doi:10.1002/2013jb010466

CrossRef Full Text | Google Scholar

Gibert, B., Seipold, U., Tommasi, A., and Mainprice, D. (2003). Thermal Diffusivity of Upper Mantle Rocks: Influence of Temperature, Pressure, and the Deformation Fabric. J. Geophys. Res. 108, 63–69. doi:10.1029/2002JB002108

CrossRef Full Text | Google Scholar

Hu, J., Liu, L. J., Liu, L., and Zhou, Q. (2018). Reproducing Past Subduction and Mantle Flow Using High-Resolution Global Convection Models. Earth Planet. Phys. 2 (3), 189–207. doi:10.26464/epp2018019

CrossRef Full Text | Google Scholar

Katsura, T., Yokoshi, S., Kawabe, K., Shatskiy, A., Manthilake, M. A. G. M., Zhai, S., et al. (2009). P-V-T Relations of the MgSiO3perovskite Determined by In Situ X-ray Diffraction Using a Large-Volume High-Pressure Apparatus. Geophys. Res. Lett. 36, L01305. doi:10.1029/2009GL039318

CrossRef Full Text | Google Scholar

Lay, T., Hernlund, J., and Buffett, B. A. (2008). Core-mantle Boundary Heat Flow. Nat. Geosci 1, 25–32. doi:10.1038/ngeo.2007.44

CrossRef Full Text | Google Scholar

Li, S., Suo, Y., Li, X., Zhou, J., Santosh, M., Wang, P., et al. (2019). Mesozoic Tectono-Magmatic Response in the East Asian Ocean-Continent Connection Zone to Subduction of the Paleo-Pacific Plate. Earth-Science Rev. 192, 91–137. doi:10.1016/j.earscirev.2019.03.003

CrossRef Full Text | Google Scholar

Li, X.-Y., Zheng, J.-P., Ma, Q., Xiong, Q., Griffin, W. L., and Lu, J.-G. (2014a). From Enriched to Depleted Mantle: Evidence from Cretaceous Lamprophyres and Paleogene Basaltic Rocks in Eastern and central Guangxi Province, Western Cathaysia Block of South China. Lithos 184-187, 300–313. doi:10.1016/j.lithos.2013.10.039

CrossRef Full Text | Google Scholar

Li, X.-Y., Zheng, J.-P., Sun, M., Pan, S.-K., Wang, W., and Xia, Q.-K. (2014b). The Cenozoic Lithospheric Mantle beneath the interior of South China Block: Constraints from Mantle Xenoliths in Guangxi Province. Lithos 210-211, 14–26. doi:10.1016/j.lithos.2014.09.028

CrossRef Full Text | Google Scholar

Li, Z.-H., Gerya, T., and Connolly, J. A. D. (2019). Variability of Subducting Slab Morphologies in the Mantle Transition Zone: Insight from Petrological-Thermomechanical Modeling. Earth-Science Rev. 196, 102874. doi:10.1016/j.earscirev.2019.05.018

CrossRef Full Text | Google Scholar

Liu, C.-Z., Liu, Z.-C., Wu, F.-Y., and Chu, Z.-Y. (2012). Mesozoic Accretion of Juvenile Sub-continental Lithospheric Mantle beneath South China and its Implications: Geochemical and Re-os Isotopic Results from Ningyuan Mantle Xenoliths. Chem. Geology. 291, 186–198. doi:10.1016/j.chemgeo.2011.10.006

CrossRef Full Text | Google Scholar

Liu, L., Peng, D., Liu, L., Chen, L., Li, S., Wang, Y., et al. (2021). East Asian Lithospheric Evolution Dictated by Multistage Mesozoic Flat-Slab Subduction. Earth-Science Rev. 217, 103621. doi:10.1016/j.earscirev.2021.103621

CrossRef Full Text | Google Scholar

Liu, S., Ma, P., Zhang, B., and Gurnis, M. (2021). The Horizontal Slab beneath East Asia and its Subdued Surface Dynamic Response. J. Geophys. Res. Solid Earth 126 (3), e2020JB021156. doi:10.1029/2020JB021156

CrossRef Full Text | Google Scholar

Liu, X., Zhao, D., Li, S., and Wei, W. (2017). Age of the Subducting Pacific Slab beneath East Asia and its Geodynamic Implications. Earth Planet. Sci. Lett. 464, 166–174. doi:10.1016/j.epsl.2017.02.024

CrossRef Full Text | Google Scholar

Lu, C., Grand, S. P., Lai, H., and Garnero, E. J. (2019). TX2019slab: A New P and S Tomography Model Incorporating Subducting Slabs. J. Geophys. Res. Solid Earth 124 (11), 11549–11567. doi:10.1029/2019JB017448

CrossRef Full Text | Google Scholar

Ma, P., Liu, S., Gurnis, M., and Zhang, B. (2019). Slab Horizontal Subduction and Slab Tearing beneath East Asia. Geophys. Res. Lett. 46 (10), 5161–5169. doi:10.1029/2018gl081703

CrossRef Full Text | Google Scholar

Ma, Q., and Xu, Y.-G. (2021). Magmatic Perspective on Subduction of Paleo-Pacific Plate and Initiation of Big Mantle Wedge in East Asia. Earth-Science Rev. 213, 103473. doi:10.1016/j.earscirev.2020.103473

CrossRef Full Text | Google Scholar

Mao, W., and Zhong, S. (2018). Slab Stagnation Due to a Reduced Viscosity Layer beneath the Mantle Transition Zone. Nat. Geosci 11 (11), 876–881. doi:10.1038/s41561-018-0225-2

CrossRef Full Text | Google Scholar

Matthews, K. J., Maloney, K. T., Zahirovic, S., Williams, S. E., Seton, M., and Müller, R. D. (2016). Global Plate Boundary Evolution and Kinematics since the Late Paleozoic. Glob. Planet. Change 146, 226–250. doi:10.1016/j.gloplacha.2016.10.002

CrossRef Full Text | Google Scholar

McNamara, A. K., and Zhong, S. (2004). Thermochemical Structures within a Spherical Mantle: Superplumes or Piles. J. Geophys. Res. 109 (B7), B07402. doi:10.1029/2003jb002847

CrossRef Full Text | Google Scholar

Meng, Q.-R. (2003). What Drove Late Mesozoic Extension of the Northern China–Mongolia Tract? Tectonophysics 369 (3), 155–174. doi:10.1016/S0040-1951(03)00195-1

CrossRef Full Text | Google Scholar

Mitrovica, J. X., and Forte, A. M. (2004). A New Inference of Mantle Viscosity Based upon Joint Inversion of Convection and Glacial Isostatic Adjustment Data. Earth Planet. Sci. Lett. 225 (1-2), 177–189. doi:10.1016/j.epsl.2004.06.005

CrossRef Full Text | Google Scholar

Müller, R. D., Cannon, J., Qin, X., Watson, R. J., Gurnis, M., Williams, S., et al. (2018). GPlates: Building a Virtual Earth through Deep Time. Geochem. Geophys. Geosyst. 19 (7), 2243–2261. doi:10.1029/2018GC007584

CrossRef Full Text | Google Scholar

Müller, R. D., Zahirovic, S., Williams, S. E., Cannon, J., Seton, M., Bower, D. J., et al. (2019). A Global Plate Model Including Lithospheric Deformation along Major Rifts and Orogens since the Triassic. Tectonics 38 (6), 1884–1907. doi:10.1029/2018TC005462

CrossRef Full Text | Google Scholar

Peng, D., Liu, L., Hu, J., Li, S., and Liu, Y. (2021a). Formation of East Asian Stagnant Slabs Due to a Pressure-Driven Cenozoic Mantle Wind Following Mesozoic Subduction. Geophys. Res. Lett. 48, e2021GL094638. doi:10.1029/2021GL094638

CrossRef Full Text | Google Scholar

Peng, D., Liu, L., and Wang, Y. (2021b). A Newly Discovered Late-Cretaceous East Asian Flat Slab Explains its Unique Lithospheric Structure and Tectonics. JGR Solid Earth 126 (10), e2021JB022103. doi:10.1029/2021JB022103

CrossRef Full Text | Google Scholar

Seton, M., Flament, N., Whittaker, J., Müller, R. D., Gurnis, M., and Bower, D. J. (2015). Ridge Subduction Sparked Reorganization of the Pacific Plate-Mantle System 60-50 Million Years Ago. Geophys. Res. Lett. 42, 1732–1740. doi:10.1002/2015GL063057

CrossRef Full Text | Google Scholar

Sorokin, A. A., Zaika, V. A., Kovach, V. P., Kotov, A. B., Xu, W., and Yang, H. (2020). Timing of Closure of the Eastern Mongol-Okhotsk Ocean: Constraints from U-Pb and Hf Isotopic Data of Detrital Zircons from Metasediments along the Dzhagdy Transect. Gondwana Res. 81, 58–78. doi:10.1016/j.gr.2019.11.009

CrossRef Full Text | Google Scholar

Sun, W., Ding, X., Hu, Y.-H., and Li, X.-H. (2007). The golden Transformation of the Cretaceous Plate Subduction in the West Pacific. Earth Planet. Sci. Lett. 262 (3), 533–542. doi:10.1016/j.epsl.2007.08.021

CrossRef Full Text | Google Scholar

Wang, F., Zhou, X.-H., Zhang, L.-C., Ying, J.-F., Zhang, Y.-T., Wu, F.-Y., et al. (2006). Late Mesozoic Volcanism in the Great Xing'an Range (NE China): Timing and Implications for the Dynamic Setting of NE Asia. Earth Planet. Sci. Lett. 251 (1), 179–198. doi:10.1016/j.epsl.2006.09.007

CrossRef Full Text | Google Scholar

Wang, Z., Kusky, T. M., and Capitanio, F. A. (2016). Lithosphere Thinning Induced by Slab Penetration into a Hydrous Mantle Transition Zone. Geophys. Res. Lett. 43 (22), 11567–11577. doi:10.1002/2016GL071186

CrossRef Full Text | Google Scholar

Wang, Y., Zhou, L., Liu, S., Li, J., and Yang, T. (2018). Post-cratonization Deformation Processes and Tectonic Evolution of the North China Craton. Earth-Science Rev. 177, 320–365. doi:10.1016/j.earscirev.2017.11.017

CrossRef Full Text | Google Scholar

Wu, F., Lin, J., Wilde, S., Zhang, X., and Yang, J. (2005). Nature and Significance of the Early Cretaceous Giant Igneous Event in Eastern China. Earth Planet. Sci. Lett. 233 (1-2), 103–119. doi:10.1016/j.epsl.2005.02.019

CrossRef Full Text | Google Scholar

Wu, F.-Y., Yang, J.-H., Xu, Y.-G., Wilde, S. A., and Walker, R. J. (2019). Destruction of the North China Craton in the Mesozoic. Annu. Rev. Earth Planet. Sci. 47 (1), 173–195. doi:10.1146/annurev-earth-053018-060342

CrossRef Full Text | Google Scholar

Xu, Y., Li, H., Hong, L., Ma, L., Ma, Q., and Sun, M. (2018). Generation of Cenozoic Intraplate Basalts in the Big Mantle Wedge under Eastern Asia. Sci. China Earth Sci. 61 (7), 869–886. doi:10.1007/s11430-017-9192-y

CrossRef Full Text | Google Scholar

Xu, Z., and Zheng, Y.-F. (2017). Continental Basalts Record the Crust-Mantle Interaction in Oceanic Subduction Channel: A Geochemical Case Study from Eastern China. J. Asian Earth Sci. 145, 233–259. doi:10.1016/j.jseaes.2017.03.010

CrossRef Full Text | Google Scholar

Yang, J., Zhao, L., Kaus, B. J. P., Lu, G., Wang, K., and Zhu, R. (2018). Slab-triggered Wet Upwellings Produce Large Volumes of Melt: Insights into the Destruction of the North China Craton. Tectonophysics 746, 266–279. doi:10.1016/j.tecto.2017.04.009

CrossRef Full Text | Google Scholar

Yang, Y.-T. (2013). An Unrecognized Major Collision of the Okhotomorsk Block with East Asia during the Late Cretaceous, Constraints on the Plate Reorganization of the Northwest Pacific. Earth-Science Rev. 126, 96–115. doi:10.1016/j.earscirev.2013.07.010

CrossRef Full Text | Google Scholar

Young, A., Flament, N., Maloney, K., Williams, S., Matthews, K., Zahirovic, S., et al. (2019). Global Kinematics of Tectonic Plates and Subduction Zones since the Late Paleozoic Era. Geosci. Front. 10 (3), 989–1013. doi:10.1016/j.gsf.2018.05.011

CrossRef Full Text | Google Scholar

Zhai, M., Zhang, Y., Zhang, X., Wu, F., Peng, P., Li, Q., et al. (2016). Renewed Profile of the Mesozoic Magmatism in Korean Peninsula: Regional Correlation and Broader Implication for Cratonic Destruction in the North China Craton. Sci. China Earth Sci. 59 (12), 2355–2388. doi:10.1007/s11430-016-0107-0

CrossRef Full Text | Google Scholar

Zhao, D., Yu, S., and Ohtani, E. (2011). East Asia: Seismotectonics, Magmatism and Mantle Dynamics. J. Asian Earth Sci. 40 (3), 689–709. doi:10.1016/j.jseaes.2010.11.013

CrossRef Full Text | Google Scholar

Zhong, S., and Gurnis, M. (1994). Role of Plates and Temperature-dependent Viscosity in Phase Change Dynamics. J. Geophys. Res. 99 (B8), 15903–15917. doi:10.1029/94JB00545

CrossRef Full Text | Google Scholar

Zhong, S., McNamara, A., Tan, E., Moresi, L., and Gurnis, M. (2008). A Benchmark Study on Mantle Convection in a 3-D Spherical Shell Using CitcomS. Geochem. Geophys. Geosyst. 9, a–n. doi:10.1029/2008gc002048

CrossRef Full Text | Google Scholar

Zhu, R., and Xu, Y. (2019). The Subduction of the West Pacific Plate and the Destruction of the North China Craton. Sci. China Earth Sci. 62 (9), 1340–1350. doi:10.1007/s11430-018-9356-y

CrossRef Full Text | Google Scholar

Zhu, R., Xu, Y., Zhu, G., Zhang, H., Xia, Q., and Zheng, T. (2012). Destruction of the North China Craton. Sci. China Earth Sci. 55 (10), 1565–1587. doi:10.1007/s11430-012-4516-y

CrossRef Full Text | Google Scholar