Calculated surface melt and lysimeter discharge were used to establish relationships between specific surface meltwater fluxes and meltwater wavespeeds through the snowpack at a continental, alpine site in 1997. These relationships were used to determine a pore-size distribution index (e ) and to calculate equivalent permeability (k _e) for the snowpack. Wavespeeds, and equivalent permeabilities derived from them, exhibited a marked decrease over the melt season and were lower then those reported in the literature. Although not a significant trend, e values seemed to increase over the melt season. The relationship between meltwater fluxes and wavespeeds was found to be dependent on the range of meltwater fluxes used. At high values the relationship is more similar to that predicted by a gravity-flow dominant regime, but breaks down for lower meltwater fluxes when capillary-flow may become more important. These results suggest that snow hydrology models need to incorporate capillary-flow processes when melt rates are less than 1 mm hr^-1. However, modeling meltwater flow as gravity-flow at higher melt rates appears to be correct. The switch between these two flow regimes will occur diurnally, with capillary-flow during the morning, and gravity-flow dominant later in the afternoon. This would also partially explain daily patterns of chemical concentration in the meltwaters.

Attempts to model meltwater flux through a heterogeneous, layered snowpack include a variety of approaches: (i) modeling flow around individual ice layers (Colbeck, 1973); (ii) simplifying the effect of ice layers and preferential flowpaths by using multiple flowpath models (Marsh and Woo, 1985); and (iii) treating the snowpack as a homogeneous anisotropic media (Colbeck, 1975). Such models provide insight into flow processes within the snowpack but, without knowledge of the extent and thickness of snowpack horizons and layers and the rates at which their properties change over time, become impractical in nearly all circumstances. An outstanding question in snow hydrology modelling is how to determine an appropriate value for the saturated permeability (k_s) of layered snowpacks. Point measurements of k_s have little meaning if the permeability is used to predict the rate at which meltwater moves through the snowpack, encountering numerous layers with different individual permeabilities. Most studies, whether analytical or numerical, are forced to treat the snowpack as an equivalent, uniform medium in which permeability is considered as a bulk or integrated property of the entire snowpack.

Unlike soils, snow permeability is almost impossible to measure in the laboratory, because of the triple phase nature of snowpacks and changes in density, grain size, type and distribution they experience during melting. As suggested by McGurk and Kattelmann (1986), measuring an equivalent permeability (k_e) of an entire snowpack in the field overcomes these problems and allows the detection of bulk changes during the melt season. k_e would be a useful parameter for a variety of approaches to modeling meltwater flow through isothermal snowpacks, including gravity-flow models (Colbeck, 1972, 1979; Dunne et al., 1976; Marsh and Woo, 1985) and numerical models (Jordan, 1983b).

The only field method for determining k_e uses measurements of surface meltwater flux (U) and the velocity at which this flux moves through the snowpack ([dz/dt]_U) (Marsh, 1991). When these data are plotted for periods of declining flux, the intercept of the line can be used to determine a value for k_e over the entire snowpack depth, and the slope of the line can be used to determine the pore size distribution index (e) of the snow. e is the exponent in the relationship between relative permeability (k_w) and k_s and effective saturation (S*) for unsaturated snow:

Where k_w is relative, or liquid water permeability, k_e is snowpack equivalent permeability and S* is effective water saturation of snow, expressed as a fraction of pore volume occupied by mobile liquid.

e has been determined from drainage experiments, with reported values ranging between 1.5 and 4.8 (Marsh, 1991). Both e and k_e depend upon physical properties of the snowpack, and may change considerably during the melt season as snowpack metamorphism affects grain size, snow density and ice layer thickness and continuity (Marsh, 1991). However, previous work has not been able to establish any relationships between either e or k_e and snowpack properties or amounts of metamorphism (Denoth et al., 1979), although it has been suggested that k_e is dependent upon the range of U used in calculations (McGurk and Kattelmann, 1986).

