1 Introduction
Several statistical analyses performed on subduction zones have indicated that a direct correlation between the slab dip and the age of the oceanic lithosphere does not exist (e.g., Jarrard, 1986) and that other parameters may have a more significant impact on slab dip variations, including: (1) the absolute motion of the overriding plate (Heuret et al., 2007; Lallemand et al., 2005); (2) the subduction rate (Cruciani et al., 2005; Lallemand et al., 2005); (3) the nature of the overriding plate (oceanic or continental) (Lallemand et al., 2005); (4) the slab width (Guillaume et al., 2009; Schellart et al., 2007). Nevertheless, although the mutual interaction between these parameters controls the slab dip variations, the low cross-correlation coefficients obtained suggest that other parameters may have relevant roles (Cruciani et al., 2005; Lallemand et al., 2005). In light of this, some analogue and numerical models were performed to investigate the influence of other causes of subduction geometry variation, such as: (5) the physical parameters of the slab (Royden and Husson, 2006); (6) the back-arc stress variation (Clark et al., 2008; Heuret et al., 2007); (7) the slab interaction with the upper-lower mantle boundary in a long-term subduction (Christensen, 1996; Guillaume et al., 2009,Olbertz et al., 1997); (8) the trenchward motion (Espurt et al., 2008). At present, no analysis has been performed to investigate the role of the thermal state of the overriding plate on the slab dip. To contribute to this subject, we used a 2-D finite element numerical model to develop a preliminary quantitative analysis of the impact of the thermal state of an ocean-continent convergent system. This kind of analysis could be very difficult if based solely on natural data. A comparison between the results of numerical modeling and the data obtained from natural systems allows the evaluation of the robustness and fidelity of the model. A convergence rate of 3 cm/a was chosen for the simulations. The primary motivation of this choice is that the present study is part of a wider analysis devoted to unravel the precollisional evolution of the Alpine belt. Within this geodynamic context 3 cm/a reveals the most appropriate subduction rate, as suggested by both geological evidences (Agard et al., 2009) and parametric analysis (Roda et al., 2010). On the other hand, while the effects induced by a variation of the subduction rate on the slab dip have been already widely investigated (e.g., Cruciani et al., 2005; Lallemand et al., 2005; Roda et al., 2010), the goal of the present work is to enlighten the effects of a factor not yet taken into account, namely, the thermal state of the overriding plate.
2 Model setup
We used the same approach as in Roda et al. (2010). The equation of continuity Eq. (1):
| (1) |
| (2) |
| (3) |
Model setup and boundary conditions for all of the simulations. The legend illustrates the marker provenance. The different fixed depths of the 1600 K are also displayed. Distances are not in scale.
Configuration du modèle et conditions aux limites pour toutes les simulations. La légende illustre la nature des couches. Les différentes profondeurs fixées à 1600 K sont aussi présentées. Les longueurs ne sont pas à l’échelle.
The model represents a strongly coupled ocean-continent convergent system (Clift and Vannucchi, 2004), affected by slab dehydration and progressive mantle wedge serpentinisation; erosion and sedimentation processes are also accounted for, while shear heating is not taken into account (Marotta and Spalla, 2007; Meda et al., 2010; Roda et al., 2010). The initial geometry, the material parameters and the boundary conditions are the same as in Roda et al. (2010) except for the convergent velocity of 3 cm/a, which here is fixed up to a 50 km depth in the trench zone, along a 45o dip plane. The lithosphere base coincides with the 1600 K isotherm. We chose four different thicknesses of the oceanic lithosphere (60, 80, 95 and 110 km), and, to modify the thermal state of the overriding continental plate, we fixed for it five different thicknesses at a distance of 800 km from the trench, with a continent/ocean ratio ranging from 100% to 200% of the thickness of the oceanic lithosphere (Fig. 1 and Table 1).
Liste des simulations numériques exécutées (O.L., lithosphère océanique ; C.L., lithosphère continentale.
| Model | O.L. | C.L. | C.L./O.L. |
| ID | Thickness (km) | Thickness (km) | Ratio (%) |
| 01 | 60 | 60 | 100 |
| 02 | 60 | 70 | 120 |
| 03 | 60 | 85 | 140 |
| 04 | 60 | 95 | 160 |
| 05 | 60 | 120 | 200 |
| 06 | 80 | 80 | 100 |
| 07 | 80 | 95 | 120 |
| 08 | 80 | 110 | 140 |
| 09 | 80 | 130 | 160 |
| 10 | 80 | 160 | 200 |
| 11 | 95 | 95 | 100 |
| 12 | 95 | 115 | 120 |
| 13 | 95 | 130 | 140 |
| 14 | 95 | 150 | 160 |
| 15 | 95 | 190 | 200 |
| 16 | 110 | 110 | 100 |
| 17 | 110 | 130 | 120 |
| 18 | 110 | 150 | 140 |
| 19 | 110 | 175 | 160 |
| 20 | 110 | 220 | 200 |
3 Model results
In Fig. 2, the results of the numerical simulations obtained for an 80-km-thick oceanic lithosphere are shown, after 25 Ma from the subduction start. The ratio of the continental/oceanic plate thickness ranges from 100% (panel a1) to 200% (panel a5) and a decrease of the mean slab dip with the increase of the upper lithosphere thickness is clearly evident (Fig. 2, from a1 to a5). Together with the slab dip decrease, a friction increase between the two plates occurs, increasing the ablation of the continental material belonging to the overriding plate, according to Faccenda et al. (2008) and Roda et al. (2010). In addition, in a low dip configuration, the efficacy of the crustal marker exhumation, driven by channel flow process, increases (Roda et al., 2010), also producing the recycling of the upper oceanic crust and sediments (Fig. 2, panels a4 and a5). The large-scale stream-line analysis points out that the increase of the continental lithosphere thickness is accompanied by a reduction in the size of the right-side convective cell, until the development of a single-cell configuration (Fig. 2, b4 and b5). The large-scale left-side convective cell concurs with the buoyancy and the suction to the slab flattening.
Simulation results obtained for 80 km thick oceanic lithosphere and for five different thermal states of the overriding plate (from 100% (panel a1) to 200% (panel a5) of the thickness of the oceanic lithosphere), after 25 Ma since the beginning of the subduction. The white isotherm represents the base of the lithosphere (1600 K). In the right-side (panels b) the stream-lines for each simulations are shown.
Résultats de la simulation obtenus pour une lithosphère océanique de 80 km d’épaisseur et pour cinq différents régimes thermiques de la plaque chevauchante (de 100 % (panneau a1) à 200 % (panneau a5) de l’épaisseur de la lithosphère océanique), après 25 Ma depuis le commencement de la subduction. L’isotherme en blanc représente la base de la lithosphère (1600 K). Dans les panneaux b, à droite, sont présentés les lignes de courant pour chaque simulation.
The influence of the lithospheric thickness variation (either oceanic and continental) on the subduction geometry was evaluated by plotting the slab top profiles starting from a depth of 50 km, which is the lower limit of the fixed velocity vectors, and at 50 km from the trench, after 30 Ma from the subduction start (Fig. 3). Generally, the slab profiles change with changes in the lithospheric thickness of the overriding plates, and this variability increases with decreasing oceanic thickness.
Traces of the slab top obtained starting from 50 km depth for each performed simulation, after 30 Ma since the beginning of the subduction. Different colors discriminate the initial thicknesses of the overriding plate, as shown in the legend.
Traces du sommet de la plaque plongeante partant de 50 km de profondeur pour chaque simulation effectuée, après 30 Ma depuis le commencement de la subduction. Les différentes couleurs distinguent les épaisseurs initiales de la plaque chevauchante, comme le montre la légende.
For 60- and 80-km-thick oceanic lithospheres, the slab dip, evaluated to depth of 250-300 km, decreases with increasing thickness of the continental lithosphere (Fig. 3, panels a and b). Below a depth of 200 km, a gradual increase of the slab dip is predicted until reaching a near vertical setting followed, in some cases, by a dip inversion. For an oceanic lithosphere thicker than 80 km, lower dip variability is observed with the increase of the continental thickness (Fig. 3, panels c and d). This is due to the colder global thermal state of the system (because of the considerable plate thicknesses, either oceanic or continental), that drives the formation of a major single-convective cell below the subducting plate and minimizes the intensity of the mantle circulation above the slab (Fig. 4, panels a and b), that is a primary controlling factor of the slab dip variability.
Stream-lines obtained for two different oceanic thicknesses (80 km, a and 110 km, b) at the same continental/oceanic thickness ratio (100%).
Lignes de courant obtenues pour deux épaisseurs de la lithosphère océanique (80 km, a et 110 km, b), avec même rapport d’épaisseur océan/continent (100 %).
The analysis of the dip distribution from 50 to 400 km depth for oceanic lithospheres of 60 and 80 km and from 50 to 300 km depth for oceanic lithospheres of 95 and 110 km shows that a strong relationship exists between the slab dip and the thickness of the continental lithosphere (Fig. 5, panel a), with a rapid dip decrease with the continental plate thickness increase, in particular for the oceanic thicknesses of 60, 80 and 95 km. A low correlation between the slab dip and the continental thickness was observed, with a near constant dip between 27o and 32o for oceanic lithosphere thickness of 110 km. The results obtained for the 80-km-thick oceanic lithosphere show the maximum slab dip variability.
Statistical analysis of the results of the performed simulations: panel a, mean slab dip vs thickness of the continental lithosphere; panel b, mean slab dip vs continental/oceanic thickness ratio; panel c, mean slab dip vs thickness of the oceanic lithosphere.
Analyse statistique des résultats des simulations effectuées : panneau a, pendage moyen en fonction de l’épaisseur de la lithosphère continentale ; panneau b, pendage moyen en fonction du rapport d’épaisseur continent/océan ; panneau c, pendage moyen en fonction de l’épaisseur de la lithosphère océanique.
The general dip decreasing with the increase of the continental lithosphere thickness also occurs when the imposed shallow angle is changed. The Fig. 6 (panels a to d) shows the slab dip variation for an entrance dip of 37o: in this case, the variability is less evident because of the higher efficiency of the overriding plate suction; the latter becomes dominant when a 30o entrance dip is imposed, driving the underplating of the subducting plate and preventing the high dip variability with the thickness of the continental lithosphere (Fig. 6, panels e to h).
Simulation results obtained for 80 km thick oceanic lithosphere and for three different thermal states of the overriding plate (100% (a–e), 120% (b–f), and 160% (c–g) of the thickness of the oceanic) lithosphere. In panels a, b and c the results for an imposed shallow dip of 37° are shown (see legend of Fig. 2). In panel d the traces of the slab top obtained starting from 50 km depth for each performed simulation, after 30 Ma since the beginning of the subduction, are displayed. In panels e, f and g the results for an imposed shallow dip of 30° are shown (see legend of Fig. 2). In panel h the traces of the slab top obtained starting from 50 km depth for each performed simulation, after 30 Ma since the beginning of the subduction, are displayed.
Résultats de la simulation obtenus pour une lithosphère océanique de 80 km d’épaisseur et pour trois différents régimes thermiques de la plaque chevauchante (100 % (a–e) 120 % (b–f) et 160 % (c–g) d’épaisseur de la lithosphère océanique. Dans les panneaux a, b et c, sont donnés les résultats pour un pendage faible imposé de 37° (voir légende de la Fig. 2). Dans le panneau d, sont présentées les traces du sommet de la plaque, partant de 50 km de profondeur pour chaque simulation réalisée, à partir de 30 Ma depuis le commencement de la subduction. Dans les panneaux e, f, et g, sont présentés les résultats pour un pendage faible imposé de 30° (voir légende de la Fig. 2). Dans le panneau h, sont montrées les traces du sommet de la plaque, obtenues en partant de 50 km de profondeur pour chaque simulation réalisée, après 30 Ma depuis le commencement de la subduction.
Although a direct relationship also exists between the slab dip and the continent-ocean thickness ratio, a weaker correlation is observed in comparison to the absolute thickness of the continental lithosphere (Fig. 5, panels a and b). This suggests that the thermal state of the overriding plate has more influence on the control of the slab dip than the variation of the continent-ocean thickness ratio.
In agreement with Jarrard (1986), Cruciani et al. (2005) and Lallemand et al. (2005), the model predictions suggest that a direct correlation between the slab dip and the age of the oceanic lithosphere does not exist (Fig. 5, panel c).
4 Comparison with the natural data
To compare the model predictions and the natural data we considered the mean dips obtained for several worldwide ocean-continent subduction zones, extracted from the databases in Cruciani et al. (2005) and Lallemand et al. (2005). The thickness of the overriding plate at a distance of 800 km from the trench for each natural data was based on the global thermal model of Artemieva (2006). Due to the lacking of thermal predictions in the submerged areas, we excluded Japan and Aleutins subductions and some data on the Mexico, Philippines and Indonesia subductions, and we grouped the data into four classes of different oceanic lithosphere thickness: 60 km (data from 30 km to 70 km), 80 km (from 70 km to 85 km), 95 km (from 90 km to 100 km) and 110 km (from 100 km to 145 km). The thicknesses of the natural subduction plates derived as a function of the slab ages were obtained by Cruciani et al. (2005) and Lallemand et al. (2005) using the following relation (Turcotte and Schubert, 2002) (Eq. (4)):
| (4) |
Slab dip vs thickness of the continental lithosphere: comparison between natural data (empty circles) and model predictions (black symbols) grouped for four different thicknesses of the oceanic lithosphere (see the legend).
Pendage de la plaque plongeante en fonction de l’épaisseur de la lithosphère continentale : comparaison entre observations (cercles vides) et prévisions du modèle (symboles pleins), regroupées pour quatre épaisseurs différentes de la lithosphère océanique (voir légende).
Even if the convergent velocities of the natural subduction systems were not discriminated, a good agreement between the natural data and the model predictions occurs, with a similar trend of the dip related to the thickness variation of the continental lithosphere: high slab dips seem to characterize the subduction systems with thinner overriding plates, and the dips decrease with increasing lithospheric thickness (Fig. 7). Furthermore, in natural cases, the maximum variability of the slab dip is associated with oceanic lithosphere thicknesses of 60 and 80 km, while lower variability characterizes the thicker oceanic plate thickness (e.g., 110 km), in agreement with the results obtained from numerical simulations (Fig. 5, panel a, and Fig. 7).
5 Conclusions
The model predictions suggest that the slab dip is influenced by the thermal state of the overriding plate and, in particular, for the cold overriding plate configuration, the efficacy of the right-side convective cell decreases and a shallower slab dip setting occurs. In contrast, in the hot overriding plate configuration, a highly efficient convective flow is developed in the right-side of the model domain and a steeper subduction occurs. Furthermore, a thickness increase of the subducting plate decreases the variability of the slab geometry, with the development of a stable single-convective cell, disregarding the thermal state of the overriding plate.
The model results also indicate that a direct correlation between the slab dip and the age of the oceanic lithosphere does not exist, at least for subduction plates thinner that 110 km, as suggested by the statistical analysis of natural subductions. On the contrary, our comparative analysis strengths the correlation between the slab dip and the thermal state of the overriding plate.
Therefore, this study indicates that a factor not yet taken into account, namely, the thermal state of the overriding plate, exerts a relevant influence on the subduction geometry, although other simulations with different convergent rates are needed to enforce the cross-correlation between the different factors and to quantify the role played by each one.
Acknowledgements
We would like to thank the anonymous reviewers for their constructive criticism of the text. Prin 2008 “Tectonic trajectories of subducted lithosphere in the Alpine collisional orogen from structure, metamorphism and lithostratigraphy”, CNR-IDPA and PUR 2008 “La ricerca geofisica: esplorazione, monitoraggio, elaborazione e modellazione” are gratefully acknowledged.