1 Introduction

The glacial retreat observed in many mountain ranges in recent decades is now believed to have contributed to the initiation of many catastrophic events (Carey et al. 2012). This retreat leads to the formation of new glacial lakes as well as to increasing the volume of water already stored within the existing lakes (Zapata 2002). These waters are often held behind young unstable frontal moraines which are considered to be most susceptible to sudden breach and the subsequent generation of a glacial lake outburst flood (GLOF) (Clauge and Evans 2000). Some frontal moraines may contain dead ice, and this further diminishes moraine dam stability (Worni et al. 2012). The triggers that cause GLOFs may range widely, but the most significant are usually associated with slope movements (Awal et al. 2011; Emmer and Cochachin 2013). The requisite displacement wave may be generated by rock/ice falls or avalanches as well as landslides (Kershaw et al. 2005) although it is also clear that there are significant regional differences in the most important triggering factors (Emmer 2011). The increasing global interest in GLOFs has been stimulated by a greater appreciation of their immense destructive power and by the increasing number of vulnerable people who now live below glacial lakes (Kershaw et al. 2005; Bajracharya et al. 2007; Bolch et al. 2008; Wang et al. 2008; Narama et al. 2010). There have been numerous examples of this phenomenon in the Cordillera Blanca, Peru (Lliboutry et al. 1977; Zapata 2002) (Fig. 1). The event at Palcacocha Lake in 1941 was one of the most devastating as it led to the destruction of a large part of the regional capital of Huaraz and killed 4,000 people (Vilímek et al. 2005). Another event triggered by ice calving into Jancarurish Lake in 1950 caused enormous economic losses as it destroyed a hydropower plant that was close to completion (Lliboutry et al. 1977). In contrast, a rock avalanche into Safuna Alta Lake in 2002 initiated an impact wave exceeding 100 m, but resulted in only a minor outburst flood due to the fact that the moraine dam was not breached (Hubbard et al. 2005).

Fig. 1
figure 1

The study area and its location (L lake, V village, C Cordillera, R river)

Full size image

The number of inhabitants living below glacial lakes in the Cordillera Blanca is increasing, while simultaneously its inhabitants are increasingly reliant on the lakes as they provide water. This water is important not only for domestic consumption but also for electricity generation and agricultural irrigation projects along the Pacific coast. This situation demands effective hazard and risk mitigation in the inhabited areas of high mountains, and more than 50 years of glacier-related hazard mitigation and management efforts in the Cordillera Blanca were recently summarised by Carey et al. (2012). It was noted that despite the fact that experts working in Peru have acquired considerable experience in hazard mitigation using technical measures, hazard assessments have been largely neglected (e.g. the spatial or temporal prediction of potentially dangerous GLOFs). This may be due to a strong and historically demonstrable resistance from the local inhabitants to hazard zoning (Carey 2010), and it is, therefore, necessary to find alternative strategies that are able to lower the risk to local population. These strategies should be able to define hazard zones with a reliable scientific basis irrespective of the fact that this task is intrinsically difficult in high and remote mountain areas due to the lack of adequate hydrological records (Bohorquez and Darby 2008). Therefore, the accurate numerical backwater analysis of the past GLOF events may be a crucial information sources for future hazard assessment as presented in Dussaillant et al. (2010), where one-dimensional steady-state calculations of peak discharges were made. It seems that comparatively little progress has been made on this topic in the Cordillera Blanca in comparison with other mountainous regions (e.g. Cenderelli and Wohl 2001; Carrivick 2006; Bohorquez and Darby 2008). The existing studies include modelling of possible debris flows from Rajucolta Lake conducted using the software Flow2D, but these results have not yet been published (Oscar Vilca, personal communication). Furthermore, a scenario-based debris flow model from the Lake 513 was made prior to the event using the model Flow2D (Valderrama and Vilca 2010), while after the event, it was modelled using the debris flow simulation model RAMMS (Huggel et al. 2012). The debris flow models made using Flow2D were also prepared for other sites within the Peruvian Andes (Valderrama et al. 2007), while the models made using RAMMS were then used to prepare GLOF hazard maps for the Chucchun River catchment for the purpose of helping local authorities develop planning praxis.

On 11 April 2010, an ice/rockfall into Lake 513 triggered a GLOF along the Chucchun River in the Cordillera Blanca of Peru. In this study, we use a widely applied approach for reconstructing unrecorded flood events (Carrivick 2006; Bohorquez and Darby 2008), which is deemed necessary because there is almost no hydrological information about the event (e.g. unavailable flood hydrograph). The main objectives of the study are to determine the extent of the water surface, the depth of the water surface, and the peak discharge of the flood event. The latter is used as flood hazard indicator, as it represents the worst flood conditions, and this is then used to define the hazard classes (Lateltin et al. 2005). In addition, we undertook scenario-based modelling of temporary stream damming, as this may possibly increase the flood hazard along specific stream reaches.

