# Role of erosion and isostasy in the Cordillera Blanca uplift: Insights from landscape evolution modeling (northern Peru, Andes)

Audrey Margirier, Jean Braun, Xavier Robert, Laurence Audin

<figure><img src="image_1.png" /></figure>

## To cite this version:

Audrey Margirier, Jean Braun, Xavier Robert, Laurence Audin. Role of erosion and isostasy in the Cordillera Blanca uplift: Insights from landscape evolution modeling (northern Peru, Andes). Tectonophysics, 2018, 728-729, pp.119 - 129. 10.1016/j.tecto.2018.02.009. hal-03135906

HAL Id: hal-03135906

https://hal.science/hal-03135906v1

Submitted on 19 Jun 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Originally published as:

Margirier, A., Braun, J., Robert, X., Audin, L. (2018): Role of erosion and isostasy in the Cordillera Blanca uplift: Insights from landscape evolution modeling (northern Peru, Andes). - *Tectonophysics*, 728-729, pp. 119-129.

DOI: http://doi.org/10.1016/j.tecto.2018.02.009
1

# Role of erosion and isostasy in the Cordillera Blanca uplift: insights

2

## from landscape evolution modeling (northern Peru, Andes)

3

4

Audrey Margirier¹,², Jean Braun¹, Xavier Robert², Laurence Audin²

5

¹Helmholtz-Zentrum Potsdam, GeoForschungsZentrum (GFZ) Potsdam, Potsdam, Germany

6

²Université Grenoble, Alpes, CNRS, IRD, ISTerre, F-38000 Grenoble, France

7

8

## Highlights:

9

- Inversion of the landscape evolution coupled with thermochronological data provides constraints on erosion efficiency factor, uplift rates and geothermal gradient

10

11

- Isostatic effect of eroding a denser rock mass represent a not negligible contribution to the Cordillera Blanca uplift on a < 5 Ma time scale

12

13

- Cordillera Blanca drainage divide location is controlled by initial drainage network rather by maximum uplift rates and precipitation distribution

14

15

16

## Abstract

17

The processes driving uplift and exhumation of the highest Peruvian peaks (the Cordillera Blanca) are not well understood. Uplift and exhumation seem closely linked to the formation and movement on the Cordillera Blanca normal fault (CBNF) that delimits and shapes the western flank of the Cordillera Blanca. Several models have been proposed to explain the presence of this major normal fault in a compressional setting, but the CBNF and the Cordillera Blanca recent rapid uplift remain enigmatic. Whereas the Cordillera Blanca morphology demonstrates important erosion and thus a significant mass of rocks removal, the impact of erosion and isostasy on the evolution of the Cordillera Blanca uplift rates has never been explored. We address the role of erosion and associated flexural rebound in the uplift and exhumation of the Cordillera

18

19

20

21

22

23

24

25

26
Blanca with numerical modeling of landscape evolution. We perform inversions of the broad features of the present-day topography, total exhumation and thermochronological data using a landscape evolution model (FastScape) to provide constraints on the erosion efficiency factor, the uplift rate and the temperature gradient.

Our results evidence the not negligible contribution of erosion and associated flexural rebound to the uplift of the Cordillera Blanca and allow us to question the models previously proposed for the formation of the CBNF.

**Keywords:** Rock uplift; Isostatic effect of eroding a denser rock mass; Normal fault; Low-temperature thermochronology; Numerical modeling of the landscape evolution; Cordillera Blanca; Peruvian Andes

# 1. Introduction

In mountain ranges, surface uplift is usually assumed to be the result of shortening and crustal thickening. Surprisingly, in northern Peru, uplift of the footwall of an active normal fault is responsible for the formation of the highest Peruvian summits in the Cordillera Blanca (Fig. 1). Several models have been proposed [Dalmayrac and Molnar, 1981; McNulty and Farber, 2002] to explain this unusual situation, but the processes driving both the Cordillera Blanca uplift and extensional deformation along the Cordillera Blanca normal fault (CBNF) remain poorly constrained. The CBNF trends parallel to the Andean range and is the most spectacular normal fault in the Andes [Fig. 1; Margirier et al., 2017]: the CBNF is ~200 km long and shows ~7 km of vertical offset in total [Margirier et al., 2016], it has been active since ~5.4 Ma [Bonnot, 1984; Giovanni, 2007]. The CBNF is located above the Peruvian flat-slab [Barazangi and Isacks, 1976], a section of the convergent plate boundary between the Nazca Plate and the South
52 American Plate characterized today by near-horizontal subduction geometry. The
53 Cordillera Blanca and the Cordillera Negra form, respectively, the hanging wall and the
54 footwall of the CBNF (Fig. 1).

55 The Cordillera Blanca fast exhumation rate (~1 mm/yr) has been previously linked to
56 motion on the CBNF [e.g., Bonnot, 1984; McNulty and Farber, 2002; Giovanni, 2007;
57 Margirier et al., 2015]. New thermobarometry data and erosion rates reconstruction
58 based on thermochronological data indicate a recent, i.e. Early Pleistocene, increase in
59 erosion rate in the Cordillera Blanca [~2-0 Ma; Margirier et al., 2016]. Margirier et al.
60 [2016] suggested that an important isostatic contribution from glacial erosion may
61 explain the recent exhumation of the Cordillera Blanca batholith. Indeed, the removal of
62 such a mass of material represents a significant upward unloading on the lithosphere,
63 which should drive substantial flexural uplift. This unloading and flexural uplift would
64 have also generated large differential stresses in the lithosphere, which could have
65 caused the reactivation of pre-existing structures such as the CBNF. Previous studies
66 demonstrated that the flexural uplift driven by alpine-type valley incision could reach
67 rates similar to those caused by tectonic processes [Montgomery, 1994; Small and
68 Anderson, 1995; Cederbom et al., 2004; Stern et al., 2005]. Recently Braun et al. [2014]
69 proposed that erosion-driven isostatic rebound should scale with the density of surface
70 rocks: denser rocks, such as a granitic body intruded in sedimentary rocks, rebound and
71 therefore are exhumed faster than the surrounding less dense rocks.

72 The rapid uplift of the Cordillera Blanca, the large volume of eroded rocks since the
73 emplacement of the Cordillera Blanca batholith and its location in the footwall of an
74 active normal fault in a compressive plate boundary, make the Cordillera Blanca the
75 perfect place to question the nature and efficiency of potential feedbacks between
76 erosion and uplift along the CBNF. The aim of this paper is thus (i) to test whether the
77 increase of erosion rate suggested for the last 2 Ma in the Cordillera Blanca [Margirier et
al., 2016] could be due to an increase of rock uplift rates since 2 Ma rather than a change of climate and/or erosion process and (ii) to quantify the importance of isostatic rebound associated with valley incision and erosion of denser rocks to explain the uplift of the Cordillera Blanca and to test the adequacy of a flexure-driven model in such a setting. To address this issue, we have attempted to model landscape evolution in the Cordillera Blanca using a landscape evolution model (LEM), in this case based on the FastScape algorithm [Braun and Willett, 2013].

## 2. Context

### 2.1 Geologic and climatic context

