Thermodynamic and dynamic contributions to seasonal Arctic sea ice thickness distributions from airborne observations

Arctic sea ice is thinning at a rate of about 10% per decade, and its seasonal variability is increasing (e.g., Kwok, 2018; Massonnet et al., 2018; Stroeve and Notz, 2018). The difference in air temperature and solar radiation between polar night and day induces a strong, thermodynamic seasonal cycle in the mean ice thickness, which is amplified or damped by internal feedbacks (Wadhams, 2000; Massonnet et al., 2018). In addition, sea ice dynamics alter the ice thickness through sea ice deformation. Colliding floes form pressure ridges or rafts, and floes breaking apart from each other create leads, where new ice can form. The complex interaction between thermodynamics and dynamics shapes the ice thickness distribution throughout the life cycle of the Arctic sea ice cover. Dynamics create an icescape with heterogeneous ice thicknesses ranging from a few centimeters to several meters (Wadhams, 2000). As most air-sea interactions and biogeochemical processes are sensitive to only a fraction of the thickness range, e.g., the very thin ice or the thick ridges, knowledge of the relative coverage of different ice thicknesses is crucial for the parameterization of important climate processes (Holland et al., 2006; Lipscomb et al., 2007; Notz, 2009). For example, the ocean-atmosphere heat and moisture fluxes are limited by thin ice (e.g., Maykut, 1978). Thick ridges change the air and ocean drag, provide shelter for biota, and increase the likelihood of the sea ice surviving the summer melt (e.g., Castellani et al., 2014; Fernández-Méndez et al., 2018). The large range of ice thicknesses can be summarized in an ice thickness distribution (ITD), i.e., a histogram of all occurring ice thicknesses (Thorndike et al., 1975). Because thermodynamics and dynamics shape the ITD differently, this study used ITDs to distinguish them (e.g., Thorndike et al., 1975; Haas et al., 2008).

Driven by thermodynamics, the large losses of thick, multi-year ice (MYI) increase the seasonality, i.e., the annual differences, in the Arctic sea ice mass balance (Kwok, 2018). Rising air and ocean temperatures have enhanced summer MYI melt while thermodynamic growth in winter has remained slow, leading to a decrease of 50% in the MYI fraction (Kwok, 2018; Stroeve and Notz, 2018, and references therein). The MYI is replaced by seasonal ice, which grows faster in winter because the ocean heat is conducted more efficiently from the ice-water interface to the atmosphere through the reduced snow and ice thickness (ice thickness–growth feedback, e.g. Perovich et al., 2003; Notz and Bitz, 2017). Nevertheless, the thin ice will not get thick enough to survive the summer, and in some regions increased ocean heat flux may even overwhelm the ice thickness–growth feedback (Ricker et al., 2021). With ongoing MYI loss, the Arctic sea ice thickness will be controlled to an increasing extent by seasonal ice processes, which will have widespread consequences for the Arctic ecosystem (Post et al., 2013; Kwok, 2018; Massonnet et al., 2018; Stroeve and Notz, 2018; Krumpen et al., 2019).

The role of dynamics in this increasingly seasonal Arctic sea ice state is not well understood. Recent studies have quantified the dynamic contribution to the sea ice mass balance in several case studies, but there is no evidence yet for more or less dynamic thickening in a warming Arctic (Kwok, 2006; Kwok and Cunningham, 2016; Itkin et al., 2018; Kwok, 2018; von Albedyll et al., 2021b). The decrease in ice thickness and loss of MYI is suspected of changing the ice dynamics by causing more mobile ice and deformation (Rampal et al., 2009; Spreen et al., 2011; Lei et al., 2020; Krumpen et al., 2021). Because thinner ice is easier to break, the observed increase in sea ice deformation, i.e., ridging, rafting, and opening of leads, may be linked to the thinning (Rampal et al., 2009). In contrast, convergent sea ice deformation, i.e., the formation of pressure ridges and rafts, increases the ice thickness locally and may help to maintain a thick ice cover (Kwok, 2015; Itkin et al., 2018). Collectively, the net effect of increasing dynamics on the ice thickness distribution remains unclear. Hansen et al. (2015) analyzed a thickness time series between 1990 and 2011 and found different trends for mean and modal thicknesses, which may hint at profound changes in the sea ice mass balance, raising the question: how will thermodynamic and dynamic contributions change in the warming Arctic? To address such questions, we require models with realistic ice ridging parameterizations. Those parameterizations lack sufficient observational evidence to constrain them well enough (Thorndike et al., 1975; Hibler, 1979; Lipscomb et al., 2007). Previous attempts to compare ridging theory and observations were limited by the difficulties of sampling the same ice at different points in time (e.g., Amundrud et al., 2004; Ungermann and Losch, 2018; von Albedyll et al., 2021b).

Precise knowledge of the drivers of seasonal changes in sea ice thickness and extent is crucial for interpreting the current summer and winter sea ice decline (Cavalieri and Parkinson, 2012; Hansen et al., 2015; Stroeve and Notz, 2018) and for improving predictions of sea ice loss in the 21st century using climate models (Massonnet et al., 2018). Since the contribution of dynamics and thermodynamics varies over the year, an evaluation of the annual cycle can help to improve understanding of their respective changes (Hansen et al., 2015). In addition, climate models that can simulate realistic seasonal growth and melt seem to predict long-term trends better (Massonnet et al., 2018). Previously, most of what we knew about the seasonal cycle was based on observations from drifting ice mass balance buoys, drifting ice camps, moorings, and submarine surveys, which demonstrated large regional and temporal variability in thermodynamic growth and decay (Rothrock et al., 2008; Hansen et al., 2015; Lindsay and Schweiger, 2015; Perovich and Richter-Menge, 2015; Mahoney et al., 2019; Belter et al., 2020). Recently, satellite thicknesses retrievals added large-scale insights but generally lack observations in the melt season (e.g., Ricker et al., 2017; Kwok, 2018; Petty et al., 2020). Few airborne Lagrangian campaigns were carried out over the same ice, and only over periods of a few months at best (Itkin et al., 2018; Lange et al., 2019). Lagrangian studies conducted during longer drift stations were normally collected on foot and therefore the data lack information about very thin and thick ice (Untersteiner, 1961; Haas et al., 2011). Because those extreme thicknesses are the most affected by ice dynamics, we need more Lagrangian, airborne, large-scale observations that cover the full range of sea ice thickness to study the effects of ice dynamics on the ice thickness distribution.

The international Multidisciplinary drifting Observatory for the Study of Arctic Climate (MOSAiC) research project aimed to provide such long-term, Lagrangian observations (Shupe et al., 2020; Nicolaus et al., 2022). Passively drifting with the Transpolar Drift, MOSAiC enabled studying the ice thickness distribution evolution from fall to late summer on several spatial scales (Figure 1). At the same time, the seasonal evolution also represents a journey along the Transpolar Drift from the interior of the ice pack in winter towards the ice edge in summer. Krumpen et al. (2021) and Koo et al. (2021) have presented MOSAiC sea ice thickness observed by satellite sensors. CryoSat-2 radar altimeter thickness increased from 0.77 m in October to 2.40 m in April (Krumpen et al., 2021), and ICESat-2 indicated thicker ice, increasing from 1.48 m in December to 2.56 m in April (Koo et al., 2021). Both studies found evidence for the impact of several dynamic events on the ice thickness in the vicinity of the MOSAiC study region. However, the satellite observations are limited to the winter period, by coarse spatial resolution (CryoSat-2) or rather small spatial coverage (ICESat-2), and uncertainties due to unknown snow depth and snow and ice densities. Therefore, we used a high-resolution airborne electromagnetic-induction (EM) ice thickness data set to study the changes of the ITD and the contributions of dynamic processes in detail. The airborne EM measurements were carried out on scales of up to 80 km. Thus, they can bridge the scales from local floe-scale measurements to Arctic-wide satellite observations. This bridging is an essential prerequisite for up-scaling the local measurements to satellite footprints and model grid cells. With the regional-scale, airborne data set presented here, we also identified differences in the seasonal cycle between the Central Observatory, the ice camp with the main installations within 10-km distance, and the Distributed Network, the hierarchy of autonomously drifting systems within 50-km distance.

Analyzing the temporal evolution of the ice thickness distribution in the Transpolar Drift, we aimed at answering the following questions:

1. How did the ice thickness distribution change from fall to late summer? (Section 3.1)

2. Was the ice thickness distribution in the Central Observatory representative of the wider surroundings? (Section 3.2)

3. How much did thermodynamics and dynamics contribute seasonally to the mean thickness and relative abundance of thin and thick ice in the ice thickness distribution? (Section 3.3)

4. Can we use this data set to evaluate the representation of ice thickness redistribution in modeling? (Section 3.4)

Ice thickness from airborne EM sounding and surface roughness from airborne laser altimetry are the core data sets of this study. If not stated otherwise, all thickness information were retrieved by those systems. They were analyzed with previously established methods, described below, to evaluate the seasonal ice thickness evolution. Supplementary data from ice mass balance and GPS buoys provided insights into the magnitude of the thermodynamic and dynamic contribution to the ice thickness, respectively.

