Next Article in Journal
Water Resources Carrying Capacity Analysis of YarLung Tsangpo River Basin (I)
Next Article in Special Issue
New Hybrids of ANFIS with Several Optimization Algorithms for Flood Susceptibility Modeling
Previous Article in Journal
Development of Bubble Characteristics on Chute Spillway Bottom
Previous Article in Special Issue
Flood Hydrograph Prediction Using Machine Learning Methods
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Flood Routing in River Reaches Using a Three-Parameter Muskingum Model Coupled with an Improved Bat Algorithm

1
Department of Water Engineering and Hydraulic Structures, Faculty of Civil Engineering, Semnan University, Semnan 3513119111, Iran
2
Department of Biological and Agricultural Engineering, Zachry Department of Civil Engineering, Texas A&M University, 321 Scoates Hall, 2117 TAMU, College Station, TX 77843-2117, USA
3
Faculty of Natural Sciences and Engineering, Ilia State University, Tbilisi 0162, Georgia
4
Civil and Structural Engineering Department, Faculty of Engineering and Built Environment, Universiti Kebangsaan Malaysia, Bangi 43600, Malaysia
5
Civil Engineering Department, Faculty of Engineering, University of Malaya, Kuala Lumpur 50603, Malaysia
*
Author to whom correspondence should be addressed.
Water 2018, 10(9), 1130; https://doi.org/10.3390/w10091130
Submission received: 1 August 2018 / Revised: 14 August 2018 / Accepted: 17 August 2018 / Published: 24 August 2018
(This article belongs to the Special Issue Flood Forecasting Using Machine Learning Methods)

Abstract

:
Design of hydraulic structures, flood warning systems, evacuation measures, and traffic management require river flood routing. A common hydrologic method of flood routing is the Muskingum method. The present study attempted to develop a three-parameter Muskingum model considering lateral flow for flood routing, coupling with a new optimization algorithm namely, Improved Bat Algorithm (IBA). The major function of the IBA is to optimize the estimated value of the three-parameters associated with the Muskingum model. The IBA acts based on the chaos search tool, which mainly enhances the uniformity and erogidicty of the population. In addition, the current research, unlike the other existing models which consider flood routing, is based on dividing one reach to a few intervals to increase the accuracy of flood routing models. Three case studies with lateral flow were considered for this study, including the Wilson flood, Karahan flood, and Myanmar flood. Seven performance indexes were examined to evaluate the performance of the proposed Muskingum model integrated with IBA, with other models that were also based on the Muskingum Model with three-parameters but utilized different optimization algorithms. The results for the Wilson flood showed that the proposed model could reduce the Sum of Squared Deviations (SSD) value by 89%, 51%, 93%, 69%, and 88%, compared to the Genetic Algorithm (GA), Particle Swarm Optimization (PSO) algorithm, Pattern Search (PS) algorithm, Harmony Search (HS) algorithm, and Honey Bee Mating Optimization (HBMO), respectively. In addition, increasing the number of intervals for flood routing significantly improved the accuracy of the results. The results indicated that the Sum of Absolute Deviations (SAD) using IBA for the Karahan flood was 117, which had reduced by 83%, 88%, 94%, and 12%, compared to the PSO, GA, HS, and BA, respectively. Furthermore, the achieved results for the Myanmar flood showed that SSD for IBA relative to GA, BA, and PSO was reduced by 32%, 11%, and 42%, respectively. In conclusion, the proposed Muskingum Model integrated with IBA considering the existence of lateral flow, outperformed the existing applied simple Muskingum models in previous studies. In addition, the more the number of intervals used in the model, the better the accuracy of flood routing prediction achieved.

1. Introduction

Flood routing is fundamental to the design of structural, as well as nonstructural, flood control measures [1]. Routing involves the calculation of changes in the magnitude, velocity, and shape of a flood wave, as a function of time at one or several points of the river [2]. There are two types of flood routing methods: Hydraulic and hydrological. Hydraulic methods have complex computations and are more data-intensive, but describe the complete flood wave profile, whereas hydrological methods are much simpler, but yield the flood hydrograph at the end of a reach [3,4,5,6,7]. The hydrologic methods need only the inflow hydrograph for a river reach. A common hydrologic method is the Muskingum method, which has several versions with parameters ranging from two to five. The two and three parameter versions of the Muskingum methods are more popular. In recent years, optimization methods, especially evolutionary algorithms, for estimating the Muskingum parameters have been popular [3]. A brief background of such algorithms is now given.

1.1. Background

The flood routing models mainly include two different types of modeling: the hydraulic and hydrologic models. The hydraulic model is usually developed in a one or two-dimensional domain. Full three water shallow models and two diffusive models were used for an urban site, and the results had the same difference with each other because of different representation of a numerical and hydraulic method in the model algorithm process [5]. Hunter et al [6] successfully set three explicit hydraulic models based on the inertia, diffusive, and shallow water models for flood simulation. The results indicated that the models with the shallow water equation were simple and could provide good accuracy, for the prediction of depth and velocity of the flood [5]. Dottori and Todini [7] evaluated two-dimensional models based on the diffusive wave for urban floods, and the results indicated that the model could simulate the overall phenomenon well. Kim et al. [8] evaluated the different meshes in diffusive models, to investigate the effect of different meshes on the flood hydrograph. Prestininzi [9] applied the diffusive models based on the impulsive wave for inundation areas, and the results showed the model could simulate the flood conditions even in complex topography, based on a good match of simulated results with the observed data. Aricò et al. [10] applied diffusive wave equation based on 2-D numerical models for a slow varying flood, and the results showed that the simulated depth of the flood had a relatively good match with observed data. Classical, explicit finite differences in the hydraulic models and simple Muskingum model were used to investigate the flood routing [11]. It has been reported that the applied numerical methods had numerical instability, for some case studies. As a result, the Muskingum model showed superior performance compared to the same applied hydraulic models [11]. It has been reported that the 2-D models could be developed, based on the availability of enough information about topography and topology. With this information, a 2-D hydraulic model could successfully simulate the flood characteristics for different urban conditions. In fact, digital maps helped to identify all the required information about the boundary conditions, and to differentiate between numerous transitions within the urban hydraulic modes [12].
In addition, it has been reported that the 2-D models had had trouble in application to particular cases under small water depths, especially when the status comes close to wet/dry boundary conditions, so that there is a need for specific algorithms for simulation [6]. Costabile et al. [13] reported that the main advantages of the one-dimensional model over the 2-D models for flood routing, are a simpler run process and low computational time where topographic data was unnecessary.
Fassoni-Andrade et al. [14] considered the development of a one-dimensional model based on the equation of hydrological models, which include the continuity equation and mass equation, such as the equation of the Muskingum models [14]. It was observed that one-dimensional flow routing inertial models, based on the explicit solution were superior to the other models. These models simplified the Saint-Venant equation, and the main advantages of these models were good simulated results with a simple structure. Singh and Arvamuthan et al. [15] applied two hydraulic models that were developed based on the Kinematic and diffusion waves, in addition, the results were compared to the Muskingum model for the flood routing. The results showed that the simulation of hydraulic models was dependent on the kinematic wave number, so that when the value of this parameter was not considered based on accurate computation, the results for the hydraulic model could be worse than hydrologic models. Costabile et al. [16] reported that 2-D models could overcome the limitation of 1-D models, when the case study characterized as unsteady flow in irregular topography. The reports showed that if the flow was not one dimensional for the urban hydraulic, the one-dimensional channel network should be used instead of a one-dimensional model. In addition, the results showed the significant difference between 1-D and 2-D models to simulate the velocity and depths.
However, the results showed that the complex nonlinear form characteristic, numerical stability, high computational time, and complexity in the run process of hydraulic models, meant that the simpler and more accurate models have high importance [13]. In fact, the hydraulic models need to measure the flow depth and discharge based on applying stream gaging. These models are known as complex models and difficult to use, whilst the hydrologic models need only to use the discharge data. In addition, the hydrologic models can be effective for the initial planning level, where the measuring system is undeveloped for accurate measurement [13]. For example, Chatila [17] simulated flood routing based on the Muskingum model and EXTRAN hydraulic model. The hydraulic model developed was based on finite difference. Both hydrologic models and hydraulic models, were applied on simple and compound channels for flood routing. The results revealed that the Muskingum model had achieved higher accuracy compared with the hydraulic model because of its flexibility in calibration, where even the river bed geometry was not considered for this model. It has been demonstrated that the Muskingum model could simulate the peak discharge, achieving a close fit with the actual one, compared to the hydraulic model. Furthermore, it has been reported that hydraulic models are dependent too many assumptions, such as reach geometry, channel slope, and flow velocity, which causes the application of some hydraulic models to be limited to the specific case studies.
The Muskingum model is a useful and important hydrological model, due to its high accuracy and simplicity. Hydrological models could be accomplished after estimating the value of parameters, on the other hand, hydraulic models are required to simulate the complex boundary hydraulic conditions that causes an increase in the computational time [17].
Therefore, this model was used as a model with free access, fast computation, highly accuracy, and low cost. Furthermore, it can be used as a good tool, instead of complex hydraulic models, for flood simulation. Additional background of the application of the Muskingum model and its integration with an evolutionary algorithm, will be presented and discussed hereinafter [18].
Under a two parameter Muskingum method, Luo and Xie [19] applied the immune clonal selection algorithm (ICSA) for flood routing in a river in China, and found that the algorithm had faster convergence than the GA and PSO; and routed discharges had a high correlation with observed discharges.
Geem [20] obtained the two parameters of Muskingum method using a harmony search algorithm (HAS) for the Wilson flood in the USA, and obtained less root mean square between the predicted and observed discharges, than for GA and PSO, and less computational time.
Nelder-Mead simplex algorithm (NMSA) was considered for flood routing, and a case study in the USA [21]. The parameters of the Muskingum model were considered as decision variables, and the results indicated that the RMSE (root mean square error) based on NMSA decreased by 20%, compared to the genetic algorithm [21].
Karahan et al. [22] applied a hybrid of GA, HAS, and nonlinear programming to a three-parameter Muskingum method for flood routing in a river, and found the hybrid algorithm more accurate.
Orouji et al. [23] used a genetic programming algorithm (GPA) for flood routing by the Muskingum method, and showed that GPA was more accurate than GA and PSO. The Muskingum model by 4 parameters was considered under different case studies [23].
Four parameters were considered as decision variables, and the results indicated that the model based on the considered parameters with the genetic algorithm decreased the RMSE and mean absolute error (MAE) by 20% and 25%, respectively, compared to the nonlinear programming methods [24].
The hybrid PSO and harmony search algorithm was considered for flood routing [25]. The results indicated that the new hybrid method could increase the convergence velocity of the harmony algorithm, and decreased the error indexes, RMSE, and MAE, compared to the simple harmony and particle swarm algorithm.
Ouyang et al. [26] applied a hybrid of PSO and GA for Muskingum flood routing and showed that the hybrid algorithm was faster, and more accurately predicted the peak discharge and time to peak.
Under the three-parameter Muskingum method, Geem [27] found the harmony algorithm (HA) to have higher convergence than PSO and GA. Under the four-parameter Muskingum method, the Frog Leaping Algorithm (FLA) was found to have a lower computational error, than PSO and GA [28]. Niknazar and Afzali [29] used an improved Honey Bee Algorithm (IHBA) to optimize three parameters of the Muskingum method, and found it to be superior to GA and PSO.
Using an Invasive Weed Optimization Algorithm (IWOA) for parameter optimization, Hamedi et al. [2] found the five-parameter Muskingum method to be more accurate than the four-parameter version. With PSO for parameter estimation, Moghadam et al. [30] found the four-parameter Muskingum method to be more accurate, than the three-parameter Muskingum method with GA and linear programming.
Using the gravitational search algorithm for parameter optimization, Kang et al. [31] found the four-parameter Muskingum method to be accurate for flood routing. Flood routing for a case study in China based on real code genetic algorithm was considered [3], and the results indicated that the four-parameter Muskingum model based on genetic algorithm decreased the RMSE and MAE of the two- and three-parameter Muskingum models.
Barati et al. [21] applied different kinds of GPA to the four-parameter Muskingum method for flood routing, and found the fixed genetic programming to be more accurate. For flood routing using the Muskingum-Cunge method, Wang et al. [23] found PSO to yield better results than GA.
In 2018, Lee [32] developed and applied an advanced Muskingum flood routing model by considering continuous stream flow, utilizing weighted inflow. Several statistical indicators have been used to evaluate the performance of the suggested model. The proposed model provided acceptable results, compared to those obtained from previous studies. The results showed that the vision corrected algorithm (VCA) had experienced a relatively small error index compared to the GA, PSO, and nonlinear programming models, for different flood case studies.
The literature review showed that the evolutionary algorithm has a high ability for obtaining the parameter values of the Muskingum model, but some algorithms have limitations, such as trapping in local optimums, slow convergence velocity, or insufficient accuracy for the simulation [33]. For example, the GA can trap in the local optimums or the PSO may have an immature solution due to fast convergence [1,27]. Some evolutionary algorithms have many random parameters, and accurate determination of these parameters is difficult. Thus, the improvement of previous algorithms or definition of new algorithms is necessary.

