Network response to disturbances in large sand-bed braided rivers

The reach-scale effects of human-induced disturbances on the channel network in large braided rivers are a challenge to understand and to predict. In this study, we simulated different types of disturbances in a large braided river to get insight into the propagation of disturbances through a braided channel network. The results showed that the disturbances initiate an instability that propagates in the downstream direction by means of alteration of water and sediment division at bifurcations. These adjustments of the bifurcations change the migration and shape of bars, with a feedback to the upstream bifurcation and alteration of the approaching flow to the downstream bifurcation. This way, the morphological effect of a disturbance amplifies in the downstream direction. Thus, the interplay of bifurcation instability and asymmetrical reshaping of bars was found to be essential for propagation of the effects of a disturbance. The study also demonstrated that the large-scale bar statistics are hardly affected.


Bar and channel dynamics in braided rivers
The complicated and dynamic network of bars and branches in large braided rivers poses a challenge to scientists and engineers.In particular, the morphological effects of river training and other human-induced perturbations on this network are still a puzzle.Both human-induced disturbances and intrinsic instability put the enormous social, economic, and ecological values of large braided rivers under pressure.For example, fertile land along the Brahmaputra River (India) and Jamuna River (Bangladesh) has been consumed by the rivers due to severe bank erosion (e.g., Sarker et al., 2003;Baki and Gan, 2012, Fig. 1c), while navigation is hampered by large and still unpredictable channel shifts.Furthermore, engineering works within and along these rivers are, in general, susceptible to failure by the massive rates of channel erosion and deposition.Despite the existence and application of basic engineering rules, attempts to control large braiding rivers have rarely been successful (Mosselman, 2006;Rah-man et al., 2012a).A crucial reason for this is the inability to predict migration and reshaping of bars and branches within the river.Another issue is that identifying morphological effects of human-induced disturbances is difficult, and in most cases it is impossible to isolate these from the autonomous morphodynamics.Enhanced insight and prediction capability of the dynamics within large braiding rivers are required to improve the success rate of river training and other engineering works in large braiding rivers, and to reduce undesired side effects and large-scale morphological reaction.
The dynamics within large braiding rivers include the interplay among bars, branches, islands, and floodplains (Bridge, 1993;Ashworth et al., 2000).A major role is played by channel bar bifurcations that distribute discharge and sediment through the braided channel network (Bolla Pittaluga et al., 2003).Commonly, bifurcations are unstable, meaning that the distribution changes over time (Kleinhans et al., 2013).The distribution of discharge and sediment determines the migration and reshaping of bars, and it determines the initiation and closure of branches (Schuurman and Kleinhans, 2015).At the same time, discharge and sediment distribution are controlled by the bifurcation topography and local flow pattern.For example, the branch with the smallest angle to the approaching flow is likely to experience the least amount of sedimentation and to become the dominant branch (Mosselman et al., 1995;Schuurman and Kleinhans, 2015).Bar migration and reshaping might change the local flow pattern, thus affecting the nearby upstream bifurcations through the backwater effect and nearby downstream bifurcations by changing the approaching flow direction.This starts a cascade of effects that links all bars and branches together.
Hence, a single disturbance in a large braiding river could affect an entire reach, beyond the extent of individual engineering projects.
The bar and branch dynamics of braided rivers have been studied by means of flume experiments (e.g., Fujita, 1989;Ashmore, 1991a), numerical modeling (e.g., Nicholas, 2013;Schuurman et al., 2013;Yang et al., 2014), and field observations (e.g., Bristow, 1987;Klaassen and Masselink, 1992;Klaassen et al., 1993;Ashworth et al., 2000;Best et al., 2003).However, in these studies, any artificial constraints and disturbances such as non-erodible (flume) walls, engineering works (Fig. 1a), and dredging were not considered.Also, these studies did not consider natural constraints such as rock outcrops (Fig. 1b) and vegetation.Another group of studies identified and explained the hydrodynamic and morphodynamic effects of engineering works, but without placing them in the wider context of a river reach or the interconnected network of branches and bars (e.g., Uijttewaal, 2005;Mosselman, 2006;Yossef and de Vriend, 2010;Rahman et al., 2012a, b), or only considered the spatial distribution of bank erosion (e.g., Bhuiyan et al., 2010).A third group of studies applied statistical analyses and metrics based on regime theory to describe the long-term effects of river engineering and human interferences on rivers.Commonly used metrics are, for example, mean channel width and longitudinal slope (e.g., Church, 1995;Brandt, 2000;Surian and Rinaldi, 2003;Ronco et al., 2010).However, these studies have not addressed the short-term response of bars and branches on the long-term equilibrium conditions.
Physics-based numerical models provide a way to explore the morphological effects of river training, discharge regulation, and other human-induced pressures beyond the scale of pilots and flume experiments, and without risks of undesired social and environmental impacts.By comparison of scenarios, modeling can be used to isolate the effects of a change in a single boundary condition, in model schematization, or in model settings, which is impossible in the field and subject to noise in flume experiments.Application of numerical models for decision making is common in, for example, the highly regulated river Rhine in the Netherlands.Furthermore, the natural behavior and general bar dynamics in large braiding rivers were successfully modeled by Nicholas (2013), Schuurman et al. (2013), andYang et al. (2014).However, the application of numerical models in morphologically dynamic rivers for decision making, especially in large braiding rivers, is still in its infancy.

