Ylona van Dinther

The underestimation of the size of recent megathrust earthquakes illustrates our limited understanding of their spatiotemporal occurrence and governing physics. To unravel their relation to associated subduction dynamics and long-term deformation, we developed a 2D continuum viscoelastoplastic model that uses an Eulerian-Lagrangian finite difference framework with similar on- and off-fault physics. We extend the validation of this numerical tool to a realistic subduction zone setting that resembles Southern Chile. The resulting quasi-periodic pattern of quasi-characteristic M8-M9 megathrust events compares quantitatively with observed recurrence and earthquake source parameters, albeit at very slow coseismic speeds. Without any data fitting, surface displacements agree with GPS data recorded before and during the 2010 M8.8 Maule earthquake, including the presence of a second-order flexural bulge. These surface displacements show cycle-to-cycle variations of slip deficits, which overall accommodate ~5\% of permanent internal shortening. We find that thermally (and stress) driven creep governs a spontaneous conditionally stable downdip transition zone between temperatures of ~350C and ~450C. Ruptures initiate above it (and below the fore-arc Moho), propagate within it, interspersed by small intermittent events, and arrest below it as ductile shearing relaxes stresses. Ruptures typically propagate upward along lithological boundaries and widen as pressures drop. The main thrust is constrained to be weak due to fluid-induced weakening required to sustain regular subduction and to generate events with natural characteristics (fluid pressures of ~75-99% of solid pressures). The agreement with a range of seismological, geodetic and geological observations demonstrates the validity and strength of this physically consistent seismo-thermo-mechanical approach.

The physics governing the seismic cycle at seismically active subduction zones remains poorly understood due to restricted direct observations in time and space. To investigate subduction zone dynamics and associated interplate seismicity, we validate a continuum, visco-elasto-plastic numerical model with a new laboratory approach (part 1). The analogous laboratory setup includes a visco-elastic gelatin wedge underthrusted by a rigid plate with dened velocity-weakening and -strengthening regions. Our geodynamic simulation approach includes velocity-weakening friction to spontaneously generate a series of fast frictional instabilities that correspond to analog earthquakes. A match between numerical and laboratory source parameters is obtained when velocity-strengthening is applied in the aseismic regions to stabilize the rupture. Spontaneous evolution of absolute stresses leads to nucleation by coalescence of neighboring patches, mainly occurring at evolving asperities near the seismogenic zone limits. Consequently, a crack-, or occasionally even pulse-like, rupture propagates toward the opposite side of the seismogenic zone by increasing stresses ahead of its rupture front, until it arrests on a barrier. The resulting surface displacements qualitatively agree with geodetic observations and show landward and, from near the downdip limit, upward interseismic motions. These are rebound and reversed coseismically. This slip increases adjacent stresses, which are relaxed postseismically by afterslip and thereby produce persistent seaward motions. The wide range of observed physical phenomena, including back-propagation and repeated slip, and the agreement with laboratory results demonstrate that visco-elasto-plastic geodynamic models with rate-dependent friction form a new tool that can greatly contribute to our understanding of the seismic cycle at subduction zones.

"Subduction megathrust earthquakes occur at the interface between the subducting and overriding plates. These hazardous phenomena are only partially understood because of the absence of direct observations, the restriction of the instrumental seismic record to the past century, and the limited resolution/completeness of historical to geological archives. To
overcome these restrictions, modeling has become a key-tool to study megathrust earthquakes. We present a novel model to investigate the seismic cycle at subduction thrusts using complementary analog (part 1) and numerical (part 2) approaches. Here we introduce a simple scaled gelatin-on-sandpaper setup including realistic tectonic loading, spontaneous rupture nucleation, and viscoelastic response of the lithosphere. Particle image velocimetry allows to derive model deformation and earthquake source parameters. Analog earthquakes are characterized by “quasi-periodic” recurrence. Consistent with elastic theory, the interseismic stage shows rearward motion, subsidence in the outer wedge and uplift of the “coastal area” as a response of locked plate interface at shallow depth. The
coseismic stage exhibits order of magnitude higher velocities and reversal of the interseismic deformation pattern in the seaward direction, subsidence of the coastal area, and uplift in the outer wedge. Like natural earthquakes, analog earthquakes generally
nucleate in the deeper portion of the rupture area and preferentially propagate upward in a crack-like fashion. Scaled rupture width-slip proportionality and seismic moment-duration
scaling verifies dynamic similarities with earthquakes. Experimental repeatability is statistically verified. Comparing analog results with natural observations, we conclude that
this technique is suitable for investigating the parameter space influencing the subduction interplate seismic cycle."

