Morphodynamics of river bed variation with variable bedload step length

Here we consider the 1-D morphodynamics of an erodible bed subject to bedload transport. Fluvial bed elevation variation is typically modeled by the Exner equation, which, in its classical form, expresses mass conservation in terms of the divergence of the bedload sediment flux. An entrainment form of the Exner equation can be written as an alternative description of the same bedload processes, by introducing the notions of an entrainment rate into bedload and of a particle step length, and assuming a certain probability distribution for the step length. This entrainment form implies some degree of nonlocality, which is absent from the standard flux form, so that these two expressions, which are di fferent ways to look at same conservation principle (i.e., sediment continuity), may no longer become equivalent in cases when channel complexity and flow conditions allow for long particle saltation steps (including, but not limited to the case where particle step length has a heavy tailed distribution) or when the domain of interest is not long compared to the step length (e.g., laboratory scales, or saltation over relatively smooth surfaces). We perform a systematic analysis of the e ffects of the nonlocality in the entrainment form of the Exner equation on transient aggradational /degr dational bed profiles by using the flux form as a benchmark. As expected, the two forms converge to the same results as the step length converges to zero, in which case nonlocality is negligible. As step length increases relative to domain length, the mode of aggradation changes from an upward-concave form to a rotational, and then eventually a downward-concave form. Corresponding behavior is found for the case of degradation. These results may explain anomalously flat, aggradational, long profiles that have been observed in some short laboratory flume experiments.


Introduction
The Exner equation of sediment conservation, when combined with a hydrodynamic model and a sediment transport model, is a central tool to evaluate the bed evolution (e.g., aggradation and degradation) in the field of morphodynamics of Earth's surface.
The Exner equation, in its classical formulation, relates the bed evolution to the divergence of the bedload sediment flux (q), which is assumed to be a local function of the flow and the topography.However, certain sediment dynamics, such as (i) particle diffusion in river bedload (e.g., Nikora et al., 2002;Bradley et al., 2010;Ganti et al., 2010;Martin et al., 2012), (ii) bed sediment transport along bedrock channels (Stark et al., 2009) and (iii) particle displacements on hillslopes (Foufoula-Georgiou et al., 2010) may show nonlocal behavior that is not easily captured by the classical form of the Exner equation (the notation used throughout the manuscript is defined and listed after the conclusions).
The nonlocality of interest here is embedded in the step length r of a bedload particle, i.e., the distance that a particle, once entrained into motion, travels before being deposited.The existence of a finite step length r implies a Published by Copernicus Publications on behalf of the European Geosciences Union.
A. Pelosi and G. Parker: Morphodynamics of river bed variation nonlocal connection between point x (where a particle is deposited) and point x − r (where it was entrained).The degree of nonlocality can be characterized in terms of the probability density (PDF) of step lengths f s (r).This PDF can be hypothesized to be thin-tailed (e.g., exponential) or heavytailed (e.g., power).
In recent years, considerable emphasis has been placed on asymptotic nonlocality associated with heavy-tailed PDFs for step length (e.g., Schumer et al., 2009;Bradley et al., 2010;Ganti et al., 2010).This is motivated by the desire to preserve nonlocality in the limit of long time, thus leading to fractional advective-diffusive equations (fADE) for pebble tracer dispersion corresponding to the now-classical fADE model (e.g., Schumer et al., 2009).Here we consider nonlocality in a more general sense, as outlined below.
Experiments conducted under the simplest possible conditions (including steady, uniform flow, single-sized sediment and the absence of bedforms) yield thin-tailed and, more specifically, exponential distributions for step length PDF (Nakagawa and Tsujimoto, 1980;Hill et al., 2010).Ganti et al. (2010), however, showed (a) the bed to consist of a range of sizes, (b) the PDF of size distribution to obey a gamma distribution and (c) the PDF of for step length of each grain size to be exponential, the resulting PDF for step length would be heavy-tailed.Hassan et al. (2013) analyzed 64 sets of field data on pebble tracer dispersion in mountain rivers (which by nature contain a range of sizes).They found that all but 5 cases either showed thin-tailed PDFs, or could be rescaled as thin-tailed PDFs.Their results, combined with those of Ganti et al. (2010), however, do suggest that the gradual incorporation of the many factors in nature that lead to complexity can also lead to nonlocal behavior mediated by heavy-tailed PDFs.
Here, however, we focus on the case of nonlocality mediated by thin-tailed (exponential) PDFs for step length.Regardless of the thin tail of the PDF, the degree of nonlocality nevertheless increases with increasing mean step length r.This nonlocality may become dominant when r approaches the same order of magnitude as the domain length L d under consideration.We show that patterns of bed aggradation and degradation are strongly dependent on the ratio r/L d , a parameter that may be surprisingly large in some smallscale experiments.Our results may explain anomalously flat, aggradational, long profiles that have been observed in some short laboratory flume experiments, without relying on either of the fractional partial differential equations or heavy-tailed distributions invoked or implied by Voller and Paola (2010).We use our framework to explore the consequences of heavytailed PDFs for step lengths as well.

