4.1. Mass-balance and climate modeling

The actualistic principle (assumption that the modern processes were the same and active also in the past) is used as a basis and prerequisite for mass-balance modeling as a function of climate. The respective empirical-statistical models for this area were established by Reference VuilleVuille (1996) and Reference KullKull (1999), and are listed in Table 2A and B.

Table 2.

(A) Parameterization for the daily, annual cycles, amplitudes of the climate and correction factors for temperature and global radiation. (B) Empirical–statistical models developed in the research area for the mass-balance modeling. (C)Equations for the ice-flow and mass-flux calculation for the considered cross-sections in the two case-study glaciers

The modeling scheme is shown in Figure 2 (for a detailed discussion see Reference KullKull, 1999). Annual and daily amplitudes, mean values and lapse rates were calculated from station data for the parameters temperature (T), relative humidity (RH), wind velocity (W), global radiation (G), cloudiness (C) and precipitation (P). Correction terms for ΔT as a function of C and P, and ΔG (as a function of C and “topography” (shielding)) were determined individually with climate data, and allowed us to parameterize the local climate on an hourly basis for the summer (October–March), and on a daily basis for the winter (Table 2A; Reference KullKull, 1999). In the Late-glacial climate scenario, G is corrected for the changing Earth orbital parameters (Reference Blatter, Funk and OhmuraBlatter and others, 1984). Reference AmmannAmmann (1996) established a correlation model for this area which describes the correlation between “cloudiness” and “precipitation”. As shown in Figure 3, the model climate averaged for the winter and summer compares well with the station data. Melt, sublimation and accumulation rates are quantified as a function of the respective climate conditions and the degree-day factor (DDF) which is variable and depends on snow characteristics, climate and season (Reference VuilleVuille, 1996; Reference KullKull, 1999; Fig. 2; Tables 1 and 2A). DDF is calibrated for the dry cold winter conditions in the research area (Reference VuilleVuille, 1996) and for humid warm summer conditions on Glaciar Zongo (Reference SicartSicart and others, 1998; Reference KullKull, 1999; Reference Wagnon, Ribstein, Francou and PouyaudWagnon and others, 1999). In order to account for the influence of the former glacier on the local climate, the air temperatures in our test areas were corrected using climatic data from Glaciar Zongo, where Reference SicartSicart and others (1998) and Reference Wagnon, Ribstein, Francou and PouyaudWagnon and others (1999) compared air temperatures on the glacier and on adjacent lateral moraines. Thus a value of dT = −2°C is applied to calculate local air temperatures on glaciated surfaces from station data in our currently ice-free test areas (Reference KullKull, 1999). Subsequently, the mass-balance terms are calculated for individual altitudinal segments of the glacier, and the mass flow is calculated after Reference OerlemansOerlemans (1997; Table 2C) for given cross-sections within the reconstructed glacier bed. Finally, the paleoclimate scenario is tuned (iteration, Fig. 2) such that the model glacier is in steady state, in equilibrium with the paleoclimate, matches the observed glacier geometry and thus fulfills two conditions:

1. (1) total overall mass balance of the glacier (annual average) = 0.

   

   (1)

   where Δh are the elevation zones (m a.s.l.) of the reconstructed glacier surface, b_h is the specific annual balance (kg m⁻²) of a given altitudinal zone and A_h is the glacier surface (m²) of a given altitudinal zone.

2. (2) difference between mass influx into a cross-section and the total mass balance below the cross-section (DMM) = 0 for individual cross-sections in the glacier valley.

   

   (2)

   where DMM_cr1/2 is the difference (kg a⁻¹) between annual mass influx in a cross-section and mass balance below the cross-section, h _cr1/2 is the altitudinal zone (m a.s.l.) of a cross-section, and M_cr1/2 is the mass influx (kg a⁻¹) in a cross-section (Table 2C).

Fig. 2.

Flow chart of the climate–glacier model (Reference KullKull, 1999) with cloudiness C, temperature T, precipitation P, relative humidity RH, global radiation energy G and wind velocity W. The degree-day factor (DDF) is variable and depends on the respective seasonal climate conditions. C regulates RH, P, T and G to a considerable degree, whereas the accumulation and ablation modules are largely controlled by T, P, G, W and RH. Thus the principal climate parameters which govern the climate scenario and control the model are C and T.

Fig. 3.

Comparison of precipitation and temperature values (monthly mean) between station Chungara Ayata (18°14′ S, 69°07′ W; 4570 m a.s.l.; Table 1) and the model climate for the same location. The model climate data were calculated using the station data (Table 1) and the parameterization shown in Table 2A.