"At oceanic margins, syn-convergent exhumation, subduction erosion, and inter-plate coupling are intimately related, but ample questions remain concerning their interaction and individual mechanisms. To analyze these interactions for a thick-skinned, visco-elastic wedge, we focus on properly modeling stresses, energies, and topographies at the inter-plate and wedge bounding interfaces using a Coulomb frictional contact
algorithm. In this innovative plane-strain, free surface, Lagrangian finite element model, fault dynamics is modulated by retreating subduction. Subduction is dynamically driven by slab-pull due to a slab sinking in a semi-analytic, computationally favorable approximation of three-dimensional induced mantle flow. Nodal trajectories show that continuous underthrusting of a slab induces a steady state corner flow through forced underplating and subsequent trenchward extrusion due to gravitational spreading. This flow pattern
confirms early-proposed models of syn-orogenic deep-seated rock exhumation propelled by coexisting extension and continuous shortening at depth. A distinct reduction in upward flowing material and accompanying decrease of exhumation velocities, to millimeters per year as observed in nature, is induced by a diversion of orogenic wedge material toward the mantle once a subduction channel is formed. The key
parameter affecting model evolution and spontaneous formation of a subduction channel is basal friction, which modulates the amount of erosion. However, formation of a subduction channel entrance needs to be ensured through the deformability of the overriding plate, which is influenced by applied pressure at
the overriding plate tip and material properties. The down dragging of the overriding plate is sufficient above a threshold inter-plate shear stress of about 2–7 MPa."

Active convergentmargins are primarily shaped by the interplay among the subducting plate, overriding plate, and mantle. The effect of important forces, like far-field mantle flow, overriding plate motion, and inter-plate coupling, however, remains partially ambiguous. In a preliminary attempt to clarify their role, a self-consistent, viscoelastic, plane-strain, mechanical finite element model, in which subducting plate, overriding plate and mantle interact dynamically, is developed. In this quasi-static framework with a freelymoving slab, trench, and inter-plate fault, the role of a compressive overriding plate on subduction zone kinematics, morphology and stress-state is characterized. A slab interacting solely with a semi-analytical three-dimensional mantle flow formulation shows that local non-induced mantle flow influences slab geometry and kinematics, adding an important dynamic term to the system. The impact of an overriding plate on this system is determined completely by overriding plate trench-ward motions and is only pertinent if the overriding plate actively advances the trench. A trench-ward moving overriding plate indents the slab and thereby enforces trench retreat and decreases slab dip. It also stimulates over-thrusting of the overriding plate onto the slab, and thereby permits mountain building within the overriding plate. Frictional resistance is observed to have a dominant local effect within the overriding plate as it is increasingly dragged down, thereby inhibiting the growth of overriding plate topography. A distinguishable effect on large-scale trench motions and deep slab dip is, however, absent for re-normalized friction coefficients ranging up to about 0.2. Minor additional effects include a decrease in plate motions of about 15% and slab bending stresses of about 10%.

Earthquakes and their faulting dynamics evolve over a broad range of space and time-scales. In this study, we model the long-term earthquake cycle on a convergent margin, starting from geodynamic large space-time scales, sequentially progressing to small space-time scales to cover the range of interest. Our numerical simulations link kinematic observables and long-term deformation phenomena to earthquake occurrences. Observed outer-rise localization features, influenced by both regional and bending stresses, are characterized in terms of seismic properties, like stress drop, absolute and relative motion, and the earthquake recurrence interval. Subsequently, we attempt to quantitatively link these properties to the cycle of large under-thrusting events at the plate interface. To fully quantify earthquake cycle characteristics in a subduction setting we use a plane-strain conservative finite-difference scheme, combined with a characteristics-based marker-in-cell technique (code I2ELVIS). This visco-elasto-plastic, coupled petrological and thermo-mechanical code solves for the conservation of momentum, mass, and energy, and allows for the spontaneous subduction of an oceanic plate below a continental overriding plate into the mantle, while sedimentation, erosion, partial melting, and slab dehydration occur. Plasticity is implemented in a Drucker-Prager-like fashion through a viscosity-like parameter, while the second invariant of the stress tensor is corrected towards the pressure-dependent yield stress. In particular, we investigate slab age, convergence velocity, slab dehydration, and inter-plate coupling to constrain their role in the occurrence of outer-rise events and their link to under-thrusting events. Initial results show the existence of spatially limited clusters of outer-rise intra-slab localizations of plastic strain that relate to slab bending and whose dimensions depend on e.g. slab age and convergence velocity. Stress-strain history of markers over different depths shows that high stress drops roughly follow the occurrence of regions of high strain-rates, confirming their seismic character. In general, this outer-rise seismicity is observed to anti-correlate in time to movements on the inter-plate thrust interface, suggesting a seismic cycle useful for assessing future seismic potential.