2 Study site

Lake 513 is located beneath Hualcán Mountain (6,122 m asl), and its water level is found 4,431 m asl. The lake is a high-altitude glacier lake dammed by a combination of bedrock and frontal moraine, while there are also another two glacial lakes located nearby (Fig. 1). It is one of the sources of the Chucchun River, which runs for 14.4 km before entering the Santa River. The distinct alluvial fan formed at this confluence provides the location for the city, and local administrative centre, of Carhuaz. Its highly hazardous position has been recognised for a long time, and considerable efforts to apply structural mitigation measures in order to lower the GLOF hazard have been made during the second half of the twentieth century. The first scientific recognition of the hazard posed by GLOFs triggered by ice/rock avalanches from Hualcán Mountain was presented by Lliboutry et al. (1977), while the first mention of the formation of Lake 513 appears in a routine report by the Unidad de Glaciology y Recursos Hídricos (in Reynolds et al. 1998). However, only 8 years later, the lake contained some 1,500,000 m3 of water, enough to inundate the city of Carhuaz (Carey et al. 2012). The lake was, at that time, considered to be very dangerous, and a scheme began in order to gradually lower its water level—the final technical mitigation measures were completed in 1993 (Carey et al. 2012). It is presently contained behind a bedrock dam with a 19-m freeboard in addition to the moraine, which originally dammed the lake which may be up to 11 m high. The water exits through an artificial tunnel whose entrance is fixed at 4,431 m asl. The bedrock dam comprises the coarse-grained granodiorites and tonalites (Wilson et al. 1995) that occur across the entire upper part of the Chucchun Basin.

The Cordillera Blanca normal fault (Wilson et al. 1995) runs along the foot of the range and divides the Cordillera Blanca batholith from the various Mesozoic and Cenozoic rocks that are found on the graben-like valley of the Santa River. These rocks are partly covered by fluvioglacial deposits. The sediments from the Mesozoic era include limestones, claystones, sandstones, and quartzites. Locally, erosional remnants of Neogene pyroclastic flows are preserved above the Mesozoic rocks. The basic morphological characteristics of the lake basin predominantly reflect glacial activity and variations in the bedrock lithology. The mountain peaks and intervening ridges, which dip between 40° and 55°, are glaciated, while the exposed bedrock continues downslope until it is covered by predominately glacial deposits below 4,240 m asl and on down to 3,400 m asl. There are distinct moraines formed here, and these bound a former palaeo-lake basin known as Pampa Shonquil. It is a broadly horizontal surface between 3° and 10°, which has been filled with lacustrine deposits, while the drinking water facility for the city of Carhuaz has been constructed here. The valley narrows farther downstream as it cuts through the moraines and before widening again as it enters the Neogene and Mesozoic rock formations. In this part, the main right-hand (Upacoto Basin) and left-hand (Cachipampa Basin) tributaries enter the Chucchun River before cutting through the more resistant members of the Mesozoic formations. The final reach crosses an alluvial fan that comprises flood and debris flow deposits before its confluence with the Santa River, located at 2,650 m asl. The river maintains a single channel for the whole of its course; no braiding is found at all, not even on the alluvial fan.

3 Methods

3.1 Field investigations

The intensive field investigations aimed to collect all the data needed for the model input as well as the information necessary for model validation. It comprised the following: cross-sectional geodetic surveying of the channel and adjacent floodplain; describing and photographically documenting the surface characteristics of these cross-sections; mapping the flood sediments; defining bed gradation; identifying flood marks; describing the damage; and interviewing local inhabitants. These data were collected during three separate campaigns. The first was conducted by experts from the UGRH in Huaraz within the hours of the ice/rock avalanche. The second was conducted by the authors in June 2010, just 2 months after the GLOF, and focused on mapping the extent of the water surface during the flood and to identify flood marks; these were then used as the boundary conditions for the model. The third was conducted by the authors in September 2011 and focused on additional cross-section measurements that were then used to improve the mathematical model. A total of 120 cross-sections along the 14.4-km-long Chucchun River were measured using a hand-held laser range finder localised with Garmin GPSMAP 60CS. The LaserAce 1000 laser range finder measures distances up to 150 m with precision of 0.05 m and angles with precision of 0.2° and, thereby, provides submeter accuracy for measurements of distances up to 100 m. It can obtain measurements through direct terrain reflection, which is imperative when surveying in difficult terrain, such as that of the study area. The steep slopes would make it extremely difficult and time-consuming to move a reflector along the selected cross-sections, while the acquisition of each individual point along the cross-profiles enabled all of important slope changes to be captured. The location of the cross-sections was selected so that all the important changes that could affect the modelling results in both the channel (i.e. sudden changes in cross-profiles and longitudinal profiles as well as around confluences) and overbank conditions were described (Brunner 2010, Fig. 2).

