Change in Subduction Dip Angle of the Indian Continental Lithosphere Inferred From the Western Himalayan Eclogites

The occurrence of ultrahigh-pressure (UHP) and high-pressure (HP) rocks in the Himalayan orogen has been conventionally attributed to the different subduction dip angles along the strike. The western Himalayan UHP eclogites point to a steep continental subduction in the Eocene. The present-day geophysical data show low subduction dip angles of the Indian lithosphere beneath southern Tibet and Karakoram, implying that a shift from steep to low-angle subduction probably happened in the western Himalaya. However, the timing and mechanism of such a subduction-angle change are still unknown. Here we present a combined analysis of zircon geochronology and geochemistry of eclogites and gneiss in the Stak massif, western Himalaya. Metamorphic zircons equilibrated with garnet and omphacite show flat heavy rare earth element patterns without Eu anomalies and, thus, yield similar eclogite-facie ages of ca. 31 Ma. The Stak HP eclogite-facie metamorphism is at least 15 Ma younger than those measured in the western Himalayan UHP eclogites, but broadly contemporaneous with other Himalayan HP rocks. Therefore, all the Himalayan HP rocks record higher peak geothermal gradients and younger ages than those of the UHP rocks. Our new data, combined with the magmatic lull observed in the Kohistan–Ladakh–Gangdese arc and with the convergent rate of the Indian plate, suggest a change in subduction dip angle over time. Consequently, we suggest that the entire Indian continental lithosphere experienced an approximately coherent shift from steep to low-angle subduction after the breakoff of the Neo-Tethyan slab since the middle Eocene. This critical change in subduction geometry is interpreted to be responsible for the transition from continental subduction to collision dynamics.


INTRODUCTION
The Indo-Asia collision after the closure of the Neo-Tethys produced the present-day largest ongoing continent-continent collisional orogen of the Earth, the Himalayan orogen (Yin, 2006). Because subduction angle, convergence rate, and slab rollback/breakoff strongly influence the spatial and temporal distribution of metamorphic and magmatic rocks (Kay and Coira, 2009;Paterson and Ducea, 2015), the subduction style of the orogen can be revealed by studying these aspects.
High-pressure (HP) metamorphic rocks (blueschist, eclogite, HP granulite facies rocks, and HP garnet amphibolite facies rocks) are exposed along the Himalayan orogen ( Figure 1). However, ultrahigh-pressure (UHP) eclogites only occur at Kaghan and Tso Morari in the western syntaxis (O'Brien et al., 2001;Sachan et al., 2004). The UHP eclogites in the western Himalaya were subducted earlier (53-46 Ma) and exhumed faster (45-40 Ma) than the HP rocks (buried at 38-15 Ma and exhumed at 25-13 Ma) in the central and eastern Himalaya (Parrish et al., 2006;Zhang et al., 2010;Donaldson et al., 2013;Rehman, 2019;Wang et al., 2021), possibly reflecting different subduction dip angles or depths of the Indian continental slab along strike . Therefore, it is proposed that the Indo-Asia plates in the western syntaxis collided earlier followed by steep subduction, and that those in the central and eastern Himalaya collided later with low-angle subduction (Chemenda et al., 2000;Guillot et al., 2008;Zhang et al., 2015). Numerical models imply that a subduction angle decreases over time due to continental lithosphere buoyancy (Duretz and Gerya, 2013;Magni et al., 2017). Geophysical data point to the present-day low subduction angles (or underthrusting) of the Indian continental lithosphere beneath east and west Tibet and Karakoram (Hazarika et al., 2017;Parsons et al., 2020). If correct, the Indian continental slab underwent a shift from steep subduction to low-angle subduction or underthrusting in the western Himalaya. However, such a subduction-angle change has not been evidenced by the rocks in this region, making the timing and mechanism of such a process poorly resolved.
In addition to UHP eclogites, HP eclogites were discovered in the Stak massif in the western syntaxis (Le Fort et al., 1997). These HP-UHP eclogites record clockwise P-T paths (Wilke et al., 2010;Lanari et al., 2013;St-Onge et al., 2013) and, thus, define a large HP-UHP province Lanari et al., 2013). However, the Stak eclogite suffered a greater degree of overprint and recorded younger metamorphic zircon ages of 32 Ma (Kouketsu et al., 2016) than the UHP eclogites with peak metamorphic ages of 53-46 Ma (e.g., Rehman, 2019). Is this a reflection of zircon recrystallization during exhumation in Stak, or a real age of eclogite formation? It remains uncertain whether the western Himalayan HP-UHP eclogites were formed synchronously at different subducting depths and then experienced different exhumation histories or were formed by diachronous subduction of various continental slices.
Here we show, based on petrology, zircon geochronology, and geochemistry, that both the Stak eclogites and their country gneiss share the same peak metamorphic age of ca. 31 Ma. This age is younger than that of the UHP eclogites, but coeval to a magmatic lull in the entire Himalayan orogen. Steep subduction of the continental lithosphere, driven by downgoing oceanic lithosphere, would generate UHP metamorphism in the subducting crust ( Figure 2A). However,  Frontiers in Earth Science | www.frontiersin.org February 2022 | Volume 9 | Article 790999 3 the continental lithosphere, without driven force from subducted oceanic lithosphere, would undergo low-angle subduction or underthrusting  and, thus, more likely produce HP metamorphism ( Figure 2B). Based on geological and geophysical records, we further infer a uniform change in the subduction dip angle over 2,500 km throughout the Himalayan belt since the middle Eocene.