Theoretical framework
1-D riverbed elevation variation is classically described by the 1-D Exner equation of sediment conservation in flux form (or equivalently in the 2-D case, divergence form): where η (L) denotes the bed elevation, t (T) denotes the time, x (L) denotes the streamwise distance and q (L 2 T −1 ) is the volume bedload transport rate per unit width.(Here, the porosity of the bed sediment is set to 0 and bedload only is considered, both for the sake of simplicity.)There is, however, a completely equivalent entrainment form of sediment conservation (e.g., Tsujimoto, 1978): where E (L T −1 ) denotes the volume rate of entrainment of bed particles into bedload per unit area per unit time and D (L T −1 ) denotes the volume rate of deposition of bedload material onto the bed per unit area per unit time.
The deposition rate can be related to the entrainment rate by means of the probability density of the step length f s (r) (L −1 ); that is, the probability density of the distance that an entrained particle moves before being re-deposited.Assuming that, once entrained, a particle undergoes a step with length r before depositing, and that this step length has the probability density f s (r) (PDF of step length), the volume deposition rate D can be specified as follows in terms of entrainment rate upstream and travel distance (e.g., Parker et al., 2000;Ganti et al. 2010), so that the entrainment form of sediment mass conservation can be written as As has been shown by Tsujimoto (1978), the two Eqs.( 1) and (4), are in principle equivalent in so far as the following equation precisely describes the bedload transport rate: Yet in any given implementation, Eqs. ( 1) and ( 4) are rarely equivalent.More specifically, in most implementations of the flux Eq. ( 1), q is taken to be a local function of the flow (e.g., Earth Surf.Dynam., 2, 243-253, 2014 www.earth-surf-dynam.net/2/243/2014/bed shear stress), whereas in most implementations of the entrainment Eq. ( 4), E is taken to be a local function of the flow (again, e.g., bed shear stress).The presence of the spatial convolution term in the entrainment Eqs. ( 3) and (4) ensures nonlocality in the entrainment form as compared to the flux form.This nonlocality is present regardless of whether the PDF of step length f s (r) is thin-tailed or heavy-tailed, and vanishes only when f s (r) becomes proportional to δ(r), where δ denotes the Dirac function.In the present implementation of Eq. ( 4), then, we take E to be a local function of flow conditions, so that q is nonlocal according to Eq. ( 5).
It should be pointed out that the formulation of Eq. ( 4) involves a purely kinematic description of particle step length, with the trajectory of the particle unmodified by intervening flow conditions.This is in line with the work of Einstein (1950), Nakagawa and Tsujimoto (1980) and Ganti et al. (2010).In a more detailed analysis, particle momentum balance, and in particular relaxation effects involving, e.g., particle inertia (Parker, 1975;Charru, 2006), should be included.In so far as a step length generally consists of many individual particle saltations (Nino et al., 1994), however, the present kinematic formulation may be sufficient for a firstorder analysis.One other way to formulate the problem is in terms of mass and momentum conservation of two sediment phases; i.e., a static bed phase and a moving bedload phase above it, with exchange between the two (e.g., Charru, 2006).An Eulerian version of such a model would, however, preclude the analysis of nonlocality associated with varying step length, which is a purely Lagrangian parameter.A Lagrangian analysis that includes the dynamics of a particle as it saltates its way through one step length could lead to an improved formulation.
Before continuing, it is of value to specifically indicate what we mean by nonlocality.Equations (3) and (4) are nonlocal in so far as the deposition rate is not determined at a point, but is instead determined from a convolution involving the entrainment rate at every point upstream.The problem is thus nonlocal in the sense of Du et al. (2012).The problem becomes nonlocal in the asymptotic sense only when the PDF of step length f s (r) is heavy-tailed, such that moments beyond a specific value fail to exist (e.g., Schumer et al., 2009).Both cases are considered here; we specifically address the problem of asymptotic nonlocality in Sect. 4.
Here we explore the consequences of nonlocality by comparing the local and nonlocal Eqs.(1) and (4) for Exner over a range of conditions.To do this, we assume that the PDF f s (r) has a mean, and consider the dimensionless parameter ε: where r (L) denotes the mean particle step length and L d (L) denotes the length of the domain of interest (e.g., flume length or length of river reach).The flux and entrainment forms become strictly equivalent only under the constraint: Here we demonstrate that this equivalence for ε 1 breaks down with increasing ε.This is because a finite mean step length r in and of itself implies nonlocality, regardless of whether or not the probabilistic distribution of particle step length f s (r) is thin-or heavy-tailed.A further degree of nonlocality can be introduced by adopting a heavy-tailed distribution for f s (r).
The standard thin-tailed form for the particle step length probability density function is the exponential distribution (e.g., Nakagawa and Tsujimoto, 1980;Hill et al., 2010): The heavy-tailed Pareto distribution with a shift, which ensures that the maximum value of the distribution is realized at r = 0, can be considered as an alternative: where α is the shape parameter and r 0 (L) is the scale parameter.The mean value r of the distribution of Eq. ( 9) can be written as

