During the Arctic winter, the conductive heat flux through the sea ice and snow balances the radiative and turbulent heat fluxes at the surface. Snow on sea ice is a thermal insulator that reduces the magnitude of the conductive flux. The thermal conductivity of snow, that is, how readily energy is conducted, is known to vary significantly in time and space from observations, but most forecast and climate models use a constant value. This work begins with a demonstration of the importance of snow thermal conductivity in a regional coupled forecast model. Varying snow thermal conductivity impacts the magnitudes of all surface fluxes, not just conduction, and their responses to atmospheric forcing. Given the importance of snow thermal conductivity in models, we use observations from sea ice mass balance buoys installed during the Multidisciplinary drifting Observatory for the Study of Arctic Climate expedition to derive the profiles of thermal conductivity, density, and conductive flux. From 13 sites, median snow thermal conductivity ranges from 0.33 W m−1 K−1 to 0.47 W m−1 K−1 with a median from all data of 0.39 W m−1 K−1 from October to February. In terms of surface energy budget closure, estimated conductive fluxes are generally smaller than the net atmospheric flux by as much as 20 W m−2, but the average residual during winter is −6 W m−2, which is within the uncertainties. The spatial variability of conductive heat flux is highest during clear and cold time periods. Higher surface temperature, which often occurs during cloudy conditions, and thicker snowpacks reduce temporal and spatial variability. These relationships are compared between observations and the coupled forecast model, emphasizing both the importance and challenge of describing thermodynamic parameters of snow cover for modeling the Arctic as a coupled system.

Snow and sea ice are the important components of the Arctic system, regulating the exchange of heat and moisture between the ocean and atmosphere. In order for sea ice to grow, energy must be conducted away from the bottom of the ice for the underlying sea water to cool and freeze. However, snow is a thermal insulator, and its presence on sea ice reduces sea ice growth by slowing energy transport at the surface with an impact that increases as sea ice thickness decreases (e.g., Maykut, 1978; Sturm et al., 2002).

Quantifying the exchanges of energy at the surface is crucial not only for understanding how the Arctic evolves on seasonal to interannual timescales but also for correctly representing the Arctic system in models. In winter over sea ice, net longwave radiation (LWNET), and turbulent heat fluxes are balanced at the surface by heat conduction through the sea ice and snow. This conductive flux represents the loss of energy at the sea ice bottom in order for sea ice to grow. Conduction is proportional to the temperature gradient through ice and snow, which is a function of sea ice and snow thicknesses and the temperature difference between the surface and underlying ocean. Because conduction depends on the surface temperature, it couples atmospheric variability to sea ice growth. As the Arctic continues to warm and sea ice continues to thin, the magnitude and variability of conductive fluxes are likely to change (e.g., Bigdeli et al., 2020; Landrum and Holland, 2022).

Conductive flux also depends on the thermal conductivity of the medium, or how readily energy is transported via conduction. The measurements of snow thermal conductivity range from less than 0.1 W m−1 K−1 for fresh snow to 0.5 W m−1 K−1 or higher for wind slabs from field observations and laboratory measurements (e.g., Sturm et al., 2002; Calonne et al., 2011). This range reflects differences in snow properties, including density, grain structure, and bonding (e.g., Brandt and Warren, 1997). To further complicate matters, these snow properties are known to vary vertically, spatially, and temporally as the snowpack undergoes metamorphism (e.g., Colbeck, 1982). Due to density considerations, all values for snow are less than the thermal conductivity of sea ice, which is around 2 W m−1 K−1, although this too can vary depending on sea ice properties, such as density, porosity, and salinity (e.g., Pringle et al., 2007; Petrich and Eicken, 2010). The spatial variability of snow thermal properties and depth can lead to large variability in conductive flux. At the Surface Heat Budget of the Arctic Ocean project, Sturm et al. (2001) found 4-fold variations in conductive flux across distances of tens of meters due to differences in snow and ice thicknesses and properties.

Despite the large and influential variability of snow properties, the representation of snow on sea ice is relatively simple in most forecast and global climate models as well as reanalyses. Typically, these models treat snow thermal conductivity as a constant for snow on sea ice, commonly 0.3 W m−1 K−1 (e.g., Hunke et al., 2017; Ridley et al., 2018). Not all forecast models even include snow on sea ice (e.g., Bauer et al., 2016; Bengtsson et al., 2017). This exclusion has been related to surface temperature biases, which is important to consider given the common practice of assessing models and forcing models with reanalyses (Batrak and Müller, 2019). The complex modeling of snow is possible (e.g., Liston et al., 2020), but it is typically only done on land because the transient nature of sea ice makes it computationally more expensive to trace snow properties on sea ice that is advected and undergoes strong seasonal evolution (Hunke et al., 2010).

Correctly modeling the variability of conductive fluxes is important for simulating meaningful representations of sea ice thickness and, therefore, sea ice variability. Urrego-Blanco et al. (2016) found that in the Los Alamos Sea Ice Model (CICE), sea ice was more sensitive to the representation of snow than any other parameter, specifically snow thermal conductivity and grain size. Lecomte et al. (2013) tested several parameterizations of snow thermal conductivity in a coupled sea ice-ocean model, and both the mean sea ice state and variability were impacted by the representation of snow thermal conductivity. Wu et al. (1999) decreased snow thermal conductivity from 0.3 to 0.16 W m−1 K−1 in a “zero-layer” coupled atmosphere-sea ice model. This change led to thinner sea ice distributions for both poles and unexpectedly warmer surface temperatures in the Arctic. They explained the warmer surface temperatures as being influenced by changes in sea ice thickness more so than snow thermal conductivity. Blazey et al. (2013) found that increasing snow thermal conductivity in a fully coupled global climate model increased sea ice extent in the Arctic and decreased pan-arctic surface temperatures. Low-level cloud cover also decreased, consistent with the fall sea ice cloud feedback (Kay and Gettelman, 2009), and suggesting that the representation of snow on sea ice is also important for modeling the atmosphere. Yet despite the known importance of snow on sea ice, and in particular its thermal conductivity, better representation is still needed in models (e.g., Batrak and Müller, 2019).

Observations from the Multidisciplinary drifting Observatory for the Study of Arctic Climate (MOSAiC) expedition provide a unique opportunity to study snow as a part of the coupled ice-snow-atmosphere Arctic system. During MOSAiC, the German icebreaker Polarstern was frozen into the sea ice and observations were taken across a regional domain beginning in October 2019 and continuing to the following September (Nicolaus et al., 2022; Shupe et al., 2022). MOSAiC provides some of the most comprehensive observations from the Central Arctic during winter, including coordinated observations of sea ice, snow, and the surface energy budget (SEB).

The goal of this work is to study the variability of snow thermal conductivity on sea ice and conductive flux during the MOSAiC winter. For context, we first present a sensitivity study using a regional coupled forecast model to demonstrate the impact of varying snow thermal conductivity on the SEB. Next, we introduce the observations and present a method for retrieving profiles of snow thermal conductivity and density that are used to calculate conductive flux. Results are then discussed for spatial and temporal variability of snow thermal conductivity and conductive flux through the snowpack. Finally, relationships between temperature, sea ice, and conductive flux from observations are used to further assess the coupled forecast model.

Most modeling studies have focused on the sensitivity of sea ice to changing snow thermal conductivity, but we also expect the adjacent atmosphere to respond since conduction couples the sea ice to the overlying air via the surface skin temperature. Here, we use the Coupled Arctic Forecast System (CAFS; Solomon et al., 2023), a fully coupled regional forecast model, to demonstrate the impact of snow thermal conductivity on the winter SEB over sea ice. This example simply serves as a context for the importance of snow thermal conductivity for modeling the near-surface boundary layer, while more comprehensive simulations and sensitivity studies are warranted but beyond the scope of this study. CAFS uses the following components: CICE version 5 with 3 snow layers and 7 ice layers, WRF ARW version 3.6.1, POP ocean model version 2, and the NCAR CLM land model version 4.5. Global 0.5° GFS forecast fields are used for lateral forcing. Output data are 6-h averages with 10-km resolution over the Arctic basin domain, shown in Figure S1. In addition to a control simulation with a snow thermal conductivity of 0.3 W m−1 K−1, 2 additional simulations were conducted with the values of 0.15 W m−1 K−1 and 0.6 W m−1 K−1 changed over the 16 grid cells poleward of 89.78°N. For each snow thermal conductivity value, a single 10-day simulation is analyzed for January 30 through February 8, 2020. After February 8, the large-scale impact of the varied snow thermal conductivity, for example, shifts in cyclone tracks, makes isolating the impacts on the local SEB difficult. However, the chosen week has both clear and cloudy conditions that uniquely respond to changes in snow thermal conductivity. Average snow depth increases from 17.5 to 20.5 cm during the week, and sea ice thickness is 2.5 m. Only grid cells from the domain where snow thermal conductivity is changed are analyzed.

Snow thermal conductivity determines how efficiently energy is transported through the snow to the snow–atmosphere interface. The domain-average conductive flux increases by 7.3 W m−2, or about 34%, when snow thermal conductivity is doubled to 0.6 W m−1 K−1, and it decreases by 6.6 W m−2, or about 30%, when snow thermal conductivity is halved to 0.15 W m−1 K−1 (Figure 1a and c). Although conductive flux is linearly proportional to thermal conductivity, the conductive flux is not twice as large after doubling snow thermal conductivity. A larger conductive flux means more energy flows to the surface and the surface temperature is higher, on average by 0.85°C when snow thermal conductivity is doubled from 0.3 to 0.6 W m−1 K−1 (Figure 1b and d). This increase in surface temperature acts as a negative feedback since warmer surface temperatures decrease conduction, mitigating some of the impact of larger snow thermal conductivity. The opposite is also true for halving snow thermal conductivity. The surface temperature decreases by an average of 0.8°C, and the conductive flux decreases by less than half.

Figure 1.

Conductive flux (a) and surface temperature (b) time series for 3 simulations with different values of snow thermal conductivity (ks) in the Coupled Arctic Forecast System (CAFS). Panels (c) and (d) show the changes of the given parameter for one of the perturbation simulations relative to the control (ks = 0.3 W m−1 K−1) simulation. Each curve represents one of the 16 grid cells, a 30-km square region poleward of 89.78°N, where snow thermal conductivity was changed. Conductive flux is positive downward.

Figure 1.

Conductive flux (a) and surface temperature (b) time series for 3 simulations with different values of snow thermal conductivity (ks) in the Coupled Arctic Forecast System (CAFS). Panels (c) and (d) show the changes of the given parameter for one of the perturbation simulations relative to the control (ks = 0.3 W m−1 K−1) simulation. Each curve represents one of the 16 grid cells, a 30-km square region poleward of 89.78°N, where snow thermal conductivity was changed. Conductive flux is positive downward.

Close modal

The impacts of changing snow thermal conductivity depend on atmospheric conditions. When the surface is warm, the temperature difference between the surface and ocean is smaller, so less conduction occurs. When the surface is coldest, the temperature gradient is relatively large and conduction is greater. This basic relationship is how the atmosphere can impact sea ice growth. When the atmosphere is cloudy, there is more downwelling longwave radiation (LWD) to insulate the surface, such as February 1–2, 2020. Accordingly, the snow thermal conductivity has comparatively less impact on the conductive flux during this time than when the surface is colder. The largest impacts of changing snow thermal conductivity occur when the atmosphere is clear and the surface is coldest, for example, January 30, 2020. Because of the different responses of conduction when the surface is warm or cold, that is, the atmosphere is cloudy or clear, changing snow thermal conductivity also alters the relationship of conduction to LWD (Figure 2a). Note that there are no changes to LWD between simulations. The response, or sensitivity, of conduction to LWD, given by the slope of the best fit line in Figure 2a, increases from 0.21 to 0.38 (an 81% increase) when snow thermal conductivity increases from 0.15 W m−1 K−1 to 0.6 W m−1 K−1. This difference means that the coupling between the surface and atmosphere changes as a result of snow properties.

Figure 2.

Relationships between downwelling longwave (LWD) and conductive flux (a), upwelling longwave (LWU) (b), and sensible heat flux (SH) (c) for 3 runs with different values of snow thermal conductivity (ks) in Coupled Arctic Forecast System (CAFS). Panel (d) shows the altered relationship between net longwave (LWNET) and SH. Fluxes are positive downward. Lines of best fit are calculated using linear least-squares regressions for ks = 0.15 W m−1 K−1 (dashed line), ks = 0.3 W m−1 K−1 (solid line), and ks = 0.6 W m−1 K−1 (dashed-dotted line).

