Lithospheric structure of the Pampean flat slab region from double-difference tomography

Journal Pre-proof Lithospheric structure of the Pampean flat slab region from double-difference tomography Lepolt Linkimer, Susan Beck, George Zandt, Patricia Alvarado, Megan Anderson, Hersh Gilbert, Haijiang Zhang PII:

S0895-9811(19)30367-0

DOI:

https://doi.org/10.1016/j.jsames.2019.102417

Reference:

SAMES 102417

To appear in:

Journal of South American Earth Sciences

Received Date: 23 August 2019 Revised Date:

8 November 2019

Accepted Date: 12 November 2019

Please cite this article as: Linkimer, L., Beck, S., Zandt, G., Alvarado, P., Anderson, M., Gilbert, H., Zhang, H., Lithospheric structure of the Pampean flat slab region from double-difference tomography, Journal of South American Earth Sciences (2019), doi: https://doi.org/10.1016/j.jsames.2019.102417. 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.

Submitted to Journal of South American Earth Sciences.

Lithospheric Structure of the Pampean Flat Slab Region from Double-Difference Tomography Lepolt Linkimer1, Susan Beck2, George Zandt2, Patricia Alvarado3, Megan Anderson4, Hersh Gilbert5, Haijiang Zhang6 1

Escuela Centroamericana de Geología, Universidad de Costa Rica, San José, Costa Rica. E-mail:

[email protected] 2

Department of Geosciences, University of Arizona, Tucson, AZ, USA

3

Centro de Investigaciones de la Geósfera y Biósfera (CONICET-Universidad Nacional de San

Juan). Departamento de Geofísica y Astronomía, Facultad de Ciencias Exactas, Física y Naturales, UNSJ, San Juan, Argentina 4

Washington Geological Survey, Olympia, WA, USA

5

Department of Geosciences, University of Calgary, Calgary, AB, Canada

6

Laboratory of Seismology and Physics of Earth’s Interior, School of Earth and Space Sciences,

University of Science and Technology of China, Hefei, Anhui, China

SUMMARY We obtain earthquake locations and a detailed three-dimensional velocity model of the flat slab subduction zone in west-central Argentina (latitudes: 32-30ºS and longitudes: 70-66ºW) using a regional-scale double-difference tomography algorithm with earthquake data recorded by the SIEMBRA (2007-2009) and ESP (2008-2010) broadband seismic networks. In this region, the flat subduction of the Nazca Plate including the Juan Fernandez Ridge is spatially correlated with a volcanic gap and the basement-cored uplifts of the Sierras Pampeanas in the overriding South American Plate. Our results show the subducting Nazca Plate as a continuous band of mostly increased P-wave velocities coinciding with the Wadati-

2 Benioff Zone. In the overriding South American Plate, the lithospheric mantle appears to be heterogeneous but mostly characterized by a ratio between P- and S-wave velocities (Vp/Vs) of 1.75-1.77, which is consistent with depleted peridotites. Two Vp/Vs anomalies deviate from this mantle with lower (1.70-1.73) and higher (1.78-1.82) Vp/Vs, which are interpreted as localized dry and hydrated regions, respectively. The lower Vp/Vs is consistent with an enrichment of 40-80 % of orthopyroxene and the higher Vp/Vs with up to 5% mantle hydration. The size, orientation, and location of these seismic anomalies suggest the progressive eastward dehydration of the subducting slab and the presence of an east-dipping large-scale lithospheric suture, which is interpreted as evidence of an ancient subduction zone and also as a weak zone that facilitates the hydration of the upper plate. Our inversion results suggest a thicker South American crust in the Western Sierras Pampeanas and the partial eclogitization of the lower crust beneath that region where velocities match three types of eclogites at depths of 40-60 km. In the middle-to-upper crust, velocities are reduced in the Precordillera and Vp/Vs is higher in the Cuyania and Chilenia terranes (>1.75) than in the Pampia terrane (1.67-1.75). These observations are consistent with the presence of a thick carbonate sequence in the Precordillera, mafic-ultramafic rocks in Cuyania and Chilenia, and felsic rocks in Pampia. The higher variability in Vp/Vs and strong velocity changes at crustal depths within the Precordillera and the Cuyania Terrane agree with more complexity in crustal structure for these regions and reveal two mid-crustal discontinuities as well as the Chilenia-Precordillera suture zone. Finally, the relocated slab earthquakes refine the slab geometry and suggest that at depths of ~100 km, the flat slab segment is ~240 km wide and has a slight westward dip (~2º) before it resumes its descent into the mantle with a steep angle (~25º). The observation of a wider flat slab segment than the width of the Juan Fernandez Ridge offshore (~100 km) implies that there might be additional contributing factors for the flattening besides the subduction of the overthickened oceanic crust.

3 Key words: Seismic tomography, Subduction zone processes, Crustal structure, Composition and structure of the mantle, South America.

1

INTRODUCTION

In the south-central Andes (latitude 32-30ºS), flat subduction of the Nazca Plate, including the Juan Fernandez Ridge (JFR), coincides with the cessation of arc volcanism and the basement-cored uplifts of the Sierras Pampeanas within the South American Plate (e.g. Pilger 1981; Ramos et al. 2002) (Fig. 1). This geographic correlation makes the Pampean flat slab region an ideal place to investigate the role of water in subduction zone processes and the effects of flat slab subduction on the overriding plates. The origin of the Pampean flat slab has been traditionally attributed to the subduction of the JFR, which marks a hotspot track and may have provided the additional buoyancy that prohibits the slab from sinking into the mantle (Yáñez et al. 2001, 2002; van Hunen et al. 2002; Espurt et al. 2008). A high degree of hydration (e.g. Porter et al. 2012) and a delay in the basalt-eclogite transition during subduction (e.g. Gutscher et al. 2000b) may facilitate the extra buoyancy of a hotspot-generated, overthickened oceanic crust. However, the association between flat slab regions and bathymetric anomalies has also been debated by Skinner & Clayton (2013) who showed little spatial correlation in current examples globally and in plate reconstructions for subduction zones in South America. A comparison of several flat slab regions worldwide shows an absence of a unique set of conditions that provide an explanation for the flat slabs (Manea et al. 2017). Instead, there are many potential contributing factors including: the hydrous weakening of the mantle wedge (van Hunen et al. 2002), the trenchward motion of a thickened overriding lithosphere and suction forces in the mantle wedge (van Hunen et al. 2004; Manea et al. 2011; Antonijevic et al. 2015), a combination of seafloor age variation, hydrodynamic suction from the overriding plate and subduction of buoyant features (Hu et al. 2016; Huangfu et al. 2016), a source of buoyancy in the subducting mantle

4 (Bishop et al. 2017), a faster trench retreat than the slab rollback (Schepers et al. 2017), and an over-pressured sub-slab mantle that lifts the slab (Boutelier & Cruden, 2008). Volcanic activity in the main Andean Cordillera above the flat slab region terminated about 9 Ma as the flattening slab presumably squeezed out the mantle wedge (Kay & Mpodozis 2002; Ramos et al. 2002). Besides the cessation of volcanism (e.g. McGeary et al. 1985) other significant effects attributed globally to the slab flattening include: the cooling of both the overriding and subducting plates (Gutscher et al. 2000a), the decrease in subducting-plate velocity (Martinod et al. 2005; Espurt et al. 2008), and the removal of a significant part of the upper-plate mantle lithosphere (Bishop et al. 2017). Some of these effects favor a strong tectonic uplift and shortening (e.g. Allmendinger et al. 1990; Martinod et al. 2013), and an increase in plate coupling, which enables an increase in seismic energy release and deformation to be transferred far into the upper plate (Scholz & Small 1997; Gutscher, et al. 2000b; Alvarado et al. 2007; Richardson et al. 2012, 2013; Bishop et al. 2017). Understanding the lithospheric structure of the current Pampean flat slab segment provides a framework for studying previous flat slab dynamics, such in western North America, where many have hypothesized that the Farallon Plate subducted beneath the North American Plate at a shallow angle during the Laramide orogeny (e.g., Henderson et al. 1984; Jordan & Allmendinger 1986). Determining whether the flat slab retains or releases fluids during subduction is also key to the study of the deformation in the overriding plate as fluids decrease the bulk rock density, increase buoyancy, and reduce rock strength (e.g. Reynard 2013). Additionally, studying the upper-plate structure is important for understanding the regional seismic hazard, as large and damaging crustal earthquakes (e.g. 1894, 1944, 1977, 1985, 2002, 2004) have occurred in the Precordillera and Sierras Pampeanas (e.g. Alvarado & Beck 2006; Alvarado et al. 2009; Alvarado & Ramos 2011; Rivas et al., 2019). Previous tomography studies in the Pampean flat slab region have outlined the first-order larger-scale structure (e.g. Pardo et al. 2002; Wagner et al. 2005, 2006, 2008; Deshayes 2008;

5 Marot et al. 2014). In particular, Wagner et al. (2008) identified a large zone (300 x 120 km) of low (1.65-1.72) P- to S-wave velocity ratio (Vp/Vs) in the mantle above the flat slab, suggesting dry conditions despite the high rate of slab seismicity (e.g. Anderson et al. 2007). More recently, Marot et al. (2014) also found a low Vp/Vs region but the extent was smaller (~40 x 70 km). This is important to understand well because seismicity may be linked to the amount of crystallographically bound water in hydrated slabs (Hacker et al. 2003b). If a substantial amount of water is present in the upper plate above where it was released by the slab, seismic tomography should exhibit characteristics consistent with hydrous mantle rocks, such as serpentinite, which may have higher Vp/Vs (e.g. Carlson & Miller 2003; Hyndman & Peacock 2003; Christensen 2004; Bezacier et al. 2010). The unexpected low Vp/Vs found above the flat slab seismicity provided the motivation for a higher-density seismic deployment to produce high-resolution tomographic images and study the hydration state in the Pampean flat slab region. In this context, our goal is to present new, higher-resolution images of the lithosphere above the eastern portion of where the Nazca plate subducts nearly horizontally into the mantle. To this end, we apply the double-difference (DD) tomography method to seismic data acquired by two temporary broadband seismic networks: SIEMBRA (December 2007-November 2009) and ESP (August 2008-August 2010) (Fig. 1). The number and distribution of these seismic stations make our temporary seismic array the densest and longest ever used for a tomography study in the Pampean flat slab. As a starting point for the three-dimensional (3D) modeling, we derived a onedimensional (1D) velocity model. In addition, we performed a series of tests to investigate the resolution and sensitivity to the model parametrization. We present and interpret our results in three sections: the subducting Nazca Plate, the South American mantle lithosphere, and the South American crust. For the discussion, we compare our absolute seismic velocities with those predicted for subduction-type rocks from Hacker et al. (2003a) and Abers & Hacker (2016) at the appropriate pressure (P) and temperature (T) conditions. We also compare our results with those of Wagner et al. (2008) and Marot et al. (2014) who both performed a similar analysis but using a smaller data

6 set inverted over a larger region, as well as other previous studies using different geophysical techniques. Our final velocity model shows the finest-scale detail of the eastern segment of the Pampean flat slab region to date.

2

TECTONIC SETTING

The western edge of the South American Plate along the Andean Cordillera displays significant along- and across-strike variations in topography, crustal thickness, crustal shortening, WadatiBenioff Zone (WBZ), and magmatism (e.g. Barazangi & Isacks 1976; Allmendinger et al. 1990; Cahill & Isacks 1992; Kay & Abruzzi 1996). In the south-central Andes, the continental margin is divided into three along-strike segments. North of 28ºS and south of 33ºS the subducting Nazca Plate dips ~30º beneath the active volcanic arc. In contrast, between 33ºS and 28ºS, the Nazca Plate subducts nearly horizontally at depths of ~100 km and the upper-plate volcanism is absent (e.g. Smalley & Isacks 1987; Cahill & Isacks 1992; Ramos et al. 1996). This segment is spatially correlated offshore with the subduction of the JFR, a ~100-km-wide hotspot chain on the Nazca Plate (von Huene et al. 1997; Yáñez et al. 2001; 2002), and onshore with the overthickened lithosphere of the Rio de la Plata Craton (Rapela et al., 2007). A westward shift of the highest peaks of the Andes and the thick-skinned, basement-cored uplifts of the Sierras Pampeanas coincide with the subduction of the JFR (e.g. Pilger, 1981; Ramos et al. 2002). Different narrow, north-trending morphotectonic units extend along the Andes. From west to east, they are: the Principal and Frontal Cordillera, Precordillera, and Sierras Pampeanas (Fig. 1b). The Principal and Frontal Cordillera, which correspond to the main part of the Andes, are mostly composed of Mesozoic sedimentary rocks and Paleozoic-Triassic magmatic rocks (Ramos et al. 1996). The Precordillera is a thin-skinned, fold and thrust belt developed on Paleozoic carbonate and younger rocks (Baldis et al. 1984) and is bound to the west by a belt of mafic-ultramafic rocks (e.g. Ramos et al. 1986; Pérez-Luján et al. 2015). The Sierras Pampeanas correspond to thickskinned, basement-cored uplifts mostly composed of crystalline Precambrian-Early Paleozoic rocks

7 bound by reverse faults and separated by relatively undeformed broad basins (e.g. GonzálezBonorino 1950; Jordan & Allmendinger 1986; Ramos et al. 2002). The western and eastern Sierras Pampeanas have different rock compositions. In the western Sierras Pampeanas (e.g. Pie de Palo, Fig. 1b), metamorphic rocks predominate over granites (e.g. Vujovich et al. 1994). In contrast, granitic rocks dominate the Eastern Sierras Pampeanas (e.g. Sierra de Chepes, San Luis, and Cordoba; e.g. Baldo et al. 1996; von Gosen 1998; Leal et al. 2003; Dahlquist & Galindo 2004). A significant east to west variation is also observed in crustal thickness: the Moho depth is at: ~70 km beneath the Principal and Frontal Cordilleras, ~50 km in the western Sierras Pampeanas, and ~35 km in the eastern Sierras Pampeanas (Gilbert et al. 2006; Gans et al. 2011; Perarnau et al. 2012). In addition, more complexity in crustal structure is observed beneath the Cuyania and Chilenia terranes than in the Pampia Terrane (e.g. Gilbert et al. 2006; Gans et al. 2011; Ammirati et al. 2013, 2015, 2016, 2018). The western part of the South American Plate is also a mosaic of accreted terranes, with generally north-south trending boundaries accreted to the much older Rio de la Plata craton (e.g. Ramos et al. 1986; Astini et al. 1995). From west to east, this includes the Chilenia, Cuyania, and Pampia terranes, which are bound to the east by the Rio de la Plata craton (Fig. 1b). Some studies interpret Cuyania as a composite terrane formed by the amalgamation of the Precordillera and Pie de Palo terranes (Ramos 2004). All these terranes are remnants of the tectonic evolution of the proto-Andean margin of Gondwana during the Paleozoic (e.g. Ramos et al. 1986; Rapela et al. 1998). The Pampean flat slab has been spatially associated to the subduction of the JFR (e.g. Pilger 1981). The geometry of the JFR with respect to the direction of plate convergence has led to the southward migration of the ridge-trench intersection over time, decreasing from a rapid (~20 cm/yr) southward migration during the Miocene to ~3.5 cm/yr during the last 11 Ma (Yáñez et al. 2001). The arrival of the JFR at the current latitude occurred 10-12 Ma (Yáñez et al. 2001). In the upper plate, arc magmatism progressively migrated eastward, and shut off completely in the Main

