Evaluation of 1D hydraulic models for the simulation of mountain fluvial floods: a case study of the Santa Bárbara River in Ecuador

River flood events are the most costly and frequent natural hazard (Timbe & Willems 2011). Although floods are a natural part of the hydrological cycle, they cause damage to the environment, and can lead to human and financial casualties (Vojtek & Vojteková 2016). They affect not only the local population but also infrastructure located on floodplains (Stoica & Iancu 2011; Jongman et al. 2012; Turner et al. 2013). Among the tools used to prevent the impact of this phenomena, a flood risk assessment – commonly based on flood modeling – is the preferred approach to obtain flood maps, which in turn permit the delineation of vulnerable areas. Thereby, development of more accurate and efficient flood modeling techniques are necessary. On the other hand, since coastal areas are more prone to inundations than highlands, flood modeling studies have been mostly performed in areas below 500 m above sea level (asl) (e.g., Casas et al. 2006; Timbe & Willems 2011). This trend has led to there being scarce information about the behavior of mountain rivers when flooding, as although they have steep slopes (typical slopes from 0.1 to 5%), in their course they pass by flat areas such as valleys where settlements and infrastructure exist (Timbe & Timbe 2012). Mountains rivers are characterised by their vast spatial complexity and dynamics, and generally are poorly equipped with hydrologic measuring networks (Weingartner et al. 2003). Therefore, frequently there is a lack of flood information for these regions.

Flood modeling can be performed using numerical hydraulic modeling (Alho & Aaltonen 2008; Timbe & Timbe 2012; Paredes et al. 2014), mostly based on solving the hydrodynamic equations for one (1D models) or two (2D models) spatial dimensions (Bladé et al. 2014). Since 1D models require less input data than 2D models (e.g., cross sections or roughness coefficient), they are easier to set up and calibrate (Vojinovic et al. 2011). It is generally acknowledged that 1D model results are robust when the main flow occurs along the direction of the main river channel (e.g., Timbe & Timbe 2012); however, for the case of over bank flows, which are generally prone to multiple flow directions, those results should be taken with care (Bates et al. 1997). In contrast, although the increase of computing power has multiplied the use of 2D (Vojinovic et al. 2011; Bladé et al. 2014), their use are still limited to smaller spatial scale studies (Di Baldassarre et al. 2010; Timbe & Willems 2011; Ali et al. 2015). As compared to 1D models, 2D models permit a more accurate simulation of flow processes. However, for practical applications where large flooded areas are to be involved, 2D or even 3D models do not necessarily perform better (Di Baldassarre et al. 2010; Timbe & Willems 2011). For flood modeling purposes, Di Baldessarre (2012) presents a summary of numerical tools and their potential application based on their dimensionality. For example, 1D and 2D models can be used for design scale modeling, the former of the order of the tens to hundreds of km, and the latter of the order of tens of km. On the other hand, 3D models only can be applied at the local scale (e.g., predictions of velocity in the main channel and floodplains) or specific scale (e.g., simulation of flood-induced bridge scour), due to its high computational requirement. Besides, models are only good as the quality of the input data for their parameterization, calibration and validation (Paredes et al. 2014). Consequently, although more complex approaches exist, when dealing with simulation of flood events, the use of 1D models is still common. The latter is especially true for cases where steady state conditions of flows are considered (Vojtek & Vojteková 2016).

The application of hydraulic models for simulation of floods has considerably grown in recent years. This in part is due to the fact that data acquisition methods for Digital Elevation Models (DEMs) have passed from conventional topographic surveys to sophisticated remote sensing techniques (Marks & Bates 2000). This shift has resulted in easier accessibility to high spatial resolutions, e.g., DEM of riverbanks and floodplains obtained through LiDAR (Light Detection and Ranging) techniques (Zerger & Wealands 2004). For hydraulic modeling, this denser and higher resolution data brought to light differences in the performance of hydraulic models (Ali et al. 2015). Uncertainties of hydraulic model results, linked to topographic spatial resolution data, have been evaluated in multitude of studies (e.g., Casas et al. 2006; Ali et al. 2015). In any case, higher spatial resolution data is preferable, however, its availability has been commonly restricted by economic factors or geographical conditions (Casas et al. 2006).

