Geodynamics of the Broad Triangle Area (active seismic zone) in Asia: Stress field modeling

Journal Pre-proofs Geodynamics of the Broad Triangle Area (active seismic zone) in Asia: Stress field modeling Shangchang Duan, Guiting Hou, Lihui Yang PII: DOI: Reference:

S1367-9120(19)30425-0 https://doi.org/10.1016/j.jseaes.2019.104073 JAES 104073

To appear in:

Journal of Asian Earth Sciences

Received Date: Revised Date: Accepted Date:

23 May 2019 4 October 2019 10 October 2019

Please cite this article as: Duan, S., Hou, G., Yang, L., Geodynamics of the Broad Triangle Area (active seismic zone) in Asia: Stress field modeling, Journal of Asian Earth Sciences (2019), doi: https://doi.org/10.1016/j.jseaes. 2019.104073

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier Ltd.

Geodynamics of the Broad Triangle Area (active seismic zone) in Asia: stress field modeling

Shangchang Duan, Guiting Hou*, Lihui Yang The Key Laboratory of Orogenic Belts and Crustal Evolution, School of Earth and Space Sciences, Peking University, Beijing 100871, China * Corresponding author, E-mail: [email protected] (G. Hou).

Abstract The Broad Triangle Area (BTA) is a significant active seismic zone located mainly in East Asia. The Pamir-Baikal tectonic belt, North-South tectonic belt and Himalayan tectonic belt compose the three boundaries of the BTA. Previous studies suggested that the geodynamics of the BTA may be related to coupling between the Indian and Eurasian plates, yet quantitative analysis is lacking. In this study, we use the elastic finite element method in the ANSYS software to establish a spherical shell model and analyze the stress field of the model. We obtain a best-fit model that matches the geological evidence and seismicity by considering the real geological background; a series of models with different initial settings are also discussed to analyze the boundary influence on the different models. Our modeling results indicate that the India-Eurasia collision, the stable Siberia craton and the blocking effects of the Yangtze and North China cratons play the most important roles in the formation of the BTA. Crustal heterogeneity has little influence on the distribution of seismicity within the BTA. Further, the effects of subduction of the Pacific and Philippine plates may not reach the BTA area, which is different from the results of other studies.

Key words: Active Seismic Zone, Elastic Spherical Shell Modeling, Indian-Eurasian Collision, Strain Energy Density

1. Introduction In addition to the circum-Pacific seismic belt and the Alpine seismic belt, Gutenberg and Richter (1954) indicated that there is another significant seismic activity zone called the Broad Triangle Area in East Asia (Fig. 1). Ma and Zheng (1981) more precisely defined the Broad Triangle Area as delimited by the Pamir-Baikal tectonic belt in the northwest margin of the area, the North-South tectonic belt in the east margin of the area and the Himalayan tectonic belt in the southwest margin of the area (Fig. 1). Fig. 1 Here In several decades, few studies have been devoted specifically to the Broad Triangle Area. However, there are still large numbers of significant related regional background studies that provide some effective ways to understand the formation of the Broad Triangle Area and basic data for further studies. Some studies have revealed that the spatial distribution of epicenters is consistent with other geological features (Fu, 2012). For investigating the indentation tectonics in Central Asia, Cobbold and Davy (1988) found that surface topography exceeds 1 km and crustal thickness exceeds 40 km throughout almost all of an approximately triangular area that is highly coincident with the Broad Triangle Area, and this feature was interpreted to result from the Cenozoic northward motion of India. Feng et al. (2007) found a spatial pattern in which the shape of the gravity low resembles a huge triangle, including the Tibetan Plateau with three vertexes in the Hindu Kush, Baikal and Yunnan of China, which is also coincident with the Broad Triangle Area, and they proposed that the NW margin of this area divides the Eurasian Plate into the stable Russia-Siberia subplate to the north and

the active China-Southeast Asia subplate to the south. It seems that these similar geological features should be related to the tectonic setting of the East Asian continent (Fu, 2012). It has been generally admitted that many of the tectonic features and active structures and the stress field of East Asia could be related to the India-Eurasia collision (e.g., Molnar and Tapponnier, 1975, 1977; Le Pichon et al., 1992; Zhao et al., 1990; Xu, 2001; Zhang et al., 2003), the far-field effect of which has been assumed to reach the Baikal rift (Molnar and Tapponnier, 1975; Jolivet et al., 2013). In addition, Li et al. (2010) suggested that Cenozoic tectonism in China is associated with eastern jumping and rollback of Pacific Plate subduction zones. The dynamic model of Asian continental deformation has long been the object of intense debate. For the ‘edge-driven block models’ (Tapponnier et al., 1979; Peltzer and Saucier, 1996), the deformation is directly controlled by the forces along the boundaries of rigid lithospheric blocks. For the ‘buoyancy-driven continuous deformation models’ (England and Houseman, 1986; Houseman and England, 1993), when considering the viscosity of the material, the crustal thickening and the buoyancy should not be ignored. As numerical techniques have made it possible to explore large-scale crustal deformation, some scientists have achieved quantitative results for East Asian tectonics (e.g., Zhang and Gao, 1989; Zhu et al., 2010; He et al., 2013; Shen et al., 2007; Lei et al., 2009), especially the impact of the India-Eurasia collision (Chen et al., 2011). Further quantitative studies have relied on refined data about the present-day kinematics by Global Positioning System (GPS) investigation (Calais et al., 2003, 2006; Ren and Ma, 2003). To date, previous studies in related fields have led to an important number of publications and have shown that the geodynamics of the Broad Triangle Area may be complex and related to the India-Eurasia collision, the subduction of the Pacific Plate

and Philippine Plate and deep mantle processes (Gao et al., 2010). Although early studies have provided only some primary interpretations regarding the spatial distribution of epicenters in the Broad Triangle Area, there has been no specific effort to construct a numerical model for the mechanism of the Broad Triangle Area. In this study, we focus on the geodynamic mechanism of the Broad Triangle Area, consider a finite element numerical shell model, calculate the strain fields under the constraints of different boundary conditions, and compare the model results to geological evidence to analyze the mechanism of the Broad Triangle Area.

2. Geological background The Broad Triangle Area shows the distribution of historical epicenters in East Asia and lies between the Himalaya and Lake Baikal (Fig. 1). It is also one of the most seismically active zones in the world, including a great number of destructive earthquakes, e.g., the Wenchuan earthquake (Ms=8, 2008.5.12). Although the distribution of epicenters is dispersed on the whole, most of the large earthquakes (Ms>7) have been concentrated in the pre-existing zones of weakness, especially along the margins of the blocks (Zhang et al., 2003; Zhang et al., 2004). As mentioned above, the geodynamic factors giving rise to the Broad Triangle Area may include the India-Eurasia collision, the subduction of the Pacific Plate and Philippine Plate and deep mantle processes. For this reason, it is necessary to expand the scope of the investigation; more specifically, we define the study area to include Central Asia, East Asia and northern Asia, based on several tectonic studies (Fig. 2) (Ren et al., 2013; Wan, 2013; Ren et al., 2003; Pospelov et al., 2016). Fig. 2 Here 2.1. Basement framework of the study area

