Performances of the New HEC-RAS Version 5 for 2-D Hydrodynamic-Based Rainfall-Runoff Simulations at Basin Scale: Comparison with a State-of-the Art Model

: The Hydrologic Engineering Centre-River Analysis System (HEC-RAS), developed by the US Army Corps of Engineers, is one of the most known, analyzed and used model for flood mapping both in the scientific literature and in practice. In the recently released version (release 5.0.7), the HEC-RAS model has been enriched with novel modules, performing fully 2-D computations based on the 2-D fully dynamic equations as well as the 2-D diffusion wave equations; moreover the application of rainfall to each cell of the two-dimensional domain is now possible. Contrarily to the common applications for flood propagation in river reach, this specific module has never been analyzed in the literature. Therefore, the main purpose of this work is to assess the potential and the capabilities of the 2-D HEC-RAS model in rainfall-runoff simulations at the basin scale, comparing the results obtained using both the options (fully dynamic equations and diffusion wave equations) to the simulations obtained by using a 2-D fully dynamic model developed by the authors for research purposes. Both models have been tested in a small basin in Northern Italy to analyze the differences in terms of discharge hydrographs and flooded areas. The application of a criterion for hazard class mapping has shown significant variations between the two models. These results provide practical indications for the water engineering community in the innovative research field related to the use of 2-D SWEs at the basin scale.


