Articles | Volume 11, issue 4
Research article
02 Aug 2023
Research article |  | 02 Aug 2023

Stream laws in analog tectonic-landscape models

Riccardo Reitano, Romano Clementucci, Ethan M. Conrad, Fabio Corbi, Riccardo Lanari, Claudio Faccenna, and Chiara Bazzucchi

The interplay between tectonics and surface processes defines the evolution of mountain belts. However, correlating these processes through the evolution of natural orogens represents a scientific challenge. Analog models can be used for analyzing and interpreting the effect of such interaction. To fulfill this purpose it is necessary to understand how the imposed boundary conditions affect analog models' evolution in time and space. We use nine analog models characterized by different combinations of imposed regional slope and rainfall rates to investigate how surface processes respond to the presence of tectonically built topography (imposed slope) under different climatic conditions (rainfall rate). We show how the combination of these parameters controls the development of drainage networks and erosional processes. We quantify the morphological differences between experimental landscapes in terms of a proposed Se/R ratio, accounting for both observables and boundary conditions. We find few differences between analog models and natural prototypes in terms of parametrization of the detachment-limited stream power law. We observe a threshold in the development of channelization, modulated by a tradeoff between applied boundary conditions.

1 Introduction

Accurately interpreting the continuous interaction between tectonics and surface processes in mountain belts is one of the main challenges that geologists have faced in the last century. Many limitations exist due to the different depths and spatial and temporal scale at which crustal and mantle processes impact the surface. Thus, it is difficult to univocally interpret how different factors interact to build the present-day landscape. Analog models allow for a useful direct control of the evolution of the studied physical process (e.g., Reber et al.2020), overcoming many of these limitations. Tectonics, erosion, and sedimentation play an integrated role in the evolution of mountain belts with complex and poorly constrained feedbacks. During the last decades modelers analyzed these feedbacks, from the rejuvenation of streams (e.g., Schumm and Parker1973; Schumm and Rea1995) to the more complex evolution of whole orogenic systems (e.g., Bonnet2009; Graveleau and Dominguez2008; Guerit et al.2016; Lague et al.2003; Tejedor et al.2017; Viaplana-Muzas et al.2019; Reitano et al.2022b). However, all the previous analog modeling efforts are based on the robustness of the characterization of material used in the experiments (e.g., Graveleau et al.2011; Reitano et al.2020) and on the scaling to natural prototypes (e.g., Graveleau et al.2011; Paola et al.2009). While the mechanical, frictional, and erosional properties have been characterized empirically and analytically by different authors (Bonnet and Crave2003; Lague et al.2003; Graveleau et al.2011; Reitano et al.2020), the application of erosion laws to analog systems is still a matter of debate. In particular, a definition of the response of the analog materials to the applied boundary conditions and an understanding of the variables that modify this response are still missing. These characterizations become even more important considering that a perfect scaling between natural- and experimental-flow laws is currently missing, limiting the reach of analog studies to insights derived from qualitative process similarity (Paola et al.2009).

In this work we analyze how different boundary conditions affect the evolution of analog landscapes, in terms of tectonically built topography (imposed regional slope) and climate (rainfall rate). The methodologies implemented here allow us to isolate how these different a priori conditions control features like channelization, morphometrics, and incision rates. We then define the ranges of boundary conditions in which different morphological features develop (e.g., channels or diffusive processes). Finally, we discuss if and how the erosional parametrization implemented in nature (stream power law) applies to analog models. Analog models run in this study are representative of the slow-tectonic regions (e.g., very low uplift/erosion rate such as the Anti-Atlas of Morocco; Lanari et al.2022; Clementucci et al.2022) or passive margins such as the Western Ghats escarpment in the eastern margin of peninsular India (Mandal et al.2015; Wang and Willett2021) and the southern Australian escarpment (Godard et al.2019), where the uplift variation in space is negligible.

2 Experimental setup, scaling, and erosional laws

The analog material is a granular, water-saturated mixture made of 40 wt % of silica powder, 40 wt % of glass microbeads, and 20 wt % of PVC powder (Reitano et al.2020). We fill a 35×30×5 cm3 box with the material and place it on a reclinable table (Fig. S1 in the Supplement). The rainfall system is made of nozzles that provide a dense fog to trigger surface processes. The droplet size in the fog is lower than 100 µm to avoid rain splash erosion (e.g., Graveleau et al.2012). We apply three different rainfall rates in the models (9, 22, and 70 mm h−1) by controlling the number of implemented nozzles. The homogeneity of rainfall distribution has been tested as described in detail in the Supplement of Reitano et al. (2022b). The angle between the table and the horizontal (i.e., imposed regional slope) is also fixed at three different values, 10, 15, and 20. Both rainfall rate and imposed regional slope are kept constant throughout the experiments. In total, we investigate nine different imposed regional slope–rainfall rate configurations (Table S1 in the Supplement). We use a camera to capture top-view digital images of the experiment evolution and a high-resolution laser scanner to acquire DEMs at defined time steps. The vertical and horizontal resolutions of the laser are 0.07 and 0.05 mm, respectively.

