Ice volumetric changes on active volcanoes in southern Chile

Abstract

Most of the glaciers in southern Chile have been retreating and shrinking during recent decades in response to atmospheric warming and decrease in precipitation. However, some glacier fluctuations are directly associated with the effusive and geothermal activity of ice-covered active volcanoes widely distributed in the region. The aim of this paper is to study the ice volumetric changes by comparing several topographic datasets. A maximum mean ice thinning rate of 0.81 ± 0.45 m a⁻¹ was observed on the ash/debris-covered ablation area of Volcan Villarrica between 1961 and 2004, whilst on Volcan Mocho the signal-to-noise ratio was too small to yield any conclusion. An area reduction of 0.036 ±0.019 km² a⁻¹ since 1976 was obtained on Glaciar Mocho, while on Volcan Villarrica the area change was −0.090 ± 0.034 km² a⁻¹ between 1976 and 2005. Glaciers on active volcanoes are therefore shrinking, mainly in response to climatic driving factors. However, volcanic activity is affecting glaciers in two opposite ways: ash/debris advection is helping to reduce surface ablation at lower reaches by insulating the ice from solar radiation, while geothermal activity is probably enhancing melting and water production at the bedrock, resulting in negative ice-elevation changes.

Introduction

Glaciers in the Chilean Lake District (38–41˚ S) have been retreating and shrinking during recent decades, presumably in response to reduction in precipitation and upper-atmosphere warming observed during the second half of the 20th century (Rivera and others, 2002). In the Chilean Lake District, surface atmosphere temperatures below 850 hPa geopotential height (approximately 1500 m a.s.l.: the minimum altitude of regional glaciers) showed a decreasing trend between 1933 and 1992, which was mainly explained by a strong atmospheric cooling between the 1950s and 1970s (Rosenblüth and others, 1997), while upper-atmosphere temperatures (850–300 hPa geopotential heights) showed a significant increase since the beginning of radiosonde measurements in 1958 (Aceituno and others, 1993).

Significant negative trends in precipitation were detected in the Chilean Lake District between 1930 and 2000, with a maximum reduction at around 39˚S, where a decrease of 450 mm in 70 years was measured (Quintana, 2004). These trends could be partially explained by recurrent rainfall deficits observed during the maximum intensity of El Niño events in the South Pacific Ocean (Montecinos and Aceituno, 2003). El Niño events have been more frequent and intense since the climatic shift of 1976 (Giese and others, 2002), but in contrast to the positive correlation with precipitation in central Chile (Rütllant and Fuenzalida, 1991) resulting in positive glacier mass balances (Escobar and others, 1995), a negative correlation in the Chilean Lake District, reinforced by upper-atmospheric warming, is generating negative glacier mass balances (Rivera and others, 2005).

Many of the glaciers in the Chilean Lake District are located over active volcanic edifices, where few glaciological researches have been done, especially on the relationship between glaciers and volcanoes. On active volcanoes, for instance, frequent ash deposition and lava flows may be responsible for ice area shrinkage, but if ash layers are thick enough, they may result in glacier growth, by insulating the ice from direct solar radiation, reducing ablation (Adhikary and others, 2002). In order to distinguish between climatic- and volcanic-related glacier retreat, it is necessary to analyze separated volcanoes, with distinctive and different recent eruptive histories and behaviours.

Volcanological studies in southern Chile have mainly focused on geological and geochemical characteristics of main active volcanoes, their eruptive history, effusion rates, lava and pyroclastic flow generation and dynamics, Holocene explosive eruptions, lahar deposits and volcanic hazards (e.g. Clavero and Moreno, 2004; Lara, 2004). However, most of the human volcanic-related casualties in the region have been caused by lahars rather than other volcanic processes such as lava flows (e.g. Naranjo and Moreno, 2004). Therefore, the study of the volume of water equivalent storage on active volcanoes is highly necessary. In this sense, ice and snow covering active volcanoes should be an important issue in terms of risk assessment, as during eruption events large volumes of snow/ice can potentially melt and flow downstream as lahars (e.g. Jόhannesson, 2002).