Fulfillment of both equations results in steady-state conditions of a glacier. Figure 4a and d show the temperature–precipitation diagrams for both case-study areas with numerical solution for Equation (1) (overall mass balance = 0), which is a prerequisite for the subsequent calculations. Figure 4b, c, e and f show the respective differences between the mass influx and the mass balance below each of the two considered cross-sections (DMM values Equation (2), fulfillment of Equation (1) as calculated using the climate scenarios where the temperature–precipitation combinations fulfill condition 1 (Fig. 4a and d)). For steady-state conditions of the glacier, DMM values must be zero (Equation (2)). In general, DMM decreases with increasing precipitation and temperature, because higher precipitation requires higher temperatures in order to fulfill Equation (1). Also mass-balance gradients increase, resulting in more negative mass-balance values below the considered cross-sections.

Fig. 4.

(a, d) Temperature–precipitation diagrams for both case-study glaciers, showing the numerical solutions for total mass balance = 0 (Equation (1)). P and T refer to 4000 m altitude. (b, c, e, f) Differences between the mass influx and the total mass balance below the respective cross-sections (DMM) as a function of temperature (b, e) or as a function of precipitation (c, f). DMM = 0 in both cross-sections fulfills Equation (2) for steady-state conditions of a glacier and provides the numerical solution for the P climate scenario which matches the observed glacier geometry.

Because the bed geometry at the considered cross-sections is field-measured, the only parameters not defined for calculating the mass influx into a cross-section are the ice-flow parameters. Using standard ice-flow parameters (Table 2c), it appears from Figure 4 that DMM_cr1/2 = 0. Figure 4b, c, e and f also show that there is only one numerical solution (DMM_cr1 = DMM_cr2 = 0) with one particular combination of “temperature” and “precipitation” (and all the other dependent climatic variables (Fig. 2)) which fulfills Equations (1) and (2), matches the observed paleoglacier geometry and is thus the most likely paleoclimate reconstruction (for discussion see section 5).

4.2. Sensitivity of the model

It appears from Figure 4b, c, e and f that the uncertainty of the paleoclimate scenario is ultimately determined by the DMM values of the cross-sections (Equation (2)), which in turn are controlled by two groups of data: (i) geometrical data based on field measurements and (ii) ice-flow parameters from experimental data (standard values, Table 3a).

Table 3.

(A) Geometry of both case-study glaciers with the altitudes of cross-sections 1 and 2cr _1,2, the lengths of the tongue below the cross-section X, ice height at the cross-section H, cross-section area A, mean valley slope β, ice-flow parameters for basal sliding f _s , and internal deformation f _d . Using the best-fit paleoclimate scenario (Table 4), the combination of the parameters yields the exact numerical solution for Sensitivity of to ±5% change in the geometrical parameters (ice height at cross-section H, bed width at cross-section B, mean valley slope β) and ice-flow parameters (basal sliding f _s and internal deformation f _d )

We allowed an error of ±5% for the mapped geometrical data (ice height H and bed width B at two cross-sections, valley slope β) and for the standard ice-flow parameters (f _d, f _s). As shown in Table 3b, ±5% changes in geometrical or flow parameters result in ΔDMM_cr1/cr2 (DMM_cr1/cr2 ≠ 0) variations between +0.5/−0.6 × 10⁶ kg a⁻¹. Thus we allow this range of ΔDMM_cr1/cr2 as a possible numerical solution for a paleoclimate scenario (Table 4; shaded range in Fig. 7a and b).

Table 4.

Best-fit late Pleistocene climate scenario and glaciological conditions for both case-study glaciers calculated for 4000 m a.s.l.

Whereas the geometry of the glacier bed at the cross-sections is well defined, the ice-flow parameters remain debatable to some extent. However, using extreme changes in the order of ±50% from standard ice-flow parameters results in an increase of ΔDMM_cr1/cr2 values to +18 × 10⁶ kg a⁻¹ and −10 × 10⁶ kg a⁻¹ for both test areas, and no numerical solution for Equations (1) and (2) is found within the observed glacier geometry. Only ice-flow parameters within the range of ±5% from standard values yield numerical results for Equations (1) and (2) which are compatible with the observed geometry (range ±5%).

In order to validate our empirical–statistical approach, we compared our model (mass-balance modeling) with a physical energy-balance model. H. Oerlemans (personal communication, 1998) used our modeled best-fit paleoclimate scenario for El Tatio as input climate scenario for his physical mass-balance model and showed that the resulting mean annual mass balance of the whole glacier is close to zero (fulfilling Equation (1)) and thus in agreement with our results.

