Modeling sediment transport in river networks

Physica A 387 (2008) 6421–6430

Contents lists available at ScienceDirect

Physica A journal homepage: www.elsevier.com/locate/physa

Modeling sediment transport in river networks Xu-Ming Wang ∗ , Rui Hao, Jie Huo, Jin-Feng Zhang School of Physics and Electric Information, NingXia University, Yinchuan, 750021, PR China

article

info

Article history: Received 31 December 2007 Received in revised form 11 July 2008 Available online 30 July 2008 Keywords: River network Sediment transport Sediment-carrying capability Erosion Sedimentation

a b s t r a c t A dynamical model is proposed to study sediment transport in river networks. A river can be divided into segments by the injection of branch streams of higher rank. The model is based on the fact that in a real river, the sediment-carrying capability of the stream in the ith segment may be modulated by the undergone state, which may be erosion or sedimentation, of the i − 1th and ith segments, and also influenced by that of the ith injecting branch of higher rank. We select a database about the upper-middle reach of the Yellow River in the lower-water season to test the model. The result shows that the data, produced by averaging the erosion or sedimentation over the preceding transient process, are in good agreement with the observed average in a month. With this model, the steady state after transience can be predicted, and it indicates a scaling law that the quantity of erosion or sedimentation exponentially depends on the number of the segments along the reach of the channel. Our investigation suggests that fluctuation of the stream flow due to random rainfall will prevent this steady state from occurring. This is owing to the phenomenon that the varying trend of the quantity of erosion or sedimentation is opposite to that of sediment-carrying capability of the stream. © 2008 Elsevier B.V. All rights reserved.

1. Introduction Natural river networks have attracted a good deal of attention in physics and geophysics communities [1–12]. As early as 1953, empirical observations on natural rivers by Leopold revealed many power-law relationships between some parameters of river networks [1]. A pioneering and systematic investigation on such relationships was made by RodriguesIturbe and Rinaldo [2]. For instance, in their paper [3], the theoretical prediction conducted by a simple finite scaling ansatz indicates the scaling properties, which are in good agreement with those of natural river basins. And then a large number of models have appeared for the purpose of getting a deeper understanding of how these natural events occur in evolutive processes [4–12]. The majority of them is based on local erosion rules to model the course of the initial earth surface gradually evolving into a real landscape. Some of these simulations can depict the landscape very vividly [4,13], show the famous power laws of the slope-area, show Horton’s laws that are relations between the number and length of different parts of the network [13], as well as other power laws of the river parameter distributions in the drainage area [14]. The main spirit of these models embodies the nature of water: flows downhill and changes the local height of the landscape due to erosion. Another explanation for the formation of river networks and their scaling properties are the theoretical studies based on the minimum energy dissipation [15–17], which introduces the assumption that an open system with fixed energy input tends to form a structure that minimizes the total energy dissipation in the process of precipitated water flowing downhill. The realistic looking drainage patterns and the scaling laws given by models suggest that these theories have captured the physical essential, at least, in the aspect of potential energy transfer in the fluvial erosion process. Indeed, erosion is fundamental and intrinsic for the development of the landscape and river patterns [18–20]. It is known that

∗

Corresponding author. Tel.: +86 951 2061182. E-mail address: [email protected] (X.-M. Wang).

0378-4371/$ – see front matter © 2008 Elsevier B.V. All rights reserved. doi:10.1016/j.physa.2008.07.027

6422

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

river patterns bear striking fractal forms and self-organized characteristics [3]. As stated by Dhar [17], this self-organized structure is the result of feedback between water flow and the landscape itself via the variation of the erosion. This implies that sedimentation cannot be neglected, at least, for the purpose of investigating the evolution of the river patterns [18,20]. In fact, sedimentation, together with erosion, is also essential for the transport of sediment in mature river networks, since the channel of a river will alter its pattern if sedimentation occurs. Therefore, a consideration of the aforementioned feedback in which sedimentation participates, is necessary to get an understanding of sediment transport in river networks. This paper will discuss the transport of sediment in river networks by a dynamical model that involves such a feedback mechanism. It is known that the confluence of streams is one of the basic processes of a river network, which may cause a change of runoff in sediment-carrying capability (SCC), and further result in the variation of the quantity of erosion or sedimentation (which can be abbreviated to QES, and defined as the local difference between the sediment carried by the inflow and that carried by the outflow). The variation of the QES inversely causes the SCC to modify. This co-adjustment may play a fundamental role in forming a relatively steadier configuration of the channels. We define here the rank of a stream as the times of the confluence, while it drains to the sea. Please note that the confluence refers to not the lateral injection of the other into this stream, but that of this stream into the other. It is, without exception, a fact that one stream laterally injects into the other when confluence takes place in a river network, so the rank of a given stream can be uniquely determined. If one goes against a stream of rank h and distinguishes between the lateral injecting streams and the others at each of the confluences, all branches of rank h + 1 on this branch can be found. The injecting streams are of rank h + 1 and the others are actually the segments (will be defined later on) of the same stream of rank h. In this way, all of the branches can be ranked. This definition differs from the conventional one, Horton ordering, and may be convenient for us to study the sediment transport along a main channel of a river. Obviously, a main stream is of rank 0. To formulate the self-adjusting ability, the channel of a branch (including the main stream), for instance, of rank i, can be naturally divided into segments by injections of branches of rank i + 1. The first segment of a branch is at the upper end, which is actually taken as the source. And then the segment number increases one by one as branches of rank i + 1 inject into the stream in turn. For convenience, we denote the stream flowing on a channel segment by the segment itself. 2. The model For one thing, we consider the confluence of the segments in a time interval, t → t + δ t, which can be simply denoted by t (one of the discrete time variables, i.e., t = 1, 2, . . .). So the confluence, in the tth time step, can be directly expressed as Qhl (i + 1, t ) = Qhl (i, t ) + Qhi +1 (ilast , t )

