Finite-volume model

The finite-volume model is based on the differential equation for heat conduction and advection:

where T is temperature, b is the advection rate (in this application it is the snow accumulation rate), k is the thermal conductivity, ρ is the snow density and C _p is the heat capacity. k/ρC _p is the snow thermal diffusivity. S is the source term, to be described below.

This formulation comes from an advective–diffusive equation for enthalpy (Reference PatankarPatankar, 1980). The equation is applicable to any quantity that is subject to diffusion under a gradient or advection under fluid flow. Each term in Equation (1) has units K s⁻¹. The boundaries of the 1-D solution domain are fixed at the snow surface and at a specified depth.

In our case, the uppermost boundary is the transient surface (i.e. it is fixed with respect to the snow surface regardless of the accumulation rate b). When significant accumulation is incorporated into the model, snow of a given T is pushed (advected) downward relative to the upper boundary of the model. On the righthand side of Equation (1), the first term represents the diffusion of heat down a temperature gradient. The second term on the righthand side, the source term, accounts for any process that may change the temperature of a control volume but is not related to the processes of thermal diffusion or advection. One potential source-term process is absorption of solar radiation by near-surface snow. Heating of a control volume by absorption of solar radiation depends on the amount of the incident solar radiation, and on the microstructure and composition of the snow. In the model, it does not depend on the thermal diffusivity of the snow nor on the accumulation rate, but in reality both solar absorption and thermal diffusivity depend on snow grain size.

The model is illustrated in Figure 1. The solution for temperature at point P is fully implicit, meaning it depends upon the most immediate previous values of temperature at point P as well as on the present unknown values at points U and D. Other numerical solutions to Equation (1) are possible, but are often not optimal either in interpretation (e.g. Taylor-series expansions) or in general applicability (e.g. variational formulation) (Reference PatankarPatankar, 1980).

Fig. 1. The box represents the 1-D finite-volume model used to calculate snow temperatures, heating rates and vapour pressures. The thick arrows indicate vertical advection (snow accumulation). The black dots are the points for which implicit solutions are determined at each time-step. The hatched area is one control volume. The sun represents a physical example of a source term for Equation (1), where the solar heating rate is measured in J m⁻³ s⁻¹. The nodes at which we solve for temperature are represented by the solid horizontal lines. The upper boundary condition is skin-surface temperature measured at 9 min intervals, and the lower boundary condition (at 6.5 m depth) is a measured seasonally varying temperature gradient.

Application of the finite-volume model to near-surface snow at the South Pole

Our goal is to understand the heat and vapour content of the near-surface snow at the South Pole on timescales from minutes to years. The high degree of horizontal homogeneity of East Antarctica and the slow horizontal movement of the ice at the South Pole (approximately 10 m a⁻¹; Reference Bingham, Siegert, Young and BlankenshipBingham and others, 2007) eliminates the need for horizontal conduction and advection of temperature in the model, allowing us to use a simplified 1-D finite-volume model. There is horizontal inhomogeneity in the snow surface in the form of sastrugi, which are wind-carved snow structures of height 5–20 cm at the South Pole. They form primarily in the winter at the South Pole, but become flattened by differential solar heating of their sides during summer (Reference GowGow, 1965). Such small-scale horizontal temperature gradients are not represented in this model. Sensitivity tests with the model indicated that the accumulation rate at the South Pole (∼80 mm w.e. a⁻¹; Reference Mosley-Thompson, Paskievitch, Gow and ThompsonMosley-Thompson and others, 1999) does not advect temperatures downward into the snow significantly relative to the rate at which temperatures were conducted vertically during the 9 year time period, so advection is not included below. Similar results were found in a modelling study of stable isotopes in firn at Taylor Mouth, Antarctica, above the McMurdo Dry Valleys (Reference Neumann, Waddington, Steig and GrootesNeumann and others, 2005).

We present results primarily from the top metre of the snow because deeper snow retains little memory of synoptic variability. The model depth was set to 6.5 m. This is a compromise between computation time, available seasonal data to constrain the lower boundary, a lower boundary that is far enough away from the near-surface snow so as not to substantially affect the results in the near-surface snow and efficient use of the available upper boundary condition data.

