Articles | Volume 11, issue 2
Research article
21 Mar 2023
Research article |  | 21 Mar 2023

River incision, 10Be production and transport in a source-to-sink sediment system (Var catchment, SW Alps)

Carole Petit, Tristan Salles, Vincent Godard, Yann Rolland, and Laurence Audin

Detrital 10Be from continental river sands or submarine sediments has been extensively used to determine the average long-term denudation rates of terrestrial catchments, based on the assumption that the rate of cosmogenic nuclide production by the interaction of source rocks with cosmic radiation balances out the loss of these nuclides by surface denudation. However, the 10Be signal recorded in sediments may be affected at the source by the response time of mountainous catchments to high-frequency forcings. In addition, transient sediment storage in piedmonts, alluvial plains and lakes or near the coast may also induce a difference between the erosive signal and its record in the sedimentary sink. Consequently, a significant part of the signal recorded in shallow-water sediments can be lost, as deep marine sediments may simultaneously record a signal coming from newly eroded source rocks along with one coming from the destabilization of previously deposited sediments.

In this paper, we use the landscape evolution model Badlands to simulate erosion, deposition and detrital 10Be transfer from a source-to-sink sedimentary system (the Var River catchment, southern French Alps) over the last 100 kyr. We first compare model-based denudation rates with the ones that would be extracted from the 10Be record of local continental sediments (equivalent to river sands) and from sediments deposited offshore over time in order to examine if this record provides an accurate estimate of continental denudation rates. Then, we examine which conditions (precipitation rate, flexure, ice cover) satisfy published measured river incision rates and 10Be concentration in submarine sediments.

Our results, based on the Var catchment cosmic ray exposure dating and modelling indicate that, while river sands do accurately estimate the average denudation rate of continental catchments, this is much less the case for deep submarine sediments. We find that deep-sea sediments have a different and often much smoother 10Be signature than continental ones and record a significant time lag with respect to imposed precipitation rate changes, representing the geomorphological response of the margin. A model which allows us to fit both measured 10Be concentration in marine sediments and river incision rates on land involves an increase in precipitation rates from 0.3 to 0.7 m yr−1 after 20 ka, suggesting more intense precipitation starting at the end of the Last Glacial Maximum.

1 Introduction

Sedimentary deposits are important archives of the tectonic and climatic history of continents: for instance, the geometry, grain size, mineralogy and geochemical signature of deposits are impacted by changes in environmental conditions (e.g. relative sea level changes and precipitation over geologic times as well as human activities during the Anthropocene) (Syvitski et al.2022). Provided there are good enough estimates of the transfer function between these sedimentary records and their external forcing, they can be reliable tools for reconstructing climatic cycles, subsidence curves or monsoon onset, for instance (Bentley et al.2016; Li et al.2016; Liu et al.2016; Wan et al.2006). However, depending on the timescales considered, the signatures in submarine sediments of some of these external (i.e. climatic or tectonic) forcings affecting aerial catchments depend on a myriad of processes which still remain difficult to extract from the deep-sea record. To this end, one would need to evaluate not only how the eroded source responds to specific forcing but also for how long and where detrital sediments are temporarily stored and when they are re-injected into the system and eventually reach their sink.

Concerning the sediment source, mountainous catchments may not be very sensitive to high-frequency forcing, and the response time of these catchments may already affect the signal recorded in locally produced sediments (e.g. Armitage et al.2013; Godard and Tucker2021; Goren2016; Jerolmack and Paola2010). Second, transient sediment storage in piedmonts, alluvial plains and lakes or near the coast may induce a large time lag between the external signal and its record in the sedimentary sink (e.g. Blöthe and Korup2013; Clift and Giosan2014; Malatesta et al.2018; Phillips and Slattery2006; Romans et al.2016). Depending on considered timescales, the erosive signal itself can be completely buffered by this process (see a complete review in Romans et al.2016). Finally, submarine sediments can be reworked by gravitational processes, especially during sea level falls (Phillips and Slattery2006). As a consequence, a significant part of the signal recorded in shallow-water sediments can be lost, whereas deep marine sediments may simultaneously record a signal coming from newly eroded source rocks and another one coming from the destabilization of previously deposited sediments. In addition, relative sea level variations may affect the connectivity between aerial rivers and submarine canyons, therefore, from time to time, limiting the efficiency of sediment transport in the offshore domain (Fryirs et al.2007).

Detrital terrestrial cosmogenic nuclide (TCN; mostly 10Be) concentrations from continental river sands or submarine sediments have been extensively used to determine the average long-term denudation rates of aerial catchments, provided enough quartz-bearing rocks outcrop at the surface to give a representative sampling of the whole catchment denudation (e.g. Bierman and Steig1996; von Blanckenburg2005; Lupker et al.2012; Mandal et al.2015; Siame et al.2011; Vanacker et al.2007). Denudation rate estimates from 10Be concentration in quartz-rich sediments are often based on the assumption that the rate of TCN production by the interaction of source rocks with cosmic radiation balances out the loss of these elements by surface denudation (Lal1991). Denudation rates can vary in time and space, which puts this steady-state assumption into question and may lead to under- or overestimates of the true denudation rates (Bierman and Steig1996). The abundance of the target mineral (i.e. quartz in the case of 10Be) in surface rocks may also vary and has to be taken into account in order to correctly estimate the total production rate of a given catchment (Bierman and Steig1996; Safran et al.2005; Carretier et al.2015).

Moreover, sediments can be seen as an amalgamation of individual grains of different sizes and with different histories: alluvial terraces contain various clasts, the TCN concentration of which is linked to the erosion rate of their catchment (inheritance) and to posterior TCN production (Repka et al.1997). The distribution of TCN concentration in individual grains depends on the geomorphic processes acting at the source and on the post-erosion TCN production (Codilean et al.2010). Numerical models of individual grains' journeys have shown that some grains may have very long residence times in the piedmont (i.e. of the order of 100 kyr) and may therefore be exported to their final depositional area long after they have been produced by bedrock erosion (Carretier et al.2020). Landslides may also significantly affect denudation rate estimates from TCN if the catchment area is small (Yanites et al.2009). Hence, understanding how detrital 10Be concentrations recorded in submarine sedimentological archives reflect denudation rates at the time of their deposition requires quantifying the following: (i) how, how fast and where 10Be is produced; (ii) how 10Be concentration in the sediments produced is representative of average catchment denudation rates at any given spatiotemporal scale (Zerathe et al.2022); and (iii) how long it has taken for sediments, once they are produced, to reach the sedimentary sink where they have potentially been sampled.

In this paper, we adapt the surface process model Badlands (Salles2016) to simulate erosion, deposition and detrital 10Be transfer from a source-to-sink sediment system (the Var River catchment, southern French Alps, and its marine depositional system in the Mediterranean Sea) over the last 100 kyr. We first compare simulated denudation rates with the ones that would be inferred from the 10Be concentration of continental sediments (equivalent to river sands) or sediments deposited offshore at each time step in order to assess at which timescales the steady-state assumption is valid and whether 10Be record in detrital sediments provides an accurate estimate of continental denudation rates. Then, we examine which conditions (precipitation rate, flexure, ice cover) satisfy published river incision rates (Cardinal et al.2022; Petit et al.2019; Rolland et al.2017, 2020; Saillard et al.2014) and 10Be concentration in marine sediments (named hereafter 10BeMS) in this particular, small-scale source-to-sink system (Mariotti et al.2021).

2 Geomorphological and geological setting

The Var catchment in the southern French Alps is ideally suited to constrain source-to-sink processes: it is a relatively small catchment (∼2800 km2), which encompasses some of the high-altitude (∼3000 m) summits of the Alpine Mercantour crystalline massif (Fig. 1). The Var River has three main large tributaries: the Tinée and Vésubie rivers, the headwaters of which are in the crystalline massif, and the Estéron River, which flows only across the overlying Mesozoic to Cenozoic sedimentary sequence. The hydrologic regime of the Var River, dominated by flash floods, is responsible for frequent hyperpycnal flows in the submarine domain (Mulder et al.1998). The continental floodplain and shelf at the mouth of the Var River are very narrow, and most detrital sediments are deposited offshore at the foot of the Ligurian margin (Mediterranean Sea) at depths below −2000 m. The average annual discharge at the mouth of the Var River is ∼50 m3 s−1, but it can reach more than 500 m3 s−1 during floods (Mulder et al.1998). According to Mulder et al. (1998) and Syvitski et al. (2000), the suspended sediment load Cs can be obtained from the following expression:

(1) C s = a Q w b ,

where parameters a [M L−3] and b relate the sediment concentration Cs [M L−3] to mean discharge Qw [L3 T−1] at the river mouth, and have been estimated to be 7.7×10-4 and ∼1.65 for the Var River, respectively, from direct measurements of suspended sediment concentration and water discharge (Mulder et al.1998). Modern average denudation rate estimates of the Var River catchment from 10Be measurements in fluvial sediments range between 0.1 and 0.8 mm yr−1 (Mariotti et al.2019). Two other coastal rivers flow east of the Var catchment: the Paillon and Roya rivers, with much smaller drainage areas (258 and 601 km2, respectively).

Figure 1Topographic and bathymetric map and schematic outline of main quartz-bearing rock outcrops of the study area (transparent red and yellow patches). Dashed white lines indicate the contours of the Var, Paillon and Roya river catchments. Stars, open dots and triangles indicate the location of data constraints provided by the following: surface exposure dating of river-polished surfaces (after Saillard et al.2014 – S2014; Rolland et al.2017, 2020 – R2017 and R2020; Petit et al.2019 – P2019; Cardinal et al.2022 – C2022), geochemical analyses of submarine sediments (Bonneau et al.2017 – B2017) and detrital 10Be in river sands (Mariotti et al.2019 – M2019). For river-polished surfaces, the numbers refer to the following sites: 1 – Salso Moreno; 2 – Isola; 3 – Saint-Sauveur; 4 – lower Tinée; 5 – Estéron; 6 – Vésubie; 7 – Paillon; 8 – Bevera; 9 – Roya. Italic expressions refer to sample names in Mariotti et al. (2019).

A previously published paper presented detailed sedimentological and geochemical analyses of sediment cores in the sedimentary ridge located at the outlet of the Var submarine canyon (Bonneau et al.2017) (Fig. 1). These analyses revealed a larger frequency of turbidite flows and slightly larger epsilon-Nd (ϵNd) values during the Last Glacial Maximum (LGM), which is interpreted as reflecting more intense erosion, especially in the crystalline massif, and larger sediment production during glacial periods. In these cores, 10BeMS varies between 2×104 and 7×104 at. g−1 (atoms per gram of quartz), which corresponds to average denudation rates of 0.2 to 0.5 mm yr−1 between 70 and 4 ka (Mariotti et al.2021). On land, cosmic ray exposure (CRE) ages of polished river cliffs have revealed fast incision rates of ∼0.5 to 2 mm yr−1 during the Late Pleistocene at most sites of the Var catchment and in the Bevera River (Cardinal et al.2022; Petit et al.2019; Rolland et al.2017; Saillard et al.2014). River gorges located at high altitudes in the Mercantour massif (Sites 1 and 2 in Fig. 1; red and brown dots in Fig. 2) show very fast incision (up to 4 mm yr−1) starting after the Younger Dryas (YD), which can be ascribed to a transient response of formerly glaciated valleys. Most of the lower-altitude river gorges (Sites 3, 4, 5, 6 and 8 in Fig. 1; orange, yellow and green dots in Fig. 2) start to be incised around 20 ka (i.e. close to the LGM). Two other sites (7 and 9 in Fig. 1; cyan and blue dots in Fig. 2) show much lower incision rates (<0.5 mm yr−1) extending from 0 to 80 ka.

Figure 2Altitude above the riverbed vs. cosmic ray exposure (CRE) ages and corresponding error bars in samples of river-polished surfaces of the Var, Roya and Paillon catchments (Cardinal et al.2022; Rolland et al.2017, 2020; Petit et al.2019; Saillard et al.2014). See location of the sites in Fig. 1. Dashed lines indicate constant incision rate curves for comparison.


Previously published interpretations of these data have suggested that rivers with sources in previously glaciated areas (Tinée, Vésubie) incised faster during and after the LGM not only because of increasing precipitation but also because of the massive release of glacier meltwaters and stored sediments occurring at that time (e.g. Saillard et al.2014; Rolland et al.2020). More recently, Cardinal et al. (2022) have pointed out a complex response of river systems of the SW French Alps to deglaciations, depending on their connection with glaciated areas and on the presence of lithological knickpoints.

To summarize, river incision data suggest the following: (1) transient and rapid incision following the Younger Dryas in high-altitude areas; (2) steady, fast incision rates of ∼1 mm yr−1 for the last 15–20 kyr at almost all other points; and (3) lower incision rates of ∼0.2 to 0.5 mm yr−1 in the Paillon and Roya rivers, east of the Var catchment. In all catchments, it appears that data points in the last ∼20 kyr range along slopes that define larger incision rates (sometimes by 1 order of magnitude) than average catchment denudation rates estimated from detrital 10Be in river sands or marine sediments (Mariotti et al.2019, 2021). This is not necessarily contradictory, as incision is a local phenomenon compared to the average catchment surface denudation. In particular, gorges where fast incision occurs can typically induce a transient decoupling from the catchment base level for the surrounding hillslopes (Reinhardt et al.2007). However, there is a fundamental difference in the interpretation of river incision rates and 10BeMS from deep marine sediments: while gorge bedrock surface exposure ages in the last 20 kyr define high incision rates, suggesting that the post-LGM period was characterized by enhanced incision and gorge entrenchment, the geochemical signature of marine sediments is interpreted (Bonneau et al.2017; Mariotti et al.2021) as reflecting more intense erosion during the glacial episodes than during deglaciations.

3 Methods

3.1 Landscape evolution model

We use the landscape evolution model (LEM) Badlands (Salles2016) to simulate erosion and deposition on an TIN topographic grid under various tectonic (uplift) and climatic (precipitation, hillslope processes) time and space parameters. River incision is simulated using the stream power law (SPL) (Whipple and Tucker1999):

(2) ε = K P A m S n ,

where ε is the erosion rate [L T−1], K is a dimensional coefficient that describes the erosional efficiency [L1–2m T−1], A is the drainage area [L2], P is the spatial and temporal variation in precipitation rate P relative to a mean precipitation rate P0 (1 m yr−1; the product being used as a proxy for the discharge), S is the channel slope, and m and n are positive exponents. An m/n value close to 0.5 was previously estimated in the Var and Vésubie catchments (Saillard et al.2014; Petit et al.2017), so in all models, we chose to set  m and n to the classically used values of 0.5 and 1, respectively, and to tune the K value in order to fit measured incision rates. Hillslope processes on land and short-distance sediment transport at sea can be simulated by a linear diffusion law:

(3) d h d t = K D 2 h ,

where dh/dt is the altitude change due to diffusive processes [L T−1], h is the altitude [L] and KD is the hillslope diffusion coefficient [L2 T−1], which can vary between the continental and marine domains. We fix the diffusion coefficient to low values (0 to 0.025 m2 yr−1) both on land and at sea in order to ensure that any observed smoothing effect on the 10Be record in deposited sediments is not due to diffusive processes (Table 2). For the same reason, we do not consider non-linear diffusion components. For river systems, we use the detachment-limited law (Eq. 2), but we impose sediment deposition either when the channel slope falls below a given threshold (alluvial plain deposition; see Table 3) or when the rivers reach their base level or an endorheic depression. In these cases, available sediment fluxes carried by rivers are used to compute the volume of sediments to deposit. If transported sediment fluxes, when deposited, are insufficient to fill the depression or to reach the prescribed channel slope threshold, all sediments would be deposited and the outgoing river sediment concentration would be zero. If, on the other hand, the available sediment flux exceeds the required deposition volume, the excess flux will be carried out to the downstream nodes. A low critical slope of 0.5 % is applied, except for two models where a very low threshold of 0.01 % is used in order to drastically limit alluvial plain deposition. In addition, we take into account submarine sediment transport in order to simulate the occurrence of hyperpycnal flows. Following an approach similar to Petit et al. (2015) and Thran et al. (2020), we assume that hyperpycnal flows occur when the sediment load at the river mouth is larger than a given threshold. If the flow density exceeds this threshold, instead of being deposited near the base level, sediments continue their route along the submarine slope. In addition, we assume that the flow does not incorporate water along its path. The flow density at the river mouth ρf [M L−3] can be computed as follows: (i) either with a mass estimate from water and sediment discharge (Qw and Qs, respectively [L3 T−1]) and densities (ρw and ρs, respectively),