Considering a length scaling factor L*=10-5/10-6 (1 cm = 1–10 km), a gravitation acceleration scaling factor g*=1, and a velocity scaling factor v*=104-105 (1 cm h−1= 0.1–1 cm yr−1), the time scaling factor t* is computed by (Reitano et al.2020, 2022b)

(1) t * = 4 L * v * ,

such that 1 min in the models corresponds to 5 to 50 kyr in nature.

The change in elevation (dz) over time (dt) of a point on the surface results from the competition between the rock uplift rate U and the erosion or sedimentation rate E:

(2) d z d t = U - E .

For fluvial networks that show exposed bedrock at the base, the erosion rate is typically expressed as a function of a “detachment-limited stream power law” (e.g., Goren et al.2014b, a; Howard1994; Howard and Kerby1983; Tucker and Whipple2002):

(3) E = κ Q m S n ,

where κ is the bedrock erodibility; Q is the channel discharge; S is the channel slope; and m and n are exponents accounting for channel geometry, basin hydrology, and erosion processes. Since Q is a function of drainage area, A, and the rainfall rate, P (Q=AP), Eq. (3) can be rewritten as

(4) E = K A m S n ,

where K=κPm is controlled by bedrock lithology, incision process, and climate (here rainfall rate). In a stream channel, Flint's law describes the power relationship between channel slope S and drainage area A (Flint1974):

(5) S = k s A - θ ,

where ks and θ are the channel steepness and the concavity indexes, respectively (e.g., Wobus et al.2006a; Kirby and Whipple2012; Hack1957, 1960). Merging Eqs. (4) and (5) and considering that θ=m/n, we can express the result in a logarithmic form:

(6) log 10 E = n log 10 k s + log 10 K .

Equation (6) is thus the equation of a line (y=mx+q), where n and log 10K are the slope and the intercept, respectively.

Figure 1Geomorphological and channel metrics of the performed experiments. The solid black lines indicate the median, while the bottom and top edges of the colored boxes indicate the 25th and 75th percentiles, respectively. The black whiskers outside the boxes cover the data points at the < 25th and > 75th percentiles that are not considered outliers, here indicated by black crosses. The color saturation in the boxes is related to the applied rainfall rate (less saturated, lower rainfall rate). The blue, red, and green boxes refer to models at imposed slopes of 10, 15, and 20, respectively.


