1. Introduction

1.1. General Information

Cyclone GAMANE hit Madagascar between 15 and 29 March 2024, and the heavy rains that accompanied it caused a great deal of damage, loss of life, destruction of infrastructure, etc. as a result of overflowing rivers and the flooding that followed. Among the vital infrastructure destroyed was the bridge over the Mahavavy River at Ambilobe on Route Nationale 6 (RN6), built in 1951, as shown in Figures 1 and 2. Also as a result of GAMANE, many other bridges over other rivers and parts of roads have also failed on the RN6.

Figure 1. View of the bridge over the Mahavavy river and the partially destroyed section

Figure 2. The left-bank abutment and part of the bridge were swept away by the flooding of the Mahavavy river

A prolonged closure of the RN6 would have a particular impact, however, as it is the only road linking the DIANA Region to the rest of the country, and is the route used by all supplies to this northern part of Madagascar, including fuel, medicines and food products. The bridge over the Mahavavy was destroyed on 27 March 2024, and immediately afterwards, the cost of foodstuffs soared, such as rice, a staple food, the price of which increased three or fourfold. In addition, as a result of the bridge being cut, the inhabitants of the town of Ambilobe could only link the two parts of the town by paying for canoes, which was a serious obstacle to travel, particularly for schoolchildren whose schools had initially been closed for more than two weeks.

The only flow data available on the Mahavavy River is contained in a 1993 book [1] and concerns the monthly flows measured between 1967 and 1969 and between 1980 and 1983. Exceptional flood flows were evaluated on the basis of a 28-year sample (1949-1984) and the flow corresponding to a return period of 100 years was evaluated at 8250 m3/s [1]. Nor is there a calibration curve (height-flow curve), but only the depth of water measured for different values of annual flood flow between 1928 and 1983. In short, the Mahavavy river catchment is not gauged, which is why hydrological modelling is needed to assess the flood hydrograph corresponding to the passage of cyclone GAMANE. One of the most widely used models in the world for simulating such extreme events under a variety of conditions is the HEC-HMS (Hydrologic Engineering Center - Hydrologic Model System) model [2,3]. The HEC-HMS model is a mature model that has been used for several decades on several continents, for example in Indonesia [4], southern Africa [5], Romania [6], Jordan [7], India [8], Morocco [9], China [10], Madagascar [11], etc.

1.2. Location of the Study Area

The Mahavavy River is the most important river in the DIANA Region and its catchment area is located in the south-eastern part of the DIANA Region (Figure 3).

Figure 3. Catchment area of the Mahavavy River in the DIANA Region. Inset: location of the DIANA Region on the island of Madagascar

1.3. Outlet of the Catchment Studied

Although the main target of the study concerns the bridge over the Mahavavy river, the outlet of the catchment area was not defined directly on the bridge but approximately 1800 m upstream of this bridge. This outlet was chosen for two reasons:

1. There is still a diversion on the left bank of the Mahavavy river (Canal SIRAMA) between the chosen outlet and the bridge (Figure 4).

2. Figure 2 shows that the destruction of the bridge was probably caused by excessive water velocities eroding the left bank and associated thrust forces on the deck. Therefore, in order to properly assess these velocities and forces in hydrodynamic modelling, the flow upstream of the bridge must be studied and the resulting flood hydrograph used as the boundary conditions for hydrodynamic modelling.

Figure 4. Location of the outlet chosen to delimit the catchment area (Google Satellite image)

2. Materials

2.1. Rainfall Data

In the study area, there are no weather stations capable of providing the rainfall data required as input to the hydrological model. Consequently, gridded daily rainfall data (GPM IMERG), downloaded from and reduced to 21 stations distributed in and around the catchment, were used. For each of these stations, the data is in the form of a time series from 15 March to 30 March 2024.

2.2. Raster Data

There are three types of raster data required for hydrological modelling: the DEM raster for catchment extraction and the hydrographic network, the land use and land cover raster and the raster containing hydrological soil classes.

The DEM was obtained by merging tiles covering the study area and downloaded from . These rasters had a resolution of 12.5 m. For the Land Use and Land Cover (LULC) raster, it was downloaded from . This raster is derived from ESA-Sentinel 2 satellite images, which have a resolution of 10m, and comprises 10 land cover and land use classes. The raster containing the hydrological soil groups was downloaded from the database available at . This data was derived from the soil texture and bedrock depth classes provided by the Food and Agriculture Organization of the United Nations (FAO) and is consistent with the Soil Conservation Service (SCS) Curve Number method [12]. The resulting raster had a resolution of 225 m at the latitudes of the study area.