1.2. Problem Statement

Flood routing is nonlinear and multimodal, with noise. The most widely used model for flood routing is the Muskingum model. One of the main components of the Muskingum model procedure, is the inclusion of a particular number of parameters ranging between two and five, which have to be estimated based on the case study characteristics to be able to accurately route the flood. There are different procedures to estimate the Muskingum model parameters. Actually, the better the estimation of the parameters’ value, the better the prediction of the flood routing characteristics achieved [33,34,35,36,37]. Having an optimization algorithm that is able to optimize the value of these parameters, whether two or three, is needed to enhance the ability of the Muskingum model to identify the flood routing characteristics. Thus, a robust algorithm that leads to the global optimal solution is needed. Bats exhibit a mysterious behavior that has long been attractive. They are able to orientate themselves to their surroundings and food acquisition, without depending on their eyesight. Bats consistently emit echolocation signals. Through analyzing the returning echoes in the auditory system, bats can distinguish their environment and find preys. By continuously watching and concentrating on the abilities of bats, scientists have suggested different bat-inspired algorithms (i.e., bat algorithm (BA) and bat intelligence (BI) algorithm) for the solution of optimization problems [33,34]. Bozorg-Hadad et al. [35] used BA to optimize the use of a repository to reduce hydroelectricity energy shortage. This algorithm has been shown to have a faster convergence, than GA and PSO. Ahmadinafar et al. [36] used a hybrid BA to exploit a 10-repository system to increase energy production, which reduced computational time compared to other evolutionary algorithms. Studies have shown that BA is a powerful method but it has weaknesses, such as trapped in the local optimum or premature convergence [37,38]. Thus, it needs to be improved.

1.3. Objective

In the light of the above, the use of Muskingum model showed its success when applied in flood routing prediction, but it has a few limitations. For example, it cannot be used for the complex boundary condition. In addition, one of the weaknesses of the model is the consideration of the lateral flow. In fact, most of the previous researches ignore the lateral flow while using the simple Muskingum model. Therefore, if lateral flow exists and has high volume, the Muskingum model simulates the flood without consideration of lateral flow [23,24,25,26,27,28,29,30,31], and hence, the achieved accuracy is relatively low. In fact, simplification of the Muskingum equation was a reason to omit the effect of lateral flow for flood routing, whilst there is lateral flow in the reach when the flood happens in nature. Although, there are a few number of references considering the lateral flow, some of them are limited to specific case studies with the low flow lateral condition and using Muskingum models with a simple structure [39,40,41,42].
Limited studies showed that the consideration of the lateral flow, with the help of hydraulic models, could simulate the actual situations that close to the real conditions [12]. Osolivan et al. [17] reported the disadvantages of the hydraulic models for flood routing in the floodplain. Their reports showed that the initial and boundary conditions, and the resistance characteristics of main channels in the floodplain, are neccessary for the hydraulic modeling.
The present study develops the new Muskingum model for flood routing considering lateral flow. In addition, the Muskingum Model has been coupled with a new Improved Bat Algorithm (IBA), to optimize the estimation of the three-parameter Muskingum model.
The objective of this study, was to couple the three-parameter Muskingum method (TPMM) with an improved bat algorithm (IBA) for flood routing, with a multi-reach method and consideration of lateral flow. Citing easy trapping in the extremum for bat algorithm, an improved bat method is suggested and used to simulate flood routing. Chaos search tool is defined to enhance the uniformity and ergodicity of the population. Adapting weight is defined to balance the local and global search tools, for the bat algorithm. One of the advantages of the new bat algorithm is related to decreasing the search range, based on dynamic contraction.
The innovation of the study is the use of a new bat algorithm for flood routing lateral flow. Moreover, whilst previous studies usually considered one reach for flood routing without lateral flow, the present study compared the effect of dividing a river into different reaches on flood routing, and the prediction of peak discharge.

2. Materials and Methods

2.1. Flood Routing

