Effects of uncertainty in determining the parameters of the linear Muskingum method using the particle swarm optimization (PSO) algorithm

Floods are natural disasters, and preventing economic, social, socio-economic and other damage from floods has long been a concern for human beings (Dutta et al. 2010; Vafaei & Harati 2010; Kadam & Sen 2012; Farzin et al. 2018; Fotovatikhah et al. 2018; Vatankhah 2018). Flood routing is a process in which the flood hydrograph is determined at a section of the river path using known or assumed data available in one or more points upstream (Raghunath 1997). Overall, flood routing methods are classified into two groups: hydraulic and hydrologic. If the water flow is routed just as a time-dependent function in a known place, the method is called hydrologic routing (Weinmann & Laurenson 1979) in which the continuity equation and discharge–storage relation is used (assuming the steady flow). On the other hand, if the routing is carried out along the system as a function of time and space, it is called hydraulic routing (Chow 1959) in which the continuity equation and momentum equation are used focusing on unsteady flows (Saint-Venant equation) (Shaw 1994). Mohan (1997) developed a model based on a genetic algorithm (GA) to optimize parameters of nonlinear Muskingum. Chu (2009) developed a new model called Muskingum FIS to predict the outflow associated with different inflow patterns, incorporating fuzzy inference system (FIS) rules and the Muskingum relation. Vatankhah (2014) applied an improved Euler's method, fourth-order Runge–Kutta method and Runge–Kutta–Fehlberg method to solve nonlinear Muskingum, stating that the explicit Euler's method lacks acceptable accuracy in calculation of storage time variation and proposed the fourth-order Runge–Kutta method as the most appropriate solution due to its simplicity and acceptable accuracy. Hirpurkar & Ghare (2014) analyzed three different nonlinear forms of parameter (m) in the nonlinear Muskingum method using Microsoft Excel. Their results confirm higher accuracy of nonlinear Muskingum when using the model developed by Chow (1959) as S = K[XI + (1–X)O]^m. Hamedi et al. (2016) improved the accuracy of the proposed models for the nonlinear Muskingum method; instead of using a constant value for the initial storage volume (S₀), they optimized this parameter as a variable parameter in the weed optimization algorithm (WOA). Niazkar & Afzali (2017) investigated a nonlinear Muskingum method using 14 new models and optimization performed through an MHBMO algorithm and the results suggest that using a three variable–parameter model is more accurate (Yadav & Mathur 2018). The extended VPMM method and SVM and WASVM has been used for flood routing and the results suggest that to predict flood wave movement, the extended VPMM is accurate compared to the SVM and WASVM models. The particle swarm optimization (PSO) algorithm is a population-based evolutionary algorithm and is applicable to civil engineering and water resource optimization issues. These include reservoir operations by Kumar & Reddy (2007), water quality management by Afshar et al. (2011), Lu et al. (2002) and Chau (2005), and water resources management in the basin by Shourian et al. (2008). Chu & Chang (2009) optimized the parameters in the nonlinear Muskingum method using the PSO algorithm. The comparison between this method and previous ones, including harmony search (HS), linear regression (LR) and genetic algorithm (GA), indicates the higher accuracy and speed of the PSO algorithm in estimating the parameters of the nonlinear Muskingum. Moghaddam et al. (2016) proposed a new four-parameter model for the nonlinear Muskingum method that was used for four flood routings. The results indicated that the PSO algorithm optimized the four parameters of the presented model with higher accuracy and a faster convergence rate. Also, to enhance the linearization process and to increase accuracy of the linear Muskingum method, Bazargan & Norouzi (2018) divided the inflow hydrograph into the start, peak and end sections and optimized the parameters (X, K, ⁠) of each section using the PSO algorithm. Norouzi & Bazargan (2020) used the PSO algorithm and two basic floods to optimize the parameters of the linear Muskingum method. Okkan & Kirdemir (2020) used a combination of PSO with the Levenberg–Marquardt (LM) algorithms to optimize the parameters of the Muskingum method. Owing to low computational time, algorithms are very capable of optimizing the parameters of the Muskingum method. Increasing the number of parameters of the Muskingum method will lead to the increased calculation time of the algorithm, while the accuracy of results does not change significantly (Farahani et al. 2018).