A deeper lower boundary condition is possible, which would allow the use of a constant lower boundary condition rather than the seasonally varying one employed. However, a deeper lower boundary condition model with a constant temperature would substantially shorten the duration of our simulations because it would require many more surface data for model initialization. A deeper lower boundary condition could also induce a temperature bias in the model because the high interannual variability of mean annual surface temperature (standard deviation of about 1 K at the South Pole) is still felt at 10 m depth. Therefore, setting a constant lower boundary condition at 10 m depth, for example, would require assuming a temperature at that depth that is likely to be biased by ±0.3 K due to interannual variability of surface temperature.

The model snow properties are taken from Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others (1966) (Fig. 2b–d); they vary with depth but are kept constant in time. The heat capacity of the snow is 1710 J kg⁻¹ K⁻¹. It varies insignificantly with depth (Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others, 1966), so is not shown. The vertical resolution is 1 cm in the top 30 cm, then becomes step-wise coarser to 50 cm at 2 m depth and uniform at 50 cm from 2 to 6.5 m, as shown schematically in Figure 1. The model is constrained at its lower boundary by a climatological temperature gradient taken from the temperature-profile annual cycle of Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others (1966).

Fig. 2. Near-surface snow properties of the South Pole during the International Geophysical Year (IGY) (Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others, 1966) used in the finite-volume model of subsurface temperatures: (a) model initialization temperature profile (31 December 1957); (b) mean annual snow density profile; (c) mean annual thermal conductivity profile; and (d) mean annual thermal diffusivity (κ = k/ρC _p).

The upper boundary condition is a continuous time series of skin-surface temperature from 1994 to 2003 at 9 min resolution derived from routine measurements of upwelling infrared radiation collected by NOAA ESRL-GMD with a downward-looking broadband pyrgeometer (sensitive to wavelengths of 4–50 μm) with a ventilated dome, deployed 1 m above the snow surface. The infrared data were converted to skin-surface temperatures using the Stefan–Boltzmann law and a snow emissivity value of 0.98 (Reference WarrenWarren, 1982). The longwave emission from snow comes from the top millimetre of snow. Therefore, the skin-surface temperature retrieved from longwave data represents the temperature of the top millimetre. Comparison of skin-surface temperatures retrieved from longwave data to snow-surface temperatures measured by thermistors resting on the snow surface from the 2001 South Pole winter (Reference Hudson and BrandtHudson and Brandt, 2005) shows a high correlation (r ² > 0.99) and mean residual of less than 1 K. The longwave infrared dataset has 3 min resolution from 1994 to 1997, and 1 min resolution from 1998 to 2003. A 9 min running mean was applied to the entire dataset, and missing data were acquired using linear interpolation prior to the conversion to temperature.

The only energy input to the model comes from the upper and lower boundary conditions. Solar heating of snow occurs during summer, but we do not include it explicitly. Most of the absorption is of near-infrared radiation, which occurs in the top few millimetres of snow (fig. 4 of Reference Brandt and WarrenBrandt and Warren, 1993), within the topmost half-volume (5 mm) of the model. This absorption is of the same order as the effective infrared emission depth of snow (<1 mm). Therefore, it is largely included by using the skin-surface temperature from the upwelling infrared data. Inclusion of solar heating in the source term is not physically compatible with a skin-surface temperature upper boundary condition. Solar absorption of snow results in a small temperature maximum, as much as 0.2 K, a few millimetres below the snow surface during December (Reference Brandt and WarrenBrandt and Warren, 1993).

If our model were instead constrained at the surface by radiative and turbulent fluxes, then absorption of solar radiation in the snow as a function of depth could be included. We decided against this approach because of the large potential uncertainty associated with estimating sensible heat fluxes on short timescales in extremely stable boundary layers, as mentioned previously.

Heat can also be forced into the snow by wind pumping (Reference ColbeckColbeck, 1989). The extent to which this effect is significant depends on the square of the surface wind speeds and the height of the surface topography, both of which are greater on average during the winter. This process will likely serve to increase the effective thermal conductivity and pore-space water-vapour transport within the snow. Windpumping was found to be insignificant to heat transport below approximately 20 cm at the South Pole (Reference Brandt and WarrenBrandt and Warren, 1997), but may have a significant effect on snow temperatures closer to the surface under appropriate conditions (i.e. lowtemperature gradients or high wind speeds) (Reference Albert and McGilvaryAlbert and McGilvary, 1992).

Surface and subsurface melting has been modelled in Antarctic snow, but it has been found to be insignificant above 2500 m a.s.l. (Reference Liston and WintherListon and Winther, 2005). The elevation of the South Pole is 2835 m a.s.l. Snowmelt has never been reported at the South Pole; we therefore do not consider its effects here in our model.