The Earth's active convergent margins are characterized by dynamic feedback mechanisms that interact to form an intricate system in which a crustal wedge is shaped and metamorphosed at the will of two large, converging plates. This framework is accompanied by complicated processes, such as seismogenesis and the exhumation of high pressure rocks. To honor the dynamic interaction between different entities and advance on these persisting issues, we model the interaction between the subducting and overriding lithospheres, the mantle and the crustal wedge explicitly, and observe how a crustal wedge evolves in detail within a set of rigid, internally-driven boundary conditions. We model crustal wedge evolution in an intra-oceanic subduction setting by using a plane-strain implicit solid-mechanical Finite Element Model, in which the mechanical conservation equations are solved using the software package ABAQUS. The crustal wedge is modeled as a thick-skinned accretionary wedge of inter-mediate thickness with a linear visco-elastic bulk rheology. The dynamic interaction between the subducting plate, the overriding plate, and crustal wedge is implemented using a Coulomb frictional algorithm. The interaction with the mantle is incorporated using a computationally favorable mantle drag formulation that simulates induced three-dimensional mantle flow. This results in a quasi-static framework with a freely moving slab, trench, and fault, where a weaker wedge deforms in response to self-regulating, rigid boundary conditions formed by single, frictional bounding faults. The self-regulating evolution of crustal wedge architecture follows three phases; 1) initial vertical growth, 2) coeval compression and extension leading to internal corner flow, and 3) a steady-state taper with continuous corner flow. Particle trajectories show that, as shortening continues throughout the second phase, wedge material is constantly forced upward against the backstop, while extension and ocean-ward extrusion occur on top due to gravitational spreading. This dynamically forced corner flow provides a mechanism for the syn-orogenic exhumation of high-pressure oceanic blueschists. This deformation pattern with a vertical reversal in shear transport direction is similar to that observed in analogue viscous paraffin wedges. Resulting exhumation rates are in the order of millimeters a year, regulated by the frictional parameters and wedge strength. A change in fault frictions lead to subduction erosion instead of accretion, as crustal wedge material is dragged down in a spontaneously formed subduction channel. Subduction erosion is essentially promoted by inter-plate friction, since this determines whether a subduction channel entrance can be formed through down dragging of the overriding plate. This down dragging is also promoted by an increase in relative strength of the wedge compared to the overriding plate. Increased basal wedge drag only slightly increases the amount of subducted material. The type of wedge behavior, accretion or subduction erosion, and its amount considerably influences inter-plate stresses and seismic coupling at the inter-plate fault.

Subduction thrust earthquakes occur over different spatial and temporal scales. In this study, the long-term earthquake cycle is investigated in a 4000x200 km ocean-continent subduction setting. Events or periods of rapid deformation, present within spontaneously formed localizations of plastic strain, are collected and quantified for comparison to observations and scaling to nature. The resulting database gives insight into the cyclic nature and timing of large earthquakes. The characteristics of the subduction zone seismic cycle also differ for different convergent margins. Subduction zones that are weakly coupled and receive a lot of sediment are likely to produce large mega-thrust earthquakes, while more and smaller earthquakes are observed next to compressive overriding plates, where the coupled fault receives less sediment. This distinction is tested in our model by producing different databases with events for different overriding plate tectonic regimes, different inter-plate frictions and different amounts of sediments. To model the long-term seismic cycle in a subduction zone we use a plane-strain finite-difference scheme with marker-in-cell technique to solve the conservation of momentum, mass, and energy for a visco-elasto-plastic rheology (code I2ELVIS). In a generic, petrologically complex continent-ocean subduction zone, localizations of plastic strain are formed when the second invariant of the deviatoric stress tensor exceeds the Drucker-Prager yield criterion, leading to a correction to the pressure-dependent yield stress by decreasing a viscosity-like parameter. We assume a seismic event occurs when a sudden considerable stress drop occurs simultaneously with an immediate strong increase in strain rate. For these events we collect estimates of depth, stress drop, and slip area, assuming equi-dimensional slip areas, and thereby acquire a very rough indication of earthquake size. Previous results have shown the existence of several clusters of plastic strain localizations, at the thrust interface, bending outer-rise, back-arc, and accretionary wedge, whose cycle is related to the by far most energetic events at the thrust interface. These thrust events spread over depths of about 15 km and occur about every 10.000 years within a seismogenic zone extending from about 20 to 35 km depth. They also showed the presence of a long-term cycle in dissipated strain energy rate of the order of 30.000 years.