We analyzed 11 surveys along the MOSAiC drift that represent varying ice conditions from the interior of the ice pack to the ice edge. Between October 2019 and July 2020, 10 helicopter-borne surveys were conducted from the Polarstern (Alfred-Wegener-Institut Helmholtz-Zentrum für Polar- und Meeresforschung, 2017) and in September 2020, the IceBird airplane campaign with Polar6 (Alfred-Wegener-Institut Helmholtz-Zentrum für Polar- und Meeresforschung, 2016) surveyed the remaining parts of the Distributed Network (Herber et al., 2021). Due to logistical reasons, the airborne EM system was only operated during the polar day and when Polarstern was close by, which resulted in large data gaps from October 2019 to March 2020 and May 2020 to June 2020. Ground-based thickness surveys are available during winter, but their limited spatial extent of a few kilometers is insufficient to observe the impact of sea ice deformation on the regional ITD when no airborne data are available. The airborne surveys were conducted over three different spatial scales covering (1) the main measurement hub of MOSAiC, called the Central Observatory, and ice within 10-km distance, (2) the extensive network of autonomous systems, called the Distributed Network, within 50-km distance, and (3) beyond 50-km distance, called the Extended Distributed Network (Figure 1; Krumpen and Sokolov, 2020). The majority of this study focuses on the Distributed Network, but in Section 3.2, we compare all three spatial scales. Table 1 contains an overview of the surveys.

One central aspect of MOSAiC was that the surveys were conducted in a Lagrangian sense, i.e., carried out over approximately the same collection of sea ice floes while drifting with the Transpolar Drift. To keep track of the relative motion of the ice and to measure approximately the same ice in every survey, we used drifting buoys (see data accessibility for IDs) as corner points of the Distributed Network surveys (Figures 1a and S1). Changes in the survey area, which covered on average 2557 km², were linked directly to ice dynamics. When the study area shrank, the ice in the study area experienced overall more convergence, and when the study area expanded, divergence dominated within it. Thus, using buoys for marking the study area could eliminate the influence of sea ice advection on the local ITD to a great extent. Figure 1a displays the drift track of the four buoys marking the survey area and Figure 1b shows the survey extents. Figure S1 displays the survey area at 12 different snapshots. The complete time series can be seen as video available from von Albedyll (2021a). There were two exceptions to this survey pattern: (1) In October 2019, the buoy setup was not yet completed, and thus the pattern of the October 14, 2019, flight deviated but covered approximately the same survey area (Figure S1). (2) In September 2020, three of the four buoys had stopped working or drifted out of the ice pack (Figure 1a). We tracked the ice surveyed on September 2, 2020, with five additional buoys located within the survey area on June 21, 2020 (Figure S1). We concluded from the tracks that the ice surveyed in September represented most likely a subset of the June survey area.

While the survey lines were anchored at the buoy positions, the data along each line were not an exact repeat of the same sea ice surfaces of a previous survey due to small-scale ice dynamics. Because the surveys covered only a finite number of ridges, leads, and different ice types, we needed to ensure that the surveys were sufficiently long to represent the surrounding ice as well as possible. To test for adequate representation, we used the established criterion of the standard error that serves as a measure for the expected variability in a statistically homogeneous field (Wadhams, 1997; Rabenstein et al., 2010; Lange et al., 2019). We calculated the subsection length at which the standard error dropped below the 0.1-m measurement uncertainty of the EM instrument (Section 2.1.1) and compared it to the total survey length (Rabenstein et al., 2010). Following that approach, all of the surveys were long enough to be representative. On average, after 38.5 ± 22.9 km (22.0 ± 15.5%) of the survey length, the threshold was reached. We concluded that each survey represents a representative snapshot of the ice conditions in the study area.

The benefit of airborne EM thickness measurements is the combination of large-scale total (snow and sea ice) thickness observation with a high spatial resolution. Suspended below the helicopter or aircraft at 10–20 m above the sea ice surface, the instrument’s footprint and point spacing was about 40–50 m (diameter), and 5 m, respectively. The towed EM instrument retrieved the distance to the ocean-ice interface based on the large conductivity difference between seawater and sea ice. A laser altimeter (Jenoptik LDM301.120) on the instrument provided the distance to the air-snow interface, or, in the case of no snow, air-ice interface. Subtracting these two distances gave the combined snow and ice thickness (Haas et al., 1997; Haas et al., 2009). An anomaly in the sensor causing interference in the EM receiver coil after April 2020 required additional manual post-processing and filtering. Details on the data processing and quality control of the results are provided in Haas et al. (2009) and von Albedyll et al. (2021d). Here we only used data that fulfilled the highest quality criterion marked by a quality flag of 1 (excellent) or 1.5 (very good; details in von Albedyll et al., 2021d).

Over undeformed, level ice, the EM thickness has an accuracy of 0.1 m which was quantified by comparing the EM thicknesses to manual drill-hole measurements (Pfaffling et al., 2007; Haas et al., 2009). However, due to the large footprint of 50 m (diameter) and the detection of seawater in ridge pores, the EM measurement system tends to underestimate the maximal thickness of ridges by up to 50% (Pfaffling et al., 2007). Therefore, we additionally analyzed the surface elevation derived from the laser altimeter located on the EM instrument as described below. They provided complementary information on ridges not affected by the EM-related underestimation of the ridge thickness.

The Jenoptik LDM301.120 laser altimeter (905-nm wavelength) measured the altitude variation at 100 Hz with a spatial resolution of 0.3–0.4 m and an uncertainty of ±2 cm. Altitude variations were caused by the surface roughness and helicopter movements. Influences of the helicopter were removed using a filter method presented by Hibler (1972). It consists of three steps and models the helicopter motion by first applying a high-pass filter, then finding the lowest points in the ice thickness survey and connecting them with straight lines, and finally smoothing the curve with a low-pass filter. The curve can then be subtracted from the original laser measurement to give the surface elevation above level ice. The result was the surface elevation above the level ice. In the filtered data set, ridges were detected when the surface elevation was higher than 0.8 m. This height was chosen to detect as many ridges as possible while avoiding confusion with other topographic features such as snow dunes or sastrugi. Independent ridges were identified if they met the Rayleigh criterion (Wadhams and Davy, 1986). Results from the linear profiling with the laser altimeter likely differ from those detected with 2D scanning lasers as the latter are more likely to pick up a peak in a region where surface features are otherwise sparse (Petty et al., 2016, and references therein). However, Beckers et al. (2015) showed that for survey lengths of several kilometers, the roughness parameters derived from 1D and 2D sampling converged.

The location, spacing, width, and height of the ridges identified by the laser altimeter were summarized as mean surface elevation and mean ridge spacing. They were used to describe the evolution of the ice thickness as complementary parameters to the EM-based parameters as described in the next section.

In this study, we used a set of well-established ice thickness and surface roughness measures to characterize the effects of thermodynamics and dynamics on the ice thickness. Below, we describe how those measures can be interpreted in the context of convergence and divergence.

The presence of pressure ridges is clear evidence of convergence. To detect them in the EM thickness data, we quantified the mean thickness and volume and the area fraction that was covered by ridges thicker than 3 m (see Rack et al., 2021 for details). Ridges also increase the surface roughness, a property that we quantified by the mean surface elevation and mean ridge spacing (e.g., von Saldern et al., 2006; Beckers et al., 2015; Petty et al., 2016). From the laser altimeter data, we calculated mean surface elevation that describes the deviation from the relative level ice height, and mean ridge spacing that is the average distance between the ridges (Section 2.1.1).

Convergence and divergence alike broaden the width of the ITD that the interdecile range quantifies as the difference between the 90% and 10% percentile. In freezing conditions, the open water fraction is a sign of divergence and was defined as the area fraction covered by ice thinner than 0.1 m. Open water was included in calculating the mean thickness to account for the effects of divergence on the mean thickness.

Level ice thickness describes the thermodynamic ice growth, and changes in the area fraction of level ice are indicators of dynamics. Increased level ice fractions point to divergence and new ice formation, whereas ridging often consumes level ice. We classified level ice based on a low, vertical thickness gradient along with the thickness measurements of less than 0.04 that extends continuously for at least 40 m in flight direction (Rabenstein et al., 2010). We excluded the data from September 2020 from our analysis because a large number of small floes made the characterization inaccurate.

The full range of different ice thicknesses is well illustrated in an ITD. We display the distribution of the thicknesses as probability density function (PDF), where the y-axis indicates the probability of occurrence per thickness increment.

The EM device measured combined snow and ice thickness, and thus knowledge of snow coverage is important to understand its contribution to the EM total thickness. Our conclusions on the snow thickness are based on three data sets. First, Krumpen et al. (2020) published snow observations from October 5, 7, 9, and 11, 2019, from several floes in the Distributed Network with a mean snow thickness of 0.10 m and maximum spatial variability of 0.05 m (standard deviation). Second, thermistor-string based sea ice mass balance arrays (SIMBAs) from mid-October until the beginning of July estimated daily snow cover changes on level ice (Section 2.2; Koo et al., 2021). Third, weekly snow thickness measurements on foot in the Central Observatory along a loop, 1 km or longer called transects, covered level and ridged ice from October 31, 2019, to July 26, 2020. They were conducted with a Magnaprobe, an automated snow depth probe equipped with a GPS logger (Sturm and Holmgren, 2018).

