Observations of capillary barriers and preferential flow in layered snow during cold laboratory experiments

We carried out observations of dyed water infiltration in layered snow during cold laboratory experiments. We considered three different finer-over-coarser textures and three different water input rates. Nine samples of layered snow were sieved and subjected to constant supply of tracer at 0◦C. By means of visual inspection, horizontal sectioning, and measurements of liquid water content, capil5 lary barriers and associated preferential flow were characterized. The dynamics of each sample were also simulated solving Richards equation within the 1D multi-layer physically-based SNOWPACK model. Results reveal that capillary barriers and preferential flow are relevant processes ruling the speed of liquid water in stratified snow. Both are marked by a high degree of spatial variability at cm scale and complex 3D patterns. During unsteady percolation of liquid water, observed peaks in bulk 10 volumetric liquid water content (LWC) at the interface reach ∼ 33 36 vol% when the upper layer is composed by fine snow (grain size smaller than 0.5 mm). A comparison with expected wetness at water entry suction suggests that LWC might locally reach saturated conditions. Spatial variability in water transmission increases with grain size, whereas we did not observe a systematic dependency on water input rate for samples containing fine snow. The comparison between observed and simu15 lated LWC profiles reveals that the implementation of Richards equation reproduces the existence of a capillary barrier for all observed cases and yields a good agreement with observed peaks in LWC at the interface between layers.

Water flow in snow emerges as a complex, 3-D process when observed in the field or in the laboratory.The coexistence of water and ice grains causes fast metamorphism (Brun, 1989) and phase change, hence melt-freeze and the possible development of ice lenses when water infiltrates in subfreezing snow (Pfeffer and Humphrey, 1996).The interaction with topography can redistribute water at slope scale (Eiriksson et al., 2013).Moreover, water movement is usually marked by high spatial variability due to the occurrence of preferential flow or fingering (Marsh andWoo, 1984a, 1985;McGurk and Marsh, 1995;Schneebeli, 1995;Waldner et al., 2004;Williams et al., 2010;Katsushima et al., 2013).

F. Avanzi et al.: Capillary barriers in snow
These processes complicate the modeling of liquid water in snow, which has been often simplified in the past by using a simple Darcian theory that neglects effects of capillary gradients on the flow (Colbeck, 1972;Wankiewicz, 1978).In particular, fingering is deemed to play a key role in ruling water arrival time at snow base, hence runoff (Wever et al., 2014(Wever et al., , 2015) ) and snowpack stability (Techel et al., 2011;Mitterer et al., 2011;Mitterer and Schweizer, 2013).
The exact physics of preferential flow in snow is still not known (Katsushima et al., 2013) and modeling strategies are therefore still preliminary.For instance, Marsh and Woo (1985) propose an explicit definition of multiple-path routes, whereas Katsushima et al. (2009a) introduce a threshold value for θ triggering preferential runoff.Wever et al. (2014) report that solving the Richards equation accounting for suction (ψ) gradients improves runoff estimations at different temporal resolutions, but it also accelerates meltwater front progress when compared with data from an upGPR and simulations by a bucket scheme (Wever et al., 2015).This result has been attributed to an unexpected simulation of some effects of preferential flow by the water scheme used.
The use of the Richards equation for modeling wetting front instability in porous media is still a matter of debate (Egorov et al., 2003;Waldner et al., 2004;DiCarlo, 2013) due to the occurrence of peculiar pore-scale processes when water infiltrates as fingers (Egorov et al., 2003;DiCarlo, 2013;Katsushima et al., 2013;Baver et al., 2014).Observations reveal that in soils an unstable infiltration profile may be marked by an overshoot profile in terms of θ (saturation overshoot) or ψ (capillary pressure overshoot; see DiCarlo, 2004DiCarlo, , 2007DiCarlo, , 2013;;Baver et al., 2014).Examples of pressure overshoots have been observed in homogeneous snow samples during preferential infiltration by Katsushima et al. (2013).In addition, Hirashima et al. (2014) report promising attempts to reproduce similar dynamics using a 3-D model; they show that solving the Richards equation by including spatial heterogeneity of snow properties and water entry suction (ψ WE ) enables to simulate preferential flow effects.These results suggest that preferential flow in isothermal snow at 0 • C might be explained (and modeled) with a similar approach to the theory of gravity-driven instability of fingers in soils (Katsushima et al., 2009b(Katsushima et al., , 2013)).
Unsaturated hydraulic properties of snow may impede water infiltration in a finer-over-coarser profile.This impedance is usually referred to as a capillary barrier (Khire et al., 2000;Waldner et al., 2004) and is due to the infiltrating water being generally marked by a very high suction when it initially moves in the finer layer.This prevents water from entering the lower layer, thus causing local accumulation of water at the interface (henceforth, simply ponding), a deceleration in the undisturbed advancement of the wetting front, horizontal diversion of water, and a delay in the expected travel time of water.Hill and Parlange (1972), Baker andHillel (1990), Hillel andBaker (1988), Stormont and Morris (1998), Stormont and Anderson (1999), and Khire et al. (2000) discuss this process in soils, whereas Wakahama (1963), Jordan (1995), Waldner et al. (2004), Peitzsch et al. (2008), andMitterer et al. (2011) report some examples for layered snowpack.According to the results by Stormont and Anderson (1999) in soils, water will enter the underlying coarser layer when ψ at the interface decreases to ψ WE ; at this suction, the coarser soil layer firstly becomes conductive (Stormont and Morris, 1998;Khire et al., 2000).A decrease in ψ during ponding is caused by the fact that θ and ψ are related by a hysteretic relation called the water retention curve (WRC) (Daanen and Nieber, 2009;Yamaguchi et al., 2010;Adachi et al., 2012;Yamaguchi et al., 2012).In soils, Hillel and Baker (1988) and Baker and Hillel (1990) note also that, after reaching ψ WE , subsequent flow in the coarser layer will be marked by fingers if, in steady conditions, the hydraulic conductivity of the lower layer at ψ WE is greater than the flux through the top layer q (due to mass conservation).Thus, ponding of water above a capillary barrier is prone to subsequent flow instability, namely, to the development of preferential channels.
Understanding water flow around capillary barriers may be an important step toward efficiently modeling liquid water flow in snow.Furthermore, capillary barriers can play an important role for triggering wet snow avalanches (Mitterer et al., 2011;Wever et al., 2016).For example, Wever et al. (2016) report that predicted local accumulations of water like those expected during ponding at capillary barriers can be used to separate avalanche from non-avalanche days.The position of peak LWC within the snow cover correlates with avalanche size.These processes also affect the timing of snowmelt runoff (Wever et al., 2014), especially during initial infiltration in dry snow.However, their characterization in the literature is still very limited (Eiriksson et al., 2013).Indeed, existing real-time observations in the laboratory or in the field consider a restricted variety of grain size (g S ) combinations (Jordan, 1995;Waldner et al., 2004).Under field conditions, LWC profiles are usually measured using destructive manual methods that strongly limit the temporal and spatial resolution of profiles.Thus, the evaluation of promising results from physically based models (Hirashima et al., 2010;Mitterer et al., 2011;Wever et al., 2014Wever et al., , 2015) ) are often hampered by a lack of a proper high-resolution experimental database.
Here, we collected quantitative information about the liquid water flow around a capillary barrier in snow using laboratory experiments.We considered nine layered snow samples with different grain size combinations and different water input rates.We measured, for each sample, the thickness of the volume of the upper layer affected by ponding of water at the textural boundary, LWC profiles, wet snow fraction at different depths, and the arrival time of water at the sample base.These experiments were performed choosing a quite high vertical resolution of measurements (2 cm) and a broad set of input rates and textures.All the laboratory experiments are compared with numerical simulations of the Richards equation in snow by the 1-D multi-layer physically based snow cover model SNOWPACK, in order to investigate how well 1-D snowpack models are able to capture the behavior of water flow over capillary barriers.
We consider isothermal conditions at 0 • C, thus avoiding any investigation about wetting front advancement in subfreezing snow, which presents additional challenges.Indeed, Marsh and Woo (1984a, b) report that water infiltration in initially subfreezing snow is marked by an alternation of wet snow at 0 • C and dry snow in subfreezing conditions (see also Pfeffer et al., 1990;Pfeffer and Humphrey, 1996).The impact of these processes on runoff response time of snow in sub-freezing conditions is still a matter of debate (Marsh and Woo, 1984b), especially on seasonal timescales (Wever et al., 2014).