Figure 2.

Relationships between downwelling longwave (LWD) and conductive flux (a), upwelling longwave (LWU) (b), and sensible heat flux (SH) (c) for 3 runs with different values of snow thermal conductivity (ks) in Coupled Arctic Forecast System (CAFS). Panel (d) shows the altered relationship between net longwave (LWNET) and SH. Fluxes are positive downward. Lines of best fit are calculated using linear least-squares regressions for ks = 0.15 W m−1 K−1 (dashed line), ks = 0.3 W m−1 K−1 (solid line), and ks = 0.6 W m−1 K−1 (dashed-dotted line).

Close modal

Because changing snow thermal conductivity alters the surface temperature, upwelling longwave radiation (LWU) and turbulent fluxes also change between simulations. Larger snow thermal conductivity and surface temperature lead to larger LWU, by up to 9 W m−2, and larger upward sensible heat (SH) flux, by up to 15 W m−2. Changes to latent heat flux are relatively modest given their already small values (not shown). Again, the sensitivities of LWU and SH to LWD differ between the simulations. As conduction becomes more sensitive to the atmosphere with larger snow thermal conductivity, LWU and SH become less sensitive to LWD (Figure 2b and c). The largest changes in LWU and SH occur when the atmosphere is clear and the surface is colder. With larger snow thermal conductivity, the surface is warmer than in simulations with smaller snow thermal conductivity, so there is more LWU (i.e., Boltzmann’s law). The surface temperatures change less in response to increased LWD when the surface is already warm, so the slope between LWU and LWD decreases from −0.48 to −0.41 when snow thermal conductivity is increased from 0.15 W m−1 K−1 to 0.6 W m−1 K−1 (blue circles to red diamonds in Figure 2b). At the same time, because the surface is warmer with higher snow thermal conductivity, the near surface temperature gradient is reduced during clear and cold periods. Therefore, SH is smaller and its sensitivity to LWD is reduced from −0.29 to −0.19 when changing snow thermal conductivity from 0.15 W m−1 K−1 to 0.6 W m−1 K−1 (Figure 2c). Overall, with higher snow thermal conductivity, conduction becomes more responsive to atmospheric radiative forcing, while SH and LWU become less responsive.

The relationship between LWNET and SH also changes between these simulations, notably when the atmosphere is radiatively clear (LWNET ≤ −20 W m−2; Figure 2d). With the lowest snow thermal conductivity (blue circles), LWNET is less negative because the colder surface results in less LWU. SH is more positive because of the larger temperature difference between the atmosphere and surface. Thus, the thermal properties of snow in a coupled model are important not only for representing sea ice but also for the lower atmospheric structure and the interactions among different heat transfer processes.

This experiment presents a practical demonstration for the need to constrain the mean and variability of the thermal properties of snow in coupled sea ice models. This will be the objective of the remainder of this study using observations made across the MOSAiC distributed network (DN).

3.1. Observations

For the MOSAiC expedition, Polarstern was frozen into an ice floe on October 4, 2019, in the Central Arctic Ocean. We use observations during winter from October 2019 until March 1, 2020, during which time the Polarstern passively drifted with the surrounding sea ice. The Central Observatory (CO) consisted of a suite of instruments on the ice floe within 2 km of the Polarstern. We use broadband LWU and LWD measured at 3 m and 1 m above the surface, respectively, from a station located within a few tens of meters from the 10-m tower. Additional instruments were deployed as part of the DN (Krumpen and Sokolov, 2020) up to tens of kilometers from the CO. At the large sensor sites (L-sites), Atmospheric Surface Flux Stations (ASFSs) were deployed to measure the components of the SEB and meteorological variables (Cox et al., 2023). Broadband LWU and LWD were measured using pyrgeometers mounted at 2 m above the surface, with estimated uncertainty of 2.6 W m−2 and 1 W m−2 for LWD and LWU, respectively. Surface skin temperature was derived from longwave fluxes using an assumed surface emissivity of 0.985, which is appropriate for winter Arctic snow conditions. Two conductive heat flux plates were installed within the snowpack within several meters of each ASFS, providing estimates of the conductive flux at a single depth near the sea ice–snow interface. Three-dimensional wind speed was measured by ultrasonic anemometers, while turbulent fluxes were calculated from 10-Hz anemometer data using the eddy covariance technique. SH estimates have an estimated uncertainty of 4 W m−2. All data used here from the ASFSs and tower have 10-min resolution.

As context for the analyses of snow and ice thermal conductivity and its interaction with the SEB, cloud measurements made from onboard Polarstern are also included. Radar reflectivity is derived from the vertically pointing Ka-band Atmospheric Radiation Measurement (ARM) Zenith Cloud Radar, operated by the U.S. Department of Energy ARM Program (Lindenmaier et al., n.d.).

We use data from 3 seasonal ice mass balance (SIMB) buoys deployed at the DN’s large L-sites (Perovich et al., 2022; Perovich et al., 2023). These SIMBs are spar buoys that included a temperature chain with 192 thermistors at 2-cm spacing and acoustic rangefinders above and below the ice that measure the top and bottom of the snow-ice system, used to determine the ice thickness (Polashenski et al., 2011; Planck et al., 2019). Data were collected at 4-h intervals. We use temperature and sea ice thickness measurements from the SIMB. At installation in October, sea ice thickness at the 3 SIMB sites ranged from 0.30 m to 1.20 m (Table 1), and snow depths are in the range of approximately 0.06–0.12 m. At installation, 0 cm refers to the nominal sea ice–snow interface.

Table 1.

Ice mass balance buoys (IMBs) used in the study

IMB NumberStart DateEnd DateBeginning Ice Thickness [m]Ending Ice Thickness [m]Nearest Site
SIMB L1a November 1, 2019 March 1, 2020 1.17 1.67 L1 
SIMB L2a November 1, 2019 March 1, 2020 0.85 1.59 L2 
SIMB L3a December 15, 2019 January 31, 2020 0.92 1.32 L3 
SIMBA 2019T47 October 23, 2019 February 7, 2020 0.42 1.24 P64 
SIMBA 2019T58 October 17, 2019 March 1, 2020 0.82 1.48 M4 
SIMBA 2019T62 October 30, 2019 March 1, 2020 0.80 1.48 Second year ice coring site 
SIMBA 2019T63a October 12, 2019 March 1, 2020 1.12 1.56 L2 
SIMBA 2019T64 October 12, 2019 March 1, 2020 1.72 2.08 M6 
SIMBA 2019T65a October 15, 2019 March 1, 2020 1.24 1.54 L2 
SIMBA 2019T67a October 22, 2019 February 28, 2020 1.36 1.78 L1 
SIMBA 2019T68 October 26, 2019 March 1, 2020 1.71 2.01 M1 
SIMBA 2019T70a October 11, 2019 March 1, 2020 0.48 1.52 L3 
SIMBA 2019T72 October 10, 2019 March 1, 2020 1.00 1.74 M5 
IMB NumberStart DateEnd DateBeginning Ice Thickness [m]Ending Ice Thickness [m]Nearest Site
SIMB L1a November 1, 2019 March 1, 2020 1.17 1.67 L1 
SIMB L2a November 1, 2019 March 1, 2020 0.85 1.59 L2 
SIMB L3a December 15, 2019 January 31, 2020 0.92 1.32 L3 
SIMBA 2019T47 October 23, 2019 February 7, 2020 0.42 1.24 P64 
SIMBA 2019T58 October 17, 2019 March 1, 2020 0.82 1.48 M4 
SIMBA 2019T62 October 30, 2019 March 1, 2020 0.80 1.48 Second year ice coring site 
SIMBA 2019T63a October 12, 2019 March 1, 2020 1.12 1.56 L2 
SIMBA 2019T64 October 12, 2019 March 1, 2020 1.72 2.08 M6 
SIMBA 2019T65a October 15, 2019 March 1, 2020 1.24 1.54 L2 
SIMBA 2019T67a October 22, 2019 February 28, 2020 1.36 1.78 L1 
SIMBA 2019T68 October 26, 2019 March 1, 2020 1.71 2.01 M1 
SIMBA 2019T70a October 11, 2019 March 1, 2020 0.48 1.52 L3 
SIMBA 2019T72 October 10, 2019 March 1, 2020 1.00 1.74 M5 

aIndicates IMBs that are included in Section 4.2.

Additional snow and sea ice temperatures are obtained from Snow and Ice Mass Balance Arrays (SIMBA), which are thermistor string-type IMBs (Jackson et al., 2013). We use 10 SIMBAs installed during October at the CO and as part of the DN at the L-sites, medium instrumentation sites (M-sites), and at position nodes (P-sites; Lei et al., 2022). SIMBAs consist of a 4.8-m-long thermistor chain with sensors at 2-cm intervals that are sampled at 6-h intervals. SIMBAs were installed vertically through the sea ice and snow, extending from the ocean to the atmosphere. In addition to temperature profiles, SIMBAs measure temperature differences after pulsative heating over 120 s. The ratio of temperature differences between 30 and 120 s after the heating is used to distinguish the sea ice–water interface and calculate sea ice thickness. At installation, sea ice thicknesses ranged from 0.34 m to 1.81 m, and snow depths ranged from 0.06 m to 0.18 m (Table 1). As with the SIMB, 0 cm refers to the sea ice–snow interface at installation. Locations of IMBs (generally referring to all SIMB and SIMBA), the CO, and ASFSs are shown in Figure 3.

Figure 3.

Locations of ice mass balance buoys (IMBs) and flux measurements on January 1, 2020. Shapes correspond to snow depth (s.d.), and shading represents sea ice thickness. Solid black markers are atmospheric flux measurements.

Figure 3.

Locations of ice mass balance buoys (IMBs) and flux measurements on January 1, 2020. Shapes correspond to snow depth (s.d.), and shading represents sea ice thickness. Solid black markers are atmospheric flux measurements.

Close modal

To prevent the artificial gradient produced by the calibration offsets between thermistors from propagating into the calculations of thermal conductivity, temperature profiles are smoothed at each time step using a running mean window of 5 levels below 0 cm and 3 levels above. This smoothing is done for both SIMBA and SIMB data used in deriving thermal conductivity and density profiles, discussed next.

3.2. Deriving thermal conductivity

To derive thermal conductivity and conductive flux from the temperature profiles measured by IMBs, we solve the one-dimensional heat diffusion equation largely based on the work of Lipscomb (1998). Our derivation further provides the estimates of density profiles. Temperature profiles have been previously combined with forward models and additional constraints to solve for thermal conductivity. Brandt and Warren (1997) and Oldroyd et al. (2013) solved the inverse problem to derive thermal diffusivity in a seasonal snowpack and thermal conductivity at the South Pole. Marchenko et al. (2019) used a forward-model constrained by density measurements to derive the profiles of thermal conductivity on Svalbard.

For this derivation, we assume energy is only transported via vertical conduction through the sea ice-snow column, described by a one-dimensional heat diffusion equation:

tcpT)=z(kTz),
1

where ρ is the density, cp is the heat capacity, k is the thermal conductivity, z is the vertical spacing, t is the time, and T is the temperature. The left-hand side of this equation represents storage, and the right-hand side is the vertical divergence of conduction. Equation 1 conserves energy and explicitly states that any divergence of the vertical conductive heat flux must manifest itself as a temperature change at that depth. It is appropriate for nonmelting conditions. Equation 1 can be approximated with finite differencing at time step m and level n (n increasing downward, i.e., at the top level n = 0) as

ρcpΔt[Tnm+1Tnm]=1Δz(kn+1Δz[Tn+1m+1Tnm+1]knΔz[Tnm+1Tn1m+1]).
2

This approximation assumes that density can change with height but does not vary between individual adjacent time steps. Heat capacity is assumed to only depend on temperature (Cuffey and Paterson, 2010). Equation 2 is written for each level between the top of the snowpack (n = 0) to within the sea ice and rearranged in matrix form to solve for k. For each time step, the lowest level N, where Equation 2 is written, is varied from −100 cm to −30 cm, with respect to the nominal sea ice–snow interface from installation, at every 2-cm interval, that is, 36 profiles if sea ice is thicker than 120 cm. At the lowest level N, kice is assumed to be 2 W m−1 K−1. We vary the depth of the lowest level because the thermal conductivity at a particular level may not be exactly 2 W m−1 K−1, but averaged over the bulk of the ice, it should be representative; 2 W m−1 K−1 is within the measured range of core sea ice thermal conductivity from Trodahl et al. (2001) and Pringle et al. (2007) for multiyear ice (about 1.8–2.3 W m−1 K−1). Measured thermal conductivity for first year ice is larger, with profile averages of 2.06–2.26 W m−1 K−1 (Pringle et al., 2007) and subsurface averages of 2.06–2.19 W m−1 K−1 (Pringle et al., 2006). If the sea ice is thinner than 120 cm, the lowest depth is set to 20 cm less than the sea ice thickness, so that knowledge of fluxes at the ocean–ice interface is not needed. We also exclude the bottom of the sea ice because it can have different properties than the bulk interior, such as a weak skeletal layer (Weeks, 2010) or brine overturning (Notz et al., 2005).