8 Cordillera and Precordillera between 7 and 11 Ma (Kay et al. 1988). Further east, in the Sierras Pampeanas, the most recent volcanism occurred in the Pocho volcanic field at ~7.9-4.5 Ma (Gordillo & Linares 1982; Kay et al. 1988; Kay & Gordillo, 1994) and in the Sierra del Morro at ~9.5-1.9 Ma (Urbina et al. 1997).

3

DATA

We used data from two temporary broadband deployments: SIEMBRA and ESP (Fig. 1 and Table S1). The SIerras pampeanas Experiment using a Multicomponent BRoadband Array (SIEMBRA) consisted of 43 broadband stations (Beck & Zandt 2007) and the Eastern Sierras Pampeanas (ESP) experiment included 12 (Gilbert 2008). We built an initial earthquake catalogue using both the SIEMBRA and ESP experiments. By using a short-term average over a long-term-average trigger algorithm within the Antelope software package developed by Boulder Real Time Technologies, we captured 69,691 events detected by > 8 stations. To build the earthquake catalogue for the tomography, we then selected only the events with the largest number of station detections using a trigger threshold value of > 30 stations and completed this catalogue with all reported felt events from the catalogues of the Argentinean Instituto Nacional de Prevención Sísmica (INPRES) and the Preliminary Determination Epicenter (PDE-NEIC) of the U.S. Geological Survey National Earthquake Information Center. In order to achieve good ray coverage distribution for the tomography in areas with sparse seismicity, we completed this catalogue by adding the events with > 10 stations in those areas. The earthquake catalogue built for the tomography contains 1,562 events. We handpicked Pand S-wave phases on all these earthquakes, weighted the wave arrivals based on the uncertainty caused by noise (Fig. S1), and used these weights in the inversion. The average P- and S-wave reading uncertainties of our data set were ±0.06 s and ±0.10 s, respectively. More details about the earthquake data and field deployments are available in the Supporting Information.

9 A total of 54,414 P-wave and 30,576 S-wave observations constrain the initial earthquake locations, computed with the programs MULPLT and HYP (Lienert & Havskov 1995) included in the software package SEISAN (Havskov & Ottemöller 1999). For the velocity inversions, we selected a subset of 1,157 well-located earthquakes (Fig. 2). These events have an azimuthal gap of observations (GAP) ≤ 180º, have at least 10 P-wave and 6 S-wave observations, and include at least one observation within a distance of 1.5 times their focal depth. This data set has mostly earthquakes with depths from 0 to 40 km and from 100 to 200 km (Fig. 2d).

4

MINIMUM 1D MODEL

We used the concept of the Minimum 1D Model (M1DM) introduced by Kissling (1988) to generate the initial velocity model and the high precision hypocenter locations required as input for the 3D earthquake inversion (Kissling, 1988; Kissling et al. 1994, 1995). The M1DM incorporates the iterative simultaneous inversion of velocity and hypocenter parameters, with constraints from station corrections that account for lateral velocity heterogeneity. Station CONE (Fig. 3b, Table S1) is the reference station for relative station corrections, because it lies toward the center of the array, it operated during the whole period of the data used here, has a high signal-to-noise ratio, and it is located on hard rock. To find the M1DM, we used the program VELEST (Kissling et al. 1995) using a total of 40,729 P-wave observations from 1,157 well-located earthquakes. We first performed inversions using starting model based on previous studies (e.g. Smalley et al. 1993; Alvarado et al. 2005; Gilbert et al. 2006; Anderson et al. 2007), and later tested a range of layer geometries following a trial-and-error approach. The model that resulted in the minimum root mean square (RMS) has layer interfaces at depths of 0, 20, 40, 60, and 100 km with velocities of 6.25, 6.46, 7.19, 8.10, and 8.32 km/s, respectively (Fig. 3a). The final RMS and variance for the data set were 0.186 s and 0.038 s respectively, contrasting with the initial values of 0.430 s and 0.202 s, respectively.

10 In order to study the dependence of the solution on the initial model, we conducted additional inversions both increasing and decreasing the initial velocities by 0.5 and 1.0 km/s for each layer. All our results show a strong convergence to the velocity model noted above, which we considered the M1DM for this data set (Fig. 3a red line). Following Husen et al. (1999) we conducted several tests by randomly shifting hypocenters by a maximum of 15 km in any direction in order to assess the stability of the final M1DM. Using these randomly shifted starting locations we still retrieved the original hypocenter locations, indicating that they are not systematically biased (Fig. S2).

5

DOUBLE-DIFFERENCE TOMOGRAPHY

We used a regional-scale double-difference DD tomography algorithm (tomoFDD) (Zhang et al. 2004; Zhang & Thurber 2006), which is based on the DD location method of Waldhauser & Ellsworth (2000). This algorithm jointly solves for 3D velocity structure and earthquake locations using a combination of absolute arrival times and differential times, leading to a substantially higher resolution in the source region compared to that obtained via conventional tomography (Zhang & Thurber 2007). The absolute arrival times correspond to the ones directly picked on the seismograms. Our absolute arrival-time catalogue contains 40,729 P- and 25,183 S-wave observations from 1,157 earthquakes recorded at 55 seismic stations. Subtracting absolute travel times from event pairs with inter-event distances of < 20 km at common stations located within 400 km of the centroid of the earthquake cluster creates a differential travel time catalogue. In this catalogue, each event is linked to a maximum of 10 neighboring events by at least 8 pairwise observations, which results in an average separation between strongly linked events of ~9 km. Our resulting differential time catalogue comprises 170,241 observations for P waves and 92,325 for S waves. We produced a model for P-wave velocity (Vp) and S-wave velocity (Vs) using our complete catalogue of absolute and differential times. Given that computing Vp/Vs directly from P- and Swave tomography models is only meaningful if both models are equally well constrained (e.g.

11 Eberhart-Phillips 1990; Zhang et al. 2009), we repeated the set of inversions for both Vp and Vs by using a data set that contains the same number of P and S phases for the same events. This subset is composed of 1,092 events, which corresponds to 24,555 absolute arrival times and 90,132 differential times for both P and S waves. A priori weights were used based on the uncertainty of the arrival readings (Fig. S1). In addition, for the first few set of iterations, 20 times greater weight was applied to the absolute data to establish the large-scale structures and for the final iterations, 20 times greater weight was applied for the differential data to refine the event locations and the velocity structure near the source regions. TomoFDD uses a finite-difference method for calculating travel times and ray paths, and it maps the desired regional portion of the spherical Earth into a uniform Cartesian grid (Zhang & Thurber 2006). The origin of our inverted region is 30.8ºS and 67.1ºW and grid nodes are uniformly spaced by 20 km horizontally and extend for 640 km in the east-west direction and 450 km in the north-south direction. With depth, the nodes are positioned every 10 km from 0 to 210 km. The map space for computing travel times is a 3-km uniform cartesian grid. As a starting model for the tomography inversion, we used our M1DM (Fig. 3a) and a Vp/Vs of 1.74. This Vp/Vs corresponds to the regional average calculated using the Wadati diagrams of 1,092 events with more than six observations of both P- and S-waves per station. In order to account for the known crustal thickness variations, we followed the approach of Lin et al. (2010), in which the inversion uses a reference model that incorporates a rough 3D-Moho interface. For doing that, we used Moho depths obtained by Gans et al. (2011), which progressively increases from 40 km in the east to 60 km in the west (Fig. 2b dashed line). To make the inversion more stable, the algorithm includes smoothing and damping parameters as regularization methods. Optimum values for these parameters were obtained via trade-off curves constructed by running a series of inversions with a large range of damping and smoothing values and plotting the data variance versus model variance (Fig. 4a). An analysis of the

12 “L-shape” curve defines the values that reduce most of the data variance without causing large increase in the solution variance (e.g. Eberhart-Phillips 1986; Zhang & Thurber 2007; Thurber et al. 2009). We chose a smoothing weight of 100 and damping of 150 to 80 (varying with iterations), which lead to condition numbers of 100-150. Since slowness models converge faster than earthquake locations (e.g. Thurber 1992), we performed eight sets of iterations alternating the simultaneous inversion with event relocation. The final RMS and variance for the data set after the 3D inversion were 0.138 s and 0.019 s, respectively, contrasting with the initial values of 0.526 s and 0.3057 s. After the 3D inversion, travel-time residuals also drastically reduced and sharply peaked symmetrically around zero (Fig. 4b).

6

RESOLUTION AND SENSITIVITY ANALYSIS

We utilized the Derivative Weighted Sum (DWS; e.g. Thurber & Eberhart-Phillips 1999) to qualitatively characterize the model resolution. DWS is proportional to the number of rays passing through a grid; hence, areas of good resolution will have high DWS values (Zhang & Thurber 2007) (Figs 5a and b). DWS values of 100 and 2500 are significant references: the first bounds the area where resolution is too low so these results are not presented due to lack of confidence, and the second delimits the area of high resolution where we focus our interpretations. We carried out two types of synthetic tests to examine the resolution of our study: checkerboard (Figs 5 and 6) and restoration (Fig. S3b) (e.g. Zhao et al. 1992; Haslinger 1998). For all these tests, inversions included the same number of P and S phases. For the checkerboard test, we calculated synthetic travel times through a checkerboard model. In this test, the source-receiver distribution matches the actual data set, but with Gaussian random noise added to the arrival times whose level is comparable to the picking uncertainty of the real data (0.06 s and 0.10 s, for P and S waves, respectively). We also added a constant noise term to the arrivals at each station, based on the final RMS of travel-time residuals (0.16 s and 0.13 for absolute and differential time data,

13 respectively). This constant noise term simulates systematic errors (model errors and pick bias) associated with the arrival times that are larger than the random ones (Zhang & Thurber 2003). Inversion of these travel times for Vp and Vs, using the same procedure as for the real data produces a checkerboard pattern that we compare to the original. Each checkerboard model consists of alternating high- and low-velocity anomalies with one grid node left open in between to test for horizontal smearing; at the same time, every other layer is free of anomalies to test for vertical smearing. We repeated the exercise shifting the checkerboard pattern for the initially free layers, to explore the smearing in all layers. We designed four sets of checkerboard tests (Figs S4-S13). The first one consisted of high- and low-velocity anomalies (±10%) in both Vp and Vs, which resulted in a Vp/Vs model without anomalies. The second and third sets consisted of either the Vp or Vs model without anomalies and either the Vs or Vp model with high- and low-velocity anomalies (±10%), which resulted in a Vp/Vs model with high and low anomalies (±10%). The fourth set consisted of high- and low-velocity anomalies of ± 5% with opposite sign between Vp and Vs, which resulted in a Vp/Vs model with high and low anomalies of ±10%. For the fourth set of tests, we present an example of a cross section at latitude 31.4ºS (Figs 5c and d) and depth slices at 20, 50, 70, and 90 km (Fig. 6). More examples of the aforementioned tests are in the Supporting Information (Figs S4-S13). Areas with DWS above 2,500 show the best recovery of the size and shape of the checkerboard pattern (Figs 5 and 6 black line), which in general correspond to the region between latitudes 32ºS and 30ºS and longitudes 69.5ºW and 67ºW at depths from 10 to 100 km (Figs 5 and 6). Amplitude recovery is especially good for the upper mantle at depths of 70-90 km (Figs 6, S10, S11, and S13) and somewhat underestimated between 0 and 60 km. Vertical smearing is low at all depths where DWS > 2,500 (e.g., Fig. S12), but it is present in some degree between 40 and 70 km, which corresponds to layers without seismicity. For the sets of checkerboard tests with Vp, Vs, or Vp/Vs models without anomalies, the recovered model does not “bleed over” between models (i.e. the recovered blank models do not

14 show checkerboard anomalies included in the initial models with anomalies) (Figs S4-S13). This was not the case when repeating this set of tests but using a higher number of P waves than S waves, giving a special relevance to the use of exactly the same amount of P- and S- wave data for deriving the Vp/Vs model, as we do for our results. For the restoration test (e.g. Zhao et al. 1992), we took the velocity model obtained by inverting the real data set (Fig. S3a) and used it as an initial model to calculate synthetic arrival times. We added the same random and constant noise described above for the checkerboard tests. By following the same inversion strategies as those for the real data, we obtained the restoring image (Fig. S3b). We found that the main anomalies and trends from the initial model are restored over well-resolved areas with DWS > 2500 (Fig. S3b). Following Zhang et al. (2004), we also tested the inter-event distances of less than 10 and 20 km for calculating the differential times (Figs 5a and b). Our results are insensitive to the choice of these cut-off distances for both P and S waves, and consequently both values can be used without impacting the well-resolved area. We present results for the distance cut-off of 20 km, which has a slightly higher DWS value. Additional tests explored the sensitivity of our results to the starting model, starting locations, parametrization, and the influence of random error in the data and are available in the Supporting Information (Fig. S3).

7

RESULTS

In this section, we present the tomographic results for Vp, Vs, and Vp/Vs using our final data set with a common set of P- and S-wave observations. Several cross sections and depth slices facilitate the description and interpretation (Figs 7, 8, 9, and 10). For each tomographic image, the velocity model is displayed for areas with DWS above 100 and interpreted for regions with DWS above 2,500, where we observed the best resolution on synthetic tests (Figs 5 and 6).