Model behaviour, initialization and validation

The finite-volume model was tested in several ways to assess its behaviour, precision and accuracy. A simulation of a simple scenario is shown in the inset in Figure 3. The snow-pack is initially set to be isothermal at −40°C, with a constant surface-temperature forcing of −30°C and a uniform thermal diffusivity. Figure 3 shows the error in the numerical results after 2 days by comparison to the analytical solution for this scenario (Reference Carslaw and JaegerCarslaw and Jaeger, 1959, p.62–63). The figure shows that the errors in the numerical finite-volume model drop to the order of mK for a resolution of 1 cm and 2 min. These results indicate that the spatial resolution of the model in the top 30 cm is adequate to determine snow temperatures on short timescales under a steep temperature gradient, which often occurs in the top 30 cm of snow at the South Pole.

Fig. 3. Accuracy of the numerical solution with respect to the analytic solution as a function of depth and time resolution (K). The error is shown for day 2 of the simulation, shown in the inset. The conditions for this scenario were a constant upper boundary condition of −30°C and an initial isothermal temperature profile of −40°C.

We also tested whether the model leaked energy and, if not, how long the model would take to converge on a repeatable annual cycle. Energy leakage is not likely in a finite-volume model of this construction. Initialized in this test as isothermal at −49.5°C, the mean annual surface temperature at the South Pole, the model took 1 year to converge. Our definition of convergence here is that the mean absolute value of the temperature profile residual between one day during the last (tenth) year of the simulation and the same day during a prior year is less than 0.5 K. For the simulations presented below, the model was initialized using the 31 December temperature profile from Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others (1966) taken during their 1957/58 field season. Using a realistic profile of snow temperature facilitates much faster model spin-up.

The finite-volume model results are sensitive to the shape and magnitude of the seasonally varying temperature gradient used as the model’s lower boundary condition (not shown). The seasonally varying temperature gradient at 6.5 m from Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others (1966) is not a simple sinusoid or annually symmetric; it ranges from 0.95 K m⁻¹ in November to −0.5 K m⁻¹ in March. It does not integrate to zero. Thus, there is a bias of 0.2 K m⁻¹ in the lower boundary condition, which results in a positive absolute bias less than 0.5 K in the top 1 m of the model, but the bias is as large as 1.5 K at a depth of 5 m. Given the larger potential biases in the lowest 2.5 m of the model, we generally present data only for the top 4 m of the snow. Rather than adjust the observed lower boundary condition to average to zero annually, which would impose our own expectations on the length and/or magnitude of the cooling and warming seasons at 6.5 m, we accept that the measurements may cause a small bias in our results.

A further positive 1 K bias may exist in our results due to the use of a low snow emissivity in the skin-surface temperature retrieval from longwave upwelling fluxes (LUF). We used an emissivity of 0.98 (Reference WarrenWarren, 1982); however, the effective emissivity may actually be 0.99 or greater, including the reflection of atmospheric longwave emission (Reference HoriHori and others, 2006). A final bias, as stated earlier, exists due to the exclusion of solar heating of the snow. This bias is small and negative, −0.2 K (Reference Brandt and WarrenBrandt and Warren, 1993), in the top few millimetres during summer.

Uncertainties in temperature, heating rate and net heat flux into snow were determined through propagation of errors in the Stefan–Boltzmann law, heat-flux and heating-rate equations. Results are given for two different periods, as the long-wave upwelling data were available at 3 min frequency for 1994–97 and at 1 min frequency for 1998–2003. Uncertainties in the upper boundary condition and snow properties used in these calculations are: LUF = ±4 W m⁻² (personal communication from E. Dutton, 2004); k = ±0.1 W m⁻¹ K⁻¹ (Reference Brandt and WarrenBrandt and Warren, 1993); C _p = ±6 J kg⁻¹ (Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others, 1966); and ρ = ±35 kg m⁻³ (personal communication from H. Conway, 2007). For the period 1994–97, the 1σ errors in temperature, heating rate and net heat flux are ±0.8 K, ±6.4 K d⁻¹ and ±3.8 W m⁻², respectively. For the period 1998–2003 the 1σ errors in temperature, heating rate and net heat flux are ±0.4 K, ±3.6 K d⁻¹ and ±2.6 W m⁻², respectively. The values are for 9 min averages; random errors decrease substantially when averaged over hours or days.

The climate of the South Pole