(1)

where h = 0, 1, 2, . . ., i, l = 1, 2, . . ., and i = 1 denotes the source of the stream. The character Q denotes the stream flow. The left side of Eq. (1) represents the outflow of the i + 1th segment of lth branch of rank h, the right side is the inflows respectively coming from the ith segment of this branch and the ith branch of rank h + 1 (which is in fact denoted by its last, ilast , segment). Obviously, the outflow balances the inflows. It implies that in a segment, the loss of water via evaporation, irrigation, and infiltration into the earth is neglected, and the increase of water can only be attributed to confluence, in the case of steady flow. Since the SCC of the segments changes due to their different channel configurations and hydraulic factors, the sediment carried by the outflow may not always be equal to that by the inflows. So erosion or sedimentation will occur in these segments. Based on the confluence expressed by Eq. (1) and the definition presented above, the QES can be determined by equation,

∆Shl (i + 1, t ) = Shl (i + 1, t ) − Shl (i, t ) − Shi +1 (ilast , t ).

(2)

On the right hand side of the equation, the first term denotes the sediment carried by the outflow, the second and third denote that carried by the inflows, respectively from the upper segment and the injecting branch of higher rank. ∆Shl (i + 1, t ) > 0 indicates that the i + 1th segment is in a scouring state in the tth time interval, while ∆Shl (i + 1, t ) < 0 implies that it is in a deposition state. A special, but nontrivial, case is ∆Shl (i + 1, t ) = 0, which denotes a balanceable state. The quantity of sediment carried by a stream may depend on the landform, the scale of the basin, the configuration and slope of the channel and the velocity of the runoff, etc. It may be very difficult to explicitly give a formula to express the dependance of the transporting sediment on so many factors. However, we take an empirical or statistical relation between the sediment and the stream flow respectively obtained on the Yellow River (in China) [21], the Arbúcies River (in Spain) [22] and the Mississippi River (in America) [23] to reveal the dependence of the sediment erosion or deposition on the water discharge, which maybe help us to understand the main dynamics in the sediment transport process. That is, S = kQ n

(3)

where n is a constant which partially integrates the above-mentioned factors, and k is a coefficient. Many similar statistical relations have been obtained in rivers elsewhere, for instance, in catchments on the tableland of New South Wales [24] or in laboratory sediment transport experiments [25]. These indicate that the relationship may be a general description, and

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

6423

show the dependence of sediment transport capability mainly on the water discharge. Such an equation can also be derived from a well-known formula that was conventionally called sediment-carrying capacity S = k(u3 /(hω))n

(4)