The Cordillera Blanca hosts the highest Peruvian summits with a cluster of 6000 m peaks (Fig. 1). It hosts a large 14–5 Ma granitic pluton [zircon U-Pb; Mukasa, 1984; Giovanni, 2007] emplaced at ~3 km depth into deformed Jurassic sediments [Margirier et al., 2016]. The Cordillera batholith is elongated (150 × 15 km) and trends parallel to the Andean range (Fig. 1A). Based on apatite fission-tracks and (U-Th/He) dating several studies gave estimations of exhumation rates ranging between 1 and 2 mm/yr in the central part of the Cordillera Blanca for the last 3-4 Myr [Giovanni, 2007; Hodson, 2012; Margirier et al., 2015]. This exhumation phase is likely associated to rock uplift and CBNF activity. On a shorter time scale (30-0 ka), the uplift rates have been constrained by ¹⁰Be dating of scarps along the CBNF [Siame et al., 2006] and geomorphic features (moraines) displaced by the fault [Schwartz, 1988; Gérard et al., in prep]. The vertical slip rates decrease from north to south ranging from 5.1 ± 0.8 mm/yr to 0.6 ± 0.2 mm/yr [Schwartz, 1988; Siame et al., 2006; Margirier et al., 2017; Gérard et al., in prep]. Whereas the higher peaks in the Cordillera Blanca are located close to the CBNF and should therefore correspond to the region of maximum tectonic uplift, the drainage
divide is located in the eastern part of the Cordillera Blanca, ~15 km away from the
103
CBNF (Fig. 1B). The batholith is deeply incised by deep U-shaped valleys resulting from
104
recent glacial erosion. Several other glacial landforms (moraines, roches moutonnées)
105
are evidence of the large extension of the glaciations that shaped the Cordillera Blanca
106
morphology [Farber et al., 2005]. Today, the Cordillera Blanca summits form the Andean
107
drainage divide [Wise and Noble, 2003] and act as an orographic barrier to moisture
108
coming from the Amazon basin [Fig. 1B; Bookhagen and Strecker, 2008]. Even if the most
109
important orographic effect of the Andes is controlled by the lower relief of the sub-
110
Andes, the Cordillera Blanca high elevation prevents moisture from reaching the
111
western flank of the Andes, resulting in wetter climatic conditions in the Cordillera
112
Blanca (mean rainfall ~1.5 m/yr) than in the Cordillera Negra (mean rainfall ~0.5
113
m/yr), farther to the west [Bookhagen and Strecker, 2008].
114

The morphology of the Cordillera Negra is that of a 4500 m high plateau incised by 1-2
115
km deep valleys along its western flank (Fig. 1). The range hosts Cretaceous and
116
Paleogene plutons [73-48 Ma; Beckinsale et al., 1985] intruded into Jurassic sediments.
117
Neogene volcano-sedimentary deposits cap the Cordillera Negra [54-15 Ma Calipuy
118
Formation; Cobbing et al., 1981]. A regional surface uplift associated with the subduction
119
geometry change and dynamic topography process has been evidenced in the Cordillera
120
Negra from 15 Ma [Eakin et al., 2014; Margirier et al., 2015]. Rare and discrete moraines,
121
only seen above ~4200 m, indicate limited ice cover and resulting glacial erosion in the
122
Cordillera Negra [Bonnot, 1984].
123

124

## 2.2 Paleogeography

Several studies have documented the Late Miocene paleogeography of the Cordillera
126
Blanca region [Wise and Noble, 2003; Giovanni, 2007; Hoorn et al., 2010]. Based on
127
pollen analyses, Hoorn et al. [2010] constrained the elevation in the Cordillera Blanca
128
region to be ~4 km in the Middle Miocene. However, at that time, the Cordillera Negra 129 topography formed the main drainage divide [Wise and Noble, 2003]. Based on the 130 location, age and stratigraphy of volcanic deposits in the Callejón de Huaylas Basin that 131 separates the two ranges, Wise and Noble [2003] and Giovanni [2007] suggested that a 132 depression already existed between the Cordillera Blanca and the Cordillera Negra at 133 the end of the Late Miocene. In addition, based on δ¹⁸O analyses of paleolake deposits, 134 Giovanni et al. [2010] showed that highest elevations in the Callejón de Huaylas Basin 135 were attained in the latest Miocene. Little is known, however, about the topography of 136 the Cordillera Blanca before the emplacement of the batholith. 137

138

## 3. Landscape evolution model: FastScape
139

### 3.1 Model
140

We used the FastScape algorithm [Braun and Willett, 2013] to solve the stream power 141 law to predict landscape evolution following a set of tectonic forcing (uplift) and initial 142 topography (geomorphic setting). Because of the optimum ordering of the nodes, the 143 algorithm is implicit in time and computationally very efficient, requiring only $O(n)$ 144 operations where $n$ is the number of points used to discretize the topography. 145 Consequently, FastScape can be used repetitively, even if using a very high spatial 146 discretization, to attempt to reproduce the main features of present-day topography 147 and, using a simple 1D thermal model, predicted cooling ages (apatite fission-tracks and 148 (U-Th/He) ages which can be related to erosion) in the hope of deducing a plausible 149 first-order topographic and uplift history. To do this, we coupled FastScape to the 150 neighborhood algorithm (NA) [Sambridge, 1999a; Sambridge, 1999b] to “invert” the 151 thermochronological ages, the barometric constraints and the known, final topography 152 in order to find the best-fitting values of several unknown parameters, including the 153
erosion efficiency ($K_f$), the elastic thickness ($T_e$), the uplift rates ($U_1$) at several periods in the past and the temperature at the base of the model ($T_{max}$). Using NA, we first carried out a large number of FastScape runs with parameter values randomly selected between fixed limits. For each run the model predictions were compared to the data and a misfit function was estimated, which was used to select a new set of model parameters to be tested. This procedure was repeated several times until an optimum combination of parameter values was found and, more importantly, the shape of the misfit function could be mapped. This allows us to analyze and discuss the range of model parameters where the model prediction fits the observations, as well as the sensitivity of the model predictions to the model parameters.

## 3.2 Erosion and isostasy

To represent surface erosion, we solve the stream power law that is commonly used to parameterize bedrock incision by rivers in steep mountainous terrains:

$$ (1) \quad \frac{\partial h}{\partial t} = U - K A^m S^n $$

where $h$ is topographic height, $t$ is time, $U$ is uplift rate, $A$ is drainage area and $S$ is slope in the direction of water flow, taken here to be the steepest path between a point on the landscape and any among its eight neighbors. Parameters $K$, $m$ and $n$ are poorly constrained constants that mostly depend on lithology and climate. We will use $m = 0.4$, $n = 1$ and vary $K$ around a mean value of $1.5 \times 10^{-5} m^{1-2m}yr^{-1}$ (see Croissant and Braun [2013] for a discussion on the value of these parameters). We also assume that hillslopes are subject to mass transport at a rate that is simply linearly proportional to topographic slope. Assuming mass conservation, this leads to a law for the rate of
change of topography on hillslopes that is proportional to the curvature of topography.

Combining this with the stream power law leads to the following evolution equation for the surface topography:

182

183

$$ (2) \quad \frac{\partial h}{\partial t} = U - K A^m S^n + K_D \left( \frac{\partial h^2}{\partial x^2} + \frac{\partial^2 h}{\partial y^2} \right) $$

184

where $K_D$ is a transport coefficient or diffusivity. We will assume an arbitrary value of $K_D$ = 0.3 m²/yr. We realize that using the stream power law and a simple diffusion-like representation of hillslope processes may seem inappropriate to model topographic evolution in a glaciated landscape. We justify our choice to use the stream power law by the fact that (i) no algorithm exists to solve the glacial erosion equations in an efficient manner that would allow for the type of inversion we have conducted here, (ii) we are interested in the effect of denser rocks erosion and the resulting isostatic rebound on the range-scale shape of the Cordillera Blanca and not in the details of the landform, and (iii) for both fluvial and glacial erosion processes, erosion rate is primarily controlled by slope and the geometry of a first order drainage system.

195

To account for the isostatic rebound associated with surface erosion, we also solve the following bi-harmonic equation representing the flexure of a thin elastic plate subject to surface loading/unloading:

199

200

$$ (3) \quad D \left( \frac{\partial^4 w}{\partial x^4} + 2 \frac{\partial^4 w}{\partial x^2 \partial x^2} + \frac{\partial^4 w}{\partial y^4} \right) = \Delta \rho g w + \rho_s g \Delta h $$

201

where $w$ is the surface displacement associated with the isostatic adjustment of an increment in erosion $\Delta h$, $D$ is flexural, $\rho_s$ is the density of surface rocks and $\Delta \rho$ is the