According to the 1:5,000,000 International Geological Map of Asia project led by Ren et al. (2013), our study area is a composite continent composed of the Siberian craton, North China craton, Yangtze craton, Tarim craton, and two orogenic belts with a number of small blocks (Fig. 2). The Siberian craton occupies the northern part of the main body of the continent and has been assumed to have been one part of the MesoNeoproterozoic supercontinent Rodinia, which was subsequently reassembled into Laurasia (Gladkochub et al., 2006). In contrast, the reported geological and paleomagnetic data demonstrate that the North China craton, Yangtze craton and Tarim craton were part of the northern margin of Gondwana (Veevers, 2004; Cawood et al., 2013; Han et al., 2016; Ren et al.,1999). In the late Paleozoic, the Siberian craton was assembled with the North China and Tarim cratons to form the Central Asian Orogenic Belt (CAOB) (Windley et al., 2007; Xiao et al., 2003; Jahn, 2004). At the end of the Triassic, the North China craton and Yangtze craton were assembled (Zhang et al., 1997; Gao et al., 1995). In the MesoCenozoic, the East Asian continent was assembled with the Indian and Arabian cratons (Ren et al., 2013; Guiraud and Bellion, 1995). The two Phanerozoic orogenic belts immediately adjacent to the rigid cratons, including the Central Asian Orogenic Belt and the Tethyan orogenic belt, were dominated by accretionary tectonics and properly belong to two global tectonic domains: the Paleo-Asian domain that evolved from the Paleo-Asian Ocean dynamic system and the Tethyan domain derived from the Tethyan and Indian Ocean dynamic systems (Ren et al., 2013; Windley et al., 2007; Şengör et al., 1988). Within the “soft” orogenic belts, there are still some small “rigid” blocks, including the Turan-Karakum, Junggar, Qaidam, Qiangtang-Qamdo and Lhasa, which were associated with the closures of the Paleo-Asian Ocean and the Tethys Ocean (Ren et al., 2013; Wan, 2013;

Ren et al., 2003). The cratons, blocks and Phanerozoic orogenic belts formed the typical dual structure of the basement framework of the study area. In addition, our study area includes a small part of the western Pacific active continental margin in the east, which has been thought to belong to the Pacific domain originating from the Paleo-Pacific and Pacific dynamic systems (Ren et al., 2013). 2.2. Boundaries of the study area The Uralide orogen (Fig. 2), which is the geographic and geologic divide between Europe and Asia, marks the Paleozoic collision zone of the East European craton with the Asian collage of terranes and acts as the west boundary of Asia in our study (Berzin et al., 1996). The Urals were also one of the major zones of continental convergence that contributed to the assembly of the late Paleozoic Pangea supercontinent (Hamilton, 1970; Şengör et al., 1993). To the north, the segments of the Gakkel Ridge and the Chersky Ridge construct the northern boundaries of the study area (Fig. 2). The Gakkel Ridge is connected to the global mid-ocean ridge system (Jokat and Schmidt‐Aursch, 2007). The Gakkel Ridge is anomalously deep and is the slowest spreading ridge in the global ridge system; its spreading rates decrease from 13 mm yr-1 at the western end to 6 mm yr-1 in the Laptev Sea (Edmonds et al., 2003). With combined GPS measurements and seismological data, the Chersky Ridge is considered the modern boundary between the Eurasian and North American Plates in northeastern Russia (Timofeev et al., 2012). The eastern and southeastern boundaries of our Asian study area mainly consist of Sakhalin Island, the Japanese islands and the Ryukyu arc (Fig. 2). The boundary between the Okhotsk Plate and the Eurasian Plate stretches along Sakhalin Island and Hokkaido, and its convergence rate increases from north to south (Jin et al., 2007). Following Sakhalin Island, the Japanese islands act as the easternmost boundary of Asia

and the boundary between the Pacific and Eurasian plates. It has generally been suggested that the subduction of the Pacific Plate beneath the Eurasian Plate strongly influences the tectonic evolution of easternmost China (Kusky et al., 2007). The boundary between the Philippine and Eurasian plates is approximately near the Nankai trough and the Ryukyu arc. Since there is no subduction zone beneath the Nankai trough, Seno and Maruyama (1984) concluded that the subduction of the Philippine Plate in this area had ceased. It is hard to estimate its influence on the stress field of the East Asia continent because the Nankai trough is not a very long boundary (Zang and Ning, 2002). On the other hand, the Philippine and Eurasian plates are weakly coupled along the Ryukyu arc (Zang et al., 1990), and the Okinawa trough behind the Ryukyu arc is still opening, where the southern part opened earlier and gradually extended toward the north (Eguchi, 1983; Letouzey and Kimura, 1986). These studies imply that the Philippine Plate has little influence on the Chinese continent along the Ryukyu trench (Zang et al., 2002). We draw the southern boundary similar to the definitions by Kreemer et al. (2000) and Bird (2013) that extend from Taiwan Island to Hainan Island and then to the southern end of the Indoburman Ranges, separating the study area from the Sunda Plate to the south (Fig. 2). It is important to note that this boundary is constrained by the geodetic velocity anomalies among different GPS stations (Bird, 2013). The following boundary, the Indoburman Ranges in western and northwestern Burma, is the eastern boundary between the Indian Plate and Eurasia Plate (Brunnschweiler, 1966). The southwestern boundary of Asia can be separated into two segments in the Tethyan tectonic domain, which are associated with the assembly processes of the Indian Plate and Arabian Plate (Fig. 2). The India-Eurasia collision is considered a conspicuous Cenozoic tectonic event, and the indentation of India has resulted in

deformation distributed over a vast area of the Asian continent, directly contributing to the Tibetan Plateau uplift with crustal thickening and the eastward and southeastward escape of crustal materials (e.g., Patriat and Achache, 1984; Le Pichon et al., 1992; Tapponnier et al., 1982). The load applied by the Indian Plate, especially along the arcuate orogen between the eastern and western Himalayan syntaxes, has received considerable attention since the highest convergence velocity and stress occur here (e.g., Banerjee and Bürgmann, 2002; Prasath et al., 2017; Bollinger et al., 2004). The Cenozoic collision between the Arabian and Eurasian plates generated the Zagros orogen (Smith, 1971). The Kopeh Dag along the southwestern-most boundary of the study area is the north border of the Zagros orogen in Iran (Fig. 2) (Reilinger et al., 2006)