Glacier studies in Chile have predominantly focused on frontal variations, ice area changes and surface mass balance during recent decades (Rignot and others, 2003). Nevertheless, the effect of geothermal activity on ice-capped volcanoes is still unknown, and furthermore studies of the relationship between hot-spots underneath the glacier and ice volumetric changes are missing.

To answer some of these questions, two ice-covered active volcanoes of the Chilean Lake District were selected: Volcán Villarrica (39˚25^’12^" S, 71˚56^’27^" W; 2847 ma.s.l.) and Volcán Mocho (39˚55^’48^" S, 72˚01^’46^" W; 2422ma.s.l.) (Fig. 1), both located within 60 km in the same climatic region, but exhibiting contrasting geological behaviour during historical times (Clavero and Moreno, 2004; Echegaray, 2005; Pérez, 2005).

Fig. 1. Map of the Chilean Lake District showing location of Volcan Villarrica and Volcan Mocho.

Volcán Villarrica is considered highly active and characterized in historical times mainly by mild strombolian activity (González-Ferrán, 1995; Lara, 2004), permanent degassing and periodic explosions, with the lava lake remaining at a high level (90–180m below surface) at least since 1984 and very sensitive to the magmatic conduit activity (Calder and others, 2004). Concentrations of acid gases measured in the summit of the crater have been recognized as a hazard to climbers ascending the volcano, who may be exposed to concentrations above limits defined by the US National Institute of Occupational Safety and Health (Witter and Delmelle, 2004). Its eruptive history indicates (Lara, 2004) low frequency of large explosive eruptions (volcanic explosivity index (VEI): 3–4), but more than 50 eruptive events have been documented since 1558 (Petit-Breuihl and Lobato, 1994). The most recent violent eruption took place in 1971/72 when lava flows were generated, as well as 30–40 km h⁻¹ laharic flows (Naranjo and Moreno, 2004) descending toward Villarrica and Calafquén lakes (Fig. 1). Lahars produced by eruptions of Volcán Villarrica in 1948/49, 1963/64 and 1971/72 resulted in the deaths of more than 75 people (Stern, 2004), and they are considered the main hazard factor of the volcano. The volcano is covered by a glacier of 30.3 km² measured during this work in 2005, mainly distributed toward the south and east where the main glacier basin (Glaciar Pichillancahue-Turbio, 17.3 km²) composed of partially ash/debris-covered ice is located. The energy balance of this glacier has been monitored since 2003 (Brock and others, 2005), and global positioning system (GPS) as well as radio-echo sounding (RES) measurements were carried out in January 2005.

In contrast, the Volcán Mocho–Volcán Choshuenco system (Fig. 1) has not erupted since 1864 (Rodríguez and others, 1999; Echegaray, 2005), with the latest observations of activity located on the western side of the main cone of Volcán Mocho (González-Ferran, 1995) showing no signs of fumarolic activity at present. However, important pyroclastic and laharic fans of Holocene eruptive activity can be recognized on major flanks, especially to the north, east and west (Pérez, 2005). The volcano was almost totally covered by a glacier of 17.6 km² in 2004, mainly distributed toward the south and east, where the main glacier (Glaciar Mocho, 4.8 km²) has been monitored monthly since 2003 (Rivera and others, 2005).

The scope of this paper is the determination of recent glacier volumetric changes at both volcanoes and the analysis of geothermal properties that may be related to glacier responses. For that purpose, several satellite images and topographic datasets were collected, as well as several field campaigns carried out on both volcanoes, where RES and GPS systems were used to survey the glaciers.

Methods

Glacier extent

The glacier extensions at both volcanoes were digitally mapped from satellite images (Table 1) which were originally corrected based upon the orbital parameters provided by NASA. These images were later orthorectified using digital elevation models (DEMs) available for each study area yielding horizontal errors smaller than the pixel size (Table 1) on low-altitude and relatively flat areas. In spite of the above corrections, the images were still deformed in the vicinity of steep volcanic cones, requiring further resampling procedures based upon GPS data collected in the field. For that purpose, a GPS kinematic survey allowed mapping of the main crater areas of Villarrica and Mocho volcanoes, resulting in a better rectification of the images, yielding a few tens of metres of horizontal errors on the steep cones.