(4) ρ f = ρ w Q w + ρ s Q s Q w + Q s ,

or (ii) using the rating parameters a and b (Syvitski et al.2000), which can be determined for each river system from discharge and sediment load measurements. While Eq. (4) considers the total mass of sediment transported by the river (i.e. bedload and suspended load), the rating parameters allow an empirical estimation of the suspended load, the one that effectively contributes to the increase in flow density. We consider that submarine flow can trigger bedrock erosion as for aerial channels, but the parameters of the stream power law are adjusted in order to account for (1) constant drainage area along the channel length in the submarine domain and (2) lower shear stress on the submarine channel bed compared to aerial rivers, which is simplified assuming an effective slope Seff such that

(5) S eff = S ρ s - ρ w ρ s .

Deposition occurs in the submarine domain either for gentle slopes (similar to alluvial plains on land) and beneath a certain depth (−2300 m in most models) corresponding to the depth of the abyssal plain. Flexural isostasy can be incorporated with a constant or space-variable effective elastic thickness (EET) used to compute the vertical motion resulting from the response of the lithosphere to loading (by ice, sedimentation or sea level rise) or unloading (deglaciation, erosion or sea level drop). A flexural isostatic response of the lithosphere is computed using the flexure equation:

(6) D 4 w ( x , y ) + Δ ρ g w ( x , y ) = L ( x , y ) ,

where w is the vertical deflection, Δρ is the density contrast between the mantle and the filling material, L is the load (N) and D is the flexural rigidity of the lithosphere (N m):

(7) D = E EET 3 / 12 1 - ν 2 ,

with E, Young's modulus, and ν, Poisson's ratio, equal to 1011 Pa and 0.25, respectively. The flexure module in Badlands uses the gFlex package (Wickert2016). Apart from flexural response to erosion/sedimentation and ice and seawater loading and unloading during the model run, no vertical motions are applied to the topography.

3.2 Ice cover and sea level changes

As the Mercantour massif was periodically covered by glaciers during the Quaternary, we simulate the ice thickness and extent at every time step assuming that the LGM corresponds to the maximum ice extent map (Brisset et al.2015). We consider that ice thickness varies with Mediterranean sea surface temperatures (SSTs), which ranged between ∼5C (LGM) and ∼15C during the considered time period (Hayes et al.2015; Rodrigo-Gamiz et al.2013). Glacial periods with full ice extent are imposed for SST lower than 6.5 C and complete deglaciation for SST above 11 C. Between these thresholds, the ice thickness is assumed to vary linearly with SST. In order to avoid rapid variations in the ice cover, the SST curve is smoothed using a 5 kyr sampling interval, which is then resampled at a 1 kyr step using cubic interpolation (Fig. 3). When ice thaws, an equivalent amount of water (assuming a ratio between the ice and water heights of 0.93) is released as runoff. Sea level variations can be imposed according to the Mediterranean eustatic sea level curve published in Waelbroeck et al. (2002). Badlands does not simulate glacial erosion; however, we must consider a non-zero erosion rate beneath glaciated areas in order to avoid overestimation of the 10Be concentration in sediments produced by basement erosion after glaciers retreat. For this purpose, we simulate the in situ erosion and sediment production due to glacial processes by increasing the local hillslope diffusion coefficient proportionally to the ice thickness, which will locally enhance denudation beneath glaciated areas, while river discharge is set to zero. This simplified representation is based on the assumption that glacial sediment fluxes are proportional to the topographic slope (hence to ice velocity) and that glacial erosion is related to the shear stress exerted by the glacier on the bedrock (Boulton1996). As glacial erosion is simulated by local processes (diffusion) and river discharge is set to zero beneath glaciers, the export of glacial sediments is very low during glaciations. Besides erosion, the effect of ice coverage is twofold: it blocks cosmic radiation so that 10Be production is zero beneath areas covered by glaciers, and it creates a positive vertical load and downward flexure of the lithosphere.

Figure 3(a) Topography, bathymetry and maximum ice thickness (after Brisset et al.2015) and contour of the main catchments (dashed line); (b) smoothed variations in the sea surface temperature (SST; after Hayes et al.2015, and Rodrigo-Gamiz et al.2013) in blue and ice thickness ratio (1 is for full ice, 0 for no ice) in red; (c) sea level variations (after Waelbroeck et al.2002).

3.310Be production and transport

10Be production rates by neutron spallation and muon capture are the same as in Mariotti et al. (2019) and computed according to Braucher et al. (2011) and Martin et al. (2017) for a latitude of 40 (Table 1) using the scaling parameters by Stone (2000). Earth magnetic field variations are not considered in this study. The topographic shielding is computed from the TIN topographic grid. A shielding correction can be applied to account for the topographic smoothing due to the digital elevation model (DEM) resolution, which tends to underestimate the actual shielding (Norton and Vanacker2009).

Table 1Sea level high-latitude parameters for 10Be production (after Braucher et al.2011).

Download Print Version | Download XLSX

A map of quartz-bearing rocks is defined according to the geological map of Nice and its hinterland in Rouire et al. (1980); most quartz-bearing rocks correspond either to granitic and metamorphic Palaeozoic basement rocks in the Mercantour massif or to Cenozoic sandstones in the sedimentary cover (Fig. 1). The average quartz concentration in source rocks is fixed at 50 %. The initial 10Be concentration in quartz-bearing rocks is computed assuming a steady-state average denudation rate for the whole grid. 10Be concentration N(z,t) varies with time and depth, and we simply compute it at the surface (z=0) of eroded domains (Lal1991):

(8) d N ( 0 , t ) d t = P ( 0 , t ) - λ + ρ ε ( t ) Λ N ( 0 , t ) ,

where N is the 10Be concentration (at. g−1), P is the production rate (at. g−1 yr−1), λ is the 10Be radioactive decay constant (yr−1), ρ is the rock density (g cm−3), Λ is the attenuation length (g cm−2) and ε is the erosion rate (cm yr−1). At each time step, the production rate is computed taking into account the quartz abundance of the source rock and the potential shielding of cosmic rays by the surrounding topography and/or by the ice or sea cover.

10Be production results primarily from neutron spallation and fast and slow muon capture with different production rates and attenuation lengths (Braucher et al.2011) (Table 1). Assuming that erosion and production rates are constant during a given time step of the model, we can explicitly compute the 10Be concentration in each eroded node of the source rocks at each time step, without any a priori steady-state assumption (Knudsen et al.2019):

(9) N ( 0 , t ) = N ( 0 , t - Δ t ) exp - λ + ρ ε ( t ) Λ Δ t + P ( 0 , t ) λ + ρ ε ( t ) Λ 1 - exp - λ + ρ ε ( t ) Λ Δ t .

This Eulerian formulation, where the erosion rate ε is equal to the vertical velocity at which the rock material is vertically advected up to the surface is chosen because it does not necessitate computing vertical concentration profiles for each grid point.

For the initial 10Be concentration N(0, 0) we assume a steady state between 10Be production and erosion by imposing a mean long-term erosion rate ε0:

(10) N ( 0 , 0 ) = P ( 0 , 0 ) ρ ε 0 Λ + λ .