3. Modeling The finite element method (FEM) is widely used to analyze kinematics, deformation, and seismological and geodynamic problems (Richardson, 1978; IsmailZadeh and Tackley, 2010; Ju et al., 2013; Dai et al., 2014; Zhu et al., 2006; Hou et al., 2010). We construct spherical shell models to explore the geodynamics of the Broad Triangle Area using the ANSYS 8.0 (University version) finite element software package. In this model, the element that we choose is a 3D 20-node homogeneous structural solid element (SOLID186 type in the ANSYS finite element software package), and the model is divided into 253,843 elements. Prior to calculation, material properties are assigned to the elements representing the various lithospheric characteristics. The finite element models describe elastic rock deformation. The mechanical behavior in the elastic domain is described using Hooke’s law. The elastic portion of this material model requires the specification of two elastic constants (Young's modulus and Poisson's ratio) that relate the stress (σ) and strain (ε).

For the two-dimensional case discussed here, the relationship is a simplified form of Hooke’s law:

{}

]{ }

𝜎𝑥 𝜀𝑥 1―𝑣 𝑣 0 𝐸 𝜎𝑦 = 𝜀 𝑣 1―𝑣 0 𝑦 (1 + 𝑣)(1 ― 2𝑣) 0 0 1 ― 2𝑣 𝜀𝑧 𝜎𝑧