203
204 density contrast between asthenospheric density, $\rho_a$ and surface rock density [Turcotte and Shubert, 1982].
205

206

207 A short word of explanation is necessary to understand how we have coupled Equation
208 (2) and Equation (3). As explained above, we assume that the imposed uplift rate does
209 not contribute to the load applied to the thin elastic plate representing the lithosphere;
210 only the eroded material does. This means that we can rewrite Equation (2) in the
211 following form:

$$ \frac{\partial h}{\partial t} = U - F(e) $$

212 where e is erosion rate:

$$ e = KA^m S^n - K_D \left( \frac{\partial^2 h}{\partial x^2} + \frac{\partial^2 h}{\partial y^2} \right) $$

213 and $F(e)$ is the rate of isostatic response to this erosion rate. $U$ is the uplift rate that
214 would result in the absence of erosion. We can make this equation more readily
215 understandable, if we assume local isostasy, in which case, it becomes:

$$ \frac{\partial h}{\partial t} = U - \left(1 - \frac{\rho_s}{\rho_a}\right) e $$

216 We see that when erosion rate is nil, the rate of surface uplift is $U$ but that the steady
217 state erosion rate, $e_{ss}$, i.e. corresponding to $\partial h / \partial t = 0$:

$$ e_{ss} = U / \left(1 - \frac{\rho_s}{\rho_a}\right) $$

218 is larger than $U$. Consequently, during a numerical experiment, we expect to see the
219 erosion rate progressively increase to values that can reach 5-6 times the imposed uplift
220 rate. Flexure should dampen this progressive increase of the erosion rate as a function
221 of time.

222

223 **3.3 Age predictions**
As stated earlier, we will constrain our inversion by testing model predictions against observed cooling ages obtained for a series of samples collected in the Cordillera Blanca [Montario, 2001; Giovanni, 2007; Hodson, 2012; Margirier et al., 2015]. At each point of the model where we wish to predict a cooling age, we solve the 1D heat transport equation by conduction and advection:

229

230

231

where $κ$ is thermal diffusivity and $v$ is the erosion rate predicted by the landscape evolution model. We neglect the effect of radiogenic heat production and assume that the top and base of the crust/lithosphere are held at constant temperature: $T(z=0)=0$ and $T(z=L)=TL$. The thermal history obtained by solving this equation through time is then used to compute cooling ages by using standard methods developed to simulate the annealing of fission tracks and the diffusion of He in apatite [Braun et al., 2006].

238

The thermochronological dataset consists of apatite fission-track (AFT) ages and (U- Th)/He (AHe) mean ages [Montario, 2001; Giovanni, 2007; Hodson, 2012; Margirier et al., 2015]. The AFT and AHe thermochronological systems record the cooling histories of rocks below 120°C and 80°C, respectively; at a given location subject to steady exhumation, AFT age should be older than AHe age. For our inversions we exclude AFT ages that are older than the age of the emplacement of the Cordillera Blanca batholith (~7 Ma; Mukasa, 1984; Giovanni, 2007) and the AHe ages that are older than AFT age for a same sample. Indeed, the AHe ages are scattered, raising the question of their reliability [Margirier et al., 2015].

248

$$ (6) \quad \frac{\partial T}{\partial t} + v \frac{\partial T}{\partial z} = \kappa \frac{\partial^2 T}{\partial z^2} $$
249 The misfit function used to guide the inversion procedure is defined by:

$$ (7) \quad \mu = \frac{1}{N} \sqrt{ \sum \frac{(Age_o - Age_p)^2}{\partial \Delta Age^2} + \frac{|Hmax_o - Hmax_p|}{\Delta Hmax} + \frac{|Vmax_o - Vmax_p|}{\Delta Vmax} } $$

251 where N is the number of thermochronological ages, Ageoi and Agepi the observed and predicted ages and ΔAgei the uncertainty on the observed ages. Hmaxo and Hmaxp are the observed (present-day) and predicted maximum topography of the Cordillera Blanca (6600 m) and ΔHmax its regional variability (200 m). Vmaxo and Vmaxp are the observed (present-day) and total exhumed volume of the Cordillera Blanca (1 × 10¹³ m³), which we derive from amphibole barometry data [Margirier et al., 2016], and ΔVmax its assumed uncertainty (5 × 10¹² m³).

258

259 **3.4 Parameterization**

260 For each inversion, we ran 15000 models; during the first iteration, we sampled 3000 models to explore the parameter space, then we performed 15 iterations of 800 models with 400 cells resampled. The runs were performed on a 200 × 200 km square domain (Fig. 2), discretized by 500 × 500 regularly space nodes. All model runs last 7 Ma, which is approximately equivalent to the age of the Cordillera batholith emplacement (see geochronological data compilation in Margirier et al., 2016). We chose to take the geological ages as references for the modelisations, thus the model runs begin at 7 Ma and end at 0 Ma. We fixed the time step length to 1 ka (other model parameters are detailed in the Table A1). We used an arbitrary initial topography including both a "proto Cordillera Negra" (3500 m) and, for inversions 1b and 2, a "proto Cordillera Blanca" (3000 m), following estimates from Wise and Noble [2003] and Hoorn et al. [2010]. We imposed a uniform precipitation rate of ~0.5 m/yr in the Cordillera Negra and ~1.0 m/yr in the Cordillera Blanca in accordance with present-day estimates from Bookhagen and Strecker [2008]. According to the Cordillera Blanca batholith lithology
and shape [granite/granodiorite with a typical density ~ 2800 kg/m³; Petford and Atherton, 1992; Sharma, 1997] and the nature of the surroundings sedimentary rocks [alternating sandstone and marl with sporadic limestone with an estimated density ~2400 kg/m³; Bonnot, 1984; Sharma, 1997], we included a 150 x 20 km ellipsoidal intrusion that is characterized by an anomalously high density (400 kg/m³ heavier than the surrounding rocks). The top of the intrusion is initially at 1.5 km beneath sea level (~3 km depth). This is constrained by amphibole barometry data obtained from the roof of the Cordillera Blanca batholith [Margirier et al., 2016]. We fixed a regional uplift rate $U_0$ = 0.2 mm/yr in agreement with the thermal history provided by Margirier et al. [2015] for the Cordillera Negra. In addition to the widespread uplift ($U_0$) we simulated the presence of an active fault by imposing a steep gradient in vertical displacement rate along a north-south boundary located a distance of 100 km from the left side of the model (Fig. 2). Along the boundary the uplift rate is equal to $U_1$ and it decreases linearly to reach a value of 0 at 25 km away from the boundary. The uplift history along this boundary is divided into 3 successive episodes based on observations concerning the initiation of the CNBF and an apparent increase in erosion rate evidenced by Margirier et al. [2016]. From 7 to 5.4 Ma (before the CBNF initiation) the uplift rate ($U_1$) is equal to 0; between 5.4 and 2 Ma, the uplift rate is constant (and equal to $U_1$); between 2 Ma and the present, the uplift rate is increased by a factor $f_{U1}$. Perez-Gussinye et al. [2009] estimated elastic thickness ($T_e$) to be in the range [0-10] km, which are relatively low values for continental regions; we use a value of 3 km in inversions where the elastic thickness is fixed. The four parameters that we want to constrain by inversion of the thermochronological ages and the range-scale present-day topography and eroded volume in our inversions, are the erosion efficiency ($K_f$), the mean uplift rate ($U_1$), the magnitude of the recent increase of uplift rate in the Cordillera Blanca ($f_{U1}$), the
299

temperature at the base of the model ($T_{max}$) and, in a third set of inversions, the elastic thickness of the lithosphere ($T_e$).

300

301

302

# 4. Results

303

## 4.1 Role of tectonics, erosion and initial topography

304

