Carbon dioxide flux in the ablation area of Koxkar glacier, western Tien Shan, China

Abstract

Accelerating loss of glacial mass caused by rising global temperatures has significant implications. The global cycle of the greenhouse gas CO₂ is also associated with mineral weathering and glaciation. In glaciated areas, most estimates of atmospheric CO₂ consumption are confined to chemical ionic mass balance or analog modeling methods. We applied the gradient method to the surface of Koxkar glacier, western Tien Shan, China, for a 5 month period during 2012. The overall net glacier-system CO₂ exchange (NGE) rate was measured as -0.05 and -0.07mmol m^–2 s^–1 for regions of exposed ice and supraglacial moraine, respectively. This suggests that atmospheric CO₂ drawdown may occur during ice melting because of consumption of H⁺ by CO₂ hydrolysis that occurs in solutions. Using the degree-day model to calculate glacier ablation in bare-ice regions and considering characteristics of the NGE rate in the supraglacial debris region with the support of GIS, the daily NGE rate was estimated to be -1.23 ± 0.17μmmol m–² d–¹ between Julian days 125 and 268 of 2012. These findings present a new approach for modeling the dynamics of glacial CO2 sinks undergoing melting, and develop an understanding of the mechanism of atmospheric CO₂ exchange.

Introduction

Global mean surface air temperature has increased by ~0.74°C over the past 100 years, and is predicted to rise an additional 1.1-6.4°C by the end of the present century (Solomon and others, 2007). This warming has influenced alpine glacier fluctuation and runoff, resulting in annual mass losses of up to 3 ma–¹ (e.g. the European Alps; Zemp and others, 2005). Glacier mass loss is expected to continue, and simulations based on climate warming from greenhouse gases predict a 20% increase in mean annual discharge of the great Arctic rivers (Macdonald and others, 2005). In China, 82.2% of all monitored glaciers have retreated, with a total area loss of 4.5% from the 1950s to the late 1990s (Liu and others, 2006).

Because warming-induced alpine glacier melting is likely to intensify, water-rock interactions such as sulfide oxidation and carbonate dissolution are also expected to increase. Glaciation is responsible for a significant drawdown of atmospheric CO₂, generating a negative feedback on climate (Sharp and others, 1995; Tranter, 1996; Hodson and others, 2002; Yde and others, 2005; Bartoszewski, 2008). It has long been recognized that rates of terrestrial chemical erosion may influence the atmospheric concentration of CO₂, a radiatively important gas contributing to mean global temperature regulation (Raymo and Ruddiman, 1992). Because of its potential impact on climate at glacial-interglacial timescales, CO₂ sequestration during episodes of deglaciation has attracted considerable attention (Brown, 2002).

Research on carbon dioxide drawdown in glaciated catchments has been launched in recent years (Hodson and others, 2002; Yde and others, 2005; Krawczyk and others, 2008; Roberts and others, 2008). However, these studies have been almost entirely reliant on the ionic mass-balance method of carbon dioxide flux calculations, without direct CO₂ concentration and flux monitoring systems. Moreover, eddy covariance data of CO₂ flux have been collected over the horizontal surfaces of large ice bodies (Mölg and others, 2003; Cullen and others, 2007; Mac-Donell and others, 2010; Jarosch and others, 2011) and in complex terrain (Rotach and others, 2003; Hiller and others, 2008; Guo and others, 2011; Jocher and others, 2012). Such an approach may be unsuitable for the study of alpine glaciers, given complex terrain, frequent precipitation and fog, and other factors. The gradient method for collecting and calculating CO₂ flux has been used to study net CO₂ exchange with the atmosphere in a wide variety of terrains (snow, grassland, forest, bare land, and alpine) (Steffen and others, 2008; Bowling and Massman, 2011; Pumpanen and others, 2011), including complex terrains. Thus, the gradient method may be useful for the study of processes underlying CO₂ exchange between the atmosphere and a glacial surface, providing continuous data without interference from precipitation, fog and refrozen condensed water (Bowling and Massman, 2011; Pumpanen and others, 2011).