[

where E is Young's modulus and ν is Poisson's ratio. More detailed descriptions of the FEM technique can be found elsewhere (Richardson et al., 1979; Logan, 2011; Fluent, 2012). We have built a great number of models to explore the geodynamics of the Broad Triangle Area and to find a best-fit model that can explain the spatial distribution of historical earthquakes in East Asia. We first introduce the best-fit model and then discuss the other typical comparison models. 3.1. Construction of a spherical shell FEM for the Broad Triangle Area It is necessary to mention that the observable earthquakes in the Broad Triangle Area mostly took place from 1900 to the present; therefore, we build a spherical shell FEM model based on the present-day plate boundary conditions. Seton et al. (2012) constructed a new model of global plate motion consisting of a set of continuously closing topological plate polygons with associated plate boundaries since the breakup of the supercontinent Pangea, which made our study possible. In general, thin elastic plates are used to simulate small tectonic plates; however, our study area is too large to ignore its spherical curvature. Meanwhile, Shi and Zhu (2006) suggested that it is necessary to use a spherical coordinate system to avoid systematic errors when using GPS displacements to calculate strain. For these reasons, it is reasonable to generate an elastic spherical shell model to simulate the lithosphere in this study area.

Based on the abovementioned studies concerning the geological background of the study area, the model can be simplified to a spherical shell that is subdivided into four types of geological units with different mechanical properties. The Siberian craton, Tarim craton, Qaidam block, Junggar block, western North China craton and western Yangtze craton are assumed to be the most rigid parts with Young's modulus of 80 GPa and Poisson's ratio of 0.1. The eastern Yangtze craton, Turan-Karakum, Lhasa and Qiangtang-Qamdo blocks are less rigid with Young's modulus of 80 GPa and Poisson's ratio of 0.2. The eastern North China craton is assumed to be less rigid due to the craton destruction process. The orogenic belts cover the largest area and have Young's modulus of 70 GPa and Poisson's ratio of 0.3 (Xie et al., 2007; Wang et al., 2012; Gu et al., 2014). Considering that seismicity mainly occurs in the shallow crust, it is rational to set the thickness of the spherical shell in our model to the crustal thickness. The thickness of the crust of our study area is assumed to decrease from approximately 70 km in the west to 30 km in the east (Yang et al., 2013). Fig. 3 Here 3.2. Constraints on the model The GPS technique has provided an effective and unique tool to precisely measure large-scale deformation of continents (Dixon, 1991). During the last decade, steadily growing numbers of GPS studies have started to provide key constraints on the mode of deformation in East Asia (e.g., Calais et al., 2003; Calais et al., 2006; Ren and Ma, 2003). Moreover, the observable historical earthquakes in the Broad Triangle Area occurred in the last century; therefore, we define the present-day boundary conditions of our model based on the GPS measurements (Fig. 3a). The western boundary is the Uralide orogen between the Siberia craton and Eastern European platform (Fig. 2 and Fig. 3a). The final tectonic evolution of the

Uralide orogen took place in the late Carboniferous and Permian (Zonenshain et al., 1984), and the Urals have remained relatively stable since then; thus, it is reasonable to assume that the Uralide orogen is fixed in our models now (Fig. 3a). Another fixed boundary in our models is the Chersky Ridge in northeast Russia since different studies have found that the rotation poles for the Eurasian and North American plates lie along the Chersky Ridge (Altamimi et al., 2002; Steblov et al., 2003) (Fig. 2 and Fig. 3a). According to GPS data in Eurasia, the displacement of the Siberian craton could be negligible relative to the adjacent regions; therefore, it is reasonable to assume that it is fixed in our models (Fig. 2 and Fig. 3a). Inside the continent, we also consider two large thrust belts in the North-South structural belt: the Longmenshan thrust belt (LMB) on the southwest margin of the North China craton and the Liupanshan thrust belt (LPB) on the west margin of the North China craton (yellow lines in Fig. 2). Because the North China Plate and the Yangtze Plate block the crustal material of the Tibetan Plateau from escaping eastward and southward, resulting in significant attenuation of the GPS displacement on both sides of the thrust belt, we apply corresponding displacement constraints on the two thrust belts according to the measured GPS data. (Fig. 3b). In addition, the distance from the geocenter to the bottom of the spherical shell is constrained to be constant to avoid the effect in which the spherical shell model is subjected to unreasonable radial deflection under lateral compression (Fig. 3b). The displacements of other boundaries are constrained based on different studies that all use Eurasia as the reference frame (Fig. 4). The displacement along the Gakkel Ridge (L1 in Fig. 3a) decreases from 10 mm/yr in the north to 0 mm/yr in the south and is perpendicular to the ridge (Edmonds et al., 2003). As mentioned above, the eastern boundaries (L2, L3, and L4 in Fig. 3a) include the convergent boundaries where prominent subduction takes place in the western Pacific Ocean. These boundaries are

assumed to be free boundaries that allow East Asian blocks to move eastward (Molnar and Tapponnier, 1975, 1977; Jolivet et al., 1990). The displacements on L5 and L6 (Fig. 3a) are 5 mm/yr trending approximately S45°E and 5 mm/yr trending approximately S65°E, respectively (Simons et al., 2007; Ren and Ma, 2003). The displacements on L7 along the Indoburman Ranges (Fig. 3a) can be divided into three parts: 20 mm/yr trending approximately N45°E, 15 mm/yr trending approximately N45°E and 15 mm/yr trending approximately S50°E from south to north (Sahu et al., 2006). The displacements on L8 vary continuously from 25 mm/yr trending eastward at the eastern Himalayan syntaxis to 30 mm/yr trending approximately N10°E at the western Himalayan syntaxis and peak at 40 mm/yr trending approximately N25°E in the middle of the Himalaya Ranges (Banerjee, 2002; Calais et al., 2003; Calais et al., 2006) (Fig. 3a). The displacements on L9 are the same as those on L10, which are 5 mm/yr trending northward (Karakhanyan et al., 2013; Vernant P and J. Chéry, 2006; Mousavi, 2017) (Fig. 3a). To explore the geodynamics of the Broad Triangle Area, the trajectories of the maximum horizontal compressional stress, the displacement vector, the distribution and the magnitude of the strain energy density in the study area are calculated.

4. Results and analysis 4.1. Trajectories of the maximum horizontal compressional stress The World Stress Map (WSM) released in 2016 displayed the mean contemporary orientations of the maximum horizontal compressional stress (blue arrows in Fig. 4) and provided a reasonable reference for the fitness of our results. Zhang et al. (2004) also reported the mean orientations of the maximum compressional stress (red solid arrows in Fig. 4) and the maximum compressional strain (red arrows in Fig. 4) based on focal mechanisms and GPS data in China.

Fig. 4 Here The Central-East Asian present-day stress field of our best-fit model refers to a visual general similarity between the calculated orientations of the maximum horizontal compressional stress (𝜎𝐻𝑚𝑎𝑥) and the observed WSM and GPS data as a whole (Fig. 4). Our calculated orientations of 𝜎𝐻𝑚𝑎𝑥 (blue short bars in Fig. 4) in our best-fit model fit well with the observed data (red arrows in Fig. 4) in Qiangtang-Qamdo, Tarim, North China, and South China (the deviations are <10°). To the south and north of Astana in Central Asia, the orientations of 𝜎𝐻𝑚𝑎𝑥 range from N10°W to N45°W, generally deviating slightly from the in situ observations of the WSM, which approximately range from N25°W to N54°W (Fig. 4). Nonetheless, our model correctly predicts that from south to north in this region, the orientations of 𝜎𝐻𝑚𝑎𝑥 gradually change from northward to northwestward. The calculated NNEtrending trajectories of 𝜎𝐻𝑚𝑎𝑥 also deviate from the observations in the Siberian craton (Fig. 4). In the southernmost part of the model, the average orientation of calculated 𝜎𝐻𝑚𝑎𝑥 is in general accord with the observed mean orientations of 𝜎𝐻𝑚𝑎𝑥 by Zhang et al. (2004). However, the calculated data are discrete (Fig. 4). Wang et al. (2006) developed a 2D model and found that taking the GPS data surrounding the Tibetan Plateau as boundary conditions is not sufficient to fit the observations within this area due to the complex geological environment. 4.2. The displacement vector map The displacement vector map (Fig. 5) shows the possible motions of the crustal material in the study area. This map indicates that in the Cenozoic, the crustal material in the Turan-Karakum blocks moves northward; the crustal material in the Tibetan Plateau moves northward and northeastward; the crustal materials in the eastern and

southeastern regions of Asia move eastward and southeastward, respectively. In addition, the results show that the displacement of the crustal material decays rapidly along the boundary of the Tibetan Plateau. Our results (red arrows in Fig. 5) fit well with the GPS observations (black arrows in Fig. 5) of Calais et al. (2006). Fig. 5 Here 4.3. The distribution and magnitude of the strain energy density The strain energy density is defined as the strain energy stored per unit volume and can be determined as follows: 𝑈 = [(𝜎21 + 𝜎22 + 𝜎23) ― 2ν(𝜎1𝜎2 + 𝜎1𝜎3 + 𝜎2𝜎3)]/2𝐸

where U represents the strain energy density; ν represents Poisson's ratio; E represents Young's modulus and 𝜎1, 𝜎2 and 𝜎3 represent the stress state. It has been generally considered that the strain energy density is one of the most important factors related to seismicity (e.g., Du and Aydin, 1996; Xu et al., 2002; Zhang et al., 2008). By calculating the strain energy density, we can explore the distribution of seismic activity and assess the possibility of future earthquakes. In the map showing magnitude contours of the strain energy density for the study area (Fig. 6), four regions with high strain energy density can be seen to occur in the western syntaxis, eastern syntaxis, North-South tectonic belt and southwest of Lake Baikal at the southern vertex of the Siberian craton. These regions with high strain energy density form the Broad Triangle Area (Fig. 6). Along the northwestern and southwestern boundaries of the BTA, the locations of seismic activity are also distributed within the high energy density regions. Meanwhile, the localization of high strain energy density forms a distinct triangular area, which is consistent with the location of the BTA. The variance of the strain energy density inside this area is also generally in accordance with the distribution of seismic activity among the blocks in

the Tibetan Plateau. For instance, the higher energy density region in the Tianshan orogen experiences more earthquakes than the lower energy density areas in the Junggar block, which have fewer earthquakes. Fig. 6 Here 4.4. Discussion In summary, the results given by the best-fit model are highly consistent with observed seismic activity and other geological evidence. However, there is a discrepancy in the results: the poor fit between the calculated results and observed data in the southern and northern parts of the study area. The factors affecting the fitness of the model can be identified as the boundary conditions and the differences in rock mechanical properties (Bott, 1992; Hou et a., 2006). To test the reasonability of our best-fit model, we also developed other models with different initial conditions to improve simulated results and investigate the reasons for the discrepancies. Fig. 7 Here In Model I, we ignore the heterogeneity of the crustal material in the study area: the cratons, blocks and orogenic belt units share the same mechanical parameters, including Young's modulus and Poisson's ratio. Meanwhile, we hold the boundary constraints invariant (Fig. 7). In other words, the continent is simulated as a homogeneous spherical shell. It is shown that the calculated maximum horizontal compressional stress trajectories and displacement vector map for model I are nearly the same as the results in our best-fit model (Fig. 7), which indicates that the anisotropy of the crust makes little difference to the orientation of the principal stress and the motion of crustal material. Although the geometry of the high strain energy density area seems to be a triangle, the internal details of the results present obvious deviations from the distribution of seismicity (Fig. 7). Therefore, Model I is not an appropriate

explanation for the dynamics of the BTA, and the anisotropy of the crust is an important factor for the formation of the BTA. In model II, the boundaries of the Siberian craton are not fixed to examine the influence of the Siberian craton on the Broad Triangle Area. Beyond this change, other model settings are the same as those in the best-fit model (Fig. 7). The magnitude of the calculated displacement vector is obviously higher than that from the GPS measurements, especially in the Siberian craton and the eastern part of the study area (Fig. 7). The orientations of maximum horizontal compressional stress in Model II fit better than those in the best-fit model compared with the observed orientations in the northern part of the study area (Fig. 7). However, the triangular area with high strain energy density cannot be formed in Model II. These comparative results between Model II and the best-fit model indicate that the stable Siberian craton (moving very slowly in the Eurasian Plate in relation to other parts but not absolutely fixed) plays a vital role in the formation of the BTA. It is generally suggested that the accelerating subduction of the Pacific Plate and Philippine Plate beneath the Eurasian Plate from 10 Ma cannot be ignored when exploring the tectonic setting in East Asia (Zhang et al., 2014). Therefore, in model III, we discuss the influence of the subduction of the Pacific Plate and the Philippine Plate beneath the Eurasian Plate on the Broad Triangle Area. For this purpose, we apply reasonable displacement constraints on boundaries L2, L3 and L4 (Fig. 7) based on the GPS measurements (Iaffaldano, 2012; Zang and Ning, 2002; Lin et al., 2016; Sagiya et al., 2000) without changing other model settings of the best-fit model. There are few differences in the calculated maximum horizontal compressional stress trajectories, displacement vector and strain energy density between Model III and the best-fit model. However, the best-fit model does better than Model III in simulating the correct

orientation of the principal stress in the southeastern part of the study area. These results lead to the conclusion that the direct effects of the subduction of the Pacific and Philippine plates beneath the Eurasian Plate cannot reach the North-South tectonic belt (the eastern boundary of the BTA), which is different from the results of other studies (Gao et al., 2010). In other words, subduction has little influence on the formation of the BTA. In model IV, we remove the displacement constraints on the western North China craton and the western Yangtze craton to explore the reasons for the formation of the eastern boundary of the BTA. Although the change in the maximum horizontal principal compressive stress trajectory is not large, the GPS displacement in the east is significantly increased, and the strain energy density is lower in the area of the NorthSouth tectonic belt (Fig. 7). These results contradict the actual observations. Therefore, we can determine that the blocking effects of the North China craton and the Yangtze craton play decisive roles in the formation of the eastern boundary of the BTA. To explore the impact of the collision between the Indian Plate and the Eurasian Plate, we construct model V by removing the displacement constraint on L8 and setting it as a free boundary (Fig. 7). The three maps (Fig. 7) show that the degree of deviation of the calculated results in model IV from the observed geological and seismic evidence around the BTA is obviously greater than those of the other three models (Fig. 7). Therefore, model V is not a viable option to explain the formation of the BTA. It is also suggested that the collision between the Indian Plate and the Eurasian Plate plays the most important role in the formation of the BTA and has a significant influence on the displacement field, stress field and localization of the strain energy density in our study area.

5. Conclusions

The Broad Triangle Area (BTA) is one of the most seismically active zones and a type locale for intraplate earthquakes. The geodynamics of the BTA are discussed here by constructing a series of elastic finite element method spherical shell models. By comparing the calculated results with the observed tectonic evidence and seismicity, the following conclusions can be reached: 1) The collision between the Indian Plate and Eurasian Plate plays the most important role in the geodynamics of the BTA and forms the Himalaya Range, which constructs the bottom edge of the BTA with the eastern and western Himalayan syntaxes as two vertexes of the BTA. Meanwhile, the southern end of the Siberia craton forms the northern vertex of the BTA because the stable Siberia craton is “fixed” relative to surrounding areas. 2) The heterogeneity of crustal material in the study area mainly controls the internal seismic distribution in the BTA but has little influence on the general geometry of the BTA. 3) After comparing Model IV and Model V with the best-fit model, we determine that the effects of the subduction of the Pacific and Philippine plates cannot reach the BTA and that the main reason for the formation of the eastern boundary of the BTA is the blocking effects of the North China craton and Yangtze craton.

