Predicting flood plain inundation for natural channels having no upstream gauged stations

Since the earliest recorded civilizations, humans have tended to settle near streams because of the proximity to water supplies, ecological conditions, biodiversity and advantages for favourable conditions for agricultural activities (Baldassarre et al. 2013). As a result, a considerable percentage of the global population live in the floodplain. This clarifies why millions of people are affected by floods. Individuals exposed to a disaster not only suffer physical injuries but also long lasting psychological damage.

Flood forecasting is one of the fundamental non-structural measures for mitigating the economic and social damage of flooding (Mapiam & Sriwongsitanon 2009). Recent research has seen significant advances in flood forecasting models (Reed et al. 2007; Tayfur & Moramarco 2008; Mapiam & Sriwongsitanon 2009; Tayfur et al. 2009). Wang et al. (2017) presented the annual maximum flood peak discharge forecasting method by using the hermite-PPR model with SSO and LS algorithms. They used the annual maximum flood peak discharge data from 1982 to 2004 in Yichang station. The data set from 1995 to 2004 is used for validating performances of the model whilst the rest is used for training purposes. Meesuk et al. (2014) developed an urban flood modelling approach combining top-view LiDAR data with ground-view SfM (Structure from Motion) observations. An urban area of Kuala Lumpur, Malaysia was chosen as the study area to simulate the extreme flood event that occurred in 2003. The study presented that the technique based on fusion of LiDAR data and Structure from Motion observations can be very beneficial for flood modelling applications. Costabile et al. (2015) generated flood mapping using LIDAR Digital Elevation Model (DEM) with 1D and 2D approaches. The main aim of their study was to determine how the improvements in the topographic description could affect the performance of 1D and 2D models. They also highlighted the critical aspects and the limitations of 1D approach in the hydraulic simulation. For both 1D and 2D approaches, discharge hydrographs, having a 20-min time step, were computed for the main river and its major tributaries to consider the runoff volumes coming from slopes. Timbadiya et al. (2014) carried out flood simulation for the years 2003 and 2006 using MIKE11 hydrodynamic model. Tam (2014) developed a framework by using integration of geospatial technology and hydrodynamic modelling for a local-based flood risk map based on an optical L, hydrological data (water level and discharge), LiDAR DEM, river networks and cross-sections, cadastral data and real estate value information in Kota Tinggi town in Johor state, Malaysia. Wortmann et al. (2014), using a SWIM hydrological model, investigated the influence of Merzbacher Lake outburst floods on discharge in Aksu headwaters in Kyrgyzstan in northwest China. In all these studies, upstream hydrographs were available. However, in most developing countries, due to lack of gauging stations, upstream hydrograph prediction becomes very crucial.