In most of the previous studies, the parameters of the linear Muskingum method were obtained using single flood data occurring in a river reach and then, for verification of the estimated parameters, the outflow (downstream) hydrograph of the same flood was calculated using the inflow (upstream) hydrograph of the same flood and the above parameters. However, the main advantage of the linear Muskingum method is that it can use the obtained parameters by using the basic flood (the flood that has already occurred in the study river reach with recorded inflow and outflow hydrographs at the hydrometric stations), and outflow hydrograph of other floods occurring in the same river reach can be obtained, provided that the morphology of the river had not changed. In the present study, using the PSO algorithm, the outflow hydrograph of four floods was first calculated using the parameters obtained from the single basic flood. Then, in the second part of the study, considering the uncertainty of the flood occurring in the desired river reach, instead of using a single basic flood, three basic floods were used simultaneously to calculate the parameters of the linear Muskingum method. In addition, the outflow discharge at the beginning of the flood was also optimized as a percentage of the inflow discharge at the beginning of the flood. The proposed approach in this study increased the accuracy of the linear Muskingum method in estimating the outflow hydrograph.

The Muskingum method was introduced by McCarthy (1938) in studies of floods of the Ohio river in the United States.

In the present study, four floods’ hydrographs occurring upstream (hydrometric station of Mollasani, station no.21-307, 48°53′ E, 31°35′N) and downstream (hydrometric station of Ahwaz, station no. 21-309, 48°40′E, 31°20′N) (distance 60.5 km) of a river reach located in the Karun River, Iran, were used (Figure 1). The inflow and outflow hydrographs of the mentioned four floods are given in Figure 2. In this research we attempted to first apply both inflow and outflow hydrographs of all four floods in the desired river reach in order to calculate the parameters of the linear Muskingum method (X, K and ⁠) and estimate the downstream hydrograph (outflow hydrograph). Second, considering the uncertainty and simultaneous use of the three floods to optimize the parameters of the linear Muskingum method (X, K and ⁠) and taking into account the inequality of the outflow discharge at the beginning of the flood (O₁) with the inflow discharge at the beginning of the flood (I₁), flood in the Karun River was routed using the linear Muskingum method and PSO algorithm so that the calculation error for different floods would be reduced.

It is worth noting that recorded data by the research department of Iran Water Resources Management Co. was used in this study.

The first flood occurred on 2012/02/26 to 2012/03/01 with an inflow discharge range of 221–565 m³/s. The second flood occurred on 2008/11/30 to 2008/12/03 with an inflow discharge range of 105–1,154 m³/s. The third flood occurred on 2012/02/02 to 2012/02/05 with an inflow discharge range of 222–494 m³/s. The fourth flood occurred on 2011/11/21 to 2011/11/25 with an inflow discharge range of 349–651 m³/s.

The basic idea in PSO algorithm is based on the assumption that potential solutions are flown through hyperspace with acceleration towards more optimum solutions. Each particle adjusts its flying according to the experiences of both itself and its companions. During the process; the overall best value attained by all the particles within the group and the coordinates of each element in hyperspace associated with its previous best fitness solution are recorded in the memory (Chau 2007; Kumar & Reddy 2007). Details of the PSO algorithm can be obtained elsewhere (Shi & Eberhart 1998; Clerc & Kennedy 2002; Chau 2007; Gurarslan & Karahan 2011; Karahan 2012).

Various algorithms have been used to optimize the parameters of the Muskingum method and other methods and other issues that need to be optimized. One of the fastest, most acceptable and most widely used is the PSO algorithm, which has been confirmed by previous researchers.

In previous studies, efficiency, high convergence speed and proper accuracy of the PSO algorithm have been examined and approved. Therefore, in the present study, the PSO algorithm was used to optimize the linear Muskingum method (X, K, and ⁠).