Numerical model
Here we solve the flux and entrainment formulations under parallel conditions, the only exception being the formulation for step length.To simplify the problem and focus on this point, we approximate the flow as obeying the normal (steady, uniform) approximation.Momentum conservation then dictates that bed shear stress τ b (M L −1 T −2 ) can be represented as proportional to the product of depth H (L) and slope S (1): where u * (L T −1 ) is the shear velocity.
The dimensionless Shields number governing particle mobility is defined as where ρ (M L −3 ) is water density, D c (L) is characteristic bed grain size (here taken to be uniform for simplicity) and R denotes the submerged specific gravity of the sediment (∼ 1.65 for quartz).
A. Pelosi and G. Parker: Morphodynamics of river bed variation The flow can be computed by introducing the Manning-Strickler resistance relation: where U (L T −1 ) is the depth-averaged flow velocity, α r is a dimensionless coefficient between 8 and 9 (Chaudhry, 1993), and k c (L) denotes a composite roughness height.In absence of bedforms, k c is equivalent to the roughness height k s (L), which is proportional to grain size D c by means of a dimensionless coefficient with typical values between 2 and 5 (Parker, 2004).Here, α r is set equal to 8.1, as suggested by Parker (1991) for gravel-bed streams, while k c , in absence of bedforms, is taken to be 2.5 times the grain size D c (Parker, 2004).
The equation for water conservation for quasi-steady flow is where Q w (L 3 T −1 ) is the water discharge and B (L) denotes the channel width.
Combining Eqs.(11-14), we relate the dimensionless Shields number to the flow properties: The basis for our morphodynamic calculations is the form of Meyer-Peter and Müller (1948), as modified by Wong and Parker (2006).It takes the form where g (L T −2 ) denotes the gravitational acceleration.The parameter τ c denotes the threshold Shields number and γ is a coefficient of proportionality; these parameters take the respective values 0.0495 and 3.97 (as specified by Wong and Parker, 2006).The volume bedload transport rate per unit width q at equilibrium can also be written as (Einstein, 1950) so that the entrainment rate takes the form Here β is a dimensionless parameter.Einstein (1950), suggested, based on a simple flume-like configuration, that r/D c takes a value on the order of 100-1000, so that a step length is about 100-1000 grain sizes.This order of magnitude has been confirmed by the experiments of Nakagawa and Tsujimoto (1980), Wong et al. (2007) and Hill et al. (2010).
In systems with higher degrees of complexity, however, β is likely to vary over a wide range.Combinations of multiple grain sizes, bedforms, scour and fill and partially exposed bedrock are likely to give rise to connected pathways along which particles may travel for an extended distance, so giving rise to larger values of r (e.g., Parker, 2008).In order to capture this effect in a simplified 1D model, we allow the ratio r, and thus β = r/D c to vary freely, so that the ratio r / L d of step length to domain length can vary from 0 (in which case the flux and entrainment formulations become equivalent) to unity (in which a particle starting at the upstream end of the domain reaches the downstream end in a single step).
Linking Eq. ( 18a), the following relation arises at equilibrium conditions: Our formulation is such that increased step length is adjusted against reduced entrainment, so that the equilibrium bedload transport rate is the same whether the flux or entrainment formulation is used.A difference, however, arises under disequilibrium conditions, in which case Eq. ( 16) is solved in conjunction with Eq. ( 1) in the flux case, and Eq. ( 18a) is solved in conjunction with Eq. ( 4) in the entrainment case.This allows us to capture the difference between the two formulations in a comparable way.The flux formulation, Eq. ( 1), corresponds to a nonlinear diffusion equation, i.e., where according to Eqs. ( 11), ( 15) and ( 16), the kinematic diffusivity ν (L 2 T −1 ) is a function of bed slope S = −∂η / ∂x: The governing equation is second order in x, and thus requires two boundary conditions.Here we require that the bed elevation at the downstream end is zero, and that the sediment transport rate at the upstream end is given as a constant, specified feed rate: The entrainment formulation of Eq. ( 4), however, is only first order in x, in so far as the entrainment rate E is a specified function of bed slope S = −∂η/∂x according to Eqs. ( 4) and (18a).Thus there can be only one boundary condition in x; here we use Eq.(22a) for this, so that both the flux and entrainment formulations satisfy the condition of vanishing bed elevation (corresponding to set base level) at the downstream end.
Although no boundary condition can be set at the upstream end for the entrainment formulation, it is still possible to choose conditions so that the sediment transport rate at the upstream equals the feed value under equilibrium conditions.
To do this, we assume that the entrainment rate everywhere upstream of x = 0 equals a specified value E f , specified as follows: The deposition rate D(x) of Eq. ( 3) can then be rewritten in terms of the sum of particles that originate within the domain (x − r ≥ 0) and those that originate upstream of the domain (x − r < 0): where is the probability (L −1 ), that a particle travels at least a distance x.
The entrainment form of sediment mass conservation thus takes the ultimate form For the numerical computation, we nondimensionalize Eqs. ( 1) and ( 26).We assume that the computation begins from some equilibrium initial condition with spatially constant slope S in , bedload transport rate and entrainment rate q in = r E in .At t = 0, however, the supply of sediment is impulsively altered, causing subsequent bed aggradation or degradation, but with an altered sediment feed rate for t > 0.
We normalize against initial equilibrium conditions using the following definitions: In addition, we nondimensionalize the entrainment rate (for the entrainment formulation) and the bedload transport rate (for the flux formulation) as Then, the nondimensional flux and entrainment forms of the sediment mass conservation, Eqs. ( 1) and ( 26) take the respective forms where is the dimensionless step length PDF for the exponential distribution, and is the corresponding form for the Pareto distribution, where r0 is the dimensionless scale parameter equal to r 0 /L d .These are the upstream conditions, for the entrainment formulation and for the flux formulation The downstream boundary condition is the same for both Here Ê f is an imposed upstream entrainment rate, and ε Ê f is an imposed upstream bedload feed rate, chosen to be different from the initial equilibrium values so  that the bed is forced to aggrade (or degrade) toward a new equilibrium state.
Manipulating the relations of Eqs. ( 15) and (18a), with the definitions of Eqs. ( 27), Ê can be at any given time as where τ * in is the dimensionless Shields number, calculated from Eq. ( 15) with the initial flow and bed conditions and ŝ is the local dimensionless slope.
The key parameter of interest here in describing the difference between the entrainment and flux formulations is ε.In the case ε 1, both formulations become identical.We show below, however, that as ε increases, the response to change in sediment supply differs between the two cases.
We discretize the relation between dimensionless slope and dimensionless bed elevation as follows: The discretization of the domain is schematized in Fig. 1: a central finite-difference scheme is used to solve Eqs. ( 28) and ( 29).