Here we aim to test (i) if the increase of erosion rate suggested for the last 2 Ma in the Cordillera Blanca [Margirier et al., 2016] could be due to an increase in uplift rates since 2 Ma rather than a change of climate and/or erosion process (glacial erosion vs fluvial erosion) and (ii) the role of initial topography in the present day Cordillera Blanca drainage divide location. In these inversions, the four parameters that we wanted to constrain by inversion of the thermochronological ages, the total exhumation and the range-scale present-day topography are the erosion efficiency ($K_f$), the mean uplift rate ($U_1$), the magnitude of the recent increase of uplift rate in the Cordillera Blanca ($f_{U1}$) and the temperature at the base of the model ($T_{max}$).

305

306

307

308

309

310

311

312

313

We performed two inversions, one without an initial topography in the Cordillera Blanca and a second one with a "proto Cordillera Blanca", in order to evaluate the paleogeography and inherited drainage network influence on the location of the present day drainage divide in the Cordillera Blanca.

314

315

316

317

318

### 4.1.1 Inversion 1a (without a "proto Cordillera Blanca")

319

Figure 3 presents the results of inversion 1a (without a "proto Cordillera Blanca") as values of the misfit function, $\mu$, displayed in parameter space. In each of the two panels, each dot represents a single model run and the color scale refers to the value of the misfit of the model (from blue for high misfit value or poor fit to the data, to red for low

320

321

322
323
324

misfit value or good fit to the data). In each of the panels, the red star indicates the
position of the best-fitting model in parameter space.

325
326
327
328
329
330
331
332

The mean uplift rate ($U_1$), the amplitude of its recent increase ($f_{U1}$) and the erodibility ($K_f$) are all well constrained by the thermochronological data (Fig. 3) although a trade-off exists (as shown in figure 3B) between $U_1$ and $f_{U1}$. These two parameters are clearly correlated with higher mean uplift rate requiring a smaller recent increase in uplift rate and vice-versa. The arcuate shape of the region of minimum misfit in Figure 3B can also be interpreted by stating that the product $U_1 \times f_{U1}$ (or in other words the uplift rate between 2 and 0 Ma) is very well constrained by the data at a value of approximately 1.8 mm/yr (the best fitting parameters are indicated by a red star on figures 3A and 3B).

333
334
335
336
337
338
339
340
341
342

The temperature at the base of the model is not as well constrained (Fig. 3A) although values between 650 and 850°C yield the lowest misfit values. These values should be interpreted as suggesting best fitting surface geothermal gradient values of 21 to 28°C/km. The inversion results suggest that the data is best explained by values of $K_f$ = 1.9 × 10⁻⁵ m¹⁻²m·yr⁻¹, $U_1$ = 0.80 mm/yr, $f_{U1}$ = 1.8 and $T_{max}$ = 677°C. The set of the good fitting value of $f_{U1}$ (misfit < 0.6) indicates that a substantial increase in erosion rate in the recent past is necessary to explain the age dataset and the total exhumation. The best value for the basal temperature ($T_{max}$ = 677°C) corresponds to a geothermal gradient of ~23°C/km, which is consistent with Henry and Pollack [1988] terrestrial heat flow measurements.

343
344
345
346
347
348

The misfit is 0.53 for the best fitting parameters model, suggesting that the model can reproduce the data (present-day topography, total exhumation and thermochronological ages, Fig. 3). In addition, the thickness of rock eroded in the model after 7 Ma reproduces the total eroded volume we fixed to correspond with field observations and granite emplacement depth obtained by Margirier et al. [2016]. However, the apatite fission track ages are not well reproduced by the best fitting parameters model, the
349

apatite (U-Th)/He ages predicted by the model are too young (Fig. 3C, D). Moreover, the best fitting parameters model does not reproduce all characteristic features of the present-day Cordillera Blanca topography. Notably, the drainage divide of the predicted topography is located close to the CBNF, which is not the case in the Cordillera Blanca (Figs. 1, 4).

350

351

352

353

354

355

### 4.1.2 Inversion 1b (with a "proto Cordillera Blanca")

356

We performed a second inversion (inversion 1b) using the same constraints as for inversion 1a but with a finite initial topography in the Cordillera Blanca (see initial topography on Figure 2). Figure 5A, B presents the results of the inversion as values of the misfit function, $\mu$, displayed in parameter space. The misfit is 0.44 for the best fitting parameters model, indicating that the model can reproduce the present-day topography, total exhumation and thermochronological ages (Fig. 5C, D) as well as the best model of inversion 1a.

357

358

359

360

361

362

363

The mean uplift rate ($U_1$), the amplitude of its recent increase ($f_{U1}$) and the erosion efficiency ($K_f$) are all well constrained by our data. Similarly to inversion 1a, $U_1$ and $f_{U1}$ are correlated with higher mean uplift rate requiring a smaller recent increase in uplift rate and vice-versa. The arcuate shape of the region of minimum misfit in Figure 5B can also be interpreted by stating that the product $U_1 \times f_{U1}$ is very well constrained by the data at a value of approximately 1.0 mm/yr. The temperature at the base of the model ($T_{max}$) is not as well constrained (Fig. 5A, B) although values between 850 and 1050°C yield the lowest misfit values. These values should be interpreted as suggesting best fitting surface geothermal gradient values of 28 to 35°C/km. The inversion result suggests that the data is best explained by values of $K_f = 1.3 \times 10^{-5} m^{1-2m}.yr^{-1}$, $U_1 = 0.52$ mm/yr, $f_{U1} = 2.0$ mm/yr and $T_{max} = 945°C$. The best value for the basal temperature ($T_{max} = 945°C$) corresponds to a geothermal gradient of ~32°C/km.

364

365

366

367

368

369

370

371

372

373

374
Figure 6 presents the topographic evolution and the uplift rates predicted in the Cordillera Blanca using the best fitting parameters. Predicted uplift rates obtained at the end of the model are higher than the imposed uplift rate due to isostasy (Figs. 5, 6). Considering the low elastic thickness of the lithosphere used in this model ($T_e = 3$ km), according to existing estimate of the elastic thickness of the lithosphere in the Cordillera Blanca region [0-10 km; Perez-Gussinye et al., 2009], the flexural rebound is substantial. As explained in the method section, if we chose to ignore the flexural rebound due to erosion, the uplift rates would correspond to the vertical displacement (i.e. tectonic forcing) imposed on the fault ($U_1 = 0.52$ mm/yr and $U_1 \times f_{U1} = 1.0$ mm/yr) but if we consider the flexural rebound associated with the erosion, the uplift rates significantly increase with time (up to ~2 mm/yr; Fig. 6). The predicted uplift rates at the end of the model are consistent with the Quaternary vertical slip rates on the CBNF ranging from 5.1 ± 0.8 mm/yr to 0.6 ± 0.2 mm/yr [Schwartz, 1988; Siame et al., 2006; Margirier et al., 2017] and with $^{10}$Be catchment-wide erosion rates in the Cordillera Blanca [mean of 0.7 mm/yr; Hodson, 2012].

### 4.1.3 Comparison between the results of the two inversions

Best fitting parameters for inversion 1b are in the same order of magnitude as for inversion 1a. However, we obtained a larger value for the product $U_1 \times f_{U1}$ for the inversion 1a, indicating larger total uplift in the Cordillera Blanca since 5.4 Ma. This can be easily explained by the lower initial elevation in the inversion without an initial topography; thus to reach the present day elevation of the Cordillera Blanca, the total uplift has to be larger.

The most striking difference between these two inversions is that without a "proto Cordillera Blanca" (inversion 1a), even if the predicted elevations are well correlated with the observed topographic bulge in the Cordillera Blanca, the model still does not
reproduce the present-day location of the main drainage divide (Fig. 7). Notably, for the model without initial relief in the Cordillera Blanca the drainage divide is located close to the CBNF, which is not the case for the present-day topography. The existence of an early relief, the “proto Cordillera Blanca”, probably controls the early drainage network and localizes the drainage divide farther east in the present day Cordillera Blanca. Our results suggest that the presence of an initial relief is needed to localize the drainage divide at its present position (eastern part of the Cordillera Blanca), although we cannot preclude a different scenario.