3. Methodology

3.1. Preliminary GIS Processing

As the rasters containing the different data described above have different resolutions, they had to be resampled to a common resolution of 30 m and then projected into the same Reference Coordinate System, which is the UTM / WGS 84 Zone 39S.

3.2. HEC-HMS Model

In the present study, the approach is conceptual and the global model is semi-distributed, which led us to choose the HEC-HMS (Hydrologic Engineering Center - Hydrologic Modeling System) software developed by the US Army Corps of Engineers (USACE). HMS has been designed to simulate precipitation and runoff processes in dendritic catchment systems and contains integrated lumped and distributed model tools for modelling hydrological processes. It consists of several components for calculating precipitation losses, direct runoff and routing [13,14]. In this process, the contribution of the various rainfall stations on each sub-catchment was carried out using Thiessen polygons.

3.3. Direct Runoff Model

The direct runoff model is a model that describes the movement of water that has not been infiltrated or stored in the catchment. To apply this direct runoff model, it is necessary to have a method for assessing rainfall losses (the part of rainfall that is not transformed into infiltration or storage) and a method for transforming excess rainfall into runoff. The methods used in this study are described in the following two subsections.

3.3.1. Losses Method

The SCS CN (Soil Conservation Service - Curve Number) method was used to estimate losses and, consequently, excess rainfall. The CN is a function of land use and land cover, antecedent humidity and cumulative rainfall. According to this method, excess rainfall is given by the following equation [15,16,17]:

(1)

where = excess rainfall accumulated at time [mm]; = depth of rainfall accumulated at time t [mm]; = initial abstraction [mm]; = maximum potential retention after runoff begins [mm]. As equation (1) shows, there is no runoff when .

The relationship between the maximum potential retention [mm] and the curve number CN is given by the following empirical equation [15]:

(2)

The USDA [15] had produced tables giving CN values according to land use and soil type for average antecedent soil moisture conditions (CN(II)). For dry initial conditions and wet initial conditions, the CN(II) values are corrected according to [18]:

(3)

(4)

In addition, for a catchment made up of different land uses and soil types, the method recommends using a composite CN value, i.e. [13]:

(5)

where = composite CN used to calculate runoff volume; = index of catchment subdivisions with uniform land use and soil type; = CN for subdivision i; and = catchment area of subdivision .

3.3.2. Transformation Method

Consistent with the loss method used, the transformation method chosen was the SCS synthetic unit hydrograph, which is a parametric model for calculating the transformation of excess precipitation into direct runoff. The advantage of this hydrograph is that it requires only one parameter, which is the lag time for each sub-catchment, given by [19,20]:

(6)

in which is the time of concentration.

There are many empirical formulae for calculating , but Watt and Chow's formula has been chosen as the most appropriate. Indeed, of the many formulas giving the time of concentration present in the literature, it is one of the few formulas that can be used for large rural catchments of up to 5840 km2 [21,22,23]:

(7)

: time of concentration [hours]

: hydraulic length of main watercourse [km]

: average slope of main watercourse [m/m]

3.4. Routing Model

On reaching a watercourse (or thalweg) within a particular sub-catchment, surface runoff combines with that watercourse and the resulting flows are routed to the outlet of that sub-catchment. The same routing phenomenon occurs at the scale of the entire catchment. Unlike hydraulic methods, which focus on resolving the Saint-Venant system, hydrological methods are based on the principle of continuity and the temporary storage of excess volumes in the section. One such method is the Muskingum method used here [24-26], whose equations are:

(8)

(9)

where is water storage, is time, is inflow and is outflow.

In equations (8) and (9), and are the two parameters of the model, which are respectively the storage time constant and a weighting factor. Based on the topographic and physiographic data of the different river sections identified in the catchment, the parameter was evaluated as follows [24]:

(10)

= wave travel time [s]; = reach length [m]; = celerity [m/s].

The weighting factor was set at = 0.25 in the absence of field data.

3.5. Analysis of the Effect of Initial Abstraction

In the SCS-CN loss method, initial abstraction is related to maximum storage according to the relationship [27]:

(11)

where is the initial abstraction rate [-].