Results
Here we compare the results for aggradation and degradation for the entrainment formulation with varying values of   8).In the case of the flux form and the entrainment form with ε = 0.01, slope increase is first realized upstream and then propagates downstream in time.For the case ε = 0.5, slope more or less increases simultaneously everywhere, corresponding to the rotational evolution in Fig. 2c.In the case ε = 1, slope first increases downstream, the effect then gradually propagating upstream in time.
ε against those for the flux formulation.In Fig. 2, bed elevation profiles are shown, having set as an upstream boundary condition Ê f = 2, so forcing the bed to aggrade.Case (a) is the solution for the flux form of Eq. ( 28), while cases (b), (c)  8).As ε increases, the differences between the results from the two forms increase because of the nonlocality of particle movement, with an evolution from downward-concave transient profiles to upward-concave ones.
and (d) are the solutions for the entrainment form of Eq. ( 29) solved, respectively, for ε = 0.01, 0.5, and 1.As expected, the solutions of Eqs. ( 28) and ( 29) collapse to nearly the same results in the case ε = 0.01; i.e., when the mean particle step length is short compared to the length of the domain.Under this condition the local (flux) form, essentially coincides with the entrainment form.For higher values of ε, however, the differences between the results increase because the entrainment form is able to capture the nonlocal feature of the particle movement.For the flux form and the case ε = 0.01, the aggradational profile is strongly upward concave, with bed slope declining downstream.The transient aggradational bed profiles tend to assume a nearly linear profile, and thus the bed rotates upward for values of ε close to 0.5.For higher values a downward-concave form profile is realized.
To highlight and quantify this change in shape, we introduce a concavity parameter δ, which measures the deviation, in the center of the profile at x = 0.5 relative to the constant initial slope: where η| x=0 denotes the dimensionless bed elevation at x = 0 and η| x=0.5 denotes the same quantity in the center of the profile ( x = 0.5).Positive δ indicates upward concavity, while negative δ indicates downward concavity.In Fig. 3, the variation in time of δ is shown for the flux case, and different values of ε for the entrainment case.It is seen that δ is positive for smaller ε and but becomes negative for ε   8).The observed behavior corresponds to that of Fig. 4. In the case of the flux form and the entrainment form with ε = 0.01, slope decrease is first realized upstream, and then propagates downstream in time.For the case ε = 0.5, slope more or less decreases simultaneously everywhere, corresponding to the rotational evolution in Fig. 5c.In the case ε = 1, slope first decreases downstream, the effect then gradually propagating upstream in time.
greater than 0.5.The results for the flux form overlap with the form for ε = 0.01.In Fig. 4, the slope evolution is plotted: the typical upward concave shape for the flux case and ε = 0.01 is due to the preferential proximal deposition of sediment, which causes the sediment load and thus the Shields number τ* to decrease downstream (Parker, 2004).Thus, according to Eq. ( 15), a downstream decreasing slope is realized (Fig. 4a, b).However, a downward concave shape for ε = 1 is characterized by an increasing slope downstream (Fig. 4d).This corresponds to bedload particles that can jump from the upstream end of the domain to the downstream end in one step.
For completeness, the case of degradation, due to an imposed entrainment and feed rate upstream Ê f = 1/2, is described by Figs. 5, 6 and 7.The results show a congruent behavior with the aggradation case.In Fig. 5, for ε = 0.01 and Ê f = 1/2, it is seen that the two profiles more or less agree.In Fig. 6, the concavity parameters δ also more or less agree for this case.When ε increases to 1, the concavity of the transient degradational profiles changes from downward to upward.In Fig. 7, slope changes from increasing downstream to decreasing upstream.When ε = 0.5, it is shown in Fig. 7 that the transient profiles tend to keep a straight shape, and the evolution of the bed is essentially rotational about the downstream end.
Summarizing, (i) the flux model and the entrainment model yield essentially the same results for ε = 0.01; (ii) for ε = 0.5, nearly rotational aggradation and degradation are ob- tained; and (iii) for ε = 1, the pattern of concavity is reversed compared to the flux case.
Then, a Pareto distribution with a shift, i.e., Eq. ( 9) for particle step length distribution, is considered as well so as to compare the case of heavy tail of the PDF of step length with the thin-tail exponential form.In the calculations for the entrainment rate with Ê f = 2, two cases are evaluated: (a) ε = 0.015 and (b) ε = 1.It is seen that the two profiles more or less agree for case (a).A more substantial difference is seen for case (b), but the concavity is quite small for both the cases of thin-tailed and heavy-tailed PDFs for step length.Assuming L = 200 m, with a thin-tailed PDF the value ε = 0.015 corresponds to a mean step length equal to 3 m, and the value ε = 1 corresponds to 200 m.We have set the shape parameter α in the Pareto PDF equal to 1.5, and the scale parameter r 0 equal to 1.5 m for case (a), and t at 100 m for case (b).This yields values of r from Eq. ( 10) that are respectively equal to 3 and 200 m; i.e., the same values as for the thin-tailed case.
The analysis shows that the shape of the tail of the step length PDF does not significantly change the results for ε = 0.015 but does result in some change compared to the thintailed case ε = 1. Figure 8 shows the long profiles resulting from both the thin-tailed and heavy-tailed cases, and Fig. 9 shows the corresponding evolution of concavity.As seen in Fig. 9c and d corresponding to the case of aggradation with ε = 1, the profiles are downward-concave for the thin-tailed PDF of step length, and upward-concave for the heavy-tailed case.The concavity in both cases, however, is so small that the same rotational behavior for profile adjustment is seen, as documented in Fig. 8c and d.
In interpreting the results regarding the thin-tailed and heavy-tailed cases, it should be recalled that the problem is solved numerically only over the domain 0 < x ≤ 1, with the further constraint 0 < ε = r/L d ≤ 1.This constraint prevents attainment of an asymptotic nonlocal state.An example of an asymptotical form is given below.