409

## 4.2 Role of the granite

In the Cordillera Blanca the granite reached the surface and starts to be eroded at ~3 Ma [Bonnot, 1984]. This enables us to explore the influence of eroding denser rocks on the evolution of topography and uplift rates on a Ma time scale. We performed additional inversions (2a and 2b) of the topography, the total exhumation and thermochronological ages to assess the role of the intrusion in the recent increase of uplift rates in the Cordillera Blanca and in high topography building. In these two inversions the initial topography includes a proto Cordillera Blanca and we fixed $f_{U1} = 1$ (no uplift increase after 2 Ma) in order to test if the increase of the uplift rate due to the batholith erosion could explain by itself the increase of erosion rates observed by [Margirier et al., 2016]. The inversion 2a didn't include an intrusion. Then, the inversion 2b included an intrusion that is characterized by an anomalously high density (400 kg/m³ heavier than the surrounding rocks). In these inversions the parameters that we aimed to constrain were $K_f$, $T_e$, $U_1$, and $T_{max}$. Figure 8 presents the results of the inversions 2a and 2b as values of the misfit function, $\mu$, displayed in parameter space.

411

412

413

414

415

416

417

418

419

420

421

422

423

424

425

### 4.2.1 Inversion 2a (without contrasting density intrusion)
The misfit of inversion 2a is 0.42 for the best fitting parameters model, indicating that the model can reproduce the present-day topography, total exhumation and thermochronological ages (Fig. 8D). In this inversion without intrusion, the mean uplift rate ($U_1$) and the erodability ($K_f$) are also well constrained whereas the temperature at the base of the model ($T_{max}$) and the elastic thickness of the crust ($T_e$) are not as well constrained (Fig. 8A, B, C), although values between 750 and 1050°C and 3 and 7 km respectively yield the lowest misfit values. These values should be interpreted as suggesting best fitting surface geothermal gradient values of 25 to 35 °C/km. The data is best explained by very similar values of $K_f = 1.2 \times 10^{-5} m^{-2m}.yr^{-1}$, $T_e = 5.6$ km, $U_1 = 0.86$ mm/yr, and $T_{max} = 1001°C$. The best value for the basal temperature ($T_{max} = 1001°C$) corresponds to a geothermal gradient of ~33°C/km.

### 4.2.2 Inversion 2b (with contrasting density intrusion)

The lowest misfit is 0.42 for inversion 2b, suggesting that the model can also reproduce the present-day topography, total exhumation and thermochronological ages (Fig. 8). In this inversion that includes a dense intrusion, the mean uplift rate ($U_1$) and the erosion efficiency ($K_f$) are all well constrained (Fig. 8E, F, G). The temperature at the base of the model ($T_{max}$) and the elastic thickness of the crust ($T_e$) are not as well constrained although values between 750 and 1050°C and 3 and 5 km respectively yield the lowest misfit values (Fig. 8E, F). These $T_{max}$ values should be interpreted as suggesting best fitting surface geothermal gradient values of 25 to 35°C/km. $U_1$ and $T_e$ are clearly correlated with smaller elastic thickness of the crust requiring smaller uplift rate and vice-versa. The inversion result suggests that the data is best explained by values of $K_f = 1.2 \times 10^{-5} m^{-2m}.yr^{-1}$, $T_e = 3.3$ km, $U_1 = 0.69$ mm/yr, and $T_{max} = 954°C$. The best value for the basal temperature ($T_{max} = 954°C$) corresponds to a geothermal gradient of ~32°C/km.
453

454

### 4.2.3 Comparison between the results of the two inversions

455

456

457

458

459

460

461

462

463

464

465

466

467

468

469

470

471

472

473

474

475

## 5. Discussion

476

### 5.1 Paleogeography
The substantial surface uplift and resulting erosion in the Cordillera Blanca make it difficult to study the paleogeography in this area based on remnants geomorphological features. However, the drainage network geometry and the drainage divide location provide information on the topography of the Cordillera Blanca in the past and its evolution. Notably, the Cordillera Blanca drainage divide is located in the eastern part of the range whereas both the higher peaks and higher uplift rates are located in the western flank of the Cordillera Blanca [Margirier et al., 2016]. But the precipitation gradient, with higher mean precipitation on the eastern flank of the Cordillera Blanca, should favor the location of the drainage divide further west. Our results suggest that the location of the drainage divide in the eastern part of the Cordillera Blanca is controlled by an initial topography in the Cordillera Blanca area and an inherited drainage network. Based on these results we propose that a proto-Cordillera Blanca existed before the CBNF initiation ~5 Myr ago. The presence of such a relief is in agreement with Giovanni et al. [2010] and Wise and Noble [2003], which already proposed that the Callejón de Huaylas was a topographic depression during the Late Miocene. Wise and Noble [2003] suggested that at that time, before the Cordillera Blanca exhumation, the Cordillera Negra corresponded to the drainage divide. However, the presence of a proto-Cordillera Blanca calls into question the position of the Andean drainage divide during the Miocene.

496

## 5.2 The Cordillera Blanca uplift: role of granite erosion and flexural rebound

Our results show a difference in imposed uplift rates for the best fitting models of the inversions 2a and 2b that do not or do include a granite ($U_1=0.86$ mm/yr versus $U_1= 0.69$ mm/yr; Fig. 8). In addition, our models do not take into account the presence of the CBNF for the isostatic response calculation despite the existence of a 200 km-long crustal fault [Dalmayrac and Molnar, 1981; McNulty and Farber, 2002]. If, to reflect the
regional geology, we had used a broken plate approximation to calculate the isostasy
contribution to uplift and erosion, the isostatic response would, potentially, have been
more important. In addition, even when considering the best fitting value for T<sub>e</sub> the role
of the granite seems to be small; if we consider the range of acceptable values for T<sub>e</sub>, we
notice that smaller values of T<sub>e</sub> would increase the importance of isostasy, which, in
turn, would require smaller valued of imposed uplift, U<sub>1</sub> (Fig. 8E). For example, if we
consider the best fitting models or the models involving a 3 km elastic thickness which
also show low misfit values, the isostasy contributes c.a. 20% of the total rock uplift (Fig.
8B, F). Considering that for the best fitting parameters model of the inversion 2b, the
granite reaches the surface approximately 3 Ma before the end of the model run, our
results indicate that the isostatic effect of eroding a denser rock mass is not negligible
even on a < 5 Ma time scale. Finally, we note that the inversion that includes a denser
granite (inversion 2b) converges toward smaller values for the imposed uplift rate, U<sub>1</sub>,
than the inversion which doesn't include an intrusion (inversion 2a) implying that
erosion of the Cordillera Blanca dense intrusion may substantially contribute to the
present-day higher uplift rates as proposed by Braun et al. [2014].

Interestingly, our results suggest that the recent increase in exhumation rate
documented by *Margirier et al.* [2016] might not be the result of a change in tectonic
since 2 Ma or climate forcing, but the result of the progressive flexural isostatic response
of the area to erosion. At the onset of uplift, slopes are relatively low, erosion is limited
to the sides of the uplifting region and the weight of uplifted material acts as a negative
feedback to further uplift; as the uplifting region becomes affected by surface processes
(fluvial erosion in our model, but more likely glacial erosion in the Cordillera Blanca) the
resulting erosional unloading reduces the negative feedback and uplift and erosion rates
increase until a steady-state situation is reached between uplift, erosion and isostasy.
528
529

The time scale over which this steady-state balance is reached must be of the order of a
few million years, at least [Whipple and Meade, 2006].

530
531
532
533
534
535
536
537
538
539
540
541
542

