Our results suggest that the timing of the shortening events is a direct consequence of the interaction between the buckling subducting plate and the weakened overriding plate. We distinguish four notable deformation phases that correspond in amplitude, timing and space to the shortening rate from the geological compilation (Onken et al., 2012). Overall, deformation migrates across the orogenic domain to the eastern foreland.
Phase I : Central domain deformation (~6.5 to ~11 My, Fig. 4): Plastic strain is localized in the Central domain due to flat slab steepening and the partial removal of the lithosphere.
Phase II : Eastern Cordillera domain deformation (~11 to ~20 My, Fig. 4): Distributed plastic strain slowly accumulates in the east. No significant deformation is observed in the continent due to efficient trench retreat.
Phase III : Deformation migrates from the Eastern Cordillera to the foreland domain (~20 to ~29 My, Fig. 4): Strain intensifies in the Eastern Cordillera domain and migrates to the foreland, where the Brazilian Cratonic shield starts to underthrust below the orogen. The delamination follows this migration.
Phase IV : Foreland domain deformation (~29 to ~38 My, Fig. 4): Underthrusting of the shield slows down. At ~33.5 My, it re-accelerates until ~35 My before decelerating until 38 My.
The compressive stress generated by the difference of velocity between the trench and the overriding plate is accommodated in one of two ways: 1) orogenic shortening, 2) underthrusting of the foreland. The effectiveness of deformation localization depends on the strength of the overriding plate and the interplate coupling. Here, we discuss the key processes that affect the strength of the overriding plate, the subduction and deformation dynamics of the slab, and, finally, the interaction between the two plates.
Delamination. Extensive lithospheric delamination is known to have taken place under the Altiplano-Puna plateau (Kay & Kay, 1993; Beck & Zandt, 2002; Beck et al., 2015) and contributed to present-day elevations (Garzione et al., 2006, 2008, 2017; Wang et al., 2021). This process is thought to be the result of the eclogitization of the mafic lower crust and lithospheric mantle, which is likely facilitated by the hydration of the sub-lithosphere from the ~200 Ma subduction history (Babeyko et al., 2002, 2006), and thick (~45 km) initial crust at ~30 Ma (Hindle et al., 2005; Sobolev et al., 2006; Armijo et al., 2015). Model S3 demonstrates that without eclogitization delamination and shortening are inhibited (supplementary Fig. 2). Moreover, the lithospheric removal due to eclogitization leads to a localization of deformation and subsequent weakening in the overriding plate. Nevertheless, model S4 shows that a very weak orogenic domain localizes too much deformation in the orogen and does not guarantee the migration of the deformation to the foreland (supplementary Fig.5).
Due to flat slab steepening in Phase I, we observe two delamination stages after the first lithospheric removal of the overriding plate. First, the initial removal exposes the crust at the western edge that is directly in contact with the asthenosphere, thereby increasing its temperature and decreasing the viscosity at its base. As a result, the lower crust delaminates faster in the west, causing it to asymmetrically drip to the east (i.e., Stage 1 in Fig. 5a). The pure shear deformation localizes in the orogenic domain until delamination is complete. Second, when the viscous deformation of the orogen connects with the plastic deformation of its foreland at 26 My, the foreland underthrusts beneath the orogen due to weak sediments. This results in orogenic thickening and a switch from pure shear to simple shear shortening. Consequently, deformation migrates to the east causing delamination to accelerate (Stage 2 in Fig. 5b).
Mechanical weakening of the foreland sediments. In the Altiplano, the presence of weak sediments is the key factor in switching deformation from pure to simple shear at ~10 Ma. Simple-shear shortening is associated with higher strain localization over fewer faults and the formation of deep low-angle detachments. In the foreland, these faults are situated at the base of the sediment cover and are characteristic of the thin-skin deformation style. Porous sediment layers, in particular the paleozoic layers (Allmendinger and Gubbels, 1996), may have accumulated enough fluids at the front of the orogen to reduce their frictional strength to ~0.05 or less and initiate the underthrusting of the Brazilian cratonic shield (Babeyko et al., 2006).
This thin-skin style of deformation is often opposed to the thick-skin style, where strain is more distributed throughout the domain and may involve basement rock. At the latitude of the Puna, thick-skin deformation resulted in a final shortening amount much lower than in Altiplano (~150 km versus ~300 km; Kley et al., 1999; Sobolev and Babeyko, 2005; Babeyko & Sobolev, 2005). This shortening difference suggests that forces were accommodated elsewhere, which we suggest to be the retreating trench.
Commonly, thick-skin deformation is thought to result from the reactivation of pre-existing normal faults that formed in past extensional events (Carrera and Muñoz, 2013). The weak faults localize strain faster and enhance the shortening magnitude. However, their reactivation could also compete against an efficient switch from pure to simple shear deformation.
In the reference model, underthrusting takes place in two stages. The first stage happens during trench blockage at ~20.5 My, causing the deformation to migrate to the foreland. When the active brittle shear zone, from the failure of the foreland sediments, connects to the ductile shear zone accommodating the on-going delamination underthrusting becomes more efficient. The delamination also facilitates the underthrusting of the Brazilian cratonic shield that meets less resistive forces. Underthrusting of the shield forces the upper crust to viscously flow and thicken. The topography uplifts, reaching present-day elevations (~4 km) at ~31 My (~7 Ma ago). A second stage of underthrusting occurs in the last ~4 Ma when the trench again becomes blocked, but this event does not significantly change the topography (Fig. 4).
While the absolute motion of the South American plate provides the main force (Martinod et al., 2010; Husson et al., 2012) for the tectonic shortening, the magnitude of the compressive stress in the South American plate margin is determined by the the resistance of the Nazca plate (i.e., by the ability of the trench to retreat Lallemand et al., 2005; Funiciello et al., 2008; Lallemand et al., 2008; Holt et al., 2015). In the central Andes, the trench has migrated west over the last ~40 Ma as a result of the rollback and subsequent sinking of the bending slab in the asthenosphere, as well as the forced trench retreat from the excess velocity of the overriding plate (Schepers et al., 2017). Recent studies have proposed that the trench velocity can also be affected by deep subduction dynamics (Faccenna et al., 2017; Briaud et al., 2020, Boutoux etal., 2021). In this section, we discuss the implications of these subduction dynamics.
Flat slab steepening. The cause of flat subduction is still debated. It likely results from the shallowing of the slab from long lasting subduction, as well as greater buoyancies related to the Juan Fernandez ridge (Schellart, 2020; Schellart & Strak, 2021) that has migrated to the south in the last ~35 Ma (Fig. 1; Yáñez et al., 2001; Bello-González et al., 2018). In this study, we are interested in the consequence of slab steepening after the passage of the ridge. Our models suggest that a flat slab at ~100 km depth, analogous to the Pampean flat slab, could scrape the base of the lithosphere. Eventually at ~7 My, the slab steepens and accelerates as the trench becomes blocked (Fig. 6a). The continental mantle coupled to the flat slab segment is pulled down and viscously accommodates the deformation. When the lower crust eclogitizes, plastic strain localizes in the top portion of the crust, slab steepening then accelerates due to the eclogitization until, eventually, parts of the lithosphere are removed. This flat subduction plays a key role in triggering the initial weakening of the overriding plate, and is facilitated by lower-crustal eclogitization .
Buckling instability cycles. We identified two buckling cycles, at ~20 My and at ~30 My. Within each cycle, three main events are distinguished that may affect the trench migration rate:
(1) Slab impediment (Fig. 6b) takes place when the slab meets viscous resistance. This is the case when the slab is impeded by the viscous lower mantle at the beginning of a buckling cycle (~17 My and ~29.5 My), or before steepening. For instance, when the slab reaches the viscous lower mantle it does not immediately penetrate it. The first slab segment in contact with the lower mantle slows down and viscously resists the new, still sinking, segment. This difference of velocity between the two segments is accommodated through bending in the slab. During these slab impediment events the dip of the slab becomes shallower and the trench continues retreating. This mechanism differs from slab anchoring (Faccenna et al., 2017), in which the difference of velocity between the two segments is too small to cause the folding of the slab.
(2) Slab folding (Fig. 6c) events occur when, after slab impediment, the slab dip flips in the transition zone. The now shallower slab dip enable the trench retreat, though no significant deformation is observed. Each buckling cycle consists of two folding events, the first to the west and the second to the east at ~20, 21 My and ~30, 33 My.
(3) Slab steepening (Fig. 6d) is a drastic event that occurs at the end of a buckling cycle after the second folding event, (~23.5 My and ~33.5 My). Chronologically, the sinking slab meets resistance from the last fold to the east (i.e., Impediment) and bends to the west as for the first folding event. However, the overriding plate has forced the trench to retreat during the previous events, which, prevents the slab from piling up. The slab continues to sink in the transition zone, steepens and accelerates. The trench slows down and blocks the overriding plate that shortens to accommodate the ongoing deformation. When the trench is blocked the horizontal stress in the overriding plate can reach values of ~350 MPa (supplementary Fig. 1, movie S1b). Overall, slab shallowing is associated to periods of trench retreat related to the folding events, whereas slab steepening is associated to periods of trench blockage following folding events folding events.
Interaction between overriding and subducting plates
Interplate coupling. Our models predict that an effective friction of 0.35 to 0.05 is required in the Central Andes to obtain significant deformation that is consistent with previous estimates (Sobolev and Babeyko, 2005; Sobolev et al., 2006). Higher friction values result in lower oceanic-plate velocities. The effective friction is dependent on the sediment thickness at the trench, which at present day may vary from ~0.5 km to ~2 km in the Central and Southern Andes, respectively (Lamb and Davis, 2003). This latitudinal variation results from the efficiency at which the surface processes supply sediments to the trench. In the last ~6 Ma, glacial erosion supplied a large amount of sediments to Southern Andes trench. Whereas in the Central Andes, the internal drainage of the Altiplano-Puna plateau is related to low erosional rates that have contributed to sediment starvation at the trench (Lamb & Davis, 2003).
Slab buckling and overriding plate interaction. The unusual timing of the growth of the Andes results from a sequence of events generated by plate interactions. While subduction dynamics exert a major control on the deformation of the sinking plate by blocking trench migration, the strength of the overriding plate ultimately controls where strain localizes and forces the trench to retreat when it is not blocked. This plate strength is evolving, first, with the passage of the flat slab that may have initially weakened the lithosphere through partial removal of the mantle lithosphere, and through thermal weakening related to crustal exposure near the hotter asthenosphere (Isacks, 1988), and second, by triggering the subsequent delamination (see previous section).
Pulsatile behavior in the deformation of the Nazca plate is observed in paleoelevation reconstructions (Boschman, 2021 ; Garzione et al., 2008), the magmatic activity (Decelles et al., 2009), and from stable isotope data (Leier et al., 2013), We suggest that buckling instabilities in a subducting plate offer a plausible explanation in the variability and timing of the Nazca plate deformation during the last ~20 Ma as well as the present-day deep seismicity distribution (supplementary Fig.7b). We find that shortening rate pulses occur at the end of each buckling cycle when slab steepening inhibits trench retreat (Fig. 7cd), and that these pulses reproduce similar signals to what is seen in the geological data. Additionally, in the last ~2 Ma the geological data shows a decrease in the shortening rate, which is also predicted by our model through underthrusting. At later stages, the trench retreat resumes and underthrusting loses its efficiency, which could indicate the beginning of a new buckling cycle.
Previous studies have suggested that the lower mantle viscosity and the dip, age, thickness and strength of the oceanic plate may affect the buckling periodicity and timing of slab stagnation in the transition zone (Ribe et al., 2007; Lee & King, 2011; Quinteros et al., 2010; Capitanio et al., 2011; Quinteros & Sobolev, 2013; Cerpa et al., 2014; Marquardt & Miyagi, 2015; Briaud et al., 2020, Boutoux et al., 2021). Analyzing the variety of interchangeable parameters affecting the buckling process exceeds the scope of this study.
Previous seismic tomography studies indicate two large negative seismic anomalies near the transition zone ( at depths of 600 km and 900 km) that are attributed to slab accumulations (Widiyantoro, 1997; Liu, 2003; Chen et al., 2019). The deeper accumulation may relate to a slab anchoring (Faccenna et al., 2017, Supplementary Fig.7), suggesting that previous accumulation cycles could have occurred before and have “avalanched” in the lower mantle (Briaud et al., 2020; Hu & Gurnis, 2020), wherein they may have become detached from the shallower slab. Indeed, over the last ~200 Ma quick alternations between compressive and extensive phases (e.g., the compressive peruvian phase or extensive Salta rift between ~120 Ma and ~60 Ma; Faccenna et al., 2017) may indicate that slab buckling events have happened earlier in the subduction history. However, because of the absence of an efficient weakening mechanism to trigger delamination and too thin crust to facilitate eclogitization, the orogen experienced no significant deformation. Potentially, we suggest that these avalanche events may have repeated at least 3 times over the last ~90 Ma, as suggested by the 3 cycles of convergence rate recognized in Martinod et al., (2010).