For each time step and lowest level, initially density is assumed to be constant in the sea ice and snow with values of 900 kg m−1 and 300 kg m−3, respectively. The generalized inverse of the thermal conductivity matrix is calculated using its singular-value decomposition (SVD) and including all large singular values with the Python function numpy.linalg.pinv(). If there is no solution to the inverse problem, that is, the determinant is zero, the profile is flagged for removal later in the derivation.

After a first profile of k is derived, density is calculated using Equation 5 from Calonne et al. (2019) for densities below 450 kg m−3. For densities above 800 kg m−3, Equation 4 from Pringle et al. (2007) is used assuming salinity is 5 ppt. Varying salinity had no measurable impact on snow thermal conductivity. For densities between 450 and 800 kg m−3, the following weighted average is used:

k=(1ϕ)*kCalonne[ρ,T]+ϕ*kPringle[ρ,T,S],
3

where kCalonne, thermal conductivity from Calonne et al. (2019), and kPringle, thermal conductivity from Pringle et al. (2007), are both functions of density. In Equation 3, ϕ=(1+exp(2*a*(ρ450)))1, and a = 0.02 m3 kg. Alternative density-conductivity relationships, for example, Yen (1981) and Macfarlane et al. (2023b), lead to differences in snow density estimates ≤30 kg m−3 or 10%.

Using the new density profile calculated from k, the generalized inverse of the thermal conductivity matrix is again solved. A density profile is again calculated from the new k profile, and this process is iterated 5 times in total. Profiles typically converge after 1 or 2 iterations, and profiles that do not converge tend to be excluded from the final results because at least 1 profile in the timestep has a determinant equal to zero. The final profiles of thermal conductivity and density for each time step are the averages over the 36 profiles calculated by starting at different depths. For time steps when there is at least 1 profile with no solution for the inverse problem, that is, the determinant is zero, thermal conductivity and density are linearly interpolated temporally over each level. A schematic of these steps can be found in Figure S2.

The same inversion method is used to derive thermal conductivity and density profiles in the bottom part of the sea ice column using an upper boundary condition. Again, the upper boundary is set to kice = 2 W m−1 K−1 for each of the multiple calculations with starting depths ranging from −30 to −100 cm. Thermal conductivity profiles are solved from the upper boundary to the first level above the sea ice bottom, and these profiles are appended below the profiles for the ice and snow that were solved using a lower boundary condition. At some time steps, thermal conductivity increases within several centimeters of the ice bottom above realistic values (>3 W m−1 K−1) (Figure S3). These large values could be due to the assumption of constant salinity, the sea ice density-thermal conductivity equation not being representative, or other fluxes (e.g., ocean heat flux) that are not included in Equation 1. Because we are primarily interested in snow and near surface fluxes, further refinement of this method near the sea ice bottom is left for future work.

Finally, conductive flux, FC, is calculated with unsmoothed temperature profiles:

FC=kTz.
4

In the end, each profile of temperature measurements results in profiles of effective thermal conductivity, density, and conductive heat flux, all of which are related to each other via the described relationships.

Key interfaces within the vertical profile, that is, the atmosphere–snow and snow–ice interfaces, can manifest in the measurements in different ways. In particular, due to differences in thermal conductivity and density, temporal and vertical variations in temperature are markedly different for thermistors that are in different media. Additionally, the SIMB systems use sonic rangers to identify the top of the snow, but the SIMBA systems do not have these ranger measurements. To maintain consistent definitions of interface depths across different measurement systems, we identify the top of the snowpack by where the derived thermal conductivity increases beyond 0.6 W m−1 K−1. This threshold generally agrees well with the change in vertical temperature gradient between snow and air. Determining the interface between snow and sea ice is complicated by the presence of refrozen surface scattering layers or transition layers that have density and thermal conductivity between those expected of snow and sea ice. Because the snow–ice interface is not clearly distinguished for many IMBs, we determine the bottom of the snowpack as where the mean vertical profile of thermal conductivity is less than or equal to 0.6 W m−1 K−1, discussed further in Section 4.1.

3.3. Sources of uncertainty

There are various sources of uncertainty in these derivations of snow thermal properties, ranging from measurement uncertainties to the assumptions within the derivation. Individual thermistor measurements in thermistor strings have an uncertainty of ±0.5°C to −40°C ambient and ±0.25°C to −20°C ambient for SIMB (Planck et al., 2019) and ±0.1°C applicable to the full temperature range for SIMBA (Jackson et al., 2013; Lei et al., 2022). These uncertainties partly manifest as offsets for a given thermistor relative to others in the string. These offsets are of a similar magnitude as the true vertical temperature gradient and its second derivative through the sea ice; hence why we apply the vertical smoothing of temperature profiles when deriving thermal conductivity. Temperature gradients within snow are typically an order of magnitude larger than those in sea ice when the snow is not isothermal, so possible offsets are generally smaller than real temperature differences.

Due to variations in salinity and, to a lesser extent, temperature, the assumption of kice = 2 W m−1 K−1 at the lower boundary may not be true at any particular level or time step, and the relative error in this assumption propagates linearly into the derived k values throughout the profile. This assumption impacts the magnitude of derived thermal conductivities in snow but not the relative vertical variability. Typically, derived thermal conductivity profiles for different starting depths at the same time step have the same vertical profile but can be shifted relative to each other (Figure 4a). The standard deviations between profiles derived with different starting depths have median values of 0.021–0.052 W m−1 K−1 in the snowpack (Figure 4c). Varying the level of the lower boundary condition and taking the average at each time step is meant to minimize the impact of this issue, in the sense that if the bulk thermal conductivity of sea ice is 2 W m−1 K−1, the thermal conductivity in the sea ice from the different solved profiles averages to approximately 2 W m−1 K−1.

Figure 4.

Variability of thermal conductivity from varying lower boundary condition depth (a and c) and temporal variability in 1-day windows (b and d). Half-violin plots show distributions of standard deviations from varying the lower boundary (c) and subdaily variability (d) for only levels above the snow–ice interface.

Figure 4.

Variability of thermal conductivity from varying lower boundary condition depth (a and c) and temporal variability in 1-day windows (b and d). Half-violin plots show distributions of standard deviations from varying the lower boundary (c) and subdaily variability (d) for only levels above the snow–ice interface.

Close modal

Subdaily variability of thermal conductivity in the snowpack is shown in Figure 4b. Standard deviations are calculated over each day, that is, across 6 or 4 consecutive time steps for SIMB and SIMBA, respectively, for 6 of the IMBs. Standard deviations range from near zero up to 0.18 W m−1 K−1 with medians around 0.02–0.06 W m−1 K−1 (Figure 4d). This variability is larger from SIMB instruments than SIMBA instruments (Figure S4). We believe this variability is primarily due to different scales between the temporal and vertical resolutions, meaning that there could be significant temperature changes between the time steps that are not adequately sampled by IMBs causing unrealistic variability. For example, several IMBs show temporary changes in thermal conductivity on the order of ±0.7 W m−1 K−1 around November 17, 2019, during a cyclone (Rinke et al., 2021). During this period, ASFS skin temperatures fluctuated by roughly ±10°C within hours. Neither 4 or 6 hourly observations fully capture these temperature fluctuations and their thermal impacts on the snow and sea ice.

It is difficult, maybe impossible, to confidently determine where exactly the snow surface is based on temperature profiles alone (Jackson et al., 2013), although various algorithms try (e.g., Zuo et al., 2018). For this reason, we focus on averages in the snowpack and “near-surface” values rather than values directly at the snow–atmosphere interface. This issue is a significant problem for estimating the surface conductive flux for SEB closure since it can impact the effective surface skin temperature that is consistent with the IMB measurements. Additionally, any potential biases in the IMBs, including those that are temperature dependent, can impact the surface temperature estimated from IMBs compared to that measured by other means. Comparing 7 IMBs closest to ASFS at the L-sites, temperatures at the highest level of resolved profiles (i.e., at the determined snow–atmosphere interface) are typically warmer than ASFS-measured skin temperatures (Figure 5a). Mean differences are +0.7°C–1.2°C; 13% of the data have differences larger than +2°C, and 1% of the data have differences larger than +4°C. Only 7% of data have colder temperatures than ASFS skin temperatures. While horizontal variability in snow-ice properties between the specific sampling sites certainly leads to some of the spread in differences, the statistics of these differences indicate that the IMB temperatures at the implied snow–atmosphere interface are on average warm relative to independent measures of the skin temperature. This apparent warm bias could be due to misidentification of the snow–atmosphere interface as well as a systematic bias in the IMB measurements that appears to be temperature dependent (Figure 5).

Figure 5.

The distribution (a) and temperature dependence (b) of the difference between the skin temperature measured by Atmospheric Surface Flux Stations (ASFSs) and the best-estimated surface temperature from 7 ice mass balance buoys (IMBs) located near the ASFS. The impact of these different surface temperatures on a bulk calculation of surface conductive flux is also shown as a distribution (c) and as a function of temperature (d). Best fit linear relations are shown in (b) and (d).

Figure 5.

The distribution (a) and temperature dependence (b) of the difference between the skin temperature measured by Atmospheric Surface Flux Stations (ASFSs) and the best-estimated surface temperature from 7 ice mass balance buoys (IMBs) located near the ASFS. The impact of these different surface temperatures on a bulk calculation of surface conductive flux is also shown as a distribution (c) and as a function of temperature (d). Best fit linear relations are shown in (b) and (d).

Close modal

The uncertainty of the snow–atmosphere interface creates additional uncertainty for the near-surface conductive flux calculations. To assess this impact, we examine bulk estimates of conductive flux using the temperature difference between the snow–atmosphere interface, Tsnow−atm, and the snow–ice interface, Tsnow−ice:

FCbulk=kTsnowatmTsnowiceheff,snow,
5

where heff,snow is the effective snow depth:

heff,snow=n=snowicesnowatmΔz*krefkn,
6

with bounds at the snow–ice interface (n = snow-ice) and snow–atmosphere interface (snow-atm). Effective snow depth is the equivalent thickness of snow if it had the same thermal conductivity as a reference value (kref = 2 W m−1 K−1). Δz is the spacing between individual thermistors, 2 cm. The snow–ice interface level is defined as where mean thermal conductivity profiles are first less than or equal to 0.6 W m−1 K−1. Equation 5 is calculated twice using ASFS measured skin temperatures and IMB surface temperatures for Tsnow-atm. The mean difference between bulk conductive fluxes calculated with these different temperatures is approximately 3 W m−2, with larger differences at colder temperatures (Figure 5c and d). Thus, while not all details of this apparent IMB measurement bias are clear, this independent assessment suggests that it leads to a negative bias ranging from near 0 W m−2 at relatively warm temperatures to an average of about 5 W m−2 at the coldest observed temperatures.

We have quantified various sources of uncertainties in this section; however, it is not a simple task to combine them into a single estimate of uncertainty. Because the individual uncertainty estimates are not independent of each other, they cannot be added in quadrature. The different uncertainties in snow thermal conductivity we have discussed each have similar magnitudes of approximately 0.02–0.06 W m−1 K−1 or approximately 10%–15%. A 10%–15% uncertainty in thermal conductivity would mean an uncertainty in conductive flux around 2–3 W m−2, on average. The uncertainty due to temperature measurements in the thermal conductivity estimates also contributes to uncertainty in conductive flux. The uncertainty associated with the surface location and SIMB temperature bias is on average 3 W m−2. Overall, we estimate the conductive flux uncertainty to be around 5 W m−2 and the thermal conductivity uncertainty to be 0.05 W m−1 K−1.

In this section, we first quantify thermal conductivity in snow and its variability vertically, spatially, and temporally. Next, we assess the conductive flux near the surface as it relates to other, independently measured, surface fluxes. We further determine how environmental factors drive variability in conductive flux through the snowpack, for example, sea ice thickness, snow depth, and atmospheric state. These observation-based relationships between the environment and conductive fluxes are finally used to assess output from CAFS.