From the time when the Cordillera Blanca summits reached high elevation, the
Cordillera Blanca itself has acted as an orographic barrier for moisture carried from the
Amazon basin [Montgomery et al., 2001; Bookhagen and Strecker, 2008]. The induced
rain shadow included in our model prevents erosion of the Cordillera Negra, and results
in higher erosion rates in the Cordillera Blanca [i.e. Montgomery and Brandon, 2002] in
the recent past. It is also likely that Quaternary glaciations [~1.5-0 Ma; Farber et al.,
2005; Smith et al., 2005] have increased erosion rates through the formation of deeply
incised U-shaped valleys, as evidenced by Montgomery [2002] for the Olympic
Mountains on the Pacific coast of the North America and suggested by Margirier et al.
[2016] for the Cordillera Blanca. Finally, the recent glacial retreat [since ~21 ka; Seltzer
et al., 2002; Farber et al., 2005; Smith et al., 2005] could have induced a flexural rebound
that both increased uplift rates in the Cordillera Blanca and slip-rates on the CBNF [e.g.,
Hetzel and Hampel, 2005; Hampel et al., 2007].

543
544

5.3 *Models for the CBNF*

545
546
547
548
549
550
551
552
553

Two models have been proposed to explained extension on the CBNF. *Dalmayrac and Molnar* [1981] proposed that the fault is the result of a gravitational collapse of the thickened crust, whereas *McNulty and Farber* [2002] suggested that the subduction of the buoyant Nazca Ridge below the Cordillera Blanca drove the footwall uplift. However, there is still no consensus on a model to explain the CBNF and the processes driving the Cordillera Blanca uplift. *Margirier et al.* [2015] already provided evidence for exhumation in the Cordillera Negra since 15 Ma, which is likely to be related to a regional surface uplift of the Western Andes [*McLaughlin, 1924; Farrar and Noble, 1976, Myers, 1976; Wipf, 2006; Hoorn et al., 2010*]. The 1300 m-thick sedimentary filling of the
Callejón de Huaylas suggests only low subsidence of the basin since ~5.4 Ma [Bonnot, 1984; Giovanni et al., 2010] even if the vertical displacement on the CBNF is estimated to be ~7 km [Margirier et al., 2016]. In addition, our models of the landscape evolution provide new constraints on uplift rates in the Cordillera Blanca region. We suggest that the Cordillera Blanca has been uplifted in relation to the Cordillera Negra during the past 7 Ma. Observations clearly suggest that the displacement on the CBNF is mostly the result of the Cordillera Blanca uplift rather than subsidence of the hanging wall. The rock uplift evidenced in this study for the Cordillera Blanca (7-0 Ma) and proposed in the Cordillera Negra [15-0 Ma; Margirier et al., 2015] is not compatible with the collapse model proposed by Dalmayrac and Molnar [1981]. Moreover, new reconstructions of the timing (15-11 Ma) and location of the initial Nazca Ridge subduction [Hampel, 2002; Rosenbaum et al., 2005; Antonijevic et al., 2015] are not compatible with the timing and location of normal faulting in the Cordillera Blanca [5.4 Ma; Bonnot, 1984; Giovanni et al., 2010]. Indeed, even if the Nazca Ridge subduction initiate at 11°S, the latitude of the Cordillera Blanca between 15 and 11 Ma, at ~5 Ma, the time of initiation of normal faulting the ridge is farther south (~13°S) and at the present day the ridge is at 15°S whereas the CBNF is still active. This suggests that the Nazca Ridge did not trigger the Cordillera Blanca uplift as proposed by McNulty and Farber [2002]. Our landscape evolution models do not permit to test the role of geodynamic processes such as slab flattening on the CBNF initiation and Cordillera Blanca uplift. However, the subduction geometry and change of the mechanical coupling on the subduction interface may have had an influence on the tectonic regime of the upper plate as suggested by Margirier et al. [2017] for the Cordillera Blanca region. After the initiation of normal faulting, when the Cordillera Blanca batholith started to be eroded, since ~3 Ma, the erosion of denser rocks is likely to have trigger an increase of the uplift rates along the CBNF.

579
580 **6. Conclusions**

581 Our study provides new constraints on the erosion efficiency, elastic thickness of the
582 lithosphere, temperature gradient in the crust and uplift rates in the Andes of northern
583 Peru. The absolute rock uplift rates obtained at the end of the models for the Cordillera
584 Blanca (ranging from 1.5 to 2.5 mm/yr) are coherent with Quaternary slip rates
585 documented on the CBNF [5.1 ± 0.8 mm/yr to 0.6 ± 0.2 mm/yr; Schwartz, 1988; Siame et
586 al., 2006; Margirier et al., 2017; Gérard et al., in prep]. Our results suggest an
587 acceleration of rock uplift in the Cordillera Blanca at 2 Ma in agreement with the
588 increase in erosion rate evidenced by Margirier et al. [2016]. We show, however, that
589 such an increase in erosion rate may be the result of erosional unloading and isostasy
590 rather than being related to a change in either tectonic or climatic forcing. We also show
591 that this acceleration may have been amplified by the unroofing of a dense granitic
592 intrusion.

593 Based on the present-day drainage divide location in the Cordillera Blanca and the
594 results of our modeling, we propose that the Callejón de Huaylas was already a
595 depression before the Cordillera Blanca batholith emplacement and that a proto
596 Cordillera Blanca already existed at that time.

597 Finally, in the light of our modeling of the landscape evolution, we show that the two
598 models previously published for the CBNF [e.g., gravitational collapse of the thickened
599 crust / footwall uplift due to the subduction of the Nazca Ridge; Dalmayrac and Molnar,
600 1981; McNulty and Farber, 2002] are inconsistent with the thermochronological data,
601 geodynamic and geologic context. Further investigations are needed to assess if the
602 Quaternary extensional tectonics have been triggered by flat slab subduction, which may
603 also have triggered at the same time the uplift of the Cordillera Blanca. We suggest,
604 however, that, in this context, the construction of high-relief topography was strongly
influenced by the isostatic flexural rebound and, to some degree, by the exhumation of
the dense Cordillera Blanca granite.

607

608 Acknowledgements

609 We thank Emily Richards for proofreading the manuscript and English improvements,
610 Jessica Stanley and Benoît Bovy for their help with the inversions on the cluster. We
611 acknowledge the work of the editor and of our two anonymous reviewers for their
612 critical and helpful reviews. Last but not least, Audrey Margirier acknowledges
613 Blockzone team for the nice boulders in Potsdam.

614

615 References

616 Antonijevic, S. K., L. S. Wagner, A. Kumar, S. L. Beck, M. D. Long, G. Zandt, H. Tavera, and
617 C. Condori (2015), The role of ridges in the formation and longevity of flat slabs,
618 *Nature*, **524**(7564), 212–215, doi:10.1038/nature14648.

619 Barazangi, M., and B. L. Isacks (1976), Spatial distribution of earthquakes and
620 subduction of the Nazca plate beneath South America, *Geology*, **4**(11), 686–692.

621 Beckinsale, R. D., A. W. Sanchez-Fernandez, M. Brook, E. J. Cobbing, W. P. Taylor, and N. B.
622 Moore (1985), Rb-Sr whole rock isochron and K-Ar determination for the Coastal
623 Batholith of Peru, in *Magmatism at a Plate Edge: The Peruvian Andes*, edited by W. S.
624 Pitcher, M. P. Atherton, E. J. Cobbing, and R. D. Beckinsale, pp. 177–202, Blackie
625 Halstead press, Glasgow.

626 Bonnot, D. (1984), Néotectonique et tectonique active de la Cordillère Blanche et du
627 Callejon de Huaylas (Andes nord-péruviennes), *Thèse présentée pour obtenir le grade
628 de docteur, Université de Paris-Sud, Centre d'Orsay*, 1-202.