From the DEMs we compute (i) topographic metrics such as basin slope (max and mode) and local relief (max and mode); (ii) channel metrics (ks and θ) and channel profiles (concavity index, θ, as well as steepness index, ks, has been extracted by power law regression between local slope and area – Flint's law); and (iii) eroded volumes, incision rate, and erosion maps. We describe the methodologies used for extracting the eroded volumes and incision rates in the Supplement. The abovementioned analyses are performed using the TopoToolbox package (Schwanghart and Scherler2014).

3 Results

Table S1 shows a list of the performed experiments, where the first and last two digits of the models' name represent the imposed regional slope and the imposed rainfall rate, respectively (e.g., mod1522: imposed regional slope of 15 and rainfall rate of 22 mm h−1). These model runs are interpreted highlighting the effect of boundary conditions on the evolution of geomorphic models. Results refer to 10 time steps that highlight key stages in the evolution of each model (Supplement) or to the entire experimental run (300 min). For every model at every time step, we defined and extracted all the basins forming and developing into the landscape, together with the main river trunk (max length) for every basin. For these basins we calculated the maximum surface slope (SSmax), the surface slope mode (SSmode), the maximum local relief (LRmax), the local relief mode (LRmode), and the concavity index (the last one from the main river in the basin). Increasing the imposed regional slope (from mod1009 to mod2070), the SSmax increases systematically from 30 to 80 (between the 25th and 75th percentile; Fig. 1a). The general increase in SSmax between models 1009, 1509, and 2009 (and the same for different rainfall rates) reflects the imposed slope increase since values are not normalized for the imposed slope. Still, increasing the imposed rainfall rate also results in an increase in SSmax. Conversely, SSmode decreases for models with the same regional slope but increasing rainfall rates (e.g., mod1509, mod1522, mod1570), except for mod1009 and mod1022 (Fig. 1b), whereas the mode increases for models with the same rainfall rate but different imposed regional slopes (e.g., mod1009, mod1509, mod2009). Models with a rainfall rate of 70 mm h−1 show the lowest values in SSmode yet the broadest range of values (between 5 and 15 from mod1070 to mod2070). The decreasing trend in the SSmode is also observed in the morphological expression of the model's surface (Fig. 2a). The models with a rainfall rate up to 22 mm h−1 show clear channelization during the final stage of the evolution, while models with rainfall rates equal to 70 mm h−1 shown little to no channelization in the final stage. Only mod2070 develops well-defined channels (branching, narrow channels, incision concentrated in valleys and channels).

Figure 2(a) DEMs and (b) erosion DEMs of the performed experiments. The erosion DEMs are obtained by computing the difference in elevation (Δz) of the same cell at consecutive times.


The concavity index ranges from 0.1 to 0.4 (between the 25th and 75th percentile; Fig. 1c). Models with the highest rainfall rate typically show the highest concavity values, correspondingly with the highest variability (mod1070, mod1570, mod2070; Fig. 1c). The local relief was extracted using a moving window with a radius of 10 mm. The maximum local relief LRmax and local relief mode LRmode show a similar pattern with respect to the surface slope (Fig. 1d and e). LRmax increases from mod1009 to mod2070, while the LRmode decreases for models with rainfall rates equal 70 mm h−1 with respect to models with the same imposed regional slope but different rainfall rates. The same information can be deduced by visual inspection of DEMs (Fig. 2a).

Figure 3(a) Amount of eroded material for the performed experiments (the color coding follows Fig. 1). (b) Maximum incision rate as a function of the imposed boundary conditions.


The amount of eroded material ranges from 0.2×106 (mod1009) to 2.3×106 mm3 (mod2070) in the last stage of the experiments' evolution (Fig. 3a). Both increased rainfall rate (e.g., mod1009, mod1022, mod1070) and increased imposed regional slope (e.g., mod1022, mod1522, mod2022) result in higher amounts of eroded material. All models show an initial phase where the material flux is highest with a later phase of decline (Reitano et al.2020), tending toward stability. For example, in mod2070, 60 min of experimental time is enough to erode 1.4×106 mm3 of material, while in the next 240 min only an additional 0.9×106 mm3 of material is eroded. This different behavior is most apparent in models with high rainfall rates (≥22 mm h−1) and slopes (≥15). The maximum incision rate increases with the imposed regional slope and rainfall rate (Fig. 3b), from <5 mm h−1 (mod1009) to ca. 55 mm h−1 (mod2070).

We extract values for E and ks of four main channels for each time step and for every model (40 channels per model, total = 360). From these values, we obtain n and K by linear regression (Sect. 2). Despite the low R2 values (0.01–0.28; Fig. 4), n ranges between −0.18 and 0.14 with K values between 0.77 and 23.51 mm1−2m h−1. Values of K increase as a function of the slope and rainfall rate, while n does not show a clear trend. On the right-hand side of Fig. 4, we plot data for models with the same imposed regional slope but different rainfall rates. Models with the same regional slope show a gradually increasing ks and incision rate in response to increased rainfall rates (Fig. 4b). Interestingly, estimates of K gradually increase in models with imposed regional slopes of 10 to 20 (Fig. 4b).

4 Discussion

4.1 Type of erosion as a function of boundary conditions

Our analog models are controlled exclusively by imposed regional slope and rainfall rate, as no other boundary conditions are applied (e.g., vertical uplift or horizontal advection of material). Higher rainfall rates (70 mm h−1 in this work) tend to inhibit the development of a channelized and branching channel network in favor of more diffusive and mass wasting processes. This trend can be deduced by analyzing the DEMs of mod1070 and mod1570 (Fig. 2a) or simply by noting the diffuse nature of erosion under high-rainfall-rate conditions (Fig. 2b). In this case, the water cannot channelize into discrete features since the amount of water is higher than the channelization capacity of the model, mainly controlled by the imposed slope. Indeed, more diffusive processes develop over a greater area, accounting for the lack of channelization. Even if these models are mainly controlled by diffusive processes, mod2070 exhibits more channelization than mod1070 and mod1570, but, similarly, mod1570 shows slightly more channelization than mod1070. These observations suggest that, in terms of channelization, the higher the slope, the more effectively a system responds to high rainfall rates. Considering the relationship between channelization and boundary conditions, the results of our experiments suggest that low imposed regional slope (mod1009, mod1022) or low rainfall rate with an average imposed regional slope (mod1509) results in a channelization characterized by low incision values (<25–30 mm in the final stage). It thus appears that a threshold exists in rainfall rate for the development of channelized networks, modulated by the slope over which erosion acts. For example, a rate of 70 mm h−1 (mod1070 and mod1570) is too high for a proper channel network to develop, while a higher imposed regional slope (mod2070) provides sufficient potential energy for the system to channelize (e.g., Burbank and Anderson2011). Thus, the tradeoff between imposed regional slope and rainfall rate controls the channelization. For higher rainfall rate (mod1070 and mod1570), a sheetlike runoff lowers the model slope homogeneously (Fig. 1b). Furthermore, both the SSmode and the LRmode drop with respect to the models with the same imposed regional slope but a lower rainfall rate. From these observations, we argue that at a high rainfall rate, channelization is subordinate to diffusive processes (controlled by ridge stability) in the final stages of model evolution. In landscapes where incision is a function of the detachment of particles from the riverbed, the erosion rate is proportional to the shear stress (e.g., Whipple and Tucker1999; Yanites et al.2010). Higher rainfall rates (i.e., higher water discharge) increase the effective shear stress on riverbeds (Thoman and Niezgoda2008). Since water discharge also increases the channel width (e.g., Shibata and Ito2014; Wobus et al.2006b), for high water discharges (i.e., rainfall rate) the shear stress can be distributed over time over wide and flat surfaces instead of concentrating in valleys (Lamb et al.2015, and references therein). At high rainfall rates our models show incipient channelization in the initial stages (Fig. 2b). Channelization is lost in favor of more diffusive processes during late stages of model evolution. High rainfall rates (i.e., high water discharge) lead to a higher sediment supply that can widen channels (Finnegan et al.2007; Johnson and Whipple2010), eventually erasing their morphological expression. Finally, we must address the erosional threshold defined in the works of Lague et al. (2003) and Graveleau et al. (2011). This threshold must be overcome before significant erosion and transport occurs and specifically apply to models with low imposed regional slopes, which may thus lead to even less channelization.

Figure 4(a) Logarithm of the incision rate over the logarithm of the steepness index ks for all the models. The imposed regional slope increases from top to bottom. Rainfall rate increases from left to right. Every plot shows 4 channels at every time step (40 channels in total; colored dots). The linear regression is shown by the red line. Values related to the linear regression (n and K) are shown in the bottom right corner of every plot, together with R2 (units for K are mm1−2m h−1). (b) Same data as (a), but plotted for every slope without differentiating for rainfall rate.


Interestingly, effective channelization does not affect the volume of eroded material, which increases with the imposed regional slope and rainfall rate (e.g., the erosion flux is higher in mod1570 than mod1522 and mod1070).

We plotted the volume of eroded material at every time step normalized by the imposed slope and the rainfall rate (Vnorm; Fig. 5a). We observe a transient behavior for all the models during the first 100–150 min, indicating the adjustment of the model to the applied boundary conditions. After this transient phase, a steady phase develops, where the system tends toward equilibrium in terms of eroded material. Ideally, due to how Vnorm is computed, all lines should collapse into one. This is not the case because analog models are intrinsically stochastic; thus unconstrained variability between the models causes deviation in the calculation of Vnorm.

Figure 5(a) Evolution of the normalized volume (Vnorm) by imposed slope and rainfall rate for all the models. (b) Se/R in controlling the development of channelization (No chan.: no channelization develops; Chan. low inc.: development of channelization with low incision; Chan.: channelization develops).


Figure 5b clearly shows the relationship between the imposed boundary conditions – imposed slope (Se) and rainfall rate (R) – and the evolution of the models. We group the models into three main landscapes, following the rainfall rate.

  • i.

    When Se/R is lower than 0.5 (dark gray, Fig. 5b), channelization is less important than diffusive processes in controlling the erosion. Models with high rainfall rate (mod1070, mod1570, mod2070) or low rainfall rate and low imposed regional slope (mod1009) show low channelization with respect to the other models. Even if mod1009 shows channelization in its DEMs (Fig. 2a), the level of incision (Fig. 2b) shows that erosion is broadly distributed and not concentrated in channels or valleys so that detachment-limited conditions for channel incision do not apply.

  • ii.

    When Se/R is comprised between 0.5 and 1 (medium light gray, Fig. 5b), channelization processes are present and are responsible for the erosion. Still, channels and valleys are wide (ca. 6 cm) with respect to models where Se/R is higher than 1 (mod1009, mod1509, mod2009; light gray in Fig. 5b).

  • iii.

    In these latter cases, channelization is the main process in controlling erosion. The maximum incision rate is lower than 20 mm h−1 (Fig. 3b), and the eroded volumes range between 0.1 and 0.75×106 mm3.

The only exception for this trend is mod2070, which should fall into the “no channelization” landscape. Mod2070 is classified as “channelization with low incision” landscapes instead. We previously described how the extremely high R is partly balanced by a very high Se so that even if the rainfall rate guides the model toward the “no channelization” landscape, the imposed slope promotes enough potential energy to allow the formation and the development of channels, even with low incision (Fig. 2b). It is then clear that a threshold for channelization exists. This threshold is modulated mainly by the rainfall rate, but also by the imposed regional slope. We propose that Se/R can be used as a threshold parameter to define if and how channelization processes may develop in analog models as an indicator for the applied boundary conditions.

The applied boundary conditions are needed to constrain the reliability of analog models in understanding the evolution of natural prototypes. The possibility of applying scaled data from analog models to natural landscapes is reliant on the boundary conditions. For example, the extremely high rainfall rate used in this work is difficult to compare with nature.

4.2 Geomorphological metrics in analog models and natural prototypes

The concavity index θ of the selected channels (for every model) is usually lower than 0.4 (values between the 25th and 75th percentile; Fig. 1c). The lower values of θ are related to straight longitudinal channel profiles (Fig. S2) that are established over the course of a model run (Duvall et al.2004; Reitano et al.2020; Whipple and Tucker1999). Still, θ is comparable between analog models and natural prototypes (within 0.1; Reitano et al.2020). However, we do not normalize the steepness index ks by a reference concavity index ksn, as is usually done in the literature (e.g., Cyr et al.2010; DiBiase et al.2010; Kirby and Whipple2012; Lanari et al.2020; Tucker and Whipple2002). This approach is used so that the steepness index directly reflects channel characteristics, without modifying the data with a normalization parameter. If n=1, there exists a linear relationship between the incision rate and ks; i.e., an increase in the steepness index increases the incision rate and vice versa (Eq. 6). Previous works that computed values of n for different natural landscapes show that n is typically greater than or equals 1 (e.g., DiBiase et al.2010; Harel et al.2016; Ouimet et al.2009), meaning that incision rates are extremely sensitive to variation in ks of rivers in active-tectonic settings (Kirby and Whipple2012). In slow-tectonic settings, the incision rates show a lower sensitivity to the steepness index (Clementucci et al.2022; Olivetti et al.2016). For this case, we show that n has values generally lower than 1, and even negative values (Fig. 4a). Since n is dimensionless, a one-to-one comparison between models and nature is possible. The n values from our analog models are closer to estimates of n values from slow-tectonic settings compared to active domains, where n is greater than 1 (Figs. 4b and 6Clementucci et al.2022; Kirby and Whipple2012). Similar to slow-tectonic landscapes in nature, incision rates of analog models are less sensitive to variations in ks. Moreover, the low R2 (Fig. 4) testifies to a poor relationship between incision rate and ks (n<1). Nevertheless, we observe a trend in the relationship between incision rate and ks as a function of the applied rainfall rate (Fig. 4b). For higher rainfall rates, the incision rate increases as expected, but only with a slight increase in ks. This consideration supports the fact that ks plays a minor role in controlling the incision rate, like slow-tectonic natural domains.

Figure 6Logarithm of the incision rate over the logarithm of the steepness index ksn for the selected natural prototypes (θref=0.45). The linear best fit is shown by the red line. Values related to the linear regression (n and K) are shown in every plot, together with R2.


Many works in the literature compute values of incision rate and ks (or ksn) for tectonically active regions. The models presented here are similar to natural prototypes where tectonics are absent or subordinate with respect to surface processes. To show this, we collected values of basin-wide denudation rate and ksn (ks normalized by the same concavity index θref=0.45) for natural prototypes that show no or almost no active tectonics, divided according to the main rock types. By using this approach, we are quantifying ksn over the denudation rate relationship and filtering it by the erodibility or rock types, as in the presented analog landscapes there is no spatial variation in bedrock erodibility. Then, we computed values of n and K (Fig. 6). Estimates of n from the analog models present a good match with data from the selected natural prototypes with slow tectonics (Fig. 6), showing values between −0.06 and 0.32 for analog models (Fig. 4) and −0.12 and 0.62 for natural systems (Fig. 6). This relationship is apparent despite the low R2 value from linear regression of data from analog and natural cases. Importantly, the difference in the range of n values can likely be decreased by filtering out the basins affected by higher uplift rates within the slow-tectonic settings (e.g., granite-dominated basins; Clementucci et al.2022). Values of K are more difficult to compare between models and nature, due to the complex dimensions that characterize the parameter (mm1−2m h−1), which require not only the application of a scaling factor but also the definition and comparison of the m variable. Since K is mainly a function of the lithology and climate (by considering uniform uplift rates; Sect. 2), and we use the same homogenous material for all the models, K only highlights the effect of the rainfall rate, and there is a nearly perfect linearity between K and the rainfall rate. For example, from 9 to 22 mm h−1 and from 22 to 70 mm h−1, the rainfall rate increases by a factor of 2.4 and 3.18, respectively. This is also apparent for K, with the only exception being the experiment with both moderate imposed slope and rainfall rate (mod1522). Interestingly, there is also a linear relationship between K and the slope (for example, in Fig. 4 when moving from high to low imposed slope). Although this relationship is still poorly constrained (Clementucci et al.2022; Harel et al.2016; Stock and Montgomery1999), we speculate that K is a function of not only lithologies and climate but also of the regional slope and, consequently, of the integrated topographic response to tectonic rates, as has been observed for the natural prototypes (Peifer et al.2021).

5 Conclusions

We investigate the role of boundary conditions (imposed regional slope and rainfall rate) on the morphological evolution of nine different analog models. Some of the applied boundary conditions (i.e., the very high rainfall rate at 70 mm h−1) should be considered as end-members used to fully cover the behavior of analog models under different conditions. The models systematically test various configurations of boundary conditions. We analyze how the stream power law used in natural landscapes apply to analog models, and we find that a threshold exists for the development of channelization in terms of boundary conditions (Se/R ratio). This threshold is the result of a tradeoff between the imposed regional slope (tectonically driven topography) and the rainfall rate (climate). We find that combining imposed regional slope and rainfall rate results in three possible changes to scenarios.

  • Se/R lower than 0.5. Diffusive processes are dominant with respect to channelization so that detachment-limited conditions for channel incision do not apply (mod1070, mod1570, mod2070).

  • Se/R greater than 0.5 and lower than 1. Channelization controls the erosion of the surface, with wide channels and valleys (ca. 6 cm).

  • Se/R higher than 1 (mod1009, mod1509, mod2009). Channelization processes are the main factors controlling the morphologies of the analog landscapes.

In summary, high rainfall rates (70 mm h−1) inhibit the development of channels, but not if coupled with high imposed regional slope (20). Similarly, low imposed regional slope and low rainfall rate (10 and 9 to 22 mm h−1, respectively) develop channelization characterized by low incision. When looking at the parameterization of the stream power law for erosion, we find that analog models are described by parameters whose values slightly differ from the natural prototypes (e.g., n exponent) and by parameters that are tuned by applied boundary conditions, imposed regional slope, and rainfall rate (e.g., K coefficient). Given these findings, we propose that the insights, limitations, and range of applied boundary conditions discussed here should be considered in the future interpretation of interactions between tectonics and surface process.

Data availability

Data are open-access and have been uploaded at the GFZ Data Service (, Reitano et al.2022a).


The supplement related to this article is available online at:

Author contributions

RR and RC proposed the original idea. RR, CF, and EMC designed the experiments, and RR and CB carried them out. RR and FC extracted the metrics used for the model analysis, performed by all authors. RC and RL collected data from natural prototypes. Interpretation of results, writing, reviewing, and editing were performed by all authors.

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.


In Figs. 2 and 3 we used the perceptually uniform colormap, Roma and Davos, by Fabio Crameri. This work was funded by the “Fund for University Departments of Excellence 2023–2027”. We thank the associate editor Greg Hancock and the anonymous reviewer for the comments and suggestions, which increased the quality of the work. We also thank Francesca Funiciello for her constant support.

Review statement

This paper was edited by Greg Hancock and reviewed by one anonymous referee.


Bonnet, S.: Shrinking and splitting of drainage basins in orogenic landscapes from the migration of the main drainage divide, Nat. Geosci., 2, 766–771,, 2009. a

Bonnet, S. and Crave, A.: Landscape response to climate change: Insights from experimental modeling and implications for tectonic versus climatic uplift of topography, Geology, 31, 123–126,<0123:LRTCCI>2.0.CO;2, 2003. a

Burbank, D. W. and Anderson, R. S.: Tectonic Geomorphology, John Wiley & Sons, Ltd,, 2011. a

Clementucci, R., Ballato, P., Siame, L. L., Faccenna, C., Yaaqoub, A., Essaifi, A., Leanni, L., and Guillou, V.: Lithological control on topographic relief evolution in a slow tectonic setting (Anti-Atlas, Morocco), Earth Planet. Sc. Lett., 596, 117788,, 2022. a, b, c, d, e

Cyr, A. J., Granger, D. E., Olivetti, V., and Molin, P.: Quantifying rock uplift rates using channel steepness and cosmogenic nuclide–determined erosion rates: Examples from northern and southern Italy, Lithosphere, 2, 188–198,, 2010. a

DiBiase, R. A., Whipple, K. X., Heimsath, A. M., and Ouimet, W. B.: Landscape form and millennial erosion rates in the San Gabriel Mountains, CA, Earth Planet. Sc. Lett., 289, 134–144,, 2010. a, b

Duvall, A., Kirby, E., and Burbank, D.: Tectonic and lithologic controls on bedrock channel profiles and processes in coastal California, J. Geophys. Res.-Earth, 109, F03002,, 2004. a

Finnegan, N. J., Sklar, L. S., and Fuller, T. K.: Interplay of sediment supply, river incision, and channel morphology revealed by the transient evolution of an experimental bedrock channel, J. Geophys. Res.-Earth, 112, F03S11,, 2007. a

Flint, J.-J.: Stream gradient as a function of order, magnitude, and discharge, Water Resour. Res., 10, 969–973,, 1974. a

Godard, V., Dosseto, A., Fleury, J., Bellier, O., Siame, L., and ASTER Team: Transient landscape dynamics across the Southeastern Australian Escarpment, Earth Planet. Sc. Lett., 506, 397–406,, 2019. a

Goren, L., Fox, M., and Willett, S. D.: Tectonics from fluvial topography using formal linear inversion: Theory and applications to the Inyo Mountains, California, J. Geophys. Res.-Earth, 119, 1651–1681,, 2014a. a

Goren, L., Willett, S. D., Herman, F., and Braun, J.: Coupled numerical–analytical approach to landscape evolution modeling, Earth Surf. Proc. Land., 39, 522–545,, 2014b. a

Graveleau, F. and Dominguez, S.: Analogue modelling of the interaction between tectonics, erosion and sedimentation in foreland thrust belts, Comptes Rendus Geoscience, 340, 324–333,, 2008. a

Graveleau, F., Hurtrez, J.-E., Dominguez, S., and Malavieille, J.: A new experimental material for modeling relief dynamics and interactions between tectonics and surface processes, Tectonophysics, 513, 68–87,, 2011. a, b, c, d

Graveleau, F., Malavieille, J., and Dominguez, S.: Experimental modelling of orogenic wedges: A review, Tectonophysics, 538-540, 1–66,, 2012. a

Guerit, L., Dominguez, S., Malavieille, J., and Castelltort, S.: Deformation of an experimental drainage network in oblique collision, Tectonophysics, 693, 210–222,, 2016. a

Hack, J. T.: Studies of longitudinal stream profiles in Virginia and Maryland, USGS,, 1957. a

Hack, J. T.: Interpretation of erosional topography in humid temperate regions, Am. J. Science, 258, 80–97, 1960. a

Harel, M.-A., Mudd, S., and Attal, M.: Global analysis of the stream power law parameters based on worldwide 10Be denudation rates, Geomorphology, 268, 184–196,, 2016. a, b

Howard, A. D.: A detachment-limited model of drainage basin evolution, Water Resour. Res., 30, 2261–2285,, 1994. a

Howard, A. D. and Kerby, G.: Channel changes in badlands, Geol. Soc. Am. Bull., 94, 739–752,<739:CCIB>2.0.CO;2, 1983. a

Johnson, J. P. and Whipple, K. X.: Evaluating the controls of shear stress, sediment supply, alluvial cover, and channel morphology on experimental bedrock incision rate, J. Geophys. Res.-Earth, 115, F02018,, 2010. a

Kirby, E. and Whipple, K. X.: Expression of active tectonics in erosional landscapes, J. Struct. Geol., 44, 54–75,, 2012. a, b, c, d

Lague, D., Crave, A., and Davy, P.: Laboratory experiments simulating the geomorphic response to tectonic uplift, J. Geophys. Res.-Solid, 108, ETG 3-1–ETG 3-20,, 2003. a, b, c

Lamb, M. P., Finnegan, N. J., Scheingross, J. S., and Sklar, L. S.: New insights into the mechanics of fluvial bedrock erosion through flume experiments and theory, Geomorphology, 244, 33–55,, 2015. a

Lanari, R., Fellin, M. G., Faccenna, C., Balestrieri, M. L., Pazzaglia, F. J., Youbi, N., and Maden, C.: Exhumation and surface evolution of the western high atlas and surrounding regions as constrained by low-temperature thermochronology, Tectonics, 39, e2019TC005562,, 2020. a

Lanari, R., Reitano, R., Giachetta, E., Pazzaglia, F., Clementucci, R., Faccenna, C., and Fellin, M.: Is the Anti-Atlas of Morocco still uplifting?, J. Afr. Earth Sci., 188, 104481,, 2022. a

Mandal, S. K., Lupker, M., Burg, J.-P., Valla, P. G., 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

Olivetti, V., Godard, V., Bellier, O., Team, A., et al.: Cenozoic rejuvenation events of Massif Central topography (France): Insights from cosmogenic denudation rates and river profiles, Earth Planet. Sc. Lett., 444, 179–191,, 2016. a

Ouimet, W. B., Whipple, K. X., and Granger, D. E.: Beyond threshold hillslopes: Channel adjustment to base-level fall in tectonically active mountain ranges, Geology, 37, 579–582,, 2009. a

Paola, C., Straub, K., Mohrig, D., and Reinhardt, L.: The “unreasonable effectiveness” of stratigraphic and geomorphic experiments, Earth-Sci. Rev., 97, 1–43,, 2009. a, b

Peifer, D., Persano, C., Hurst, M. D., Bishop, P., and Fabel, D.: Growing topography due to contrasting rock types in a tectonically dead landscape, Earth Surf. Dynam., 9, 167–181,, 2021. a

Reber, J. E., Cooke, M. L., and Dooley, T. P.: What model material to use? A Review on rock analogs for structural geology and tectonics, Earth-Sci. Rev., 202, 103107,, 2020. a

Reitano, R., Faccenna, C., Funiciello, F., Corbi, F., and Willett, S. D.: Erosional response of granular material in landscape models, Earth Surf. Dynam., 8, 973–993,, 2020. a, b, c, d, e, f, g

Reitano, R., Clementucci, R., Conrad, E. M., Corbi, F., Lanari, R., Faccenna, C., and Bazzucchi, C.: Raw data (pictures, DEMs, .mat files) about analogue landscapes evolution, GFZ Data Services [data set],, 2022a. a

Reitano, R., Faccenna, C., Funiciello, F., Corbi, F., Sternai, P., Willett, S. D., Sembroni, A., and Lanari, R.: Sediment recycling and the evolution of analog orogenic wedges, Tectonics, 41, e2021TC006951,, 2022b. a, b, c

Schumm, S. and Rea, D. K.: Sediment yield from disturbed earth systems, Geology, 23, 391–394,<0391:SYFDES>2.3.CO;2, 1995. a

Schumm, S. A. and Parker, R. S.: Implications of complex response of drainage systems for Quaternary alluvial stratigraphy, Nat. Phys. Sci., 243, 99–100,, 1973. a

Schwanghart, W. and Scherler, D.: TopoToolbox 2 – MATLAB-based software for topographic analysis and modeling in Earth surface sciences, Earth Surf. Dynam., 2, 1–7,, 2014.  a

Shibata, K. and Ito, M.: Relationships of bankfull channel width and discharge parameters for modern fluvial systems in the Japanese Islands, Geomorphology, 214, 97–113,, 2014. a

Stock, J. D. and Montgomery, D. R.: Geologic constraints on bedrock river incision using the stream power law, J. Geophys. Res.-Solid, 104, 4983–4993,, 1999. a

Tejedor, A., Singh, A., Zaliapin, I., Densmore, A. L., and Foufoula-Georgiou, E.: Scale-dependent erosional patterns in steady-state and transient-state landscapes, Sci. Adv., 3, e1701683,, 2017. a

Thoman, R. W. and Niezgoda, S. L.: Determining erodibility, critical shear stress, and allowable discharge estimates for cohesive channels: Case study in the Powder River basin of Wyoming, J. Hydraul. Eng., 134, 1677–1687,, 2008. a

Tucker, G. and Whipple, K.: Topographic outcomes predicted by stream erosion models: Sensitivity analysis and intermodel comparison, J. Geophys. Res.-Solid, 107, ETG 1-1–ETG 1-16,, 2002. a, b

Viaplana-Muzas, M., Babault, J., Dominguez, S., Van Den Driessche, J., and Legrand, X.: Modelling of drainage dynamics influence on sediment routing system in a fold-and-thrust belt, Basin Res., 31, 290–310,, 2019. a

Wang, Y. and Willett, S. D.: Escarpment retreat rates derived from detrital cosmogenic nuclide concentrations, Earth Surf. Dynam., 9, 1301–1322,, 2021. a

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

Wobus, C., Whipple, K. X., Kirby, E., Snyder, N., Johnson, J., Spyropolou, K., Crosby, B., and Sheehan, D.: Tectonics from topography: Procedures, promise, and pitfalls, in: Tectonics, Climate, and Landscape Evolution, Geological Society of America,, 2006a. a

Wobus, C. W., Tucker, G. E., and Anderson, R. S.: Self-formed bedrock channels, Geophys. Res. Lett., 33, L18408,, 2006b. a

Yanites, B. J., Tucker, G. E., Mueller, K. J., Chen, Y.-G., Wilcox, T., Huang, S.-Y., and Shi, K.-W.: Incision and channel morphology across active structures along the Peikang River, central Taiwan: Implications for the importance of channel width, GSA Bull., 122, 1192–1208,, 2010. a

Short summary
Tectonics and surface processes work together in shaping orogens through their evolution. Laboratory models are used to overcome some limitations of direct observations since they allow for continuous and detailed analysis of analog orogens. We use a rectangular box filled with an analog material made of granular materials to study how erosional laws apply and how erosion affects the analog landscape as a function of the applied boundary conditions (regional slope and rainfall rate).