Sample asymptotic nonlocal relation for entrainment form of mass conservation
Taking the spatial Fourier transform of Eq. ( 4) results in the form where the Fourier transform of any parameter Z(x) is given as Following the analysis of e.g., Ganti et al. (2010), we assume that the PDF f s (r) has a mean r but no standard deviation, so that fs (k) can be expanded in asymptotic form where 1 < α < 2. The implication of this is that f s (r) has a power-law tail.Substituting Eq. ( 40) into Eq.( 38), inversetransforming back to real space and reducing with Eq. ( 17), it is found that Eq. ( 4) reduces to In so far as q is specified by Eq. ( 16), Eq. ( 41) takes the form of a nonlinear fractional PDE (partial differential equation).
While the asymptotic form is of interest from a theoretical point of view, numerical solutions of specific problems are more easily carried out in terms of the original convolution form of Eq. ( 4).

Discussion and conclusions
The main goal of the work is to show how the entrainment form of the Exner equation of sediment continuity diverges from the flux form of the Exner equation when nonlocal behavior in particle motion arises: (i) as the mean particle step length r increases from 0 to the order of magnitude of the domain length L d for a thin-tailed step length PDF and (ii) as a heavy-tailed PDF for particle step length is used.The dimensionless parameter ε is defined as the ratio between the mean step length r and the length of the domain of interest L d .We analyzed the effect of variation of ε on bed aggradational/degradational profiles by solving the entrainment form of the Exner equation, with the assumption of a thin-tailed PDF for particle step length.As expected, the two forms collapse in the case ε 1.
For high values of ε, however, the differences between the results from the two forms increase because of the nonlocality of particle movement, which is not captured by the classical flux form of the Exner equation: the transient aggradational (degradational) bed profiles tend to assume, for ε greater than 0.5, a downward (upward) concave shape, rather than the upward (downward) concave shape of the flux form.When the value of ε is close to 0.5, an interesting behavior for both cases of aggradation and degradation has been found: the transient profiles tend to rotate around the downstream point, keeping almost a straight shape.For a value of ε in the range [0,0.5), the concavity of the bed profiles is still upward for aggradation and downward, for degradation, but by increasing ε to 0.5, the concavity is nearly vanishing.These results may serve as an explanation for relatively flat aggradational bed profiles, which have been achieved in some short laboratory experiments (e.g., Muto, 2001;and Voller and Paola, 2010), where the value of the ratio between mean particle step length and length of the domain of interest may not be negligible.At the laboratory scale, the mean step length becomes comparable to domain length so that the inclusion of nonlocal effects in the PDF of step length, which this circumstance entails, should clearly be evaluated in order to properly model the bed evolution.
The analysis also investigates the effect of the heavy tailedness in the PDF of step length on the bed profile.For the case studied, we show that the variation of the shape of the step length distribution from thin-to heavy-tailed does not significantly influence the results when step length is small.This is probably due to the "short" domain length compared to the tail of the power law distribution.There is a somewhat larger difference in the case when step length equals domain length, but the bed elevation profiles are nearly linear for both thintailed and heavy-tailed PDFs.Voller and Paola (2010) introduced heavy-tailed behavior to explain profiles that evolve with concavity that is small compared to the standard flux case of Eq. ( 1).Here we find that a heavy-tailed behavior is not necessary to obtain this result.
Recently Falcini et al. (2013) have presented a nonlocal formulation for sediment transport and bed evolution that bears comparison to the present work.They assume a locallydetermined "reference [sediment transport rate] q L whose physical interpretation requires some care", and then integrate this with a nonlocal weighting function to determine the actual sediment transport rate q at a section.Their analysis can yield an upward concave nonlinear final equilibrium state in the absence of subsidence (which is not included in this analysis as well), whereas the present analysis predicts only equilibrium states with constant slope.That is, in our analysis, profile concavity or convexity is a transient phenomenon.We suggest that our analysis has a somewhat clearer basis than that of Falcini et al. (2013), who determine their weighting function from heuristic considerations.
Long step lengths of bedload particles in the field may result from any bed pattern that induces preferential paths for transport, including grain size mixtures (Ganti et al., 2010), bedforms, scour and fill, and intermittent bedrock exposure (Stark et al, 2009).Thus our results may be applicable to these cases.The case of sediment suspension can also be represented in entrainment form (e.g., Parker, 2004).This case is generally associated with much longer mean path lengths than the case of bedload.As a result, the suspension-dominated case may show much more nonlocal behavior than the bedload case.This case deserves further investigation.volume bedload transport rate per unit width (L 2 T −1 ) q normalized q against initial conditions q f sediment transport rate at the upstream end (L 2 T −1 ) Q w water discharge (L 3 T −1 ) r step length r of a bedload particle (L) r mean particle step length (L) r dimensionless particle step length

