First, a novel nonlinear Muskingum flood routing model with a variable exponent parameter and simultaneously considering the lateral flow along the river reach (named VEP-NLMM-L) was developed in this research. Then, an improved real-coded adaptive genetic algorithm (RAGA) with elite strategy was applied for precise parameter estimation of the proposed model. The problem was formulated as a mathematical optimization procedure to minimize the sum of the squared deviations (SSQ) between the observed and the estimated outflows. Finally, the VEP-NLMM-L was validated on three watersheds with different characteristics (Case 1 to 3). Comparisons of the optimal results for the three case studies by traditional Muskingum models and the VEP-NLMM-L show that the modified Muskingum model can produce the most accurate fit to outflow data. Application results in Case 3 also indicate that the VEP-NLMM-L may be suitable for solving river flood routing problems in both model calibration and prediction stages.

Accurate forecasting of flood wave movement in natural river channels is extremely important for the development of flood protection and alert systems. Generally, there are two categories of approaches for flood routing: hydraulic and hydrologic methods. The former routes flood by numerically solving the famous Saint-Venant equations, which have strict requirements for topographical data of the investigated river reach (such as the channel cross section and roughness) and complicated computations. In contrast, the latter is based on continuity and empirical storage equations and is more popular in engineering applications owing to its simplicity. The Muskingum model, developed by McCarthy (1938), is the most frequently used hydrologic method. In the original linear Muskingum model (LMM), the following continuity equation (Equation (1)) and storage equation (Equation (2)) which involves a storage parameter K and a weighting parameter w are used:
1
2
where is the absolute channel storage at time t; and are observed rates of inflow and outflow at time t, respectively; K is the storage-time constant, a value usually close to the flow travel time through the routing river reach; w is the dimensionless weighting factor that represents the inflow–outflow relative effects on the storage, for stream channels and for reservoir storage (Mohan 1997).
As expressed in Equation (2), the LMM is established on a basic hypothesis that the channel storage within the river reach is some weighted function of the inflow and outflow rates. However, the relationship between them is not always essentially linear in many rivers, thereby making use of LMM may be inappropriate. To deal with this characteristic, an additional exponent parameter m has been introduced to account for the effects of nonlinearity in the following two forms of three-parameter nonlinear Muskingum model (named NLMM1 and NLMM2, respectively) (Gill 1978; Singh & Scarlatos 1987):
3
4
Although the linear and nonlinear Muskingum models (Equations (1)–(4)) have been widely discussed (Tung 1985; Yoon & Padmanabhan 1993; Das 2004; Chen & Yang 2007; Geem 2011; Xu et al. 2012; Niazkar & Afzali 2014), they all neglect the lateral flow along the investigated river reach. However, the lateral flow really exists in many river reaches in actual flood events. Based on the assumptions that the lateral flow varied linearly along the river reach and could be expressed as a ratio of the inflow rate , O'Donnell (1985) first proposed a three-parameter linear Muskingum model which considers the lateral flow (LMM-L) as Equations (5) and (6):
5
6
Apart from the necessity of considering the lateral flow, the nonlinear storage–discharge relationship and its variation (Easa 2013) are also important in Muskingum model building. Thus, a modified nonlinear Muskingum model, which assumes the exponent parameter in the NLMM2 to be varying with the inflow dividing levels and simultaneously takes the lateral flow into consideration, is presented in this paper. In addition, the precise estimation of parameters is the key point in applying the Muskingum model for real-time flood forecasting. During the past decades, various techniques have been applied to solve this problem and can be classified into two types: mathematical methods and heuristic optimization algorithms. The mathematical methods include the least-squares method (LSM) (Gill 1978), the Hooke–Jeeves (HJ) pattern search in conjunction with linear regression (HJ + LR), the conjugate gradient (HJ + CG) or Davidon–Fletcher–Powell (HJ + DFP) methods (Tung 1985), nonlinear least-squares regression (NONLR) (Yoon & Padmanabhan 1993), the nonlinear programming (NLP) algorithm (Kshirsagar et al. 1995), the Broyden–Fletcher–Goldfarb–Shanno (BFGS) technique (Geem 2006), the Nelder–Mead simplex (NMS) algorithm (Barati 2011), and the generalized reduced gradient (GRG) method (Hirpurkar & Ghare 2014; Easa 2015). The above-mentioned methods inevitably have the drawbacks of some special derivation conditions, a suitable guess of initial solution or being time-consuming, and they cannot guarantee convergence to the global optima in most cases (Orouji et al. 2013; Bozorg Haddad et al. 2015). Therefore, an increasing number of researchers have focused in recent years on heuristic algorithms that exhibit strong global search ability and are easy to design, such as genetic algorithms (GAs) (Mohan 1997; Chen & Yang 2007), harmony search (HS) (Kim et al. 2001), particle swarm optimization (PSO) algorithm (Chu & Chang 2009), differential evolution (DE) (Xu et al. 2012), the shuffled frog leaping algorithm (SFLA), the simulated annealing (SA) algorithm (Orouji et al. 2013), and the honey bee mating optimization (HBMO) (Niazkar & Afzali 2014).