4.1. Variability of snow thermal conductivity

Figure 6a and c shows thermal conductivity and density at SIMBA 2019T62, located at the CO near the MOSAiC second year ice coring site. The sea ice grew from 0.80 to 1.48 m over November 1, 2019, to March 1, 2020, and the mean estimated snow depth was 19 ± 2 cm. No flooded sea ice has been identified at this site during this period (Lei et al., 2022). Above 0 cm, the estimated snow–ice interface when the SIMBA was installed, thermal conductivities range from 0.06 to 1.43 W m−1 K−1, with a corresponding density range of 142–721 kg m−3. Snow with the lowest density and thermal conductivity (<0.2 W m−1 K−1) appears at the top of the snowpack, associated with time periods where the snowpack increases in depth, for example, January 15–25, 2020, that are consistent with fresh snow (Matrosov et al., 2022). Thermal conductivity near the snow–ice interface increases from less than 0.75 W m−1 K−1 at +2 cm depth to greater than 1.5 W m−1 K−1 below −4 cm depth. Some noise is also apparent in Figure 6a. Notably, on November 17, 2019, there are sudden, temporary changes around 0–10 cm depth. This episode is during a cyclone (Rinke et al., 2021) that caused warmer surface temperatures and potentially changes to temperatures that were poorly sampled by the SIMBA, as discussed in Section 3.3.

Figure 6.

Effective thermal conductivity (keff) (a and b) and density (c and d) from 2 Snow and Ice Mass Balance Arrays (SIMBAs). SIMBA 2019T62 was located near the second year ice coring site in the Central Observatory, and SIMBA 2019T63 was located at L2. 0 cm represents the nominal sea ice–snow interface at installation. Effective thermal conductivity and density from all ice mass balance buoys are given in Figure S4.

Figure 6.

Effective thermal conductivity (keff) (a and b) and density (c and d) from 2 Snow and Ice Mass Balance Arrays (SIMBAs). SIMBA 2019T62 was located near the second year ice coring site in the Central Observatory, and SIMBA 2019T63 was located at L2. 0 cm represents the nominal sea ice–snow interface at installation. Effective thermal conductivity and density from all ice mass balance buoys are given in Figure S4.

Close modal

SIMBA 2019T63 at L2 shows a different vertical structure than SIMBA 2019T62. Sea ice grew from 1.12 to 1.56 m at L2 during the considered period, and snow depth was 10 ± 2 cm from November through January and increased to 19 cm after February 1, 2020. SIMBA 2019T63 (Figure 6b and d) has a lower range of snow thermal conductivity and density compared to SIMBA 2019T62, 0.10–0.55 W m−1 K−1 and 180–503 kg m−3 above 0 cm. Excluding time periods with noise due to cyclones (e.g., November 17, 2019), thermal conductivity varies less vertically within the snowpack than at SIMBA 2019T62, although low thermal conductivity values also appear near the top of the snow. SIMBA 2019T63 additionally shows a more gradual transition between snow and sea ice, at least in terms of density and thermal conductivity values. Mean thermal conductivity increases from 0.31 W m−1 K−1 at 0 cm to 1.5 W m−1 K−1 at −20 cm, whereas at SIMBA 2019T62, thermal conductivity increases the same amount in about half the thickness.

Different vertical structures are apparent between IMBs when comparing mean vertical profiles of thermal conductivity (Figure 7a and b). SIMB L3 has the largest thermal conductivity throughout the snowpack (approximately 0.45 W m−1 K−1). Several SIMBAs (2019T58, 2019T62, 2019T67, 2019T63, and 2019T65) all have minima in thermal conductivity within the middle of the snowpack with lowest values just above 0.2 W m−1 K−1. 2019T62 (solid blue) varies the most within the snowpack, from approximately 0.4 W m−1 K−1 at 2 cm to approximately 0.2 W m−1 K−1 at 12 cm. Three SIMBAs (2019T47, 2019T64, and 2019T62) show larger thermal conductivity in the topmost layer. These values could be due to wind slab formation or the result of how the top of the snowpack is determined using a threshold of keff = 0.6 W m−1 K−1.

Figure 7.

Vertical and temporal variability of snow thermal conductivity (keff). Panels (a) and (b) show the mean vertical profiles of thermal conductivity in and below the snow. 0 cm represents the apparent sea ice–snow interface at installation. These profiles are used to define the snowpack at each ice mass balance buoy (IMB) based on where the mean profile is less than or equal to 0.6 W m−1 K−1 (dashed vertical line), indicated with stars. Panel (c) shows 2-week running means of snow thermal conductivity. Values are medians from above the calculated snow–ice interface.

Figure 7.

Vertical and temporal variability of snow thermal conductivity (keff). Panels (a) and (b) show the mean vertical profiles of thermal conductivity in and below the snow. 0 cm represents the apparent sea ice–snow interface at installation. These profiles are used to define the snowpack at each ice mass balance buoy (IMB) based on where the mean profile is less than or equal to 0.6 W m−1 K−1 (dashed vertical line), indicated with stars. Panel (c) shows 2-week running means of snow thermal conductivity. Values are medians from above the calculated snow–ice interface.

Close modal

Thermal conductivity profiles from all IMBs have some “transition layer” at least 6 cm thick between the snow and the sea ice, where mean thermal conductivity values are between 0.7 W m−1 K−1 and 1.5 W m−1 K−1. Ten of the IMBs have transition layers at least 8 cm thick, and 3 of those 10 have a transition layer that is 14–16 cm thick. IMBs installed in first year ice (SIMBAs 2019T47 and 2019T70 and SIMB L3) have thinner transition layers (6 or 8 cm), while the thickest transition layers are all observed in second year ice. Past observations found sea ice cores with lower densities and therefore lower expected thermal conductivities, near the top of sea ice due to granular ice, transition layers, and/or the high desalination above the sea level (e.g., Eicken et al., 1995; Huang et al., 2013). Surface scattering layers were also observed at MOSAiC (Macfarlane et al., 2023a), and these transition layers could be remnants of the surface scattering layer from the prior melt season for the second year ice. Due to vapor diffusion and other processes, some sort of transition layer is expected between the snow and ice, and we have no reason to believe the observed transition layers are an artifact of the measurements or processing. Moving forward, we focus on thermal conductivity and conductive fluxes within the snowpack only, defined based on where the mean vertical profile of thermal conductivity is less than or equal to 0.6 W m−1 K−1 (Figure 7a and b).

Spatial differences in snow thermal conductivity are typically larger than temporal variability at a given location. Figure 7c shows the 2-week running means of median thermal conductivity above the snow–ice interface. For most IMBs, snow thermal conductivity varies in time by approximately 0.05–0.1 W m−1 K−1. Local minima, such as around November 25, 2019, for IMBs at L2 (orange solid and dashed lines), are due to higher derived snow top estimates with lower thermal conductivities consistent with snow accumulation. When the upper snow is redistributed, compressed, or undergoes other transformations, the median snow thermal conductivity increases again. Over the full winter, 7 IMBs show increasing trends, consistent with Macfarlane et al. (2023b). Three of the IMBs show overall decreases in snow thermal conductivity, and the remaining 3 have no clear trend and/or large variability.

Median snow thermal conductivity values can vary as much between sites as within a given site. In Figure 8, values are included for levels above the sea ice–snow interface (Figure 7) that are also less than 0.8 W m−1 K−1. This threshold is used to exclude time periods with synoptic events that temporarily and unrealistically increase the apparent thermal conductivity. How the sea ice–snow interface is defined does impact the magnitude of the following statistics but not the relative relationships between IMBs or changes over time. Over winter, median snow thermal conductivity ranges from 0.33 W m−1 K−1 (SIMBA 2019T65 at L2) to 0.47 W m−1 K−1 (SIMB at L3). At L1 (red) median, thermal conductivities are 0.39 and 0.41 W m−1 K−1. The SIMB at L2 has a similar median as IMBs at L1, 0.39 W m−1 K−1 but is significantly higher than the SIMBAs at L2 that have the medians of 0.34 and 0.33 W m−1 K−1 (orange). The SIMB at L3 median is also larger than the SIMBA at L3, 0.47 compared to 0.41 W m−1 K (yellow). At a given L-site, the individual systems are typically within tens of meters of each other. SIMBAs at 3 of the 4 M-sites (light blue) have medians that have a similar range, 0.38–0.43 W m−1 K−1, despite being kilometers apart. The largest interquartile ranges are from SIMB measurements. Over all sites, the median is 0.39 W m−1 K−1 with an interquartile range of 0.34–0.44 W m−1 K−1 (black).

Figure 8.

Mean thermal conductivity (keff) of snow. The time periods used for each ice mass balance buoy are listed in Table 1. The top and bottom of each box are the 25th and 75th percentiles, and the whiskers extend to the 25th percentile minus 1.5 × inter quartile range and the 75th percentile plus 1.5 × interquartile range. The snow–ice interface is based on mean thermal conductivity profiles (Figure 7a and b).

Figure 8.

Mean thermal conductivity (keff) of snow. The time periods used for each ice mass balance buoy are listed in Table 1. The top and bottom of each box are the 25th and 75th percentiles, and the whiskers extend to the 25th percentile minus 1.5 × inter quartile range and the 75th percentile plus 1.5 × interquartile range. The snow–ice interface is based on mean thermal conductivity profiles (Figure 7a and b).

Close modal

4.2. Variability of conductive flux

In this section, we focus on conductive flux and its temporal and spatial variability during winter and among sites. First, a case study is used to highlight the relationships between SEB terms, their variability, and the difficulty of trying to close the SEB. Next, we investigate variability due to sea ice thickness, snow depth, and temperature and quantify how differences in snow impact the coupling between conduction and the atmosphere. Finally, we use these relationships to reassess conductive flux in CAFS.

A case study from November 18 to December 13, 2019, shows the temporal variability of conductive flux and its relationships to other surface fluxes (Figure 9). The first half of this period is dominated by clouds, with some intermittent cloud-free times around November 25, 28, 30, and December 1–2, 2019 (Figure 9a). The atmosphere is mostly cloud-free after December 9. The changes in LWNET (blue) and skin temperature (red) follow cloud cover: LWNET and skin temperature are most negative during cloud-free time periods, for example, November 27 and December 2, 2019 (Figure 9a and c). When the atmosphere transitions to being cloudy, LWNET increases closer to 0 W m−2 and the surface warms by approximately 10°C. SH (orange in Figure 9c) frequently follows changes in skin temperature, although high wind speed can also increase it. SH is directed into the surface during cloud-free and cold periods (December 1–2) and out of the surface during cloudy and warm periods (December 3–5), both on the order of tens of W m−2. LH (green in Figure 9c) is less than 5 W m−2 except for December 3–5.

Figure 9.

Meteorological and energy budget terms over November 18–December 13, 2019. Ka-band radar reflectivity (a), skin temperature and mean wind speed (b), and surface fluxes (c and d) from 3 Atmospheric Surface Flux Stations (ASFSs) and 7 nearby ice mass balance buoys (IMBs) (indicated with “a” in Table 1). Net longwave (LWNET) and conductive flux are positive into the surface, and turbulent fluxes (SH and LH) and positive upward. The conductive flux from IMBs (d; pink Fc IMB) is taken 2 levels below the derived snow top. Panel (e) shows mean 6-h averages of conductive flux from IMBs (pink) and flux plates (black) and -(net LW–SH–LH) (cyan). Shading shows the maximum and minimum values within the 6-h windows across sites.

Figure 9.

Meteorological and energy budget terms over November 18–December 13, 2019. Ka-band radar reflectivity (a), skin temperature and mean wind speed (b), and surface fluxes (c and d) from 3 Atmospheric Surface Flux Stations (ASFSs) and 7 nearby ice mass balance buoys (IMBs) (indicated with “a” in Table 1). Net longwave (LWNET) and conductive flux are positive into the surface, and turbulent fluxes (SH and LH) and positive upward. The conductive flux from IMBs (d; pink Fc IMB) is taken 2 levels below the derived snow top. Panel (e) shows mean 6-h averages of conductive flux from IMBs (pink) and flux plates (black) and -(net LW–SH–LH) (cyan). Shading shows the maximum and minimum values within the 6-h windows across sites.

Close modal