Acknowledgments This work was supported by funds from the National Natural Science Foundation of China (Grant No. 41530207) and State Key Projects (Grant No. 2016ZX05051004).

References

Altamimi, Z., P. Sillard, and C. Boucher, 2002. ITRF2000: A new release of the International Terrestrial Reference Frame for earth science applications, J. Geophys. Res., 107(B10), doi:10.1029/2001JB000561, 2002. Banerjee P, Bürgmann R, 2002. Convergence across the northwest Himalaya from GPS measurements[J]. Geophysical Research Letters, 29(13): 30-1-30-4. Berzin R, Oncken O, Knapp J H, et al, 1996. Orogenic evolution of the Ural Mountains: results from an integrated seismic experiment[J]. Science, 274(5285): 220-221. Bird P, 2013. An updated digital model of plate boundaries[J]. Geochemistry, Geophysics, Geosystems, 4(3). Bollinger L, Avouac J P, Cattin R, et al., 2004. Stress buildup in the Himalaya[J]. Journal of Geophysical Research: Solid Earth, 109(B11). Bott, M.H.P., 1992. The stress regime associated with continental break-up. Geol. Soc. Lond. 68 (1), 125–136 Special Publications. Brunnschweiler R O, 1966. On the geology of the Indoburman ranges: (Arakan Coast and Yoma, Chin Hills, Naga Hills)[J]. Journal of the Geological Society of Australia, 13(1): 137-194. Calais E, Dong L, Wang M, et al, 2006. Continental deformation in Asia from a combined GPS solution[J]. Geophysical Research Letters, 33(24):194-199. Calais E, Vergnolle M, San'Kov V, et al., 2003. GPS measurements of crustal deformation in the Baikal‐Mongolia area (1994–2002): Implications for current kinematics of Asia[J]. Journal of Geophysical Research: Solid Earth, 108(B10). Cawood P A, Wang, Y., Xu, Y., et al., 2013. Locating South China in Rodinia and Gondwana: a fragment of greater India lithosphere?[J]. Geology, 41(8): 903-906. Chen, L.W., Zhang, Z.M., Ye, J.Y., et al, 2011. Numerically modeling the influence of rheological properties on tectonic deformation of Tibet plateau. Journal of Geodesy and Geodynamics (in Chinese), 31(3): 8-14. Cobbold P R, Davy P H, 1988. Indentation tectonics in nature and experiment. II: Central Asia[J]. Bulletin of the Geological Institution of the University of Upsala, 14: 143-162.

Dai, L.M., Li, S.Z., Lou, D., et al., 2014. Numerical modeling of Late Miocene tectonic inversion in the Xihu Sag, East China Sea Shelf Basin, China[J]. Journal of Asian Earth Sciences, 86(2):25-37. Dixon T H, 1991. An introduction to the Global Positioning System and some geological applications[J]. Reviews of Geophysics, 29(2): 249-276. Du Y, Aydin A, 1996. Is the San Andreas big bend responsible for the Landers earthquake and the eastern California shear zone? [J]. Geology, 24(24):219-222. Edmonds H N, Michael P J, Baker E T, et al., 2003. Discovery of abundant hydrothermal venting on the ultraslow-spreading Gakkel ridge in the Arctic Ocean[J]. Nature, 421(6920):252-6. Eguchi T, 1983. Seismotectonic of the Okinawa Trough and Ryukyu Arc[J]. Mem. Geol. Soc. China, 5: 189-210. England P, Houseman G, 1986. Finite strain calculations of continental deformation: 2. Comparison with the India‐Asia collision zone[J]. Journal of Geophysical Research: Solid Earth, 91(B3): 3664-3676. Feng, R., Ma, Z., Fang, J., et al., 2007. A Developing Plate Boundary: Tianshan-Baikal Active Tectonic Belt[J]. Earth Science Frontiers, 14(4):1-17. Fluent, A., 2012. 14.5 Theory Guide. ANSYS Inc, Canonsburg, PA, USA. Fu, Z., 2012. Summary of Seismic Activity in the Broad Triangle Area in the eastern Asia[M]. Seismological press (in Chinese). Gao, S., Zhang, B.R., Gu, X.M., et al., 1995. Silurian-Devonian provenance changes of South Qinling basins: implications for accretion of the Yangtze (South China) to the North China cratons[J]. Tectonophysics, 250(1-3): 183-197. Gao, X.L., Ma, X.J., Li, X.L., 2010. A surrounding and deep dynamic context of the great TriangleShaped seismic region in the eastern Asia continent. Earth Science Frontiers (in Chinese), 17(4): 33-42. Gladkochub D, Pisarevsky S, Donskaya T, et al., 2006. Siberian Craton and its evolution in terms of Rodinia hypothesis[J]. Episodes-Newsmagazine of the International Union of Geological Sciences, 29(3): 169-174.