Between October 2019 and July 3, 2020, the three snow thickness data sets revealed a mean snow thickness decrease of 0.2–0.3 m. The spatial variability (standard deviation) increased until late spring, peaking at ±0.18 m on May 7, 2020. The transect and SIMBA snow measurements showed strong agreement from autumn to spring. In summer, the SIMBA snow thickness values reached zero in late June, when the transect data still indicated 0.1 m of snow (Section 3.3.2; Webster et al., n.d.). There were several likely reasons for this discrepancy of which two are rooted in the measurement principles of both techniques. On the one hand, solar radiation could heat up the SIMBA chains which led to faster, predominant melt around the chain. On the other hand, the magnaprobe measurement was incapable of distinguishing melting snow from the surface scattering layer in melting conditions because both surfaces are easily penetrable by a metal probe. Here, we refer to the surface scattering layer as granular melting ice with a ‘crumbly’ texture which resembles large-grained melting snow (Light et al., 2008). Therefore, the summertime values in the transect data included the surface scattering layer thickness in addition to the melting snow thickness. In September, the Fram Strait 2020 cruise on board of the Norwegian R/V Kronprins Haakon was in the survey area. Their in situ measurements and IceWatch/ASSIST observations revealed a virtually snow-free sea ice surface for September 2 (personal communication, DV Divine, Norwegian Polar Institute; IceWatch/ASSIST, 2021).

The measurement uncertainties in summer and the large spatial variability in snow thickness made an accurate estimate challenging. Snow accumulation and ablation are complex functions of ice type, deformation history, age, surface roughness, and thickness, and require detailed knowledge about those parameters, which is beyond the scope of this study. Because the temporal change in mean snow thickness was small compared to the large changes of the ice thickness of 1.4 m, we refrained from bulk-correcting the total thickness measured by the EM system for the snow layer. Hereafter, we refer to total thickness as ‘ice thickness’ unless specified differently. We accounted for snow accumulation and ablation only in the thermodynamic contribution, and discuss the related uncertainties in Section 3.6.

The EM thickness observations documented the integrated thermodynamic and dynamic thickness growth. To estimate the thermodynamic contribution individually, we used daily SIMBAs observations (Jackson et al., 2013; Koo et al., 2021). SIMBAs measure the vertical temperature profiles through air-snow-sea ice-ocean. The differentiation of the thermal conductivities of those materials causes different vertical temperature gradients from which the air/snow, snow/ice and ice/ocean interfaces can be identified by either manual (e.g., Lei et al., 2018) or automatic SIMBA algorithms (Cheng et al., 2020) to obtain snow depth and ice thickness. We used manually processed daily snow depth and ice thickness from ten SIMBA buoys observed between October and July (Koo et al., 2021). The SIMBAs were located in the Central Observatory and the Distributed Network. The SIMBAs had a vertical resolution of 2 cm and were deployed on initial snow and ice thicknesses ranging from 0.05 to 0.30 m and 0.4 to 1.70 m, respectively. The buoys were deployed over level ice based on the visual surface reconnaissance.

The good coverage of initial thicknesses allowed us to make a 2D linear interpolation of the ice thickness for each time step. We interpolated and extrapolated the thickness changes at each step linearly with 0.01-m spacing and smoothed the growth/melt rates with a running mean of 0.1 m. We integrated the total thickness change for initial thicknesses between 0.3 and 3 m and obtained the total thickness time series. For ice thinner than the IMB observations, we modeled the ice growth with a thermodynamic model (see below). We assumed that the ice thicker than 3 m did not continue to grow thermodynamically. However, we accounted for snow thickness changes as indicated by the SIMBAs between October and April. Further, we assumed that the ice thicker than 3 m experienced the same ice melt rates as the ice with a thickness of 2.99–3.0 m. The interpolated time series is presented in Section 3.3.

We used those interpolated thickness time series to estimate the undisturbed thermodynamic ice thickness growth from October 2019 to April 2020. The interpolated thicknesses have the advantage that we can account for thickness-dependent differences in ice growth rate. We proceeded as follows. For example, on October 14, 2019, we measured a point with a thickness of 0.6 m. The interpolated SIMBA time series indicated that ice with an initial thickness of 0.6 m grew 1.3 m between October 14, 2019, and April 17, 2020. Therefore, on April 17, 2020, the ice thickness of the point adjusted for the undisturbed thermodynamic growth was 1.9 m. We applied this approach to all measurement points from the survey of October 14, 2019. Eventually, we calculated the mean thickness and an ITD from all thermodynamically adjusted thicknesses. The results are presented in Section 3.3.

We used a thermodynamic model based on Thorndike (1992) and Pfirman et al. (2004) to estimate sea ice thickness growth in four cases: (1) growth of the thin ice not covered by the IMBs from October 14, 2019, to April 14, 2020, with initial thicknesses between 0.0 and 0.3 m (see above); (2) growth of the newly formed ice from September 1, 2019, the earliest estimated time of freeze-up until the first EM survey on October 14, 2019; (3) growth of the newly formed ice from September 13, 2019, the latest estimated time of freeze-up until the first EM survey on October 14, 2019; and (4) growth of the second-year ice (SYI) from its estimated thickness on September 25, 2019, given in Krumpen et al. (2020), up to the time of the first EM survey on October 14, 2019. We applied the same model setup as in Belter et al. (2021), Krumpen et al. (2019), and Krumpen et al. (2021), i.e., with a constant ocean heat flux of 2 W m^–2. Ice growth started with an initial ice thickness of 0.1 m. We forced the model with NCEP/NCAR reanalysis 2-m air temperatures (Kalnay et al., 1996) and snow observations from Krumpen et al. (2020). The 2-m air temperature was extracted along the MOSAiC drift trajectory. For the time before Polarstern arrived in the study area, we used the Lagrangian tracking system ICETrack (Krumpen, 2018; Krumpen et al., 2020) to reconstruct the drift trajectory. We validated the model with individual measurements from the SIMBAs from October 2019 until April 2020 and found good agreement throughout the growth phase. We further used the 2-m air temperature in this study as a proxy to describe the thermodynamic growth conditions.

Marking ice floes with GPS buoys is a well-established method to analyze sea ice dynamics (e.g., Heil et al., 2008; Hutchings et al., 2011; Itkin et al., 2017; Lei et al., 2020). We have limited our analysis to sea ice area changes. Divergence causes an increase in area, while convergence causes a decrease in area. We calculated the area from the four GPS buoys marking the corner points of the EM study area using the python package GeographicLib (Figures 1a and S1; Karney, 2021). As the buoy setup was not yet completed in October 2019 (Section 2.1), we backward-extended the area time series to October 2019 with a second buoy at 14-km distance to the original one. Due to the change in area size, the first part (October 14, 2019 to April 17, 2020) and the second part (April 17 to July 31, 2020) of the area time series were not directly comparable to each other. When the buoys approached the ice edge after July 31, 2020, they drifted quickly apart from each other, resulting in a strongly distorted shape for which area estimates were too uncertain to be useful. Therefore, we restricted the area time series to October 14 to July 31. This approach also ensured that we were excluding data from freely floating GPS buoys. We calculated two metrics based on the time series: net area change and accumulated area change. Net area change was the difference between the start and end of the time series. This metric identifies whether divergence or convergence dominated. We estimated the accumulated area increase and reduction by adding up the daily (hourly) opening and closing motion, treating the two separately. Because deformation events occurred on very short time scales (several hours), the accumulated area change can describe more accurately how much divergence and convergence the ice experienced.

The main aim of this study was to describe the temporal evolution of the ITD from fall to late summer in the Distributed Network with comparisons to the Central Observatory. The contributions of ice dynamics to the ITD and mean thickness are the focus of our analyses.

We characterized the sea ice conditions in fall, late winter, early summer, and later summer based on four large-scale, Distributed Network surveys from October 14, April 17, June 21, and September 2 (bold in Table 1, Figure 2). For simplicity, we refer to winter, spring and summer as the time between the surveys, i.e., winter: October 14 to April 17, spring: April 17 to June 21, and summer: June 21 to September 2. For completeness, we display the full time series of EM surveys in Figure 3, and Table 1. An important aspect while analyzing the data is that they also present a spatial evolution along the Transpolar Drift from the Siberian Arctic into Fram Strait (Figure 1a).

In the Siberian Arctic in fall (October 14), 31% of the ice was level with a few ridges and frozen leads (Figure 3b). The EM survey revealed a rather uniform ice thickness distribution with most ice between 0.4 m and 0.8 m thick (Figure 2a). From its formation history, we knew that the ice consisted of two ice types: (1) new ice formed during fall freeze-up, and (2) SYI that had survived the summer and was advected from upstream of the Transpolar Drift (Krumpen et al., 2020).

