Nomenclature
Introduction
Snow, a sintered porous material made of ice grains, has a complex porous microstructure that continuously changes with time and external conditions. Its effective mass transport properties, strongly dependent on the complex microstructure, are relevant for investigating a wide range of environmental processes.
Permeability has a direct effect on snow–air exchange processes with an impact on atmosphere chemistry (Reference GrannasGrannas and others, 2007; Reference Clifton, Manes, Rüedi, Guala and LehningClifton and others, 2008), on snow metamorphism (Reference Albert, Shuman, Courville, Bauer, Fahnestock and ScambosAlbert and others, 2004) and on water flow through snow (Reference Waldner, Schneebeli, Schultze-Zimmermann and FlühlerWaldner and others, 2004). Reference Bader and BaderBader (1939) gave the first quantitative data of snow permeability. Reference ShimizuShimizu (1970) and Reference Sommerfeld and RocchioSommerfeld and Rocchio (1993) parameterized it in relation to density and grain size, and a permeameter for field measurements was developed by Reference Conway and AbrahamsonConway and Abrahamson (1984). The experimental characterization of the effective transport properties is difficult, in part due to the rapid change of the snow’s microstructure with temperature and temperature gradient (Reference Albert and ShultzAlbert and Schultz, 2002). Lately, permeability and specific surface area were used to develop a new textural characterization of snow (Reference Arakawa, Izumi, Kawashima and KawamuraArakawa and others, 2009). The Dupuit–Forchheimer coefficient is taken into account at a higher Reynolds number, when inertial effects become important (Reference KavianyKaviany, 1995).
Theoretical and empirical correlations for the determination of permeability have been developed for simplified two-phase media such as capillary drag and the Carman–Kozeny models (Reference ErgunErgun, 1952; Reference DullienDullien, 1979; Reference Macdonald, El-Sayed, Mow and DullienMacdonald and others, 1979; Reference KavianyKaviany, 1995) for fibrous beds (Reference DaviesDavies, 1952; Reference ChenChen, 1955) and cellular foams (Reference Moreira, Innocentini and CouryMoreira and others, 2004). No previous studies on the Dupuit–Forchheimer coefficient of snow were found; correlations were proposed for other porous materials (Reference DullienDullien, 1979; Reference KavianyKaviany, 1995).
Direct pore-level simulations (DPLS) have become a powerful tool for the characterization of a wide range of porous materials. In previous studies (Reference Fredrich, DiGiovanni and NobleFredrich and others, 2006; Reference Petrasch, Meier, Friess and SteinfeldPetrasch and others, 2008; Reference Haussener, Lipiński, Petrasch, Wyss and SteinfeldHaussener and others, 2009, Reference Haussener, Coray, Lipiński, Wyss and Steinfeld2010), micro-computed tomography (μCT) was applied to obtain the precise digital three-dimensional (3-D) geometrical representation of complex porous media, such as reticulate ceramic foams, porous rocks and packed beds of opaque or semi-transparent particles, and subsequently used in DPLS to calculate the effective transport properties. Recently, DPLS has been applied for the characterization of polar firn (Reference Courville, Hörhold, Hopkins and AlbertCourville and others, 2010) and shown to describe with good accuracy morphological properties as supported by experimental validation. A fundamental advantage of DPLS compared with direct measurements is that stratigraphically complex snow samples with thin layers can also be characterized, leading to different properties for each layer, whereas an experimental measurement would yield only an average over the whole sample. In the present paper, μCT is applied to obtain the 3-D digital geometry of seasonal snow types. The governing mass and momentum conservation equations are numerically solved at the pore scale (DPLS) by the finite-volume method, allowing for the determination of the permeability and Dupuit–Forchheimer coefficient.
Morphological Characterization
Five different snow samples are considered: decomposing snow (ds), metamorphosed I (mI), metamorphosed II (mII), depth hoar (dh) and wet snow (ws). They correspond to the grain shape classifications: DFdc, RGsr/DFdc, RGsr, DHcp and MFcl (International Classification for Seasonal Snow on the Ground (ICSSG; Reference FierzFierz and others, 2009)). Tomographic scans were carried out with a Scanco μCT 80 desktop X-ray tomographic set-up (Reference Kerbrat, Pinzer, Huthwelker, Gäggeler, Ammann and SchneebeliKerbrat and others, 2008). The voxel sizes were 10 μm for the ds, mI and mII snow samples, and 18 μm for the dh and ws samples, with a scanned volume of 600 × 600 × 400 voxels, corresponding to 144 mm3 and 840 mm3, respectively. As an example, Figure 1 depicts the 3-D surface rendering of the snow sample (ws) with fluid flow streamlines.
Table 1 summarizes the morphological characteristics by experimental methods (Reference Kerbrat, Pinzer, Huthwelker, Gäggeler, Ammann and SchneebeliKerbrat and others, 2008) and by computation of the two-point correlation function and opening size distribution with spherical structuring elements on the μCT scans (Reference HaussenerHaussener, 2010). The pore and particle sizes must be read with care as they describe the smallest dimension of the pore and particle spaces, which might not characterize well the size of complex, non-spherical pores and particles.
The representative elementary volume (REV), i.e. the smallest cubic volume that can be considered as a continuum, is determined by calculating a continuum property of the sample on subsequently growing volumes until it reaches a constant value within a band of ±ξ, with ξ < 1. Figure 2 shows an example of the convergence value of the length of cubic REV calculated based on porosity for the wet snow sample with ξ = 0.05. The lengths of cubic REV based on porosity for the five snow samples are listed in Table 1 and are 2.4–13.8 times larger than the calculated pore and particle diameter, respectively. dh needs a relatively small REV while ws needs a relatively large REV compared with their characteristic lengths.
REV calculated based on pressure drop – and consequently on permeability and Dupuit–Forchheimer coefficient – is also shown in Figure 2, which indicates that larger REVs are required. The need for larger REV based on heat transfer properties was previously discussed (Reference Haussener, Coray, Lipiński, Wyss and SteinfeldHaussener and others, 2010).
Methodology
The pressure drop over a spatially averaged isotropic porous medium is given by the extended Darcy’s law (Reference Petrasch, Meier, Friess and SteinfeldPetrasch and others, 2008; Reference HaussenerHaussener, 2010):
where K is the permeability, F DF is the Dupuit–Forchheimer coefficient, ρ is the fluid density, μ is the dynamic viscosity of the fluid and u D is its superficial velocity, , with volume V larger than or equal to REV. The first term is the result of viscous effects, predominant at low Reynolds numbers, whereas the second term describes the inertial effects, which become important at higher fluid velocities (Re>1) (Reference Petrasch, Meier, Friess and SteinfeldPetrasch and others, 2008). Non-dimensionalization of Equation (1) for the one-dimensional case yields:
where d is a characteristics length scale, p is the pressure, c 0 and c 1 are constants and the normalized pressure drop IIpg is linearly dependent on Reynolds number Re. The characteristic length scale, d, used throughout this study is the numerically calculated pore diameter (see Table 1). A linear least-square-fitting method was used to fit the numerically calculated Re-dependent IIpg, allowing for the determination of K and F DF. DPLS of fluid flow across the five characteristic snow samples was performed. An in-house tetrahedron-based mesh generator was used to create the computational grid directly on the μCT scans. A commercial CFD code (ANSYS, 2009) based on the finite volume technique was used to solve the continuity and Navier–Stokes equations. The computational domain, shown in Figure 3, consists of a square duct containing a sample of the porous material. The boundary conditions are: uniform inlet velocity and temperature and outlet pressure, no-slip and constant wall temperature at the solid–fluid interface, and symmetry at the lateral duct walls.
Preliminary calculations were carried out for various sample sizes and mesh element sizes of the subset to elucidate the trade-off between computational time and accuracy. A representative sample size of 600 × 600 × 200 voxels (10.8 × 10.8 × 3.6 mm3) with a largest mesh element size of 225 μm was chosen for ws and dh samples. For ds and mI samples, the chosen sample size was 600 × 600 × 300 voxels (6 × 6 × 3 mm3) with a largest mesh element size of 125 μm. For the mII sample, the sample size was 600 × 600 × 200 voxels (6 × 6 × 2 mm3) with a largest mesh element size of 125 μm. Convergence was achieved for a termination residual root mean square (RMS) of the iterative solution below 9 × 10−5. The sample sizes were chosen with a relative difference in which the highest possible size was 1.7–6.2%, whereas the largest mesh element sizes were chosen with a relative difference in which the smallest possible mesh element was 2.8–10.3%. These differences were obtained by calculating the pressure drop for each mesh and sample size.
Results and Discussion
The dimensionless pressure gradient IIpg is plotted as a function of Re in Figure 4 for the five snow samples. The calculated permeability K and Dupuit–Forchheimer coefficients F DF are plotted in Figure 5 versus the pore diameter; their values and the goodness of fit are listed in Table 2.
K is lowest and F DF highest for the ws sample, as its density is highest and porosity lowest. On the other hand, K increases with ds, mI, mII and dh samples because of the increasing pore size, which reduces pressure loss and leads to a higher K and a smaller F DF (Reference Haussener, Coray, Lipiński, Wyss and SteinfeldHaussener and others, 2010). The unexpected decrease of ws in K and increase in F DF highlights that one morphological characteristic (e.g. d pore) does not describe sufficiently well the microstructure and supports the importance of the CT-based determination of the effective mass transport properties.
The values of K and F DF of Table 2 are compared with theoretical and empirical models using simplified microstructure. The models are listed in Table 3.
The K and F DF of the five characteristic snow samples, calculated by the CT-based DPLS method (Table 2) and by the simplified models of Equations (3–10), are shown in Figures 6 and 7, respectively. The conduit flow model, Equation (3), compares well with DPLS for mI, mII, dh and ds, particularly for dh (relative difference from DPLS of 17%). DPLS gives results close to those of the fibrous bed model for all types of snow (relative difference from DPLS from 18% for ws to 77% for ds). Equation (7) (Reference ShimizuShimizu, 1970) gives values comparable to DPLS for all types of snow, but with a higher relative difference for most of them (from 41% for ds to 97% for ws). Equations (4) and (6) give results far from DPLS, with relative differences up to 4300% and 6000%, respectively. It has been shown that the Carman– Kozeny model does not fit to experimental data when the porous medium has high porosity, its particles are far from a spherical shape, or the porous medium is consolidated (Reference Mauran, Rigaud and CoudevylleMauran and others, 2000), which is the case for snow.
Equation (8) compares well with F DF values found with DPLS, in particular for mI, dh and ds, for which relative differences from DPLS are only 7%, 10% and 9%, respectively. The values obtained from Equation (9) compare well with that of DPLS, in particular for ws (relative difference of 14%). Finally, Equation (10) gives results comparable to DPLS for ws (relative difference from DPLS of 46%), but not for the other types of snow.
The permeability obtained by DPLS can also be compared to lattice Boltzmann modeling of permeability in firn (Reference Courville, Hörhold, Hopkins and AlbertCourville and others, 2010). However, these results have to be taken cautiously, as the samples analyzed do not come from the same snow. Firn has generally higher grain sizes than those studied here; however, ws can be similar to small-grained firn. Permeability of ws is compared with lattice Boltzmann results in Table 4. Results are similar, however, for similar grain diameter; specific surface area is much higher in firn than in ws. This could explain the higher firn permeability. Again, this difference in specific surface area shows the high differences in microstructure.
Experimental data on the permeability of snow are scarce because handling and precise measurements are difficult. This is especially true for the more permeable snow types (e.g. depth hoar). For these snow types, DPLS will possibly become the method of choice because μCT of samples casted with diethyl phthalate has become possible (Reference Heggli, Frei and SchneebeliHeggli and others, 2009). Snow is often anisotropic at different scales. Even a homogeneous sample may be anisotropic in permeability at the pore level, which is accounted for in Equation (1) when introducing a permeability tensor K (Reference KavianyKaviany, 1995). The technique applied in this study is able to capture effects as previously shown for porous ceramic structures (Reference HaussenerHaussener, 2010). At a larger scale, snow is highly layered at the mm scale (Reference Pielmeier and SchneebeliPielmeier and Schneebeli, 2003), therefore pore-level anisotropy is usually masked when measuring anisotropy in the field. DPLS is able to detect these effects and would therefore be advantageous to experimental techniques.
Summary and Conclusions
Mass transfer properties, namely permeability (K) and Dupuit–Forchheimer coefficient (F DF), were determined for five characteristic snow samples. The methodology involved first obtaining the complex 3-D geometrical representation of the snow microstructure by computer tomography. The μCT scans were digitalized and used in direct pore-level simulations (DPLS). An in-house tetrahedron-based mesh generator was used to create the computational grid directly on the μCT data. Mass and momentum conservation equations were numerically solved at the pore scale by the finite-volume method. Pressure drop over the snow sample was determined and fitted to Darcy’s law extended by the Dupuit–Forchheimer term, allowing for the determination of K and F DF. A larger pore size led to higher K, except for wet snow, for which a large pore size was compensated by a low porosity and high density.
As expected, the low K of wet snow led to a high F DF. The values of K and F DF computed by DPLS were compared with those obtained by analytical and empirical models of porous media with simplified microstructure. The conduit flow model compared particularly well with DPLS for four types of snow. Shimizu’s prediction gave reasonable agreement. The extension of the hydraulic radius theory for F DF yielded particularly good results compared with DPLS for three types of snow. The applied methodology is able to accurately account for the complex snow microstructure, which cannot be described by only a few morphological characteristics such as porosity, pore or particle size. Furthermore, it can be applied to investigate anisotropy on multiple scales. The calculated effective transport properties can be readily applied in volume-averaged (continuum) models of snow-pack for a wide range of environmental applications.