The remainder of this paper is organized as follows. First, a modified nonlinear Muskingum model called VEP-NLMM-L is described in the section below. The detailed structure and standard routing procedure of VEP-NLMM-L are also specified in this section. Second, the implementation steps of applying the improved real-coded adaptive genetic algorithm (RAGA) to estimate parameters of the VEP-NLMM-L are given in the following section. The next section defines the five evaluation criteria for assessing the quality of the calculated outflows. Then, three typical cases are employed to demonstrate the validity of VEP-NLMM-L in both model calibration and prediction stages. Finally, some conclusions about our research work are drawn.

Model structure

As can be seen in both traditional LMM, NLMM1, NLMM2, and LMM-L, the parameters are assumed to be constant during the entire flood routing periods. It is worth noting that the exponent parameter m in NLMM2 has no physical meaning but only represents the average nonlinearity in the whole flood routing procedure. For the given complex, unsteady flow conditions, Easa (2013) considered m should be varying with the inflow levels. Hence, after integrating the later flow assumptions by O'Donnell (1985) (Equation (5)) and the variable exponent parameter suggestion by Easa (2013) into the NLMM2 (Equation (4)), a new nonlinear Muskingum model (VEP-NLMM-L) was proposed. The model structure is described below.

Let the variation range of inflow be divided into L levels by (L–1) inflow dividing values , . Define the real variable , where is the maximum inflow. Thus, the lower and upper limits of are 0 and 1, respectively. For inflow level i, , denoting the exponent parameter by . Variation range of are typically between 1 and 3. The constraints for real variable are given by:
7
For a flood hydrograph that goes through time intervals T, the exponent value for flood routing in time interval t, , is determined by:
8
Then, the structure of the VEP-NLMM-L can be written as:
9

Routing procedure

The main steps of using the VEP-NLMM-L to route floods are as follows:

  • Step 1: Assume values of the parameters K, w, , , .

  • Step 2: Calculate the storage amount using Equation (9), where the value of is determined by Equation (8) beforehand.

  • For t = 1, the channel storage is calculated by:
    10
  • Step 3: Calculate the time rate of the storage change using Equation (11):
    11
  • Step 4: For subsequent time intervals, the storage is calculated by Equation (12):
    12
  • where is the rate of change of storage volume at ; is the time interval.

  • Step 5: Calculate the routed outflow at time t using modified Equation (13):
    13
  • where most previous studies suggested using the inflow at the previous time-point rather than at the current time-point in the calculation of the outflow (Geem 2011).

  • Step 6: Repeat Steps 2–5 for all time steps.

The genetic algorithm is a population-based, global search approach inspired by the evolutionary ideas of natural selection and survival of the fittest. GAs have been demonstrated to perform well in handling hydrological modeling problems (Wang 1997; Chang & Chen 1998; Kumar et al. 2012; Yang et al. 2013). A real-coded adaptive genetic algorithm with elite strategy for parameter estimation of the VEP-NLMM-L is described as follows.

Individual coding and fitness assignment