The Muskingum method of flood routing is based on the continuity equation and a storage-discharge relation [30]. The present study considers the lateral flow for the flood routing for the case studies based on a ration of inflow rate O l a t e r a l = β I t , while the other studies do not consider the effect of lateral flow for flood regime, and thus it adds one term to the equation of the Muskingum model with three parameters. If β coefficient equals to zero, the lateral flow has not been considered for the flood routing.
d s t   d t = O t ( 1 + β ) I t
S t = K [ X ( 1 + β   ) I t + ( 1 X ) O t ]
where O t is the output flow at time t, I t is the input flow at time t, S t is the storage at time t, d s t d t is the storage time variation at time t, K is the time coefficient of storage, and X is the weighting factor showing the effect of input and output flows on storage.
Equation (2) expresses a linear relation between storage, and input and output flows. However, a non-linear has also been presented as:
S t = K [ X ( 1 + β   ) I t + ( 1 X ) O t ] m
where m is an exponent. Using Equations (1) and (3), one can obtain [21]:
O t = ( 1 1 X   ) ( S t K ) 1 m ( X ( 1 + β ) 1 X ) I t
Δ S t   Δ t = ( 1 1 X ) ( S t K ) 1 m + ( 1 ( 1 + β ) 1 X ) I t
Using S t , Δ S t , the storage later can be expressed as:
S t + 1   = S t + Δ S t
Flood routing can be done using the following steps:
  • Consider initial values for parameters K, X, β , and m and enter them into the optimization algorithm, in the form of initial population.
  • Calculate the storage based on Equation (3), assuming the equality of input and output flow.
  • Calculate the change in storage relative to time, based on Equation (5).
  • Calculate the storage based on t + 1, according to Equation (6).
  • Calculate the output flow at t + 1, based on Equation (4).
  • Repeat steps 2 to 5.

2.2. Optimization of Multi-Reach Muskingum Coefficients

The multi-reach Muskingum method is introduced to enhance the accuracy of the Muskingum method. The river under study was divided into several smaller reaches, and for each reach, routing was done separately. In other words, for each reach, parameters X, K, β, and m were calculated separately, and the output hydrograph was obtained based on the input flood hydrograph and the assumed values of X, K, β, and m for the first reach [3,21,28,31]. This output hydrograph was considered as the input hydrograph for the second reach, and so on. For the second reach, the assumed values of X, K, β, and m were used, and the output hydrograph was calculated. This process was repeated for all the reaches, until the output hydrograph of the last was obtained. By comparing the computed hydrograph with the observed hydrograph, the error was calculated, and to reach the minimum error, Muskingum coefficients were optimized at all reaches. Figure 1 shows the division of the river into several reaches.

2.3. BAT Algorithm

Bat algorithm (BA) is based on bat sound reproduction and sound reflection. The difference in loudness that comes from the surrounding environment, allows the bat to identify the barrier from food. Bats produce very high sound pulses and listen to their return from the objects around. Each pulse remains only for a few milliseconds. BA is based on the following assumptions [30,39,40,41]:
  • All bats have a high ability to receive sound, so that they can detect food after producing loud sounds.
  • Bats fly randomly at a velocity at place yl, capable of producing sound with fmin frequency and λ wavelength. The sound produced by bats also has loudness A 0 .
  • The loudness of sound, of the bats ranges from A 0 to A min .
Each sound produced by the bat has a pulse rate (r), between 0 and 1. The sound frequency speed and position of the bats are updated as:
f l = f min   + ( f max f min ) × β
v l ( t )   = [ y l ( t 1 ) Y ] × f l , t = 1 , 2 , T
y l ( t )   = y l ( t 1 ) + v l ( t ) , t = 1 , 2 , , T
where f l is the frequency of sound of bats, f min is the minimum frequency, f max is the maximum frequency, β is the random coefficient between 0 and 1, v l ( t ) is the velocity of the bat, Y is the best position of the bat, y l ( t ) is the position of the bat, and T is the number of periods evaluated.
The following equation is used for local search in the bat algorithm:
y ( t )   = y ( t 1 ) + ε A ( t ) , t = 1 , 2 , , T
where ε is the random variable between −1 and 1, and A ( t ) is the loudness of sound. Loudness and pulse rates are updated, according to each stage of the iteration. For example, zero sound loudness means that the bat has found its prey and has temporarily stopped the search.
r l t + 1   = r l 0 [ 1 exp ( γ t ) ] A l t + 1 = α A l t + 1
where α and γ are fixed as constant coefficients. For any value of α between 0 and 1, and γ greater than zero, A l t 0 and r l t t l 0 are true. Figure 2 shows the mathematical procedure of the bat algorithm. In addition, it should be noted that the random walk is considered as a parameter for the local search, for the bat algorithm.

2.4. Improved Bat Algorithm (IBA)

The initial arrangement for the initial version of the BA is defined randomly, and it can be a reason for the uneven distribution that causes premature convergence. The chaos is a technique for improving different algorithms, where the basic idea is related to the exchange of members in the range of (−1,1). The logic mapping function is used for modulating the algorithms. Then individuals are inserted into the chaos sequence, so that it should satisfy the chaos variable space. Then, linear transformation is used to return the members to the corresponding position. The convert space is shown based on following mathematical equation:
L i = 2 ( x i a a   ) ( b a ) 1
where x i a : the initial position of the members.
The following equation is used to show the logic mapping function:
L i + 1   = 1 2 × L i 2
Then, the elements and values are returned to the corresponding position by the following linear transformation:
y i 0 = 1 2 ( b a   ) L i + 1 2 ( b + a )
Furthermore, adapting weight is applied to the bat algorithm to have a good balance between global search ability and local search ability.
y l t = w ( t )   · y l t 1 + v l t
The weight is computed based on the following equation:
w ( k )   = w max · ( w max w min ) · ( T max k ) T max
w max : the initial weight, w min : the final weight, and k: the current iteration number.
The final level is related to the application of dynamic contraction, for adjusting the convergence speed:
y min , i   = max { y min , i , x r a n d × ( y max , i y min , i ) } y max , i = min { y max , i , x r a n d × ( y max , i y min , i ) }
where y min , i : the lower bound position, and y max , i ; the upper bound position.
The algorithm functions using the following steps:
  • Adjust the random parameters for the algorithm, such as loudness, pulsation rate, frequency, and other parameters.
  • The individual position is computed using Equations (13)–(15), and then the objective function is computed for each member, and the best solution is considered as Y .
  • The frequency and velocity are updated using Equations (7) and (8), and the position is computed using Equation (17).
  • The randomness value is compared with rl, and if rl is less than the randomness value, the distribution of the best position is acted based on 0.01 times the random disturbance.
  • The local search is considered for this level. If the loudness is less than rand, the loudness should be updated and the pulsation rate should be improved using Equation (12).
  • Compute the objective function and change the range using Equation (16).
  • The convergence criterion is checked and if it is satisfied, the algorithm finishes or else the algorithm goes to step 2.

2.5. Genetic Algorithm (GA)

In GA, the initial version of the population is composed of different solutions [40]. During an iterative process, subsequent populations are generated to improve the objective function. At each stage, some members from the current population are selected to generate individuals or children of the next generation, based on the fact that the likelihood of selecting people with better performance than others is more likely [41]. The selected individuals produce the next population based on two genetic operators, composition and mutation. The following equations can be used for the composition operator [33].
P o p i n e w   = α P o p i o l d + ( 1 α ) P o p j o l d
P o p n e w   j = α P o p j o l d + ( 1 α ) P o p i o l d
where P o p i n e w is the i-th child, P o p i o l d is the i-th parent, P o p j o l d is the j-th parent, P o p n e w j is the j-th child, and α is a coefficient between 0 and 1. Moreover, mutation is based on the following equation:
P o p i , j   n e w = V a r i , j l o w + β ( V a r j , i h i V a r j , i l o w )
where V a r i , j l o w is the lower limit of the i-th gene in the j-th chromosome, V a r j , i h i is the upper limit of the i-th gene in the j-th chromosome, and β is a random coefficient between 0 and 1. In composition, the production of both new individuals is done by changing the gene. The mutation operator is used for the change in chromosomes and transforming their genes to create diversity in the population.

2.6. Particle Swarm Algorithm (PSO)

Initially, the process starts with a particle set. Each particle is considered as a random solution. In the next step, searches are performed sequentially to achieve the optimal answer. The i-th particle is associated with a position in an s-dimension space, where the value of s shows decision-making variables of the problem [42]. The values of s variables, which determine the positions of particles are a possible solution for the optimization problem. Each particle i is completely determined by three vectors. Vector Xi is the current position of the particle, Yi is the best position where the particle is iterated, and the vector of the particle velocity is shown by V i . Then the particle position and particle velocity vector are updated as:
V i i t e r + 1   = w V i i t e r + c 1 r a n d ( Y i i t e r X i i t e r ) + c 2 r a n d ( Y i t e r X i i t e r )
X i i t e r + 1   = X i i t e r + V i i t e r + 1
where V i i t e r + 1 is the new velocity of the particle, the personal learning coefficient c 2 is the global learning coefficient, Y i t e r is the best solution among the solutions, and X i i t e r + 1 is the new position of the particle. Moreover, w is the coefficient of inertia.

Indices of Error Measurement