Table 1. Datasets

Colour composite images were generated by combining bands 1, 2 and 3N of Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) images and bands 1 , 3 and 4 of Landsat Multispectral Scanner (MSS) images. An automatic non-supervised classification was performed in order to distinguish between snow/ice and rock areas. In spite of that, and considering the huge number of fractals generated around the main cones and the underrepresentation of ash/debris-covered ice, the glacier boundaries were finally manually digitized on screen, based upon the colour composite images overlapped with the available contour lines. The resulting glacier extents in 1976 and 2004/ 05 were compared in order to measure area change rates. An error assessment was carried out for each change in glacier area, assuming a worst-case scenario where ice-margin delineation error yielded ±1 pixel size (Table 1) which is multiplied by the perimeter of the changed area (Williams and others, 1997).

Surface topography

Based upon aerial photographs acquired on 12 December 1961, Instituto Geografico Militar (IGM) of Chile photogrammetrically derived the first regular cartography of both volcanoes at 1 : 50 000 scale. The digitized contour lines of the cartography were interpolated using an inverse distance-weight method yielding a 90m pixel size DEM. The main problem with this cartographic method is realized on Volcan Villarrica where the cone was not covered by the IGM cartography due to stereo-matching problems between the aerial photographs (Fig. 2). Volcan Mocho was much better represented by contour lines every 25 m which covered the whole glacierized area.

Fig. 2. DEMs of Volcán Villarrica based upon different datasets used in this study. Areas with lack of data are shown in white, and areas where data were acquired are in greyscale.

Shuttle Radar Topography Mission (SRTM) data collected in February 2000 by NASA and the US Department of Defense were obtained at 90 m grid size, covering almost the whole of the volcanic cones, with minor gaps located in data shadows of steep flanks. In March 2004, the Jet Propulsion Laboratory (JPL)-NASA, USA, in collaboration with Centro de Estudios Científicos (CECS), Chile, carried out an airborne synthetic aperture radar (AirSAR) topographic mission on board a DC-8 aircraft, which surveyed several volcanoes and glaciers of southern Chile. The instrument on board this aircraft operates an interferometric synthetic aperture radar at C-band (5.6cm wavelength) and L-band (24 cm wavelength) frequency, allowing generation of several products, including DEMs based upon both frequencies. C-band topography is expected to be closer to the snow/ice surface, as previous investigations comparing AirSAR C-band to a surface reference (laser altimetry) suggested penetration depths of <1 m on temperate snow/ice, and greater (a few metres) penetration in cold/dry ice (Dall and others, 2001; Rignot and others, 2001). Although C-band seems to be more accurate for measuring surface topography, both DEMs (AirSAR C based upon C-band and AirSAR L based upon L-band) were compared to GPS data collected on site.

Geodetic GPS measurements were used to provide information on the glacier’s surface topography independent of/ complementary to datasets derived by aerophotogrammetric and remote-sensing techniques. Therefore, ice elevations were surveyed in differential mode during field campaigns at both volcanoes in 2004 and 2005 using dual-frequency receivers and antennae (JNS Lexon GD with MarAnt GD). In each case, a GPS station mounted on bedrock was used as local reference for stop and go as well as kinematic observation technique. The former was used to remeasure stake networks installed on both glaciers for monitoring mass balance and ice surface velocities. The latter was employed simultaneously for ice-thickness surveys using RES. The coordinates of both reference stations were determined using 48 hours of static GPS measurements in the framework of the most recent densification of the Chilean realization of SIRGAS (Geocentric Reference System for The Americas) in November 2004. Baselines were in the range up to 5 km, allowing the reliable resolution of ambiguities at double-difference level of the carrier-phase observations. GPS data were processed using Javad’s commercial software package ENSEMBLE. The resulting precision of the vertical component is assessed to be about 20 cm or better.