Disturbances in braiding rivers
River training works and other human activities such as sand mining and discharge regulation are disturbances to a "natural" river, additional to disturbances caused by the intrinsic instability of braiding rivers described by, for example, Ashmore (1991b).If we consider a river reach on the order of 100 km and without downstream tidal influence, four groups of additional disturbances can be identified: (1) external at the upstream inflow, (2) external along the outer channel banks, (3) external at the downstream end, and (4) internal within the reach.
Discharge is one of the dominant external boundary conditions for a river, regarding the abundance of hydraulic geometry relations that use discharge as the independent parameter (e.g., Leopold and Wolman, 1957;Latrubesse, 2008).Discharge variation is attenuated by human-made hydropower dams and water storages, and affects the downstream morphodynamics (Brandt, 2000).Many river modeling studies have applied constant discharge, assuming a morphologically dominant or representative discharge exists that gives similar yearly morphodynamics to the "real" hydrograph (e.g., Nicholas, 2010;Schuurman et al., 2013).However, other studies have shown that discharge variation has a large effect on river morphology (e.g., Kiss and Sipos, 2007;Crosato and Saleh, 2010) due, among other things, to vegetation colonization on exposed bar sections (Gordon and Meentemeyer, 2006;Tealdi et al., 2011).Also, Egozi and Ashmore (2009) demonstrated that braiding intensity increased with increasing discharge, although this was temporary and braiding intensity decreased after the channel adapted to the new discharge.In the context of both river training and morphological studies, the effects of discharge variation on the braided river network are still largely unknown.
In addition, the direction of the flow pathway needs further attention.Asymmetrical inflow stimulates bar and bend formation, which has been deployed in flume experiments to generate meander bends (e.g., Peakall et al., 2007;Van Dijk et al., 2012).Inflow asymmetry enhances the initiation of steady bars and subsequent channel bending that propagates over a distance of at least several meander lengths (Van Dijk et al., 2012), but the direct effect of the disturbance damps rapidly in the downstream direction (Struiksma et al., 1985;Schuurman et al., 2016).Linear theory also explained that a disturbance in a river with, among other things, a sufficient width-depth ratio amplifies in the downstream direction (Struiksma et al., 1985;Crosato and Mosselman, 2009;Kleinhans and Van den Berg, 2011).However, this theory is based on the initial stage of bars on a flat bed, and its application to a developed braiding river is questionable.
The second group of disturbances involves bank erosion (Fig. 1c) and non-uniform channel width.Although braiding rivers are known for their dynamics of bars and branches within their braidplain (Lewin and Ashworth, 2014), channel migration and local widening of the braidplain are common (Khan and Islam, 2003;Ashworth and Lewin, 2012).Bank erosion results in local braidplain widening and thus higher braiding intensity as predicted by linearized models (e.g., Struiksma et al., 1985;Blondeaux and Seminara, 1985) and observed in nature (Xu, 1997;Ahktar et al., 2011).It could also result in fixation of bars (Wu and Yeh, 2005).At the same time, Rahman et al. (2012b) and Takagi et al. (2007) demonstrated that the bank erosion along the braidplain is linked to bar dynamics, as mid-channel bars steer the flow towards the braidplain banks.Bank erosion is also an important sediment source for mid-channel bars (Xu, 1997;Ahktar et al., 2011).Furthermore, lateral constraints by non-erodible banks can cause local bed degradation (Mosselman, 2006) and attract flow.In numerical models, both erodible floodplains and fixed walls have been applied with variable results.Erodible floodplains in the simulations of Nicholas (2013) resulted in major local braidplain widening.In contrast, Schuurman et al. (2013) reported a relatively small difference in bar pattern statistics between a braided channel with erodible floodplains and non-erodible walls.However, that model failed to produce sustained bar and branch dynamics that would cause bank erosion, because the grid resolution was too low to produce cross-bar channels and thus new bifurcations.Thus, a robust comparison of bar and channel dynamics in a braiding river with and without erodible floodplains is still lacking.
The third group of disturbances is related to engineering and training works, such as groynes (Mosselman, 2006), bridges (Bhuiyan et al., 2010, Fig. 1a) and sand mining.Although the structures are static, they introduce a disturbance to the original situation.River training is a common practice in meandering rivers to control meander migration and channel depth, but it is scarcely applied in large braiding rivers.This is due to the enormous dimensions of these rivers, the high costs, and the large uncertainties and risks of uncontrollable negative impact.Furthermore, bar and channel dynamics affect the efficiency of the training works (Nakagawa et al., 2013).Both the capability of the river to destroy the training works and the difficulty of predicting the effects of training works are problems engineers face in controlling large braiding rivers.

Research questions, hypothesis and approach
In this study, we analyzed the natural response of several simplified generic types of human-induced disturbances in large sand-bed braiding rivers using the numerical model Delft3D.
The main research question is, how do disturbances within and along a large braiding river affect the reach-scale braided channel network?Minor research questions are, what are the local effects of engineering works and how do these effects propagate through the braided channel network, and what are the effects of discharge attenuation and bank protection along the river on the bar and branch dynamics in braiding rivers?
The hypothesis is that the local bed level change due to a disturbance can be estimated using basic engineering rules.For example, a decline of channel width is expected to result in deposition upstream of the narrowing, bed degradation in the vicinity of the narrowing, and deposition further downstream.These bed level changes are expected to modify the discharge and sediment division over nearby bifurcations and the bifurcations become instable ("bifurcation instability").Then, the network aspect described by Schuurman and Kleinhans (2015) is expected to emerge: bars are reshaped and migration directions are altered by the bifurcation instability, which again affects both the upstream bifurcation through the backwater curve and the downstream biwww.earth-surf-dynam.net/4/25/2016/Earth Surf.Dynam., 4, 25-45, 2016 F. Schuurman et al.: Response to disturbances in braided rivers furcation through redirection of the approaching flow.Thus, it is hypothesized that a single disturbance within, along, or upstream of a braiding river reach causes a local bed level change that triggers a cascade of morphological changes in the network, eventually affecting the entire reach.We tested the hypotheses using a physics-based numerical model to systematically set up a "data set" of braiding rivers with different types and degrees of disturbances.These disturbances are loosely comparable with real cases, and our goal was to show the generic effects of perturbations rather than evaluating specific training works.We compared the morphodynamics in the disturbance scenarios with a reference case without the disturbance.In general, first the local morphological effects were determined, and second the larger-scale effects were identified and analyzed.