Preparation of samples
The main prerequisite to observe capillary barriers in initially dry snow is a finer-over-coarser profile in layering.For this purpose, three combinations of grain size were considered here: (1) FC, i.e., fine-over-coarse snow; (2) FM, i.e., fine-over-medium snow; (3) MC, i.e., medium-over-coarse snow.We classify snow with 0.25 mm ≤ g S ≤ 0.5 mm as fine, snow with 1 mm ≤ g S ≤ 1.4 mm as medium, and snow with 2 mm ≤ g S ≤ 2.8 mm as coarse.Note that this nomenclature is convenient for the presented work, but it is not consistent with the international classification proposed by Fierz et al. (2009), which defines for instance medium snow grain size as 0.5 mm ≤ g S ≤ 1 mm.
Experimental evidence revealed that the area occupied by fingers in snow and the value of ψ WE may be both functions of water input rate (Katsushima et al., 2013).In order for our conclusions to be more general, we carried out experiments with three different water inputs W : these are ∼ 10, 30 and 100 mm h −1 .These water input rates are a compromise between the need for exploring the properties of capillary barriers over a broad range of W , expected melt rates in natural conditions (DeWalle and Rango, 2011), and operational constraints (specifically the expected duration of the tests).Because the saturated conductivity of snow is rather high compared with the chosen input rates, most existing instability criteria (Saffman and Taylor, 1958;Baker and Hillel, 1990;de Rooij, 2000;DiCarlo, 2013) will predict unstable flow in these conditions.Accordingly, Katsushima et al. (2013) have already observed preferential infiltration in snow with average g S between 0.421 and 1.439 mm for different water input rates (21.7 mm h −1 ≤ W ≤ 205.5 mm h −1 ).Nine samples were prepared in a cold room at −20 • C using refrozen melt forms (one sample for each of the three grain size combinations and three water input rates).Henceforth, numbers 1, 2, and 3 differentiate samples with the same grain size combination but subjected to different water input rate (10, 30, and 100 mm h −1 , respectively).Fragmented snow particles were firstly partitioned in several classes of grain size.Afterwards, the three g S chosen were sieved a second time to prepare the samples.Snow was packed in a cylindrical container.The container was composed by a number of acrylic rings (height equal to 20 mm, diameter equal to 50 mm) that were previously taped on the external side.After sieving the lower layer, its dry density (ρ D,L ) was measured by gravimetry.The dry density of the upper layer (ρ D,U ) was measured by gravimetry at the end of sieving operations (by considering the difference between sample total weight and sample weight before sieving the upper layer).After preparation, each sample was moved to a second cold room at 0 • C, where it was stored for at least 12 h to reach initial conditions of dry snow at 0 • C.
We report in Table 1 the details of each experiment.Water input rates are reported both in mm h −1 and in g min −1 (samples diameter equal to 5 cm).The coefficients of variation of ρ D,U and ρ D,L read 0.06 and 0.03.We did not apply any tamping during sieving operations so we had no direct control on the values of ρ D,U and ρ D,L .Given the low variability of these two variables, we point out that this work investigates how capillary barrier effects and associated preferential flow vary with grain size only.Future investigations should focus on the generalization of this work to layers of different density.Some samples (namely, FC2, FM2, and MC2) are shorter than the others.However, the thickness of the upper layer is the same for all the samples.This is important as ponding occurs in the upper layer.