The climate of the South Pole has been studied extensively since the International Geophysical Year (IGY, 1957/58). Aspects of the climate and weather of the South Pole relevant to this work are presented or reviewed by Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others (1966), Reference Schwerdtfeger and OrvigSchwerdtfeger (1970, Reference Schwerdtfeger1984), Reference CarrollCarroll (1982), Reference King and ConnolleyKing and Connolley (1997), Reference NeffNeff (1999), Reference Hudson and BrandtHudson and Brandt (2005) and Reference Town, Walden and WarrenTown and others (2005, Reference Town, Walden and Warren2007). It has been shown that the South Pole, and most of the East Antarctic plateau, has an annual temperature cycle with a coreless winter (Fig. 4), i.e. a winter with no well-defined minimum in temperature (Reference Warren and SchneiderWarren, 1996). Figure 4 shows that the day-to-day variability in 2 m atmospheric temperature is much greater during winter than during summer. The surface winds (not shown) are stronger on average during winter (approximately 6 m s⁻¹) than summer (approximately 4 m s⁻¹). However, daily variability in wind speed varies little throughout the year.

Fig. 4. Monthly mean 2 m atmospheric temperature for the South Pole for 1994–2003. The dashed curves show the standard deviation of daily average temperatures about the monthly mean.

The association between 2 m temperature, wind speed, cloud cover and regional weather patterns at the South Pole was examined by Reference NeffNeff (1999) and Reference Town, Walden and WarrenTown and others (2007). In general, cloud cover is associated with higher temperatures and higher wind speeds throughout the year. This association is due to alternating influences of relatively strong cyclonic weather systems originating from the coast and calm geostrophic flow circulating along elevation contours of the East Antarctic plateau. In general, the cyclonic activity advects heat, moisture and cloud cover to the South Pole from the Weddell and Bellingshausen Seas. The anticyclonic, geostrophic flow aloft is set up by radiative cooling of the surface to space. This cooling generates a downslope pressure gradient that is balanced by the Coriolis force. Near the surface, the geostrophic force balance is perturbed by friction, directing some of the flow downslope. This is known as an inversion wind (Reference SchwerdtfegerSchwerdtfeger, 1984), which is the source of katabatic winds at the coast of Antarctica.

The correlation of higher temperatures with cloud cover is stronger during winter than during summer. Episodic fluctuations in 2 m atmospheric temperatures of up to 30 K are possible during winter, whereas synoptically similar episodes result in much smaller temperature fluctuations during summer. This is due primarily to three factors. The Southern Hemisphere winter is typically stormier than the Southern Hemisphere summer (Reference Simmonds, Keay and LimSimmonds and others, 2003). Over the Antarctic plateau this means that there is more cyclonic activity to advect heat and moisture towards the South Pole during winter. In addition, the temperature gradient from the Southern Ocean to the continental interior is steeper during winter. Therefore, cyclones that do reach the South Pole will be warmer relative to ambient conditions than their summertime counterparts. Finally, the lack of solar heating during winter allows strong surface-based atmospheric temperature inversions under clear skies. The temperature inversions are very sensitive to changes in wind speed and cloud cover. Therefore, 2 m atmospheric temperatures fluctuate more on short timescales during winter than summer in response to similar changes in synoptic conditions.

Subsurface temperatures

The seasonal cycle of near-surface snow temperatures has been reported or utilized by many researchers around Antarctica: in Dronning Maud Land (e.g. Reference BintanjaBintanja, 2000; Reference Reijmer and OerlemansReijmer and Oerlemans, 2002; Reference Liston and WintherListon and Winther, 2005; Reference Van As, van den Broeke and van de WalVan As and others, 2005; Reference Van den Broeke, Reijmer, van As, van de Wal and OerlemansVan den Broeke and others, 2005, Reference Van den Broeke, Reijmer, van As and Boot2006), in West Antarctica (e.g. Reference Morris and VaughanMorris and Vaughan, 1994), at the Antarctic coast (e.g. Reference King, Anderson, Smith and MobbsKing and others, 1996), the South Pole (Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others, 1966; Reference CarrollCarroll, 1982; Reference JacksonJackson, 1982; Reference Brandt and WarrenBrandt and Warren, 1993, Reference Brandt and Warren1997; Reference McConnell, Bales, Stewart, Thompson, Albert and RamosMcConnell and others, 1998) and multiple other locations (e.g. Reference King and ConnolleyKing and Connolley, 1997; Reference Bailey and LynchBailey and Lynch, 2000). The short-term variability of near-surface snow temperatures has been illustrated by a subset of those above (Reference King, Anderson, Smith and MobbsKing and others, 1996; Reference McConnell, Bales, Stewart, Thompson, Albert and RamosMcConnell and others, 1998; Reference BintanjaBintanja, 2000; Reference Reijmer and OerlemansReijmer and Oerlemans, 2002; Reference Van As, van den Broeke and van de WalVan As and others, 2005; Reference Van den Broeke, Reijmer, van As, van de Wal and OerlemansVan den Broeke and others, 2005, Reference Van den Broeke, Reijmer, van As and Boot2006). However, beyond summertime investigations of the diurnal cycle, the short-term variability and the interannual variability of near-surface snow temperatures has yet to be quantified.

