Snow energy balance
This page is a work in progress. If you have any questions or notice any errors, please raise an issue.
Overview
The prognostic snow energy is the depth-integrated (column) internal energy $\bar{U}_\text{snow}$ (J/m²). Its tendency is computed as,
\[\begin{equation} \frac{\partial \bar{U}_\text{snow}}{\partial t} = Q_\text{base} - Q_\text{top} + Q_\text{precip} + Q_\text{subl}, \end{equation}\]
where $Q_\text{top}$ is the net surface heat flux (the surface energy balance closure flux over the snow), $Q_\text{base}$ is the conductive heat flux at the snow base, $Q_\text{precip}$ is the nondirectional (positive) sensible/latent heat advected by precipitation (see Snow mass balance), and $Q_\text{subl}$ is an advective correction for energy loss due to sublimation. Meltwater drains as liquid water at 0 °C, which is the zero-enthalpy reference of the enthalpy closure below, so it carries no enthalpy out of the snowpack and thus no explicit meltwater energy flux is needed.
The sublimation correction arises from the same enthalpy reference. $Q_\text{top}$ includes the surface latent heat flux which reduces the energy content of the snowpack by $\rho_w L_{sg} E_\text{subl}$, but the the mass leaving the snowpack departs as ice, whose specific enthalpy relative to liquid water at 0 °C is $-L_{sl}$. Adding back $Q_\text{subl} = \rho_w L_{sl} E_\text{subl}$ leaves a net pack loss of $\rho_w L_{lv} E_\text{subl}$ since $L_{sg} = L_{sl} + L_{lv}$ by definition.
Energy–temperature closure
The bulk snowpack is treated as an ice-water-air mixture, and its temperature and liquid water fraction are recovered from the depth-integrated internal energy using the same medium-agnostic FreeWater enthalpy relations used for the soil (see Soil energy balance). The volumetric energy $U_\text{snow} = \bar{U}_\text{snow} / d_\text{snow}$ is computed from the snow depth $d_\text{snow}$ (see compute_volumetric_snow_energy). The FreeWater map references the internal energy to liquid water at 0 °C ($U_\text{snow} = 0$), so the phase-change range is $U_\text{snow} \in [-L\theta, 0]$ with 0 °C ice at $-L\theta$ where $\theta$ is here the total volumetric content of water/ice. The volumetric latent heat is $L\theta = \rho_\text{snow} L_{sl}$ — equivalently $\rho_w L_{sl}$ times the water-equivalent content $\theta = \rho_\text{snow}/\rho_w$ — so the density factor cancels in the depth integral, $L\theta\, d_\text{snow} = \rho_w W_\text{snow} L_{sl}$, leaving the total latent heat set by the snow water equivalent alone. The bulk heat capacity $C_\text{snow} = c_{p,i}\rho_i\theta_\text{ice} + c_{p,w}\rho_w\theta_\text{liq}$ is the ice/liquid volume-weighted sum (see compute_snow_volumetric_heat_capacity); the pore air enters only through the bulk density $\rho_\text{snow}$ and its sensible heat is neglected. Because snow temperature cannot exceed 0 °C, the recovered free-water temperature is clipped at zero.
Terrarium.SnowEnergyTemperatureClosure — Type
struct SnowEnergyTemperatureClosure{NF} <: Terrarium.AbstractEnergyClosureEnergy–temperature closure for snow volumes. For SingleLayerSnow, the depth-averaged snow temperature T_snow (°C) and liquid water fraction θ_liq are recovered from the depth-integrated internal energy Ū_snow (J/m²) using the medium-agnostic FreeWater enthalpy relations, treating the bulk snowpack as an ice-water-air mixture. The internal energy is defined as:
\[U(T) = T_{ ext{snow}} \times C(T) - \rho_{snow} L_{sl} (1 - F(T))\]
with C(T) the temperature-dependent volumetric heat capacity of the snowpack (J/m³/K), ρ_snow L_sl = ρ_w L_sl θ the volumetric latent heat of fusion (J/m³), and F(T) = θ_liq/θ the fraction of the total (liquid water + ice) volumetric water content that is liquid.
Terrarium.compute_volumetric_snow_energy — Function
compute_volumetric_snow_energy(Ū_snow, d_snow) -> Any
Volumetric snow internal energy U_snow = Ū_snow/d_snow [J/m³] from the depth-integrated energy Ū_snow [J/m²] and the snow depth d_snow [m]. The denominator is regularized with a machine-eps offset so the result stays finite (and differentiable) as W_snow → 0: with Ū_snow → 0 and d_snow → 0 together, U_snow → 0. The thin-snow indeterminacy is additionally masked downstream by f_snow → 0. A finite offset is used rather than safediv because safediv returns Inf at exactly d_snow = 0 (even for Ū_snow = 0), which would produce NaN when multiplied by f_snow = 0 in the surface energy balance blend.
Terrarium.compute_snow_volumetric_heat_capacity — Function
compute_snow_volumetric_heat_capacity(
_::Terrarium.AbstractSnow{NF},
constants::PhysicalConstants,
ρ_snow,
liq
) -> Any
Bulk volumetric heat capacity C_snow [J/m³/K] of the snowpack treated as an ice–liquid–air mixture, given the bulk snow density ρ_snow [kg/m³] and the liquid water fraction liq ∈ [0,1] of the water substance. The water-substance mass per unit snow volume is ρ_snow, of which (1 − liq) is ice and liq is liquid; the remaining void space is dry air. The corresponding constituent volume fractions are
\[\begin{aligned} θ_{ice} &= ρ_{snow}·(1 − liq)/ρ_i,\ θ_{liq} &= ρ_{snow}·liq/ρ_w,\ θ_{air} &= 1 − θ_{ice} − θ_{liq} \end{aligned}\]
and the heat capacity is the volume-weighted sum over the ice and liquid constituents, C_snow = cp_i·ρ_i·θ_ice + cp_w·ρ_w·θ_liq. The air's own sensible-heat storage (≈0.1%) is neglected; note the air is still reflected in the bulk density ρ_snow < ρ_ice.
Advected heat from precipitation
Precipitation carries sensible and latent heat into the snowpack, referenced to liquid water at 0 °C (the enthalpy reference of the closure, see Snow energy balance).
Terrarium.compute_snow_precip_heat_flux — Function
compute_snow_precip_heat_flux(
_::Terrarium.AbstractSnow,
constants::PhysicalConstants,
P_s,
R_on_snow,
T_air
) -> Any
Advected heat flux [W/m²] carried into the snowpack by precipitation, relative to liquid water at 0 °C (the U = 0 reference of the FreeWater enthalpy closure). Fresh snow P_s arrives as ice, which sits L_sl below the liquid reference, plus sensible heat for T_air < 0; rain-on-snow R_on_snow arrives as liquid carrying only its sensible heat for T_air > 0. The latent heat released when rain refreezes in a cold pack is captured implicitly by the enthalpy closure, so it is not added here (adding L_sl would double-count relative to the liquid-water reference).
Tendency
Terrarium.compute_snow_energy_tendency — Function
compute_snow_energy_tendency(
i,
j,
grid,
fields,
snow::SingleLayerSnow,
atmos::Terrarium.AbstractAtmosphere,
constants::PhysicalConstants
) -> Any
Depth-integrated snow energy tendency [W/m²] at grid cell i, j (all fluxes positive upward):
dŪ_snow/dt = Q_base − Q_top + Q_precip + Q_sublwhere Q_top/Q_base are the surface/basal heat fluxes, Q_precip the advected precipitation heat (see compute_snow_precip_heat_flux), and Q_subl an advective correction for sublimation.
The sublimation correction Q_subl = ρ_w·L_sl·E_subl is required because the latent heat flux carries the full sublimation enthalpy ρ_w·L_sg·E_subl, whereas the mass leaving the snowpack departs as ice, whose specific enthalpy relative to the liquid-water reference is −L_sl. Adding back ρ_w·L_sl·E_subl leaves the snowpack with a net loss of ρ_w·(L_sg − L_sl)·E_subl = ρ_w·L_lg·E_subl, the vaporization enthalpy carried by the departing vapor.
Note that no explicit meltwater energy term appears because meltwater drains as liquid water at 0 °C, which is the zero-enthalpy reference (U = 0) of the FreeWater closure, so it carries no enthalpy out of the snowpack.
Boundary and coupling heat fluxes
The snow→soil basal conductive flux and the blended soil-top flux, together with the snow surface sublimation rate, are diagnosed after the surface energy balance solve.
Terrarium.compute_snow_basal_heat_flux — Function
compute_snow_basal_heat_flux(
κ_snow,
T_soil,
T_snow,
d_snow
) -> Any
Snow→soil basal conductive heat flux Q_base [W/m²], positive upward (soil → snow). The snowpack is a strong insulator, so only its resistance is retained (snow-resistance-only closure): Q_base = 2·κ_snow·(T_soil − T_snow)/d_snow. The depth is regularized with a machine-eps offset so the flux stays finite as the snowpack vanishes.
Terrarium.compute_snow_soil_heat_flux — Function
compute_snow_soil_heat_flux(
i,
j,
grid,
fields,
snow::SingleLayerSnow,
constants::PhysicalConstants,
soil::Terrarium.AbstractSoil
) -> Any
Blended soil-top heat flux [W/m²] at grid cell i, j: the snow-cover-fraction-weighted combination of the snow→soil basal conductive flux Q_base and the bare-ground surface energy balance closure flux G (ground_heat_flux), f_snow·Q_base + (1 − f_snow)·G. The bulk snow thermal conductivity entering Q_base is recovered lazily from the density scheme rather than stored.
Terrarium.compute_snow_sublimation_flux — Function
compute_snow_sublimation_flux(
i,
j,
grid,
fields,
_::Nothing,
atmos,
constants,
skinT
) -> Any
Snow-surface sublimation rate [m/s SWE] at grid cell i, j, area-weighted by the snow-covered fraction f_snow. The snow surface is treated as saturated: a bulk-aerodynamic vapor flux Δq/rₐ evaluated at the skin temperature, with Δq taken over ice for a sub-freezing surface (the saturation humidity already dispatches over ice for T ≤ 0 — see saturation_specific_humidity_vapor). The water-vapor mass flux ρₐ·Δq/rₐ is converted to a snow-water-equivalent rate via ρ_w. Zero without snow (snow === nothing).
Terrarium.compute_snow_interface_fluxes! — Function
compute_snow_interface_fluxes!(
state,
grid,
_::Nothing,
args...
)
Launch compute_snow_interface_fluxes! to diagnose the snow↔surface/soil coupling fluxes from the snow state and the surface energy balance outputs: the blended soil-top heat flux (soil_heat_flux) and the snow surface sublimation rate (sublimation, the snow-fraction bulk-aerodynamic vapor flux at the converged skin temperature). Must run after the surface energy balance (which sets ground_heat_flux and the skin temperature). No-op when there is no snowpack (snow === nothing).
compute_snow_interface_fluxes!(
out,
i,
j,
grid,
fields,
snow::SingleLayerSnow,
seb::Terrarium.AbstractSurfaceEnergyBalance,
soil::Terrarium.AbstractSoil,
constants::PhysicalConstants,
atmos::Terrarium.AbstractAtmosphere
)
Diagnose the snow↔surface/soil coupling fluxes at grid cell i, j, run after the surface energy balance: the blended soil-top heat flux (see compute_snow_soil_heat_flux) and the snow surface sublimation rate (see compute_snow_sublimation_flux, the same snow-fraction vapor flux the surface energy balance uses for the latent-flux partition).