Supplied Sediment Tracking for Bridge Collapse with Large-Scale Channel Migration

: Here, we provide a numerical model that assigns an identiﬁcation number to trace sediments and also identify the source of sediment supply. We analyze the e ﬃ cacy of our model by reproducing the reach-scale ﬁeld observations from ﬂooding events in 2010 and 2016 that a ﬀ ected Kyusen Bridge over the Bebetsu River, Hokkaido, Japan. Our simulation results can successfully reproduce and trace the formation of bars caused by sediment supply in the study region. Our study also suggests a strong relationship between bank erosion rate, sediment supply and ﬂow-discharge. The bank erosion rate is higher when sediment supply increases, and bank erosion reduces as ﬂow discharge goes down. The model can also replicate the changes in a bed concerning sediment supply and was used to reproduce the bridge-abutment failure caused by the 2016 ﬂooding with large sediment supply and the bridge-pier failure caused by the 2010 ﬂooding with less nonlinear, unsteady hydrodynamics and subsequent bed and bank morphodynamics. The equations governing the flow field are the two-dimensional depth-averaged equations for mass and momentum conservation. Since these equations are the same as those on the homepage of iRIC [42], we omit their description in this manuscript. Here, we describe the sediment transport and bed evolution models that are strongly associated with tracer particle tracking. iRIC Nays 2DH uses equations in the boundary-fitted coordinate system, but for simplicity, we express the equations here in an orthogonal coordinate system. The changes in the bed and tracer concentrations are calculated using the active layer model [39,40] shown below. The active layer model is a model developed to calculate the bed deformation of a mixed grain size bed and gives the identification number for each grain size to calculate


Introduction
The stability of bridge foundations across the globe suffers from erosion problems and bridge failures during floods. Although many studies have been conducted on predicting local scour depths around piers, and theoretical and empirical equations and numerical techniques have been developed (see a review by Pizarro et al. [1]), the effects of stream channel instabilities such as lateral migration, river widening, bed degradation and aggradation have received less attention [2,3]. Channel instabilities are not due to local bridge effects but to larger scale interaction between flow and sediment transport. For example, 128 bridges were damaged in Hokkaido (Japan) during the heavy rainfalls in August 2016, 58% of which were accompanied by the erosion of road embankments connected to the bridges due to channel migrations ( Figure 1). In total, 96% of the damaged bridges were over steep slope rivers with bed slopes from 0.0025 to 0.04 (mainly, channel sections lying from lower parts of mountain valleys to upper parts of alluvial fans), and the river beds tended to rise due to large amounts of sediment transported from mountain regions [4]. These abovementioned data suggest that increased sediment supply can induce stream meandering and bank erosion in steep slope rivers. Therefore, an understanding of channel morphodynamics, including sediment supply effects, is important for the protection and maintenance of bridges. been proposed and applied to disaster research. Sensitivity analyses using the above numerical models revealed that channel shape is dependent on the fluctuation of flow discharge, erodibility of riverbanks, growth rate of bars, and expansion rate of vegetation. However, the effects of sediment supply on channel migration have not been analyzed numerically. In particular, it is difficult to separately evaluate the effects of sediments supplied from a mountain landslide and sediments supplied by riverbank erosion when both sediments are similar in size and material. Placing tracer particles on a bed and tracking them by identifying them by color, magnetic properties, and radioisotope properties is a widely used method for directly measuring sediment transport [35][36][37]. Although these studies have made a great contribution to understanding the characteristics of bed load transport such as travel distance and waiting time, detailed measurements of particle motion with sufficient spatiotemporal resolution are still limited to experimental scales because tracing particles in natural rivers is time-consuming and costly [38]. Numerical models are powerful tools that can overcome this limitation. Recently, Iwasaki et al. [38] verified the difference between an entrainment-based model that includes a probability density function and a flux-based In many rivers, sediment supply from mountain regions is decreasing due to human activities such as dam construction [5]. A decrease in sediment supply leads to riverbed degradation, including upstream migration of knickpoint [6][7][8], and suppresses the formation of gravel bars [9][10][11]. When landslides occur due to heavy rainfall, the landslide debris that cannot be captured by dams passes through a mountain river and is supplied to an alluvial fan. The increase in sediment supply develops incised meandering in a mountain bedrock channel [12][13][14] and increases avulsion frequency in an alluvial fan [15][16][17]. Church qualitatively suggests that meandering depends strongly on channel slope, grain size, bank strength, and sediment supply, based on field observations [18]. According to his channel planform classification, the sinuosity of a gravel river increases with an increase in sediment supply, and finally, a sinuous channel changes to a braided pattern [18]. Theoretical studies on channel migration show that channel curvature is a key factor for meander development and progression [19][20][21]. Although these theoretical models assume that the inner and outer banks migrate at the same rate, observation of channels meandering without changing their width has been poor [22]. In the past few decades, numerical models that include the effects of width changes [23][24][25][26], cutoffs [27,28], land accretion [27][28][29], vegetation growth [29][30][31][32] and bank strength [33,34] have been proposed and applied to disaster research. Sensitivity analyses using the above numerical models revealed that channel shape is dependent on the fluctuation of flow discharge, erodibility of riverbanks, growth rate of bars, and expansion rate of vegetation. However, the effects of sediment supply on channel migration have not been analyzed numerically. In particular, it is difficult to separately evaluate the effects of sediments supplied from a mountain landslide and sediments supplied by riverbank erosion when both sediments are similar in size and material.
Placing tracer particles on a bed and tracking them by identifying them by color, magnetic properties, and radioisotope properties is a widely used method for directly measuring sediment transport [35][36][37]. Although these studies have made a great contribution to understanding the characteristics of bed load transport such as travel distance and waiting time, detailed measurements of particle motion with sufficient spatiotemporal resolution are still limited to experimental scales because tracing particles in natural rivers is time-consuming and costly [38]. Numerical models are powerful tools Water 2020, 12, 1881 3 of 13 that can overcome this limitation. Recently, Iwasaki et al. [38] verified the difference between an entrainment-based model that includes a probability density function and a flux-based model that does not include it, and showed that both models exhibit almost the same advection-dispersion characteristics in a channel with alternate bar formation. This suggests that, in a channel with large morphological changes, tracer dispersion is dominated by the randomness inherent in bar morphodynamics, rather than the randomness of individual particle movement. Iwasaki et al. [38] traced bed load particles in a straight experimental channel without bank erosion, but tracer advection and diffusion in a channel with bank erosion have not been investigated. In this study, we use an active layer model [39,40], which is one of the flux-based models, to trace sediments supplied from upstream and riverbanks. This is the first attempt to combine basic tracer research limited to laboratory scale with disaster research in a natural river. The objective of this study is to evaluate the effect of sediment supply on channel migration and bridge damages through numerical tracer simulation.