Data collection
Before starting each experiment, we placed a thin cotton ring on the top of the sample to enable the point source of the tracer to spread over the surface of the upper layer.Then, each experiment was started by supplying dyed water into samples using a micro-tube pump.The dye used was blue ink, diluted by a factor of 10 in water at 0 • C. We monitored W during each experiment by automatically measuring the weight of the tracer reservoir (1 min resolution).Absolute relative differences between experimental (Table 1) and reference (10, 30, and 100 mm h −1 ) values of W range between 6 and 19 % as it is difficult to apply a constant, low input rate.
When the tracer reached the base of each sample, tracer supply was stopped.The arrival time of the tracer at sample base (t t ) was registered with a manual chronometer watch by visually inspecting samples during the experiments.Since samples had different heights, we define a specific travel time as with h equal to sample height.Note that τ is in min cm −1 as it is the reciprocal of velocity.After each experiment, pictures of the external sides of the sample were taken to estimate the approximate thickness of the upper layer marked by liquid water accumulation (p, in cm).Soon afterwards, we took pictures of the top section of each acrylic ring (by gradually removing them from the column, snow included).At the same time, the liquid water mass w, in grams, in each of the rings was measured using a portable calorimeter (Kawashima et al., 1998).These measurements were translated into profiles of volumetric liquid water content by converting w to θ .Fractions of wet areas over total area (f ) were also estimated for each section by manually delimiting fingers in all the pictures taken and calculating their extension using the ImageJ software (http://imagej.nih.gov/ij/,v. 1.48; see Abramoff et al., 2004).

The comparison with SNOWPACK
We simulated the dynamics of each sample using the physically based 1-D multi-layer snow cover model SNOWPACK (Bartelt and Lehning, 2002).These simulations aim at comparing observations of capillary barrier development with predictions by a physically based model, as previously done by, e.g., Hirashima et al. (2010), Mitterer et al. (2011), andWever et al. (2015) mainly using field observations.The relatively high resolution of LWC measurements (2 cm) enables a rather detailed discussion of both the physical process and its simulation by a physically based model.This comparison will not include preferential flow patterns and arrival times, as the model does not include an explicit treatment of preferential flow regimes.The model discretizes snow using a finite element grid.It simulates the evolution in time of a broad set of variables along a vertical profile of snow starting from external forcings.The original version of SNOWPACK considers a bucket-type approach to simulate water percolation in snow.Accordingly, water is retained at a given position in the profile until it exceeds a threshold (see Coleou and Lesaffre, 1998).After exceeding, excess water is transmitted downwards.Hirashima et al. (2010) have introduced in SNOW-PACK a water transport model based on the model by van Genuchten (1980) and on an equilibrium approximation of water flow to tackle numerical instability (see Hirashima et al. (2010) for details).Recently, Wever et al. (2014) have also introduced a discretization of the Richards equation that significantly improves several aspects of liquid water content simulation in snow.We used the numerical scheme by Wever et al. (2014) in this paper (WE, SNOWPACK version 3.3, https://models.slf.ch/).
The initial spatial resolution of simulations was set to 2 cm.The time step was set to 1 min, but the numerical scheme by Wever et al. (2014) reduces this initial time step basing on an iteration rule (see Wever et al. (2014) for details).Snow initial conditions were chosen to replicate the grain type, size, and density (Table 1), θ (initially dry), and temperature (0 • C) of the physical samples.Using the same sieves that we used here, Katsushima et al. (2013) obtained a median grain size (hereinafter, ḡS ) for the class 0.25-0.5 mm and the class 1-1.4 mm equal to 0.406 and 1.463 mm, respectively.ḡS for medium snow is greater than the upper boundary of the sieve probably because snow grains used were not perfectly spherical.In the simulation, we therefore set ḡS = 0.406 mm for fine snow, ḡS = 1.463 mm for medium snow, and ḡS = 2.926 mm for coarse snow (by assuming this last value as 2 times the median medium grain size).Bond size was assumed equal to one-third of grain radius.Input data were chosen to replicate experimental conditions in the cold chamber, i.e., a constant precipitation flux (equal to the measured water input flux W , see Table 1) and a fixed air temperature of 0 • C. The threshold temperature for classifying solid and liquid events was set to −0.01 • C, in order for W to be classified as liquid.Wind speed and solar radiation were set to 0, while incoming longwave radiation was calculated as σ T 4 , where σ = 5.67 × 10 −8 W m −2 K −4 and T = 273.15K. Parametrizations for snow permeability, unsaturated hydraulic conductivity, water retention curve, residual water content, and averaging method for hydraulic conductivity at the interface were all kept at the default settings discussed for SNOWPACK in Wever et al. (2014Wever et al. ( , 2015)).
Particular attention is paid to the comparison between observed and simulated LWC peak over the interface between layers, as this is an important variable involved in capillary barrier formation and in wet snow avalanche forecasting (Wever et al., 2016).Another key feature of capillary barriers is the vertical profile in LWC (Hirashima et al., 2010).However, choosing a single snapshot of simulated LWC for the comparison with our observations is problematic, as the model is 1-D and, at this stage, does not include an explicit treatment of preferential flow patterns.These are expected to play a key role in water flow around capillary barriers as water concentrates in fingers that are characterized by a higher-than-average unsaturated hydraulic conductivity (due to a higher-than-average LWC).It follows that restricting this comparison to a single profile (i.e., only one time step) may be misleading as possible differences between observations and simulations might be due to a process that is currently not treated by the model.This may limit a comparison aiming at assessing the capability of a model to reproduce capillary barriers.Thus, we will compare observed profiles of LWC with model results at two different times.The first one (WE1) is the observed arrival time of water at sample base; the second one (WE2) is the simulated arrival time of water at sample base, which is chosen by identifying the instant when θ at sample base reaches ∼ 3 vol % in the simulation (Mitterer et al., 2013).Note that WE1 and WE2 for each sample were obtained from the same simulation.
On the one hand, WE1 is advantageous as supplied mass in experiments and in simulations matches because both profiles refer to the same time step.On the other hand, flow in the lower layer will be at the beginning spatially variable, strongly accelerated, and highly fingered (Katsushima et al., 2013), which are all features that are not explicitly included in the 1-D model and that might hamper the application of the Richards equation.Thus, when considering WE1, we will focus on the profile over the interface, where the peak of LWC develops.Conversely, WE2 enables a comparison of a full profile of LWC, but the simulated mass of liquid water will be greater than observed due to the possible mismatch between observed and expected arrival time of water in simulations (Wever et al., 2014).This is particularly evident in the upper layer of FC and FM samples, as water speed in fine snow during matrix flow is slow.Thus, when considering WE2, we will focus on the profile in the lower layer.Because liquid water flow in MC samples turned out to be highly fingered (see next section), observations in these samples are compared only with WE2 (full profile), which is probably less affected by effects due to preferential flow, including the time arrival of water at a certain point of the profile.Note that both model version and the evaluation methodology are different from the preliminary results reported in the discussion paper (Avanzi et al., 2015a).3 Results and discussion

Overview
Figure 1 reports the horizontal sections of all the samples (2 cm vertical resolution) at the end of the experiments (i.e., when dyed water arrived at sample base).Dyed water is visible as blue stains.Generally, the darker the color is, the greater local LWC is (Waldner et al., 2004).We report in Fig. 2 three examples of samples at the end of the experiment.These are FC2 (as an example of FC tests), FM2 (as an example of FM experiments), and MC1 (as an example of MC experiments).Table 2 reports observations in terms of thickness of the upper layer marked by liquid water accumulation (p), arrival time of water at the base of each sample (t t ), and specific travel time of water in snow (τ ).In Figs. 3 and 4, profiles of wet snow fractions f and LWC are given.In Fig. 4, each point represents bulk LWC in the underlying 2 cm.As an example, any value reported at a depth equal to 8 cm is bulk LWC between 8 and 10 cm.This represents the LWC measured immediately over the interface between layers.
Figure 5 compares observed and SNOWPACK-based profiles of volumetric LWC for each sample.Each point represents bulk LWC in the underlying 2 cm.As described in Sect.2.3, two simulated profiles are reported for FC and FM samples (WE1 and WE2), whereas only WE2 is reported for MC samples.Dashed lines indicate the parts of the simulated profiles that are not considered in the discussion (see again Sect.2.3).transition in wetness between finer and coarser layers (see Fig. 2).In four out of six samples of these two classes, a homogeneously blue area is observed even at a depth equal to 8 cm (i.e., 2 cm above the interface).MC samples show a more variable behavior.Indeed, almost no water spreading was observed for MC1 (see also Fig. 2), whereas marked horizontal redistribution of water is visible in MC2 and MC3.MC samples show a smaller p than FC and FM samples.The difference between FC-FM and MC samples in terms of ponding behavior may be explained considering the retention properties of snow with different grain size.ψ WE for medium and coarse snow can be estimated from grain size using the relation reported in Katsushima et al. (2013) and Hirashima et al. (2014): ψ WE ∼ 0.025 m in coarse snow and ∼ 0.04 m in medium snow.In contrast, ψ in fine snow for 5 and 10 % LWC is ∼ 0.22 and 0.21 m, respectively, and ψ in medium snow for 5 and 10 % LWC is ∼ 0.09 and 0.08 m, respectively.This implies that for typical low saturation values in snow, the difference between the suction pressure in the finer snow and the water entry pressure of the coarser snow is larger for fine snow.Furthermore, unsaturated conductivity of coarser snow is likely to quickly decrease with increasing ψ, as already observed in soils (Stormont and Anderson, 1999).It follows that, at the same W , the greater the mismatch between unsaturated properties of finer and coarser snow is, the larger the mass of water accumulated and horizontally diverted until the underlying layer is able to convey the supplied flux.These approximate values of suction in snow were estimated using the WRC parametrization in snow by Yamaguchi et al. (2012), assuming a residual LWC equal to 2.4 vol % (Yamaguchi et al., 2010;Hirashima et al., 2014) and considering 5 and 10 vol % as reference values for relatively low saturation.Note that the WRC by Yamaguchi et al. (2012) refers to a drying process and this may cause some additional uncertainty when estimating ψ for a wetting process, due to hysteresis.Indeed, ψ for a primary wetting process is expected to be smaller than ψ for a primary drying process.Available data of hysteresis in snow are, however, very preliminary (Adachi et al., 2012).A specific discussion  about the role of hysteresis for interpreting these results is given in Sect.3.4.

Development of capillary barriers
All layering types are characterized by similar LWC profiles with the only difference in absolute values for LWC (Fig. 4).LWC increases with depth in the upper layer, it presents a marked peak at the textural boundary, and it decreases again below the capillary barrier.Peaks in LWC at the interface may be associated with capillary barriers as water ponds until ψ reaches ψ WE and the underlying layer becomes conductive (Stormont and Anderson, 1999): similar examples are reported or discussed in, e.g., Waldner et al. (2004)  (2015b), and Wever et al. (2015Wever et al. ( , 2016)).All FC-FM samples yield a similar LWC in the upper layer at the interface: ∼ 33 vol % in FC samples and ∼ 34-36 vol % in FM samples.Again, this may be explained by considering that the peak of LWC at the boundary is ruled by the retention properties of snow with different grain size and by their contrast (see Khire et al. (2000) for a similar discussion in soils).Note that such values of LWC are much greater than those usually reported in field profiles (Fierz et al., 2009).LWC drives, among others, snow settling (Marshall et al., 1999) and wet snow metamorphism (Brun, 1989).Both processes experience a dramatic acceleration with increasing LWC.This supports the idea that capillary barriers may play an essential role for snow stability in wet conditions (Wever et al., 2016).
In MC samples, a smaller, but distinct, peak in LWC was measured at the interface: in MC1, LWC over the boundary is ∼ 4.5 vol %, whereas LWC values immediately above and below are 2.7 and 1.7 vol %, respectively.In MC2 and MC3, the peak in LWC over the interface is ∼ 9 vol %.A smaller peak of LWC in MC layering is attributed to the small difference between unsaturated hydraulic properties in medium and coarse snow (for example, a smaller difference between ψ in medium snow and ψ WE in coarse snow).Note again that this difference could be even smaller than expected if hysteresis were explicitly taken into account (Adachi et al., 2012).The observed peaks in MC samples are generally greater than the peak value observed by Waldner et al. (2004) during snowmelt infiltration through a 1.5 over 2.5 mm transition in an artificially sieved snowpack.This difference may be due to (microstructural) heterogeneity in snow, infiltration rate, experiment durations, and the larger measurement area of the TDR system used by Waldner et al. (2004) compared to a calorimeter (see Sects. 3.4 and 3.5).
The occurrence of capillary barrier causes horizontal redistribution of water.Thus, spatial homogenization of liquid water patterns at the interface is promoted.Indeed, f increases with depth over the boundary (where f = 1 for all FC and FM samples; see Fig. 3).However, sections in Fig. 1 reveal a remarkable spatial variability of this process at centimeter scale.For example, some pockets of dry snow persist at depths equal to 8 cm in samples FC2 and FC3 (i.e., 2 cm above the interface).Other indications are the observed spatial variability of p in each sample (Table 2) and the differences of coloring in some sections at depths equal to 8 or 10 cm (e.g., FC2, FC3, and FM3), which can be linked to differences in LWC (Waldner et al., 2004).Isolated clusters of liquid water surrounded by dry snow are also visible in MC1 and MC2.All these observations show that the distribution of liquid water above a capillary barrier has a marked 2-D (or even 3-D) structure at local scale, probably due to heterogeneity in snow microstructure.The difference in wet areas in MC samples for different water input rates might be on the contrary an effect of water input rate on ψ WE , as observed in snow by Katsushima et al. (2013).

Preferential flow patterns and travel time of water in snow
Observed profiles of f suggest that water movement in samples was marked by high spatial variability and that this variability is lower in fine snow layers than in medium or coarse snow.Overall, preferential flow turns out as the predominant pattern of water infiltration in snow (Schneebeli, 1995).We observed that new fingers created during the percolation, and that sometimes fingers stopped their vertical percolation at some locations but continued to develop at others.The movement of fingers in the lower layer was very rapid and represented a small fraction of the total duration of each experiment (typically in the range of minutes).Katsushima et al. (2013) report that the total area of preferential flow (hence f ) in snow samples made by vertically homogeneous snow decreases with increasing grain size but increases with increasing input flux.We also observed a decrease in f with increasing g S , whereas a clear increase of f with increasing W was detected only for MC samples.On the one hand, the expected dependency of f on sublayer unsaturated conductivity (Hillel and Baker, 1988;Baker and Hillel, 1990) and the possible relation between ψ WE and velocity (DiCarlo, 2007(DiCarlo, , 2013;;Katsushima et al., 2013) support the existence of a relation between f and W (albeit both effects have been mainly observed only in soils).On the other hand, these experiments included a capillary barrier, contrary to experiments by Katsushima et al. (2013), and this represents a major driver of liquid water content patterns at a more local scale (see the previous Section).Accordingly, water speed in the upper layer was locally affected by ponding, whereas inflow rate in the sublayer was driven by breakthrough of water when reaching ψ WE .Both processes limit the impact of external water input rate on f , at least until steady conditions are reached.In MC samples, the difference in retention prop-erties between layers is lower; thus, the effect of a capillary barrier is spatially very localized.This may explain why observations in MC samples agree with previous observations in homogeneous snow.
Considering outcomes of different experiments in sand, DiCarlo (2013) reports that finger width might increase with both very high (i.e., close to saturated conductivity) and very low supplied fluxes, while finger width keeps constant for a broad range of supplied flux in between.Water input rates in our experiments span 11 and 113 mm h −1 ; these values are very small when compared with expected values of saturated conductivity in snow (∼ 10 4 -10 6 mm h −1 ; see Katsushima et al., 2013).We therefore suggest that future developments of this work should investigate the relation between flux and f extensively, i.e., enlarging the range of W considered during the experiments and/or reaching steady conditions.Furthermore, additional experiments should be carried out using containers of different size in order to assess whether the experimental geometry used may induce possible boundary effects (see also Katsushima et al. (2013) on this).
τ increases with decreasing W , as clearly expected (Table 2).In the case of FM2 (fine-over-medium snow, W = 27.7 mm h −1 ), we can compare the τ measured during the experiment (2.2 min cm −1 ) with the τ observed during the experiment by Katsushima et al. (2013), since this is the only g S − W combination that these two works share.The τ measured by Katsushima et al. (2013) for fine snow and W = 22.3 mm h −1 is equal to 1.7 min cm −1 , while the τ for medium snow and W = 21.7 mm h −1 is equal to 0.7 min cm −1 .These results suggest that τ in a FM sample is higher than the τ observed in a homogeneous sample composed by medium snow.This is clearly expected since permeability in medium snow is higher than in fine snow.However, this τ is even higher than the specific travel time observed in a homogeneous sample made by fine snow, which is marked by a low saturated conductivity.This comparison helps to quantify the relevance of capillary effects in ruling water speed in snow and the arrival time of meltwater at snow base.

The comparison with SNOWPACK
According to Fig. 5, SNOWPACK clearly reproduces an increasing LWC with depth in the upper layer and a peak of LWC at the interface at WE1.Furthermore, LWC profiles below the barrier are generally in good agreement, once water has reached the base in the simulation (WE2).Point differences between observed and simulated LWC at the interface at WE1 read ∼ 2-5 vol % in FC1-FC2, ∼ 3-8 vol % in FM1-FM2, and ∼ 0.1-5 vol % in MC samples.A larger difference (∼ 9-13 vol %) is found for higher input rates.However, note that FC3 and FM3 were subjected to an extremely high water input rate compared with natural conditions.
Previous evaluations of SNOWPACK already show that the inclusion of the Darcy-Buckingham equation in snow en- The Cryosphere, 10, 2013-2026, 2016 ables a correct prediction of the onset of capillary barriers at textural discontinuities (Hirashima et al., 2010;Mitterer et al., 2011;Wever et al., 2015).Here, we enlarged previous findings by considering a broad set of snow textures and in-put rates and a relatively high resolution of measurements.Note that the version of SNOWPACK used does not implement a parametrization of water entry suction (Wever et al., 2014), but the model is anyway able to provide a sufficiently good performance in reproducing the profile around a capillary barrier.Results by Stormont and Morris (1998), Stormont and Anderson (1999), Khire et al. (2000), andDiCarlo (2007) show that ψ and θ at an infiltrating front (hence, ψ WE ) may follow a wetting WRC.Both variables are also strictly coupled with unsaturated conductivity (Mualem, 1976;Stormont and Anderson, 1999) and the impedance mismatch given by the low unsaturated conductivity of the coarser layer compared to the applied flux plays a key role in delaying water on the barrier.SNOWPACK currently includes a parametrization of both a WRC and unsaturated conductivity and solves the Richards equation.It may be that implementing the Richards equation in 1-D is sufficient to mimic some essential features of capillary barriers in snow (e.g., ponding) even without an explicit calculation of ψ WE .Note that an increase in LWC with depth as well as abrupt transitions in LWC at the interface may occur even in equilibrium, i.e., when suction increases with height.This is due to the different retention properties of fine, medium, and coarse snow.Furthermore, suction profiles over the barrier may depend on applied flux too (Stormont and Morris, 1998;Stormont and Anderson, 1999).Thus, a more detailed analysis of capillary barrier dynamics in snow necessarily needs observations of suction profiles during infiltration; this represents an important step of future research.
Another important limitation for this discussion may be the present lack of an exhaustive investigation of WRC hysteresis in snow (Adachi et al., 2012).As already noted, the absence of a proper parametrization of hysteresis may hamper the estimation of expected LWC at ψ WE during ponding (i.e., wetting), if different WRCs for a wetting and a drying process are not known.Also, the hysteretic behavior of the WRC is considered an important factor in driving preferential flow in general, since for example it promotes the persistence of fingers in soils (Liu et al., 1994;DiCarlo, 2013).The magnitude of hysteresis is also related with the magnitude of capillary and saturation overshoot (DiCarlo, 2007;Katsushima et al., 2013), although it is not the prime cause of instability (DiCarlo, 2013).In this context, note that, according to Khire et al. (2000), hysteresis plays a less important role than the difference in unsaturated hydraulic properties between the finer and the coarser layer when studying the general properties of capillary barriers in soils and how they depend on layer parameters; this may again support the idea that the existing implementation of unsaturated flow in a complex 1-D model may be sufficient to mimic LWC distribution around a capillary barrier.A possible improvement may be represented by the set of parametric models proposed by Luckner et al. (1989) for porous media, which includes hysteresis.

F. Avanzi et al.: Capillary barriers in snow
Observations show that both breakthrough of liquid water below a capillary barrier and wet conditions in the upper layer or in fingers may present a high spatial variability at centimeter scale.This is because natural snow is spatially heterogeneous (Hirashima et al., 2014) and this may affect 3-D patterns of capillary barriers (e.g., see the already discussed pockets of dry snow in FM3).Alternation of dry and wet snow can sensibly decrease the bulk LWC in a ring, although local LWC can still be very high.This may partially explain some differences between observations and 1-D simulations.For example, predicted peak LWC in FC3 and FM3 is ∼ 43-46 vol %, which is close to saturated conditions (the porosity of fine snow in both samples is ∼ 0.5) but greater than observations.An approximate estimation of LWC at ψ WE in fine snow (obtained assuming continuity of suction at the interface) reads 50 vol %, which is closer to SNOWPACK simulations than data.Thus, saturated conditions might be reached at a very local scale, while bulk LWC in each ring can be lower due to heterogeneity in wetness at a larger scale.Another example is water flow below the interface, which showed a high degree of spatial variability.The good agreement between the model and the data (at WE2) might suggest that at this measurement resolution differences in LWC between a highly channeled flow and a matrix-only simulation balance; that is, fingers are usually highly saturated (Waldner et al., 2004) but occupy only a small fraction of total volume.Thus, the average LWC at ring scale is much lower than saturation and close to matrix conditions.This result suggests that an exhaustive process understanding of the physics of capillary barriers in snow and water flow instability may need that a proper measurement and/or modeling scale are established to clearly separate modeldata significant discrepancies and effects due to the sampling strategy (see Blöschl (1999) for a definition).Importantly, the spatial resolution needed to capture 3-D patterns of capillary barriers might be smaller than that usually used to sample LWC in the field (see again Fig. 1).Increasing the spatial resolution of LWC measurements is challenging as measuring LWC in snow is still marked by high uncertainties (Colbeck, 1978;Fierz and Föhn, 1995;Techel and Pielmeier, 2011;Avanzi et al., 2014).It is only recently that undisturbed, non-destructive, and repetitive measurements of LWC have been obtained (Heilig et al., 2010(Heilig et al., , 2015;;Schmid et al., 2015;Kinar and Pomeroy, 2015).A promising alternative might be given by pore-scale measurements of liquid water flow (Adachi et al., 2012;Walter et al., 2013).
This discussion also reveals the role played by heterogeneity (Hirashima et al., 2014) in introducing possible differences in LWC between 3-D (bulk) and 1-D conditions.Additional uncertainty in this comparison may be caused by instrumental precision (see next Section), ambiguity in the identification of the correct snapshot of LWC for this comparison, possible air trapped in voids at saturation (Yamaguchi et al., 2010), and possible boundary effects due to the experimental geometry.To summarize, we suggest that addi-tional investigations should be carried out to establish proper frameworks for the high-resolution comparison of complex models and laboratory (or field) observations.

The role of instrumental precision
A mass balance between supplied and measured liquid water mass reveals that the measured mass ranges between 93 and 176 % of supplied mass in eight out of nine samples, while in MC1 measured mass is 434 % of supplied mass.Note that in this last sample the total mass supplied is nonetheless very small due to the short duration of this experiment (∼ 2.88 g).
This discrepancy can be explained by instrumental noise.Melting calorimetry has been widely used to measure LWC for decades (Yosida, 1960), but Colbeck (1978) points out that this method may be inaccurate as it implies the calculation of a difference between large numbers (Stein et al., 1997).According to Kinar and Pomeroy (2015), absolute errors in measuring LWC using calorimetry span 1 and 5 %.The instrument we used here (the so-called Endo-type snowwater content meter) was proposed by Kawashima et al. (1998).They note that measured LWC span ±2 % of known LWC (by weight) in 87 % of the cases.By comparing measurements by the Endo-type calorimeter with those by a dielectric device in snow pits (see Kawashima et al. (1998) for details), they note that this device returns alternatively higher or lower LWC if compared with high and low readings by the dielectric device.
We estimated an absolute error for these experiments by comparing the height-integrated LWC measured using calorimetry within each sample with the ratio between supplied liquid water volume and total volume of samples.The absolute difference spans 0.8 and 2.97 vol %, thus it is consistent with the literature (Kinar and Pomeroy, 2015).Measured and simulated LWC by SNOWPACK are also in fair agreement (see previous section) which underlines the above mentioned range of absolute error, since SNOWPACK bases an energy conservation on mass.
Capillary barriers and associated preferential flow represent a large challenge for LWC measurements.On the one hand, peaks in LWC at the interface are rather high and this is a problem for those instruments that may lose accuracy for high LWC, such as a snow fork (Techel and Pielmeier, 2011).On the other hand, bulk LWC in fingered snow may be very low, as water accelerates and occupies a small fraction of total volume.This means that such experiments need an instrument that guarantees a comparable performance for both high and low LWC.This may represent a benefit of the Endo calorimeter (see Fig. 4 in Kawashima et al., 1998), which seems also appropriate given the small dimension of each ring.Furthermore, Fierz and Föhn (1995) report that the absolute error in measuring water content using dielectric methods spans 0.2 and 0.9 vol %, while Techel and Pielmeier (2011) note that the expected difference between measurements taken using a Denoth meter and a snow fork is ∼ 1 vol %.Thus, measuring low LWC is generally very challenging for several existing techniques.Finally, a highly fingered flow may be missed and/or disturbed by using larger instruments (Stein et al., 1997).For future experiments, we will also consider alternative portable techniques, such as a dilution method (Davis et al., 1985;Kinar and Pomeroy, 2015;Mitterer, 2016).

Conclusions
We focused on the systematic observation of capillary barriers and associated preferential flow during laboratory experiments in a cold chamber.We sieved nine samples of finerover-coarser snow.These samples were subjected to controlled supply of dyed water until water arrived at sample base.Liquid water patterns in stratified snow were characterized using visual inspection, LWC measurements, and horizontal sectioning.Results were also compared with SNOW-PACK simulations.
Overall, results confirmed that a finer-over-coarser transition in snow layering causes ponding of water when it arrives at the textural boundary.Measured peaks in LWC over the boundary are large with respect to usual measurements in the field (up to ∼ 33 vol % in fine-over-coarse samples and ∼ 34-36 vol % in fine-over-medium samples), while peaks in medium-over-coarse samples are usually ≤ 10 vol %.Differences in peak LWC between samples were explained by varying unsaturated hydraulic properties of snow with different grain sizes.A more detailed analysis of horizontal sections revealed marked variability of wetness conditions at centimeter scale, thus suggesting that local LWC might even be greater than measured.
Horizontal sectioning of samples confirmed that preferential flow seems the dominant process in water transmission in snow.The area occupied by fingers (f ) increases with grain size, while no definitive result was obtained to establish a relation between f and water input rate.This is explained by the strong perturbation introduced by the capillary barrier in liquid water content patterns when compared with previous observations in homogeneous snow.
The comparison with SNOWPACK showed that, in general terms, the implementation of the Richards equation clearly reproduces the existence of a capillary barrier and yields a good agreement with observed peaks in LWC at the interface.The marked spatial variability of liquid water content in snow represents a source of uncertainty when comparing measurements at a relatively high resolution with a 1-D model.Future steps of this work will compare these measurements with a 3-D simulation of liquid water infiltration in snow (Hirashima et al., 2014).

Data availability
The data sets plotted in Fig. 3 and 4 are available upon request to the authors.

Figure 1 Figure 1 .
Figure1confirms that liquid water movement through a finer-over-coarser snow texture is subjected to ponding and horizontal diversion of water when the wetting front comes to the textural interface.In FC and FM samples, horizontal spreading of water at the interface introduces a clear textural

Figure 2 .
Figure2.Three samples at the end of the experiments: FC2 (on the left, as an example of FC samples), FM2 (at the center, as an example of FM samples), and MC1 (on the right, as an example of MC samples).

Figure 3 .
Figure 3. Measured f profiles.f is the ratio between wet and total area for all the sections in Fig. 1.(a) FC samples; (b) FM samples; (c) MC samples.The vertical coordinate refers to the depth of the section from sample top surface.Serial numbering 1-3 represent the three different input rates.

Figure 4 .
Figure 4. Measured LWC (vol %).(a) FC samples; (b) FM samples; (c) MC samples.Each point represents bulk LWC in the underlying 2 cm.This convention is consistent with Figs. 1 and 3. Serial numbering 1-3 represent the three different input rates.

Figure 5 .
Figure 5.Comparison between observed and simulated profiles of volumetric LWC.(a), (b), and (c) refer to samples FC1, FC2, and FC3.(d), (e), and (f) refer to samples FM1, FM2, and FM3.(g), (h), and (i) refer to samples MC1, MC2, and MC3.Note that panels (g) and (h) have a different horizontal range from the others.Each point represents bulk LWC in the underlying 2 cm.As described in Sect.2.3, two simulated profiles are reported for FC and FM samples (WE1 and WE2), whereas only WE2 is reported for MC samples.Dashed lines indicate the parts of the profiles that are not considered in the discussion (see again Sect.2.3).

Table 2 .
Experimental results: observed ponding layer thickness p, experiment duration t t , and specific travel time τ .As for p, approximated lower and upper values are reported due to spatial heterogeneity in this variable.