GEOLOGICAL SETTING AND SAMPLES
The Stak massif is located in the north Indian continental margin, northeast of the Nanga Parbat-Haramosh massif (NPHM), southwest of the Ladakh arc ( Figure 1B), and close to the Main Mantle Thrust (MMT). This massif consists of gneisses, schists, metabasites, with minor marbles, and develops fold and imbricate structures . Eclogites with extensive retrogression in the Stak massif occur as boudins or dikes in the gneisses. A detailed petrological study on the Stak eclogites yields HP peak conditions at~750°C and~2.5 GPa and retrograde conditions at 650°C-700°C and 0.9-1.6 GPa (Lanari et al., 2013). This retrograde temperature is slightly higher than those (<600°C at 1.0-1.7 GPa) of the Kaghan and Tso Morari UHP eclogites (Wilke et al., 2010;Wilke et al., 2015).
The peak metamorphic age of the Stak eclogite is still unknown. Sensitive high-resolution ion microprobe (SHRIMP) zircon U-Pb data yielded scattered ages between 70 and 50 Ma (Riel et al., 2008). However, Kouketsu et al. (2016) used the same method to show that a small cluster of low Th/U (<0.03) and Yb (<10 ppm) zircon ages are concentrated between 36 and 28 Ma with a lower intercept age of~32 Ma. This age was interpreted as recrystallization after eclogite-facie metamorphism, which was possibly induced by heating from nearby NPHM at lower crustal levels (Kouketsu et al., 2016).
The present study focuses on the eclogite boudins (16PK190 and 16PK194) and country gneiss (16PK181) in the Stak massif. Both eclogite samples contain an eclogite-facies assemblage of garnet, omphacite, quartz, and rutile ( Figure 3). Garnet is surrounded by amphibole and plagioclase kelyphite and has a compositionally homogeneous core (Alm 55-57 Prp 7-13 Grs 24-28 Sps 1-2 ) (see the Supplementary Table S1) with a narrow retrograde zoning rim. Omphacite with a 23.4-30.0 mol.% jadeite component (Supplementary Table S2) is partially replaced by fine-grained symplectite of diopside + plagioclase + amphibole. The symplectite is rimmed by coarser (200-500 μm) amphibole and biotite. The eclogite sample zircons occur as inclusions within garnet and symplectite, and as an intergranular phase in the matrix ( Figure 4). The zircon grains in the symplectite are larger (~100 μm) than the plagioclase and diopside grains, and show disequilibrium texture with these phases. The gneiss (16PK181) contains garnet, muscovite, biotite, plagioclase, and quartz with minor rutile and zircon ( Figure 3). The garnet is surrounded by plagioclase and mica, and could be a relict eclogite-facie phase. The gneiss sample zircons exist as matrix phases and inclusions within garnet and muscovite ( Figure 4).