Study Site
The study site is Kyusen Bridge over the Bebetsu River (43 • 35 59.0" N 142 • 37 00.1" E). The Bebetsu River is a steep gravel bed river located in Kamikawa Basin, Hokkaido, Japan. The catchment area and the annual average rainfall of the Bebetsu River are 195.6 km 2 and 1090 mm/year [41], respectively. Kyusen Bridge is located 23 km upstream from the downstream end of the Bebetsu River (confluence with Biei River), and a two-span bridge with a total length of 58.7 m. The catchment area upstream of the bridge is 76.4 km 2 .
An aerial photograph from the year 2009 shows that the width of the Bebetsu River is about 30 m, and the channel curves upstream of Kyusen Bridge (Figure 2a). There is also another small dam (or bed sill) with a height of 1.5 m, upstream of the bridge to prevent sediment-related disasters. Airborne LiDAR data measured in 2009 indicate that the right bank height (i.e., the difference between the elevation of the flood plain and the thalweg bed) was 7.4 m near the bridge; however, it was 3.4 m nearly 600 m downstream from the bridge, suggesting that the riverbed near the bridge may have been lowered due to dam construction.  In August 2016, a huge flood in the Bebetsu River caused by heavy rainfall in the Hokkaido region strongly affected channel morphology. The flood stream amplified the bend upstream of Kyusen Bridge and shifted the bend apex downstream (Figure 2b). Since the right bank of the bridge eroded about 175 m in length and 23 m in width, the right bank abutment collapsed, and cracks from twisting were found in the main girder of the superstructure (Figure 2c,d).
Our field survey shows that a large amount of gravel was deposited on the riverbed, and that the riverbed near the bridge rose about 3.5 m before and after the flood, forming a point bar on the inner bank of the bridge. D60 (60% of the particles are finer than this size) of the gravel deposited near the bridge was 43 mm, which was almost the same as D60 of the bank materials (37 mm). Hence, it is difficult to identify the source of the sediment that forms the point bar.
Kyusen Bridge also suffered damage during the 2010 floods. The central pier was damaged due to bed degradation in the 2010 flood. The peak discharge of the 2010 flood was smaller than that of the 2016 flood ( Figure 3). Although not observed, the amount of sediment production in the mountains and the amount of sediment deposition in the dam were probably different between the 2010 and 2016 floods. Therefore, we perform numerical experiments using a morphological model to analyze the effects of flow discharge and sediment supply on bridge damages.