Guiraud R, Bellion Y, 1995. Late Carboniferous to recent, geodynamic evolution of the West Gondwanian, cratonic, Tethyan margins[M]//The Tethys Ocean. Springer, Boston, MA: 101124. Gutenberg, B. and C.F. Richter, 1954. Seismicity of the earth and associated phenomena. Princeton University Press, 16-20 Gu, S., Liu, Y., Chen, Z., 2014. Numerical study of dynamic fracture aperture during production of pressure-sensitive reservoirs. Int. J. Rock Mech. Min. Sci. 70, 229–239. Han, Y., Zhao, G., Cawood P A, et al., 2016. Tarim and North China cratons linked to northern Gondwana through switching accretionary tectonics and collisional orogenesis[J]. Geology, 44(2): 95-98. Hamilton W, 1970. The Uralides and the motion of the Russian and Siberian platforms[J]. Geological Society of America Bulletin, 81(9): 2553-2576. He, J.K., Lu, S.J., Wang, W.M., 2013. Three-dimensional mechanical modeling of the GPS velocity field around the northeastern Tibetan plateau and surrounding regions[J]. Tectonophysics, 584:257-266. Hou, G., Liu, Y., Li, J., 2006. Evidence for ∼1. 8 Ga extension of the eastern block of the North China craton from SHRIMP U–Pb dating of mafic dyke swarms in Shandong province. J. Asian Earth Sci. 27 (4), 392–401. Hou, G.T., Wang, Y.X., Hari K R, 2010. The Late Triassic and Late Jurassic stress fields and tectonic transmission of North China craton. Journal of Geodynamics, 50, 318-324. Houseman G, England P, 1993. Crustal thickening versus lateral expulsion in the Indian‐Asian continental collision[J]. Journal of Geophysical Research: Solid Earth, 98(B7): 12233-12249. Iaffaldano G, 2012. The strength of large-scale plate boundaries: Constraints from the dynamics of the Philippine Sea plate since∼ 5 Ma[J]. Earth and Planetary Science Letters, 357: 21-30. Ismail-Zadeh, A., Tackley, P., 2010. Computational Methods for Geodynamics. Cambridge University Press. Jahn B, 2004. The Central Asian Orogenic Belt and growth of the continental crust in the Phanerozoic[J]. Geological Society, London, Special Publications, 226(1): 73-100.

Jin, S.G., Park P H , Zhu, W.Y., 2007. Micro-plate tectonics and kinematics in Northeast Asia inferred from a dense set of GPS observations[J]. Earth and Planetary Science Letters, 257(34):0-496. Jolivet L, Davy P, Cobbold P, 1990. Right‐lateral shear along the northwest Pacific margin and the India‐Eurasia collision[J]. Tectonics, 9(6): 1409-1419. Jolivet M, Arzhannikov S, Chauvet A, et al, 2013. Accommodating large-scale intracontinental extension and compression in a single stress-field: A key example from the Baikal Rift System[J]. Gondwana Research, 24(3–4):918-935. Jokat W, Schmidt ‐ Aursch M C, 2007. Geophysical characteristics of the ultraslow spreading Gakkel Ridge, Arctic Ocean[J]. Geophysical Journal International, 168(3): 983-998. Ju, W., Hou, G.T., Hari, K.R., 2013. Mechanics of mafic dyke swarms in the Deccan large igneous province: palaeostress field modelling. J. Geodyn. 66, 79–91. Karakhanyan A , Vernant P , Doerflinger E , et al., 2013. GPS constraints on continental deformation in the Armenian region and Lesser Caucasus[J]. Tectonophysics, 592(Complete):39-45. Kreemer C, Haines J, Holt W E, et al., 2000. On the determination of a global strain rate model[J]. Earth, Planets and Space, 52(10): 765-770. Kusky T M , Windley B F , Zhai, M.G., 2007. Tectonic evolution of the North China Block: from orogen to craton to orogen[J]. Geological Society London Special Publications, 280(1). Lei, J.S., Zhao, D.P., Su, J.R., et al., 2009. Fine seismic structure under the longmenshan fault zone and the mechanism of large Wenchuan earthquake[J]. Chinese Journal of Geophysics (in Chinese), 52(2): 339-345. Letouzey J, Kimura M, 1986. The Okinawa Trough: genesis of a back-arc basin developing along a continental margin[J]. Tectonophysics, 125(1-3): 209-230. Le Pichon X, Fournier M, Jolivet L, 1992. Kinematics, topography, shortening, and extrusion in the India‐Eurasia collision[J]. Tectonics, 11(6): 1085-1098. Li, S.Z., Suo, Y.H., Dai, L.M., et al, 2010. Development of the Bohai Bay Basin and destruction of the North China Craton[J]. Earth Science Frontiers, 17(4): 64-89

Lin W , Byrne T B , Kinoshita M , et al., 2016. Distribution of stress state in the Nankai subduction zone, southwest Japan and a comparison with Japan Trench[J]. Tectonophysics, 692(4):120130. Logan, D.L., 2011. A First Course in the Finite Element Method. Cengage Learning. Ma, Z.J., Zheng, D.L., 1981. The central axis zone and its seismicity in the china mainland and Mongolia [J]. Earthquake Research, 4(4): 421-436 (in Chinese with English abstract). Molnar, P, and P. Tapponnier.,1975. Cenozoic Tectonics of Asia: Effects of a Continental Collision [J]. Science, 189(4201):419-426. Molnar, P, and P. Tapponnier., 1977, Relation of the tectonics of eastern China to the India-Eurasia collision: Application of slip-line field theory to large-scale continental tectonics: Geology, v. 5, p. 212–216. Mousavi S M, 2017. Spatial variation in the frequency-magnitude distribution of earthquakes under the tectonic framework in the Middle East[J]. Journal of Asian Earth Sciences, 147: 193-209. Patriat P, Achache J, 1984. India–Eurasia collision chronology has implications for crustal shortening and driving mechanism of plates[J]. Nature, 311(5987): 615-621. Peltzer G, Saucier F, 1996. Present‐day kinematics of Asia derived from geologic fault rates[J]. Journal of Geophysical Research: Solid Earth, 101(B12): 27943-27956. Pospelov I, Petrov O, Shokalsky S, et al., 2016. New tectonic map of Northern-Central-Eastern Asia: Position and evolution of Mesozoic sedimentary basins[J]. Global Geology, 19(4): 261276. Prasath R A, Paul A, Singh S, 2017. Upper crustal stress and seismotectonics of the Garhwal Himalaya using small-to-moderate earthquakes: Implications to the local structures and free fluids[J]. Journal of Asian Earth Sciences, 135: 198-211. Reilinger R , Mcclusky S , Vernant P , et al., 2006. GPS constraints on continental deformation in the Africa-Arabia-Eurasia continental collision zone and implications for the dynamics of plate interactions[J]. Journal of Geophysical Research Solid Earth, 111(B5):43-55. Ren, J.W, Ma, Z.J., 2003. Crustal movement and tectonic deformation of eastern Asia [J]. Earth Science Frontiers, 10(Suppl): 58-65(in Chinese)

