Graphically interpreting how incision thresholds influence topographic and scaling properties of modeled landscapes
We examine the influence of incision thresholds on topographic and scaling properties of landscapes that follow a landscape evolution model (LEM) with terms for stream-power incision, linear diffusion, and uniform uplift. Our analysis uses three main tools. First, we examine the graphical behavior of theoretical relationships between curvature and the steepness index (which depends on drainage area and slope). These relationships plot as straight lines for the case of steady-state landscapes that follow the LEM. These lines have slopes and intercepts that provide estimates of landscape characteristic scales. Such lines can be viewed as counterparts of slope–area relationships, which follow power laws in detachment-limited landscapes but not in landscapes with diffusion. We illustrate the response of these curvature–steepness index lines to changes in the values of parameters. Second, we define a Péclet number that quantifies the competition between incision and diffusion, while taking the incision threshold into account. We examine how this Péclet number captures the influence of the incision threshold on the degree of landscape dissection. Third, we characterize the influence of the incision threshold using a ratio between it and the steepness index. This ratio is a dimensionless number in the case of the LEM that we use and reflects the fraction by which the incision rate is reduced due to the incision threshold; in this way, it quantifies the relative influence of the incision threshold across a landscape. These three tools can be used together to graphically illustrate how topography and process competition respond to incision thresholds.
Processes that shape landscapes leave topographic signatures, which can often be visualized by plotting different topographic metrics against one another. An example is the relationship between river gradient and drainage area, which has been used to analyze landscapes and river profiles, as well as to diagnose the processes that shape them (e.g., Montgomery and Foufoula-Georgiou, 1993; Howard, 1994; Montgomery and Dietrich, 1994; Dietrich et al., 2003). For example, the stream-power incision model predicts that if tectonics, climate, and rock properties are uniform, then bedrock rivers should approach a steady state in which their gradient scales as a power law of drainage area (e.g., Tucker, 2004; Lague, 2014). This power-law scaling implies that river gradient data should plot as a straight line against drainage area data on logarithmic axes. The properties of this line can give estimates of properties of the landscape; e.g., its slope gives the concavity index (Whipple, 2004). Plotting synthetic topographic data from landscape evolution models (LEMs) helps to illustrate the effects of different model formulations or parameterizations. For example, including a threshold in the incision term of an LEM affects the resulting slope–area line (e.g., Tucker, 2004; Lague et al., 2005; Deal et al., 2018).
In the case of landscapes that are influenced by diffusion, topographic slope does not scale as a power function of drainage area (e.g., Howard, 1994). Thus, slope and area data from these landscapes do not plot as straight lines. In Theodoratos et al. (2018), we presented a counterpart relationship for the case of landscapes produced by an LEM that includes linear diffusion (along with stream-power incision and uplift). This relationship predicts that in steady state, curvature and the steepness index (which depends on drainage area and slope; e.g., Whipple, 2001) plot as a straight line against each other on linear (i.e., non-logarithmic) axes. The slope and intercept of this line depend on characteristic scales of length and height of the landscape, which in turn depend on the relative strengths of the processes that shape it. Thus, this relationship predicts a link between topographic and scaling properties of landscapes that follow the LEM.
Here, we demonstrate an example of the explanatory power of plots of the curvature–steepness index relationship. Our example shows that these plots can visualize topographic and scaling effects of incision thresholds. Incision thresholds can markedly influence erosion, as shown by numerous studies. For instance, incision thresholds can influence the relationship between river gradient and the uplift rate (e.g., Snyder et al., 2003), the dependence of long-term erosion rates on the average, variability, and duration of precipitation events (e.g., DiBiase and Whipple, 2011; Scherler et al., 2017), and the dynamics of migrating knickpoints (e.g., Lague, 2014). Here, we are not further elaborating on the insights of these studies. Instead, we focus on the effects of incision thresholds on the competition between incision and diffusion and on the topographic and scaling properties of landscapes reflecting this competition. The topographic and scaling effects that we examine have been studied before (e.g., Montgomery and Dietrich, 1992; Howard, 1994; Tucker, 2004; Perron et al., 2008). Here, however, we present a novel, purely graphical method to identify, quantify, and interpret these effects based on the relationship between curvature and the steepness index.
In Theodoratos et al. (2018), we dimensionally analyzed a frequently used LEM with terms for uplift, linear diffusion, and stream-power incision without an incision threshold. In Theodoratos and Kirchner (2020), we added an incision threshold to this LEM and dimensionally analyzed it. Here, we summarize the definitions of characteristic scales and dimensionless numbers that emerged from the dimensional analyses of these two LEMs in Sect. 2. Then, in Sect. 3, we show that these characteristic scales and dimensionless numbers have a geomorphologic meaning that can be expressed graphically using plots of curvature vs. the steepness index. The graphical explanatory power of these plots is further highlighted by comparing plots of LEMs with and without an incision threshold (Figs. 1 and 2).
2.1 Governing equations
The LEM without an incision threshold follows the governing equation (e.g., Howard, 1994; Dietrich et al., 2003):
This equation gives the rate of elevation change as the sum of three terms, namely (a) stream-power incision , where K is the incision coefficient, A is drainage area, and |∇z| is topographic slope; (b) linear diffusion D∇2z, where D is the diffusion coefficient and ∇2z is the Laplacian curvature; and (c) the uplift rate U. We assume that Eq. (1) has base dimensions of horizontal length L, height H (which we treat as dimensionally distinct from L), and time T. All quantities in Eq. (1) have dimensions that are combinations of L, H, and/or T, which we show in Table 1.
Note that the incision term is a special case of the more general incision term KAm(|∇z|)n. As we explained in Theodoratos et al. (2018), dimensional analysis of an LEM with generic exponents m and n would lead to equivalent results as the analysis of Eq. (1), but these results would be expressed with much more complicated mathematical formulas. Therefore, in Theodoratos et al. (2018) we focused on the case of exponents m=0.5 and n=1, and we presented the results for generic m and n in an appendix. Likewise, in the current study, the main presentation focuses on the case of m=0.5 and n=1, and in Appendix A we demonstrate that our graphical method is also valid for the case of generic exponents m and n.
Following Perron et al. (2008), we can add an incision threshold to the LEM by recasting the incision term as , where θ is the incision threshold. This formulation assumes that the incision rate is reduced everywhere by the constant quantity Kθ. The LEM examined here is based on the assumption that sediment transport is detachment-limited. Thus, it does not include deposition, and negative incision rates would not be meaningful. Therefore, the incision term is set to zero where the term would be negative, i.e., where , and the governing equation becomes
The incision threshold θ has the same dimensions as , i.e., dimensions of H.
Equation (2) assumes that precipitation rates are constant in time and uniform in space, and it incorporates climatic effects into the incision coefficient K. Other LEMs use stochastic precipitation to drive their incision terms (e.g., Tucker, 2004, Whipple, 2004; Lague et al., 2005; DiBiase and Whipple, 2011; Deal et al., 2018). The incision thresholds of these LEMs define limiting values of shear stress or stream power, below which no incision occurs. At any given location in the landscape, these limiting values might be exceeded during some stochastic events and not exceeded during other events, depending on their intensities. By contrast, in the case of the LEM that we examine, the assumption of constant and uniform precipitation implies that any given combination of drainage area A and slope |∇z| would lead to the same value of stream power (or shear stress) for any storm event (as all events would be equal), and this value of stream power would either be above or below the incision threshold. In this idealized case, defining a topographic threshold based on is exactly equivalent to defining a threshold of stream power (or shear stress).
We acknowledge that the LEMs with stochastic precipitation allow much more realistic integration of incision rates over time compared to the LEM that we examine here. Therefore, these LEMs are more appropriate for studying the influence of incision thresholds on erosion rates compared to the LEM that we use. However, our study has a different focus. Our study focuses on how the incision threshold θ influences topographic and scaling properties of landscapes and on how this influence can be graphically expressed with curvature–steepness index lines. For these tasks, the simplified formulation of the incision term of Eq. (2) is more practical. It may be possible in future work to extend this approach to include incision thresholds driven by stochastic precipitation.
2.2 Characteristic scales
The two governing equations (Eqs. 1 and 2) can be non-dimensionalized using characteristic scales of length, height, and time (lc, hc, and tc) defined as (Theodoratos et al., 2018; Theodoratos and Kirchner, 2020)
We summarize these and other definitions for this presentation in Table 2. The characteristic scales lc, hc, and tc can be viewed as intrinsic properties of a landscape in the sense that they depend exclusively on the values of the parameters K, D, and U and not on extensive properties of the landscape such as the size of its domain or its maximum relief. We present geomorphologic interpretations of these characteristic scales in Sect. 3.
By combining lc, hc, and tc we can define additional characteristic scales (Theodoratos et al., 2018). For example, given that drainage areas A have dimensions of L2, we can define a characteristic area Ac as the square of the characteristic length:
Likewise, we can define a characteristic gradient Gc as
and a characteristic curvature κc as
2.3 Incision-threshold number Nθ
In Theodoratos and Kirchner (2020), we derived a dimensionless number, whose definition and interpretation we summarize here. Dimensional analysis of the governing equation with an incision threshold θ (Eq. 2) yielded the dimensionless grouping of parameters . Specifically, all terms in Eq. (2) give rates of elevation change and have dimensions of H T−1. Therefore, to non-dimensionalize Eq. (2), we divided all of its terms by the uplift rate U. The quantity Kθ, which is included in the incision term in Eq. (2) and gives the reduction in the rate of incision due to the threshold, also has dimensions of H T−1. Therefore, dividing the incision term in Eq. (2) by U yielded the dimensionless ratio . We defined this dimensionless ratio as an incision-threshold number Nθ:
This analysis led to a dimensionless version of Eq. (2) that includes only one parameter: the incision-threshold number Nθ. This implies that Nθ is a control on the topography of landscapes that follow Eq. (2). Specifically, model landscapes that have equal incision-threshold numbers Nθ can be set up such that they follow geometrically similar evolutions. Model landscapes that have different Nθ cannot evolve geometrically similarly, and their topographies differ in ways that depend on their Nθ values. Simulation results illustrating these points are presented in Theodoratos and Kirchner (2020).
We proposed two interpretations of the incision-threshold number Nθ in Theodoratos and Kirchner (2020). First, Nθ is defined as the incision rate reduction Kθ relative to the uplift rate U (Eq. 9). The uplift rate U can be viewed as a characteristic rate of elevation change because it is equal to the ratio of the characteristic height to the characteristic time, i.e., (Eqs. 4 and 5). Consequently, Nθ is a normalized incision rate reduction with respect to U. Second, if we rearrange Eq. (9) as , then we can interpret Nθ as giving the magnitude of θ relative to the parameter ratio . Thus, the definition of Nθ shows that incision thresholds from different landscapes should not be compared to each other according to their own values but instead according to their values relative to the ratio of each landscape.
3.1 Defining a steady-state topographic relationship between the steepness index and curvature
In Theodoratos et al. (2018), we presented a relationship that describes the steady-state topography of landscapes that evolve according to Eq. (1). Specifically, if we set and we solve the governing equation for curvature ∇2z, we obtain
The quantity is equal to the steepness index (defined as for drainage area and slope exponents m and n; e.g., Whipple, 2001). For this reason, we refer to Eq. (10) as the curvature–steepness index relationship.
In a coordinate system in which the steepness index ( and curvature (∇2z) are plotted on the horizontal and vertical axes, respectively, Eq. (10) plots as a straight line (for example, see Fig. 1, which we describe in more detail further below). Equation (10) is a testable, quantitative prediction; if a landscape is in steady state and has evolved according to Eq. (1), then curvature should plot as a straight line against the steepness index. Furthermore, this line can give estimates of the parameters K, D, and U because its slope is , and its intercepts are and . While we have not validated this prediction with data, Eq. (10) is a rearranged version of Eq. (5) in Perron et al. (2009), which has been successfully tested with real-world landscape data and has been used to estimate model parameters. Testing Eqs. (10), (11), and (12) as well as Figs. 1 and 2, which are described further below, would be a reasonable next step after the current study.
3.2 Characteristic scales and the curvature–steepness index relationship
If we substitute the characteristic scales lc and κc for the parameter ratios and , then the curvature–steepness index relationship (Eq. 10) becomes
As this equation shows, an interpretation of lc and hc is that they control steady-state topography. Specifically, for a landscape to be in steady state, drainage area A, topographic slope |∇z|, and curvature ∇2z must obey Eq. (11), which is parameterized by the characteristic scales lc and κc or, equivalently, by lc and hc because (Eq. 8). We can graphically illustrate the control of lc and hc on the topography by plotting the curvature–steepness index line described by Eq. (11). As Fig. 1 shows, the properties of such a line are controlled by lc and hc; specifically, its slope is , and its intercepts are and . Note that the slope of this line can be represented either as units of curvature per 1 unit of steepness index or 1 unit of curvature per units of steepness index. For simplicity, we use the latter notation to express the slopes of the curvature–steepness index lines in Figs. 2–4.
Likewise, the curvature–steepness index relationship that corresponds to the LEM with an incision threshold θ is controlled by the characteristic scales lc and κc. This relationship, however, is also controlled by the incision-threshold number Nθ. To derive this relationship, we set in Eq. (2), and we solve it for ∇2z. When we do this for the second subdomain (where , we encounter the ratio . This ratio can be rewritten as (Eqs. 8 and 9). Thus, we obtain the curvature–steepness index relationship:
We plot this equation in Fig. 2 in black, and for comparison we also plot the curvature–steepness index line without an incision threshold (Eq. 11) in gray. The black line consists of two segments that correspond to the two subdomains in Eqs. (2) and (12). The first segment is horizontal and describes a uniform steady-state curvature value of for points with , where incision is fully suppressed by the threshold and only diffusion and uplift operate. The second segment is inclined and corresponds to points with where all three processes operate.
Equations (11) and (12) and Figs. 1 and 2 show that the characteristic scales lc, hc, and κc describe the steady-state topography at points of special interest (see also Theodoratos et al., 2018). Furthermore, some effects of incision thresholds on landscape properties can be visualized by comparing the curvature–steepness index lines with and without an incision threshold (black and gray lines in Fig. 2).
First, the vertical-axis intercept of the curvature–steepness index line without an incision threshold (Fig. 1, Eq. 11) corresponds to ridges and drainage divides, which have A=0 and/or |∇z|=0, i.e., . This intercept shows that the steady-state curvature of ridges and drainage divides is (see also Roering et al., 2007; Perron et al., 2009). Note that −κc is the most negative value of curvature. The horizontal segment of the black line in Fig. 2 (described by the first subdomain in Eq. 12) expresses the fact that, in landscapes with an incision threshold θ, the points with have the same steady-state curvature as ridges and drainage divides, i.e., the most negative value of curvature. This shows that adding an incision threshold to the LEM results in more convex hillslopes (e.g., Howard, 1994; Theodoratos and Kirchner, 2020).
Second, the curvature–steepness index line without an incision threshold (Fig. 1, Eq. 11) has a horizontal-axis intercept of . This intercept corresponds to points with curvature ∇2z=0, which can be viewed as defining the transition between hillslopes and valleys (e.g., Howard, 1994). Thus, points with a steepness index equal to the characteristic height hc can be used to map hillslope–valley transitions (Theodoratos et al., 2018). Adding an incision threshold θ to the LEM makes landscapes steeper and decreases the drainage density; i.e., it makes first-order basins bigger (e.g., Montgomery and Dietrich, 1992; Howard, 1994; Perron et al., 2008). These two effects lead to steeper gradients |∇z| and larger drainage areas A at hillslope–valley transitions. Specifically, as Fig. 2 shows, the horizontal-axis intercept increases from (gray line) to (black line).
3.3 Quantifying and visualizing the effect of the incision threshold on the scales of landscape dissection
In Theodoratos et al. (2018), we derived an interpretation of the characteristic length lc by analyzing the competition between the advection and diffusion of elevation perturbations (e.g., knickpoints), which gives rise to ridges and valleys and controls their characteristic sizes (e.g., Smith and Bretherton, 1972; Howard, 1994; Perron et al., 2008). Following Perron et al. (2008, 2009, 2012), we quantified the relative strength of advection vs. diffusion using a Péclet number (Pe). The definition of our Péclet number differs somewhat from that of Perron et al. Specifically, our definition includes a length scale l that we termed flow path length and that we defined as the distance along flow paths from a given point to the farthest ridge.
The Péclet number is defined (e.g., Perron et al., 2008) as the ratio of a diffusion timescale tD to an incision timescale tI, each of which gives a measure of the strength of the respective process. Specifically, a diffusion timescale can be defined as (e.g., Perron et al., 2008)
This timescale characterizes diffusive propagation over a distance l. In Theodoratos et al. (2018), to define tI, we first calculated the celerity c that corresponds to the incision term in Eq. (1), which is a kinematic wave term (e.g., Whipple and Tucker, 1999). This celerity is equal to . Perturbations can be assumed to be advected at this celerity (e.g., Berlin and Anderson, 2007; Perron et al., 2008). Lague (2014) has criticized this assumption because it does not take into account the effects of knickpoints on hydraulics (e.g., on stream width) and their feedbacks on the rate of knickpoint propagation, especially in the presence of incision thresholds. While we acknowledge this limitation, we nonetheless assume that the rate of knickpoint advection is equal to the celerity c in Eq. (17) because our current focus is on interpreting the characteristic scales lc, hc, and tc that pertain to Eqs. (1) and (2), which do not describe hydraulics explicitly. Therefore, in Theodoratos et al. (2018), we defined the incision timescale tI as the ratio of the flow path length l, which characterizes the location of points within drainage basins, to the celerity c, which characterizes the strength of advection:
Note that small values of tI and tD correspond to strong advection and diffusion, respectively.
We can quantify the relative strengths of advection and diffusion using the ratio of the respective timescales, which defines the Péclet number (Theodoratos et al., 2018):
Diffusive propagation is stronger at points with a Péclet number smaller than 1, and advective propagation is stronger where the Péclet number is larger than 1. Where the Péclet number is roughly equal to 1, diffusion and advection will be roughly equal (when measured by tD and tI). Equation (15) shows that if a point's flow path length l is roughly equal to the characteristic length lc and its drainage area A is roughly equal to the characteristic area Ac, then its Péclet number will be roughly equal to 1, i.e.,
Note that if the incision term has a slope exponent n≠1, then the condition |∇z|≈Gc must be included along with l≈lc and for the Péclet number to be Pe≈1.
The conditions and l≈lc (Eq. 16) are not the only combination of A and l that give Pe≈1, but they are significant because they lead to an interpretation of lc. Specifically, these conditions show that advective propagation, which promotes valley dissection, is dominant at points farther than lc from the ridge and with a drainage area greater than . Therefore, in Theodoratos et al. (2018) we interpreted the characteristic length lc as giving a measure of the smallest scales of dissection. This interpretation does not imply that valley heads are exactly 1 lc away from ridges or that they have drainage areas exactly equal to 1 . Rather, it implies that flow path lengths and drainage areas of valley heads are of a similar order of magnitude as lc and , respectively. Furthermore, it implies that valley heads in different landscapes have l and A that scale with lc and , respectively.
Adding the threshold θ to the incision term of the LEM changes the kinematic wave celerity c and thus the incision timescale tI and the Péclet number (Pe). Specifically, the celerity becomes
and thus the incision timescale tI becomes
Note that the diffusion timescale tD is not affected by the incision threshold. Thus, we can use Eqs. (13) and (18) to define a Péclet number (Pe) for the LEM with an incision threshold θ (Eq. 2), specifically
It can be shown that Eq. (19) can be rewritten as
where Nθ is the incision-threshold number (Eq. 9). Equation (20) shows that adding an incision threshold θ to the LEM reduces the Péclet number relative to the Péclet number for the LEM without a threshold (Eq. 15). This agrees with the fact that the threshold weakens the incision term. More specifically, the Péclet number for the LEM with θ is reduced by the quantity .
Note that the Péclet number definition by Perron et al. (2008) also includes a reduction that depends on Nθ (denoted as θ′ in Perron et al., 2008). The two definitions differ in that ours includes the product (where A is the drainage area and l is the flow path length), whereas the Perron et al. definition includes only a length scale (squared). By including , our definition can account for the scaling of A with l, which depends on the convergence or divergence of topography. The implications of this property of our Péclet number are discussed in Sect. 4.2.3 of Theodoratos et al. (2018).
Using Eq. (20) we see that the conditions l≈lc, , and |∇z|≈Gc, which lead to a Péclet number roughly equal to 1 for the case without an incision threshold (Eq. 16), will lead to when θ is included. The fact that the difference between Pe and 1 is equal to Nθ suggests that we could obtain the value Pe≈1 by adjusting the values of l, A, and |∇z| such that they depend on Nθ. Indeed, we observe that
Note that is larger than lc, which agrees with observations that incision thresholds reduce landscape dissection (e.g., Montgomery and Dietrich, 1992; Howard, 1994; Perron et al., 2008).
Equation (21) shows that, in the case of a landscape that includes an incision threshold θ, the smallest scales of dissection are not characterized by the characteristic length lc on its own, but rather jointly by lc and the incision-threshold number Nθ through the quantity . Consequently, the presence of θ changes the dependence of the scales of dissection on the LEM parameters. Without an incision threshold, the scales of landscape dissection depend on lc, which depends on the incision and diffusion coefficients K and D (Eq. 3). On the other hand, when θ is included in the LEM, the scales of dissection depend on , which depends on K and D, but also on the uplift rate U and the incision threshold θ. We illustrate an example of the dependence on U in Fig. 4b.
The length scales lc and can be expressed graphically by the horizontal- and vertical-axis intercepts of curvature–steepness index lines, specifically by the ratio of these intercepts (or, more precisely, by the ratio of their absolute values). This ratio is equal to in the case without an incision threshold (see Fig. 1) and equal to in the case that includes the incision threshold θ (see Fig. 2). Note that the first ratio is equal to the inverse of the slope of the curvature–steepness index line, which is . On the other hand, the second ratio is not the inverse of this slope, which remains when the threshold θ is included. Instead, it is the inverse of the slope of an auxiliary line connecting the two intercepts. In Fig. 2, we show this auxiliary line as a black dashed line. The effect of the incision threshold on valley dissection can be visualized graphically by comparing the slope of the curvature–steepness index line against the slope of the black dashed auxiliary line. We show this comparison as a thick white arrow.
It should be noted that the characteristic length lc depends only on K and D when the slope exponent is n=1. However, for other values of n, lc will also depend on the uplift rate U; in this more general case, (see Appendix A in Theodoratos et al., 2018). Therefore, the degree of landscape dissection in general is not independent of U. Specifically, an increase in U leads to a decrease in landscape dissection for n>1 and to an increase in landscape dissection for n<1, which agrees with previous observations of the dependence of drainage density on the uplift rate (e.g., Clubb et al., 2016). Interestingly, as revealed by the current study, if an incision threshold is included, the degree of landscape dissection depends on U even for n=1.
3.4 How the curvature–steepness index line responds to parameter value changes
As we show in Figs. 3 and 4, differences in the properties of landscapes with different parameters K, D, U, and θ can be graphically summarized by curvature–steepness index lines because the slopes and intercepts of these lines depend on the characteristic scales lc, hc, and κc, as well as on the incision-threshold number Nθ, which in turn depend on the parameters.
Figure 3 shows curvature–steepness index lines without incision thresholds. It consists of three panels, each showing how the lines respond to an increase in one of the three parameters U, K, and D. In Fig. 3a, an increase in the uplift rate U shifts the curvature–steepness index line downward and to the right without changing its slope. This illustrates that the characteristic height and curvature (hc and κc), which control the intercepts of the line, are proportional to U (Eqs. 4 and 8), while the characteristic length lc, which controls the line's slope, is independent of U (Eq. 3). The parallel shift of the line corresponds to more convex ridges (so that diffusion can keep up with uplift), steeper gradients (so that incision can keep up with uplift), and unchanged landscape dissection. Analogously, Fig. 3b shows that an increase in the incision coefficient K leads to a counterclockwise rotation of the line around the vertical-axis intercept, which corresponds to a more dissected landscape (smaller lc), milder gradients (smaller hc), and unchanged ridge convexity (unchanged κc). Finally, in Fig. 3c, an increase in the diffusion coefficient D results in a clockwise rotation of the line around the horizontal-axis intercept. This corresponds to a smoother landscape with less dissection (larger lc), less convex ridges (smaller κc), and an unchanged steepness index at hillslope–valley transitions (unchanged hc).
Figure 4 illustrates in four panels how curvature–steepness index lines respond to increases in the value of either the incision threshold θ or one of the parameters U, K, and D. As a reminder, a curvature–steepness index line with an incision threshold consists of two segments, a horizontal and an inclined. Note that, as we explain in the previous subsection (Sect. 3.3), a curvature–steepness index line that includes an incision threshold does not express landscape dissection through the slope of its inclined segment, which depends only on the characteristic length lc, but rather through the ratio of the horizontal- and vertical-axis intercepts, which is equal to . This ratio can be graphically illustrated by the slope of an auxiliary line that connects the two intercepts, such as the dashed black line in Fig. 2. In each panel of Fig. 4, we show two dashed black auxiliary lines to illustrate how the ratio of intercepts responds to the parameter changes.
In Fig. 4a we illustrate an increase in θ. The steepness index must reach a greater value before exceeding the increased θ, and thus the horizontal segment of the curvature–steepness index line becomes longer. The vertical position of this segment (along with the vertical-axis intercept) does not change because the characteristic curvature κc does not depend on θ. The slope of the curvature–steepness index line also does not change because lc does not depend on θ. Thus, the increase in θ parallel-shifts the inclined segment of the line to the right. Consequently, the horizontal-axis intercept increases, which expresses the steepening of gradients and the decrease in landscape dissection. The decrease in dissection is also expressed by the fact that the ratio of the horizontal- to the vertical-axis intercept increases, as shown by the clockwise rotation of the dashed auxiliary line.
In Fig. 4b we show that an increase in the uplift rate U parallel-shifts the curvature–steepness index line downward and to the right. Furthermore, the horizontal- and vertical-axis intercepts move to the right and downward, respectively (κc and hc are proportional to U), and the slope of the inclined segment remains unchanged (lc does not depend on U). As we explain in Sect. 3.3, the value of U affects the value of , which expresses the scales of landscape dissection. Specifically, the increase in U leads to a decrease in Nθ. This reflects the fact that θ becomes less important relative to the increased U. Thus, the decrease in dissection due to the threshold is somewhat moderated by the increase in U. This moderation is graphically illustrated by the slopes of auxiliary lines connecting the intercepts of the curvature–steepness index lines. These auxiliary lines do not intersect, and thus their slopes cannot be readily compared visually. Therefore, we plot them again in an inset such that they start from the same point. In this way, we can see that the increase in U leads to a counterclockwise rotation of the auxiliary lines, which expresses the increase in dissection.
In Fig. 4c we illustrate the response of the curvature–steepness index line to an increase in the incision coefficient K. The horizontal segment of the line remains unchanged and the inclined segment is rotated counterclockwise around the point of transition between the two segments. Likewise, the dashed auxiliary line connecting the horizontal- and vertical-axis intercepts is rotated counterclockwise. These responses express the fact that dissection is decreased and that gradients become milder when K is increased. Finally, in panel (d), we show that increasing the diffusion coefficient D leads to a clockwise rotation of the inclined segment of the line around its horizontal-axis intercept, which remains unchanged. The rotation results in moving the horizontal segment up and in rotating the dashed auxiliary line clockwise. These changes express the reduction in landscape dissection and the reduction in the convexity of ridges and hillslopes.
Thus far, we have examined how the influence of the incision threshold θ varies between different landscapes with different parameters using the incision-threshold number Nθ (Eq. 9). This number is constant for any given landscape if the parameters of the landscape are constant. Now, we turn our attention to how the influence of the incision threshold θ varies within a given landscape.
We can quantify the relative influence of the threshold θ on the rate of incision using the fraction . This fraction is equal to Kθ, which is the reduction in the incision rate due to the threshold, divided by , which is the incision rate if there were no threshold. Therefore, shows by what fraction the incision rate is reduced due to the threshold. Where , the fraction is equal to 1, which agrees with the incision rate being reduced by 100 % (i.e., being reduced to zero). At points with smaller than θ, calculating the fraction would not be meaningful; instead, because the threshold completely suppresses incision under these conditions, we assign a value of 1 to the fractional reduction in incision rate.
We can associate the fractional reduction in incision rate with the Tucker (2004) threshold factor Φ. Tucker (2004) defined Φ to quantify the fraction of precipitation events leading to shear stress above a threshold value, i.e., the fraction of events that lead to erosion. Tucker (2004) used Φ to express the incision term of his LEM as . In the case of the LEM examined here (Eq. 2), following Tucker's notation, we can express the incision term as , where the threshold factor Φ is equal to for and to 0 for . Thus, the quantity 1−Φ is equal to the fractional reduction in incision rate, i.e.,
Consequently, in what follows we denote the fractional reduction in incision rate as 1−Φ. We illustrate the properties of the quantity 1−Φ with plots and maps in Figs. 5–7.
In Fig. 5 we plot 1−Φ vs. the steepness index according to Eq. (22). The curve consists of two parts. The first is a horizontal segment that describes the value and corresponds to points with , where incision is fully suppressed by the threshold. The second part corresponds to points with , forming part of a hyperbola that asymptotically approaches 0. This asymptotic approach expresses the fact that, at points with a steepness index much larger than θ, the incision threshold has a very small relative influence on the incision rate.
To indicate how different parts of the 1−Φ curve in Fig. 5 correspond to different regimes of a landscape, we identify the point that corresponds to hillslope–valley transitions. As explained in Sect. 3.2, hillslope–valley transitions can be defined as points with zero curvature and therefore with a steady-state steepness index of . Consequently, the fractional reduction in incision rates at these points is . We can rewrite this value in terms of the incision-threshold number Nθ as . Thus, in Fig. 5, hillslope–valley transitions correspond to the point with coordinates , which we mark with a black dot. The part of the curve above and to the left of this dot corresponds to hillslopes, and the part below and to the right corresponds to the valley network.
With Fig. 6, we examine how the value of the incision-threshold number Nθ of a landscape controls the relationship between the quantity 1−Φ and the steepness index . Specifically, in Fig. 6 we show curves of 1−Φ vs. for four landscapes with incision-threshold numbers Nθ equal to 0.2, 0.4, 1, and 2. The landscapes are assumed to have equal parameters K, D, and U and therefore to have equal characteristic scales. The curves with greater values of Nθ also have greater incision thresholds θ, and thus they have longer horizontal segments. Furthermore, the curves with greater values of Nθ go towards zero more slowly. On each curve, we show the hillslope–valley transition using a black dot. The value of the quantity 1−Φ corresponding to each dot becomes larger as Nθ increases. Thus, in landscapes with smaller Nθ, the incision rate is reduced by large fractions only on the hillslopes, and in valleys it is reduced by small fractions. By contrast, in landscapes with greater Nθ, incision can be reduced by large fractions both on hillslopes and in valleys.
Figure 7 shows maps of the quantity 1−Φ across four steady-state landscapes. We simulated these landscapes with the CHILD (Channel–Hillslope Integrated Landscape Development) numerical model (Tucker et al., 2001). Details about these simulations and additional results are presented in Theodoratos and Kirchner (2020). Here, we provide brief information about the parameters and setup of these simulations in Appendix B. To illustrate how the spatial distribution of 1−Φ depends on the incision-threshold number Nθ, we ran four simulations with Nθ values of 0.2, 0.4, 1, and 2, i.e., the same Nθ values as in Fig. 6. The pixels of the four maps are colored according to their values of 1−Φ using a grayscale that ranges from white to black. Lighter colors correspond to larger values of 1−Φ, i.e., to a stronger influence of the incision threshold. As expected, lighter colors appear near ridges and on hillslopes where the incision threshold has a stronger influence.
The patterns in Fig. 7 reflect the spatial distribution of drainage area and slope because the incision threshold in Eq. (2) is defined as a topographic threshold. However, maps of the quantity 1−Φ would be useful for other formulations of the incision threshold as well. For example, the Tucker (2004) formulation of the incision threshold assumed stochastic precipitation. Tucker quantified the influence of this incision threshold using the threshold factor Φ, which ranges between 0 and 1 (and on which our quantity 1−Φ is based, as mentioned above). Therefore, the quantity 1−Φ could be calculated for the case of the Tucker (2004) LEM, and maps of this quantity would visualize how the influence of the incision threshold is spatially distributed across landscapes.
The fractional reduction in incision rate 1−Φ and the threshold factor Φ can be used to simplify the definition of the Péclet number (Pe). Specifically, we can rearrange the definition of Pe (Eq. 19) such that it includes the fraction :
If this equation is combined with the definition of 1−Φ (Eq. 22), then we can rewrite the definition of Pe in compact form as
where Peθ=0 is the Péclet number for the LEM without an incision threshold (see Eq. 15). Equations (23) and (24) reveal that the influence of the incision threshold on the Péclet number varies across the landscape. Specifically, larger values of Pe, which correspond to larger values of the steepness index , are less sensitive to the incision threshold.
We present graphical methods that summarize topographic and scaling properties of landscapes following a simple stream-power incision and linear diffusion LEM (Eq. 1) and that illustrate the effects of adding an incision threshold θ (Eq. 2). Our results referring to the LEM without an incision threshold (Eq. 1) have been presented before (Theodoratos et al., 2018), but we show them here again to contrast them against those referring to the LEM with the threshold θ (Eq. 2). The two LEMs (Eqs. 1 and 2) assume that the incision term has drainage area and slope exponents m=0.5 and n=1 because this combination significantly simplifies the mathematical derivations. However, as we show in Appendix A, our results are also valid for generic exponents m and n.
For the first graphical method, we plot steady-state relationships between curvature ∇2z and the steepness index (Eqs. 10, 11, and 12), which we obtain from the governing equations (Eqs. 1 and 2). These relationships can be viewed as counterparts of the relationship between topographic slope and drainage area, which is typically assumed to follow a power law in detachment-limited landscapes but not in landscapes that are also influenced by hillslope diffusion. These relationships plot as straight lines (Figs. 1 and 2), whose properties (slope and intercepts) depend on the incision threshold θ and on the characteristic scales of the landscape, which in turn depend on the parameters of the LEM, i.e., on the incision coefficient K, the diffusion coefficient D, and the uplift rate U. (Eqs. 3, 4, and 8). A reasonable follow-up study would be to validate these results against real-world landscapes and specifically explore whether curvature and steepness index data from real landscapes would plot against each other as straight lines.
With Fig. 2, we show that curvature–steepness index lines can graphically illustrate the effects of incision thresholds on landscapes. Specifically, the ways in which curvature–steepness index lines with and without a threshold differ from each other illustrate that thresholds make hillslopes more convex, make gradients steeper, and reduce the drainage density. These effects have been presented elsewhere (e.g., Montgomery and Dietrich, 1992; Howard, 1994; Perron et al., 2008), but the curvature–steepness index lines offer new ways to visualize them graphically. In Figs. 3 and 4, we illustrate the dependence of these properties on the parameters K, D, U, and θ by showing how curvature–steepness index lines respond to increases in these parameters, one at a time. These figures demonstrate an advantage of curvature–steepness index lines: the topographic effects of model parameter changes are expressed as simple shifts and rotations of these lines.
In Sect. 3.3, we examine in more detail the effects of the incision threshold θ on drainage density and the scales of landscape dissection, as well as how these effects can be visualized with curvature–steepness index lines. We assume that dissection is controlled by the competition between the advection and diffusion of perturbations (e.g., Smith and Bretherton, 1972; Howard, 1994; Perron et al., 2008), and thus we examine the effects of θ using a Péclet number (Pe) (Eqs. 15 and 19; see also Perron et al., 2008, 2012; Theodoratos et al., 2018). For the LEM that does not include an incision threshold, we found in Theodoratos et al. (2018) that the characteristic length lc characterizes the smallest scales of dissection. Note that the slope of curvature–steepness index lines is ; therefore, this slope graphically expresses the scales of the dissection of landscapes without incision thresholds. Adding the incision threshold θ, we find that the smallest scales of dissection are characterized by the length scale , where Nθ is a dimensionless incision-threshold number defined as (Eq. 9). This length scale is longer than lc, which expresses the fact that incision thresholds reduce the drainage density. The square of this length scale is and is equal to the ratio between the horizontal- and vertical-axis intercepts of the curvature–steepness index line. As we show in Fig. 2, an auxiliary line connecting these two intercepts would have a slope of . Thus, we can graphically visualize the effect of the incision threshold on landscape dissection by comparing the slope of this auxiliary line with the slope of the curvature–steepness index line.
The second graphical method consists of plots of the dimensionless fraction , which gives the fraction by which the incision rate is reduced due to the threshold (see the governing equation, Eq. 2). We found that this fraction is equal to 1−Φ (Eq. 22), where Φ is a threshold factor (see Tucker, 2004) that subsumes the effect of the incision threshold θ on the incision term of the LEM (see Eqs. 2 and 22). Thus, we denote the fractional reduction in the incision rate as 1−Φ. In Figs. 5–7, we present plots and maps of 1−Φ that illustrate how the relative influence of incision thresholds will vary across a given landscape and how the variation of this relative influence depends on the landscape's incision-threshold number Nθ.
The two dimensionless numbers examined here, Nθ and 1−Φ, quantify the relative influence of the incision threshold θ: the first with respect to the parameters K and U and the second with respect to the steepness index. Thus, Nθ quantifies how θ affects different landscapes with different parameters, and 1−Φ quantifies how the influence of θ varies across different points of a given landscape. We find that the definition of the Péclet number (Pe) can be rewritten in two equivalent forms (Eqs. 20 and 24), which reveal how Pe depends on Nθ and on Φ, respectively.
The three dimensionless numbers, Pe, Nθ, and Φ, along with the characteristic scales lc, hc, and tc, provide a thorough characterization of landscapes that follow the governing equation (Eq. 2). Furthermore, plots of the curvature–steepness index relationship offer a straightforward way to graphically express the geomorphologic meaning of these dimensionless numbers and characteristic scales. Even though the specific definitions of these quantities refer only to the LEMs examined here (Eqs. 1 and 2), the approach that underpins our graphical methods is more generally applicable. For example, an LEM with an incision threshold and stochastic precipitation would have a different governing equation than Eq. (2) and thus a different curvature–steepness index relationship than Eq. (12) and Fig. 2. However, curvature and the steepness index would still be reasonable axes for plotting data from such an LEM if it included diffusion. Likewise, the quantity 1−Φ would follow a different formula than Eq. (22), but maps of this quantity would be useful in visualizing spatial patterns of the influence of the incision threshold across a landscape. Consequently, our graphical methods could potentially be helpful for the analysis of a broader range of models than those examined here.
In this Appendix, we demonstrate that our graphical method remains valid for the case of LEMs with incision terms that have generic drainage area and slope exponents m and n, respectively.
For generic exponents m and n, the governing equations Eqs. (1) and (2) respectively become
Given that the steepness index is defined as (e.g., Whipple, 2001), the quantity Am(|∇z|)n in the above equations is equal to the steepness index raised to the power n, i.e., , and the incision threshold θ is a threshold of the quantity .
Setting in Eqs. (A1) and (A2), we can derive the corresponding steady-state relationships between curvature and the steepness index:
where Nθ is the incision-threshold number, still defined as .
When plotted in axes of ∇2z and , Eq. (A3) has the same basic properties as Eq. (11), which is the curvature–steepness index relationship with m=0.5, n=1, and θ=0. Specifically, Eq. (A3) plots as a straight line with a slope equal to , a vertical-axis intercept equal to , and a horizontal-axis intercept equal to . Consequently, for generic exponents m and n, changes in the values of the parameters K, D, and U are still expressed graphically as shifts and rotations of the curvature–steepness index line, as seen in Fig. 3 for the case of m=0.5 and n=1.
Note that the characteristic scales of length and height (lc and hc) are not equal to and for generic exponents m and n. Rather, they are defined by the more complicated formulas
whose derivation can be seen in Appendix A of Theodoratos et al. (2018). However, the parameter ratios and still express the relative strengths of incision vs. diffusion and incision vs. uplift. By contrast, note that the parameter ratio remains equal to the characteristic curvature κc, which expresses the relative strength of diffusion vs. uplift. Consequently, the shifts and rotations of the curvature–steepness index line still express changes in scaling and topographic properties of landscapes, such as changes in the curvature of ridges, the degree of dissection, and gradients at hillslope–valley transitions.
Likewise, when plotted in axes of ∇2z and , Eq. (A4) has the same properties as Eq. (12): the curvature–steepness index relationship with an incision threshold and with m=0.5 and n=1. Specifically, Eq. (A4) plots as a line with two segments: a horizontal segment at for and an inclined segment with slope equal to and horizontal-axis intercept equal to for . This line, too, responds to changes in the parameters θ, K, D, and U, with shifts and rotations equivalent to the shifts and rotations shown in Fig. 4 for the case of m=0.5 and n=1.
Finally, for generic exponents m and n, the fractional reduction in incision rate 1−Φ is defined as
which plots as shown in Figs. 5 and 6 but in axes of .
We prepared the maps in Fig. 7 with results from numerical simulations that we performed using the CHILD model, originally for the work discussed in Theodoratos and Kirchner (2020). In that work, we present much more information about these simulations and their results. Here, we briefly summarize the model setup and parameterization.
All four landscapes in Fig. 7 have an incision coefficient K = 2 × 10−6 a−1, diffusion coefficient D = 0.5 × 10−2 m2 a−1, and uplift rate U = 0.5 × 10−4 m a−1. Each landscape's incision threshold θ depends on the value of its incision-threshold number Nθ according to (see Eq. 9), where = 25 m for all landscapes. Therefore, the landscapes have the incision thresholds seen in Table B1.
We simulated the four landscapes on triangular irregular networks (TINs) with a total extent of 7.5 km × 11.25 km and an average TIN edge length of 20 m, which resulted in around 250 000 TIN points. Each map in Fig. 7 shows a part of the TIN, specifically a rectangular region with a size of 5 km × 4 km centered around the largest drainage basin of the corresponding landscape.
Details about the implementation of the governing equation (Eq. 2) in CHILD (Tucker et al., 2001) can be found in Theodoratos et al. (2018) and in Theodoratos and Kirchner (2020).
The data presented here were synthesized using the CHILD model (Tucker et al., 2001). The input files needed to reproduce them are available from the corresponding author upon request.
NT derived and analyzed the theoretical, numerical, and graphical results, and NT and JWK interpreted them. NT drafted the paper, and NT and JWK edited it.
The authors declare that they have no conflict of interest.
Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
This study was made possible by financial support from ETH Zurich. We thank Eric Deal for helpful discussions.
This paper was edited by Paola Passalacqua and reviewed by Fiona Clubb and Philippe Steer.
Berlin, M. M. and Anderson, R. S.: Modeling of knickpoint retreat on the Roan Plateau, western Colorado, J. Geophys. Res., 112, F03S06, https://doi.org/10.1029/2006JF000553, 2007.
Clubb, F. J., Mudd, S. M., Attal, M., Milodowski, D. T., and Grieve, S. W. D.: The relationship between drainage density, erosion rate, and hilltop curvature: Implications for sediment transport processes, J. Geophys. Res.-Earth, 121, 1724–1745, https://doi.org/10.1002/2015JF003747, 2016.
Deal, E., Braun, J., and Botter, G.: Understanding the role of rainfall and hydrology in determining fluvial erosion efficiency, J. Geophys. Res.-Earth, 123, 744–778, https://doi.org/10.1002/2017JF004393, 2018.
DiBiase, R. A. and Whipple, K. X.: The influence of erosion thresholds and runoff variability on the relationships among topography, climate, and erosion rate, J. Geophys. Res., 116, F04036, https://doi.org/10.1029/2011JF002095, 2011.
Dietrich, W. E., Bellugi, D. G., Sklar, L. S., Stock, J. D., Heimsath, A. M., and Roering, J. J.: Geomorphic Transport Laws for Predicting Landscape form and Dynamics, in: Prediction in Geomorphology, Geophysical Monograph Series, 135, edited by: Wilcock, P. R. and Iverson, R. M., American Geophysical Union, Washington, DC, USA, 103–132, https://doi.org/10.1029/135GM09, 2003.
Howard, A. D.: A detachment-limited model of drainage basin evolution, Water Resour. Res., 30, 2261–2285, 1994.
Lague, D.: The stream power river incision model: evidence, theory and beyond, Earth Surf. Proc. Land., 39, 38–61, 2014.
Lague, D., Hovius, N., and Davy, P.: Discharge, discharge variability, and the bedrock channel profile, J. Geophys. Res., 110, F04006, https://doi.org/10.1029/2004JF000259, 2005.
Montgomery, D. R. and Dietrich, W. E.: Channel Initiation and the Problem of Landscape Scale, Science, 255, 826–830, 1992.
Montgomery, D. R. and Dietrich, W. E.: Landscape Dissection and Drainage Area-Slope Thresholds, in: Process Models and Theoretical Geomorphology, British Geomorphological Research Group Symposia Series (Book 8), edited by: Kirkby, M. J., Wiley, New York, 221–246, 1994.
Montgomery, D. R. and Foufoula-Georgiou, E.: Channel network source representation using digital elevation models, Water Resour. Res., 29, 3925–3934, 1993.
Perron, J. T., Dietrich, W. E., and Kirchner, J. W.: Controls on the spacing of first-order valleys, J. Geophys. Res., 113, F04016, https://doi.org/10.1029/2007JF000977, 2008.
Perron, J. T., Kirchner, J. W., and Dietrich, W. E.: Formation of evenly spaced ridges and valleys, Nature, 460, 502–505, https://doi.org/10.1038/nature08174, 2009.
Perron, J. T., Richardson, P. W., Ferrier, K. L., and Lapôtre, M.: The root of branching river networks, Nature, 492, 100–103, https://doi.org/10.1038/nature11672, 2012.
Roering, J. J., Perron, J. T., and Kirchner, J. W.: Functional relationships between denudation and hillslope form and relief, Earth Planet. Sc. Lett., 264, 245–258, https://doi.org/10.1016/j.epsl.2007.09.035, 2007.
Scherler, D., DiBiase, R. A., Fisher, G. B., and Avouac, J.-P.: Testing monsoonal controls on bedrock river incision in the Himalaya and Eastern Tibet with a stochastic-threshold stream power model, J. Geophys. Res.-Earth, 122, 1389–1429, https://doi.org/10.1002/2016JF004011, 2017.
Smith, T. R. and Bretherton, F. P.: Stability and the conservation of mass in drainage basin evolution, Water Resour. Res., 8, 1506–1529, 1972.
Snyder, N. P., Whipple, K. X., Tucker, G. E., and Merritts, D. J.: Importance of a stochastic distribution of floods and erosion thresholds in the bedrock river incision problem, J. Geophys. Res., 108, 2117, https://doi.org/10.1029/2001JB001655, 2003.
Theodoratos, N. and Kirchner, J. W.: Dimensional analysis of a landscape evolution model with incision threshold, Earth Surf. Dynam., 8, 505–526, https://doi.org/10.5194/esurf-8-505-2020, 2020.
Theodoratos, N., Seybold, H., and Kirchner, J. W.: Scaling and similarity of a stream-power incision and linear diffusion landscape evolution model, Earth Surf. Dynam., 6, 779–808, https://doi.org/10.5194/esurf-6-779-2018, 2018.
Tucker, G. E.: Drainage basin sensitivity to tectonic and climatic forcing: Implications of a stochastic model for the role of entrainment and erosion thresholds, Earth Surf. Proc. Land., 29, 185–205, https://doi.org/10.1002/esp.1020, 2004.
Tucker, G. E., Lancaster, S., Gasparini, N., and Bras, R.: The channel-hillslope integrated landscape development model (CHILD), in: Landscape erosion and evolution modeling, edited by: Harmon, R. S. and Doe III, W. W., Kluwer Academic/Plenum Publishers, New York, USA, 349–388, 2001.
Whipple, K. X.: Fluvial landscape response time: How plausible is steady-state denudation?, Am. J. Sci., 301, 313–325, 2001.
Whipple, K. X.: Bedrock rivers and the geomorphology of active orogens, Annu. Rev. Earth Pl. Sc., 32, 151–185, 2004.
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, https://doi.org/10.1029/1999JB900120, 1999.