In the original version of the SCS-CN method, = 0.2 but this value was later challenged by many researchers who stated that = 0.05 or even a lower value was much more appropriate [28-31]. Thus, the effect of initial abstraction on peak flow was analysed in this study, varying from 0.01 to 0.2.

3.6. Analysis of the Effect of Time of Concentration

For ungauged catchments, as is the case here, determining the time of concentration remains a major challenge. It has been stated that 75% of the error in estimating the peak discharge of a flood hydrograph may be due to the error in estimating the time of concentration [32]. Several empirical formulae exist in the literature and, in this study, in addition to the Watt and Chow formula (equation (7)), five other empirical formulae were compared, namely:

• Kirpich formula [33]:

(12)

• Haktanir and Sezen formula [34]:

(13)

• Passini formula [34]:

(14)

• Johnstone and Cross formula [34]:

(15)

• Ponce (SCS) formula [35]:

(16)

In equations (12) to (15): = time of concentration [hr]; = length of main channel from outlet to watershed [km]; = slope between maximum and minimum elevations [m/m]; = mean slope of main channel [m/m].

In equation (16), = time of concentration [hr]; = hydraulic length measured along the main channel from the outlet to the watershed [m]; CN = Curve Number; = average slope of the catchment [m/m].

4. Results and Discussion

4.1. GIS Processing Results

4.1.1. Catchment Extraction and Modelling in HEC-HMS

Figure 5 shows the Mahavavy river catchment area delimited according to the chosen outlet (Figure 4) and its representation schematised in HEC-HMS, which has subdivided it into 7 sub-catchment areas SB1 to SB7.

Figure 5. Left: Mahavavy catchment area and physiographic characteristics. Right: Symbolisation in HEC-HMS

For modelling purposes, it is also necessary to know the contribution of each rainfall station for each sub-catchment identified. This is shown in Figure 6.

Figure 6. Location of rainfall stations (left) and contribution according to Thiessen polygons (right)

4.1.2. LULC and HSG

As already indicated, to obtain the CN we need to consider the land cover and vegetation cover associated with the hydrological soil groups (Figure 7).

Figure 7. Left: Land use and vegetation cover. Right: Hydrological soil groups

As can be seen in Figure 7 (left), the vast majority of the catchment is occupied by savannah/shrubland and forest in its southern part. In terms of soil types, there are only two groups, Group C (low infiltration rate) and Group D (very low infiltration rate). This indicates that the runoff rate is relatively very high in this catchment.

4.2. Results of the SCS-CN Loss Method

After applying a matrix algebra to the two rasters shown in Figure 7, the resulting CNs are shown in Figure 8. To implement the CN in the HEC-HMS software, the average composite CN was calculated for each sub-catchment according to equation (5).

Figure 8. Left: CN gridded. Right: average composite CN for each sub-catchment

It should be noted that the values shown in Figure 8 are normal CN values, i.e. CN(II). As the cyclonic event occurred more than three months after the start of the rainy season, the soil can be considered to be saturated with moisture and CN(III) should therefore be used. This was done using equation (4) and the results with the initial abstraction (λ =0.05) are shown in Table 1.

Table 1. CN and initial abstraction by sub-catchment

4.3. Results of the Transformation Method (Unit Hydrograph SCS)

By applying equation (6) of Watt & Chow for the time of concentration and equation (5) for the lag time, the results gave Table 2.

Table 2. Lag time for the various sub-catchments

4.4. Modelling Results in HEC-HMS

After implementation in HEC-HMS with the adopted sub-models, the flood hydrograph resulting from the passage of cyclone GAMANE is shown in Figure 9.

Figure 9. Flood hydrograph for the Mahavavy river during the passage of cyclone GAMANE, using the Watt & Chow formula for time of concentration

Although no field data is available for the value of the peak discharge, the hydrograph shown in Figure 9 is consistent with the statements made by local residents, who confirmed that most of the flooding occurred on 27 March 2024 but were unable to specify either the value or the exact time.

4.5. Results of the Analysis of the Effects of Initial Abstraction

By varying from 0.01 to 0.2, and considering the time of concentration according to the Watt and Chow formula (equation (7)), the effects of initial abstraction are shown in Figure 10.

Figure 10. Effects of initial abstraction on peak flow rate according to values of λ