5.1. Modern glaciological conditions

In a first step, we tested the suitability of our climate-glacier model for modern climatic and glaciological conditions in the case-study areas. The results (Fig. 5) compare favorably with the field data and show the major characteristic glaciological features of today. In the southern case-study area of El Tatio at 22° S (Fig. 5a), the calculated mass balance is negative in all the altitudinal zones up to 6000 m. This is in agreement with field observations, confirming that the modern climate is too arid to sustain glaciers even in the continuous permafrost belt above 5600 m. Peaks as high as 6176 m (Volcán Aucanquilcha (21°13′ S) and Cerro Pili (23°17′ S; 6040 m a.s.l.) are currently free of glaciers and perennial snowfields, so no ELA currently exists. In the slightly more humid and warmer northern test area at Choquielimpie (18° S, Fig. 5b), the annual mass balance changes to positive values above altitudes of 5850 m, suggesting that there is enough moisture to sustain glaciers at high altitudes. However, the humidity is still at critically low levels. Indeed, the highest peaks in this area such as Volcán Sajama (6542 m) and Parinacota (6200 m) show ice caps with modern ELAs at about 5800 m (Reference HastenrathHastenrath, 1971; Jordan, 1991), and glacier tongues as low as about 5200 m on the southwest-facing slope (comparable to the case-study glacier at Choquielimpie). Interestingly, peaks at about 6000 m (such as Guaillatiri (6024 m)) still exhibit small snow/firn patches, whereas nearby Volcán Tacora (5980 m) and all the other peaks below 5800 m are currently free of glaciers.

Fig. 5.

Model results for mass-balance terms at different altitudes for both hypothetical case-study glaciers under modern climate conditions (station data, Table 1). In the southern test area at 22° S, El Tatio, the mass balance is negative up to 6000 m a.s.l., confirming the lack of glaciers even in the continuous permafrost belt. In the northern test area at 18° S, Choquielimpie, the mass balance turns to positive values above 5850 m altitude, effectively reflecting the existence of small glaciers on peaks higher than 6000 m. Cerro Choquielimpie (5327 m) is currently ice-free.

We conclude from the good comparison between calculated and observed modern ELAs that the statistical–empirical models which link climate data with sublimation, melt and accumulation terms (Reference VuilleVuille, 1996; Reference KullKull, 1999) are suitable for both case-study areas.

5.2. Late-glacial climate conditions

Here, we use the climate–glacier model to reconstruct paleo-climatic conditions from two late Pleistocene case-study glaciers with well-known geometries of the former ice mass, assuming steady-state conditions of the glacier. Table 4 lists the results of the paleoclimate reconstruction (best numerical solution) reduced to 4000 m a.s.l., and Figure 6 shows the respective mass-balance terms in different altitudinal belts. In this first step the paleoclimate scenario is based on the assumption that the seasonality of precipitation is comparable to the modern situation (i.e. predominantly tropical summer precipitation).

Fig. 6.

Mass-balance terms for both case-study glaciers, calculated with the best-fit paleoclimate scenariofor late Pleistocene (Late-glacial) times (Table 4). The ELA (specific annual mass balance = 0) is at 4830 m (El Tatio) and at 4900 m (Choquielimpie). Comparison between summer (b, e) and winter (c, f) shows significant differences. For (sub) tropical glaciers, summer is the main period for accumulation and ablation, whereas accumulation and ablation (mainly sublimation) are strongly reduced during the cold, dry winter season.

The results from the case-study areas at El Tatio and Choquielimpie are in good agreement. Accounting for a ±5% margin of uncertainty in the geometry of the glacier, the model suggests a depression of the mean annual temperature of 4.4 (− 0.1/+ 0.2) °C, and an increase in annual precipitation rates in the order of 840 (− 50/+ 10) mm for the northern test area at 18° S. The results for the southern test area at 22° S show a mean annual temperature depression in the order of 3.2 (±0.1) °C and an annual precipitation increase of 1000 (− 10/+ 30) mm (see Table 4). The vertical distribution of the annual mass-balance terms (Fig. 6a and d) suggests a Pleistocene ELA at 4900 m (Choquielimpie) and at 4830 m (El Tatio), which is about 50 m higher than the Pleistocene ELA inferred from field data (starting-point of lateral moraines; Reference Jenny and KammerJenny and Kammer, 1996). In the northern area, the Pleistocene depression of the ELA amounts to approximately 1000 m, which is comparable with the range 800–1200 m reported for adjacent areas in Bolivia and Peru (Reference HastenrathHastenrath, 1967, Reference Hastenrath1971; Reference GrafGraf, 1991; Reference Klein, Seltzer and IsacksKlein and others, 1999). The model also suggests an accumulation-area ratio (AAR) of 0.61–0.66 for the Pleistocene glaciers, which agrees with field data (Reference Jenny and KammerJenny and Kammer, 1996).

Figure 6b, c, e and f show the mass-balance terms separated for the summer and winter seasons. Typically for (sub)-tropical glaciers, the main ablation period coincides with the main accumulation period during summer. On the other hand, this effectively reflects the cold, dry climate conditions during winter, when the major proportion of the energy available for ablation is consumed for sublimation. Sublimation is approximately eight times more energy-intensive than melting and is the dominant ablation term at high altitudes during winter in this climate. Thus Pleistocene summer mass-balance gradients in the ablation area (about 12 kg m⁻² m⁻¹) are greater than winter values (about 4 kg m⁻² m⁻¹) and compare favorably with values for (sub) tropical glaciers (Reference Kaser, Hastenrath and AmesKaser and others, 1996; Reference Kaser and GeorgesKaser and Georges, 1999). The annual mean value is 8.8–9.2 kg m⁻² m⁻¹. Interestingly, both seasonal extremes (∼12 kg m⁻² m⁻¹ during summer, ∼4 kg m⁻² m¹ during winter) result in a mean annual mass-balance gradient and, finally, in an AAR which is close to values for alpine glaciers in oceanic climates of mid-latitudes (AAR 0.66, mass-balance gradient 7.7 kg m⁻² m⁻¹; Reference KaserKaser, 1996). These modeling results are in agreement with the field data (AAR of 0.66) by Reference Jenny and KammerJenny and Kammer (1996).

The inferred Late-glacial temperature depressions of 4.4 (− 0.1/+ 0.2) °C for the northern test area at 18° S, and 3.2 (±0.1) °C for the southern test area at 22° S compare favorably with results from glaciological studies in adjacent areas of Bolivia and Peru (Reference SeltzerSeltzer, 1992, Reference Seltzer1993; Reference Clayton and ClappertonClayton and Clapperton, 1997). However, our results as well as the findings of Reference Clayton and ClappertonClayton and Clapperton (1997) suggest substantially higher precipitation (ΔP between +600 and +1000 mm) than the paleo-climate reconstruction from water-balance models of paleo-lakes (Reference Hastenrath and KutzbachHastenrath and Kutzbach, 1985; Reference KesslerKessler, 1985; Reference GrosjeanGrosjean, 1994; ΔP around 300 mm). Given broad synchronicity between the paleolakes and the maximum glacier advances (Reference Clayton and ClappertonClayton and Clapperton, 1997), we might speculate about redeposition of snow in wind-protected glaciated valleys which would add to the accumulation of direct precipitation and affect the mass balance. Such processes are not likely for paleolakes. Preliminary data for net mass accumulation measured from a new ice core of Cerro Tapado at 30° S in Chile, and modern climate conditions (personal communication from P. Ginot, 1999) suggest additions of redeposited snow in the order of 100–200 mm a⁻¹ w.e. on Glaciar Cerro Tapado, which is located in the wind-protected area of the mountain. Although nothing is known about this effect in the past, we consider such processes in the order of 100–200 mm a⁻¹ as a possible explanation for the high values of total precipitation drawn from glacier–climate models in Late-glacial northern Chile.

We also emphasize that glaciers responded likely on a time-scale of decades to centuries to climate changes. Thus our climate reconstruction may reflect a relatively short period of maximum humid conditions with maximum ice expansions (the absolute age determination is still missing), whereas the ground-water-fed lakes responded at much larger time-scales and may not have recorded this event. Therefore, differences are expected between the climate reconstructions from glaciers and lakes.

We might also argue that mountains trigger local precipitation due to their topography. Glaciers would benefit directly from such local effects, whereas the water balance of paleolakes integrates regional precipitation patterns of the lake basin where such effects are subordinate. The mean altitude of a glacier is also substantially higher than the mean altitude of a paleolake basin, which, given a regional precipitation lapse rate of +12 mm per 100 m altitude, may explain a difference in the order of 100 mm a⁻¹ between paleoclimate reconstructions from glaciers and paleolakes.

5.3. Late-glacial seasonality: winter or summer precipitation?

In section 5.1 and 5.2, the late Pleistocene climate reconstructions are based on the assumption that the seasonal distribution of precipitation remained similar to the modern situation, whereby approximately 80% of the precipitation falls during summer. This implies that the late Pleistocene increase in humidity is mainly explained by a reinforcement of the tropical precipitation belt. However, the source of the late Pleistocene humidity and hence the seasonality of precipitation remains a matter of debate (winter vs summer precipitation). Here, we use both case-study glaciers to test the hypothesis that reinforced extratropical winter precipitation may result in mass-balance terms and, ultimately, in a glacier geometry which matches the “real” geometry as observed in the field.

In our first scenario, winter precipitation increases during the 4 months June–September. Nowadays this is the typical season for episodic frontal precipitation penetrating into the case-study areas from the west wind belt (Pacific moisture source). Figure 7a (and Table 4) shows that the total amount of winter precipitation must not exceed 15% (180 mm a⁻¹) of the total annual precipitation, otherwise the geometry of the model glacier mismatches the field observation, again allowing for a measuring error of ±5%. Apparently, a winter precipitation scenario with high accumulation rates during winter and changed accumulation/ablation rates during the summer with mainly cloud-free conditions leads to different mass-balance gradients, and thus to a glacier with a geometry which differs from that observed in the field. Dry, cloud-free summer conditions favor energy-intensive sublimation compared to melt. This results in reduced mass-balance gradients and lower AARs. Such glaciers are currently found in the subtropical Andes at around 30° S in the extratropical winter precipitation regime. This finding has wide implications and suggests from the glaciological standpoint that a substantial increase in winter precipitation (June–September) must be discarded as a possible paleoclimate scenario for Late-glacial times.

Fig. 7.

ΔDMM_cr1/cr2 (DMM_cr1/cr2 ≠ 0) as a function of winter precipitation (% of total annual P). Based on the sensitivity test (±5% for geometrical parameters and flow parameters), we accept a range between −0.6 × 10⁶ kg a⁻¹ and +0.5 × 10⁶ kg a⁻¹ ( shaded area, dotted bar) as a possible numerical solution for a steady-state glacier. The scenario based on 4 months of winter precipitation (Fig. 7a) is indicative for frontal precipitation, whereas the scenario based on 6 months (Fig. 7b) includes frontal precipitation and “cut-off” events.

However, there is, today at least, a second mechanism which brings Pacific moisture from the west wind circulation to the research area: isolated drops of cold polar air which migrate to lower latitudes and trigger precipitation when colliding with continental warm humid tropical air masses (“cut-off” events; Reference Vuille and AmmannVuille and Ammann, 1997). These synoptic situations are typical for precipitation between the winter and the summer modes (i.e. in May–October) and involve the tropical and extratropical circulation belts. To account for these mechanisms, for the second scenario with increased winter precipitation we extended the period of precipitation to 6 months (May–October). Figure 7b (and Table 4) shows that the maximum proportion of winter precipitation increases markedly by 0–30% for Choquielimpie and by 25–50% for El Tatio. This result points to the importance of precipitation during interseasonal periods (May and October) and to the potential significance of “cut-off” events during Late-glacial times. It is also found that these mechanisms were probably more important in the southern case-study area, El Tatio, than in the more tropical northern test area, Choquielimpie.

However, in terms of changes in the seasonal precipitation pattern, one might argue that the lapse rates and gradients of the climate parameters (used as input data for the model) are different for summer and winter precipitation regimes, and thus the actualistic principle (one of the prerequisites of our model) may not be appropriate. In particular, the lapse rates for precipitation are expected to vary substantially between convective summer precipitation and advective winter precipitation. Modern meteorological data show that there is no significant difference in the lapse rates of winter and summer precipitation on the isolated peaks in the research area. Hence this line of argument does not seem to question our conclusions.

In summary, our results suggest that the postulated increase in Late-glacial precipitation was mainly due to reinforced summer precipitation. This compares favorably with earlier findings from long-distance transported pollen, and spatial patterns of paleolake levels and ELA depressions (Reference KesslerKessler, 1991; Reference Markgraf and WrightMarkgraf, 1993; Reference Grosjean, Geyh, Messerli and SchottererGrosjean and others, 1995, Reference Grosjean1996). A paleoclimate scenario based primarily on increased austral winter precipitation between June and September must be discarded as a possible explanation, but may, combined with precipitation of the “cut-off” type during the interseasonal months (May and October), account for up to 30% of the total annual precipitation in the northern test area, and 25–50% in the southern test area.