1. The sum of squared deviations (SSD): SSD index is used as the objective function in the present study. The index calculates the total of squared deviations between observed and real discharges [28,42,43,44,45,46]:
M i n i m i z e ( S S Q ) = t = 1   n ( O b t O s t ) 2
where O o b t is the observed discharge, O s t is the simulated discharge, and n is the number of data.
2. The sum of absolute deviations (SAD): SAD is the total sum of total deviations between observed and predicted discharges [20,30]:
M i n m i z e ( S A D   ) = i = 1 n ( O b t O s t )
3. Error of Peak discharge (EP): EP index measures the difference between predicted and observed discharges [43,44,45,46].
E Q p = | O o b s e r v e d   p e a k O r o u t e d p e a k | O o b s e r v e d p e a k
4. Error of time to peak (ETP): The ETP index measures the difference between predicted and observed time differences of discharge [24,37,38].
E T p = [ T o b s e r v e d   p e a k T r o u t e d p e a k ]
T o b s e r v e d p e a k is the observed discharge, and T r o u t e d p e a k is the time related to the routed discharge.
5. Mean absolute relative error (MARE): The mean of the relative error between observed and predicted discharges:
M A R E = 1 N i = 1   n ( Q t o b s e r v e d Q t r o u t e d ) Q t o b s e r v e d
N is the number of data.
6. Varex Q (Variance index): This indicator shows the proximity of predicted and observed hydrographs with each other.
V a r e x Q = [ 1 i = 1   N ( O t o b s e r v e d O t r o u t e d ) i = 1 N ( O t o b s e r v e d O m e a n o b s e r v e d ) ] × 100
O m e a n o b s e r v e d is the observed average discharge. The closer the coefficient is to one, the more accurate the ability to predict the flood will be.
7. The agreement index (d) based on follow equation, shows the performance of the model well, so that the value of index can change from 0 to 1 [45,46].
d = 1 i = 1   N ( O i o b s e r v e d O i r o u t e d ) 2 i = 1 N ( | O i r o u t e d O ¯ i o b s e r v e d | + | O i o b s e r v e d O ¯ i o b s e r v e d | )
O ¯ i o b s e r v e d : average of observed data.

3. Results and Discussion

This paper considered three case studies for flood routing. Two case studies were considered as bench problems, which have been used by different researchers using many methods for flood routing (Wilson and Karahan floods), and one case study was related to a river in Myanmar that had an important flood. The Wilson flood is considered as an important case study and different researchers tested different algorithms on this case study [2,3,28,29,30], and thus a comprehensive study can be considered for this case study. The Karahan flood is considered as one of the case studies that have been investigated by different researchers as a benchmark problem [28,29,30,31,35].

3.1. Wilson Flood

Wilson Flood [44] is one of the most important benchmark applications, to investigate the performance of the Muskingum model and other hydrologic models. In fact, this flood pattern was generated under experimental conditions. Different mathematical models were examined using this flood data pattern, which received great attention from researchers in examining their models.
The data was extracted from Wilson [44]. This information includes single peak inflow and outflow hydrographs with lateral flow. In addition, the applied algorithms in this research consider the lateral flow, although as shown in Table 4, previous researches have ignored the lateral flow because of low value, as shown in Table 6.
Several methods, such as the Segmented Least Square Method (SLSM), Hook and Jeeves (HJ)method, in combination with the Conjugate Gradient (HJ + CG), HJ method, in combination with Davidson Fletcher-Powell (HJ + DFP), nonlinear least squares (NONLER), Genetic Algorithm (GA), Harmony Search (HS), Particle Swarm Optimization (PSO), and Honey Beaming Optimization (HBMO), have used this flood without consideration of lateral flow. Given that evolutionary algorithms have random parameters, sensitivity analysis was used to determine the exact values of parameters. The evolutionary algorithms have random parameters, where the accurate values of these algorithms are computed based on sensitivity analysis. It means that the variation of objective function is determined versus the variation of parameter values, and when the objective function has the best value for a parameter, the value of this parameter is introduced as the optimal. The SSD was considered as an objective function for the current study. The frequency parameters were used to update the velocity and then the position of bats was computed based on velocity. When the objective function value is minimized, the value of the different parameters is considered at its optimal value.
Table 1, Table 2 and Table 3 show the sensitivity analysis of IBA parameters, BA, GA, and PSO for flood routing in a single reach. The objective function was considered SSD for this study.
The best population associated with IBA was 60, with the lowest SSD. Moreover, the maximum frequency was 5, with the objective function as 5.01. The maximum loudness of sound was 0.6, with a random walk rate of 5. Furthermore, the mutation rate for GA was 0.6, and the recombination rate for GA was 0.7. In addition, the personal and global learning coefficient in PSO was 2, and the inertia coefficient was 0.6. When the inertia coefficient was 0.6, the objective function had the least value (10.82), and thus, the best value for the inertia coefficient was selected to be equal to 0.6 for the PSO. Other parameters can be seen in the Table 1, Table 2 and Table 3.
Table 4 compares IBA with other evolutionary algorithms, for a single reach routing. The SSD value for IBA was 4.123, with SSD being reduced by 89%, 51%, 93%, 69%, 88%, and 97% compared to GA, PSO, PS, HS, HBMO, and SLSM, respectively. In addition, the SAD index was 7.112, reduced by 69%, 22%, 75%, 69%, 81%, and 84% relative to GA, PSO, PS, HS, HBMO, and SLSM, respectively. In addition, EP and MARE showed the superiority of IBA to other methods. VarexQ index for IBA, compared to other methods, showed more consistency with the predicted hydrograph. Furthermore, the performance of IBA was better than BA, so that SSD, SAD, and other error indexes for IBA were less than BA. For example, SSD and SAD for IBA were 4.123 and 7.112, whilst SSD and SAD for BA were 5.123 and 8.114. In addition, focusing on the peak discharge value, it could be depicted that the computed peak discharge based on IBA (85.11) was the most nearest estimated value to the real observed one (85), and in general, showed a worthy match with observed discharge during the whole period. In addition, for the estimation of the peak time which was predicated based on IBA, it was same to the observed one 60 h. For further assessment, Table 5 shows the performance of all the models, and it could be observed that the proposed IBA in the present study outperformed all the other models.

3.2. Multi-Interval Flood Routing (Wilson Flood)

It can be seen from Table 5 that SSD, SAD, and MARE for IBA, BA, GA, and PSO for three reaches were less than for two reaches, and their values for two reaches were less than for the single reach. The best results by IBA were obtained for three reach flood routing. For flood routing using single, two, and three reaches, all evolutionary algorithms predicted the peak discharge accurately, such that the difference with the observed value was 0. For example, SSD for IBA, for one reach was 4.123, whilst it was 3.988 for three reaches. SSD and SAD for IBA for one, two, and three reaches were 1 and 85%, compared to other algorithms, and thus, IBA better performed than other algorithms. The results based on Varex Q showed that the generated hydrograph based on IBA with consideration of more value for Varex Q had a better performance than the other algorithms, and it had a high match with the observed hydrograph. The investigation of Ep showed the value of the index had decreased from 50% to 97% based on IBA, compared to the other methods for two and three intervals. In addition, the value of MARE had decreased from 6% to 56% based on IBA, for the two and three reaches. Furthermore, the agreement index (d) showed a better performance for IBA, achieving a value closer to 1 compared to the other algorithms.
An increase in the number of reaches for evolutionary algorithms increased the accuracy of flood routing. Figure 3 shows the performance of IBA for flood routing with one, two, and three reaches, which is more consistent with observed discharges. Moreover, the computational time showed better performance for the IBA, compared to the other algorithms. For example, the computational time for the IBA based on one reach was 5s, whilst it was 7, 8, and 9s for the BA, PSO, and GA, respectively.

3.3. Karahan Flood