Perhaps most intriguing, McGurk and Kattelmann (1986) report that k_e decreased over time in a Sierra Nevada snowpack. The decrease in k_e with time and a thinning snowpack is counter-intuitive since one would expect that k_e would increase with time as grain size becomes larger and more uniform under a wet metamorphism regime (Colbeck, 1973). It is unclear whether the values for k_e and the change in k_e over time reported by McGurk and Kattelmann (1986) were site specific or whether they apply to other sites and climatic regimes.

Here we estimate values for k_e from a continental, mid-latitude site in the Colorado Rocky Mountains. Measurements of surface meltwater flux, discharge through snow lysimeters, and snowpack properties were conducted in 1997. Specific objectives of our research were to: 1) establish a range of values for e and k_e at this site; 2) investigate whether they vary consistently in a way that can be used to parameterize gravity-flow and numerical models during the melt season; 3) investigate whether k_e is dependent upon the range of U used in calculations; 4) evaluate whether this new information can provide insights into the physical processes that govern meltwater flow through snow.

DQ_S + DQ_M = Q* + Q_H + Q_E + Q_G + Q_R Eq. 2

where DQ_S is the cold content of the snowpack, DQ_M is the latent heat flux during melting or freezing, Q* is the net all wave radiation flux, Q_H is the sensible heat flux, Q_E is the latent heat flux, Q_G is the ground heat flux, and Q_R is the heat advected by precipitation. Cline (1997a) found Q_G and Q_R were not significant heat exchanges until the snowpack became very thin, and so are not considered any further here. The sign convention used was that energy transferred from the snow surface to the atmosphere was negative and energy transferred from the atmosphere to the snow surface was positive (Marks et al., 1992). The accuracy of the instruments and their proximity to the lysimeters means that surface meltwater flux could be modeled with good precision over a wide range of fluxes and short time intervals, a considerable improvement over previous studies such as Jordan (1983a) and McGurk and Kattelmann (1986).

Where e is the pore size distribution index; k_s is the saturated permeability; f* is the effective porosity, a is the flow factor (5.47 x 10⁶ m^-1 s^-1 at 0^oC).

The pore-size distribution index (e) for the snowpack was calculated for each day on which [dz/dt]_U values were calculated. The slope (b) of log-log plots of [dz/dt]_U against U gives the flux exponent (e-1)/e in equation (2), such that e = 3 when b is 0.67. Mean U for all values used in calculating the [dz/dt]_U - U regression on an individual day was substituted into the calculated regression equation to give [dz/dt]_meanU. Values of mean U and [dz/dt]_meanU were then substituted into equation (3), and given a value for e, and that a is a constant, rearranging equation (3) gives k_s^1/ef*^-1. This single value characterizes an important property of the snowpack, which Colbeck (1978) suggested should have typical values assigned, and should then be used as a parameter in operational forecasting, changing as the snowpack matures. Having calculated effective porosity (f*) for a given day using density and irreducible water content values following the procedure of McGurk and Kattelmann (1986), the saturated permeability k_s for the snowpack can be determined from k_s^1/ef*^-1. As these k_s values were determined from wavespeeds calculated for the entire snowpack depth, they represent an 'averaged' value for the snowpack, an 'equivalent permeability' (k_e) as defined in Marsh (1991).

Fifteen values of e were calculated from the 15 flow exponents b (equation 4), which is equivalent to (e - 1)/e value (equation 3). Values of e ranged from 1.16 to 5.58, with a mean of 2.24 and standard deviation of 1.10 (Table 1). Removing the outlier of 5.58 on 17 June, e values ranged from 1.16 to only 3.32. Perhaps more important, e was found to be 2.93 when all surface meltwater fluxes > 1 mm hr^-1 were used to determine a value for b (discussed above as representing an apparent inflection point in the [dz/dt]_U – U relationship, Figure 7). This is similar to the value of 3 often used as a default when modeling meltwater flow through snow. A simple linear regression shows that e increased with time especially with the outlier value of 5.58 removed, but the relationship remained non-significant ( a = 0.05).

