Surface Heat and Moisture Exchange in the Marginal Ice Zone: Observations and a New Parameterization Scheme for Weather and Climate Models

Aircraft observations from two Arctic field campaigns are used to characterize and model surface heat and moisture exchange over the marginal ice zone (MIZ). We show that the surface roughness lengths for heat and moisture over uninterrupted sea ice vary with roughness Reynolds number (  R ; itself a function of the roughness length for momentum, 0 z , and surface wind stress), with a peak at the transition between aerodynamically smooth (  R <0.135) and aerodynamically rough (  R >2.5) regimes. A pre-existing theoretical model based on surface-renewal theory accurately reproduces this peak, in contrast to the simple parameterizations currently employed in two state-of-the-art numerical weather prediction models, which are insensitive to  R . We propose a new, simple parameterization for surface exchange over the MIZ that blends this theoretical model for sea ice with surface exchange over water as a function of sea ice concentration. In offline tests, this new scheme performs much better than the existing schemes for the rough conditions observed during the ‘Iceland Greenland Seas Project ’ field campaign. The bias in total turbulent heat flux across the MIZ is reduced to only 13 W m −2 for the new scheme, from 48 and 80 W m −2 for the Met Office Unified Model and ECMWF Integrated Forecast System schemes, respectively. It also performs marginally better for the comparatively smooth conditions observed during the ‘Aerosol-Cloud Coupling and Climate Interactions in the Arctic ’ field campaign. The new surface exchange scheme has the benefit of being physically-motivated, comparatively accurate and straightforward to implement, although to reap the full benefits an improvement to the representation of sea ice topography via 0 z is required. scalar a calculations of for