Numerical three-dimensional model
We used the physics-based numerical model Delft3D to construct a braided channel morphology for different scenarios of disturbance.This approach is similar to our earlier work in Schuurman et al. (2013) and Schuurman and Kleinhans (2015).The hydrodynamics were computed in three dimensions by applying conservation of momentum (Eqs. 1 and 2) and conservation of mass (Eq.3).The hydrostatic pressure assumption was adopted (Eq.4).
where x is the downstream coordinate (m), y is the lateral coordinate (m), z is the vertical coordinate (m), z w is the water surface level (m), u is the flow velocity in the x direction (m s −1 ), v is the flow velocity in the y direction (m s −1 ), w is the flow velocity in the z direction (m s −1 ), h is the water depth (m), C is the Chézy roughness (m 1/2 s −1 ), g is the gravity acceleration constant (m s −2 ), V h is the horizontal eddy viscosity (m 2 s −1 ), V v is the vertical eddy viscosity (m 2 s −1 ), ρ w is the water density (kg m −3 ), and P is the pressure (N m −2 ).The bed friction terms in Eqs. ( 1) and (2) are only applied in the first near-bed layer.Near the bed w = 0 m s −1 , and at the water surface w = dh/dt.A detailed description of the hydrodynamics and numerical scheme of Delft3D can be found in Lesser et al. (2004), Van der Wegen andRoelvink (2008), andDeltares (2009).
The bed level change in Delft3D is the result of sediment transport, bed slope effects, bank erosion, and mass conservation in the bed.The sediment transport rate in each grid cell is equal to the sediment transport capacity calculated with Engelund and Hansen (1967): where q s is the total sediment transport per unit width (m 2 s −1 ), U is the depth-averaged flow velocity in the streamline direction (m s −1 ), is the relative mass density of sediment underwater (−), and D 50 is the median grain size (m).
The amount of upstream sediment inflow at the upstream boundary was set equal to the local sediment transport capacity, which keeps the bed level along the upstream boundary constant.The transverse bed slope effect, which is the downslope pulling of sediment by gravity and essential in morphodynamic models (e.g., Struiksma et al., 1985;Talmon et al., 1995;Schuurman et al., 2013), is computed according to Koch and Flokstra (1981).After each time step, the bed level was updated using the Exner equation.To reduce computational time, an acceleration factor of 25 was used for bed level change on the basis of spatial sediment transport gradients, which is allowed because morphology changes much slower than hydrodynamics.The chosen acceleration factor has no significant effect on morphology (Roelvink, 2006;Schuurman et al., 2013).Sediment transport was only calculated above threshold water depths of 0.1 m.Grid cells with smaller water depth were considered to be inactive.Inactive grid cells reactivated when the threshold water depth was exceeded, either by water level rise or by a simplified formulation of bank erosion.Here, "bank erosion" of a dry grid cell occurred when a neighboring wet grid cell eroded, where 50 % of the incision in the wet cell was applied to the dry cells (Van der Wegen and Roelvink, 2008).This prevents unnatural effects of accidentally emerged cells.Test runs showed that the resulting morphology is relatively insensitive to the bank erosion percentage.

Default model settings and boundary conditions
We adopted the river dimensions and conditions from Schuurman and Kleinhans (2015) for the default scenario conditions (Table 1): a straight initially plane bed with 3200 m width, 80 km length for scenario series A and 40 km for scenario series B, an initial bed slope of 9.3 × 10 −5 , uniform fine sand (D 50 = 200 µm), and a constant discharge of 40 000 m 3 s −1 (which is close to the effective discharge of the Brahmaputra River).Fixed channel walls were applied by default.The computational domain was discretized by 50 m×20 m grid cells, and the water column was divided into seven grid cells with boundaries at constant fractions of the water depth (thus σ grid).Thus vertical grid resolution was relatively high at low water depths.The length of each grid cell was 2.5 times the grid cell width in order to keep the aspect ratio around 2 and to optimize computational speed at the same time.
The hydraulic boundary conditions were as follows: the inflow condition was set at the upstream and the water level was specified for the downstream.The upstream boundary was split into 20 separate boundary sections, i.e., eight grid cells per boundary section.For each boundary section, the amount of inflow was defined in a time series.The water level at the downstream boundary was based on initial uniform flow conditions.The hydrodynamic time step was 6 s; thus the morphodynamic time step was 150 s.
Following, for example, Nicholas (2003), Lesser et al. ( 2004), Nicholas et al. (2012), and Schuurman et al. (2013), a constant uniform bed roughness was applied, assuming bed forms were subgrid and thus captured by the bed roughness parameter.We applied a uniform Nikuradse k s of 0.15 m, which is recomputed to a Chézy roughness in Delft3D.
By default, the inflow and initial bed level were perturbed in the same way as in Schuurman et al. (2013).The upstream inflow disturbance was a random time-and spatially varying noise added or subtracted to the inflow at each of the upstream boundary sections.The standard deviation of the inflow disturbance was 0.5 % of the total discharge and the discharge disturbances changed every 2.3 days.The initial bed level disturbance was spatially varying with a maximum of 1 cm added to or subtracted from the initially smooth bed.The maximum bed level disturbance was 2.4 % of the initial water depth.As these disturbances had a much shorter time and spatial scale than the bars, they were considered noise rather than forcing.

Scenarios
An overview of the model scenarios is given in Table 2 and Fig. 2. We used two series of simulations: series A to simulate a "realistic" situation starting with a "realistic" bed topography constructed by Delft3D as described by Schuurman and Kleinhans (2015), and series B, which is a simplified situation starting with a regular symmetrical bar pattern.
The focus of our study is on series A, as it provides a better representation of reality than series B. The reason for using series B is that it enables analysis of the instability of initially regular and symmetrical bifurcations.
In series A, we applied different kinds of engineering works: a non-erodible construction on top of one of the bank attached bars (Run A4) or mid-channel bars (Run A5), protection at the upstream of a large mid-channel bar along one side of the bar (Run A2) or two sides (Run A3), and closure of one branch (Run A1).As a reference situation, the simulation was also run without engineering work (Run A0).
In series B, we also applied closure of one branch (Run B1), upstream bar protection (Run B2) and a reference case without engineering work (Run B0).Furthermore, we conducted two extra situations in series B: removal of a bar, as a simplified sort of sand mining (Run B3), and asymmetrical inflow at the upstream boundary (Run B4).
The dams and bank protection works were schematized as impermeable infinitely high dams, called "thindams" in Delft3D.Thindams are infinitely thin, only block the flow in the direction perpendicular to the dam, and do not add roughness to the flow parallel to the dam.Because of these properties, the dams and bank protection in the model schematization deviated from their real-world counterparts.same height as the bank, thus allowing overflow in the case of a sufficiently high water level.