To be able to forecast potential areas which are likely to be affected by a flooding, flow hydrograph prediction becomes a crucial element, especially at upstream sections of a river. Determination of an upstream hydrograph based upon knowledge of downstream hydrograph and hydraulics characteristics of a river channel is known as ‘reverse routing process’ (Das 2009; D'Oria & Tanda 2012; Zucco et al. 2015). Eli et al. (1974) performed reverse flow routing in James River in Virginia, USA using the knowledge of a downstream hydrograph. They employed the implicit four-point finite difference scheme for the solution of the flow equations. Szymkiewicz (1993) presented the formulation of the inverse problem to solve the Saint–Venant equations backwardly in an open channel. Sui (2005) developed an approach that may be used to determine peak flows for an ungauged watershed based on the environmental features of a region and the meteorological data. Das (2009) analyzed the performance of the Muskingum method for the determination of the upstream hydrograph given the downstream hydrograph and proposed a requirement for a fresh calibration of the Muskingum model. Canovas et al. (2011) proposed a methodology to forecast realistic peak discharge for flash flood in an ungauged mountain catchment by using dendrogeomorphic evidence (i.e. scars on trees) and 2D hydraulic model. D'Oria & Tanda (2012) put forth the Bayesian Geostatistical approach, based on the concept of flow hydrograph that can be statistically described because of the nature of continuous function of the time displaying autocorrelation, for reverse flow routing in open channels. Zucco et al. (2015) developed a new methodology to predict flow discharge at an upstream site by using downstream hydrograph and channel characteristics. They applied their model to three river segments along Tiber River in central Italy which is characterized by low slope and smooth geometry. They used the basic continuity equation for flood routing and they obtained optimal values of the parameters and coefficients by using the genetic algorithm (GA) method. Their methodology is quite advantageous since it does not require rainfall time series and topographical data.

In this study, flood risk area in a mountainous region is simulated using an approach which incorporates the reverse flood routing method of Zucco et al. (2015) with the numerical model of HEC-RAS.

Guneysu Basin is situated in Rize in Eastern Black Sea Region in Turkey (Figure 1). It has an area of around 40 km². Generally, Guneysu Creek and its basin is in a mountainous area and it exhibits high topographic relief, which is elevated from 0 to 1,050 m. As such, it has typical characteristics of the Eastern Black Sea topography. The slope varies between approximately 30 and 60%. The stream length is about 10 km and conduit width of the stream is not constant, the widest section of the stream is about 70 m and the narrowest section is about 10 m, and the average stream width is almost 30 m. In the southwest, southeast, and south parts of the basin, the altitudes are high and the heights are decreasing towards the north. The highest point in the basin is 1,050 m and the region shows a typical Eastern Black Sea Climate feature, which is unique in Turkey. Annual average rainfall is around 2,000 mm, winters are cold and rainy while summers are warm and mild.

The return period of severe floods in this study area is generally two years. Long-term heavy rainfall in the region has an effect on the character of floods, as it increases the intensity of the floods in the high slope topographic structure of the region. So far, the largest events recorded had peak discharges occurring between August and October. One of the most destructive flood events in Guneysu happened in 2002, with a peak discharge of 450 m³/s corresponding to a 500-year return period flood which killed 22 people and affected residential, industrial, and agricultural areas. It destroyed the transportation infrastructure in the region. Figure 2 shows the devastation after such flood events that occurred in recent years.

A large-scale topographic map, scale 1:1000, and orthophotos (projected into UTM-WGS84) were used for creation of the DEM of about 20–30 cm spatial resolution for Guneysu Basin. All cross-sections used in HEC-GEO RAS and z values were extracted from the Triangulated Irregular Network surface. Likewise, the stream centreline, riverbanks, flow direction, bridges/culverts, ineffective flow areas, blocked obstructions, blocked positions etc., were created by using the HEC-GEO RAS based on the Digital Terrain Model and the imagery in GIS (Geographical Information Systems) environment. The extraction of land use classes was carried out by the supervised classification. The signatures were generated from the selected training samples. The land use map was produced by performing the maximum likelihood classification. Discharge data (22–26 August 2015) recorded in the downstream gauging station, flood records (between 1973 and 2015), and historical information about damage caused by previous flood disasters were used.

Manning's roughness coefficient (n) is quite crucial to simulate open channel flows (Ding & Wang 2005). Hence, choosing an appropriate value for Manning's n is very important for accuracy of calculated water profiles (HEC-RAS User's Manual 2010). In this study, the roughness coefficient for Guneysu Creek is calculated according to Cowan's method (Cowan 1956), which computes an average roughness value using the formulation of ⁠. Table 1 summarizes values of n_b, n₁, n₂, n₃, n₄, and m for different roughness cases. Using Cowan's formulation and the parameter values in Table 1, the computed roughness values for Guneysu Creek are summarized in Table 2.

The modelling procedure involves three steps: (1) the reverse flood routing model predicts upstream hydrograph at an upstream ungauged station by means of information regarding channel geometric characteristics, downstream flow stage and downstream flow hydrograph; (2) the so-obtained upstream hydrograph becomes an inflow for HEC-RAS model which computes flood wave spreading; (3) the GIS based flood map produced by HEC-RAS is overlapped to the related orthophoto images to delineate the flood-risk areas. The flow chart of the methodology is illustrated in Figure 3.

GA is a nonlinear search and optimization method inspired by biological processes of natural selection and survival of the fittest. It makes relatively few assumptions and does not rely on any mathematical properties of the functions such as differentiability and continuity and this makes it more generally applicable and robust (Liong et al. 1995; Goldberg 1999).

Basic units of GA consist of ‘bit’, ‘gene’, ‘chromosome’ and ‘gene pool’. Gene consisting of bits [0 and 1] represents a model parameter (or a decision variable) to be optimized. The combination of genes forms the chromosome, each of which is a possible solution for each variable. Finally, a set of chromosomes form the gene pool.

The main GA operations basically consist of ‘generation of initial gene pool’, ‘evaluation of fitness for each chromosome’, ‘selection’, ‘cross-over’, and ‘mutation’. The details of GA can be obtained from Goldberg (1989) and Tayfur (2012), among others.

The reverse flood routing procedure can be simply summarized as follows:

• 1.

  Assign initial values for parameters and coefficients:⁠, ⁠, γ, η, a₃, β₃.

• 2.

  Compute upstream hydrograph Q_u(t) by Equation (1).

• 3.

  Compute storage S(t) and S(t − T_L) by Equation (4).

• 4.

  Compute lateral discharge q_l(t) by Equation (5).

• 5.

  Compute downstream hydrograph Q_d(t) by Equation (3).

• 6.

  Compute the MAE for the model produced downstream hydrograph and observed hydrograph by Equation (6).

• 7.

  Update the current values of the parameters and coefficients and go to Step 2. Perform the operations from Step 2 to 5.

• 8.

  Continue minimizing the error, while trying to reach the optimal values of the parameters and coefficients in Step 1, by performing Steps 2–5.

• 9.

  Stop the iterations when the global minimum error is reached and the optimal values for the parameters and coefficients are obtained.

Note that single iteration in GA involves Steps 2–5. At the end of a certain number of iterations (satisfying predetermined tolerance limit), a downstream hydrograph is generated as close as possible to the observed downstream hydrograph. The parameter and coefficients values that result in this best solution are reserved as the optimal values. Then, so-obtained optimal values of the parameters are employed in Equation (1) to produce the upstream hydrograph.

There are many types of hydraulic–hydrodynamic modelling approaches to determine flood inundation accurately. In this study, the hydraulic–hydrodynamic model of HEC-RAS developed by the US Army Corps of Engineers was used. HEC-RAS allows users to perform one dimensional steady and unsteady river flow hydraulics calculations by using the implicit-forward finite difference scheme between consecutive sections of flexible geometry for natural streams (HEC-RAS User's Manual 2010). The HEC-RAS model can simulate unsteady flow in a complex network of open channels. The model has an ability to apply different external and internal boundary conditions. It can also consider off-channel storage and overbank storage areas in the simulations. Although HEC-RAS calculates flow propagation only in a longitudinal direction, the model provides a powerful representation of topography because it is not raster based (Dimitriadis et al. 2016). Furthermore, the model has the ability to simulate flow over several hydraulic structures such as weirs, culverts, gated and uncontrolled spillways, road overtopping, etc. (Castellarin et al. 2009; Sharkey 2014; Dimitriadis et al. 2016; Papaioannou et al. 2016). The EU Directive on floods (2007/60) also emphasizes the necessity to use the most capable and suitable tools for flood modelling, therefore, the HEC-RAS model is preferred because of its efficient results and low computational cost.

Two dimensional models can provide more detailed and accurate simulation results, provided that substantial data (soil maps, topographical maps, land cover maps, land use maps, etc.) related to river basins are available (Werner 2001). This may not be a major problem in developed countries, yet it is an important issue in developing countries, such as Turkey. Hence, in such a case, 1D models become more advantageous. 1D models can provide flood maps in two-dimensions, using less data. Another advantage of 1D models is that they require shorter CPU time while performing the simulations (Cook 2008). Hence, in this study, we employed a 1D HEC-RAS model, which has gained confidence of users over the years and has no cost since it is an open source.

Flood risk is strongly dependent on topography and land features. Therefore, to estimate probable flood damage, cartographic representation of inundation areas is carried out by using the geospatial environment and the remote sensing (aerial imageries). Figure 4 shows the DEM (Figure 4(a)), the orthophoto image (Figure 4(b)), the flood simulation (Figure 4(c)), and the schematic view of the flooded area in Güneysu, Rize (Figure 4(d)). Flooding extent and inundation depth maps (e.g. see Figure 4(d)) are superimposed as a layer over the imageries (e.g. see Figure 4(b)) and the high-resolution DEM in a geospatial environment (Figure 4(a)). Next, they were interpreted with 1, 2 and 3-dimensional visualization, such as 3D in Figure 4(c). The results of the analysis reveal the zones and features (vegetation, land use, roads, etc.) at risk.

A discharge hydrograph in the upstream section of Guneysu Creek is obtained by the reverse flood routing method using gauged station records in the downstream. The distance between the downstream and upstream is about 10 km and the conduit width of the stream varies; the narrowest part is 10 m, the widest part is 70 m and the average width is about 30 m. To generate the upstream hydrograph for the Guneysu Creek by the reverse flood routing procedure, the model employed 30 m channel width, 10 km channel length, downstream stage hydrograph and downstream discharge hydrograph (Figure 5(a)). The model employed GA to determine the optimal values of the coefficients and exponents of the model by generating the downstream hydrograph as closely as possible. The model estimated the downstream hydrograph with mean absolute error, MAE = 3.65 m³/s.

Figure 5(b) shows the reverse flood routing simulation for the August 2015 event in the Guneysu Creek reach. It is the model prediction of the upstream hydrograph which was used as input for HEC-RAS for the flood simulation in the valley.

Figure 6 shows the simulated flood area against the observed one obtained from the General Directorate of State Hydraulic Works (DSI). Figure 6(a) shows the model-produced flooding area while Figure 6(b) gives the report of flooding map by the DSI. As seen, the County centre, Pekmezli, Ortakoy, Kibledag ve Ulucami districts were especially affected by the flood event. At the same time, when inundation areas predicted by the model in Guneysu Creek and the related historical records are overlaid, one can see that approximately 78% of the areas overlapped. Hence, one can say that the model can successfully predict flooded areas with almost 80% success.

The high risk areas predicted by the model are compared in detail against the observed ones in Figures 7–10. As seen in Figure 7, there is a collapsed road due to the flood wave in Kibledag. Figure 8 shows the demolishing of the bank protection along the creek close to centre of the city, where a 550 m long stacked stone fortification and a stone brick were demolished. In addition, 100 m of stone pavement was destroyed on the left bank. Figure 9 shows the collapsed pedestrian bridge located in Ulucami district. Extreme debris was accumulated from the bridge in the County centre and some parts of roads were destroyed (see Figure 9). Figure 10 shows the demolishing of the bank protection along the creek where the Pekmezli Bridge was located. As seen in these figures, the model was able to capture these high risk flooding areas.

This study presented a flood modelling method for Guneysu Creek, having no gauged station in the upstream region. The method first predicts the upstream hydrograph by the reverse flood routing. This hydrograph in turn becomes an input for HEC-RAS model for simulating the flooding area. Results show that the model was successful in determining the flooding areas. The model successfully predicted the high-risk flood zones. This implies that this methodology can be employed for flood analysis or flood mapping in creeks (or river reaches) where there is no upstream gauging stations.