We first present the mean annual cycle of skin-surface temperature and the mean annual cycle of snow temperatures (Fig. 5a and b). These results compare well with the in situ data from 1958 (Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others, 1966) and the winter of 1992 (Reference Brandt and WarrenBrandt and Warren, 1997) within the interannual variability shown in Figure 5d.

Fig. 5. Mean seasonal temperatures, short-term variability and interannual variability for 1995–2003: (a) climatological skin-surface temperature from longwave upwelling radiation measurements, the upper boundary condition for the 9 years of simulation (°C); (b) 9 year averages of mean 2 week subsurface temperatures (°C); (c) 9 year averages of mean 2 week standard deviation of 9 min temperatures (K); and (d) 9 year standard deviation (K) of means in (b).

The mean short-term variability is shown in Figure 5c. Most of the temperature variability on short timescales is contained in the topmost metre. In general, the snow shows greatest variability during winter. The 2 week variability at 1 m depth is greatest during February–April and November–December. During these times, the snowpack is cooling or warming in response to the large seasonal change in atmospheric temperature as well as responding to synoptic variability. Interannually, greater variability is found during winter than during summer (Fig. 5d).

Figures 6 and 7 further illustrate the dynamics of subsurface temperatures on short timescales. The largest variability in temperature is at the surface, as indicated by the 5% and 95% extremes shown in Figure 6. The variability in temperature in the top 50 cm during January and July is due primarily to synoptic weather influences (i.e. changes in radiative and turbulent fluxes) at the South Pole, with changes in solar elevation playing a minor role during January. The spread in temperatures is larger in July than in January due to the sensitivity of the surface-based atmospheric temperature inversion to changes in the surface energy balance. During July, there is little variability below 1 m because there is no strong, consistent temperature gradient between the surface and the bottom boundary of the model to force significant changes in temperature. January shows some variability in temperatures between 2 and 4 m depth as the warm pulse from November and December diffuses into the snow. During March and November, however, the near-surface snow shows a broad distribution of temperatures due both to synoptic weather influences and changes in solar elevation.

Fig. 6. Thick curves show monthly snow-temperature profiles (°C) from the finite-volume model for January (a), March (b), July (c) and November (d) of 1996. The 5%, 25%, 75% and 95% distribution values are also indicated as thin solid and dashed curves (from left to right).

Fig. 7. Snow temperatures for the months of January (a), March (b), July (c) and November (d) of 1996 (°C). The time series in the panel above each contour plot is the skin-surface temperature used as forcing for that month.

Sub-daily variability in the top metre of snow is shown for these four months of 1996 in Figure 7. The panel above each contour plot shows the skin-surface temperature. The influence of synoptic variability on the temperature profiles is evident from these panels. The snow shows significant variability in the top 30 cm during summer (January) due to fluctuations of the skin-surface temperature on hourly, or longer, timescales. Below 30 cm, temperature variations with time and depth are less dramatic. Larger variability in skin-surface temperatures during winter leads to increased temperature variations in the topmost snow. These large temperature fluctuations penetrate deeper during winter due to the stronger surface forcing. The dual influences of dramatically changing solar irradiance and synoptic variability are evident in Figure 7 for March and November. The seasonal warming or cooling of the snow during these transitional months, illustrated as wide distribution boundaries in Figure 6, is shown as steadily sloping temperature contours in Figure 7b and d.

Figure 7 shows clearly that temperature gradients in the near-surface snow can change direction several times a month in response to variable synoptic conditions. This will be important in estimating short-term net heat fluxes into the snow, and may have significant implications in studies of snow metamorphism.

Net heat fluxes and heating rates