Method for analysis
We used different methods to analyze and compare the model simulations.The first method was the use of 2-D time series of the bed level detrended by the initial slope.In addition, we used 2-D depth average flow velocities to identify the dominant branches and flow diversion by the bars and islands.Furthermore, metrics were applied to describe the bar and branch morphology: the active braiding index (ABI), active channel width, and bar height, following Schuurman et al. (2013).Additionally, we used a channel width ratio, which we defined as the ratio between the width of the widest branch and the width of the second widest branch.It gives a measure of the dominance of the largest branch.The ABI indicates the reach average number of parallel active channels, using cross-sectional average sediment transport rates as a threshold to discriminate between active channels and both bars and non-active channels.The active channel width is defined as the percentage of the braidplain width occupied by the active channels, indicating the degree of flow concentration within the cross section.The bar height is defined as the difference in height between the lowest 5 % and highest 5 % of the bed level, indicating the bed level difference between the channels and the bars.
The downstream celerity or propagation rate of the effects of the disturbances was compared with the propagation rate predicted from theory (Ribberink, 1985): with q s and h being the sediment transport and water depth above the disturbance.For the initial conditions, the predicted celerity, c, is around 17 km year −1 .Sarker and Thorne (2006) estimated the celerity in the Brahmaputra at 16 to 32 km year −1 , and compared this with the propagation celerity observed in the field after an enormous earthquake.They observed different propagation celerities for different types of responses: around 50 km year −1 for bed level change, which is much faster than the bars; 10 to 37 km year −1 for width adaptation; and 13 km year −1 for braiding index adaptation.
The evolution of bifurcations was analyzed using the ratio between the discharges and sediment transport towards the left and right branches.In our definition, a symmetrical bifurcation has a ratio of 1, and a ratio < 1 indicates that more than half of the discharge or sediment transport goes through the southern branch.transport capacity and incision in the branches.Together with this, braiding intensity declined to a year-averaged ABI of 2.5.Consequently, bar and channel dynamics declined and a dynamics equilibrium was reached, demonstrated by the steady statistics with only seasonal variation (Fig. 4).

Unconfined channel with natural discharge variation
Although the bar and branch dynamics declined after reaching the dynamic equilibrium, the river remained active with new branches formed by cross-bar flow and existing branches closed.The channel network statistics varied because of the seasonal water level variation.During low discharge, the ABI was around 1.5 to 2 and only around 10 % of the total channel width was transporting significant amount of sediment.During high discharge, the ABI increased to around 3.5 and active channel width increased to 30 %.
If we take a closer look at the short-term bar dynamics within a year, then we can identify characteristic processes of bar dynamics in each stage of the hydrograph (Fig. 5).During low discharge, sediment mobility was low and bar dynamics were negligible.When the discharge increased, the unit bars reactivated (feature A in Fig. 5).Also, large bartail limbs formed along and downstream of the bars (features B and C in Fig. 5).During the peak discharge period, many bars were overtopped and aggradation of the bars occurred as the flow velocity over the bars rapidly declined (features D and E in Fig. 5).At the same time, new branches formed by cross-bar flow.During the declining discharge period, these newly formed branches incised and widened (features F to J in Fig. 5), whereas other branches were closed by bars blocking their entrance (feature K in Fig. 5).Also, the bar margins became steeper as these branches incised.

Constant discharge
With a constant discharge of 40 000 m 3 s −1 , the time to reach a dynamic equilibrium reduced to around 13 months (Fig. 4).
From that moment, the network statistics were similar to the year-average equilibrium statistics reached after 3 years in the runs with variable discharge: an ABI of around 2.5, an active channel width of around 20 to 30 %, and a bar height of around 30 m.Thus, the bar pattern statistics were similar, but the exact pattern of the bars and branches was different and bar formation occurred at a higher pace.

Non-erodible banks
Figure 4 shows the bed level after 16 months for a case with fixed walls and with erodible floodplains, both with a constant discharge.Obviously, the exact bar and branch patterns were different and difficult to compare directly.In both runs, large mid-channel bars and bank-attached bars formed, and many sections were dominated by a single branch.However, the reach-averaged bed level along the non-erodible channel walls was clearly lower than in the case of erodible floodplains (Fig. 4f).On average, the incision along the non-erodible walls was around 6 m deeper than with erodible floodplains.Despite the erodible floodplains, overall incision along the initial channel occurred in both runs.
A comparison of the bar pattern statistics is given in Fig. 4c-e.Because of the bank erosion along the floodplains and thus larger channel width, the ABI was higher in these runs, and the active channel width and bar height were lower.The spatially average floodplain erosion distance was around 300 m after 16 months, which gives an annual braidplain www.earth-surf-dynam.net/4/25/2016/Earth Surf.Dynam., 4, 25-45, 2016