At each time step, the total sediment, quartz and 10Be fluxes are computed on each grid node. In the case of deposited sediments, we then compute the mean detrital 10Be concentration and quartz content of sediments, knowing the volume contribution, quartz proportion and 10Be concentration of each eroded source to the total amount of deposited sediments, assuming a perfect mixing between all sources. As initial conditions, we use a smoothed topographic and bathymetric DEM with a spatial resolution of 500 m × 500 m.

In a first test, a simulation is run with a constant precipitation rate of 0.5 m yr−1 over 5000 years with a time step of 100 years in order to calibrate SPL parameters, mainly the erodibility coefficient K against the results of Mariotti et al. (2019). From this model (M1), we compute the 10Be production rate in the Mercantour massif, the mean 10Be concentration of continental river sediments corresponding to the catchments sampled by Mariotti et al. (2019) and the related mean catchment denudation rate deduced from these concentrations with an initial steady state 10Be concentration computed using a mean denudation rate of 0.2 mm yr−1, in agreement with Mariotti et al. (2019) (Fig. 4 and Table 2). Discrepancies between our model and the results of Mariotti et al. (2019) can arise from (1) the fact that we use a smoother topographic grid and a simpler map of quartz-bearing areas (for 10Be production rate mismatches) and (2) the fact that we compute erosion for each grid node and then average it over the catchment area instead of considering a spatially homogeneous catchment denudation rate (for erosion rate mismatches). Data from the Coulomp River, a small tributary of the Var River, could not be satisfyingly reproduced (Table 2, Cou-16.1): the predicted low 10Be concentration resulted in a denudation rate more than twice that of Mariotti et al. (2019). Similarly, one point in the Vésubie River (Ves-16.1), which has a large 10Be production rate, shows a modelled 10Be concentration slightly larger than measured by Mariotti et al. (2019). As a consequence, the erosion rate for this point is underestimated (Fig. 4). While not perfect, simulated estimates of 10Be production and erosion rates provide a good fit with observations and allow us to constrain the value of K. Several improvements might reduce the observed discrepancies, first by using a higher DEM resolution and second by accounting for underground (karstic) water circulation known to take place in the region, like for the Coulomp River (Audra et al.2009), which should modify river discharge and related incision rates.

Table 2Calibration of 10Be production and erosion rates after the results of Mariotti et al. (2019) for the present-day rates.

Download Print Version | Download XLSX

Table 3Model parameters.

Download Print Version | Download XLSX

Figure 4Comparison between the results of Mariotti et al. (2019) (x axis) and the short-term model (this study, y axis) for several rivers of the Var catchment for 10Be production rates (a), 10Be concentration in river sediments (b) and steady-state erosion rates (c). Red line indicates the 1:1 slope. See Table 2 for values.


4 Results

4.1 Record of time-variable erosion rates in river sands and at sea

We then test the response of a topographic grid representing the Var aerial and submarine systems to climatic (precipitation) variations. The model runs for 100 kyr with an adaptive time step of a maximum of 1000 years. We keep a constant rock erodibility of 5×10-6 yr−1 and a hillslope diffusion coefficient of 2.5×10-2 m2 yr−1. The river sediment load is computed from water discharge using rating parameters (Syvitski et al.2000) a and b  for the Var River equal to 1×10-3 and 1.6, respectively, and a low-threshold flow density (equal to water density), which ensures that all river sediments are exported to the deep submarine basin and that no large river delta is formed near the coastline.

In a first series of tests, we want to investigate how well the 10Be concentrations in continental (i.e. in river sand) and submarine (turbidite-like) deposits compare with the average catchment denudation rate directly output by the model. For this purpose, we select the total area of the Var catchment where 10Be is produced and compute the average 10Be concentration Nex in the volume of rock eroded from this source for each time step (hereafter called “in situ” sediments because they are not yet transported or deposited) from the contribution of each catchment node i, such that

(11) N ex = i = 1 n N i ε i Q i i = 1 n ε i Q i ,

where Qi is the quartz abundance, Ni is the 10Be concentration (at. g−1) and εi the eroded mass of sediments (in g). Then we use Eq. (10) to compute the catchment denudation rate assuming a steady-state condition and compare it with the average “actual” denudation rate directly output by the model (i.e. the average volume of eroded sediments in the same area per time step). Finally, we extract the average 10BeMS value in an 8×8 km square areas located at the mouth of the Var River for model M6, in the submarine canyon for model M7 and in the deep submarine basin for models M2 to M5. We then use the same equation to compute the average catchment denudation rate as recorded by deep-sea sediments (as shown in Fig. 5). Simulations are run with alternating low (0.25 m yr−1) and high (1 m yr−1) precipitation rate periods lasting 20 kyr each (Fig. 6). The first simulation (M2) is run with a constant sea level, no ice cover and no lithospheric flexure; then we successively implement a variable Mediterranean sea level, ice cover in the Mercantour massif (M3 and M4; see Fig. 3) and lithospheric flexure with a constant EET of 20 km (M5), which corresponds to a moderately rigid lithosphere where the crust and mantle elastic lids are decoupled (Burov and Diament1995). The last two models of this series (M6 and M7) have similar parameters to M2 except that (i) we drastically decrease the slope threshold for alluvial plain deposition, which almost impedes sediment storage on land, and (ii) we force deposition in the submarine domain at shallower depths: either in the coastal delta (−200 m) for M6 or in the submarine canyon (−800 m) for M7.

Figure 5Schematic description of the approach used to compare the apparent denudation for a given time step computed from the mean 10Be concentration of sediments produced in situ (left part of the topographic profile, apparent denudation A) and from the mean 10Be concentration of marine sediments (10BeMS; right part of the profile, apparent denudation B) with the true denudation.


Figure 6Test of the apparent denudation rates given by the 10Be concentration of in situ (red line) and deep marine (10BeMS, black line) sediments versus the actual one extracted from the model (blue dots). Periods of large precipitation rates are indicated by the transparent blue vertical strips. (a) Model with only variable precipitation rate; (b) model with variable precipitation rate and sea level variations; (c) model with variable precipitation rate, ice cover and sea level variations; (d) model with variable precipitation rate, sea level variations, ice cover and flexure; (e) model similar to (a) but with a lower critical slope for sediment deposition and forced deposition below −200 m; (f) model similar to (e) but with forced deposition below −800 m.

The first model (M2, Fig. 6a), where there are no sea level changes, no ice cover and no flexural isostatic response of the lithosphere, shows that the two periods of more intense precipitation are well-detected in in situ sediments, the 10Be signature of which gives erosion rates consistent with the actual ones. Marine sediments also seem to record two periods of larger erosion rates, but the first one (between 25 and 45 kyr) is barely visible and the peak of the second one (80–90 kyr) occurs 10 to 20 kyr after the middle of the second high-precipitation period. The second model (M3, Fig. 6b), where sea level variations are present, displays approximately the same behaviour as M2 for in situ sediments, with apparent erosion rates consistent with real ones. Erosion rates computed from marine sediments also show fluctuations consistent with two periods of more intense precipitation, but they appear much smoother and offset by ∼20 kyr, compared to the actual ones. In the two following models where ice cover is present (M4 and M5, Fig. 6c and d), the effect of partial coverage of the Mercantour massif by ice during almost all the duration of the model is twofold: (i) it reduces the effect of precipitation variations, since no run-off occurs beneath the ice cover, and (ii) it reduces the production and exportation of 10Be-poor sediments, except during the last period (80–100 kyr into the model run) when the massif is free of ice. As a consequence, while in situ sediments still record the actual erosion rate well, the 10Be signal in the submarine domain strongly underestimates it. Finally, models similar to M2 but with a lower critical slope and where sediment deposition is forced either in the delta (M6) or in the submarine canyon (M7) show that marine sediments of the delta record the onset of first precipitation pulse well, with only a small time lag (5 kyr). Then, the apparent erosion rate computed from marine sediments progressively decreases, and the second period of intense precipitation is not recorded. In this model, 10Be-poor sediments coming from the second precipitation pulse are not deposited in the sampled area but further south, due to delta progradation. Finally, if deposition is forced to occur in the middle of the submarine slope, the first precipitation peak is still detected, but it appears smoother and of a lower amplitude and it occurs ∼20 kyr later than on land (Fig. 6f).