where u denotes flow velocity, h water depth and ω sediment settling velocity [26]. This formula can be deduced from the principle of energy conservation and some properties of sediment in fluid [27,28], and implies that the suspended sediment carried by fluid decreases as h and/or ω increases. It indicates that the SCC of a stream is self-adjustable by h and ω. Many experiments have revealed that the suspension of sediment is mainly attributed to the eddy motion rather than the buoyancy force of water. And the sediment in the stream can suppress turbulence [28,29]. So the SCC is modulated by the sediment concentration [30]. If we incorporate the dependence of the SCC on h and ω into that of the coefficient k, and replace u3 by the stream flow Q , i.e., u3 → k0 Q (k0 is the transforming coefficient, and can also be incorporated into k), Eq. (4) then becomes Eq. (3). And thus coefficient k in Eq. (3) now actually represents the adjustable SCC for a given Q and the constant n. In the following text, we name it after sediment-carrying capability. Here n is set at an approximation, n ∼ 2, which was obtained empirically in the Yellow River in the lower-water seasons [21]. Similarly, the adjustable coefficient k can be labeled by four characters, i.e., klh (i, t ). According to the above discussions, the feedback between the erosion and sedimentation in the sediment transporting process can be incarnated by the self-adjustment of k. This adjustment can be conducted by the QES on a segment, and that on its neighboring segments. (1) For a natural river, sedimentation implies that both the deposited and the suspended sediment particles are relatively smaller, and the viscosity intensity of fluid is enhanced, so that the friction of river-bed against the flow will increase and the stream can carry more sediment. On the contrary, erosion implies that the deposited material over the river-bed and the suspended particles become coarser, friction will decrease, thus the stream can transport less sediment. Accordingly, for a segment in the model, e.g., the segment i + 1, the SCC (is now represented by k) in t + 1th time interval will increase if it is in a deposition state in the tth time interval (∆Shl (i + 1, t ) < 0), while the SCC will decrease if it is in a scouring state (∆Shl (i + 1, t ) > 0). These adjustments agree with that in Eq. (4) performed by ω. The sediment concentration will decrease as deposition occurs. This decrease will cause, according to the above statement, an enhancement of turbulence, the SCC will increase, while the increase of the sediment concentration as scouring occurs will result in the attenuation of turbulence, the settling motion of the sediment will increase, so the SCC will decrease. (2) Whereas, the QES on segment i will influence the SCC of segment i + 1, in comparison with what has been presented about segment i + 1, in the reversed way. If ∆Shl (i, t ) > 0, it is inevitable that segment i + 1 has to discharge more sediment that comes from segment i, while if ∆Shl (i, t ) < 0, the sediment coming from segment i decreases, and then segment i + 1 will have its SCC decreased. Naturally, the SCC of segment i + 1 can also be affected by the state of the ith branch of the higher rank (which will join together with the segment i in the t + 1th time interval) in the same way. The effect of the state of segment i + 2 on the SCC of segment i + 1 may be similar to that of segment i. For simplifying our model, we neglect this effect. In fact, the SCC of a segment should be determined by the integration of aforementioned factors. This integration may imply competition between the two opposite effects, increasing and decreasing the SCC, on one segment. And competition may lead the SCC to be a compromise, the real one. Thus the SCC incarnating the above mentioned self-adjustments is expressed as klh (i + 1, t + 1) = klh (i + 1, t ) − A1 ∆Shl (i, t )/(Qhl (i + 1, t + 1))n

+ A2 (∆Shl (i, t ) + ∆Shi +1 (ilast , t ))/(Qhl (i + 1, t + 1))n

(5)

where A1 and A2 are the so-called impact factors. The former represents the impact of the undergone state of segment i + 1 on the SCC of itself, and the latter denotes the impact of both the undergone state of segment i and the injecting stream of rank h + 1 on the SCC of segment i + 1. Obviously, the impact factors of the upriver neighboring segment i and the injecting branch stream can be equated for the sake of simpleness. In our model, any one upriver segment can be chosen as the source of the stream. And we assume that the SCC, QES and water discharge are changeless. That indicates the equations Qhl (1, t ) = Qhl 0

(6)

∆Shl (1, t ) = ∆Shl 0

(7)

and are satisfied for all possible h, l and t. So our model is constituted by Eqs. (1)–(3) and (5)–(7). The core spirit of the model consists in the description of self-adjustability of a river in the sediment transport process, which is manifested by the adjustment of SCC that is expressed by Eq. (5). Similar self-adjustability often dominates the dynamics of a completely open system. The simulation work may not only help us to understand the nature of sediment transport in river networks, but also provide us with some illumination for the general dynamics of a self-adaptive system elsewhere. The simplicity of this model may give prominence to relationships between the main parameters and/or variables of a river such as QES, SCC and water discharge. These relationships may be submerged in the outputs of the models if they include a good many of the river parameters.

6424

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

Fig. 1. Conceptual flowchart to show the simulation processes.

To let the readers know how our model works, a simple flowchart (Fig. 1) is provided to explain the evolution algorithm. Some details of the main approaches in each time step are presented as follows. (1) Confluences: according to the continuity of flow, the outflow of segment i in time step t will become the one of the inflows of segment i + 1 in time step t + 1. Therefore, for each i we set Qhl (i, t + 1) = Qhl (i, t ) to obtain one inflow of segment i + 1, and then the calculations by Eq. (1) are executed at each of all nodes to obtain the corresponding outflow. (2) SCC adjustment: based on the QES in time step t the pre-calculation of SCC for each segment in time step t + 1 is conducted by Eq. (5). (3) State adjustment (calculation of QES): the calculations by Eqs. (2) and (3) are carried out synchronously at each of all nodes to obtain the corresponding QES. Before conducting the first time step, A1 and A2 , only two free parameters of this model are evaluated in a fixed manner, the values of Qhl (i, 0) for all possible h, l are given randomly or regularly to simulate different processes of river network. From the second time step onwards, all parameters are confined. In any one time step while the program is running, the discharge of the highest rank branches (among the discussed branches) and/or some headstreams can be changed regularly or randomly to simulate the corresponding variations caused by rainfall (headstream conditions). Rainfall in the drainage basins can be assumed to be uniform or heterogeneous without large fluctuations (in the lower-water season). It may be important and interesting that one can know how the QES on a segment affects the state of the segments along the network downriver. That can be easily deduced from the aforementioned equations under the condition, Qhl (i, t + 1) = Qhl (i, t ), which indicates that the rainfall in the drainage area of the branch is changeless (in the lower-water season), that is,