15 7.1 The Subducting Nazca Plate The WBZ marks the subducted slab of the Nazca plate (Figs 7 and 8). Along east-west and northeast-southwest cross sections, the WBZ is nearly flat and continuous at a depth of ~100 km between longitudes 70.3ºW and 67.5ºW (Figs 7 and 8). East of longitude 67.5ºW, the WBZ dips ~25º to the east, the slab seismicity is discontinuous, and earthquakes have increasing depths from 100 km to 200 km (Figs 7 and 8a, c, and e). Along north-south cross sections, the WBZ is flat and continuous at a depth of ~100 km between latitudes 32ºS and 30ºS beneath the Precordillera (Fig. 8b). At latitudes 31-30.5ºS beneath the Cuyania Terrane, the WBZ is also flat but discontinuous (Figs 8d and f). South of 31ºS, beneath Cuyania, the WBZ is inclined towards the south (Figs 8d and f). Our earthquake locations imply that the flat slab segment stretches 240 km in the north-south direction and 270 km east to west. Initial single-event locations (Fig. 2) show that, within the flat slab, there is a progressive shallowing of the WBZ towards the northeast opposite to the direction of subduction. Relocated hypocenters in the tomography using a 3D velocity model flatten out, aligning more closely at a depth of 100 km but still retaining a small westward dip of 2º (Figs 7, and 8a, c, and e). This slight shallowing of the seismicity was also documented in previous studies (e.g. Anderson et al. 2007; Deshayes 2008; Scarfi et al. 2012; Marot et al. 2013). For the subducting plate, we only have good tomographic resolution in the uppermost section where the earthquakes are located. The velocity results for the subducting Nazca Plate show a continuous flat band of mostly increased Vp (i.e. an increase of 0-4% with respect to the P-wave initial velocity model) at latitudes 32-31.2ºS and longitudes 70.3-67.5ºW (Figs 7b and 8h). In this this segment, the WBZ is located at the bottom of this increased Vp zone, and also marks the edge of the well-resolved area in the tomography. East of longitude 67.5ºW, the WBZ is inclined ~25º to the east together with an inclined band of increased (2-6%) Vp (Fig. 8h), however, in some parts this inclined band lacks resolution (Fig. 7b). East of 67ºW, the recorded seismicity is low and our resolution is low (DWS values < 2,500) for tracking the slab at deeper levels. In the northern part, at

16 latitudes 31-30.3ºS, the subducting Nazca Plate contains regions of reduced velocity (-2-0%) (Fig. 8g). The resulting absolute velocities at the top of the WBZ (depth ~100 km), where we have good resolution, tend to be 8.0-8.5 km/s for Vp, 4.5-4.8 km/s for Vs, and 1.73-1.76 for Vp/Vs (Figs 7 and 9f). The Vs obtained is in agreement with the ambient noise tomography of Porter et al. (2012), who found a Vs of 4.4-4.6 km/s immediately beneath the zone of the slab seismicity of Linkimer (2011). Above the WBZ, at depths of ∼100 km, there is a region of increased Vp with no seismicity (Figs 7b and 8h) and is mostly characterized by reduced (-2-0%) Vs (Fig. 7d). Within the area of good resolution, the velocities of this layer tend to be 7.9-8.3 km/s for Vp, 4.4-4.7 km/s for Vs, and 1.75-1.79 for Vp/Vs (Fig. 7). In this layer, Vs is higher in the eastern part of the flat slab (Fig. 9e) and slightly increases along the path of the JFR from 4.4 km/s at longitude 70.3ºW to 4.7 km/s at longitude 68.5ºW (Figs 7c and 9e).

7.2 The South American Mantle Lithosphere The mantle lithosphere of the South American Plate exhibits a great deal of heterogeneity with overall Vp of 7.7-8.3 km/s, Vs of 4.4-4.7 km/s, and Vp/Vs of 1.75-1.77 at depths from 60 to 90 km in the area where our model is best resolved. These velocities are more common at depths of 70 km (Figs 7, 8, and 9 labeled as D). Above the flat slab, we observe two distinct regions that stand out from this velocity range, which we refer to as anomalies A and B (Figs 7, 8, and 9). Anomaly A is characterized by Vp of 7.5-8.1 km/s, Vs of 4.4-4.8 km/s, and Vp/Vs of 1.701.73. In map view, anomaly A has an area of ∼50 x 50 km at depths of 70-90 km (Fig. 9a). Along the path of the JFR (Fig. 7), this anomaly extends from longitude 69.6ºW and 68.8ºW, at depths between 70 to 90 km, and has an apparent dip towards the east. The anomaly also has increased Vp and Vs with respect to the initial model (Figs 7b and 7d).

17 Previously, Wagner et al. (2005; 2008) found a low Vp/Vs (1.65-1.72) characterized by Vp of 7.8-8.2 km/s and Vs of 4.6-4.8 km/s covering an area of 300 x 120 km at depths of 70-100 km above the flat slab. This tomography covered a much larger area than the current study (36-30ºS and 72-67ºW) and had less resolution due to a much larger station spacing (CHARGE network) (Fig. 1a) and the use of a coarser grid spacing (40 x 40 x 20 km) with a smaller data set (~16,000 P- and ~5,000 S- wave arrivals). More recently, another tomography combining data sets from CHARGE and other shorter deployments before 2006 (Marot et al. 2014) suggested that the extent of the low Vp/Vs region could be smaller (~40 x 70 km). In this study, the checkerboard tests presented show that resolution was not ideal in the eastern part of the flat slab (longitudes 69.5-67ºW) causing the authors to discuss their results with caution. Marot et al. (2014) also used a coarser grid (40 x 40 x 10 km) than the current study but inverted a larger earthquake data set (55,128 P- and 54,889 Swave arrivals) than Wagner et al. (2005; 2008). The pattern of the low Vp/Vs anomalies presented by Marot et al. (2014) are similar to our results, but their values are higher (1.74 as compared to 1.70-1.73 from this study). When compared to Wagner et al. (2008), our results show a smaller low-Vp/Vs volume characterized by higher values (1.70-1.73). Since our study incorporates a differential time catalogue (170,241 P and 92,325 S-wave observations), a smaller grid spacing (20 x 20 x 10 km), and an improved starting model (M1DM) derived from our data set, we provide greater ray coverage and better resolution for the eastern segment of the Pampean flat slab than any previous regional tomography, making our observations more accurate. We believe this accounts for the above discrepancies. Anomaly B corresponds to Vp of 7.5-8.4 km/s, Vs of 4.3-4.6 km/s, and Vp/Vs of 1.78-1.82. This high Vp/Vs anomaly occurs at depths from 70 to 90 km. In map view, within the area of good resolution, this anomaly is visible as two separate regions at depths of 70 km (Fig. 9a) and as a single region at depths of 80 km, where it displays its greatest extent covering an area of 200 x 150 km (Fig. 9b). Along the path of the JFR, anomaly B appears to be above anomaly A, and is also inclined towards the east (Fig. 7e). A similar pattern exists along the east-west profiles between

18 latitudes 32ºS and 31.3ºS (Fig. 8c and e). In north-south profiles, anomaly B is especially clear as a 20-30 km-thick layer and is not present in the southern part of the flat slab (Figs 8b, d, and f). Marot et al. (2014) observed localized patches of higher Vp/Vs (1.79) at mantle depths below 70 km with a location similar to our observations.

7.3 The South American Crust From west to east, our tomography model shows that the middle and upper crust (0-30 km depth) have distinct velocities for the different terranes and morphotectonic units (Figs 7, 8, and 10). The eastern edge of the Chilenia Terrane exhibits the fastest observed velocities at upper crustal levels for Vp (7.0-7.9 km/s) and Vs (3.6-4.2 km/s), although we could only sample a limited area of this unit at the western edge of our well-resolved model. The Precordillera has overall lower Vp (mostly 5.8-6.4 km/s) and Vs (mostly 3.2-3.6 km/s) compared to those of the Cuyania Terrane (Vp 6.3-7.0 km/s and Vs 3.5-3.9 km/s) (Figs 7 and 10). In general, the Chilenia and the Precordillera are highlighted as regions of increased and reduced velocities, respectively, for both Vp and Vs (Figs 7b and d, 8h, and 10e and f). The Pampia Terrane where the model is well resolved also exhibits lower Vp (mostly 6.0-6.7 km/s) and higher Vs (mostly 3.4-4.1 km/s) compared to the Cuyania Terrane. The lower Vs in the middle-to-upper crust observed in the Precordillera compared to the Cuyania and Pampia terranes was also detected by Porter et al. (2012) and Ammirati et al. (2015). In general, the Vp/Vs at depths 0-30 km is slightly higher in the area sampled for the Chilenia (mostly >1.74) and Cuyania (mostly >1.75) terranes compared to the Pampia Terrane, which has the lowest Vp/Vs imaged at crustal levels (∼1.67-1.75) (Figs 7e and 10a) The Precordillera and the Cuyania Terrane seem to have a quite variable Vp/Vs, but the values in general are higher in the Cuyania Terrane (mostly >1.75). These results agree with Gilbert et al. (2006) and Alvarado et al. (2007) who showed higher Vp/Vs (∼1.80-1.85) for Cuyania compared to Pampia (< 1.75), although our estimates for Cuyania are not greater than 1.80. In the case of Chilenia and the Principal Cordillera, Marot et al. (2014) observed Vp/Vs of 1.78-1.79.

19 In the deeper levels of the crust (40-60 km) beneath the Precordillera and the Cuyania terranes, seismic velocities are mostly in the range of 7.6-8.2 km/s for Vp, 4.4-4.7 km/s for Vs, and 1.73-1.75 for the Vp/Vs (Figs 7 and 8, labeled as C). Our Vp/Vs is similar to that observed by Marot et al. (2014) (1.74-1.76), but lower than that the suggested by Alvarado et al. (2007) ( 1.80). The resolution of our model is insufficient to reliably analyze the Chilenia and Pampia terranes at lower crustal depths. Regarding the crustal seismicity relocated in the tomography, it is especially abundant in the Cuyania Terrane beneath the Pie de Palo mountain range and the Precordillera. The deepest crustal earthquakes beneath these regions occurred at ∼36-40 km (Figs 7 and 8).

8

INTERPRETATION AND DISCUSSION

In this section, we discuss our results in the context of the classic dehydration embrittlement and water cycle models that explain patterns of water recycling in the subduction factory and the presence of common mantle rocks and mineralogies (Figs 11 and 12; Table 1). In these models, seismicity might be associated with the hydration level of the incoming plate and the dehydration of the subducting plate, which might subsequently hydrate the mantle above the zone of water release (e.g. Kirby et al. 1996; Hacker et al. 2003b; Ranero et al. 2005). Even though the water cycle is widely accepted in subduction zones (e.g., Rüpke et al. 2004; Abers et al., 2013; Halpaap et al. 2019), there is still some debate in the direct relationship between seismicity and dehydration embrittlement because of the complicated nature of the mechanical behavior and reactions involved shown in laboratory experiments (e.g. Proctor & Hirth, 2015; Ferrand et al. 2017). For the incoming Nazca plate, in the outer-rise region, earthquakes seem to occur along faults that act as conduits for hydrating the slab before subduction (e.g. Kopp et al. 2004; Ranero et al. 2005). Offshore of Chile, to the west of our study area, plate-bending faulting is pervasive (Fig. 12b number 1) with fault strikes orientated roughly parallel to either the JFR or the trench (e.g. Kopp et al. 2004; Ranero et al. 2005). Near the ridge-trench intersection, at ∼33ºS latitude, outer-rise

20 seismicity shows evidence of conjugate normal faults extending to depths of 30 km within the oceanic upper mantle (Fromm et al. 2006). Additionally, reduced velocities in tomographic inversions in the outer-rise region suggest the hydration of the oceanic upper mantle (Kopp et al. 2004). Once the hydrated plate subducts, it dehydrates through reactions thought to produce embrittlement and intra-slab earthquakes (e.g. Kirby et al. 1996; Peacock 2001; Hacker et al. 2003b; Abers et al. 2013). Thermal modeling shows that seismicity within the subducting slab is distributed in a narrow belt of P and T conditions, which coincides with the breakdown of highpressure hydrous phases (e.g. Abers et al. 2006). Intra-slab earthquakes are particularly abundant in the subducting slab at latitudes 33-30ºS both in the dipping slab down to 120 km (Fig. 12b number 2) and in the flat slab segment at depths of 100-120 km (e.g. Anderson et al. 2007; Marot et al. 2013) (Fig. 12b number 6). A comparison with outer-rise tensional events leads to the hypothesis that these intermediate-depth earthquakes are related to the reactivation of outer-rise faults that were formed prior subduction (Marot et al. 2012; Ranero et al. 2015). Petrologic models predict the hydration of the mantle above the subducting slab (e.g. Hacker et al. 2003b). In a steeply-dipping subduction zone, water is released as the downgoing slab dehydrates, and it migrates into the hot overlying mantle where it could induce melting (e.g. Peacock 1993; Hacker et al. 2003b). In the case of a flat slab subduction zone, the released water would penetrate into and presumably hydrate the overlying cold lithospheric mantle, however, in this case temperatures are too low to trigger melting (Kay et al. 1988). Based on the convergence between the Nazca and South American plates (Norabuena et al. 1999) and the age of the cessation of volcanism across the Pampean flat slab region at ~9 Ma (Kay and Mpodozis, 2002; Ramos et al. 2002), Booker et al. (2004) estimated that the entire subducting plate at depths < 400 km has never been below active volcanism. Therefore, none of the water released by this plate segment has triggered upper-plate volcanism, and it potentially may still exist in the subsurface.

21 8.1 The Subducting Nazca Plate We interpret the band of increased Vp, where the WBZ is located, as the subducting Nazca Plate. In particular, the region of increased Vp and reduced Vs above the WBZ (Figs 7b and d) may correspond to the subducting oceanic crust, following the interpretation of Gans et al. (2011), who detected a moderately overthickened crust of ∼13-19 km above the WBZ. The general trend of a slight eastward increase in Vs at longitudes 70.3-68.5ºW within the subducting oceanic crust (Figs 7c and 9e) is indicative of a progressive dehydration as suggested by Porter et al. (2012) who showed that the subducting oceanic crust seems more hydrated in the west than in the east (Fig. 12b number 9). This pattern may connect to the slab seismicity. More abundant slab earthquakes in the west (longitudes 70.3-68ºW), may signify a greater volume of water released through the dehydration embrittlement relationship (e.g., Hacker et al. 2003b). In addition, the slab seismicity diminishes at longitudes 68-67ºW which may imply less dehydration below where Porter et al. (2012) and our results suggest the oceanic crust may be drier. Several studies using different geophysical techniques (Gans et al. 2011; Porter et al. 2012; Burd et al. 2013; Ammirati et al. 2013, 2015, 2016) determined the lithospheric structure and compared it to the hypocenters of slab earthquakes presented by Linkimer (2011). For example, Gans et al. (2011) found that these slab earthquakes appear to be occurring in the subducting oceanic mantle. Porter et al. (2012) supported this conclusion showing that the P and T conditions for the dehydration of serpentine and the typical lower plane of the WBZ coincides with the conditions of the oceanic mantle at 100 km depth in the Pampean flat slab region. Ammirati et al. (2015) found that the slab earthquakes were located below a low-Vs zone that they interpreted as the subducting oceanic crust. However, some of the 18 earthquakes they relocated through moment tensor inversion moved to within the low-Vs zone, which led them to suggest there may be seismicity within the subducting oceanic crust. To the west of our study area, beneath Chile, Marot et al. (2013) recognized a double seismic zone (DSZ) with the lower seismic plane (LP) located 20 km below the upper seismic plane (UP)

