X-ray microtomography analysis of isothermal densification of new
snow under external mechanical stress
Stefan SCHLEEF, Henning LO
¨
WE
WSL Institute for Snow and Avalanche Research SLF, Davos Dorf, Switzerland
E-mail: loewe@slf.ch
ABSTRACT. We have investigated the isothermal densification of new snow under an external
mechanical stress. New snow samples that mimic natural snow were made in the laboratory by sieving
ice crystals grown in a snowmaker. This allowed us to assemble homogeneous initial samples with
reproducible values of low density and high specific surface area (SSA). Laboratory creep experiments
were conducted in an X-ray microtomograph at –208
8C for 2 days. We focused on the evolution of
density and SSA as a function of constant stress at a single temperature. External mechanical stresses
resembled natural overburden stresses of a snow sample at depths of
0–30 cm of new snow. We
demonstrate that densification increases with higher external stress and lower initial densities. We find
that the evolution of the SSA is independent of the density and follows a unique decay for all
measurements of the present type of new snow. The results suggest that details of the SSA decrease can
be investigated using carefully designed experiments of short duration which are convenient to conduct.
Additionally, we calculated the strain evolution and identify transient creep behavior that does not
follow the Andrade creep law of denser snow or polycrystalline ice.
1. INTRODUCTION
During snowfall, atmospherically grown ice crystals settle
on the ground and form a snowpack in which originally
isolated snow crystals are hard to discern, in accordance
with the common perception of a fast transition from
powder-like to a compacted and well-connected material.
The fastest changes occur within 24 hours immediately after
snowfall (Judson and Doesken, 1999). In this stage, snow is
subject to high densification rates, which are attributed to
structural changes caused by destructive snow metamorph-
ism (De Quervain, 1958; Anderson and Benson, 1963).
Visco-plastic deformation of snow under its own weight,
accompanied by surface-energy-induced metamorphism, is
relevant for all densities ranging from low-density snow
(Maeno and Ebinuma, 1983; Lo
¨
we and others, 2011) to firn
(Arnaud and others, 2000) and is therefore of fundamental
interest. The densification of new snow is also of practical
interest. The operational avalanche forecast in Switzerland
is based on the IMIS (Intercantonal Measurement and
Information System) network of automatic measurement
stations (Lehning and others, 1999). IMIS stations are
powered by solar cells and generally lack precipitation
gauges due to energy constraints, so only snow height is
measured continuously. In contrast, the prediction of
avalanche hazard usually requires estimates of new snow
mass, in order to assess, for example, the increase of stress
due to precipitation. For operational avalanche forecasts it
is therefore necessary to establish a link between snow
height changes and new snow mass via densification
modeling with operational snowpack models such as
SNOWPACK (Bartelt and Lehning, 2002), CROCUS (Brun
and others, 1989) and SNTHERM (Jordan, 1991). These rely
on a suitable constitutive equation which relates stresses to
strain rates via a ‘compactive viscosity’ (Bader, 1960) to
predict new snow densification. All processes which
contribute to the densification are lumped into parameter-
izations for the viscosity. However, these parameterizations
are based on empirical microstructural indices (Lehning and
others, 2002; Domine and others, 2008), which cannot be
measured unambiguously. To formulate new constitutive
laws in terms of objective parameters, it is necessary to
conduct measurements for the simultaneous evolution of
density and microstructural parameters, such as the specific
surface area (SSA).
Studies of snow metamorphism are commonly based on
the SSA, as the key objective quantity to avoid empirical
indices. Decrease of the SSA is the most prominent indicator
of snow metamorphism and has been studied for isothermal
metamorphism by various researchers (Flin and others,
2004; Legagneux and Domine, 2005; Kaempfer and
Schneebeli, 2007; Lo
¨
we and others, 2011), who focused
on long timescales of weeks, months or years. While long
experiments allow for robust conclusions about the SSA
decrease (Legagneux and Domine, 2005), they are very
time-consuming. To improve our general understanding of
SSA evolution it would be desirable if robust results could
also be obtained from short-term measurements by in-
creasing experimental accuracy.
Though the SSA has become a key parameter to
characterize snow metamorphism, dedicated experiments
for the simultaneous evolution of density and SSA for new
snow have not yet been conducted. This is partly because
experimental difficulties arise from the tenuous structure of
new snow. Fragile new snow samples usually prevent
techniques of casting samples in the field (Heggli and
others, 2009). In addition, microstructural processes are
sensitive to highly fluctuating meteorological conditions,
leading to great variability in snow densities obtained from
field measurements (McGurk and others, 1988; Judson and
Doesken, 1999; Roebber and others, 2003). During the past
10 years the X-ray micro-computertomograph (m-CT) has
become a fast, non-destructive standard instrument to
analyze snow metamorphism (Cole´ou and others, 2001;
Flin and others, 2004; Kaempfer and Schneebeli, 2007).
Besides the possibility of visualizing snow microstructure, it
allows multiple laboratory experiments to be conducted on
Journal of Glaciology, Vol. 59, No. 214, 2013 doi:10.3189/2013JoG12J076 233
https://doi.org/10.3189/2013JoG12J076 Published online by Cambridge University Press
the same sample, in order to monitor changes in density and
SSA of snow with high precision. However, new snow has
not yet been addressed with m-CT measurements.
Here we investigate the capability of m-CT measurements
to observe new snow densification and contribute to the
understanding of the coupled evolution of density and SSA
in general. We have restricted ourselves to isothermal
conditions at T ¼ –208C. On the one hand, the evolution is
sufficiently slow at low temperatures that the structure
remains invariant during the time (2 hours) required for a
single m-CT scan. On the other hand, the evolution is fast
enough to monitor significant changes in the structure
during 2 days at high temporal resolution. Conditions are
thus primarily chosen for experimental convenience to
ensure high-quality m-CT scans. We also avoid the
experimental difficulties of casting field samples and focus
on generating samples from ‘new snow’ crystals, grown in
the laboratory to reduce sample-to-sample variability. Note
that we use the term new snow for low-density snow with
high SSA. We outline the particularities of m-CT measure-
ments in reference to the common continuum description
of snow densification. This allows us to deduce the strain
evolution from density measurements. Thereby our results
for new snow can be compared with previous results of the
material behavior of denser snow and ice. Our simul-
taneous density and SSA measurements during creep
experiments enable us to suggest the relevance of SSA in
constitutive modeling of new snow. The impact of our
results on ‘natural’ snowfall conditions at higher tempera-
tures and the feasibility of m-CT experiments in these
conditions are discussed.
The paper is organized as follows. After a brief summary
of the theoretical background and the connection of m-CT
experiments to continuum modeling (Section 2), we de-
scribe our experimental method of preparing samples of
new snow from ice crystals grown in the laboratory
(Section 3). Section 4 is dedicated to the description of m-
CT measurements and the evaluation of the retrieved data.
The results are presented in Section 5. We show that
Lagrangian and Eulerian densification rates coincide,
enabling us to deduce the strain evolution for different
stresses. The influence on the SSA decrease is addressed
subsequently. Then the results of densification and meta-
morphism for a set of experiments with samples of varying
initial densities are presented. A summary and conclusions
are presented in Section 6.
2. THEORETICAL BACKGROUND
2.1. Snow densification
To show how we relate m-CT experiments to the material
behavior of snow, we summarize here the main aspects of
snow densification modeling. Based on the work of Bader
(1960), snow is generally treated as a quasi-one-dimen-
sional, purely viscous material in a state of uniaxial (vertical)
compression (Brun and others, 1989; Jordan, 1991; Bartelt
and Lehning, 2002). In a Eulerian frame of reference, the
evolution of the density, ðz, tÞ, at position z and time t is
governed by conservation laws. Mass conservation can be
written as
1
@
@t
þ v
@
@z
¼
_
" ð1Þ
in terms of the vertical velocity field, v, which is related
to the vertical strain rate,
_
" ¼ @v=@z. For slow snow
deformation, the momentum conservation reduces to the
force balance
@
@z
ðz, tÞ¼gðz, tÞð2Þ
for the stress ðz, tÞ, where g denotes the gravitational
constant. These equations govern the compaction of a self-
densifying snow column. The continuum starting point
(Eqns (1) and (2)) is identical for seasonal snow (Bader,
1960) and firn (Wingham, 2000). The difference between
snow types enters only via the constitutive law which is
required to close the equations. The constitutive law relates
stresses to strains and strain rates. For slow densification
under its own weight, snow is commonly regarded as a
purely viscous material with a constitutive law of the form
_
" ¼ f ðÞð3Þ
A simplified view relates the strain rate,
_
", and the stress, ,
directly by a fluid-like law
_
" ¼
ð4Þ
and lumps any dependence on microstructural details into a
compactive viscosity, . The viscosity is commonly assumed
to depend on density, (Lehning and others, 2002; Domine
and others, 2008), to capture the reduced compliance at
higher volume fractions. The microstructure is commonly
incorporated by parameterizations in terms of grain and
bond size (Lehning and others, 2002), which can be difficult
to define for a continuous microstructure or irregularly
shaped crystals. A stress dependence is sometimes con-
sidered (Scapozza and Bartelt, 2003), in order to capture the
underlying transitions in the creep mechanisms of ice. The
accompanying visco-plastic deformation of ice complicates
the situation, since it may give rise to an explicit time
dependence of the strain or strain rates, due to transient
creep, also referred to as primary creep (Petrenko and
Whitworth, 1999). Transient creep is commonly evaluated in
terms of a power law
" tðÞ¼At
ð5Þ
with a prefactor, A. If the exponent takes the value ¼ 1=3,
Eqn (5) is known as Andrade creep in dense polycrystalline
ice (Petrenko and Whitworth, 1999).
The continuum description (Eqns (1–3)) is valid on length
scales large compared with the microstructural scales. This
can be regarded as the definition of a representative
elementary volume (REV), which is in the order of
millimeters (Kaempfer and others, 2005). A cube of snow
of a few millimeters side length, which is the size typically
analyzed by m-CT experiments, represents the smallest unit
to which the above continuum description can be applied. A
m-CT scan can be thought of as taking a snapshot of the snow
structure in the control volume between the X-ray tube and a
sensor at a fixed vertical position, z, in the laboratory frame
of reference. A m-CT scan thus measures the Eulerian density
and densification rate, @=@t. According to Eqn (1) this
comprises two contributions, (1) an advective contribution,
v@=@z, due to new material entering or leaving the control
volume and (2) the strain-rate contribution,
_
", due to
intrinsic material deformation. Only the latter is related to
the constitutive equation via Eqn (3).
This has to be contrasted to a Lagrangian measurement,
where the densification is measured in a frame of reference
Schleef and Lo
¨
we: X-ray microtomography analysis of new snow densification234
https://doi.org/10.3189/2013JoG12J076 Published online by Cambridge University Press
attached to a material element. In a Lagrangian frame of
reference, mass conservation (Eqn (1)) simply reads
1
D
Dt
¼
_
" ð6Þ
in terms of the material derivative, D=Dt ¼ @=@t þ v@=@z.A
Lagrangian measurement of density, , and densification
rate, D=Dt, is directly related to the strain rate,
_
", which
allows us to address properties of the constitutive law
through Eqn (3).
2.2. Isothermal metamorphism
In parallel with visco-plastic deformation, snow undergoes
metamorphism. The driving force for isothermal meta-
morphism is the differences in surface energy between
different parts of the snow structure. This induces mass
transfer, which is commonly (Legagneux and others, 2004;
Legagneux and Domine, 2005; Kaempfer and Schneebeli,
2007; Lo
¨
we and others, 2011) interpreted in the manner of
classical Ostwald ripening, as explained by the Lifshitz–
Slyozov–Wagner (LSW) theory (Lifshitz and Slyozov, 1961;
Wagner, 1961) for an assembly of spheres. While the sphere
case can be solved analytically, this is not possible for a
bicontinuous structure such as snow. In analogy to the
increase of the average sphere radius in the LSW theory,
isothermal snow ripening is generally analyzed by the
decrease in SSA, as the simplest indicator for the
complicated processes during metamorphism. It has been
widely confirmed (Legagneux and Domine, 2005; Kaemp-
fer and Schneebeli, 2007) that the decrease of the SSA can
be well characterized phenomenologically by a decay law
of the form
SSAðtÞ¼SSAð0Þ
t þ
1=n
ð7Þ
where and n are parameters of the model. Monte Carlo
simulations that take into account pore diffusion and
surface diffusion are also highly compatible with this
relation (Vetter and others, 2010). Exponent n is expected
to signal the dominant process of underlying mass transport
(Legagneux and others, 2004; Legagneux and Domine,
2005; Kaempfer and Schneebeli, 2007; Lo
¨
we and others,
2011). This power-law regime can be observed after the
crossover time, , which is, however, difficult to reach
experimentally (Legagneux and Domine, 2005). The inverse
quantity, 1=, defines a characteristic rate of SSA decrease,
which is expected to depend on snow density (Legagneux
and Domine, 2005), due to an enhanced mass transfer for a
denser arrangement of grains. This density dependence is
analogous to the effect of volume fraction on the coarsening
rate in the LSW theory (Ratke and Voorhees, 2002).
3. SAMPLE PREPARATION
In order to investigate the impact of external mechanical
stresses on new snow settlement it is particularly important
to prepare snow samples with identical structural character-
istics. Experiments have to be reproducible and this requires
careful sample preparation, because of the fragile structure
of new snow. Since new snow has large natural variations in
the microstructure and it is inherently difficult to collect
well-defined new snow samples in the field, we used
artificial snow that mimics natural snow but is grown in the
laboratory. This was done using a snowmaker (in a similar
way to Nakamura, 1978). In a cold laboratory at an ambient
temperature of –258C, a box of water is heated to +308C. A
fan blows cold air over the warm water surface to produce a
flow of supersaturated air. The air then rises through a net of
nylon wires, promoting nucleation. Under these temperature
conditions dendritic ice crystals grow, and these are broken
off by shaking the wire and fall as snow. (For further details
see Lo
¨
we and others, 2011.) The fresh snow was collected
and put into a cooler at –608C. At this temperature
metamorphism is minimized (Kaempfer and Schneebeli,
2007) and is hardly detectable, and the microstructure of the
new snow was thus preserved until we prepared the samples
for the creep experiments.
Generally the snowmaker produces snow that also
contains large directionally grown ice crystals, which are
rather unrealistic for naturally fallen snow. They may lead to
different material properties when compared with natural
snow if sieved directly in the sample holders (Gergely and
others, 2010). To avoid this we employed a pre-sieving
procedure through a sequence of vibrating meshes with
decreasing mesh size, down to 500 mm, for 30 min. Thus,
larger ice crystals were discarded and we obtained grain
sizes typical of natural snow (Fierz and others, 2009). By
keeping only finer snow, a smaller REV was obtained. In this
way, snow densities of 100 kg m
3
were produced, which
fit well with the 10 : 1 forecasting rule of thumb, which states
that new snow height is ten times the snow water equivalent
(McGurk and others, 1988; Judson and Doesken, 1999;
Roebber and others, 2003). To mimic snowfall conditions
the pre-sieved snow was manually sieved into the m-CT
sample holders. A rotating table allowed uniform filling and
preparation of multiple samples at the same time. The m-CT
sample holders were polyetherimide (PEI) cylinders of
18 mm diameter (Fig. 1). As shown below, this procedure
gave rise to very similar values of initial densities and SSA.
In addition, we were interested in varying the initial
density without damaging the structure under otherwise
identical conditions. This was achieved by placing the
samples on a vibrating plate for 0.5–1.5 min directly after
sieving. Snow crystals are only poorly sintered at this stage
and the vibrations enabled particle rearrangements, which
Fig. 1. Schematic of the snow-filled sample holder with a weight on
the top of the snow layer.
Schleef and Lo
¨
we: X-ray microtomography analysis of new snow densification 235
https://doi.org/10.3189/2013JoG12J076 Published online by Cambridge University Press
resulted in a higher packing density, similar to dry granular
media (Richard and others, 2003). Though it was not
possible to adjust the density to prescribed values, this
procedure allowed us to study the main influence of the
initial density.
For all samples, filling levels of snow in the sample
holders were adjusted by carefully removing snow with a
center bit that fits the 18 mm inner diameter of the sample
holder. In this way identical filling levels of 15 mm height
were attained for all samples without compacting them. All
samples were sealed with a cap to disable sublimation and
stored in a box at –608C to preserve the state until
measurement. To avoid temperature fluctuations when
taking out individual samples, the storage box was insulated
by several-centimeter-thick styrofoam walls and also con-
tained a steel plate at the bottom.
4. X-RAY MICRO-COMPUTERTOMOGRAPH
MEASUREMENTS
One hour before starting a m-CT measurement the appropri-
ate sample was taken out of the –608C cooler and put in a
cold laboratory that was kept at constant temperature of
–208C, to allow thermal equilibration. To simulate varying
depths of the sample in the snowpack, different loads were
applied, by placing weights on top of the snow layer before
sealing the sample holder. The weights were steel or
aluminum cylinders of 18 mm diameter. We calculated a
nominal stress from the cross-sectional area and the applied
load. The resulting nominal stresses applied on top of the
snow were 0, 133, 215 and 318 Pa. For a density of
100 kg m
3
this corresponds to burial depths of 0, 10, 20
and 30 cm, all of which are reasonable values for common
snowfall events. Using weights instead of real snow enabled
us to keep sample heights small and thereby limit perturb-
ation from the contact area of the sample with the boundary.
In contrast to the gradual stress increase found in natural
conditions during snowfall, the stress was applied instant-
aneously. If the weight was placed carefully an immediate
collapse-like compaction of the sample was avoided. We
compared the density profile before and after inserting the
weight to confirm that local grain rearrangement only
occurred in the vicinity of the weight. Since the diameter
of the weight is slightly smaller than the inner diameter of
the sample holder, air permeability is allowed and thus the
potential influence of moist air compression is eliminated.
To capture sample fluctuations and assess reproducibility,
two experiments were conducted under identical conditions
for each stress. For the sample with a stress of 133 Pa we
conducted only one experiment due to scanning problems.
(The later experiments, with varying initial densities, were
performed on a separate set of samples, where a stress of
215 Pa was always applied.)
The measurements were made using a micro-computer-
tomograph from SCANCO Medical AG (m-CT 80). We used
45 keVp X-rays and a nominal resolution of 10 mm voxel
size. In total, 630 slices of the whole sample diameter were
scanned, yielding an observed snow height of 6.3 mm. The
scanned region was at a fixed position in the laboratory
frame, located vertically in the middle of the snow samples
(Fig. 1). A single measurement took 2 hours, and 16
measurements were taken with a time-step of 3 hours,
leading to a total observation time of 48 hours. The fixed
time-step between the measurements guaranteed identical
conditions in the m-CT to avoid measurement bias (e.g. by
differences in the X-ray tube temperature). The heat
generated by the X-ray tube also perturbed the desired
isothermal conditions for the sample. A fan was used to
compensate for this. Temperature measurements with a
sensor (iButton TMEX) placed in a sample holder confirmed
that the temperature increase is <18C during a measurement.
We assumed that these small temperature fluctuations had
no significant influence on the snow evolution. In all
samples, mass loss by evaporation was strictly avoided by
the sealing.
For evaluation of the m-CT raw data, a cubic volume of
6.3 mm edge length was extracted from the scan region. The
size of the cube is sufficiently large to meet the REV
requirements found by Kaempfer and others (2005) and
Cole´ou and others (2001). The position of the cube was
chosen between the boundary and the center to minimize
potential influences of the lateral boundary and center
artifacts in the rotation center of the m-CT.
To reduce measurement noise, grayscale images were
Gaussian-filtered. They were then segmented into binary
images of the evaluated volume containing only ice and air.
Unfortunately, for new snow it is not possible to determine
the threshold by analyzing the histograms of the grayscale
images. The ice matrix contains very fine structures of size in
the order of the resolution of the m-CT, which results in a
smooth transition between the adsorption coefficients of ice
and air. For our snow the density changed almost linearly
with variation of the threshold. Different threshold values
were used and volume fractions were derived and compared
with those obtained from mass and volume measurements of
the entire sample. Eventually we maintained a constant
threshold for all measurements. With the fitted threshold
value we obtained three-dimensional images of new snow,
as shown in Figure 2.
Figure 3 shows a close-up of the image, in which the
complex structure of new snow becomes apparent. It is
easy to identify structural features at the end of 2 days of
measurements, even though they are at a lower position in
the scanned volume, due to densification. With the naked
eye no major differences between the two images in
Fig. 2. Three-dimensional image of a new snow sample. Cube edge
length 6.3 mm.
Schleef and Lo
¨
we: X-ray microtomography analysis of new snow densification236
https://doi.org/10.3189/2013JoG12J076 Published online by Cambridge University Press
Figure 3 are visible, but with the m-CT a more precise
analysis is possible.
From the segmented binary image, ice volume fractions
were calculated by counting voxels and using triangulation
(Guilak, 1994) to obtain the snow density. Comparison of
the two methods reveals that volume fractions obtained by
triangulation are always smaller with a constant offset.
Henceforth, we use the volume fractions obtained from
triangulation. The snow density was obtained by multiplying
this with the ice density,
ice
¼ 919:7kgm
3
for
T
snow
¼208C. The SSA was calculated by dividing the
triangulated surface area of the ice matrix by the ice volume.
The time series of densities and SSA for each sample were
obtained by assigning the respective values to the beginning
of each scanning time-step of 3 hours. The initial time, t ¼ 0,
was set to the beginning of the first m-CT scan, which was
started exactly 1 hour after putting the sample into the
isothermal conditions of the cold laboratory, and 30 min
after applying the stress.
To determine the accuracy of the m-CT measurements it
would be necessary to compare multiple scans of an
invariant sample. Since this is not possible with our new
snow samples, we used an experiment with a different
snow type, in which the density of 100 kg m
3
remains
constant during 2 days of measurements. This series of
measurements has a standard deviation of 0.2 kg m
3
,
which we take as an estimate for our instrumental accuracy.
Note that there might be a much higher offset, of up to 5%,
for the absolute value of the density, due to the uncertainty
of the threshold determination. But, since we focus here on
density changes with time, an almost constant offset is not
essential for our results.
5. RESULTS AND DISCUSSION
5.1. Lagrangian and Eulerian densification
First we address the homogeneous nature of the densifi-
cation and potential differences between Lagrangian and
Eulerian densification rates.
As explained in Section 4, the standard m-CT scanning
and evaluation monitors the densification in a cube
fixed in the Eulerian (laboratory) frame of reference. The
evaluation cube contains different material at each scan-
ning time. This is displayed in Figure 4, in which
horizontally averaged density profiles at three different
time-steps are shown as a function of vertical position for an
experiment at an applied stress of 215 Pa. The displacement
of characteristic maxima and minima of the density profile
is easily identified at different times (arrows in Fig. 4),
which indicates the advective component of the densifica-
tion ‘flow’. Thereby mass enters the m-CT control volume
from the top while other mass leaves it at the bottom. The
snow inside the cube is advected as a whole; at the same
time the structure is strained during 2 days of densification.
In order to separate the two effects, we computed the
densification in a Lagrangian frame of reference by image
registration. To this end, we tagged two horizontal slices
(top and bottom of the evaluation cube) at t ¼ 0 and
followed their z position in the laboratory frame with time.
The vertical displacement was found by computing the
maximum cross-correlation of a slice with neighboring
slices at a lower position in later scans. The sub-volume
bounded by the two slices contains identical parts of the
material, i.e. it represents a Lagrangian material element.
The densification computed from the sub-volume reduction
thus yields the Lagrangian density. The resulting values are
shown in Figure 5, where they are compared with the
Eulerian measurements, and are shown to be in good
agreement. This implies that the sample is homogeneous
and the density gradient in the advection term in Eqn (1)
can be neglected.
Fig. 4. Evolution of the horizontally averaged density profile of one
sample in the evaluation window (z-axis oriented downwards).
Fig. 3. Zoomed-in image of the new snow structure. Top: structure
at the beginning of a series of measurements at t ¼ 0 with SSA
0
(top-left edge of Fig. 2). Bottom: same structure for the last scan at
t ¼ 45 hours (at a deeper position in the scanned volume).
Schleef and Lo
¨
we: X-ray microtomography analysis of new snow densification 237
https://doi.org/10.3189/2013JoG12J076 Published online by Cambridge University Press