Mineral Major Elements
Major elements of rock-forming minerals were analyzed using an electron microprobe analyzer (EPMA, JEOL-JXA8100) at the IGG-CAS. The operating conditions were a 15-kV accelerating voltage, a 20-nA beam current, and 3-μm spot diameter. The counting time was 20 s at peak and 10 s at the lower and upper background positions, respectively. All data were corrected online using a modified ZAF (atomic number, absorption, and fluorescence) correction procedure. The detection limits (1σ) were in the range of 0.008-0.02wt%. The precision of the major element analysis was better than 1.0%.

Zircon U-Pb Isotopes and Trace Elements
Zircon grains were separated from samples by using conventional heavy liquid and magnetic separation techniques. Approximately 200 zircons of each sample were handpicked out from under a binocular microscope and mounted in epoxy resin together with zircon U-Pb age reference materials. Then the epoxy mounts were polished to expose the interior of the crystals. All grains were photographed in reflected and transmitted light under petrographic microscope to avoid fissures and inclusions. Cathodoluminescence (CL) and back-scattered electron (BSE) images were obtained using field emission scanning electron microscope (Nova NanoSEM 450) equipped with Gatan MomoCL4 at the Institute of Geology and Geophysics, Chinese Academy of Sciences (IGG-CAS). Mineral inclusions in zircon were identified using Raman spectroscopy at the IGG-CAS. Homogeneous zircons without inclusions or fissures were chosen for U-Pb dating and trace element analyses. Zircon U-Pb ages were determined by using secondary ion mass spectrometry (SIMS) and laser ablation-inductively coupled plasma-mass spectrometry (LA-ICP-MS). Zircon trace elements were acquired in the same run of age dating by LA-ICP-MS.

Secondary Ion Mass Spectrometry Zircon U-Pb Dating
Zircon U, Th, and Pb isotope analyses were performed with a Cameca IMS-1280HR SIMS at the IGG-CAS. Detailed instrumental parameters and analytical procedures are described by Li et al. (2009). The analytical spot is about 20 × 30 μm in size. Plešovice zircon standard ( 206 Pb/ 238 U age of 337.3 Ma; Sláma et al., 2008) was interspersed with unknown grains. Nonradiogenic 204 Pb was used for common Pb correction of measured compositions, and the present-day crustal common Pb composition is used in the model of Stacey and Kramers (1975). An inhouse zircon standard Qinghu was alternately analyzed together with other unknown zircons for the purpose of monitoring external uncertainties of SIMS U-Pb zircon dating calibrated by the Plešovice standard. Excel and the add-in Isoplot 2.49 program (Ludwig, 2001) were used for data calculation. Uncertainties on individual analyses in Supplementary Table S3 are reported at 1σ level. The weighted mean U-Pb ages are quoted with 95% confidence interval. Seven measurements on Qinghu zircon yielded a Concordia age of 159.7 ± 1.8 Ma, which is identical within errors with the recommended value of 159.5 ± 0.2 Ma (Li et al., 2013).

