At present, the glaciers that have been studied for mass balance are not only few in number but also rather small. Reference VilesovVilesov (1969) and Reference Kotlyakov, Grosswald, Dyurgerov and MazoKotlyakov and others (1991) estimated the changes to glaciers and their runoffs in central Asia using mass-balance data from several small glaciers in the Tien Shan mountains, China. However, the derived changes may not represent the actual changes to larger glaciers. Other investigators (e.g. Reference OerlemansOerlemans,1997; Reference Schmeits and OerlemansSchmeits and Oerlemans, 1997; Reference Zuo and OerlemansZuo and Oerlemans, 1997; Reference Wallinga and van de WalWallinga and Van de Wal, 1998) have used ice-flow models to study the responses and sensitivities to climatic change of a series of specific glaciers. It is not practical, however, to simulate each glacier individually where there are a large number of glaciers in the same region. Moreover, there is significant variation in the way glaciers respond to climate change, not only over time but also in quantitatively, due to their different sizes and shapes (Reference ZongtaiWang,1991; Reference Yongjian and HaeberliDing and Haeberli,1996; Reference Chaohai, Jianming, Mingxie and GuodongLiu Chaohai and others, 1996). The objective of this paper is to simulate the response to climatic warming of various-sized glaciers in the Yili river basin using a glacier ice-flow model, and thus to reveal the impact of climatic change on glaciers of different sizes.
2. The Yili River Basin
The Yili river originates in the western Tien Shan, China (42°14′–44°51′ N, 80°15′–86°56′ E), then flows through the Yili region of Xinjiang Uigur Autonomous Region and finally discharges into Barkshen lake in Kazakstan (Fig. 1). The Yili river basin consists of a series of valleys and subbasins that are separated by several mountain ranges, oriented approximately east–west. According to data from the glacier inventory of China and the former Soviet Union (Reference VilesovVilesov, 1969; Reference Liangfu, Weirong, Chaohai, Chunzu and ZongtaiDing and others, 1986), 2414 glaciers with a total area of 2113 km2 exist in the upper basin. These are at elevations of 2640–5568 m, with the mean snowline at about 3730 m a.s.l., and the mean annual lower boundary of permafrost at approximately 3000 m a.s.l.
3. Glacier Ice-Flow Model
The introduced glacier ice-flow model comes from Oerlemans’ main glacial flowline model (Reference OerlemansOerlemans, 1988) and considers the features and shapes of glaciers as described by Reference Schmeits and OerlemansSchmeits and Oerlemans (1997). The model consists of glacier-continuity and glacier-flow equations with the main flowline being the x coordinate.
where S is the area of the cross-section, B the specific surface mass balance, W the width of the ice surface, U the depth-averaged ice velocity, U d the velocity caused by internal deformation, U S the velocity caused by sliding, τ d the shear stress at the bottom, H the glacier thickness, α the slope of the ice surface along the flowline, γ the slope angle on the sides of the valley (γ is assumed to be 45°), f 1 the glacier-flow parameter, f 2 the ice-sliding shape parameter, and g and ρ the acceleration of gravity and density of ice, respectively (g = 9.8 m s−2, ρ = 900 kg m−3).
As the glaciers in theTien Shan are categorized as subpolar and polar types, their movement, except for some very large valley glaciers, is mainly attributed to glacier ice deformation (Reference MaohuanHuang, 1995). Therefore, the ice deformation should be taken into account in the glacier flow equation (Reference Ye and KegongYe and Chen,1997), and the glacier sliding item U S in Equation (3) can be omitted. Thus, Equation (3) is modified to:
In this model, the effect of climatic change on glaciers is expressed as changes of glacier mass balance. Although climatic change may also alter glacier temperature, and thus lead to glacier change, this process is ignored in this study because of its complexity and the relatively small magnitudes of the glacier changes that result.
Substituting the expressions for U and S in Equation (1) leads to a non-linear diffusion equation for H. A forward time-differencing scheme is used, which is stable if the time-step is sufficiently small (Reference OerlemansOerlemans, 1997). Tests have been done on grids with spacing of 500, 300, 100 and 50 m with a time-step of 30 days in this study. The scheme is stable for small glaciers (map area S g < 0.5 km2) with 100 m grid spacing, and for large glaciers (S g > 40 km2) with 300 m spacing.
3.1. Characteristic values of glacier shapes
In order to parameterize the shapes of glaciers of various sizes, some characteristic values, such as glacier area, length, altitude difference (between the maximum and minimum elevation of a glacier), maximum elevation and average slope, were statistically analyzed for more than 2000 glaciers in the Yili river basin, based on data from the glacier inventory of China (Reference Liangfu, Weirong, Chaohai, Chunzu and ZongtaiDing and others, 1986). The results show good correlations between glacier area and altitude difference, as well as between glacier length and altitude difference. The glacier shape can thus be determined from the two correlations. The parameters representing the absolute glacier elevation can be expressed using the relationships of glacier area to maximum, mean and minimum elevation, where the relationship of the area to the maximum elevation is best. Other methods can be used to calculate initial glacier surface elevation, thickness and width along the main flowline, as described below.
3.1.1. Glacier area and maximum elevation
In a basin, the maximum glacier elevation (H max) has a significant impact on the simulation result. The relationship between H max and glacier area is found using the glacier inventory; it is not good, especially for the large glaciers since most of the glaciers used in the relation are small. A statistical average of H max for the different glacier areas of the Yili river basin is derived by means of glacier area classification. The relationship between H max and area (S g) is shown in Figure 2 and expressed by:
where the unit of S g is km2 and the unit of H max is meters.
3.1.2. Glacier length and glacier altitude difference
Using the glacier inventory data, the relationships among glacier area (S g), glacier length (L) and altitude difference (TZ) were derived by the regression method as:
where the unit of S g is km2 and the unit of L and TZ is meters. n is the sample size and R is the correlation coefficient.
3.1.3. Change in ice surface slope along the main flowline
Generally, the surface slope of a glacier strongly influences its flow velocity and thickness. The surface slope varies greatly, especially in places with icefalls, due to the underlying terrain, and the surface slope along the main flowline varies with glacier area and type. To estimate changes of surface slope along the main flowlines, 47 glaciers in the Yili river basin were selected for analysis from the 1: 50 000 topographic map, and 29 glaciers in the Tomor peak region from the 1: 200 000 glacier topographic map. Statistical analysis of the data yields the following result for the surface slope K along the main flowline:
where l is the relative distance (from the top of the glacier) along the main flowline [0,1].
3.1.4. Changes in glacier area and surface width with elevation
Based on the data from two topographic maps published in 1943 and 1977, Reference KuzmichenokKuzmichenok (1993) measured and statistically analyzed glacier areas for different elevation belts of 178 glaciers in the Akcerak peak region of the western Tien Shan. He found that the distribution of glacier area with elevation was determined by:
where Zj is the relative elevation, Zj [0, 1], and ai different factors:
where, when S g ≤ 0.1 km2, S g0 = 0.1, and when S g ≥ 70 km2, S g0 = 70, when 0.1 km2 ≤ S g ≤ 70 km2, S g0 = S g (km2).
The glacier width, W, at a given elevation belt can be calculated according to the glacier area in the belt (ΔSj ) and the corresponding distance which is determined by the surface slope and the elevation difference. The surface slope is determined by Equation (9). W canbe expressed as:
where Wj , ΔSj and Zj are in meters.
3.1.5. Estimation of initial glacier thickness
According to Reference PatersonPaterson (1981), initial glacier thickness H 0 along the main flowline can be estimated by:
where τ 0 is ice yield stress, varying from 80 to 150 kPa, and K h is a factor accounting for glacier cross-section. We determined glacier thickness in the main flowline from the slope in the main flowline. The goal of estimating glacier thickness is to determine glacier bed altitude from glacier surface altitude.
Based on the other studies (Reference ZhenSu, 1985; Reference Genxiang, Yafeng, Ersi and JianmingYou and others, 1988) that were conducted using topographic data and the thickness of Tailan glacier in the Tomor peak region and of Ürümqi glacier No. 1 (UG1) in the Tien Shan, τ 0/K h ρg changes from 8 to 17 m, with an average of 11.3 m for τ 0 = 100 kPa.
Based on the above analyses, the shape characteristics of glaciers of various sizes can be statistically determined by Equations (6–12). From these we can create the “representative glacier” which will have the statistical characteristics of glaciers of various sizes. The representative glacier can be regarded as the initial input parameters of the glacier ice-flow model to simulate the responses of glaciers to climatic change.
3.2. Calculation of glacier mass balance and runoff
The glacier mass balance can be calculated by the degree-day factor (Braithwaite and others, 1993; Reference Shiyin, Yongjian, Ye and GuodongLiu Shiyin and others, 1996), and the runoff can be estimated by the water balance in various elevation bands (Reference KonovalovKonovalov,1979). The degree-day factors of snow and ice ablation are 3.1 and 7.3 mm K−1 d−1, respectively, as estimated by the mass-balance data of UG1 (Reference Shiyin, Yongjian, Ye and GuodongLiu Shiyin and others, 1996). The elevation bands of UG1 cover a wide range (3740–4400 m a.s.l.), so the different precipitation types have been taken into account in calculating the glacier mass balance (Reference Ye and KegongYe and Chen,1997). Because most of the meteorological and hydrological stations are located below 2000 m a.s.l. in the basin, and data from high-elevation areas are limited, temperature and precipitation in high-elevation areas were estimated by the co-relationship between elevation and precipitation from data collected at mid- and low-elevation areas (Reference Ye, Zuming and YafengYe and others, 1997).
3.3. Determination of comprehensive glacier flow parameter (f 1)
The glacier flow parameter, f 1, varies from 10−25 to 10−24 m6 N−3 s−1, as suggested by other studies (Reference OerlemansOerlemans, 1997; Schmeits and others,1997; Reference Zuo and OerlemansZuo and Oerlemans,1997; Reference Wallinga and van de WalWallinga and Van de Wal,1998), with the range. In this study, f 1 is calculated using Equation (6) and the observed data from UG1, close to theYili river basin. The observed data include the surface movement velocity measured during the period 1980–86, the glacier thickness determined by radar and validated by borehole data, and the topographical data from UG1 (Reference Genxiang, Yafeng, Ersi and JianmingYou and others,1988). Values of f 1 are listed inTable 1.
For mass-balance calculation, the time interval was set to 30 days, and the distance interval was limited to 100–300 m according to glacier sizes.
4. Glacier Response to Climatic Change
Glaciers with areas of 0.5–45.0 km2 in the upper reaches of theYili river basin were selected to simulate glacier response to climatic change. Based on statistical analyses of glacier shape parameters, 15 different glacier schemes were calculated using Equations (6–12). The 15 schemes were regarded as initial conditions of the glacier ice-flow model, and the input data are mean monthly air-temperature and precipitation data from 1954 to 1997 in the Yili river basin.
4.1. Simulations of glacier sizes at steady state
Glaciers of different sizes were selected for simulation using the glacier ice-flow model. Glaciers with different areas, initial conditions and average temperature and precipitation during 1954–97 are the input of the model.The mass balance is calculated by the degree-day factor according to the monthly temperature and precipitation at different elevations, and then the calculated mass balance is used as input data for the glacier ice-flow model. The model indicates that the glaciers will reach steady state after a long period (2000 years) under the current climatic conditions. It can be seen from Figure 3 that the simulated length and area of glaciers are consistent with their initial values. In fact, the glaciers in the Tien Shan have been retreating since the 1950s. The model has somewhat overestimated the glacier size in modern climate conditions. The mass-balance and ice-flow models have been examined in UG1 (Reference Ye and KegongYe and Chen, 1997). The small overestimate may be caused by (1) overestimated precipitation in the alpine area, based on the precipitation gradient in the middle and low mountains, and (2) the H max which is a significant factor to steady glacier size. The simulated results consistent with the initial glacier areas and lengths show that the model is reliable for simulating glacier sizes under modern climatic conditions, at least for the mass-balance model and the representative glacier with various areas.
4.2. Schemes of climatic change
The schemes of climatic change were selected based on mean air temperature and precipitation observations from 1954 to 1997. Because future precipitation change is difficult to predict, the function of temperature change is only considered for the climate-change scenarios in the study. It was assumed that air temperature rises by 1 K, and the rates vary from 1 K per year (sudden change) to 1 K per 1000 years.
4.3. Simulated results
4.3.1. Processes of glacier change
The temporal changes in glacier area under different climatic scenarios are shown in Figure 4. It can be seen that there are significant differences in the sensitivity of various-sized glaciers to climatic warming, and that the larger the glacier, the faster the retreat (absolute change due to warming). However, the relative changes to the larger glaciers are rather small. This result is consistent with the comparison between two recent aerial photographs in the Ürümqi river basin (Reference Jianming, Chaohai and MingxieChen and others, 1996). In addition, the time required for glaciers to reach a steady state is prolonged with size increase under a given climatic disturbance.
4.3.2. Changes in glacier runoff
Change in glacier runoff under different climatic scenarios is illustrated in Figure 5 (the runoff is relative to the average runoff under the present climate for 1954–97). In each case, there are maxima in glacier runoff which are defined as runoff peaks. It is obvious that the smaller the glacier, the more sensitive the response of glacier runoff to climatic change, the higher the runoff peak, and the faster the attenuation of runoff. The timing of the runoff peak depends on the rate of air-temperature rise (RTR (K a−1)). The higher the rate, the earlier the appearance of the peak, and the higher the peak magnitude. The discrepancy between the changes in runoff and temperature shows hysteresis of the response of glaciers to climatic change.
4.3.3. Values of glacier runoff peaks
In the arid region of northwest China, the availability of water is of vital importance to society and the economy, and the glacier runoff forms an important part of the water resources. Accordingly, we focus on the glacier runoff change, especially on the volume and the timing of the glacier runoff peak.The dependence of runoff peaks on glacier size and different rates of air-temperature rise is shown in Figure 6. The runoff peak decreases progressively with decrease in the rate of air-temperature rise and with increase in glacier size (Fig. 6a and b). If the glacier area is large enough ( > 30 km2) or the rate of temperature rise is slow enough ( < 0.005 K a−1), the runoff peak tends to be a constant, and thus the difference between different sizes of glaciers is small, becoming significant only when climate change is rapid. Because of limitations of the model, the timing of the runoff peak was not reliably simulated when the rising rate was small and the peak was low (Fig.7). However, the results still show a distinct tendency: if the rate of temperature rise is large, the timing of the runoff peak will stay synchronous with the occurrence of maximum air temperature; with decreasing rising rate, however, the runoff peak appears before the temperature peak reaches its maximum for small glaciers.
4.3.4. Changes of glaciers of various sizes at equilibrium state
Figure 8 illustrates the change in area, length and size of initially stable glaciers for a warming of 1 K. It shows that the relative sensitivity to warming is larger for smaller glaciers. The relative changes in length in response to warming are smaller than the volume changes, and changes in glacier area due to warming are smaller than changes in length under the same conditions. The relative changes are similar if the glacier size is relatively small (e.g. < 0.6 km2). The sensitivity of glacier area to warming exceeds that of length when the glacier area is small. This shows that changes to small glaciers are dominated by area reduction rather than length shortening under warming. This result is validated by the 1964–92 observation data from the Ürümqi river basin: changes in glacier length, area and volume of 12.4%, 13.8% and 16.8%, respectively, were measured for glaciers with an average area of 0.31 km2 (Reference Chaohai, Jianming, Mingxie and GuodongLiu Chaohai and others, 1996). For larger glaciers, length decreases dominate.
Figure 9 compares the simulated change in glacier area and the change measured by distributions of glacier moraine since the Little Ice Age (LIA) in the Tien Shan (Reference ZongtaiWang,1991). Figure 10 shows the changes in glacier volume, length and area for different glacier sizes measured from the selected 66 glaciers in the upper Yili river basin using aerial photographs taken during the last 40 years. Although the climatic conditions applied to the two simulations and the degrees of glacier changes are different, the trends of glacier changes indicated by Figures 9 and 10 are consistent, suggesting that the glacier ice-flow model can simulate glacier changes with different glacier sizes in response to climatic warming.
4.3.5. Model validation
Validation of the glacier model, including the mass-balance model, has been performed on UG1 (Reference Shiyin, Yongjian, Ye and GuodongLiu Shiyin and others, 1996;Reference Ye and KegongYe and Chen,1997). However, because observation data are lacking for processes of glacier change affecting different areas, the simulation was validated by comparing initial glacier areas and lengths with the simulation of the glacier sizes at steady state, and with the actual change of glacier area in the Tien Shan since the LIA. To test the simulated results, we compared them with the actual glacier changes in the Tien Shan since the LIA and during the period 1962–92. The air temperature has risen by 0.5–0.8 K on the global surface over the last 150 years (Reference Jones, New, Parker, Martin and RigorJones and others,1999) and by 1.4 K in west China since about AD 1700 (Reference Yafeng and ShiyinShi and Liu, 2000). According to meteorological data from 97 stations above 200 m altitude in west China, the summer temperature rise is about 0.008 K a−1 and the annual temperature rise 0.013 K a−1 during the period 1955–95. Considering the difference between annual and summer temperature, and excluding the rise over the last 40 years, the change in summer temperature in northwest China is about 0.5 K between the Little Ice Age and the 1950s. The area changes for a selected 15 km2 glacier are listed inTable 2. The simulated change in glacier area with a warming of 0.5 K within 250 years is very close to the change in glacier area since the LIA. However, the glacier change during the period 1962–89 is not consistent with the simulated result. The reason is related to differences in initial conditions: glacier retreat in the Tien Shan has been ongoing since the 1950s, but the simulated initial condition is a stable one, so the simulated result is less than the actual change.
The following conclusions are drawn:
(1) The glacier ice-flow model can simulate glacier changes as functions of glacier size in response to climatic warming.
(2) The larger the glacier, the faster the retreat (absolute change due to warming). However, the relative changes of larger glaciers are rather small.
(3) Compared to glacier length and area, glacier volume undergoes the most pronounced change with warming. For larger glaciers, the change in length (retreat) is more significant than the change in area; for smaller glaciers ( < 0.6 km2), the opposite is the case.
(4) Under climatic warming, smaller glaciers with higher runoff peaks tend not only to vary more in runoff, but also to retreat more quickly than larger glaciers.
(5) The magnitude and appearance time of runoff peaks depend not only on glacier size but also on the rising rate of air temperature. The faster the air temperature rises, the higher the peaks, and the earlier the appearance time of peaks.
This research project was supported by grants from the Innovation project of the Chinese Academy of Sciences (KZCX1-10-06, KZCX2-301) and the Chinese National Scientific Fund (49971022). We thank T. S. Clarke, R. J. Braithwaite and W. D. Harrison for their valuable comments and suggestions, and J.W. Glen, the Scientific Editor. We also thank Zhang Xinping and Yang Daqing for suggestions.