Using the mass-balance method, total solute fluxes and transient carbon dioxide sinks (which accounted for 14.2% of total solutes in bulk river water) were estimated at 791.2 kg (km ^₂ d)^–1 and 81.0 kg (km² d)^–1 in the Koxkar glaciated basin (summer 2004), respectively (Wang and others, 2010). Using the flux gradient method, this paper aims to estimate CO₂ fluxes over bare ice and supraglacial moraine at Koxkar glacier during the 2012 melt season. These results, together with other data, may shed light on the influence of atmospheric factors on CO₂ fluxes. Specifically, these data provide new insight into CO2 concentrations and fluxes in glaciated basins, which may _be used for modeling CO₂ change on glacial–interglacial timescales (e.g. Tranter and others, 2002). A further aim is to seek correlations between CO₂ flux and other factors across different underlying surfaces on Koxkar glacier, which may facilitate conclusions regarding the feedback relationship of glacial change to regional atmospheric CO₂ circulation.

Site Description

The glaciated Koxkar basin (41842’–41853’ N, 79859’– 80810’E; Fig. 1) is on the southwest side of Toumuer mountain, northwest China. The glacier is 25.1 km long and has an area of 83.56km². Supraglacial debris covers -15.6 km², representing 83% of the total ablation area (Han and others, 2010). Mean annual air temperature observed near the glacier terminus is 0.77°C, and mean summer temperature is 7.74°C (Han and others, 2008). The main source of precipitation is water vapor from the Atlantic and Arctic Oceans (Kang and others, 1985). Annual precipitation, 81% of which occurs from May through September, is ~630 mm at the end of the glacier. Precipitation in the study area is in a solid state (snow and hail). A field investigation from 2003 to 2012 indicated annual discharge at the glacial terminus of >1.0 x 10⁸m³, mainly during the warmer half of the year (94.5% of total runoff flux was from May through October).

Fig. 1. Map of the Koxkar glacier region (northwest China) sampling sites.

Two CO₂ flux gradient observation stations (Fig. 1) sit atop bare ice (site A: 41°47ʹN, 808060E; 3730ma.s.l.) and supraglacial moraine (site B: 418430N, 80°09ʹE; 3212ma.s.l.). Debris thickness at site A is <0.01 m and discontinuous, while the thickness at site B is ~0.70 m. The debris is mainly gray or dark-gray granite stone powder, broken stone particles, and rock mass, which covers 83% of the total melt area. Therefore, moraine is the principal physical source of crustal dissoluble ions. The nearly horizontal surfaces at sites A and B have an average surface slope less than 28, and a fetch of several hundred meters in the direction of the prevailing wind, which most often should be a valley-mountain breeze.

Methods

During the observation period, we established two synchronous gradient observational systems to measure two-tiered CO₂ concentrations, wind speed, wind direction, temperature and humidity (2.0 and 1.0 m above the surface), and a series of radiation sensors that included total solar radiation, net radiation and reflection radiation in bare-ice and supraglacial moraine regions. Supraglacial two-tiered CO₂ concentrations (2.0 and 1.0 m above the surface) were determined using infrared CO₂ sensors (ES-D type, Sense Air, Sweden), with a variance of ±5.4% for CO₂ (Rinne and others, 2000; Lu and others, 2010; Fang and others, 2011). Data (e.g. CO₂ concentration, wind speed, wind direction, temperature, humidity, total solar radiation, net radiation, reflection radiation, and precipitation) were obtained using a PC-3 data acquisition system, and were automatically updated hourly throughout the observation period. Generally, the effective distance was less than ten times the height of the equipment probes. Our gradient observational stations were located in the central part of the glacier, away from nearby grasslands (-1.43 km away), to reduce the possibility of the grassland affecting the CO₂ measurement. Unfortunately, we could not obtain synchronous data because of complications with the solar power supply.

The net glacier-system CO₂ exchange (NGE) rate between the glacier surface and atmosphere for the temporally and spatially integrated net glacier system is defined as

(1)

where F_c is its vertical flux (upward flux being positive) and S is the change of CO₂ storage below the CO₂ infrared sensors.

We applied the gradient method, in which vertical CO₂ fluxes at the glacier surface are estimated by the aerodynamic gradient method (e.g. Sutton and others, 1993). The flux F_c of a component with concentration c is calculated from the product of friction velocity u_* , and a concentration scaling parameter c _* as (Rinne and others, 2000; Loubet and others, 2013)

(2)