22 (Fig. 12b number 2). The LP and the UP merge at a depth of ∼100 km, which is approximately where the slab becomes horizontal (Fig. 12b number 3). Petrological models predict that the UP is located within the oceanic crust and the LP occurs within the oceanic uppermost mantle (e.g. Hacker et al. 2003b). Globally, the DSZ tends to occur in the shallower parts of the WBZ at depths from 40 and 180 km (Marot et al. 2013). Once the LP and UP merge at greater depths, the seismicity appears to occur mostly within the oceanic upper mantle (e.g. Hacker et al. 2003b). In our 3D velocity model, we do not observe slab earthquakes in the layer we correlate with the oceanic crust. In addition, we do not observe a DSZ where the slab is flat. Therefore, we infer that the seismicity within the subducting oceanic crust ceased in the shallower levels above where the LP and UP merged and we suggest that the slab seismicity in flat slab segment occurs in the oceanic uppermost mantle as presented by Gans et al. (2011) and Porter et al. (2012).

8.1.1 Slab geometry The Pampean flat slab has been traditionally linked to the subduction of the JFR, which may be buoyant because of its hydrated and moderately overthickened (~20-22 km) oceanic crust (Yáñez et al. 2002; Gans et al. 2011; Porter et al. 2012). However, if slab earthquakes illuminate the slab geometry, our earthquake locations imply that the lateral extent of the flat slab region is at least 240-km wide, which is broader than the width (~100 km) of JFR offshore (von Huene et al. 1997, Yáñez et al. 2001). In addition, our interpretation that the slab earthquakes occur in the oceanic uppermost mantle may imply that the oceanic crust is mostly dry and partially eclogitized making more difficult to explain the flat geometry. These aspects may indicate that more factors contribute to the flattening, besides the subduction of the JFR. Experiments suggest that a modest ridge perpendicular to the trench such as the present-day JFR is not buoyant enough to modify the slab geometry (Martinod et al. 2005). Therefore, the flat slab segment probably results from both the subduction of the buoyant anomaly and other key contributing factors such as the trenchward motion of the thick South American lithosphere and

23 suction forces in the mantle wedge as proposed by van Hunen et al. (2004), Manea et al. (2011) and Antonijevic et al. (2015). Another potentially influential factor is the migration of the ridge over time. The orientation of the JFR and the direction of plate convergence led to the southward migration (~3.5 cm/yr) of the ridge-trench intersection during the last 12 Ma (Yáñez et al. 2001). This migration indicates that the region north of the current inland-projected path of the ridge was also affected by its subduction and that the slab geometry may still be disturbed due to the ridge passage. Additionally, the amount of crustal thickening along hotspot tracks likely varies along the length of the track, with the potential implication that the volume of overthickened crust currently offshore does not reflect what has already been subducted. We interpret that depth variations in the WBZ reflect true variations in the slab geometry. For example, there is a westward dip of ∼2º still retained within our final tomography earthquake locations at longitudes 70.3-68ºW, indicating that the slab follows that geometry before it bends to the east at longitude ∼67.5ºW to resume its descent into the mantle with a steep angle (~25º) until depths of ~200 km (Fig. 12b number 10). To the east of the well-resolved tomography area, between longitudes 65ºW and 63ºW, a tear in the Nazca slab has been proposed (Fig. 12b number 11) based on earthquake locations (Cahill & Isacks 1992) and observations in the continental mantle lithosphere, which include a low-velocity anomaly in a teleseismic P-wave tomography (e.g. Portner et al. 2017) and a high-conductivity feature visible from magnetotelluric experiments (Booker et al. 2004; Burd et al. 2013). These observations were interpreted as evidence of entrained subslab material or asthenospheric material rising through a possible slab tear to the base of the lithosphere (Fig. 12b number 12). We cannot confirm the existence of these features, even with using seismic stations located above the possible tear, because the scarce slab seismicity captured at deeper levels (> 200 km) limited the resolution on this area.

24 8.2 The South American Mantle Lithosphere We investigated the seismic properties of common upper-mantle minerals (Fig. 11a) and rocks (Fig. 11b) at the appropriate P-T conditions, using Hacker et al. (2003a) and Abers & Hacker (2016), including lherzolite and harzburgite, the common enriched and depleted rock types of the upper mantle, respectively. Based on synthetic models of Marot et al. (2014), the P and T conditions for the upper mantle above the flat slab at depths of 60-100 km are 2-3 GPa and 600-700ºC, respectively. According to Gutscher et al. (2000b) and Gutscher et al. (2002) temperatures could be higher (600-900ºC), but still cooler than in a normal-dipping subduction zone, where temperatures may be 900-1200ºC above the subducting slab at depths 60-100 km (Marot et al. 2014). In general, the Vs and Vp/Vs above the flat slab tend to be 4.4-4.7 km/s and 1.75-1.77, respectively (Figs 11a and b region D). These values are consistent with either a chlorite lherzolite (ol66opx16chl15cpx3) or a chlorite harzburgite (ol73opx16chl11), which are both stable at 3 GPa and temperatures of 600-800ºC (Fig. 11b numbers 7 and 11) and have a low water content (< 1.9%) (Hacker et al. 2003a). The seismic properties of different types of lherzolites, either dry or hydrated, are outside the range of our seismic observations, therefore, we rule out the presence of an enriched (lherzolite) upper mantle (Fig. 11b numbers 1, 2, 3, and 4). Two velocity anomalies deviate from depleted lherzolite and harzburgite pointing to a heterogeneous composition for the upper mantle. We have referred to these anomalies as A and B. Anomaly A (Vs 4.4-4.8 km/s and Vp/Vs 1.70-1.73) is consistent with either orthopyroxenites or peridotites with a high content of enstatite, the Mg end-member of orthopyroxene (Figs 11a and b). The orthopyroxene enrichment needed to fit the elastic properties observed at 3 GPa, could be at least 70% (e.g. ol20opx70chl10) at temperatures of 500-700ºC (Fig. 11b number 13) and 80% (e.g. ol20opx80) at temperatures of 500-600ºC (Fig. 11b number 14). Comparing anomaly A with the elastic properties calculated more recently by Qian et al. (2018) at temperatures of 800-1450ºC and also at 600-800ºC (Zhongqing Wu, private communication, 2018), the orthopyroxene enrichment

25 could be lower. In this case, our results are consistent with peridotites with orthopyroxene content from 40 to 80% (Fig. 11b numbers 16, 17, and 18). According to Soustelle & Tommasi (2010), reactive percolation of Si rich fluids or hydrous melts in the upper mantle results in crystallization of orthopyroxene, which induces a decrease in Vp/Vs. Previously, Wagner et al. (2005; 2006; 2008) interpreted these seismic velocities as orthopyroxene enrichment. In that study, the Vp/Vs observed was low (1.65-1.72) (Figs 11a and b, labeled as W) and extending over almost the entire flat slab. Many authors (e.g. Soustelle & Tommasi 2010; Hacker & Abers 2012; Worthington et al. 2013) have not found such a low Vp/Vs (<1.70) in experimental measurements of samples with orthopyroxene enrichment and alternatively explain it as evidence of anisotropy in metasomatized peridotites sampled with vertically incident rays. In contrast, Qian et al. (2018) obtained a low Vp/Vs (<1.70) by calculating the elasticity of orthoenstatite at mantle conditions (2-3 GPa and 800-1450ºC) (Fig. 11b numbers 15-19) using first principle calculations based on density functional theory. According to these authors, the discrepancy may originate from the extrapolation to high P and T conditions in other studies (e.g. Hacker & Abers 2012). Qian et al. (2018) concluded that even a low Vp/Vs, like the values found by Wagner et al. (2008), could be explained as local enrichment of orthoenstatite in the mantle. The extent of the low Vp/Vs anomaly imaged with higher resolution in our study is smaller than that detected by Wagner et al. (2005; 2006; 2008). Recent geophysical studies (e.g. Porter et al. 2012; Marot et al. 2014; Ammirati et al. 2015) also interpreted a small region of cold, dry continental lithosphere in the upper mantle above the flat slab. Since our Vp/Vs values are slightly higher than 1.70, our interpretation does not require mantle anisotropy, and instead we link anomaly A to composition only. Thus, our tighter constraints on the velocity values make the interpretation of local orthopyroxene enrichment more plausible. Above the central part of the flat slab we also identify anomaly B, which is characterized by Vs of 4.3-4.6 km/s and Vp/Vs of 1.78-1.82 (Figs 7e and 8). These velocities are consistent with up

26 to 5% hydration for either a depleted lherzolite or harzburgite at temperatures of 500-700ºC (Fig. 11b numbers 6 and 10). These rocks could be a chlorite-serpentine wehrlite (ol57atg27chl13cpx3) (Fig. 11b number 6) or a serpentine-chlorite dunite (ol60atg30chl10) (Fig. 11b number 10). A larger amount of water for these rocks does not match our seismic observations because it would lower Vs (e.g. Fig. 11b numbers 5 and 9). Similarly, a hydrated lherzolite does not fit the observed velocities (Fig. 11b numbers 1, 2, and 3). Hydration of peridotite will stabilize a variety of hydrous minerals, in particular, antigorite (atg), which is the main serpentine mineral in ultramafic rocks and is stable to temperatures of 620720ºC and at depths 30-150 km (Ulmer & Trommsdorff 1995). At the Pampian flat slab region, temperatures at depths of 60-100 km are cooler (600-700ºC) than in a normal-dipping subduction zone (Marot et al. 2014), therefore antigorite may be stable. We rule out the possibility that the seismic signature of anomaly B corresponds to melt because temperatures are too low to trigger melting.

8.2.1 Water circulation in the upper mantle There is a consistent relationship between anomalies A and B with respect to the Vp/Vs detected for the subducting Nazca Plate. The general Vp/Vs pattern can be best noticed in map view at depths of 80 and 100 km (Fig. 9). At the depth of the WBZ (100 km), the subducting slab near 69.5ºW and 32ºS has higher Vp/Vs (1.76-1.80) and transitions to lower values in the northeast (< 1.76) (Fig. 9d). At a depth of 80 km, the trend seems to reverse: near 69.5ºW and 32ºS the upper mantle displays a low Vp/Vs (1.70-1.73) and transitions to be higher in the northeast (1.78-1.82) (Fig. 9b). Along the path of the JFR, this pattern in Vp/Vs is also visible in cross section (Fig. 7e). At 69.5ºW, the high Vp/Vs (1.76-1.80) exists in the slab at a depth of 100 km while the mantle above it exhibits a lower Vp/Vs (1.70-1.73). To the northeast, at 68.5ºW, the lower Vp/Vs exists in the slab while the mantle above it shows a higher Vp/Vs (Fig. 7e).

27 For the flat slab segment along east-west profiles, the WBZ in the southern region has a higher Vp/Vs (Figs 8c and e) compared to the northern part (Fig. 8a). In contrast, in the overlying continental mantle, the lower Vp/Vs (anomaly A) is only present in the southern region (Fig. 8e) while in the north there is a large extent of high Vp/Vs (Fig. 8a). These patterns support the hydration and dehydration processes across this part of the flat slab region. We interpret the lower and higher Vp/Vs as localized dry and hydrated areas, respectively. The subducting slab is more hydrated in the southwest (longitudes 70-69ºW) and above it, the upper mantle appears to be drier (Anomaly A). As water leaves the oceanic mantle of the slab, further to the northeast (longitudes 68.5-69ºW), the subducting slab seems to be drier, while the Vp/Vs in the overlying mantle increases (Anomaly B) as the mantle absorbs the water released from the slab. The largest extent of hydration also appears to occur in the northern part of the flat slab (Fig. 9b). Between longitudes 67ºW and 68.5ºW, where the slab resumes its descent into the mantle, both the upper mantle and subducting slab seem to be less hydrated. These seismic signatures agree with interpretations made by Porter et al. (2012) and Ammirati et al. (2015) of an eastward progressive dehydration of the slab during flat subduction. These authors observed an eastward increase in Vs in the subducting slab (Porter et al. 2012) and pockets of dry mantle observed mainly above the western part of the hydrated flat slab (Ammirati et al. 2015). We noticed that anomaly B is not present above the entire area of abundant WBZ seismicity (Figs 7 and 8). This suggests that there may be regions across the subducting slab where fluids may escape easily and that there may be pathways to circulate and concentrate fluids in certain areas of the continental mantle. Porter et al. (2012) suggested that a shift in stress orientation would likely affect the permeability of the subducting plate, impacting dehydration. In their interpretation, the slab bending associated with the transition from a normal angle of subduction to a flat geometry places the upper part of the slab under compression, closing fractures and decreasing the permeability of the top of

28 the slab, which would inhibit flow out of the slab. This idea implies that fluid flow paths are not vertical and that the places where fluid actually leaves the slab may differ significantly from the places where fluids were generated, as proposed by Hacker et al. (2003b) and van Keken et al. (2011). Another possibility is that much of the water, which is considered expelled from the slab, actually migrates to structurally deeper levels in the slab and it is subducted to much greater depths than expected, as also proposed by those authors. Our seismic velocity observations show for the first time that water is released from the slab in larger amounts near the center of the flat slab segment (longitudes 69-68ºW), where we observed the largest extent of anomaly B (Figs 7e, 9b). The orientation of the anomalies A and B, between 70ºW and 68.5ºW, at depths from 60 to 100 km, may suggest a large-scale structure that could relate to the larger water release near the center of the flat slab. Velocity anomalies have an inclination of ∼30º towards the east, from latitude 31.3ºS to 32ºS (Figs 8c and e). Both the Vp and Vs models show a sharp inclined zone bounding regions of increased and decreased velocities (Figs 7b and d and 8h green line). An inclined zone of reduced velocities, spanning the entire width of the continental lithosphere and dipping ∼30º to the east, is also a feature of a previous tomography (Marot et al. 2014) coinciding with the location of our observations, but extending further west to the Chilean forearc. A possible explanation for this velocity configuration is a lithospheric-scale suture, probably related to an ancient subduction zone. The tectonic evolution of the proto-Andean margin of Gondwana during the Paleozoic includes the collision of the Pampia (Early Cambrian), Cuyania and Precordillera (Middle-Late Ordovician), and Chilenia (Middle-Late Devonian) terranes to the Rio de la Plata craton. It also includes several subduction zones beneath Gondwana that led all the aforementioned terranes to accrete to its margin (e.g. Ramos et. al. 1986; Ramos 2004; Rapela et al. 1998; 2007). During the Silurian-Lower Devonian, a subducting slab was dipping beneath the Cuyania and Precordillera terranes, which were already part of Gondwana, and led the collision of Chilenia to Cuyania-Precordillera during the Middle-Late Devonian (Ramos et. al. 1986). Our