Ice thickness and internal structure of glaciers

Ice thickness, the internal structure of the ice and the bed topography characteristics were obtained at each volcano by means of a low-frequency impulse radio-echo sounding (RES) system (Plewes and Hubbard, 2001). The RES system consists of the transmission of a signal into the glacier being reflected at the glacier bedrock and returning to the ice surface where the signal is captured by a receiver. The ice thickness is determined based upon the two-way travel time (Rivera and others, 2001; Rivera and Casassa, 2002). Considering that the glaciers of the area are temperate (Casassa and others, 2004), a propagation velocity of 0.161 m ms⁻¹ (Macheret and others, 1993) was used.

Two RES configurations were applied: On Volcán Mocho, a profiling system was used, pulled by a snowmobile with ropes connecting the two sledges carrying the transmitter and receiver respectively. The transmitter used here was designed by The Ohio State University (OSU; F.E. Huffman, unpublished information), having an impulse voltage amplitude of 1600V peak to peak (Vpp). Resistively loaded antennae with a half-dipole length of 20m were used for propagation, yielding a central frequency of 2 MHz. A GPS antenna was mounted on the snowmobile collecting data each 5s for positioning RES measurements. On Volcán Villarrica a portable ground system was used, with at least three persons carrying the radar transmitter, receiver and GPS. Resistively loaded antennae were used with a half-dipole of 8 m, encased in tubular webbing, yielding a central frequency of 5 MHz. At Villarrica a transmitter designed by the University of Bristol, UK, was used, having a pulse amplitude of 800 Vpp (Gilbert and others, 1996). All antennae were designed based upon a model by Watts and Wright (1981) with the same gain setting (0.045 dB) and the same radiation resistance (400 Ω).

In order to analyze subglacial characteristics of each volcano, the bedrock reflection power (BRP) of RES data was determined based upon an equation derived by Gades and others (2000), where differences in the amplitude of the return could be related to basal properties when the transmitted energy now has variations, being the main factors affecting the RES signals, the geological variations of the bedrock (e.g. sediments vs hard bedrock) and the amount of water in the ice–bed interface (Gades, 1998).

Because the RES systems employed on the two volcanoes differed (central frequency and transmitter out power), we estimated the signal attenuation with depth in both systems using an equation derived by Bogorodsky and others (1985, p. 12). Based upon this equation and assuming no difference in the propagation media, ice thickness or BRP, the mean difference in attenuation between the two systems resulted in 14% higher values in the RES system used on Volcán Villarrica (5MHz and 800 Vpp) than in the system used on Volcán Mocho (2MHz and 1600 Vpp). That means the BRP values obtained in Volcán Villarrica were underestimated by 14% compared to those from Volcán Mocho.

Results

Areal changes

The oldest available satellite images of the studied volcanoes were acquired in 1976 (Table 1), presenting low-resolution but cloud-free conditions, allowing mapping of the main glacier of each volcano, namely Glaciar Pichillancahue-Turbio of Volcán Villarrica (Fig. 3) and Glaciar Mocho of Volcán Mocho (Fig. 4). These images were compared to the most recent and higher-resolution ASTER satellite images, yielding area change rates much larger on Villarrica than on Mocho (Table 2). The percentage of area change is nevertheless more significant on Volcán Mocho, where the studied glacier is smaller than the glacier on Volcán Villarrica.

Fig. 3. Ice-thickness measurements on Volcán Villarrica in metres. AirSAR C contour lines, each 100 m, are shown in light grey, and the boundaries of Glaciar Pichillancahue-Turbio are shown in black (2005) and dotted line (1976). Ice radar tracks are shown in grey, with tones representing depth values as shown in the scale.

Table 2. Recent glacier area changes, 1976–2004/05 ^*

* Dates are based upon acquisition times in Table 1.

Ice-elevation changes