Karahan et al. [22] routed a flood using various algorithms. This study considered the 1960 flood on the River Wye, Uk. River Wye is 69.75 km from Erwood to Belmont with consideration of lateral inflow, which was ignored in the previous studies, as shown in Table 7. In this study, the proposed model was evaluated considering the flood routing based on lateral inflow, which occurred in this Karahan flood. The data was extracted from Karahan [22].
This flood was also used in the present study. The population used for IBA was 50, the maximum frequency was 5, the minimum frequency was zero, and the maximum sound loudness was 95%. In addition, the number of chromosomes for GA was 50, the probability of mutation was 0.6, and the recombination rate was 0.7. Furthermore, the number of particles used in the particle swarm algorithm was 50, the inertia coefficient was 0.7, and the personal and global learning coefficients were 2. Table 6 shows a comparison of algorithms used in flood routing in a single reach. The value of SSD for IBA was 17,120.21, which had reduced by 81%, 87%, 84%, and 10% compared with PSO, GA, HS, BA, respectively. In addition, SAD for IBA was 117, which had reduced by 83%, 88%, 94%, and 12% compared to PSO, GA, HS, and BA, respectively. The other error indices also showed a more favorable performance of IBA, compared to other algorithms. The predicted peak discharge difference with observed discharge was 0.002 cm, which was less than the other algorithms (Table 6). The time difference between predicted and observed discharge peak time for IBA was one hour, whilst this time for other algorithms was 6 hours, so IBA had a better performance. In addition, VarexQ for IBA had a larger value than other methods, which indicated a better performance. The difference of peak discharge based on IBA with observed peak was 69 cms (the closest value to the observed one), whilst it was 135 cms, 134 cms, and 70 cms for HS, PSO, and BA, respectively.
Table 7 compares IBA, BA, GA, and PSO in flood routing (Karahan flood) for single, two, and three reaches. SSD for IBA for three reaches was 16,098.21, which was 6 percent lower than for a single reach with IBA. Moreover, SAD for IBA for three reaches was 102, which had decreased by 12% compared to a single reach. Thus, IBA had an improved performance in routing with three reaches relative to single and two intervals. This was also true for the other algorithms, as shown in Table 5. Figure 4 also shows the superior performance of IBA based on three reaches. SSD using IBA for three reaches was reduced by 46%, 51%, and 5.3% compared to GA, PSO, and BA, respectively. Furthermore, SAD using IBA for three, two, and one reaches was reduced by 23–89%, compared to the other algorithms. As a result, IBA had a better performance than BA, GA, and PSO because the error indexes had the lowest value using the IBA, compared to the other algorithms. Examining the computation time showed that IBA achieved the optimal value of the objective function “under any number of multi-reach interval” faster than the other algorithms. The lowest value of the MARE index was the one associated with the IBA, on the other hand, this value decreased by 12% and 91% when using two and three reaches, respectively. Furthermore, the d index showed the IBA algorithm and simulated hydrographs using IBA had a better match with the observed hydrograph. Many studies did not consider investigate whether parameters generated by calibration could show different floods on the same river reach. In this article, this issue was considered. The optimal parameters based on IBA and three intervals for the flood event in December 1960, were used to obtain the output hydrograph for a flood on the same river on the another flood event in January 1969, based on the inflow hydrograph on January 1969, and it was considered as the first scenario. Then, the output hydrograph was extracted based on application of obtained parameters of Muskingum models for January 1969. It meant that the model was calibrated with the 1960 storm “first scenario”.
The value of SSD for the first scenario was 878 and the peak value was 280 cms, whilst the observed peak value was 285 and there was a small difference between the observed and the simulated peak value (Figure 5). The SSD for the second scenario was 869 and the peak value was 282 cms. Although, the second scenario acted better than the first scenario, the results for the first scenario were acceptable. The results were so similar because of the accurate sensitivity analysis considered for the Muskingum model for the objective function and different parameters, such as in pervious sections. Although the proposed model structure successfully predicted the flood routing, further research could consider different Muskingum models based on more parameters to investigate the performance of IBA. In addition, application of these models helps the hydraulic designers to have the accurate value for the peak discharge and the flood characteristics, to be able to optimally design the target structure.

3.4. Chindwin River

One of the most important Myanmar Rivers is the Chindwin River, shown in Figure 6. This river, over the past twenty years, has experienced many floods. Due to severe flood conditions, the construction of hydraulic structures and living conditions in the downstream areas are difficult and there is high lateral inflow for this basin. Therefore, it is important to predict river flood conditions. The length of the river is 114 km and the surface of the basin is 7485 km2. A historical flood considered in this study occurred in July 2004, and the Mawliak Station recording the flood is shown in Figure 5. Table 8 shows a comparison of different methods for the first flood with flood routing using a single reach. Comparing IBA, BA, GA, and PSO based on routing with a single reach, SSD for IBA relative to GA, BA, and PSO was reduced by 32%, 11%, and 42%, respectively. For routing with two reaches, SSD for IBA was reduced by 33%, 34%, and 12% compared to PSO, GA, and BA, respectively. For three reaches, SSD for IBA was reduced by 44%, 50%, and 37% compared to PSO, GA and BA, respectively. The SAD index for the three routings indicated the superiority of IBA. Comparison of routing with three, two, and one reaches indicated that three reaches improved routing. SSD in relation to routing with two and one reaches had reduced by 39% and 20%, respectively. ETp showed that there was no time lag in forecasting of peak flow by different algorithms. VarexQ showed that three-reach routing was better than two and one reach routings. Additionally, MARE for IBA-based routing had reduced by 64% and 62%, relative to one and two reach routing, respectively. Figure 7 shows that the three reach flood routing was superior to that based on two and single reaches. The results indicated that the IBA has the less computational time compared to the other algorithms, as shown in Table 6. In addition, considering the attained value of the d index, it was obvious that IBA had the best performance compared to the other algorithms.
To enhance further, the used structure of the proposed Muskingum model coupled with IBA showed the highest abilities, for studying the flood plain with back water. Furthermore, such enhancement for the Muskingum model could be successfully applied as an easy and inexpensive method for computing the time and shape, for an overbank flood, when there is back water and inertia influences along a river channel. In fact, with the new nature-inspired optimization algorithms, the traditional Muskingum method could be integrated with the hydrodynamic software packages, such as HEC-RAS model. The development procedure could be carried out by utilizing HEC-RAS software as the model input for hydrographs and geometric characteristic, to estimate the travel times and attenuations peak, whilst the weighted coefficient (X) value could be achieved based on the optimization model and Muskingum model, as shown in Figure 8. As a result, the proposed model in this research showed the potential to be suitable and appropriate for studying and analyzing flood propagation and flood mapping.

4. Conclusions

The present study investigated the potential of utilizing a three-parameter Muskingum model coupled with Improved Bat Algorithm (IBA) to accurately predict flood routing. Three different case studies have been used in this study, to evaluate the performance of the proposed model. These three case studies were the Wilson flood, Karahan flood, and Myanmar River. Seven different performance indices were used to examine and compare the performance of the proposed model over other algorithms. In addition, discretization of the river stream was considered to improve model accuracy. The results showed that IBA outperformed all other algorithms and was able to reduce the SSD by percentages ranging between 20% and 84%, compared with the other algorithms. In addition, the achieved results using the IBA could predict the peak discharge accurately, with a value very close to the observed one. Under the Karahan flood, IBA considerably achieved the minimum level of error indices for a single reach, compared to other algorithms. Finally, IBA in flood routing with three intervals had a better performance than with single and two reaches. The division of the river into different reaches increased the accuracy of flood routing. The performance of IBA for the river in Myanmar, also showed that the simulated hydrograph with three reaches was more accurate. For example, the computational time for the IBA based on three intervals was 2, 4, and 6 s less than, BA, PSO, and GA, respectively. Furthermore, the EP for the IBA was 33%, 96%, and 97% less than BA, PSO, and GA, respectively. As a result, the proposed Muskingum model coupled with the IBA, could be considered as a strong alternative method for predicting flood routing characteristics.

Author Contributions

S.F., H.K., A.E.S, O.K. and V.P.S. initiated the research point of the hydrological problem and supervised the study. N.F. collected the data and performed the modeling. M.E. and M.F.A. handled the writing up of the introduction and methodology. The analysis of the modeling outcomes have been handled by N.S.M, N.F. organized the whole manuscript and managed the paper submission and revision.

Funding

This research was funded by the University of Malaya Research Grant (UMRG) coded RP025A-18SUS University of Malaya, Malaysia.

Acknowledgments

The authors would like to thank so much the data supplier. We also thank all reviewers and the Editor-In-Chief for their insightful comments that have improved the quality of the final manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

IBAImproved Bat Algorithm
GAGenetic Algorithm
PSOParticle Swarm Optimization
PSPattern Search
HSHarmony Search
HBMOHoney Bee Mating Optimization
NMSANelder-Mead Simplex Algorithm
GPAGenetic Programming Algorithm
RMSERoot Mean Square Error
MAEMean Absolute Error
HAHarmony Algorithm
FLAFrog Leaping Algorithm
(IHBA)Improved Honey Bee Algorithm
(IWOA)Invasive Weed Optimization Algorithm
(BA)Bat Algorithm
BIBat Intelligence
TPMMThree-Parameter Muskingum method
SSDSum of Squares
SADSum of Absolute Deviations
EPError of Peak
ETPError of Time to Peak