In simulations without ice cover (Fig. 6a, b, e and f), denudation rates estimated from in situ sediments do record the succession of periodic pulses but generally overestimate the actual denudation rate by a few percent. This overestimation reflects a more intense erosion in low-altitude river thalwegs than on the hillslopes and high-altitude interfluves, which slightly promotes the export of 10Be-poor sediments with respect to 10Be-rich ones, hence increasing the apparent denudation rate.

For all simulations, both in situ and actual denudation rate evolutions depict more complex patterns than the imposed climatic (precipitation) forcing. Short-period, low-amplitude variations are visible, which are related to local and internal adjustments of the modelled topography and not to the external forcing, although sea level variations seem to slightly modify their amplitude (Fig. 6a and b). There is no clear time lag between the onset of higher-precipitation periods and their record in in situ sediments, but the apparent denudation rate from in situ river sediments displays a curved shape at the transition between low- and high-precipitation rate periods (well visible, for instance, in Fig. 6a between 60 and 80 kyr of the model run), which corresponds to well-known analytical solutions for periodic changes in erosion rates (Bierman and Steig1996) and is consistent with the detachment-limited hypothesis. This reflects the time needed for 10Be concentration to reach a steady state relative to the massif denudation rate (∼5 kyr); this effect is also visible at the transition from high to low precipitation rates (Fig. 6e and f).

Figure 7Difference in 10BeMS between two models with the same parameters (Table 3), one of them having imposed 10Be-rich (2×106 at. g−1) rocks above 1800 m between 30 and 35 kyr (transparent red vertical strip in panel b). Panel (a) shows the location of the area of 10Be enrichment in the Mercantour massif. The time lag (thick dashes) is defined as the interval between the middle of the 10Be enrichment period in the massif and the age of the largest 10Be peak in marine sediments. The relaxation time is defined as the interval between the end of the 10Be enrichment period in the massif and the end of excess 10Be recorded in marine sediments.

This series of tests seem to indicate that (i) the alluvial sediments produced in situ record the variations in the rate of denudation well, although they may slightly overestimate them; (ii) 10Be of submarine sediments 10BeMS in the deepest part of the basin does not allow us to retrieve sharp climate variations, especially if glaciers cover the area where 10Be is produced for a long period of time, but correctly estimate the average long-term erosion rate, with only smooth variations; and (iii) finally, there can be a significant time lag between the middle of the high-precipitation periods and the peak in denudation rates recorded in submarine sediments.

In order to refine our understanding of the smoothing and time lag effects, we aimed at tracing artificially enriched 10Be-rich sediments in the mountain range down to the submarine basin. To do so, we imposed a constant and artificial high 10Be concentration in high-altitude areas (above 1800 m) over a 5 kyr period from 30 to 35 kyr into the model run. Then, we compare the 10Be signature of deep-sea sediments between two identical models: one run with the 10Be enrichment and one without (model M8, Fig. 7). We do not show all the tests here for the sake of simplicity but rather illustrate a typical signature of this transient 10Be enrichment as visible in submarine sediments. This result is obtained accounting for a variable sea level, ice cover and lithospheric flexure and a constant precipitation rate of 0.5 m yr−1 (Table 3). If we consider the earliest arrival of 10Be-rich sediments in the submarine basin, the time lag from source to sink appears to be relatively small (1–3 kyr) since 10Be-rich sediments arrive in the basin shortly after they begin to be produced in the massif (Fig. 7); however, the time lag between the middle of the 10Be enrichment period on land and the largest 10Be peak recorded in marine sediments is rather large (∼10 kyr). Moreover, this signal takes a long time to relax: sediments that are richer in 10Be than the reference simulation still reach the basin ∼25 kyr after the end of the enrichment period. These tests suggest that, although the initial time lag is not necessarily very important, the relaxation time, i.e. the time it takes for the 10Be concentration to return to its background value, can be large. It suggests that, in the case of a succession of climatic events with a frequency greater than this relaxation time, it will therefore be impossible to differentiate between them in the 10Be signal. We want to stress here that this delay is different from the time needed to reach cosmogenic radionuclide steady state, as computed in Bierman and Steig (1996). Here, we do not force high 10Be concentration in the Mercantour massif by drastically reducing the erosion rate; it would indeed take some time for the 10Be signal to adapt and reach a new steady state. Instead, we instantaneously apply a large, artificial 10Be value in the mountain peaks, which should be immediately visible in deposited sediments if the time lag was zero.

Figure 8Stages of evolution of the reference model from 0 to 100 kyr (model time counted forward, corresponding to 100 to 0 ka in reality and on the figure time stamps). Upper panels (a) present the topography and ice thickness, middle panels (b) present the surface 10Be concentration and lower panels (c) present a blow-up of the surface 10Be concentration in submarine sediments in the last 30 kyr, with an adjusted colour map. The solid rectangle in the lower panels indicates the area where average 10BeMS has been computed (see results in Fig. 9).

Figure 9Comparison between model outputs and data from CRE ages in river-polished cliffs (ad; see Fig. 2) and 10BeMS (eg; after Mariotti et al.2021). Vertical bars in the right-hand graphs of each panel correspond to measured standard deviations in 10BeMS in the sampled area (Fig. 8). See Table 3 for model parameters.


4.2 Reference model for the Var catchment over the last 100 kyr

Finally, we try to determine if it is possible to find model parameters for which the result matches both the previously determined river incision rates during the last 30–40 kyr in the Nice hinterland (Fig. 2) and measured 10BeMS for the same period (Mariotti et al.2021). We chose not to consider the two river gorge sites that have been dated in the upper Tinée valley (Salso Moreno and Isola; see Fig. 1) because they likely correspond to very transient post-glacial sediment wash-out that it is not possible to model (Rolland et al.2020). This “reference” model has been determined by a trial-and-error method, with initial and boundary conditions similar to the previous series of runs (Figs. 8 and 9a, model M9). Sea level variations and ice cover are imposed, and simulations are run for 100 kyr with a time step of 1 kyr. We assume that for the last period of the model, which corresponds to the present-day situation, the precipitation rate should be compatible with current average values (i.e. 700–800 mm yr−1;, last access: 25 February 2023). Here again, the flexural isostatic response is computed using a constant effective elastic thickness of 20 km. We also tested different EET values of 10 and 30 km as well as a space-variable EET, with larger values on land than at sea. Since variations in EET are difficult to constrain and only cause minor differences in the final result, we chose to impose this constant value of 20 km on the entire grid for all models. We first present a model which satisfyingly reproduces measured river incision rates and 10Be measurements; then we discuss the implications of each parameter (flexure, sea level, ice cover) in the final result. This simulation involves a precipitation rate of 0.3 m yr−1 from 0 to 80 ka increasing to 0.7 m yr−1 during the last 20 kyr, an initial denudation rate of 0.2 mm yr−1, a diffusion coefficient of 2.5×10-2 m2 yr−1 in submarine and river sediments and of 0.1 m2 yr−1 in bedrock areas, and a constant erodibility coefficient of 3.5×10-6 yr−1 (Figs. 8 and 9a and Table 3). This model satisfyingly reproduces measured incision rates in river channels and yields a slight increase in 10BeMS in the Var deep-sea fan after ∼40 ka (Fig. 9a). Given the large variability of 10BeMS values (as depicted by the standard error bars), the modelled 10BeMS variation is compatible with almost all measurements published in Mariotti et al. (2021).

From the simulation outputs, we find that the increase in 10BeMS is not due to a change in erosion rate on land: indeed, river incision rates tend to increase after 20 ka due to the release of glacier meltwaters and to increased precipitation rates, while the 10BeMS also increases 2-fold. The observed increase in 10BeMS around 30–40 ka is partly due to the presence of patches of 10Be-rich sediments deposited in the upper part of the basin after about 30 ka, which slowly feed lower areas of the basin in 10Be-rich sediments during the run (Fig. 8b and c). The simulation fit with the measured 10BeMS may somehow be fortuitous, since the surface concentration in 10Be varies at lot, both locally and vertically. This simulation thus shows that (1) 10Be concentration at the surface of the submarine basin at a given time can be highly variable (hence, the vertical variation in 10BeMS can also be variable from one place to another) and that (2) it is quite possible to get both increased incision rates on land and apparent decreased denudation rates in submarine cores for the same external (climatic) forcing.