29 observations of inclined velocity anomalies across the lithosphere could correspond to the fossilized dipping interface of this relict subducting slab below Cuyania-Precordillera and therefore may correspond to the first documented tomographic evidence of the amalgamation of terranes to Gondwana. Similar to the outer-rise region where large-scale faults cut the entire oceanic crust and facilitate the hydration of the incoming plate (Ranero et al. 2005), we propose that this lithosphericscale suture may be a weak zone that allows slab fluids to migrate into the continental upper mantle and facilitates the upper-plate hydration (Fig. 12b number 8). Marot et al. (2014) explored a similar idea in the Chilean forearc region, where high Vp/Vs areas were linked to a ramp-detachment structure documented by Farías et al. (2010) at latitude 33.5ºS. In this case, the structure is dipping to the west from the upper crust at a depth of 10 km to the subduction interface at 60 km (Fig. 12b number 4). Marot et al. (2014) suggested this fault might channel slab fluids into the continental lower crust. Our finding of moderately hydrated regions is also consistent with previous studies on paleo-flat slabs in the western U.S. and beneath the Bolivian Altiplano that have suggested that flat slabs serve to hydrate the lithosphere above them (Humphreys et al. 2003; James & Sacks 1999). The current state of hydration in the Pampean flat slab may have implications for future volcanism. According to Booker et al. (2004) massive melting of hydrated lithosphere could explain the commonly observed flare-up in volcanism associated with the reestablishment of the asthenospheric wedge when the flat slab returns to a normal dip. In addition, our finding of regions where water appears to be released in larger amounts could be a target for future research of non-volcanic tremor (NVT), a subject that has not been explored yet for the Pampean flat slab region. In other areas, like the Mexican subduction zone, the NVT has been linked to slab dehydration processes in Oaxaca (Brudzinski et al., 2010) and above the flat slab of the Cocos Plate in Guerrero (e.g., Payero et al., 2008).

30 8.3 The South American Crust Our results reflect the differences in composition for the Sierras Pampeanas basement. The higher Vp/Vs for the Cuyania Terrane (mostly >1.75) agrees with a mafic-ultramafic composition for that region (e.g. Ramos 2004). Similarly, the small portion sampled for Chilenia (mostly >1.74) agrees with a belt of mafic-ultramafic outcrops near the suture zone between Cuyania and Chilenia (e.g. Pérez-Luján et al. 2015). The lower Vp/Vs for the Pampia Terrane (mostly 1.67-1.75) matches more felsic rocks (e.g. Rapela et al. 1998) and the lower Vp and Vs in the middle-to-upper crust of the Precordillera is consistent with the presence of Paleozoic carbonate rocks (e.g. Baldis et al. 1984). Seismic velocities calculated from different methodologies by Gilbert et al. (2006), Alvarado et al. (2007), Porter et al. (2012), Marot et al. (2014), and Ammirati et al. (2015; 2018) led to the same interpretation. In the lower crust (~40-60 km) beneath the Precordillera and the Cuyania Terrane seismic velocities (Vs 4.4-4.7 km/s and Vp/Vs 1.73-1.75) are consistent with the presence of eclogitized or partially eclogitized lower crust. The P and T conditions for the lower crust are 1-2 GPa and 400700ºC, respectively (Marot et al. 2004). Based on the petrologic models from Hacker et al. (2003a) and Abers & Hacker (2016), within these P-T conditions, our seismic velocities match three types of eclogites : amphibole, zoisite-amphibole, and zoisite (Fig. 11c), which have a low water content of 0.2-0.7% (Hacker et al. 2003a). This interpretation agrees with geophysical and petrological evidence presented in several studies (e.g. Calkins et al. 2006; Gilbert et al. 2006; Alvarado et al. 2007; Marot et al. 2014; Pérez-Luján et al. 2015; Ammirati et al. 2013, 2015, 2016, 2018; Venerdini et al. 2016). Our results from both the 1D and 3D inversions are also consistent with a thicker South American crust in the Western Sierras Pampeanas than in the Eastern Sierras Pampeanas. The regional crustal thickness variation is clearly marked by the P-wave station corrections of the M1DM. The station corrections show a strong pattern of positive values in the west and negative values in the east (Fig. 3b, Table S1). Since the Moho depth was set at 60 km for calculating

31 M1DM throughout the entire region, the station-correction pattern indicates that arrivals on the west side of the network are delayed with respect to those on the east, reflecting a large scale velocity variation due to the changes in crustal thickness. A similar pattern of a profound variation in station corrections with Moho depth, was also found by Pujol et al. (1991), Anderson et al. (2007), and Scarfi et al. (2012). Our 3D inversion also shows evidence of the crustal thickness variation. Vp contours in the 7.0-7.9-km/s range tend to get deeper in the Precordillera region, west of longitude 68.4ºW (Fig. 7a). For example, the 7.8 km/s contour appears at depths of 60 km beneath the Precordillera and shallows progressively to the east to depths of 40 km beneath Pampia. This trend also appeared when performing the tomography inversion using an initial flat Moho (Fig. S3). Other studies also imaged a thicker South American crust in the Western Sierras Pampeanas than in the Eastern Sierras Pampeanas using receiver functions (Calkins et al. 2006; Gilbert et al. 2006; Heit et al. 2008; Gans et al. 2011; Ammirati et al. 2013), Pn-phase velocity (Fromm et al. 2004), and S-to-P converted phases (Regnier et al. 1994).

8.3.1 Crustal structure The observed variation in crustal character, for both thickness and seismic velocities, corresponds to sutures as well as mid-crustal discontinuities, which are probably related to the terrane accretion history. The heterogeneity in the Vp/Vs and higher seismicity (Figs 7e, 8, and 10) within Cuyania and Precordillera agrees with greater complexity in crustal structure for these regions compared to other imaged regions, as was also shown by Gilbert et al. (2006), Alvarado et al. (2007), Gans et al. (2011), and Ammirati et al. (2013; 2015; 2016). There is a strong change in both Vp and Vs models in the western part of the study region (Figs 7b, d, 8a, h, and 10c, d, e, and f), which we correlate with the Chilenia-Precordillera (Cuyania) suture zone (Ramos et al. 1986; Pérez-Luján et al. 2015; Ammirati et al. 2016). Ammirati et al. (2016) also observed this strong change and interpreted it as a major west-dipping

32 Ordovician master fault affecting the entire crust (Fig. 12b number 5), allowing the Cuyania basement to slide under Chilenia during terrane accretion continental collision. This boundary spatially correlates with a belt of mafic-ultramafic rocks that crop out in the Central and Western Precordillera (Pérez-Luján et al. 2015; Ammirati et al. 2016). Our tomography model shows an east-dipping velocity feature crossing the lower crust, which also appeared in the tomography presented by Marot et al. (2014) and coincides with the location where Ammirati et al. (2016) identified a double Moho in receiver functions at longitudes 69.2-69ºW (Fig. 12 number 7). Ammirati et al. (2016) interpreted the double Moho as the displacement of the continental Moho controlled by an east-dipping shear zone that extends across the lower crust connecting to the aforementioned west-dipping master fault (Fig. 12 numbers 5 and 7). The east-dipping crustal structure may be the shallow continuation of the suture we proposed in the previous section for the upper mantle (Fig. 12 number 8). If this is a suture zone related to an ancient subduction zone, as interpreted above, it may represent a lithospheric-scale structure extending from the lower crust and crossing the continental Moho down to the flat slab (Figs 7b and d and 10h green line). This structure may promote the hydration of the upper mantle and may serve as a way to transport fluids up through the base of the South American crust. Ammirati et al. (2013; 2015; 2016; 2018) identified several mid-crustal discontinuities for the Pampean region. Beneath the Precordillera, Ammirati et al. (2013) mapped two mid-crustal discontinuities at depths of 21 and 36 km, beneath station DOCA (Figs 8a and g). Our results show sharp velocity contrasts at depths of 20 and 40 km beneath that station, which may represent those discontinuities. From 0 to 20 km, the upper crust exhibits increased velocities and from 40 to 60 km it is reduced (Figs 8g and h). The resulting Vp/Vs is high (1.77-1.81) for almost the entire crust from 0 to 40 km, which agrees with mafic and ultramafic compositions for this region located near the Chilenia-Precordillera suture zone (Fig. 8a) (e.g. Pérez-Luján et al. 2015).

33 9

CONCLUSIONS

Our regional-scale tomography reveals compositional and structural aspects of the subducting Nazca Plate and the lithosphere of the overriding South American Plate using the densest and longest temporary seismic array ever used for a tomography study in the Pampean flat slab region of Argentina. In our detailed 3D velocity model, the subducting Nazca Plate is imaged as a mostly continuous band of increased (0-6%) Vp. Above the WBZ a region of increased Vp and slightly reduced Vs is associated with the subducting oceanic crust with no related seismicity. Our relocated slab earthquakes illuminate the slab geometry in great detail highlighting a flat segment at a depth of 100 km (longitudes 70-67.5ºW), where it is at least 240-km wide and 270-km long, and has a slight westward dip of 2º. Further east, the slab resumes its descent into the mantle with a steep angle (~25º) and it cannot be traced deeper than 200 km east of longitude 65ºW. We found that the flat slab segment is wider than the width of the JFR offshore (~100 km). This implies that there might be additional contributing factors for the flattening besides the subduction of the moderately overthickened oceanic crust of the ridge. Our results show that the upper mantle of the South American Plate above the flat slab is laterally heterogeneous, but it is mostly characterized by Vp/Vs of 1.75-1.77, which is consistent with either a chlorite lherzolite or a chlorite harzburgite. We found two seismic anomalies with lower and higher Vp/Vs, which we interpret to be localized dry and hydrated areas, respectively. The area of lower Vp/Vs (1.70-1.73) is smaller than that imaged in previous studies, and may correspond to either orthopyroxenites or peridotites with a 40-80% enrichment of enstatite. Our tighter constraints on these velocity values compared to previous tomographies, make the interpretation of local orthopyroxene enrichment more plausible. On the other hand, the higher Vp/Vs (1.78-1.82) is consistent with up to 5% hydration of mantle peridotite and may be a chloriteserpentine wehrlite or a serpentine-chlorite dunite. The size, orientation, and location of the seismic anomalies found offer insights about the hydration and dehydration processes across the flat slab region. We found seismic velocity trends

34 that indicate the progressive eastward dehydration of the slab and the hydration of the overriding mantle lithosphere. Our velocity observations also suggest that the slab releases water in greater amounts near the center of the flat slab (longitudes 69-68ºW) and that there may be regions across the subducting slab where fluids escape easily. The velocity anomalies at longitudes 70-68.5ºW displays for the first time in a tomography, a clear dip of ∼30º towards the east across the entire lithosphere from depths of 40 to 100 km. We interpret this velocity configuration as a large lithospheric suture perhaps associated with an ancient subduction zone prior to the collision between the Chilenia and Cuyania terranes. We propose that this structure might facilitate the migration of water from the subducting slab to the base of the South American crust, and connecting to shallower structures, into the crust itself. Our results show lateral differences in composition and thickness for the Sierras Pampeanas basement. The middle-to-upper crust in the Precordillera has a lower Vp and Vs compared to the Cuyania and Pampia terranes. In addition, the Vp/Vs is higher in the Chilenia and Cuyania terranes (mostly >1.75) compared to the Pampia Terrane (mostly 1.67-1.75). These observations are consistent with a Paleozoic carbonate rocks in the Precordillera, a more mafic-ultramafic composition for portions of the Chilenia and Cuyania terranes and a more felsic composition for the Pampia Terrane. The greater variability in Vp/Vs and higher seismicity found within the Cuyania Terrane and the Precordillera also agrees with greater complexity in crustal structure for the western basement regions. Our results from both the 1D and 3D inversions are consistent with a thicker South American crust in the Western Sierras Pampeanas than in the Eastern Sierras Pampeanas. P-wave station corrections in the M1DM clearly mark this thickness variation with positive values in the west and negative values in the east. The 7.8-km/s- Vp contour in our 3D model, which appears at depths of 60 km beneath the Precordillera and shallows progressively to the east to depths of 40 km beneath Pampia also reflects this thickness variation.

35 The change in crustal character for both thickness and seismic velocities corresponds to sutures between terranes as well as mid-crustal discontinuities, which are probably related to the terrane accretion history. We found a strong velocity change near the contact between Chilenia and Precordillera (Cuyania) in both Vp and Vs models, which we correlate with the suture between those terranes. Additionally, mid-crustal discontinuities exist beneath the Precordillera at depths of 20 and 40 km. Our seismic observations also confirm the interpretation of a partially eclogitized lower crust (~40-60 km) below the Cuyania terrane and Precordillera.

ACKNOWLEDGMENTS. This research was supported by the National Science Foundation (EAR-0510966, EAR-0738935, EAR-0739001, and EAR-1415914). Seismic instruments were provided by IRIS through the PASSCAL Instrument Center (NSF Cooperative Agreement EAR-0552316). We sincerely thank the Instituto Nacional de Prevención Sísmica (INPRES) and Mauro Sáez in Argentina as well as Noel Barstow (PASSCAL) for their invaluable assistance in the field for the SIEMBRA deployment. Special thanks to C. Berk Biryol, Ryan Porter, Christine Gans, Ivonne G. Arroyo, and Bradley Hacker for helping with codes, discussions, and reviewing the final manuscript and Dr. Zhongqing Wu for providing the elasticity properties in the 600-800ºC range for orthoenstatite. Maps were created using the Generic Mapping Tools (GMT) software (Wessel & Smith 1998).