We explored how the two ice types differed in their thermodynamic ice thickness by initializing the simple thermodynamic model with the observations of Krumpen et al. (2020; Section 2.2). The model showed that by the time of the survey, the SYI and new ice forming after September 1 or 13, 2019, had reached very similar thermodynamic total thicknesses of 0.5–0.6 m. This thickness fits very well to the observed modal thickness of 0.5 ± 0.1 m on October 14, 2019 (Table 1). The similar thicknesses imply that we cannot separate new ice and SYI on the MOSAiC floe based on their thermodynamic thickness alone. However, other structural properties may help to distinguish between the ice types. For example, Krumpen et al. (2020) described the ice as severely weathered, with only the upper 0.3 m being solid and abundant frozen-over melt ponds. In addition, the new ice and SYI likely differed in their mean thickness. The thick (>3 m) ridges in the tail of the ITD are most likely SYI and increase the mean substantially. Those ridges account for a quarter of the ice volume and doubled the overall observed mean compared to the observed mode (Table 1).

Crossing the Central Arctic, by late winter (April 17), the ice had turned into a heterogeneous ice cover with pronounced ridges and active cracks. Over the winter, the mean total thickness doubled and became more heterogeneous, as seen by the increase of the interdecile range (Figures 2b and 3b). The EM survey showed that ridges became much more abundant, but not necessarily higher (see mean ridge spacing, area fraction, and thickness of ice thicker than 3 m in Figure 3d). The change in level ice thickness is interesting to note: 16% of the level ice from fall was consumed during ridging and 6% new, thin level ice had formed in leads (Figure 3b). Most of this new level ice was 0.9 m ± 0.1 m thick and created a secondary mode in the ITD (Figure 2b).

In early summer (June 21), the MOSAiC expedition had reached the western Nansen Basin, located to the north of Fram Strait. The floes were surrounded by open water patches and shrunk in size to 2–3 km in diameter. As the melt season progressed, melt ponds formed. The EM survey revealed heterogeneous and heavily deformed ice, and reached a maximum thickness in early summer (see mean thickness, interdecile range, level ice occurrence in Figure 3a, b, and d). The increased area and thickness of ridges account for more than a third of the total ice volume (Table 1). The volume of ice thicker than 3 m must be considered as a lower estimate of the ridge volume, because the level ice classification revealed that ridges thinner than 3 m were also abundant. The substantial fraction of deformed ice between 2 and 3 m thick can be seen in Figure 2c.

In Fram Strait in late summer (September 2), small patches of ice were surrounded by open water and some thin, new ice. The preceding five days with air temperatures below zero most likely facilitated new ice formation. The remaining snow was completely melted (personal communication, DV Divine, Norwegian Polar Institute). The airborne EM survey was consistent with those visual observations and revealed a mixture of open water (mode at 0.1 ± 0.1 m) and rotten, often strongly deformed ice (second mode at 0.9 ± 0.1 m; Figure 2d). The interaction of melting and divergent motion created more open water and caused a large drop in the mean ice thickness (Figure 3a). Interestingly, during the melting period, the fraction of very thick ridges (>5 m) stayed almost constant while the thinner ice melted more quickly which explains the increase in mean thickness of ice thicker than 3 m seen in Figure 3d compared to the overall decrease in mean ice thickness (Figure 3a). However, the causes of the increase in ridge spacing (Figure 3d) still need to be explored in detail which might require additional auxiliary data because the late summer survey is only a subset of the study region in June.

We discuss the seasonal cycle of ice thickness in the context of the satellite-derived thickness time series covering a 50-km radius around Polarstern. The EM surveys indicated less thickening over the winter than the CryoSat-2 time series (Hendricks and Ricker, 2020; Krumpen et al., 2021). This difference is likely linked to the uncertainty in the snow load parametrization of the CryoSat-2 retrieval. The CryoSat-2 and EM retrievals agree within their uncertainties in fall, but over the winter, the CryoSat-2 total thicknesses increased more strongly than the EM surveys (Figure 3a), which resulted in temporary differences of up to 0.6 m on April 4, 2020. In mid-April, the CryoSat-2 total thickness retrievals were strongly variable (1.5–3.2 m between April 9 and 19, 2020), but the average thickness of 2.5 m over these ten days agreed with the EM total ice thickness on April 17, 2020. The effects of ice dynamics were seen in the increase of interdecile and interquartile range in both the EM and CryoSat-2 ITDs (Krumpen et al., 2021). The ICESat-2 time series from December 2, 2019, to April 11, 2020, and the EM-Bird time series revealed similar growth rates of 0.779 cm/day and 0.73 cm/day, respectively (Koo et al., 2021). However, ICESat-2 systematically indicates 0.2–0.3 m thicker ice than the EM surveys (Koo et al., 2021).

The seasonal thickness cycle observed here is consistent in the timing and magnitude with those from several other Arctic studies. The peak-to-trough amplitude of 1.4 ± 0.2 m is slightly higher, yet it is within the uncertainty compared to other estimates, which range between 0.6 m and 1.20 m (Hansen et al., 2015; Perovich and Richter-Menge, 2015; Kwok, 2018). Compared to the year-long field experiment ‘Surface Heat Budget of the Arctic Ocean (SHEBA)’ that took place in 1997–1998 on MYI in the Beaufort Sea, our winter growth rates are twice as large, consistent with the initially thinner ice in the MOSAiC region, but summer thinning rates are similar (Perovich et al., 2003). Like Hansen et al. (2015) but in contrast to Rösel et al. (2018), we found the maximum mean thickness in June, although the thermodynamic growth ended in April. For a detailed discussion on the dynamic contribution that caused this discrepancy, we direct readers to Section 3.3. The strong seasonal cycle of thin and thick ice fractions emphasizes how crucial the knowledge of seasonal ITDs is for various seasonal processes, e.g., heat and gas exchange in winter and melt pond formation in summer.

In the previous section, we described the temporal evolution of the ice thickness in the Distributed Network. Many of the spatially limited observations made during MOSAiC need to be tested for their spatial representativeness on larger, pan-Arctic scales. Because ice thickness is crucial for many of the observations, we investigate if the thickness observations carried out in the vicinity of the Central Observatory (<10 km) were representative of the Distributed Network (<50 km) and the Extended Distributed Network (>50 km) by comparing observations from all three spatial scales (Figure 4).

First, we compare the Central Observatory and the Distributed Network in fall and in late winter (Figure 4a and b). In fall, the mean and modal thicknesses were alike, but a detailed look at the ITDs revealed some interesting differences. The Central Observatory consisted of above-average thin ice (0.2–0.8 m) and above-average thick, ridged ice >4 m. The differences were rather small for the thin ice (26% difference in aerial coverage, 56.0% vs. 44.3%) while the thick ice was three times more common in the Central Observatory than in the Distributed Network (3.7% vs. 1.2%). The intermediate ice thicknesses between 0.8 m and 4 m were similar or underrepresented in the Central Observatory compared to sea ice at larger distances (Figure 4a and d). The larger fraction of thicker ice (>4 m) is due to the presence of an area of heavily deformed SYI zone in the center of the Central Observatory, called ‘the fortress’. The differences in ice properties can be seen well in SAR images (Figure 5) and the ice thickness parameters (Figure 3). The fortress favored the logistical choice for the location of the Central Observatory as it promised stable ice for the whole drift experiment.

In late winter, the ice regime in the Central Observatory was characterized by an above-average mean ice thickness of 2.5 m (April 30, 2020) compared to the Distributed Network (2.3 m, April 26, 2020). Strong gains at the lower end of the ITD caused the observed difference combined with the thick deformed sea ice of the fortress.

Second, we compare the Distributed Network with the Extended Distributed Network. To elaborate, we split the large-scale flight on April 4, 2020, into a part covering the Distributed Network and the Extended Distributed Network. We found the Extended Distributed Network to be only 0.03 m thicker, which is below the uncertainty of our measurements. The comparison revealed, however, that the Extended Distributed Network had recently experienced a large divergence event that created thin ice (0.2–0.6 m), which was not found in the Distributed Network survey (Figure 4c and e).

This example of the divergence event in the Extended Distributed Network illustrates well that our comparison is limited by the small number of observations and the continuous spatial and temporal variability induced by the dynamics. Keeping this in mind, we conclude that the differences in mean thickness between the three scales were small but the relative abundance of thin and thick ice varied. We hypothesize that two processes may explain the faster growth of sea ice in the Central Observatory over the winter. (1) Faster thermodynamic ice growth of the originally thinner and likely less snow-covered ice of the Central Observatory resulted in faster ice growth. The simultaneous increase in modal and mean thickness supports this thermodynamic hypothesis. (2) The initially thinner ice experienced more ridging and dynamic thickening than the surroundings. For an enhanced dynamic thickening, we expected the large differences in ridge abundance (e.g., mean ridge spacing, ice thicker than 3 m) to decrease between the Central Observatory and the Distributed Network, which is only partly the case (Figure 3d). This points towards faster thermodynamic growth as the reason for the difference in thickness increase between the Central Observatory and Distributed Network. Based on the limited data, we did not find further evidence for a thin sea ice anomaly in the Distributed Network compared to the wider surroundings as suggested by Krumpen et al. (2020; Krumpen et al., 2021).