Starting from this reference simulation (Fig. 9a), we then evaluate the individual effects of lithospheric flexure, ice cover and sea level changes. A similar simulation without lithospheric flexure predicts slightly lower incision rates for all rivers, with a more dramatic change for the Vésubie River where it drastically decreases (Fig. 9b, model M10); meanwhile, both the absolute values of 10BeMS and their variations are still consistent with measurements. A simulation similar to the reference one but without ice cover (Fig. 9c, model M11) predicts slightly too low incision rates for some rivers (Estéron, lower Tinée) but a good fit to 10BeMS and to some other rivers (Paillon and Roya) and rather high and constant 10BeMS. Removing sea level changes does not significantly affect the simulated river incision rates but gives a lower 10BeMS signal, still with an increase after 30 ka (Fig. 9d, model M12). A simulation with a constant precipitation rate gives satisfying results for the river incision rates but a constant and too low 10BeMS, especially for the most recent (20 to 0 ka) period (Fig. 9, model M13). Finally, removing both the flexural isostatic response of the lithosphere and the main sources of load (i.e. ice cover and sea level variations) has large effects both on the river incision rates (especially for the Vésubie River) and on the 10BeMS signal (Fig. 9f and g, models M14 and M15). Interestingly, the last three models without precipitation changes or without flexure ± ice cover or sea level changes produce rather stable 10BeMS values throughout the model. In fact, the main difference between the reference model M9 and for instance models M13, M14 and M15 is that the first 50 kyr of model M9 are characterized by a lower 10BeMS concentration. This is possibly due to a relatively lower contribution of high-altitude rocks to the 10BeMS signal when the massif is partly covered by ice.

Based on these models, we show the following:

  • Simulations with constant precipitation or without flexure and without either ice cover or sea level changes do not produce the observed increase in 10BeMS after 30–40 ka.

  • The large incision rate in the Vésubie River (compared to Bevera, Roya, Estéron and Paillon) can only be explained by the effect of isostatic rebound, the latter being smaller in the simulations with no glaciers.

  • Cosmic ray exposure (CRE) data on river-polished walls are not dense enough to discriminate between constant (0.5 mm yr−1) or variable (0.3 then 0.7 mm yr−1) precipitation rates.

  • However, a constant precipitation induces more important denudation in the earliest stages of the predicted evolution (compared to the reference simulation where precipitation rate is low), hence producing on average a lower 10BeMS at the end of the model.

  • A combination of sea level changes, ice cover and lithospheric flexure provides a reasonable fit both to river incision rates and to measured 10BeMS.

5 Discussion and conclusions

These series of models, where continental sediment deposition is controlled almost exclusively by one parameter (the critical slope for alluvial deposition), show that it is sufficient to induce significant differences in the 10Be signal between its production in quartz-bearing rocks of high-altitude massifs and its record in the submarine domain. However, here, we do not capture all the details of erosion, transport, sediment mixing and sedimentation, as we do not consider different lithologies, grain sizes, or any complex mixing law between the various sediment sources. In contrast to the study by Carretier et al. (2020), the lack of grain tracking prevents us from being able to measure particle residence times and compare them to the average geochemical signature of the sediments.

Submarine sediment transport and deposition modelling is almost as simple as on land and does not consider ocean dynamics (turbulence of turbiditic flows, coastal currents), which are responsible for characteristic sedimentary features such as the Var sedimentary ridge at the bottom of the Var Canyon (Migeon et al.2001). These phenomena can also contribute to sediment dispersion and smoothing of the 10Be signal offshore.

Similarly, it could be interesting to couple our models with glacial erosion models (e.g. Seguinot and Delanay2021) in order to better control the variability of glacial erosion and sediment production during the last cycles and better quantify how it affects 10Be production and exportation.

Finally, initial model conditions are another possible source of uncertainty, especially when it comes to our initial topography, which is derived from present-day elevation and therefore limits our ability to run the model over a longer period of times as it would produce a final topography too different from the current one. A better estimate of initial quartz abundances in source rocks (i.e. crystalline lithologies of the Mercantour massif and Cenozoic sandstones of the foreland) would also help to reduce the uncertainties associated with quartz and 10Be transfer from source to sink.

Our results based on dating and modelling indicate that, while river sands do accurately estimate the average denudation rate of continental catchments for the Var region (provided the latter does not vary at high frequency, i.e. with periods smaller than the time needed to reach a steady state), this is much less the case for submarine deep-sea sediments. These sediments have a different and often much smoother signature than continental ones and record significant lag and/or relaxation times with respect to external forcing, probably due to the geomorphological response of the continental margin both on land and at sea. With this area being prone to strong and rapid geomorphological modifications (i.e. transition from narrow bedrock channels to wide, braided rivers) during violent flood events (like during the Storm Alex, which took place in October 2020), it could be of primordial importance to estimate the relaxation time of such events and their role in the long-term landscape evolution and geochemical signature of the sedimentary archives.

On the long-term, the presence of ice in the massifs where 10Be production occurs together with the reworking of alluvial and deltaic sediments during low sea level periods and vertical motions due to lithospheric flexure largely modify the signal coming from precipitation variations and can lead to poor estimates of the actual denudation rate variations from 10BeMS. All these effects have been exemplified in this study of the Var catchment where the distance from source to sink is short, precipitation rates are large and the mouths of the main rivers are devoid of any large deltas. Hence, in regions with very large catchments, alluvial plains and deltas, it could be even more difficult to reconstruct past denudation rates from deep-sea sediments.

Our reference simulation highlights the complex interactions between river incision, sea level variations, ice coverage and the resulting isostatic response of the lithosphere. It seems impossible to disentangle the respective role of any of these forcings specifically, as most of them are interdependent. However, further tests (Fig. 9) show that, depending on their location, rivers have a different sensitivity to these parameters: the Estéron, Roya, Paillon and Bevera are less affected by the parametric changes applied to the reference forcing conditions; their incision rate is only significantly reduced when considering no sea level changes, ice cover or flexure (not shown here).

In contrast, the rivers with the largest incision rates (the Vésubie and, to a lesser extent, the Tinée River) are also the ones which seem more sensitive to the effect of ice cover and flexural isostatic response of the lithosphere (Fig. 9). This is possibly because a significant part of their length (20–40 km) runs over the Mercantour crystalline massif, i.e. over the area where post-glacial isostatic rebound is important. The sampling sites of these rivers being rather close to their headwaters, they are not sensitive to sea level variations. Quantification of each river sensitivity to local or regional processes like isostatic uplift, sea level changes, precipitation or ice cover, depending on the dating site location, could therefore be useful to better understand the respective importance of these external forcings on this southern Alpine margin.

Code and data availability

The relevant version of Badlands with all the data used for this paper can be found here: (Salles2023).

Author contributions

CP developed the implementation for Badlands, performed model runs and wrote the first draft of the paper. TS checked and released the new version of the code and contributed to the paper. YR, VG and LA contributed to the paper.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


Fruitful discussions with Guillaume Duclaux (Geoazur) were greatly appreciated. The authors thank Sebastien Carretier and Yanyan Wang for their insightful comments on the first version of this paper and associate editor Simon Mudd for handling the review process and for the final review of this paper.

Financial support

This study is part of a project that has been funded by the French Geological Survey (Bureau de Recherches Géologiques et Minières; BRGM) through the national programme “Référentiel Géologique de France” (RGF-Alpes). This work has been supported by the French government, through the UCA-JEDI Investments in the Future project managed by the National Research Agency (ANR) with the reference number ANR-15-IDEX-01.

Review statement

This paper was edited by Simon Mudd and reviewed by Yanyan Wang and Sebastien Carretier.


Armitage, J., Dunkley Jones, T., Duller, R., Whittaker, A., and Allen, P.: Temporal buffering of climate-driven sediments flux cycles by transient catchment response, Earth Planet. Sc. Lett., 369–370, 200–210,, 2013. a