where u, and c. are derived from the stability-corrected gradients of wind speed (u) and CO₂ concentration (c) vs height (z):

(3)

where d is the zero plane displacement height, which represents the shift in aerodynamic ‘ground’ because of the presence of the glacier surface, z is the mounting height of the sensors (2 m), k is the von Kármán constant with values between 0.35 and 0.43 (usually 0.4), and _m is the integrated stability correction function for momentum.

Harman and others (2008) developed an analytical transfer model for the roughness layer that includes two displacement heights; the model was extended by Siqueira and Katul (2010) to account for stomatal regulation and soil respiration. Harman and others (2008) showed that these additional functions would lead to an additional integral term _H in the concentration profile, which in turn would lead to a new approximation of the concentration scaling parameter c_* :

(4)

Siqueira and Katul (2010) give the following expression for 𝜓_H:

(5)

where S_n is the Schmidt number, (β = c_*/u_h, d_c and d_m are the scalar and momentum displacement heights, and c _2c is a constant that can be taken as 1 in the first approximation. More details of 𝜓 _m, d_c and d_m are provided by Benjamin and others (2013).

The value for S was calculated from vertical CO₂ concentration profiles (Aubinet and others, 2001; Araújo and others, 2010

(6)

where P _a is atmospheric pressure, R is the molar gas constant, T _a is the air temperature (K), h _m is the maximum measurement height above ground level (m), c is the CO₂ concentration, t is time (s), and z is the height about ground level (m). Details of the S calculation are provided by Valentini and others (2000) and Kowalski (2008).

Results

CO₂ fluctuations

Typical fluctuations in hourly CO₂ concentration, atmospheric temperature and wind speed measured 2m above bare ice and supraglacial moraine are shown in Figure 2. Mean wind velocity during the observation periods was 2.3 ms–₁ at both stations. Mean temperatures were 3.4°C and 7.6°C, and the amplitudes of dominant fluctuation in CO₂ concentration were ~439.3 and 472.1 mgm–³ for sites A and B, respectively. Variations in CO₂ concentration gradually decreased, followed by increasing atmospheric temperature, which otherwise only gradually increased. For example, there is strong negative correlation between CO₂ concentration and atmospheric temperature at site A. The associated equation is

Fig. 2. Fluctuations of hourly CO₂ concentration, atmospheric temperature and wind speed above the surface of Koxkar glacier. Data are from bare ice (site A, right of the vertical dotted line) and above supraglacial moraine (site B, left of the vertical dotted line).

(7)

where C _CO2 is the CO₂ concentration (mgm–³), and T_atm is atmospheric temperature (°C) measured 2 m above the ice surface.

Based on Eqns (7) and (8) for site B, the CO₂ concentration was less influenced by atmospheric temperature. Equation (8) shows the relationship between CO₂ concentration and atmospheric temperature at site B:

(8)

Unfortunately, data from sites A and B were not collected simultaneously during the observation period. However, we conclude that the CO2 concentration at site A was dramatically lower than at site B (Fig. 2). This can be explained by the different altitudes and underlying surface conditions at the two sites.

Net glacier-system CO₂ exchange (NGE) rate

During the observation period, hourly NGE ranged from –0.74 to +0.48 mmol m^–2 s^–1 at site A and from –1.81 to +0.66 mmol m^–2 s^–1 at site B (Table 1). Amplitudes of dominant hourly NGE fluctuation were approximately –0.05 and –0.07 mmol m^–2 s^–1 for the bare-ice and supraglacial moraine observation points, respectively. These variance results suggest that the fluctuation range of NGE at the supraglacial moraine (site A variance = 60.70) was greater than at bare ice (site B variance = 33.07). Moreover, the maximum and minimum values of hourly NGE were both recorded at site B. Under normal conditions, upward flux (i.e. CO₂ release) occurs with precipitation or at night, when ice melt decreases because of reduced air temperature.

Table 1. Hourly NGE (mmol m^–2 s^–1) at sites A and B on the surface of Koxkar glacier in 2012

The negative mean values of hourly NGE at the two sites imply that the glacial surface was a sink for atmospheric CO₂ during the melt season because of the hydrochemical conditions and mineralogical reactions existing under the ice-melt water (Hodson and others, 2002 ; Tranter and others, 2002; Yde and others, 2005; Krawczyk and others, 2007, 2008; Lerman and others, 2007). In glaciated regions, these chemical reactions are primarily carbonate and silicate dissolution, for example:

(9)

(10)

Although there was a little cryoconite at the bare-ice site, there were limited amounts of carbonate, silicate and other soluble materials, which could inhibit chemical reactions such as those described by Eqns (9) and (10). Substantial glacial debris at the supraglacial moraine site, including abundant carbonate, silicate and potash/soda feldspar, generated ample chemical hydrolysis under precipitation or ice/snow meltwater conditions. Such mineral differences between sites can account for the differences in hourly NGE rates.

Discussion

Errors in the NGE calculation

In the past, errors in the aerodynamic gradient method used to calculate the atmospheric H₂O/CO₂ flux have been extensively studied for sites in such locations as grassland, forest, tundra and farmland (Wenzel and others, 1997; De Ridder, 2010). In the Koxkar glaciated area, the main factors that influenced NGE accuracy were as follows:

1. 1. The variance, ±5.4%, of CO₂ measured by infrared sensors was an objective value, which might cause 4.8-7.2% error in the NGE result, with a mean value of 5.3%. The sequential sampling of gas concentrations, sourced from ice melting, at various heights above the ice showed a variability in the gradients due to the non-stationarity of the concentrations. These were obviously not an issue for u _*, but might be an issue for CO₂ flux (Benjamin and others, 2013).

2. 2. A roughness sub-layer correction has been proposed previously (Wenzel and others, 1997; Siqueira and others, 2010; Benjamin and others, 2013). In this work, the roughness sub-layer correction was evaluated based mainly on the approaches of Garratt (1978), Cellier and Brunet (1992) and De Ridder (2010). In general, the error most attributed to the geometric correction for the measurement height of the monitoring equipment was <15%.

3. 3. The concentration and flux footprint errors are mainly the result of local advection error (Marcolla and others, 2005; Feigenwinter and others, 2010; Novick and others, 2014). Usually, the values of horizontal advection and vertical advection were opposite, and their sums were approximately zero (Feigenwinter and others, 2010; Novick and others, 2014). Calculating the horizontal advection should include at least two flux towers. Hence, the flux footprint error is ignored in this paper.

4. 4. In recent years, englacial gas release has been widely mentioned (Jaworowski and others, 1992; Smellie, 2006; Ryu and Jacobson, 2012). Many studies have suggested that characteristic gases (including CO₂) from various historical periods may be released from air bubbles within ice (Jaworowski and others, 1992; Anderson and others, 2000; Ryu and Jacobson, 2012). This raises the possibility that future atmospheric CO₂ levels will increase as glaciers retreat or disappear. Although this result appears to conflict with the CO₂ consumption of negative total NGE as a result of the chemical reactions described here, the total volume will actually decrease by -10% during ice melt. If we suppose that the CO₂ concentration of glacial air bubbles maintains the higher value of AD 1100 (295 ppmv; Barnola and others, 1995) and that gas pressure inside the bubbles was that of the standard atmosphere, CO₂ release was 2.30 mmol m–² at the bare-ice observation site. This is only ~1.16% of the NGE rate between Julian days 223 and 268 in 2012. As a result, CO₂ release from glacial ice was ignored, whereas the NGE in the glacier region was analyzed. However, released CO₂ dissolves in water and participates in chemical reactions between ice-melt water and rocks (Brown, 2002; Jacob and others, 2005).

CO₂ storage term

The CO₂ storage term (S) is an important component in the statistical inference of net ecosystem CO₂ exchange (NEE) rates (Finnigan, 2006, 2009; Kowalski, 2008; Araujo and others, 2010). There has been little study of S in glacial regions. During the observation period, daily S values ranged from -1.74 to +1.43 and -1.19 to +1.55 mmol m–² d–¹ at sites A and B, respectively (Fig. 3), with mean values of +0.06 and -0.03 mmol m^–2d–¹ (about 1.36% and 0.51% of total NGE rates). These data demonstrate that CO₂ storage (S) rates had little effect on NGE rates in this glaciated area. Given the rapidly increasing glacial surface temperature during daylight hours, S may be positive, resulting in atmospheric CO₂ releases below the CO₂ detection limit of the infrared sensors. Glacial meltwater was shown to be a fundamental carrier of chemical reactions, producing transient CO₂ drawdown and decreasing air temperature caused by precipitation. The maximum S (1.55 mmol m^–2 d^–1) was caused by sustained heavy rainfall between Julian days 194 and 196 (total precipitation 44.6 mm), leading to a decrease in air temperature and a reduction of ice melt, thereby inhibiting chemical reactions.

Fig. 3. Study area variations of daily F_c and S at sites A and B (right and left sides of vertical dotted line, respectively).

NGE rates affected by ice-melt intensity (IMI)

The most obvious characteristic of the NGE rate was that glacial meltwater was the fundamental carrier of carbonate and silicate hydrolysis but was not a determining factor in the daily NGE rates at site B (Fig. 4). Moreover, there were no significant changes in IMI, with daily NGE rates fluctuating around -6.26 mmol m–² d–¹ at that site (Table 1). For example, NGE rates at the site were -179.17 and -30.40 mmol m–² between Julian days 125 (4 May) and 153 (1 June), and 215 (2 August) and 220 (7 August), respectively. The IMI was consistent across the same periods (Table 2). This suggests that the daily NGE rates are relatively stable with abundant soluble substances, and therefore not remarkably affected by IMI. At site A, however, there was correlation between IMI and NGE rates (Eqn (11)), while daily precipitation was <2 mm, demonstrating that the daily NGE rates, including those caused by the dissolution of carbonate and silicate, depended on ice melting. A mass of ice-melt water accelerated the strong chemical reactions, leading to a decrease in NGE rate:

(11)

where IMI is in mm and n is the number of days.

Fig. 4. Relationships of NGE rates, IMI and precipitation at (a) site B and (b) site A.

Table 2. Relationship between IMI and NGE rates at site B in 2012

Based on these results, it would not be possible to use the ice-melting and NGE rates to calculate the yearly flux of atmospheric CO₂ consumption in the Koxkar glaciated region. To estimate the NGE rate for the whole glacial area, we used the degree-day model (Singh and others, 2000; Sicart and others, 2008; Mourshed, 2012; Tennant and others, 2012; Mayr and others, 2013) to calculate the IMI in the bare-ice region, and considered the characteristics of the NGE rate in the supraglacial debris region with the support of GIS, resulting in a daily NGE rate of about -1.23 ± 0.17 mmol m–² d–¹ between Julian days 125 (4 May) and 268 (24 September) in 2012. This estimate is less than the 1.84 mmol m–² d–¹ that was calculated in terms of ionic mass balance by Wang and others (2010) in 2004, likely by ignoring the chemical reactions that follow glacier surface runoff in englacial conduit geometry. The NGE rate in the Koxkar glaciated area is greater than the CO₂ drawdown of the Rhone and Oberaar glacial catchment in the Swiss Alps (Hosein and others, 2004), and also higher than that of the Scottbreen glaciated basin in the Arctic (Krawczyk and others, 2008). Substantial glacial debris might be the main reason for the higher chemical erosion and NGE rates in the Koxkar glacier area (Giles, 2002; Hosein and others, 2004; Krawczyk and others, 2008; Zhang and others, 2013).

Conclusions

At Koxkar glacier, atmospheric CO₂ concentrations have gradually decreased as atmospheric temperatures have increased. The average values of hourly NGE were approximately –0.05 and –0.07 mmol m^–2 s^–1 for bare ice and supraglacial moraine, respectively. This implies that the basin is a sink for atmospheric CO₂ above the glacial surface during melting seasons because of existing hydrochemical conditions and mineralogical reactions under ice-melting water. The CO₂ storage (S) rates had little effect on NGE rates in the glaciated area, accounting for only 1.36% and 0.51% of the total NGE rate at sites A and B, respectively. A NGE rate was analyzed in the glacier region. Using the degree-day model to calculate glacial ablation in a bare-ice region, and considering the characteristics of the NGE rate in a region of supraglacial debris with the support of GIS, the daily NGE rate was calculated as approximately -1.23±0.17mmolm^–2d-¹ between Julian days 125 and 268 in 2012. This could be lower than the actual value, because chemical reactions that follow the process of glacier surface runoff in englacial conduit geometry are not considered.