Introduction
In polar regions, where seasonal melt is virtually nonexistent, atmospheric chemical species become incorporated into the snow and firn, and decades of accumulation compact the firn into glacial ice. Ice cores drilled through glaciers in polar regions can then provide a record of changes in concentrations of chemical species over time-scales ranging from seasonal to glacial–interglacial transitions (Oeschger, 1985). For some chemical species, analysis of high-resolution ice cores (high-accumulation sites with no melt) makes interpretation on a yearly scale possible. For example, accumulation-rate inferences from ice-core data taken from Summit, Greenland, indicate that regional climate can change modes (from glacial to interglacial conditions) in as little as 1–3 years (Alley and others, 1993).
The process of air-to-snow transfer can filter and potentially distort atmospheric signals before they can be preserved in the glacial record (Gjessing, 1977). In order to interpret fully the chemical signals recorded in the ice cores, it is necessary to understand the processes by which those signals are transmitted and altered from the atmosphere, through the snow and firn finally to become incorporated into the ice. Currently, the main limitation on the use of many chemical records in the ice is a lack of understanding of the transfer functions among the air, snow, firn and ice. The purposes of this paper are to present progress on a mathematical/numerical model developed to simulate transfer through the snow and firn, and to use the model to investigate layering effects on ventilation. The model is unique in its ability to simulate the multi-dimensional effects of both diffusional and advective (air-flow) processes through the firn with heat, vapor and chemical transport.
While interpretation of the ice-core records has traditionally been based on the assumption that the process of diffusion controls gas movement in the snow and firn (Schwander, 1989; Schwander and others, 1993), there is evidence that convective air flow within the snow, induced by winds, can have significant influence on the aerosol record in ice cores (Cunningham and Waddington, 1993). Although Schwander and others (1993) employed a diffusion model, they also noted that the data cannot exclude the possibility of a convective zone below the surface. Early measurements by Benson (1962) in Greenland, and Dubrovin (1960) and Dubrovin and Petrov (1962) in Antarctica, yielded firn-temperature profiles that are not attributable to heat conduction alone, and led the investigators to attribute the profiles to advection. It has been shown that both turbulent winds over a flat surface and winds over surface relief can cause pressure perturbations to propagate vertically into the firn (Colbeck, 1989; Clarke and Waddington, 1991); these mechanisms may also cause significant lateral flow. Albert and McGilvary (1992b) demonstrated that temperature profiles resulting from windpumping arise from a balance between diffusive and advective heat-transfer processes and pointed out that, when the overall temperature gradient in the firn is strong, there could be significant ventilation that is not evident from temperature measurements.
Air flow due to ventilation affects chemical species transport. Cunningham and Waddington (1993) concluded that ventilation due to windpumping can account for the estimated present dry-deposition flux of non-sea-salt sulfate (NSS) at the South Pole. They found that changes in surface topography and wind speed that may have been possible in an ice-age climate could explain the increase of NSS in the Wisconsin age of the Vostok Station core without a change in the atmospheric concentration of NSS. Other species may also be affected by ventilation. For example, loss of hydrogen peroxide in the Bra after the first year of deposition may be attributable to ventilation of the firn at the South Pole (Neftel and others, in press). If, in fact, air flow in firn (ventilation) has a significant influence on chemical-species transport and reaction rates in firn then the magnitude and extent of ventilation effects need to be determined, and new theories will be needed for polar ice-core interpretation that include advective processes in addition to diffusion.
Model Description
The mathematical modeling described here builds upon an existing ventilation and transport model developed at CRREL (Albert and McGilvary, 1992a, b). Currently, the model is capable of simulating multi-dimensional air Slow through porous media with advective and diffusive heat transfer, including vapor transport and sublimation, and chemical-species transport. Algorithms for firn densification and grain growth, and chemical-species reaction will be added in the near future. Thermal-boundary conditions on the surface are driven either by a surface-energy budget or specified surface temperatures. Chemical-transport boundary conditions are driven either by known surface composition or known flux, The model solves the system of coupled partial differential equations that describe conservation of mass, momentum and energy in the firn, as follows.
The air flow through the snow or firn is described as Darcian flow and is driven by boundary conditions on pressure that may vary in space and time.
Here. P is pressure, k jk is permeability, t is lime, μ is air viscosity, υ j is the air-flow velocity and repeated indices of tensor notation for j and k imply summation (for dimensions 1, 2, 3).
The flow of heat is described by the time-dependent advection–diffusion equation with a source term to account for latent-heat effects of sublimation:
In this equation, T is temperature, ρ is density, ϕ is porosity, C is specific heal, λ jk is thermal conductivity and Q is the source term, which includes latent-heat terms and radiational heating effects. Subscripts a and s indicate air and snow, respectively.
Vapor movement through the snow is described by the advective–diffusive transport equation with a source term for sublimation or condensation:
Here, ρ v is vapor density in the air. D s is the diffusion coefficient for vapor in snow and S v is a source term for sublimation/condensation. The sublimation or condensation that can occur in snow and firn under given temperature and air-flow conditions is dependent on quantities that include grain-size and porosity. A unique aspect of the model is that it does not assume that the snow is saturated everywhere with water vapor and thus has the flexibility of predicting effects of unsaturated air flow should they occur, for example, in the upper parts of the pack. Undersaturated conditions induce sublimation, while supersaturated conditions are necessary for the grow of crystals.
The finite-element numerical technique is used in either two-dimensional or three-dimensional applications to solve the coupled partial differentia] equations that describe the flow of air, heat, vapor and chemical species. Finite elements are used so that the geometry of the modeled physical space can be arbitrary, allowing for simulations of irregular surface features such as sastrugi and non-uniform interior features such as snow and firn layers.
It is well known that mass and chemical transfer rates are greatly affected by the velocity of flow of the surrounding medium, in this case air. The equations for gas transport, including chemical reactions, are also described by advective–diffusive relationships (Conklin and others, 1993):
where
where C i and C a are the chemical concentrations in ice and air, respectively, θ a and θ i are the volume fractions of air and ice, x k is a coordinate direction (where k = 1, 2 or 3), D k is the dispersion coefficient in direction k, k f is the mass-transfer coefficient and K is the air-to-ice equilibrium partition coefficient. Repealed subscripts imply the summation convention. The dominant balance of terms in Equation (4b) indicates that the air-to-ice transfer process is slow relative to species movement in the air.
The meteorological influences of the air temperature, species concentration and relative humidity, wind speed, solar radiation and longwave radiation drive the heat-transfer solution at and near the surface, and measured meteorological data is usually used to drive the model (Albert and McGilvary, 1992a). The model can also be driven by specified temperatures.
Application: is layering important?
The transport of scalar quantities (heat, water vapor and chemical species) all follow advection–diffusion equations, as described above. The conventional “transfer function” for snow–air exchange in ice-core interpretation assumes that advective processes like ventilation do not exist, that exchange occurs because of diffusion alone. In contrast, theoretical results on ventilation (Colbeck, 1989; Albert, 1993; Cunningham and Waddington, 1993) have all made the simplifying assumption that the snow and firn form a homogeneous half-space. This assumption was used, for example, in an analysis that concluded that ventilation could completely explain the increase of aerosols in the Vostok Station ice core without any change in the atmospheric concentration (Cunningham and Waddington, 1993), a striking result. Both conventional and contemporary analysis of transport processes have largely ignored the layered nature of the snow and firn, yet it is well known that layering is a characteristic feature of the firn (e.g. Alley, 1988). In this section, measured firn properties are employed to investigate the relative importance of layering to ventilation processes. In particular; because wind-packed snow is commonly observed in polar regions, the effects of a less-permeable surface layer on ventilation are investigated. Measurements of the key firn property controlling ventilation, permeability, are discussed first.
Current field studies at Summit. Greenland, are focusing on detailed measurements of firn properties relevant to transport processes. We have found that the transport properties of the firn are highly inhomogeneous due to layering, and, especially in the case of depth hoar, have directional components (Albert and others, 1995). The key parameter that controls air flow through porous media is the permeability. Permeabilities were measured by cutting a sample of snow or firn in a cylindrical sampler (Shimizu, 1970), then pumping air through the sample while measuring flow rates and pressure drops using the method described by Chacho and Johnson (1987) and Hardy and Albert (1993). Darcy’s law is then used to infer permeability from measured air-flow rates and pressure drops. Permeability measurements taken in the top 3 m of firn in the summer of 1995 at Summit, Greenland, are shown in Figure 1. The measurements were taken on vertical firn samples and the lengths of the lines in the figure represent the thickness of the samples. In some cases, the samples encompassed more than one layer. Density measurements and preserved samples for quantitative microscopy were also taken and will be described in detail elsewhere. It is evident that the permeability varies by at least a factor of 10 over the top 3 m, with the surface wind-pack having much lower permeability in general than the underlying firn. In general, the permeability does not follow the firn-density profile (Albert and others, 1995) but is closely tied to the geometry of the pore space. For this investigation, the permeabilities assigned to the model vary with depth as indicated by the thin lines in Figure 1. A “mesoscale” layering is considered here. The current year’s accumulation (top 0.6 m) is considered as one layer with permeability 8 × 10−10 m2. Permeabilities for depths 0.6–1.1 m, 1.1–1.6 m, 1.6–2.0 m and below 2.0 m are assigned model values of 20 × 10−10 m2, 49 × 10−10 m2, 32 × 10−10 m2 and 27 × 10−10 m2, respectively. This mesoscale layering was observed, from several pit observations, to exist at Summit over a scale of at least kilometers. The finer-scale layering evident in Figure 1 has horizontal variability on the scale of meters, and so is not addressed in these model calculations, which cover distances of 5–10 m. Model results will be compared to results using the homogeneous half-space approach, where “representative” permeabilities of 30 × 10−10 m2 (representative of the top 3 m of firn) and 8 × 10−10 m2 (representative of the lower-permeability surface layer) are considered. In these simulations, the permeability within any one layer is considered to be isotropic (i.e. not dependent upon direction within the layer); future work will address directionality in local permeabilities as well as modeling layering on smaller and larger scales.
Fig. 1. Snow and firn permeability at Summit, Greenland. Thick lines depict measured values; thin lines represent values used in the numerical model.
For these simulations, a steady sinusoidal surface-pressure forcing with amplitude 5 Pa is imposed in simulations for two wavelengths, 1.7 and 3.3 m, for the surface-boundary condition. The amplitude and wave-lengths fall within the range of those measured over sastrugi under windy (10 m s−1 at 2 m height) conditions at Summit this summer. A 41 × 48 finite-element grid representing a depth of 3 m and surface distance of 5 or 10 m (corresponding to the 1.7 or 3.3 m wavelengths) is used. The element size increases linearly with depth, in order to give more numerical resolution in the near-surface region of higher gradients. A vector plot of a typical velocity field is shown in Figure 2; for clarity in the plot, only approximately one-third of the calculated velocities are shown. Because of higher pressure gradients in the near-surface firn, the air-flow velocities are higher there and generally decrease with depth. Although not depicted here, contour plots of the velocity magnitudes over the two-dimensional region reveal that there is little variation in magnitude over distance at any given depth: most of the variation in velocity magnitude comes with depth in the firn.
Fig. 2. A typical calculated air-flow field in response to sinusoidal pressure forcing on the snow surface.
The first set of .simulations investigates the effect of mesoscale layering on calculated ventilation velocities. In order to compare velocities in the firn, the velocity magnitudes, defined as the square root of the sum of the squares of the directional components for each vector, are plotted as a function of depth for layered firn in Figure 3. Because, in general, longer surface-pressure forcing wavelengths are known to induce lower near-surface flow but larger flow at depth (Colbeck, 1989), two wavelengths are investigated. In both cases, it can be seen that velocities in a more-permeable layer can be higher than the adjacent velocities in a less-permeable layer that lies above it, as evidenced by rises in velocity magnitudes at depths between 0.6–1.1 m and between 1.1–1.6 m. Thus, channeling, or increased horizontal and net air flow within a layer, can occur even if the microstructure of that layer is isotropic. This occurs because the pressure gradients that drive the air flow penetrate the firn and although the pressure gradients generally decrease with depth, variations in permeability allow variations in flow. Higher horizontal fluxes may then occur in layers of higher permeability. This effect is most pronounced for a steady-state surface-pressure forcing and so will most likely occur in the field when winds are sustained. At Summit and many other polar sites, once windy conditions develop, they can last for days or longer. Although the surface features such as sastrugi migrate over periods of days, they generally comprise the top 10 cm or so of snow; the more-permeable layers involved with channeling discussed here do not move in time.
Fig. 3. Calculated velocity magnitudes for the layered firn cases.
The next set of simulations consider the firn as a homogeneous hall-space using a “representative” permeability throughout the region. The first uses 30 × 10−10 m2, a value that could be judged to be representative of the top 3 m of the firn (Fig. 1) and the second uses 8 × 10−10 m2, a value that characterizes the firn as homogeneous using the surface permeability; both assumed a wavelength of 3.3 m. The results are shown in Figure 4, along with results from the layered case for reference. In this figure, the open symbols represent results from the homogeneous firn, while the solid squares depict results from the layered case. It can be seen that, in this case, the use of the larger permeability, although representative of most of the lop 3 m of firn, greatly overpredicts the air-flow velocities at all depths. The effect of the lower-permeability layer at the surface serves to decrease the air flow in the underlying firn.
Fig. 4. Calculated velocity magnitudes vs depth. Open circles are from uniform firn calculation using a permeability of 30 × 10−10 m2; open squares used a uniform permeability of 8 × 10−10 m2. Solid squares are calculations of layered firn with permeabilities shown as thin lines in Figure 1. Channeling (increased air flow) occurs in the more-permeable buried layers.
A better agreement with layered results can be had for this case by employing the properties of the surface layer (a permeability of 8 × 10−10 m2) for homogeneous simulation, as is evident from Figure 4. However, the homogeneous case still overpredicts the velocities in the surface hiver and underpredicts the velocities in the underlying layers. This is shown more clearly in Figure 5, an illustration of the normalized velocity difference, which is calculated as the quotient of the difference between the velocity magnitudes in the uniform and layered firn divided by the magnitude in the layered firn for each depth. In this case, the homogeneous model overestimates velocities in the surface snow by as much as 80% but underestimates most of the velocities at depths approaching the wavelength of the surface sastrugi by 20–60%.
Fig. 5. Normalized velocity difference between the layered case and the case of uniform firn, all having properties, of the surface windpack. The uniform firn assumption overestimates the flow in the surface windpack but underestimates the flow at depths approaching the wavelength of the surface sastrugi.
This result indicates that detailed modeling of firn properties in models of ventilation with chemical transport is necessary for accurate assessment of the firn–air transfer function for polar ice-core assessments. Estimates based on uniform firn properties (the homogeneous half-space approach) could be significantly in error. In addition, it does not necessarily follow that higher wind speeds (either at a different location or at the same location at a different time) will always induce more ventilation. A complicating factor is that higher winds induce more saltation of snow across the surface and the creation of a lower-permeability surface wind-pack which reduces the permeability, However, even with the reduced air flow due to the moderating effects of the surface wind-pack, the ventilation velocities in this case are on the order of several mm s−1 at a depth of 1 m which is still significant in terms of accelerating transport processes over that due to diffusion alone (e.g. Bales and others, 1995). Therefore, conventional estimates of chemical transport due to diffusion alone are likely to underestimate transport, while estimates of ventilation that consider the firn as a homogeneous half-space may overestimate ventilation effects, especially in the near-surface. The firn layering is important for ventilation.
Conclusion
A finite-element model for simulating multi-dimensional air flow with heat, mass and chemical-species transport through firn has been presented. It is driven by either surface energy-balance measurements or prescribed boundary conditions, and is able to represent arbitrary geometries with spatially variable material parameters such as permeability and thermal conductivity. The model is used here to investigate the importance of layering on air lion through the snow and firn.
Although layering has been neglected in all previous estimates of air movement through snow and firn layering is important. Field measurements of permeability at Summit, Greenland, show variations by at least a factor of 10 over the top 3 m, with the surface wind-pack having much lower permeability in general than the underlying firn. The effect of a lower-permeability surface layer is to decrease the air flow in the underlying firn (i.e. decreased over that which could be expected in the absence of the wind-pack), yet there is still sufficient air flow in the fop meters of the firn so that ventilation must be considered in chemical transport modeling. In addition, channeling, or increased horizontal and net air flow in a deeper layer overlain by a less permeable layer, can occur even if the microstructure of each layer is isotropic. Thus, while air-flow velocities in a low-permeability surface wind-pack may be small, under conditions of sustained pressure forcing it is possible to have higher velocity air flow in a more permeable deeper layer than in the wind-pack.
Conventional estimates of chemical transport due to diffusion alone are likely to underestimate transport, while estimates of ventilation that consider the firn as a homogeneous half-space (without layering effects) may overestimate ventilation effects at the near-surface. Effects of firn layering are important for ventilation and must be considered for accurate assessment of firn–air transport mechanisms.
Acknowledgements
I thank E. Arons, S. Colbcck and A. Gow for their helpful reviews of this paper. I am indebted to S. Decato for his assistance with the field measurements, and I thank F. Perron for help with the instrumentation. This research was made possible by U.S. National Science Foundation gram NSF-OPP # 9321420 and by Department of the Army basic research work unit AT24-SSE09.