∆Shl (i + 1, t + 1) = (1 − A1 )∆Shl (i + 1, t ) + (A1 + A2 )∆Shl (i, t ) − A2 ∆Shl (i − 1, t ) 0

+ A1 ∆Shi +1 (ilast , t ) − A2 ∆Shi +1 (ilast − 1, t ) − A2 ∆Shi +2 (i0last , t )

(8)

0 Shi +2

where i0 = ilast − 1, ∆ (i0last , t ) denotes the QES on the last segment of the i0 th branch of rank h + 2 that is about to inject into the last segment of the ith branch of rank h + 1. It follows that the QES on a segment can influence its two downriver neighboring segments in one time step. Especially, for a last segment in the branch of rank i, it is possible that such influence can reach the branch of rank i-2 in one time step. Therefore, the QES on a segment will finally influence the state of all of its downriver segments through which the stream passes when it heads to the sea. 3. Testing the model To test the model, the simulated QES, theoretically calculated with some real data, is needed for comparison with the really evolutive ones. However, any process on a natural river network is transient, because the water flow and the sediment concentration may change at any moment. The precise comparison seems very difficult to execute. But then, we suggest an approximate test method: use the transient result generated by the model to compare with the result observed in the real transient process in the lower-water season. Many databases provide the data of the stream flow and the corresponding sediment concentration. These are often obtained by taking an average over a month. It should be, in passing, pointed out that the conservation law for runoff is not likely to hold. We select the data to calculate the real QES in each segment. This QES can be roughly regarded as the evolutive result, since there exists no large fluctuation in the stream flow in this period. The instantaneous stream flow is close to the average, and thus the average can be approximatively taken as the initial and changeless condition. The simulated QES can be obtained by directly calculating with the model under the same initial values, and averaging the result in the preceding transient process. As this partial transience is short-term, for instance, several minutes, the simulated QES may be used to compare with the observed average one. So the data of the Yellow River in December, January, and February are, in the nature of things, the preferred ones for us. Because the data in the

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

6425

Table 1 The comparisons of the calculated QES, ∆S0C (i), with the observed QES, ∆S0O (i) and the steady QES, ∆S0S (i) in Feb, 1982 Segment TNH GD XH LZ XHY QTX SZS DK SHHK TDG HQ

∆S0O

∆S0C

10.898 84.069 −49.338 −38.879 21.800 −31.317 781.441 −577.783 −36.279 157.041 929.286

10.898 83.243 −47.944 −39.107 21.205 −30.732 773.308 −563.504 −42.927 155.609 921.792

∆S0S

k0 (Initial)

10.898 6.538 3.923 2.354 1.412 0.847 0.508 0.305 0.183 0.110 0.066

2.114 × 10 1.513 × 10−3 7.664 × 10−4 3.908 × 10−5 1.009 × 10−4 2.079 × 10−5 1.643 × 10−3 4.221 × 10−4 3.852 × 10−4 7.699 × 10−4 2.752 × 10−3

k0 (Steady) −4

2.114 × 10−4 2.774 × 10−4 3.584 × 10−4 8.047 × 10−5 9.202 × 10−5 9.422 × 10−5 7.774 × 10−5 7.507 × 10−5 8.305 × 10−5 8.788 × 10−5 1.188 × 10−4

The units of them are kg/s. The impact factors A and B are chosen as 1.0 × 10−4 and 6.0 × 10−5 , respectively. To be compared in detail, the initial and steady values of k0 are also listed in the 5th and 6th columns.

database [31] are only about the main stream and its first rank branches, we have to change our model by neglecting the adjustment of the SCC from the first rank branches. Then Eq. (5) becomes k0 (i + 1, t + 1) = k0 (i + 1, t ) − A1 ∆S0 (i + 1, t )/(Q0 (i + 1, t + 1))n + A2 ∆S0 (i, t )/(Q0 (i + 1, t + 1))n .

(9)