629 Bookhagen, B., and M. R. Strecker (2008), Orographic barriers, high-resolution TRMM
630 rainfall, and relief variations along the eastern Andes, *Geophys. Res. Lett.*, **35**(6),
631 L06403, doi:10.1029/2007GL032011.

632 Braun, J., P. Van der Beek, and G. Batt (2006), Quantitative thermochronology, numerical
633 methods for the interpretation of thermochronological data, *Cambridge University
634 Press*, Cambridge, UK, 258 p.

635 Braun, J., and S. D. Willett (2013), A very efficient O(n), implicit and parallel method to
636 solve the stream power equation governing fluvial incision and landscape evolution,
637 *Geomorphology*, **180-181**, 170–179, doi:10.1016/j geomorph.2012.10.008.

638 Braun, J., T. Simon-Labric, K. E. Murray, and P. W. Reiners (2014), Topographic relief
639 driven by variations in surface rock density, *Nature Geoscience*, 7, 534–540,
640 doi:10.1038/ngeo2171.

641 Cederbom, C. E., H. D. Sinclair, and F. Schlunegger (2004), Climate-induced rebound and
642 exhumation of the European Alps, *Geology*, 32(8), 709–712, doi:10.1130/G20491.1.

643 Cobbing, J., W. Pitcher, J. Baldock, W. Taylor, W. McCourt, and N. J. Snelling (1981),
644 Estudio geológico de la Cordillera Occidental del norte del Perú, *Instituto Geológico
645 Minero y Metalurgico, Serie D. Estudios Especiales*, 10(D), 1–252.

646 Croissant, T., and J. Braun (2014), Constraining the stream power law: a novel approach
647 combining a landscape evolution model and an inversion method, *Earth Surf.
648 Dynam.*, 2(1), 155–166, doi:10.5194/esurf-2-155-2014-supplement.

649 Dalmayrac, B., and P. Molnar (1981), Parallel thrust and normal faulting in peru and
650 constraints on the state of stress, *Earth and Planetary Sciences Letters*, 55, 473–481.

651 Eakin, C. M., C. Lithgow-Bertelloni, and F. M. Dávila (2014), Influence of Peruvian flat-
652 subduction dynamics on the evolution of western Amazonia, *Earth and Planetary
653 Science letters*, 404(C), 250–260, doi:10.1016/j.epsl.2014.07.027.

654 Farrar, E., and D. C. Noble (1976), Timing of late Tertiary deformation in the Andes of
655 Peru, *Geological Society of America Bulletin*, 87(9), 1247–1250.

656 Farber, D. L., G. S. Hancock, R. C. Finkel, and D. T. Rodbell (2005), The age and extent of
657 tropical alpine glaciation in the Cordillera Blanca, Peru, *J. Quaternary Sci.*, 20(7-8),
658 759–776, doi:10.1002/jqs.994.

659 Gérard, B., L. Audin, D. Farber, S. Hall, S. Zérathe, X. Robert, A. Margirier, J. Carcaillet, S.
660 Schwartz, Northward increase of Quaternary slip rates along the Cordillera Blanca
661 normal fault, Peruvian Andes, in prep.

662 Giovanni, M. K. (2007), Tectonic and Thermal Evolution of the Cordillera Blanca
663 Detachment System, Peruvian Andes: Implication for Normal Faulting in a
664 Contractionnal Orogen, Ph.D. Thesis, 1–255 pp. University of California, Los Angeles.

665 Giovanni, M. K., B. K. Horton, C. N. Garzione, B. McNulty, and M. Grove (2010),
666 Extensional basin evolution in the Cordillera Blanca, Peru: Stratigraphic and isotopic
667 records of detachment faulting and orogenic collapse in the Andean hinterland,
668 *Tectonics*, 29(6), TC6007, doi:10.1029/2010TC002666.

669 Hampel, A. (2002), The migration history of the Nazca Ridge along the Peruvian active
670 margin: a re-evaluation, *Earth and Planetary Science letters*, 203(2), 665–679.

671 Hampel, A., R. Hetzel, and A. L. Densmore (2007), Postglacial slip-rate increase on the
672 Teton normal fault, northern Basin and Range Province, caused by melting of the
673 Yellowstone ice cap and deglaciation of the Teton Range?, *Geology*, 35(12), 1107-
674 1110, doi:10.1130/G24093A.1.

675 Henry, S. G., and H. N. Pollack (1988), Terrestrial heat flow above the Andean subduction
676 zone in Bolivia and Peru, *Journal of Geophysical Research: Solid Earth*, 93(B12),
677 15153–15162.
678 Hetzel, R., and A. Hampel (2005), Slip rate variations on normal faults during glacial-interglacial changes in surface loads, *Nature*, **435**, 81-84.

680 Hodson, K. R. (2012), Morphology, exhumation, and Holocene erosion rates from a tropical glaciated mountain range: the Cordillera Blanca, Peru, Department of Earth and Planetary Sciences, Unpublished Thesis, McGill University, Montreal, 1-94.

683 Hoorn, C. et al. (2010), *Amazonia Through Time: Andean Uplift, Climate Change, Landscape Evolution, and Biodiversity*, *Science*, **330**(6006), 927-931, doi:10.1126/science.1194585.

686 Margirier, A., X. Robert, L. Audin, C. Gautheron, M. Bernet, S. Hall, and T. Simon-Labric (2015), Slab flattening, magmatism and surface uplift in the Cordillera Occidental (northern Peru), *Geology*, **1**-16, doi:10.1130/G37061.1.

689 Margirier, A., L. Audin, X. Robert, F. Herman, J. Ganne and S. Schwartz (2016), Time and mode of exhumation of the Cordillera Blanca batholith (Peruvian Andes), *Journal of Geophysical Research: Solid Earth*, **121**, 6235-6249, doi:10.1002/2016JB013055.

692 Margirier, A., L. Audin, X. Robert, A. Pêcher, and S. Schwartz (2017), Stress field evolution above the Peruvian flat-slab (Cordillera Blanca, northern Peru), *Journal of South American Earth Sciences*, **77**, 58-69, doi:10.1016/j.jsames.2017.04.015.

695 McLaughlin, D. H. (1924), Geology and physiography of the Peruvian Cordillera, Departments of Junin and Lima, *Bulletin of the geological Society of America*, **35**, 591-632.

698 McNulty, B. A., and D. L. Farber (2002), Active detachment faulting above the Peruvian flat slab, *Geology*, **30**(6), 567-570.

700 Montario, M. J. (2001), Exhumation of the Cordillera Blanca, Northern Peru, based on apatite fission track analysis, Department of Geology, Unpublished Thesis, Union College, Schenectady, New York, 1-12.

703 Montgomery, D. R. (1994), Valley incision and the uplift of mountain peaks, *Journal of Geophysical Research: Solid Earth*, **99**(B7), 13913-13921.

705 Montgomery, D. R. (2002), Valley formation by fluvial and glacial erosion, *Geology*, **30**(11), 1047-1050.

707 Montgomery, D. R., and M. T. Brandon (2002), Topographic controls on erosion rates in tectonically active mountain ranges, *Earth and Planetary Science letters*, **201**, 481-489.

710 Montgomery, D. R., G. Balco, and S. D. Willett (2001), Climate, tectonics, and the morphology of the Andes, *Geology*, **29**(7), 579-582.

712 Mukasa, S. B. (1984), Comparative Pb isotope systematics and zircon U-Pb geochronology for the Coastal, San Nicolás and Cordillera Blanca Batholiths, Peru, Ph.D. Thesis, 1-362 pp., University of California, Santa Barbara.

715 Myers, J. S. (1976), Erosion surfaces and ignimbrite eruption, measures of Andean uplift in northern Peru, *Geol. J.*, **11**(1), 29-44.
717 Perez-Gussinye, M., C. J. Swain, J. F. Kirby, and A. R. Lowry (2009), Spatial variations of the effective elastic thickness, Te, using multitaper spectral estimation and wavelet methods: Examples from synthetic data and application to South America, *Geochem. Geophys. Geosyst.*, 10(4), 1-10, doi:10.1029/2008GC002229.