Fifteen values of k_e were then calculated from the 15 values of e . A mean wave wavespeed ([dz/dt]_{mean U}) was calculated for each day that wave speed analysis was conducted, using equation 4 and mean U for that day. The effective porosity (f*^-1) was calculated from measured snowpack density, assuming an irreducible water content of 0.07, following McGurk and Kattelmann (1986). These values were then substituted into equation (3), which was then solved for k_s. As k_s was calculated from wavespeeds over the full snowpack depth, it represented an integrated property of the entire snowpack, and was thus an 'equivalent' permeability (k_e). The calculated k_e ranged over two orders of magnitude, from 3.23 x 10 ^-12 to 6.57 x 10 ^–10 m², with a mean of 6.14 x 10 ^–11 m² (SD = 1.66 x 10 ^–12 m²) (Table 1). Such a large range in k_e was in part caused by the large range in e values calculated previously. Because of the importance of e in equation (3), any variation in e was magnified greatly when k_e was calculated using these values. We also calculated k_e when e was held constant at 3 in order to allow comparison with k_e reported in the literature (Table 1). Holding e constant greatly reduces the variability in k_e. With e = 3, the range of k_e was reduced by an order of magnitude to 1.08 x 10 ^–11 m² to 2.98 x 10 ^–10 m², with a mean value of 1.11 x 10 ^–10 m² (SD = 8.22 x 10 ^–11 m²).

Because snowpack properties such as grain size change over the snow melt season, k_e may be expected to change as well. A regression of k_e versus time using the measured values of e showed no significant change over time (R² = 0.03; p > 0.05). However, when the pore size distribution index was held constant at e = 3, the calculated k_e decrease over time was significant at the a = 0.05 significance level (R² = 0.40; p = 0.02). With e = 3, k_e declined from 3 x 10 ^–10 m² to 1 x 10 ^–11 m², a 96% reduction in 30 days (Figure 8).

Our values of [dz/dt]_{2mm hr}^-1 range from 61 to 176 mm hr^-1 and are low compared to those reported in the literature. For example, McGurk and Kattelmann (1986) reported [dz/dt]_{2mm hr}^-1 values of about 360 mm hr^-1in a Sierra Nevada snowpack in March 1984, declining to approximately 180mm hr^-1 by May. They report even higher values in March 1985 (approximately 540 mm hr^-1) with a decline by May to values similar to those of the previous year. Jordan (1983a) found [dz/dt]_{2mm hr}^-1 values between approximately 250 and 300 mm hr^-1in late June/early July, at an alpine site in British Columbia. Colbeck and Davidson (1973) calculated [dz/dt]_{2mm hr}^-1 rates of approximately 288 and 576 mm hr^-1 from drainage experiments with 0.6 and 1.8 m columns of homogeneous snow and with densities between 600 and 650 kg m^-3. Colbeck (1977) and Colbeck and Anderson (1982), using data from a 55.7 m² snow lysimeter at the Central Sierra Snow Laboratory collected in April 1954 (US Army, Corp of Engineers, 1956), calculated [dz/dt]_{2mm hr}^-1 values of approximately 685 and 755 mm hr^-1. Similar values were reported in March 1973 from a snow lysimeter 0.7 m² in area located in the Sleepers River Research Watershed, Vermont (Anderson et al., 1977). All of these studies were conducted in predominately maritime climates near either the west or east coast of North America.

Similarly, despite the large range in k_e calculated in this study, there was only limited overlap between them and values calculated at other field sites. Many daily values reported here were over an order of magnitude smaller than the minimum values reported elsewhere. Denoth et al. (1979) report some of the highest k_s values in the literature, ranging between 1.0 x 10^-9 and 1.0 x 10^-8 m² for disaggregated wet snow with a porosity of approximately 0.48-0.49 and grain size > 1 mm. These values increased to 2.5 x 10^–8 m² for disaggregated, wet snow with a porosity of 0.52 and grain size of approximately 3 mm. As k _s values represent a rather different snowpack property than k _e values, as already discussed, the limited overlap in ranges, and lower k _e values might be expected. The other reported permeability values used here were calculated in a way as to be more similar to k _e as we have defined it, so direct comparison and altered notation is more appropriate. Somewhat lower values of k_e equal to 1.0 x 10^-10 m² were reported by Colbeck and Davidson (1973) for wet snow with a porosity of 0.32 and density of 650 kg m³, with a slightly higher value of 3.0 x 10^-10 m² for a higher porosity of 0.36 and lower density of 620 kg m³. For a spring snowpack in the Sierra Nevada with a porosity of 0.47 and density of 480 kg m³, more similar to the snowpack properties found in the Colorado Front Range, Colbeck and Anderson (1982) found k_e was 2.6 x 10^-9 m², and a higher value of 3.0 x 10^-9 m² for the Vermont snowpack, with a porosity of 0.55 and density of 410 kg m³. Jordan (1983a) does not determine a value of k_e from eak_e^1/ef*^-1, which he found to approximately vary between 70 and 250, the equivalent of over a range of magnitude in k_e. Calculating k_e from the data presented by Jordan (1983b) and using the standard values adopted by Jordan (1983b) of e = 3 and f* = 0.43 (from density measurements), for eak_e^1/ef*^-1 = 130, k_e = 3.29 x 10^-10 m². Over their two year study period, McGurk and Kattelmann (1986) report k_e ranging from 3 x 10^-10 m2 to 1.8 x 10^-8 m².