For the main stream, h = 0, the superscripts of Q and k are canceled. Before conducting the calculation, the data from the database are necessarily preprocessed. The database includes the stream flow Q (m3 /s) and sediment concentration Qs (kg/m3 ), obtained by averaging over a month in the mainstream and its first rank branches. Two situations about such data should be noted. One is the presence of the data of the branches with the absence of the data in some segments of the mainstream between these branches, another is the reverse. For the former case, the branches can be equivalently regarded as one branch, then the quantity of stream flow can be obtained by summing up that of the branches. For the latter case, the data of the branch can be roughly equal to zero. This processing may be acceptable, since the neglected stream flow is often small (it is the cause why the hydrological station is absent), the quantity of QES may be neglected. On the other hand, the stream flow is actually merged into that of the main stream. Moreover, to be taken for simulations, the QES should be calculated by Eq. (2) and the relationship, S = QS Q . We select the upper-middle reach, Tangnaihai(TNH)-Hequ(HQ), of the Yellow River as the studied reach for the relatively sufficient data in the database. Table 1 (only lists the results in the main stream of Feb, 1982) shows the observed QES in the segments that are denoted by the names of the hydrological stations on them, and the calculated ones yielded via averaging over 180 transient steps (3 min). It is remarkable that the latter is in good agreement with the former. The calculations performed with other data in this database support the same conclusion. For instance, we calculate the QES in three lowerwater periods including Nov and Dec of a year and Jan, Feb, Mar and Apr of the next year from 1980 to 1983. The comparison of these results with the observed ones show that in each lower-water period, the coincidence in Feb (the center of the lowerwater period) is the best. The deviation of the simulated results to the practical data get larger as it departs from Feb, in one period, in the two reversed directions. These indicate that our model can describe sediment transport in river networks in the lower-water season. The model can be employed to mimic the process in a natural river network. 4. Predicting the steady state Furthermore, it may be most interesting that after the transient process, the steady QES on each segment of this reach of the Yellow River can be predicted by this model. That is, the model can tell us what state, which may be not observed in a real network, will be if we keep the condition (water discharge) unchanged. After a super-long transient process, 5.62 × 105 time steps (seconds) or so, approximately 6.5 days, the QES on each segment converges to the steady state. The steady state in different segments are qualitatively the same — all the segments are in a scouring state. The average erosion decreases downriver, and tends to be in an erosion-sedimentation equilibrium in some segments below the studied reach. This may be the result of the competition stated above, and manifested by the scaling law,

∆S0 (i) ∝ eτ i ,

τ < 0,

(10)

as shown by Fig. 2. This dependence of QES exponentially on the number of the segment, i, can also be analytically deduced from Eqs. (8) and (9) under the simplifications: Q1i (ilast ) = Q1l (llast ) and ki1 (ilast ) = kl1 (llast ) for i 6= l. That is, ln |∆S0 (i)| = i ln

A2 A1

+ ln

|∆S0 | A1

.

(11)

The details are presented in Table 1. The SCC, k0 (i), in the segments have been adjusted to those as shown by the inset of Fig. 2 to the effect that the distribution of QES obeys the scaling law. To show this adjustment in detail, the steady and initial values of k0 (i) are also listed in Table 1.

6426

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

Fig. 2. The scaling law dominating the distribution of QES on different segments after the transience. The inset shows the SCC of the segments in Feb 1982 after the transience.

Fig. 3. The varying trends of (a) the QES and (b) the SCC on the main stream as we only consider the confluence of the first rank branches. The initial water flux is randomly chosen from 0.5–5.0. The calculation is conducted by increasing the water flux of the first branches with 0.1 and equably increasing the sediment concentration from 4.727 × 10−4 to 4.727 × 10−1 every 200 time steps in the total time steps, 60 000.

5. Understanding the sediment transport in the random rainfall processes As is well-known, sediment transport in a natural river is actually a transient process, due to changes in the runoff, at any moment, caused by random rain. As discussed in the last section, the scaling law describing the steady state implies a great deal of erosion or sedimentation, which predicts a catastrophic event might occur in the lower reaches of a river if the stream keeps unchanged. However, there is no such observation in natural rivers. Obviously, stochastic changes of runoff play a key role in preventing this trend from appearing. In this section, we will present a detailed understanding of this phenomenon. Let us first see about varying trends of the QES and SCC as the runoff changes. Fig. 3 (only considering the injections of the first rank branches into the main stream) shows the simulated results as the stream flow and the sediment concentration increase. These may mimic the transition of it from the lower-water season to the higher-water season. A strong impression of the figures on one may be that the varying trend of the SCC is almost opposite to that of the QES. It indicates the selforganized characteristics of the dynamics: the self-adjustment of the SCC conducted according to the previous QES tends to largely decrease the quantity of erosion or sedimentation, since the transient process towards the steady state is frequently interrupted by random changes of the runoff. In fact, this interruption will cause a rapid change of the QES. So the transient process can be naturally fallen into two successive phases. One is the just mentioned that we name the faster transience and followed by the other that we call the slower transience. For a given segment, the former will last for the stream flow getting a new stationary value, and then the latter starts and leads the dynamics to the steady state.

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

6427

Fig. 4. (a) The QES on the main stream as the confluence of the second rank branches is considered. The water flux and the SCC are the same for each of the second branches. The former varies with the increment of −4.9 × 10−3 every 200 time steps, the latter has a fixed value. The other parameters are set with fixed values: k0 (i, t ) = 2.217 × 10−4 , k1 (i, t ) = 4.70 × 10−4 , k2 (i, t ) = 4.54 × 10−4 , Q0 (1, t ) = 100.0 and ∆S0 (1, t ) = 5.0. The number of the segments in the main stream and each of branches of rank 1 is 30 and 50, respectively. (b) The partial magnification of (a) for 28th segment of the main stream from time step 1998 to 2086.