The conventional binary-coded mode for individuals involves two complex steps: encoding each variable as a bit string and a reverse process of decoding (Goldberg 1989), which often results in data redundancy with long string and low search efficiency. Both practical experience and theoretical research show that the real-coded (floating-point) pattern can certainly avoid these problems and perform better in GAs (Chang & Chen 1998; Chang et al. 2005). Herein, each individual is coded as a vector of floating point numbers with the same length as the number of parameters in VEP-NLMM-L, i.e., individual . The total number of decision variables is , the real-value coding of variable is expressed as:
14
where and are specified lower and upper bounds of variable ; rand(0,1) denotes a random number between 0 and 1.
To estimate the optimal parameters for the VEP-NLMM-L, minimizing the sum of the squared deviations (SSQ) between the observed and calculated outflows is taken as the objective function, given by:
15
where SSQ is the sum of the squared deviations between the observed and calculated outflows; is the observed outflow at time t; is the calculated outflow at time t; T is the total number of time intervals.
In minimization problems, the smaller the objective obj(x) is, the higher the fitness is. Thus, for each individual x, its fitness function fit(x) is defined as Equation (16):
16

Adaptive probabilities of crossover and mutation

The tournament selection operator which integrates the idea of ranking is applied for selection operation. To match with the real-coded pattern, the simulated binary crossover (SBX) operator is used for crossover operation (Deb 2000). The polynomial mutation operator is adopted for maintaining diversity in the population during the evolutionary process (Deb et al. 2002).

The crossover probability and mutation probability are vitally important in controlling GAs' performance. Typical ranges of and are 0.5–1.0 and 0.005–0.05, respectively. The higher of , the greater degree of damage for the genetic pattern of parent population inherited by the next new offspring population, which may slow down the speed of GA search at the same time. Too large will lead GA to a purely random search algorithm, while properly small can prevent the premature convergence of GA to local optimal solutions. Conventional GAs usually set the values of and to be constant. Srinivas & Patnaik (1994) first proposed the concept of adaptive probabilities of crossover and mutation that are varied depending on the fitness values of the individuals. However, Yang et al. (2013) consider that the adaptive and are not only depending on the fitness value of the individual, but are also related to the dispersion degree of the population. They can be expressed as follows:
17
18
where is the larger fitness value between the two selected individuals to be crossed; is the average fitness value of the population; is the maximum fitness value of the population; , , , are adaptive parameters.

Constraint handling

In most applications of GAs, the penalty function method has been the most popular approach in handling inequality and equality constraints. In this approach, a penalty term that is proportional to the constraint violation of a solution is added to the objective function obj(x) to form the penalized function F(x), as follows:
19
where J is the number of constraints; the term is zero if and otherwise; is the penalty parameter of the jth constraint; hence, if no constraint violation occurs, individual x is feasible and the overall constraint violation , otherwise, x is infeasible and .

As can be seen in Equation (19), the appropriate setting for values of penalty parameters is the major difficulty in using the penalty function method. This usually requires users to experiment with different values. To overcome this drawback, the following method proposed by Deb (2000) is employed.

The overall constraint violation v(x) is coupled with the tournament selection operator to implement the penalty function method according to the following three criteria when comparing two individuals or solutions:

  1. Any feasible solution (v(x) = 0) is preferred to any infeasible solution (v(x) > 0).

  2. Between two feasible solutions, the one having higher fit(x) is preferred.

  3. Between two infeasible solutions, the one having smaller v(x) is preferred.

Through this method, individuals are never compared in terms of both objective function value and constraint violation information. Thus, penalty parameters are not needed.

The inequality constraints that need to be considered in the parameter optimization problem of VEP-NLMM-L are summarized as follows. On the one hand, as defined in the VEP-NLMM-L, the real variables should satisfy the constraints for :
20
On the other hand, the calculated values of and in the routing procedure are possible to be negative if the decision variables K, w, and of some individuals are infeasible. To avoid this, they should be bounded and penalized by Equations (21) and (22):
21
22
where , are the penalty constants that have very small value (for example, ) to avoid negative and ; and are the next penalized storage and outflow, respectively.

Elite strategy

The elite strategy is the key for guaranteeing GAs converge to the global optimal solution. The idea is to directly preserve the best individuals in the current population to the next generation without crossover and mutation operations and replace the corresponding number of worst individuals of the new generated offspring generation. Here, the top 5% of individuals are preserved in each generation.