The DEMs generated by JPL/NASA based upon bands C and L (AirSAR C and AirSAR L respectively) were acquired in March 2004 at 10 m pixel size resolution, covering only Volcán Villarrica (Fig. 2). Therefore, for comparison purposes, GPS data acquired in January 2005 were rasterized to this spatial resolution, yielding several thousands of common pixels. In spite of the difference in acquisition time, both surveys were obtained during the dry season and included the accumulation, ablation and debris-covered zones of the glacier, as well as some bedrock areas surrounding the glacier. The best comparison with the GPS data was obtained by means of the AirSAR C DEM, confirming previous estimations obtained in Greenland where this band presented lower penetration on temperate snow and ice (Forsberg and others, 2000).

This AirSAR C DEM was resampled to a lower resolution (90 × 90 m), in order to compare it to the other datasets. At this lower resolution, the mean difference between GPS and and AirSAR C data yielded a small bias (0.2 m) and a random error of 6.23 m (Table 3). The resulting bias could be related to the time separation between the acquisition date of GPS and AirSAR data, while the random error is comprised by minor inaccuracies of GPS data (0.2 m), and the errors associated with AirSAR C (6.2 m).

Table 3. Height differences on Volcan Villarrica

Rock areas were compared between datasets in order to estimate vertical rms errors (Table 4), yielding low biases and a maximum rms error of 18m for the DEM generated by IGM in 1961 (Table 5). This result is consistent with previous estimations carried out in Patagonia, where flat snow-covered zones lacking stereoscopic vision on the aerial photographs prevented stereo matching, generating holes in the cartography and high vertical errors in the surrounding areas (Rignot and others, 2003). The combined error between AirSAR C and IGM data resulted in smaller random errors, probably because the compared area did not include the main cone, as in the IGM data this part of the volcano has no topographic control.

Table 4. Height differences between datasets on rock areas

Table 5. Vertical rms error for each dataset

Ice-elevation changes were determined by subtracting topographic datasets, yielding signals higher than noise levels only on Volcán Villarrica, where the lower part of the glacier showed a thinning rate of 0.69–0.81 ma⁻¹ since 1961 (Table 6). Most of the compared area on Volcan Villarrica lies within the ash/debris-covered ablation area, and very few points were located above 2000 m a.s.l., which was the approximate location of the equilibrium-line altitude (ELA) in March 2005. On Volcan Mocho the available datasets covered both zones of the volcano, with an ELA at 1956 m a.s.l. in 2003/04 (Rivera and others, 2005), but no meaningful elevation changes were obtained for the accumulation or the ablation zone (Table 7).

Table 6. Ice comparison of IGM map and the other datasets at Glaciar Pichillancahue-Turbio^*

Table 7. Ice comparison of IGM and the other datasets at Volcán Mocho^*

Subglacial topography

On Volcan Villarrica, several RES profiles were measured in January 2005 (Fig. 3) with a total length of the radar system of 40m being carried on backpacks. These profiles surveyed both the ash/debris-covered area of the glacier and the snow-covered surfaces, from the margins of the glacier up to 2436 m a.s.l., yielding a mean thickness of 7 5 ± 4 m , the error being the mean difference between 663 crossing points. Most of the profiles were noisy due to the presence of internal layers, hyperbolic features and complex internal structure of the ice. Some of the most prominent internal layers were followed up to the surface of the ice (Fig. 5), where the boundary between snow and ash/debris-covered ice was located, confirming the internal structure detected in the RES profiles consisting of pyroclastic deposits originating from the volcano, being advected by ice flow and emerging on the ablation area of the glacier. These deposits are probably related to the large Pucon Ignimbrite eruption that occurred at 3700 BP (Clavero and Moreno, 2004), as evidenced by its characterictic juvenile material, formed by phenocrysts-rich basaltic-andesite cauliform and breadcrusted bombs (Clavero and Moreno, 1994; Clavero, 1996). In many surveyed areas, the subglacial topography was not visible or was confused by internal layers, requiring new and denser data to detect deeper ice. One of these areas is illustrated in Figure 5 where the subglacial topography is interrupted. These features were visible in all the records obtained from this sector of the glacier, suggesting that large crevasses obscure the bedrock returns. These crevasses could be related to a break or a large crater structure in the subglacial topography.