Introduction
In the polar regions, sea ice provides a topographic obstacle to airflow, resulting in the generation of turbulence in the atmospheric boundary layer (Notz, 2012). It is also associated with strong gradients in surface properties such as temperature and humidity where it borders water, in what is known as the Marginal Ice Zone (MIZ). The turbulence and strong vertical gradients result in elevated surface fluxes of momentum, heat and moisture, which exert a significant influence on the weather and climate of the mid-to high-latitudes (e.g., Rae et al., 2014;Renfrew, Elvidge, & Edwards, 2019;Pope et al., 2020); on ocean circulation and deep water formation (e.g., Renfrew et al., , 2021Stössel et al., 2008); and on the extent, thickness and transport of sea ice (e.g., Fichefet & Maqueda, 1997;Tsamados et al., 2014;Vavrus & Harrison, 2003). In numerical weather prediction (NWP), climate, and sea-ice models, these surface fluxes must be parameterized.
Momentum flux (or drag) over a rough surface can be partitioned into two components: a skin component due to surface friction, and a form component due to the pressure force exerted on the airflow by roughness elements (Arya, 1973(Arya, , 1975. Over sea ice, roughness elements include the vertical edges of ice floes and melt ponds, pressure ridges, and other topographic features resulting from sea-ice dynamics. In the MIZ, peak values typically occur at ice fractions of 0.5-0.8, where the surface is roughest owing largely to the ice floe edges Elvidge et al., 2016;Held et al., 2011;Lüpkes & Birnbaum, 2005). A series of physically based parameterization schemes have been developed over several decades that provide the neutral drag coefficient, DN C , as the sum of three components: drag from ice-free water ( DNw C ), drag from uninterrupted sea ice ( DNi C ), and the form drag against ice floe edges across the MIZ ( DNf C ; Arya, 1973Arya, , 1975 Abstract Aircraft observations from two Arctic field campaigns are used to characterize and model surface heat and moisture exchange over the marginal ice zone (MIZ). We show that the surface roughness lengths for heat and moisture over uninterrupted sea ice vary with roughness Reynolds number (  R ; itself a function of the roughness length for momentum, 0 z , and surface wind stress), with a peak at the transition between aerodynamically smooth (  R <0.135) and aerodynamically rough (  R >2.5) regimes.
A pre-existing theoretical model based on surface-renewal theory accurately reproduces this peak, in contrast to the simple parameterizations currently employed in two state-of-the-art numerical weather prediction models, which are insensitive to  R . We propose a new, simple parameterization for surface exchange over the MIZ that blends this theoretical model for sea ice with surface exchange over water as a function of sea ice concentration. In offline tests, this new scheme performs much better than the existing schemes for the rough conditions observed during the 'Iceland Greenland Seas Project' field campaign. The bias in total turbulent heat flux across the MIZ is reduced to only 13 W m −2 for the new scheme, from 48 and 80 W m −2 for the Met Office Unified Model and ECMWF Integrated Forecast System schemes, respectively. It also performs marginally better for the comparatively smooth conditions observed during the 'Aerosol-Cloud Coupling and Climate Interactions in the Arctic' field campaign. The new surface exchange scheme has the benefit of being physically-motivated, comparatively accurate and straightforward to implement, although to reap the full benefits an improvement to the representation of sea ice topography via 0 z is required.
Schemes using this framework have been validated and modified using aircraft and tower observations (e.g., Elvidge et al., 2016;Lüpkes & Birnbaum, 2005) and have recently been implemented in a range of numerical models, including the Met Office Unified Model (MetUM; Renfrew, Elvidge, & Edwards, 2019), the European Centre for Medium-Range Weather Forecasts Integrated Forecast System (IFS; Roberts et al., 2018), the Los Alamos Sea Ice Model (CICE; Tsamados et al., 2014), and a regional coupled climate model (HIR-HAM-NAOSIM 2.0; Yu et al., 2020). In the parameterization schemes of NWP models such as the MetUM and IFS, DNf C is typically provided simply as a function of ice fraction, DNi C is fixed, and DNw C is provided by a separate air-sea exchange scheme that is a function of wind speed. Elvidge et al. (2016) demonstrated that this configuration is capable of accurately representing drag across the marginal ice zone, but that variability in DNi C remains a major source of uncertainty which is currently unaccounted for in parameterization schemes. The large variability in DNi C (i.e., sea ice topography) has also been highlighted in other observational studies (e.g., Castellani et al., 2014;Petty et al., 2017).
The theoretical framework for scalar (heat and moisture) exchange over the MIZ has received less attention than that for momentum exchange. On the basis of a surface-renewal model, Andreas (1987; hereafter A87) derived a theoretical conceptual model for scalar exchange over ice and snow surfaces for both aerodynamically rough and smooth regimes. Here, aerodynamic roughness is defined by the roughness Reynolds number-the ratio of characteristic roughness height to the viscous length-scale of the flow: where 0 z is the surface roughness length for momentum,  u is the friction velocity (which is positively correlated with windspeed), and  is the kinematic viscosity of the surface layer flow. 0 z is related to DN C as follows: where  is the von Karman constant (0.4) and z is the reference height at which the exchange coefficient is evaluated (typically 10 m). In a rough regime (   2.5 R ; relatively strong winds over a rough surface), turbulent eddies sweep into the interfacial sublayer, causing air to become trapped between roughness elements for a short period of time. In a smooth regime (   0.135 R ; relatively weak winds over a smooth surface), impinging eddies remain in motion. Accordingly, A87 employs a molecular diffusion model to represent scalar exchange in the rough regime; whilst an advective diffusion model is employed in the smooth regime. A87 also identifies a third, transitional regime. The A87 model provides the ratios of each of the scalar roughness lengths to the momentum roughness length (z z T 0 0 / for heat; z z q 0 0 / for moisture) as functions of  R (see Figure 8 in A87). The scalar exchange coefficients ( HN C for heat and EN C for moisture) and consequently the sensible and latent heat fluxes (SH and LH, respectively) are then provided by Monin Obukhov similarity theory as: where  is air density; U is wind speed;    s and  s q q are the surface-air differences in potential temperature and specific humidity, respectively; p c is the specific heat capacity of air; and V L is the latent heat of vaporization. The skill of the A87 scheme has been demonstrated previously using tower observations located over consolidated sea ice and snow (A87; . Note that the A87 theoretical model was developed for exchange over uninterrupted ice and snow, as opposed to the mixture of ice floes and open water of the MIZ. To represent scalar exchange over sea ice and the MIZ, numerical models have tended to adopt very simple approaches. For instance, contrary to the findings of A87, the MetUM and IFS schemes simply prescribe 0T z and 0q z as proportional to 0 z over uninterrupted sea ice, then use a mosaic technique similar to that expressed in Equation 1 to derive the exchange coefficients across the MIZ. See Lock and Edwards (2013), Renfrew, Elvidge and Edwards (2019) and Roberts et al. (2018) for details of the MIZ surface exchange schemes used in the MetUM and IFS, respectively. Lüpkes and Gryanik (2015) have presented a more complex parameterization framework for HN C across the MIZ using a similar approach to that described above for DN C ; that is, including a component that accounts for scalar exchange due to ice floe edges. In their study, it is noted that 0 0 T z z could be derived using existing parameterizations (including A87), though they do not implement this approach.
In this study, we investigate scalar exchange over sea ice using new aircraft-derived observations focused on the MIZ. At over 230 data points, we believe this is the largest scalar exchange data set currently available for the MIZ. We then develop a novel parameterization scheme to represent scalar exchange over the MIZ, given the momentum exchange, that could be easily implemented in NWP, climate or sea-ice models. Finally, we demonstrate that this new scheme is significantly more accurate than existing schemes.

Acquisition and Processing of Aircraft Observations From Over Arctic Marginal Sea Ice
The observations used for this study come from two Arctic field campaigns: The Aerosol-Cloud Coupling and Climate Interactions in the Arctic (ACCACIA) field campaign of 2013, with observations over the MIZ of the Barents Sea and Fram Strait around Svalbard (Elvidge et al., 2016); and the Iceland Greenland Seas Project (IGP) field campaign of 2018, with observations over the MIZ to the south-east of Greenland (Renfrew et al., , 2021. All observations have been obtained using the British Antarctic Survey's instrumented DH6 Twin Otter research aircraft. For details of the aircraft instrumentation see Fiedler et al. (2010) and Renfrew et al. ( , 2021. In-situ fluxes of momentum, heat, and moisture have been derived from the aircraft observations along low-level (typically ~40 m) flux-runs, that are assumed to be within the atmospheric surface layer, using the eddy covariance technique. Surface exchange coefficients are then derived from these fluxes using the conventional atmospheric surface layer assumptions of constant fluxes and logarithmic vertical profiles of wind, temperature and specific humidity. For model parameterization development, these exchange coefficients have been derived at the standard 10-m reference height and corrected to neutral stratification using the Monin Obukhov length following Elvidge et al. (2016).
The data have undergone careful quality control to detect and discard those flux-runs for which the sampled turbulence is not homogenous. For details of the assumptions and procedures employed in deriving these quantities-including quality control-see Petersen and Renfrew (2009). An estimate of sea-ice concentration (or fraction) along each flux-run has been derived from surface albedo measurements following Elvidge et al. (2016). We have used a flux-run length of ∼9 km, following Elvidge et al. (2016), that is similar to values used in other previous studies. In total, 303 flux runs were available for this study, of which 196 are from ACCACIA and 107 are from IGP. However, quality control reduces this number to 168 + 99, 154 + 88 and 159 + 74, for "good" quality ACCACIA + IGP data for momentum, heat, and moisture exchange respectively. All observations are located over, just upwind of, or just downwind of, the MIZ ( Figure 1).

Surface Fluxes of Heat and Moisture Under Cold Air Outbreak Conditions Across the Arctic MIZ
All observations used in this study were obtained under cold air outbreak conditions-that is, broadly northerly winds associated with predominantly upward and often large (up to 315 W m −2 ) turbulent heat fluxes over and downwind of the MIZ as cool air sourced from high latitudes passes southwards over a relatively warm ocean ( Figure 1). Across both field campaigns, for the observations over the sea just downwind of the MIZ, the interquartile ranges of SH and LH are 60-139 W m −2 and 34-64 W m −2 , respectively, both in the upward direction (not shown).
In Figure 2, SH and LH magnitudes ( SH and LH respectively) and exchange coefficients for heat and moisture from all available data from each field campaign are presented, binned by sea-ice fraction (these bins are denoted hereafter as bin A = x , where x is the ice fraction range). Generally, SH and LH are greatest over the sea and lower sea-ice fractions, and relatively small over higher ice fractions (Figures 2a-2d). This is due to the vertical gradients in temperature and humidity being much greater over the sea than over ice (see Equations 6 and 7). In contrast, the bin-median scalar exchange coefficients are smallest over the sea and greatest for the high ice fractions. For sensible heat exchange during ACCACIA, bin-median HN C values increase monotonically with increasing ice fraction ( Figure 2e

The Parameterization of Surface Heat and Moisture Exchange Over Uninterrupted Sea Ice
The parameterization of surface heat and moisture exchange is dependent on that of surface momentum exchange, since HN C and EN C are functions of 0 z (see Equations 4 and 5). Figures 3a and 3b show observed DN C across the MIZ (in the same format as Figure 2). In the ACCACIA data, median DN C exhibits a peak over the MIZ, at bin A = 0.75-1 , in agreement with previous studies (see Section 1). In the IGP data however, median DN C rises monotonically with ice fraction to peak over uninterrupted sea ice. Furthermore, the median DN C (and consequently 0 z ) is systematically larger in the IGP data than in the ACCACIA data. For ACCACIA and IGP, the median observed 0i z are 4.3 × 10 −4 m and 10 −2 m, respectively. To put these values in context, and to highlight the relative roughness of the ice observed during the IGP campaign, a one-year record of tower observations from the Arctic yielded a median 0i z of 3-6 × 10 −4 m (the precise value being dependent on the sensor height) (Persson et al., 2002). Owing to the higher 0 z and also stronger windspeeds (see Figure 1), the roughness Reynolds number  R is consistently greater in the IGP data set than in the AC-CACIA data set for all ice-fraction bins (cf. Figures 3c and 3d). In other words, the IGP surface flow regime ELVIDGE ET AL.
10.1029/2021JD034827 6 of 16 , R as a function of ice fraction for observations from (a and c) the ACCACIA data set and (b and d) the IGP data set. All panels show observations as box and whiskers plots with the median, interquartile range and 10th/90th percentiles (following Figure 2). Algorithms from the Met Office Unified Model (blue) and Integrated Forecast System (green) are shown in (a and b), with DNw C and DNi C anchored to observations. For each panel, the total number of data points plotted (N) are given. In (c and d), gray background shading marks the transitional aerodynamic roughness regime, below which (   0.135 R ) is the smooth regime, and above which (   2.5 R ) is the rough regime.
is aerodynamically rougher. Note that in contrast to DN C , there is no such systematic difference between the two field campaigns in either HN C or EN C (cf. Figures 2e-2h).  / being broadly insensitive to  R in the smooth regime and decreasing exponentially with increasing  R in the rough regime. Figures 4c-4e show 0T z , 0q z and 0 z , respectively, as functions of  R . Our observations show that 0 z increases monotonically with  R , whilst 0T z and 0q z are broadly positively correlated with  R in the smooth regime, and negatively correlated with  R in the rough regime (i.e., with a peak around the transitional regime). This result is predicted by the A87 theoretical model, but, as far as we are aware, it has never been demonstrated in observations. Our observations also show that  R is a much stronger function of 0 z than it is of  u , as 0 z varies by many orders of magnitude with  R , while  u only varies by about 1 order of magnitude (Figure 4e).  R is a comparatively weak function of  .
The differing sensitivities to  R of the roughness lengths for momentum and for the scalars can be explained physically within the framework of the A87 model. In the smooth regime, roughness elements are embedded in the viscous sublayer and only skin friction is present. Turbulent eddies remain in motion as they scour the ice/snow surface, transferring scalar constituents across the interface by advective diffusion. The absence of form drag means that the Reynolds Analogy holds: that is, both momentum and heat transfer depend on the same turbulent eddies, and so z z T 0 0 / and z z q 0 0 / are constant. In the rough regime, eddies in contact with the surface become stagnant-trapped by roughness elements that extend above the viscous sublayer-and scalar transfer is strictly a diffusion process. Here, momentum exchange is due to both skin friction within the viscous sublayer and form drag above. The addition of pressure forces means the Reynolds Analogy is no longer valid, and momentum exchange is enhanced relative to heat exchange. Consequently, z z T 0 0 / and z z q 0 0 / decrease as aerodynamic roughness-and consequently the influence of pressure forces-increases. Previous studies have highlighted that self-correlation due to the same values of 0 z and  u appearing in both x and y axes of plots such as those in Figure 4 can generate a fictitious correlation (e.g., Andreas, 2002;. To avoid this issue, those studies have derived  R using a bulk parameterization. In Appendix A, we explain why this approach is not appropriate here, and show that self-correlation, whilst present, does not significantly affect our results. Note that the relationships shown in Figures 4a-4d are not significantly affected by varying the ice fraction threshold used to select the data from 0.5 to 0.75. However, as the threshold increases, the number of data points reduces and the scatter increases.

Validation of Existing Model Parameterizations of Scalar Exchange in the MIZ
In Figures 2 and 3, the MIZ surface exchange parameterization schemes used in the MetUM and IFS are plotted over the observations. Since our focus here is exchange over sea ice and the MIZ, for each field campaign we have fixed the exchange coefficients over open water ( DNw C , HNw C and ENw C ) and the equivalent roughness lengths ( 0w z , 0Tw z , and 0qw z ) in the parameterization schemes to the median observed values. Likewise, since our concern is scalar as opposed to momentum transfer, we have also fixed DN C , 0 z to the median observed values at an ice fraction of 1. That is, we prescribe an accurate representation of DNi C (or 0i z ), where subscript i denotes uninterrupted sea ice. The MetUM and IFS schemes use 0i z , 0Ti z , 0qi z and Equations 4 and 5 to derive HNi C and ENi C . The scalar roughness lengths are defined as   Here, these coefficients are used to derive the offline parameterized heat fluxes shown in Figures 2a-2d using observed , U ,  ,  s , q, and s q , and Equations 6 and 7. Note that these offline fluxes are derived using flux-run averaged surface temperatures-not using the mosaic approach applied when these schemes are implemented in a NWP model. This is necessary since extracting separate sea ice and open water values for surface temperature and humidity from our observations over the MIZ is not generally possible, due to differences in the measurement area (or "footprint") of the aircraft-mounted instruments used to derive ELVIDGE ET AL. 10.1029/2021JD034827 8 of 16 ice fraction and surface temperature. Our approach allows an evaluation of the relative performance of the different schemes, given the differences in the exchange coefficients they predict.
For both field campaigns, the MetUM and IFS schemes accurately reproduce the variability in DN C with ice fraction across the MIZ when DNw C and DNi C are fixed to the observations (Figures 3a and 3b). For the ACCACIA data, they also provide a relatively good representation of scalar exchange, with the MetUM performing better than the IFS: the MetUM (IFS) scheme exhibits positive biases in bin-median SH and LH , peaking at 13 (31) and 4 (11) W m −2 , respectively (Figures 2a, 2c, 2e, and 2g). For the IGP data, the schemes significantly overestimate scalar exchange across the MIZ (Figures 2b, 2d, 2f, and 2h). Biases in the exchange coefficients increase with ice fraction, peaking for bin A = 1 where bin-median values are several times greater than observed. Biases in the heat fluxes are greatest at intermediate ice fractions, with bin-median biases for the MetUM (IFS) in SH and LH reaching 41 (88) and 20 (43) W m −2 , respectively, for bin A = 0.25-0.5 . Flux biases are smaller at the higher ice fractions, despite larger biases in the exchange coefficients, because the fluxes are much smaller.
The consistently larger bias in the IFS scheme relative to the MetUM scheme reflects the fact that ), the A87 model predicts a value of k that is in fact closer to that used in the IFS than the MetUM scheme. So in certain conditions, the IFS scheme may be expected to perform better than the Me-tUM scheme.
The reason for the disparity in performance of the existing model parameterization schemes between the two field campaigns is apparent from Figure 4. Over uninterrupted sea ice, these schemes have 0T z and 0q z proportional to 0 z . This is broadly consistent with our observations and the theoretical model of A87 in the smooth regime, but fails to capture the observed decrease in log R of ~2.5; Figure 3d). Our observations confirm A87's theoretical model, demonstrating the dependence of the relationship between momentum exchange and scalar exchange on aerodynamic roughness. This suggests that surface exchange parameterization schemes for the MIZ could be improved if they were to represent this dependency. In Section 6 we derive and test such a scheme.

Derivation and Validation of a New, Aerodynamic-Roughness-Dependent Scalar Exchange Parameterization for the MIZ
Since the A87 theoretical model is only appropriate over sea ice, to provide a useful prediction of scalar exchange across the MIZ we propose a new parameterization approach, based on blending the A87 scheme with an ocean exchange scheme. Our scheme provides the exchange coefficient for heat, blended across the MIZ, as follows:  Figure 4a). An analogous set of equations can be defined for moisture exchange. Note that in the A87 theoretical model, the derivations for heat and moisture differ slightly the derivations for heat and moisture differ slightly (i.e., the relationship given by the pink line in Figure 4a differs slightly from that in Figure 4b). Note that i R is distinct from  R as presented in Figures 3c and 3d, since the observations used for these figures are for a range of ice fractions.
We validate this Blended A87 scheme against our observations using a similar approach to that used for validating the MetUM and IFS schemes in Section 5. Again, 0i z , 0w z , 0Tw z , and 0qw z are set to the median observed values for each field campaign. We also require i R , which is derived from turbulence observations for each flux-run using an adapted version of Equation 2, as follows: where  is kinematic viscosity, derived from observed air temperature and density using Sutherland's formula (see for example Burnett, 1965);  is a stability correction function (see for example Stull, 1988); and  10 z m. In this context, for each flux-run, i R denotes  R over any sea ice and, for each field campaign, i R varies predominantly with U (but also with  and, negligibly, with  ). U , , and  are all taken from our observations. For heat exchange, Equations 10 and 11 are then used to determine HNi C , which is combined with HNw C using Equation 12 to provide a Blended A87 value of HN C for each flux-run. Again, using analogous equations, the same procedure applies for deriving EN C . To compare against observations, the Note that implementation of the Blended A87 scheme in an NWP model would be straightforward, since 0w z , 0Tw z , and 0qw z would be provided by the model's ocean surface exchange scheme; and A, 0i z and i u would be provided by the model's MIZ surface momentum exchange scheme. In this context, for each grid box, i R would denote  R over any sea ice and would vary with U , ,  and, if variable in the model, 0i z . Recall however that 0i z is currently fixed in most weather and climate models. Unlike the MetUM and IFS schemes, the Blended A87 scheme, run offline, reproduces scalar exchange reasonably accurately for the IGP field campaign (Figure 2), except for an overestimate for bin A = 1 (though this overestimate is much smaller than those in the MetUM and IFS schemes, and recall from Section 3 that there is some uncertainty in the observed exchange coefficients in bin A = 1 ). With respect to key statistics, the Blended A87 scheme performs much better in the reproduction of SH and LH for IGP (Table 1; Figure 5), whilst more modest improvements are seen for ACCACIA (Table 1). For IGP, the bias error (root mean square error; R.M.S. error) in total turbulent heat flux for the Blended A87 scheme is 13 (25) W m −2 ; which is considerably lower than the equivalent errors of 48 (63) W m −2 and 80 (103) W m −2 found for the MetUM and IFS schemes, respectively. For ACCACIA, the bias error (R.M.S. error) in total turbulent heat flux for the Blended A87 scheme is 0 (29) W m −2 ; marginally lower than the equivalent error of 3 (33) W m −2 found for the MetUM scheme and significantly lower than the 17 (49) W m −2 found for the IFS scheme.
The large improvement seen for the IGP data is due to the new scheme being sensitive to aerodynamic roughness via i R (i.e., predominantly 0i z and U ), whilst the MetUM and IFS schemes are not. Note that in the ACCACIA data there is little difference in the representation of and ENi C (and consequently also 0Ti z and 0qi z ) in each of the schemes (see Figures 2e and 2g). The more modest improvement in the Blended A87 scheme here is predominantly due to its sensitivity to variability in U across the MIZ, and also possibly to more accurate values of z z

Conclusions
This study has presented aircraft observations from two field campaigns in the Arctic, amounting to an observational data set of unprecedented spatial coverage, for the study of scalar exchange over sea ice and the MIZ. Our observations show large turbulent heat and moisture fluxes over and downwind of the MIZ during cold air outbreak conditions. They corroborate previous studies in showing that scalar exchange varies with ice fraction, surface roughness and wind speed. Where the surface roughness is small relative to the viscous length-scale of the flow (relatively weak winds over smooth sea ice), the roughness Reynolds number (  R ) is small, that is, conditions are aerodynamically smooth. In such conditions over uninterrupted sea ice, the surface roughness lengths for heat and moisture can be accurately described as being proportional to the surface roughness length for momentum, as is prescribed in the parameterization schemes currently used in two state-of-the-art NWP models-the MetUM and IFS. These schemes therefore have the potential to perform relatively well for aerodynamically smooth conditions such as those observed during the ACCACIA field campaign.
However, our observations have revealed that in an aerodynamically rough regime (relatively strong winds over rough sea ice), the scalar roughness lengths are inversely correlated with the roughness length for momentum. Consequently, scalar roughness lengths peak in the transition between the smooth and rough regimes-a finding that, as far as we are aware, has not previously been demonstrated. In the aerodynamically rough sea ice regime-such as observed during the IGP field campaign-the MetUM and IFS parameterization schemes are not suitable and become ever less so with increasing  R . For the IGP data, this leads to overestimates in the median sensible and latent heat fluxes of up to 40 and 22 W m −2 , respectively, in offline tests.  In contrast, the theoretical model of A87 for scalar exchange over uninterrupted sea ice does predict the key relationships between the scalar roughness lengths and  R from our observations over high ice concentrations. Consequently, we have derived a simple new parameterization scheme for scalar exchange across the MIZ by blending output from the A87 model with air-sea exchange as a function of ice fraction. Due to the sensitivity of this new scheme to surface roughness and wind speed, it performs significantly better than the MetUM and IFS schemes for the IGP field campaign, more than halving the bias and R.M.S. errors in total turbulent heat flux in offline tests. The new scheme also performs slightly better for the ACCACIA field campaign, due to its sensitivity to wind speed and possibly more accurate values of z z T 0 0 / and z z q 0 0 / at average wind speeds. Our Blended A87 scheme is more consistent with surface-renewal theory and, according to our observations, provides a more physically realistic treatment of scalar exchange over the MIZ than the existing, operational NWP schemes.
There is a large amount of scatter in our observations, as is expected for measurements of turbulence from a limited aircraft-derived data set. Furthermore, our results are limited to off-ice flow in predominantly statically unstable conditions. Consequently, additional observations in different regions and under different meteorological conditions would be useful in corroborating our findings and testing the veracity of our proposed parameterization. Further testing of the Blended A87 scheme is required using different observational data sets. For example, measurements from a ship would complement our aircraft measurements and, owing to the higher spatial resolution of such measurements, could enable the use of the mosaic approach for the offline derivation of parameterized surface fluxes. Recall that this mosaic approach would be applied if the scheme were implemented in an NWP model, but is inappropriate in offline experiments forced by aircraft observations such as those of the present study (see Section 5).
While our new scheme demonstrates the importance of incorporating sensitivity to aerodynamic roughness over uninterrupted sea ice (via i R ) in the parameterization of scalar exchange over the MIZ, it does not incorporate the effects of ice floe edges. The scheme proposed in Lüpkes and Gryanik (2015) does incorporate these effects via a scalar exchange component similar to the form drag component used for momentum exchange (see Section 1). Their scheme predicts peaks in the scalar exchange coefficients over the MIZ. Our observations support such peaks in three of four cases. However these observed peaks are a result of exchange coefficient values for an ice fraction of 1 that are much smaller than predicted by the A87 model and may be unphysical (see Section 3). Further work is required to verify the physical basis for the incorporation of ice floe edge effects in scalar exchange parameterization. ELVIDGE ET AL.

10.1029/2021JD034827
12 of 16 Note. The data used to derive these statistics are from aircraft measurements and offline experiments in the MetUM, IFS and Blended A87 schemes. The bias, root-mean-square (R.M.S.) error and correlation coefficient from a linear regression are shown. Blue shading indicates which of the three schemes performs best in each statistic for each field campaign. For ACCACIA (IGP), there are 103 (57), 104 (46), and 90 (38) data points used respectively for sensible, latent and turbulent (sensible + latent) heat flux statistics. Data points over open water (ice fraction of 0) and where the observed surface exchange coefficients are anomalously small (below the 10th percentile) are excluded.

Table 1 Parameterisation Performance Statistics for Heat Fluxes Using Data From the ACCACIA and IGP Field Campaigns
It is important to note the skill of the new scheme is dependent on an accurate representation of the sea ice topography (i.e., 0i z ), since i R is a strong function of 0i z . This study corroborates previous work in showing that 0i z varies considerably with region. This variability is a function of several factors including sea-ice type, age, and thickness; and deformation and erosion by wind, waves, and water drainage. It is currently not clear how to represent this variability in models, so existing parameterization schemes usually prescribe a fixed 0i z . This is a major limitation to progress in this field that will require significant further research to address. For example, Elvidge et al. (2016) discuss the possibilities of using sea-ice model output, satellite-derived roughness estimates (Petty et al., 2017) or a stochastic element to sea-ice drag parameterization. The representation of scalar exchange over the MIZ would also benefit from improved model representation of surface roughness over water surfaces within the MIZ (i.e., z0w C ). This is currently typically provided by a model's ocean exchange scheme via some form of Charnock's relation (Charnock, 1955), which is most appropriate for open ocean conditions and less applicable in the MIZ, where the presence of sea ice means that waves are generally dampened relative to in the open ocean, suggesting a fetch dependence may be worth investigating in future work.
Surface fluxes of heat and moisture in the MIZ are known to impact the weather and climate of the high latitudes and the development of weather systems downwind, for example, over Northern Europe (Renfrew, Elvidge, & Edwards, 2019). The volume and spatial coverage of sea ice and marginal ice in the Polar oceans are changing; the MIZ is widening whilst total sea-ice coverage is reducing (Strong & Rigor, 2013). These trends are expected to continue in the future (Sigmond et al., 2018). The potential impacts of these changes on the atmosphere are poorly understood as our current capacity to model the variability in airsea-ice exchanges in the polar regions is lacking. New shipping routes resulting from these changes will also require improved weather forecasts across the Arctic Ocean. Consequently, there is strong motivation to improve the representation of surface fluxes of heat and moisture over sea ice in weather and climate models. The scheme we have developed provides a straightforward framework for doing so, though improved representation of 0i z will be necessary to reap its full benefit.