The previous sections described the strong seasonal cycle of the mean ice thickness parameters in the Transpolar Drift and addressed spatial variability. Below, we analyze the processes responsible for those changes and separate contributions of dynamic and thermodynamic growth. These aspects are summarized in Table 2 and Figure 6. We used air temperature to compute the thermodynamic growth conditions and present this time series in Figure 7b. To infer the prevailing dynamic conditions, we analyzed the time series of area change as presented in Figure 8. A series of snapshots of the area change is displayed in the supplement (Figure S1) and the video available from von Albedyll (2021a).

In winter, the average air temperature was–25°C and hence strong thermodynamic growth occurred. Figure 7b shows the time series of thermodynamic growth as observed by the SIMBAs (black lines). To estimate the thermodynamic growth of each initial ice thickness on October 14, 2019, we interpolated and extrapolated the SIMBA observations (see Section 2.2, green to blue lines in Figure 7b). For ice thinner than 0.3 m not covered by the SIMBAs, we used a thermodynamic model (gray lines in Figure 7b). Results show that, on average, thermodynamic growth contributed 1.0 m to the observed thickness change. To calculate this number, we estimated the thermodynamic thickness growth for each initial thickness observation between fall and late winter based on the interpolated SIMBAs and calculated a late winter ITD when only thermodynamic growth of the different thickness classes is considered (Section 2.2). The procedure is also shown in Figure 9a which displays the three ITDs: (1) the fall ITD (dark blue), (2) the late winter ITD resulting from the fall ITD when thermodynamic growth of the different thickness classes is considered (light blue), and (3) the observed ITD in late winter (red). Thus, 1.0 m is the mean growth weighted by the observed relative abundance of thin and thick ice. As the thinner ice started to grow earlier and faster than the thicker ice, thermodynamics reduced the differences between thinner and thicker ice (Figure 7b), resulting in a narrower ITD. This effect can be seen with the dark blue (observed, fall) and narrower light blue (thermodynamically adjusted, late winter) ITDs in Figure 9a.

With a mean thickness increase of 1.1 m and undisturbed thermodynamic growth of 1.0 m, we attributed the remaining 0.1 m to deformation-related ice thickening over the winter. Additionally, dynamics caused a more heterogeneous ice cover with a broader ITD with reduced skewness and kurtosis (Figures 6 and 9a). As divergence and convergence have opposing effects on the mean thickness, we considered their contributions to the net mean growth separately. Overall, the time series of area change shows that area opening and closing balanced each other over the five months with a zero net area change (Figure 8).

To separate between the effects of divergence and convergence, we distinguished ice that formed before and after fall based on its thickness. By late winter, undisturbed thermodynamic growth had thickened all ice from fall to at least 1.6 m. In turn, all thinner ice was formed by divergence after fall (the part of the red ITD not cross-hatched in Figure 9a). This thin ice covered 12% of the area and had a mean thickness of 1.1 m. Thus, divergence-induced new ice formation added at least 0.1 m, i.e., 0.12 × 1.1 m, to the mean ice thickness. A significant portion of the new, thin ice had a rather uniform thickness of 0.9 ± 0.1 m forming a secondary mode (Figure 2b). Based on the thermodynamic model (Section 2.2), we estimated that this ice formed between March 14 and April 1, 2020, in one large or several small divergence events in quick succession.

In a Lagrangian reference system, the mean ice thickness increases when the survey area shrinks. While the net area change was close to zero, we found an accumulated area decrease of 59% that describes the isolated effect of convergence (Figure 8, Table 2). To calculate the thickening of the ice that already existed in fall, we compared the thermodynamic equivalent (light blue ITD in Figure 9a) with the part of the observed ITD in late winter that formed before fall, i.e., had a minimum thickness of 1.6 m (cross-hatched part of the red ITD in Figure 9a). Redistribution of the 2.1-m thick thermodynamic equivalent resulted in 2.4-m thick observed ice in late winter, i.e., an increase of 14% relative to 2.1 m (Figure 9a). Therefore, we concluded that 59% of accumulated convergence led to 14% ice thickening. This apparent disagreement indicated that the majority of the accumulated area decrease must have gone into closing of open water leads instead of ridging. This agrees well with the zero net area change over this period.

The accuracy of the two estimates for the effects of divergence and convergence is dependent on the ice that formed in winter but was ridged into ice classes thicker than 1.6 m. To assess this limitation, we evaluated whether the thickness increase of 14% can be explained by the area decrease alone or whether additional ice was necessary. We found that a mean thickness increase from 2.1 m (thermodynamic thickness in late winter) to 2.4 m (thickness in late winter >1.6 m) requires an area reduction by 12% to 88% of its original area, which is identical to the area occupied by the thickness classes >1.6 m (cross-hatched part of red ITD in Figure 9a) in the observed late winter ITD. This comparison leads to two important conclusions: (1) our estimates of the effects of divergence and convergence are robust, and (2) over the winter, ridging predominantly redistributed the ice pre-existing in fall.

In contrast to the winter, the SIMBAs and snow transect measurements indicated that the thermodynamic net contribution in spring was between–0.1 and 0 m. The SIMBAs showed thermodynamic growth of 0.1 m from mid-April until mid-May with constantly cold air temperatures. In mid-May, when air temperatures eventually rose to and stayed near 0°C, thermodynamics quickly turned into melting, and thickness decreased by 0.2 m. The thickness increase in the first half of this period was primarily ice growth at the bottom while melting affected mainly the snow. As described in Section 2.1.3, the SIMBAs likely overestimated the snow melting and the transects suggested that snow decreased only by 0.1 m between May 7 and June 17, 2020 (Figure 7a). The very similar EM modal thicknesses agree with the insignificant thermodynamic change. Although the modal thickness of all ice increased by 0.2 m, the modal thickness of the level ice did not (Figure 3a). Figure 2c shows that the all-ice and level-ice modes differ because deformed ice contributed to about 50% of the mode.

In the absence of any significant thermodynamic change, we attributed the increase in mean thickness of 0.3 m to deformation-related ice production. At first sight, the large dynamic contribution is surprising because divergence dominated the sea ice deformation in spring. In total, the survey area increased by 23% (Figure 8). At the same time, the study region also experienced shear that distorted the original shape of the buoy array (Figure 1b). The effect of shear can also be seen in the distance time series shown in Figure 8 where one buoy shortened its distance to Polarstern while the others moved away.

Quantifying the individual effects of divergence and convergence is very challenging for spring. We note that despite the predominant divergence, the open-water fraction was still below 1%. The slight increase in abundance of ice thinner than 0.6 m (Figure 2c) indicates some new ice formation after divergence (Figure 8).

In contrast, there are several indications of dynamic thickening of the sea ice cover (Figure 3a and d). For example, the mean ridge spacing decreased, and the areal fraction of ice thicker than 3 m increased significantly (Figure 3d). Those changes were created by 11% of accumulated convergence (daily) and most likely shear that acted on the ice (Figure 8).

How was so much dynamic thickening possible despite the predominantly divergent conditions? Based on the ITDs, we suggest a redistribution mechanism with a combination of divergence and convergence, also referred to as the ‘ridging pump.’ Thin ice was formed between April and mid-May in leads and was ridged and rafted shortly after its formation into 2-3 m thick ridges. To evaluate this hypothesis, we compared the loss of ice volume thinner than 2 m and the increase in ice volume thicker than 2 m. Indeed, the observed loss of ice thinner than 2 m can explain only 60% of the gain in ice thicker than 2 m which means that new ice formation must have taken place to explain the difference. This new ice, together with the pre-existing thin ice, was then effectively ridged (Section 3.4). These results indicate that in spring, formation and deformation of thin lead ice was the dominant process shaping the ITD.

From June 21 to September 2, thermodynamic melting was responsible for the majority of the observed thinning. During this time, average air temperatures were positive (0.85°C) and varied between–4 and 2°C. Only five days before the survey on September 2, air temperatures had dropped below 0°C. In these Arctic summer conditions with warm air temperatures, solar radiation, ocean heat flux, and strong melt, little to no new ice formation was expected. In the absence of valid SIMBA observations, we followed the approach of Lange et al. (2019) and estimated the total (ice and snow) melt by subtracting the late summer modal ice thickness from the early summer modal thickness. For the late summer modal thickness, we picked the non-open water mode, i.e., ice mode, at 0.9 ± 0.1 m (Figure 2d). Thus, a rough estimate of the thermodynamic melt is 1.2 m or 0.02 m/day. Mean ice thickness decreased by the same amount. Surprisingly, the thickness of the thickest 10% dropped more slowly from 4.8 m to 4.2 m with 0.01 m/day and the ice thickness of the ice thicker than 3 m even rose (Figure 3d), possibly indicating preferential melting of thin ice. As the late summer survey covered only a part of the early summer survey (Section 2.1), our estimates may still underestimate the overall thinning. For example, ice that had drifted eastwards towards the ice edge and could not be surveyed (Figure 1a) had most likely experienced greater ocean heat flux and will thus have melted faster.

Divergence and ocean waves dominated the last phase of the drift. The loose ice pack with ice concentrations below 60% in July (Krumpen et al., 2021) caused the survey area to spread out. Since convergence only closes open water in such low ice concentrations, we did not expect any significant dynamic thickening during this period.