The mean ice thickness obtained in Volcan Villarrica represents areas where signals were sufficiently clear to be distinguished from internal layers. Unfortunately, we have not been able to measure the total thickness of the ice in many sectors below the ELA of the glacier (Fig. 3).

On Volcan Mocho, several RES profiles were measured in September 2004 (Fig. 4), using a snowmobile pulling a 120 m long system while moving on snow surfaces between 1900 and 2180 m a.s.l. Neither the steep flanks of the volcano nor the steep cliffs in the ablation area were covered. The maximum ice thickness was obtained in one tongue of Glaciar Mocho (Fig. 6) yielding 270 m (Table 8), and the mean thickness yielded 138 ± 9 m, being the error of the mean difference between 111 points where thickness was measured in crossing points of different RES profiles. Two main subglacial U-shaped valleys were distinguished (Fig. 6), which are associated with the main tongues of Glaciar Mocho. RES records were generally clear with few internal layers, but several hyperbolic features were detected near the surface of the ice, most of them presumably related to ice crevasses (Fig. 6).

Fig. 4. Ice-thickness measurements on Volcán Mocho in metres. SRTM contour lines, each 100 m, are shown in light grey, and the boundaries of Glaciar Mocho are shown in black (2004) and dotted line (1976). Ice radar tracks are shown in grey, with tones representing depth values as shown in the scale.

Table 8. Ice-thickness measurements on Villarrica and Mocho volcanoes

Bedrock reflection power (BRP)

The BRP values obtained on Volcan Villarrica were logarithmically distributed, similar to the results obtained by Copland and Sharp (2001), with attenuation increasing with depth. In spite of this attenuation, some areas showed higher BRP values possibly in connection with enhanced presence of meltwater around hot-spots or areas where geothermal activity is presumably higher (Fig. 5). On Volcán Mocho (Fig. 6), BRP data were also logarithmically distributed, but when they were compared to the values obtained on Volcán Villarrica (approximately at the same thicknesses), BRP values were smaller, in spite of being overestimated in Mocho due to the lower central frequency used in this volcano. The smaller BRP values at Mocho can be related to smaller amounts of meltwater at the glacier bed, presumably due to less significant volcanic activity.

Fig. 5. Topographic profile A–A⁰ (see location in Fig. 3) showing surface (light grey indicates snow-covered area, light black indicates ash/debris-covered area) and subglacial topography (black line) of Glaciar Pichillancahue-Turbio of Volcán Villarrica. In the middle is the radar non-migrated corresponding profile with subglacial returns in white. At the bottom are BRP values obtained along this profile. The arrow indicates appearance of ash/debris-covered layer on top of the glacier.

Fig. 6. Topographic profile B–B⁰ (see Fig. 4 for location) showing surface (light grey) and subglacial topography (black line) of Glaciar Mocho on Volcán Mocho. In the middle is shown the corresponding non-migrated radar profile, including in white the subglacial topography. At the bottom are BRP values calculated along track.

In order to compare BRP values obtained at Volcán Villarrica independently of the ice thickness, a minimum-square polynomial curve was applied, yielding a normalized distribution where much higher BRP values are detected in several spots along a line connecting the main crater of the volcano and Pichillancahue valley (Fig. 7). This line is coincident with a regional geological lineament connecting Villarrica, Quetrupillán and Lanín volcanoes (Fig. 1), which has been interpreted as an old, deep fracture in the continental crust used by Pleistocene–Holocene magmas to reach the surface (López-Escobar and others, 1995).

Fig. 7. BRP data calculated along tracks acquired on Volcán Villarrica in January 2004. For comparison purposes, BRP values have been normalized based upon a minimum-square polynomial curve. Units are standard deviations. Dotted line shows Glaciar Pichillancahue-Turbio basin. Notice the higher density of higher BRP values which are presumably hot-spots along a line connecting the main crater and Pichillancahue valley.

Discussion