Implementation of RAGA on VEP-NLMM-L

The flowchart of applying the RAGA to estimate optimal parameters for the VEP-NLMM-L is shown in Figure 1, in which iter and MaxIter denote the iteration number and the maximum number of iterations, respectively; pop denotes the population size; and are the parent and offspring populations in each iteration, respectively.
Figure 1

Flowchart for parameter estimation of the VEP-NLMM-L using the RAGA.

Figure 1

Flowchart for parameter estimation of the VEP-NLMM-L using the RAGA.

Close modal

To comprehensively assess the accuracy of the calculated (estimated or predicted) outflow hydrographs by different Muskingum models, five relevant evaluation criteria are used in this study. They are: (1) the absolute deviation of peak outflow (DPO); (2) the deviation of peak time (DPOT) (Yoon & Padmanabhan 1993); (3) the Nash–Sutcliffe criterion (Barati 2013); (4) the mean absolute error (MAE); and (5) the mean absolute relative error (MARE) (Niazkar & Afzali 2014). The criteria DPO and the DPOT check the accuracy of magnitude and occurrence time of the calculated outflow peak, respectively. The criterion (%) represents the correctness of the shape and size of the calculated outflow hydrograph, compared with the observed outflow hydrograph. The criteria MAE and MARE measure the degree of how close the calculated outflows are to the observed outflows. These five evaluation criteria are defined as follows.

Evaluation of magnitude and occurrence time of outflow peak

For the calculated and observed peak outflows, the absolute deviation between their magnitude and the deviation between their occurrence time, are defined as Equations (23) and (24), respectively:
23
24
where and are peaks of the calculated and observed outflows, respectively; , are the peak time of the calculated and observed outflows, respectively.

Accuracy of procedure consideration

The Nash–Sutcliffe criterion in percentage is used for evaluating the goodness of fit of the calculated outflow hydrograph, defined as:
25
where = average value of the observed outflows.

Mean absolute error or relative error

Other criteria in mathematical statistics include the mean absolute error (MAE) or the mean absolute relative error (MARE), defined as Equations (26) and (27):
26
27