The main purpose of this study is to compare the performance of estimations of one public domain model (HEC-RAS) and other two commercial models (Flood Modeller and MIKE 11). These models will be applied to a case study of an Ecuadorian mountain river, using 1D modeling approach. Main goals of the present study are to quantify: (1) the sensitivity of water surface elevation (WSE) predictions, depending on the selected model; (2) the relevance of input data resolution, particularly cross sections, on the model results; and (3) show the variability of the water levels results using information of a validated model without carrying out a calibration process. The last objective was due to the lack of flood information records in most of the mountain regions, especially for the Ecuadorian case. The aim is to compare the model performance for fast flood assessment situations. In turn, it will help to improve flood management strategies developed by decision makers.

This research was conducted in the Santa Bárbara mountain river, with a basin area of 953 km². The river passes through the city of Gualaceo (Figure 1), in the province of Azuay, in Southern Ecuador. The river reach under study is 5 km, with an average slope of 0.25%. The average elevation of the study area is 2,330 m asl, mean annual temperature is 17.6 °C, and the annual rainfall is around 960 mm (INAMHI 2015).

For the same river, a 1D hydraulic model using HEC-RAS version 4.1 was implemented for a 10 km river reach (SENAGUA 2014), which used a high resolution cross section survey every 25 m along the river. The model was validated based on historical inundation areas for flood events ranging from 2 to 10 years return periods (exceedance probabilities 50% and 10%), because there is no available time series measurements of water levels and discharges. This study uses the same surveyed information but focuses only on the most critical zone, i.e., the most prone to inundation: the 5 km river reach. The HEC-RAS results for different return periods implemented with detailed cross sections (field survey complemented with information extracted from DEM) is used as a reference model to compare with the other two models. Due to the lack of flood information records in the study area, we use the reference model results and data to compare the performance of the other modeling packages.

In the former equations, Q = discharge (m³·s⁻¹), ω = vertical velocity distribution coefficient, A = cross section area (m²), g = gravitational acceleration (m·s⁻²), x = river distance in the downstream direction (m), h = water height above datum (m), t = time (s), n = the Manning coefficient (s·m^−1/3), q = lateral inflow (m²·s⁻¹), and R = the hydraulic radius (m).

Equations (2) and (3) are based on the following assumptions: (i) the main flow occurs in just one direction (longitudinal velocity) after averaging all the hydraulic parameters in the cross-section area; (ii) hydrostatic pressure prevails and vertical accelerations are negligible; (iii) the longitudinal slope of the channel bed is small; (iv) the water is incompressible and homogeneous; (v) Manning's equation can be used to describe resistance effects; and (vi) the streamline curvature is small (gradually-varying flow). Although models are based on the same equations, the main difference between them is related to the numerical method used in each model for solving the differential equations of flow.

Besides, each model has its own particular features in terms of user interface. For instance, for model implementation, MIKE 11 uses separate software modules; in contrast, Flood Modeller has one main window where all tools and options for the modeling process are available. On the other hand, HEC-RAS has independent software modules, but they can be accessed through one main interface. In this regard, as compared with MIKE 11, Flood Modeller and HEC-RAS have friendlier user interfaces and more detailed user and technical manuals. Another difference between the evaluated models is that while for HEC-RAS and Flood Modeller most of hydrodynamic parameters (e.g., solver, iteration number, flow tolerance factor, critical depth, relaxation coefficient) are pre-configured, MIKE 11 mandatorily requires the setting of several of these parameters, such as the wave approximation method (solver), the maximum number of iterations, and the computational scheme coefficients (e.g., relaxation).