In the first part of the present study, using a single basic flood data and equating the outflow discharge at the beginning of the flood (O₁) with the inflow discharge at the beginning of the flood (I₁), the parameters of the linear Muskingum method (X, K, and Δt) were optimized using the PSO algorithm. Then, using the optimized parameters of a single flood, the outflow hydrographs of all four floods were calculated and compared with the observed outflow hydrographs (recorded at the Ahwaz hydrometric station). The optimized parameters for all four floods are listed in Table 1. For example, when the first flood was considered as the basic flood, the parameters of the linear Muskingum method were optimized using the inflow and outflow hydrographs of the first flood and then, using the parameters obtained, the outflow hydrographs of the first, second, third and fourth floods were obtained.

In order to calculate X, K and Δt, if the first flood is considered as the basic (observed) flood and all four floods as the computed floods, then the mean relative error (MRE) of the total flood, the MRE of the flood peak section and the absolute value of the deviations of peak of observed and routed outflows (DPO) value is as described in Table 2. Moreover, in order to calculate X, K and Δt, if the second, third and fourth flood is considered as the basic (observed) flood and all four floods as the computed floods, then the MRE of the total flood, the MRE of the flood peak section and the DPO value is as described in Table 2. It is worth noting that the time of the flood peak section for the first, second, third and fourth flood is 33–48 hours, 32–55 hours, 27–45 hours and 27–41 hours, respectively.

The average values of errors considering all four floods as the basic flood in the case of considering the first, second, third and fourth floods as the computed flood is listed in Table 3.

Moreover, the values of inflow hydrograph, observational outflow hydrograph and computational outflow hydrographs for all four floods, in the case of using a single basic flood, are shown in Figure 4.

As can be seen from Table 2, if the inflow discharge of the basic (observed) flood is closer to the inflow discharge of the computed flood, the above-mentioned error values (the MRE of the total flood, the MRE of the flood peak section and DPO) will decrease. According to Figure 1, the inflow discharge range of the second flood is different from the other three floods. According to Table 3, the average values of MRE of the total flood, the average values of MRE of the flood peak section and the average values of DPO for the second flood were 12.42, 14.19 and 6.16%, respectively, according to the results presented in Table 2. This indicates that the average of errors of the second flood is higher than the other three floods. The average error results are shown in Table 3 for all four floods.

In other words, in order to calculate the downstream (outflow) hydrograph of a new flood occurring in the desired river reach, it is better to use the X, K and parameters obtained from the basic flood whose inflow discharge range would be closer to the inflow discharge range of the new computed flood.

In the second part of the present study, considering the flood uncertainty, instead of using a single basic flood, three basic floods were used simultaneously to calculate the parameters of the linear Muskingum method (X, K, and Δt) using the PSO algorithm. In addition, the outflow discharge at the beginning of the flood (O₁), which was considered equal to the inflow discharge at the beginning of the flood (I₁) in the previous studies, was optimized as a percentage of I₁ of all three floods so that the linear Muskingum method would result in maximum accuracy.

As explained above, the main advantage of the linear Muskingum method is that using the basic flood parameters (a flood that has already occurred in the study river reach with the inflow and outflow hydrographs recorded at the hydrometric stations) any other flood occurring at the same river reach could be obtained, provided that the morphology of the river reach had not changed. On the other hand, according to the results presented in the first part, the calculations would be more accurate if the inflow discharge range of the basic flood would be closer to that of the computational flood. However, in the case of a flood with an inflow discharge range very different from that of all other floods that have already occurred, the basic flood parameters would not yield accurate results and a more appropriate solution should be taken to calculate the outflow hydrograph of the flood.

In other words, considering the uncertainty of the inflow hydrograph of the flood occurring in the study river reach, three basic floods were used simultaneously in the present study to optimize the parameters of the linear Muskingum method.

In addition, when a new flood occurs in the river, the outflow discharge at the start of the flood is unknown and, inevitably, it was considered equal to the inflow discharge at the beginning of the flood in the past studies. However, considering the uncertainty of the outflow discharge at the beginning of the flood (O₁), this value was optimized as a percentage of the inflow discharge at the beginning of the flood (I₁) in the present study.

According to Table 4, for example, if the third flood is considered as the first basic flood, the first flood is considered as the second basic flood and the fourth flood is considered as the third basic flood, then the optimum values of X, K and would be equal to 0.378, 13.229 (h) and 1.009 (h), respectively. In addition, the ratio of O₁ to I₁ for the first, second and third basic floods was optimized as 81.93, 87.65 and 85.46%, respectively, with an average value of 85.02%.

