Introduction
Volcano impacts on glaciers include rapid thermal melting, ice and snow avalanches, surficial abrasion, gullying or mechanical scouring and basal melting. Pyroclastic flows and surges, blasts, hot avalanches and the melting of basal ice through eruptive and hydrothermal activity are typically the agents of these impacts (Reference Major and NewhallMajor and Newhall, 1989; Reference Guðmundsson, Sigmundsson, Björnsson and HögnadóttirGuðmundsson and others, 2004). The total volume of lahars, floods and/or jökhulhaups generated by heat transfer from pyroclastic density currents and/or basal melting is of the order 107–109m3. In contrast, lava flows, water ejected from crater lakes and tephra fallout are less efficient at melting snow and ice.
Volcano impacts on glaciers have been observed at several composite volcanoes, such as Nevado del Ruiz (Colombia) in 1985, Mt Redoubt and Mt Spurr (Alaska) in 1990 and 1993, Ruapehu (New Zealand) in 1995, Tokachi-Dake (Japan) in 1926 and Mt St Helens (USA) in 1980. Eruption-induced multiple-phase gravity flows that result from hot pyroclasts interacting with snow, firn and ice on composite cones can be summarized in a ternary diagram (Fig. 1). Primary lahars and floods can result from eight interaction processes on ice-clad composite volcanoes (Figs 1 to 4): (1) snow avalanches, snowmelt and slushflows triggered by dry high-energy pyroclastic surges (Reference Waitt, Pierson, MacLeod, Janda, Voight and HolcombWaitt and others, 1983; Reference WaittWaitt, 1989); (2) wet-surge transformation to lahars (Reference Janda, Scott, Nolan, Martinson, Lipman and MullineauxJanda and others, 1981; Reference PiersonPierson, 1985; Reference ScottScott, 1988; Reference WaittWaitt, 1995); (3) melting and incorporation of metre-thick snow-slab avalanches (Reference Fairchild, Costa and WieczorekFairchild, 1987); (4) scouring and gullying effects of pyroclastic flows in steep glaciers, transforming to ice-rich mixed diamicts (Reference Waitt, Gardner, Pierson, Major and NeallWaitt and others, 1994); (5) turbulent pyroclastic flows and dilute pyroclastic surges that trigger dynamic mixing and fluid drag within a thick cover of loose snow and low-density firn (Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990; Reference ThouretThouret, 1990, Reference Thouret1993; Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995); (6) mass failure of glacier ice caused by eruptive activity and seismic shaking (Reference Thouret, Cantagrel, Salinas and MurciaThouret and others, 1990a, Reference Thouret, Salinas and Murciab); (7) breakouts of glacier-dammed lakes by volcanic debris avalanches (Reference WaythomasWaythomas, 2001); and (8) tephra-fall deposits interacting with snow, ice and liquid water (Reference Manville, Hodgson, Houghton, Keys and WhiteManville and others, 2000).
Fig. 1. Multiphase gravity flows and deposits linked to volcano–glacier interactions on composite cones. This ternary diagram (modified from Reference Manville, Hodgson, Houghton, Keys and WhiteManville and others, 2000) is based on the Mt St Helens (1980), Nevado del Ruiz (1985), Mt Redoubt (1990), Mt Spurr (1992) and Ruapehu (1995) case studies.
Fig. 2. Sketch map of Nevado del Ruiz ice cap and summit region. Extent and effects of pyroclastic flows, surges and mixed avalanches on and beyond the ice cap (modified from Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990). The five drainages that conveyed channeled lahars are also shown. A–E locate the sites of the photographs in Figure 4. Inset: the Ruiz-Tolima massif in Colombia (left) and its volcanoes (right).
Fig. 3. Effects of the 13 November 1985 eruption (modified from Reference ThouretThouret, 1990). Pre- and post-eruption boundaries of ice cap. 1: ice cap affected by tephra fall and pyroclastic surges and flows; 2: missing ice resulting from avalanching, melting and erosion; 3: ice covered by thin tephra; 4: glacial marshes, ponded meltwater and temporary fumarolic activity; 5: deposits of sediment-laden snow and ice avalanches; 6: Arenas crater with fissures and fumaroles; 7: crevasses and fractures caused by the 1985 eruption; 8: gullies and grooves in the firn and ice formed by pyroclastic flows and surges; 9: rockslides, ice avalanches, and small debris-avalanche deposits; and 10: major channelled lahars.
Fig. 4. Degradations recorded on the ice cap after the 13 November 1985 eruption. (a) Facies of the 1 m thick mixed avalanche deposit on ice, Farallón-Lisa glacier catchment (hammer in box for scale). (b) Lava rock septum (19 November 1985) between the Arenas crater and the Azufrado headvalley wall showing scars of ice and snow-slab avalanches and abraded glacier surfaces. (c) Tephra-covered englacial tunnel about 2 m in diameter at the Gualí ice margin. (d) Mass-flow (laharic) deposit of glacial moraine mixed with ice blocks >1m across on 15 January 1995 in the Rio Lagunillas valley (photograph courtesy of L.F. Guarnizo). (e) Upper Nereidas glacier showing levee deposits of pyroclastic flows and surges that lowered the snow and firn cover by 3 to 6 m. (f) Mixed pyroclastic-surge deposit, 1 km southwest of the Arenas crater, showing metre-scale, cross-bedded layers of pulverized ice and lapilli ‘a’, covered by a mud layer ‘b’ 2 cm thick, which is the base of the 1 m thick massive pumice-rich pyroclastic-flow deposit ‘c’.
We examine how hot pyroclasts interact with surficial snow, firn and ice on active composite cones, release meltwater and eventually trigger lahars. On the basis of the Nevado del Ruiz case study in Colombia (Thouret and others, 1999), this report pursues the following objectives: record the degradation of the ice cap after the 1985 eruption, using remote sensing and field survey to identify preserved features and deposits; compute surface area and volume and estimate the geometry of the ice cap and of the volcano summit bedrock; reassess lahar hazards around Nevado del Ruiz; and examine the process of heat transfer from hot erupted pyroclasts to snow with lahar-generating potential and discuss how far Reference WalderWalder’s (2000a, Reference Walderb) ‘thermal scour’ model can be applied to the Ruiz data set.
Methods
Changes in the ice cap due to the 13 November 1985 eruption have been recorded by aerial and ground photographs captured in 1981 and 1985, satellite remote sensing and ground mapping (Figs 2 and 3) (Reference ThouretThouret, 1990).We used the 1:12 500 scale orthophotomap from Reference FinsterwalderFinsterwalder (1991) and a DEM digitized by us at a similar scale (Fig. 5). A 1986 SPOT scene was geocoded and draped on this DEM and compared to an aerial photograph taken in 1981, four years before the eruption, with the aim of estimating the changes due to the eruption (Figs 2, 3 and 5). The volume of ice loss was computed on the basis of ground and aerial photographs taken before and after the 1985 eruption, a 1986 SPOT XS satellite scene, 3D-orthophotographs draped on a DEM (Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995) and on thickness and specific weight measured along cores.
Fig. 5. Observed changes in ice cap of Nevado del Ruiz between 1987 and 1997 (based on the 1:12 500 scale DEM after Reference FinsterwalderFinsterwalder, 1991). 1: thick (>60 m) glacier ice covered by 6–9m thick snow and firn; 2: thin, hard and bare glacier ice with crevasses; 3: dead or dirty glacier ice, mantled by pyroclastic debris, mixed avalanche and tephra; 4: tephra >1m thick produced by the 1985 and post-1985 eruptions; 5: crevasses and ice kettles related to glacier flowage and break-in-slopes on the bedrock surface; 6: probable depression (large craters or summit caldera) filled by 100–190m thick ice; 7: glacier valley tongues affected by pyroclastic-flow and surge deposits; 8: path of the 13 November 1985 pyroclastic flows and surges that abraded and scoured the surface of glaciers; 9: bedrock showing scars of glacier-ice avalanching in hydrothermally altered rocks; 10: lava flows and scarps; 11: scar of collapse outlining a horseshoe-shaped amphitheatre; and 12: lava dome probably hidden beneath the ice cap.
Ice cap changes between 1986 and 2003 were determined from aerial and ground photographs, remote sensing by SPOT and ASTER satellites and field surveys. We used the 1:12 500 scale orthophotomap from Fin-sterwalder (1991), our DEM (Fig. 5), vertical aerial photographs (1991, 1995) and oblique photographs (1995–99) from Instituto Geográfico Agustín Codazzi (IGAC), GPS measurements (Reference Ramírez and GuarnizoRamírez and Guarnizo, 1994; Reference Vargas, Ramírez, Ruiz and NaranjoVargas and others, 2002) and ground photographs, captured before and after the eruption. These were used as references for ice thicknesses, which were measured by hand core drillings and portable radar. Analysis of deposit spectral signatures on a SPOT XS scene (8 September 1986) enabled us to outline the deposit extent using a maximum likelihood criterion (Reference Vandemeulebrouck, Thouret and DedieuVandemeulebrouck and others, 1993; Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995). Using the method of Principal Component Analysis, the ice cap image displays three different textures of ice and snow, as shown in Figure 6 (Aster images of March 2003 and March 2005): (1) a white area indicating thick ice covered by thick snow and firn; (2) dark grey/green areas corresponding to thin, often fractured or crevassed ice, which thins toward the southwest and is mantled by tephras toward the north, west and northeast; and (3) the darkest margins of the ice cap, representing dead or dirty ice mixed with tephra and debris apart from the ice-free crater area that is covered by thick tephra. The 1997 boundary of the ice cap (Figs 5 and 6) roughly corresponds to the contact between the white and grey areas as defined in the PCA image obtained from the 1986 SPOT scene (Reference Vandemeulebrouck, Thouret and DedieuVandemeulebrouck and others, 1993).
Fig. 6. Ice cap and margins of Nevado del Ruiz as observed by ASTER satellite on 7 March 2003 and 21 March 2005. The asymmetrical shape of the ice cap reveals the long-term effects of the eruption on the northwest and western catchments (cf. Fig. 5). Widespread dead or dirty ice and snowfield form the west, northwest, north, and east ice cap margins.
Portable impulse radar-based profiles were measured along eight lines totalling a distance of 6434m (Fig. 7a) across the ice cap between 1996 and 1999 to determine the ice cap geometry and volume (radar measurements are described in Appendix A). Data were analysed using GIS ILWIS 2.23 software. Where profiles were not obtained, ice thickness was interpolated from nearby profiles (Fig. 7b). Additional physical parameters were acquired from cores obtained by hand drilling in snow, firn and ice between 5100m and 5200m from 1986 to 1989 (personal communication from Reynaud, Laj and Boutron, 1990). Measured and interpolated data were laid over a 12 500 scale DEM (Fig. 7a). Annual snowfall averaged 1.5–2m (water equivalent). The snow cover was typically 3–6m thick, with a density of 300–400 kgm–3. Average firn density was 560 kgm–3, with the transition between firn and ice located ∽10m deep.
Fig. 7. Topographic DEM of the Ruiz summit. (a) Boundaries of the ice cap in 1987 and 2000, lines along which ice thicknesses were measured with a portable impulse radar, and GPS interpolated points for cross sections. (b) Two ice thickness profiles show ellipses measured from electromagnetic waves across the Nereidas and Alfombrales glaciers.
Degradation of the Nevado Del Ruiz Ice Cap by the 1985 Eruption
The Nevado del Ruiz case study helps to clarify how a large volume of meltwater can be rapidly released from a summit ice cap during a short and moderate eruption. The ice cap may have retreated, albeit slowly, between 1981 and the 1985 eruption. The 1985 eruption lasted only 20 to 90 min but reduced the area of the ice cap by 16% from ∽25 km2 to ∽20.8km2 (Reference ThouretThouret, 1990). The corresponding volume of snow and ice loss was estimated to be 6×107m3, i.e. about 9% of the pre-eruption total. An additional 25% of the ice cap area was fractured and/or destabilized. The most significant loss of ice occurred where high-energy pyroclastic flows and surges scoured or removed the snow and fractured ice mass that failed from steep-sided glaciers (fig. 3 in Reference ThouretThouret, 1990). In contrast, the ice and snow loss was lower on glaciers with gentle slopes, where tephra was passively deposited, even though the temperature of the deposit exceeded 5008C upon emplacement (Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990; Reference ThouretThouret, 1990). Other eruption-induced effects included tunnel formation at the ice margin (Fig. 4c). The pattern of parallel gullies emerging from these tunnels at the ice margin reflects englacial drainage within the glacier, which, once established, enhanced heat and meltwater transfer. During the eruption, the englacial drainage passages collected pyroclastic material and increased in size, allowing hotter meltwater to be drained to the margins. This process enhanced ice melting and probably detached glacier ice from the bedrock, favouring ice mass instability.
Mixed avalanches on the northeast ice cap (Figs 2 and 3) have long been recognized as forming unique primary deposits at Nevado del Ruiz (Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990; Reference ThouretThouret, 1990; Reference Pierson and JandaPierson and Janda, 1994). These are poorly sorted and comprise ungraded mixtures of pumice and scoria lapilli, coarse lithics and coarse fragments of slush snow and ice (Fig. 4a). Such deposits are coarser-grained than pyroclastic surge, mudflow and hyperconcentrated stream-flow deposits, but finer-grained and better sorted than many pyroclastic flow and debris flow deposits. Mixed avalanches with a volume of 6–7×106m3 travelled up to 6 km from the ice margin, forming lobate deposits with steep fronts. An H/L ratio of 0.2 (where H is height and L is length) and the corresponding energy line (arc tanH/L) of 128 suggest an apparent basal friction coefficient of 0.2, similar to that for debris avalanches. The calculation of flow velocity is based on the equation derived by Reference Johnson, Brunsden and PriceJohnson (1984):
where R is the radius of the curvature of the bend, g is the acceleration due to gravity, S is the channel slope and B is the flow surface tilt angle perpendicualt to flow direction. Velocity reached >15ms–1 on the steep slope between the Arenas crater and the Azufrado head-valley (Fig. 4b), but decreased to 4–6ms–1 on the gently sloping lava flow forming the Azufrado-Plazuela middle ridge (Reference Pierson and JandaPierson and Janda, 1994). The mixed avalanche deposits were probably transported by unchannelized, dense gravity currents similar to high-concentration pyroclastic flows (Froude number > 1; Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995; Pierson and Reference Janda, Scott, Nolan, Martinson, Lipman and MullineauxJanda, 1994).
Post-1985 Evolution of the Nevado Del Ruiz Ice Cap
We have computed the post-1985 variations in the surface of the ice cap in order to re-assess lahar hazards and to distinguish the short-term effects of the eruption from the long-term impact of climate change. The Nevado del Ruiz ice cap, with surface area 25 km2 before and 21.3 km2 after the eruption, encompasses seven glaciers (Reference Thouret, Hammen, Salomons and JuvignéThouret and others, 1997). By 2003, the Nevado del Ruiz ice cap had lost as much as 52% in area and roughly 30% in volume (Table 1). Ice cap boundaries (Figs 5, 6 and 8) show that the rate of ice loss has been steady since 1986, with 13 m of vertical loss in 1987–88 and 8.8m in 1990–91, a rapid retreat also computed by Reference Linder, Jordán and ChristkeLinder and others (1994) between 1986 and 1991. Dead or regenerated ice bodies at the foot of lava scarps to the north and northeast, which existed until 1991, have disappeared since (Figs 5 and 6; Reference Ramírez and GuarnizoRamírez and Guarnizo, 1994). From 1986 to 1997, the average annual areal retreat of the ice cap reached 4.5%, i.e. 20–30m in altitude per year. The retreat has likely accelerated between 1997 and 2003. In addition to climate changes, the decrease in albedo caused by the November 1985 pyroclastics that accumulated until 1995 on the north, northeast, northwest, west and southwest ice cap induced a gain in solar energy and provoked a rapid retreat (Fig. 5).
Table 1. Distribution of the surface areas of each of the glacier watersheds of the Nevado del Ruiz ice cap before and after the 1985 eruption (Reference ThouretThouret, 1990), and comparison with the 1997 distribution (see Figs 5 and 6)
Based on the measured ice thicknesses and interpolated cross-sections, we estimated the ice cap volume to be 0.57± 0.2 km3 in 2000. Despite several sources of error attributed to the signal reception, location of antennae, interpolation of thickness ellipses and small-scale irregularities in bedrock topography, the accuracy of the ice thickness ranges between 2–5% (Reference Bauder, Funk and GuðmundssonBauder and others, 2003). There is a striking contrast between the ≥60m thick ice field, forming and blanketing the summit plateau, and the steep northern and eastern plateau margins covered by ≤30m of hard ice devoid of snow. Along the Nereidas Glacier–Arenas crater profile (Fig. 7b), the ice thickness, which is 49 m on average, increases up to 135 m towards the Arenas crater and to as much as 190m beneath the Nereidas Glacier to the southwest. The summit–Alfombrales glacier profile shows that the maximum ice thickness is 60 m south of the crater (Fig. 7b). The ice-free Arenas crater lies at the northern rim of a 2 km wide depression filled with 50–190m of ice, suggesting the existence of a summit caldera or a crater complex. For example, the bedrock topography reveals a 700m wide depression, probably an old crater, west of the Arenas crater.
Long-term effects due to eruption aftermath and climate change
By observing the southern and southwest glaciers, which were unaffected by the 1985 eruption, long-term effects attributed to climate change can be distinguished from impacts of the 1985 eruption. On the basis of Figures 5 and 6, greater ice loss occurred on glaciers that were more affected by pyroclastic flows and the tephra of the 1985 eruption (cf. fig. 3 in Reference ThouretThouret, 1990). Pyroclastic flows and surges left the ice surface abraded, scoured and gullied on the Nereidas glacier to the west, the Gualí and Azufrado glaciers to the north, the Molinos glacier to the northwest and the Lagunillas glacier to the east (Fig. 5). Glacier ice of the Lisa and Farallón to the northwest, the Azufrado to the north and the Plazuela catchments to the northeast was partly removed by additional ice mass failure, and snow-slab and mixed avalanches (Thouret, 1990). However, the loss rate has changed in time and space. The rate of retreat from 1987 to 1991 increased between 1991 and 1997. The 1997 boundary (Fig. 5) likely corresponds to the equilibrium line altitude (ELA) of the northern and eastern ice cap (between 5200 and 5100m a.s.l.), which is lower in the southern and western ice cap (4900–5000m a.s.l.). In contrast, no change is perceptible on the central flat area of the thick ice cap.
The southern and southwest glaciers are retreating at a slower rate than the western, northern and eastern glaciers. The extensive Recio glacier tongue in 1981 was already retreating in 1985. After 1985, it turned into a dirty glacier surrounded by dead ice at the foot of a U-shaped trough and a thin fractured head-valley glacier retreating toward the edge of the summit plateau. The glacier tongues, restricted to troughs (Nereidas, Recio) and to horseshoe-shaped amphitheatres (Azufrado, Lagunillas), eventually disappeared because they were steep (>358), thin (<30 m) and were fractured and destabilized during the 1985 eruption. In addition, these glaciers were covered by thin tephra layers, which decreased the albedo therefore favouring ice melting.
The March 2005 Aster image (Fig. 6) helps to confirm that the 1985 eruption caused long-term effects on the glaciers. The northeast–southwest aligned body of ice covered by thick snow southeast of Arenas crater contrasts with the retreating, discontinuous ice fields on the north, northwest and west ice cap margins. This post-1985 evolution has enhanced the asymmetric distribution of the ice cap that existed prior to the 1985 eruption (Reference ThouretThouret, 1990). The asymmetric distribution also reflects a structural asymmetry of the volcano that results from its bisection by the N208E strike slip Palestina fault (inset, Fig. 2).
Topography of the volcano summit unveiled by ice retreat and hazard assessment
Back-stripping of the ice cap on the DEM allows determination of the summit topography of Nevado del Ruiz (Fig. 8). In contrast with simple composite cones, Nevado del Ruiz is a compound edifice, which encompasses a cluster of lava domes, two small parasitic edifices (La Piraña and La Olleta; Fig. 2), and four U-shaped amphitheatres related to flank collapse and to faults (Figs 5, 6 and 8). We hypothesized previously that five lava domes encircled the summit plateau, whose flat surface was thought to be the top of a thick ice fill (Reference Thouret, Cantagrel, Salinas and MurciaThouret and others, 1990a). Although mantled by thin ice, the presence of lava domes was inferred from lava flows emanating out of prominently convex-shaped areas of the ice cap.
Fig. 8. Topographic DEM of the Ruiz summit bedrock without the ice cap. Volcano summit structures and landforms unveiled by the ice cap retreat can be compared with those inferred by Reference Thouret, Salinas and MurciaThouret and others (1990b, fig. 9) and Thouret (1999, fig. 4).
Ice retreat has since uncovered four domes surrounding the plateau, where 60–190m thick ice probably fills a summit caldera south and west of the active Arenas crater (Figs 5, 6 and 7). The central part of the summit plateau between the Arenas crater and the Recio glacier probably holds the largest volume of ice. A lava dome and a vent appeared in the area between Farallón-Lisa and Nereidas where the glacier disappeared following the 1985 eruption. The Arenas crater, now 750m across, shows asymmetrical rims: the northern wall is almost vertical whereas the southern rim shows a collapsed terrace and several arcuate fissures. The Arenas crater is enlarged towards the southwest, showing a narrow collapse structure that parallels the N10-trending strike slip Palestina fault zone, coinciding with the west rim of the Azufrado amphitheatre (Figs 5, 6 and 8).
Glacier retreat has also uncovered three horseshoe-shaped scars that cut back into edges of the summit plateau (Figs 5, 6 and 8). The most extensive Azufrado amphitheatre (5×1 km) is related to a sector collapse prior to the 1985 eruption (Reference Thouret, Salinas and MurciaThouret and others, 1990b). The deep, elongated and narrow Lagunillas and Recio U-shaped troughs with vertical scarps may correspond to two additional sector collapses. Alternatively, the Recio trough may be the southern end of a northnorthwest–southsoutheast aligned graben, which branches out from the Palestina fault in the area of the Arenas crater. The head-valley scar and hummocky topography in the trough floor suggests that repeated retrogressive landslides have contributed to amphitheatre formation in addition to glacier erosion. Landslides are favoured by hydrothermal alteration of lava flows in the Azufrado, Recio and Lagunillas head-valleys. Reference López and WilliamsLópez and Williams (1993) have stressed the role of hydrothermal weathering and fluids in destabilizing volcanic edifices.
Re-Assessing Lahar Hazards at Nevado Del Ruiz
Reassessment of lahar hazard is necessary because an area of 10 km2 of snow and ice still caps the active volcano, and debris flows occurred without an associated eruption ten years after the 1985 event. Analyses of post-1985 ice cap evolution and geometry enable us to reassess lahar hazards posed by ice and snow interactions with pyroclastic debris. In 1990, we identified the northwest and northeast areas of the ice cap as hosting the most unstable and rapidly retreating glaciers (Lisa, Farallón, Molinos, Gualí and Azufrado-Plazuela), which would eventually disappear. In 1997, the thick and less disturbed ice on the flattish summit of the volcano became separated from the thin, bare and fractured glaciers hanging on the steep-sided west, north and northeast edges (Figs 5, 6 and 8). Two groups of glaciers, which can fail or be swept by pyroclastic flows and surges, would likely feed future mixed avalanches and lahars: (1) the Azufrado and Lagunillas catchments, north and east of Arenas crater, where fractured glaciers have retreated on unstable cliffs cut in hydrothermally altered rocks; and (2) the Nereidas, Farallón, Gualí and Molinos ice fields west and northwest of Arenas, whose long gentle (<108) slopes can be scoured by pyroclastic flows and surges. Assuming that future pyroclastic flows and surges would efficiently interact with snow, firn and ice, lahars could be triggered from the Nereidas, Lagunillas and Recio glaciers.
However, the situation in 2003 was quite different from that of 1985, because the ice cap had lost about 52% in area and probably 30% in volume (Fig. 5, Table 1). On the one hand, the volume of snow and ice decreased to a point where the volume of meltwater would be smaller than that of 1985, unless the magnitude of the eruption increased or interactions were long lasting on widespread snow cover. On the other hand, small eruptions and seismic activity can still generate ice and snow avalanches on the northwest, east and northeast fractured ice cap margins.
Due to the increased retreat rate since 1997, dead or dirty ice masses (e.g. Recio, Gualí-Farallón, Nereidas and Lagunillas) are and will be abandoned near the ice cap margin (Fig. 5). Despite the ice retreat, this situation is still hazardous for four reasons: (1) dirty or dead ice represents a source of material and water for a mixed avalanche in case of eruption; (2) dead ice masses or moraines can dam water in front of rapidly retreating glaciers: subsequent lake breakouts may trigger mudflows and debris flows; (3) solar radiation on thin ice hanging on steep rocky slopes enhances ablation; and (4) seismic shaking and small eruptions can induce ice and rock avalanches. For example, a mass flow of ∽550 000m3 including ice blocks (Fig. 4d) travelled 6 km down the Rio Lagunillas valley on 15 January 1995 (Reference Guarnizo, Ramírez, Haeberli, Hoelzle, Suter and Frauen-felderGuarnizo and Ramírez, 1996). Ice- and rockslides were observed in the head-valley cliffs in April 1994 but the triggering factor may be attributed to an earthquake as the flow was first reported 15 min after a long tremor. The Lagunillas, Nereidas and Recio head catchments remain the areas most prone to lahar generation.
Discussion
Heat transfer rate for available snow, firn and ice surface
The case study of the 1985 Nevado del Ruiz eruption provides critical parameters (duration and area of contact between snowand hot pyroclasts, volume and characteristics of pyroclastic density currents and volume of meltwater) to estimate the efficiency of heat transfer from hot pyroclasts to snow, firn and ice (Table 2). Derived from the volume of lost snow, firn and ice at Ruiz, a volume of 22–43×106m3 of meltwater was released from 20–90 min (Reference ThouretThouret, 1990); the first lahars were heard about 20 min after the eruption started. Assuming a steady heat transfer with time and discarding the thermal effect of hot water that results from the melting (and external heat losses), the computed volume of meltwater implies a melting rate of 0.5–1.6mms–1 during 20–90 min (Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995). Such a melting rate is compatible with the removal of 2.3–8.2m on average (Table 2; Fig. 4e) of snow and firn on the most affected 18km2 of the ice cap (based on isopach maps: Reference CalvacheCalvache, 1990 and Thouret, 1990) (Tables 1 and 2; Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995; Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990).
Table 2. Distribution of the glacial drainage basins, of the pyroclastic impact (type and area) and calculation of meltwater production
Laboratory experiments have been carried out to reproduce the melting processes induced by a turbulent fluid (i.e. hot water) and a snow/firn ‘bedrock’, represented by polyethylene glycol (PEG) with a density of 1.13 g cm–3 (Reference DonnadieuDonnadieu, 1995; Fig. 9). Hot water contained in a bucket was conveyed at a constant rate through a pipe towards a sheet feed at the top of a flume (with a surface of 1400 cm2 and slope between 1.88 and 358) covered by PEG. The mixture of melted PEG and water (hereafter referred to as melt) was collected in a bucket at the bottom of the flume and measured. Computed parameters were water and ‘bedrock’ temperature, discharge, surface of erosion, slope of the flume and volume of melt. Discharge and temperature of water were kept constant for 20 experiments. Melt rate was controlled by temperature (51–778C) and flow duration (1.25–5 min), temperature of ‘bedrock’ (22–238; PEG melts at 458C), and in part by discharge. It was determined that slope did not play a major role in melt production. The surface area of PEG covered by hot water is more important than the volume of PEG to be melted. This confirms that the larger the contact area with the snow, firn and ice, the more efficient the mechanical and thermal scouring effect of the pyroclastic density currents, provided they are neither too thick nor too thin.
Fig. 9. Comparison between volumes of meltwater as computed from experiments by Reference DonnadieuDonnadieu (1995). (a) Conduction model: differences between curves are small for low temperatures but increase for high temperatures. (b) Turbulence model (see Reference Huppert, Sparks, Turner and ArndtHuppert and others, 1984): differences between curves are higher than in (a) except for the highest tested temperatures. Volumes of meltwater, either measured or computed, show a similar range of magnitude.
Thermal exchanges between the hot water flow and PEG have been measured according to two distinct models. For both pairs of curves depicted in Figure 9, volumes of melted mixture are normalized with respect to volumes of water (V PEG /V water)×100 during the experiment so that they remain independent of time. Experiments are ranked on the basis of the role of increasing temperature defined by (T water – T melt)=(T melt – T PEG).
The first conductive heat transfer model, based on the passive, slow and non-turbulent deposition of hot tephra-fall on snow and ice, results in a low melting rate of 0.015– 0.03 mms–1 (Reference MeierMeier, 1980; Reference Huppert, Sparks, Turner and ArndtHuppert and others, 1984; Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990; Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995). If a constant melting rate is applied to the 18km2 Ruiz snowpack (assuming an average snow/firn density of 0.6 cm3 g–1; Thouret, 1990), the released volume of meltwater would be far less (0.66–1.33 to ∽3–6×106m3 during 20 to 90 min) than the volume of liquid water estimated from the missing ice and snow (22–43×106m3; Thouret, 1990). The passive deposition of pyroclastic material would release a volume of meltwater even smaller than the volume of water (11– 22×106m3) inferred from the total lahar volume (Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990). Thus, the passive emplacement of hot tephra on snow cannot produce enough meltwater to feed large-scale lahars. To reach the computed volume of meltwater included in the 1985 lahars with a constant rate of 0.015–0.03mms–1, tephra fall should be applied to an ice cap of 234–470km2 or to the Ruiz ice cap during a period of 572–1144 min. As an example, the southern and southeast glacial areas of the Ruiz ice cap mantled by a tephra-fall deposit <10cm thick, yielded no significant meltwater and no lahar fed the south and southeast drainages (Figs 2 and 3).
Thermal exchanges between the hot water flow and PEG have also been measured according to the model of turbulence (Fig. 9b) used by Reference Huppert, Sparks, Turner and ArndtHuppert and others (1984). Results show that melt volumes produced by turbulent hot water (of Reynolds number 131–793) on the PEG flume are much larger than melt volumes calculated with the conduction model (Fig. 9a). Rates of PEG melt calculated with the turbulent model range from 0.3 to 0.9 mms–1, i.e. ten to sixty times faster than the ‘conductive’ melt rate. Firstly, turbulent hot water flow enhances the capacity of the flow to remove the melted PEG because hot water drags the surficial PEG layer as soon as melt is initiated. The uncovered underlying PEG layer is re-incorporated, as we envisage the motion of a convective pyroclastic flow over snow. Secondly, turbulence helps to keep temperature constant within the water flow, particularly at the contact with the PEG. Melt efficiency is enhanced by turbulence. The implications of a model of rapid heat transfer more appropriate to turbulent pyroclastic density currents over snow and ice are similar to those drawn from Walder’s experimental results (2000a). The turbulence increases thermal exchanges, maintaining a steady temperature within the flow, and enhances the surface of contact between the hot pyroclasts and snow, firn and ice. Thermal erosion generated by turbulent pyroclastic density currents is efficient enough to produce tens of millions of m3 of meltwater in a few tens of minutes, at a rate compatible with observations carried out at Nevado del Ruiz and elsewhere (Reference Thouret, Vandemeulebrouck, Komorowski, Valla and SlaymakerThouret and others, 1995; Reference PiersonPierson, 1999).
Applicability of Walder’s ‘thermal scour’ model to the Ruiz case study
The complex glacier–volcano interactions encompass mechanical and thermal effects (Fig. 10), accounting for efficient and rapid heat transfer of pyroclastic density currents (PDCs) to snow, firn and ice (Reference PiersonPierson, 1995, Reference Pierson1999; Reference WalderWalder, 2000a, Reference Walderb). All processes that generate meltwater are thermally driven. Snow is melted and vaporized in part; the vapour flux fluidizes the pyroclastic flow and makes it convect. Bubbling causes convection when a pyroclastic flow runs over snow (Walder, 2000, fig. 10), vapour and ingested air fluidize the flow, which engulfs all the meltwater produced if this flow is not too thin or too thick. Vigorous bubbling and convection at the interface between PDCs and snow are favoured by high temperature, small grain size and small layer thickness (Reference WalderWalder, 2000b). The unstable vapour-driven fluidization and convective motion scours and mixes snow into the PDC, resulting in a slurry. According to the fluidization domain (Walder, 2000, figs. 6 and 9), fluidization occurs only for a certain combination of sand grain size D s, emplacement temperature T e, and layer overburden pressure P s (proportional to layer thickness). Certain categories of PDCs scour and efficiently melt snow and firn if they are neither thin nor thick, and if they become inflated and fluidized as air and gas exsolved from pyroclasts, they are ingested in the basal flow layer. These PDCs encompass hot (>4008C) and relatively thin (<2 m) flows such as pyroclastic flows at Ruiz, as well as high-energy ‘dry’ surges that generate slush snow and slurries (e.g. Mount St Helens, 1980, 1983 in Pierson, 1999; Ruapehu, 1995–1996 in Reference Cronin, Neall, Lecointre and PalmerCronin and others, 1996), but also wet surges that transform into lahars owing to the segregation of their dense and wet part (Mount St Helens, 1980).
Fig. 10. Dynamic transformations between multiple gravity flows and generation of eruption-triggered lahars on snow-capped composite cones: a summary (‘thermal scour’, ‘vapour-induced particle convection’; after Reference WalderWalder, 2000b).
Analogue experiments carried out by Reference WalderWalder (2000a) have reproduced the effects of thermally driven erosion processes generated by hot pyroclastics dropped on shaved ice taken as snow substrate. His ‘thermal scour’ model (2000b) has been extended to polydisperse pyroclastic flows. Walder has drawn a series of fluidization curves, defining a domain of efficient heat transfer by analysing distinct grain sizes, overburden pressure and temperature upon emplacement. Melt rate is proportional to (P s =4) (overburden pressure approached by pyroclastic layer thickness)=D s (grain size given by M z and σϕ where M z is mean and σϕ is the sorting index).
We have reconsidered our Ruiz data in order to test the applicability of Walder’s experiments. This model emphasizes three parameters: flow temperature on emplacement T e, overburden load P s (approached by the flow layer thickness) and grain size D s. At Ruiz, flow temperature and velocity are unknown and grain size has been measured (M z and σϕ ) on two surge deposits and three flow deposits (Reference Pierson, Janda, Thouret and BorreroPierson and others, 1990); however, we lack the grain size parameters of additional categories of flows and of the key welded tuff around the crater (Table 3). We cannot apply Walder’s model fully to the Ruiz case study for the following reasons.
Table 3. Dataset relevant to the pyroclastic deposits on the ice cap and parameters relevant to the applicability of Walder’s thermal scour model
The analogue experiments deal with ‘static’ hot debris dropped on shaved ice. No relative shear motion is taken into account, and no purely mechanical scour has been modelled. At Ruiz, the widespread mixed avalanches and the ice removed from the northern, western and eastern ice edges (Figs 2 and 3) indicate that shear stress has played a major role on these areas most affected by pyroclasts and snow interactions. Basal shearing imparted by fluid drag and particle collisions, within the high concentration base of PDCs in rapid motion at the contact with the snow and firn cover, has a role in transferring kinetic energy from the grains to snow and firn. Entrained snow slab (Fig. 4b) and failed ice masses have fed mixed avalanches on the surface of the destabilized snowpack. The mechanical interactions have not been taken into account in Walder’s experiments (2000a) as there was no relative shear motion between sand and snow and no purely mechanical scour. Furthermore, heat transfer was limited at this stage as snow and ice blocks were included in mixed avalanches as observed three weeks after the event.
Reference WalderWalder (2000a) claims that ‘contact geometry of transitory deposits laid down on snow or ice cannot be used to make inferences about processes because an ice or snow substrate will be differentially melted by an overlying pyroclastic layer, even after that layer has come to rest’. At Ruiz, however, evidence for mechanical scour is inferred from erosional features carved in the ice nearby the Arenas crater in 1985 (Fig. 4f) and preserved beneath an unaffected mud layer. A rather unusual mixed surge deposit (Fig. 4f) consists of alternate layers of pulverized ice and lapilli layers, labelled as ‘a’, overlain by a layer of mud ‘b’ at the base of massive pyrolastic-flow deposits ‘c’.
Overburden pressure exerted by PDCs does a better job if the pyroclastic flow is not too thin or too thick. Pyroclastic surges have affected the largest part of the ice cap (at least 17 km2) and their extent encompasses that of mixed avalanches on the northwest, north, northeast and east flanks. These were too thin (1–20 cm) and not hot enough (as shown by mud, accretionary lapilli and the uncharred alpine hut) to have exerted any major interaction with snow if we consider figure 4a–d of Reference WalderWalder (2000b). Pyroclastic surges probably exerted a more mechanical and gravitational role on the steep-sided ice margins.
Temperatures must be in the range of 250–8008C to favour fluidization (Reference WalderWalder, 2000b, fig. 4). At Ruiz, the welded airfall tuff was the hottest pyroclastic deposit (probably above the temperatures indicated by Reference WalderWalder, 2000b, fig. 4) but this deposit is the coarsest and the thickest (up to 560 cm) and occupies the smallest area around the crater (see Table 3). The coarse grain size range and the overburden pressure must exceed the efficient range (5– 15 kPa) indicated by Reference WalderWalder (2000b, fig. 9).
The pyroclastic flows (50–200 cm thick each) have affected 9 km2 of the ice cap. They correspond to the most efficient agents in figure 8 of Reference WalderWalder (2000b). However, flow temperature and velocity are unknown at Ruiz. The grain size of these pyroclastic flows (1.91–11.88 mm) is much larger than the efficient range (0.25–0.75 mm) and more poorly sorted (2.33–3.65) than the best values (2) for fluidization obtained in the numerical tests of Walder. From this, only the pyroclastic surges at Ruiz (grain size M z = 0.05–0.5mm and sorting σϕ 1.39–3.62) correspond closely to the most efficient data used by Walder but they are much thinner than the threshold shown in his figure 8 (M d = 0.1–1.1 and σϕ = 2–4).
The thermal regime depends on specific parameters of the pyroclasts but also on the characteristics of the snow and firn cover (Reference Arons and ColbeckArons and Colbeck, 1995). Contributing factors including fluid drag and percolating meltwater in the snowpack depend on its permeability and density, which is variable in space and depth over the Ruiz ice cap. The processes that destabilize snow at the contact of the pyroclastic materials were analysed by Keys (1996, cited in Reference Manville, Hodgson, Houghton, Keys and WhiteManville and others, 2000); Reference Conway and RaymondConway and Raymond (1993); and Reference Manville, Hodgson, Houghton, Keys and WhiteManville and others (2000).
Conclusions
Several lessons have been learnt from the Nevado del Ruiz case study: (1) A moderate-sized VEI 2 eruption (10– 20×106m3 of pyroclastics) can trigger large-volume (108m3) lahars on snow-capped composite volcanoes. (2) Hot pyroclasts interacting with snow, firn and ice can rapidly transfer heat to produce an equal or greater volume of meltwater within a few tens of minutes. (3) Such meltwater volumes can generate catastrophic lahars and floods that can travel as far as 100 km from the ice cap. (4) On snow-capped composite volcanoes, the contact area between hot pyroclasts and snow or firn is more critical for heat transfer than ice volume. (5) Turbulent pyroclastic flows and surges are critical in producing meltwater as they erode and incorporate snow and ice. (6) Mechanical processes (fluid drag, mass failure, shear stress imparted by pyroclastic flows and particle collisions within the base of the flows) also play a role in destabilizing snowpack and ice mass. (7) The rapid retreat of the north, northwest and northeast glacier catchments of the Nevado del Ruiz ice cap suggests that eruption-induced effects have been superimposed on the glaciological consequences of climate change. (8) Multi-temporal remote sensing of an ice cap can help to diagnose potentially unstable slopes, as well as potential source areas for lahar generation.
In summary, the complex interactions of hot tephra with snow and ice on composite volcanoes include (Fig. 10): the thermal erosion of snow, firn and ice; the mechanical abrasion and incorporation of snow and ice particles into pyroclastic density currents; the turbulent regime of the pyroclastic density currents; and the englacial drainage producing tunnels or increasing their dimensions and scouring the ice margins at their outlets. All these processes contribute to the production of meltwater that eventually feeds lahars.
Acknowledgements
This is an outcome of the Ice and Magma Interaction Processes IMIP research project funded by the European Union (Contract ERB-CII*CT94–0139). We thank our Colombian colleagues of Ingeominas (OVS) and Universidad de Caldas (Manizales), as well as our colleagues F. Donnadieu and J.-L. Froger at Laboratoire Magmas et Volcans, Université Blaise Pascal. We thank V. Manville, V.E. Neall, and T.C. Pierson for their reviews, which led to substantial improvements of this paper.
Appendix A Radar Measurements
The radar system used is a standard low frequency type pulse system suitable for radar applications in temperate ice conditions. It consists of a transmitter and a receiver (Reference Funk, Guðmundsson and HermannFunk and others, 1994). The transmitter generates an electromagnetic impulse with a continuous frequency range between 0.1 and 3 kHz. The receiver is a digital oscilloscope that acquires and transmits the signal to a laptop. The two antennae are resistively loaded dipoles of 10 and 20m half length depending on ice thickness. Both antennae are symmetrical around a free point. The transmitter A and the receiver B are located on the ice surface at a typical distance of 30 m, where the signal can be reflected. The EM wave travels in the air between the transmitter toward the receiver as well as through the ice towards the glacier bed. The reflection of the EM wave at the glacier bed follows Snell’s law. The location of the reflection points is initially unknown. An ellipsoid of rotation gives all the possible reflection positions to a certain travel time between sending and receiving antenna. The antennas are deployed to receive the reflections of a vertical plane, perpendicular to the ice surface and to the longitudinal axis of the antennas. This allows reduction of an ellipsoid to an ellipse, which is located in this vertical plane (Reference Funk, Guðmundsson and HermannFunk and others, 1994). The distance which separates the transmitter A and the receiver B is AB. The path that the EM wave travels in ice after exiting the transmitter and hitting the rock surface at R is AR þ RB.
The time interval Δt measured by the oscilloscope is the time difference between the EM wave in air and its respective echo in ice, which is given by:
where v E = 168mms–1 and v L = 300mms–1 are the propagation velocities of the EM wave in the ice and in the air respectively (Reference Bauder, Funk and GuðmundssonBauder and others, 2003).
To calculate ice thicknesses, the semi-major axis a and the semi-minor axis b of the ellipse are then obtained with the following equations:
The known coordinates of the A and B points and the time interval Δt allow equations (A1)–(A3) to be solved for each pair of measurement points. The a and b parameters are computed for each ellipse (Reference Funk, Guðmundsson and HermannFunk and others, 1994). According to Reference Narod and ClarkeNarod and Clarke (1994), if the length between antennas is d (m), the ice thickness D (m) can be calculated using the following equation:
To determine the line which defines the bedrock, successive measurements of ice thicknesses were carried out and ellipses of each pair of points were calculated. The most probable glacier bed is the envelope of all ellipse functions. An error analysis which accounts for the uncertainties suggests that the accuracy of the ice thicknesses ranges between 2 and 5% (Reference Bauder, Funk and GuðmundssonBauder and others, 2003).
