Atmospheric sensitivity to marginal‐ice‐zone drag: Local and global responses

The impact of a physically based parametrization of atmospheric drag over the marginal ice zone (MIZ) is evaluated through a series of regional and global atmospheric model simulations. The sea‐ice drag parametrization has recently been validated and tuned based on a large set of observations of surface momentum flux from the Barents Sea and Fram Strait. The regional simulations are from March 2013 and make use of a collection of cold‐air outbreak observations in the vicinity of the MIZ for validation. The global model analysis uses multiple 48 h forecasts taken from a standard test suite of simulations. Our focus is on the response of the modelled atmosphere to changes in the drag coefficient over the MIZ. We find that the parametrization of drag has a significant impact on the simulated atmospheric boundary layer; for example, changing the surface momentum flux by typically 0.1–0.2 N m−2 (comparable to the mean) and low‐level temperatures by 2–3 K in the vicinity of the MIZ. Comparisons against aircraft observations over and downwind of the MIZ show that simulations with the new sea‐ice drag scheme generally have the lowest bias and lowest root‐mean‐square errors. The wind speed and temperature biases are reduced by up to 0.5 m s−1 and 2 K respectively, compared to simulations with two settings of the previous drag scheme. In the global simulations the atmospheric response is widespread – impacting most of the Arctic and Antarctic sea‐ice areas – with the largest changes in the vicinity of the MIZ and affecting the entire atmospheric boundary layer.

The impact of a physically based parametrization of atmospheric drag over the marginal ice zone (MIZ) is evaluated through a series of regional and global atmospheric model simulations. The sea-ice drag parametrization has recently been validated and tuned based on a large set of observations of surface momentum flux from the Barents Sea and Fram Strait. The regional simulations are from March 2013 and make use of a collection of cold-air outbreak observations in the vicinity of the MIZ for validation. The global model analysis uses multiple 48 h forecasts taken from a standard test suite of simulations. Our focus is on the response of the modelled atmosphere to changes in the drag coefficient over the MIZ. We find that the parametrization of drag has a significant impact on the simulated atmospheric boundary layer; for example, changing the surface momentum flux by typically 0.1-0.2 N m −2 (comparable to the mean) and low-level temperatures by 2-3 K in the vicinity of the MIZ. Comparisons against aircraft observations over and downwind of the MIZ show that simulations with the new sea-ice drag scheme generally have the lowest bias and lowest root-mean-square errors. The wind speed and temperature biases are reduced by up to 0.5 m s −1 and 2 K respectively, compared to simulations with two settings of the previous drag scheme. In the global simulations the atmospheric response is widespread -impacting most of the Arctic and Antarctic sea-ice areas -with the largest changes in the vicinity of the MIZ and affecting the entire atmospheric boundary layer.