To show the details of the transience, we calculated the QES by considering the influences of the second order branches on the dynamics of the mainstream, and showed the results in Fig. 4(a). For the convenience of analytic work, we only considered the most simple case, that is, the stream flow and sediment concentration of each second rank branches are the same, and each of all the first rank branches has the same number of segments. Fig. 4(b), the partial magnification of Fig. 4(a) for segment 28, shows that the faster transience can be divided into three sub-phases. We call them the sub-faster transiencies (SFTs). It is not difficult to deduce a conclusion about what creates these SFTs. SFT (1) is caused by the changing branches (in the water flux) of first rank injecting, one by one, into this segment. It will last for all of the changing streams from its upriver branches of the first rank reach the segment; SFT (2) is induced by the more and more changed streams from the branches of rank 2 reach the segment. This SFT will go on up to its nearest upriver stream of rank 1 becoming stationary; then the transience enters into SFT (3), and comes to the end time step, at which the uppermost stream of rank 1 also becomes stationary. To demonstrate the validity of the above mentioned deductions, we define two relative variables, δ ∆S0 (i, t ) and δ A0 (i, t ). The former denotes the variation of ∆S0 from time step t to time step t + 1, the latter denotes that of the SCC. They are expressed as

δ ∆S0 (i, t ) = ∆S0 (i, t + 1) − ∆S0 (i, t ),

(12)

δ k0 (i, t ) = k0 (i, t + 1) − k0 (i, t ).

(13)

and

δ ∆S0 (i, t ) can be generally deduced from the expressions of the model and the aforementioned understandings, and will appear in the Appendix. The result is shown in Fig. 5(a), which is in good agreement with that obtained numerically and shown by Fig. 5(b), so the deductions prove reasonable. To show the details of the self-adjustment in the faster transience, we give the corresponding δ k0 (i, t ) in Fig. 5(c) that are calculated numerically. Since the SCC is adjusted by the QES, the response of the QES to the abrupt change of the stream flow lags behind that of the SCC for one time step. For a given segment, at beginning of SFT (1), the sediment carried by the inflow drops with the decrease of the stream flow, while the sediment carried by the outflow changes only according to the slighter adjustment of the QES. It leads to an increase in erosion or a decrease in sedimentation. The situation is opposite when the stream flow increases. These imply the self-organized characteristics of the sediment transport process: river trends to depress a mass of erosion when the stream flow increases, while it controls a vast amount of sedimentation when the stream flow decreases. At the beginning of SFT (2) the decrease of inflow becomes smaller by comparing with that at the end of SFT (1), which results in sudden adjustments of SCC and QES. From Fig. 5(c) one can see in this SFT the decrease of the SCC becomes smaller; however the erosion will go up or the sedimentation will go down for the nonlinear dependence of the sediment carried by the flow on the water flux. The outset of SFT (3) corresponds to the nearest injecting branch of rank 1 having a stationary water flux, and then all of the upriver branches become stationary one by one. In this process, the stream flow gradually trends to changeless, the increment of SCC also trends to be the smallest one, and therefore the transported sediment decreases. Of course, the natural fluctuation of water flux is attributed to all of the upriver branches of higher rank with random precipitation. This may weaken or erase the divisions of the three SFTs. However, our simulations suggest that the response

6428

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

Fig. 5. δ ∆S0 (28, t ) (a) obtained analytically and (b) numerically, and (c) δ k0 (28, t ) to show the three SFTs. Both insets indicated by the arrows of (a) and (b) are the magnifications of SFT (2), and show the absolute value of δ ∆S0 (28, t ) decease linearly.

of QES to the random changes of the discharge can hold, in the main, similar features of that to the regular changes. As the stream flow randomly increases, the QES manifests the opposite trend of faster transience, that appears as the water flux decreases. And that the frequent change of the stream will prevent slower transience occurring. One may now draw the conclusion that these understandings have captured the essentials of the transient process: a random change in stream flow leads to a relatively balanced sediment transport. 6. Conclusions and discussions In this article, we suggest a dynamical model for sediment transport in river networks in the lower-water season. The model is based on feedback between erosion and sedimentation occurring in the channel. The core of the feedback embodies that the sediment-carrying capability of the stream is adjusted by the undergone state, erosion or sedimentation. The model is tested by the observed data of the Yellow River in China via a comparison of these data with those produced by the preceding short-term transient process in the lower-water season. With this model, we can predict the steady state after transience. The steady state of erosion or sedimentation is shown to be dominated by a scaling law that the QES distributes exponentially along the channel in the downriver direction. For the river with A1 > A2 , the scaling law (the scaling constant A is τ = ln A2 ) implies a decay of the QES to zero, which corresponds to the transformation from erosion to sedimentation 1 at some segment of the channel. For the river with A1 < A2 , the QES, whether it denotes erosion or sedimentation, will be exponentially enlarged along the channel in the downriver direction. This may predicate a catastrophe due to a great deal of sedimentation or erosion in the downriver reach of the channel. However, the way into the steady state will be interrupted by frequent fluctuations in runoff, caused by random rainfall, which can be simulated by introducing random changes of water flux in the model. Our numerical and analytical investigations demonstrate that QES will rapidly respond to abrupt changes of the runoff in the so-called faster transience. Since the SCC is modulated by the QES, the response of the QES to abrupt changes in the stream flow lags behind that of the SCC for one time step. This results in the varying trend of QES to be opposite to that of water flux in the faster transience, the period that is from the beginning of the change of the stream flow to the time at which the stream becomes stationary. It will then be followed by the so-called slower transience, that leads the dynamics to the steady state. The basic characteristic

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