721 Petford, N., and M. P. Atherton (1992), Granitoid emplacement and deformation along a major crustal lineament: the Cordillera Blanca, Peru, *Tectonophysics*, 205(1), 171-185.

724 Rosenbaum, G., D. Giles, M. Saxon, P. G. Betts, R. F. Weinberg, and C. Duboz (2005), Subduction of the Nazca Ridge and the Inca Plateau: Insights into the formation of ore deposits in Peru, *Earth and Planetary Science letters*, 239(1-2), 18-32, doi:10.1016/j.epsl.2005.08.003.

728 Sambridge, M., (1999a), Geophysical inversion with a neighbourhood algorithm — I. Searching a parameter space, *Geophys. J. Int.*, 138, 479-494.

730 Sambridge, M., (1999b), Geophysical inversion with a neighbourhood algorithm — II. Appraising the ensemble, *Geophys. J. Int.*, 138, 727-746.

732 Schwartz, D. P. (1988), Paleoseismicity and neotectonics of the Cordillera Blanca fault zone, northern Peruvian Andes, *Journal of Geophysical Research, Solid Earth*, 93(B5), 4712-4730.

735 Seltzer, G. O., D. T. Rodbell, P. A. Baker, S. C. Fritz, P. M. Tapia, H. D. Rowe, and R. B. Dunbar (2002), Early warming of tropical South America at the last glacial-interglacial transition, *Science*, 296, 1685-1686.

738 Sharma, P. (1997), Environmental and Engineering Geophysics, *Cambridge Univ. Press*, New York.

740 Siame, L. L., M. Sébrier, O. Bellier, and D. Bourles (2006), Can cosmic ray exposure dating reveal the normal faulting activity of the Cordillera Blanca Fault, Peru? *Revista de la Asociación Geológica Argentina*, 61(4), 536-544.

743 Small, E. E., and R. S. Anderson (1995), Geomorphically driven late Cenozoic rock uplift in the Sierra Nevada, California, *Science*, 270, 277-278.

745 Smith, J. A., R. C. Finkel, D. L. Farber, D. T. Rodbell, and G. O. Seltzer (2005), Moraine preservation and boulder erosion in the tropical Andes: interpreting old surface exposure ages in glaciated valleys, *J. Quaternary Sci.*, 20(7-8), 735-758, doi:10.1002/jqs.981.

749 Stern, T. A., A. K. Baxter, and P. J. Barrett (2005), Isostatic rebound due to glacial erosion within the Transantarctic Mountains, *Geology*, 33(3), 221-224, doi:10.1130/G21068.1.

752 Turcotte, D. L., and G. Schubert (2002), *Geodynamics*, Cambridge University Press, 456 pp.

754 Whipple, K. and B. Meade (2006), Orogen response to changes in climatic and tectonic forcing, *Earth and Planetary Science letters*, 243, 218-228.
756

Wipf, M. (2006), Evolution of the Western Cordillera and Coastal Margin of Peru:

757

Evidence from low-temperature Thermochronology and Geomorphology,

758

Department of Earth Sciences, PhD Thesis, Swiss Federal Institute of Technology

759

Zürich, Zürich, 1–163.

760

Wise, J. M., and D. C. Noble (2003), Geomorphic evolution of the Cordillera Blanca,

761

Northern Peru, *Boletin de la sociedad Geologica del Peru*, 96, 1–21.

762

763

764

765

766

767

<figure><img src="image_3.png" /><figcaption>Figure 1: A) Topographic map of the northern Peru showing active tectonic features (Neotectonic Open Database, neotec-opendata.com), location of the Cordillera Blanca (CB, the batholith is highlighted in purple), the Cordillera Negra (CN) and the Callejón de Huaylas (yellow). Inset shows map location within the South America. B) E-W cross-section of the Andes at the latitude of the Cordillera Blanca both the topography (black line) and the rainfall (blue line) are represented [*modified from Bookhagen and Strecker, 2008*]. On the topographic cross-section the Cordillera Blanca batholith is represented in purple and a black arrow points the drainage divide location. C) Photograph of the CBNF showing the 1 km high triangular facets along the active fault scarp.</figcaption></figure>

768

769

770

771

772

773

774

775

776
<figure><img src="image_4.png" /><figcaption>Figure 2: Parameterization of the models.</figcaption></figure>

777

778

779

780

781
**Inversion 1a** *without a «proto Cordillera Blanca»*

<figure><img src="image_5.png" /></figure>

Figure 3: Results of the NA inversion 1a. A, B) Scatter plots colored by likelihood values for erosion efficiency ($K_f$), basal temperature ($T_{max}$), uplift rate ($U_1$) and factor of uplift rate increase for the last 2 Ma ($f_{U1}$). The most likelihood model is indicated with a red star. C) Comparison between observed ages and synthetic ages predicted by the best-fit model for apatite fission-track ages and apatite (U-Th)/He ages. D) Observed ages and synthetic ages vs. distance to the Cordillera Blanca normal fault.

**Inversion 1a**

<figure><img src="image_6.png" /></figure>

Figure 4: Evolution of the topography and rock uplift rate for the best fitting parameter model of inversion 1a. A-D) Topography, E-H) Rock uplift rate.
**Inversion 1b** with a «proto Cordillera Blanca»

<figure><img src="image_7.png" /></figure>

Figure 5: Results of the NA inversion 1b. A, B) Scatter plots colored by likelihood values for erosion efficiency ($K_f$), basal temperature ($T_{max}$), uplift rate ($U_1$) and factor of uplift rate increase for the last 2 Ma ($f_{U1}$). The most likelihood model is indicated with a red star. C) Comparison between observed ages and synthetic ages predicted by the best-fit model for apatite fission-track ages and apatite (U-Th)/He ages. D) Observed ages and synthetic ages vs. distance to the Cordillera Blanca normal fault.

**Inversion 1b**

<figure><img src="image_8.png" /></figure>

Figure 6: Evolution of the topography and rock uplift rate for the best fitting parameter model of inversion 1b. A-D) Predicted topography, E-H) Map of the predicted rock uplift rate.
<figure><img src="image_9.png" /><figcaption>Figure 7: DEM and swath profiles through the Cordillera Blanca region and modeled topography. A) DEM of the Cordillera Blanca region. B) DEM of the best model of inversion 1a, which involves a fault, a granitic batholith, a flexural rebound and a "proto" Cordillera Negra. C) DEM of the best model of inversion 1b with both a "proto" Cordillera Negra and a "proto" Cordillera Blanca. D-F) Swath profiles through the Cordillera Blanca and from the two models showing elevation (mean, maximum and minimum), the drainage divide position (green arrow), the CBNF location (black arrow) and the maximum elevation (vertical dot line).</figcaption></figure>
<figure><img src="image_10.png" /></figure>

Figure 8: Results of the NA inversions 2a and 2b as scatter plots colored by likelihood values (from blue for high misfit value or poor fit to the data, to red for low misfit value or good fit to the data) for erosion efficiency ($K_f$), basal temperature ($T_{max}$), uplift rate ($U_1$) and elastic thickness of the crust ($T_e$). The most likelihood models are indicated with purple stars (inversion 2a) and red stars (inversion 2b). A-C) Scatter plots and D) comparison between observed thermochronological ages and synthetic ages predicted by the best-fit model for the inversion 2a (without intrusion). E-G) Scatter plots for inversion 2b (with an intrusion). H) Comparison between observed thermochronological ages and synthetic ages predicted by the best-fit model for the inversion 2b.

<figure><img src="image_11.png" /></figure>

Figure 9: Map of the topography at the end of the model A) without granite and B) with a dense granite. Map of the uplift rates at the end of the model C) without granite and D) with a dense granite. The oval dotted line indicates the location of the granite.