Introduction
Flood hazard and flood risk maps are key elements for flood risk management around the world.They are an essential tool for flood warming, mitigation of property damage and loss of life and for flood risk communication to the stakeholders.In Europe, the EU Flood Directive (2007/60/EC) requires all Member States to assess flood hazard and risk producing, for example, flood extent, flow depth and flow velocity maps [1].In the United States, the Federal Emergency Management Administration (FEMA) produced the Flood Insurance Rate Maps (FIRMs) that identify flood insurance rate zones.HEC-RAS [35], and cellular automata approaches, that were not able to estimate inundation depth in areas with higher velocity as well predicted by HEC-RAS [36].Besides these studies, more focused on the purposes of this work, several recent applications of HEC-RAS model can be found in the literature, confirming the interest of the international water research for this software [37][38][39].
It should be observed that the new HEC-RAS 5.0.7 has been developed also to deal with a challenging issue for both scientist and technicians, related to the use of the SWEs for the description of the hydrodynamic-based surface runoff computations in rainfall-runoff simulations at the catchment scale.This approach, sometimes referred as direct rainfall method [36,40], consists in the application of rainfall to each cell of the two-dimensional domain.The consequent runoff is then routed by the 2D hydrodynamic model.This approach seems to be attractive for modelers and scientists because it significantly reduces the need for hydrologic modelling for the runoff estimation and, consequently, the number of models for generation and propagation of surface runoff at a basin scale [41].Moreover, it has been recently observed that the hydrodynamic-based rainfall-runoff simulations can shed new light to the analysis of basin drainage network, introducing novel scaling laws and an innovative method for channel heads detection [42].This kind of approach is more and more analyzed in the literature, as testified by the plethora of the papers in this field (see among others, the most recent ones [43][44][45][46][47][48][49].
Contrarily to the common applications for flood propagation in river reach, this specific module implemented in the new HEC-RAS 5.0.7 has never been analyzed in the literature.The authors believe that, due to the popularity of this software and the innovative research related to the use of 2-D SWEs at the basin scale, there is the need to analyze the performances of HEC-RAS also for this specific kind of simulations, because of the relevance for the water engineering community.
Therefore, the main purpose of this work is to assess the potential and the capabilities of the 2-D HEC-RAS model in rainfall-runoff simulations at the basin scale, comparing its results to the simulations obtained by using a 2-D fully dynamic model developed by the authors already used in the literature.

Case Study: Scuropasso Basin
The case study analyzed here is a part of the Scuropasso Basin, a tributary of the Po River, flowing in Lombardy Region.Scuropasso is a small basin (43 km 2 ) with an average slope of 12.95%.The sub-catchment area is equal to 13.75 km 2 with an elevation at the basin closure section of 187 m a.s.l., a maximum elevation of 537 m a.s.l. and a mean elevation of 187 m above sea level.The subcatchment mean slope is 24.57%.Agricultural areas represent the dominant land use while urbanization is composed by very small isolated urban agglomeration.A high resolution airborne LiDAR covers the sub-catchment area, providing a sampling density of four points over 1 m 2 (see Figure 1).

Hydrodynamic Models for Surface Runoff
3.1.1.HEC-RAS 5.0.7 HEC-RAS is a complex modelling tool which include a large number of applications for 1D model and recently it solves the fully 2D shallow water equations written as follows: where H(x,y,t) = z(x,y) + h(x,y,t) is the surface elevation (m), z is the cell elevation in the cartesian coordinates x, y, h is the water depth (m), p = hu and q = hv are the specific flow in the x and y directions (m 2 s −1 ), u and v are the velocities in x and y respectively, r is the net rain (m), g is the gravity acceleration (ms −2 ), n is the Manning's roughness coefficient (s m −1/3 ), ρ is the water density (kg m −3 ), τxx , τyy and τxy are the components of the stress tensor and f is the Coriolis parameter (s −1 ).When the diffusive wave is selected, the inertial terms in Equations ( 2) and (3) are neglected.The above equations are solved with an implicit finite-volume scheme.HEC-RAS 5.0.7 introduced the capability to compute 2D fully momentum shallow water equation based on unstructured grid joint to a sub-grid based on a high-resolution DEMs, (e.g., LiDAR) which may drastically affect models results as well as performance [33,50].Hence, the subgrid approach provides information given by the high resolution DEMs bathymetry.Thus, Equation (1) is presented in the integral form: where Ω is the volumetric three-dimensional space occupied by the fluid, R is the source term (net rainfall), Vk is the average velocity at the k face with the normal vector nk, Ak(H) the averaged area at the k face.The superscript n represents the time step and ∆t the integration time.The integral form of Equation ( 4) is needed to compute extra information (e.g., hydraulic radius, volume, crosssectional data) in a preprocessing stage for every computational cell face [27,50].The sub-grid approach can be used to increase the cell dimensions preserving the details of the bathymetry [50].
Finally, in the HEC-RAS 5.0.7 the Courant-Friedrichs-Lewy (CFL) condition is implemented into the model.

Shallow Water Equations-Finite Volume Code
The implemented code (hereinafter called SWE-FVM) for surface runoff simulation is based on the 2-D fully-dynamic SWE that can be expressed in the following form: where: in which: t is time; x, y are the horizontal coordinates; h is the water depth; u, v are the depth-averaged flow velocities in x-and y-directions; g is the gravitational acceleration; r is the rain intensity and f are the infiltration losses; S0x, S0y are the bed slopes in x-and y-directions; Sfx, Sfy are the friction slopes in x-and y-directions.For the numerical integration of Equation ( 5), the finite volume method (FVM) approach has been considered.In particular, the Roe scheme [51], first order accurate in time and space, has been used for the computation of the numerical fluxes.Moreover, suitable numerical treatments of the source terms and wet/dry fronts have been applied [52,53].The mathematical aspects and performances in real cases have been already presented in previous works [42,54,55] and, therefore, are not reported here.The code is parallelized using MPI directives.

Computational Grids and Boundary Conditions
The generation of computational grid is one of the key elements for the application of 2-D models because of its influence on the accuracy and on the computational times [56].Though some methods have been proposed to generate a priori grids, having different resolution throughout the basin [57], the evaluation of the performances of the SWE-FVM and HEC-RAS models will be shown on two different computational grids with uniform resolution.The first one is a fine triangular grid, of 5 m side.This grid represents the most accurate computational domain used in this work (fine grid, hereafter).On the contrary, the second one has an average triangular side of 10 m, being much coarser (coarse grid hereafter).The coarse grid dimension choice is in accordance to literature which assert that a resolution up to 10 m does not influence geomorphologic indices like flow patterns and basin hydrologic responses [58].The fine grid dimension is based on the availability of the HEC-RAS geometry editor to manage a large number of elements.In order to run properly the two overland flow models, the use of the same grid in both the models is needed.Namely, the computational grid has been created through the external mesh generator SMS 13.0 (Aquaveo-Surface-water Modeling System, https://www.aquaveo.com/software/sms-surface-water-modeling-system-introduction).The SMS software is well-known in the scientific context (see for example [18]).It is not an open source package but is able to export the computational grid in a 0.2 dm format readable in the SWE model as well as in a format compatible with HEC-RAS.This process ensures the uses of the same computational domain by the two models.Finally, the fine grid is composed by 1,293,762 unstructured triangular elements, whereas in the coarse mesh number of elements is 323,161.
Whatever is the grid used by HEC-RAS, the software applied the so-called high-resolution subgrid model [59], which allows coarse computational meshes to be used without losing the detail of the underlying terrain [50].
As regards the boundary conditions, we have not specified the upstream boundary condition due to the nature of the simulations.Specifically, the approach used here is devoted to generate the surface runoff and, as consequence, the flood wave using only one model (the 2-D SWEs) without using the traditional modelling chain in which the flood hydrograph is computed using a hydrological model and, then, a hydrodynamic model is used for flood propagation.Therefore, there is no need to compute an upstream boundary condition, so the use of the net rainfall on grid avoids the computation of synthetic hydrograph with a given return period.The method used for the computation of the net rainfall is presented in the next section.
As regards the downstream boundary condition, transmissive boundary conditions have been set for the SWE-FVM model while the "normal depth" option has been selected for the HEC-RAS model, that means uniform flow condition.

Net Rainfall Computation and Roughness Estimation
The intent of the comparative analysis shown in this work is to analyze the discrepancy between the models results in cases for which no observed data, related to rainfall, flooded areas and discharge, are available.This is the typical situation in small basins that, very often, are ungauged catchments for which is not possible to have data for model calibration and validation.As a consequence, "blind" simulations should be performed.In our view, in situations like these, is important to assess the performance of a commercial software using a state-of-the art model, developed for research purposes, as the benchmark solution.For these reasons, we have not considered an observed event for model comparisons.Instead, following a technical approach, two rainfall events were considered to compare the two models, having a return period of T1 = 30 and T2 = 200 years respectively.The Chicago hyetograph was adopted using the rainfall intensity-duration curve taken from the Regional Agency for the Protection of the Environment (ARPA) website, which is evaluated taking into accounts data from 1987 to 2011 (http://idro.arpalombardia.it).The intensityduration curve used is the following: h = 24.53t 0.307 (7) The net rainfall was evaluated with the well-known SCS-CN method [60], considering an average Curve Number of 86 for the sub-basin and adopting the initial loss as the 20% of the maximum infiltration losses.The basin Runoff Curve Number was determined for soil in average moisture conditions.Figure 2 shows the net hyetographs for the two adopted return periods.
Therefore, a total of 4 simulations were performed for each model, considering rainfall events having 30 years and 200 years of return period and for the two grids described in the previous section.As regards the Manning's roughness coefficients, we set the following values according to the available map related to the use of soil: 0.035 s m −1/3 , 0.06 s m −1/3 and 0.1 s m −1/3 for agricultural areas forest and urbanized areas respectively.

Flood Hazard Estimation
Several methods are available in the literature for flood hazard estimation.In this work, we refer to the method proposed by [61], recommended by [62] and recently analyzed in [13,[63][64][65].This criterion, also known as Australian Institute for Disaster Resilience (AIDR) hazard classification, combines information related to water depth h and velocity V, defining six hazard classes.We have selected this criterion essentially due to its simplicity and because it considers the main targets exposed to flood risk (see Table 1).hV > 4.0 1 A minimum limit (0.1) has been assumed to restrict the analysis only to the flooded areas.

Performance Measures
The overlapping of the flooded areas belonging to the different AIDR classes can be performed through several indexes, widely used in the literature.Specifically, in this work we use the Hit Rate, False Alarm Ratio, Critical Success Index [66].
The Hit Rate (HR), sometimes referred to as the probability of detection, is a simple measure that indicates how well the model replicates the benchmark data without penalizing for overprediction.It is defined as: and it is equal to 0 when all wet areas of the two models are totally separated whereas is 1 when all wet cells in the benchmark data are wet in the model x, being Ax the modeled inundated area and AR is the benchmark inundated area.The False Alarm Ratio (FAR) is a measure of model overprediction.It is defined as: The values of FAR ranges between 0, when no false alarms are present, and 1 when there are all false alarms.
The Critical Success Index (CSI) is defined as: CSI = 0 means no match between model x and the benchmark, CSI = 1 implies perfect match.Finally, a simple measure of error bias B can be defined as: When B ranges between 0 and 1, the model x shows a tendency to underpredict the benchmark, whereas B between 1 and ∞ indicated the tendency of the model to overpredict.

Results and Discussion
The comparison between the performances of HEC-RAS and the SWE model is presented here.Attention will be first focused on the simulated discharges hydrographs at the basin outlet, then will be oriented to the flooded areas mapping and finally to the flood hazard assessment according to the predictions of the two models.

Discharge Hydrographs at the Basin Outlet
HEC-RAS has been used according the two options: full momentum and the diffusion wave approximation, hereafter identified as FDW (Fully Dynamic Wave) and DSW (DiffuSive Wave) respectively.Figures 3 and 4 show the comparison of the discharge hydrographs at the basin outlet simulated by the two models for a return period of 200 years and 30 years respectively, using two different grid resolutions.The FDW version of HEC-RAS gives very similar results in terms of the flood wave shape, peak discharge and time to peak in respect of SWE-FVM (see Table 2).In particular, for the smaller grid resolution, HEC-RAS results are practically superimposed to those ones obtained by SWE-FVM model (see Figures 3a and 4a).However, it is interesting to stress that the grid resolution has a very small influence in terms of flood hydrographs.Specifically, the SWE-FVM model results are quite insensitive to the analyzed grid resolution while, considering HEC-RAS, the peak discharge for the 10-m resolution is slightly lower (approximately 4%) than the homologous one related to the 5-m resolution (Figures 3b and 4b).This behavior can be observed for both the analyzed return periods.The flood hydrographs simulated by the DSW version of HEC-RAS are slightly different from the other two simulated waves.In particular, the peak discharges values are overestimated (about 12%) if compared to the other two simulations for both the return periods.The DSW model shows a tendency to simulate faster flood waves (see Table 2).These differences seem to decrease with the increase of the grid size.The sensitivity of the DSW to the grid resolution might be explained considering that, in this approach, the velocity is depending almost entirely on the slope of the terrain.
Finally, it is interesting to discuss the HEC-RAS performance in terms of overall computational times.HEC-RAS simulations were computed on a AMD Ryzen 2700X processor Eight-Core Processor 3.70 GHz (AMD, Santa Clara, CA, United States), with 32 GB RAM.
Table 2 shows that the simulations carried out using the finest grid resolution are not particularly feasible for practical study, resulting more than four times greater the duration of the real flood event.The reduction of the computational times obtained using the DSW option is not significant (about 5-6%).Conversely, the performances related to the coarsest grid are interesting.The FDW option has provided results in a shorter time respect to real event duration (see Table 2) for both the return periods, while the DSW option allowed a further decrease of the runtimes in respect of the FDW ones (about 40%).Obviously, these performances cannot be compared to the SWE-FVM, because the latter model has been parallelized using MPI directives and computed on a cluster using 192 threads.This justifies the very good performances of the SWE-FVM, which is faster than HEC-RAS (about oneorder of magnitude).
Table 2. Information related to discharge hydrographs computed at the outlet and overall computational times of the models as a function of the return period and grid element side.The analysis that follows will be focused only on the comparison between the FDW option of HEC-RAS and SWE-FVM, due to the very similar performances of the two models in terms of simulated discharge hydrographs.

Simulated Flooded Area Extent
Beside the discharge hydrograph at the outlet, one of most important information for flood hazard assessment is the flooded areas delimitation.For rainfall-runoff simulations at the basin scale, since the rainfall input is applied to each computational cell, it is necessary to define a water depth threshold for the identification of the flooded cells.
Though some methods have been proposed in the literature (see [42]), in this work we have analyzed the flooded areas for different, very small water levels.Table 3 shows the extent of the flooded areas, computed by the two models, equal to or greater than pre-fixed water depth thresholds (5, 10, 15 and 20 cm).In Figure 5, the comparison between the flooded areas computed by the two models for each return period and grid size is shown.For the fine grid, the overall results seem to be comparable while, for the coarse one, some significant differences appear.Figure 6 shows the relative variations in the flooded area, computed by Equation ( 12), taking SWE-FVM results as reference solution, for both the different grid resolutions and the return periods: where AHR(h) and ASF(h) represent the flooded area related to the threshold h as simulated by HEC-RAS and SWE-FVM respectively.Greater variations are found for the bigger cell dimension.Neglecting the first threshold value, HEC-RAS always overestimates SWE-FVM flooded areas with variations approximately equal to 5-10% and 20-30% for the fine grid and coarse grid, respectively.For some water depth thresholds, the variation due to the 10-m resolution grid is up to three times the one computed for the 5-m resolution grid (see Figure 6a).The decrease in the difference between model results, as the grid size increases, might be explained considering that, at lower resolutions, the small-scale characteristics, which may be handled differently in the models, become negligible and thus the results converge.If a more severe event is considered (Tr 200 year vs. 30 year), the variations between the two models generally decrease for the coarse grid while are substantially constant for the fine one (Figure 6b).Therefore, the variations in flooded areas increase as the cell size increases and as the return period decreases, at least for the coarse grid.
The smaller water depth threshold is the only value for which HEC-RAS provided an underestimation of the flooded area.In order to better understand what happens in this situation, the extent of the flooded area computed by the two models in two particular zones of the sub-basin, one upstream and one close to the closure section, have been analyzed.Figures 7 and 8 show the flood extent related to a return period of 200 years, using the fine grid resolution.In those figures, the river network simulated by HEC-RAS is characterized by quite disconnected flooded areas.Looking at the Figure 7, it is clear that the tributary highlighted by the red arrow does not flow into the main channel or, in other words, the tributary is not connected to the main channel.Moreover, the part of the network having smaller Horton orders is quite fragmented and, in general, underestimated.This is enhanced in Figure 8 where the network described by HEC RAS gets to a cloud of small areas disconnected one from the other.This observation can be considered a generalization of what has been already noted by [33], who observed, in a traditional flood propagation study, a certain amount of hydraulically disconnected flooded areas, that is areas simulated as flooded but disconnected from the main inundated area.This could be induced by the way in which HEC-RAS distributes the water within a computational cell, based also on the volume-elevation curve drawn for each cell-face in pre-processing stage.However, the current study did not look into the specific features of this property.

Flood Hazard Classification
The water depths and velocity fields simulated by the two models have been used for flood hazard mapping purposes.As mentioned in Section 3.4., the criterion proposed in [62] have been considered here according to which six hazard classes can be identified as a function of h and V.The results obtained using this approach are reported in Figure 9 for each return period and grid resolution, showing the distribution of the flooded areas belonging to each class.
Looking at the H2, H3 and H4 classes, the flooded areas simulated by HEC-RAS are generally comparable to those computed by SWE-FVM.However, the flooded area falling in class H1 is significantly overestimated.The latter discrepancy may be ascribed to the way HEC-RAS redistributes the water depth in the computational cell, defined as the sub-grid approach.Moreover the two models give different flooded area for smaller water depths.These variations seem to decrease with the return period (Figure 9a,b, or Figure 9c,d) and to increase with the cell size (Figure 9a,c, or Figure 9b,d).Another significant difference between the two models can be found in the class H5, in which the HEC-RAS results underestimate the results obtained by SWE-FVM model.This variation seems to increase as the return period decreases (Figure 9a,b, or Figure 9c,d).Figure 9 gives an overall description of the flooded areas extent simulated by the two models but it cannot provide any information about their overlapping from a topographic point of view.To that purpose, it is necessary to use performance indexes mentioned in Section 3.4.For the sake of brevity, only the finer grid results (cell size equal to 5 m) are reported in Tables 4 and 5, related to the return period 30 years and 200 years, respectively.
The values of the CSI index highlight that the areas flooded by the two models are far to be coincident, at least for the first five classes.This is further confirmed by the HR values that show a value close to 1 only for the H6 class.The values assumed by B seem to suggest the HEC-RAS simulation underestimate the flooded area belonging to the classes H4 and H5 while it is overestimated especially in the first three classes Generally, these measures seem not to be significantly affected by the flood magnitude and therefore by the return periods considered here.Figures 10 and 11 show the comparison between the hazard maps computed by the two models in the two selected areas discussed in Figures 7 and 8.While in the river network the hazard classes computed by the two models seem to be close together, it is interesting to note that HEC-RAS defines zones far from the network in H1 and H2 classes.In these areas, water depths are smaller than the selected threshold (0.05 m), meaning that very high velocities are found.

Conclusions
This work has been focused on the performance of the new HEC-RAS 5.0.7 for the hydrodynamic-based surface runoff computations in rainfall-runoff simulations at the catchment scale.Contrarily to the common applications for flood propagation in river reach, this work shows for the first time a critical discussion of this specific new module, never analyzed in the literature.
Specifically, in order to assess its potential and the capabilities, HEC-RAS has been applied to a small basin in Italy and compared to the results obtained by a 2-D fully dynamic model developed by the authors (SWE-FVM).The main conclusions of this paper are listed below:


Using the "full momentum" options, HEC-RAS has provided very similar results in terms of the flood wave shape, peak discharge and time to peak in respect of SWE-FVM.The diffusive wave version (DSW) of HEC-RAS overestimates the peak discharges values computed by the other two models, showing a tendency to simulate faster flood waves;  The reduction of the computational times obtained using the DSW option is not significant for finer grid but it becomes important using a coarser one, with time savings up to 40%;  The analysis of flood extent, carried out for different pre-fixed water depth thresholds, highlights that, except for the smaller threshold, HEC-RAS has always overestimated the SWE-FVM flooded areas with variations approximately equal to 5-10% and 20-30% for the fine grid and coarse grid respectively.The variations in terms of flooded areas increase as the cell size increases and the return period decreases, at least for the coarse grid;  The flooded area related to the first threshold (0.05 m) is underestimated in respect to the SWE-FVM one.In particular, the part of the network having lower Horton-Strahler order is described by HEC-RAS as disconnected, while HEC-RAS and the SWE-FVM models describe the main channel in a similar way;  The application of the AIDR criterion for hazard class mapping has highlighted significant variations between the two models.For example, in respect to the SWE-FDW results, HEC-RAS has significantly underestimated the flooded areas related to the H5 class while an important overestimation has been observed for the H1 class;  The application of the performance measures has highlighted that the areas flooded by the two models belonging to a specific hazard class are far to be coincident, at least for the first five classes.
These results highlight the need of further investigations and a larger number of case studies in order to generalize these conclusions.However, on the basis of this first study, it seems that HEC-RAS can be considered as a reliable model for discharge hydrograph computation in rainfall-runoff simulations.

Figure 3 .
Figure 3. Discharge hydrograph simulated with HEC-RAS and SWE-FVM models for a return period of 200 years: element size 5 m (a) and 10 m (b).

Figure 4 .
Figure 4. Discharge hydrograph simulated with HEC-RAS and SWE-FVM models for a return period of 30 years: element size 5 m (a) and 10 m (b).

Figure 5 .
Figure 5. Flooded areas computed by HEC-RAS and SWE-FVM for pre-fixed water depth thresholds considering: T = 30 years, grid size 5 m (a) and 10 m (b), and T = 200 years, grid size 5 m (c) and 10 m (d).

Figure 6 .
Figure 6.Relative variations in terms of flooded areas extent for given water depths: T = 30 years (a) and T = 200 years (b).

Figure 7 .
Figure 7. Detail of the flooded areas computed by HEC-RAS and SWE-FVM for water depth threshold equal to 0.05 m, considering T = 200 years and grid size 5 m.

Figure 8 .
Figure 8. Detail of the flooded areas computed by HEC-RAS and SWE-FVM for water depth threshold equal to 0.05 m, considering T = 200 years and grid size 5 m.

Figure 9 .
Figure 9. Flooded areas belonging to specific AIDR hazard classes as simulated by the two model: T = 200 years and grid resolution equal to 5 m (a) and 10 m (c), T = 30 years and grid resolution equal to 5 m (b) and 10 m (d).

Figure 10 .
Figure 10.Hazard maps computed by HEC-RAS and SWE-FVM for water depth threshold equal to 0.05 m, considering T = 200 years and grid size 5 m.

Figure 11 .
Figure 11.Hazard maps computed by HEC-RAS and SWE-FVM for water depth threshold equal to 0.05 m, considering T = 200 years and grid size 5 m.

Table 3 .
Flooded areas for a given water depth: comparison between the models.

Table 4 .
Value of the performance indexes related to flooded areas extent of each hazard class (T = 30 years).

Table 5 .
Value of the performance indexes related to flooded areas extent of each hazard class (T = 200 years).