These results show that, within the limits of the value of , the peak flow varies from 8376.8 m3/s (λ = 0.01) to 8371.0 m3/s ( = 0.2). This represents an absolute variation of 5.8 m3/s or, compared with the adopted peak flow of 8375.6 m3/s, a relative variation of 0.07%, which is quite negligible.

It can therefore be concluded that initial abstraction has a negligible effect on peak flow, at least for the case studied.

4.6. Results of the Effects of Time of Concentration

For the time-of-concentration formulae used, the results for peak flows are summarised in Table 3.

Table 3. Peak flows, runoff volumes and time of onset of peak flow according to the formula used to calculate the time of concentration

In Table 3, the variations are calculated with reference to the Watt & Chow formula. These results are also shown in Figure 11.

Figure 11. Flood hydrographs according to the time-of-concentration formula

It can be seen from Table 3 and Figure 11 that the values of the peak flow rates and the times at which the peak appears are different according to the time-of-concentration formula used. The minimum peak flow rate is obtained with the Ponce (SCS) formula, whereas the Kirpich formula gives the maximum peak flow rate. Table 3 also shows that the volume of runoff remains practically the same (0.5%) but that it is the dynamics of the entire flood hydrograph which is variable according to the time-of-concentration formula used, and the variation in peak discharge reaches 53.5%.

To better understand these differences, we need to go back to the origin of these different time-of-concentration formulae [34,36] given in Table 4.

Table 4. Characteristics of time-of-concentration formulae

It can be seen from Table 4 that the empirical formulas used were developed in temperate zones, which is not at all the case in the zone in which the present study was carried out (tropical zone). Despite this, and despite the fact that most of these empirical formulas are fairly old, they continue to be used - successfully - on all continents.

On the other hand, the definition of concentration time is not identical for all authors: in fact, the literature reports at least eight different definitions of concentration time, including six so-called "computer" definitions and two theoretical definitions [36]. This makes comparison difficult because the definition used to establish the different formulas used is not known. This situation was also mentioned by Grimaldi et al [37] who described the time of concentration as a "paradox" of modern contemporary hydrology because the peak flows found could show a variation of up to 500%.

The abnormally low value of the peak flow of the Ponce-SCS formula (4531.1 m3/s) is not surprising as numerous studies have shown that empirical equations involving rainfall attributes and coefficients related to land use generally have a high uncertainty and this is the case for the Ponce-SCS formula across CN [34,38]. On the other hand, the other five formulae, which only involve geomorphological parameters, have relatively low uncertainty.

In the end, therefore, the Ponce formula can be discarded and it can be stated that, during the passage of cyclone GAMANE, the peak flow occurred on 27 March 2024 and its value was reasonably between 8376 m3/s and 9009 m3/s.

5. Conclusions

The aim of this work was to produce a reliable flood hydrograph that could be used as a boundary condition for hydrodynamic modelling in order to correctly identify the velocities and forces involved during the passage of cyclone GAMANE in March 2024, which destroyed part of the bridge over the Mahavavy river at Ambilobe. Indeed, identical reconstruction is not an option, especially as extreme cyclonic events of this kind are becoming increasingly frequent and intense.

As no field data was available for the Mahavavy catchment, the runoff and routing processes were reconstructed using the SCS-CN method for estimating losses, the SCS synthetic unit hydrograph for transforming excess rainfall into runoff and the Muskingum method for routing.

A sensitivity analysis was then carried out on initial abstraction and time of concentration using six different empirical formulae. By varying the initial abstraction rate λ from 0.01 to 0.2, it was found that this had little influence, as the variation in peak flow was only of the order of 0.07%. On the other hand, according to the time-of-concentration formula used, the results showed a significant variation (53.5%) in the value of the peak flows found and a smaller difference in the time of appearance of these peak flows (from 00:10 to 10:20 on 27 March 2024).

If the Ponce (SCS) formula, which involves the CN, is discarded, it can finally be concluded that the peak flood discharge caused by cyclone GAMANE has a value reasonably ranging between 7232 m3/s and 9009 m3/s. Compared with the historical value of 8250 m3/s with a return period of 100 years, these values seem to justify the failure of the bridge without prejudging its structural condition, as it was built in 1951.

The next logical step in this work will therefore be to carry out hydrodynamic modelling of the flows (solving the Saint-Venant equations) using one or other of the selected flood hydrographs as the upstream boundary condition (Figure 11).