Methods
To simulate the observed channel migration in the Bebetsu River, we used free software called iRIC Nays 2DH [42], which can capture nonlinear, unsteady hydrodynamics and subsequent bed and bank morphodynamics. The equations governing the flow field are the two-dimensional depthaveraged equations for mass and momentum conservation. Since these equations are the same as those on the homepage of iRIC [42], we omit their description in this manuscript. Here, we describe the sediment transport and bed evolution models that are strongly associated with tracer particle tracking. iRIC Nays 2DH uses equations in the boundary-fitted coordinate system, but for simplicity, we express the equations here in an orthogonal coordinate system.
The changes in the bed and tracer concentrations are calculated using the active layer model [39,40] shown below. The active layer model is a model developed to calculate the bed deformation of a mixed grain size bed and gives the identification number for each grain size to calculate

Methods
To simulate the observed channel migration in the Bebetsu River, we used free software called iRIC Nays 2DH [42], which can capture nonlinear, unsteady hydrodynamics and subsequent bed and bank morphodynamics. The equations governing the flow field are the two-dimensional depth-averaged equations for mass and momentum conservation. Since these equations are the same as those on the homepage of iRIC [42], we omit their description in this manuscript. Here, we describe the sediment transport and bed evolution models that are strongly associated with tracer particle tracking. iRIC Nays 2DH uses equations in the boundary-fitted coordinate system, but for simplicity, we express the equations here in an orthogonal coordinate system.
The changes in the bed and tracer concentrations are calculated using the active layer model [39,40] shown below. The active layer model is a model developed to calculate the bed deformation of a mixed grain size bed and gives the identification number for each grain size to calculate volumetric concentration. In this study, the grain size is treated as a uniform grain size, and instead, the identification numbers are assigned to different production sources to trace the sediments: where, t is the time, λ is the porosity of the alluvium, η is the elevation of the alluvial bed, q bx and q by are the sediment transport rates per unit width in the x and y direction, respectively, e a is the thickness of the active layer (i.e., bed load layer) and p ak is the volumetric concentration of tracer sediments in the active layer, and the subscript k is the identification number of tracer sediment. In this study, sediments supplied from upstream are set to k = 1, sediment on the right bank near Kyusen bridge are set to k = 2 and sediments on the banks and bed other than the above are set to k = 3. p Ik is the concentration of sediment tracer between the active layer and the lower layer (transition layer). The concentration in the active layer p ak is used as p Ik when ∂η/∂t > 0, and the concentration in the transition layer p tk is used as p Ik when ∂η/∂t < 0. The thickness of the transition layer e t changes with the riverbed deformation, and p tk is stored as the concentration of the nth deposited layer p dk,n , when the thickness exceeds the specified value e d . Hence, the total thickness of each tracer sediment T k is calculated by e a p ak + e t p tk + n e d p dk,n /(1 − λ).
The components of the bed load transport rate are calculated using the relations: where U x and U y are flow velocities in the x and y directions obtained from the simulation of the two-dimensional flow field, and q bs and q bn represent the bed load transport rate per unit width in the streamwise and transverse directions, respectively. q bs and q bn can be estimated using the equations proposed by Meyer-Peter and Müller [43] and Hasegawa [44]: q bn = q bs u bn u bs − τ * c µ s µ k τ * ∂η ∂n (6) where, τ * is the Shields number (dimensionless shear stress) obtained from the calculation of the two-dimensional flow field, τ * c is the critical Shields number (here, set to 0.05 [45,46]), s g is the specific gravity of submerged sediment, and d is a grain diameter. (u bs , u bn ) are the along-streamwise and transverse components of the vector of near-bed velocity, ∂η/∂n is the local bed slope in the transverse direction, µ s is a static friction coefficient (here, set to 0.7) and µ k is a kinematic friction coefficient (here, set to 0.7). The ratio of u bs to u bn is calculated as N * h/r, where N * is a parameter dependent on the strength of the secondary flow (here, set to 7 [47]) and r s is the radius of streamline curvature. The bank erosion rate is simulated as the sum of two methods: one is a method of calculating the sediment transport rate on the riverbank from the Shields number, and the other is a method of calculating the slope failure rate from the angle of the riverbank slope. The sediment transport rate on the riverbank is estimated by Equations (5) and (6). The slope failure rate is estimated by a simple slope failure model [30,42]. If the angle of bed slope θ exceeds a critical angle θ c (which is approximately equal to the angle of repose of submerged gravel in general; here, set to 35 degrees), the bed elevation is corrected to satisfy the local slope of θ c . Since this simple slope failure model does not include the influence of cohesive materials and pore water pressure, the model is limited to non-cohesive banks, such as this study site.
The computational domain is 1200 × 100 m in the streamwise and transverse directions, respectively, constructed on the boundary-fitted coordinate system with 150 × 40 grids. The initial bed topography is made from the LiDAR data measured in 2009. The grain size of the bed and banks is set as 40 mm, based on the survey result after the flood. The dam is treated as a fixed bed. Due to the unavailability of the observed sediment supply from upstream, we need to assume the sediment supply rate. Sediment produced on the mountain slope passes through the upstream mountain channel and is supplied to the study site. Therefore, the amount of sediment that cannot be transported by the upstream mountain Water 2020, 12, 1881 6 of 13 channel will not be supplied to the study site. We assume that the sediment supply rate is equal to the sediment transport capacity calculated from Equations (6) and (7) using the Shields number at the upstream end of the computational region ( Figure 4). The initial thickness of k = 1 and 2 (i.e., tracer sediments) is 5 m.
the unavailability of the observed sediment supply from upstream, we need to assume the sediment supply rate. Sediment produced on the mountain slope passes through the upstream mountain channel and is supplied to the study site. Therefore, the amount of sediment that cannot be transported by the upstream mountain channel will not be supplied to the study site. We assume that the sediment supply rate is equal to the sediment transport capacity calculated from Equations (6) and (7) using the Shields number at the upstream end of the computational region ( Figure 4). The initial thickness of = 1 and 2 (i.e., tracer sediments) is 5 m.
In Run 1, we use the 2016 flood discharges and compare the calculated and observed bank erosion lengths and widths to validate the numerical model. In Runs 2 to 5, we simply use two constant discharges and two supply rates shown in Table 1 to investigate the model sensitivity.