The studied glaciers have been affected by similar climatic conditions within the Chilean Lake District. However, due to the differences between the volcanoes, glacier surfaces in the ablation areas are characterized by bare ice on Volcán Mocho and ash/debris-covered ice on Volcán Villarrica, resulting in different glacier responses to the same climatic driving factors. Considering that Volcán Mocho is less affected by volcanic activity and exhibited fewer areas with high BRP, it was possible to presume that its glaciers were responding more directly to climate changes with area shrinkage and steepest mass-balance gradient. This argument is supported by the behaviour exhibited by Glaciar Mocho, where area reduction has taken place without significant ice thinning. In contrast, on Glaciar Pichillancahue-Turbio, the glacier is experiencing both thinning and area reduction, in spite of presenting a thick layer of ash and debris (thicker than 1m in places) covering most of the ablation area (Fig. 7), which insulates the ice from direct solar radiation. In spite of this thick ash layer, in some areas where the glacier is more crevassed, backwasting seems to be an important process on steep ice walls.

Apart from the possible volcanic component discussed above, which is affecting Glaciar Pichillancahue-Turbio, the main glacier variations are driven by climate change or decadal atmosphere/ocean oscillations (i.e. the 1976 shift; Giese and others, 2002). However, not all changes are a direct response to warmer/drier conditions, as an important role is played by feedbacks triggered by climatic changes. Among these feedbacks, the ice surface elevation and glacier length responses are the most important.

In order to test these feedbacks, Equation (1) below, modified from Raymond and others (2005), is used to account for the contribution from climate change described in terms of ELA changes (first term on the righthand side), the elevation feedback caused by the ice-elevation changes (second term) and the glacier area feedback caused by the area changes (third term). These three contributions are a standard consideration for most non-calving glaciers (Els-berg and others, 2001).

(1)

where Δb) is the average change in mass-balance rate per unit area of glacier, A_g is the glacier area measured relative to a reference climate and geometry (A_g and h in 1976), G is a scaling factor with time associated with effective vertical gradients in ice equivalent thickness balance rate, Δh is the average ice-elevation change experienced by the glacier relative to a reference geometry (1976), H is the difference in elevation between the ELA and terminus of the glacier and Δ A _g is the area change experienced by the glacier.

Ice-elevation changes have been calculated since 1961, but area changes only since 1976. Due to the lack of more detailed data between 1961 and 1976, we have assumed a similar geometry in those years in order to calculate the main feedbacks. In order to estimate the ELA migration experienced by the region between 1976 and 2004/05 (29years), we have considered an ELA at 2000 m a.s.l. in 2004/05, an upper-atmosphere warming of 0.023˚Ca⁻¹ at 2000ma.s.l. since 1976 (Bown and Rivera, in press) and a regional temperature lapse rate of 0.0065˚Cm⁻¹ (DGA, 1987), resulting in a vertical rise of ∼100 m in the ELA since 1976 (0.023 x 29/0.0065).

Ice thinning was detected only on the ablation area of Glaciar Pichillancahue-Turbio (Volcán Villarrica). Based upon direct observations and a satellite image acquired in 2005 (Table 1), it was possible to estimate the ELA of the glacier at 2000 m a.s.l., defining an accumulation-area ratio of 0.38, yielding a mean thinning of 0.5 m a⁻¹ for the whole glacier.

The scaling factor (G), representing the glacier mass-balance gradient, has not been measured in the region, but a value of 0.015 has been obtained on nearby Volcan Mocho (Rivera and others, 2005), which was also applied to Volcan Villarrica. Similar values were calculated on temperate glaciers located in Patagonia (Naruse and others, 1995).

The above feedbacks have different significances: the elevation feedback is traditionally important (Raymond and others, 2005), but area shrinkage and ELA migration seem to be the main factors in the volcanoes studied here (Table 9). The average change in mass-balance rate relative to the 1976 reference climate and geometry resulted in the same value for both glaciers (−0.9 m a⁻¹), suggesting a climate-driven response, but the thinning on the ash/debris-covered ablation area of Volcán Villarrica is unlikely to be totally related to these climatic conditions.

Table 9. Main parameters used for Equation (1)