Ren, J.S., Niu, B.G., Wang, J., et al., 2013. Advances in research of Asian geology—A summary of 1: 5M International Geological Map of Asia project[J]. Journal of Asian Earth Sciences, 72: 311. Ren, J.S., Wang, Z.X., Chen, B.W., 2003. Tectonics of China from a global viewing -- brief description of the tectonics map of China and adjacent regions. Acta Geosicientia Sinica, 24(1): 1-2(in Chinese) Ren, J., Wang, Z., Chen, B., Jiang, C., Niu, B., Li, J., Xie, G., He, Z., Liu, Z., 1999. The Tectonics of China from a Global View—A Guide to the Tectonic Map of China and Adjacent Regions. Geological Publishing House, Beijing, p. 32. Richardson, R.M., 1978. Finite element modeling of stress in the Nazca plate: driving forces and plate boundary earthquakes. Tectonophysics 50 (2-3), 223–248. Richardson, R.M., Solomon, S.C., Sleep, N.H., 1979. Tectonic stress in the plates. Rev.Geophys. 17 (5), 981–1019. Sagiya T, Miyazaki S, Tada T, 2000. Continuous GPS array and present-day crustal deformation of Japan[J]. Pure and applied Geophysics, 157(11-12): 2303-2322. Sahu V K, Gahalaut V K, Rajput S, et al., 2006. Crustal deformation in the Indo-Burmese arc region: implications from the Myanmar and Southeast Asia GPS measurements[J]. Current Science: 1688-1693. Şengör A M C, Altıner D, Cin A, et al., 1988. Origin and assembly of the Tethyside orogenic collage at the expense of Gondwana Land[J]. Geological Society, London, Special Publications, 37(1): 119-181. Şengör A M C, Natal'In B A, Burtman V S, 1993. Evolution of the Altaid tectonic collage and Palaeozoic crustal growth in Eurasia[J]. Nature, 364(6435): 299. Seno T, Maruyama S, 1984. Paleogeographic reconstruction and origin of the Philippine Sea[J]. Tectonophysics, 102(1-4): 53-84. Seton, M., Müller, R.D., Zahirovic, S., Gaina, C., Torsvik, T., Shephard, G., Talsma, A., Gurnis, M., Turner, M., Maus, S., Chandler, M., 2012. Global continental and ocean basin reconstructions since 200 Ma. Earth Sci. Rev. 113 (3), 212–270.

Simons W J F, Socquet A, Vigny C, et al., 2007. A decade of GPS in Southeast Asia: Resolving Sundaland motion and boundaries[J]. Journal of Geophysical Research: Solid Earth, 112(B6). Shen, Z.K., Jackson D D , Kagan Y Y, 2007. Implications of Geodetic Strain Rate for Future Earthquakes, with a Five-Year Forecast of M5 Earthquakes in Southern California[J]. Seismological, Research Letters, 78(1):116-120. Shi, Y.L., Zhu, S.B., 2006. Discussion on method of calculating strain with GPS displacement data[J]. Journal of Geodesy and Geodynamics, 26(1): 1-8.(in Chinese) Smith A G, 1971. Alpine deformation and the oceanic areas of the Tethys, Mediterranean, and Atlantic[J]. Geological Society of America Bulletin, 82(8): 2039-2070. Steblov G M, Kogan M G, King R W, et al., 2003. Imprint of the North American plate in Siberia revealed by GPS[J]. Geophysical Research Letters, 30(18). Tapponnier P, Molnar P, 1979. Active faulting and Cenozoic tectonics of the Tien Shan, Mongolia and Baykal regions. Journal of Geophysical Research. 84(B7): 3425-3459． Tapponnier P, Peltzer G, Le Dain A Y, et al., 1982. Propagating extrusion tectonics in Asia: New insights from simple experiments with plasticine[J]. Geology, 10(12): 611-616. Timofeev V Y, Ardyukov D G, Solov'ev V M, et al., 2012. Plate boundaries in the Far East region of Russia (from GPS measurement, seismic-prospecting, and seismological data)[J]. Russian Geology and Geophysics, 53(4): 376-391. Veevers J J, 2004. Gondwanaland from 650–500 Ma assembly through 320 Ma merger in Pangea to 185–100 Ma breakup: supercontinental tectonics via stratigraphy and radiometric dating[J]. Earth-Science Reviews, 68(1-2): 1-132. Vergnolle M , Calais E , Dong L, 2007. Dynamics of continental deformation in Asia[J]. Journal of Geophysical Research, 112(B11):B11403. Vernant P , J. Chéry, 2006. Low fault friction in Iran implies localized deformation for the Arabia– Eurasia collision zone[J]. Earth and Planetary Science Letters, 246(3-4):0-206. Wan T.F. A new Asian tectonic unit map. Geol. China, 2013, 40(5): 1351-1365(in Chinese). Wang, H., Zhang, G.M., Shi, Y.L., et al., 2006. Numerical simulation of movement and deformation of Qinghai-Tibet Plateau [J]. Journal of Geodesy and Geodynamics, 26(2): 15-23(in Chinese).