The contributions of thermodynamics and divergent motion to the observed ITDs are shown in Figure 9c. The figure shows the early summer ITD and an ITD for late summer which was obtained by subtracting uniform melting of 1.2 m (see thermodynamics paragraph) from the early summer ITD. The calculated late summer ITD agrees reasonably well with the observed late summer ITD for all ice thicker than the ice mode of 0.9 ± 0.1 m, but underestimates the presence of ice thinner than the ice mode.

We suggest that mutually reinforcing processes have most likely caused ice thinner than the ice mode to melt faster. One possible process could be dynamics that have fractured the thinner ice along with weaknesses, thus exposing more area to lateral melt, and releasing heat stored in the upper ocean (Arntsen et al., 2015). Melt pond formation has likely contributed as well. Besides the summer connection between dynamics and melting, several studies have also shown that winter dynamics precondition the summer melt by creating zones of enhanced melt at thin lead ice and ridge subduction zones (Perovich et al., 2001; Arntsen et al., 2015; Hwang et al., 2017). Given the location at the end of the Transpolar Drift, the summer decay continued for the MOSAiC ice floe until the complete melt.

We conclude that sea ice dynamics impact the ice thickness from fall to early summer and potentially also throughout the summer decay. They require a comprehensive treatment in sea ice models to gain realistic predictions of the future ice thickness development. The following section is concerned with how this data set could be used to improve the representation of dynamic thickening in sea ice models.

State-of-the-art sea ice models use mainly two parameterizations to describe the effects of sea ice deformation on the ice thickness based on Thorndike et al. (1975; hereafter TH75) whose theory was refined by Lipscomb et al. (2007; hereafter LI07). There have been several suggestions based on theory and observations for the parameterizations (Thorndike et al., 1975; Babko et al., 2002; Amundrud et al., 2004; Lipscomb et al., 2007); however, thorough validation with a sequence of ITDs of the same ice is still missing. We demonstrate how to use this data set to evaluate the redistribution parameterizations in a sea ice model. First, we summarize the redistribution theory and then highlight some aspects of the comparison of theory and observations.

Ridging is parameterized in sea ice models by a set of equations that contain tuning parameters and loosely constrained constants (e.g., Thorndike et al., 1975; Lipscomb et al., 2007). A summary of the core equations and parameters is given in the  appendix. The central element of the redistribution is the redistribution function ψ that describes how ice of a certain thickness category is redistributed into others during ridging and shearing (Thorndike et al., 1975). The redistribution function is the sum of the ice that is lost in ridging and the ice that is redistributed, i.e., forms ridges. The redistribution function increases proportionally to the deformation magnitude.

We calculated the ice redistribution following TH75 and LI07 for late winter to early summer and compared the result to the early summer observations (equations in the  appendix). We initialized the redistribution with the early summer ITD. Figure 10 compares the observed ITD to the model output for TH75 and LI07 integrated between April 17, 2020, and June 21, 2020. We forced the redistribution with the closing and opening rates from the time series of area change (Figure 8) at daily and hourly resolution. Opening and closing amounted to a net deformation rate of 23% where the positive number indicated area increase. Because we considered the whole ITD, our model setup is analogous to one model grid cell with an extent of about 50 × 50 km. We discretized the ITD into regular bins of size 0.2 m from 0 to 8 m where the last bin contained all ice thicknesses between 8 and 999 m. As our analysis showed that thermodynamic growth and melt of the existing ice were negligible for this period, we did not consider them in our model (Section 3.3).

Figure 10a displays the comparison between observations and the model output at the two different temporal resolutions. We note that the model underestimates the redistribution and suggests a much higher open-water fraction. Theory and observations agree particularly well on the fact that only ice up to a thickness of 2 m is subject to ridging (Figure 10d). We found only small differences between hourly and daily resolution. For the hourly resolution, ice was less often redistributed into thicker categories. We hypothesize that either unaccounted divergence-induced new ice formation or underestimating the deformation magnitude could explain the discrepancy between model and observations.

First, we roughly estimated the effect of new ice formation on the redistribution combining the ridging with the thermodynamic growth model (Section 2.2). Substantial thermodynamic growth was found between April 14 and May 15 and in the thinnest ice categories (0–0.4 m). For daily time steps, we calculated new ice formation for the thickness categories 0–0.4 m and added the ice in the ITD. Figure 10b shows that new ice formation in this simple implementation cannot explain the discrepancy alone.

Second, we explored the impact of the deformation magnitude by increasing/decreasing the buoy area change time series in 10% increments from net closing rates of –60% to 78%. We chose this range because it covers the most extreme cases of accumulated divergence (hourly) and accumulated convergence (hourly; Table 2). We found the best fit for a net closing rate of –30% (Figure 10c). The fact that a larger closing rate than observed from the buoy area change is needed to obtain realistic model results may hint at spatial limitations as deformation rates are larger for smaller spatial scales. The buoy area change time series cannot resolve convergence and divergence that balances each other within the study area and thus remains unaccounted for. SAR satellite-based derived deformation may provide higher spatial resolution (1.4 km), but at the downside of lower, daily temporal resolution and a gap between mid-January and mid-March 2020 due to poor satellite coverage (Krumpen et al., 2021). A merged, buoy-satellite product may help overcome this shortage to gain realistic forcing data.

For a net closing rate of –30%, we also compared the TH75 and LI07 parameterizations (Figure 10c and d). We note that LI07 created more ridges between 3-5 m than TH75, which fits well with the observations. We noted that the tail of the observed ITDs, characterized by the e-folding of an exponential fit to all ice thicker than the mode, changed only slightly. As TH75 predicted a smaller e-folding change, the TH75 suggestion for the e-folding agrees better with the observations than LI07.

Our comparison is limited by the one-dimensional character of the EM surveys and the footprint smoothing that underestimates the ridge peak height. Also, in our calculations, thermodynamic growth and decay were treated very simply. A resolution of 50 × 50 km may be too coarse, and as a next step, subsets of the EM surveys could be created. A Lagrangian model (e.g., ICEPACK, Hunke et al., 2021), initialized with the conditions and forced with the in-situ observations, will overcome those limitations and help explore the nature of the redistribution functions in much more detail. Potential questions to explore include the following. Can a more complex treatment of new ice formation explain the pronounced differences between 2 and 3 m? Can deformation on a higher spatial resolution produce realistic results? Could scaling the e-folding in the LI07 parameterization with the deformation magnitude as suggested by von Albedyll et al. (2021b) improve the agreement with the observations?

Figure 6 and Table 2 summarize the thermodynamic and dynamic contributions over the study period. One unexpected finding was the large, positive dynamic net contribution during the divergence-dominated drift regime in spring. In the following, we explore the prevailing conditions and consequences of this phenomenon.

First, we ruled out that more deformation was responsible for the dynamic thickening. Even though the mean area changes were larger in spring than in winter, the accumulated area change was not because we examined a much shorter time. In addition, we noted that the mean thickness was largest (Table 1).

Second, dynamic thickening seems to work more efficiently when paired with sufficient divergence that creates thin, new ice that is easy to deform. Krumpen et al. (2021) confirmed the higher abundance of leads with satellite imagery for this phase of the drift. We suggest that dynamic thickening works more efficiently in a mobile, divergent drift regime as long as new ice formation occurs.

Although we have discussed the dynamic contribution more in a temporal sense, the Lagrangian drift also allows a spatial interpretation. The MOSAiC floe was approaching Fram Strait between late winter and early summer. The distance to the ice edge defined as 15% ice concentration reduced from on average 801 km (winter) to 374 km (spring; Krumpen et al., 2021) and the MYI fraction in the wider surroundings fell from on average 55% to 0% beginning of March (Ye et al., 2016a; Ye et al., 2016b; data not shown; University of Bremen, 2021). The region in which the MOSAiC floe was located during this time showed higher lead fractions in almost all seasons (Krumpen et al., 2021; their figure 15a). Taken together, we suggest that a mobile, divergence-dominated drift regime with tidal motions favors dynamic thickening.

Third, the spring drift differed from the winter in the enhanced sub-daily motion. Because reoccurring opening and closing of leads can increase the sea ice mass balance, it may also contribute to the enhanced dynamic thickening (e.g., Heil and Hibler, 2002; Kwok et al., 2003; Hutchings and Hibler, 2008). To quantify the presence of sub-daily motion, we calculated the accumulated area expansions on an hourly and daily basis and compared them (Table 2). For winter, the hourly accumulated area expansion and reduction is 2.2 times larger than when calculated daily. The larger ratio for spring, i.e., 2.7 (area expansion) and 5.2 (area convergence), points to more sub-daily opening and closing. Stronger tides and inertial motions in spring may have a role in the enhanced dynamic thickening.

What are the consequences of dynamic thickening under those conditions? First, the maximum mean thickness was reached later in the season after the melt onset. This finding is in agreement with Hansen et al. (2015) who reported the same for their time series of mean and modal sea ice thickness in Fram Strait. Second, these results open up interesting perspectives on how dynamic thickening could develop if the warming Arctic implied a more divergent, mobile ice regime. Modeling studies are required to explore such scenarios.