Between 2004 and 2005, ice melt beneath ash depths of 0.16–0.32m was recorded at four stakes located in the ash/ debris-covered ablation area of Glaciar Pichillancahue-Turbio, resulting in superficial thinning in the range 0.42– 0.56 ma⁻¹. Monitoring of ice melt beneath ash depths in the range 0.005–0.13 m, and on bare ice, during the mid-ablation season identified the critical ash depth between melt enhancement and melt reduction as <0.01m, demonstrating the extremely low thermal conductivity of the supraglacial material. For comparison purposes, on the lower bare-ice area of Glaciar Mocho an annual ablation of 4–8 m was measured in 2003/04 (Rivera and others, 2005). In this context, most of the surface ablation taking place on Glaciar Pichillancahue-Turbio must be related to backwasting of steep slopes not able to support ash or debris.

At Volcán Villarrica the amount of annual ablation of the ash/debris-covered ice is approximately two times smaller than the long-term ice-elevation change rate experienced by the glacier, suggesting that subglacial melting due to volcanic activity is responsible for at least one-third of the thinning, assuming no changes in ice velocities.

The ice-thickness measurements carried out on both volcanoes support the above suggestion, as very low BRP values were obtained on Volcán Mocho, whilst higher BRP values were obtained in several places on Volcán Villarrica (Fig. 7), which are presumably hot-spots related to high geothermal activity along a geological lineament connecting three active volcanoes in the area (Lόpez-Escobar and others, 1995; Lara, 2004). If these spots enhance subglacial melting, they could also be responsible for the observed thinning on the ash/debris-covered area.

Another possible explanation is related to a dynamic component (e.g. changes in ice velocity), but the grounded and non-calving condition of the glacier front suggests minor contributions from this factor. Another factor that could be responsible for higher thinning on this volcano is the edifice deformation commonly associated with magma rise from a deep magmatic chamber towards shallower levels (i.e. Voight and others, 1981). This process affects glaciers on volcanoes by increasing crevassing, which could result in enhanced melting (Fountain and others, 2005). Before, during and after the last main eruption of Volcán Villarrica in 1984/85, several new ice crevasses were detected on the upper flanks of the main cone (Fuentealba and others, 1984). Volcán Villarrica has an almost permanent lava lake within its crater at least since 1984 (Clavero and Moreno, 2004). The activity of this lava lake commonly increases in springtime, showing more explosions due to magma rising towards the upper parts of the volcanic system. This activity is usually accompanied by an increase in ice crevasse formation on the upper flanks of the volcanic edifice (Clavero and Moreno, 2004).

Conclusions

Based upon all available measurements, and considering a mean ice density of 900 kgm⁻³, the total volume of water equivalent storage on the two volcanoes was 4.2 ± 1.8 km³ w.e., the error being comprised by the uncertainties in area delineation and ice-thickness measurements. This volume is much smaller than the previous estimation for both volcanoes, 10.7 km³ (Rivera, 1989), which was calculated using a mean ice thickness based upon the area rank of each glacier because of the lack of field measurements. The previous estimations were based upon aerial photographs from 1961, so part of the smaller volume detected in 2004/05 is due to area shrinkage. This more accurate determination of volume storage on these active volcanoes will certainly improve volcanic risk assessments.

In spite of recent upper-atmosphere warming, Glaciar Mocho has shown no significant ice-elevation change, the glacier responses being restricted to area shrinkage presumably in response to the ELA migration. On the other hand, Glaciar Pichillancahue-Turbio has suffered frontal retreats and ice thinning two times higher than the annual surface ablation experienced by the ash/debris-covered area of the glacier, suggesting an important contribution to melting from geothermal activity at the bedrock of the glacier. In this sense, the volcanic activity on Volcan Villarrica is affecting the glacier in two opposed ways: by insulating the ice with ash and debris, resulting in reduced surface ablation, and by enhancing subglacial melting due to geothermal activity, resulting in higher thinning than in non-active volcanic environments. The effect of geothermal activity plays a more important role since Glaciar Pichillancahue-Turbio is experiencing thinning.