Using specified profiles of thermal conductivity, heat capacity and temperature, we computed snow heat fluxes and heating rates from the computed snow temperatures. The snow heat fluxes were calculated as linear gradients in snow temperature multiplied by the thermal conductivity of the snow at that depth. Heating rates were computed as the divergence of the snow heat fluxes into a volume. The net heat fluxes into the snow (G) are shown in Figure 8. They were computed from the topmost two temperatures in the finite-volume model and the thermal conductivity of the snow at that level. A positive G indicates a downward-directed heat flux into the snow (i.e. heating the snow).

Fig. 8. Monthly mean net heat fluxes into the snow (G) for 1995–2003 are shown by the thin black curves (W m⁻²). The 9 year mean of monthly mean G is shown by the thick black curve. The case-study months used in this paper are shaded. Positive G is directed downward into the snow. The month of January is repeated.

The monthly mean net heat fluxes into the snow shown here compare well with other fluxes modelled or observed at the South Pole (Reference Dalrymple, Lettau, Wollaston and RubinDalrymple and others, 1966; Reference CarrollCarroll, 1982; Reference King and ConnolleyKing and Connolley, 1997) and are similar in magnitude and timing to other sites around the continent, despite large differences in latitude, longitude, altitude and continentality (Reference King and ConnolleyKing and Connolley, 1997; Reference BintanjaBintanja, 2000; Reference Reijmer and OerlemansReijmer and Oerlemans, 2002; Reference Van den Broeke, Reijmer, van As, van de Wal and OerlemansVan den Broeke and others, 2005, Reference Van den Broeke, Reijmer, van As and Boot2006).

The seasonal cycle of net heat flux into snow has an unusual shape; it is not sinusoidal or symmetric in amplitude. The rate of energy exchange between snow and atmosphere is greatest during the months when the temperature difference between snow and atmosphere is greatest. The heat content of the snow is greatest during January, but the largest downward temperature gradient (upward flux of energy) occurs in March as the sun sets. In the multi-year mean, cooling continues throughout the winter. Heating of the snow during November and December is more dramatic and shorter-lived because the sun must rise to elevations large enough to significantly heat the snow and overlying atmosphere. The net heat flux into snow integrates to 0.07 W m⁻², which is zero within the standard error of the 9 year mean. The snowpack is therefore not changing temperature in the annual mean over this time period.

The monthly mean winter net heat flux into snow is negative (upward), and progresses steadily to zero from April to September in the 9 year mean (thick line in Fig. 8). However, it is clear that this steady increase is a result of multiyear averaging. The net heat flux into snow is consistently directed upward from February to April. After April, the near-surface snow has lost most of the memory of the previous summer, so the influence of heat advected in from the coast can dominate the monthly mean net heat flux into snow. In general, energy is still being drawn from the snow surface by turbulent fluxes and radiation loss to space throughout the winter. Any given month between May and September will therefore likely have a negative net heat flux into snow, but its magnitude depends more on the recent history of energy advection to the Antarctic plateau from the coast than on the remaining summer heat content of the snow.

In the monthly mean, turbulent and radiative energy transfer are larger components of the surface energy balance than the net heat flux into snow shown here, often an order of magnitude larger (Reference TownTown, 2007). The seasonal heat stored in the snow from the summer and released to the atmosphere during autumn is limited by the diffusivity of the snow. There is a significant amount of energy storage and energy reflux on shorter timescales, which consistently averages out on monthly timescales. The effect of the net heat flux into snow is shown as snow heating rates in Figure 9.

Fig. 9. Snow heating-rate climatology for 1995–2003: (a) net heat flux into snow (G), as shown in Figure 8 (W m⁻²); (b) 2 week mean heating rates (K d⁻¹); (c) 2 week standard deviation of heating rates (K d⁻¹); and (d) 1σ interannual variability about the mean shown in (b) (K d⁻¹). Note the different scales on the vertical axes of (b), (c) and (d).

Prior work on short-term variability of net heat flux into snow in Antarctica is limited primarily to summertime investigations of diurnal variations in the surface energy balance. In Dronning Maud Land, the diurnal cycle can range from ±10 W m⁻² to ±20 W m⁻², proceeding inland from the coast (Reference BintanjaBintanja, 2000; Reference Reijmer and OerlemansReijmer and Oerlemans, 2002; Reference Van den Broeke, Reijmer, van As and BootVan den Broeke and others, 2006). The diurnal range is largest under clear skies. It is attenuated under overcast skies, due to the cloud’s absorption of downwelling shortwave radiation and emission of longwave radiation (Reference Van den Broeke, Reijmer, van As and BootVan den Broeke and others, 2006). Variability in mean daily net heat flux into snow during summer is on the order of ±8 W m⁻² on the high plateau in Dronning Maud Land (Reference Van As, van den Broeke and van de WalVan As and others, 2005).