Figure 1 .
Figure 1.Discretization of the domain.

Figure 2 .
Figure 2. Bed profile evolution for the case Ê f = 2: (a) flux form; (b) entrainment form for ε = r/L d = 0.01, (c) entrainment form for ε = 0.5 and (d) entrainment form for ε = 1, using the thin-tailed exponential step length function of Eq. (8).As ε increases, it is clearly seen that the differences between the results for the two formulations increases.More specifically, as ε increases, nonlocality effects mediate a transition from upward concave transient profiles to downward concave transient profiles.

Figure 3 .
Figure 3. Aggradation case: variation in time of the dimensionless concavity parameter δ in the case of the flux formulation, and in the cases of the entrainment formulation for different values of ε ranging from 0.01 to 1.The result for the flux form overlaps with the result for the entrainment form with ε = 0.01.Note the reversal in behavior as ε increases beyond about 0.5.

Figure 4 .
Figure 4. Slope profile evolution for the case Ê f = 2: (a) flux form; (b) entrainment form for ε = r/L d = 0.01; (c) entrainment form for ε = 0.5; and (d) entrainment form for ε = 1, using the thin-tailed exponential step length function of Eq. (8).In the case of the flux form and the entrainment form with ε = 0.01, slope increase is first realized upstream and then propagates downstream in time.For the case ε = 0.5, slope more or less increases simultaneously everywhere, corresponding to the rotational evolution in Fig.2c.In the case ε = 1, slope first increases downstream, the effect then gradually propagating upstream in time.