Runs
Flow

Results
In Run 1, where the 2016 flood hydrograph was given, the channel migrated with time, and the flood stream reached the right abutment 4 h after the start of the simulation (Figures 5 and 6). In the early stage of the hydrograph, the water flowed along the thalweg of the initial channel (Figure 5a, b). Because the bed elevation decreased and the channel width narrowed due to the dam construction before the flood, the flow velocity was very fast and exceeded 7 m/s (Figure 6b). The high-velocity zone moved toward the outer bank at the channel bend upstream of Kyusen Bridge, from where the bank erosion started (Figures 5c and 6c). After that, the bend apex gradually migrated downstream, with its amplitude increasing (Figures 5d, e, and 6d, e). The increase in the channel width and the decrease in the flow discharge reduced the flow depth and velocity, and consequently, the bank erosion rate (Figures 5f and 6f). The calculated bank erosion width and length on the upstream bend were 22 and 190 m, respectively, which were in good agreement with the observed results (23 and 175 m). A point bar formed near the left bank of the bridge 2 h after the start of the simulation, and the width of the bar increased until after 6 h, and then, almost reached equilibrium (Figure 7). The initial thalweg bed increased by up to 4 m, which is also in agreement with our field survey observations (3.5 m).
A meandering thread of deep water depth and high-velocity flow shown in the downstream of Kyusen Bridge (Figures 5c and 6c) suggests the formation of alternate (single-row) bars. Previous In Run 1, we use the 2016 flood discharges and compare the calculated and observed bank erosion lengths and widths to validate the numerical model. In Runs 2 to 5, we simply use two constant discharges and two supply rates shown in Table 1 to investigate the model sensitivity.