The South Pole has no 24 hour cycle of solar elevation, so all the sub-daily variability in the net heat flux into snow comes from variations in synoptic conditions such as cloud cover, wind speed and temperature. The effect of synoptic variability can clearly be large, sometimes as large as the diurnal solar variation at those sites. Mean variability (1σ) in 9 min net heat flux into snow ranges from a minimum in December of ±5 W m⁻² to a mean winter (April–September) variability of ±10 W m⁻². Histograms of 9 min net heat flux into snow for January, March, July and November of 1996 are shown in Figure 10. Figure 11 shows histograms of the resulting heating rates in the topmost 10 cm of snow.

Fig. 10. Histograms of net heat flux into snow (G) for four months of 1996 (W m⁻²). The data have 9 min time resolution and the bin width is 1 W m⁻². These distributions are representative of 1995–2003. The means and 1σ standard deviations are (a) 2.1 ± 5.6 W m⁻² for January; (b) −3.9 ± 8.1 W m⁻² for March; (c) 0.3 ± 10.7 W m⁻² for July; and (d) 4.0 ±6.1 W m⁻² for November.

Fig. 11. Histograms of 9 min heating rates averaged over the top 10 cm of snow for four months of 1996 (K d⁻¹). The data have 9 min time resolution, and bin widths are 2 K d⁻¹. These distributions are representative of 1995–2003.

The histograms in Figures 10 and 11 are asymmetric. The distribution is always skewed to larger net heat fluxes and heating rates, regardless of season, a phenomenon also noticed in the distribution of winter surface temperatures at Plateau Station (Reference Kuhn, Riordan, Wagner, Weller and BowlingKuhn and others, 1975). This is largely a result of synoptic activity, which at the ‘pole of cold’ can only heat the snow surface. The limit on the rate of downward heat flux is related to the rate that heat can be advected in from the coast. The upward (negative) heat fluxes are limited by the thermal conductivity of the snow, and the rate of radiative cooling at the surface. The radiative cooling rate drops quickly with decreasing temperature; it is proportional to the fourth power of temperature. Winter months show a wider spread than summer months due to the surface-based atmospheric temperature inversion.

Sub-daily (9 min) heating rates can vary by as much as 40 K d⁻¹ (1σ) at the snow surface (Figs 9c and 11). The sub-monthly variability in heating rates extends deeper in winter than in summer because the atmospheric forcings are larger and more variable during winter. The interannual variability of heating rates is largest close to the surface, and also larger during winter than summer (Fig. 9d). Below a depth of approximately 60 cm, the interannual variations have diffused into a mean, climatological heating rate.

The heating rates in Figure 11 are more symmetric about their means than the corresponding net heat flux into snow. The shape of the distribution changes with season as the atmospheric and radiative influences at the surface change. The narrow distributions for January and November in Figure 11 result from weaker synoptic activity and solar forcing of snow-surface temperatures during spring and summer. In terms of the variability of temperatures and heating rates, 1996 is representative of the 9 years modelled here. Despite variability within each month, the net heat flux into snow consistently averages out to very similar monthly mean values each year. This is unexpected because it requires the net synoptic influence on the surface to be approximately the same from year to year in a given month.

In Figure 12, we examine the same months from 1996 that were shown in Figure 7. Heating rates can be as large as 10 K d⁻¹ or greater in the bulk of the near-surface snow, but such heating rates are often followed immediately by significant cooling. The narrow panel above each contour plot shows net heat flux into snow. Whereas the monthly mean net heat flux into snow never exceeded ±6 W m⁻² in the 9 years simulated here, the 9 min values of net heat flux can be as large as ±20 W m⁻². The amount of energy refluxed each month between the snow and the atmosphere is greater during the winter than during the summer. Even though the net heat flux into snow on short timescales is large, the peaks are smaller than the largest mean daily excursions observed by Reference CarrollCarroll (1982). However, our results generally confirm the short-term behaviour of the net heat flux into snow estimated by Reference CarrollCarroll (1982).

