Flood simulation to determine flood hazard susceptibility of downstream Singkil watershed in Aceh Province

. The damage to the Singkil watershed, one of the largest watersheds in Aceh Province, has become a concern for many parties, primarily due to increased flood events in the downstream area of the watershed. Based on previous studies, the problems which cause flooding are watershed damage due to illegal logging and high rainfall intensity reaching 3000-4500 mm/year, and increased erosion of 0.887 tons/ha/year, which causes a decrease in river capacity due to sedimentation. For this reason, this study aims to do 2D flood hydraulic modeling for five return periods of 2, 5, 10, 25, and 50 years which can be a reference for flood management in Aceh Singkil District. The hydrological analysis of the design flood discharge for several return period was carried out using Nakayasu and SCS-CN method. 1D-2D flood simulation generated using GeoHECRAS software to study flood hazard characteristic at downstream by combining tidal effects. The flood event in September 2012, a 10-year return period flood, became the basis for validating and calibrating the simulation model. Based on flood simulation result, it is known the downstream area of the Singkil watershed, Aceh Singkil districts is an area with a high flood potential and frequently occurs.

Floods are frequent and destructive natural disasters, posing significant threats to lives and economies [1], [2].Data from the Indonesian Disaster Information List (DIBI) reveals that flooding accounts for 33.53% of all recorded disasters in Indonesia until June 2023 (Source: https://dibi.bnpb.go.id/kbencana2).Floods typically fall into two categories: natural and human-induced.Natural factors include heavy rainfall, temperature changes, dam breaches, snowmelt, tides, and obstructions to water flow [3].Human activities like improper land use, deforestation, waste disposal in rivers, and construction in flood-prone zones also contribute to flooding.Mismanagement in upstream areas can harm downstream regions in watersheds.Increased resource exploitation, such as logging for plantations and mining in upstream forests, accelerates watershed degradation [4].
Flood control is a complex endeavour involving various engineering disciplines, including hydrology, hydraulics, river engineering, river morphology, sedimentation, flood control system engineering, urban drainage systems, water structures, and more [5].To effectively address flooding issues, especially in areas prone to annual flooding, flood analysis is an absolute necessity.It serves as the foundation for planning flood management and mitigation activities.Precisely determining and mapping flood hazard zones, along with implementing appropriate mitigation measures, can substantially minimize flood damage [6]- [8].Furthermore, flood hazard mapping plays a vital role in land-use planning, early warning systems, emergency response design, and the implementation of measures to reduce flood risks [7], [9].
To effectively manage flood occurrences and their associated risks, comprehensive flood control system planning is crucial.This planning process involves a thorough evaluation, including assessing the extent of flood inundation within a specific area [5].It requires detailed mapping of flood profiles and characteristics, such as inundation area, depth, and water flow velocity, all of which are essential for flood control planning.Increased housing demand due to population growth and climate change-induced surface runoff has led to more frequent floods.Consequently, flood modeling has become vital for future flood management, emergency response planning, and spatial planning efforts [10]- [12].For precise micro-scale assessments, 1D/2D E3S Web of Conferences 447, 01007 (2023) https://doi.org/10.1051/e3sconf/202344701007The 15 th  AIWEST-DR 2023 hydrodynamic models are commonly used [11].Over the past decade, the integration of Geographic Information Systems (GIS) and hydraulic modeling systems has gained popularity in floodplain modeling.The model developed by the Centre for Hydraulic Engineering is widely employed for studying inundation and mapping flood-prone areas [13].
This paper centers on dynamic flood modeling using GeoHECRAS 3.00.078 to map flood inundation.The study zone, the downstream area of Aceh Province's Singkil watershed, endures recurrent flooding due to environmental degradation like logging, illegal mining, and non-conservation-oriented plantations.These activities have escalated erosion and sedimentation, diminishing the watershed's capacity to handle water flow.With an average potential water discharge of 17-18 lt/dt/km 2 and annual rainfall reaching 3,000-4,500 mm/year, the region experiences frequent and severe flooding, with inundation heights of 2-3 meters [14].Indonesian Disaster Information Data (DIBI) recorded 25 floods from 2009-2021, indicating over one flood annually.Eight out of 11 sub-districts, including Gunung Meriah, Simpang Kanan, North Singkil, Suro, Singkohor, Kota Baharu, Singkil, and Danau Paris, are consistently affected.In September 2012, a 7-day flood incurred losses of 27 billion Indonesian Rupiah (Head of BPBD Aceh Singkil District, 2012).These frequent floods caused extensive damage to infrastructure and agricultural lands, disrupting daily life.This research aims to deliver valuable flood simulation outcomes to aid stakeholders in managing and mitigating flood risks in the downstream Singkil watershed area.

Flood mapping
Flood assessment and mapping can be approached through four methods: detailed studies, limited detailed studies, approximate studies, and redelineation [14].In detailed assessments, hydrological analysis relies on river measurement data or rainfall-runoff models, while hydraulic analysis involves comprehensive flow modeling with detailed river survey data.Limited detail assessments use river measurement data for hydrological analysis and steady flow modeling for hydraulics without the need for detailed river data.Approximate studies estimate inundation boundaries through simpler techniques like topographic map reading and field surveys.The redelineation method, the most straightforward, recreates flood maps by overlaying historical flood elevations onto new topographic maps.
To collect comprehensive flood inundation data, including discharge, flow velocity, inundation depth, and extent, hydrodynamic modeling and simulation are essential [15].In a detailed flood mapping study, both hydrological and hydraulic analyses are conducted sequentially.Initially, hydrological analyses are performed, utilizing hydrometeorological data like rainfall and river flow parameters (velocity, depth, discharge) as inputs.These inputs undergo processing using methods such as rainfall-runoff models, frequency analysis, or regression to generate flood hydrographs and peak discharges (Qp) [16].The output from the hydrological analysis serves as the basis for establishing boundary conditions in the hydraulic analysis.The hydraulic analysis employs a hydrodynamics model to calculate floodwater elevation.Incorporating topographic data in the form of a digital elevation model (DEM) and a base map enables accurate simulation of flood inundation characteristics, resulting in a comprehensive flood map.

Flood dynamic simulation modeling
A dynamic flood model simulates floodwater movement through various elements, including waterways, storage components, and hydraulic structures.These models calculate flood levels, flow patterns, and account for factors like backwater effects, levee overtopping, confluence at waterways, and bridge behaviour [17].Representing river and floodplain topography is crucial in hydraulic flood modeling.Flood inundation can be simulated using one-dimensional (1D), two-dimensional (2D), or three-dimensional (3D) spatial flow models.For urban areas, 2D or hybrid 1D/2D models offer more realistic calculations and visualizations compared to 3D models, which are rarely used for large-scale flood modelling [18].Numerous 2D flow hydrodynamics models, like Hecras 5 developed by The Centre for Hydraulic Engineering, have gained recognition from FEMA [19].In 2020, GeoHECRAS software was introduced, streamlining the process into a single interface combining GIS and HEC-RAS for efficiency.It simplifies validation, error reduction in floodplain mapping, and intuitive graphic interactions.GeoHECRAS builds upon HEC-RAS 6 but is underutilized, with limited studies conducted [20], [21].
GeoHECRAS 2D modeling is hydraulic flow modelling with a flow that has two directions, the first is the flow that leads from upstream to downstream of the river and the second direction is the flow that leads outside the river flow (inundation area) [22].
HECRAS 2D modeling can model variability, particularly along the river channel, where the model area is divided into topography-based mesh-shaped cells with the river geometry represented by break lines with smaller cell values, resulting in increased computation time [23].Flow is then controlled by Manning's number based on the land cover used [24].This advanced model facilitates 2D dynamic flow flood inundation simulations uses a shallow water equation approach that describe water flow in terms of the average depth of 2dimensional flow velocity and water level in response to the effects of gravity and friction.Hydraulic simulation in GeoHECRAS uses the conservation of mass and momentum field approach.The finite volume method is applied in the programme because it is conservative, uncomplicated, and the geometry is easy to change or flexible [22].Historical inundation height data, topography and flood inundation maps are used to analyse the inundation that has occurred so that the modeling is able to represent the situation in the field [25].
In addition to obtaining flood depths, dynamic flood simulations are also carried out to obtain flow velocities that are useful in flood management [26].Flow velocity is usually divided into 5 based on the category of water velocity or flow, namely very slow flow <10 cm/sec, slow flow 10-25 cm/sec, moderately fast river flow 25-50 cm/sec, fast river flow 50-100 cm/sec, and very fast river flow >100 cm/sec [27].

Study area
The research location taken place at lower areas of DAS Singkil which located at Aceh Singkil District, Aceh Province, Figure 1.The main river is Alas-Singkil river which connect into the Indonesian ocean so its flow is influenced by the tides.

2D flood hydraulics simulation
The 2D modeling in this research is 2 Dimensional Horizontal (2DH) hydraulic flow modelling will be carried out using GeoHECRAS software which will be Fig. 1.Research Location calibrated and validated based on flood events that have occurred.The flood that will be the reference is the flood with the most recent event, namely the September flood in 2012.The initial modelling schematisation that has been calibrated and validated will become the modelling schematisation to obtain flood inundation areas with return periods of 2, 5, 10, 20 and 50 years.The data used in the development of 2D Flood hydraulics simulation using GeoHECRAS software consists of secondary and primary data presented in Table 1.The series of analyses conducted in this study are sequentially detailed as follows:

A. Catchment area delineation
Watershed boundaries are determined using Digital Elevation Model (DEM) satellite data downloaded from http://tides.big.go.id/DEMNAS/.Catchment boundaries are established with the hydrology tool in ArcGIS 10.5 software.This catchment area serves as the basis for flow hydrograph calculations, which, in turn, are used as input for the 2D hydraulic flood modeling for upstream boundary conditions.

B. Calculation of design rainfall
The calculation conducted through frequency analysis to determine T-year return period rainfall using normal, log-normal, log-Pearson III, and Gumbel distribution methods.The suitability tests used are the chi-square and Smirnov Kolmogorov methods.

C. Flow hydrograph analysis
The input discharge in the 2D hydraulic flood simulation in this study is the flow hydrograph from the calculation of surface runoff using the Nakayasu method for main river and SCS-CN for tributaries.The initial step in calculation the volume of surface runoff is determined by identifying the observed flood rainfall, calculating the return period rainfall based on the hydrological method of rainfall data.
The step for Nakayasu method is analysing the characteristics as follow : a.The time from the onset of rain to the peak of the hydrograph (time to peak magnitude).b.The time from the rainfall centre to the hydrograph centre (time log).c.Time base of hydrograph.d.Area of watershed.e. lenght of the longest channel.The equation is using to calculate as follow : where : Qp = peak flood discharge (m3/sec) Ro = rainfall unit (mm) Tp = grace time (time log) from the beginning of the rain to the peak of the flood (hour).T0.3 = time required by the decrease in discharge, from peak discharge to 30% of peak discharge (hours).
The SCS-CN method is empirical, employing a curve number (CN) approach based on land cover, soil type, and prior soil moisture conditions.CN values are derived from the analyzed land cover and soil type data sourced from DSMW (Digital Soil Map of the World).These CN values are determined using the CN number table, considering the catchment's average Antecedent Moisture Content (AMC) conditions, ranging from 0 to 100 (Table 2).Additionally, the difference between rainfall and surface runoff (S) is computed using equation ( 4), while the effective rainfall depth is determined using equation ( 3).
where : CN = curve number; A = total area catchment Runoff volume (flow hydrograph) can be obtained by multiplying the thickness of surface runoff with watershed area [28].Peak discharge can be obtained using the following equation [29]: =   2 +   (7)   = 0,6    (8) where : qp = peak discharge (m3/det); C = constant (2.08);A = watershed area (km2); Tp = rise time or the time required between the onset of rain until it reaches the peak of the hydograph (hours); tr = effective rainfall duration (hours); tp = lag time (hours); Tc = concentration time (hours).The concentration time can be calculated by Kirpich's formula as follows [29]: where : Tc = time of concentration (minutes); L = stream length (m); S = slope

E. Schematisation of 2D Flood Hydraulics Simulation
Flood simulation schematisation consists of several stages.The schematisation carried out includes: -River geometry: The channel geometry parameters required in simulation modelling are flow path, cross section, input river junction, manning coefficient, and energy loss values (contraction and expansion coefficients).The river geometry set up consists of input boundary, bank station, flow path, and cross section.-2D geometry set up Schematization begins with GeoHECRAS 2D flow area tool in the geometric data menu, creating the river's geometry.This links the simulation from the river's flow to the overflow area.GeoHECRAS forms an adaptive mesh, adjusting its size and shape using input data like DEMNAS.This mesh incorporates elevation values governing flow direction.For this study, the mesh will align with the Singkil watershed's boundaries, extending from the Singkil river boundary downstream.

-Set up 2D Flow Areas
In GeoHECRAS, the mesh generated is an adaptive mesh, which can vary in size and shape, according to the shape and origin of the data to be made into mesh data, in this case DEMNAS.The mesh spacing value (x and y) will be 8.5 m, which means that every 1 mesh/cell represents 8.5 metres in the field.The mesh contains elevation values that will determine the flow direction and will be given a roughness value (manning) to make the flow velocity match the natural conditions.-Set up flow Boundary In this study, we'll employ unsteady flow hydraulic modeling.The upstream boundary is defined by input flood hydrographs and the river's normal water level.In the downstream section, which is also part of the 2D flow area, we'll use fluctuating tidal values via the Time Series Gate Openings Boundary.This approach enables the integration of tidal elevations' impact at the river mouth on flood behaviour.As for the 2D flow area boundaries, we will adhere to the Singkil watershed boundary, extending from the upper river boundary to the downstream watershed boundary.

F. Validation and calibration of flood simulation
To validate the flood hydraulics model, actual flood depth data from a prior flood event will serve as a reference in the simulation.In this study, flood height data collected through field measurements, community interviews, and discussions with stakeholders will be used for validation.Equations will be employed to calculate and validate the flood simulation results as follows: Calibration occurs when there's a notable gap between flood model results and validated observations.In this research, we'll employ the surface roughness (Manning) value as the calibration parameter.The goal is to pinpoint the optimal parameter values that minimize the disparity between simulated and observed flood data, achieving a precise alignment.
For more detail, the development flow of 2D flood hydrodynamic simulation can be seen in Figure 3.

Delineation catchment area
Based on the morphometry of the Singkil watershed which has many tributaries that affect flood conditions, especially in the downstream part of the watershed, the catchment area will be divided according to the influential river network, Singkil river and Lae Cinedang River.

Fig. 3. Catchment Area and Domain 2D Simulation
The catchment area as shown in Figure 4 was determined based on DEMNAS data with a resolution of 8.6 m, using ArcGIS application with ArcHydro tool.From the delineation results, the total catchment area for runoff discharge calculation for 2D hydraulics flood simulation is 9,869.27km2 .Table 4 provides details of the area for each sub-catchment.

Design rainfall/return period
Design rainfall was computed from a decade's worth of maximum daily rainfall data (2013 -2022) obtained from the Blangkejren rainfall station managed by Meteorological, Climatological, and Geophysical Agency (BMKG) Aceh Province, as depicted in Figure 5. Four probability distributions -normal, log-normal, log Pearson type III, and Gumbel -were applied in this study, with their equations detailed in [25].Two probability distribution testing methods, chi-square and Smirnov Kolmogorov, were employed.Both tests identified log Pearson III as the most suitable probability distribution for the rainfall data series used as shown in Table 5, contain analysis results for probability distribution requirements base on chi-square (C s ) and Smirnov Kolmogorov (C k ) requirement.The planned rainfall is calculated based on the selected distribution which is log Pearson III distribution as distribution requirement.The results of the log Pearson III design rainfall calculation for the catchment area are shown in Table 6.The planned rainfall was transformed into planned flood discharge using hourly rainfall.The percentage distribution of rainfall that occurs in the study area for 24 hours is obtained from the IDF (Intensity-Duration-Frequency) curve.The results of rainfall intensity calculations using the Mononobe method for return periods of 2, 5, 10, 25, and 50 years are presented in Table 7.The calculation of the plan flood discharge is carried out using the Nakayashu Synthetic Unit Hydrograph (HSS) method for Singkil river and SCS-CN method for Lae Cinedang river catchment.The complete plan discharge calculation is presented in Table 8 and the plan flood hydrograph of the Cathment Area of the Singkil River flood and Lae Cinedang River presented in Figure 5 and 6.
In the Singkil river catchment, the duration of direct runoff reaches its peak at 37 hours after rainfall onset and persists for up to 150 hours.This extended duration is attributed to the watershed's large size, elongated shape, numerous river bends, and relatively gentle slopes across most of its terrain.In contrast, for the Lae Cinedang river, direct runoff lasts for 43 hours, with the peak occurring at 29 hours.This shorter duration and lower runoff volume in Lae Cinedang are a result of the catchment's smaller size and considerably shorter river length  In the first step, the initial schematisation conducted as the base of flood model for different return periods 2, 5, 10, 25 and 50 years.In this study, flood simulation using unsteady flow data with direct runoff input was carried out to see the characteristics of floods that occur such as duration, depth, water velocity and flood distribution.

Fig. 8. Tidal graph data
The Singkil and Lae Cinedang rivers serve as breaklines, representing river geometry in the 2D hydraulic Flood Simulation schematization, as depicted in Figure 7.This simulation encompasses a total 2D flow area covering 1,072.21km2 of land.For the upstream boundary of the 2D flow domain, the direct runoff flow hydrograph from a 10-year return period was used, with added rainfall intensity in the 2D flow areas.This setup aligns with field-collected flood depth data from a November 2012 event, identified as a 10-year return period flood based on rain intensity and event recurrence.Rain intensity is introduced to the 2D flow areas because the flood event in the Singkil district area was consistently accompanied by rain for approximately 4 days.The downstream boundary accounts for the maximum sea level tide to accommodate the influence from the sea under severe conditions, as seen in Figure 8.
This study considers the impact of Manning's surface roughness values in overland flow calculations, which are determined based on land-use types (Figure 9).For validation accuracy and model calibration purposes, error results were calculated by comparing the simulated flood height with the flood event (November 2012) using equation (10) described in the previous section.From result of calculation which given the average error value 0.301 m, we can conclude that the 2D flood model able to represent the flood characteristic and behaviour for the study area.In this case, the calibration not necessary to be carries out.For more detail the result calculation can be seen in Table 9.    Flooding starts with a 2-year return period, reaching heights of 2 meters.As return periods increase in the simulation, flood conditions worsen, with heights exceeding 5 meters and inundation durations lasting up to 10 days for a 50-year return period (Figure 11).A summary of the flood simulation results is provided in Table 10.It can be said that the flood simulation results confirm that the geomorphological conditions, elevation and slope of the downstream area of the Singkil watershed are dominated by swamp, river meander and terraces with fairly low elevations and very gentle slope which play important role in flood event [30].For more details on the effect of geomorphological conditions on flood extent can be seen from Table 11.From Table 7, we can see how the geomorphological conditions and slope of the downstream area of the Singkil watershed affect the flood extent, both slopes below 2% and swampy marshes have a flood extent of more than 50% of the total area for all return periods.Based on velocity classifications, reveal that the highest current speeds, exceeding 50 cm/sec (considered fast), occur at the outer bends of the river for all simulated return periods (Figure 12).Conversely, in the river's straight sections and flood-prone areas, flow velocity only reaches 25 cm/sec, considered slow.This reduced current velocity results from the river's width and gradient, contributing to overflow and inundation.In flood-affected regions, water flow velocity decreases due to the gentle plain slope and dense vegetation, particularly in swampy areas.Based on the relationship between return period discharge and simulated flood extent (Figure 13), it is known that the change in flood extent is in line with the increase in flood discharge rate.This is in line with the typical coastal area or downstream watershed which is low-lying and flat with a very gentle slope like the condition of the study area.Additionally, the flow velocity simulation results will offer insights for reinforcing riverbanks in high-velocity areas (>50 cm/sec), preventing erosion and potential collapses that could worsen flooding.These simulations are also envisioned to support the district's master planning, aiding in the selection of appropriate residential, office, and commercial areas, especially in light of extended flood durations.
The coastal area or the downstream part of the Singkil watershed is a flood-prone area with a high potential.From the simulation results it is known that the height of the flood can reach 5 m with the flood area reaching 45% of the total area of Aceh Singkil district.The duration of the flood also occurred for quite a long time, up to days.
2D flood hydraulics modeling with the GeoHECRAS application is able to represent flood conditions if the DEM or topography data for building a 2D simulation domain and break line to accommodate river geometry are reliable.
rainfall depth (mm); P = rainfall (mm); S = differentiate between rainfall and runoff (mm); CN = curve number Based on data on soil type, AMC and land use, the CN in the watershed can be calculated with the following equation:

Figure 7 .Fig. 5 .Fig. 6 .
Figure 7. Singkil river flood hydrograph Initial simulation results reveal a flood distribution covering an area of 12 km2, with flood heights exceeding 3 meters.The flood inundation encompasses 8 out of 9 sub-districts in the Aceh Singkil district (Figures 10), closely resembling the November 2012 flood event used for model validation.

Table 1 .
2D flood hydraulics simulation development data and their sources

Table 2 .
Soil hydrological groups

Table 3 .
Manning's roughness coefficient used to model floodplain surface roughness

Table 4 .
Area of per catchment

Table 5 .
Determination of distribution type

Table 6 .
Log Pearson III design rainfall

Table 8 .
Design flood discharge

Table 9 .
Calculation of model error validation

Table 10 .
Summary flood simulation

Table 11 .
Correlation flood inundation, geomorphology unit and slope