Inflow asymmetry
The effect of an asymmetrical inflow is illustrated in Fig. 6.In Run B4 (high discharge inflow along the south), the asymmetrical inflow caused asymmetrical reshaping of the most upstream bars.Bar-tail limbs along the dominating branches grew parallel to the prevailing flow and faster than along recessive branches.The expansion of the bar-tail limbs resulted in merging of bars, forming long and slim compound bars.After 6 months, the river reach was dominated by two parallel branches.
The many examples of asymmetrical deformation of bars indicate instability -resulting in bifurcation asymmetry -of the directly upstream-located bifurcations, which provoked instability of the directly downstream-located bifurcations.This generated a cascade of bifurcation instability, bifurcation asymmetry, asymmetrical bar deformation, and instability of the next bifurcation.This cascade started at the upstream boundary and propagated in the downstream direction with a fairly constant celerity of 0.3 km day −1 (Fig. 6d).
The development in Run B4 is in contrast to the development in Run B0 with uniform inflow.The upstream bars in Run B0 remained almost symmetrical, and the asymmetrical bar deformation started at the downstream end of the reach (Fig. 6d).This deformation propagated slowly in the upstream direction, which could only be the result of backwater effects that cause a reduction in upstream water level gradient and thus deposition.This backwater effect also occurred in Run B4 and interfered with the downstream propagating cascade triggered by the inflow disturbance.Although the downstream bars were more complex than the upstream bars -which is an indication of bifurcation instability as the instability results in closure of branches and merging of bars -the effect of the inflow disturbance in Run B4 on the bar shapes and network morphology throughout the entire reach was clear.

Branch closure
Disturbance by an engineering construction in one of the branches also affected the network morphology with its bars, bifurcations, and branches.The bed level evolution after building a dam at x = 20 km in the idealized situation is shown in Fig. 7a.As expected, the flow through the closed branch steered around the dam, causing major outflanking scour at the dam heads.The scour depth along the southern edge of the dam was 50 m, and along the northern edge around 30 m.The deeper scour hole in the southern branch could be explained by the non-erodible channel wall along the southern branch.At the same time, sediment was deposited in front of the dam.Also, deposition occurred downstream of the dam along the dam heads where flow decelerated and lost sediment transport capacity.These deposits formed long bar-tail limbs that reached the mid-channel bar downstream of the dam.The scour holes and deposition were a temporal response to the dam, because after 2 months, large bar-tail limbs of the mid-channel bar upstream of the dam extended and reached the dam (Fig. 7a).The flow was then guided by the upstream bar and diverted from the dam, instead of being blocked by the dam directly.
If we look at the network on the reach scale, we can clearly see an insignificant upstream impact and a major downstream impact of the dam.The bifurcations upstream of the dam remained stable and symmetrical (features A and B in Fig. 7), similar to the reference B0.Downstream of the dam, we could see, besides the first-order morphological response to the dam, a sequence that is the same as in Run B4 (Fig. 6a).The flow asymmetry caused by the dam propagated in the downstream direction, which we could see in month 2 at the long bar that starts from x = 22 km and extends to downstream of x = 35 km.This long bar was formed by merg-ing of bars.However, we could also see a difference with Run B4: the long bar was deformed significantly in month 6 by annexation of a bank-attached bar, which almost doubled the surface area of the bar complex.Eventually, the dam ended in the middle of a large bar complex.From the perspective of river management, this would be a disappointing development if the dam was a hydropower dam but would be a good development if the purpose of the dam was enhancement of river navigability.
Figure 8 shows the morphological development after dam construction in one of the branches of Run A1.When the northern channel was closed by the dam, the first effect was a 1 m water level impoundment upstream of the dam.This impoundment, and thus reduction of longitudinal water surface slope, did not result in a clear deposition upstream of the dam.Along the dam head, the southern channel incised several meters to compensate for the channel width loss (Fig. 8b).The local incision was around 6 m after 2 months.Further downstream, at around x = 35 km, the eroded sediment was deposited.Incision around the dam and downstream deposition occurred immediately after dam construction.For example, after only 2 days, a layer of 2 m was deposited (Fig. 8b).
The water level impoundment caused by the dam enabled flow from the northern branch to cross the large compound bar and drain into the southern branch.It was remarkable that a large bar (D2) blocked more than half of the remaining channel width, leaving the southern branch with a width of around 800 m.An explanation for this is the relict of bar D (D1) that redirects the flow towards the south and protects bar D2.
The incision and partial erosion of bar D provided sediment that was deposited downstream at x = 38 to 41 km, forming a large bar over the entire channel width.Furthermore, it resulted in a long bar-tail limb extending from x = 33 to 37 km parallel to the prevailing flow.As shown in Fig. 8b, the bed level change due to the dam propagated in the downstream direction with a celerity of around 0.5 km day −1 , which was on the same order of magnitude as predicted with Eq. ( 6).However, it was much faster than the migration rate of the bars themselves.Thus, the effect of the dam propagated through a change in flow field and a sediment wave initiated by the local incision.For example, the flow direction at x = 35 km changed around 30 • towards the north, favoring the northern branch around bar F in Run A1, instead of the southern branch in Run A0.Subsequently, this changed the flow at the confluence downstream of bar F and the approaching flow of bars G and H.
In addition to the local incision and deposition, as well as adjustment of the location, shape, and dimensions of individual bars, the dam also affected the bar pattern statistics (Fig. 9).For example, the active channel width near the dam reduced by around 50 % and local incision increased the bar height.This reduction in width-depth ratio decreased the ABI in the vicinity of the dam.The statistics near the dam in the section x = 25 to 40 km adapted in the first month to the new situation and remained constant afterwards.Downstream of x = 40 km, however, the statistics changed in the downstream direction and fluctuated around the statistics of the reference Run A0.

