1. Introduction
When a snow avalanche is flowing in a steady state, the dissipated work Φ is in balance with the gravitational work rate W g,
This is a simple statement of the principle of energy conservation for non-conservative systems in which the change in mechanical energy is equal to the work done by the frictional or dissipative processes. In our steady-state system, the change in mechanical energy is the work done by gravity, the change in kinetic energy being zero. The two primary processes to dissipate work in granular snow avalanches are (1) the churning and collisional movement between the snow granules within the avalanche body (Fig. 1a) and (2) the abrasive sliding of the granular core on the running surface, usually a hard, densified layer of old snow. Chute experiments with snow (Reference Kern, Tiefenbacher and McElwaineKern and others, 2004) show that the granular churning and collisional movement within the avalanche body is concentrated in a highly sheared fluidized layer, located at the bottom surface of the avalanche. The bulk of the avalanche flows with uniform velocity in a flow plug where no work is dissipated. This flow behaviour has also been observed in field experiments (Reference GublerGubler, 1987; Reference Dent, Burrell, Schmidt, Louge, Adams and JazbutisDent and others, 1998). The dissipated work is transferred entirely into heat and represents the change in internal energy of the granular flow system. Recent temperature measurements in snow avalanches have shown that temperature increases of a few degrees are possible (Reference Miller, Adams, Schmidt and BrownMiller and others, 2003).
Fig. 1. (a) The granular deposits of a large dry mixed flowing/powder avalanche released at the Vallée de la Sionne test site, Switzerland. The granules are 5_10cm in size. (b) The experimental snow chute at the Weissfluhjoch, Davos, Switzerland. The picture shows devices to measure the sliding friction and velocity profile in the fluidized layer. For more information see Reference Tiefenbacher and KernTiefenbacher and Kern (2004) or Reference Kern, Tiefenbacher and McElwaineKern and others (2004).
The purpose of this paper is to demonstrate four further properties of the dissipation function Φ. Namely,
In steady state, the dissipated work is as small as possible within the given system constraints, δΦ = 0.
The magnitude and form of the steady-state velocity profile can be found either by balancing the momentum, Λ = 0, or, equivalently, finding the stationary state of the dissipated work functional, δΦ = 0.
The dissipated work functional Φ describes the stability of the flow system with respect to fluctuations in velocity. Unstable regimes lead to flow transitions.
The statement δΦ = 0 is a special case of a general theorem of irreversible thermodynamics: the principle of minimum entropy production.
These four statements will be demonstrated using large- scale, steady-state chute experiments with snow (Fig. 1b). The experimental set-up allows the investigation of the two primary dissipational mechanisms because both the internal deformations in the granular core and the basal slip velocity are measured (Reference Tiefenbacher and KernTiefenbacher and Kern, 2004).
Granular snow avalanches are complex dissipative systems containing a vast number of interacting snow particles. From a practical standpoint it is both unreasonable and unnecessary to track the motion of each particle, especially when the initial conditions (the starting mass) or the boundary conditions (the snow cover and track terrain) are not exactly known. Macro-continuum models (Reference Bartelt, Salm and GruberBartelt and others, 1999), which treat the complex collisional and frictional interactions in bulk, are now used extensively to predict avalanche run-out distances and flow velocities in complex terrain. Although these models work well (provided one has a clear idea of the empirical parameters), they cannot treat important physical phenomena arising from the complex interactions, usually when the avalanche is far from steady state. Such phenomena include fluctuations in velocity and pressure at the head of the avalanche – an important factor for structural engineers charged with the task of designing structures in avalanche run-out zones – and flow regime bifurcations leading to different avalanche types. Whether an object is struck by a heavy flowing avalanche, a ‘mixed’ avalanche or a powder avalanche plays an important role in developing hazard mitigation strategies.
The statement that the dissipated work is as small as possible within the given system constraints, δΦ = 0, is a statement of tendency which differs significantly from the use of mass-continuity and momentum-balance relationships that are presently the basis of the macro-continuum avalanche models. Mass and momentum conservation describe the unique state of equilibrium exactly. The principle δΦ = 0 allows us to treat states of flow outside but on the way to equilibrium (Reference Glansdorf and PrigogineGlansdorf and Prigogine, 1974). In fact, it allows us to pose the question whether, for a specific set of flow constraints (terrain roughness, slope, snow characteristics), it is even possible to reach equilibrium. Thus, if we can show the applicability of the statement δΦ = 0, we might have an instrument to treat important practical problems such as the question of whether flow fluctuations will regress on the way to steady-state equilibrium, which avalanches do in fact exhibit (Reference Vallet, Turnbull, Joly and DufourVallet and others, 2004).
This idea is not new. Helmholtz and Kroteweg showed that δΦ = 0 holds for a viscous Newtonian fluid (Reference LambLamb, 1945). Later the principle was shown by Reference PrigoginePrigogine (1955) to be a general principle of macroscopic irreversible thermodynamics, that is, the thermodynamics of dissipative systems such as avalanches. Presently, it is the starting point for new theories concerning the ‘form’ of the natural world (constructal theory) (Reference BejanBejan, 2000) or self-organization in nature (Reference Nicolis and PrigogineNicolis and Prigogine, 1989). Our interest, however, is much more practical: we seek a tool to help us quantify the complex, seemingly random behaviour of avalanches that we can use in conjunction with our numerical macrocontinuum models.
The statement δΦ = 0 is a variational statement. Thus, to see if δΦ = 0 holds for snow flows – non-Newtonian fluids with complex slip boundary conditions – means using variational calculus on non-trivial systems. For this reason, we do not begin immediately with the avalanche problem. Instead, we solve two simple problems (Couette flow between two moving plates and Hagen-Poiseuille flow in a pipe) before turning our attention to flowing avalanches. This serves to introduce the reader to constrained variational problems. We will show that without the mass- continuity constraint, only trivial solutions are possible. Afterwards, the avalanche problem can be more efficiently treated.
2. Two Hydrodynamic Systems With Λ = 0 and δΦ = 0
In this section, we apply the variational principle to two simple hydrodynamic flow systems. Our first goal is to establish the equivalence of the momentum and stationary dissipation solutions. Secondly, we want to demonstrate the importance of the mass-continuity constraint for the variational solution (Reference LiuLiu, 1972). This constraint leads to variational problems of the isoperimetric kind (Reference Gelfand and FominGelfand and Fomin, 2000). Without it, the equivalence between Λ = 0 and δΦ = 0 cannot be established.
The first example problem considers an incompressible, viscous Newtonian fluid (viscosity μ) placed between two moving rigid plates of infinite extent in the x direction (Fig. 2a). The plates are located a distance h apart. The upper plate is moving with constant velocity u p; the lower plate is moving with constant velocity u 0. Furthermore, there is no movement of the fluid in the z direction (parallel flow). Solution of the momentum equation for this system (Reference Incropera and DeWittIncropera and DeWitt, 2002)
leads to
since there is no pressure gradient ∂p/∂x along x, ∂p/∂x = 0. The constants of integration C 1 and C 2 are found by applying the boundary conditions u(0) = u 0 and u(h) = u p:
Fig. 2. (a) A rigid plate moves with velocity u p and shears a Newtonian fluid of height h. The velocity profile that minimizes the work dissipation is linear in z. (b) Hagen-Poiseuille flow. A Newtonian fluid moves with mean velocity Um in a tube of height h. Non-slip boundary conditions are imposed at the upper and lower surfaces. The velocity profile that minimizes the work dissipation is parabolic. However, this solution is found only after respecting the mass-flow constraint.
Hence, the solution u(z) to this problem is the well-known linear Couette velocity profile.
Now consider the minimum dissipation solution. The dissipated work per unit volume is (Glansdorf and Prigogine,
since the shear stress pzx is related to the velocity gradient by the phenomenological relation
for a Newtonian fluid. We are seeking the function u(z) which renders the integral
a minimum, or δΦ″ = 0. The double-prime superscript denotes a quantity per unit area; the triple-prime superscript a quantity per unit volume. The stationary solution δΦ″ = 0 is most easily found by solving the associated Euler- Lagrange equation for u(z) (Reference Gelfand and FominGelfand and Fomin, 2000):
The subscripts denote partial derivatives with respect to the corresponding arguments. The function u(z) which satisfies the Euler-Lagrange differential equation is of the form
where the constants of integration satisfy the boundary conditions at the top and bottom plates. Because the solutions, Equations (3) and (9), are exactly the same, we conclude that the linear velocity profile not only satisfies the momentum balance (Λ = 0) but is also the form that minimizes the work dissipation (δΦ″ = 0). Previously, Reference BejanBejan (1996) determined the viscous dissipation for Couette flow assuming a linear velocity profile. The approach here differs somewhat since we make no assumption regarding the form of the profile, but show that the linear profile is the form, of all possible forms, that dissipates the least work possible for the prescribed boundary conditions.
Of importance in the Couette flow case is the fact that ∂p/∂x = 0; that is, the motion of the fluid between the plates is not sustained by a pressure gradient, but rather by the shearing motion of the plates and the non-slip boundary conditions. For our second example, Hagen-Poiseuille flow (Fig. 2b), this is not the case since there is a non-zero, constant pressure gradient ∂p/∂x driving the flow. Again, our goal is to find the form of the velocity profile u(z) such that the work dissipation is a minimum. The solution of the momentum-balance equation, Equation (2), assuming a constant pressure gradient, leads to a parabolic velocity profile,
where C 1, C 2 and C 3 are the constants of integration,
Two constants are found from the boundary conditions u(0) = 0 and (du/dz)(h/2) = 0 (or u(h) = 0 and (du/dz)(h/2) = 0 from symmetry). Non-slip boundaries are imposed at the upper and lower surfaces of the flow which are not moving. The third constant is found by substituting Equation (10) into the momentum equation, Equation (2). The mean velocity of the flow is calculated according to
The variational δΦ″ = 0 solution of the Hagen-Poiseuille problem appears at first inspection to provide the wrong answer. Since the dissipation is still given by an incompressible Newtonian fluid, the Euler-Lagrange equation of the Couette flow case (Equation (8)) is still applicable. Imposing the boundary conditions, we find
That is, no flow over the entire domain. Not surprisingly, the minimum work (in this case no work) is dissipated when there is no flow. This is a correct, but trivial solution. More useful solutions can be found by imposing an integral constraint, written generally as:
where l is some non-zero value and G is a constraint functional. If G = u(z)/h, then we recover the equation for the mean or depth-averaged velocity U m:
In the following (also in the upcoming avalanche case) we will use this specific functional, which we term a mass- continuity or mass-flow constraint because the density p is constant. To find δΦQ″ = 0 with the constraint functional, we must find a velocity profile u(z) such that the functional
is a minimum where λ is the Lagrangian multiplier for the constraint condition, Equation (14). The profile u(z) for which δΦ″ = 0 must satisfy the differential equation is (Reference Gelfand and FominGelfand and Fomin, 2000)
In variational calculus this is an isoperimetric problem, discussed in many standard textbooks (e.g. Reference SmithSmith, 1998; Reference Gelfand and FominGelfand and Fomin, 2000). The solution to Equation (17) is a parabolic function, similar to Equation (10). The two boundary conditions u(0) = 0 and u′(h/2) = 0 can be used to evaluate two of the three integration constants,
The integration constant C 1 is still undetermined; this can be found using the mass-continuity constraint, Equation (15):
Thus
Finally, the Lagrangian multiplier can be found from the Euler-Lagrange equation, Equation (17)
For the Λ = 0 and δΦ = 0 solutions to be equivalent, the constants of integration (Equations (11) and (20)) must be equal. This is the case if, and only if,
but this is true by definition of the mean velocity, Equation (12). Therefore, we conclude once more that the Λ = 0 and δΦ = 0 solutions are equivalent. For Hagen- Poiseuille flow the parabolic velocity profile satisfies the momentum balance and is the only form that dissipates the least work.
Differences between the two solutions, however, exist. In the Λ = 0 solution, the pressure gradient – or, applied external forces – is assumed to be known, whereas in the δΦ = 0 solution the mass flux, represented by the mean velocity U m, is known. The pressure gradient and mass flux cannot be chosen independent of each other in steady state. This fact was well known to Lamb in 1879 (Reference LambLamb, 1945). Note that the Lagrangian multiplier λ has dimensions of pressure and acts as a force in the variational solution imposing the flow condition, but also restoring equilibrium. Nonetheless, a fundamental difference exists between the two solution treatments. In the Λ = 0 solution, forces are balanced which are vectors and therefore have a directional component. In the δΦ = 0 solution, we consider a scalar – the dissipated work – which is constrained by another scalar – the mass continuity. As we shall see in the next section, this is advantageous when considering multidirectional, complicated flow phenomena, such as snow avalanches.
3. Avalanches: Simple Shear Flow Systems with Λ = 0 and δΦ = 0
The methodology discussed in the preceding section will now be applied to granular avalanches flowing in simple shear. This requires a slight modification to our variational procedure. In the previous two examples, we stated the dissipation functional Φ‴ and then found and solved the Euler-Lagrange equation for a continuous velocity function u(z). Since the non-slip boundary conditions were known, both problems could be solved; we obtained the form (linear or parabolic) of the velocity profile.
In the avalanche problem we are confronted with slip boundary conditions, governed by some sliding friction law, usually expressed as a function of the sliding velocity u 0. This changes the nature of the variational problem, because work is now dissipated at the sliding boundary, unlike the previous two cases. In fact, most flowing-avalanche models, including the well-known Voellmy models used for run-out calculations in Switzerland (Reference Bartelt, Salm and GruberBartelt and others, 1999), assume that all the work is dissipated at this sliding boundary. This means that in the avalanche problem the slip velocity u 0 has an increased importance because it is at the interface of two dissipational processes: it determines the degree of shear deformation in the core as well as the degree of sliding friction. In the avalanche problem, we are no longer primarily interested in the form of the velocity profile, which, as we have shown, can be determined equally well from the momentum Λ = 0 solution, but are interested in the magnitude of the slip and core velocities, which are discrete values.
This problem is resolved by formulating a discrete variational problem in terms of the slip (u 0) and plug (u p) velocities of the avalanche. The mass-continuity constraint will remain the same as before. That is, in order to obtain the correct solution, we must prescribe a flow rate that is in agreement with the given external gravity force. The first step is to state the dissipation functional Φ‴ in terms of u 0 and u p, which, in turn, requires constitutive formulations for flowing snow and basal friction. The assumption of a Newtonian fluid cannot be applied with conviction in the face of the many experimental investigations (Reference Dent and LangDent and Lang, 1983; Reference NishimuraNishimura, 1990; Reference Kern, Tiefenbacher and McElwaineKern and others, 2004; Reference Tiefenbacher and KernTiefenbacher and Kern, 2004) which we briefly review in section 3.1. Therefore, unlike the simple problems discussed above, we apply two non-Newtonian constitutive models that have been proposed by Reference Dent and LangDent and Lang (1983) and Reference Norem, Irgens and SchieldropNorem and others (1987) to describe steady-state avalanche motion. Although we would like to find the best constitutive model for flowing snow, we are primarily interested in the general properties of the functional for snow avalanches.
3.1. Chute experiments and field observations of granular snow avalanches
Chute experiments with snow show that the velocity profile above the basal layer can be divided into two distinct flow layers: a flow plug moves with a uniform velocity above a highly sheared fluidized layer with a non-zero basal slip velocity (Reference Nishimura and MaenoNishimura and Maeno, 1987;Kern and others, 2003) (Fig. 3). The measurements show that the shear strain rates 
in the fluidized layer are of the order of 50 s-1 whereas the shear strain rates in the flow plug are 
The few velocity profile measurements in real avalanches appear to corroborate the chute experiments. The velocity profile of the Aulta avalanche in Switzerland, captured with ground radar, exhibits the two layers (Reference GublerGubler, 1987). The fluidized layer of this avalanche had a depth of 0.2-0.4m (the total height of the avalanche was between 0.6-1.0 m). The mean shear strain rate in the fluidized layer was of the order of 65 s-1. Velocity profile measurements in the steady-state region of a small dry-snow avalanche in 1994 at the ‘Revolving Door’ (Montana, USA) test site again show a 0.15 m high passive flow plug 
above a highly sheared, thin (0.01-0.05 m high) layer with strain rates 
(Reference Dent, Burrell, Schmidt, Louge, Adams and JazbutisDent and others, 1998). Recent velocity measurements of a mixed flowing/powder avalanche with optical velocity sensors mounted on the pylon of the Swiss Vallee de la Sionne test site (Reference TiefenbacherTiefenbacher, 2003) reveal small shear deformations in the flow body, 
indicating again the existence of a near-uniform flow regime. The measurements were captured in the lower 2.5 m of the avalanche where a granular core is expected (but not certain). Because the sensors are spaced 1 m apart, a fluidized layer could not be identified.
Fig. 3. Velocity profile of a two-layer flow model defining the internal flow parameters: u0 slip velocity, up plug velocity, hs fluidized or shear layer height, h p plug layer height. The avalanche is flowing in steady state on a slope with inclination angle ϕ with mean velocity U m and height h. Variational solution. Discrete variational solution Ṡ‴ per unit volume. Work is dissipated within the fluidized layer and at the sliding surface. For a Voellmy fluid, work is not dissipated in the avalanche flow body.
In field experiments, it is difficult, perhaps impossible, to determine the location of the slip surface and therefore the nature of the slip velocity. A sensor must be placed at the slip surface which is usually not known in advance. Both Reference GublerGubler (1987) and Reference Dent, Burrell, Schmidt, Louge, Adams and JazbutisDent and others (1998) assume zero slip velocities when analyzing their experimental results. Likewise, Reference Nishimura and MaenoNishimura and Maeno (1987) fitted the results of their chute experiments assuming the non-slip condition, u 0 = 0. Nishimura and Maeno, however, did acknowledge, and their results also show, definite non-zero slip velocities, u 0 ≠ 0. Chute experiments with ice spheres showed significant amounts of slip (Reference Nishimura, Kosugi and NakagawaNishimura and others, 1993). Alternatively, Reference Kern, Tiefenbacher and McElwaineKern and others (2004) assumed a non-zero u 0 because they were able to resolve the shear layer near to the slip surface. The experimental results indicate that the slip velocity is a function of the surface roughness of the chutes. Reference Nishimura and MaenoNishimura and Maeno (1987) covered their chute with a polyethylene film in order to minimize the basal friction. However, at the locations where the velocity profile was measured, the surface was covered with a rough cloth or sifted snow and they therefore obtained small slip velocities relative to the plug velocities. It should be noted that Kern and others placed rubber mats on the surface of their chute in order to get realistic chute avalanches (Reference Tiefenbacher and KernTiefenbacher and Kern, 2004) and found large slip velocities.
In the following we will use four chute experiments reported in Reference Kern, Tiefenbacher and McElwaineKern and others (2004). A description of the chute is contained in Reference Tiefenbacher and KernTiefenbacher and Kern (2004). The results of the chute experiments where mean velocities, basal slip velocity and velocity profiles were measured are summarized in Table 1.
Table 1. Summary of four fluidized-layer experiments from the Swiss Weissfluhjoch snow chute. For more details see Reference Kern, Tiefenbacher and McElwaineKern and others (2004). The measured velocity profiles are also depicted in Figures 5, 7, 8 and 9
The observed velocity profiles are given by the general formula:
where u p is the plug velocity, u 0 is the slip velocity at the basal layer, h is the avalanche flow height, h s is the fluidized layer height and α is a parameter which determines how the velocity in the fluidized layer varies with z, the height coordinate measured perpendicular to the basal surface (z = 0). The parameters are depicted in Figure 3 and specific values for these parameters are listed in Table 1.
The shear strain rate is
If further experiments do not show the bi-layer behaviour, the formulas can still be applied by setting h s = h. Of course, pure plug-flow behaviour is also possible by setting h s = 0.
3.2. Work dissipation in the avalanche core Φc
The experimental results show that work is dissipated by two primary mechanisms. The first is the internal work dissipation in the fluidized layer, Φc; the second is the sliding work dissipated at the base of the avalanche, Φ0. No work is dissipated in the flow plug because the shear strain rates are zero (Fig. 3). Additionally, we assume that the longitudinal strain rates are zero in steady-state flow du/dx ≈ 0, although small fluctuations have been measured in the chute experiments (Reference Kern, Tiefenbacher and McElwaineKern and others, 2004). This assumption avoids the inclusion of an additional material parameter, either an active or passive pressure term (Reference SalmSalm, 1993) or a normal viscosity (Reference Norem, Irgens and SchieldropNorem and others, 1987). These can be included, if necessary, in a next step.
The shear stress in the avalanche core is given by
where p yield is the yield stress, p is the flow density and m the shear viscosity. The parameter n defines the constitutive model: n = 1 for a Bingham fluid and n = 2 for a dilatant fluid. The dimension of the shear viscosity m is a function of the choice of n. If n = 1, then m has dimensions m2s-1; if n = 2 then m has dimensions m2. The n = 1 Bingham model has been investigated by Reference Dent and LangDent and Lang (1983); n = 2 models have been proposed by Reference Norem, Irgens and SchieldropNorem and others (1987). Reference NishimuraNishimura (1990) found that the shear stress is independent of the shear rate at low shear rates (<10-1 s-1). Above 20 s-1 the shear stress depends on the square of the shear rate (n = 2). For shear rates between these two limits, Reference NishimuraNishimura (1990) proposed a linear relationship, i.e. n = 1. Thus, Nishimura proposed a constitutive model similar to Equation (27), but containing both the n = 1 and n = 2 terms simultaneously.
The yield stress p yield accounts for the plug behaviour observed in the experiments. When the shear stress is less than p yield the core can sustain an applied stress without shearing. A Mohr-Coulomb yield equation is used to define the yield limit:
where c is the flowing-snow cohesion and b is the internal friction parameter (related to the internal friction angle δ according to b = tan δ). The normal stress pzz is given by
where ϕ is the angle of the slope and g is the gravitational acceleration. The normal pressure varies linearly in height. Thus, the constitutive model assumes that the shear stress pxz is a function of the overburden stress. Finally, because the top surface of the avalanche is stress-free, we set c = 0 in order to fulfil this boundary condition.
Experimental measurements which are analyzed assuming a linear relationship between shear stress and shear strain rates have found m values between 10—5m2s-1 (Reference Nishimura and MaenoNishimura and Maeno, 1987) and 10—3m2s-1 (Reference Dent and LangDent and Lang, 1983). In fluidized snow experiments, Reference Maeno, Nishimura and KanedaMaeno and others (1980) found constant viscosity values of the order of 10—4 m2 s-1, but a function of the particle size and density of the snow. Reference NishimuraNishimura (1990) analyzed his snow chute experiments assuming a dilatant law; he found dilatant viscosities of the order of m ≈ 10—5 m2.
Not directly expressed in the shear stress relationship is the fact that the velocity profile exponent is α = 2 for a Bingham fluid and α = 1.5 for a dilatant fluid. These values are found by noting that in steady state the shear stress must be linear in z. (For a derivation see Reference Nishimura and MaenoNishimura and Maeno, 1987.)
The work dissipated per unit volume in the core is
The work dissipated per unit area of the avalanche is found by integrating the flow height h:
Evaluation of this integral leads to:
For clarity, we state the work dissipated in the core for a Bingham fluid (n = 1, α = 2)
and a dilatant fluid (n = 2, α = 1.5)
3.3. Work dissipation at the basal sliding layer Φ0
At the bottom (z = 0) flow surface of the avalanche, friction is given by:
Basal friction is thus composed of a dry Coulomb-like friction, dependent on the normal pressure pzz and a velocity-dependent friction (Reference Norem, Irgens and SchieldropNorem and others, 1987; Reference SalmSalm, 1993). The friction parameter b is a dry snow-to-snow friction and is assumed to be equal to the snow-to-snow friction within the fluidized layer. The friction parameter s is dependent on the roughness of the sliding surface. The velocity-dependent friction term is introduced to ensure states of steady flow. Without this term, avalanches reach unrealistic terminal velocities on long avalanche tracks. This law has been used successfully to predict avalanche run-out distances in Switzerland and is part of the Swiss guidelines on avalanche run-out calculation (Reference Salm, Burkard and GublerSalm and others, 1990); b and s values for avalanches are known (Reference Buser and FrutigerBuser and Frutiger, 1980).The work dissipated per unit flow area on the sliding surface is
3.4. Total dissipation and gravitational work rate
The total dissipated work per unit area for our steady-state flowing avalanche is given by the sum of the work dissipated in the core (Equation (32)) and at the sliding surface (Equation (36)):
As stated in the introduction, conservation of energy demands that the total dissipation be in balance with the gravitational work rate,
The above equation assumes that the change in potential energy affects only the mean translational velocity U m. There is no change in the angular velocity of the granules and subsequently there is no change in the rotational energy of the system.
The total dissipated work is a function of the internal flow variables, u p, u 0 and h s; α and n are given by the constitutive formulation; b, m and s are the constitutive parameters given by the materials involved in the flow; the flow height h and mean velocity U m are constrained variables in the sense that they are given by the mass continuity.
The total dissipated work as a function of the slip velocity u 0 is shown in Figure 4 for (a) a Bingham fluid (n = 1, α = 2) and (b) a dilatant fluid (n = 2, α = 1.5). The mean flow height (h = 0.4 m) and shear layer height (h s = 0.05 m) are taken from chute experiment A (Table 1). The mean velocity of the event is U m = 7.30 m s-1.
Fig. 4. The total dissipated work Φ″ as a function of the slip velocity u
0 in chute experiment A. (a) Bingham fluid with viscous sliding law. (b) Dilatant fluid with sliding friction proportional to the velocity squared. The total dissipated work is the sum of the dissipated work in the core 
and the basal sliding surface 
. At steady-state equilibrium the total dissipated work is in balance with the gravitational work rate, 
The graphs show that the sliding dissipation Φ0 increases with increasing sliding velocity u
0: the faster the avalanche flows, the higher the sliding dissipation. 
thus 
and thus the curve represents the total dissipation for fluids of the Voellmy type. Note that for such fluids no dissipation minima exist when the avalanche is flowing. For the Bingham and dilatant models, the internal dissipation function
decreases with increasing sliding velocity. When 
is zero since the avalanche core moves as a plug; the core contains no shear deformations. For u
0 > u
p (slip velocity larger than the plug velocity) the dilatant internal work dissipation function becomes negative, implying negative dissipation, a violation of the second law of thermodynamics.
Because the sliding dissipation function increases with u
0 and the internal dissipation function decreases with u
0, the two mechanisms can be viewed as competing processes. Slight variations in u
0 will increase the work dissipation of one mechanism and decrease the work dissipation of the other. This means that in most cases there exists a velocity u
0 such that the total dissipation Φ″ is a minimum. We denote this velocity (u
0)min and the corresponding dissipation 
If the external constraints on the system allow the avalanche to flow with the sliding velocity (u
0)min, the avalanche will dissipate as little work as possible. For other values of u
0 the avalanche will be flowing with the same mean height h and mean velocity U
m, but will be dissipating more work for the same mass flow.
3.5. Dissipation minima
The graphs of the work dissipation function indicate that dissipation minima for u 0 exist for particular combinations of the constitutive parameters b, s and m. To find the minima we will take the variation of the dissipation function with respect to the discrete internal velocities, u 0 and u p.
The constraint condition is introduced into the variational method by noting that
The constraint condition G is thus
This constraint ensures that the internal velocities u 0 and u p are such that the mean velocity of flow is U m. The mean velocity is known from the experiments. Minimization of the dissipation function without this constraint would give internal velocities that are not in agreement with the observed (or calculated) mass flux.
The variation of Φ″ with respect to u 0 and u p is (λ is the Lagrangian multiplier)
Since the variations δu 0 and δu p are independent, the constrained variational problem consists of the two equations
and the constraint condition (Equation (40)). These three equations must be solved for u 0, u p and λ.
Figure 5 depicts the solution of the variational equations (Equation (42)) for chute experiment A (Table 1) using both constitutive relations. Figure 5a shows the dissipated work using the dilatant model as a function of the sliding velocity u
0. The dissipation functional has a well-defined minimum that corresponds exactly with the gravitational work rate, 
as theoretically required. The minimum value of the slip velocity (u
0)min satisfies the steady-state momentum solution. (Note the asymmetric form of the dissipation functional. At low slip velocities, large shear deformations exist in the flowing core of the avalanche. Since all points on the dissipation curve have the same mass flux, this means that non-equilibrium systems, induced say by a fluctuation in slip velocity, would dissipate over twice as much work as the steady-state solution.) Figure 5b compares the calculated velocity profile in the shear layer with the measurements (Table 1). The constitutive parameters are listed in Table 2. The agreement can be considered good, although the strongly linear character of the solution is disturbing. The linearity arises from the fact that the overburden stresses, which vary linearly in z, are larger than the viscous stresses. Non-linear solutions are possible, but would require reducing the influence of the overburden term in the constitutive formulation.
Fig. 5. Comparison between predicted velocity profiles and measurements for chute experiment A. (a) Dilatant fluid dissipation function showing location of dissipation minimum. (b) Predicted velocity profile using dilatant fluid, comparison to measurements. Constitutive parameters: b = 0.565, m = 0.00014m2, s = 2.5 kgm-3. (c) Comparison between Bingham and dilatant fluid dissipation functions. Both satisfy A = 0 solution. Bingham solution uses viscous sliding law. (d) Predicted velocity profile using Bingham fluid, comparison to measurements. Constitutive parameters: b = 0.57, m = 0.005 m2s-1, s = 12.8kgm-2s-1. Flow parameters: U m = 7.3ms-1, h = 0.4m, hs = 0.05 m, Wg = 6.1 kWm-2.
Table 2. Summary of the parameters found such that 
for the Weissfluhjoch chute experiments. The dissipation functionals for the dilatant fluids are depicted in Figures 5-9
The Bingham solutions presented immediate problems. We found that the momentum and dissipation solutions, Λ = 0 and δΦ = 0, did not match. Although the mass- continuity constraint is exactly fulfilled in both solutions, different velocity profiles are obtained, as depicted in Figure 6. The dissipation solution predicts a much smaller slip velocity and a slightly higher plug velocity than the momentum solution.
Fig. 6. When a Bingham model is used with a velocity-squared sliding friction law, the variational δΦ = 0 solution diverges from the momentum solution Λ = 0. Depicted above is the calculated velocity profile in the fluidized layer of chute experiment A. The mean velocity of both systems is U m = 7.3 ms_1. Thus, the variational solution satisfies the mass-continuity constraint, but not the momentum balance at every height z.
Since the two solutions must provide the same result, we concluded that our Bingham constitutive formulation was in error. Since the simple example problems discussed earlier employed Newtonian laws (n = 1), we believed that the Newtonian part of the formulation was correct. Subsequently, we changed the nature of the sliding friction law from Equation (35) to
Thus, the dependence of the sliding shear stress on a velocity-squared term was removed; a viscous sliding law was employed. This had the effect that the Λ = 0 and δΦ = 0 solutions no longer diverged. Apparently, all possible combinations of sliding friction and constitutive laws do not meet the requirement that Λ = 0 and δΦ = 0, a necessary condition for viscous fluids, according to the theorem of Helmholtz and Kroteweg (Reference LambLamb, 1945). Since this theorem is equally applicable to granular flows, we discarded the Bingham formulation with a velocity-squared sliding friction law from further consideration. The Bingham formulation with the viscous sliding law (Equation (43)) is very capable of reproducing the experimental results. This is shown in Figure 5c and d.
In Figure 7 the dissipation functional using the dilatant law is plotted for chute experiment B in terms of both the slip velocity (Fig. 7a) and plug velocity (Fig. 7c). A comparison between the calculated and measured velocity profiles is provided in Figure 7b. The agreement can be considered good, but the predicted velocity profile appears to be almost linear, once again indicating that the overburden stress influence is too great. Note that the measured and predicted plug velocities are in good agreement.
Fig. 7. Comparison between predicted velocity profiles and measurements for chute experiment B. (a) Dissipation function showing location of dissipation minimum. (b) Velocity profile, comparison with measurements. Constitutive parameters: b = 0.54, m = 0.0002 m2, s = 2.4 kgm_3. Flow parameters: U m = 8.5 ms_1, h = 0.4 m, h s = 0.05 m, W g ″ = 7:1k Wm_2.
The dissipation functional can be drawn in terms of the slip velocity (Fig. 7a) or the plug velocity (Fig. 7c). In both cases, a well-defined minimum value exists, a consequence of the discrete variational formulation which contains two stationary equations for u 0 and u p. Each velocity pair (u 0, u p) satisfies the mass-continuity equation.
For completeness, we present the results of the remaining two experiments, C and D, in Figures 8 and 9. These two experiments are significant since they were performed with lower-density snow. The approximate density of p = 250 kgm-3 is near the flow density of real granular avalanches (usually assumed to be p = 300 kg m-3). In all the experiments the values of viscosity (see Table 2), whether for the Bingham or dilatant model, are in good agreement with the values reported by Reference Dent and LangDent and Lang (1983), as well as Reference Nishimura and MaenoNishimura and Maeno (1987) and Reference NishimuraNishimura (1990). Note that the viscosity values decrease with decreasing snow density, also in agreement with other work (Reference Maeno, Nishimura and KanedaMaeno and others, 1980).
Fig. 8. Comparison between predicted velocity profiles and measurements for chute experiment C. (a) Dissipation function showing location of dissipation minimum. (b) Velocity profile, comparison with measurements. Constitutive parameters: b = 0.565, m = 0.00005 m2, s = 2.5 kgm_3. Flow parameters: U m = 7.3 ms_1, h = 0.4 m, h s = 0.05 m, W g ″ = 3.8 kWm_2.
Fig. 9. Comparison between predicted velocity profiles and measurements for chute experiment D. (a) Dissipation function showing location of dissipation minimum. (b) Velocity profile, comparison with measurements. Constitutive parameters: b = 0.60, m = 0.00002 m2, s = 2.25 kgm_3. Flow parameters: U m = 5.9 ms_1, h = 0.4 m, h s = 0.05 m, W g ″ = 3.1 kWm_2.
4. Fluctuations and Stability
The flowing avalanches we are considering consist of a large number of interacting snow granules. The calculated velocity profiles statistically represent the average motion of the granules over a long time interval. Instantaneous pressure and velocity measurements of real snow avalanches (Reference Dent, Burrell, Schmidt, Louge, Adams and JazbutisDent and others, 1998) show strong deviations from these average or reference values. These deviations can be loosely termed ‘chaotic’;however, we would like to show that these are fluctuations around the average state exhibiting order and structure. For this purpose the properties of the dissipation $ must be examined in more detail. When
the system is in the uninteresting state of static equilibrium where nothing happens. For the case
the dissipative system is in steady state. For brevity, we shall use the term equilibrium to mean dynamic equilibrium.
Two questions must be asked. Firstly, how do the fluctuations originate and secondly, will the fluctuations regress or will the flow system evolve into an alternative state? The latter question is one of stability.
In the case of flowing avalanches, fluctuations arise due to random perturbations of velocity (or, alternatively, kinetic energy). These can be induced by small-scale variations in boundary conditions (slope angle and roughness) or internal material properties (size, strength and density of the granules) or, especially in the case of the chute experiments, from the initial conditions. The manifold ways that deviations from the equilibrium state can be created are opposed to the singularity of the equilibrium state where the fluctuations diminish, given enough time without additional perturbations. In this sense the state of equilibrium, Λ = 0, or equivalently, the state of minimum dissipation, δΦ = 0, acts as an attractor. Consider Figure 10 which depicts the velocity phase space of the flowing avalanches in terms of the slip and plug velocities u 0 and u p, respectively (chute experiment A). There exists only one pair (u 0, u p) that satisfies steady-state equilibrium. The non-equilibrium values lie on a straight line that passes through equilibrium. The values on the line outside equilibrium satisfy the mass- continuity condition (a linear relation), but dissipate more work. Fluctuations in velocity are represented as points on the line; if the fluctuations regress they will follow this linear path to equilibrium.
Fig. 10. Phase space (u 0, u p) of the granular flow system found for non-equilibrium values of the total dissipation Φ″. Equilibrium is an attractor at which ΛΦ = 0. Because of the mass-continuity constraint, the phase space is linear around equilibrium, as well as far from equilibrium. This particular example is calculated for chute experiment A.
The linearity of the phase space (Fig. 10) both near and far from equilibrium has significance for avalanche flow. The linearity implies that the combined action of several deviations from equilibrium can be treated individually since the principle of superposition holds. Even if these deviations are large, we would expect our system to return to steady-state equilibrium, so long as the system constraints (mass flow, slope angle, surface roughness) remain unchanged. We do not want to extrapolate the linearity of the avalanche phase space to cases outside the investigated chute flows, since this would mean that flow state transitions are impossible. Snow avalanches exhibit flow regime transitions, for example, the as yet unexplained transition from a flowing to powder avalanche. The calculated dissipation functions predict that the chute flows are stable and will remain in a single flow regime. No flow bifurcations are possible. This theoretical prediction is very much in agreement with the chute experiments, where no flow transitions were observed. The results, however, cannot be extrapolated outside the chute experiments.
Inspection of all the dissipation diagrams presented in this work (e.g. Fig. 5) reveals that Φ″ has a particular form: a valley bounded by two sides. As stated previously, the lefthand side is formed by work dissipation in the fluidized layer which decreases with increasing slip velocity. The righthand side is formed by work dissipation at the basal layer which increases with increasing slip velocity. The dissipation function will have the same form if it is plotted as a function of the plug velocity (see Fig. 7c). Dissipation landscapes of this form mean that perturbations around the equilibrium state, which rests at the valley bottom, will remain bounded for all values of time (Reference Nicolis and PrigogineNicolis and Prigogine, 1989). Thus, we have found that dissipation in granular avalanches is stable in the sense of Lyapunov (Reference Glansdorf and PrigogineGlansdorf and Prigogine, 1974) with respect to the slip and plug velocities.
Are unstable flows possible? Consider Figure 11 which shows the dissipated work Φ″ as a function of the slip velocity for different material parameters b, m and s. For this particular example we will only consider a dilatant fluid.
Fig. 11. Dilatant fluid. Influence of constitutive parameters b, m and s on the dissipation function Φ″. The total dissipated work Φ″ as a function of the slip velocity u 0. (a) Three different b values with p = 300 kgm_3, m = 0.0001 m2, s = 2.5kgm_3. (b) Three different s values with p = 300 kgm_3, b = 0.4, m = 0.0001 m2. (c) Three different m values with p = 300 kgm_3, s = 2.5 kgm_3, b = 0.4.
Figure 11a depicts the influence of dry Coulomb friction on the dissipation function. We note that the material parameter b does not influence the form of the curve, only the magnitude of the dissipated work. Moreover, dry Coulomb friction will not influence the stability of the system. This is not the case for the remaining two parameters s, representing the surface roughness, and m, the dilatant viscosity. For decreasing s values the righthand side of the dissipation function flattens (Fig. 11b). We expect unbounded oscillations in the slip velocity on smooth surfaces. The onset of this instability is characterized by slip velocities that are higher than the plug velocities. Alternatively, for decreasing m values, the lefthand side of the dissipation function flattens. This is shown in Figure 11c. Note the existence of very flat dissipation functions where both the right- and lefthand valley sides have disappeared. In such cases, any perturbation in dissipation, caused by slight fluctuations in velocity, will have a significant effect on the granular flow system. Physically, this means that the viscous properties of the flow are not strong enough to maintain or recover the steady-state equilibrium of the system.
The fact that unstable systems are possible, depending on the nature of the flow surface and material viscosity, leads us to believe that the equilibrium state we have described for granular snow avalanches is only locally stable. Stability will not prevail for all deviations from equilibrium; flow transitions are destined. Interestingly, the dependence of stability on surface roughness and viscosity is already in good agreement with our empirical understanding of powder avalanche formation: powder avalanches form on long smooth slopes with a low-density (i.e. low-viscosity) snow cover.
5. Temperature and Entropy
So far we have assumed that the snow temperature T is constant. This cannot be strictly true since the dissipated potential energy is directly converted into heat, resulting in a temperature increase ΔT. We shall see that the relative temperature change ΔT/T is small: take an avalanche that descends vertically H = 1000 m and consider 1 kg of snow or ice, M = 1kg. The potential energy is MgH = 104J; the specific heat of ice c v being 2000 J kg K-1, the avalanche cannot be heated up by more than 5 K, or ΔT = T ≈ 2%. This shows that the assumption of a constant temperature in our consideration of minimum dissipation is acceptable. (Thermal effects, such as the local melting of the basal layer, however, could be important.) Since the dissipation is proportional to the entropy production Ṡ by the factor T
we find that this principle fits as a special case the more general statement of minimum entropy production (Reference Glansdorf and PrigogineGlans-dorf and Prigogine, 1974). We assume that the mechanical energy is irreversibly transformed into heat.
The application of the principle of minimum entropy production is subject to severe restrictions. Firstly, the relations between currents (or fluxes) Jk and the generalized forces Xk , defined according to
for m dissipative (entropy-producing) processes, is linear in the forces Xk (Reference PrigoginePrigogine, 1955):
Secondly, the phenomenological coefficients L must obey Onsager’s reciprocity relation (Reference OnsagerOnsager, 1931):
One possible formulation is to divide the entropy production into two parts corresponding to the two frictional layers, the sliding surface (subscript 0) and fluidized layer (subscript f). For a dilatant model, we have (again, the double-prime subscript denotes a quantity per unit area)
with
The phenomenological coefficients then take on the remarkable (and dimensionless) form:
where Fr is the Froude number 
and Saf is the Savage number 
is the diameter of a snow granule, and 
The advantage of this formulation is that it underscores the importance of the internal friction angle b since both coefficients L 00 and L ff are linear in b. Thus, the internal friction angle is similar to the conductivity coefficient in Fourier’s law of heat conduction or the diffusion coefficient in Fick’s law, in that it is a homogeneous, empirical value relating the generalized flows Jk to the forces Xk of our system. This definition of b, as a proper phenomenological coefficient in an irreversible system, is more attractive than the standard definition of the angle of repose, usually defined for a static system at rest. It has long been known that the b values used to match the run-out distances of real avalanches, b ≈ 0.15, have little or no connection with the angle of repose for snow (Reference Buser and FrutigerBuser and Frutiger, 1980). We regard b as having a bulk phenomenological coefficient controlling the work dissipation in granular snow avalanches.
Another feature of this formulation is that it provides information about what limits where the variational principle can be applied. Note that b ≫ (m/d
2)Saf
2 for snow avalanches (assuming 
so the coefficient 
This fact indicates that neglecting the internal deformation, represented by Saf, is not too bad an approximation for run-out distances or flow velocities. A similar argument does not hold for L
00. At high Froude numbers on rough surfaces the magnitude of b is of the same order as 
Unlike the internal work dissipation, sliding friction s cannot be neglected in a first-order model. For the chute experiments of Reference Kern, Tiefenbacher and McElwaineKern and others (2004) (Table 1), we note that 
and we are well within the application range. This indicates that surface roughness might produce instabilities since we are no longer in a regime where superposition is valid.
The phenomenological coefficients b, s and m used might be temperature-dependent, especially near the melting point of snow. We cannot exclude the possibility that in the shear and frictional layers the temperature rise could be such that the values of the coefficients might change. However, this would simply increase the fluctuations, the steady state being unaffected.
In a steady state there is no time dependence. This applies to the minimum entropy production as well. However, considering changes from one steady state to another steady state, both producing the least entropy, it has been conjectured (Reference JaynesJaynes, 1980) that the changes will occur as fast as possible. Therefore, the entropy production per unit time will be a maximum.
6. Conclusions and Outlook
We began our investigations in the hope of finding an additional thermodynamic constraint for studying flow avalanches. We found that a constrained principle of minimum entropy production is not in contradiction with statements of mass and momentum conservation that are usually applied to study the dynamics of snow avalanches. However, the principle additionally provides us with a tool to study flow stability, transition phenomena and the relationship between constitutive equations and irreversible processes.
The state of minimum entropy production is a steady state in dynamic equilibrium because the entropy is steadily increasing in time. This minimum state is a result of two competing dissipational processes: the viscous shear deformation in the fluidized layer and the sliding friction at the basal surface. These two processes also stabilize velocity fluctuations within the flow. We have found that the investigated dissipation functionals for the snow chute avalanches are stable with respect to the slip and plug velocities in the sense of Lyapunov. Thus, we expect stable (non-increasing) fluctuations in velocity at steady state for these particular flows.
From the dissipation functionals, we derived the two generalized forces of our dissipative system (Equation (51)). These are the generalized forces that will restore steady state. The phenomenological coefficients have been written in dimensionless form as a function of the Froude and Savage numbers. The generalized forces and conjugate fluxes are linear with respect to the parameter b which we consider as the primary phenomenological dissipation coefficient of granular avalanches.
We believe that these results do not hold for all avalanche flows. Because snow avalanches exhibit a wide variety of forms, the stable state we have found can only be a local one. However, we now have the means to consider the stability of all flows and the fluctuations within the flow, such that we might be able to find the point when the flowing avalanche ignites into a powder snow avalanche. In future, constitutive models should predict a transition point between flowing and powder-snow avalanches. Perhaps a constitutive model will be found where the generalized forces include the interaction between the competing processes, such that L 12 = L 21 ≠ 0.
The treatment of flowing avalanches by thermodynamical methods like minimum entropy production is not simply another way to come to the same results as with conventional mechanical methods such as mass and momentum balance. It opens the possibility of investigating changes from one steady state to another; that is, when the entropy production is no longer constant in time. If the change of entropy production occurs as fast as possible (under the prescribed boundary conditions), it must be a maximum, namely
Whether avalanches are dissipative systems that obey such a postulate, further investigations will show.