DATA AVAILABILITY The SIEMBRA (Beck & Zandt 2007) and ESP (Gilbert 2008) data sets are available via IRIS DMC under ZL (https://doi.org/10.7914/SN/ZL_2007) and XH (https://doi.org/10.7914/SN/XH_2008) FDSN

codes.

The

arrival-time

catalogue

is

available

upon

request

via

email

to

[email protected]. Figures were made using the Generic Mapping Tools version 4.5.13 (www.soest.hawaii.edu/gmt; Wessel &Smith, 1998).

36 DECLARATION OF COMPETING INTEREST There are no conflicts of interest to declare.

REFERENCES Abers, G. A., Nakajima, J., van Keken, P. E., Kita, S. & Hacker, B. R., 2013. Thermal–petrological controls on the location of earthquakes within subducting plates. Earth Planet. Sci. Lett., 369, 178-187. Abers, G. A., & Hacker B. R., 2016. A MATLAB toolbox and Excel workbook for calculating the densities, seismic wave speeds, and major element composition of minerals and rocks at pressure and temperature, Geochem. Geophys. Geosyst., 17, 616–624. Abers, G.A., Keken, P.E. van, Kneller, E.A., Ferris, A. & Stachnik, J.C., 2006. The thermal structure of subduction zones constrained by seismic imaging: Implications for slab dehydration and wedge flow, Earth Planet. Sci. Lett., 241, 387–397. Allmendinger, R.W., Figueroa, D., Snyder, D., Beer, J., Mpodozis, C. & Isacks, B.L., 1990. Foreland shortening and crustal balancing in the Andes at 30°S latitude, Tectonics, 9, 789– 809. Alvarado, P. & Beck, S., 2006. Source characterization of the San Juan (Argentina) crustal earthquakes of 15 January 1944 (Mw 7.0) and 11 June 1952 (Mw 6.8), Earth Planet. Sci. Lett., 243, 615–631. Alvarado, P. & Ramos, V., 2011. Earthquake deformation in the northwestern Sierras Pampeanas of Argentina based on seismic waveform modelling., J. Geodyn., 51, 205–218. Alvarado, P., Beck, S., Zandt, G., Araujo, M. & Triep, E., 2005. Crustal deformation in the southcentral Andes backarc terranes as viewed from regional broad-band seismic waveform modelling, Geophys. J. Int., 163, 580–598. Alvarado, P., Beck, S. & Zandt, G., 2007. Crustal structure of the south-central Andes Cordillera and backarc region from regional waveform modelling, Geophys. J. Int., 170, 858–875.

37 Alvarado, P., Pardo, M., Gilbert, H., Miranda, S., Anderson, M., Saez, M. & Beck, S., 2009. Flatslab subduction and crustal models for the seismically active Sierras Pampeanas region of Argentina, Geol. Soc. Am. Mem., 204, 261–278. Ammirati, J.-B., Alvarado, P., Perarnau, M., Saez, M. & Monsalvo, G., 2013. Crustal structure of the Central Precordillera of San Juan, Argentina (31°S) using teleseismic receiver functions, J. South Am. Earth Sci., 46, 100–109. Ammirati, J.-B., Alvarado, P. & Beck, S., 2015. A lithospheric velocity model for the flat slab region of Argentina from joint inversion of Rayleigh wave phase velocity dispersion and teleseismic receiver functions, Geophys. J. Int., 202, 224–241. Ammirati, J.-B., Pérez Luján, S., Alvarado, P., Beck, S., Rocher, S. & Zandt, G., 2016. Highresolution images above the Pampean flat slab of Argentina (31–32°S) from local receiver functions: Implications on regional tectonics, Earth Planet. Sci. Lett., 450, 29–39. Ammirati, J.-B., Venerdini, A., Alcacer, J.M., Alvarado, P., Miranda, S. & Gilbert, H., 2018. New insights on regional tectonics and basement composition beneath the eastern Sierras Pampeanas (Argentine back-arc region) from seismological and gravity data, Tectonophysics, 740–741, 42–52. Anderson, M., Alvarado, P., Zandt, G. & Beck, S., 2007. Geometry and brittle deformation of the subducting Nazca Plate, Central Chile and Argentina, Geophys. J. Int., 171, 419–434. Antonijevic, S.K., Wagner, L.S., Kumar, A., Beck, S.L., Long, M.D., Zandt, G., Tavera, H. & Condori, C., 2015. The role of ridges in the formation and longevity of flat slabs, Nature, 524, 212–215. Astini, R.A., Benedetto, J.L. & Vaccari, N.E., 1995. The early Paleozoic evolution of the Argentine Precordillera as a Laurentian rifted, drifted, and collided terrane: a geodynamic model, Geol. Soc. Am. Bull., 107, 253–273. Baldis, B.A., Beresi, M., Bordonaro, O. & Vaca, A., 1984. The Argentine Precordillera as a key to Andean structure, Episodes, 7, 14–19.

38 Baldo, E.G.A., Demange, M. & Martino, R.D., 1996. Evolution of the Sierras de Cordoba, Argentina, Tectonophysics, 267, 1–4. Barazangi, M. & Isacks, B.L., 1976. Spatial distribution of earthquakes and subduction of the Nazca plate beneath South America, Geology, 4(11), 686–692. Beck, S. & Zandt, G., 2007. Lithospheric Structure and Deformation of the Flat Slab Region of Argentina. International Federation of Digital Seismograph Networks. Other/Seismic Network. DOI: 10.7914/SN/ZL_2007. Bezacier, L., Reynard, B., Bass, J.D., Sanchez-Valle, C. & Van de Moortèle, B., 2010. Elasticity of antigorite, seismic detection of serpentinites, and anisotropy in subduction zones, Earth Planet. Sci. Lett., 289, 198–208. Bishop, B.T., Beck, S.L., Zandt, G., Wagner, L., Long, M., Antonijevic, S.K., Kumar, A. & Tavera, H., 2017. Causes and consequences of flat-slab subduction in southern Peru, Geosphere, 13, 1392–1407. Booker, J.R., Favetto, A. & Pomposiello, M.C., 2004. Low electrical resistivity associated with plunging of the Nazca flat slab beneath Argentina, Nature, 429, 399–403. Boutelier, D.A. & Cruden, A.R., 2008. Impact of regional mantle flow on subducting plate geometry and interplate stress: insights from physical modelling, Geophys. J. Int., 174, 719– 732. Brudzinski, M. R., Hinojosa‐Prieto, H. R., Schlanser, K. M., Cabral‐Cano, E., Arciniega‐Ceballos, A., Diaz‐Molina, O., and DeMets, C. (2010), Nonvolcanic tremor along the Oaxaca segment of the Middle America subduction zone, J. Geophys. Res., 115, B00A23 Burd, A.I., Booker, J.R. Mackie, R. Pomposiello, M.C. & Favetto, A., 2013. Electrical conductivity of the Pampean shallow subduction region of Argentina near 33S: Evidence for a slab window, Geochem. Geophys. Geosyst., 14, 3192–3209. Cahill, T. & Isacks, B.L., 1992. Seismicity and shape of the subducted Nazca Plate, J. Geophys. Res., 97, 17503–17529.

39 Calkins, J.A., Zandt, G., Gilbert, H.J. & Beck, S.L., 2006. Crustal images from San Juan, Argentina, obtained using high frequency local event receiver functions, Geophys. Res. Lett., 33, 1–4. Carlson, R.L. & Miller, D.J., 2003. Mantle wedge water contents estimated from seismic velocities in partially serpentinized peridotites, Geophys. Res. Lett., 30, 1250. Christensen, N., 2004. Serpentinites, Peridotites, and Seismology, Int. Geol. Rev., 46, 795–816. Dahlquist, J.A. & Galindo, C., 2004. Geoquímica isotópica de los granitoides de la sierra de Chepes: un modelo geotectónico y termal, implicancias para el orógeno famatiniano, Rev. Asoc. Geológica Argentina, 59, 57–69. Deshayes, P., 2008. Tomographie en vitesse et en atténuation de la zone de subduction au Chili central-ouest de l’Argentine (29S–34S) à partir de données sismologiques locales: apport à l’étude de la composition minéralogique, PhD dissertation, Université de Nice, SophiaAntipolis, France. Eberhart-Phillips, D., 1986. Three-dimensional velocity structure in Northern California Coast Ranges from inversion of local earthquake arrival times, Bull. Seismol. Soc. Am., 76, 1025– 1052. Eberhart-Phillips, D., 1990. Three‐dimensional P and S velocity structure in the Coalinga Region, California, J. Geophys. Res. Solid Earth, 95, 15343–15363. Espurt, N., Funiciello, F., Martinod, J., Guillaume, B., Regard, V., Faccenna, C. & Brusset, S., 2008. Flat subduction dynamics and deformation of the South American plate: Insights from analog modeling, Tectonics, 27, TC3011. Ferrand, T.P., Hilairet, N., Incel, S., Deldicque, D., Labrousse, L., Gasc, J., Renner, J., Wang, Y., Green, H. W., Schubnel, A., 2017. Dehydration-driven stress transfer triggers intermediatedepth earthquakes, Nat. Commun., 8, 15247. Fromm, R., Zandt, G. & Beck, S.L., 2004. Crustal thickness beneath the Andes and Sierras Pampeanas at 30°S inferred from Pn apparent phase velocities, Geophys. Res. Lett., 31, L06625.

40 Fromm, R., Alvarado, P., Beck, S.L. & Zandt, G., 2006. The April 9, 2001 Juan Fernández Ridge (Mw 6.7) Tensional Outer-Rise Earthquake and its Aftershock Sequence, J. Seismol., 10, 163– 170. Gans, C.R., Beck, S.L., Zandt, G., Gilbert, H., Alvarado, P., Anderson, M. & Linkimer, L., 2011. Continental and oceanic crustal structure of the Pampean flat slab region, western Argentina, using receiver function analysis: new high-resolution results, Geophys. J. Int., 186, 45–58. Gilbert, H., 2008. Lithospheric Structure above the variably dipping Nazca Slab. International Federation

of

Digital

Seismograph

Networks.

Other/Seismic

Network.

DOI:

10.7914/SN/XH_2008 Gilbert, H., Beck, S. & Zandt, G., 2006. Lithospheric and upper mantle structure of central Chile and Argentina, Geophys. J. Int., 165, 383–398. González-Bonorino, F., 1950. Algunos problemas geológicos de las Sierras Pampeanas. Rev. Asoc. Geológica Argentina, 5, 81–110. Gordillo, C.E. & Linares, A., 1982. Geocronología y petrografía de las vulcanitas terciarias del Departamento de Pocho, Provincia de Córdoba, Rev. Asoc. Geológica Argentina, 36, 380–388. Gutscher, M.-A., 2002. Andean subduction styles and their effect on thermal structure and interplate coupling, J. South Am. Earth Sci., 15, 3–10. Gutscher, M.-A., Maury, R., Eissen, J.-P. & Bourdon, E., 2000a. Can slab melting be caused by flat subduction?, Geology, 28(6), 535–538. Gutscher, M.-A., Spakman, W., Bijwaard, H. & Engdahl, E.R., 2000b. Geodynamics of flat subduction: Seismicity and tomographic constraints from the Andean margin, Tectonics, 19, 814–833. Hacker, B. R., & Abers, G. A., 2012. Subduction Factory 5: Unusually low Poisson’s ratios in subduction zones from elastic anisotropy of peridotite, J. Geophys. Res., 117, B06308. Hacker, B.R., Abers, G.A. & Peacock, S.M., 2003a. Subduction factory 1. Theoretical mineralogy, densities, seismic wave speeds, and H2O contents, J. Geophys. Res. Solid Earth, 108, 2029.

41 Hacker, B.R., Peacock, S.M., Abers, G.A. & Holloway, S.D., 2003b. Subduction factory 2. Are intermediate-depth earthquakes in subducting slabs linked to metamorphic dehydration reactions?, J. Geophys. Res. Solid Earth, 108, 2030. Halpaap, F., Rondenay, S., Perrin, A., Goes, S., Ottemöller, L., Austrheim, H., Shaw, R., Eeken, T., 2019. Earthquakes track subduction fluids from slab source to mantle wedge sink, Sci. Adv., 5, eaav7369. Haslinger, F., 1998. Velocity structure, seismicity and seismotectonics of Northwestern Greece between the Gulf of Arta and Zakynthos, PhD dissertation, Swiss Federal Institute of Technology, Zuerich, Switzerland. Havskov, J. & Ottemoller, L., 1999. SeisAn Earthquake Analysis Software, Seismol. Res. Lett., 70, 532–534. Heit, B., Yuan, X., Bianchi, M., Sodoudi, F. & Kind, R., 2008. Crustal thickness estimation beneath the southern central Andes at 30°S and 36°S from S wave receiver function analysis, Geophys. J. Int., 174, 249–254. Henderson, L.J., Gordon, R.G. & Engebretson, D.C., 1984. Mesozoic aseismic ridges on the Farallon Plate and southward migration of shallow subduction during the Laramide Orogeny, Tectonics, 3, 121–132. Hu, J., Liu, L., Hermosillo, A. & Zhou, Q., 2016. Simulation of late Cenozoic South American flatslab subduction using geodynamic models with data assimilation., Earth Planet. Sci. Lett., 438, 1–13. Huangfu, P., Wang, Y., Cawood, P.A., Li, Z.H., Fan, W. & Gerya, T. V., 2016. Thermo-mechanical controls of flat subduction: Insights from numerical modeling, Gondwana Res., 40, 170–183. Husen, S., Kissling, E., Flueh, E. & Asch, G., 1999. Accurate hypocentre determination in the seismogenic zone of the subducting Nazca Plate in northern Chile using a combined on/offshore network, Geophys. J. Int., 138, 687–701.

42 Humphreys, E., Hessler, E., Dueker, K., Farmer, G.L., Erslev, E., & Atwater, T., 2003. How Laramide-age hydration of North American lithosphere by the Farralon slab controlled subsequent activity in the western United States, International Geology Review, 45, 575-595. Hyndman, R.D. & Peacock, S.M., 2003. Serpentinization of the forearc mantle, Earth Planet. Sci. Lett., 212, 417–432. James, D. E., & Sacks, I. S., 1999. Cenozoic formation of the Central Andes: A Geophysical Perspective, in Skinner, B. J., ed., Geology and ore deposits of the Central Andes: Society of Economic Geologists Special Publication, 7, p. 1–25. Jordan, T.E. & Allmendinger, R.W., 1986. The Sierras Pampeanas of Argentina; a modern analogue of Rocky Mountain foreland deformation, Am. J. Sci., 286, 737–764. Kay, S.M. & Abbruzzi, J.M., 1996. Magmatic evidence for Neogene lithospheric evolution of the central Andean “flat-slab” between 30°S and 32°S, Tectonophysics, 259, 15–28. Kay, S.M. & Gordillo, C.E., 1994. Pocho volcanic rocks and the melting of depleted continental lithosphere above a shallowly dipping subduction zone in the central Andes, Contrib. to Mineral. Petrol., 117, 25–44. Kay, S.M. & Mpodozis, C,. 2002. Magmatism as a probe to the Neogene shallowing of the Nazca plate beneath the modern Chilean flat-slab., J. South Am. Earth Sci., 15, 39–57. Kay, S.M., Maksaev, V., Moscoso, R., Mpodozis, C., Nasi, C. & Gordillo, C.E., 1988. Tertiary andean magmatism in Chile and Argentina between 28°S and 33°S: Correlation of magmatic chemistry with a changing Benioff zone, J. South Am. Earth Sci., 1, 21–38. Kendrick, E., Bevis, M., Smalley, R., Brooks, B., Vargas, R.B., Laurıá , E. & Fortes, L.P.S., 2003. The Nazca–South America Euler vector and its rate of change, J. South Am. Earth Sci., 16, 125–131. Kirby, S., Engdahl, R.E. & Denlinger, R., 1996. Intermediate-depth intraslab earthquakes and arc volcanism as physical expressions of crustal and uppermost mantle metamorphism in

43 subducting slabs, in Subduction Top to Bottom, Geophysical Monograph Series, 96, pp. 195– 214, American Geophysical Union. Kissling, E., 1988. Geotomography with local earthquake data, Rev. Geophys., 26, 659–698. Kissling, E., Ellsworth, W.L., Eberhart-Phillips, D. & Kradolfer, U., 1994. Initial reference models in local earthquake tomography, J. Geophys. Res. Solid Earth, 99, 19635–19646. Kissling E. Kradolfer U. & Maurer H., 1995. Program VELEST User's Guide-short Introduction, Internal Report, Institute of Geophysics, Swiss Federal Institute of Technology, Zuerich, Switzerland, pp 26. Kopp, H., Flueh, E.R., Papenberg, C. & Klaeschen, D., 2004. Seismic investigations of the O’Higgins Seamount Group and Juan Fernández Ridge: Aseismic ridge emplacement and lithosphere hydration, Tectonics, 23, TC2009. Leal, P.R., Hartmann, L.A., Santos, J.O.S., Ramos, V.A. & Miró, R.C., 2003. Volcanismo postorogénico en el extremo norte de las Sierras Pampeanas Orientales: Nuevos datos geocronológicos y sus implicancias tectónicas, Rev. Asoc. Geológica Argentina, 58, 593–607. Lienert, B.R. & Havskov, J., 1995. A Computer Program for Locating Earthquakes Both Locally and Globally, Seismol. Res. Lett., 66, 26–36. Lin, G., Thurber, C.H., Zhang, H., Hauksson, E., Shearer, P.M., Waldhauser, F., Brocher, T.M. & Hardebeck, J.,, 2010. A California Statewide Three-Dimensional Seismic Velocity Model from Both Absolute and Differential Times, Bull. Seismol. Soc. Am., 100, 225–240. Linkimer, L., 2011. Lithospheric Structure of the Pampean Flat Slab (Latitude 30–33S) and Northern Costa Rica (Latitude 9–11N) Subduction Zones, PhD dissertation, University of Arizona, Tucson, AZ, USA. Manea, V.C., Manea, M., Ferrari, L., Orozco, T. & Valenzuela, R.W., 2017. A review of the geodynamic evolution of flat slab subduction in Mexico, Peru, and Chile, Tectonophysics, 695, 27–52.

44 Manea, V. C., Pérez-Gussinyé, M. & Manea, M., 2011. Chilean flat slab subduction controlled by overriding plate thickness and trench rollback. Geology, 40(1), 35–38. Marot, M., Monfret, T., Pardo, M., Ranalli, G. & Nolet, G., 2012. An intermediate-depth tensional earthquake (MW 5.7) and its aftershocks within the Nazca slab, central Chile: A reactivated outer rise fault?, Earth Planet. Sci. Lett., 327–328, 9–16. Marot, M., Monfret, T., Pardo, M., Ranalli, G. & Nolet, G., 2013. A double seismic zone in the subducting Juan Fernandez Ridge of the Nazca Plate (32°S), central Chile, J. Geophys. Res. Solid Earth, 118, 3462–3475. Marot, M., Monfret, T., Gerbault, M., Nolet, G., Ranalli, G. & Pardo, M., 2014. Flat versus normal subduction zones: a comparison based on 3-D regional traveltime tomography and petrological modelling of central Chile and western Argentina (29o–35oS), Geophys. J. Int., 199, 1633–1654. Martinod, J., Funiciello, F., Faccenna, C., Labanieh, S. & Regard, V., 2005. Dynamical effects of subducting ridges: insights from 3-D laboratory models, Geophys. J. Int., 163, 1137–1150. Martinod, J., Guillaume, B., Espurt, N., Faccenna, C., Funiciello F., & Regard, V., 2013. Effect of aseismic ridge subduction on slab geometry and overriding plate deformation: Insights from analogue modeling. Tectonophysics, 588, 39–55. McGeary, S., Nur, A. & Ben-Avraham, Z., 1985. Spatial gaps in arc volcanism: The effect of collision or subduction of oceanic plateaus, Tectonophysics, 119, 195–221. Norabuena, E.O., Dixon, T.H., Stein, S. & Harrison, C.G.A., 1999. Decelerating Nazca-South America and Nazca-Pacific Plate motions, Geophys. Res. Lett., 26, 3405–3408. Pardo, M., Monfret, T., Vera, E., Eisenberg, A., Gaffet, S., Lorca, E. & Perez, A., 2002. Flat-slab subduction zone in central Chile–Argentina: seismotectonic and body-wave tomography from local data. 5th International Symposium on Andean Geodynamics (ISAG 2002), Toulousse, France, Extended Abstracts, Institut de recherche pour le développement (IRD), pp. 469–472. Payero, J. S., Kostoglodov, V., Shapiro, N., Mikumo, T., Iglesias, A., Pérez‐Campos, X., and Clayton, R.

45 W. (2008), Nonvolcanic tremor observed in the Mexican subduction zone, Geophys. Res. Lett., 35, L07305. Peacock, S.M., 1993. Large-scale hydration of the lithosphere above subducting slabs, Chem. Geol, 108, 49–59. Peacock, S.M., 2001. Are the lower planes of double seismic zones caused by serpentine dehydration in subducting oceanic mantle?, Geology, 29(4), 299–302. Perarnau, M., Gilbert, H., Alvarado, P., Martino, R. & Anderson, M., 2012. Crustal structure of the Eastern Sierras Pampeanas of Argentina using high frequency receiver functions, Tectonophysics, Vol. 580, p. 208-217, doi: 10.1016/j.tecto.2012.09.02.1 Pérez Luján, S.B., Ammirati, J.B., Alvarado, P. & Vujovich, G.I., 2015. Constraining a mafic thick crust model in the Andean Precordillera of the Pampean flat slab subduction region, J. South Am. Earth Sci., 64, 325–338. Pilger, R.H., 1981. Plate reconstructions, aseismic ridges, and low-angle subduction beneath the Andes, Geol. Soc. Am. Bull., 92, 448–456. Porter, R., Gilbert, H., Zandt, G., Beck, S., Warren, L., Calkins, J., Alvarado, P. & Anderson, M., 2012. Shear wave velocities in the Pampean flat-slab region from Rayleigh wave tomography: Implications for slab and upper mantle hydration, J. Geophys. Res. Solid Earth, 117, B11301. Portner, D. E., Beck, S., Zandt, G. & Scire, A., 2017. The nature of subslab slow velocity anomalies beneath South America, Geophys. Res. Lett., 44, 4747–4755. Proctor, B. & Hirth, G., 2015. Role of pore fluid pressure on transient strength changes and fabric development during serpentine dehydration at mantle conditions: Implications for subductionzone seismicity, Earth Planet. Sci. Lett., 421, 1–12. Pujol, J., Chiu, J.M., Smalley, R., Regnier, M., Isacks, B., Chatelain, J.L., Vlasity, J., Vlasity, D., Castano, J. & Puebla, N., 1991. Lateral velocity variations in the Andean foreland in Argentina determined with the JHD method, Bull. Seismol. Soc. Am., 81, 2441–2457.

46 Qian, W., Wang, W., Zou, F. & Wu, Z., 2018. Elasticity of orthoenstatite at high pressure and temperature: Implications for the origin of low VP/VS zones in the mantle wedge., Geophys. Res. Lett., 45, 665–673. Ramos, V.A., 2004. Cuyania, an Exotic Block to Gondwana: Review of a Historical Success and the Present Problems, Gondwana Res., 7, 1009–1026. Ramos, V.A., Jordan, T.E., Allmendinger, R.W., Mpodozis, C., Kay, S.M., Cortés, J.M. & Palma, M., 1986. Paleozoic terranes of the central Argentine-Chilean Andes, Tectonics, 5, 855–880. Ramos, V.A., Cegarra, M. & Cristallini, E., 1996. Cenozoic tectonics of the High Andes of westcentral Argentina (30–36°S latitude), Tectonophysics, 259, 185–200. Ramos, V.A., Cristallini, E.O. & Pérez, D.J., 2002. The Pampean flat-slab of the Central Andes, J. South Am. Earth Sci., 15, 59–78. Ranero, C.R., Villaseñor, A., Phipps Morgan, J. & Weinrebe, W., 2005. Relationship between bendfaulting at trenches and intermediate-depth seismicity, Geochem. Geophys. Geosys., 6, Q12002. Rapela, C.W., Pankhurst, R.J., Casquet, C., Baldo, E., Saavedra, J. & Galindo, C., 1998. Early evolution of the Proto-Andean margin of South America, Geology, 26(8), 707. Rapela, C.W., Pankhurst, R.J., Casquet, C., Fanning, C.M., Baldo, E.G., González-Casado, J.M., Galindo, C. & Dahlquist, J., 2007. The Río de la Plata craton and the assembly of SW Gondwana, Earth-Science Rev., 83, 49–82. Regnier, M., Chiu, J., Smalley, R., Isacks, B.L. & Araujo, M., 1994. Crustal thickness variation in the Andean foreland, Argentina, from converted waves, Bull. Seismol. Soc. Am., 84, 1097– 1111. Reynard, B., 2013. Serpentine in active subduction zones, Lithos, 178, 171–185. Richardson, T., Gilbert, H., Anderson, M. & Ridgway, K.D., 2012. Seismicity within the actively deforming Eastern Sierras Pampeanas, Argentina, Geophys. J. Int., 188, 408–420.

47 Richardson, T., Ridgway, K. D., Gilbert, H., Martino, R., Enkelmann, E. Anderson, M. & Alvarado, P., 2013. Neogene and Quaternary tectonics of the Eastern Sierras Pampeanas, Argentina: Active intraplate deformation inboard of flat-slab subduction, Tectonics, 32, 780–796. Rivas, C., Ortiz, G., Alvarado, P., Podesta, M., Martin, A., 2019. Modern crustal seismicity in the northern

Andean

Precordillera,

Argentina.

Tectonophysics,

762,

144-158.

doi.org/10.1016/j.tecto.2019.04.019. Rüpke, L.H., Morgan, J.P., Hort, M. & Connolly, J.A.D., 2004. Serpentine and the subduction zone water cycle, Earth Planet. Sci. Lett., 223, 17–34. Scarfi, L., Raffaele, R., Badi, G., Ibanez, J.M., Imposa, S., Araujo, M. & Sabbione, N., 2012. Seismotectonic features from accurate hypocentre locations in southern central Andes (western Argentina), Tectonophysics, 518–521, 44–54. Schepers, G., van Hinsbergen, D.J.J., Spakman, W., Kosters, M.E., Boschman, L.M. & McQuarrie, N., 2017. South-American plate advance and forced Andean trench retreat as drivers for transient flat subduction episodes, Nat. Commun., 8, 15249. Scholz, C.H. & Small, C., 1997. The effect of seamount subduction on seismic coupling, Geology, 25(6), 487–490. Skinner, S.M. & Clayton, R.W., 2013. The lack of correlation between flat slabs and bathymetric impactors in South America, Earth Planet. Sci. Lett., 371–372, 1–5. Smalley R. & Isacks B.L., 1987. A high-resolution local network study of the Nazca Plate WadatiBenioff Zone under western Argentina, J. Geophys. Res., 92(B13), 13903–13912. Smalley, R., Pujol, J., Regnier, M., Chiu, J.-M., Chatelain, J.-L., Isacks, B.L., Araujo, M. & Puebla, N., 1993. Basement seismicity beneath the Andean precordillera thin-skinned thrust belt and implications for crustal and lithospheric behavior, Tectonics, 12, 63–76. Soustelle, V. & Tommasi, A., 2010. Seismic properties of the supra‐subduction mantle: Constraints from peridotite xenoliths from the Avacha volcano, southern Kamchatka, Geophys. Res. Lett., 37, L13307.

48 Thurber, C.H., 1992. Hypocenter-velocity structure coupling in local earthquake tomography, Phys. Earth Planet. Inter., 75, 55–62. Thurber, C. & Eberhart-Phillips, D., 1999. Local earthquake tomography with flexible gridding, Comput. Geosci., 25, 809–818. Thurber, C., Zhang, H., Brocher, T. & Langenheim, V., 2009. Regional three‐dimensional seismic velocity model of the crust and uppermost mantle of northern California, J. Geophys. Res. Solid Earth, 114, B01304. Ulmer, P. & Trommsdorff, V., 1995. Serpentine stability to mantle depths and subduction-related magmatism., Science, 268, 858–861. Urbina, N.P., Sruoga, P. & Malvicini, L., 1997. Late Tertiary Gold-Bearing Volcanic Belt in the Sierras Pampeanas of San Luis, Argentina, Int. Geol. Rev., 39, 287–306. van Hunen, J., van den Berg, A.P. & Vlaar, N.J., 2002. On the role of subducting oceanic plateaus in the development of shallow flat subduction, Tectonophysics, 352, 317–333. van Hunen, J., van den Berg, A.P. & Vlaar, N.J., 2004. Various mechanisms to induce present-day shallow flat subduction and implications for the younger Earth: a numerical parameter study, Phys. Earth Planet. Inter., 146, 179–194. van Keken, P.E., Hacker, B.R., Syracuse, E.M. & Abers, G.A., 2011. Subduction factory: 4. Depthdependent flux of H2O from subducting slabs worldwide, J. Geophys. Res., 116, B01401. Venerdini, A., Sánchez, G., Alvarado, P., Bilbao, I. & Ammirati, J-B., 2016. Nuevas determinaciones de velocidades de ondas P y ondas S para la corteza sísmica del terreno Cuyania en el retroarco andino. Rev. Mexicana de Ciencias Geol., 33(1), 59-71. von Gosen, W., 1998. Transpressive deformation in the southwestern part of the Sierra de San Luis (Sierras Pampeanas, Argentina), J. South Am. Earth Sci., 11, 233–264. von Huene, R., Corvalán, J., Flueh, E.R., Hinz, K., Korstgard, J., Ranero, C.R. & Weinrebe, W., 1997. Tectonic control of the subducting Juan Fernández Ridge on the Andean margin near Valparaiso, Chile, Tectonics, 16, 474–488.

49 Vujovich, G., Miller, H. & Ramos, V.A., 1994. Proterozoic metavolcanics from western Sierras Pampeanas terrane, Argentina, J. South Am. Earth Sci., 7, 309–323. Wagner, L.S., Beck, S. & Zandt, G., 2005. Upper mantle structure in the south central Chilean subduction zone (30° to 36°S), J. Geophys. Res., 110, B01308. Wagner, L., Beck, S., Zandt, G. & Ducea, M., 2006. Depleted lithosphere, cold, trapped asthenosphere, and frozen melt puddles above the flat slab in central Chile and Argentina, Earth Planet. Sci. Lett., 245, 289–301. Wagner, L.S., Anderson, M.L., Jackson, J.M., Beck, S.L. & Zandt, G., 2008. Seismic evidence for orthopyroxene enrichment in the continental lithosphere, Geology, 36(12), 935-938. Waldhauser, F. & Ellsworth, W.L., 2000. A Double-Difference Earthquake Location Algorithm: Method and Application to the Northern Hayward Fault, California, Bull. Seismol. Soc. Am., 90, 1353–1368. Wessel, P. & Smith, W.H.F., 1998. New, improved version of generic mapping tools released, EOS, Trans. Am. geophys. Un., 79, 579. Worthington, J.R., Hacker, B.R. & Zandt, G., 2013. Distinguishing eclogite from peridotite: EBSDbased calculations of seismic velocities, Geophys. J. Int., 193, 489–505. Yáñez, G.A., Ranero, C.R., von Huene, R. & Díaz, J., 2001. Magnetic anomaly interpretation across the southern central Andes (32°–34°S): The role of the Juan Fernández Ridge in the late Tertiary evolution of the margin, J. Geophys. Res. Solid Earth, 106, 6325–6345. Yáñez, G., Cembrano, J., Pardo, M., Ranero, C. & Selles, D., 2002. The Challenger–Juan Fernández–Maipo major tectonic transition of the Nazca–Andean subduction system at 33– 34°S: geodynamic evidence and implications, J. South Am. Earth Sci., 15, 23–38. Zhang, H. & Thurber, C.H., 2003. Double-difference tomography: The method and its application to the Hayward Fault, California, Bull. Seismol. Soc. Am., 93, 1875–1889. Zhang, H. & Thurber, C.H., 2006. Development and Applications of Double-difference Seismic Tomography, Pure Appl. Geophys., 163(B4), 6325–6345.

50 Zhang, H. & Thurber, C.H., 2007. Estimating the model resolution matrix for large seismic tomography problems based on Lanczos bidiagonalization with partial reorthogonalization, Geophys. J. Int., 170, 337–345. Zhang, H., Thurber, C.H., Shelly, D., Ide, S., Beroza, G.C. & Hasegawa, A., 2004. High-resolution subducting-slab structure beneath northern Honshu, Japan, revealed by double-difference tomography, Geology, 32(4), 361–364. Zhang, H., Thurber, C. & Bedrosian, P., 2009. Joint inversion for Vp, Vs, and Vp/Vs at SAFOD, Parkfield, California, Geochem. Geophys. Geosys., 10, Q11002. Zhao, D., Hasegawa, A. & Horiuchi, S., 1992. Tomographic imaging of P and S wave velocity structure beneath northeastern Japan, J. Geophys. Res. Solid Earth, 97(B13), 19909–19928.

SUPPORTING INFORMATION Additional supporting information may be found in the online version of this article.

51 Table 1: Summary of results and interpretations Plate Depth Region Trends in seismic (km) South

0-30

American

Interpretations

velocities (km/s) Chilenia

Vp 7.0-7.9

Mafic-ultramafic composition near

Terrane

Vs 3.6-4.2

Chilenia-Precordillera suture.

Vp/Vs > 1.74 South

0-30

Precordillera

American

Vp 5.8-6.4

Carbonate rocks and complexity in

Vs 3.2-3.6

crustal structure.

Vp/Vs 1.71-1.76 South

0-30

American

Cuyania

Vp 6.3-7.0

Terrane

Vs 3.5-3.9

Mafic-ultramafic composition.

Vp/Vs > 1.75 South

0-30

American

Pampia

Vp 6.0-6.7

Terrane

Vs 3.4-4.1

Felsic composition.

Vp/Vs 1.67-1.75 South

40-60

American

South

60-90

American

Lower crust

Vp 7.6-8.2

Partial eclogitization of lower

Cuyania

Vs 4.4-4.7

crust.

(Anomaly C)

Vp/Vs 1.73-1.75

Upper mantle

Vp 7.7-8.3

(Region D)

Vs 4.4-4.7

Depleted lherzolite or harzburgite.

Vp/Vs 1.75-1.77 South

70-90

American

Upper mantle

Vp 7.5-8.1

Orthopyroxenites or peridotites

(Anomaly A)

Vs 4.4-4.8

with 40-80% of enstatite.

Vp/Vs 1.70-1.73 South

70-90

American

Upper mantle

Vp 7.5-8.4

Up to 5% of hydration of depleted

(Anomaly B)

Vs 4.3-4.6

lherzolite or harzburgite.

Vp/Vs 1.78-1.82 Nazca

90-

Subducting

Vp 7.9 to 8.3

Progressive eastward dehydration.

105

oceanic crust

Vs 4.4 to 4.7

Without seismicity.

Vp/Vs 1.75-1.79 Nazca

105-

Subducting

Vp 8.0 to 8.5

120

oceanic

Vs 4.5-4.8

mantle

Vp/Vs 1.73-1.76

Increased Vp. Contains the WBZ.

52

Figure 1. Location map. (a) Tectonic context of the Pampean flat slab region. Plate-convergence velocity from Kendrick et al. (2003). The inland projection of the subducted Juan Fernandez Ridge (JFR) is based on Yánez et al. (2001). Yellow lines mark the slab contours every 10 km according to Anderson et al. (2007). (b) Morphotectonic units and terrane boundaries based on Ramos et al. (2002) and Rapela et al. (2007). Red lines mark the profiles shown in Figs 5, 7, 8, and S3. P1 means profile 1 and PP Pie de Palo range.

53

Figure 2. A. Distribution of the 1,157 earthquakes selected for the 1D and 3D inversions shown in (a) map view, (b) cross section in E-W direction, and (c) cross section in N-S direction. (d). Histogram of the number of events by depth. Gray ticked line represent the inland projection of the subducted Juan Fernandez Ridge (JFR) based on Yánez et al. (2001). Dashed lines in cross sections represent the Moho used in the initial velocity model.

54

Figure 3. Results of the M1DM for the P wave. (a) Initial and final velocity models. The gray area shows the velocity space explored from different initial velocities (gray lines). All initial models converge into a M1DM (shown in red). (b) P-wave station corrections for the M1DM. A yellow star marks the reference station CONE and near zero station corrections are indicated by squares. Crosses and circles denote positive and negative corrections, respectively, and are scaled by the magnitude of the correction (see values in Table S1). Minimum and maximum corrections are indicated below the map. Gray line represents the inland projection of the subducted Juan Fernandez Ridge (JFR).

55

Figure 4. (a) Trade-off curves between model variance and data variance for selecting damping and smoothing weight parameters. Each trade-off curve is indicated in a different color for a particular smoothing weight. Selected damping values are indicated by symbols. Our preferred smoothing value is 100 (red curve). (b) Residual misfit for both P and S waves in the initial locations and after the 1D and 3D inversions.

56

Figure 5. Solution quality assessment along latitude 31.4°S (Fig 1b). (a) The square root of the DWS for P waves when differential times are calculated by using event pairs with inter-event distance < 10 km. Contour spacing is every 10. (b) Same as (a) but for inter-event distance < 20 km. (c) Checkerboard initial model for Vp/Vs with high- and low-velocity anomalies of 10%. Only regions with DWS > 100 are shown. Black lines encompass regions with DWS > 2500 where checkerboard patterns are recovered. Small crosses denote the grid nodes. Topography is exaggerated by a factor of 8. Earthquakes (circles) and stations (triangles) within 0.2° of the profile are plotted. Contour spacing is every 2% of deviation. (d) Recovered model for Vp/Vs. Contour spacing is shown every 1% of deviation.

57

Figure 6. Synthetic checkerboard models for the case of high- and low-velocity anomalies (± 5 %) of opposite sign for Vp and Vs, so that the resulted Vp/Vs model has high and low anomalies of ±10 % (see more examples in Supporting Information). Contour spacing is shown every 2 % for the initial model and 1% for the recovered model. Only regions with DWS > 100 are shown. Black lines encompass regions with DWS > 2500 where checkerboard patterns are recovered. Small crosses denote the grid nodes. (a) Initial Vp/Vs model at 20 km depth. (b) Recovered Vp/Vs. (c) Initial Vp/Vs model at 50 km depth. (d) Recovered Vp/Vs. (e) Initial Vp/Vs model at 70 km depth. (f) Recovered Vp/Vs. (g) Initial Vp/Vs model at 90 km depth. (h) Recovered Vp/Vs.

58

Figure 7. Cross section along profile 1 (Fig 1b). Contour spacing is 0.2 km/s for Vp, 0.1 km/s for Vs, 0.02 for Vp/Vs, and 2% for deviation. Black lines encompass regions with DWS > 2500 and areas with DWS < 2500 are masked with lighter colors for absolute velocities. Earthquakes (circles) within 0.2° of the profile are shown for the deviation plots. Anomalies “A”, “B”, “C”, and “D” are marked by dotted red lines for absolute velocities and gray lines for deviation. Possible tectonic sutures are represented by green lines. CPFC means Chilenia and Principal and Frontal Cordilleras, PP Pie de Palo range, and PR Precordillera. Other symbols as in Fig. 5. (a) Absolute P-wave velocity. (b) P-wave velocity as a percentage change with respect to the initial model. (c) Absolute S-wave velocity. (d) S-wave velocity as a percentage change with respect to the initial model. (e) Vp/Vs. (f) Vp/Vs as a percentage change with respect to the initial model.

59

Figure 8. Cross sections through the velocity models. Latitude and longitude of each section is indicated on top of the profiles (see location in Fig. 1b). Contour spacing is 0.02 for Vp/Vs and 2% for deviation. Black lines encompass regions with DWS > 2500. Areas with DWS < 2500 are masked with lighter colors for the Vp/Vs results. Anomalies “A”, “B”, “C”, and “D” are marked by dotted red lines and possible tectonic sutures by green lines. Other symbols as in Fig. 7. Cuy means Cuyania, and RPC Rio de la Plata Craton.

60

Figure 9. Slices at a fixed depth through the velocity model in the upper mantle. Depth is indicated above each plot. Triangles represent the stations and the gray line the inland projection of the subducted Juan Fernandez Ridge (JFR). Contour spacing is 0.1 km/s for Vp and Vs, and 0.01 for Vp/Vs. Black lines encompass regions with DWS > 2500. Areas with DWS < 2500 are masked with lighter colors. Anomalies labeled “A”, “B”, and “D” are discussed in the text.

61

Figure 10. Slices at a fixed depth through the velocity model for the continental crust. Depth is indicated above each plot. Contour spacing is 0.1 km/s for Vp and Vs, 0.01 for Vp/Vs, and 2% for deviation. Terrane boundaries as in Fig. 1b, except for the Cuyania-Chilenia, which is shown as a dotted line in (c), (d), (e), and (f) to match the velocities changes at a depth of 30 km. Black lines encompass regions with DWS > 2500. Areas with DWS < 2500 are masked with lighter colors for absolute velocities and Vp/Vs. Other symbols as in Fig 9.

62

Figure 11. Location of the Vp/Vs anomalies discussed in the text (A, B, C, D, and W) in the context of predicted velocities of (a) common mantle mineralogies at 3GPa and temperatures from 500 to 1000ºC, (b) peridotites stable at 3GPa and temperatures from 500 to 900ºC (c) a range of eclogites at 1-3 GPa and temperatures from 400 to 700ºC. All seismic velocities were calculated using Abers and Hacker (2016) with the exception of the results taken from Qian et al. (2018). Abbreviations for each mineral are shown in (a). Types of peridotites and eclogites are based on Hacker et al. (2003a). Mineralogy of each rock is shown in parenthesis as percentage of volume. More specifically, we use iron-free olivine, chlorite (clin90daph10), orthopyroxene (en90fs10), clinopyroxene (di90hed10), and garnet (alm24gr10py66). Water content is also indicated as a percentage of weight (wt%) for the hydrated rocks. See Table 3 in Hacker et al. (2003a) for the exact composition of the eclogites mentioned.

63

Figure 12. Schematic conceptual model based on the interpretation of Vp/Vs results. (a) Perspective map with the location of profile and stations and terrane boundaries as in Fig. 1b. (b) Schematic profile along latitude 31.8ºS integrating relevant aspects in the vicinity of our study area from Booker et al. (2004), Ranero et al. (2005), Burd et al. (2013), and Marot et al. (2013). Main continental crustal structures (red lines) are from Ramos et al. (2002), Farías et al. (2010), Richardson et al. (2013), and Ammirati et al. (2016) and earthquakes (circles) from this study. ESP means Eastern Sierras Pampeanas, PFC Principal and Frontal Cordilleras, PR Precordillera, RPC Rio de la Plata Craton, and WSP Western Sierras Pampeanas.

Highlights Article Title: Lithospheric Structure of the Pampean Flat Slab Region from Double-Difference Tomography Authors: Lepolt Linkimer, Susan Beck, George Zandt, Patricia Alvarado, Megan Anderson, Hersh Gilbert, Haijiang Zhang •

Double-difference tomography reveals lithospheric composition and structure in Argentina

•

Relocated earthquakes refine slab geometry of Nazca Plate

•

Flat slab geometry suggest other factors for the flattening besides ridge subduction

•

Lithospheric suture may facilitate migration of water across the Pampean flat slab region

•

Mantle velocities in South American Plate suggest local enstatite enrichment and peridotite hydration