Bar protection
Figure 10a shows the bed level change in the idealized situation of Run B2, with protection of the bar head of one of the bars.After 2 months, the effect of the protection works on the bar shapes is still small, although the discharge distributions at the bifurcations downstream of the protection works became increasingly asymmetric (Fig. 10b).This started in cross section D after around 40 days, and in cross section E after around 50 days.Interestingly, we could see a periodic behavior in the discharge asymmetry.For example, after around 60 days, the discharge towards the southern branch (Q2) in cross section E was slightly higher than to the northern branch.However, after 80 days, Q1 was nearly twice as large as Q2, but after 110 days, Q2 became 4 times larger than Q1.It should be noted that this behavior was not caused by migration of the bars, as the bars hardly migrated within the 6 months.
During the first 4 months, the bar protection favored discharge through the northern branch, with Q1 being around 60 % of the total discharge (panel C, Fig. 10b).This is attributed to a lack of bar-tail limbs along the protected bar, while the other bars formed bar-tail limbs that diverge the flow.Between months 2 and 6, elongation of the protected bar and a lack of resupply from erosion of the bar head re-  sulted in severe trimming of the bar flanks.This increased the local branch widths and attracted discharge in expense of discharge through the northern branch in cross section C.
As the dominant branch has a tendency to meander through the channel -which we saw in the other runs -the shift in discharge from the northern branch to the southern branch also changed the discharge distribution in cross section D, and later in E. Through closer inspection of the network-aspect effect of bar protection in the simulation with regular bars, we can see that the bar protection affected the discharge and sediment distribution at the downstream bifurcation (panel C, Fig. 10b).Also, the bar along the right bank at x = 20 km was largely removed and pushed downstream, as the bar erosion attracted discharge towards the right branch.This had major consequences for the bars and branches further downstream, as much larger compound bars were formed compared to the scenario where bar protection was absent.
The effect of bar protection in Runs A2 and A3 is demonstrated in Fig. 11a.The effect of the bar protection after 3 months can be split into three sections: (1) local effects mainly covering the protected bar itself, (2) medium-distance effects at around x = 35 to 50 km with hardly any morphological effect, and (3) long-distance effects downstream of x = 50 km.If we compare the exact bar shapes and locations, the long-distance effects exceeded the medium-distance ef-fects.This might be partly because the bars downstream of x = 45 km were not yet developed at the moment we built the dam and thus might be more susceptible to small flow adjustments from upstream, although the bars at x = 40 to 45 km were also not yet completely developed at the moment of dam building and still have almost similar positions and shapes.Thus, this is another example of increasing morphological response in the downstream direction.
The local effects of bar protection works are illustrated in Fig. 11b.Without bar protection, the upstream bar head migrated around 1.5 km in 3 months in the downstream direction and the bifurcation angle slightly increased.At the same time, at the downstream side, a left bar-tail limb formed and extended around 2 km, almost reaching the next downstreamlocated bar.Also, the bar-tail limb along the right-hand side expanded by around 3 km in 3 months, superpositioned on relict inactive bar-tail limbs.
In the case of bar protection along the left-hand side, the protected bar-head side remained fixed, minor erosion occurred along the protected bar side, and a slim bar-tail limb formed.The left bar-tail limb had the same length as the left bar-tail limb in the case without bar protection, but with only half of the width.Because the bar head did not migrate in the downstream direction and the upstream bars still migrated in the downstream direction, the entrance of the left branch narrowed, reducing discharge towards the left branch.Conse- quently, discharge towards the right branch increased, causing deepening of the right branch.At the same time, the righthand side of the bar, including the right-hand side of the bar head, migrated in the downstream direction, similar to the case without bar protection.This bar-head erosion was accompanied with expansion of the upstream bar.Although the right branch became more dominant, the bar-tail limb along the right-hand side did not extent as far as without bar protection.
With bar protection along both sides of the bar head, the local effect along the left-hand side was similar to the single bar protection.At the right-hand side, the branch entrance deepened.Interestingly, the downstream bar-tail limb along the right-hand side had more resemblance to the case without bar protection than with the case of bar protection along the left-hand side.

Structures on bar
The effect of structures on a bar is demonstrated in Fig. 12.The structures blocked flow over the bar, thus diverting the flow around the bar.It is remarkable that the bar morphology near the structureswas hardly affected by the structures.For example, the large mid-channel bar at x = 30 to 35 km was almost the same for both runs.Nevertheless, the bar morphology further downstream was clearly different.The bar at x = 45 km showed a minor difference, the bar at x = 50 km showed more difference, and the bars downstream of x = 55 km were completely different.
If we compare the bar morphology downstream of x = 55 km with the bar morphology of Runs A0, A2, and A3 (Fig. 11a), then Run A4 had many similarities to Runs A0, A2, and A3, whereas Run A5 had a clearly different bar morphology, with one large bar between x = 55 km and x = 60 km.An explanation for this is that the structure in Run A4 was built on a relatively high bar which already had minor overflow, and thus the effect of the structure was relatively small.The structure in Run A5, however, was built in the middle of the river on a relatively low bar.

Sand mining
The simulations show that after removal of a complete sand bar, a new bar emerged on the empty spot (Fig. 13).The bar was formed by aggradation of the unit bar upstream of the gap.While the other unit bars migrated downstream and wrapped around the droplet bars, the unit bar forming the new bar was free to migrate downstream.Sediment deposition on top of the unit bar diverted the flow, stimulating further deposition on top of the unit bar.The new bar was shorter, but with longer bar-tail limbs, and lower than the original bar.The bar width was similar to the original bar.
Despite the appearance of a new bar, the sand mining significantly affected the bar and channel morphology further downstream.For example, the bank-attached bar downstream of the empty spot completely disappeared.The empty spot also attracted flow, resulting in enhanced channelization.This channelization stimulated elongation of the bars by bar-tail limb formation and bar flank trimming, resulting in merging of the droplet-shape bars into tall compound bars.Furthermore, an enormous bar complex was formed downstream of x = 21 km along the south by merging of bars, with around 75 % of the discharge flowing through the northern part (panel D, Fig. 13b).This merging of bars was much more pronounced than in the case without sand mining (Fig. 6).Upstream of the empty spot, however, there was no significant effect from the sand mining.