Audra, P., Mocochain, L., Bigot, J.-Y., and D'Antoni-Noblecourt, J.-C.: The Grand Coyer karst, exploration at the Coulomp spring (Alpes-de-Haute-Provence, France, in: Proceedings of the 15th International Congress of Speleology, Kerrville, Texas, USA, 1755–1759, 2009. a

Bentley, S., Blum, M., Maloney, J., Pond, L., and Paulsell, R.: The Mississippi River source-to-sink system: Perspectives on tectonic, climatic, and anthropogenic influences, Miocene to Anthropocene, Earth-Sci. Rev., 153, 139–174,, 2016. a

Bierman, P. and Steig, E.: Estimating Rates of Denudation Using Cosmogenic Isotope Abundances in Sediment, Earth Surf. Proc. Land., 21, 125–139, 1996. a, b, c, d, e

Blöthe, J. and Korup, O.: Millennial lag times in the Himalayan sediment routing system, Earth Planet. Sc. Lett., 382, 38–46,, 2013. a

Bonneau, L., Toucanne, S., Bayon, G., Jorry, S., Emmanuel, L., and Silva Jacinto, R.: Glacial erosion dynamics in a small mountainous watershed (Southern French Alps): A source-to-sink approach, Earth Planet. Sc. Lett., 458, 366–379,, 2017. a, b, c

Boulton, G.: Theory of glacial erosion, transport and deposition as a consequence of subglacial sediment deformation, J. Glaciol., 42, 43–62,, 1996. a

Braucher, R., Merchel, S., Borgomano, J., and Bourlès, D.: Production of cosmogenic radionuclides at great depth : A multi element approach, Earth Planet. Sc. Lett., 309, 1–9,, 2011. a, b, c

Brisset, E., Guiter, F., Miramont, C., Revel, M., Anthony, E. J., Delhon, C., Arnaud, F., Malet, E., and de Beaulieu, J.-L. Lateglacial/Holocene environmental changes in the Mediterranean Alps inferred from lacustrine sediments, Quaternary Sci. Rev., 110, 49–71,, 2015. a, b

Burov, E. and Diament, M.: The effective elastic thickness Te of continental lithosphere What does it really mean, J. Geophys. Res., 100, 3905–3297, 1995. a

Cardinal, T., Rolland, Y., Petit, C., Audin, L., Zerathe, S., Schwartz, S., and the ASTER Team: Fluvial bedrock gorges as markers for Late-Quaternary tectonic and1 climatic forcing in the French Southwestern Alps, Geomorphology, 418, 108476,, 2022.  a, b, c, d, e

Carretier, S., Regard, V., Vassallo, R., Martinod, J., Christophoul, F., Gayer, E., Audin, L., and Lagane, C.: A note on 10Be-derived mean erosion rates in catchments with heterogeneous lithology: examples from the western Central Andes, Earth Surf. Proc. Land., 40, 1719–1729,, 2015. a

Carretier, S., Guerit, L., Harries, L., Regard, V., Maffre, P., and Bonnet, S.: The distribution of sediment residence times at the foot of mountains and its implications for proxies recorded in sedimentary basins, Earth Planet. Sc. Lett., 546, 116448,, 2020. a, b

Clift, P. and Giosan, L.: Sediment fluxes and buffering in the post-glacial Indus Basin, Basin Res., 26, 369–386,, 2014. a

Codilean, A., Bishop, P., Hoey, T., Stuart, F., and Fabel, D.: Cosmogenic 21Ne analysis of individual detrital grains: opportunities and limitations, Earth Surf. Proc. Land., 35, 16–27,, 2010. a

Fryirs, K., Brierley, G., Preston, N., and Kasai, M.: Buffers, barriers and blankets: The (dis)connectivity of catchment-scale sediment cascades, Catena, 70, 49–67,, 2007. a

Godard, V. and Tucker, G.: Influence of Climate-Forcing Frequency on Hillslope Response, Geophys. Res. Let., 48, e2021GL094305,, 2021. a

Goren, L.: A theoretical model for fluvial channel response time during time-dependent climatic and tectonic forcing and its inverse applications, Geophys. Res. Lett., 43, 10,753–10,763,, 2016. a

Hayes, A., Kucera, M., Kallel, N., Sbaffi, L., and Rohling, E.: Lateglacial/Holocene environmental changes in the Mediterranean Alps inferred from lacustrine sediments, Quaternary Sci. Rev., 110, 49–71,, 2015. a, b

Jerolmack, D. and Paola, C.: Shredding of environmental signals by sediment transport, Geophys. Res. Lett., 37, L19401,, 2010. a

Knudsen, M., Egholm, D., and Jansen, J.: Time-integrating cosmogenic nuclide inventories under the influence of variable erosion, exposure, and sediment mixing, Quatern. Geochronol., 51, 110–119,, 2019. a

Lal, D.: Cosmic ray labeling of erosion surfaces: in situ nuclide production rates and erosion models, Earth Planet. Sc. Lett., 104, 424–439, 1991. a, b

Li, C., Yang, S., Zhao, J., Xin, J., Dosseto, A., Bi, L., and Clark, T.: The time scale of river sediment source-to-sink processes in East Asia, Chem. Geol., 446, 138–146,, 2016. a

Liu, Z., Zhao, Y., Colin, C., Stattegger, K., Wiesner, M., Huh, C., and Li, Y.: Source-to-sink transport processes of fluvial sediments in the South China Sea, Earth-Sci. Rev., 153, 238–273,, 2016. a

Lupker, M., Blard, P.-H., Lavé, J., France-Lanord, C., Leanni, L., Puchol, N., Charreau, J., and Bourlès, D.: 10Be-derived Himalayan denudation rates and sediment budgets in the Ganga basin, Earth Planet. Sc. Lett., 333–334, 146–156,, 2012. a

Malatesta, L., Avouac, J.-P., Brown, N., Breitenbach, S., Pan, J., Chevalier, M., Rhodes, E., Saint-Carlier, D., Zhang, W., Charreau, J., Lavé, J., and Blard, P.-H.: Lag and mixing during sediment transfer across the Tian Shan piedmont caused by climate-driven aggradation–incision cycles, Basin Res., 30, 613–635,, 2018. a

Mandal, S., Lupker, M., Burg, J.-P., Valla, P., Haghipour, N., and Christl, M.: Spatial variability of 10Be-derived erosion rates across the southern Peninsular Indian escarpment: A key to landscape evolution across passive margins, Earth Planet. Sc. Lett., 425, 154–167,, 2015. a

Mariotti, A., Blard, P.-H., Charreau, J., Petit, C., and Molliex, S.: Denudation systematics inferred from in situ cosmogenic 10Be concentrations in fine (50–100 µm) and medium (100–250 µm) sediments of the Var River basin, southern French Alps, Earth Surf. Dynam., 7, 1059–1074,, 2019. a, b, c, d, e, f, g, h, i, j, k, l, m

Mariotti, A., Blard, P.-H., Charreau, J., Toucanne, S., Jorry, S., Molliex, S., and Keddadouche, K.: Nonlinear forcing of climate on mountain denudation during glaciations, Nat. Geosci., 14, 16–22,, 2021. a, b, c, d, e, f, g

Martin, L., Blard, P.-H., Balco, G., Lavé, J., Delunel, R., Lifton, N., and Laurent, V.: The CREp program and the ICE-D production rate calibration database: A fully parameterizable and updated online tool to compute cosmic-ray exposure ages, Quatern. Geochronol., 38, 25–49,, 2017. a

Migeon, S., Savoye, B., Zanella, E., Mulder, T., Faugères, J.-C., and Weber, O.: Detailed seismic-reflection and sedimentary study of turbidite sediment waves on the Var Sedimentary Ridge (SE France): significance for sediment transport and deposition and for the mechanisms of sediment-wave construction, Mar. Petrol. Geol., 18, 179–208, 2001. a