There is no obvious single reason for the generally lower values of [dz/dt]_U and k_e reported here compared with elsewhere. The lower values may be explained by a combination of factors, including physical properties of the snowpack, and characteristics of the continental, alpine climatology of the site. Compared with some investigations, the level of heterogeneity in our alpine snowpacks may have been much higher. For example, Colbeck and Davidson (1973) used repacked, homogenized snow columns in their experiments. Denoth et al. (1979) used glacier firn. Both of those experiments were likely to involve lower grain size variations and fewer ice lenses, resulting in higher wavespeeds and higher k_e.

The 50% reduction in [dz/dt]_{2mm hr}^-1 values over the 1997 season observed at this continental, alpine site was similar to the decrease reported for two seasons by McGurk and Kattelmann (1986). The decrease in wavespeed (Figure 6) and the marked reduction in k_e (Figure 8) over time may be caused by the break down of ice layers in the snowpack. Furbish (1988) suggested that ice layers in a snowpack can reduce transit times. Ice layers may act to increase wavespeeds because saturated flow over the impermeable ice layers has a much greater velocity than unsaturated flow. This leads to a rapid concentration of meltwater in a limited area of the snowpack, and once this water breaks through the ice layer, k _w of the snowpack below will be markedly higher because k _w increases exponentially with increasing effective saturation (Colbeck, 1978). This will result in preferential flowpaths, with local meltwater travel times much faster than in the surrounding snowpack. Preferential flowpaths are known to occur frequently in the snowpack on Niwot Ridge (Williams et al., 1999a), where their hydrological effects are a subject of on going investigation.

Wankiewicz (1979) attempted to quantify the effect that concentration of meltwaters caused by initial ponding over ice layers may have on flow. He suggested that preferential flowpaths would yield travel times faster by a factor F than those in the uniformly wetted case. He further proposed that F would be related to e as a power function, such that meltwater would move faster by F ⁽e-1)/e for a snowpack with an initial water excess and and F if there is a water deficit. Once the snowpack becomes "ripe", the increase in vertical flow velocities which resulted along these preferential pathways may be sufficient to compensate for the delay resulting from horizontal flow above the ice layers (which would be relatively short if saturated flow occurred). The ice layers would have the overall effect of reducing transit times for a meltwater flux moving from the surface through the entire snowpack depth to the lysimeters at the snowpack base. Consequently, we would have calculated higher wavespeeds and k_e values in our study.

As ice layers breakdown during the melt season, these preferential pathways are lost, resulting in more uniform meltwater movement through the snowpack, with decreased relative permeability due to reduced levels of saturation. This would cause longer transit times for the same snowpack depth, resulting in decreased k_e. Our measurements of snow stratigraphy lend qualitative support to this argument. On 13 May the snowpack was characterized by four distinct ice layers. Over time these ice layers decayed and disappeared. This corresponds with the observed decrease in effective permeability.