Figure 5 .
Figure 5. Bed profile evolution for the case Ê f = 1/2: (a) flux form; (b) entrainment form for ε = r/L d = 0.01; (c) entrainment form for ε = 0.5; and (d) entrainment form for ε = 1, using the thin-tailed exponential step length function of Eq. (8).As ε increases, the differences between the results from the two forms increase because of the nonlocality of particle movement, with an evolution from downward-concave transient profiles to upward-concave ones.

Figure 6 .
Figure 6.Degradation case: variation in time of the dimensionless concavity parameter δ in the case of the flux formulation and in the cases of the entrainment formulation for different values of ε ranging from 0.01 to 1.The result for the flux form overlaps with the result for the entrainment form with ε = 0.01.Note the reversal in transient behavior as ε increases beyond about 0.5.

Figure 7 .
Figure 7. Slope profile evolution for the case Ê f = 1/2: (a) flux form; (b) entrainment form for ε = r/L d = 0.01; (c) entrainment form for ε = 0.5; and (d) entrainment form for ε = 1, using the thintailed exponential step length function of Eq. (8).The observed behavior corresponds to that of Fig.4.In the case of the flux form and the entrainment form with ε = 0.01, slope decrease is first realized upstream, and then propagates downstream in time.For the case ε = 0.5, slope more or less decreases simultaneously everywhere, corresponding to the rotational evolution in Fig.5c.In the case ε = 1, slope first decreases downstream, the effect then gradually propagating upstream in time.