Conductive flux from 7 IMBs that are within tens of meters of the ASFS are plotted in Figure 9d (pink). As previously discussed, the precise location of the top of the snowpack is not known from the IMBs and is difficult to extrapolate from the ASFS snow depth measurements due to strong spatial variability in snow depth. Therefore to be certain the flux is within the snow, we plot conductive flux from 2 levels (4 cm) below the estimated snow top from the thermal conductivity derivation. Conductive flux from the 2 flux plates installed just above the snow–ice interface at each ASFS is also shown in black. Comparing conductive flux estimates from IMBs and flux plates is complicated because of the constant height of the flux plates relative to the sea ice–snow interface and the spatial heterogeneity of the snow. As snow accumulates over the flux plates, the measured conductive flux is reduced, while a storage term should increase. Relating the change in conductive flux to snow depth is complicated, though, because there was only 1 snow depth measurement at each site, while the 2 flux plates can vary by a factor of two. Despite these caveats, conductive flux estimates from IMBs and flux plates in Figure 9d have similar magnitudes and temporal variability, for example, November 18–24. By the end of winter, conductive flux measurements from the different flux plates have smaller magnitudes and temporal variability compared to the IMB estimates (not shown).

The conductive flux responds to changes in the atmosphere, both in terms of magnitude and spatial variability, since it partly depends on the surface temperature. Higher surface temperatures, coincident with cloud occurrence, reduce the temperature gradient through the snowpack and therefore reduce the conductive flux. For example, during the transition from cloudy to cloud-free conditions around November 20, the conductive flux from IMBs increases from near 0 W m−2 to approximately 15–40 W m−2. Variability among flux plates and IMBs also tends to be larger when the surface is cold (November 26–December 1) and smallest when the surface is warm and the atmosphere is cloudy (November 19 and 20). The flux plates at ASFS 50 differ by upward of 40 W m−2 during colder time periods (November 26–December 1) despite being within meters of each other. The differences between flux plates that are meters apart are often larger than differences between IMBs that are kilometers apart, consistent with Sturm et al. (2001). This local variability also demonstrates the modulating effects of sea ice and snow on conduction, although it complicates comparisons between instruments.

SEB closure from surface-based observations, such as those employed here, is a well-recognized challenge (i.e., Foken, 2008). There are multiple reasons why the energy budget from these observations may not fully close, including different scales of spatial variability for each term, different footprints of influence for each measurement, different time resolutions for the measurements, and different locations and snow-ice properties between IMBs and ASFSs. Unsurprisingly, individual IMBs rarely balance the net atmospheric flux (Figure 9d). Before November 23, conductive flux estimates from IMBs are within a few W m−2 of the net atmospheric flux, but the conductive flux is generally smaller than the net atmospheric flux. The significant heterogeneity of snow makes it difficult to compare point measurements of conductive flux, either from IMBs or flux plates, to the footprint of the atmospheric fluxes, which range from a few meters to a few hundred meters depending on parameters and conditions. Therefore, we also compare averages. The pink shading in Figure 9e shows the range of conductive flux from IMBs within 6-h windows, and the solid pink line gives the mean. The cyan shading and line are ranges and mean from ASFS net atmospheric flux, multiplied by −1, during the same 6-h windows; 6-h averages are also calculated from the flux plates, shown in black, to show what another instrument measures even though we know a storage term should be included for them to balance the SEB. At the beginning and end of the time period, the range of conduction overlaps the net atmospheric flux, but the mean conductive flux is still biased low by 10 W m−2 at most. Means from the IMBs and flux plates agree fairly well during these times. Over November 23 through November 30, mean conductive flux is approximately 20–30 W m−2 smaller than the net sum of atmospheric fluxes. This time period has largest difference between IMBs and ASFS fluxes during the winter. For most of the winter at least, the ranges from IMBs and ASFS net atmospheric flux estimates overlap.

Overall, the residual from 6-h mean IMB conduction and ASFS net atmospheric flux has a mean of −6 W m−2 and standard deviation of 7 W m−2 over the entire winter (Figure 10a). It is unsurprising that there is a bias considering the imprecise knowledge from IMBs, including uncertainty in the thermistor measurements, the uncertainty associated with the location of the surface, and the general aggregate uncertainty in deriving thermal conductivity. Additionally, the atmospheric fluxes carry significant uncertainty themselves (Cox et al., 2023), which, when combined in quadrature, result in an aggregated uncertainty for the net atmospheric flux of approximately 5 W m−2. The residuals tend to be larger during radiatively clear time periods, with large negative LWNET, and during periods of high wind (Figure 10b and c). Over the entire winter, the average residual during radiatively clear times (LWNET < −20 W m−2) is −8 W m−2 but only −3 W m−2 during radiatively cloudy times (LWNET > −20 W m−2). Although there is uncertainty in IMB and ASFS measurements, the conductive flux being too small during cold time periods would be consistent with the potential temperature-dependent bias in thermistors discussed in Section 3.3. Some disagreement between IMBs and ASFS could be due to spatial variability, which is generally greater during cold time periods, as discussed more in the following sections. The association between larger wind speeds and residuals could mean that energy is transported at the surface by means other than vertical conduction and/or that the atmospheric turbulent heat flux estimate has a larger uncertainty. However, further investigation into closing the SEB is beyond the scope of this work. While the average residual is nonzero, the fact that it is within the combined uncertainties from atmospheric and conductive fluxes provide some confidence that the IMB estimates are physically reasonable.

Figure 10.

Distribution of the surface energy budget residual (a) and its relationships with net longwave (LWNET) (b) and mean wind speed (c). The residual is calculated from mean net atmospheric heat flux and ice mass balance buoy conductive flux in 6-h windows (Figure 9e). Radiatively clear and cloudy time periods are defined as having LWNET values less than or greater than −20 W m−2, respectively. Means in panel (a) are shown in parentheses.

Figure 10.

Distribution of the surface energy budget residual (a) and its relationships with net longwave (LWNET) (b) and mean wind speed (c). The residual is calculated from mean net atmospheric heat flux and ice mass balance buoy conductive flux in 6-h windows (Figure 9e). Radiatively clear and cloudy time periods are defined as having LWNET values less than or greater than −20 W m−2, respectively. Means in panel (a) are shown in parentheses.

Close modal

4.2.1. Drivers of variability in conduction

Given the large variability, both spatially and temporally, seen in Figure 9d, we now turn to exploring what environmental variables drive that variability. Because snow is a stronger thermal insulator than sea ice but physically thinner, we look at the relationships between conduction and sea ice and snow thicknesses as well as effective depth (heff). Effective depth is the equivalent thickness of snow and sea ice scaled by a reference value (kref = 2 W m−1 K−1):

heff=n=30cmsnowatm[Δz*krefkn]+krefk¯ice*(hice30),
7

where kn are the thermal conductivities over each 2-cm layer Δz that are summed over the portion of the column for which thermal conductivity is derived above −30 cm. The upper bound is the snow–atmosphere interface “snow–atm.” Depths are relative to the 0 cm defined at IMB installation. The right-hand term represents the depth of sea ice below the snowpack and possible transition layer (hice) using the average thermal conductivity for the bulk sea ice at and below −30 cm (k¯ice), generally close to 2 W m−1 K−1. Effective snow depth is calculated in the same manner as effective depth but only for levels within the snowpack (Equation 6).

Effective depth is useful for considering the thickness impacts of snow compared to sea ice. For example, 10 cm of snow with ksnow = 0.3 W m−1 K−1 is equivalent to 67 cm of sea ice with kice = 2.0 W m−1 K−1, in terms of how readily each medium conducts energy. We further use effective snow depth on its own to compare the sensitivity of conduction to sea ice and snow depths. In Figure 11, snowpack-average conductive fluxes are plotted against snow depth, effective snow depth, sea ice thickness, and total effective depth, and colored by surface temperature. Best fit lines are calculated for 3 temperature ranges shown in the legends of each subplot; 27% of the data falls within the coldest temperature range, 17% is in the warmest temperature range, and 55% is in the middle temperature range. Using temperature ranges that evenly split the percentage of data in each temperature range results in the same basic relationships but different magnitudes.

Figure 11.

Snowpack averaged conductive fluxes regressed against physical snow depth (a–c), effective snow depth (d–f), sea ice thickness (g–i), and effective depth (j–l). Columns are separated into temperature ranges, and shading represents surface temperature. Solid black lines show linear fits to data within temperature ranges (slopes (r2) values provided in legends), and cyan dashed lines show exponential fits. Values are taken at ice mass balance buoy time resolutions, either 4 or 6 hourly.

Figure 11.

Snowpack averaged conductive fluxes regressed against physical snow depth (a–c), effective snow depth (d–f), sea ice thickness (g–i), and effective depth (j–l). Columns are separated into temperature ranges, and shading represents surface temperature. Solid black lines show linear fits to data within temperature ranges (slopes (r2) values provided in legends), and cyan dashed lines show exponential fits. Values are taken at ice mass balance buoy time resolutions, either 4 or 6 hourly.

Close modal

Thicker snowpack and warmer surface temperature reduce the variability of conductive flux. Differences in snow depth, both physical and effective, become less important at warmer temperatures because the temperature gradient is reduced. From linear regressions (solid lines in Figure 11a–c), the sensitivity of conductive flux to physical snow depth decreases from −0.81 W m−2/cm for the coldest temperatures (−40°C to −28°C) to −0.47 W m−2/cm or the warmest temperatures (−17°C to −4°C). The sensitivity to effective snow depth shows a similar decreases, by over a factor of two, from −0.11 to −0.04 W m−2/cm. Put another way, differences in snowpack thickness result in more spatial variability in conductive flux when the surface is coldest. In addition to linear regressions, exponential fits are shown with cyan dashed lines in Figure 11a–f. Exponential fits may be more representative of the relationship between conduction and snow depth than a linear regression, given that the analytical solution of the heat equation using Fourier series is exponential. Exponential fits are also more physically intuitive since conduction will not become increasingly negative for thicker snowpacks. This reasoning appears true for the middle and warmest temperature ranges; however, the exponential fit for the coldest temperature range still appears linear, reiterating the point that conductive fluxes in the coldest temperature range have higher sensitivity, that is, a steeper curve, to effective snow depth than the warmest temperatures. The exponential fits for the medium and warmest temperature ranges approach asymptotes around 5–10 W m−2, whereas in the case of infinitely thick snow, we might expect the conductive flux to approach 0 W m−2. The nonzero asymptote could be due to the relatively thin sea ice below that still allows some conduction to occur. We can also consider the inverse relationship: how snow thickness influences the sensitivity of conduction to surface temperature. As snow depth increases, the range of conductive fluxes decreases (Figure 11a). Similarly, for effective snow depth less than 100 cm, the standard deviation of conductive flux is 10.4 W m−2, whereas for effective snow depth thicker than 200 cm, the standard deviation is only 3.7 W m−2, less than half the variability (Figure 11d–f).

In contrast to the decreasing sensitivity of conduction to snow depth with temperature, conductive flux in snow is about equally sensitive to sea ice thickness in each temperature range, −0.09 to −0.11 W m−2/cm (Figure 11g–i). This consistency could be due to the fact that sea ice grows as the surface temperature overall decreases throughout the winter. For the coldest temperature range, conductive flux is about as sensitive to sea ice thickness as effective snow depth, but at warmer surface temperatures, conductive flux is more sensitive to sea ice thickness. This sensitivity could drive the variability between IMBs during warm periods in Figure 9d. Because conduction is still sensitive to sea ice thickness when the surface temperature is warm, the sensitivity of conduction to effective depth is not quite as varied between temperature ranges as the sensitivity to effective snow depth. Conduction is less than half as sensitive to effective depth over the warmest temperatures (−0.03 W m−2/cm) as compared to the colder temperatures (−0.08 W m−2/cm; Figure 11j–l).

Conductive flux is coupled to atmospheric forcing through the surface temperature, which is driven by the balance of longwave radiation, turbulent heat fluxes, and conduction. The slope from the line of best fit from regressing conductive flux against LWD (Figure 12a) represents the sensitivity of conductive flux to winter radiative forcing. For each additional W m−2 of LWD, conductive flux decreases by 0.22 W m−2. For comparison, an additional W m−2 of LWD increases LWU by approximately 0.5 W m−2 (not shown). The sensitivity of conduction to LWD is actually the same regardless of effective snow depth, although the range is larger for thinner effective snow depths than thicker effective snow depths.

Figure 12.

Average conductive fluxes through the snowpack versus downwelling longwave (LWD) (a and c) and net longwave (LWNET) (b and d) colored by frequency of occurrence. Longwave fluxes are the mean from Atmospheric Surface Flux Stations (ASFS) and Central Observatory (CO) observations interpolated from 10-min observations to ice mass balance buoy (IMB) time steps.