Overview of channel pattern statistics
Linear analyses (e.g., Crosato and Mosselman, 2009) suggest that the region near a disturbance should have a different channel pattern than regions further away.For example, the braiding intensity should be lower in the case of a dam due to reduction of the effective channel width.Indeed, this difference can be seen in Fig. 14, in which the bar pattern statistics are given for Run A1 to Run A5 relative to the reference Run A0.The most pronounced differences occurred, as expected, in Run A1 with branch closure: the active chan-nel width and ABI near the dam reduced, and the dominance of the dominant branch and bar height increased.This pattern extended to around 2 km downstream of the dam, after which the channel compensated and showed an opposite behavior: slightly higher ABI, higher active channel width, and less dominance of the dominant branch.Considering the five runs together, an increase in effect emerged with time and a large effect arose near the disturbance in regions c and d (Fig. 14).Further away from the disturbance, the effects on the channel statistics were relatively small.Figure 9 shows that the effects on the statistics in region f of Fig. 14 highly fluctuated, both along the river and with time.Although the specific location, shape, and planimetry of the bars and branches were clearly affected by the disturbances, the reach-average statistics were insensitive to the disturbances.Only the statistics of the region near the disturbance were affected.

Discussion
In this study, we conducted computer simulations of a large hypothetical sand-bed braiding river and perturbed the river in different ways.First, at the upstream inflow, the discharge was varied: from the simplest inflow condition of uniform and steady inflow to a steady asymmetrical inflow and a hydrograph.Second, along the channel we applied fixed walls and erodible floodplains, both perturbing the river in their own way.Finally, we perturbed the river internally by adding dams and bar protection works or by mining a bar.We analyzed the effects on a local scale, which was either near the upstream boundary, along the channel walls, or in the vicinity of the construction/mining, and on the reach scale.On the reach scale, the propagation of both the local morphological effects and the bifurcation instability was found to be important.
Figure 15 shows a conceptual model, inspired by the model results, of how disturbances affect the local bed level and how this effect propagates through the channel network by means of bifurcation instability and asymmetrical reshaping of bars.An adjustment to one bar, bifurcation, or branch initiates a sequence of adjustments in the downstream direction through (1) asymmetrical division of discharge and sediment transport over bifurcation branches, (2) elongation of the bar along the dominant branch, and (3) change in approaching flow towards the successive bifurcation.The celerity of this propagating wave was several orders of magnitude larger than the migration rates of the bars themselves, which is in agreement with the observations of Sarker and Thorne (2006) and with theory.A crucial driver behind the propagation was found to be the asymmetrical reshaping of mid-channel bars in response to an unequal division of discharge and sediment over the directly upstream-located bifurcations.The importance of bifurcations for the evolution of rivers, as well as the link between bifurcation asymmetry and bar asymmetry, has already been demonstrated by Schuurman and Kleinhans (2015).The novelty in this study is the downstream propagation of disturbances by means of bifurcation asymmetry, caused by bifurcation instability and bar reshaping.
With downstream propagation of the effect of a disturbance, the effect amplifies each time it destabilizes a bifurcation (Fig. 16).This way, even small disturbances, for example a relatively small dam on top of a bar, may cause a major impact on the bar and branch planimetry and dynamics, with closure and initiation of branches.In addition to the destabilization of bifurcations and asymmetrical bar growth, a change in bifurcation division almost instantaneously affects the division of downstream bifurcations, before the morphology responses.With a change in bifurcation division, the flow conditions automatically changed at the first downstream confluence and the next bifurcation.Such a purely hydrodynamic response is expected to decay with distance and to shift downstream simultaneously with the morphodynamic response (Fig. 16a).
Besides the propagating wave, we could identify different regions of morphological effects of disturbances (Fig. 17), starting with the morphological effect in the vicinity of the disturbance.This local effect is incision in the case of a flowblocking structure or deposition in the case of sand mining.The next region is the compensation region.Here, the local bed level change compensates for the upstream bed level change: deposition in case of upstream incision, or incision in case of upstream deposition.Further downstream, the effect of the compensation region is still present, caused by flow steering and thus alteration of the bifurcation stability.
Discharge variation had a relatively small effect on the long-term bar pattern, demonstrated by the bar pattern statistics that fluctuated around the steady statistics of the constant discharge runs.However, it affected the short-term bar dynamics and bifurcation stability, with the dominance of processes depending on discharge stage.It also doubled the time required to reach an equilibrium state, because for a large part of the year the discharge and water level were too low for significant bar dynamics.Based on these results, we conclude that it is correct to use a single representative discharge for long-term bar pattern analyses.For short-term modeling, on the order of months to a few years, it is preferable to use a hydrograph.The argumentation for this is based on a distinctive bar and branch dynamics within each stage of the hydrograph.This said, differences in bar and branch dynamics between the discharge stages were relatively small, and it was the sequence of importance of the processes that differed between discharge stages.For example, bar trimming and incision of the branches dominated during the declining limb of the hydrograph, whereas bar migration and formation of bar-tail limbs dominated during the rising limb of the hydrograph.
Erodible floodplains along large braiding rivers had a small effect on the bar and branch dynamics and statistics.
As predicted by theory of Struiksma et al. (1985), Blondeaux and Seminara (1985), and Crosato and Mosselman (2009), the braiding index increased with widening of the channel by bank erosion.The widening of the channel had a similar rate to that observed along the Brahmaputra.The small difference between fixed walls and erodible floodplains can be explained by the large initial channel width and the simulation time, considering that the simulation conducted only covered a couple of years.In the long term, erosion of the floodplains may have a major impact.Despite the similarity in bank erosion rates with the Brahmaputra, the highly simplified bank erosion procedure in Delft3D needs to be improved and more physics-based, with, among other things, accounting for bank instability and failure.Furthermore, the opposite process of bank erosion, which is bar-floodplain conversion by, for example, vegetation encroachment, is not considered in Delft3D.The necessity of this bar-floodplain conversion for channel migration was demonstrated by Schuurman et al. (2016), and the large effect of riparian vegetation on braiding river morphology has previously been demonstrated by, for example, Murray and Paola (2003) and Crosato and Saleh (2010).This missing mechanism must be included to fully understand the contribution of floodplain-channel interaction on the morphodynamics in braiding rivers.
This study shows that disturbances in large braided sandbed rivers affect the bar pattern -described by statisticsas well as the location, reshaping, and migration of individual bars and branches throughout the entire downstream river.This finding has implications for river training works and other interferences, as these may affect the river over a large distance, far downstream of the project area.Con-versely, this mechanism gives the opportunity to adjust the river over a long distance by means of a simple and lowcost disturbance.However, more research is necessary to develop quantitative predictors of reach-scale morphological responses to these types of disturbances.For this, it is necessary that fluvial morphologists, river engineers, and river managers join forces and collaborate more extensively.