Laser Ablation-Inductively Coupled Plasma-Mass Spectrometry Zircon U-Pb Dating and Trace Element Analyses
U-Pb dating combined with in situ trace element analyses of zircons were carried out in a single run by Agilent 7500a ICP-MS instrument equipped with Geolas-193 UV laser ablation system at the State Key Laboratory of Continental Dynamics, Northwest University, China. Operating conditions and data processing are described by Liu et al. (2007) in detail. The spot diameter is 32 μm for two eclogites (16PK190 and 16PK196) and 44 μm for gneiss (16PK181) with a laser repetition rate of 6 Hz. Helium was used as the carrier gas. Laboratory standards (GJ-1, 91500, NIST 610) were interspersed with unknown grains. U-Th-Pb isotope ratios and trace element contents were calculated using the GLITTER 4.0 program (Macquarie University). Harvard zircon 91500 was used as external standard. To monitor the external uncertainties of U-Pb zircon dating calibrated against 91500, a second zircon standard GJ-1 was analyzed as an unknown together with other unknown zircons. Trace element concentrations were calibrated by using 29 Si as internal standard and NIST 610 as external standard. Data reduction was carried out using Isoplot/Excel version 2.49 (Ludwig, 2001). Fifteen measurements on GJ-1 and 32 measurements on 91500 yielded weighted mean 206 Pb/ 238 U ages of 602 ± 3.7 Ma (1σ) and 1,062.5 ± 2.8 Ma (1σ), respectively,

RESULTS
The zircons from the two eclogites are mostly subhedral and consist of two major parts under cathodoluminescence (CL) images ( Figure 5): 1) a core with relatively bright and slightly oscillatory zoning luminescence, and 2) a rim with relatively dark, homogeneous luminescence. Inclusions of diopside, amphibole, plagioclase, quartz, and apatite are common in the zircon cores, while the rims have garnet, rutile, and quartz inclusions without plagioclase and coesite (Supplementary Figure S1). Fourtheen zircons from sample 16PK190 and 10 zircons from sample 16PK194 were dated by the SIMS method. The Th/U ratios of the zircon cores and rims are 2.05-3.17 and ≤0.05, respectively (Supplementary Table S3). The SIMS dating results show that the two eclogites, 16PK190 and 16PK194, have similar 206 Pb/ 238 U weighted mean ages of 284.3 ± 5.2 Ma (MSWD = 1.7, n = 7) and 280.4 ± 7.4 Ma (MSWD = 2.0, n = 5) for the zircon cores and 31.0 ± 0.5 Ma (MSWD = 0.55, n = 6), 31.9 ± 1.4 Ma (MSWD = 2.8, n = 5) for the zircon rims ( Figures 6A, B). The LA-ICP-MS results (see the Supplementary Figure S2) are consistent with those of SIMS. The zircon cores have higher REE concentrations than the rims and are characterized by positive Ce anomalies, negative Eu anomalies, and steep heavy REE patterns ( Figure 6C). The zircon rims for both samples show flat heavy REE patterns with a distinct lack of Eu anomalies ( Figure 6C).
Zircons from the country gneiss 16PK181 are oval or prismatic crystals characterized by weak luminescence and patchy zoning under CL, with rare bright cores ( Figure 5). All the analyses were conducted on the core-free zircons. They show low Th/U ratios (<0.01) (Supplementary Table S3) and REE concentrations ( Supplementary Table S5

Timing of Eclogite-Facies Metamorphism
The 206 Pb/ 238 U weighted mean ages of ca. 284-280 Ma measured in the zircon cores from both eclogites represent the crystallization age of the magmatic protolith, which is demonstrated by their oscillatory zoning, high Th/U ratios, and steep REE patterns. This age is contemporaneous with the protolith formation of the Kaghan and Tso Morari UHP eclogites Frontiers in Earth Science | www.frontiersin.org February 2022 | Volume 9 | Article 790999 6 (Rajkumar, 2015;Rehman et al., 2016), which is ascribed to the eruption of the Panjal Traps in the north Indian continental margin (Shellnutt, 2018). Our results, coupled with the similar whole-rock chemistry and Nd isotopic compositions between the Stak eclogites and the Panjal Traps volcanic rocks (Kouketsu et al., 2017), indicate that the western Himalayan eclogites likely shared the same protoliths during the Permian magmatism of the Panjal Traps.
The rims of the eclogite zircons exhibit metamorphic characteristics and equilibrium texture with HP phases. The date of ca. 31 Ma most likely reflects eclogite-facies zircon recrystallization because, 1) the zircon rims have inclusions of garnet and rutile without plagioclase, 2) the disequilibrium texture of zircon with symplectitic diopside and plagioclase implies zircon inclusions presented in preexisting omphacite (Figure 4), and 3) the REE data indicate recrystallization in the presence of garnet and absence of plagioclase ( Figure 6C). A similar metamorphic age (ca. 32 Ma) reported by Koukestsu et al. (2016) was interpreted as a recrystallization time after eclogitefacies metamorphism. However, the low Yb concentrations in those zircons only reflect equilibration with garnet, which could occur both under eclogite-and granulite-facies conditions. Moreover, metamorphic zircons from the country gneiss show REE characteristics and ages similar to those from the eclogites ( Figures 6C,D). Therefore, the Stak massif was coherently buried beneath the Asian plate and underwent eclogite-facies metamorphism at ca. 31 Ma.
The HP metamorphic time of the Stak massif is at least 15 Ma later than the UHP metamorphic time (53-46 Ma) in the Kaghan and Tso Morari massifs (de Sigoyer et al., 2000;Kaneko et al., 2003;Donaldson et al., 2013). Therefore, the HP-UHP massifs in the western Himalaya were asynchronously buried tõ 2.5-3.8 GPa, indicating a continuous subduction/collision process spanning more than 15 Ma. After the Kaghan and Tso Morari massifs experienced UHP metamorphism (~46 Ma) and rapid exhumation to the middle to lower crustal levels (~44-40 Ma) (de Sigoyer et al., 2000;Parrish et al., 2006;Wilke et al., 2010), the Stak massif was buried to HP eclogitefacie conditions at ca. 31 Ma. Such a continuous burial process recorded by HP-UHP rocks is common in other collisional belts such as Dabie-Sulu (e.g., Zhang et al., 2009), Alps (e.g., Berger and Bousquet, 2008), and the Western Gneiss Region (e.g., Kylander-Clark et al., 2009). In the following section, we will focus on how this process proceeded through time.