Figure 12.

Average conductive fluxes through the snowpack versus downwelling longwave (LWD) (a and c) and net longwave (LWNET) (b and d) colored by frequency of occurrence. Longwave fluxes are the mean from Atmospheric Surface Flux Stations (ASFS) and Central Observatory (CO) observations interpolated from 10-min observations to ice mass balance buoy (IMB) time steps.

Close modal

Previous work has pointed to the bimodal winter radiative states of the Arctic sea ice system, which are largely defined by the presence or absence of clouds (e.g., Shupe and Intrieri, 2004; Stramler et al., 2011). This bimodality is evident in the distribution of LWNET (Figure 12b), with −25 W m−2 serving as the approximate threshold between radiatively clear and cloudy conditions. The joint histogram between LWNET and conductive fluxes also has distinct clear and cloudy clusters (Figure 12d). Despite the distribution of conductive fluxes appearing unimodal (Figure 12e), mean conductive fluxes in clear versus cloudy states are statistically different: 25 W m−2 (10 W m−2 standard deviation) for radiatively clear conditions and 14 W m−2 (10 W m−2 standard deviation) for radiatively cloudy times. The difference in conduction between clear and cloudy states has a slight dependence on effective snow depth. The difference between radiatively clear and cloudy states is slightly larger for thinner snowpacks than thicker snowpacks. For thinner effective snow depths (<92 cm), the means are 29 W m−2 versus 16 W m−2 for radiatively clear versus cloudy. The means drop to 20 W m−2 versus 11 W m−2 for thicker effective snow depths (>92 cm). These differences imply that thicker snow modestly dampens the coupling between the atmosphere and sea ice.

4.2.2. Assessment of a coupled model

Given these relationships between conductive flux and the environment, we now return to CAFS to reassess how it compares to the observations. CAFS data are analyzed from a simulation over the MOSAiC year with output from grid cells nearest the Polarstern and ASFSs. For this run, the first 48 h of each daily forecast are used. Fluxes are positive toward the surface. For comparisons, near-surface conductive flux estimates are taken from IMBs, and data from CAFS and observations are averaged to daily timescales.

The relationships between sea ice thickness, snow depth, and conductive flux differ between CAFS and observations. In CAFS, sea ice thickness ranges from 161 to 254 cm over the winter (blue in Figure 13b). These values are generally larger with less variability over the winter compared to observations that vary from 50 to 208 cm (orange in Figure 13b). Snow depths in CAFS are within the observed range (4–29 cm) but have a narrower spread (7–20 cm). Conductive flux in CAFS depends less on both sea ice thickness and effective snow depth than in observations (Figure 13a). In observations, the thicker the sea ice or effective snow is, the less conduction. However, in CAFS even the thickest sea ice (approximately 250 cm) experiences the full range of conductive flux (about −10 to 60 W m−2 in Figure 13d). Sea ice and snow do not appear to modulate conductive flux in CAFS in the same way that they do in observations.

Figure 13.

Comparison of relationships between conductive flux and snow depth (a and c), sea ice thickness (b and d), downwelling longwave (LWD) (f and h), and surface temperature (g and i) from the Coupled Arctic Forecast System (CAFS) and observations. Daily averages from observations (seasonal ice mass buoys [IMBs] and ASFS) and CAFS are from October through February. Conductive fluxes are positive upward. Conductive fluxes from IMBs are near-surface estimates. Normalized distributions are shown for snow depth (a), sea ice thickness (b), LWD (f), surface temperature (g), and conductive flux (e and j).

Figure 13.

Comparison of relationships between conductive flux and snow depth (a and c), sea ice thickness (b and d), downwelling longwave (LWD) (f and h), and surface temperature (g and i) from the Coupled Arctic Forecast System (CAFS) and observations. Daily averages from observations (seasonal ice mass buoys [IMBs] and ASFS) and CAFS are from October through February. Conductive fluxes are positive upward. Conductive fluxes from IMBs are near-surface estimates. Normalized distributions are shown for snow depth (a), sea ice thickness (b), LWD (f), surface temperature (g), and conductive flux (e and j).

Close modal

Differences between observed and modeled conductive fluxes are also tied to atmospheric state. In CAFS, daily average conductive fluxes span a larger range (−10.9 to 62.4 W m−2) than observations (−3.7 to 48.0 W m−2). The larger conductive fluxes in CAFS are associated with colder surface temperatures. In Figure 13g, CAFS has more occurrences of surface temperatures below −30°C compared to surface temperatures from IMBs and skin temperatures from ASFS. These colder surface temperatures occur during cloud-free periods, which CAFS has slightly more occurrences of compared to observations (Figure 13f). This difference suggests that CAFS simulates an optically thinner atmosphere during clear-sky conditions. CAFS also has more instances of conduction going downward from the surface, that is, negative conductive fluxes, than in the observations; 4% of conductive fluxes in CAFS are less than 0 W m−2 compared to 1% from IMBs and less than 1% from the flux plates. In observations, downward conduction at the surface typically only occurred during cyclones when the atmosphere was warmer than the surface. If the surface is more often colder than the atmosphere, negative conduction is more likely to occur. These connections between conduction, surface temperature, and the atmosphere underscore the importance and complications of modeling the Arctic as a coupled system.

Our thermal conductivity values in snow are within the range of past observations and other measurements from MOSAiC, albeit on the high end. Including all measurements leads to a median snow thermal conductivity of 0.39 W m−1 K−1, which is higher than the value used by large scale models, although it is sensitive to how the sea ice–snow interface is defined. Macfarlane et al. (2023b) measured a median snow thermal conductivity of 0.26 ± 0.05 W m−1 K−1 using a finite element method and microCT scans of snow samples at MOSAiC. They defined snow as samples with density less than 550 kg m−3. From our estimates, none of the 13 IMBs have median values in snow less than 0.3 W m−1 K−1, and our range extends above the maximum of 0.5 W m−1 K−1 in Macfarlane et al. (2023b). Our larger estimates could, in part, be due to differences in defining snow and processes are included in the “effective” definition of thermal conductivity. Macfarlane et al. (2023b) consider processes that contribute to vertical conduction, while our methods assume any measured temperature changes are due to vertical conduction. Perovich et al. (2023) also derived snow thermal conductivity from the SIMBs at MOSAiC around 0.39–0.48 W m−1 K−1 by assuming flux continuity at the sea ice–snow interface. Our median estimates from the same SIMB instruments are similar; however, including results from additional IMBs increases the spread of thermal conductivity. Future work should be done to synthesize the different thermal conductivity estimates from MOSAiC.

A key point of this work is that the representation of snow thermal conductivity is important for modeling the SEB and the lower atmosphere. In agreement with past work (Wu et al., 1999; Blazey et al., 2013), changing snow thermal conductivity in a coupled model also changes the surface temperature and balance of surface fluxes. While we are not the first to acknowledge this, the representation of snow has received far less consideration in atmospheric modeling as compared to sea ice modeling. To wit, not all forecast models even include snow on sea ice (e.g., Bauer et al., 2016; Bengtsson et al., 2017). Given the observed variability of snow thermal conductivity, future studies should test the extent to which spatial and temporal variability of snow thermal conductivity is important in coupled models.

Although conductive flux estimates are often within measurement uncertainties for closing the SEB, future work can improve upon SEB closure. The estimates of conductive flux at the snow–atmosphere interface from thermistor string measurements are uncertain due to both the challenge in identifying the surface location as well as potential temperature-dependent biases in the thermistors. If the thermistors have larger positive biases at colder temperatures (as suggested by Figure 5), the temperature gradients during these time periods would be biased low. Underestimated temperature gradients would cause underestimated conductive flux, consistent with our conductive flux estimates being typically smaller than the net atmospheric flux. Future studies could further assess, quantify, and correct this potential bias in the temperature profile data.

We have estimated effective snow thermal conductivity during October–February over the Central Arctic sea ice during the MOSAiC expedition using temperature profiles from sea ice mass balance buoys. Our method provides vertical profiles of thermal conductivity, density, and conductive flux that are continuous in time. Median values from 13 IMBs range from 0.33 to 0.47 W m−1 K−1, with spatial variability greater than the estimated 0.05 W m−1 K−1 uncertainty. The median from all available data during October through February is 0.39 W m−1 K−1, and the interquartile range is 0.34–0.44 W m−1 K−1. The temporal and vertical variability of snow thermal conductivity and density are consistent with environmental changes, such as new snowfall.

Estimated thermal conductivity was used to calculate conductive flux through the snowpack. Conductive fluxes are larger and more spatially variable during cold and radiatively clear periods as compared to warmer and cloudier times. Typically, the derived conductive flux is smaller than the observationally estimated net atmospheric flux, leading to a mean residual of −6 W m−2 over the winter, with larger residuals during cold, cloud-free periods and diminished residuals during warm, cloudy periods. The mean residual is on the same order of magnitude as uncertainties in the conductive flux and the net atmospheric flux (each approximately 5 W m−2), but there can still be considerable residuals at any specific time.

Using a coupled regional forecast model, we found that varying snow thermal conductivity impacts the balance of all surface heat fluxes. Importantly, the snow thermal conductivity modulates the interplay between atmospheric radiation and turbulent heat fluxes. Detailed and time-variable thermal conductivity and conductive flux measurements can be coupled with the time-varying estimates of atmospheric surface flux terms to better explore the processes that help to drive the surface temperature and the overall SEB. In both observations and models, coupling the atmosphere and surface is essential for linking atmospheric variability to the thermal state of the snow and sea ice as well as to sea ice growth in winter. Future work should further explore these processes for model assessment and improvement.

Coupled Arctic Forecast System output is being archived on PANGAEA. Atmospheric Surface Flux Stations and Met City measurements are available from https://doi.org/10.18739/A2PV6B83F, https://doi.org/10.18739/A2XD0R00S, https://doi.org/10.18739/A25X25F0P, and https://doi.org/10.18739/A2FF3M18K. SIMB data are available at https://doi.org/10.18739/A2SB3X05G. SIMBA data are available from PANGAEA at the following: https://doi.org/10.1594/PANGAEA.940393 (2019T58), https://doi.org/10.1594/PANGAEA.940387 (2019T47), https://doi.org/10.1594/PANGAEA.940231 (2019T62), https://doi.org/10.1594/PANGAEA.940593 (2019T63), https://doi.org/10.1594/PANGAEA.940617 (2019T64), https://doi.org/10.1594/PANGAEA.940634 (2019T65), https://doi.org/10.1594/PANGAEA.938128 (2019T67), https://doi.org/10.1594/PANGAEA.940650 (2019T68), https://doi.org/10.1594/PANGAEA.940659 (2019T70), and https://doi.org/10.1594/PANGAEA.940668 (2019T72). Derived snow and ice thermal properties are available from https://doi.org/10.18739/A2CF9J80P.

The supplemental files for this article can be found as follows:

Supplemental material is available in sledd_elementa_supplemental_20231207.pdf.

Data sets used in this article were produced as part of the International Multidisciplinary drifting Observatory for the Study of the Arctic Climate (MOSAiC) with the tag MOSAiC20192020 and the Project_ID: AWI_PS122_00. Some data used here were collected by the Atmospheric Radiation Measurement User Facility, a DOE Office of Science user facility operated by the Biological and Environmental Research program. The authors thank all people involved in the expedition of the Research Vessel Polarstern (Alfred-Wegener-Instut, Helmholtz-Zentrum für Polar-und Meeresforschung, 2017) during MOSAiC in 2019–2020 as listed in Nixdorf et al. (2021). They additionally thank Darren Jackson for feedback on the manuscript.