Results
In Run 1, where the 2016 flood hydrograph was given, the channel migrated with time, and the flood stream reached the right abutment 4 h after the start of the simulation (Figures 5 and 6). In the early stage of the hydrograph, the water flowed along the thalweg of the initial channel (Figure 5a,b). Because the bed elevation decreased and the channel width narrowed due to the dam construction before the flood, the flow velocity was very fast and exceeded 7 m/s (Figure 6b). The high-velocity zone moved toward the outer bank at the channel bend upstream of Kyusen Bridge, from where the bank erosion started (Figures 5c and 6c). After that, the bend apex gradually migrated downstream, with its amplitude increasing (Figure 5d,e, and Figure 6d,e). The increase in the channel width and the decrease in the flow discharge reduced the flow depth and velocity, and consequently, the bank erosion rate (Figures 5f and 6f). The calculated bank erosion width and length on the upstream bend were 22 and 190 m, respectively, which were in good agreement with the observed results (23 and 175 m). A point bar formed near the left bank of the bridge 2 h after the start of the simulation, and the width of the bar increased until after 6 h, and then, almost reached equilibrium (Figure 7). The initial thalweg bed increased by up to 4 m, which is also in agreement with our field survey observations (3.5 m).
A meandering thread of deep water depth and high-velocity flow shown in the downstream of Kyusen Bridge (Figures 5c and 6c) suggests the formation of alternate (single-row) bars. Previous studies have shown that, unlike point bars formed by the channel curvature, alternate free bars are formed due to the instability between the bed and flow, and migrate downstream e.g., in reference [48]. The meandering thread moved with the downstream migration of the alternate free bars (Figure 5c-e), making relatively uniform incisions in both banks in the longitudinal direction. Section 3. Figure 8 reveals that the sediment particles released near the upstream end passed over the small dam and deposited near the inner bank of the channel bend because of the effect of the secondary flow included in Equation (6), and caused the formation of the point bar. The deposited area of tracer particles gradually expanded toward the outer bank, and the flow path also migrated outward so as to be pushed out by the development of the point bar (Figures 5 and 6). These findings clearly show that the sediment supply induced the development of a point bar and the lateral migration of the sinuous channel.     In order to investigate the influence of the supplied sediments on the bed aggradation, we plotted the total thickness of the supplied sediments T k=1 using the tracer technique mentioned in Section 3. Figure 8 reveals that the sediment particles released near the upstream end passed over the small dam and deposited near the inner bank of the channel bend because of the effect of the secondary flow included in Equation (6), and caused the formation of the point bar. The deposited area of tracer particles gradually expanded toward the outer bank, and the flow path also migrated outward so as to be pushed out by the development of the point bar ( Figures 5 and 6). These findings clearly show that the sediment supply induced the development of a point bar and the lateral migration of the sinuous channel. Next, we plot the total thickness of the sediments on the right bank near Kyusen Bridge (Figure 9). Sediments produced by bank erosion deposited on the alternate bars downstream of Kyusen Bridge, but hardly deposited on the point bar. This suggests that the sediments from the right bank affected the erosion of riverbanks downstream of Kyusen Bridge.
The results of the sensitivity analysis (i.e., Runs 2 to 5) show that the bank erosion width and Next, we plot the total thickness of the sediments on the right bank near Kyusen Bridge T k=2 (Figure 9). Sediments produced by bank erosion deposited on the alternate bars downstream of Kyusen Bridge, but hardly deposited on the point bar. This suggests that the sediments from the right bank affected the erosion of riverbanks downstream of Kyusen Bridge.

Discussion
Rainfall intensity mainly controls landslides in mountain regions [49]. However, upstream channel topography and artificial structures such as dams affect the arrival time of landslide debris [50][51][52]. Sediment supply to an alluvial fan does not necessarily depend on rainfall intensity or flow The results of the sensitivity analysis (i.e., Runs 2 to 5) show that the bank erosion width and riverbed variation in the cross-section of Kyusen Bridge strongly depends on the flow discharge and sediment supply rate ( Figure 10). The bank erosion width increased over time, but the increased rate of the erosion width gradually decreased in the runs with sediment supply (i.e., Runs 2 and 3 in Figure 9a) because the point bar in the upstream bed approached its equilibrium. The bed elevation increased and the erosion width at 7 h was larger when the flow discharge was larger in the runs with sediment supply (Runs 2 and 3 in Figure 10b). The channel incised diagonally downward and the erosion width at 7 h was larger when the flow discharge was larger in the runs without sediment supply (Runs 4 and 5 in Figure 10b). The results of Runs 4 and 5 resembled the damage from the 2010 flood (peak discharge is 237 m 3 /s shown in Figure 3), in which the central bridge pier was damaged from riverbed scouring. Note that our model underestimates the local scour depth around the pier, as the 3D vortexes around piers are not included in our model. In addition, comparing the runs performed under the same flow discharge (i.e., Run 2 vs. Run 4 and Run 3 vs. Run 5), we can find that the deposition of supplied sediments not only decreased the bank erosion depth, but also significantly increased the bank erosion width ( Figure 10).