Implications for change in subduction dip angle through time
The HP metamorphism of the Stak massif is broadly contemporaneous with that of quartz eclogites, HP granulites, and garnet amphibolites (mostly 40-25 Ma) in the central and eastern Himalayas ( Figure 7A) (Corrie et al., 2010;Zhang et al., 2010;Rubatto et al., 2013;Wang et al., 2021). This indicates that the subducted Indian continental slab probably underwent a coherent process along the Indus-Yarlung Suture Zone since 40 Ma. However, UHP metamorphism spatially restricted to the western Himalaya (Kaghan and Tso Morari) was only recognized before ca. 46 Ma (Kaneko et al., 2003;Donaldson et al., 2013). The discrepancy of HP-UHP timing in the western and central Himalayas could be induced by continental subduction of a jagged Indian margin , by different exhumation styles (residence times at crustal level) of eclogites (O'Brien, 2019;Rehman, 2019), or by changes in subduction dip angle ( Figure 2). However, it is unlikely that an irregular northern margin of India would trigger the long discrepancy (>15 Ma) between UHP and HP peak metamorphic ages. The different P-T-t paths observed between UHP and HP eclogites were interpreted to result from the prolonged residence times in the central and eastern Himalayas but rapid exhumation of UHP rocks in the west (Rehman, 2019). However, the different exhumation styles are not the first-order mechanism responsible for the distinct peak metamorphic ages of   Vance and Harris, 1999) and the Ayilari garnet amphibolite (~32 Ma; Chen et al., 2021) in the western Himalaya ( Figure 1B), indicating that these HP rocks were synchronously buried to shallow depths much later than the UHP rocks. Notably, the old UHP eclogites and young HP rocks record different peak-pressure thermal gradients (temperature change with depth, referred to as metamorphic T/P) (Figure 8), implying a change in the subduction geometry of the Indian continental slab during 46-40 Ma. The subduction of oceanic lithosphere (Neo-Tethys) after the initial India-Asia collision (~60 Ma; e.g., Parsons et al., 2020) would continuously release fluids or melts, generate convective corner flow in the mantle wedge, and, thus, induce magmatism at convergent plate margins ( Figure 7B). The successive steep subduction of Indian continental crust driven by slab pull of the Neo-Tethys underwent UHP metamorphism (low T/P) at 53-46 Ma (Figures 2A and 9A). The steep subduction of the Indian continental slab is also supported by the presence of the UHP eclogite units in the western Himalaya now directly adjacent to the Indus-Tsangpo suture zone (O'Brien, 2019). However, after~40 Ma, the sole occurrence of HP metamorphism (high T/P) along the strike reflects an increase of geothermal gradient that requires a change in subduction dip angle ( Figures 2B and 9B). These young HP rocks most likely resulted from low-angle underthrusting during continental collision Soret et al., 2021). Therefore, the change in subduction dip angle of the Indian continental slab during 46-40 Ma is critical for the transition from continental subduction to collision dynamics.
The change in subduction dip angle can be induced by the continuous subduction of buoyant continental lithosphere (Parsons et al., 2021) or by the breakoff of the oceanic (Neo-Tethyan) lithosphere (Davies and von Blanckenburg, 1995;Magni et al., 2017). Breakoff of the Neo-Tethyan lithosphere would lead to fast exhumation of UHP rocks (Kohn and Parkinson, 2002), uplifting of the buoyant continental lithosphere, and shifting to gently dipping subduction or underthrusting (Davies and von Blanckenburg, 1995; FIGURE 8 | Integrated peak P-T conditions and ages of eclogites and HP granulites in the Himalayan orogen. The age data sources are listed in Supplementary Table S6. Source references for P-T conditions: Annapurna (Ap) (Kohn and Corrie, 2011), Ama Drime (AD) (Li et al., 2018), Arun (Ar) (Corrie et al., 2010), Ayilari (Ay) (Chen et al., 2021), Bhutan (Bh) (Grujic et al., 2011), Dinggye (Dg) (Wang et al., 2017), Kaghan (Kaneko et al., 2003), Sikkim (Sk) (Faak et al., 2012), Stak (Lanari et al., 2013), Tso Morari (Mukherjee et al., 2005, Yadong (Ya) , Namche Barwa Syntaxis (NBS) (Zhang et al., 2015), and Zanskar (Zk) (Vance and Harris, 1999). Frontiers in Earth Science | www.frontiersin.org February 2022 | Volume 9 | Article 790999 8 Chemenda et al., 2000). In addition, slab breakoff mostly occurs at greater depths (~130-240 km) than the base of the overriding lithosphere and would not trigger significant magmatism at convergent margins (e.g., Freeburn et al., 2017), which is consistent with the weak magmatism within the Asian plate at 46-40 Ma ( Figure 7B). There are many lines of evidence supporting the breakoff of the Neo-Tethyan slab around 46-40 Ma after the initial India-Asia collision. 1) According to paleomagnetic reconstructions and tomographical observations, Negredo et al. (2007) estimated the breakoff of oceanic slab occurring at 48-44 Ma. 2) The rapid exhumation of the Kaghan UHP eclogite from~100 to~35 km during 46-44 Ma (Parrish et al., 2006) also agrees with the breakoff model (Kohn and Parkinson, 2002). 3) Based on tomography, the Asian tectonics reconstructions and the Indian plate kinematics, Replumaz et al. (2014) proposed that the major breakoff between India and the Tethys ocean occurred at~45 Ma. However, Parsons et al. (2021) proposed a different timing for slab breakoff based on interpretations of the same tomographic anomalies. 4) The occurrence of 42-40 Ma intraplate-type mafic dykes in eastern Tibet supports the breakoff of subducting Neo-Tethyan slab from the Indian continental slab during the middle Eocene (Xu et al., 2008). 5) The 45 Ma oceanic island basalt-type gabbro in southern Tibet was used to constrain the breakoff time of the Neo-Tethyan slab (Ji et al., 2016). Numerical studies suggest that the timing of slab breakoff after the onset of continental collision varies from 10 to 25 Ma, which is largely affected by the strength and age of the subducting oceanic slab (e.g., van Hunen and Allen, 2011). Our proposed slab breakoff timing (46-40 Ma) for the old Neo-Tethyan slab after the onset of continental collision (~60 Ma) broadly fits this modeled result. Therefore, we suggest that the different T/P and ages between the Himalayan UHP and HP rocks were induced by the breakoff of the Neo-Tethyan slab and subsequent transition from steep continental subduction to low-angle underthrusting (collision) at 46-40 Ma ( Figure 9B).
The synchronous HP metamorphism induced by laterally largescale low-angle underthrusting is also supported by palaeomagnetic and magmatic data. Since the Late Cretaceous, the India-Asia convergence rates in the eastern and western Himalayas show high consistency ( Figure 7C), decreasing from 140-160 to 80-100 mm/year at 52-50 Ma and then decreasing to 40-60 mm/year at ca. 45 Ma (van Hinsbergen et al., 2011). The first deceleration is attributed to the initial India-Asia collision (van Hinsbergen et al., 2011), whereas the second most likely resulted from the loss of oceanic slab pull after breakoff (e.g., Bercovici et al., 2015;Ji et al., 2016) or the India-Asia second collision (45-40 Ma) after an earlier collision (~60 Ma) of either Indian continent + intraoceanic arc or greater India microcontinent + Asia (e.g., Parsons et al., 2020Parsons et al., , 2021. A low convergence rate benefits sufficient heat exchange between the subducting plate and the overlying mantle, thus inducing HP metamorphism with high T/P . In addition, the Himalayan magmatic activities from east to west were systematically complementary to the HP metamorphic events after 40 Ma. The Kohistan-Ladakh-Gangdese arc system and south Karakoram exhibit a coherent magmatic lull at 40-25 Ma ( Figure 7B). The Stak eclogite in the western Himalaya and most of the other HP rocks in the central and eastern Himalaya also formed at this stage ( Figures 7A and 9C). The low-angle underthrusting after slab breakoff would have driven the asthenosphere beneath the Asian plate away and then shielded the active continental margin from convective heat, which would have led to the magmatic lull in the convergent margin (Chung et al., 2005;Ji et al., 2016;Ji et al., 2020). The driving force for the continuous convergence and underthrusting after the breakoff of the Neo-Tethyan slab could be subduction of Australian oceanic lithosphere (e.g., Parsons et al., 2021) or subduction of Indian continental lithosphere (e.g., Capitanio et al., 2010). In addition, the continuous low-angle underthrusting during continental collision at 40-25 Ma potentially enhanced crustal thickening (Soret et al., 2021) and elevated the thermal structure of orogenic wedge due to radiogenic heat production in the thickened crust (Berger et al., 2011), contributing to the high-or ultrahigh-temperature overprints in HP eclogites during exhumation (Wang et al., 2021). In this regard, the change in subduction dip angle of continental slab would have influenced both magmatic activities and peak geothermal gradients of HP-UHP rocks. To summarize, a lateral (>2,500 km) change in the subduction geometry of the Indian continental slab during the middle Eocene would be broadly consistent with the igneous and metamorphic records.

CONCLUSION
The Stak eclogites and their country gneiss in the western Himalayan syntexis provide insight into the change in subduction dip angle of the Indian continental lithosphere. Precise geochronological data show that the Stak massif underwent HP eclogite-facie metamorphism at~31 Ma, which is at least~15 Ma later than the Himalayan UHP metamorphism. Considering the discrepancies in peak ages and geothermal gradients between the Himalayan HP and UHP metamorphic rocks, we suggest that the Indian continental lithosphere underwent a coherent change in subduction dip angle during 46-40 Ma. Based on geological and geophysical evidence, this change in subduction geometry is likely induced by the breakoff of the Neo-Tethyan slab. We suggest that the change in subduction dip angle is a critical process controlling the transition from continental subduction to collision dynamics.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Materials. Further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
SC, YC, and QL performed the analyses. SC contributed to the data sorting and compilation. YC interpreted the data. YC, SC, and SG wrote the manuscript.