Using 2D HEC-RAS Modeling and Embankment Dam Break Scenario for Assessing the Flood Control Capacity of a Multi-Reservoir System (NE Romania)

Abstract

Using hydraulic modeling techniques (e.g., one-dimensional/two-dimensional (1D/2D) hydraulic modeling, dam break scenarios) for extracting the flood settings is an important aspect of any action plan for dam failure (APDF) and flood mitigation strategy. For example, the flood hydraulic models and dam break scenario generated based on light detection and ranging (LiDAR)-derived digital elevation models (DEMs) and processed in the dedicated geographic information systems (GIS) and hydraulic modeling software (e.g., HEC-RAS—Hydrologic Engineering Center River Analysis System, developed by USACE HEC, Davis, CA, USA) can improve the flood hazard maps in case of potentially embankment dam failure. In this study, we develop a small-scale conceptual approach using 2D HEC-RAS software according to the three embankment dam break scenarios, LiDAR data (0.5 m spatial resolution), and 2D hydraulic modeling for the Başeu multi-reservoir system which belongs to the Başeu River (NE Romania) including R1—Cal Alb reservoir, R2—Movileni reservoirs, R3—Tătărăşeni reservoirs, R4—Negreni reservoirs, and R5—Hăneşti reservoirs. In order to test the flood control capacity of the Bașeu multi-reservoir system, the Cal Alb (R1) dam break scenario (piping failure) was taken into account. Three 2D stream flow modeling configurations based on R1 inflow rate with a 1% (100 year), 0.5% (500 year), and 0.1% (1000 year) recurrence interval and the water volume which can be accumulated with that specific inflow rate (1% = 10.19 × 10⁶ m³; 0.5% = 12.39 × 10⁶ m³; 0.1% = 17.35 × 10⁶ m³) were computed. The potential flood wave impact was achieved on the basis of different flood severity maps (e.g., flood extent, flood depth, flood velocity, flood hazard) generated for each recurrence interval scenario and highlighted within the built-up area of 27 settlements (S1–S27) located downstream of R1. The results showed that the multi-reservoir system of Bașeu River has an important role in flood mitigation and contributes to the APDF in the context of climate change and the intensification of hydrological hazard manifestation in northeastern Romania.

1. Introduction

The location of human settlements near water sources (e.g., lotic ecosystems—rivers, streams; lentic ecosystems—lakes, ponds, marshes) has proven to be a very important aspect of habitation practice since ancient times [1,2,3]. In this context, when the water supply was limited during drought periods, the water reserve necessity led to the construction of the first earth dams and storage lakes on Earth (e.g., Jawa Dam, Eastern Jordan—ca. 4000 BCE; Jaddel-Kafara Dam, Egypt—ca. 2700 BCE; ca. 10,000 rice-paddy irrigation dams, Japan—ca. 200–300 CE) [4,5,6]. Initially used for seasonal crops irrigation (e.g., rice, vegetables) or for domestic water supply, currently, earth-filled dam lakes are common around the world and bring tremendous benefits to society such as flood control and management, fish farming and aquaculture, hydropower, irrigation, water supply, navigation, recreation, tailing industry, and ecological management [7]. According to the International Commission on Large Dams (ICOLD) worldwide database which monitors dams with heights ≥15 m from the lowest foundation to crest or with heights between 5 m to 15 m and ≥3 million m³ volume of water stored, 65% of all dams are earth-filled (ca. 37,500) of which more than 10,000 dams were built-up for flood mitigation [8]. For example, the Nurek embankment dam on the Vakhsh River, Tajikistan, constructed primarily for hydroelectric power generation, is the world’s second highest embankment dam (300 m high) [9]. Currently, the Nurek reservoir has multiple functions including the flood control in the context of climate change and the intensification of hydrological hazards manifestation [9].

Even if the embankment dams bring tremendous benefits to society, the risk of various types of failure exists, and the flood waves produced by dam breaks are more destructive than the natural floods [10,11,12]. The causes that can lead to a dam’s failure are generally classified as natural causes (e.g., heavy rainfall, storms, tsunami, landslides, and earthquakes) or anthropogenic causes (e.g., dam design errors, mismanagement in water storage, maintenance errors, deliberate dam destruction, terrorism) [13]. In this context, there are many examples of dam failure around the world [5]. One of the most catastrophic events took place in 1938 due to the destruction of Henan Dam on the Yangtze River, China. The flood wave affected 44 provinces, 800,000 inhabitants were killed, and the number of refugees reached about 4 million [5]. In Europe, the most destructive hydrological event due to dam failure occurred on 9 October 1963, on the Erto e Casso River, Italy, when a large landslide produced a 250 m high flood wave which destroyed the Vajont Dam and killed 2600 inhabitants [14]. In Romania, due to the heavy rainfall which occurred in 1991 in the northeastern part of the country, the Belci Dam on the Tazlău River was destroyed by floods that killed 25 people [15]. Overall, according to the last studies related to this topic [16], no matter what triggering agents exist, the embankment dam failure modes can be categorized as either overtopping [17] or piping [18]. Therefore, according to the well-documented dam failure reports provided by the United States Geological Survey (USGS), approximately 34% of dam failures were caused by overtopping, 30% by foundation defects, and 28% by the piping process [5].

