T54B-01 INVITED
A Subduction Scissor and Development of the Intervening Plate Collision at South Island, New Zealand
The continental plate collision across South Island, New Zealand is bounded to the north by west-dipping Hikurangi subduction of the Pacific plate and to the south by east-dipping Fiordland-Puysegur subduction of the Australian plate. As such, South Island is in the midst of a type of 'subduction scissor'. The tectonic behaviour of the deep plate portion of South Island (viz., the sub-crustal portion of the lithosphere) is unresolved. Previous studies suggest it may be undergoing subduction-like consumption or Rayleigh-Taylor-type viscous instability. We consider these possibilities but in the three-dimensional tectonic context of the subduction scissor. In particular, using 3D computational geodynamic models we explore how the evolution of the continental plate collision is controlled by the adjacent and opposite-facing subduction zones. In the absence of the subduction forcing, the subduction/drip structure develops as a near-vertical sheet. With imposed subduction on the edges, the plate consumption is still steep, but there is appreciable flow along the strike of the collisional plate boundary that results in localized accumulations and gaps in the consumed plate. The response of the model South Island to the subduction scissor depends largely on the rates of subduction and the continental lithospheric rheology. Recent interpretations of the morphology of deep lithospheric structure beneath South Island from geophysical observables indicate that there may be similar along-strike variations in the consumed plate/blob.
T54B-02
Strong plate coupling along the Nazca/South America convergent margin
The force balance in plate tectonics is fundamentally important, but poorly known. Much information on the dynamics is embedded in the record of past and present plate velocities, featured with long- as well as short- term variations, but a precise budget, in particular of resistive coupling forces along convergent margins, is hard to come by. Building on substantial, yet separate progress in modeling lithosphere dynamics and mantle convection, we couple global lithosphere models with high-resolution (more than 100 million grid points) 3-D circulation models of Earth's mantle and demonstrate that an accurate budget of plate boundary forces can be obtained. We prove the effectiveness of our approach by computing a detailed force budget along the Nazca/South America subduction zone, showing that a large portion of it comes from the recent uplift of the Andes. We find that forces computed with our global, coupled models provide simultaneous explanations for three seemingly unrelated key observations along the South American margin: (1) trench parallel gravity anomalies, (2) pronounced bathymetry variations, as well as (3) a substantial reduction in Nazca/South America plate convergence recorded over the past 10 million years. All these observations can be explained from along- trench, lateral and temporal variations in plate coupling forces that are predicted from our simulations. Interestingly enough, the distribution of great earthquakes such as the recent M 8.0 event of Peru coincides with moderate to low coupling between subducting and overriding plates. For the same convergent margin we also show that frictional forcing due to trench sediment infill is, by comparison, of minor importance. Finally, we provide an intriguing explanation for the peculiar convex shape of the South American margin. Paleomagnetic and geodetic data indicate substantial rotation over the past m.y. and continuing at present day. We tie the bend of the margin to variations in plate coupling associated with Andean growth. Specifically, the symmetric distribution of plate coupling forces along the margin provides strong torques that could have caused the bend in a feedback mechanism between mountain belt growth and plate coupling.
T54B-03
Continental subduction, surface processes, P-T-t-z conditions and unstable vs. stable plate dynamics: Insights from thermo-mechanical modelling
We analyze major mechanisms of shortening of continental lithosphere (simple shear subduction, pure shear collision, folding, Rayleigh-Taylor instabilities). We use a thermo- dynamically coupled thermo-mechanical numerical model that accounts for brittle-elastic-ductile rheology, surface processes and metamorphic phase changes. The model also traces P-T-t-z paths of metamorphic facies that can be compared with petrology data. The experiments suggest that continental subduction occurs in case of relatively strong lithospheres with a competent mantle part (TMoho < 550° C) , at relatively high initial convergence rates (> 1.5-5 cm/yr). Depending on the lower-crustal rheology (strong or weak), either the entire (upper and lower) crust or only the lower crust can be involved in subduction. Pure shear collision is dominant when TMoho > 550° C or convergence rates are lower than 1.5-3 cm/yr (subduction number, S > 0.5). Large-scale folding is favored in case of TMoho=500-650° C and is more effective in case of mechanical coupling between the crust and mantle (e.g., strong diabase lower crust). Gravitational R-T instabilities overcome other mechanisms for very high values of TMoho (>800° C) and lead to the development of subvertical "cold spots." In case of weak metamorphic rheologies, phase changes improve chances for stable subduction. In general, exhumation of UHP-HP rocks to the surface is favored if the crustal rheological profile is characterized by two internal ductile decolement levels (between the upper and lower or intermediate crust and the lower crust and mantle lithosphere). Finally, we investigate the impact of surface processes (erosion/sedimentation) on the amount of continental subduction. The maximal amount of subduction is achieved for intermediate erosion rates when tectonic uplift rates are fine-balanced by denudation rates. In case of India-Asia collision-like scenario (fast convergence > 5 cm/y, stiff lower plate), the optimal balance is achieved for k ~ 3000 m2/yr. The experiments suggest that both extra slow (k < 50-100 m2/yr) and extra rapid erosion (k > 6000-8000 m2/yr) limit, by up to 50%, the total amount of subduction, if not totally prevent it. We suggest that most orogenic belts could have started their formation from continental subduction, yet, in case of slow convergence (< 3 cm/yr) or weak lithosphere the subduction channel locks up after about few Myr, and subduction is then relayed by a different deformation mode such as pure shear collision. In case of India-Asia-like convergence settings, continental subduction may continue for tens of Myr allowing for subduction of 600-800 km continental "slab".
T54B-04 INVITED
Common Observables of Trench Migration and Plate Motion in Different Global Reference Frames
Plate velocities and trench migration velocities are commonly described in some sort of global "absolute" reference frame. From calculating such motions for all plates and subduction zones on Earth, one might obtain insight into the importance of various driving and resistive forces of plate tectonics and plate boundary migration. Trench migration velocities and plate velocities have been calculated for all subduction zones on Earth in eight global reference frames. The calculations show that such velocities can differ substantially between different global reference frames (up to 4 cm/yr), in particular between one Pacific hotspot reference frame (HS3- NUVEL1A) and all the others. In addition, this reference frame shows a bimodal distribution of trench velocities, while all the others show a Gaussian distribution. Nevertheless, some common features are observed irrespective of the reference frame. First, trench retreat always dominates over trench advance, with 62-78% of the trench segments retreating, while the mean and median trench velocities are always positive (retreating). Second, trench retreat is always slow in the middle of wide subduction zones, i.e. far (>2000 km) from lateral slab edges (<2 cm/yr in seven reference frames). Third, fast trench retreat (>6 cm/yr) is only found close (<1500 km) to lateral slab edges. Fourth, plates with a substantial percentage of their circumference attached to a subducting slab (Pacific, Nazca, Cocos, Philippine, Australia) move trenchward. These calculations are predicted by three-dimensional geodynamic models of free subduction with a variable slab width (300-7000 km), in which the slab to upper mantle viscosity ratio is low (100-200). This suggests that trench velocities and plate velocities are indeed primarily controlled by the negative buoyancy and width of subducting slabs. It further suggests that slab/upper mantle viscosity ratios in nature are 100-200, as the models show trench motion dominated by retreat, and forward dipping "slab-draping" geometries (e.g. Tonga) or steep slab geometries (e.g. Kermadec) as observed in nature. Overturned "roll-over" slab geometries are not observed.
T54B-05
Dynamics of plate bending at the trench and slab-plate coupling
The bending strength of subducting lithosphere plays a critical role in the Earth's plate tectonics and mantle convection. It is the subject of a lively debate how much of the potential energy of the downgoing plate is consumed in bending the plate, and how much of the lithospheric strength is lost in the process. This modulates the amount of slab pull transmitted to the surface, thus setting the boundary conditions under which plates move and deform. We model the subduction of a viscoelastic lithosphere, driven solely by the downgoing plate's buoyancy, freely sinking in the passive mantle, represented by drag forces. To investigate the dynamics of bending, we vary the viscosity profile within the plate, from isoviscous, where the strength is distributed, to strongly layered, where the strength in concentrated in a thin core, and map stress, strain and forces along the downgoing plate. We found that stress and strain distributions and thus the force transmitted to the surface depend strongly on the viscosity profile. While isoviscous plates deform throughout and bending is accompanied by significant thinning, layered plates keep their integrity by localizing strain in the weak outer layers, and stress in the strong core, allowing stress propagation through the bend far into the plate. The low dissipation in both rheological models, the stress magnitude largely below the plastic limit, as inferred from rock mechanics, and the continuity of stress through the bend into the plate suggest strongly that layered lithosphere keeps its integrity while bending. The rheological model presented offers a consistent description of the dynamics of convergent margins and an its relevant properties to global plate tectonics modeling.
T54B-06
Numerical Modeling of Free Subduction: a Boundary-Element Approach
Many important aspects of subduction can be understood by studying an idealized model of purely buoyancy-driven ("free") subduction of an isolated thin fluid sheet with constant viscosity in an ambiant fluid with a lower viscosity. The problem can be efficiently solved using the boundary-element method, wherein a Stokes flow in a given region of space is represented in terms of the velocities and tractions at its boundaries. Essential advantages of the method include reduction of the dimensionality of the problem (from 3D to 2D or 2D to 1D); the capacity to implement a true free surface; accurate interface tracking; and the complete absence of wall effects (unless they are desired.) Moreover, the two-fluid configuration implemented by the method closely resembles that of many recent laboratory experiments on subduction. The study I will describe proceeds hierarchically from simple models to more complex ones. The simplest reference case is free subduction of a 2D sheet in an infinite half-space, with surface topography included to ensure the resolution of Stokes's paradox. Next, an impermeable boundary at 660 km depth is included, and then a viscosity jump. Finally, the same sequence is followed for a 3D sheet of finite width, for which toroidal flow around the sides of the slab is important. I will present a number of quantitative scaling laws for key subduction parameters such as the surface plate speed, the rate of trench rollback, and the slab geometry, and compare their predictions with relevant geophysical observations and with laboratory experiments (e.g., those of the group at the University of Rome-III.)
T54B-07
The Role of Viscous Dissipation on the Thermal Structure of Subduction Zones
The thermal structure of subduction zones is related to mechanisms of magma generation and surface heat flow. For active arc volcanism and high heat flow in the back arc, the presence of high temperature magma (> 1200 °C) and mantle wedge are essential. To explain the thermal structure, frictional heating, exothermic metamorphism, and radiogenic heat production have been suggested but, they are insufficient to increase the temperature of the mantle wedge. Experiments using compressible fluid approximations show significant heat generation in the regions of downwelling (subduction) from viscous dissipation; however, most subduction experiments use incompressible fluid approximations (e.g. Boussinesq Approximation, BA hereafter) so that the effects of viscous dissipation on the thermal structure can not be assessed. In this study, we quantitatively evaluate the effects of viscous dissipation on the thermal structure of subduction zones using a compressible fluid approximation (Truncated Anelastic Liquid Approximation, TALA hereafter). We use a 50 km-depth overriding plate over the mantle wedge and a kinematically driven subducting slab with the dip of 45 degrees, the age of 130 Ma and the subduction velocity of 5 cm/year. We use the rheology of dry and wet olivine for the viscous mantle wedge by assuming that the mantle wedge is wet down to 200 km. At the corner of the mantle wedge, low viscosity of serpentinite and partial coupling between the subducting slab and the mantle wedge are assumed down to 70 km and 120 km, respectively. We vary the mantle rheology by considering diffusion creep, combined creep of diffusion and dislocation, and constant viscosity. Uniform radiogenic heat production of 3×10- 10 J/kg·s is assumed in the upper 20 km of the lithosphere. For comparison, we conduct the corresponding experiments using BA with the same parameters. Slight increases in temperature of the mantle wedge and the subducting slab as well as increases in surface heat flow are observed in the experiments using the rheology of diffusion creep and constant viscosity. However, there are significant differences between the thermal structures in the experiments using the rheology of the combined creep of diffusion and dislocation. The TALA experiments develop high temperature cornerflow but the BA experiments do not. The BA experiments produce low heat flow values in the arc. This suggests that viscous dissipation plays an important role in the thermal structure of subduction zones because dislocation creep is thought to be a major flow mechanism in the upper mantle. Like previous studies, we do not observe high heat flow in the back arc and thus require further study.
T54B-08
A new Remesh-Lagrange technique for advecting temperature that minimizes numerical diffusion
The proper treatment of heat-advection is a generally underappreciated problem within CFD, yet particularly critical for calculating physically sound erosion in plume-lithosphere interactions and temperature sensitive melting processes. Typically, Eulerian (fixed-mesh) codes have been preferred to solve for fluid flow and they are almost essential for finite-difference-based algorithms. Unfortunately, the Eulerian approach introduces numerical artifacts into the solution of the advection-diffusion heat transport problem that can only be suppressed by adding 'too-diffusive' artificial diffusion to the equations, as for example in the Smolarkiewicz formulation for heat advection. We have developed a 'Remesh-Lagrange' method using a partly deforming finite element mesh and find it to be significantly more accurate than our previous methods. In several test scenarios we show the large improvement in accuracy that can be obtained by using a Lagrangian approach for 10-30 time steps (depending upon the distortion of the finite elements in the deformed Lagrangian mesh) and then regridding to the initial mesh. When an element becomes too distorted the nodes connected to it become fixed and we switch from Lagrange to a Semi-Lagrange formulation for these nodes. Instead of the standard 'linear backward' Semi-Lagrange we are also experimenting with a more accurate interpolation scheme for an unstructured mesh that additionally includes the nodal derivatives of the temperature field when calculating the value at the Semi-Lagrange traceback point. The same bicubic interpolation method for an unstructured grid is used to remesh the 'too-distorted' Lagrange grid back to the initial undistorted mesh. We compare the Remesh-Lagrange technique against the following Eulerian methods in a series of 2-D numerical experiments advecting stripes and Gaussian peaks in steady circulating flow: linear back-interpolation Semi-Lagrange method; bicubic back-interpolation Semi-Lagrange method; SUPG; and the Smolarkiewicz (flux- limiting diffusion) method, which is arguably the 'best-in-practice' of the Eulerian finite-difference methods. We also compare these methods for a test-case of diapiric upwelling that is forced by a pure chemical (non-diffusive) buoyancy associated with tracer particles that are advected by the plume upwelling. In this case, we find that the Eulerian methods typically underpredict the temperature of the plume head by about 20C, while overheating a ~50km region outside the plume-stem by about 10C in comparison to the Remesh-Lagrange technique.