INTRODUCTION
The marginal ice zone (MIZ) is the band of partially ice-covered water that separates the ice-free ocean and the main Arctic or Antarctic sea-ice pack. Morphologically it is heterogeneous, consisting of a variety of sea-ice types with broken or rafted floes on a range of spatial scales, interspersed with leads and larger patches of open water. The MIZ is relatively dynamic, responding to the wind, ocean currents and internal dynamics (Notz, 2012), as well as ocean waves (Kohout et al., 2014). In recent years the summertime Arctic MIZ has widened, from ∼100 km to almost 150 km (Strong and Rigor, 2013). Furthermore, this trend is projected to continue over coming decades with most climate models suggesting that a transition to a "new Arctic" is underway; where, in summer, sea-ice conditions for the whole Arctic region will be analogous to an MIZ (e.g. Sigmond et al., 2018). For this reason, determining the sensitivity of the atmosphere to drag imparted by the MIZ is of increasing importance for improving our understanding and predictability of weather and climate in the Arctic.
The atmospheric boundary layer (BL) usually changes rapidly over the MIZ. During off-ice flow there is typically a rapid increase in BL temperature and humidity concomitant with increases in surface turbulent heat fluxes, and a transition from stable to saturated neutral or unstable stratification (e.g. Brümmer, 1996;1997;Renfrew and Moore, 1999). BL cloud development downstream is common (e.g. Young et al., 2016), often in the form of cloud streets or other forms of shallow convection. Numerical simulations show this BL development is sensitive to the distribution of sea ice, with significant differences in temperature, humidity, wind, cloud and surface fluxes over the MIZ and downstream depending on the ice fraction and distribution (e.g. Liu et al., 2006;Gryschka et al., 2008;Chechin et al., 2013). During on-ice flow an internal stable boundary layer will typically develop as the BL is cooled through downward surface heat fluxes (e.g. Overland, 1985). Cloudy boundary-layers are common as the moisture-laden on-ice flow cools. In short, the MIZ marks a zone of change for the BL, as well as the ocean surface; for a broader context see Vihma et al. (2014).
Atmospheric boundary-layer changes across the MIZ are a result of changes in the surface exchange of momentum, heat and moisture. These are challenging quantities to measure directly, but observations are vital to understand BL processes and represent these fluxes in numerical models. Over water, momentum exchange can be represented by the Charnock formula, which provides an aerodynamic roughness length (z 0 ) that is dependent on the friction velocity (u * ) (Charnock, 1955). Over solid sea-ice the momentum flux has often been represented by a constant z 0 or constant neutral drag coefficient (C DN ); observations have suggested C DN between 1.2 and 3.7 × 10 −3 depending on ice morphology (e.g. Overland, 1985;Castellani et al., 2014;Petty et al., 2016). Over the MIZ, momentum flux observations have been scarce, but limited observations have found relatively large C DN values and a peak in C DN over the MIZ -as reviewed in Overland (1985), Lüpkes and Birnbaum (2005), and Elvidge et al. (2016b). The recent study of Elvidge et al. (2016b used low-level aircraft observations to generate over 200 estimates of aerodynamic drag over the MIZ -doubling the number of observations available in the literature at that time. E2016 confirmed there was a peak in drag coefficient for ice fractions of 0.6-0.8. They also illustrated that the morphology of the ice was critical in dictating its roughness; for example, for the relatively small unconsolidated floes of the Barents Sea the median C DN was 2.5 × 10 −3 , while for the larger smoother floes of Fram Strait the median C DN was 1.2 × 10 −3 . The aerodynamic drag over the MIZ is composed of two components: a skin drag and a form drag (Figure 1), the skin drag being due to friction and the form drag due to pressure forces on vertical faces such as floe edges or ridges (Arya, 1973). In numerical weather and climate prediction models this surface drag must be parametrized via a surface exchange scheme. Over sea ice this has traditionally been treated rather crudely, often with one value of C DN for 100% sea ice and so no accounting for ice morphology (i.e. no accounting for differences in form drag). Over partially ice-covered grid boxes a "mosaic method" is commonly employed, which takes a weighted average of the fluxes over water and over ice in proportion to their surface areas (e.g. Vihma, 1995). However, this approach is not appropriate over the MIZ. It results in a linear function of C DN with ice area (A), in contrast to the peak in C DN for A = 0.6-0.8 seen in observations (E2016). In a meticulous series of articles, Lüpkes and colleagues have developed a hierarchy of physically based parametrizations for momentum exchange over sea ice which aim to capture the effects of form drag. Here the total surface exchange is parametrized by where C DN_water and C DN_ice, are the neutral drag coefficients over water and ice and C DNf is the neutral form drag coefficient representing the drag due to ice floe edges. This approach follows from Arya (1973) and is developed and discussed most recently in Birnbaum and Lüpkes (2002), Lüpkes and Birnbaum (2005), Lüpkes et al. (2012), and Lüpkes and Gryanik (2015). Lüpkes et al. (2012) provide details of a hierarchy of schemes, prescribing C DN10 as a function of various sea-ice properties, such as ice fraction and (optionally) freeboard height and floe size. Note C DN10 is the neutral drag coefficient referenced to a height of 10 m. One Lüpkes et al. (2012) scheme is becoming widely adopted (Table 1; Figure 2); for example, it is available as an option in the CICE sea-ice model (Tsamados et al., 2014;CICE Consortium, 2019) and has recently been implemented by us in the Met Office Unified Model (MetUM). This scheme is summarized by E2016, who show that the functional form of this parametrization is a good fit to their large set of aircraft-based observations. E2016 provide recommended values for key parameters for this scheme in their Table 1 and it is the impacts of this new sea-ice drag scheme (hereafter "L12E16") that we investigate here. The surface exchanges of heat and moisture over the MIZ are less well observed than momentum exchange, as often similar ocean and atmospheric temperatures lead to poorly constrained flux estimates (e.g. E2016). In the absence of more observations and improved understanding, most numerical models (including the MetUM) use constant scalar roughness lengths. In this study we leave the scalar roughness lengths at their default value (z 0t = z 0q = 0.2 × z 0_skin , where z 0_skin is a prescribed interfacial roughness length, that is, it does not include any form drag contribution; see Table 1 for After Lüpkes et al. (2012), where parameters are set as the E2016A values in Elvidge et al. (2016b) and C DN10_ice = 1.6 × 10 −3 . Scalar roughness lengths: z 0t_ice = z 0q_ice = 0.1 × 10 −3 m.
GL3 (Walters et al., 2011) GSI4 (Rae et al., 2015) The drag settings are those of the new L12E16 scheme and the previous drag scheme which prescribes a roughness length for the MIZ (z 0_MIZ ) and another for solid sea ice (z 0_ice ), where two prescriptions have often been used, referred to here as Rough and Smooth. Met Office Global Land (GL) configurations that use these setting are noted in column 3. In atmosphere-only configurations surface exchange over the ocean/sea ice is defined in these GL components. In coupled atmosphere-ice-ocean configurations surface exchange over the ocean/sea-ice is defined in the Global Sea Ice (GSI) component - Rae et al. (2015). Note the operational global forecast model has used the Rough settings from prior to 1993 to 2018 (Smith (1996) Neutral drag coefficient (C DN10 ) as a function of ice concentration for the three parametrizations tested: the L12E16 scheme (following Lüpkes et al., 2012;Elvidge et al., 2016b) and the Rough and Smooth roughness lengths (see Table 1). Note C DN_ice is set to 1.6 x 10 −3 . Aircraft-based observations of C DN10 from Fig values). This prescription has some justification from observations (Andreas, 2002;Andreas et al., 2010), but larger MIZ datasets are required for a more sophisticated approach. Consequently, the differences in the heat and moisture exchange come from any differences in z 0_skin and -noting that the scalar exchange coefficients have a dependency on z 0 -the inclusion of form drag for the L12E16 scheme. Note Lüpkes and Gryanik (2015) make some theoretical suggestions for how z 0t and z 0q could vary over the MIZ, but without comprehensive observations to verify their suggestions we have chosen the simpler approach of using default values.
Here we investigate the impact of the L12E16 scheme on the atmosphere through a series of numerical model experiments, that is, we evaluate atmospheric model sensitivity to a physically based sea-ice drag scheme. The model is run in an atmosphere-only configuration for both regional and global domains, so assessing both local and global impacts. We find the impact on the BL is significant, especially on near-surface winds and temperatures, and beneficial. In section 2 we present the details of our approach; sections 3, 4 and 5 describe the surface-layer, BL and global responses to the new scheme, including a quantitative assessment of accuracy; section 6 provides a discussion and conclusions.

The model
The MetUM is a state-of-the-art atmospheric model solving the deep non-hydrostatic compressible equations for a fluid and incorporating numerous parametrizations for physical processes such as radiation, microphysics, boundary-layer turbulence and surface exchange. It is used by the Met Office for operational weather forecasting and as a component in all their climate models (e.g. Walters et al., 2017). Here we use both regional (limited-area) and global atmosphere-only configurations of version 10.6 of the MetUM. Scientific configurations (which are typically available at several successive versions) are also named, numbered and documented according to their function. For example, Walters et al. (2017) describe the Global Atmosphere "GA6" configuration and the Global Land "GL6" configuration that defines the treatment of the land surface and also, in uncoupled simulations without interactive ocean and sea-ice models, surface exchange over open sea and sea ice. Here, global simulations will be based on GA6/GL6. The regional domain covers most of the Greenland, northern Norwegian and western Barents Seas and is centred at about 77 • N, 10 • E, just south of Svalbard; it uses a horizontal grid spacing of 4 km with 70 vertical levels.  Generally physical parametrizations follow those used operationally for the "UKV" configuration (Kendon et al., 2012;Tang et al., 2013). This set-up of the MetUM has proven reasonably accurate at simulating cases of cold-air outbreaks and polar lows in this area (e.g. Sergeev et al., 2017), as well as cases of orographic flows in the Antarctic (Orr et al., 2014;Elvidge et al., 2015;2016a;Elvidge and Renfrew, 2016). The global configuration uses a horizontal grid spacing of ∼40 km (N320) with 70 vertical levels and also generally follows operational settings.

Observational data
The regional simulations are for a 10-day period from 21 to 31 March 2013 during the Aerosol Cloud Coupling and Climate Interactions in the Arctic (ACCACIA) field campaign. This period was dominated by off-ice flow, that is, northerly cold-air outbreak conditions. Low-level research aircraft observations over the MIZ in the Barents Sea and Fram Strait were made on eight separate days (see E2016) and these are used as validation data for our numerical experiments. The low-level (typically 30-40 m) aircraft legs have been divided into 209 runs of 9 km in length. This ACCACIA "surface-layer dataset" includes standard meteorological variables (e.g. temperature, humidity, wind, pressure, etc.), turbulent and radiative fluxes, and infrared-derived sea-surface temperature. The meteorological variables have been interpolated to standard heights using the Coupled Ocean-Atmosphere Response Experiment (COARE) stability-dependent bulk flux algorithm. The turbulent fluxes are estimated from covariances and have been carefully quality controlled (following Petersen and Renfrew, 2009; Cook and Renfrew, 2015); 195 runs have good quality momentum flux values. Figure 3 shows where these observations were obtained, colour-coded by sea-ice fraction (as determined from the aircraft's sea-surface temperature -see E2016, appendix B). The surface-layer dataset covers a full range of sea-ice fractions over two separate locations: in the Barents Sea, to the southeast of Svalbard, where sea ice was characterized by relatively small unconsolidated floes; and in Fram Strait, to the north of Svalbard, where sea ice was characterized by larger smoother floes. In addition to these surface-layer observations, meteorological data from dropsondes provides some observations of BL development over the MIZ for one case.

Experimental design
The regional experiments are for the 21-31 March 2013 ACCACIA period. Simulations of 24 h were run, initialised at 0000 UTC every day from the Met Office's operational global analyses. Hourly model data have been interpolated in space and time to match that of the surface-layer dataset. These observations typically correspond to a lead time of 12-15 h. The global model experiments incorporate 24 case-study hindcasts of 5-day duration, with 12 cases across the extended-winter months (November to April inclusive) and 12 cases across the extended-summer months (May-October) between 2011 and 2014. This experimental design is standard practice for sensitivity testing and trialling parametrization schemes in the MetUM (e.g. Elvidge, 2016).
Three model experiments have been carried out with different specifications of surface roughness over the MIZ and sea-ice. The existing MetUM scheme uses two values of z 0 over ice: one for the MIZ (z 0_MIZ ) set for A = 0.7 and one for solid sea ice (z 0_ice ) set for A = 1 (see Table 1 and Figure 2). The C DN for each grid point is calculated by interpolation between these values and that of open water, that is, a mosaic method -see UM Documentation (2016). Two settings of the existing scheme are used -Smooth and Rough -while the third experiment uses the new sea-ice drag scheme "L12E16". The Smooth z 0 values have often been used in climate model configurations, while the Rough z 0 values have been used in the global operational forecast model for many years and a few climate model configurations (see Table 1). The L12E16 parametrization became part of the global operational forecasting system on 25 September 2018 and will also be used in future configurations of the Met Office's climate models via the GL8 (Global Land 8) configuration. Figure 2 shows the neutral drag coefficient, C DN , as a function of ice fraction for the three roughness length settings. C DN is very large for intermediate ice fractions in the Rough setting, due to a very large value of z 0_MIZ (originally motivated by some indirect Antarctic MIZ estimates of roughness by Andreas et al. (1984)). It has been known for some time that this very large z 0_MIZ value is not representative of many MIZ environments and lower values have been substituted where appropriate (e.g. Outten et al., 2009;. However, it has remained in the operational configuration of the MetUM partly as the alternative (the Smooth setting) appeared too smooth and partly because it was not clear how to replace it. Figure 2 illustrates that L12E16 has a C DN that is up to 35% higher than the Smooth setting, and much lower than the Rough setting. By design, L12E16 fits the observations of E2016 very well: approximately equal to the median values of the 0.4, 0.6 and 0.8 ice fraction bins and well within the interquartile range for all bins. The Smooth C DN is also within the observed interquartile ranges, whilst the Rough C DN is well outside them (the Barents Sea and Fram Strait do not have substantial sea-ice ridges).
It is worth noting that although the ACCACIA aircraft observations are used as validation data for the regional model experiments, the C DN values in the model are derived from the Rough, Smooth and L12E16 parametrizations with ice fraction taken from the model's operational surface boundary conditions, that is, from the Operational Sea-surface Temperature and sea Ice Analysis (OSTIA) dataset. In other words, the ACCACIA period experiments are a test of these three schemes' ability to reproduce the observed atmosphere without in situ observations of the surface.
In both the regional and global model experiments we are testing the model's sensitivity to MIZ drag, both near the surface and in the BL. We provide a quantification of this sensitivity both locally and globally which will likely be similar in other atmospheric models that employ comparable BL parametrizations. Note in comparing the simulations against meteorological observations we are investigating which drag scheme performs the best, with all other factors held the same, that is, all other parametrizations and the initialisation data are identical for each set of simulations.

SENSITIVITY OF THE ATMOSPHERIC SURFACE LAYER TO MIZ DRAG
A comparison of the three regional MetUM simulations (Rough, L12E16 and Smooth) against the ACCACIA surface-layer (SL) observations is summarized in Table 2. Generally, the MetUM simulates the atmospheric SL reasonably well for this collection of cold-air outbreak cases. Overall the L12E16 simulations are the most accurate, most frequently having the highest correlation coefficient, lowest bias and lowest RMS (root-mean-square) error of the three experiments. The correlation coefficients for the L12E16 and Smooth experiments are similar, with both noticeably higher than for the Rough experiment. In addition, the RMS errors for the L12E16 and Smooth experiments are often similar. The clearest statistical difference is in a significantly lower bias for the L12E16 simulations for wind speed, temperature Comparison variables are wind speed (U, units of m s −1 ), temperature (T, K), potential temperature ( , K), specific humidity (q, g kg −1 ), relative humidity w.r.t. water (RH, %), surface temperature (T sfc , K) and pressure (p, hPa); as well as covariance momentum flux ( , N m −2 ), sensible heat flux (SHF, W m −2 ) and latent heat flux (LHF, W m −2 ). The highest correlation coefficient, the lowest bias and the lowest root-mean-square errors are in bold; the bias being defined as model -observations and potential temperature, compared to the other experiments. The biases in wind speed and temperature reduce to only −0.06 m s −1 and 0.02 K, which are remarkably low for this sort of observation versus grid-point comparison in a remote polar region. For example, in Renfrew et al. (2009) wind speed and temperature biases of −0.7 m s −1 and −1.3 K were found for the MetUM compared to a set of aircraft-based SL observations during cold-air outbreak conditions. Although the L12E16 bias in relative humidity (RH) is the lowest, in the Rough simulations a negative RH bias combined with a positive temperature (T) bias result in a very low specific humidity bias through a compensation of errors.
The SL comparisons are illustrated as scatter plots for wind speed and T in Figure 4. For wind speed the L12E16 experiment illustrates the large scatter in simulated winds that is common to all the model experiments (the Rough and Smooth comparisons appear similar -not shown). Even for the best performing drag setting (L12E16) the correlation coefficient is only 0.53 and the RMS errors are >2 m s −1 for wind speeds of 4-10 m s −1 , which compares poorly to similar comparisons for more homogeneous surface conditions, such as over the Irminger or Iceland Seas Harden et al., 2015). The relatively large scatter in wind is also evident from the measured sampling variance in momentum flux which is approximately three times the expected sampling variance (cf. Drennan et al., 2007). We suggest the large scatter in winds and momentum flux is a result of the heterogeneity of the surface roughness elements of the MIZ which adds a "stochastic element" that is simply not accounted for in the model.
The differences in temperature between the three experiments are clear in the scatter plots (Figure 4). The relatively large T bias and scatter of the Rough simulation data are evident. While there are obvious improvements for the L12E16 and Smooth data, the lower bias and RMS errors are clear, and the least-squares regression lines are closer to the 1:1 line. Both the L12E16 and Smooth T comparisons have a low linear regression slope: low temperatures are biased warm, higher temperatures are biased cold, by about 1 K.
Also included in Table 2 is a comparison of model output against eddy covariance turbulent flux estimates based on the aircraft observations (see E2016 for details). Comparing model output directly to such flux estimates is problematic, as covariance flux estimates have a relatively large random error associated with the finite sampling of a turbulent process (Donelan, 1990;). This sampling error confers an uncertainty (of typically 30%) to every flux estimate and leads to relatively low correlation coefficients and relatively high RMS errors, compared to   (Table 2). There is a small warm bias in surface temperature in all the experiments, but this is nowhere near large enough to account for the SHF or LHF biases. Instead we believe there to be an inadequacy with the prescription of the model's roughness lengths for heat and moisture (set here as z 0t = z 0q = 0.2 × z 0_skin , where z 0_skin is the interfacial roughness length). Our comparison suggests z 0t and z 0q may be too high. But the lack of reliable observations in the literature has restricted further investigation here. In future work we plan to analyse new turbulent heat flux datasets, then develop and test an improved parametrization for z 0t and z 0q and evaluate the new surface heat fluxes. Even so, it does not seem plausible that improving the scalar roughness lengths could fully account for the very large heat flux biases seen here. It seems that there are other model inadequacies at play, which we also aim to investigate in future work.
An advantage of validating against a set of aircraft observations is we can evaluate the model's ability to capture spatial distributions. Figure 5 Figure 5a shows 10 m wind speed observations of between 2 and 12 m s −1 from a north to east-northeasterly direction and generally slightly higher wind speeds over the MIZ than over the sea ice to the north (i.e. an ice breeze: Chechin et al., 2013). Figure 5b shows flight-level (∼35 m) temperature generally increasing to the south for each flight and noticeably higher for the most southeasterly flight when winds had a more easterly component. Observations of specific humidity echo the spatial changes seen in the temperature (not shown). The models all generally capture the observed spatial distributions in surface-layer meteorology. But there are some systematic differences; for example, the wind speed is generally too high for the two flights to the northwest and the upwind air temperature is too low for the most northerly flight (not shown). These systematic differences are common to all the model experiments -and are likely due to inadequate data for initialisation or problems simulating the large-scale flow evolution -so are not discussed further here. Instead we focus on how each experiment simulates the spatial changes across the MIZ. Figure 5c-f show the change in temperature (ΔT) from the most northerly point in each flight (i.e. the coldest points for these flights). The observed ΔT is an increase of about 4-6 K across the MIZ with some spatial differences between each flight. All three model experiments capture the spatial distribution of ΔT but overestimate the magnitude. In the Rough simulations ΔT reaches 12 K in the two most northerly flights, an overestimate of more than 6 K for many points. In the L12E16 and Smooth simulations ΔT is closer to that observed, but still an overestimate by typically 2-4 K. The model overestimates of ΔT are primarily due to overestimates of the SHF for these cases (not shown), especially for the two most northerly flights where the wind speeds are systematically too high. These flights are illustrative of the large biases in SHF found over all the flights (Table 2). We have also examined spatial distributions of the wind speed for these three flights, but no clear spatial differences are obvious and there is too much scatter in the observed winds to judge which experiment performs best at simulating the spatial distribution.
The comparison of the three simulations against surface-layer observations over the MIZ is encouraging as the L12E16 parametrization represents a significant improvement in accuracy over the existing parametrization -especially for the Rough setting. The Rough experiment has too much momentum and heat being exchanged with the surface, leading (on average) to significant biases in wind speed and temperature. The sensitivity of the atmospheric SL to the parametrization changes is surprisingly large -there are mean differences of 0.8 m s −1 and 2.5 K.

SENSITIVITY OF THE ATMOSPHERIC BOUNDARY LAYER TO MIZ DRAG
We now investigate whether the simulation differences seen in the atmospheric surface layer are carried into the boundary layer. Fortunately, we have an excellent set of observations of BL evolution across the MIZ from a set of dropsondes during a cold-air outbreak on 23 March 2013 during the ACCACIA field campaign. Figure 6 shows cross-sections of potential temperature ( ), specific humidity (q) and wind speed from the dropsonde observations. It illustrates an archetypal internal boundary layer, which is deepening with distance downstream where it is close to neutral static stability (a cross-section of e illustrating the stability of the BL was, largely, saturated neutral -not shown). This BL structure is typical for off-ice flows (e.g. Brümmer, 1996;1997;Renfrew and Moore, 1999). The BL top (from airborne lidar) rises from around 1,000 to 2,000 m; although note that over the sea ice (north of 76 • N) there is a cold pool of air with a secondary inversion at 500 m situated within the larger-scale BL. The BL and q increase steadily with distance south, while the BL wind speed increases between 74.5 and 72.5 • N perhaps due to an ice-breeze effect or a decrease in the surface roughness over the water.
The model cross-sections (not shown) are qualitatively similar to each other and, in general, to the observations. There are some differences: the simulations do not have a cold pool (or secondary inversion) over the sea ice; the simulated BL depth is lower than observed, rising from around 500 to 1,500 m, seemingly due to a lower BL over the sea ice; and the BL top inversions are too diffuse (a common feature in such simulations, e.g. . The  Figure 6. Panel (a) shows sea-ice fraction from OSTIA and from the MetUM; the remaining panels are mixed-layer averages of (b) potential temperature, (c) specific humidity, and (d) wind speed. The black triangles mark dropsonde release locations (as in Figure 6) [Colour figure can be viewed at wileyonlinelibrary.com]. differences between the simulations are only quantitative, for example the BL is slightly deeper (∼100 m) in the Rough simulation than in the L12E16 and Smooth simulations. Note that above the BL there is almost no difference between the simulations.
The BL is well-mixed, so looking at differences in mixed-layer averaged variables is representative of the BL response. Figure 7 shows observed and simulated mixed-layer average potential temperature ( ml ), specific humidity (q ml ) and wind speed (U ml ), where the mixed-layer averages are between 50 and 400 m. All of the simulations have significant discrepancies in ml at the most northerly dropsonde location. What is interesting is that the downstream change in ml is simulated reasonably well in the L12E16 and Smooth experiments, but less well in the Rough experiment which has too great an increase in ml over the MIZ, and then not enough warming further downstream. This is consistent with the average overestimated SHF in the Rough simulations (Table 2). There are qualitatively similar downstream changes in q ml : the L12E16 simulation matches the observations reasonably well, whereas the Rough simulation has too high a gradient in q ml across the MIZ. All the simulations show an increase in BL wind speed in the centre of the cross-section. The Rough simulation has this maximum around 74.5 • N, while the L12E16 and Smooth simulations have it around 74 • N in better agreement with the observations. The largest difference between the three simulations is to the north. The quantities tend to converge with downstream distance, although the differences persist for hundreds of kilometres.
In short, the L12E16 simulation reproduces the observed spatial gradients in BL temperature, humidity and wind speed more accurately than the Rough simulation for this cross-section. This case is an archetypal internal BL development and there is every reason to suppose these results would be generic for such cold-air outbreak cases. Indeed, examining our other ACCACIA observations does provide some support here (although this corroboration is limited to the surface layer). It should be noted that we have only examined off-ice flows (i.e. cold-air outbreaks) where the BL rapidly transitions from statically stable to neutral conditions across the MIZ. The neutral stability reflects the convectively driven mixing which links the SL and BL and means their responses are similar in the model experiments. During on-ice flows, where the BL is likely to be stably stratified with shear-driven mixing, it is possible the link between the SL and BL will be weaker. At present we do not have access to an appropriate dataset for on-ice flow, so testing this drag parametrization's impact on the BL during such situations is reserved for future work.

GLOBAL IMPACTS
To examine the broader-scale impacts of sea-ice drag on the atmosphere we have run a standard test suite of cases in a global configuration of the MetUM with a grid size of 40 km -so a typical MIZ would be 5-10 grid points wide. Figures 8 and 9 provide results for the Northern Hemisphere and Southern Hemisphere respectively, showing mean sea-ice fraction and then the differences between the L12E16 and Rough output for 10 m wind velocities, near-surface air temperature and mean-sea-level (m.s.l.) pressure; where each panel represents the mean difference over 12 hindcasts, each at a lead time of 48 h. The impact in the Northern Hemisphere differs greatly depending on the time of year, essentially because the Arctic MIZ is poleward of the major land masses, so both extended wintertime (November-April) and summertime (May-October) results are shown. For the Southern Hemisphere, impacts are qualitatively similar between the seasons, essentially because the Antarctic MIZ is equatorward of the land mass so cold-air outbreaks are still common in summer (Papritz et al., 2015). Consequently, only wintertime impacts, which are generally higher amplitude, are shown. The impact of the L12E16 scheme is widespread, affecting much of the Arctic and Antarctic regions. The differences in wind speed and air temperature are concentrated over the MIZs of both hemispheres and reflect reductions in surface momentum flux and upward heat fluxes, generally leading to stronger winds and lower temperatures. The exception to this is during summer in the Northern Hemisphere, where the mean heat flux is predominantly downward (not shown) due to relatively warm air and weaker off-ice flow; the L12E16 scheme reduces this downward flux, resulting in a weak positive impact on air temperature (Figure 8f). The areas that are most sensitive to the new scheme include large parts of the northern subpolar seas -such as the Greenland Sea, Barents Sea, Sea of Okhotsk, Bering Sea, Labrador Sea and Hudson Bay -where sea-ice cover advances and retreats seasonally and ice areas can have large variability (e.g. Cavalieri and Parkinson, 2012). In the Southern Hemisphere the impact is almost circumpolar, reflecting where ice concentrations are often below 100% (Comiso and Steffen, 2001) or are highly variable . There is almost no part of the Southern Ocean that is unaffected. The impact in the Southern Hemisphere summer is spatially similar for all variables (not shown), but it is lower in magnitude; for example, temperature differences of <1 K in summer, compared to more than 2 K in winter.
The impact of the L12E16 scheme is also significant. Compared to the Rough scheme, it reduces the wintertime MIZ surface drag by typically 0.1-0.3 N m −2 (i.e. comparable to the mean -not shown), and resultant 10 m wind speeds and air temperatures by 1-3 m s −1 and 1-3 K. These differences are consistent with the regional experiments (Table 2) and are substantial compared to typical model errors. For example, compared to model biases in wind speed and air temperature of 0.5 m s −1 and 2.1 K for our Rough simulations against our aircraft observations over the MIZ (Table 2); compared to biases of 0.12 m s −1 and 0.43 K for ERA-Interim reanalyses over the central Iceland Sea (Harden et al., 2015); and compared to biases of 0.7 m s −1 and −1.3 K for MetUM analyses off southeast Greenland during cold-air outbreaks . Note the shading in Figures 8 and 9 has been scaled to emphasise significant differences -bold colours illustrate differences that are greater than typical model biases.
The impact on the m.s.l. pressure field is also significant but more variable spatially and seasonally (Figures 8g-h and  9d). The greatest differences are up to 0.7 hPa in magnitude. These are substantial differences compared to typical model biases, for example 0.3 hPa for our Rough simulations over the MIZ (Table 2), 0.1 hPa over the central Iceland Sea (Harden et al., 2015), and −1.0 hPa off southeast Greenland .
In the Northern Hemisphere during winter, m.s.l. pressure generally increases over the MIZ and decreases elsewhere. The greater pressure over the MIZ is primarily due to colder temperatures (Figure 8e), but there are clearly some compensations in the mass field more widely, leading to decreases in pressure. These broader pressure changes are commensurate with the positive impact on the low-level wind field, which results in a strengthening of circulation patterns. The impact is one of widespread polar divergence -or an increase in the mean northerly, off-ice winds in particular over the Greenland Sea and the Bering Strait (Figure 8c) -and consequently a decrease in m.s.l. pressure over the highest latitudes and over some land masses, for example, Alaska. During the Northern Hemisphere summer, anticyclonic circulation is strengthened (Figure 8d), leading to a broad amplification of the North Polar high-pressure centre (Figure 8h). The impacts on Southern Hemisphere m.s.l. pressure are more complicated, but also correspond with the changes in temperature and circulation ( Figure 9).
The impact of the L12E16 scheme on simulated m.s.l. pressure and temperature is generally to reduce known biases in the operational MetUM (i.e. with the Rough settings) when compared to ERA Interim (taken as truth). In particular a long-standing positive MetUM bias in wintertime Northern Hemisphere polar m.s.l. pressure is significantly improved.

CONCLUSIONS AND DISCUSSION
The impact of MIZ surface drag on the atmosphere has been evaluated through three sets of numerical weather prediction model simulations. The "L12E16" simulations make use of a physically based drag scheme (Lüpkes et al., 2012) that has recently been tuned with a relatively large set of observations (Elvidge et al., 2016b); while the other two simulations make use of the drag scheme used in the former operational forecast configuration of the Met Office Unified Model (Rough), and that commonly used in climate model configurations of the MetUM (Smooth). The impact on the atmosphere is significant and surprisingly widespread. The differences in simulated surface fluxes are typically 0.1 N m −2 , 70 and 40 W m −2 for momentum, sensible and latent heat respectively; leading to differences of typically 0.5 m s −1 , 2 K and 0.2 g kg −1 for surface-layer wind speed, temperature and specific humidity respectively. Comparing regional simulations against a collection of aircraft-based observations of the atmospheric surface layer over and downstream of the MIZ, we find that the simulation employing the L12E16 scheme performs the best overall. It has the lowest bias in momentum flux, wind speed and temperature and significantly higher correlation coefficients and lower RMS errors than the Rough scheme. The improvement in the simulated momentum flux is comparable to the mean. Based on the above we would recommend that a physically motivated and validated parametrization scheme for atmospheric drag over marginal sea ice is used in weather and climate prediction models. The L12E16 scheme became operational in the Met Office's global forecasting system on 25 September 2018 1 , whilst a similar scheme was implemented on 12 May 2015 (cycle 41) in the ECMWF Integrated Forecast System. This study has limitations. Firstly, the prescription of roughness lengths for heat and moisture has not been properly addressed. Here their prescription is left unchanged from the Smooth settings. But this decision is currently a pragmatic one: we do not yet have confidence in how to prescribe the scalar roughness lengths. Further research is required to define scalar roughness-length parametrizations, validate these against large and reliable observation datasets and then test the parametrizations in models. Secondly, the validation dataset used here is from one region of the Arctic and for off-ice flows resulting in an atmosphere transitioning from being statically stable to neutral or unstable. A subsequent study would ideally make use of observations from other regions of the Arctic or Antarctic and would include some observations from on-ice flows. Thirdly, the comparison of the three sets of simulations against observations may be affected by factors that are beyond our control, such as changes in the way other parametrization schemes react to the changes in surface forcing or the accuracy of the surface representation (e.g. the sea-ice concentration). We have done our best to control for this by using the same initialisation fields and leaving all other parametrizations the same. The consistency of our findings over the regional and global domains does provide some evidence that our results are robust, that is, that we have managed to isolate the impacts of sea-ice drag on the atmosphere.
Although the results presented here are all for one specific model (the MetUM) we would assert that the sensitivity to surface drag that is found is likely to be qualitatively and quantitatively similar in other numerical weather and climate prediction models. At least for the fundamental atmospheric properties that are exchanged with the surface -namely momentum, heat and moisture -and the surface-layer meteorological variables these directly impact -wind speed, temperature and humidity.
We would contend that the relationship between atmospheric surface drag and ice concentration is now reasonably well constrained. However, this relationship is an indirect one. The surface drag is actually a function of ice morphology (ridge heights, melt pond depths, floe heights, etc.). This is demonstrated by the fact that the drag coefficient for 100% sea ice (C DN_ice , an "end point" in the L12E16 parametrization) varies dramatically with ice morphology; for example, Elvidge et al. (2016b) found median values of C DN_ice = 2.5 × 10 −3 for the Barents Sea (smaller unconsolidated floes), compared to C DN_ice = 1.2 × 10 −3 for Fram Strait (larger smoother floes) -see also Castellani et al. (2014). These are substantial differences and in some regions such differences may occur over relatively small distances, even subgrid-scale. The sea-ice pack can be as heterogeneous as the land surface, but at present is all treated as one surface type. One idea to address this would be to incorporate a stochastic element to surface exchange over sea ice and the MIZ. Stochastic parametrization is now widely used in forecast models and, in particular, in ensemble prediction systems where a greater ensemble spread is often sought (e.g. Palmer et al., 2005). Its efficacy in high-latitude locations has rarely been examined (Jung et al., 2016), one example being the impact on sea-ice prediction of a stochastic element in the parametrization of sea-ice strength (Juricke et al., 2014). A stochastic surface exchange scheme for sea-ice areas may be worth investigating as one way to increase ensemble spread of the high latitude atmosphere.
Surface exchange can be dramatically different due to contrasting sea-ice morphologies and it is not currently clear how to account for this in models. Where models include a sophisticated sea-ice component then characteristics of the simulated sea-ice field could be used to derive C DN_ice directly, and indeed C DN in the MIZ directly too, so dispensing with the need for the sort of MIZ scheme investigated here. Developing such a surface exchange scheme that is reliable and well-constrained will require considerable effort. Where models do not include a sea-ice component then assimilating sea-ice morphology properties directly may be appropriate, again with considerable research required. For the time being we believe the surface exchange scheme for the MIZ investigated here will remain appropriate for many years and particularly for uncoupled models.