We compared our results to Koo et al. (2021) who analyzed the ICESat-2 mean thickness and the same SIMBAs as here along the MOSAiC track from December 2019 to mid-April 2020. Both studies agreed that the dynamics increased the mean thickness. Koo et al. (2021) found a larger dynamic contribution of 34% (individual SIMBAs) and 43% (spatial average) while this study indicated 10% for roughly the same time (October 2019 to April 2020). We attribute the discrepancy mainly to the different thermodynamic growth rates. The rates of Koo et al. (2021) are lower because they did not explicitly consider the ITD while calculating the average of the thermodynamic contribution, which favors thicker ice with lower growth rates. In addition, their study period starts in December, resulting in a different average growth rate.

Our estimates for the dynamic thickening in the Transpolar Drift support previous studies highlighting the importance of dynamic thickening elsewhere in the Arctic. Deformation-related ice production of 0.4 m or 30% of the total thickness change of 1.4 m corroborates the satellite-based estimates (0.3–0.5 m) from Kwok (2006) that were derived from RGPS data covering the Pacific part of the Arctic. They are higher than the extrapolated estimates from storm events north of Svalbard (7% increase in ice volume; Itkin et al., 2018), where a deformation event after a single minor storm was sampled over a confined area. Although those studies examined different locations and seasons, they add to the growing body of research that indicates that the divergence-induced new ice production contributes on the order of 30% to the sea ice mass balance. What is unique about this research is that it quantifies dynamic thickening along the life-cycle of ice in the Transpolar Drift. Thus, we argue that our results are representative of a large portion of the drifting Arctic pack ice. Satellite-based studies like, e.g., Koo et al. (2021) using ICESat-2 will help to extend our knowledge on dynamic thickening to other dynamic regimes, as, e.g., the Beaufort Gyre and over multiple years.

Below we discuss the main sources of uncertainty to address the limitations of our results.

The redistribution and spatial variability of snow is the largest source of uncertainty for the dynamic contribution. Unaccounted snow may have led to an overestimation in winter and an underestimation in spring.

In winter, the largest source of uncertainty is the differing and irregular snow accumulation in space and time over level ice and ridges. Strong winds during storm events can cause abrupt wind erosion followed by rapid snow deposition close to ridges, inducing a strong temporal and spatial variability. Wagner et al. (2021) reported that in the Central Observatory more than 50% of the precipitation was eroded and mostly deposited close to ridges or in the depressions of frozen leads. The transects showed snow accumulation of about 0.2 m while the SIMBAs deployed on level ice indicated only 0.1 m. Thus, we assumed that the difference between SIMBAs and transects reflects the different conditions over level ice (SIMBAs) and all ice (transects). Hence, for the majority of the ice which is level ice with the modal thickness, snow thickness changes over the growth season were as small as 0.1 m. Close to ridges, the accumulation was higher, which is also reflected in the large standard variability of 0.16 m in the transect measurements in this period (April 16, 2020). What are the implications for the estimated dynamic contribution of 0.1 m? Since the SIMBAs most likely underestimated snow accumulation close to ridges, our estimated total ice growth including snow accumulation in the thicker ice categories may be too small. For example, the more frequent occurrence of ice with a total thickness larger than 3 m in spring might not only be related to ridging, but also to snow accumulation in or adjacent to ridges. Therefore, in turn, we may have overestimated the dynamic contribution. The difference between the SIMBA and transect snow accumulation of 0.1 m gives a first rough estimate of the magnitude of this overestimation. Against this estimate we observed that the fraction of newly formed ice (thinner than 1.6 m) adds up nicely to a contribution of 0.1 m to the total mean, which supports our former estimate of a dynamic contribution of 0.1 m (Section 3.3).

In spring, from late winter to early summer, the transects indicate that about 0.1 m of snow was melted. Potentially, this loss was even up to 0.1 m larger, as indicated by the SIMBAs. The spatial variability on the transects was still large (0.18 m, May 7, 2020). Because we assumed a snow accumulation and melt of ± 0.1, any larger snow melt that actually took place reduced our estimate of the dynamic contribution. Thus, for spring, the snow thickness uncertainty may imply an underestimation of the dynamic contribution by up to 0.1 m.

Taken together, the spatial variability of snow has introduced an uncertainty of ± 0.1 m to the ice thickness. Seasonal effects probably offset each other at least in part between winter and spring. The consistent message from parameters dependent and independent of snow, e.g., interdecile range and the ridge spacing, clearly strengthens the validity of our key messages. Snow ablation and accumulation due to wind with respect to surface roughness require a sophisticated model approach (e.g., Liston et al., 2007; Liston et al., 2020) that is beyond the scope of this study but should be an aim of future analyses.

Besides the uncertainties related to the EM footprint (Section 2.1), care should also be taken when interpreting the thickness of the porous SYI in October 2019 (Krumpen et al., 2020) and likely the ice in September 2020 during the break-up. The seawater-filled pores in the ice and the brackish water in the broken-through melt ponds pick up the electromagnetic currents and likely underestimate the true ice thickness, including the rotten part. Thus, the FYI and SYI possibly may have differed in ice thickness in October 2019, but the EM survey could not resolve the difference. Platelet ice that was observed mid-winter beneath the MOSAiC floe and its potential contribution to ridge consolidation (Katlein et al., 2020) could have biased the EM-thickness as well; however, no direct observations exist for the April flights. We also note that on the MOSAiC floe at the end of June, refreezing freshwater in porous ridges and a freshwater lens below parts of the MOSAiC floe were observed with yet unclear consequences for the EM retrieved ice thickness. The absence of conductive seawater would likely lead to thicker ice retrievals and introduce a bias for non-consolidated ridges.

Dynamic and thermodynamic growth depends on the ice thickness and thus influences each other. A simple separation, as presented here, identifies the dominant processes but is not completely accurate. However, from the good agreement between the thermodynamically estimated modal ice thicknesses and the observed level ice thickness (Figure 9), we are confident that errors are smaller than the bin width of 0.2 m. A higher temporal resolution of ice thickness observations could overcome this limitation partly, and we expect a further in-depth analysis of thermodynamic and dynamic contributions with, e.g., ICESAT-2.

For simplicity, we have chosen air temperature and area change as proxies to describe thermodynamic and dynamic processes, deliberately neglecting the influence of humidity, wind, incoming radiation, and ocean heat flux, and shear. As a next step, we suggest using a full numerical sea ice model to explore the full range of interactions of dynamics and thermodynamics.

In the warming Arctic, sea ice is thinning and will likely experience a stronger seasonality. Understanding the consequences of the extensive changes requires knowledge of the relative, seasonal contributions of thermodynamics and dynamics to the sea ice thickness distribution. The MOSAiC research project provided a unique opportunity to collect Lagrangian large-scale, airborne sea ice thickness and surface roughness data. Drifting with the ice, we surveyed approximately the same ice between October 2019 and September 2020 in the Transpolar Drift to analyze the temporal evolution of the sea ice thickness distribution.

The main goal of this study was to describe and quantify the seasonal contributions from thermodynamics and dynamics to the mean sea ice thickness. We have identified a seasonal cycle in mean thickness from 1.1 m (October 2019) to 2.5 m (June 2020) to 1.4 m (September 2020) followed by the complete melt at the end of the Transpolar Drift. The peak-to-trough amplitude and its timing are in good agreement with previous Arctic studies (e.g., Hansen et al., 2015; Perovich and Richter-Menge, 2015). From fall to late winter, we observed the strongest thermodynamic growth. Ridging redistributed the ice already existing in fall, and deformation-related ice production increased the mean thickness by 0.1 m compared to an undisturbed thermodynamic growth. In spring, when the MOSAiC Central Observatory was in the western Nansen Basin, approaching Fram Strait from the north, dynamic and thermodynamic growth was significantly different from the winter period. There was a negligible thermodynamic contribution but 0.3 m of dynamic thickening. The effective ridging of thin lead ice led to this large contribution. We identified the mobile, divergence-dominated drift regime that experienced substantial sub-daily, periodic motion as the most likely mechanism of the larger dynamic contribution. Melting and divergent motion in the marginal ice zone contributed interactively to the decay in thickness and the final break-up of the sea ice cover in summer. Taken together, thermodynamics dominated the Lagrangian sea ice thickness change. Ice dynamics increased the mean thickness from fall to early summer by 0.4 m, i.e., 30% of the total change and enhanced the heterogeneity of the ice thicknesses by creating thin lead ice and thick ridges. The temporal evolution of the ITDs presented here may be of great use for interpreting other MOSAiC time series such as heat, moisture and momentum flux, ocean- and air drag, and biological activity. Extending the analysis to the laser scanner time series on surface elevation, freeboard, and roughness is a crucial next step to complement our results with additional information, especially during the long gap between October 2019 and April 2020.

Many of the spatially limited observations made during MOSAiC need to be scaled up to provide a more complete picture of the Arctic. Therefore, the second aim of this study was to investigate whether the sea ice thickness observations carried out in the vicinity of the Central Observatory (<10 km) were representative of the regional (<50 km, Distributed Network) and large-scale (>50 km, Extended Distributed Network) surroundings. We showed that the differences in mean thickness between the three scales were generally small, but we noted some significant variations in the abundance of thin and thick ice. The Central Observatory consisted of above-average thin ice and above-average thick ridged ice, while the intermediate thicknesses were underrepresented. We speculated that the above-average thin ice grew thermodynamically faster over the winter resulting in a thickness difference of 0.2 m between Central Observatory and Distributed Network in late winter. As possible from the limited data, we did not find a significant difference in mean thickness between the Distributed Network and the Extended Distributed Network. The implications of those observations for up-scaling heat and moisture flux estimates, the internal ice strength, and small-scale sea ice dynamics are subject to further research.

One promising approach for investigating those implications are sea ice models. Because those processes depend heavily on the thin and thick ice fraction, accurate modeling of the ITD is crucial. The temporal evolution of the ITDs presented here may be of great use for sea ice modeling to develop and validate ice thickness redistribution parameterizations. We demonstrated the potential of this data set for future modeling efforts by comparing two commonly used parameterizations of ridging to the observations and identified the e-folding as one potential parameter to improve.

This comprehensive description of the seasonal sea ice thickness properties will serve as a base for several key MOSAiC time series from different disciplines and, thus, advance our knowledge of the Arctic climate and ecosystem. The study also extends our understanding of the magnitude and favorable conditions of the dynamic contribution to the sea ice mass balance. As the presented data set is suitable for modeling those processes, we expect it to provide insights on the ambiguous question of how dynamic thickening may change in the warming Arctic.

We present the set of equations that we used to calculate the redistribution functions, following TH75 and LI07. More details, e.g., normalizing factors, are found in the respective publications.

Changes due to ice redistribution are given by the redistribution function ψ, which describes how ice of a certain category is redistributed into another one during ridging and shearing.

Following LI07, the redistribution function (ψ) is proportional to the deformation, expressed by the strain rate magnitude |ε̇|:

where δ(h) is the rate of opening (delta-pulse at zero thickness) and w_r(h,g) the ridging mode that describes closing and ridging. The parameters α and β are defined by the rheology of the ice and are derived from the respective yield curve (Lipscomb et al., 2007). They distribute ε̇ into contributions from closing and opening.

We used a simplified setup in our experiment where we considered only pure closing and opening. For those extreme cases, the time series of area change (Δ(t)) is a good proxy for the strain rate magnitude.

Thus, for area expansion (pure opening) and α(0) = 1 and β(0) = 0, the equation simplifies to:

and for area decrease (pure closing) and α(π) = 0, and β(π) = 1:

The ridging mode w_r is a function of the ice participating in ridging a(h) and of the ridged ice n(h). A normalization factor N serves for area conservation:

The participating function a(h) = b(h)g(h) describes which ice is removed from each thickness category h during ridging and is thus a function of the original ITD (g(h)). The weighting function b(h) ensures that thinner ice is more likely to participate in ridging than thicker ice. TH75 suggested a linear weighting function b(h) that ensures that leads are closed before thin ice is ridged.

TH75 proposed that only 15% of the ice participates in ridging, which is why b(h) is defined up to a cutoff thickness that is reached if the cumulative thickness distribution exceeds G* = 0.15. LI07 introduced a smooth approximation of the TH75 function to suppress numerical noise which allows all ice to participate while still clearly favoring the thin ice.

with the empirical parameter a* = 0.05.

The ridged ice function n(h) describes the ice that is gained in each thickness category h during ridging and it is defined by:

where h_in denotes the ice thickness of the ice before participating in ridging while h is the ice thickness of the ice after having participated. There are again two choices for the redistribution function γ.

TH75 proposed a uniform redistribution function γ.

with the empirical thickness H_* = 100. This means ridging creates thicknesses between h_min = 2h_in (rafting) and $hmax=20hin$

LI07 developed an exponential redistribution function γ:

They limited the maximum raft thickness and defined the minimum thickness as $hmin=min(2hin,hin+hraft)$ where h_raft = 1 m is the thickest allowed raft. Here, γ₀ is a normalization factor. The e-folding scale depends on the thickness of the ice participating in ridging and μ, which can be used to tune the redistribution function. Recently, von Albedyll et al. (2021b) suggested that the e-folding scale also depends on the deformation magnitude. We chose μ = 2 following Ungermann and Losch (2018).

The processed airborne ice thickness data collected during MOSAiC are available under: https://doi.org/10.1594/PANGAEA.934578 (von Albedyll et al., 2021c). The raw airborne ice thickness and laser scanner data collected during MOSAiC are available under: https://doi.pangaea.de/10.1594/PANGAEA.934814 (Hendricks et al., 2021). The raw airborne ice thickness and laser scanner data collected during IceBird 2020 are available under: https://doi.pangaea.de/10.1594/PANGAEA.924916 (Belter et al., 2020). The snow thickness data collected during MOSAiC are available under: https://doi.pangaea.de/10.1594/PANGAEA.937781 (Itkin et al., 2021). Sea ice extent and autonomous sea ice measurements (drift) from 01-10-2019 to 2-09-2020 were obtained from https://www.meereisportal.de, grant: REKLIM-2013-04 (Grosfeld et al., 2016). We used drift from the following buoys that marked the corner points of the study region: 2020P217_300234068281900, 2019P188_300234068121990, 2019P124_300234066089220, 2019P92_300234067700760, 2019P90_300234067608220. A time series of buoys is available as video from: https://doi.org/10.5446/55704 (von Albedyll, 2021a). SIMBA ice and snow thickness is available under: https://doi.pangaea.de/10.1594/PANGAEA.938244 (Lei et al., 2021). The IDs of the used SIMBAs are: PRIC0901_2019T62, PRIC0902_2019T63, PRIC0903_2019T64, PRIC0904_2019T65, PRIC0905_2019T66, PRIC0906_2019T67, FMI0506_2019T56, FMI0509_2019T58, FMI0601_2019T68, FMI0603_2019T70. TerraSAR-X images shown in Figure 5 were provided by the German Aerospace Center (DLR) and acquired using the TerraSAR- X AO OCE3562_4 (PI: Suman Singha). Multi-year ice fraction was downloaded from the University of Bremen (Ye et al., 2016a; Ye et al., 2016b; last access: 13 August 2021, University of Bremen, 2021).

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

Text S1. Buoy area change.

Figure S1. Changes in study area size and shape marked by buoys.

We thank all MOSAiC participants and people supporting on land for their dedicated work, especially those from the transects (Nixdorf et al., 2021). A special thank you goes to Heli Service, the Polar6, and the Polarstern crew for supporting the airborne measurements. In situ sea ice observations in the survey area conducted during Fram Strait 2020 cruise onboard Norwegian R/V Kronprins Haakon were kindly provided by DV Divine, Norwegian Polar Institute. We thank the editor, the associate editor, and two reviewers for their detailed input that helped us improve the style and clarity of this manuscript.

Data used in this manuscript were produced as part of the international MOSAiC project with the tag MOSAiC20192020, Project ID: AWI_PS122_00. This work was mainly funded by the German Federal Ministry of Education and Research (BMBF) through financing the Alfred-Wegener-Institut Helmholtz-Zentrum für Polar- und Meeresforschung (AWI) and the Polarstern expedition PS122 under the grant N-2014-H-060_Dethloff and the AWI through its projects: AWI_ICE and AWI_SNOW. RL was supported by the National Natural Science Foundation of China (41976219). BC was supported by the European Union’s Horizon 2020 research and innovation program (727890-INTAROS). SA was supported by the German Research Council (DFG) in the framework of the priority program ‘Antarctic Research with comparative investigations in the Arctic ice areas’ (grant nos. SPP1158 and AR1236/1). MW was supported by NASA’s New Investigator Program in Earth Science (80NSSC20K0658). We acknowledge support by the Open Access Publication Funds of Alfred-Wegener-Institut Helmholtz-Zentrum für Polar- und Meeresforschung.

The authors declare that they have no conflict of interest.

Conceptualized the study: LvA, CH, SH.

Carried out the analysis and prepared the original draft: LvA.

Processed the EM data and carried out data quality control: CH, LvA, RG.

Calculated the surface roughness from the laser-altimeter: MS.

Provided ideas in an enriching exchange: LvA, CH, SH, TK, JH, DW, HJB.

Contributed to the data curation: LvA, HJB, TK, SA, GB, MH, RR, JR, AT, MW, PI, JH, RL, MN, BC, DW.

Reviewed and contributed to the final draft: All authors.

How to cite this article: von Albedyll, L, Hendricks, S, Grodofzig, R, Krumpen, T, Arndt, S, Belter, HJ, Cheng, B, Birnbaum, G, Hoppmann, M, Hutchings, J, Itkin, P, Lei, R, Nicolaus, M, Ricker, R, Rohde, J, Suhrhoff, M, Timofeeva, A, Watkins, D, Webster, M, Haas, C. 2022. Thermodynamic and dynamic contributions to seasonal Arctic sea ice thickness distributions from airborne observations. Elementa: Science of the Anthropocene 10(1). DOI: https://doi.org/10.1525/elementa.2021.00074

Domain Editor-in-Chief: Jody W. Deming, University of Washington, Seattle, WA, USA

Associate Editor: Lisa A. Miller, Institute of Ocean Sciences, Fisheries and Oceans, Sidney, BC, Canada

Knowledge Domain: Ocean Science

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