6429

of the dynamics, in the process that it varies from the lower-water season to the higher-water season, is that the varying trend of the QES is contrary to that of SCC. It embodies the self-organized mechanism of river networks: the self-adjustment of SCC trends to depress the rapid and vast increase of the erosion or sedimentation, and is important for the channel to hold a relatively steadier configuration. In spite of its simplicity, our model can provide us with a basic understanding of the dynamical process of the sediment transport in river networks. More simulated investigations and understanding of the natural sediment transport process are expected. Especially, work that can give some suggestions as to relevant hydraulic engineering is expected. Acknowledgements This study is supported by the National Natural Science Foundation of China under Grant No. 10565002 and the Program for New Century Excellent Talents in University of China under Grant No. NCET-06-0914. Appendix. The expressions of δ∆S0 (i, t ) for the three SFTs

δ ∆S0 (i, t ) caused by the change, δ Q = Q2 (t 0 ) − Q2 (t 0 − 1), of the stream flow of each branch of rank 2 in time step t 0 (Q2 (t + 1) = Q2 (t ) if t > t 0 ), can be obtained analytically via the so-called linearization approximation of Eq. (2) for the main stream (h = 0) in time step t + 1 by carrying out a Taylor expansion of its right-hand side in the neighborhood of Q0 (i, t ), and neglecting the higher-order terms in comparison with the linear terms. The linearized equations for the three SFTs differ mainly in the change of the stream flow of the mainstream in two neighboring time step, i.e., δ Q0 (i, t ) = Q0 (i, t +1)−Q0 (i, t ). For the lengthy derivation, we only present the final expressions without providing the details. For SFT (1): after the change of the stream flow, δ Q , of the second rank branches occurs in time step t 0 , the stream flow of the mainstream will follow

δ Q0 (i + 1, t ) = (t + 1 − t 0 )δ Q ,

t = t 0 , t 0 + 1, t 0 + 2, . . . , t 0 + i − 1.

(14)

This is due to the fact that the number of changing (in water flux) branches of rank 1 that visit the i + 1th segment of the mainstream increases one in every time step of this SFT. Therefore, the increment of the number of the changed streams of rank 2 that injected into this mainstream segment in every time step increases linearly with time. Considering this relation and conducting a linearization approximation, we have

δ ∆S0 (i + 1, t ) = 2δ Q (t (k0 (i + 1, t + 1))Q0 (i + 1, t ) − (t − 1)k0 (i, t + 1)Q0 (i, t ) − ki1 (ilast , t + 1)Q1i (ilast , t )) +

2 X (−1)p Ap Rp .

(15)

p=1

For SFT (2): from time step t 0 + i, more and more changed streams of rank 2 reach segment i + 1 of the main stream, and the number of such streams of rank 2 increases i in every time step. Therefore, the change of the stream flow of the main stream caused by δ Q can be expressed as

δ Q0 (i + 1, t ) = iδ Q ,

t = t 0 + i, t 0 + i + 1, t 0 + i + 2, . . . , t 0 + N − 1

(16)

where N is the segment number of each branch of rank 1. Obviously, the range of time here implies that the discussions are valid for i < N. So a linearization approximation gives δ ∆S0 as follows:

δ ∆S0 (i + 1, t ) = 2δ Q ((i + 1)(k0 (i + 1, t + 1)Q0 (i + 1, t ) − ik0 (i, t + 1)Q0 (i, t ) − ki1 (ilast , t + 1)Q1i (ilast , t ))) +

2 X (−1)p Ap Rp .

(17)

p=1

For SFT (3): the faster transience enters SFT (3) from time step t 0 + N, in which the upriver stream of rank 1 that is nearest to segment i + 1 of the mainstream becomes stationary. From then on, the number of the streams of rank 1 those injected into the discussed mainstream segment and evolved into stationary ones increases one in every time step, that is, the increment of the number of the changed streams of rank 2 that injected into this mainstream segment decreases linearly with time. So δ Q0 (i + 1, t ) becomes

δ Q0 (i + 1, t ) = (i − (t + 1 − (t 0 + N )))δ Q ,

(18)

t = t 0 + N , t 0 + N + 1, t 0 + N + 2, . . . , t 0 + N + i − 1. From the above discussions, one may know that the faster transience can last for N + i time steps. Similarly, we have