This work was supported by DOE Atmospheric System Research Program (DE-SC0021341), the National Science Foundation (OPP-1724551), NOAA’s Physical Sciences Laboratory (PSL; NOAA Cooperative Agreement NA22OAR4320151), and by NOAA’s Global Ocean Monitoring and Observing Program (GOMO; FundRef https://doi.org/10.13039/100018302). RL is funded by the Program of Shanghai Academic/Technology Research Leader (22XD1403600). DP was supported by NSF-OPP 1724424, NSF–OPP–2034919, and NSF–OPP–2138785.

The authors declare that they have no conflict of interest.

Contributed to conception and design: AS, MDS, AS, CJC.

Contributed to acquisition of data: MDS, AS, CJC, DP, RL.

Contributed to analysis and interpretation of data: AS, MDS, AS, CJC, DP.

Drafted and/or revised the article: All authors.

Approved the submitted version for publication: All authors.

Alfred-Wegener-Instut, Helmholtz-Zentrum für Polar-und Meeresforschung
.
2017
.
Polar research and supply vessel POLARSTERN operated by the Alfred-Wegener-Institute
.
Journal of Large-Scale Research Facilities JLSRF
3
:
A119
.
Batrak
,
Y
,
Müller
,
M.
2019
.
On the warm bias in atmospheric reanalyses induced by the missing snow over Arctic sea-ice
.
Nature Communications
10
(
1
):
4170
.
Bauer
,
P
,
Magnusson
,
L
,
Thépaut
,
JN
,
Hamill
,
TM.
2016
.
Aspects of ECMWF model performance in polar areas
.
Quarterly Journal of the Royal Meteorological Society
142
(
695
):
583
596
.
Bengtsson
,
L
,
Andrae
,
U
,
Aspelien
,
T
,
Batrak
,
Y
,
Calvo
,
J
,
de Rooy
,
W
,
Gleeson
,
E
,
Hansen-Sass
,
B
,
Homleid
,
M
,
Hortal
,
M
,
Ivarsson
,
K-I
,
Lenderink
,
G
,
Niemelä
,
S
,
Nielsen
,
KP
,
Onvlee
,
J
,
Rontu
,
L
,
Samuelsson
,
P
,
Muñoz
,
DS
,
Subias
,
A
,
Tijm
,
S
,
Toll
,
V
,
Yang
,
X
,
Køltzow
,
.
2017
.
The HARMONIE–AROME model configuration in the ALADIN–HIRLAM NWP system
.
Monthly Weather Review
145
(
5
):
1919
1935
.
Bigdeli
,
A
,
Nguyen
,
A
,
Pillar
,
H
,
Ocaña
,
V
,
Heimbach
,
P.
2020
.
Atmospheric warming drives growth in Arctic Sea ice: A key role for snow
.
Geophysical Research Letters
47
(
20
):
e2020GL090236
.
Blazey
,
B
,
Holland
,
M
,
Hunke
,
E.
2013
.
Arctic Ocean sea ice snow depth evaluation and bias sensitivity in CCSM
.
The Cryosphere
7
(
6
):
1887
1900
.
Brandt
,
RE
,
Warren
,
SG.
1997
.
Temperature measurements and heat transfer in near-surface snow at the South Pole
.
Journal of Glaciology
43
(
144
):
339
351
.
Calonne
,
N
,
Flin
,
F
,
Morin
,
S
,
Lesaffre
,
B
,
du Roscoat
,
SR
,
Geindreau
,
C.
2011
.
Numerical and experimental investigations of the effective thermal conductivity of snow
.
Geophysical Research Letters
38
(
23
).
Calonne
,
N
,
Milliancourt
,
L
,
Burr
,
A
,
Philip
,
A
,
Martin
,
CL
,
Flin
,
F
,
Geindreau
,
C.
2019
.
Thermal conductivity of snow, firn, and porous ice from 3-D image-based computations
.
Geophysical Research Letters
46
(
22
):
13079
13089
.
Colbeck
,
SC.
1982
.
An overview of seasonal snow metamorphism
.
Reviews of Geophysics
20
(
1
):
45
61
.
Cox
,
CJ
,
Gallagher
,
MR
,
Shupe
,
MD
,
Persson
,
POG
,
Solomon
,
A
,
Fairall
,
CW
,
Ayers
,
T
,
Blomquist
,
B
,
Brooks
,
IM
,
Costa
,
D
,
Grachev
,
A
,
Gottas
,
D
,
Hutchings
,
JK
,
Kutchenreiter
,
M
,
Leach
,
J
,
Morris
,
SM
,
Morris
,
V
,
Osborn
,
J
,
Pezoa
,
S
,
Preußer
,
A
,
Riihimaki
,
LD
,
Uttal
,
T.
2023
.
Continuous observations of the surface energy budget and meteorology over the Arctic sea ice during MOSAiC
.
Scientific Data
10
(
1
):
519
.
Cuffey
,
KM
,
Paterson
,
WSB
.
2010
.
The physics of glaciers
.
Amsterdam, the Netherlands
:
Academic Press
.
Eicken
,
H
,
Lensu
,
M
,
Leppäranta
,
M
,
Tucker
,
W
III
,
Gow
,
A
,
Salmela
,
O.
1995
.
Thickness, structure, and properties of level summer multiyear ice in the Eurasian sector of the Arctic Ocean
.
Journal of Geophysical Research: Oceans
100
(
C11
):
22697
22710
.
Foken
,
T.
2008
.
The energy balance closure problem: An overview
.
Ecological Applications
18
(
6
):
1351
1367
.
Huang
,
W
,
Lei
,
R
,
Ilkka
,
M
,
Li
,
Q
,
Wang
,
Y
,
Li
,
Z.
2013
.
The physical structures of snow and sea ice in the Arctic section of 150°–180° W during the summer of 2010
.
Acta Oceanologica Sinica
32
:
57
67
.
Hunke
,
EC
,
Lipscomb
,
WH
,
Jones
,
P
,
Turner
,
A
,
Jeffery
,
N
,
Elliott
,
S.
2017
.
CICE, the Los Alamos sea ice model
.
Los Alamos, NM
:
Los Alamos National Lab (LANL
).
Hunke
,
EC
,
Lipscomb
,
WH
,
Turner
,
AK.
2010
.
Sea-ice models for climate study: Retrospective and new directions
.
Journal of Glaciology
56
(
200
):
1162
1172
.
Jackson
,
K
,
Wilkinson
,
J
,
Maksym
,
T
,
Meldrum
,
D
,
Beckers
,
J
,
Haas
,
C
,
Mackenzie
,
D.
2013
.
A novel and low-cost sea ice mass balance buoy
.
Journal of Atmospheric and Oceanic Technology
30
(
11
):
2676
2688
.
Kay
,
JE
,
Gettelman
,
A.
2009
.
Cloud influence on and response to seasonal Arctic sea ice loss
.
Journal of Geophysical Research: Atmospheres
114
(
D18
).
Krumpen
,
T
,
Sokolov
,
V.
2020
.
The Expedition AF122/1: Setting up the MOSAiC Distributed Network in October 2019 with Research Vessel AKADEMIK FEDOROV
.
Berichte zur Polar-und Meeresforschung (Reports on Polar and Marine Research) 744
.
Landrum
,
LL
,
Holland
,
MM.
2022
.
Influences of changing sea ice and snow thicknesses on simulated Arctic winter heat fluxes
.
The Cryosphere
16
(
4
):
1483
1495
.
Lecomte
,
O
,
Fichefet
,
T
,
Vancoppenolle
,
M
,
Dominé
,
F
,
Massonnet
,
F
,
Mathiot
,
P
,
Morin
,
S
,
Barriat
,
PY.
2013
.
On the formulation of snow thermal conductivity in large-scale sea ice models
.
Journal of Advances in Modeling Earth Systems
5
(
3
):
542
557
.
Lei
,
R
,
Cheng
,
B
,
Hoppmann
,
M
,
Zhang
,
F
,
Zuo
,
G
,
Hutchings
,
JK
,
Lin
,
L
,
Lan
,
M
,
Wang
,
H
,
Regnery
,
J
,
Krumpen
,
T
,
Haapala
,
J
,
Rabe
,
B
,
Perovich
,
DK
,
Nicolaus
,
M.
2022
.
Seasonality and timing of sea ice mass balance and heat fluxes in the Arctic transpolar drift during 2019–2020
.
Elementa: Science of the Anthropocene
10
(
1
):
000089
.
Lindenmaier
,
I
,
Feng
,
YC
,
Johnson
,
K
,
Nelson
,
D
,
Isom
,
B
,
Hardin
,
J
,
Matthews
,
A
,
Wendler
,
T
,
Castro
,
V
,
Deng
,
M
.
n.d.
Ka ARM Zenith Radar (KAZRCFRGE)
. DOI: http://dx.doi.org/10.5439/1498936.
Lipscomb
,
WH.
1998
.
Modeling the thickness distribution of Arctic sea ice [PhD thesis]
.
Seattle, WA
:
University of Washington
.
Liston
,
GE
,
Itkin
,
P
,
Stroeve
,
J
,
Tschudi
,
M
,
Stewart
,
JS
,
Pedersen
,
SH
,
Reinking
,
AK
,
Elder
,
K.
2020
.
A Lagrangian snow-evolution system for sea-ice applications (SnowModel-LG): Part I—Model description
.
Journal of Geophysical Research: Oceans
125
(
10
):
e2019JC015913
.
Macfarlane
,
AR
,
Dadic
,
R
,
Smith
,
MM
,
Light
,
B
,
Nicolaus
,
M
,
Henna-Reetta
,
H
,
Webster
,
M
,
Linhardt
,
F
,
Hämmerle
,
S
,
Schneebeli
,
M.
2023
a.
Evolution of the microstructure and reflectance of the surface scattering layer on melting, level Arctic sea ice
.
Elementa: Science of the Anthropocene
11
(
1
). DOI: https://doi.org/10.1525/elementa.2022.00103.
Macfarlane
,
AR
,
Löwe
,
H
,
Gimenes
,
L
,
Wagner
,
DN
,
Dadic
,
R
,
Ottersberg
,
R
,
Hämmerle
,
S
,
Schneebeli
,
M
.
2023
b.
Temporospatial variability of snow's thermal conductivity on Arctic sea ice
.
The Cryosphere
17
:
5417
5434
.
Marchenko
,
S
,
Cheng
,
G
,
Lötstedt
,
P
,
Pohjola
,
V
,
Pettersson
,
R
,
Van Pelt
,
W
,
Reijmer
,
C.
2019
.
Thermal conductivity of firn at Lomonosovfonna, Svalbard, derived from subsurface temperature measurements
.
The Cryosphere
13
(
7
):
1843
1859
.
Matrosov
,
SY
,
Shupe
,
MD
,
Uttal
,
T.
2022
.
High temporal resolution estimates of Arctic snowfall rates emphasizing gauge and radar-based retrievals from the MOSAiC expedition
.
Elementa: Science of the Anthropocene
10
(
1
):
00101
.
Maykut
,
GA.
1978
.
Energy exchange over young sea ice in the central Arctic
.
Journal of Geophysical Research: Oceans
83
(
C7
):
3646
3658
.
Nicolaus
,
M
,
Perovich
,
DK
,
Spreen
,
G
,
Granskog
,
MA
,
von Albedyll
,
L
,
Angelopoulos
,
M
,
Anhaus
,
P
,
Arndt
,
S
,
Belter
,
HJ
,
Bessonov
,
V
,
Birnbaum
,
G
,
Brauchle
,
J
,
Calmer
,
R
,
Cardellach
,
E
,
Cheng
,
B
,
Clemens-Sewall
,
D
,
Dadic
,
R
,
Damm
,
E
,
de Boer
,
G
,
Demir
,
O
,
Dethloff
,
K
,
Divine
,
DV
,
Fong
,
AA
,
Fons
,
S
,
Frey
,
MM
,
Fuchs
,
N
,
Gabarró
,
C
,
Gerland
,
S
,
Goessling
,
HF
,
Gradinger
,
R
,
Haapala
,
J
,
Haas
,
C
,
Hamilton
,
J
,
Hannula
,
H-R
,
Hendricks
,
S
,
Herber
,
A
,
Heuzé
,
C
,
Hoppmann
,
M
,
Høyland
,
KV
,
Huntemann
,
M
,
Hutchings
,
JK
,
Hwang
,
B
,
Itkin
,
P
,
Jacobi
,
H-W
,
Jaggi
,
M
,
Jutila
,
A
,
Kaleschke
,
L
,
Katlein
,
C
,
Kolabutin
,
N
,
Krampe
,
D
,
Kristensen
,
SS
,
Krumpen
,
T
,
Kurtz
,
N
,
Lampert
,
A
,
Lange
,
BA
,
Lei
,
R
,
Light
,
B
,
Linhardt
,
F
,
Liston
,
GE
,
Loose
,
B
,
Macfarlane
,
AR
,
Mahmud
,
M
,
Matero
,
IO
,
Maus
,
S
,
Morgenstern
,
A
,
Naderpour
,
R
,
Nandan
,
V
,
Niubom
,
A
,
Oggier
,
M
,
Oppelt
,
N
,
Pätzold
,
F
,
Perron
,
C
,
Petrovsky
,
T
,
Pirazzini
,
R
,
Polashenski
,
C
,
Rabe
,
B
,
Raphael
,
IA
,
Regnery
,
J
,
Rex
,
M
,
Ricker
,
R
,
Riemann-Campe
,
K
,
Rinke
,
A
,
Rohde
,
J
,
Salganik
,
E
,
Scharien
,
RK
,
Schiller
,
M
,
Schneebeli
,
M
,
Semmling
,
M
,
Shimanchuk
,
E
,
Shupe
,
MD
,
Smith
,
MM
,
Smolyanitsky
,
V
,
Sokolov
,
V
,
Stanton
,
T
,
Stroeve
,
J
,
Thielke
,
L
,
Timofeeva
,
A
,
Tonboe
,
RT
,
Tavri
,
A
,
Tsamados
,
M
,
Wagner
,
DN
,
Watkins
,
D
,
Webster
,
M
,
Wendisch
,
M.
2022
.
Overview of the MOSAiC expedition: Snow and sea ice
.
Elementa: Science of the Anthropocene
10
(
1
):
000046
.
Nixdorf
,
U
,
Dethloff
,
K
,
Rex
,
M
,
Shupe
,
M
,
Sommerfeld
,
A
,
Perovich
,
DK
,
Nicolaus
,
M
,
Heuzé
,
C
,
Rabe
,
B
,
Loose
,
B
,
Damm
,
E
,
Gradinger
,
R
,
Fong
,
A
,
Maslowski
,
W
,
Rinke
,
A
,
Kwok
,
R
,
Spreen
,
G
,
Wendisch
,
M
,
Herber
,
A
,
Hirsekorn
,
M
,
Mohaupt
,
V
,
Frickenhaus
,
S
,
Immerz
,
A
,
Weiss-Tuider
,
K
,
König
,
B
,
Mengedoht
,
D
,
Regnery
,
J
,
Gerchow
,
P
,
Ransby
,
D
,
Krumpen
,
T
,
Morgenstern
,
A
,
Haas
,
C
,
Kanzow
,
T
,
Rack
,
FR
,
Saitzev
,
V
,
Sokolov
,
V
,
Makarov
,
A
,
Schwarze
,
S
,
Wunderlich
,
T
,
Wurr
,
K
,
Boetius
,
A.
2021
. MOSAiC extended acknowledgement.
Zenodo
.
Notz
,
D
,
Wettlaufer
,
JS
,
Worster
,
MG.
2005
.
A non-destructive method for measuring the salinity and solid fraction of growing sea ice in situ
.
Journal of Glaciology
51
(
172
):
159
166
.
Oldroyd
,
HJ
,
Higgins
,
C
,
Huwald
,
H
,
Selker
,
JS
,
Parlange
,
M.
2013
.
Thermal diffusivity of seasonal snow determined from temperature profiles
.
Advances in Water Resources
55
:
121
130
.
Perovich
,
D
,
Raphael
,
I
,
Moore
,
R
,
Clemens-Sewall
,
D
,
Lei
,
R
,
Sledd
,
A
,
Polashenski
,
C.
2023
.
Sea ice temperature and mass balance measurements during the MOSAiC drift campaign
.
Elementa: Science of the Anthropocene
11
(
1
):
00017
.
Perovich
,
D
,
Raphael
,
I
,
Moore
,
R
,
Clemens-Sewall
,
D
,
Polashenski
,
C
,
Planck
,
C.
2022
.
Measurements of ice mass balance and temperature from autonomous seasonal ice mass balance buoys in the Arctic Ocean, 2019-2020
.
Arctic Data Center
. DOI: http://dx.doi.org/10.18739/A2SB3X05G.
Petrich
,
C
,
Eicken
,
H.
2010
.
Growth, structure and properties of sea ice
.
Sea Ice
2
:
23
77
.
Planck
,
CJ
,
Whitlock
,
J
,
Polashenski
,
C
,
Perovich
,
D.
2019
.
The evolution of the seasonal ice mass balance buoy
.
Cold Regions Science and Technology
165
:
102792
.
Polashenski
,
C
,
Perovich
,
D
,
Richter-Menge
,
J
,
Elder
,
B.
2011
.
Seasonal ice mass-balance buoys: Adapting tools to the changing Arctic
.
Annals of Glaciology
52
(
57
):
18
26
.
Pringle
,
D
,
Eicken
,
H
,
Trodahl
,
H
,
Backstrom
,
L.
2007
.
Thermal conductivity of landfast Antarctic and Arctic sea ice
.
Journal of Geophysical Research: Oceans
112
(
C4
).
Pringle
,
D
,
Trodahl
,
H
,
Haskell
,
T.
2006
.
Direct measurement of sea ice thermal conductivity: No surface reduction
.
Journal of Geophysical Research: Oceans
111
(
C5
).
Ridley
,
JK
,
Blockley
,
EW
,
Keen
,
AB
,
Rae
,
JG
,
West
,
AE
,
Schroeder
,
D.
2018
.
The sea ice model component of HadGEM3-GC3. 1
.
Geoscientific Model Development
11
(
2
):
713
723
.
Rinke
,
A
,
Cassano
,
JJ
,
Cassano
,
EN
,
Jaiser
,
R
,
Handorf
,
D.
2021
.
Meteorological conditions during the MOSAiC expedition: Normal or anomalous?
Elementa: Science of the Anthropocene
9
(
1
):
00023
.
Shupe
,
MD
,
Intrieri
,
JM.
2004
.
Cloud radiative forcing of the Arctic surface: The influence of cloud properties, surface albedo, and solar zenith angle
.
Journal of Climate
17
(
3
):
616
628
.
Shupe
,
MD
,
Rex
,
M
,
Blomquist
,
B
,
Persson
,
POG
,
Schmale
,
J
,
Uttal
,
T
,
Althausen
,
D
,
Angot
,
H
,
Archer
,
S
,
Bariteau
,
L
,
Beck
,
I
,
Bilberry
,
J
,
Bucci
,
S
,
Buck
,
C
,
Boyer
,
M
,
Brasseur
,
Z
,
Brooks
,
IM
,
Calmer
,
R
,
Cassano
,
J
,
Castro
,
V
,
Chu
,
D
,
Costa
,
D
,
Cox
,
CJ
,
Creamean
,
J
,
Crewell
,
S
,
Dahlke
,
S
,
Damm
,
E
,
de Boer
,
G
,
Deckelmann
,
H
,
Dethloff
,
K
,
Dütsch
,
M
,
Ebell
,
K
,
Ehrlich
,
A
,
Ellis
,
J
,
Engelmann
,
R
,
Fong
,
AA
,
Frey
,
MM
,
Gallagher
,
MR
,
Ganzeveld
,
L
,
Gradinger
,
R
,
Graeser
,
J
,
Greenamyer
,
V
,
Griesche
,
H
,
Griffiths
,
S
,
Hamilton
,
J
,
Heinemann
,
G
,
Helmig
,
D
,
Herber
,
A
,
Heuzé
,
C
,
Hofer
,
J
,
Houchens
,
T
,
Howard
,
D
,
Inoue
,
J
,
Jacobi
,
H-W
,
Jaiser
,
R
,
Jokinen
,
T
,
Jourdan
,
O
,
Jozef
,
G
,
King
,
W
,
Kirchgaessner
,
A
,
Klingebiel
,
M
,
Krassovski
,
M
,
Krumpen
,
T
,
Lampert
,
A
,
Landing
,
W
,
Laurila
,
T
,
Lawrence
,
D
,
Lonardi
,
M
,
Loose
,
B
,
Lüpkes
,
C
,
Maahn
,
M
,
Macke
,
A
,
Maslowski
,
W
,
Marsay
,
C
,
Maturilli
,
M
,
Mech
,
M
,
Morris
,
S
,
Moser
,
M
,
Nicolaus
,
M
,
Ortega
,
P
,
Osborn
,
J
,
Pätzold
,
F
,
Perovich
,
DK
,
Petäjä
,
T
,
Pilz
,
C
,
Pirazzini
,
R
,
Posman
,
K
,
Powers
,
H
,
Pratt
,
KA
,
Preußer
,
A
,
Quéléver
,
L
,
Radenz
,
M
,
Rabe
,
B
,
Rinke
,
A
,
Sachs
,
T
,
Schulz
,
A
,
Siebert
,
H
,
Silva
,
T
,
Solomon
,
A
,
Sommerfeld
,
A
,
Spreen
,
G
,
Stephens
,
M
,
Stohl
,
A
,
Svensson
,
G
,
Uin
,
J
,
Viegas
,
J
,
Voigt
,
C
,
von der Gathen
,
P
,
Wehner
,
B
,
Welker
,
JM
,
Wendisch
,
M
,
Werner
,
M
,
Xie
,
ZQ
,
Yue
,
F.
2022
.
Overview of the MOSAiC expedition: Atmosphere
.
Elementa: Science of the Anthropocene
10
(
1
):
00060
.
Solomon
,
A
,
Shupe
,
MD
,
Svensson
,
G
,
Barton
,
NP
,
Batrak
,
Y
,
Bazile
,
E
,
Day
,
JJ
,
Doyle
,
JD
,
Frank
,
HP
,
Keeley
,
S
,
Remes
,
T
,
Tolstykh
,
M.
2023
.
The winter central Arctic surface energy budget: A model evaluation using observations from the MOSAiC campaign
.
Elementa: Science of the Anthropocene
11
(
1
):
00104
.
Stramler
,
K
,
Del Genio
,
AD
,
Rossow
,
WB.
2011
.
Synoptically driven Arctic winter states
.
Journal of Climate
24
(
6
):
1747
1762
.
Sturm
,
M
,
Holmgren
,
J
,
Perovich
,
DK.
2001
.
Spatial variations in the winter heat flux at SHEBA: Estimates from snow-ice interface temperatures
.
Annals of Glaciology
33
:
213
220
.
Sturm
,
M
,
Perovich
,
DK
,
Holmgren
,
J.
2002
.
Thermal conductivity and heat transfer through the snow on the ice of the Beaufort Sea
.
Journal of Geophysical Research: Oceans
107
(
C10
):
SHE 19-1
SHE 19-17
.
Trodahl
,
H
,
Wilkinson
,
S
,
McGuinness
,
M
,
Haskell
,
T.
2001
.
Thermal conductivity of sea ice; dependence on temperature and depth
.
Geophysical Research letters
28
(
7
):
1279
1282
.
Urrego-Blanco
,
JR
,
Urban
,
NM
,
Hunke
,
EC
,
Turner
,
AK
,
Jeffery
,
N.
2016
.
Uncertainty quantification and global sensitivity analysis of the Los Alamos sea ice model
.
Journal of Geophysical Research: Oceans
121
(
4
):
2709
2732
.
Weeks
,
W.
2010
.
On sea ice
.
Fairbanks, AK
:
University of Alaska Press
.
Wu
,
X
,
Budd
,
W
,
Lytle
,
V
,
Massom
,
R.
1999
.
The effect of snow on Antarctic sea ice simulations in a coupled atmosphere-sea ice model
.
Climate Dynamics
15
:
127
143
.
Yen
,
YC.
1981
. Review of thermal properties of snow, ice, and sea ice.
Vol. 81
.
U.S. Army Corps of Engineers, Cold Regions Research and Engineering Laboratory
.
Zuo
,
G
,
Dou
,
Y
,
Lei
,
R.
2018
.
Discrimination algorithm and procedure of snow depth and sea ice thickness determination using measurements of the vertical ice temperature profile by the ice-tethered buoys
.
Sensors
18
(
12
):
4162
.

How to cite this article: Sledd, A, Shupe, MD, Solomon, A, Cox, CJ, Perovich, D, Lei, R. 2024. Snow thermal conductivity and conductive flux in the Central Arctic: Estimates from observations and implications for models. Elementa: Science of the Anthropocene 12(1). DOI: https://doi.org/10.1525/elementa.2023.00086

Domain Editor-in-Chief: Detlev Helmig, Boulder AIR LLC, Boulder, CO, USA

Guest Editor: Zoe Courville, US Army Engineer Research and Development Center Cold Regions Research and Engineering Laboratory, Hanover, NH, USA

Knowledge Domain: Atmospheric Science

Part of an Elementa Special Feature: The Multidisciplinary Drifting Observatory for the Study of Arctic Climate (MOSAiC)

This is an open-access article distributed under the terms of the Creative Commons Attribution 4.0 International License (CC-BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. See http://creativecommons.org/licenses/by/4.0/.

Supplementary data