Different modes for considering the four floods as the basic flood, with the above-mentioned conditions, are given in Table 4.

According to the calculations listed in Table 4, the values of X, K and were considered as 0.378, 13.229 (h) and 1.002 (h) respectively to calculate the errors listed in Table 5 (the first mode of the three basic floods) and in the case in which the third, first, and fourth flood is considered as the computed flood, the computed was obtained as 81.93, 87.65 and 85.46%, respectively. However, the second flood in the first order of the floods (according to Table 4) has not been used as the basic flood to calculate the X, K and parameters and ⁠. Therefore, in order to calculate the outflow discharge at the beginning of the second flood (O₁) the average value of the three other floods (listed in the last column of Table 4), which is equal to 85.02%, has been used. The calculations of the second, third, and fourth flood modes are similar to the above mentioned process and are shown in Table 5.

The average values of errors in the case where the first to fourth floods are considered as the computational flood are listed in Table 6. As can be seen in this table, if the second flood is considered as the computational flood, the error will be maximized. It is due to a wider range of inflow discharge of the second flood in comparison to the other three floods.

Moreover, values of inflow hydrograph, observational outflow hydrograph and computational outflow hydrograph for all four floods, in the case of using each mode of three basic floods at the same time, and optimization of O₁ as a percentage of I₁, is provided in Figure 6.

By comparing the error values listed in Table 6 with the error values in Table 3, it was found that the errors in the case of using the proposed method in the second part of the study were smaller than the corresponding values in the first part of the study. In other words, considering the uncertainty, when using the three basic floods simultaneously to calculate the X, K and parameters and considering the outflow discharge at the beginning of the flood as a fraction of the inflow discharge at the same time (⁠ %), the calculation accuracy of all the four floods has increased.

The present study is organized in two parts. In the first part, the parameters of the linear Muskingum method (X, K and ⁠) are optimized using the PSO algorithm as well as data of a single basic flood, and then the outflow hydrograph of all four floods is calculated by applying those parameters. In the second part, the PSO algorithm is applied to optimize the said parameters, considering the uncertainty of inflow discharge variations of floods, the data obtained from three basic floods are used simultaneously to optimize X, K and _, rather than data of a single flood. In addition, outflow discharge at the beginning of flood (O₁) is optimized as a percentage of inflow discharge at the beginning of flood (I₁), and then the outflow hydrograph for all four floods is calculated using those optimized values.

Overall, the results of the present study include the following:

• 1.

  When a single basic flood is used to optimize the parameters of the linear Muskingum method, the average values of mean relative error (MRE) of total flood, the average values of MRE of peak section of flood and average value of DPO are calculated as follows: for the first flood as 3.47, 2.02 and 2.25%; for the second flood as 12.42, 14.19 and 6.16%; for the third flood as 4.93, 5.84 and 4.15% and for the fourth flood as 3.62, 2.28 and 2.34%, respectively.

As there are significant differences between inflow discharge variations of the second flood (105–1,154 m³/s) with inflow discharge variations in the first (221–565 m³/s), third (222–494 m³/s) and fourth (349–651 m³/s) floods, the computations of the outflow hydrograph of the second flood are made with lower accuracy. In other words, the closer the range of inflow discharge variation of basic flood is to the computational flood, the accuracy of computations in estimating the outflow hydrograph will be increased.

• 2.

  When the proposed approach in this study (simultaneous use of three basic floods in order to optimize the parameters of the linear Muskingum method and optimization of outflow discharge at the beginning of flood) is used, the above mentioned values are estimated as follows: for the first flood as 2.41, 0.99 and 2.13%; for the second flood as 10.83, 11.49 and 3.69%, for the third flood as 3.03, 5.54 and 3.84%, and for the fourth flood as 2.45, 2.08 and 1.72%, respectively. As such, considering the uncertainty of inflow discharge variations of the flood, the proposed approach will increase the accuracy of the linear Muskingum method in estimating outflow discharge of all four floods.

The authors gratefully acknowledge the research department of Iran Water Resources Management Co. for the support in collecting and providing the data required in this project.

All relevant data are included in the paper or its Supplementary Information.