δ ∆S0 (i + 1, t ) = 2δ Q (t (k0 (i + 1, t + 1)))Q0 (i + 1, t ) − (t − 1)(k0 (i, t )Q0 (i, t )) +

2 X (−1)p Ap Rp . p=1

(19)

6430

X.-M. Wang et al. / Physica A 387 (2008) 6421–6430

In the above equations, R1 and R2 are introduced only for shortening the length of the equations. On account of the same reason, they are simply expressed as: R1 =

3 X

Rn1

(20)

Rn2 .

(21)

n=1

and R2 =

3 X n=1

The terms included in these two equations are: R11 = ∆S0 (i + 1, t )(Q0 (i + 1, t )/Q0 (i + 1, t + 1))2 , R21

(22)

= −∆S0 (i, t )(Q0 (i, t + 1)/Q0 (i, t + 1)) , 2

∆S1i

(23)

(ilast , t )( (ilast , t )/ (ilast , t + 1)) ,

R31

=−

R12

= (∆S0 (i, t ) +

Q1i

∆S1i

Q1i

2

(24)

(ilast , t ))(Q0 (i + 1, t )/Q0 (i + 1, t + 1)) , 2

(25)

R22 = −(∆S0 (i − 1, t ) + ∆S1i ((i − 1)last , t ))(Q0 (i, t + 1)/Q0 (i, t + 1))2

(26)

R32 = −(∆S1i (ilast − 1, t ) + ∆S2N (Nlast − 1, t ))(Q1i (ilast , t )/Q2N (Nlast , t + 1))2 .

(27)

and

Now δ ∆S0 (i + 1, t ) is completely expressed. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31]

L. Leopold, Am. J. Sci. 251 (1953) 606. I. Rodriguez-Iturbe, A. Rinaldo, Fractal River Basins: Chance and Self-Organization, Cambridge University Press, Cambridge, 1997. A. Maritan, A. Rinaldo, R. Rigon, A. Giacometti, I. Rodriguez-Iturbe, Phys. Rev. E 53 (1996) 1510. J.R. Banavar, F. Colaiori, A. Flammini, A. Giacometti, A. Maritan, A. Rinaldo, Phys. Rev. Lett. 78 (1997) 4522. S.S. Manna, B. Subramanian, Phys. Rev. Lett. 76 (1996) 3460. S.S. Manna, Physica A 254 (1998) 190. K. Sinclair, R.C. Ball, Phys. Rev. Lett. 76 (1996) 3360. S. Kramer, M. Marder, Phys. Rev. Lett. 68 (1992) 205. H. Takayasu, H. Inaoka, Phys. Rev. Lett. 68 (1992) 966. A. Maritan, F. Colaiori, A. Flammini, M. Cieplak, J.R. Banavar, Science 272 (1996) 984. G. Caldarelli, A. Giacometti, A. Maritan, I. Rodriguez-Iturbe, A. Rinaldo, Phys. Rev. E 55 (1997) R4865. A. Giacometti, Phys. Rev. E 62 (2000) 6042. E. Somfai, L.M. Sander, Phys. Rev. E 56 (1997) R5. A. Rinaldo, A. Maritan, F. Colaiori, A. Flammini, R. Rigon, I. Rodriguez-Iturbe, J.R. Banavar, Phys. Rev. Lett. 76 (1996) 3364. T. Sun, P. Meakin, T. Jøssang, Physica A 210 (1994) 24. T. Sun, P. Meakin, T. Jøssang, Phys. Rev. E 49 (1994) 4865. D. Dhar, Physica A 369 (2006) 29. A.B. Murray, C. Paola, Nature (London) 371 (1994) 54. T.B. Liverpool, R.C. Ball, S.F. Edwards, Europhys. Lett. 30 (1995) 181. P. Meakin, T. Sun, T. Jøssang, K. Schwarz, Physica A 233 (1996) 606. W. Zhou, Q. Zeng, H. Zhao, Y. Yan, J. Sediment Research (3) (1994) 1 (in Chinese). R.J. Batalla, Earth Surf. Process. Landf. 22 (1997) 121. B. Sivakumar, J. Chen, Earth Surf. Process. Landf. 32 (2007) 269. P. Rustomji, I. Prosser, Hydrol. Process. 15 (2001) 1003. X. Tang, D.W. Knight, J. Hydraul. Eng. 130 (2004) 305. Q.-C. Guo, Y.-C. Jin, J. Hydraul. Eng. ASCE. 125 (12) (1999) 1262. R.A. Bagnold, Proc. R. Soc. Lond. A 332 (1973) 473. R.J. Zhang, Sediment Dynamics in Rivers, China Water Power Press, Beijing, 1997 (in Chinese). Z. Cao, S. Egashira, P.A. Carling, Water Resour. Res. 39 (2003) 1029 and references therein. S. Jung, S.E. Spagnolie, K. Parikh, M. Shelley, A.-K. Tornberg, Phys. Rev. E 74 (2006) 035302 and references therein. http://www.loess.csdb.cn/hyd/user/index.jsp.