We explore several statistical properties of the observed and simulated Arctic sea-ice lead fraction, as well as the statistics of simulated Arctic ocean–atmosphere heat fluxes. First we show that the observed lead fraction in the Central Arctic has a monofractal spatial scaling, which we relate to the multifractal spatial scaling present in sea-ice deformation rates. We then show that the relevant statistics of the observed lead fraction in the Central Arctic are well represented by our model, neXtSIM. Given that the heat flux through leads may be up to 2 orders of magnitude larger than that through unbroken ice, we then explore the statistical properties (probability distribution function – PDF – and spatial scaling) of the heat fluxes simulated by neXtSIM. We demonstrate that the modelled heat fluxes present a multifractal scaling in the Central Arctic, where heat fluxes through leads dominate the high-flux tail of the PDF. This multifractal character relates to the multi- and monofractal character of deformation rates and the lead fraction. In the wider Arctic, the high-flux tail of the PDF is dominated by an exponential decay, which we attribute to the presence of coastal polynyas. Finally, we show that the scaling of the simulated lead fraction and heat fluxes depends weakly on the model resolution and discuss the role sub-grid-scale parameterisations of the ice heterogeneity may have in improving this result.

The author's copyright for this publication is transferred to the Nansen Environmental Remote Sensing Center, Bergen, Norway.

Sea ice is well known to be an excellent insulator.
In the Arctic, it reduces the potential flux of heat from the ocean to the atmosphere by 2 orders of magnitude in the winter

When a lead opens in the ice during winter, relatively warm ocean waters become exposed to the cold atmosphere, resulting in heat fluxes of up to 600 W m

The transport of sensible heat flux, moreover, is significant over leads and was shown to be even more efficient over smaller than over larger leads

On the ocean side, the formation of new ice in leads removes fresh water and releases brine.
Data from various field campaigns

Leads, therefore, have a potentially significant role to play in the Arctic, through both their impact on the local weather conditions and their long-term influence on the state of the atmosphere and ocean.
Even though their role in the ocean–atmosphere interaction is being actively researched

The fractal characteristics of deformation rates and other quantities are of particular interest because their presence may provide interesting information about the underlying mechanisms at play in the physical system. In the case of deformation of geophysical solids, such as sea ice, a fractal characteristic comes about, at least in part, because of a propagation of fracturing events, which can be modelled by multiplicative cascades

In Sect.

We use the latest version of the next-generation sea-ice model, neXtSIM, presented in

The model set-up covers the central Arctic Ocean, with open boundaries at the Bering Strait, the Canadian Arctic Archipelago, Greenland, and the Barents and Kara seas (see Fig.

We run three simulations, in addition to the control simulation, covering the same space domain and time period.
The forcings are the same in the three cases.
Two of these simulations investigate the effect of changing the model's spatial resolution, and one investigates the effect of changing the model's rheological framework.
We run the simulations related to model resolution at 12.5 and 25 km resolution, with all model parameters kept the same as in the control run.
An exception to this is the cohesion parameter of the Mohr–Coulomb criterion,

In the third simulation, the MEB rheology is replaced with a linear viscous rheology.
The deficiencies of the linear viscous model are well known, and it is suited for neither a detailed study of the model physics nor model evaluation

We analyse the observed lead fraction from

The data show the area fraction of each grid cell that is covered by leads filled with open water, thin ice, or a mixture thereof.
The observations are filtered for feature orientation, and the product, therefore, shows only the fraction of leads, excluding larger, non-linear features such as seasonally thin ice and polynyas.
Although the thickness threshold for thin-ice detection in this product is not known precisely, it is unlikely that this approach classifies ice thicker than about 0.1 m as thin ice

We briefly outline the methodology for investigating the statistical properties of the lead fraction and heat fluxes here. More details can be found in

The second step is to investigate changes in the PDF of both the lead fraction and heat-flux magnitudes (hereafter referred to as heat fluxes) with respect to the scale of observation.
In this study, we focus on the space domain and, therefore, fix the temporal scale of analysis to 1 d.
We choose the daily timescale to consistently compare the simulated and observed lead fraction, and we retain the daily timescale for the sake of simplicity when analysing the heat fluxes.
The scaling analysis consists of evaluating the different moments of the distribution at different spatial scales. The moments of the distribution are calculated as

In the coarse-graining analysis, the mean and higher-order moments are calculated across a range of spatial scales,

Using the spatial coarse-graining method, we derive the moment values as a function of the scale. When presented in log–log space, the moment values should decrease linearly with increasing

For a quantity related to a so-called scale-invariant process, there is generally a monotonic change in

In this section, we demonstrate the capability of neXtSIM to reproduce lead-fraction statistics by comparing the statistical properties of simulated and observed lead fractions.
As the observed lead fraction corresponds to the fraction of open water as well as thin ice

It is important to note that the simulated lead fraction is not strictly a lead fraction as it includes all open-water areas, including polynyas (see Fig.

Observed and simulated daily mean lead fraction on 14 March 2007. The figure shows the entire model domain, and the red lines indicate the boundaries of the Arctic (outer region) and Central Arctic regions used in the study.

Figure

The probability density function of observed and simulated lead fraction in the Central Arctic areas over JFM 2007. The dashed black lines are linear fits discussed in the text.

We now consider the spatial scaling of the lead fraction for the first four moments of the distribution, as the absolute value of the slope of the PDF in the lead-fraction range of

The structure function underlines the good agreement between the simulation and observations: within the estimated uncertainties, the slopes of the observed and simulated structure functions are

We now explore the statistical properties of the heat fluxes simulated by neXtSIM.
Figure

The probability density function of modelled atmospheric heat flux in the Arctic and Central Arctic regions over JFM. The dashed lines are linear fits discussed in the text.

Given the results of the analysis of the PDF, we consider the first three moments for the spatial scaling of heat fluxes, confining the analysis to the Central Arctic region.
The results of this analysis are plotted in Fig.

To demonstrate the role of leads in the scaling obtained for the simulated heat fluxes, we ran the model with a linear viscous rheology, as described in Sect.

In the Arctic region we still obtain an exponential decay with the viscous model. This result supports the idea that this exponential function is the expression of the presence of polynyas in that region, since polynyas are present in both the MEB and viscous simulations.

The probability density functions of modelled atmospheric heat flux in the Arctic and Central Arctic regions over JFM, using the viscous (coloured lines) and MEB (grey lines) rheological models.

To briefly explore the effects of the model resolution on the statistics of the simulated lead fraction and heat fluxes, we ran the model in the same configuration but at 12.5 and 25 km resolution (see Sect.

This result means that running the model at different resolutions will give similar values for the mean, while for the higher-order moments, model runs at different resolutions will give different values. The differences between the runs at 6.25 and 12.5 km resolution are only modest, while going to a 25 km resolution results in significant differences for the third and fourth moments, with a much less clean scaling as well. These differences also result in a change in the structure function, which could be multifractal, rather than monofractal, for the 25 km resolution. The significance of this is limited by the large uncertainties associated with the scaling at 25 km resolution.

The spatial scaling and structure function of modelled lead fraction in the Central Arctic region over JFM 2007.
The two panels are the same as in Fig.

For the simulated heat fluxes there are not only similarities but also clear differences between the results of simulations at different model resolutions (Fig.

The spatial scaling of modelled heat fluxes in the Central Arctic region over JFM 2007.
The two panels are the same as in Fig.

In addition to these differences in the scaling, there also seems to be a difference in the nature of the structure function, depending on the model resolution. The structure function at 6.25 km resolution clearly indicates a multifractal scaling. This is also the case at 12.5 km model resolution, even though the uncertainty is higher and the structure function is smaller. Going from the 12.5 to the 25 km model resolution continues the shift in the structure function going from 6.25 to 12.5 km of decreasing slopes and increasing uncertainty. Consequently, there is low confidence that the structure function indicates a multifractal scaling at the 25 km resolution. This change mirrors the change in fractality for the lead fraction going to the 25 from the 12.5 km resolution.

We have shown that the observed and modelled lead-fraction statistics in the Arctic have a monofractal character, as revealed by a spatial scaling analysis. The fact that the lead-fraction statistics are fractal means that the event propagation mechanism that generates fractal deformation rates also generates a fractal lead fraction – which is to be expected. We note that the lead fraction has a monofractal character, while the deformation rates are multifractal. This means that while large deformations are more localised than small ones, all leads are localised equally. We did not investigate the difference in particular, but a likely explanation is that the lead fraction has too little variation to give multifractal statistics – it is essentially in either an on or an off state. The modelled heat fluxes, on the other hand, do show a multifractal scaling, and this is probably because the heat flux is much more sensitive to the thickness of newly formed ice then the lead fraction.

We compare the fraction of the pixel area covered with open water or thin ice within a lead as estimated from passive microwave satellite data to the modelled fraction of open water and newly formed ice thinner than 0.098 m in neXtSIM. Our analysis shows that the agreement between the two is very good and improves when we take into account the fact that small leads are under-represented in the observations. This under-representation of small leads is most likely the largest source of uncertainty in the present comparison. It is not clear whether the model needs further tuning or enhancement at this point.

The next step would be to compare neXtSIM with other lead-fraction datasets at higher resolution, e.g. from synthetic aperture radar (SAR) images or the Moderate Resolution Imaging Spectroradiometer (MODIS; see

Considering the PDF of simulated heat fluxes, we note that we obtain a linear tail on a log–log scale in the Central Arctic region and a linear tail on a semi-log scale in the Arctic region. Using the results of the viscous model simulations, we show that the log–log tail in the Central Arctic disappears when no leads are present, while much of the semi-log tail in the Arctic region remains. This suggests that the log–log tail in the Central Arctic and the observed scaling of heat fluxes is directly related to the formation of leads there. The semi-log tail in the Arctic can then be related to polynya formation, which occurs in the viscous model.

To explain how polynya formation causes the semi-log tail in the Arctic, we use a simple model of ice growth in a polynya.
This model assumes the polynya to be instantaneously opened by the wind and then closed due to ice growth where all ice that forms inside the polynya is herded by the wind to its downwind edge.
We can then assume that the closing rate is directly proportional to the area of the polynya since the total heat loss and total ice formation can be assumed to be proportional to the area of open water.
An equation for the area evolution of the polynya is then

Since we can assume the heat flux is proportional to the polynya area, to first order, we can expect the PDF of the heat flux to follow the same basic behaviour.

The linear tail on a log–log plot observed in the Central Arctic in Fig.

Having analysed the PDFs of heat fluxes for both the Arctic and Central Arctic regions, we focus on the Central Arctic, where we found a clear multifractal scaling of the simulated heat flux.
It is interesting that the heat-flux scaling is multifractal while the lead-fraction scaling is monofractal.
This difference shows that, in the case of heat fluxes, the high values are also more localised in space, as is the case for sea-ice deformation

The current model does not take lead width into account when calculating the heat flux, as suggested by

We also briefly explored the impact of model resolution on the simulated lead fraction and heat fluxes and found that the statistics of the simulated lead fraction and heat fluxes are affected by the model resolution.
This can, at least partially, be traced to the fact that neXtSIM does not exactly preserve deformation-rate (and in particular the opening-rate) scaling for the second and third moments and the structure function across different model resolutions

The sub-grid-scale thickness distribution used in neXtSIM is a very crude approximation of the highly heterogeneous thickness distribution present in reality. Furthermore, running the model at a high resolution should capture heterogeneity of the ice that would otherwise be parameterised in a coarse-resolution model. As such, even if the deformation scaling were preserved across different model resolutions, one would not necessarily expect the lead-fraction or heat-flux scaling to be preserved. Realistically, the only way to achieve this would be to use a much better parameterisation of the sub-grid-scale thickness distribution.

The main implications of these model shortcomings can be expected to be related to the influence an underestimation of the higher-order moments of the heat fluxes has on the atmosphere and ocean in a coupled set-up. It means that while different model resolutions will deliver the same mean heat fluxes and the same heat-flux statistics on the largest scales, we can expect extreme events to be underestimated in a coarse-resolution model. This in turn may lead to an underestimation of local effects, such as mixing, moisture transport, and brine release – the large-scale impact of which cannot be estimated here.

We have analysed the observed lead fraction in the Arctic using scaling analysis. We then evaluated the simulated lead fraction in neXtSIM against the observed lead fraction.
Following this we investigated the scaling of heat fluxes in neXtSIM and the conservation of lead-fraction and heat-flux scaling across different model resolutions.
The main results of this work are as follows:

The observed lead-fraction statistics in the Arctic have a monofractal scaling in space.

The model reproduces the PDF and scaling of the observed lead fraction in the Central Arctic reasonably well, especially after taking the limits of the observations into account to the extent possible. It is not clear to what extent the differences between simulated and observed lead fractions are due to model or observation deficiencies.

The model shows a clear multifractal scaling in space of the simulated heat fluxes in the Central Arctic. This scaling was shown to originate in the formation of leads.

The mean heat flux is preserved across different model resolutions.

The simulated lead-fraction and heat-flux scaling is not preserved across different model resolutions. This is most likely due in part to the misrepresentation of sub-grid-scale heterogeneity in the current model ice thickness distribution.

These results indicate that the multifractal fracturing mechanism, already identified in sea-ice deformation rates

The lead-fraction dataset was not made publicly available when the original paper was published. The authors will provide copies of these data on request. The neXtSIM source code is not yet publicly available, but the code will be made available upon request, pending usage restrictions.

EO and PR developed the scientific questions and analysis approach. EO performed the model runs and analysis and wrote the bulk of the text with PR and VD giving substantial input on the interpretation of results and the final text.

The authors declare that they have no conflict of interest.

The authors would like to thank Nils Hutter and the anonymous second reviewer for reviewing the paper and making helpful suggestions and comments, as well as the editor Petra Heil and the editorial staff at

This research has been supported by Norges Forskningsråd (grant nos. 263044 and 276730) and the Bjerknes Centre for Climate Research (through the AOI project).

This paper was edited by Petra Heil and reviewed by Nils Hutter and one anonymous referee.