Discussion
Rainfall intensity mainly controls landslides in mountain regions [49]. However, upstream channel topography and artificial structures such as dams affect the arrival time of landslide debris [50][51][52]. Sediment supply to an alluvial fan does not necessarily depend on rainfall intensity or flow discharge. It is possible that during the 2010 floods, the landslide deposits remained in mountain

Discussion
Rainfall intensity mainly controls landslides in mountain regions [49]. However, upstream channel topography and artificial structures such as dams affect the arrival time of landslide debris [50][51][52]. Sediment supply to an alluvial fan does not necessarily depend on rainfall intensity or flow discharge. It is possible that during the 2010 floods, the landslide deposits remained in mountain valleys and upstream dams and did not reach the study site. On the other hand, in the 2016 flood, landslide sediments might have been transported to the study site as a 'large sediment wave' because the upstream of the dams was filled with sediments discharged in 2010. Although we only analyzed the extreme scenarios of sediment supply, i.e., supply capacity or no supply, channel instability against various sediment supply scenarios, such as the effects of climate change and dam removal, will be an exciting challenge in the future. The model proposed in this study is capable of tracking sediments from different production sources, and hence, can be a powerful tool in carrying out the above study.
Eke et al. [29] numerically analyzed channel migration and width adjustment under bankfull conditions, and suggested that there are two patterns of meandering: "bank-pull" and "bar-push." In the bank-pull pattern of meander formation, riverbank erosion occurs first, deposition on the opposite bank accompanies this erosion and the channel migrates laterally. In the bar-push pattern, the water flow is pushed out toward the outer bank by the growing point bar, and the channel migrates laterally. Inoue et al. [13] and Mishra et al. [14] showed that large sediment supply caused "bar-push" channel migration, increasing the lateral erosion rate in a mixed alluvial-bedrock valley. The case of Kyusen Bridge can be interpreted as a bar-push channel migration induced by the large sediment supply.
Ashworth et al. [15] and Bryant et al. [53] experimentally showed that channel avulsion occurs in rivers with various slopes (bed slope = 0.008-0.13), including the slope in our study site (bed slope = 0.016). They also showed that larger sediment supply increases avulsion frequency in alluvial fans [15,53]. In addition, Church suggested that an increase in sediment supply causes channel braiding [18]. However, neither channel avulsion nor channel braiding occur despite a large sediment supply, and the meandering amplitude simply increased in our study site. Perhaps avulsion frequency depends not only on supply but also on bank height, that is, avulsion is unlikely to occur even if the riverbed rises a little in a deep valley. In the study site, because the bank height had increased due to riverbed degradation before the 2016 flood, channel avulsion may not have occurred. This suggests that long-term sediment supply also affects channel migration and bridge damages during floods.

Conclusions
Previous studies exploring the detailed tracer-dispersal have mostly been performed in controlled spatial and temporal environments in labs. In this study, we present a numerical model that can trace the sediment supplied from upstream, thereby identifying the sediment supplied by upstream landslides or riverbanks. We implemented an active-layer model and used identification number k to identify sediment sources if the sediment is supplied from upstream (k = 1) or riverbank (k = 2). The model estimated the effect of sediment and flow field in both vertical and lateral directions.
The proposed model successfully reproduced the reach-scale field observations of Kyusen Bridge over the Bebetsu river, Japan. Our simulation results suggest a strong correlation between bank erosion and the rates of sediment supply and flow discharge. An increase in the sediment supply rate contributes to the formation of a point bar, shifting sediment towards the outer bank and thereby determining the sinuosity of the channel. The numerical results show that bank erosion increases with the sediment supply rate and decreases with lower flow discharge without sediment supply variation.