Fig. 2
figure 2

The longitudinal profile of Chucchun River with the locations of the constructed cross-profiles used for modelling; the position of each profile is indicated by a thin line below the longitudinal profile

Full size image

It is assumed that the major channel changes occurred during the modelled peak flow, and therefore, the channel topography measured after the flood closely reflects the conditions found during peak discharge (Bohorquez and Darby 2008). The minimum water level during the flood was determined by the location of fine-grained flood sediments deposited during the flood. Unfortunately, only two reliable flood marks could be identified on buildings, and these were used to determine the water surface elevation during the peak flow. This problem is compounded by the fact that they are both located close to one another on the alluvial fan at Carhuaz (Fig. 3). The surface of each cross-section was described and documented, while six bridges were identified along the length of the river (Online Resource 1). These were each measured and then added to the mathematical model. Finally, at two sites (Fig. 3), the conditions that caused temporary damming of the river were described and modelled. The standard granulometric analyses of the flood sediments, disregarding large boulders, were undertaken on five samples (Novotný and Klimeš, submitted). The Wolman method was applied at another eleven sites to determine their grain size distribution (Matoušek et al. 2011). One hundred pebbles, chosen at random from the actual river bed affected by the GLOF, were measured at each site for their median size. These measurements were processed to standard granulometric curves. The results of the geomorphological field mapping were used to determine which reaches were characterised by prevailing accumulation or erosion processes as well as for identifying sites where temporary damming may occur during future events as a result of landslides.

Fig. 3
figure 3

The location of the identified flood marks, sites of flood sediment sampling, and sites where the stream may become blocked during a flood; the numbers indicate river stations

Full size image

3.2 Hydraulic modelling

The hydraulic calculations were all made using the one-dimensional mathematical model HEC-RAS version 4.1.0. This software allows the user to undertake steady-flow river hydraulic calculations (i.e. water flow is considered to be constant in time), unsteady-flow river hydraulic calculations (i.e. water flow varies in time), sediment transport/mobile bed modelling, and water temperature analysis. It is also possible to simulate water flow through bridges, culverts, and other hydraulic structures or flow obstructions. There is a digital supplement to this article in which the topographical input information used in the presented HEC-RAS model is available and may be used by any researchers who wish to use their experience to contribute further to this modelling. The geometry of the stream is specified by the ground surface cross-sections of the channel and floodplain for particular measured intervals (reach lengths), and therefore, the flow conditions (cross-section shape, surface roughness, and slope) are modified from one section to the next. The six bridges at Stations 0.235, 0.881, 1.878, 3.292 km (with destroyed deck) 4.376, and 5.061 km were modelled, while the landslide was modelled as an inline structure in the cross-section (see also Online Resource 1). A mix flow regime was used for the steady-flow analysis with a supercritical flow regime modelled in the steep upper part of the river and a subcritical regime modelled in its shallower lower part. This was imposed by the Froude number calculation made by the model. Due to the steepness of the river bed in the upper part of the model (approximately 30 %), it is necessary to emphasise that the steady-flow analysis in this reach represents a rather rough estimate (due to aeration of the stream, different velocity distribution, etc.). We tested the possibility of unsteady-flow analysis, but the model was unable to process the calculation, especially in the steep upper part of the channel (see also Chow 1959), and provided unrealistic results.