In this context, due to the significant damage caused by dam failure during the time and the high recurrence of these events determined, the stakeholders (e.g., engineers, scientists, authorities involved in emergency management) aimed to better understand the breach process of earthen dams. Accordingly, the first studies in this regard took place in France at the end of the 19th century, when some solutions were theorized on the problems caused by the earthen dam failure. Subsequently, according to [19], there are two important periods in which significant progress was made: (a) the first period from 1950 to 1990 when many laboratory experiments were performed and the first computerized hydraulic modeling programs were implemented (e.g., DAMBRK ’88, BREACH 7/88, DWOPER 8/89, SMPDBK 9/91, MIKE SHE) [20,21,22]; (b) the second period from 1990 to present, when more hydraulic modeling software programs were upgraded (e.g., one-dimensional/two-dimensional (1D/2D) HEC-RAS (Hydrologic Engineering Center River Analysis System, developed by USACE HEC, Davis, CA, USA) and determined the improvement of dams failure forecasting [23,24,25,26]. Currently, floods caused by dam-break events still force society to create more upgraded software applications and geographic information system (GIS) tools used for forecasting, prevention, and disaster risk reduction caused by dam failure. For example, through the 2007/60/CE Directive, each member state of the European Union (EU) has the obligation to generate, for each large river basin under administration, individual flood risk management plans (FRMPs) [27,28].

In Romania, according to 1422/192 Ordinance from 16 May 2012, which is additional to 2007/60/CE Directive, the owners or the administrators of a reservoir with a dam height over 10 m and/or with water volume more than 10 million m³ have the obligation to make an APDF [29]. Furthermore, in the same ordinance, it is specified that an APDF can be made for reservoirs with a dam height lower than 10 m and/or a water volume > 10 million m³ only if a hypothetical dam failure can cause significant material damage or human casualties. However, with more than 1000 registered dams in Romania (e.g., 250 dams classified as big dams, 85 classified as embankment big dams), only few of them have an APDF. In this case is also the Cal Alb reservoir (R1) located in the upper sector of the Başeu multi-reservoir system (Figure 1a). Similar to other embankment dams in the country which were built up during the communist period (before 1989), the Cal Alb reservoir has not yet received an APDF, even if it accomplishes all the 1422/192 Ordinance requirements (e.g., dam height more than 10 m, water volume greater than 10 million m³, frequent dam faced infiltrations). Thus, given the existence of a precedent on the integrity of the dam body (e.g., piping process), the fact that the Cal Alb reservoir provides the water resource for an important community in the Botoșani County, the fact that downstream of Cal Alb reservoir are located 25 rural settlements (S1–S10; S12–S23; S25–S27) and two urban settlements (S11—Săveni and S24—Ștefănești), we consider that the Cal Alb reservoir provides the optimal conditions to test the integrity of the dam body and the capacity of the multi-reservoir system (R2, R3, R4 and R5) located downstream to take over and reduce the flood wave produced by a potential dam break scenario [30] (Figure 1b).

In this study, in order to test the flood control capacity of the Bașeu multi-reservoir system, we generated the first dam break scenario for the Cal Alb (R1) reservoir using a LiDAR-derived DEM with 0.5 m spatial resolution and 2D HEC-RAS hydraulic modeling [30,31,32,33,34,35]. The different scenarios were based on three R1 inflow rates with different recurrence intervals (1% or 100 years; 0.5% or 500 years; 0.1% or 1000 years) and the water volume which can be accumulated with that specific inflow rate (1% = 10.19 × 10⁶ m³; 0.5% = 12.39 × 10⁶ m³; 0.1% = 17.35 × 10⁶ m³) [7]. We used 2D HEC-RAS hydraulic modeling instead of 1D HEC-RAS hydraulic modeling because the 1D method cannot calculate the propagation of the multidirectional flood wave [31,36,37,38]. On the other hand, the 2D hydraulic models are the most appropriate solution in order to simulate the flood wave lateral diffusion [39,40,41,42,43]. Therefore, in order to accomplish the main objective of the study, a 2D hydraulic break scenario was created using HEC-RAS v.5.0.7 software [23,24,25]. Considering the recent infiltration history of Cal Alb (R1) dam, a piping failure mode was used [7,44]. Additionally, on the basis of the adapted Federal Emergency Management Agency (FEMA) methodology on hazard classification [45], after the Australian Institute for Disaster Resilience (AIDR) hazard classification [46], the flood severity assessment was done for all 27 settlements potentially affected by R1 dam failure. The results provide for the first time the flood control capacity of the Bașeu multi-reservoir system in case of Cal Alb (R1) failure and contribute to the APDF in the northeastern region of Romania.

Study Area

Başeu river basin (96,770 ha) is located in northeastern Romania, in the administrative territory of Botoșani County. With a total length of 118 km, the Başeu River is one of the most important tributaries of Prut River (natural border between Romania and Republic of Moldova) (Figure 1a) [47,48,49]. The lowest altitude is found on the common floodplain of Bașeu and Prut River (54 m above sea level (a.s.l.)) and the highest altitude is on the top of Pădurea Cristesti Hill (318 m a.s.l.) [7,30]. The geological conditions are characterized by the presence of recent alluvial deposits (upper Holocene) accumulated over 200–300 m sandy clay deposits (Miocene). Climate conditions control 60% of the water flow rate, often overtaken in transition seasons and in periods with maximum rainfall. The average air temperature values are between 7.7 °C and 9.6 °C. The annual average precipitation is between 450 and 624 mm, with peak values occurring in the high basin areas. The groundwater contributes between 30% and 40% of the average annual flow. The multi-annual average flow rate at Ștefănești station (south of the basin) is 1.6 m³/s, with fluctuations between 0.1 m³/s (2016) and 6.9 m³/s (1969). The historical flow was recorded in July 1969 and was 330 m³/s. Negative hydrological events associated with the maximum discharge occurred in 1973 (100 m³/s) and 2005 (124 m³/s) [30].

The Başeu multi-reservoir system consists of five permanent and nonpermanent multi-reservoirs areas (total water surface 953.3 ha) from upstream to downstream as follows: permanent multi-reservoirs—Cal Alb (R1) and Negreni (R4); nonpermanent multi-reservoirs—Movileni (R2), Tătărăşeni (R3), and Hănești (R5) (

Table 1. Characteristics of the main dams within Bașeu river basin according to Romanian National Classification (RNC) [51].

Reservoir Dam	Dam Height (m)	Dam Type	Reservoir Volume (Million m3)	1 (%) of Total Storage Capacity	Class of Importance	Category
R1	14.5	Earth dam	16.3	100	III	B
R2	4	Earth dam	2.5	28.8	IV	C
R3	4.5	Earth dam	3.3	63.3	IV	C
R4	11	Earth dam	25	41.2	III	B
R5	5	Earth dam	5.95	37	IV	C

¹ The percentage of the total storage capacity used as initial condition within two-dimensional (2D) HEC-RAS (Hydrologic Engineering Center River Analysis System, developed by USACE HEC, Davis, CA, USA) hydraulic model for each reservoir.

; Figure 1b). According to Romanian legislation regarding the management and operating rules of dam reservoirs, permanent reservoirs are filled up with water throughout the operation period, and the nonpermanent reservoirs are maintained empty for flood mitigation. In the last decade, within Başeu river basin, the nonpermanent reservoirs were used more for fish farming instead of flood control. Moreover, more than 50 polders or lateral structures (total surface of 241.8 ha), which ensure the protection of 1071 ha (690 ha of agricultural land, 320 ha of fish farming, 61 ha of urban area), were built up in the Başeu floodplain over time. Downstream of Cal Alb (R1) reservoir are located 27 settlements (S1–S27) with a total number of 33,869 inhabitants (S1–S10; S12–S23; S25–S27—rural settlements with 19,094 inhabitants; S11 and S24—urban settlements with 14,775 inhabitants) (

Table 2. Settlement codes and habitation features within Başeu multi-reservoir system according to the 2016 census [30,47].

Code	Settlement Name	Settlement Type	Built-Up Area Surface (ha)	Number of Inhabitants	Population Density (Inhabitants/km2)
S1	Havârna	Village	418	2823	675
S2	Gârbeni	Village	64	420	656
S3	Tătărășeni	Village	114	797	699
S4	Balinți	Village	67	396	591
S5	Niculcea	Village	13	41	315
S6	Negreni	Village	109	909	834
S7	Știubieni	Village	214	1727	807
S8	Chișcăreni	Village	61	495	811
S9	Sat Nou	Village	28	124	443
S10	Petricani	Village	97	645	665
S11	Săveni	City	229	8145	3557
S12	Vlăsinești	Village	156	1596	1023
S13	Bozieni	Village	34	255	750
S14	Sârbi	Village	104	1064	1023
S15	Miron C.	Village	83	472	569
S16	Slobozia	Village	18	86	478
S17	Hănești	Village	128	1127	880
S18	Moara J.	Village	22	89	405
S19	Mihălășeni	Village	107	743	694
S20	Negrești	Village	32	300	938
S21	Păun	Village	41	289	705
S22	Năstase	Village	31	183	590
S23	Stânca	Village	105	812	773
S24	Ștefănești	City	242	6630	2740
S25	Bădiuți	Village	108	985	912
S26	Bobulești	Village	206	1322	642
S27	Românești	Village	183	1394	762

; Figure 1b). The total built-up area is 3014 ha (rural built-up area—2613 ha; urban built-up area—401 ha) with an average population density of 1123 inhabitants/km² (rural population density—731 inhabitants/km²; urban population density—3685 inhabitants/km²) [30] (Figure 1b).

Within the Bașeu hydrographic basin, the Cal Alb reservoir (R1) (water surface—174.5 ha at normal water level (NWL)) is placed in the upper sector of Başeu River (northeastern part of the basin) (Figure 2a), at 900 m downstream of the confluence with the Ciolac Brook. The river basin of the R1 reservoir has a total surface of 18,300 ha, with an altitude between 318 m (Pădurea Cristesti Hill) and 123 m (R1 water mirror level). The average altitudinal class is between 125 and 200 m (>50%). The embankment dam was built between 1971 and 1973 through decree No. 118 of 28 May 1971 [50]. According to the Romanian National Classification (RNC) of reservoirs and lakes, the Cal Alb reservoir is classified as category B (special importance) and as the third class of importance [51]. The embankment dam of R1 reservoir (Figure 2b) was built up transversely on the Başeu riverbed axis. The main feature characteristics of the R1 dam are as follows: the longitudinal section is trapezoidal where the base has a length of 295.07 m and the dam crest width does not exceed 5 m; the maximum height of the dam is 14.5 m (10.22 m without dam crest); the upstream slope report coefficient is 1:3.5 due to the concrete slab wall for protection, and the downstream slope is 1:3 being equipped with a grass layer for protection; the dam crest elevation is 134.13 m a.s.l. and the dam base elevation is 123.91 m a.s.l. (Figure 2c). The dam body was filled with local materials (>95% clay) and isolated with a 10 cm concrete layer on the water side and with a 20 cm grass layer on the opposite side. Based on the lithological profiles, two types of geological formations are highlighted: the geological foundation, made up of gray clays deposits, and the superficial layer (Başeu alluvial deposits) with a thickness between 6.5 m and 8.5 m

This status was acquired because the R1 reservoir ensures the water supply for irrigation of 1060 ha of agricultural land and also ensures 88 ha of water surface for fish farming. With a dam height of 11 m and a total storage volume of 25 million m³, Negreni reservoir (R4) is the second permanent reservoir downstream of Cal Alb reservoir. The nonpermanent reservoirs are R2, R3, and R5. These reservoirs have a dam height between 4 and 5 m and a storage water volume between 2.5 and 6 million m³ (Table 1).

Currently, the R1 reservoir and dam are under the management of the Prut–Bîrlad Water Basin Administration (PBWBA) [50,51].

In 1973, less than a year after the construction of the R1 reservoir, there were first infiltrations occurring through the dam foundation and embankment body with immediate effects downstream where the entire area became unstable due to the swamping process. According to the technical reports made at that time by the authorities, due to the presence of the aquifer layer at just 1 m below the dam foundation combined with the high permeability of the lithological deposits, which allowed a good circulation of groundwater, all of these affected the dam body integrity. However, the dam operated for over 35 years with these issues until the period between 2009 and 2013 when a dam stabilization project was carried out.

2. Materials and Methods

Figure 3 summarizes the workflow diagram followed in this study with LiDAR, R1 bathymetric model (Figure 3a), hydrological and built-up data obtaining processes (Figure 3b), and the key steps for HEC-RAS 2D stream flow simulation and flood modeling caused by a potential R1 dam break scenario (Figure 3c,d). The flood control capacity of the Başeu multi-reservoir system in case of R1 dam failure (Figure 3e) was evaluated on the basis of the flood wave impact (flood extent, flood depth, flood velocity) and flood hazard within 27 settlements located downstream of R1 (Figure 3f,g).

2.1. Data Acquisition

2.1.1. LiDAR Data

With the implementation of the 2007/60/EC Directive and SMIS-CSNR No.17945 project: Works to reduce the flood risk in the Prut–Bîrlad River Basin (PBRB) [51] in 2015, the northeastern territory of Romania was airborne scanned using LiDAR technology (Figure 3a). Within this project, the detailed flight level on the Başeu multi-reservoir system (study area) was classified as the first (A) level of importance—flights in the area adjacent to watercourses and urban areas. Accordingly, using the Leica ALS60 Airborne Laser Scanner (ALS) [52], a detailed elevation point cloud (.las files) of the topographic surface was achieved. After this stage, the obtained .las files were classified and vegetation areas were extracted [53], and the individual raster files (.tiff files) with an area of 0.25 km² and a spatial resolution of 0.5 m/pixel were generated. The LiDAR-derived DEMs obtained via the ArcGIS processing of more than 2000 tiles (.tiff files) were unified to a spatial resolution of 0.5 m/pixel. In this way, we obtained a high-accuracy DTM which included detailed topographic information related to hydro-technical structures (e.g., dams, reservoirs, polders, lateral structures) which considerably improved the 2D hydraulic model. Furthermore, to capture all the topographic details given by the presence of new constructions (post 2015) within the built-up area such as houses or other types of buildings (e.g., attachment buildings, administrative buildings, industrial buildings), they were integrated with the DTM for higher accuracy of the 2D hydraulic model [30]. The building database was obtained by on-screen digitizing techniques using orthophotos collected in 2018. To each building, an average height of 6 m was assigned, in order to rasterize and join to the final LiDAR-derived DTM [54].

2.1.2. Development of Bathymetric Model

In order to generate the 2D hydraulic models, we decided to integrate the R1 bathymetric model within the LiDAR-derived DTM. The bathymetric data for the R1 reservoir were obtained through the SMIS-CSNR No. 17945 project results [51]. The bathymetric measurements were made in the summer of 2015 using a motorboat equipped with a single-beam echo-sounder. More than 40,000 depth point networks were acquired in two different coordinate systems: WGS 84—Global Datum, EPSG: 4326 and Stereo 1970—National Datum, EPSG: 31700. After processing the depth points, a bathymetric model with a spatial resolution of 0.5 m/pixel was obtained [55,56]. Unfortunately, from all reservoirs within the study area, only the R1 reservoir has bathymetric measurements, whereas, for R2–R5 reservoirs, we used the normal water level extracted from LiDAR-derived DTM. However, for the final DTM used in the 2D hydraulic modeling of the R1 dam failure, the bathymetric data were integrated. In this context, a raster dataset was generated in order to join to the LiDAR-derived DTM and the R1 bathymetric model.

2.1.3. Hydrological Data

The Başeu river sector between the source and the confluence with the Cal Alb (R1) reservoir has a length of 17 km and ensures an average multiannual inflow rate of the 0.365 m³/s. In 2014, the official operating rules (OORs) for R1 reservoir were drawn up by PBWBA. According to the OORs of R1 reservoir, the inflow rate and the water volume for the three scenarios intervals (0.1% or 1000 year; 0.5% or 500 year; 1% or 100 year) were extracted (

Table 3. The Q_max (m³/s) inflow rates and volume of water contained (million m³) estimated by PBWBA for 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals for Cal Alb (R1) reservoir.

Recurrence Interval	1Qmax (m3/s)	Volume of Water Contained (Million m3)
1% (100-year)	185	10.19
0.5% (500-year)	225	12.39
0.1% (1000-year)	315	17.35

¹Q_max (m³/s)—maximum inflow rate (m³/s) estimated for 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals according to the hydrological data recorded between 1966 and 2019 at two gauging stations located downstream of Cal Alb (R1) reservoir.

). The recurrence intervals were calculated by PBWBA on the basis of hydrological data obtained from the two hydrometric stations located downstream of the R1 reservoir. The hydrological data series were between 31 years (1969–2000, Havârna gauging station) and 53 years (1966–2019 Ștefănești gauging station). The flow sizing rate Q2% (50-year recurrence interval) was 145 m³/s. The flow sizing rate Q0.5% + 20% (500-year recurrence interval + 20%) was 270 m³/s. Thus, to ensure the life of the aquatic ecosystem for R1 reservoir, an inflow rate of 0.015 m³/s is needed.

2.1.4. Land-Use Data

The land-use data were obtained via on-screen digitizing techniques on the basis of orthophotos collected in 2018 by the National Agency for Cadastre and Land Registration of Romania (NACLRR) (Figure 3b). In this way, we obtained 20 land-use categories: houses (0.52%), attachments buildings (0.36%), industrial buildings (0.11%), national road (0.11%), county road (0.11%), local road (0.84%), secondary streets (0.08%), exploitation roads (0.24%), yards (3.17%), orchard (0.08%), vineyard (0.57%), shrubbery (0.20%), forest vegetation (3.65%), grassland (28.7%), arable land (53.3%), unproductive land (0.10%), degraded land (0.11%), lakes and reservoirs (5.52%), streams (0.40%), and wetlands (1.65%). After the final validation, we obtained more than 30,000 polygons with a minimum mapping unit of 20 m² [30].

2.2. 2D HEC-RAS Modeling

The 2D hydraulic modeling requires an impressive and accurate volume of information regarding the spatial data (e.g., terrain model, bathymetric model, roughness coefficient) and hydrological data. According to these requirements, the HEC-RAS software [23,24,25] has the ability to simulate an inflow flood using three different methods: one-dimensional (1D) unsteady flow model [57], two-dimensional (2D) unsteady flow model [58], and level pool flow model [59]. Therefore, due to the high density of human settlements downstream of the R1 reservoir in conjunction with the availability of all necessary data (e.g., LiDAR-derived DTM, bathymetry data, hydrological data, and roughness information), we considered the 2D hydraulic modeling as the most appropriate method for testing the flood control capacity of Başeu multi-reservoir system. Furthermore, due to the high density of settlements (rural and urban areas) downstream of R1 reservoir, we chose to represent the built-up area by integrating each building in the LiDAR-derived DTM for a high accuracy of the 2D HEC-RAS hydraulic model. This process was detailed in Section 2.1.1 [51,54]. The advantage of using 2D modeling is that this method can capture the changes in water elevation when a flood wave arrives, as well as when the water elevation starts decreasing due to the dam break [60,61].

Overall, the HEC-RAS software can perform a 2D unsteady flow routing with two different equations: Equation (1), 2D full Saint Venant (also called full momentum equation) [23,25,31,62]; Equations (2) and (3), 2D Diffusion wave equation (which is set as the default in HEC-RAS software) [23,25]. In the process of choosing the proper equation for the breach scenario, two different hydraulic plans with Equation (1) and with Equations (2) and (3) were created, and the results were compared. If there are no significant differences between the two models, the user can proceed with the diffusion wave equation [23]. In our case, the two models were similar, and we decided to complete the analysis using the default equation (2D diffusion wave equation) which assures a greater stability property and a faster computational time.

$$\frac{\partial \zeta }{\partial t}+\frac{\partial p}{\partial x}+\frac{\partial q}{\partial y}=0,$$

(1)

$$\frac{\partial p}{\partial t}+\frac{\partial }{\partial x}\left(\frac{p^2}{h}\right)+\frac{\partial }{\partial y}\left(\frac{pq}{h}\right)=-\frac{n^{2}pg\sqrt{p^2+q^2}}{h^2}-gh\frac{\partial \zeta }{\partial x}+pf+\frac{\partial }{\rho \partial x}\left(h{\tau }_{xx}\right)+\frac{\partial }{\rho \partial y}\left(h{\tau }_{xy}\right),$$

(2)

$$\frac{\partial q}{\partial t}+\frac{\partial }{\partial y}\left(\frac{q^2}{h}\right)+\frac{\partial }{\partial x}\left(\frac{pq}{h}\right)=-\frac{n^{2}qg\sqrt{p^2+q^2}}{h^2}-gh\frac{\partial \zeta }{\partial y}+qf+\frac{\partial }{\rho \partial y}\left(h{\tau }_{yy}\right)+\frac{\partial }{\rho \partial x}\left(h{\tau }_{xy}\right),$$

(3)

where h is the water depth (m), p and q are the specific flow in the x- and y-directions (m²·s⁻¹), $\zeta$ is the surface elevation (m), g is the acceleration due to gravity (m·s⁻²), n is the Manning resistance, $\rho$ is the water density (kg·m⁻³), ${\tau }_{xx}$, ${\tau }_{yy}$, and ${\tau }_{xy}$ are the components of the effective shear stress, and f is the Coriolis (s⁻¹). When the diffusive wave is selected, the inertial terms of the momentum equations are neglected (Equations (2) and (3)).

2.3. Dam Break Scenario

On the basis of historical data obtained from the dam breaks around the world, scientists identified two main modes of dam failure: piping and overtopping [63,64,65]. Although information on dam breaks is numerous, no pattern has yet been identified regarding the shape of the breach, its location within the dam, or the time formation. To estimate breach characteristics, four methods can be used: (i) comparative analysis of the analyzed dam with another similar historical dam failure [23,66]; (ii) regression equations developed from historical dam failures [23,44,67]; (iii) using the velocity vs. erosion rates [23]; (iv) using computer models and principles of hydraulics [23,59,68,69,70,71]. In this case, we used regression equations to estimate the dam breach characteristics [44].

For this study, we assumed that a historical flood can occur in the river basin of the Cal Alb reservoir. The inflow for each scenario would create a dam breach due to the piping process. The flood wave created by the dam break would lead to the filling of the downstream reservoirs. Due to the location of the reservoirs, polders, and lateral structures downstream of the R1, the flood wave severity would be diminished. The flood wave would overlap with normal level water for each reservoir located downstream of the R1. In order to run the unsteady flow analysis with the breach scenario [15,23,25,31,72], six steps were followed: (a) the storage area extent was imported, (b) the dam characteristics were introduced in order to create the dam, (c) the breach parameters were calculated, (d) the subgrid model (mesh) was created, (e) the roughness data were introduced, and (f) the inflow data were set and the computational time was calculated [23,25].

The first step was to import the storage area extent in the HEC-RAS software. The following characteristics of the dam were introduced: dam length—295.07 m; dam crest width—5 m; dam height—14.5 m. Once the dam was created, the following parameters were introduced to breach the dam: center station, final bottom width, final bottom elevation, left- and right-side slope, breach weir coefficient, and piping coefficient. According to the equations developed by Froehlich [44], who studied 74 dam breaks, we estimated two parameters: average breach width using Equation (4) and breach formation time according to Equation (5). On the basis of the average breach width, average slide slope, and height of the final breach (Figure 3c), we calculated the breach bottom width using Equation (6).

$$B_{ave}=0.27K_0{V}_w^{0.32}{h}_b^{0.04},$$

(4)

$$t_f=63.2\sqrt{\frac{V_w}{g{h}_b^2}},$$

(5)

$$W_b=B_{ave}-\left(\frac{H}{V}{h}_b\right),$$

(6)

where B_ave is the average breach width (m), K₀ is constant (1.0 for piping failure), V_w is the reservoir volume at time of failure (m³), h_b is the height of the final breach (m), g is the gravitational acceleration (9.80665 m·s⁻²), t_f is the breach formation time (s), W_b is the final bottom width, H is the left slope side (0.7 for piping failure), and V is the right slope side (1.0 for piping failure). According to Equation (4), the average breach width is 55.6 m; according to Equation (5), the breach formation time is 1.5 h; according to Equation (6), the breach bottom width is 42.6 m. The connection between the storage area and the dam was realized in the moment when the break line for the dam was drawn [23,73].

The next step was to create the subgrid model (mesh). The subgrid models represent a polygonal network where each cell has information about the underlying terrain. In the 2D model, each computational cell is similar to each cross-section from the 1D model, and a property table is created for each cell [25,31]. The table contains information about the roughness, the elevation–volume relationship, area, maximum elevation, minimum elevation, average elevation, etc. On the basis of these data obtained from the underlying terrain, the flow between the cells is computed. To obtain an accurate representation of the flood wave, a proper cell size must be chosen [74,75]. On the one hand, a large cell may lead to loss of details about the flow moving; however, on the other hand if the cell size is too small, this could cause computational instability [23] (Figure 3d). For the current study, the 2D flow area was created on the basis of the floodplain extent with a surface of 160 km² and more than 2,500,000 computational cells with an average size of 64 m². The 2D flow area was directly connected with the dam. Regarding the Manning’s roughness values, i.e., the study area being a typical plain with a low urbanization process, the roughness of the land-use categories was as follows: degraded land—0.03; lakes and reservoirs—0.04; shrubbery—0.04; streams—0.04; wetlands—0.07; exploitation roads—0.08; railroad—0.08; forest vegetation—0.10; orchard—0.01; unproductive land—0.01; attachment buildings—0.015; county roads—0.015; houses—0.015; local roads—0.015; main streets—0.015; national roads—0.015; secondary streets—0.015; arable land—0.035; grassland—0.035; vineyard—0.035; yards—0.035 [23,76,77,78].

For the unsteady flow data, boundary conditions and initial conditions were set. As initial conditions, the elevation versus volume curve was used for the storage area as input. For each scenario, we assumed that the reservoir has a water volume of 17.35 million m³ for the 0.1% recurrence interval (which corresponds to 315 m³/s flow rate), 12.39 million m³ for the 0.5% recurrence interval (which corresponds to 225 m³/s flow rate), and 10.16 million m³ for the 1% recurrence interval (which corresponds to 185 m³/s flow rate) (Figure 3e,f). For each reservoir located downstream of the R1, the initial conditions were set as a partial fill with water because these reservoirs are used for fish farming. A percentage of the total storage capacity of 28.8% was used for R2, 63.3% for R3, 41.2% for R4, and 37% for R5. As boundary conditions, a normal depth (the energy slope of the riverbed) was set for the 2D flow area [31].

The last step within the dam break modeling process was to calculate the computational time to ensure the stability of the hydraulic model. The main effect of a wrong computational time is the instability of the model [23]. Explained in more detail, if the time step is too big, in addition to the instability of the model, attenuation of the peak can occur, and, if the time step is too small, the computational period will be very long. However, in order to estimate the computational time, the Courant condition was used (Equation (7)) [23,31].

$$C=\frac{V_w\Delta T}{\Delta X}\le 1, \Delta T=\frac{\Delta x}{V_w}, V_w=\frac{dQ}{dA},$$

(7)

where C is the Courant number, ∆T is the time step (s), ∆x is the distance step in m (average two-dimensional cell size), V_w is the flood wave speed (m/s), dQ is the change in discharge over a short time interval (Q₂ − Q₁), dA is the change in cross-section area over a short time interval (A₂ − A₁). According to Equation (7), a time step of 10 s was used to run the model.

2.4. Flood Hazard Assessment

Worldwide, flood hazard classification studies have been performed [31,45,46,62,79], where the most common hazard classification is based on the flood depth and flood extent [79]. In this study, we generated the flood hazard maps on the basis of 2D HEC-RAS modeling and Cal Alb embankment dam break multi-scenarios via the depth × velocity (D × V) raster exported from the RAS Mapper module. In order to generate the hazard severity classes, the flood D × V was classified according to AIDR [46] in six categories: H1 (D × V ≤ 0.3 m²/s), H2 (D × V range between >0.3 m²/s and ≤0.6 m²/s), H3 (D × V range between >0.6 m²/s and ≤1.2 m²/s), H4 (D × V range between >1.2 m²/s and ≤2 m²/s), H5 (D × V range between >2 m²/s and ≤4 m²/s), and H6 (D × V >4 m²/s) (

Table 4. Flood hazard classification based on the flood depth and flood velocity according to the AIDR [46].

Flood Hazard	D × V (m2/s)	Hazard Description
H1	≤0.3	Generally sage vehicles, people, and buildings
H2	≤0.6	Unsafe for small vehicles
H3	≤1.2	Unsafe for vehicles, children, and the elderly
H4	≤2	Unsafe for vehicles and people
H5	≤4	Unsafe for vehicles and people. All the building types vulnerable to structural damage. Some less robust building types vulnerable to failure.
H6	>4	Unsafe for vehicles and people. All building types considered vulnerable to failure.

) [46].

3. Results

3.1. Flood Pattern

On the basis of the 24 h breach simulation period, using the RAS Mapper interface [25], we managed to export, as raster and vector, two types of data for each break scenario (Figure 3g). To analyze the impact of the flood created by the hypothetical dam breach, depth and velocity (raster files) were needed. To analyze the extent of the flood, we exported the inundation boundary (vector file) from RAS Mapper [25,80].

3.1.1. Flood Extent

The hypothetical dam break in the case of the three simulated scenarios can cause the flooding of 1191 ha (20.1 ha within built-up area) in the case of the 1% (100-year) scenario (Figure S1, Supplementary Materials), 1229 ha (21.74 ha within built-up area) in the case of the 0.5% (500-year) scenario (Figure S2, Supplementary Materials), and 2850 ha (87.2 ha within built-up area) in the case of the 0.1% (1000-year) scenario (Figure 4, Figure S3, Supplementary Materials). Considering the flood extent in each dam break scenario, we can say that the multi-reservoir system proved to be effective in diminishing the flood extent in the case of two scenarios (1% and 0.5%). In the case of the third scenario (0.1%), the multi-reservoir system did not cope with the flow rate discharged upstream, and the flood extent was over 58.2% (with 1659 ha more) higher than the 1% scenario and 56.8% (with 1621 ha more) higher than the 0.5% scenario. Regarding the built-up area affected by floods, the multi-reservoir system led to the protection of 22 settlements (of which two are urban settlements) in the case of the 0.5% and 1% scenarios and just four rural settlements in the case of the 0.1% scenario.

For the first scenario (1%), the built-up area of five settlements is potentially affected. A total surface of 20.1 ha is potentially affected. In terms of land-use categories, the most potentially affected categories (within the built-up area) are the arable land (12.8 ha), wetlands (0.75 ha), grassland (1.99 ha), forest vegetation (1.78 ha, yards (1.58 ha), and buildings (0.45 ha—25 houses and 28 attachment buildings). Other categories have a potentially affected surface below 0.2 ha (e.g., streams, local roads, orchards, national roads, degraded land, shrubbery, vineyard). Referring to the affected area of the settlements, Havârna (S1) is the most affected settlement with 10.44 ha (

Table 5. Flood extent within built-up area of 27 settlements located downstream of the Cal Alb reservoir and number of inhabitants potentially affected by flood waves computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

Settlement Code	Built-Up Area (ha)	Number of Inhabitants	1% (100-year)		0.5% (500-year)		0.1% (1000-year)
			1 ha	2 inhab.	1 ha	2 inhab.	1 ha	2 inhab.
(S1) Havârna	418	2823	10.44	25–50	11.20	50–75	25.78	150–200
(S2) Gârbeni	64	420	2.67	5–25	3.33	25–50	15.4	150–200
(S3) Tătărășeni	114	797	0.44	<5	0.54	5–25	5.6	25–50
(S4) Balinți	67	396	5.79	25–50	5.87	25–50	9.45	50–75
(S5) Niculcea	13	41	-	-	-	-	0.08	<5
(S6) Negreni	109	909	0.79	5–25	0.80	5–25	1.28	5–25
(S7) Știubieni	214	1727	-	-	-	-	4.72	25–50
(S8) Chișcăreni	61	495	-	-	-	-	1.75	5–25
(S9) Sat Nou	28	124	-	-	-	-	1.18	5–25
(S10) Petricani	97	645	-	-	-	-	3	5–25
(S11) Săveni	229	8145	-	-	-	-	0.41	5–25
(S12) Vlăsinești	156	1596	-	-	-	-	0.84	5–25
(S13) Bozieni	34	255	-	-	-	-	2.1	5–25
(S14) Sârbi	104	1064	-	-	-	-	0.7	5–25
(S15) Miron C.	83	472	-	-	-	-	1.63	5–25
(S16) Slobozia	18	86	-	-	-	-	-	-
(S17) Hănești	128	1127	-	-	-	-	1.4	5–25
(S18) Moara J.	22	89	-	-	-	-	0.01	<5
(S19) Mihălășeni	107	743	-	-	-	-	-	-
(S20) Negrești	32	300	-	-	-	-	-	-
(S21) Păun	41	289	-	-	-	-	0.76	5–25
(S22) Năstase	31	183	-	-	-	-	0.01	<5
(S23) Stânca	105	812	-	-	-	-	0.71	5–25
(S24) Ștefănești	242	6630	-	-	-	-	4.1	75–100
(S25) Bădiuți	108	985	-	-	-	-	1.43	5–25
(S26) Bobulești	206	1322	-	-	-	-	1.92	5–25
(S27) Românești	183	1394	-	-	-	-	-	-

¹ S: Built-up area (ha) potentially affected by the flood waves computed for each dam breaks scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals; ² inhab.: interval number of inhabitants potentially affected by the flood waves computed for each dam breaks scenarios with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

), which is also the first settlement situated downstream of Cal Alb reservoir. A number of 25–50 inhabitants can be potentially affected by the flood wave. The second most affected settlement is Balinți (S4) with 5.79 ha and 25–50 inhabitants potentially affected. The least affected settlement is Tătărășeni (S3) with 0.44 ha and fewer than five inhabitants potentially affected. Other affected built-up areas are Gârbeni (S2–2.67 ha) and Negreni (S6–0.79 ha). The potential number of affected inhabitants for these two settlements is between 5 and 50.

According to the second scenario (0.5%), a built-up area of 21.74 ha (five settlements) is potentially affected and the result is similar to the first scenario (0.1%). The most affected land-use categories (within the built-up area) are arable land (13.6 ha), grassland (2.29 ha), forest vegetation (1.9 ha), yards (1.8 ha), wetlands (0.78 ha), and buildings (0.5 ha—28 houses and 32 attachment buildings). Other categories have a potentially affected surface below 0.2 ha (streams, lakes and reservoirs, local roads, county roads, orchards, degraded land, shrubbery, vineyard). Havârna (S1) remains the most affected settlement with 11.2 ha and 50–75 inhabitants potentially affected (

Table 6. Potentially affected land-use categories within built-up area by flood waves computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

Land-Use Category	Recurrence Interval
	1% (100-year)	0.5% (500-year)	0.1% (1000-year)
Houses (no.)	25	28	179
Attachments buildings (no.)	28	32	192
Industrial buildings (no.)	-	-	2
Buildings (ha)	0.45	0.53	3.45
Streams (ha)	0.01	0.01	0.08
Lakes and reservoirs (ha)	0.06	0.06	0.12
Yards (ha)	1.588	1.84	10.55
Exploitation roads (ha)	-	-	0.02
Local road (ha)	0.27	0.32	2.28
County road (ha)	0.01	0.01	0.04
National Road (ha)	-	-	0.001
Orchard (ha)	0.01	0.01	0.26
Grassland (ha)	1.9	2.29	12.31
Forest vegetation (ha)	1.78	1.9	7.85
Secondary streets (ha)	-	-	0.01
Arable land (ha)	12.8	13.6	44.64
Unproductive land (ha)	-	-	0.2
Degraded land (ha)	0.09	0.09	0.17
Shrubbery (ha)	0.27	0.27	0.36
Vineyard (ha)	0.01	0.01	0.05
Wetlands (ha)	0.75	0.78	1.82

). Balinți (S4) is the second most affected built-up area with a number of 25–50 inhabitants potentially affected. The least most affected settlements are Gârbeni (S2–3.33 ha), Tătărășeni (S3–0.54 ha), and Negreni (S6–0.80 ha), and all three settlements have a potentially number of affected inhabitants of 5–25.

In the third scenario (0.1%), a built-up area of 87.2 ha is potentially affected. A total of 23 settlements are at risk to be affected by the potential flood wave (Table 6). The most affected land-use categories (within the built-up area) over 10 ha are arable land (44.6 ha), grassland (12.3 ha), and yards (10.5 ha). Land-use categories with a potential affected surface between 1 and 10 ha are forest vegetation (7.8 ha), buildings (3.4 ha—179 houses, 194 attachment buildings, two industrial buildings), local roads (2.3 ha), and wetlands (1.82 ha). Land-use categories with above 1 ha of potentially affected surface are streams, lakes and reservoirs, exploitation roads, county roads, national roads, orchards, secondary streets, unproductive land, degraded land, shrubbery, and vineyard (Table 6). The most potentially affected settlements remain those located in the vicinity of the dam. Havârna (S1—25.78 ha) and Balinți (S3—15.4 ha) settlements are the most affected with a potential number of affected inhabitants of 150–200. Other settlements that have a significant surface with high potential to be affected are Tătărășeni (S3—5.6 ha), Balinți (S4—9.45 ha), Știubieni (S7—4.72 ha), and Ștefănești (S24—4.1), with a number of potentially affected inhabitants between 5 and 100. The least affected settlements (less than 1 ha) are Niculcea (S5—0.08 ha), Sârbi (S14—0.7 ha), Moara Jorii (S18—0.01 ha), Păun (S21—0.76 ha), Năstase (S22—0.01 ha), and Stânca ((S23—0.71 ha), with a number of potentially affected inhabitants between 1 and 25. Due to the topographical location on the first terrasse of Bașeu River and the high flood control capacity of Hănești reservoir, four settlements (Slobozia—S16, Mihălășeni—S19, Negrești—S20, and Românești—S27) are fully protected from the flood wave. According to all three flood scenarios developed in the R1 dam break scenario, the first four settlements (S1—Havârna; S2—Gârbeni; S3—Tătărășeni; S4—Balinți) have the highest potential of being affected by the flood wave. The high flood control capacity of the multi-reservoir system is proven by the small number of settlements (four from 27 settlements) that are subject to a potential flood event.

3.1.2. Flood Depth

The second parameter that was calculated and generated by HEC-RAS software was flood depth, which was calculated for each cell on the basis of the underlying terrain, taking into consideration the maximum flood depth of the computational cell no matter the time when the maximum values were registered [25]. The maximum flood depth for the first scenario is 7.4 m and corresponds to the main channel of the Bașeu River (Figure S4, Supplementary Materials). In terms of affected buildings, 75.5% (40 buildings) of them can be potentially affected by a flood depth that does not exceed 1 m, and 24.5% (13 buildings) are potentially affected by a flood depth between 1 and 2 m (

Table 7. Number of buildings with a high potential to be affected and the flood wave depth computed for each dam break scenarios with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

Flood Depth (m)	Recurrence Interval
	1% (100-year)		0.5% (500-year)		0.1% (1000-year)
	1 build	2 rel.freq.	1 build	2 rel.freq.	1 build	2 rel.freq.
<1	40	75.5	43	71.7	175	46.9
1–2	13	24.5	15	25.0	85	22.8
2–3	-	-	2	3.3	56	15.0
3–4	-	-	-	-	35	9.4
4–5	-	-	-	-	10	2.7
>5	-	-	-	-	12	3.2

¹ build: the number of potentially affected buildings by the flood depth waves computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals; ² rel.freq.: the relative frequency of the potentially affected buildings by the flood depth waves computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

). In the case of the second scenario, the maximum depth is 7.6 m (Figure S5, Supplementary Materials). The result is very similar to the first scenario, even in flood depth classes and affected buildings. A percentage of 71.7% of the buildings are potentially affected by a flood depth that does not exceed 1 m, 25% by a flood depth between 1 and 2 m, and, unlike the first scenario, two buildings (3.3%) are affected by a flood depth between 2 and 3 m (Table 7). In the case of the third scenario, a maximum depth of 9.95 m (Figure 5, Figure S6, Supplementary Materials) was registered, and 46.9% (175 buildings) of buildings are potentially affected by a flood depth that does not exceed 1 m, 22.8% (85 buildings) by a flood depth between 1 and 2 m, 15% (56 buildings) by a flood depth between 2 and 3 m, 9.4% (35 buildings) by a flood depth between 3 and 4 m, 2.7% (10 buildings) by a flood depth between 4 and 5 m, and only 3.2% (12 buildings) are potentially affected by a flood that exceed 5 m (S1–S5 settlements) (Table 7).

3.1.3. Flood Velocity

The velocity parameter was calculated in the same way as the flood depth [25]. For each computational cell, the maximum value of velocity from the 24 h simulation time was registered. For all the computed scenarios, the flood velocity exceeds 7 m/s: 7.3 m/s in the case of the first and second scenario and 9.1 m/s in the case of the third scenario (Figure 6). The high values of the flood velocity can be attributed to the geomorphologic and vegetation conditions which lead to low roughness values, accentuated by the narrow areas of the drainage channel. In the case of the first scenario, 88.7% (47 buildings) are potentially affected by a flood velocity that does not exceed 1 m/s, and 11.3 (six buildings) are potentially affected by a flood velocity between 1 and 2 m/s (

Table 8. Number of buildings with a high potential to be affected and the flood wave velocity computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

Flood Velocity (m/s)	Recurrence Interval
	1% (100-year)		0.5% (500-year)		0.1% (1000-year)
	1 build	2 rel.freq.	1 build	2 rel.freq.	1 build	2 rel.freq.
<1	47	88.7	52	86.7	257	68.9
1–2	6	11.3	8	13.3	105	28.2
2–3	-	-			11	2.9
3–4	-	-	-	-	-	-
4–5	-	-	-	-	-	-
>5	-	-	-	-	-	-

¹ build: the number of potentially affected buildings by the flood depth waves computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals; ² rel.freq.: the relative frequency of the potentially affected buildings by the flood depth waves computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

, Figure S7, Supplementary Materials). In the case of the second scenario, 86.7% (52 buildings) are potentially affected by a flood velocity that does exceed 1 m/s, and 13.3% (eight buildings) are potentially affected by a flood velocity between 1 and 2 m/s (Table 8, Figure S8, Supplementary Materials). In the case of the third scenario, 69.9% of the total buildings (257 buildings) are potentially affected by a flood velocity that does not exceed 1 m/s, 28.2% (105 buildings) are potentially affected by a flood velocity between 1 and 2 m/s, and only 2.9% (11 buildings) are potentially affected by a flood velocity between 2 and 3 m/s (Table 8, Figure S9, Supplementary Materials).

3.2. Flood Hazard Assessment

According to the 1% (100-year) scenario, the maximum value of D × V registered along the Başeu floodplain was 16.8 m²/s (Figure S10, Supplementary Materials). In the built-up area, 53 buildings are potentially affected by floods, and the D × V values do not exceed 4 m²/s, of which 28.3% (15 buildings) are situated in the H1 class (0.3 m²/s), 20.8% (11 buildings) are in the H2 class (≤0.6 m²/s), 26.4% (14 buildings) in the H3 (1.2 m²/s), 9.4% (five buildings) in the H4 (≤2 m²/s), and 15.1% (eight buildings) in the H5 (≤4 m²/s) (

Table 9. Number of potentially affected buildings by flood vs. D × V hazard classes computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

Flood Hazard Classes	D × V (m2/s)	Recurrence Interval
		1% (100-year)		0.5% (500-year)		0.1% (1000-year)
		1 build	2 rel.freq.	1 build	2 rel.freq.	1 build	2 rel.freq.
H1	≤ 0.3	15	28.3	31	51.7	125	33.5
H2	≤ 0.6	11	20.8	11	18.3	56	15.0
H3	≤ 1.2	14	26.4	8	13.3	55	14.7
H4	≤ 2	5	9.4	4	6.7	48	12.9
H5	≤ 4	8	15.1	6	10.0	67	18.0
H6	> 4	-	-	-	-	22	5.9

¹ build: the number of potentially affected buildings by flood based on D × V hazard classes computed for each dam break scenario with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals; ² rel.freq.: the relative frequency of the potentially affected buildings by flood based on D × V hazard classes computed for each dam break scenarios with 1% (100-year), 0.5% (500-year), and 0.1% (1000-year) recurrence intervals.

). In the case of 0.5% (500-year) scenario, the maximum value of D × V was 16.9 m²/s, similar to the maximum value of D × V from the 1% (100 year) scenario (16.8 m²/s) (Figure S11, Supplementary Materials). However, the differences between the 0.5% (500-year) and 1% (100-year) scenarios can be observed on the affected buildings and hazard classes. Therefore, more than 50% (31 buildings) of the total potentially affected buildings (60 buildings) identified according to the 0.5% (500-year) scenario are situated in the H1 class of hazard (≤0.3 m²/s), 18.3% (11 buildings) in the H2 class (0.6 m²/s), 13.3% (eight buildings) in the H3 class (≤1.2 m²/s), 6.7% (four buildings) in the H4 class (≤2 m²/s), and the remaining 10% (six buildings) in the H5 class (≤4 m²/s) (Table 9). According to the 0.1% (1000-year) scenario, the maximum value of D × V was 21.6 m²/s (Figure 7, Figure S12, Supplementary Materials). Related to the distribution of potentially affected buildings on flood hazard classes, 33.5% (125 buildings) of the total potentially affected buildings (373 buildings) according to the 0.1% (1000-year) scenario are situated in the H1(≤0.3 m²/s) class of hazard, 15% (56 buildings) in the H2 (≤0.6 m²/s), 14.7% (55 buildings) in the H3 (≤1.2 m²/s), 12.9% (48 buildings) in the H4 (≤2 m²/s), 18% (67 buildings) in the H5 (≤4 m²/s), and the remaining 5.9% (22 buildings) in the H6 (>4 m²/s) (Table 9).

4. Discussion

Permanent and nonpermanent reservoirs have played an important role in protecting human settlements which are located in floodplains, especially rural settlements where subsistence farming is still practiced [47,48,49]. The intensification of extreme hydrological events (e.g., floods, flash floods, overtopping) and uncontrolled urbanization within flood hazard areas led to the point where the reservoirs became a key factor in defense against flood hazard [31]. The newly identified purpose of the reservoirs, in addition to those attributed to agricultural, social, and economic activities, makes reservoirs important strategic structures. In addition to the many benefits that dams bring to human communities located in the proximity of the reservoirs, their management also poses the risk that these water storages will yield to certain special conditions. The failure of an embankment dam, depending on its characteristics (e.g., surface area of the reservoir, water storage volume, dam characteristics) can generate flood waves that are significantly higher in terms of water volume release and water velocity than floods waves triggered by natural conditions [4] and, consequently, they can lead to more significant economic damages and human loses [4,5]. However, the flood and flood wave generated by the dam failure are categorized as accidental floods [2,5].

Therefore, complex and in-depth studies on accidental floods, especially those caused by embankment dam failure, are now mandatory due to the dynamics of hydroclimatic conditions and the destructive impact that such events can have on areas downstream of the reservoir. In this context, the awareness and perception of the risk that these hydrological phenomena induce to human society have led to the adoption of strategies to reduce and prevent the negative effects that may occur through their manifestation [48]. In Romania, these strategies should be implemented at national, regional, and even local level for the most important river basins and reservoirs, but they are not implemented in all the cases (e.g., R1—Cal Alb embankment dam and reservoir).

The present study complemented the national and local flood defense plans, given that the two permanent reservoirs on the Bașeu River (R1—Cal Alb reservoir and R4—Negreni reservoir) did not have an APDF. At the same time, the importance of conducting such a study is given by the fact that, at national and regional level, no unanimously accepted methodology is developed to assess the risk of accidental floods caused by the failure of a dam. Also, it is the first study of this type applied to the Başeu river basin, i.e., the R1–R5 multi-reservoir system.

For high accuracy of the hydraulic model, a digital surface model (DSM) would have been needed. The advantage of using a DSM is that the houses, attachment buildings, road network, or other artificial structures are integrated in the hydraulic modeling. The presence of these solid objects in the DSM can cause a different flood behavior, as well as change the flood parameters (flood extent, flood depth, flood velocity) and the number of the potential affected buildings. Due to the unavailability of this dataset, the buildings and the road network were used to increase the accuracy level of the DTM. The availability of a LiDAR model with a resolution of 0.5 m made it possible to create a complex hydraulic model that simulated the maximum flood wave and its overlap with the local topographic (e.g., old drain channels, abandoned meanders). On the basis of LiDAR-derived DTM (0.5 m), a bathymetrical model (0.5 m), hydrological data (correlation between R1 inflow with water volume and water level for three different recurrence intervals), and built-up area and land use data from Bașeu floodplain, we developed and applied a methodology which can be easily adapted also by regional authorities (e.g., PBWBA) and used within flood mitigation strategy.

In order to create the hydraulic model, a polygonal network with a computation point spacing of 7.5 m was created. For each computational cell, the LiDAR terrain characteristics were captured in the hydraulic model and a detailed elevation–volume relationship was developed by the HEC-RAS software. According to the underlying terrain process used in the HEC-RAS modeling, a large computation cell can be used without decreasing the accuracy of the model [25]. An important role in hydraulic modeling is played by the roughness coefficient. The flood parameters, especially flood velocity, depend on the roughness coefficient. Between roughness coefficient and flood velocity exists an inversely proportional relationship. For example, an increase in the roughness coefficient can cause a decrease in the flood velocity. Thus, for this study, we updated and optimized the roughness coefficient according to the spatial distribution of land-use categories, which is imperative and recommended for high-accuracy hydraulic modeling [30,31]. Moreover, by creating a complex hydraulic model and simulating a piping breach in the body of the R1—Cal Alb reservoir, the flood control capacity of the R1–R5 multi-reservoir system was assessed. Water depth levels, flood wave velocities, and flood extent were extracted and calculated on the basis of the results obtained from the RAS Mapper module. In the future, we want to redo the 0.1% (1000 year) recurrence interval scenario, while replacing the LiDAR DTM with LiDAR DSM in order to simulate the presence of artificial structures (e.g., roads, houses, attachment buildings).

The efficiency of the R1–R5 multi-reservoir system in the present study was also proven by the study [30] where the flood risk on the Bașeu River was analyzed by creating a 1D steady flow hydraulic model. A part of the previous study was a scenario generated by the introduction of hydrotechnical works in the 1D hydraulic model, but without the dam failure. Comparing the results of the two studies, the neighboring settlements (S1) Havârna, (S2) Gârbeni, (S3) Tătărășeni, (S4) Balinți, and (S6) Negreni of the R1—Cal Alb reservoir, R2—Movileni reservoir, R3—Tătărăşeni reservoir, and R4—Negreni reservoir showed the most important damages regarding the potentially flooded surface. Starting with R4—Negreni reservoir, the flood wave generated by the failure of the Cal Alb dam is diminished, and the drainage is made through the riverbed and floodplain. No built-up area represents a potentially flooded area. Downstream of R5—Hăneşti reservoir, the drainage is made exclusively through the riverbed, only in the area of (S23) Stânca, (S24) Ștefănești, (S25) Bădiuți, (S26) Bobulești, and (S27) Românești settlements (where the altitude drops to 84 m), and the flood wave covers the floodplain to a width of about 1–2 km. None of the settlements at the confluence with the Prut River (S23—Stânca, S24—Ștefănești, S25—Bădiuți, S26—Bobulești, and S27—Românești) present potentially flooded areas.

5. Conclusions

The 2D hydraulic modeling using HEC-RAS (v 5.0.7) and high-density LiDAR data (0.5 m spatial resolution) produced sufficiently accurate information regarding the flood control capacity of the Başeu multi-reservoir system (R1—Cal Alb reservoir; R2—Movileni reservoir; R3—Tătărăşeni reservoir; R4—Negreni reservoir; R5—Hăneşti reservoir) in case of R1 dam failure. On the basis of the three 2D stream flow hydraulic modeling configurations using the R1 inflow rate with 1% (100 year), 0.5% (500 year), and 0.1% (1000 year) recurrence intervals and the water volume which can be accumulated with that specific inflow rate (1% = 10.19 × 10⁶ m³; 0.5% = 12.39 × 10⁶ m³; 0.1% = 17.35 × 10⁶ m³), the following conclusions can be summarized:

• Combining 2D hydraulic modeling of the R1 dam break scenario (piping failure) with high-density LiDAR data (0.5 m spatial resolution) and local hydrological parameters (correlation between R1 inflow (m³/s) with water volume (m³) and water level (m)) proved to be an efficient method to improve the action plan for dam failure (APDF) and flood mitigation strategy within the Başeu multi-reservoir system.
• The multi-scenario approach using the inflow rate with 1% (100 year), 0.5% (500 year), and 0.1% (1000 year) recurrence intervals allowed the testing of the flood control capacity of Başeu multi-reservoir system according to R1 water volume which can cause a flood event in the case of dam failure. Accordingly, the first two scenarios (100 year and 500 year) indicate that only first four settlements (S1—Havârna, S2—Gârbeni, S3—Tătărășeni, S4—Balinți) located downstream of R1 are potentially affected by floods due to the location between R1 and R2 reservoirs, and only in case of 1000-year scenario are all 27 settlements potentially affected.
• The 2D hydraulic models were exported into a set of flood hazard parameters (e.g., flood extent, flood depth, flood velocity) and can be used for improve the flood hazard maps and answer real questions regarding the flood hazard threat at the local level in case of a dam failure.