References

  1. Barati, R.; Badfar, M.; Azizyan, G.; Akbari, G.H. Discussion of “Parameter Estimation of Extended Nonlinear Muskingum Models with the Weed Optimization Algorithm” by Farzan Hamedi, Omid Bozorg-Haddad, Maryam Pazoki, Hamid-Reza Asgari, Mehran Parsa, and hugo a. Loáiciga. J. Irrig. Drain. Eng. 2018, 144, 7017021. [Google Scholar] [CrossRef]
  2. Hamedi, F.; Bozorg-Haddad, O.; Pazoki, M.; Asgari, H.R.; Parsa, M.; Loáiciga, H.A. Parameter estimation of extended nonlinear muskingum models with the weed optimization algorithm. J. Irrig. Drain. Eng. 2016, 142, 4016059. [Google Scholar] [CrossRef]
  3. Zhang, S.; Kang, L.; Conservation, B.Z. Parameter estimation of nonlinear muskingum model with variable exponent using adaptive genetic algorithm. In Environmental Conservation, Clean Water, Air & Soil (CleanWAS); IWA Publishing: London, UK, 2017; pp. 231–234. [Google Scholar]
  4. Bradbrook, K.F.; Lane, S.N.; Waller, S.G.; Bates, P.D. Two dimensional diffusion wave modelling of flood inundation using a simplified channel representation. Int. J. River Basin Manag. 2004, 2, 211–223. [Google Scholar] [CrossRef]
  5. Neal, J.; Villanueva, I.; Wright, N.; Willis, T.; Fewtrell, T.; Bates, P. How much physical complexity is needed to model flood inundation? Hydrol. Process. 2012, 26, 2264–2282. [Google Scholar] [CrossRef]
  6. Hunter, N.M.; Bates, P.D.; Neelz, S.; Pender, G.; Villanueva, I.; Wright, N.G.; Liang, D.; Falconer, R.A.; Lin, B.; Waller, S.; et al. Benchmarking 2D hydraulic models for urban flooding simulations. Proc. Inst. Civ. Eng. Water Manag. 2008, 161, 13–30. [Google Scholar] [CrossRef] [Green Version]
  7. Dottori, F.; Todini, E. Testing a simple 2D hydraulic model in an urban flood experiment. Hydrol. Process. 2013, 27, 1301–1320. [Google Scholar] [CrossRef]
  8. Kim, B.; Sanders, B.F.; Schubert, J.E.; Famiglietti, J.S. Mesh type tradeoffs in 2D hydrodynamic modeling of flooding with a godunov-based flow solver. Adv. Water Resour. 2014, 68, 42–61. [Google Scholar] [CrossRef]
  9. Prestininzi, P. Suitability of the diffusive model for dam break simulation: Application to a cadam experiment. J. Hydrol. 2008, 361, 172–185. [Google Scholar] [CrossRef]
  10. Aricò, C.; Sinagra, M.; Begnudelli, L.; Tucciarelli, T. MAST-2D diffusive model for flood prediction on domains with triangular Delaunay unstructured meshes. Adv. Water Resour. 2011, 34, 1427–1449. [Google Scholar] [CrossRef] [Green Version]
  11. Schubert, J.E.; Sanders, B.F.; Smith, M.J.; Wright, N.G. Unstructured mesh generation and landcover-based resistance for hydrodynamic modeling of urban flooding. Adv. Water Resour. 2008, 31, 1603–1621. [Google Scholar] [CrossRef]
  12. Bates, P.D.; Horritt, M.S.; Fewtrell, T.J. A simple inertial formulation of the shallow water equations for efficient two-dimensional flood inundation modelling. J. Hydrol. 2010, 387, 33–45. [Google Scholar] [CrossRef]
  13. Costabile, P.; Costanzo, C.; Macchione, F. Performances and limitations of the diffusive approximation of the 2-D shallow water equations for flood simulation in urban and rural areas. Appl. Numer. Math. 2017, 116, 141–156. [Google Scholar] [CrossRef]
  14. Fassoni-Andrade, A.C.; Fan, F.M.; Collischonn, W.; Fassoni, A.C.; Paiva, R.C. Comparison of numerical schemes of river flood routing with an inertial approximation of the Saint Venant equations. RBRH 2018, 23, e10. [Google Scholar] [CrossRef]
  15. Singh, V.P.; Aravamuthan, V. Errors of kinematic-wave and diffusion-wave approximations for steady-state overland flows. Catena 1996, 27, 209–227. [Google Scholar] [CrossRef]
  16. Costabile, P.; Macchione, F.; Natale, L.; Petaccia, G. Flood mapping using LIDAR DEM. Limitations of the 1-D modeling highlighted by the 2-D approach. Nat. Hazards 2015, 77, 181–204. [Google Scholar] [CrossRef]
  17. Chatila, J.G. Muskingum Method, EXTRAN and ONE-D for routing unsteady flows in open channels. Can Water Resour. J. 2003, 28, 481–498. [Google Scholar] [CrossRef]
  18. O’Sullivan, J.J.; Ahilan, S.; Bruen, M. A modified Muskingum routing approach for floodplain flows: Theory and practice. J. Hydrol. 2012, 470, 239–254. [Google Scholar] [Green Version]
  19. Yoo, C.; Lee, J.; Lee, M. Parameter Estimation of the Muskingum Channel Flood-Routing Model in Ungauged Channel Reaches. J. Hydrol. Eng. 2017, 22, 5017005. [Google Scholar] [CrossRef]
  20. Geem, Z.W. Parameter Estimation of the Nonlinear Muskingum Model Using Parameter-Setting-Free Harmony Search. J. Hydrol. Eng. 2011, 16, 684–688. [Google Scholar] [CrossRef]
  21. Barati, R. Parameter estimation of nonlinear Muskingum models using Nelder-Mead simplex algorithm. J. Hydrol. Eng. 2011, 16, 946–954. [Google Scholar] [CrossRef]
  22. Karahan, H.; Gurarslan, G.; Geem, Z.W. Parameter estimation of the nonlinear Muskingum flood-routing model using a hybrid harmony search algorithm. J. Hydrol. Eng. 2013, 18, 352–360. [Google Scholar] [CrossRef]
  23. Orouji, H.; Bozorg Haddad, O.; Fallah-Mehdipour, E.; Mariño, M.A. Estimation of Muskingum parameter by meta-heuristic algorithms. Proc. Inst. Civ. Eng. Water Manag. 2013, 166, 315–324. [Google Scholar] [CrossRef]
  24. Easa, S.M. New and improved four-parameter non-linear Muskingum model. Proc. Inst. Civ. Eng. Water Manag. 2014, 167, 288–298. [Google Scholar] [CrossRef]
  25. Talatahari, S.; Sheikholeslami, R.; Farahmand Azar, B.; Daneshpajouh, H. Optimal Parameter Estimation for Muskingum Model Using a CSS-PSO Method. Adv. Mech. Eng. 2013, 5, 480954. [Google Scholar] [CrossRef]
  26. Ouyang, A.; Li, K.; Truong, T.K.; Sallam, A.; Sha, E.H.-M. Hybrid particle swarm optimization for parameter estimation of Muskingum model. Neural Comput. Appl. 2014, 25, 1785–1799. [Google Scholar] [CrossRef]
  27. Geem, Z.W. Issues in optimal parameter estimation for the nonlinear Muskingum flood routing model. Eng. Optim. 2014, 46, 328–339. [Google Scholar] [CrossRef]
  28. Haddad, O.B.; Hamedi, F.; Fallah-Mehdipour, E.; Orouji, H.; Mariño, M.A. Application of a hybrid optimization method in Muskingum parameter estimation. J. Irrig. Drain. Eng. 2015, 141, 4015026. [Google Scholar] [CrossRef]
  29. Niazkar, M.; Afzali, S.H. Assessment of modified honey bee mating optimization for parameter estimation of nonlinear Muskingum models. J. Hydrol. Eng. 2015, 20, 4014055. [Google Scholar] [CrossRef]
  30. Moghaddam, A.; Behmanesh, J.; Farsijani, A. Parameters estimation for the new four-parameter nonlinear muskingum model using the particle swarm optimization. Water Resour. Manag. 2016, 30, 2143–2160. [Google Scholar] [CrossRef]
  31. Kang, L.; Zhang, S. Application of the elitist-mutated PSO and an improved GSA to estimate parameters of linear and nonlinear Muskingum flood routing models. PLoS ONE 2016, 11, e0147338. [Google Scholar] [CrossRef] [PubMed]
  32. Lee, E.; Lee, H.; Kim, J. Development and Application of Advanced Muskingum Flood Routing Model Considering Continuous Flow. Water 2018, 10, 760. [Google Scholar] [CrossRef]
  33. Wang, L.; Lapin, S.; Wu, J.Q.; Elliot, W.J.; Fiedler, F.R. Accuracy of the Muskingum-Cunge method for constant-parameter diffusion-wave channel routing with lateral inflow. arXiv, 2018; arXiv:1802.04429. [Google Scholar]
  34. Jamil, M.; Zepernick, H.J.; Yang, X.S. Multimodal Function Optimization Using an Improved Bat Algorithm in Noise-Free and Noisy Environments. In Nature-Inspired Computing and Optimization; Patnaik, S., Yang, X.-S., Nakamatsu, K., Eds.; Modeling and Optimization in Science and Technologies; Springer International Publishing: Basel, Switzerland, 2017; Volume 10, pp. 29–49. [Google Scholar]
  35. Bozorg-Haddad, O.; Karimirad, I.; Seifollahi-Aghmiuni, S.; Loáiciga, H.A. Development and application of the bat algorithm for optimizing the operation of reservoir systems. J. Water Resour. Plan. Manag. 2015, 141, 4014097. [Google Scholar] [CrossRef]
  36. Ahmadianfar, I.; Adib, A.; Salarijazi, M. Optimizing multireservoir operation: hybrid of bat algorithm and differential evolution. J. Water Resour. Plan. Manag. 2016, 142, 5015010. [Google Scholar] [CrossRef]
  37. Cai, X.; Wang, H.; Cui, Z.; Cai, J.; Xue, Y.; Wang, L. Bat algorithm with triangle-flipping strategy for numerical optimization. Int. J. Mach. Learn. Cybern. 2018, 9, 199–215. [Google Scholar] [CrossRef]
  38. Chakri, A.; Khelif, R.; Benouaret, M.; Yang, X.S. New directional bat algorithm for continuous optimization problems. Expert Syst. Appl. 2017, 69, 159–175. [Google Scholar] [CrossRef] [Green Version]
  39. Ehteram, M.; Allawi, M.; Karami, H.; Mousavi, S. Optimization of Chain-Reservoirs’ Operation with a New Approach in Artificial Intelligence. Water Resour. 2017, 31, 2085–2104. [Google Scholar] [CrossRef]
  40. Allawi, M.F.; Jaafar, O.; Ehteram, M.; Mohamad Hamzah, F.; El-Shafie, A. Synchronizing Artificial Intelligence Models for Operating the Dam and Reservoir System. Water Resour. Manag. 2018, 32, 3373–3389. [Google Scholar] [CrossRef]
  41. Allawi, M.F.; Jaafar, O.; Mohamad Hamzah, F.; Abdullah, S.M.S.; El-shafie, A. Review on applications of artificial intelligence methods for dam and reservoir-hydro-environment models. Environ. Sci. Pollut. Res. 2018, 25, 13446–13469. [Google Scholar] [CrossRef] [PubMed]
  42. Allawi, M.F.; Jaafar, O.; Mohamad Hamzah, F.; Ehteram, M.; Hossain, M.S.; El-Shafie, A. Operating a reservoir system based on the shark machine learning algorithm. Environ. Earth Sci. 2018, 77, 366. [Google Scholar] [CrossRef]
  43. El-shafie, A.; Mukhlisin, M.; Najah, A.A.; Taha, M.R. Performance of artificial neural network and regression techniques for rainfall-runoff prediction. Int. J. Phys. Sci. 2011, 6, 1997–2003. [Google Scholar]
  44. Wilson, E.M. Engineering Hydrology; Palgrave: London, UK, 1990; pp. 1–49. [Google Scholar]
  45. Legates, D.R.; McCabe, G.J. Evaluating the use of “goodness-of-fit” Measures in hydrologic and hydroclimatic model validation. Water Resour. Res. 1999, 35, 233–241. [Google Scholar] [CrossRef] [Green Version]
  46. Ehteram, M.; Binti Othman, F.; Mundher Yaseen, Z.; Abdulmohsin Afan, H.; Falah Allawi, M.; Najah Ahmed, A.; Shahid, S.; Singh, V.P.; El-Shafie, A. Improving the Muskingum flood routing method using a hybrid of particle swarm optimization and bat algorithm. Water 2018, 10, 807. [Google Scholar] [CrossRef]