With regard to the geomorphological conditions, both upstream and downstream, boundary conditions have to be set in HEC-RAS. The boundary condition in the uppermost cross-section was set as critical depth because this section has a considerable stream gradient and is located at the edge of Lake 513, while the boundary condition in the lowermost cross-section was set as normal depth (this section is located in the vicinity of Santa River, but at a sufficient distance, 820 m, so as not to be affected by the trunk stream). The peak discharge represents another necessary boundary condition for steady-flow analysis, and this was not known prior to modelling. It had to be estimated repeatedly until there was a correspondence between the entered peak discharge and the level of two measured flood marks. These flood marks, one on a building and one on a fence, were found approximately 825 and 885 m from the confluence with the Santa River, respectively (Fig. 3). Manning’s n value is a type of loss coefficient used for evaluating friction loss between two adjacent cross-sections: the value is highly variable and depends on a number of factors including surface roughness, size and shape of channel, channel irregularities, channel alignment, scour and deposition, vegetation, obstructions, seasonal changes, stage and discharge, and suspended material and bed load. The specific values for each cross-section were determined by comparing field photographs with values in hydraulic tables (Chow 1959). Further adjustments were made while considering the modelled depth of the flow and bed slope as the higher the flow depth, the lower the applied friction coefficient in the channel. On reaches with steep channels, the coefficient values were adjusted using the formula of Jarrett (i.e. the values were increased compared to those suggested by Chow 1959, see also http://hecrasmodel.blogspot.cz). The comparison between our field photographs and those taken from catalogues of rivers with known roughness values (e.g. http://manningsn.sdsu.edu/) was used as a reference only since they do not define coefficients during floods and, therefore, cannot be applied without further adjustments to the performed hydraulic modelling.

4 Results

4.1 The GLOF event based on eyewitness accounts and fieldwork

The personnel responsible for maintaining the intake of drinking water for the city of Carhuaz reported that the avalanche began at 07:40 on the 11 April 2010. The avalanche detachment area was located at around 5,200 m asl (Fig. 4), and its estimated volume was between 300,000 and 1,000,000 m3 (Schneider et al. 2012) or 200,000–400,000 m2 (Carey et al. 2012). It moved across steep glaciers and down a 190-m high cliff into Lake 513, whose constant water level is maintained by the outflow at 4,431 m asl. The displacement wave is thought to have attained a total height of around 24 m and overtopped the bedrock dam by around 5 m. It partly eroded the overlying moraine dam and continued in a bedrock stream channel until it entered a steep deep channel incised by up to 15 m into the glacial and colluvial sediments below Rajupaquinan Lake (‘B’ in Fig. 5). From this point, intensive erosion and material entrainment continued for about 1 km until energy loss led to deposition as the steep stream gradient shallowed to the broadly horizontal surface of the Pampa Shonquil (‘C’ in Fig. 5). It has been determined that between 50,000 and 80,000 m3 of mainly coarse sediment (11.8 in Figs. 2, 5) was deposited over an area of about 100,000 m2. Further downstream, although still within the palaeo-lake basin, the river eroded into the sediments during the highest water discharge. Thereafter, the stream slowed and accumulation was dominant. The accumulated material was comparatively fine-grained (10 in Figs. 3, 6) with more than 23 % of the particles finer than 10 mm and 50 % finer than 80 mm. This contrasts with the accumulated material at the onset of the palaeo-lake basin in which only 19 % of the particles were finer than 100 mm and 50 % were greater than 270 mm—this was the coarsest sediment sampled during the fieldwork (11.8 in Fig. 6).

Fig. 4
figure 4

The source area of the ice/rockfall that occurred on the south-western facing wall of the mountain of Hualcán in April 2010 (indicated by red circle)

Full size image
Fig. 5
figure 5

The prevailing types of fluvial processes that occurred along Chucchun River during the flood event in April 2010

Full size image
Fig. 6
figure 6

The granulometric curves obtained from the flood sediments using the method of Wolman. The table shows the diameter of the particles that comprise 50 % of the sample (the samples are named according to the river station—for locations, see Fig. 3), while the lower part of the figure shows the sample location within each cross-profile (indicated by asterisks)

Full size image

The stream was then predominantly erosive for the next 2.5 km after leaving the Pampa Shonquil (‘F’ in Fig. 5). The erosion first affected moraine sediments, but in the final part of this reach, the stream passes through a bedrock canyon with widths as narrow as 7 m and walls as high as 60 m. The bedrock hereabouts suggests that the erosion was minimal, while accumulation occurred only very locally. The flood then left the canyon some 300 m above Paricacaca Village and continued to erode laterally. There was considerable lateral erosion which affected the left bend around Paricacaca Village with up to 10 m of old alluvial sediment was washed away which caused damage to agricultural fields and the loss of livestock. In certain restricted parts of this reach, there was also deposition in places where the sediment transport capacity was low. The sample from site 5.8 in Fig. 6 comprises 50 % cobbles and boulders with diameters of more than 120 mm, while only 6 % of the particles were finer than 10 mm. The next reach is characterised by both accumulation and erosion processes (‘G’ in Fig. 5). The final reach is characterised by accumulation over the last 5 km (‘H’ in Fig. 5), starting with the deposition of more coarse material (3.3, 2.3, and 1.4 in Figs. 3, 6) and finishing with the deposition of the finest material near the confluence of the Chucchun and Santa rivers (0.4 in Figs. 3, 6). This sediment comprises 50 % of particles smaller than 8 mm and was deposited on a floodplain with a width of up to 280 m.

The direct damage caused by the flood comprised aesthetic or functional damage (i.e. a destroyed fence) to several residential houses located either within the river bed itself or on the floodplain in the lower part of the course (‘I’ in Fig. 5). In some places, the original surface was covered by 0.6 m of flood sediments, which damaged crops and houses. In Pariacaca Village, several cattle were lost, and in one instance, a bridge deck was destroyed. The intake of the drinking water treatment plant at Pampa Shonquil was blocked by fine-grained sediment (10 in Figs. 3, 6), and consequently, several thousand inhabitants in Carhuaz were left without drinking water for 5 days (Valderrama and Vilca 2010). The event also deeply affected local population from a psychological point of view as it rekindled memories of previous catastrophic events in the Cordillera Blanca.

4.2 The GLOF event based on the modelling

As mentioned in the methods, unsteady-flow simulations were attempted considering a very rapid increase in the flood wave (compared with Fig. 1 in Walder and Costa 1996) within a hypothetical hydrograph, which requires a very dense set of temporal calculations. The time step can be calculated from, for example www.nws.noaa.gov. We ran several experiments using the smallest available time step in HEC-RAS, while the high channel slope in the upper reach of the stream was also lowered to ensure better modelling conditions, and all flow obstacles were eliminated. Despite these simplifications and the use of mix flow rather than reaches with supercritical and subcritical flows, and the interpolation of cross-sections in the upper reach, the model was unable to process the calculation and provided unrealistic results. Therefore, only results of the steady-state modelling are discussed further.

The results of the HEC-RAS model show that the peak discharge of the event was 580 m3 s−1 and the total volume of water to leave the lake is estimated to have been about 840,000 m3. It is estimated that the peak discharge took 40 min to travel from Lake 513 to the Santa River, with an average flow velocity of 4.5 m s−1. This is in general agreement with eyewitness accounts, which suggest that the whole event spanning from rock/ice avalanche detachment until the peak discharge reached the alluvial cone took about 60 min. The results of a selected cross-section, the reach situated in the lower part of the model between Stations 0.454 and 3.055, are shown in Table 1. The selection of an appropriate value for Manning’s n is important for the accuracy of the computed water surface profiles. Figure 7 shows an example of a cross-section with its Manning’s n values at Station 3.797 as well as a field photograph of the site. The Manning’s n values used in the model for the channel ranged from 0.035 m1/6 (for bare cross-sections with almost no vegetation) to 0.075 m1/6 (for cross-sections with large stones and vegetation) with typical values of 0.045 or 0.05 m1/6 and correspond well with, for example, the values suggested by Dussaillant et al. (2010). The values used in the model for the floodplain ranged from 0.035 m1/6 (for bare cross-sections without any vegetation or other obstacles to water flow) to 0.13 m1/6 (for profiles with dense vegetation and other obstacles to water flow).

Fig. 7
figure 7

An example cross-section from Station 3.797 with Manning’s n values (numbers at the top of a) together with a field photograph of the site: the view downstream corresponds to the orientation of the cross-section

Full size image

The HEC-RAS software calculates the water surface elevation during the peak discharge of the flood at each of the cross-sections used for the model. Figure 8 shows the level of the flood at the site where the bridge crossing the river was destroyed (Station 3.292 km) and the level of the flood at the bridge to Hualcán Village (Station 5.061 km), which survived with only aesthetic damage. The modelling results suggest that the destroyed bridge was about 3 m below the water level during the flood event, while the bridge to Hualcán Village was only 1.2 m below it. The extent of the water surface during peak flow at each cross-section was then used to define flooding lines along the whole length of the stream (Fig. 9), which were adjusted to the actual topography. This allows the total flooded area to be estimated as 700,000 m2. In most parts of the basin, the flood was confined to the river channel, but in some parts of the alluvial fan at Carhuaz, it extended up to 70 m onto the floodplain on each side of the river (Fig. 9; Table 1).

Fig. 8
figure 8

A, A′ The water surface elevation during peak flow at the cross-section where the bridge was destroyed along with a photograph of the site; B, B′ the water surface elevation during peak flow at the cross-section adjacent to the bridge to the Hualcán Village which survived with only aesthetic damage (the bridge foundations are shown in light grey)

Full size image
Fig. 9
figure 9

The water surface extents were used to delimit flooding along the stream (blue line), while the calculated water depths were used to define the intensity of the GLOF at each cross-profile; intensity classes are shown for the lowest river reach on the alluvial fan in Carhuaz (inset of the figure)

Full size image
Table 1 The profile output table presenting the fundamental characteristics of the water flow at each cross-section
Full size table

4.3 Verification of the model

The obtained hydraulic calculations do not represent scenario-based modelling undertaken on the same basin (e.g. Valderrama et al. 2007, Valderrama and Vilca 2010, Haeberli 2013); it is intrinsically difficult to verify the latter as they do not represent past events. In our work, the hydraulic characteristics of the flood event have to be compared with field evidence to assess its reliability. The uppermost cross-section in the model is that located at the end of the Lake 513, where the water surface elevation is maintained at a constant level. Thus, using the computed water surface elevation in the uppermost cross-section and comparing it to the known water surface elevation at the beginning of the flood event, it was possible to calculate the height of the wave in the lake during the peak flow. This height was determined to be 26.5 m, whereas field observations after the flood event suggested a height of 24 m. These values correspond well especially when considering that the peak of the overtopping wave is likely to have been higher than the indicative erosion of the moraine sediment measured in the field.

In addition, there is good agreement between the computed travel time needed for the displacement wave to reach Carhuaz and that of eyewitness accounts recorded during fieldwork in 2010. The extent of the flood mapped by the HEC-RAS model was verified through comparison with field photographs taken shortly after the event (Fig. 10) as well as with the field notes of the experienced scientists who conducted fieldwork immediately after the event. The flood extent, as defined by fine-grained sediments overlying the preflood surface, was mapped in 25 cross-sections spaced along the length of the stream. In only four cases, the calculated water level was below the mapped position, but this model error was never larger than 1 m, and in most cases, the discrepancy was only 0.5 m. Therefore, there are no notable discrepancies between the model and the field observations. In addition, the flood marks found in the lowermost stream reach also accord with the modelled water surface elevations (Fig. 11).

Fig. 10
figure 10

The flood extent as evidenced by a photograph taken within few hours of the event during the first field investigation conducted by experts from ANA-UGRH in Huaraz

Full size image
Fig. 11
figure 11

The flood mark on Station 0.885 km is indicated by grey mud on the partly destroyed adobe wall. The photograph was taken by the authors during fieldwork in 2010

Full size image

It is also possible to compare the presented model with other independent models in order to indirectly assess the reliability of the presented hydraulic calculations and geomorphological mapping. It is, however, important to note that there are no other models that have tried to replicate the event described here and this makes the process of validation more difficult. Valderrama and Vilca (2010) presented a Flow2D model, prepared before the event, which modelled a possible flood from Lake 513. The characteristics of their model results are similar to those seen during the GLOF. It was predicted that the total area flooded by the debris flow would be about 1,500,000 m2, which is roughly twice the flooded area found in our model. It was also predicted that the total volume of water and debris mobilised during the event would be about 909,562 m3, which is close to our estimation of 840,000 m3. This volume is also in general agreement with the expert estimation published in Carey et al. (2012) who suggested the total released water volume was <1,000,000 m3.

The calculations of Valderrama and Vilca (2010) modelled turbulence within the flow, which distinguishes the point at which the hyperconcentrated flow transforms into a flood. The results of the model accord with the field data. The geomorphological evidence gathered during the fieldwork identified two reaches in which this transition may have occurred. The first is on Pampa Shonquil at the point where the flood sediments change from coarse- to fine-grained, while the second is at the apex of the Carhuaz alluvial fan where the transition between coarse- and fine-grained sedimentation was described (Fig. 5). The geomorphological evidence and results of modelling (Valderrama and Vilca 2010) suggest that a complex flood/debris flow model is needed to precisely describe conditions during the event. The results of the scenario-based model made using RAMMS suggest maximum discharge about 600 m3 for the river reach below Pampa Shonquil (Haeberli 2013). This value is somewhat larger than our results for the peak discharge.

4.4 The model results applied to GLOF hazard assessment

The depth of the flood calculated at each cross-section was used to define the flood intensity as classified by Loat and Petrascheck (1997). This states that high-intensity flooding is characterised by depths greater than 2 m, while low-intensity flooding is characterised by depths of less than 0.5 m. These data, together with the extent of the water during peak flow, provide information about the severity of the flood (Fig. 9). The houses damaged by the flood were all, with only one exception, located in areas where the modelled flood intensity was medium. The destroyed bridge at Station 3.292 km was located in an area where the modelled flood intensity was high as the deck was more than 3 m below the water surface. In contrast, other bridges along the river were located in places where the modelled flood intensity was medium (e.g. Hualcán Bridge at Station 5.061 km was only 1.2 m below the water level).

There is a Peruvian regulation (Law No. 29338, Water Sources Law and Regulations) that suggests a safety zone of 50 m along streams, in where no development should occur. The water surface extent of 82 % of the cross-profiles will fit within this theoretical safety zone. The regulation, however, is largely ignored, and many houses have been constructed close to the river. There were, despite this, only a limited number of buildings that were damaged by the flooding; for example, while houses in Pariacaca and Carhuaz are often only 20 m from the banks of the river, but only seven instances of damage were recorded (Fig. 9). The HEC-RAS software enables the modelling of temporary stream damming by an obstacle, which may represent a possible scenario during future flood events. The geomorphological field mapping identified three sites in which the stream could be blocked by landslides (Fig. 3). The water levels in the respective cross-sections increased by maximum of 7.7 m (Station 7.478), 6.1 m (Station 5.536), and 8.7 m (Station 5.258). The temporary damming at Station 5.258 caused an upstream increase in the water level for as much as 140 m. There are, however, no houses that could possibly be affected by increased water levels along each of those river reaches.

5 Discussion

The Cordillera Blanca has long-term historical record of GLOFs (Zapata 2002), but estimates of their volumes are available in only a few cases. Such estimates exist for floods from Jancarurish Lake (10,000,000 m3 Costa and Schuster 1988), Palcacocha Lake (8,000,000–10,000,000 m3 Evans and Clague 1994; Vilímek et al. 2005), and Artizon Alto (200,000 m3 (Emmer et al. 2014). The shorter historical record available for the British Columbia shows that GLOF events usually exceed volumes of 1,000,000 m3 although smaller volume events also occur (e.g. 300,000–400,000 m3) (Clauge and Evans 2000).The 2010 GLOF from Lake 513 probably attained a total volume of around 1,000,000 m3 when considering that the calculated water volume of 800,000 m3 was enhanced by an unknown amount of transported debris. Therefore, it can be considered as a comparatively small magnitude event, compared to other recorded GLOFs in the Cordillera Blanca. In addition, it did not cause extensive losses although it attracted large public attention, which in turn inspired specific responses from administrative and investigative bodies, not only in Peru (Carey et al. 2012) but also abroad (Huggel et al. 2012).

Our field investigations corroborate studies (Valderrama and Vilca 2010; Worni et al. 2012) which show that GLOFs reflect complex processes with potentially multiple transitions between flood and debris flow behaviour as the sediment was repeatedly mobilised and deposited according to local channel properties (e.g. available material for transport and sufficient carrying capacity of the stream). This complexity is reflected in the ambiguity of outflow flood modelling as some models consider the events as floods (Jonathan and Carrivick 2006; Bajracharya et al. 2007; Bohorquez and Darby 2008), while others consider them as debris flows (Huggel et al. 2003, 2012). It was decided to approach this event as though it were a flood as we think, this behaviour better characterises the studied event and imposes less uncertainties with respect to the modelling calibration.

We also considered the possibility of calculating an input hydrograph for the HEC-RAS model using empirical equations for discharge calculations based on a number of recorded events (e.g. ten events in Clague and Mathews (1973) or seventeen events in Walder and Costa (1996)). The common drawback of this approach, summarised in, for example, Dussaillant et al. (2010), is that the formulas are affected by the site specifics of the recorded events and these may differ significantly from place to place. Perhaps the most important site-specific condition is the process of lake water release, and this may occur, for example, as a result of glacier dam outburst following ice melt, release through a subglacial channel, or by an overtopping wave as described in this article; each of these processes generates a significantly different flood hydrograph at the lake outflow point. It is argued that the use of empirical equations is reliable in cases where the site-specific conditions are comparable, but when applied in a case where significantly different conditions occur, such as in the described GLOF, an unknown error is introduced into the modelling results. As there is no rigorous way to qualify or quantify this error, the application of this type of model for flood hazard assessment, as intended in this work, is highly unsuitable.

In order to obtain reliable model outputs regarding the peak discharge, water surface extent, and flood depths, it is important to collect accurate detailed information about the topography of the river and its floodplain because this input information is known to significantly influence the quality of the model (Worni et al. 2012). The topographical data used in this study fulfil this requirement although it is acknowledged that the model in the lowermost reach of the river, around the broadly flat alluvial fan, could be refined with more detailed topographical information. The prepared model is based on reliable topographical data that characterise the channel and adjacent floodplain, while the cross-sections are located so that they reflect all the important variations that occur along the stream (the average cross-sectional spacing is 120 m). In addition, the technique used to measure the cross-sections ensures submeter accuracy for distances of up to 100 m. Such detailed topographical input data are not typical for this type of study nor for studies in such remote areas. Somos-Valenzuela and McKinney (2011) undertook HEC-RAS modelling of the nearby Cojup Valley using a DEM with a resolution of 30 m, while other hydrological GLOF models (Bajracharya et al. 2007; Worni et al. 2012) or debris flow models (Huggel et al. 2003) use DEMs with pixel sizes of 7–10 m. There is an absence of reliable high-resolution topographical data for the studied catchment such as that provided by LiDAR. The only readily available topographical data are the 1:25,000 topographical map sheets with contour intervals of 20 m or the commercially available WorldView DEM with a pixel resolution of 8 m. Furthermore, considerable attention was given to assigning representative values to the surface characteristics (e.g. Manning’s n) at each cross-section. These values were assigned after comparing the field photographs with all other available information and cross-referencing with published work (e.g. Dussaillant et al. 2010). The detailed investigation of all known obstacles to the flow (e.g. bridges) has further improved the reliability of the flood modelling, especially in those reaches characterised by lower gradients, where it is more likely that houses may be located close enough to the stream to be affected by flooding.

The possibility of incorrectly assigning Manning’s n values and the limited availability of flood marks are sources of possible error in this model. The Manning’s n values have a significant effect on GLOF modelling in HEC-RAS as shown by the sensitivity analyses undertaken Carrivick (2006). It was found that a 10 % increase in Manning’s n values triggered a decrease in the peak discharge by up to 9 %. Other work suggests that the most commonly applied resistance coefficients usually underestimate hydraulic resistance in high-gradient channels but typically overestimate flow velocities (Yochum et al. 2012). The next set of uncertainties relates to the complexity of the modelling process and the use of a 1D approach. It is an approach commonly used for GLOF modelling (Bajracharya et al. 2007), and a detailed comparison with a two-dimensional model led to the conclusion that it is suitable for reconstructing the discharge for GLOFs (Bohorquez and Darby 2008). 1D simplification can be made as the magnitude of the flow perpendicular or vertical to the slope is smaller than the downslope component as a result of the high channel gradient, and therefore, those components may be neglected (Chow 1959).

The observed flood marks are located only in the lowermost reach of the stream on the rather flat alluvial cone, and this leaves considerable uncertainties about the behaviour of the flood in the far steeper upper reaches of the stream. The lack of flood marks thereabouts is especially problematic as the applied model tends to give less reliable results in higher-gradient channels due to air admixture into the flow, which has not yet been modelled (e.g. Jarrett 1992). It is, therefore, suggested that the modelled water extent and flood depths reflect the minimum values and perhaps an area with greater depths was affected by the flood. However, in such a scenario, it is only areas with low flood intensities that are likely to be missing and it is suggested that majority of the areas with higher intensities have been modelled accurately. It is thought that the less precise GPS elevation measurements that were used to construct the longitudinal profile of the stream may offer a potential source of inaccuracy. There are also other uncertainties inherent to any flood model that relate to changes in the channel topography as a result of anthropogenic activities. These are impossible to incorporate into any model but are important to stress to the end-users of the modelling results—new bridges, levees, or other mitigation measures may greatly influence the propagation of the flooding at specific sites and thus alter the potential losses caused by future floods.

There are other uncertainties associated with the flood event relating to the preparatory and triggering conditions of the ice/rockfall from Hualcán Mountain. It is only possible to speculate about the conditions that led to the avalanche (e.g. Carey et al. 2012). It is pertinent to note that the origin of the Huascarán ice/rock avalanche in 1962 has yet to be plausibly explained (Plafker and Ericksen 1978). In both instances, it is possible to exclude a seismic trigger for the avalanches. However, interesting information about the complexity of the avalanche behaviour during the April 2010 event has come from photographs of the source area as these clearly show that the main ice/rockfall was followed by second icefall (Valderrama and Vilca 2010). The timing of this avalanche is not known with precision nor is its effect on the lake.

6 Conclusions

A GLOF in April 2010 along the Chucchun River in the Cordillera Blanca of Peru was triggered by ice/rockfall into Lake 513. The basic hydrological characteristics of the flow were not measured at the time of the event, and therefore, the only evidence for them comes from eyewitness accounts and limited field evidence (e.g. flood marks). This information was used in conjunction with the results of detailed field surveys of the topography and the surface characteristics of the river channel and adjacent floodplain in order to reconstruct the basic hydrological characteristics of the flood using the software HEC-RAS. The one-dimensional hydrological calculations were used because of the absence of input data (e.g. detailed digital elevation data and a measured flood hydrograph) and information regarding the steepness of the channel (in places up to 30 %). The reconstructed hydrological characteristics constrain the extent of the water surface and its elevation at each of the 120 cross-sections in addition to constraining the peak discharge (580 m3 s−1). These data have been used to produce a map of the flooded area and to classify flood intensity based on calculated depth. This then provides retrospective hazard information and shows that the majority of houses exposed to the flood were subjected to medium-intensity flooding. It is suggested that the calculated flood characteristics represent rather minimum values, possibly underestimating peak discharge and flood depths while overestimating flow velocities, given the possible inaccuracies of the calculations caused by the limited availability of flood marks; uncertainties in the characterisation of surface roughness; and general difficulties associated with undertaking hydrological calculations along high-gradient streams. Despite these potential problems, there is reasonably good agreement between the obtained results and the evidence obtained during fieldwork, eyewitness accounts, and the work of independent research teams. It is proposed that the presented calculations may be used to calibrate possible future flood scenarios along the Chucchun River as well as for neighbouring streams with similar hydrological and geomorphological characteristics.