Mulder, T., Savoye, B., Piper, D., and Syvitski, J.: The Var submarine sedimentary system: understanding Holocene sediment delivery processes and their importance to the geological record, Geol. Soc. Lond. Spec. Publ., 129, 146–166, 1998. a, b, c, d

Norton, K. and Vanacker, V.: Effects of terrain smoothing on topographic shielding correction factors for cosmogenic nuclide-derived estimates of basin-averaged denudation rates, Earth Surf. Proc. Land., 34, 145–154,, 2009. a

Petit, C., Migeon, S., and Coste, M.: Numerical models of continental and submarine erosion: Application to the northern Ligurian Margin, Earth Surf. Proc. Land., 425, 681–695, 2015. a

Petit, C., Goren, L., Rolland, Y., Bourlès, D., Braucher, R., Saillard, M., and Cassol, D.: Recent, climate-driven river incision rate fluctuations in the Mercantour crystalline massif, southern French Alps, Quaternary Sci. Rev., 165, 73–87,, 2017. a

Petit, C., Rolland, Y., Braucher, R., Bourlès, D., Guillou, V., and PetitPerrin, V.: River incision and migration deduced from 36Cl cosmic-ray exposure durations: The Clue de la Cerise gorge in southern French Alps, Geomorphology, 330, 81–88,, 2019. a, b, c, d

Phillips, J. and Slattery, M.: Sediment storage, sea level, and sediment delivery to the ocean by coastal plain rivers, Prog. Phys. Geogr., 30, 513–530,, 2006. a, b

Reinhardt, L. J., Bishop, P., Hoey, T. B., Dempster, T. J., and Sanderson, D. C. W.: Quantification of the transient response to base-level fall in a small mountain catchment: Sierra Nevada, southern Spain, J. Geophys. Res.-Earth, 112, F03S05,, 2007. a

Repka, J., Anderson, R., and Finkel, R.: Cosmogenic dating of fluvial terraces, Fremont River, Utah, Earth Planet. Sc. Lett., 152, 59–73, 1997. a

Rodrigo-Gamiz, M., Martinez-Ruiz, F., Rampen, S., Schouten, S., and Sinninghe Damsté, J.: Sea surface temperature variations in the western Mediterranean Sea over the last 20 kyr: A dual-organic proxy, Paleoceanography, 29, 87–98,, 2013. a, b

Rolland, T., Darnault, R., Braucher, R., Bourlès, D., Petit, C., and Bouissou, S.: Deglaciation history at the Alpine-Mediterranean transition (Argentera-Mercantour, SW Alps) from 10Be dating of moraines and glacially polished bedrock, Earth Surf. Proc. Land., 45, 393–410,, 2020. a, b, c, d, e

Rolland, Y., Petit, C., Saillard, M., Braucher, R., Bourlès, D., Darnault, R., and Cassol, D.: Inner gorges incision history: A proxy for deglaciation? Insights from Cosmic Ray Exposure dating (10Be and 36Cl) of river-polished surfaces (Tinée River, SW Alps, France), Earth Planet. Sc. Lett., 457, 271–281,, 2017. a, b, c, d

Romans, B., Castelltort, S., Covault, J., Fildani, A., and Walsh, J.: Environmental signal propagation in sedimentary systems across timescales, Earth-Sci. Rev., 153, 7–29,, 2016. a, b

Rouire, J., Bodelle, J., and Autran, A.: Carte Géologique de la France à 1/250000, Feuille no. 40 (Nice), BRGM Editions, Orléans, ISBN 036073598, 1980. a

Safran, E. B., Bierman, P. R., Aalto, R., Dunne, T., Whipple, K. X., and Caffee, M.: Erosion rates driven by channel network incision in the Bolivian Andes, Earth Surf. Proc. Land., 30, 1007–1024,, 2005. a

Saillard, M., Petit, C., Rolland, Y., Braucher, R., Bourlès, D., Zerathe, S., and Jourdon, A.: Late Quaternary incision rates in the Vésubie catchment area (Southern French Alps) from in situ-produced 36Cl cosmogenic nuclide dating: Tectonic and climatic implications, J. Geophys. Res.-Earth, 119, 1121–1135,, 2014. a, b, c, d, e, f

Salles, T.: Badlands: A parallel basin and landscape dynamics model, Software X, 5, 195–202,, 2016. a, b

Salles, T.: Badlands-model/badlands-Be: 10Be Production & Transport from Source-to-Sink (eSurf), Zenodo [data set],, 2023. a

Seguinot, J. and Delanay, I.: Last-glacial-cycle glacier erosion potential in the Alps, Earth Surf. Dynam., 9, 923–935,, 2021. a

Siame, L., Angelier, J., Godard, R. C. V., Derrieux, F., Bourlès, D., Braucher, R., Chang, K.-J., Chu, H.-T., and Lee, J.-C.: Erosion rates in an active orogen (NE-Taiwan): A confrontation of cosmogenic measurements with river suspended loads, Quatern. Geochronol., 6, 246–260,, 2011. a

Stone, J.: Air pressure and cosmogenic isotope production, J. Geophys. Res., 105, 753–759, 2000. a

Syvitski, J., Morehead, M., Bahr, D., and Mulder, T.: Estimating fluvial sediment transport: The rating parameters, Water Ressour., 36, 2747–2760, 2000. a, b, c

Syvitski, J., Angel, J., Saito, Y., Overeem, I., Vörösmarty, C., Wang, H., and Olago, D.: Earth's sediment cycle during the Anthropocene, Nat. Rev. Earth Environ., 3, 179–196,, 2022. a

Thran, A., East, M., Webster, J., Salles, T., and Petit, C.: The influence of carbonate platforms on thegeomorphological development of amixed carbonate-siliciclastic margin (Great Barrier Reef, Australia), Geochem. Geophy., Geosy., 21, e2020GC008915,, 2020. a

Vanacker, V., von Blanckenburg, F., Hewawasam, T., and Kubik, P.: Constraining landscape development of the Sri Lankan escarpment with cosmogenic nuclides in river sediment, Earth Planet. Sc. Lett., 253, 402–414,, 2007. a

von Blanckenburg, F.: The control mechanisms of erosion and weathering at basin scale from cosmogenic nuclides in river sediment, Earth Planet. Sc. Lett., 237, 462–479,, 2005. a

Waelbroeck, C., Labeyrie, L., Michel, L., Duplessy, J., McManus, J., Lambeck, K., Balbon, E., and Labracherie, M.: Sea-level and deep water temperature changes derived from benthic foraminifera isotopic records, Quaternary Sci. Rev., 21, 295–305,, 2002.  a, b

Wan, S., Li, A., Clift, P., and Jiang, H.: Development of the East Asian summer monsoon: Evidence from the sediment record in the South China Sea since 8.5 Ma, Palaeogeogr. Palaeocli. Palaeoecol., 241, 139–159,, 2006. a

Whipple, K. and Tucker, G.: Dynamics of the stream-power river incision model: Implications for height limits of mountain ranges, landscape response timescales, and research needs, J. Geophys. Res., 104, 17661–17674, 1999. a

Wickert, A.: Open-source modular solutions for flexural isostasy: gFlex v1.0, Geosci. Model Dev., 9, 997–1017,, 2016. a

Yanites, B., Tucker, G., and Anderson, R.: Numerical and analytical models of cosmogenic radionuclide dynamics in landslide-dominated drainage basins, J. Geophys. Res., 114, F01007,, 2009. a

Zerathe, S., Litty, C., Blard, P.-H., Delgado, F., Audin, L., and Carcaillet, J.: Cosmogenic 3He and 10Be denudation rates in the Central Andes: Comparison with a natural sediment trap over the last 18 ka, Earth Planet. Sc. Lett., 599, 117869,, 2022. a

Short summary
We present new tools in the landscape evolution model Badlands to simulate 10Be production, erosion and transport. These tools are applied to a source-to-sink system in the SW French Alps, where the model is calibrated. We propose a model that fits river incision rates and 10Be concentrations in sediments, and we show that 10Be in deep marine sediments is a signal with multiple contributions that cannot be easily interpreted in terms of climate forcing.