The trend of increasing e over time is also consistent with accelerated movement of meltwater through preferential flow channels early in the melt season followed by a more uniform and slower wave front later in the melt season. Increasing e over time indicates that the distribution of pore-sizes increases over the same time period. This is opposite of what may be expected given that it is generally accepted that wet snow metamorphism leads to increasingly homogeneous, large, and rounded snow grains (e.g. Colbeck, 1986). However, if meltwater flow was predominately through preferential flow channels early in the melt season, wet snow metamorphism in these channels may have resulted in a relatively small range of pore sizes. The change to matrix flow later in the melt season may result in meltwater coming into contact with a greater range of pore sizes. Thus, the range of effective pore sizes that liquid water moving through the snowpack comes into contact with may increase even while the range of total pore sizes in the snowpack is decreasing, as the two would only be the same if the snowpack was saturated and so a two phase system.

An alternative suggestion is that the apparent trend in e values is an artifact of our experimental methodology. Peak melt rates, and thus the range of U values used in analysis, were slightly lower earlier in the melt season. As the relationship between [dz/dt]_U - U is dependent on the range of U (Figure 7), and e values were calculated from these relationships, it is possible that this small change in the overall nature of the relationship results in an increase in e over the course of the melt season.

Other studies have reported the dependence of the [dz/dt]_U - U relationship on the range of U values (e.g. McGurk and Kattelmann (1986)). However, there has been no attempt at explanation, or discussion of its implication for the use of simple gravity-flow models for describing meltwater movement through snowpacks. In our study there was no significant relationship between [dz/dt]_U - U for surface meltwater fluxes below 1 mm hr^-1. However, there was a clear relationship between these two variables when higher fluxes were included. For the [dz/dt]_U - U relationships at higher ranges of U, the slope of the relationship tends to be higher, and when used to calculate e values, gives results more similar to those reported in other snow studies and values reported for soil textures with similar grain sizes as snow (Muelem, 1978). For example, in our study the regression line of the [dz/dt]_U - U relationship for all wavespeeds calculated in 1997 from U values greater than 1 mm hr^-1 was y = 68.457 x ^0.608, giving e = ~ 3, the value used in most snow hydrology studies. As these e values derived from calculation and laboratory experimentation in other snow and soil hydrology studies generally assume gravity-flow, it appears that only at relatively high surface meltwater fluxes, >1 mm hr^-1, can meltwater movement be accurately described by this simple model.

This change between flow regimes may also occur diurnally. There may be a switch from capillary-flow and low lysimeter discharge early in the day, to gravity-flow and high lysimeter discharge when meltwater fluxes exceed 1 mm hr^-1. These diurnal changes in flow regime may explain daily patterns of chemical concentration in meltwaters. For example, Bales et al. (1993) and Harrington and Bales (1998) reported that the conductance of snowpack meltwaters increased on the recession limb of the diurnal hydrograph and then decreased sharply on the rising limb. This linkage of chemistry and hydrology is consistent with capillary water at low U rates flowing around grain boundaries and removing impurities, while at high U rates there is gravity-flow which has been moving through snowpack pore-space with limited interaction with individual grains and hence lower conductance.

Surface meltwater fluxes and lysimeter discharge were used to calculate wavespeeds of vertical meltwater movement at a continental, alpine site. These wavespeeds, and equivalent permeabililties for the entire snowpack depth derived from them, were found to be considerably lower then those reported in the literature. Similar to other studies, they exhibited a marked decrease over the melt season, possibly due to the declining influence of ice layers. Although not a significant trend, e values seemed to increase over the melt season. This is contrary to the change expected given simple wet snow metamorphism, and is interpreted to suggest a change from meltwaters moving through predominantly preferential flow channels early in the melt season to a matrix-type flow as the melt season progresses. The nature of the relationship between U and [dz/dt]_U , and thus k_e and e values, were found to be dependent on the range of [dz/dt]_U values used. At high U values the relationship is more similar to that predicted by a gravity-flow dominant regime, but breaks down for lower U values. At these lower values capillary-flow may become more important.

a flow factor constant = 5.47 x 10⁶ m^-1 s^-1 at 0^oC

• pore-size distribution index in the k-S* relationship

k_w relative, or liquid water permeability (m²)

k_e snowpack equivalent permeability (m²)

k_s saturated permeability (m²)

f* effective porosity