Figure 1. Discretization of the river stream.
Figure 1. Discretization of the river stream.
Water 10 01130 g001
Figure 2. The flow chart of the bat algorithm.
Figure 2. The flow chart of the bat algorithm.
Water 10 01130 g002
Figure 3. The simulated Hydrographs for Wilson flood using different methods.
Figure 3. The simulated Hydrographs for Wilson flood using different methods.
Water 10 01130 g003
Figure 4. Simulated hydrograph for Karahan flood.
Figure 4. Simulated hydrograph for Karahan flood.
Water 10 01130 g004
Figure 5. The first and second scenario for calibration of model.
Figure 5. The first and second scenario for calibration of model.
Water 10 01130 g005
Figure 6. The location of the River in Myanmar.
Figure 6. The location of the River in Myanmar.
Water 10 01130 g006
Figure 7. Simulated hydrograph for Myanmar River using IBA.
Figure 7. Simulated hydrograph for Myanmar River using IBA.
Water 10 01130 g007
Figure 8. Advanced Muskingum model for flood plain with complex boundary condition.
Figure 8. Advanced Muskingum model for flood plain with complex boundary condition.
Water 10 01130 g008
Table 1. Analysis for the performance of Improved Bat Algorithm (IBA) algorithm (Wilson flood).
Table 1. Analysis for the performance of Improved Bat Algorithm (IBA) algorithm (Wilson flood).
SSD
Objective Function (cms)Random Walk RateObjective Function (cms)Maximum LoudnessObjective Function (cms)Maximum FrequencyObjective Function (cms)Population Size
6.2316.010.26.1216.2320
5.6635.890.45.7835.8940
4.1254.120.64.1254.1260
5.1475.240.805.7675.1580
Table 2. Analysis for the performance of the Genetic Algorithm (GA) (Wilson flood).
Table 2. Analysis for the performance of the Genetic Algorithm (GA) (Wilson flood).
SSD
Objective FunctionCrossover RateObjective Function (cms)Mutation RateObjective Function (cms)Population Size
46.120.1047.120.2045.3920
43.210.3042.240.4038.9440
39.190.5039.240.6039.2360
40.120.7040.230.8040.1280
Table 3. Analysis for the performance of the Particle Swarm Optimization (PSO) (Wilson Flood).
Table 3. Analysis for the performance of the Particle Swarm Optimization (PSO) (Wilson Flood).
SSD
Objective Function (cms)wObjective Function (cms)c2Objective Function (cms)c1Objective Function (cms)Population Size
12.220.211.211.612.111.612.2410
10.900.410.891.811.891.810.4530
10.820.610.802.010.822.010.8050
11.320.811.122.211.242.211.2370
Table 4. Computed error indexes for Wilson flood.
Table 4. Computed error indexes for Wilson flood.
MethodKXmSSDSADEPETPMAREVarexQ
SLSM0.00100.25002.3470143.60046.400.021600.056198.33
HJ + CG0.00690.26851.929149.64025.200.005900.030199.59
HJ + DFP0.07640.26771.898745.64024.800.003500.033199.63
NONLR0.06000.27002.360041.28025.200.008310.025199.60
GA0.10330.28731.828239.23023.800.008200.031199.70
HS0.08330.28731.863036.78023.400.010700.031299.63
PSO0.07550.29813.6818.8209.7710.000500.026199.93
PS0.48910.27141.828162.6529.480.290100.034599.25
HMBO0.63040.33991.853336.24237.4510.700100.028199.69
BA0.03110.29340.82355.1238.1120.000400.031299.96
Present study IBA0.03120.29971.86784.1237.1120.000200.024599.98
Table 5. Flood routing for Wilson flood according to three different intervals.
Table 5. Flood routing for Wilson flood according to three different intervals.
MethodKXm β SSDSADEPMAREVarexQTime (s)d
IBA0.03120.29971.86780.02124.1237.1120.00020.024599.9850.96
BA0.03140.29961.89230.02105.1238.1250.00040.031299.9670.87
PSO0.07550.29813.6810.01998.8209.7710.00050.026199.9380.76
GA0.10330.28731.82820.011139.23023.800.00820.031199.7090.65
IBAK1 = 0.0378
K2 = 0.0345
X1 = 0.2267
X2 = 0.2245
m1 = 1.9435
m2 = 1.8912
β 1 = 0.0124
β 2 = 0.0114
4.0117.0110.00020.023199.9880.95
BAK1 = 0.0368
K2 = 0.0355
X1 = 0.2167
X1 = 0.2457
m1 = 1.9735
m2 = 1.8812
β 1 = 0.0021
β 2 = 0.054
4.0217.1050.00030.024199.97100.89
PSOK1 = 0.0871
K2 = 0.0881
X1 = 0.2676
X2 = 0.2512
m1 = 1.2311
m2 = 1.2212
β 1 = 0.0034
β 2 = 0.045
8.1239.1230.00040.025199.93120.84
GAK1 = 0.0861
K2 = 0.0882
X1 = 0.2214
X2 = 0.2312
m1 = 1.1211
m2 = 1.1112
β 1 = 0.0074
β 2 = 0.054
38.1122.1210.00820.028199.70140.82
IBAK1 = 0.0871
K2 = 0.0851
K3 = 0.0812
X1 = 0.2876
X2 = 0.2745
X3 = 0.2212
m1 = 2.0121
m2 = 2.111
m3 = 2.123
β 1 = 0.0174
β 1 = 0.0162
β 1 = 0.0151
3.9886.9890.00010.022199.98150.90
BAK1 = 0.0841
K2 = 0.0852
K3 = 0.0822
X1 = 0.2976
X2 = 0.2641
X3 = 0.2314
m1 = 2.1122
m2 = 2.221
m3 = 2.2231
β 1 = 0.0139
β 1 = 0.0111
β 1 = 0.0167
4.0016.9990.00020.023199.97160.87
PSOK1 = 0.078
K2 = 0.0812
K3 = 0.0816
X1 = 0.4567
X2 = 0.4569
X3 = 0.4745
m1 = 5.112
m2 = 5.114
m3 = 5.116
β 1 = 0.0094
β 1 = 0.0082
β 1 = 0.0065
7.1268.9890.00030.024199.94190.86
GAK1 = 0.0651
K2 = 0.0612
K3 = 0.0724
X1 = 0.3212
K2 = 0.3414
K3 = 0.3515
m1 = 6.123
m2 = 6.178
m3 = 6.115
β 1 = 0.0056
β 1 = 0.0065
β 1 = 0.0033
37.12321.1230.00720.027199.72220.89
Table 6. Inflow and outflow for Karahan flood.
Table 6. Inflow and outflow for Karahan flood.
Time (h)Inflow (cms)Observed Outflow (cms)HS
[2,25]
GAPSOBAPresent Study IBA
0154102154132102102102
6150140154152.21154137.89137.24
12219169152153.44152.1165.78166.12
18182190181178.11179.4185.43186.11
24182209191190.45190.9209.01207.12
30192218185185.1185.4212.32214.33
36165210187188.21186.9204.45205.24
42150194179179.45180.20191.32192.12
48128172162163.11164.1010.45171.25
54168149141142.11143.70141.44141.38
60260136154151.12152.8132.22133.56
66471228198197.11196.3221.14222.21
72717303264265.21267.3299.12301.12
781092366344349.10351.4387.12385.21
841145456416423.11431.8451.22453.12
90600615599600.12617.4610.34611.21
96365830871872.32881.5826.34827.12
102277969834835.11836.6899.34900.12
108277665689690.11696.2667.24665.21
114187519535534.11549.2522.34520.21
120161444397400.1416.8455.67453.11
126143321283287.10305.10314.32316.11
132126208202203.11221.4212.22210.25
138115176152155.21164.9177.54170.10
144102148124131.10131.20151.23145.11
15093125106108.12110.0127.34119.14
1568811494106.2196.04116.34112.10
162821068888.2389.20107.21105.10
16876978281.2182.7092.1293.43
17473897576.1176.3091.2388.11
18070817373.1073.1082.3480.21
1866776696969.8078.1275.10
1926371666666.772.3469.21
1985966626262.4065.2164
SSD--37,944.1432,944.1431,099.5219,122.2317,120.21
SAD 21621012695134117
EP 0.2780.0780.0900.0680.002
ETP 66611
MARE 0.330.100.090.020.01
VarexQ 83.2984.7898.0598.1299.15
Table 7. Results of four models (i.e., IBA, BA, PSO and GA) for Karahan flood according to three different intervals.
Table 7. Results of four models (i.e., IBA, BA, PSO and GA) for Karahan flood according to three different intervals.
MethodKXm β SSDSADEPETPMAREVarexQTime (s)d
One section
IBA0.6120.4011.6330.014217,120.211170.00210.0199.1560.94
BA0.6160.4221.6540.013619,122.231340.06810.0298.1290.93
PSO0.5860.3651.5450.013831,099.516950.09060.0998.05100.90
GA0.4560.3221.8240.013732,944.1410120.07860.1094.078120.89
Two sections
IBAK1 = 0.672
K2 = 0.524
X1 = 0.382
X2 = 0.375
m1 = 1.723
m2 = 1.645
β 1 = 0.0174
β 2 = 0.0154
17,091.201150.00210.0199.2580.93
BAK1 = 0.652
K2 = 0.521
X1 = 0.352
X2 = 0.355
m1 = 1.623
m2 = 1.642
β 1 = 0.0214
β 2 = 0.0124
17,114.251220.05810.0199.15110.90
PSOK1 = 0.112
K2 = 0.110
X1 = 0.289
X2 = 0.284
m1 = 1.623
m2 = 1.545
β 2 = 0.0064
β 2 = 0.0039
30,235.456870.08860.0898.12140.89
GAK1 = 0.78
K2 = 0.689
X1 = 0.244
X2 = 0.232
m1 = 1.611
m2 = 1.811
β 1 = 0.0027
β 2 = 0.0029
31,231.2310090.06860.1094.79160.88
Three sections
IBAK1 = 0.692
K2 = 0.690
K3 = 0.612
X1 = 0.392
X2 = 0.391
X3 = 0.394
m1 = 1.112
m2 = 1.114
m3 = 1.116
β 1 = 0.0178
β 2 = 0.0162
β 3 = 0.0144
16,098.211020.00210.00899.56100.91
BAK1 = 0.694
K2 = 0.696
K3 = 0.622
X1 = 0.394
X2 = 0.399
X3 = 0.394
m1 = 1.115
m2 = 1.117
m3 = 1.116
β 1 = 0.0074
β 2 = 0.0124
β 3 = 0.0134
16,999.211080.003810.00999.54150.86
PSOK1 = 0.237
K2 = 0.312
K3 = 0.321
X1 = 0.298
X2 = 0.321
X3 = 0.312
m1 = 1.311
m2 = 1.411
m3 = 1.512
β 1 = 0.0044
β 1 = 0.0029
β 3 = 0.0020
30,230.216670.08160.0699.11170.84
GAK1 = 0.900
K2 = 0.878
K3 = 0.815
X1 = 0.296
X2 = 0.294
X3 = 0.224
m1 = 1.655
m2 = 1.652
m3 = 1.651
β 1 = 0.0014
β 2 = 0.0012
β 3 = 0.0010
30,298.119870.05760.0996.12190.82
Table 8. The computed error indexes for Myanmar River.
Table 8. The computed error indexes for Myanmar River.
MethodKXm β SSDSADEPETPMAREVarexQTime (s)d
One section
IBA0.4110.3011.6110.03448.242.120.00200.01599.2270.93
BA0.4100.3041.6120.3219.223.140.00400.01799.1690.92
PSO0.3860.2551.5450.026512.223.240.07800.09598.15100.89
GA0.2560.3121.5240.022214.254.250.08900.10294.78110.87
Two sections
IBAK1 = 0.472
K2 = 0.424
X1 = 0.392
X2 = 0.365
m1 = 1.721
m2 = 1.635
β 1 = 0.0212
β 2 = 0.0139
7.992.100.00200.01499.25100.91
BAK1 = 0.479
K2 = 0.428
X1 = 0.394
X2 = 0.368
m1 = 1.821
m2 = 1.433
β 1 = 0.0210
β 2 = 0.0131
9.112.890.00300.01699.18120.90
PSOK1 = 0.102
K2 = 0.110
X1 = 0.279
X2 = 0.224
m1 = 1.622
m2 = 1.542
β 1 = 0.0197
β 1 = 0.0190
11.953.090.08800.08898.22140.89
GAK1 = 0.78
K2 = 0.789
X1 = 0.244
X2 = 0.212
m1 = 1.610
m2 = 1.611
β 1 = 0.0196
β 2 = 0.0188
12.244.550.06800.10094.89160.87
Three sections
IBAK1 = 0.492
K2 = 0.491
K3 = 0.512
X1 = 0.292
X2 = 0.261
X3 = 0.294
m1 = 1.110
m2 = 1.112
m3 = 1.114
β 1 = 0.0342
β 2 = 0.0288
β 3 = 0.0488
5.121.980.00200.00599.56180.91
BAK1 = 0.412
K2 = 0.471
K3 = 0.514
X1 = 0.291
X2 = 0.254
X3 = 0.292
m1 = 1.112
m2 = 1.114
m3 = 1.116
β 1 = 0.0342
β 2 = 0.0148
β 3 = 0.0138
8.112.100.00500.00399.41200.90
PSOK1 = 0.236
K2 = 0.311
K3 = 0.319
X1 = 0.295
X2 = 0.320
X3 = 0.310
m1 = 0.211
m2 = 0.221
m3 = 0.212
β 1 = 0.0221
β 2 = 0.0212
β 2 = 0.0217
9.272.120.08100.061299.31220.87
GAK1 = 0.910
K2 = 0.876
K3 = 0.815
X1 = 0.396
X2 = 0.396
X3 = 0.324
m1 = 1.655
m2 = 1.652
m3 = 1.651
β 1 = 0.0251
β 2 = 0.0218
β 2 = 0.0217
10.123.250.05700.09096.18240.85