To reasonably evaluate the performance of the VEP-NLMM-L, the data set from Wilson (1974), a flood in the River Wyre in October 1982 (O'Donnell 1985), and two different flood events that occurred in the same river reach of the River Wye (Natural Environment Research Council (NERC) 1975) were selected as three illustrative case studies (Cases 1–3). Search spaces for the 2L + 2 parameters of the VEP-NLMM-L are specified in Table 1. The parameter represents the contribution of lateral flow to the outflow and is related to the area under the inflow–outflow hydrograph. As is shown in Table 1, the lateral flow is very small in Case 1 and Case 3, but it is large and makes a significant contribution to the outflow in Case 2.

Table 1

Parameter search spaces for each case in calculations

 Kαwmi (i = 1,…, L)ri (i = 1,…, L-1)
Case 1 [0.01, 1.00] [−0.10, 0.10] [0.01, 0.50] (1.00, 3.00) (0.00, 1.00) 
Case 2 [0.01, 6.00] [−3.00, 3.00] 
Case 3 [0.01, 1.00] [−0.10, 0.10] 
 Kαwmi (i = 1,…, L)ri (i = 1,…, L-1)
Case 1 [0.01, 1.00] [−0.10, 0.10] [0.01, 0.50] (1.00, 3.00) (0.00, 1.00) 
Case 2 [0.01, 6.00] [−3.00, 3.00] 
Case 3 [0.01, 1.00] [−0.10, 0.10] 

In all cases, the VEP-NLMM-L was solved over 50 runs using the RAGA to search for the optimal solution and the parameter settings of the RAGA were as follows: population size pop = 200; distribution indexes for the SBX and the polynomial mutation are 1.0 and 100.0, respectively; adaptive controlling parameters , , , ; the maximum number of iterations MaxIter is 1,000.

Case 1: Data set of Wilson (1974) 

The data set from Wilson (1974), which has been demonstrated to have an obvious nonlinear relationship between the storage and the weighted-flow (Yoon & Padmanabhan 1993; Mohan 1997), was selected as the first case. These flood data are a single peak hydrograph and have been previously studied by many researchers (O'Donnell 1985; Yoon & Padmanabhan 1993; Al-Humoud & Esen 2006; Niazkar & Afzali 2014), where and T = 21.

For this case, optimal solutions for the VEP-NLMM-L with different numbers of inflow dividing levels (L = 1–5) are listed in Table 2, and corresponding results by the LMM-L (O'Donnell 1985) and the NLMM2 (Niazkar & Afzali 2014) are also given for comparisons. As is shown in Table 2, the SSQ values by LMM-L and NLMM2 are 815.680 and 36.242, respectively. The value 36.242 is the best existing value for the NLMM2 so far in the literature. The SSQs by estimating the VEP-NLMM-L of different L have obvious reductions, at least (L = 1) 98.71% and 70.91% lower than those obtained by the LMM-L and NLMM2, respectively. Comparing the VEP-NLMM-L with the classical NLMM2, the variation range of the exponent parameter is small and the optimal values of other parameters are close, but the improvement on accuracy of flood routing is excellent. For example, the optimal inflow dividing values are , , , when L = 5.

Table 2

Optimal results by estimating different Muskingum models for the data set from Wilson (1974) 

ModelKwαm/mi (i = 1,…,L)SSQ
LMM-L (O'Donnell 19855.4300 0.1393 0.0235 – 815.680 
NLMM2 (Niazkar & Afzali 20140.6589 0.3399 – 1.8456 36.242 
VEP-NLMM-L (this study) L = 1 0.5362 0.3005 −0.0215 1.8634 10.541 
L = 2 0.4921 0.2978 −0.0195 1.8996, 1.8833 7.475 
L = 3 0.5875 0.3048 −0.0199 1.8583, 1.8446, 1.8361 5.730 
L = 4 0.6176 0.3115 −0.0199 1.8474, 1.8344, 1.8299, 1.8231 4.911 
L = 5 0.5565 0.2938 −0.0192 1.8889, 1.8697, 1.8540, 1.8573, 1.8505 4.535 
ModelKwαm/mi (i = 1,…,L)SSQ
LMM-L (O'Donnell 19855.4300 0.1393 0.0235 – 815.680 
NLMM2 (Niazkar & Afzali 20140.6589 0.3399 – 1.8456 36.242 
VEP-NLMM-L (this study) L = 1 0.5362 0.3005 −0.0215 1.8634 10.541 
L = 2 0.4921 0.2978 −0.0195 1.8996, 1.8833 7.475 
L = 3 0.5875 0.3048 −0.0199 1.8583, 1.8446, 1.8361 5.730 
L = 4 0.6176 0.3115 −0.0199 1.8474, 1.8344, 1.8299, 1.8231 4.911 
L = 5 0.5565 0.2938 −0.0192 1.8889, 1.8697, 1.8540, 1.8573, 1.8505 4.535 

The five prescribed evaluation criteria in the section ‘Evaluation criteria for calculated outflows' are calculated according to the estimated outflows and shown in Table 3. In general, the VEP-NLMM-L produced the most accurate routing outflows in terms of all criteria and the forecasting accuracy increased with the number of inflow dividing levels. All the models except the LMM-L estimate the peak outflow on the correct time interval (DPOT = 0). The VEP-NLMM-L with L = 4 estimates the closest peak outflow to the observed value among all models.

Table 3

Evaluation criteria of the estimated outflows for the data set of Wilson (1974) 

ModelDPODPOTηMAREMAE
LMM-L (O'Donnell 19855.80 −1 93.33 0.137 4.918 
NLMM2 (Niazkar & Afzali 20140.70 99.70 0.028 0.994 
VEP-NLMM-L (this study) L = 1 0.20 99.91 0.017 0.593 
L = 2 0.31 99.94 0.014 0.481 
L = 3 0.10 99.95 0.012 0.390 
L = 4 0.03 99.96 0.012 0.354 
L = 5 0.26 99.96 0.010 0.334 
ModelDPODPOTηMAREMAE
LMM-L (O'Donnell 19855.80 −1 93.33 0.137 4.918 
NLMM2 (Niazkar & Afzali 20140.70 99.70 0.028 0.994 
VEP-NLMM-L (this study) L = 1 0.20 99.91 0.017 0.593 
L = 2 0.31 99.94 0.014 0.481 
L = 3 0.10 99.95 0.012 0.390 
L = 4 0.03 99.96 0.012 0.354 
L = 5 0.26 99.96 0.010 0.334 

The estimated outflow hydrograph by optimizing the VEP-NLMM-L with L = 4 for Wilson's data is shown in Figure 2.
Figure 2

Observed and estimated outflows for Wilson's data.

Figure 2

Observed and estimated outflows for Wilson's data.

Close modal

Case 2: River Wyre October 1982 flood

This case study employed the River Wyre flood event in October 1982 (O'Donnell 1985), which exhibits a considerable increase in the flood volume (lateral flow) between the inflow and outflow sections (some 25 km apart). The flood data have a multi-peaked inflow, and T = 31.

Table 4 lists the optimal solution vectors for the VEP-NLMM-L of different inflow dividing levels (L = 1–5). The SSQs calculated by the VEP-NLMM-L of different L have obvious reductions, at least 88.23% (L = 1) lower than the value obtained by the LMM-L. For VEP-NLMM-L with L = 5, the optimal inflow dividing values are , , , .

Table 4

Optimal parameter vectors by estimating different Muskingum models for River Wyre flood, October 1982

ModelKwαm/mi (i = 1,…,L)SSQ
LMM-L (O'Donnell 19855.5400 0.1370 2.6200 – 468.840 
VEP-NLMM-L (this study) L = 1 5.1485 0.2282 2.5285 1.0000 55.193 
L = 2 4.9473 0.2524 2.5212 1.0120, 1.0059 28.462 
L = 3 4.7858 0.2518 2.5225 1.0188, 1.0102, 1.0134 26.239 
L = 4 4.2020 0.2730 2.5142 1.0135, 1.0498, 1.0458, 1.0398 20.490 
L = 5 3.7303 0.2856 2.5166 1.0477, 1.0841, 1.0772, 1.0722, 1.0658 17.574 
ModelKwαm/mi (i = 1,…,L)SSQ
LMM-L (O'Donnell 19855.5400 0.1370 2.6200 – 468.840 
VEP-NLMM-L (this study) L = 1 5.1485 0.2282 2.5285 1.0000 55.193 
L = 2 4.9473 0.2524 2.5212 1.0120, 1.0059 28.462 
L = 3 4.7858 0.2518 2.5225 1.0188, 1.0102, 1.0134 26.239 
L = 4 4.2020 0.2730 2.5142 1.0135, 1.0498, 1.0458, 1.0398 20.490 
L = 5 3.7303 0.2856 2.5166 1.0477, 1.0841, 1.0772, 1.0722, 1.0658 17.574 

Table 5 shows the evaluation criteria of the estimated outflow hydrographs by the LMM-L and VEP-NLMM-L. Similarly to Case 1, the VEP-NLMM-L still produces more accurate routing outflows in terms of all criteria than the LMM-L and the accuracy increases with the number of inflow dividing levels L. All models can estimate the peak outflow on the correct time interval (DPOT = 0). The VEP-NLMM-L with L = 4 also gives the closest peak outflow to the observed value.

Table 5

Evaluation criteria of estimated outflows for River Wyre flood, October 1982

ModelDPODPOTηMAREMAE
LMM-L (O'Donnell 19851.10 98.39 0.101 3.231 
VEP-NLMM-L (this study) L = 1 1.88 99.81 0.030 0.950 
L = 2 0.29 99.90 0.024 0.665 
L = 3 0.17 99.91 0.025 0.689 
L = 4 0.07 99.93 0.023 0.662 
L = 5 0.16 99.94 0.021 0.624 
ModelDPODPOTηMAREMAE
LMM-L (O'Donnell 19851.10 98.39 0.101 3.231 
VEP-NLMM-L (this study) L = 1 1.88 99.81 0.030 0.950 
L = 2 0.29 99.90 0.024 0.665 
L = 3 0.17 99.91 0.025 0.689 
L = 4 0.07 99.93 0.023 0.662 
L = 5 0.16 99.94 0.021 0.624 

The estimated outflow hydrograph by optimizing the VEP-NLMM-L with L = 4 for the River Wyre flood, October 1982, is shown in Figure 3.
Figure 3

Observed and estimated outflows for River Wyre flood, October 1982.

Figure 3

Observed and estimated outflows for River Wyre flood, October 1982.

Close modal

Case 3: River Wye December 1960 and January 1969 floods

The above two case studies demonstrated how well the VEP-NLMM-L performs in the parameter calibration stage. However, the optimal parameters obtained by the RAGA have not been examined as to whether they are suitable for routing other different floods occurring in the same river reach. Thus, two flood events that occurred in the River Wye in December 1960 (Scenario I) and January 1969 (Scenario II) were used in this case to further validate the applicability and reliability of the VEP-NLMM-L. First, parameters of the VEP-NLMM-L (L = 4) for Scenario I were estimated. Then, the obtained parameters were applied to predict the outflow hydrograph for Scenario II. The optimal parameter vectors by optimizing the LMM and VEP-NLMM-L (L = 4) for Scenario I are listed in Table 6. The corresponding estimated outflows of Scenario I are shown in Figure 4.
Table 6

Optimal parameters by estimating the LMM and VEP-NLMM-L (L = 4) for the 1960 flood

ModelKwαm/mi (i = 1,…,L)
LMM-L 5.0500 0.1730 0.0720 – 
VEP-NLMM-L (L = 4) 0.0981 0.3168 0.0617 1.8190, 1.7742, 1.8131, 1.8169 
ModelKwαm/mi (i = 1,…,L)
LMM-L 5.0500 0.1730 0.0720 – 
VEP-NLMM-L (L = 4) 0.0981 0.3168 0.0617 1.8190, 1.7742, 1.8131, 1.8169 
Figure 4

Observed and estimated outflows for the River Wye, December 1960 flood.

Figure 4

Observed and estimated outflows for the River Wye, December 1960 flood.

Close modal
Predicted outflow hydrographs for Scenario II, by making use of the optimal parameters from estimating Scenario I, are shown in Figure 5. As expected, it can be clearly seen from Figure 5 that the predicted SSQ value (12,006) by VEP-NLMM-L (L = 4) is much lower than the value obtained by LMM-L (SSQ = 26,113). The peak of predicted outflows by the proposed model could be controlled with a 2.3% error (DPO = 7.46 and ) and the peak occurrence time was correct (DPOT = 0), and also performed much better than the LMM-L.
Figure 5

Observed and predicted outflows for the River Wye, January 1969 flood.

Figure 5

Observed and predicted outflows for the River Wye, January 1969 flood.

Close modal

This paper proposed a new multi-parameter nonlinear Muskingum model (VEP-NLMM-L). The VEP-NLMM-L includes a variable exponent parameter accounting for the nonlinearity of flood wave and a flow parameter that represents the lateral flow along the investigated river reach. The variable exponent parameter varies depending on the inflow level. An improved real-coded adaptive genetic algorithm (RAGA) with elite strategy was designed and applied for the parameter calibration of the VEP-NLMM-L. From the comparative results of different Muskingum models in three typical case studies, several conclusions can be drawn as below:

  1. The VEP-NLMM-L was successfully modeled and applied in three watersheds (Cases 1–3) with different characteristics. The lateral flow contribution in both Case 1 and 3, i.e., the data set of Wilson (1974) and the two floods in the River Wye, are minor. However, an obvious nonlinear storage–discharge relationship exists in Case 1. The River Wyre, October 1982 flood (Case 3) has a major lateral flow contribution.

  2. For Case 1 and Case 2, optimal results by estimating the VEP-NLMM-L indicate that the SSQ values have obvious reductions, at least 70.91% and 88.23%, compared to the best SSQs by traditional models (NLMM2 and LMM) reported in the literature, respectively. Although the improvement on accuracy of flood routing is excellent, the variation range of the optimal variable exponent parameter is small and the optimal values of all parameters in the VEP-NLMM-L are not substantially different from the NLMM2 constant exponent parameter. Comparisons of the five evaluation criteria based on the estimated outflows also show that the VEP-NLMM-L can produce the most accurate flood routing hydrograph than the other methods in the model calibration stage and its accuracy will increase with the number of inflow dividing levels.

  3. Furthermore, in Case 3, the River Wye, December 1960 flood was used for estimating model parameters (model calibration stage) and the River Wye, January 1969 flood for model forecasting validation (model prediction stage); the results show that the proposed VEP-NLMM-L has excellent applicability and reliability in solving flood routing problems.

This study was financially supported by the Hubei Support Plan of Science and Technology (No. 2015BCA291), the Wuhan Planning Project of Science and Technology (No. 2014060101010062) and the Natural Science Foundation of China (Grant No. 51509099).

Al-Humoud
J. M.
Esen
I. I.
2006
Approximate methods for the estimation of Muskingum flood routing parameters
.
Water Resour. Manage.
20
(
6
),
979
990
.
Bozorg Haddad
O.
Hamedi
F.
Orouji
H.
Pazoki
M.
Loáiciga
H. A.
2015
A re-parameterized and improved nonlinear Muskingum model for flood routing
.
Water Resour. Manage.
29
(
9
),
3419
3440
.
Chang
F.-J.
Chen
L.
Chang
L.-C.
2005
Optimizing the reservoir operating rule curves by genetic algorithms
.
Hydrol. Process.
19
(
11
),
2277
2289
.
Das
A.
2004
Parameter estimation for Muskingum models
.
J. Irrig. Drain. Eng.
130
(
2
),
140
147
.
Deb
K.
2000
An efficient constraint handling method for genetic algorithms
.
Comput. Meth. Appl. Mech. Eng.
186
(
2–4
),
311
338
.
Gill
M. A.
1978
Flood routing by the Muskingum method
.
J. Hydrol.
36
(
3–4
),
353
363
.
Goldberg
D. E.
1989
Genetic Algorithms in Search, Optimization and Machine Learning
.
Addison-Wesley Longman Publishing Co., Inc.
,
Boston, MA
,
USA
.
Kim
J. H.
Geem
Z. W.
Kim
E. S.
2001
Parameter estimation of the nonlinear Muskingum model using Harmony Search
.
J. Am. Water Resour. Assoc.
37
(
5
),
1131
1138
.
Kshirsagar
M. M.
Rajagopalan
B.
Lal
U.
1995
Optimal parameter estimation for Muskingum routing with ungauged lateral inflow
.
J. Hydrol.
169
(
1–4
),
25
35
.
Kumar
K.
Hari Prasad
K. S.
Arora
M. K.
2012
Estimation of water cloud model vegetation parameters using a genetic algorithm
.
Hydrol. Sci. J.
57
(
4
),
776
789
.
McCarthy
G. T.
1938
The unit hydrograph and flood routing
. In:
Conference of North Atlantic Division
.
US Army Corps of Engineers, Washington, DC
.
Natural Environment Research Council (NERC)
1975
Flood studies report
.
Report
,
Institute of Hydrology
,
Wallingford
,
UK
.
Orouji
H.
Haddad
O. B.
Fallah-Mehdipour
E.
Mariño
M. A.
2013
Estimation of Muskingum parameter by meta-heuristic algorithms
.
Proc. Inst. Civil Engrs Water Manage.
166
(
6
),
315
324
.
Singh
V.
Scarlatos
P.
1987
Analysis of nonlinear Muskingum flood routing
.
J. Hydraul. Eng.
113
(
1
),
61
79
.
Srinivas
M.
Patnaik
L. M.
1994
Adaptive probabilities of crossover and mutation in genetic algorithms
.
IEEE Trans. Syst. Man Cyber.
24
(
4
),
656
667
.
Tung
Y.
1985
River flood routing by nonlinear Muskingum method
.
J. Hydraul. Eng.
111
(
12
),
1447
1460
.
Wang
Q. J.
1997
Using genetic algorithms to optimise model parameters
.
Environ. Modell. Softw.
12
(
1
),
27
34
.
Wilson
E. M.
1974
Engineering Hydrology
.
Macmillan Book Company
,
London
,
UK
.
Yoon
J.
Padmanabhan
G.
1993
Parameter estimation of linear and nonlinear Muskingum models
.
J. Water Resour. Plann. Manage.
119
(
5
),
600
610
.