Fig. 12. Snow heating rates for months of 1996 (K d⁻¹). The time series in the panel above each contour plot is net heat flux into snow (G) for that month (W m⁻²).

Interannual variability of heating rates (Fig. 9d) is of the same order as that of the climatological heating rates in the topmost few centimetres during March and November. The heating rates are much more variable interannually and at greater depth during winter without the dominating influence of the sun.

Subsurface vapour pressures

Understanding the seasonal cycle and variability of subsurface vapour pressures is important for a number of applications: snow metamorphism, snow microstructure and analysis of paleoclimate records. Here we present vapour pressures in the near-surface snow at the South Pole calculated from the subsurface temperatures by assuming saturation with respect to ice in the snow pore spaces (Goff and Gratch formula from Reference ListList, 1949).

Figure 13 shows statistics of pore-space vapour pressure for the near-surface snow. The maximum climatological vapour pressure is 60–70 Pa, which occurs at the surface during December and January. The amplitude of the seasonal cycle of vapour pressure drops off much more steeply than the subsurface temperatures in Figure 5 due to the non-linearity of the Clausius–Clapeyron relationship.

Fig. 13. Climatology of pore-space vapour pressure for 1995–2003: (a) climatological skin-surface temperature from longwave upwelling measurements (°C); (b) 9 year average of 2 week mean vapour pressures (Pa) in the top 2 m; (c) 9 year average of 2 week standard deviation of vapour pressures (Pa) in the top 30 cm; and (d) 1σ interannual variability (Pa) in the top 30 cm about the mean shown in (b).

It can be seen from Figure 13c that synoptically driven changes in vapour pressures and vapour-pressure gradients occur primarily in the top 10 cm of the snow. The relatively higher temperatures of summer and the changing vapour-pressure gradients during summer can cause significant snow metamorphism, consistent with the idea that most snow metamorphism occurs during summer.

Water vapour can be transported primarily by three mechanisms in cold snow: diffusion across temperature and snow-grain radius gradients; forced ventilation by surface winds (Reference ColbeckColbeck, 1989); and convection within the snow. Reference Brandt and WarrenBrandt and Warren (1997) found that forced ventilation (i.e. wind-pumping) is not significant at the South Pole below the top 20 cm of the snow. It may turn out to be significant at shallower depths. Even if significant ventilation is found at the South Pole, the snow may still be saturated with water vapour. However, forced ventilation does have implications for paleoclimate records in that the isotopic composition of water vapour in the pore spaces in the top 10 cm is likely a mixture of the atmospheric δ¹⁸O and the δ¹⁸O of the surrounding snow. The isotopic signature of the snow may therefore change after deposition due to ventilation and subsequent vapour transport within the snow (Reference Neumann, Waddington, Steig and GrootesNeumann and others, 2005; Reference TownTown, 2007), which is likely greater than changes in isotopic content due only to pore-space diffusion down isotopic gradients (e.g. Reference Johnsen, Clausen, Cuffey, Hoffmann, Schwander, Creyts and HondohJohnsen and others, 2000; Reference Helsen, van de Wal, van den Broeke, van As, Meijer and ReijmerHelsen and others, 2005, Reference Helsen2006).

Figure 13d shows significant interannual variability in subsurface vapour pressure due to interannual variability in synoptic forcing of snow temperatures. Figure 14 shows that synoptic variability affects vapour pressure in snow pore spaces. Vapour pressures on short timescales are extremely sensitive to synoptic activity at the surface during summer due to the relatively high temperatures. January 1996 experienced a change in vapour pressure at the surface of more than 30 Pa due primarily to variability in synoptic activity. The drop in temperature from January to March produces the low mean vapour pressures for March. Despite the large temperature swings during late autumn, winter and early spring, there is not much vapour activity in the snow in March or July because of the low temperatures. The latter half of November begins to show some significant vapour pressure, greater than 1% of the surface pressure (approximately 600–700 mbar), as the skin-surface temperatures rise with increasing solar elevation. The snow shows relatively high vapour pressures down to approximately 150 cm from December to February. Again, the variability in the monthly means is confined to the upper 10 or 20 cm of the snow during summer and adjacent months. Vapour pressures in the snow during the rest of the year are depressed by the low temperatures.

Fig. 14. Pore-space vapour pressures for: (a) January; (b) March; (c) July; and (d) November during 1996. The time series in the panel above each contour plot is the time series of skin-surface temperature for that month.