Data for detailed cross sections (XS) of the main river channel (riverbed and banks) and part of the floodplains were obtained with a topographic survey (186 XS and three bridges) using a total station (SOKKIA SET 650X, precision [5 ± 2 ppm x D]mm). To complete the floodplain areas, coarser data were extracted from a DEM (3 × 3 m spatial resolution), obtained from a LiDAR survey performed by SIGTIERRAS in the year 2012 (http://www.sigtierras.gob.ec) during the drought period. Therefore, the riverbed and banks have suitable representation for that spatial resolution. For each modeling software, two scenarios were evaluated: (a) XS with high spatial resolution for the main river channel, plus extensions of these XS for floodplain areas with data extracted from the available DEM (hereafter so-called type I data); and (b) cross sections entirely extracted from a DEM (i.e., main river channel and floodplain areas, from now on named type II data) (see Figure 2).

Model implementation was based on the input data used by SENAGUA (2014). For our case, the three software packages were fed with the same geometric information (according to scenarios described in the previous section), roughness coefficients and boundary conditions (BC). Model implementation was straightforward, although the main disadvantage when using Flood Modeller and MIKE 11 was the time-consuming task of inputting the geometry data. Five return periods of discharges, as estimated through a flood frequency analysis (SENAGUA 2014), were selected to perform steady flow simulations (Table 2). HEC-RAS and Flood Modeller allow the analysis of steady and unsteady flow, while MIKE 11 only works with unsteady flow. For the latter case the solution was to use a constant discharge hydrograph for a period of time long enough in order to have a constant flow along the entire reach. For the steady flow computations in the gradually-varying flow regime (assumption iv), the upstream boundary condition is a constant discharge, which in time plays the role of the Bernoulli equation parameter. Contrasting, for unsteady flow, each discharge hydrograph was set up as an upstream boundary condition. For HEC-RAS and Flood Modeller, normal flow depth (slope) was set up as the downstream boundary. The normal flow depth (slope) value is different because its computation is based on the two bathymetric information sources (type I and type II). The BC types in MIKE 11 are: inflow, Q-h relation and water level. This package does not allow explicitly for the use of normal depth as downstream BC. Therefore, the water level option was chosen. For this purpose, the water level results from the reference model (HEC-RAS) were used. The objective of this process was to have the same BC in the models as the reference. Also, we want to determine if there is a significant impact on model results by using different downstream BC. Nonetheless, if normal flow depth condition is required at downstream boundary in MIKE 11, a Q-h curve must be constructed. This curve can be developed using the Manning equation based on the bed slope of the reach at the outlet. This Q-h curve is equivalent to normal depth condition in HEC-RAS.

Some highlights in BC were: HEC-RAS and MIKE 11 allow the user to increase the total discharges (tributary discharges) along the river; these values are set in the steady flow by directly selecting a cross section, while when using Flood Modeller, a tool called ‘junction’ is needed to include point discharges (tributaries) in the river network.

Usually, river hydraulic modeling is performed through a calibration process, adjusting the roughness values to fit the model results with observations. This study set out to investigate the sensitivity of the modeling packages on water levels. The analysis and comparison was based on a previous validated model developed in HEC-RAS, therefore every package was implemented with the same input data: (i) cross sections, (ii) calibrated roughness Manning's coefficients, Table 3 shows the calibrated values of reference model, and (iii) return period discharges. No calibration process of roughness values was carried out, because the specific goal of the authors was to know to what extent the numerical scheme used to solve the Bernoulli and Saint-Venant equations can impact the computation of water levels.

Figure 3 shows the comparison between the WSE of the reference model (HEC-RAS type I) versus those of the other five model scenarios (HEC-RAS type II, MIKE 11 type I, MIKE 11 type II, Flood Modeller type I, Flood Modeller type II), and for two return periods (5 and 20 yr). The water surface profiles were obtained by simulating flow discharges for return periods under steady flow conditions. Despite the differences between the solution scheme methods used by the models, the simulated WSEs for type I scenarios were all similar to the reference model. Among results, taking the 20 years event as a sample, Flood Modeller reached the closest fit to HEC-RAS type I, with an average difference of 0.30 m; while MIKE 11 showed the highest difference (0.56 m). On the other hand, except for the HEC-RAS type II model scenario, which showed a good match for all discharge return periods, when compared to the reference model, the remaining type II model scenarios presented bigger differences among them and also when compared to the reference. In general, type II models underestimate the water levels compared with the type I models. Except for HEC-RAS, for which the water levels of both type I and type II models are similar. Specifically, WSE between the type I and type II for each model show approximately a constant variance, on the average 0.26 m (HEC-RAS), 0.68 m (Flood Modeller) and 0.54 m (MIKE 11) for discharges with return periods ranging from 10 to 50 yr. Moreover, for lower events, 2 and 5 yr, the difference is bigger because the flow occurs mainly on the river channel, where there is the most significant gap among type I and type II data. Differences up to 1.3 m between type I and type II model results were found at the downstream boundary (Figure 3(a) and 3(b)). These findings suggest that the variances in water levels between type I and type II models could be explained mainly due to the differences in geometry data (see Figure 2). Considering the lowest point in each cross section (river thalweg), an average deviation of 0.62 m and a maximum difference of 2.46 m are found between the two data sources.

Table 4 presents the comparison of WSEs results with the reference model (type I) through evaluation of statistical metrics. For this, the reference model has been taken as the ‘observed data’. In terms of elevations, the MAE shows that type I models have differences between 0.23 and 0.42 m for a return period of 2 yr. In contrast, type II model scenarios were between 0.36 and 0.84 m for the same return period. In general, the variance increases as the return period increases (e.g., type I models are between 0.30 and 0.56 m, and type II models between 0.26 and 1.08 m). Although for a 20 yr discharge return period, the HEC-RAS type II model has a better performance than type I models. Models that include surveyed cross sections (type I models) have, in general, a better match to the reference model.

Table 4 shows that although all type I models fit very well to the reference model in terms of EF and dr, when compared to other type II derived models, only HEC-RAS type II performed as good as type I model scenarios. This interesting finding shows that, for our study case, HEC-RAS model performance was similar regardless of the geo-spatial information sources (type I and type II). On the other hand, simulated results from Flood Modeller type I yielded EF statistics of 0.99 for return periods of 2 and 5 yr, and 0.98 for return periods of 10 and 20 yr. These results show the higher accuracy when compared to MIKE 11 type II, for which we obtained lower values (e.g., 0.84 for a return period of 10 yr). The same interpretation can be drawn if dr is considered.

In hydraulic modeling, bridges are considered important benchmarks. During flood events, bridges generally reduce the width of the free surface of the water, increasing the water level upstream of the bridge (e.g., backwater effect). According to Figures 1 and 3, for our study, this effect is more evident for bridge number 1, where there was a noticeable and significant difference between upstream and downstream water levels. Besides, for this effect, HEC-RAS is more sensitive. The water levels upstream and downstream of bridge 1 are shown in more detail in Figure 4. If we compare the WSE deviations of the models with the reference, the results upstream and downstream of the bridge are similar for type I models; with a maximum average deviation of 0.41 m (upstream) and 0.51 m (downstream) for MIKE 11. The same behavior is observed for type II models at the downstream position, with a maximum average deviation of 0.45 m for MIKE 11. In contrast, at the upstream position, the maximum average deviation of type II models is twice that of type I models (0.87 m for Flood Modeller). In summary, Figure 4 reveals that the numerical models are stable, with a clear and constant trend for the five events analyzed. Only Flood Modeller is more sensitive at the upstream position, as the deviation increases with the magnitude of the flood event.

As representative examples, for the evaluated return periods, water levels at four cross sections are shown in Figure 5. It can be noticed that, excluding cross section 315.33, for which MIKE 11 type I depicts considerably underestimations for return periods of 10, 20, 50 yr (Figure 5(a)), results of type I models are in the same range than those of the reference model (average deviations ranging from 0.22 m to 0.38 m for Flood Modeller, and from 0.32 m to 0.59 m for MIKE 11). In general, when comparing the WSE between the type I models (Figure 5(a), 5(c), 5(e) and 5(g)) and type II models (Figure 5(b), 5(d), 5(f) and 5(h)) a bigger difference is noticeable; with average deviations between 0.50–1.36 m for Flood Modeller and 0.37–1.31 m for MIKE 11. Among them, Flood Modeller and MIKE 11 type II had the largest difference while HEC-RAS type II the closest match (average differences of 0.17–0.29 m).

We compared three 1D hydraulic modeling software packages (HEC-RAS, MIKE 11, Flood Modeller) for simulation of flood events of the Santa Bárbara River. In terms of flood modeling of mountain rivers framework, and as guidance for implementations of more appropriate 1D hydraulic models, the present model intercomparison provided valuable insights. In addition, the usefulness of a DEM for the implementation of 1D flood models was evaluated. This research provides general guidelines on how to select the best hydraulic model in terms of model (i) implementation, (ii) effect of BC, and (iii) resolution scheme.

According to the findings of this study the following conclusions can be derived:

• In general, based on the results we can infer that the equations and the numerical schemes have minor effects on the WSE computations (type I models, have similar performance). Larger differences were found for type II models, using cross sections extracted from a DEM. Thus, the geometric data is the most important factor on the water levels estimation.

• Flood Modeller set up with type I data has the better fit with the reference model in terms of water surface elevations. The similar behaviour observed between Flood Modeller and the reference model is related to the similarity in the numerical methods that both models (HEC-RAS and Flood Modeller) use to solve the equations.

• When compared with type I models, those models implemented with the type II cross sections showed underestimations of water levels. The systematic water levels underestimation is of the same order of the difference in bed levels between the two bathymetric datasets, which is expected.

• HEC-RAS developed with the type II cross sections has a good agreement with the reference model. This suggest that the HEC-RAS is less sensitive to the elevations along the river thalweg, and has almost the same model performance as the reference model.

• MIKE 11 type II presents the lowest performance, with the highest underestimation in water levels. The results show that MIKE 11 type II is governed by the downstream BC (using the water levels instead of normal flow depth).

• When compared to MIKE 11 and Flood Modeller, for the model implementation, the HEC-RAS interface is easier and less time consuming. In addition, HEC-RAS is free software in comparison to MIKE 11 and Flood Modeller which are commercial software.

• The three models present not major inconvenience in computational time, all of them have similar modeling times.

• The use of DEM elevations to implement 1D models shows underestimation of WSE higher than 0.5 m for Flood Modeller and MIKE 11 compared with the reference model. This variance can have a considerable impact on the flooded area.

• Flood Modeller and MIKE 11 type I models have a similar level of performance than the reference model. Therefore, the three modeling packages could be used for evaluated conditions, when detailed field survey (cross section data) is available.

• Floods in high mountain topography can be simulated correctly with one-dimensional models using a detailed XS survey. Model results using DEM elevations should be taken with caution because the systematic underestimation of WSE. This drawback could be mainly attributed to the high variability of the riverbed elevations due to the DEM resolution (3 × 3 m).

• HEC-RAS is the most appropriate model for modeling mountain river floods using different geometric input data (type I and type II) with more accurate results in comparison to the other two models.

• Finally, the lack of flood data records is a serious problem which limit a proper and traditional procedure for the calibration and validation of flood hydraulic models. However, this study presents the sensitivity of different hydraulic models on the computation of inundation levels. Future work must focus on the monitoring of actual inundation events to improve the prediction accuracy in the mountain range of Ecuador and evaluate the effect of sediment or large wood transport on flood hazard assessment.

The authors thank the DIUC (Research Department of the University of Cuenca) through the project ‘Evaluación del riesgo de inundación del río Santa Bárbara’. We also thank the SENAGUA (Demarcación Hidrográfica del Santiago) and the Municipality of Gualaceo for their valuable contribution with information from the study area. We also thank the anonymous reviewers for their insightful comments and suggestions.