Share and Cite

MDPI and ACS Style

Farzin, S.; Singh, V.P.; Karami, H.; Farahani, N.; Ehteram, M.; Kisi, O.; Allawi, M.F.; Mohd, N.S.; El-Shafie, A. Flood Routing in River Reaches Using a Three-Parameter Muskingum Model Coupled with an Improved Bat Algorithm. Water 2018, 10, 1130. https://doi.org/10.3390/w10091130

AMA Style

Farzin S, Singh VP, Karami H, Farahani N, Ehteram M, Kisi O, Allawi MF, Mohd NS, El-Shafie A. Flood Routing in River Reaches Using a Three-Parameter Muskingum Model Coupled with an Improved Bat Algorithm. Water. 2018; 10(9):1130. https://doi.org/10.3390/w10091130

Chicago/Turabian Style

Farzin, Saeed, Vijay P. Singh, Hojat Karami, Nazanin Farahani, Mohammad Ehteram, Ozgur Kisi, Mohammed Falah Allawi, Nuruol Syuhadaa Mohd, and Ahmed El-Shafie. 2018. "Flood Routing in River Reaches Using a Three-Parameter Muskingum Model Coupled with an Improved Bat Algorithm" Water 10, no. 9: 1130. https://doi.org/10.3390/w10091130

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop