Interactive comment on “ Technical Note : On the use of the mushy-layer Rayleigh number for the interpretation of sea-ice-core data ”

General description: This is a short notice paper which discusses the possible effects of the uncertainties in physical parameter on the mushy-layer Rayleigh number. As the uncertainties, the formulation of permeability of sea ice, brine losses during core extraction, and lack of sampling frequency were mainly taken into account. The formulation of permeability was checked by comparing Freitag et al.(1999) and Eicken et al.(2004). The amount of brine losses during core extraction were assumed to depend on brine volume fraction. Regarding the error caused by lack of sampling frequency, the observational results obtained in the Ross Sea were examined. As a result, it was found that significant error can be caused in thin growing ice and spring ice when sea ice temperature is relatively high. Since the mushy-layer Rayleigh number is closely re-


Introduction
Brine convection desalinates sea ice (Untersteiner, 1968;Notz and Worster, 2009) and drives the exchange of biogeochemical tracers between sea ice and the underlying sea water (e.g.Thomas and Dieckmann, 2010).The convection is driven by the density difference between high-salinity brine in the ice and the sea water underneath.Its onset and strength can be described by a mushy-layer Rayleigh number of the form convection: a diffusive time scale z 2 /κ and an advective time scale zµ/(g∆ρΠ) as determined from Darcy's law with g∆ρ as the driving pressure gradient.If the diffusive time scale is shorter than the advective time scale, brine temperature, and hence, because of the required phase equilibrium, its salt content, will adjust to the surrounding ice; convection stops.The Rayleigh number can also be interpreted as the ratio of the potential energy density g∆ρz available for convection and the energy density κµ/Π that is dissipated during convection through thermal diffusion and internal friction because of the brine's viscosity.
In situ measurements of the salinity evolution of sea ice reveal that indeed, this evolution can be understood by an analysis of the according mush-Rayleigh number: as long as that number is super-critical, the bulk salinity of the ice decreases, once the number becomes sub-critical, the bulk salinity no longer changes.The Rayleigh number is particularly large in the highly permeable region close to the ice interface during ice growth, causing efficient desalination there (Notz and Worster, 2008).Recently the relationship of the mushy-layer Rayleigh number to brine flux has been theoretically studied, using idealized 2-D numerical simulations (Wells et al., 2010) as well as 2-D and 3-D analytical approaches (Rees Jones and Worster, 2013b).Results from these studies were successfully applied to explain salt loss and the convection depth in synthetic sea ice (Wells et al., 2011;Rees Jones and Worster, 2013a).
During the spring transition, the warming induced increase in sea-ice permeability can render the Rayleigh number super-critical throughout the entire thickness of the ice, causing full-depth desalination of the ice and its replenishment with nutrients as suggested by model studies (Vancoppenolle et al., 2010;Jardon et al., 2013;Griewank and Notz, 2013).Hence, Ra is a powerful indicator of brine convection that could in principle be very helpful for the interpretation of ice-core data: the sample-depending quantities ∆ρ and Π depend directly on temperature (T ) and bulk salinity (S), which are the two standard variables that are commonly obtained in ice-core studies.Indeed, ice core derived Ra has been provided as a component of the physical context in recent biogeochemical sea ice studies (see, e.g.Carnat et al., 2013;Zhou et al., 2013).Figures

Back Close
Full However, the interpretation of Ra derived in this manner is not straightforward, because the physical quantities involved in the computation have several possible expressions and because of observational errors in ice core S and T .Hence, conclusions based on the comparison of Ra to a specific threshold value can be misleading.
In the following section, we will describe some of the uncertainties regarding the calculation of the mush Rayleigh number from ice-core data.In Sect.3, we detail how these uncertainties affect the Rayleigh number computations.In Sect.4, we briefly discuss the possible impact of sampling errors.Based on these analyses, we provide some guidance on a possible best-practice for the interpretation of ice-core data based on a Rayleigh number in the concluding Sect. 5.

Uncertainties in the calculation of the Rayleigh number
Apart from the gravitational acceleration g = 9.81 m s −2 , different approaches have been used for all input parameters that go into the calculation of Ra according to Eq. (1).We discuss each of these individually.

Density difference ∆ρ
For the calculation of the density difference between the brine and the underlying sea water obviously the density of both must be known.For the density of sea water, the most precise results are obtained by using the full equation of state of sea water in its most modern formulation TEOS-10 (IOC, SCOR and IAPSO, 2010).For brine, no equally precise method for calculation of its density currently exists for the entire temperature range that is experienced in natural sea ice, primarily due to a lack of sufficiently precise measurements at low temperatures and the accordingly high brine salinity.A Gibbs function, which forms the foundation of modern methods to calculate the equation of state of sea water, was defined down to a temperature of −6 • C and up to a salinity of 120 g kg differences can be carried out with a high degree of accuracy, using the numerical packages that can be downloaded from the TEOS-10 website (e.g.McDougall and Barker, 2011).
For simplicity, the density difference between brine and sea water ∆ρ = ρ br −ρ oc often assumes that brine density linearly depends on the brine-sea water salinity difference, e.g., (2) with β the haline contraction coefficient and ρ oc the sea water density with a salinity S oc and a temperature equal to the freezing point.Defining B = ρ oc β, we get ∆ρ = B(S br − S oc ). (3) Sea water with S oc = 34 g kg −1 and at its freezing point (−1.85 • C) is, following TEOS-10, characterized by β = 0.79 × 10 −4 kg g −1 , ρ oc = 1027 kg m −3 and B = 0.81 kg m −3 (g kg −1 ) −1 .The actual TEOS-10 density difference for brine at −6 • C with a salinity of ∼ 100 g kg −1 compared to such sea water is 52.3 kg m −3 , whereas the linear relation gives 53.4 kg m −3 , e.g., a 2 % difference.Since uncertainties in other parameters are usually much larger, for practical purposes, using this constant B gives sufficiently reliable results.At lower temperatures and higher salinities, the calculations cannot be performed, because the TEOS-10 fit is not well-behaved.Formulation (3) requires the brine salinity to be known.This can be derived from measurements of ice-core temperature, since the brine salinity is practically in phase equilibrium with the in-situ temperature.A fit to experimental data by Assur (1958) is given by Notz (2005): with α = 0.054 • C (g kg −1 ) −1 .This gives significantly different results compared to the third-order fit from −5 • C towards low T 's, reaching nearly a 100 g kg −1 difference at −15 • C (see Fig. 1a).Hence, especially for cold sea ice, the full third order polynomial should be used.

Permeability Π
For the calculation of permeability Π, the porosity must be known.Neglecting the gas content of sea ice, the porosity reduces to brine volume fraction e (or liquid fraction) and can be calculated from known bulk salinity and brine salinity according to e = rS S br + S(r − 1) (6) (Notz, 2005).Here, r = ρ i /ρ br is the pure ice-brine density ratio.Pure ice density can be approximated as ρ i = 916.8− 0.1403T , where T is in degree centigrade (Pounder, 1965), and brine density is estimated using Eq. ( 2).Assuming a linear brine freezing point and r = 1, liquid fraction reduces to: which only slightly overestimates Eq. ( 6) (see This gives values of Π that are up to one order of magnitude larger than the estimate by Freitag (1999) (see Fig. 1c) for brine volumes between 5 and 20 %.For both formulations, a marked permeability increase occurs near e = 5 %, representing the increase in brine connectivity (Golden et al., 1998, see dotted line in Fig. 1c).
To calculate the Rayleigh number for a certain level in the ice, ideally the harmonic mean of the permeability of all ice between that level and the ice-ocean interface is calculated.However, for practical purposes, usually simply the minimum permeability of the ice between a specific level and the ice-ocean interface is used as the deciding value of Π.

Vertical distance z
Early formulations of the mushy-Rayleigh number simply used the ice thickness as the vertical distance z (e.g.Wettlaufer et al., 1997).To account for vertical variability of the mush Rayleigh number, Notz and Worster (2008) introduced a depth dependent Rayleigh number, where z is given as the distance from a specific level in the ice to the ice-ocean interface.

Thermal diffusivity κ
The choice of the correct value for the thermal diffusivity κ is subject to physical interpretation.Notz and Worster (2008)  −7 m 2 s −1 is obtained for sea water at 0 • C and 100 g kg −1 .For brine at lower temperature, only very little data is available to allow for a calculation of κ.In particular, precise measurements of heat conductivity are largely lacking.The choice of a thermal diffusivity as representative for interstitial brine is based on the assumption that the controlling process for the exchange of heat between brine and the surrounding ice is the temperature adjustment of the small mass of the moving brine as it moves through the ice (see Griewank and Notz, 2013, for discussion).
Based on the interpretation that Ra compares the time scales of vertical brine motion to those associated with heat diffusion in sea ice, where sea ice always maintains phase equilibrium with the moving brine, the thermal diffusivity of sea ice is used in the LIM1D model (Vancoppenolle et al., 2010).This thermal diffusivity includes a constant sea ice density (ρ = 917 kg m −3 ), whereas specific heat c (see Bitz and Lipscomb, 1999, Eq. 1) and heat conductivity k (Pringle et al., 2007) depend on S and T ; both vary significantly over the typical temperature range.In turn, the sea ice heat diffusivity decreases by an order of magnitude from −15 • C to the vicinity of the freezing point (see Fig. 1d).One major issue with this formulation of κ is that c and k have uncontrolled asymptotic properties.In particular, k is expected to decrease rapidly near the ice bulk melting temperature, which is not the case in the observation-derived Pringle et al. fit.
In turn, whereas κ should converge to the thermal diffusivity of salt water, the computed κ instead sharply decreases, approaching 10 −9 m 2 s −1 .
A third option for the choice of the representative thermal diffusivity is that of pure ice, based on the assumption that the surrounding sea ice does not maintain phase equilibrium with the moving brine.Finally, one could also use a hybrid diffusivity, based on the heat conductivity of pure ice and the specific heat and density of the moving brine.Such relationship is obtained from a direct scaling of the advective and diffusive terms in the mushy-layer heat equation.For our analysis, we here only retain the first two options, which span the entire range of possible choices of thermal diffusivity Introduction

Conclusions References
Tables Figures

Back Close
Full
At lower temperatures and higher salinity, this value increases, and an extrapolation of existing equations to brine at T = −6 • C and S = 100 g kg −1 gives a value of 2.5 × 10 −3 kg (m s −1 ) −1 , which is the value used in this note.Note that sometimes, kinematic viscosity ν is used in formulations of the Rayleigh number, which is simply µ × ρ.

Impact of the different physical parameterizations on the Rayleigh number
As detailed in the previous section, there are several possibilities regarding the choice of physical parameters used to compute Ra.We recommend the use of the original set of parameters proposed by Notz and Worster, as (i) this set was used to interpret the results of their experiments, (ii) using another set of parameters does not change the conclusions that can be drawn, as illustrated in the forthcoming paragraphs.
To illustrate the impact of physical uncertainties, we computed Ra using three sets of physical parameters (see Table 1).The NW08 set correspond to the original Notz and Worster ( 2008 We first note that the different formulations used give significantly different results.Ra computed using NW08 and LIM are generally close, with small differences essentially due to different κ.The large S br differences in NW08 and LIM have a small impact on Ra, because they occur below T = −5 • C, which also drastically reduces Π.In contrast, using the permeability of Eicken et al. induces a 2-fold mean increase in Ra compared to NW08, especially strong at high T , inducing large differences in Π.Despite these differences, the shape of the seasonal cycle of the mean is similar for all formulations, with maxima in fall and spring, essentially driven by relatively high temperatures.The vertical shape of Ra profiles is also similar among different formulations for every of the three typical ice regimes found in the data set (fall, winter-early spring, late spring) as illustrated in Fig. 2 (see also Carnat et al., 2013).

Observation and sampling errors
Measurements of ice-core temperature and bulk salinity include observation errors that affect the calculation of the Rayleigh number.In particular, the bulk salinity is underestimated due to brine losses during core extraction in warm and/or high salinity sections, where liquid fraction, permeability and Rayleigh number are highest (Eicken et al., 1991;Notz et al., 2005).For cold ice with low permeability, the sampling error in bulk salinity should be of a few tenths of g kg −1 at most.However, near the ice base, where the temperature is close to the freezing point and permeability is high, brine loss implies an underestimation of the bulk salinity ranging up to 20 g kg −1 compared to non-destructive measurements (Notz et al., 2005).In order to illustrate the impact of salt loss during coring on the Rayleigh number, we added a synthetic correction ∆S to the measured bulk salinity profile, designed for illustrative purposes.lower panels by the black (no correction) and gray (correction) curves.For the typical winter profile shown in Fig. 2b, e, an average salt loss of about 0.3 g kg −1 induces a 0.5 underestimation of the Rayleigh number.At the ice base, the salt loss introduced by the correction is maximum (4.8 g kg −1 ), and accounting for it increases the Rayleigh number by 2.7.This type of error virtually removes any evidence for the actively convective layer that is expected near the ice base.
A representative estimate of Rayleigh number profiles in sea ice is also hindered by the usually low ice core sampling rates, which only gives a snapshot of variables that are highly space and time dependent.Temperature strongly varies at sub-daily and synoptic time scales.For illustration, a second series of Ra computations were done based on high-frequency T and 5 ice-core S profiles from the last 25 days of the Year-Round Ocean-Sea Ice Atmosphere Exchanges (YROSIAE) study on landfast sea ice in the Ross Sea (December 2011-December 2012).More information on the physical settings at YROSIAE will be given in forthcoming papers.The specificity of this data set compared to the previous one is the presence of Ice-T, an automated sea ice buoy recording T at high frequency (15 min) every 5-12.5 cm depth.Ice-T records indicate ample T diurnal differences up to > 5 • C near the ice surface from early November (day 300) onwards (Fig. 4, upper panel).The large temperature range is due to thin snow (a few cm, typical of coastal wind-blown Antarctic fast ice), e.g., too small to provide efficient thermal insulation and protection from sunlight penetration.Subsequent Ra diurnal oscillations, reaching ∼ 3 with the NW-08 formulation (Fig. 4, lower panel), follow temperature forcing once the ice is warm and permeable enough (daily mean T −5 • C, from day 320).Ra diurnal cycle is strongest near the ice surface and is likely larger than estimated with NW-08 due to salt loss during coring (see NW-∆S curves in Fig. 4).Possible high-frequency S variations have not yet been measured.1-D model simulations of sea ice desalination indicate that super-critical Ra values, built up over several days, quickly (< 1-3 days) vanish during full-depth convective events rapidly expelling salt out of sea ice (Vancoppenolle et al., 2010;Jardon et al., 2013).Introduction

Conclusions References
Tables Figures

Back Close
Full Large spatial variations in all sea ice variables from centimetric to kilometric scales are also important, and this is particularly true for S and T , and hence, Ra.Horizontal variations in bulk salinity, related to the presence or absence of brine channels occur at cm scales (Cottier et al., 1999;Gough et al., 2012).At larger-scale, snow-depth (0-100 cm) and ice thickness (0-20 m) can widely vary mostly due to blowing snow (Massom et al., 2001) and ice deformation (Thorndike et al., 1975), respectively, which clearly affects T (Perovich et al., 2003) and S (Eicken et al., 1991).But these spatial variations rather affect the large-scale representativeness of Ra estimations rather than their actual intrinsic validity.

Discussion
We have shown that the various uncertainties in the physical parameters, observational errors, and sampling issues make comparisons of Ra numbers derived from ice core data against any observational threshold imprecise.Errors in Ra are the largest for warm and permeable ice.Uncertainties in permeability, in particular, have a large impact on Ra.Weekly ice core sampling does not resolve significant diurnal variations in Ra that could only be identified from high-frequency vertical T .In this context, it seems difficult to draw firm conclusions on brine convection from the comparison of Ra with any observational threshold (e.g.7; Notz and Worster, 2008), a threshold that is currently not well constrained either.
In contrast, the shape of the vertical Ra profiles, as well as the seasonality of the mean, are only slightly affected by errors and parameterization uncertainties.Hence, relative comparisons can be made, such as brine convective activity is likely more intense in spring than in deep winter, or in winter, the ice is more prone to convection at the base than at the top.While firm statements are not solid at this stage, Ra still provides new information nicely complementing other variables (e.g., T and S, biogeochemical tracer concentrations in sea ice, upper ocean salinity, ...).An important caution is the relatively weak evidence for the basal winter permeable and convective Introduction

Conclusions References
Tables Figures

Back Close
Full layer for all parameterizations, which can be related to salt loss during ice coring.As all parameterizations lead to similar conclusions, we advise to use the original Notz and Worster Ra formulation, for consistency with previous work (see Table 2).
A qualitative evaluation of Ra profiles derived from the CFL ice core data indicates that although convection mostly occurs near the ice base during winter, full depth brine circulation is possible in thin growing ice and in spring, when rising temperatures make the ice more permeable (see Figs. 3 and 4), as argued by Jardon et al. (2013).The latter type of brine convection could exhibit diurnal variations.
A more quantitative analysis of some of the uncertainties described in this paper is currently hindered by a lack of sufficient experimental data.In particular, a better description of brine properties at low temperatures, a better understanding of the evolution of sea-ice permeability, and a more robust quantification of systematic biases related to ice-core studies are very desirable.We hope that the review of the current state of our understanding as summarized by this study will be helpful in guiding such future work.Introduction

Conclusions References
Tables Figures

Conclusions References
Tables Figures

Back Close
Full (bottom) Ra computed with each parameter set (see Table 1).Data from the CFL study in the Beaufort Sea, 2007-2008(Carnat et al., 2013).Introduction

Conclusions References
Tables Figures

Back Close
Full Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | −1 .Hence, down to this temperature, the calculation of density Discussion Paper | Discussion Paper | Discussion Paper | ) However, often a simpler, linear relationship between brine salinity and in-situ temperature is used based on the linearised freezing point relationship of sea water, namely Discussion Paper | Discussion Paper | Discussion Paper | S br = −T/α, Fig. 1b) by a mean absolute difference of ∆e = 0.77 % over the −15 to 0 • C range.The relation between Π and porosity typically neglects gas content, thus the total porosity reduces to e.The dependence of permeability on porosity is not well constrained; two parameterizations are most widely used.Based on laboratory and field experiments, Freitag (1999) proposes Π = 1.995 × 10 −8 • e 3.1 .Discussion Paper | Discussion Paper | Discussion Paper |Based on measurements in first-year Arctic sea ice,Eicken et al. (2004) suggest Π = 4.708 × 10 −14 • exp(76.90e),e ≤ 0.096 3.738 × 10 −11 • exp(7.265e),e > 0.096.(9) use the thermal diffusivity of cold sea water κ = 1.2 × 10 −7 m 2 s −1 , calculated from a heat conductivity k = 0.5 W (K m) −1 , a heat capacity of 4000 J (kg K) −1 and a density of 1027 kg m −3 .A more precise calculation based on the equations given in Sharqawy et al. (2010) gives k = 0.55 W (K m) −1 and hence κ = 1.3 × 10 −7 m 2 s −1 for sea water at 0 • C and 34 g kg −1 .A somewhat higher value of Introduction Discussion Paper | Discussion Paper | Discussion Paper | κ = 1.4×10 Discussion Paper | Discussion Paper | Discussion Paper | ) computation.The second set (NW-E) is the same as NW08 except that the Eicken et al.Π is used instead of Freitag's.The LIM set of parameters is the one used in the LIM1D model (Vancoppenolle et al., 2010) and includes simplifications for S br , e, Freitag's Π and the sea ice diffusivity.These computations were based on vertical profiles of S and T from 31 ice cores from the Circumpolar Flaw Lead (CFL) study in the Beaufort Sea, covering an almost entire sea ice season, from December 2007 to June 2008.More information on the physical settings is given in Carnat et al. (2013).The computed vertical profiles of Ra are shown for three cores, characteristic of fall, winter and spring seasons (Fig. 2), as well as the full time series of the vertically-averaged Ra (Fig. 3).Discussion Paper | Discussion Paper | Discussion Paper | ∆S depends on porosity.It starts from zero at e = 0, quickly increases from e = 4 % on, reaches ∆S/S = 10 % for e = 5 %, and then steadily increases towards higher e values, reaching 40 % at e = 20 %.The corrected bulk-salinity profiles are shown in gray in Fig. 2 (upper panels).The impact of salt loss on the Rayleigh number is illustrated in the Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Table 1 .
The sets of physical parameterizations used to compute Ra. ∆S refers to the correction on ice salinity.