Conclusions
The model simulations carried out in this study showed how the morphological effects of disturbances in and along a large braided sand-bed river propagate through the network of bars, branches, and bifurcations.The interplay between bifurcations and bars was found to be the essential mechanism driving the propagation.Different steps and zones of disturbance propagation can be recognized.First, a disturbance changes the local bed level and flow pattern over a relatively short distance.Second, these local effects destabilize nearby bifurcations, resulting in asymmetrical division of discharge and sediment transport at the bifurcations.Third, the asymmetrical division of discharge and sediment transport cause asymmetrical reshaping and migration of the bars, which in turn destabilize bifurcations further downstream.This cascade of bifurcation instability and asymmetrical bar dynamics amplifies in the downstream direction.In addition to the downward-amplifying morphological propagation, there is an instantaneous disturbance of cross-channel flow distribution along a reach that likely fades away from the disturbance location.However, morphological effects of this hydraulic disturbance are small.Furthermore, the effects of studied disturbances in the upstream regions are minor and only occur through the backwater effect.Furthermore, the channel pattern statistics only changed in the vicinity of the disturbance, remaining unchanged further upstream and downstream.
The study also showed that discharge variation in the form of an annual hydrograph affects short-term and barscale morphodynamics but hardly affects the longer-term and reach-scale morphology.In addition, using a highly simplified bank erosion method, the study demonstrated that floodplain interaction along large braiding rivers only causes minor effects on the bar and branch morphology within the river.
Furthermore, this study illustrated that physics-based models are useful tools for fluvial morphologists and engineers to explore not only the morphodynamic effects in the direct vicinity of disturbances such as training works but also the propagation of these effects on the reach-scale braided channel network.

Figure 1 .
Figure 1.Examples of disturbances in and along large braiding sand-bed rivers: (a) engineering works in and along the Jamuna River in Bangladesh: Jamuna Bridge; (b) geological constraint by a non-erodible bank along the Irrawaddy River in Myanmar; and (c) bank erosion along the Jamuna River in Bangladesh (courtesy of Royal HaskoningDHV).

Figure 2 .
Figure 2. Model scenarios with different types of disturbance: series A with training works in the channel starting with bars from Run A0; and series B with training works and sand mining starting with droplet bars.Runs A0 and B0 are undesturbed and act as reference cases.

AFigure 3 .
Figure 3. Evolution of a braided channel network and bar pattern in the case of erodible floodplains and yearly hydrograph.

Figure 4 .
Figure 4. Evolution of the braided channel network and bar pattern in case of (1) constant Q and fixed walls, (2) constant Q and erodible floodplains, (3) hydrograph and fixed walls, and (4) hydrograph and erodible floodplains.Panels (a)-(b) show the bed level after 16 months for cases (1) and (2).Panels (c)-(e) show reach-average channel statistics, and (f) shows average cross-sectional bed level profiles after 16 months for cases (1) and (2).

Figure 5 .
Figure 5. Detail of short-term bar and branch dynamics in the case of a hydrograph and fixed channel walls, starting from low discharge in month 36 and continuing to the peak discharge in month 40 and a declining discharge thereafter.

Figure 6 .Figure 7 .
Figure 6.Evolution of the braided channel network and bar pattern in Run B4 with asymmetrical inflow: (a) time series of the bed level in Run B4, (b) bed level in Run B0 after 6 months for comparison with last time step of (a), (c) depth-average flow velocity after 6 months in Run B4 and Run B0, and (d) discharge and sediment distribution over the branches with Q1 and Qs1 for the northern branches.The black arrows indicate the position of the propagating front.

Figure 8 .
Figure 8. Evolution of the braided channel network and bar pattern in Run A1 with branch closure: (a) time series of the bed level for Run A1 (with dam) and Run A0 (no dam) and (b) width-average bed level difference between Run A1 and Run A0, with negative values indicating more incision in Run A1 than in Run A0.

FFigure 9 .
Figure9.Effect of the disturbance in Run A1 on bar pattern statistics along the river, given as a ratio of Run A1 to the reference Run A0.The dam is located at x = 32 km.A moving average filter of 2.5 km was used.

Figure 10 .
Figure 10.Evolution of the braided channel network and bar pattern in Run B2 with bar protection: (a) time series of the bed level and (b) discharge distribution with Q1 the discharge through the northern branches.The black arrows indicate the position of the propagating front.

FFigure 11 .
Figure 11.Evolution of the braided channel network and bar patter in Run A2 and Run A3 with bar protections: (a) reach-scale development and (b) detail of the development of the protected bar with the bar protection (red line) and initial bar perimeter (black line).

Figure 12 .
Figure 12.Bed levels in month 15, 3 months after building of the structures in Run A4 at x = 23 km (top) and in Run A5 at x = 22 km (below).

Figure 13 .
Figure 13.Evolution of the braided channel network and bar pattern in Run B3 with sand mining: (a) time series of the bed level and (b) discharge distribution with Q1 the discharge through the northern branches.The black arrows indicate the position of the propagating front.

Figure 14 .Figure 15 .
Figure14.Difference in bar pattern statistics between Run A0 and Runs A1 to A5.The metrics are averages over time: averaged over first month (yellow), second month (blue), and third month (red), and averages over regions.The boundaries of these regions are defined at specific distances from the disturbance ("pert.").

Figure 16 .Figure 17 .
Figure 16.Responses to a disturbance in a braiding river: (a) hydrodynamic response by means of a change in approaching flow direction and discharge division over bifurcations and (b) morphological response by means of bar shape adjustment.

Table 1 .
Default initial and boundary conditions.
a Additional 1400 m in runs with floodplains.b 40 000 m for Runs B0 to B4.