Figure 8 .
Figure 8. Bed profile evolution for the case Ê f = 2. (i) ε = 0.015: (a) thin-tailed exponential step length PDF; (b) heavy-tailed Pareto step length PDF (α = 1.5, r 0 = 1.5 m).(ii) ε = 1 (c) thin-tailed exponential step length PDF; (d) heavy-tailed Pareto step length PDF (α = 1.5, r 0 = 100 m).The shape of the tail of the step length PDF does not significantly change the results for ε = 0.015, but does result in some change compared to the thin-tailed case ε = 1.It should be realized that the numerical calculation has been carried out under the constraints 0 < ε = r/L d ≤ 1 and 0 ≤ x = x/L d ≤ 1, constraints that preclude the evolution of asymptotic behavior.

Figure 9 .
Figure 9. Variation in time of the concavity parameter δ for the case of the thin-tailed exponential distribution for step length, and the case of heavy-tailed Pareto distribution for step length.The parameter ε = r/L d takes the value 0.015 in (a) and 1.0 in (b).
of x = 0 (L T −1 ) Ê f normalized E f against initial conditions f ls (x)probability (L −1 ) that a particle travels at least a distance x f s (r) probability density (PDF) of step lengths (L α shape parameter for the Pareto step length PDF α r dimensionless coefficient in the Manning-Strickler resistance relation β ratio between particle mean step length and grain size γ coefficient of proportionality in the bedload transport formulae δ concavity parameter ε ratio between particle mean step length and domain length η