Wang, X.Q., Schubnel, A., Fortin, J., David, E.C., Guéguen, Y., Ge, H.K., 2012. High Vp/Vs ratio: saturated cracks or anisotropy effects? Geophys. Res. Lett. 39 (11). Windley B F, Alexeiev D, Xiao W, et al., 2007. Tectonic models for accretion of the Central Asian Orogenic Belt[J]. Journal of the Geological Society, 164(1): 31-47. Xiao, W., Windley B F, Hao J, et al.,2003. Accretion leading to collision and the Permian Solonker suture, Inner Mongolia, China: termination of the central Asian orogenic belt[J]. Tectonics, 22(6). Xie, R., Qu, T., Qian, G., 2007. Structural geology. China University of Mining and Technology Press, pp. 25–37. Xu, C.J., Dong, L.X., Shi, C., Li, Y.X., Hu, X.K., 2002. A study on annual accumulation of strain energy density significance by using GPS measurement in North China. Chinese[J]. Chinese Journal of Geophysics (in Chinese), 45(4): 184-506. Xu, Z.H. A present-day tectonic stress map for eastern Asia region. Acta Seismologica Sinica (in Chinese), 2001, 23(5): 497-501. Yang, T., Fu, R.S., Huang, J.S., 2013. On the inversion of effective elastic thickness of the lithosphere with Moho relief and topography data. Chinese Journal of Geophysics (in Chinese), 55(11): 3671-3680. Zang, S.X., Ning, J.Y., 2002. Interaction between the Philippine Sea Plate and the Eurasia Plate and its Influence on the Tectonics of Eastern Asia[J]. Chinese Journal of Geophysics, 45(2): 184194. Zang, S.X., Ning, J.Y., Xu, L., 1990. Distribution of earthquakes, figuration of the Benioff zone and stress state in the Ryukyu island arc[J]. Acta Seismologica Sinica, 3(2): 137-148. Zhang, D.N., Gao, L.S., 1989. Three dimensional numerical simulation of the stress field in East Asia. Earth Research in China (in Chinese), 5(4): 24-33. Zhang, G.M., Ma, H.S., Wang, H., Li, L., 2004. The relation between active blocks and strong earthquakes in China. Science in China (Seri. D) (in Chinese), 34(7): 591∼599 Zhang, H.F., Gao, S., Zhang, B.R., et al., 1997. Pb isotopes of granitoids suggest Devonian accretion of Yangtze (South China) craton to North China craton[J]. Geology, 25(11): 1015-1018.

Zhang, L.F., Lü, Y.J., Peng, Y.J., Xie, Z.J., 2008. Analysis of seismological energy density of minormoderate earthquakes[J]. earthquake research in China (in Chinese). 24(4): 407-414 Zhang, P.Z., Deng, Q.D., Zhang, G.M., et al., 2003. Active tectonic blocks and strong earthquakes in the continent China. Science in China (Series D) (in Chinese), 46(Suppl): 13-24. Zhang, P.Z., Zhang, H.P., Zheng, W.J., et al. 2014. Cenozoic tectonic evolution of continental eastern Asia［J］. Seismology and Geology(in Chinese), 36(3):574-585. Zhao Z, Oike K, Matsumura K, et al., 1990. Stress field in the continental part of China derived from temporal variations of seismic activity[J]. Tectonophysics, 178(2-4): 357-372. Zhu, S.B., Cai, Y.E., Shi, Y.L., 2006. The Contemporary Tectonic Strain Rate Field of Continental China Predicted from GPS Measurements and its Geodynamic Implications[J]. Pure and Applied Geophysics, 163(8):1477-1493. Zhu, S.B., Zhang, P.Z., Shi, Y.L., 2010 A study on the mechanisms of strong earthquake occurrence in the North China Basin. Chinese J. Geophys. (in Chinese), 53(6): 1409-1417. Zonenshain L P, Korinevsky V G, Kazmin V G, et al., 1984. Plate tectonic model of the South Urals development[J]. Tectonophysics, 109(1-2): 95-135

Figure captions Fig. 1 Topography, main active faults (black lines), seismicity (United States Geological Survey, 1900 to present, M≥6) in eastern Asia and the location of the Broad Triangle Area (BTA) (within the red dotted triangle) (modified from Vergnolle et al. (2007)). The seismicity within the BTA is represented by red circles, and the magnitude of an earthquake is indicated by the size of the circle. Fig. 2 Simplified tectonic map of Asia (modified from Ren et al. (2013), Wan (2013) and Pospelov et al. (2016)). Our study area is delimited by red and blue lines. LMB: Longmenshan belt; LPM: Liupanshan belt. The distribution of the plate boundaries is based on the digital model published by Bird (2013) and represented by dotted blue lines.

Fig. 3 (a) The boundary conditions and geological unit divisions of the spherical shell model for the study area. (b) diagram of the cross section AA’. (a) The white arrows represent the approximate orientations and magnitudes of the displacement constraints. The yellow triangles represent fixed boundaries. The yellow lines represent the boundaries of cratons, blocks and the study area. The red and white lines represent the boundaries of crustal thickness variations. 1-Siberian craton, 2-North China craton, 3-Yangtze craton, 4-Turan-Karakum, 5-Lhasa, 6-Qiangtang-Qamdo, 7Qaidam, 8-Tarim craton, 9-Junggar, 10-Orogenic belt. (b) Diagram of the cross section AA’. The crustal thickness decreases from 70 km to 30 km, and the crustal material is divided into two types: orogenic belts (gray color) and cratons and blocks (yellow color). The circles represent the fixed distance from the geocenter to the bottom of the crust while sliding motion is allowed. Fig. 4 The maximum horizontal compressional stress (𝜎𝐻𝑚𝑎𝑥) trajectory map of our best-fit model and the reference results obtained from the WSM released in 2018 and Zhang et al. (2004). The short blue bars represent our results. The blue arrows represent the results obtained from the WSM released in 2018. The red arrows represent the results obtained by Zhang et al. (2004). 1-Siberian craton, 2-North China craton, 3-Yangtze craton, 4-TuranKarakum, 5-Lhasa, 6-Qiangtang-Qamdo, 7-Qaidam, 8-Tarim craton, 9-Junggar, 10Orogenic belt. Fig. 5 The calculated displacement vector map of our best-fit model. The gray polygons represent the cratons and blocks in our study area. The red arrows represent the motions of crustal materials based on our best-fit model. The black arrows represent the GPS observations (Calais et al., 2006). Fig. 6 The magnitude contours of the strain energy density in our best-fit model

The red circles represent historical earthquakes from 1900 to the present, and their sizes indicate the earthquake magnitudes. The black lines represent the boundaries of the cratons and blocks. The white line marks the boundaries of the BTA. 1-Siberian craton, 2-North China craton, 3-Yangtze craton, 4-Turan-Karakum, 5-Lhasa, 6-QiangtangQamdo, 7-Qaidam, 8-Tarim craton, 9-Junggar, 10-Orogenic belt. Fig. 7 Comparison of different model results for the calculated displacement vectors, strain energy density and stress trajectories. Model I ignores the anisotropy of the crust. Model II removes the constraints on the Siberian craton. Model III applies appropriate displacement constraints on the eastern boundaries to simulate the effects of subduction. Model IV ignores the blocking effects of the North China craton and Yangtze craton. Model V removes the displacement on L8. The legends are the same as those in previous figures.

Fig. 1

Fig. 2

Fig. 3

Fig. 4

Fig. 5

Fig. 6

Fig. 7

Highlights 1.It is the first attempt at establishing a numerical model for the mechanism of the BTA. 2.We establish a spherical shell model to avoid system error when using GPS displacements. 3.The main factors include Indian-Eurasian collision, stable Siberian craton and blocking effect of Yangtze and North China craton.

Declaration of interest statement

We declare that we have no financial and personal relationships with other people or organizations that can inappropriately influence our work, there is no professional or other personal interest of any nature or kind in any product, service and/or company that could be construed as influencing the position presented in, or the review of, the manuscript entitled “Geodynamics of the Broad Triangle Area (active seismic zone) in Asia: stress field modeling”.