A theory for the mean ice thickness and the Transpolar Drift in the Arctic Ocean is developed. Asymptotic expansions of the ice momentum and thickness equations are used to derive analytic expressions for the leading-order ice thickness and velocity fields subject to wind stress forcing and heat loss to the atmosphere. The theory is most appropriate for the eastern and central Arctic, but not for the region of the Beaufort Gyre subject to anticyclonic wind stress curl. The scale analysis reveals two distinct regimes: a thin ice regime in the eastern Arctic and a thick ice regime in the western Arctic. In the eastern Arctic, the ice drift is controlled by a balance between wind and ocean drag, while the ice thickness is controlled by heat loss to the atmosphere. In contrast, in the western Arctic, the ice thickness is determined by a balance between wind and internal ice stress, while the drift is indirectly controlled by heat loss to the atmosphere. The southward flow toward Fram Strait is forced by the across-wind gradient in ice thickness. The basic predictions for ice thickness, heat loss, ice volume, and ice export from the theory compare well with an idealized, coupled ocean–ice numerical model over a wide range of parameter space. The theory indicates that increasing atmospheric temperatures or wind speed result in a decrease in maximum ice thickness and ice volume. Increasing temperatures also result in a decrease in heat loss to the atmosphere and ice export through Fram Strait, while increasing winds drive increased heat loss and ice export.
Arctic sea ice plays an essential role in the global climate system. Over much of the Arctic and over most of the year, the surface of the ocean is covered by ice. This provides a strong filter on how the atmosphere interacts with the ocean. Ice (or snow-covered ice) has a much higher albedo than open water. As a result, much more incoming solar radiation is absorbed in regions of open ocean than in regions covered by ice. Thus, a change in ice cover alters the heat balance for both the ocean and atmosphere. This dependence provides for an ice albedo feedback where regions of open ocean absorb more heat so that in the subsequent winter less ice is formed, leading to more open ocean and a further increase in heat absorption. The ice formed in the Arctic Ocean is exported to lower latitudes, primarily through Fram Strait, into the Nordic seas, and then into the North Atlantic subpolar gyre. The freshwater transported in solid form is of comparable magnitude to that exported in liquid form, and so represents an important component of the freshwater system (Haine et al. 2015). Once at lower latitudes, this freshwater has the potential to inhibit deep convection in the basin interior. Several such disruptions, often called Great Salinity Anomalies, have been observed (Dickson et al. 1988; Belkin et al. 1998).
Regions of ice formation not only transfer water from the liquid to the solid phase, they are also sources of salinity to the ocean through brine rejection. This results in convection and water mass transformation in the ocean. Ice cover also acts as a filter on the momentum flux from the atmosphere to the ocean. In regions of open water, the atmosphere transmits momentum directly to the ocean while in regions of ice cover the momentum is transferred to the ice. Some of this momentum will be transmitted to the ocean and some will be distributed laterally within the ice. Recent climate model studies have also implicated changes in ice cover in the Arctic in changes to storm tracks, temperature, and precipitation patterns at midlatitudes [see the recent review by Vihma (2014)]. Ice cover in the Arctic is clearly an important player in the global climate and ocean circulation.
The mean ice motion in the Arctic Ocean is dominated by two features: the Beaufort Gyre and the Transpolar Drift (TPD; Fig. 1). They are both thought to be primarily wind driven. The Beaufort Gyre lies in the western Arctic in a region of anticyclonic wind stress curl. This stress puts the ice in motion while the ice in turn drives the anticyclonic Beaufort Gyre below. While there is considerable seasonal and low-frequency variability in the wind stress curl pattern and magnitude (Proshutinsky et al. 2002), there is, in the long-term mean, a downward Ekman pumping in the ocean. This is thought to be modulated by ocean currents and balanced by lateral eddy fluxes arising from instability of the Beaufort Gyre (Manucharyan et al. 2016, 2017; Meneghello et al. 2018; Dodridge et al. 2019).
The Transpolar Drift transports ice from the eastern Arctic shelves across the North Pole and to the south through Fram Strait. The volume of ice exported from the Arctic via the TPD is on the order of 2000 km3 yr−1, although there is strong seasonal and interannual variability. This represents a significant fraction of the total freshwater export from the Arctic (Haine et al. 2015). The TPD is generally moving faster than the underlying ocean and so imparts stress to, and acts as a driving agent for, the ocean. While the Transpolar Drift is also thought to be primarily wind driven, and is readily represented in general circulation models, there is no simple theory to explain its location, strength, or dependence on forcing or internal parameters. Ice in free drift tends to be directed about 25° to the right of the wind, but the TPD crosses the surface isobars from high to low pressure (Fig. 1). While it is possible that this is due to ocean advection, the model study of Steele et al. (1997) indicates that the ice motion is driving the ocean, not the other way around. This suggests that the dynamics that determine the large-scale mean ice drift differ significantly from free drift.
Early efforts to understand the drift of ice in the Arctic focused on the momentum balance. Nansen (1902), Sverdrup (1928), and Rossby and Montgomery (1935) considered various combinations of wind stress, ocean stress, Coriolis force, and horizontal viscosity in estimating the relatively short-term drift patterns of ice floes. Their models agreed well with some observations but there was a lack of general agreement between the models and a range of observations. These models lacked a dependence on ice thickness and ice velocity shear. Campbell (1965) considered ice as a Newtonian viscous fluid by adding a Laplacian viscosity for the internal stress term, which brought into consideration lateral shear in the ice velocity, but still lacked a dependence on ice thickness.
The next step in the evolution of ice models combined the momentum budget with more complex representations of the internal stress tensor and a distribution of ice thickness. Coon et al. (1974) and Thorndike et al. (1975) formulated the thickness equation in terms of a thickness distribution function that was subject to deformation and thermodynamic forcing. They represented the ice stress tensor by a nonlinear elastic–plastic parameterization. Because the rheology depended on both ice velocity and the ice thickness, this coupled the momentum equations, thickness equation, and the stress tensor.
Hibler (1979) proposed a viscous–plastic rheology that relates the ice stress to the rate of strain and the internal ice pressure. Since the internal ice pressure depends on the ice thickness, this couples the dynamical balance (momentum and rheology) with the thermodynamic balance (advection and heat exchange with the atmosphere and ocean). The continuous thickness distribution used by Thorndike et al. (1975) was replaced with two levels: thick and thin ice. This model was specifically designed for long integrations spanning seasonal to interannual time scales. The resulting equations are sufficiently complex that direct numerical integrations are generally required. This dynamic/thermodynamic ice model forms the essential basis for ice models used in most current climate studies. A computationally simpler “cavitating fluid” form of this ice rheology that does not resist lateral shear or horizontal divergence was proposed by Flato and Hibler (1992) and will be used in the theory developed below.
The increasing complexity of ice models has provided the ability for models to fairly accurately represent ice thickness, formation, and drift over a wide range of space and time scales. Basin-scale models and model/data synthesis are now capable of reproducing the mean and low-frequency variability of ice in the Arctic Ocean (Schweiger et al. 2011; Mu et al. 2018). In addition, models can be used to identify the leading-order terms in the ice momentum budgets (Steele et al. 1997) and to calculate how the ice responds to changes in atmospheric forcing representative of low-frequency climate variability and climate change (Zhang et al. 2000; Roy et al. 2015).
While physically more complete and computationally more accurate than the early dynamical models, these latest models are sufficiently complex that they can only be solved for most problems using direct numerical integration. This makes it difficult to discern how and why ice drift and thickness depend on external forcing or internal ice parameters. The objective of the present study is to develop a simple analytic model for the lowest order characterization of the ice drift and thickness in the central and eastern Arctic Ocean.
While it is likely that seasonal and low-frequency variability in forcing is important, even for the mean ice thickness and drift, the approach taken here assumes that the basic ice thickness and drift patterns are not due to time dependence in the forcing but instead result from the mean atmospheric forcing. The starting point is the steady, linear ice momentum equations (Hibler 1979), written in tensor form as
The wind stress vector is given by τi, ui is the horizontal ice velocity vector, σij is the internal ice stress tensor, m = ρ0h is the mass per unit area, and ϵij3 = 0 if any two indices are the same, 1 if the indices are in order, and −1 otherwise. The density of sea ice is ρ0, g is the gravitational acceleration, and H is the sea surface height. The ice–ocean drag has been taken to be linear with coefficient Cd.
The mean surface pressure field in Fig. 1 indicates that the mean surface winds in the Arctic, away from the Beaufort high, are nearly unidirectional from Russia toward Greenland/Canada. Based on this, the coordinate system for the theory is defined such that x1 and u1 are parallel to the wind and x2 and u2 are perpendicular to the wind. The wind stress tensor then reduces to τ1 = τ and τ2 = 0.
Following Hibler (1979), the full stress tensor for the viscous plastic model is written as
where ζ and η are the nonlinear bulk and shear viscosities, P is the internal ice pressure, and δ is the Kronecker delta. The term is the strain rate tensor, defined as
The internal ice pressure is expressed in terms of the ice strength P*, the ice thickness h, the concentration , and an empirical constant c* as
Two major simplifications will now be made. Motivated by the findings of Flato and Hibler (1992), the stress tensor will be represented by the cavitating fluid model [retaining only the last term on the right-hand side of (2)] so that there is no resistance to shear or divergence. The second simplification is to assume that each model grid cell has complete ice cover so that in (4). This provides a great simplification for analytical purposes while still retaining the dominant term for the ice stress subject to full ice cover:
The ice thickness is determined by a continuity equation that balances the horizontal flux divergence of the ice with the formation rate controlled by the net heat loss Q, written as
where λ is the heat of fusion for sea ice. For simplicity, it will be assumed that the heat loss at the top of the ice is determined by the thermal conductivity through the ice and heat exchange with the ocean will be neglected. Both of these assumptions will be checked with a more complete numerical model in the following section. The heat loss is then related to the surface temperature of the ice as
where Δ = Tf − Ts is the temperature change across the ice, Tf is the freezing temperature, Ts is the surface temperature of the ice, and κ is the thermal conductivity of ice. The inherent time scale defined by this advective/heat loss balance in (6) is given by h2ρ0λ/κΔ = O(107) s. So the present theory is most applicable for roughly seasonal and longer time scales.
Anticipating solutions provided in the next section, and motivated by the observed ice drift and thickness distribution, it will be assumed that the ice thickness changes much more rapidly in the direction of the wind stress than in the direction orthogonal to the wind stress. In the present coordinates, that is equivalent to assuming that ∂/∂x2 = O(ϵ∂/∂x1), where ϵ ≪ 1 (defined below). The momentum and continuity equations are nondimensionalized using the following scaling:
The horizontal scale of the basin is L, the internal ice strength is P*, and a typical velocity scale is U. The thickness scale was chosen in anticipation of a leading-order balance between internal stress and wind stress.
The stretched meridional length scale Y = y/ϵ has been introduced so that ∂/∂Y = O(1).
Now the dependent variables u, υ, and h are nondimensional and subscripts refer to partial differentiation. There results four nondimensional numbers:
The importance of the Coriolis force compared to the internal stress gradient is given by ϵ, the sea surface height tilt term compared to the internal stress is given by δ; the importance of ocean drag compared to wind is given by c; and γ measures the rate of cooling by heat loss to the atmosphere compared to the advective flux divergence of the ice.
Typical values for the scaling parameters are τx = 0.03 kg m−1 s−2, ρ0 = 103 kg m−3, f = 10−4 s−1, L = 106 m, P* = 3 × 104 kg m−1 s−2, U = 3 × 10−2 m s−1, g = 10 m s−2, Cd = 10−4 m s−1, , Δ = 10°C, λ = 3 × 105 J kg−1, and κ = 2 J m−1 s−1 °C−1. Using these values, the nondimensional numbers scale as
The sea surface tilt terms (those multiplied by δ) are an order of magnitude smaller than the next largest terms and so will be neglected going forward. The next largest terms are the Coriolis terms and the drag terms (those multiplied by ϵ and c), which are of similar order and will be retained in the analysis.
The ice thickness was scaled assuming a leading-order balance between the internal stress and the wind stress. This will be shown to be a good approximation in the interior and western Arctic and so the nondimensional h = O(1) over that part of the domain. However, near the Fram Strait, where North Atlantic water enters the basin, the ice is very thin and so the nondimensional h ≪ 1 near the right-side boundary. So it is anticipated that there will be two regimes, one where h = O(1) in the interior, which will be called the “thick ice” regime, and one where h ≪ 1 near the upwind boundary, which will be called the “thin ice” regime. A similar scaling approach can be applied to the “thin ice” regime if one assumes a leading-order balance is between wind and ocean drag. This results in different nondimensional constants, although the results will be the same as that presented below. For the sake of clarity only one nondimensionalization will be used for both regimes.
a. Thick ice regime
The regime where h = O(1) will be considered first. If we expand the dependent variables u, υ, and h in terms of ϵ, for example, as h = h0 +ϵh1 + O(ϵ2), and substitute into (9), the leading-order terms in the momentum and continuity equations are
For simplicity, and consistent with the mean surface pressure field over the central and eastern Arctic in Fig. 1, the wind stress will be taken to be spatially uniform, τ = τ0. The along-wind momentum equation is then readily solved for the leading-order ice thickness as
where hc is the ice thickness at the upwind edge of the thick ice regime and the transition between the thin ice and thick ice regimes is at x = xc.
The leading-order along-wind velocity can be solved from the continuity equation (12c) to be
An integration constant has been chosen such that u0 and h0 are continuous with the thin ice regime at x = xc.
The leading-order across-wind velocity can be calculated from the ice thickness and along-wind velocity through the across-wind momentum equation (12b).
Equation (13) indicates that the ice thickness in this thick ice regime is, to leading order, controlled by the wind. If we consider a domain with x = 0 on the upwind boundary, a wind blowing toward negative x will produce an ice thickness that increases linearly with distance from the boundary. The along-wind ice velocity is, to leading order, controlled by the thermodynamics, the γ terms in (12c) and (14). The heat loss to the atmosphere is balanced by an increase in ice thickness. Because the ice thickness increases more rapidly in x (the direction of the wind) than in y, the heat loss is balanced by the along-wind advective flux divergence, and so the along-wind velocity is that required to balance the thermodynamic forcing. The across-wind velocity is derived as the residual between the ice thickness gradient, which pushes the ice toward negative y, and the Coriolis term, which, for negative ice motion, opposes the pressure gradient. The ice thickness and drift are shown graphically for both the thick and thin ice regimes in Fig. 2 and are discussed further below.
b. Thin ice regime
If it is assumed that the water flows into the Arctic from the south along the right-side boundary, where x = 0, with no ice cover (as is the case just north of Fram Strait), then close to the boundary h ≪ 1. The nondimensional continuity equation (9c) indicates that if h ≪ 1 then u ≫ 1. This changes the leading-order balance in the along-wind momentum equation, compared to the thick ice regime, so that
The continuity Eq. (12c) is now used to derive the leading-order ice thickness as
A boundary condition of h = hbdy has been applied at x = 0.
The thin ice regime is fundamentally different from the thick ice regime. The along-wind ice drift is now controlled by the wind and the ice thickness is controlled by the thermodynamics, opposite to the controls in the thick ice regime. The thickness increases along the flow direction in response to the accumulated heat loss to the atmosphere. For thin ice, the ice velocity increases, and the ice thickness decreases with increasing wind speed. However, for thick ice, the velocity decreases and the thickness increases with increasing wind speed. It is anticipated that there will be a region of transition between the thin and thick regimes in which the stress gradient and the drag terms are both important, but this regime is difficult to treat analytically because the momentum and continuity equations become coupled at leading order.
If we define the transition from thin ice to thick ice to take place for a given ice thickness hc, then the width of the thin ice regime is calculated from (16) to be
It is useful to calculate the ice drift and thickness for the simple case of a uniform wind in the negative x direction. In this calculation the eastern boundary is not perpendicular to the wind, roughly in accord with the orientation of the boundary and wind in the eastern Arctic (Fig. 1). For simplicity, the ice thickness is taken to be zero at the boundary (hbdy = 0). Since the wind stress is independent of y, y only enters the solution parametrically and h0y in this case is determined by the shape of the eastern boundary. Other configurations could also give rise to gradients in ice thickness normal to the wind, such as making hbdy or the wind stress a function of y, but the basic controlling dynamics and thermodynamics are readily revealed with the uniform wind, curved boundary case.
Equations (13), (14), (12b), (15), and (16) were solved with τ0 = −1, ϵ = 0.1, c = 0.2, hbdy = 0.0, and γ = 0.25. The equations are solved over a domain −1.5 < x < 0, 0 < y < 1. The transition between the thin regime and the thick regimes is defined by hc = 0.25. Different choices for hc will shift the transition between regimes but the qualitative behavior remains unchanged.
The ice thickness and horizontal velocities are shown in Fig. 2. The transition between thin and thick ice regimes is indicated by the red line. In the thin ice regime, the ice thickness grows slowly away from the boundary and the ice drift is primarily in the direction of the wind with magnitude 1/c. Once in the thick ice regime the ice thickness grows linearly with distance with a much stronger gradient than in the thin ice regime. The along-wind velocity weakens as the ice thickens. This is required because there is less heat loss to the atmosphere, and less ice formation, yet the along-wind gradient of ice thickness is the same (because the wind is constant), so the velocity must decrease in order to satisfy the continuity equation. As a result, the Coriolis term in the across-wind momentum Eq. (12b) decreases. Because the internal stress gradient (hy) does not decrease, the across-wind velocity gets stronger as the ice gets thicker. This gives rise to a cyclonic circulation with the ice slowing downwind and turning to the left under the wind stress vectors. Note that, although the wind has no curl, the stress imparted on the ocean by the ice has a cyclonic curl. Near y = 1, hy → 0, because the eastern boundary is perpendicular to the wind stress. In this region the Coriolis force drives the ice toward positive y.
This mean ice thickness and drift are qualitatively similar to the observed ice thickness and drift pattern in the Transpolar Drift, where ice is advected off the Siberian shelf, across the North Pole, and most turns to the south to flow out Fram Strait. A section of the spring mean ice thickness along the red line in Fig. 1 is shown by the red line in Fig. 3. In the eastern Arctic the ice grows slowly toward the pole, but at about 1850 km from the boundary the gradient abruptly increases by more than a factor of 5. The ice thickness predicted from (16) is shown by the black line for 0 < x < 1850 km from the boundary. This calculation took the observed ice thickness as the boundary condition at x = 0, and assumed Δ = 10°C, κ = 2.16 W m−1 °C−1, τ0 = −0.04 N m−2, λ = 3 × 105 J kg−1, and Cd = 0.003 m s−1. Further than 1850 km from the boundary, the black line is extended with (13) using the same τ0 and P* = 2.75 × 104 N m−2. There is quite close agreement between the theory and the observed ice thickness. While there is some uncertainty in the choice of parameters, the robust result is that the observed ice distribution shows an upwind region with thin ice and weak gradients and a downwind region with thick ice and strong gradients, consistent with the thin and thick ice regimes. Steele et al. (1997) found suggestions of a similar bimodal distribution in the dominant balances in the ice momentum equation in a general circulation model of the Arctic such that the balance was primarily either between wind and ocean drag or wind and the internal stress gradient.
It is noteworthy that the southward velocity in the theory is not a result of a no normal flow condition at Greenland since this model has no solid lateral boundaries. Consistent with this, the observed ice drift crosses under the surface isobars in the central Arctic, well away from Greenland. It is, however, implicit in the theory that there is some means to support the ice thickness gradient, such as a solid boundary somewhere at large negative x, although this boundary condition does not appear explicitly in the theory. This is akin to the Sverdrup balance for midlatitude wind-driven gyres where there is a zonal pressure gradient in the interior of the basin. There would also be a boundary layer required to satisfy the no normal flow condition at the downwind coast, which will result in a narrow boundary current, but this is not considered here, nor is it required to drive flow to the left of the wind stress vector in the interior.
3. Numerical model
The theory presented in the previous section provides simple solutions and interpretations for the leading-order ice thickness distribution and drift in two different regimes. However, in order to derive the analytic solutions several important assumptions were required. The internal ice stress tensor was taken to be linearly proportional to the ice thickness and independent of ice motion. The spatial derivatives orthogonal to the wind stress were assumed much smaller than those in the direction of the wind stress. The heat loss was assumed to be inversely proportional to the ice thickness. The theory also does not impose any no-normal flow conditions at the boundaries. Finally, contributions from sea surface height tilt and ocean-induced melt were ignored. A numerical model is now used to test the basic predictions provided by the theory in the context of more complete ice, ocean, and surface forcing physics and thermodynamics.
The model is the MITgcm (Marshall et al. 1997), which solves the primitive equations on a fixed Cartesian grid. The ocean is coupled to a dynamic/thermodynamic sea ice model and forced by momentum and heat fluxes from the atmosphere. The sea ice model uses the elastic–viscous–plastic (EVP) rheology of Hunke and Dukowicz (1997) and the two-layer thermodynamic model of Zhang et al. (1998). The EVP rheology solves the stress tensor relation proposed by Hibler (1979), (2)–(4), subject to some numerical techniques and approximations appropriate for synoptic and longer time scales. The thermodynamic model predicts ice thickness and concentration.
The model is forced by heat and momentum fluxes from the atmosphere and by restoring temperature and salinity within the ocean. The atmospheric fluxes are calculated using the bulk formulas of Large and Pond (1981). The downward longwave and shortwave radiation are specified, and the outgoing longwave radiation is calculated from the surface temperature. The latent and sensible turbulent heat fluxes are calculated from the bulk formulas from the wind speed and an atmospheric temperature. The wind stress imparted on the sea ice is calculated using a quadratic drag law, as is the stress between the ice and ocean. There is a heat flux from the ocean into the ice that is linearly dependent on the ice–ocean temperature difference, and the heat flux to the atmosphere is dependent on the heat conduction through the ice, as in (7). The surface temperature of the ice is determined by the surface heat budget between incoming longwave and shortwave radiation, outgoing longwave radiation, and the sensible and latent heat fluxes.
The model domain is an idealized representation of an Arctic basin, a Nordic seas, and a source region for warm, salty water (Fig. 4). This is similar to the model used to study the maintenance of the halocline and cyclonic circulation of Atlantic Water in the Arctic by Spall (2013). The horizontal grid spacing is 10 km and there are 30 levels in the vertical ranging from 25 m in the upper 500 m to 50 m over the lower 500 m. The model depth is 1000 m with a continental slope 100 km wide encircling the basin. There is a narrow strait connecting the northern basin (the Arctic) to the middle basin (the Nordic seas). The model is configured on an f plane with rotation rate 1.2 × 10−4 s−1.
The model temperature and salinity are restored toward values of 6°C and 35, respectively, with a time scale of 30 days for y < 350 km. This provides a source of warm, salty water to balance the ice formation and heat loss at higher latitudes. The model salinity is also restored toward a value of 31 with a time scale of 30 days over the upper 50 m within 100 km of the perimeter of the Arctic basin. This represents the source of low salinity water that is provided to the Arctic by river runoff and flow through Bering Strait. This is of course a great simplification, but the main purpose here is to maintain a halocline to insulate the ice cover from the warm Atlantic water below.
The model is initialized at rest with a temperature of 6°C and a salinity of 35 everywhere except for a 100-m-thick halocline at the freezing point with a salinity of 31 within the Arctic basin. The model is run for a period of 5 years and the diagnostics below are based on the mean over the final two years. The ice thickness and velocities spin up over the first two years, so these diagnostics are representative of the steady balance. Most calculations use steady forcing, but one calculation reported here has been carried out with a seasonal cycle in the shortwave downward radiation, atmospheric temperature, and wind stress.
The atmospheric temperature used to calculate the sensible and latent heat fluxes is shown in Fig. 5a. The atmospheric temperature is a minimum of −15°C at the center of the Arctic basin, increases toward the perimeter of the basin, and continues to increase toward the southern limit of the domain to a maximum of 6°C. Calculations with a uniform atmospheric temperature over the Arctic produce very similar results, so the main conclusions are not sensitive to the detailed specification of Ta(x, y). The choice of this symmetric pattern shows that any variation in the along-wind direction is the result of ice dynamics, not asymmetries in atmospheric forcing. The net downward radiative flux is a constant 325 W m−2 over the Arctic basin and increases linearly to a maximum of 450 W m−2 at the southern limit of the domain. The 10-m wind speed is constant at 4 m s−1 in the negative x direction. This is not very realistic away from the central Arctic, in particular there is no anticyclonic Beaufort high in the western Arctic. Calculations with a region of anticyclonic winds do produce an anticyclonic Beaufort Gyre in the model, but do not significantly influence the central and eastern Arctic, which is the primary region of interest here.
The purpose of these model calculations is not to reproduce the observed ice thickness and drift in detail, but to represent the primary features (thin ice upwind, thick ice downwind, cyclonic drift from shelf through Fram Strait) while being simple enough to test the response to parameters as predicted by the theory in section 2. This central calculation uses typical parameters for the ice model, as summarized by the top row in Table 1. The model is forced by the 10-m wind, but the theory is forced by the wind stress. The wind speed in Table 1 is the model forcing while the stress in Table 1 is calculated in the model using a quadratic law with a specified drag coefficient. The stress value shown in the table is diagnosed from the model output and used in the theory for comparison with the model results.
While the theory assumes steady forcing, there is in reality a large seasonal cycle in ice cover and drift, resulting primarily from seasonal variability in shortwave radiation, atmospheric temperature, and wind speed. Most of the model calculations discussed below and compared with the theory use steady forcing; however, for one calculation a sinusoidal seasonal cycle was added to the shortwave radiation (±60 W m−2), atmospheric temperature (±35°C), and along-wind wind speed (±0.5 m s−1) at the center of the Arctic basin. The amplitude of the seasonal cycle in shortwave radiation and atmospheric temperature were linearly reduced to zero at the southern limit of the domain. The phase is such that the atmosphere is warmest, downward radiation is largest, and winds are weakest at 25% through the model year.
This is of course very simple but does provide a seasonal cycle in the ice cover, export, and volume. The integrated ice volume and ice export are shown in Fig. 6 for both the steady forcing and the seasonal forcing calculations. The spinup time scale is about 1–2 years, after which the steady forcing case is nearly steady and the seasonal forcing case shows a clear seasonal cycle about the mean. The seasonal cycle in ice volume is about 20% of the mean and the seasonal cycle in ice export is about 50% of the mean. These are somewhat lower than observational estimates, but still provide a modest seasonal variability. The mean of the ice volume for the seasonal forcing is very similar to the volume under steady forcing, so the steady theory is a reasonable approximation for the mean of the time-dependent model. The mean ice export with seasonal forcing is 20% less than that for steady forcing. This results from a correlation between the ice thickness and ice drift in the model Fram Strait. The mean ice thickness distribution and drift across the Arctic basin are very similar between the steady and seasonal forcing (not shown). Recognizing that there are some differences due to seasonal forcing, the steady theory nonetheless represents the primary features of the mean of the seasonally forced calculation, supporting the neglect of the seasonal cycle in deriving the theory for the mean ice properties.
The mean ice thickness and drift calculated over the last two years of the steady-forcing run are shown in Fig. 7. The ice has zero thickness in the southeast part of the domain, just north of the model Fram Strait. The ice forms along the right-side boundary and is advected toward the central Arctic. In the central interior of the basin the ice turns toward the south and flows out of the Arctic basin. The ice is very thin near the upwind boundary, initially grows slowly away from the boundary, and then grows more quickly in the along-wind direction. This is broadly consistent with the observed ice thickness and mean drift in the Arctic Ocean, although the width of the observed thin ice region is larger than that produced by the model.
Before carrying out a comparison between the model and the theory, several of the main assumptions of the theory can be checked against the model. It was assumed that the heat loss to the atmosphere was inversely proportional to the ice thickness and that the ocean to ice heat flux was negligible. The heat loss from the ice to the model atmosphere is plotted in Fig. 8a as a function of ice thickness. The colors represent the area within the Arctic domain that has the given heat flux and ice thickness relation. The heat loss dependence on ice thickness is consistent with (7), indicated by the red line with Δ = 13°C (diagnosed from the model) and κ = 2.16 W m−1 °C−1. While it slightly underpredicts the heat flux for thick ice and overpredicts the heat flux for thin ice, the overall magnitude and parameter dependency is well reproduced by (7). The ocean to ice heat flux is generally less than 2 W m−2, much smaller than the ice to atmosphere heat flux (Fig. 8b), as assumed in the theory.
The linear relationship assumed between the internal ice stress and the ice thickness is supported by the scatterplot of the model ice stress σ11 and the ice thickness scaled with 0.5P* shown in Fig. 9. In the model, σ12 and σ21 are much less than σ11, as assumed in the theory.
The expansion of h, u, and υ in terms of the small parameter ϵ resulted in a great simplification of the governing equations and allowed for analytic solutions for the leading-order ice thickness and drift. The leading-order balance predicted in the along-wind momentum equation was between wind and drag in the thin ice regime and wind and internal ice stress in the thick ice regime. The mean along-wind momentum balance is shown in Fig. 10 as a function of x at y = 2900 km in the Arctic basin. Near the upwind boundary the balance is between wind stress and drag with the ocean, as assumed in the thin ice limit. In the central and western basins, the balance is primarily between wind and internal stress, as assumed in the thick ice limit. There is a transitional region, roughly within 200 < x < 600, where both drag and stress are important. The Coriolis term is much smaller and always negative (resulting from the southward flow). The sea surface tilt term is very small away from the downwind boundary. These diagnostics support the scaling (11) and the balances derived for the thin and thick ice regimes.
The theory (17) predicts the transition from the thin ice regime to the thick ice regime, taken here to be where the ice is 1 m thick, to be at a distance of 450 km from the upwind boundary, while in the model, at the midlatitude of the Arctic basin, the ice reaches 1-m thickness at 460 km from the boundary. The theory predicts a maximum ice thickness of 5.1 m while the model produces a maximum thickness of 6.5 m. These quantities are reasonably well predicted by the theory, but this is only one set of parameters and thus is not a very strict test.
The theory predicts various dependencies on the ice strength P*, the thermal conductivity κ, the wind speed τ0, and the domain radius R. While the discussion will be based on changes to κ in the theory, this is the same as changes to the difference between the surface temperature of the ice, which might be expected to change in a changing climate, and the freezing point (Δ) because κ and Δ always appear together. The model was run for a series of 13 additional calculations in which each of these parameters was varied, as summarized in Table 1. The point here is not that these parameters are poorly known, it is to test the parameter dependency encapsulated in the theory. If the behavior predicted by the theory is reproduced in the numerical model, this gives confidence that the theory is representing the essential dynamics and thermodynamics properly, or at least is consistent with the numerical model. The calculation with the seasonal cycle is also included, indicated by the green stars on Fig. 11 and run 14 in Table 1. There is generally good agreement with the steady theory, so the seasonal case will not be discussed in detail. The ice response predicted by the theory over a wider range of variations in wind stress and surface ice temperature, as might be expected in a warming climate, is discussed in the following section.
The maximum mean ice thickness in the model was diagnosed from the final two years for each of the model runs in Table 1 and plotted against the theoretical estimate from (13) with (17), in dimensional units, in Fig. 11a. There is generally good agreement between the model and the theory. The maximum thickness ranges from about 2 m for very weak winds to over 15 m for the lowest value of the ice strength P*. In general, the ice thickness increases for increasing values of wind speed and domain size, decreases with increasing ice strength, and is only weakly dependent on the thermal conductivity. A similar decreasing mean ice thickness for increasing ice strength was found in a general circulation model by Steele et al. (1997), although the mean ice thickness in that model remained constant for very large values of P*. The ice thickness is too large in the model for very strong winds. This is related to the shape of the model domain. If the western boundary is made straight the maximum thickness for the strongest wind speed is 6.4 m, very close to the theoretical prediction. The increased thickness with the curved western boundary is caused by a blocking of southward flow by the boundary, a geometric complexity not included in the theory. However, the otherwise close agreement in Fig. 11a indicates that the narrow outflow through the model Fram Strait is not a significant source of error over most of the parameter space tested here.
The transition between the thin ice regime and the thick ice regime is defined (somewhat arbitrarily) as the distance from the upwind boundary where the ice thickness at the central latitude of the Arctic domain is 1 m. Other choices have been tried and the parameter dependency is qualitatively similar. The agreement between theory and model is not quite as good as for the maximum ice thickness (Fig. 11b). Stronger winds clearly shift the transition region to the west, consistent with the dependence predicted by (17). The very strong wind case does not agree because of the blocking from the southern boundary in the western basin. This expands the thick ice regime, which is why the model ice is too thick. For the model domain with the straight western boundary the transition occurs at 1230 km from the eastern boundary, in much closer agreement with the theory. The domain size in the model also influences the extent of the thin ice regime but the theory does not predict such a dependence.
Another quantity of interest for climate is the total heat loss to the atmosphere. This also provides a direct test of the assumption that the heat exchange with the atmosphere is given by (7) where h comes from (13) and (16). The total heat loss from the ice to the atmosphere, averaged over the Arctic basin, is plotted in Fig. 11c for each of the calculations. There is increasing heat loss for increasing wind speed, increasing thermal conductivity, and increasing ice strength. The average heat loss decreases with increasing basin size (not the total heat loss) because a smaller fraction of the domain is in the thin ice regime, where the heat loss is largest. The fraction of the total heat loss that occurs in the thin ice regime is shown in Fig. 11d. In general, the relative importance of the thin ice and thick ice regimes depends strongly on the model parameters, especially the wind speed. However, the main point is that, in general, both regimes are important for the total heat loss.
The total volume of ice in the Arctic basin represents an important storage of freshwater in the climate system. The predicted ice distribution (13) and (16), with (17), can be integrated over the model domain to predict the ice volume. Here it is assumed that the ice thickness is a function of x only, but the integral takes into account the circular shape of the domain. The total volume diagnosed from the model generally compares well with the theory (Fig. 11e), with the exception of the very strong wind case where the model ice thickness is too large in the western basin. There is increasing ice volume for increasing wind speed, basin size, and thermal conductivity, and for decreasing ice strength.
One of the most important quantities related to ice in the Arctic is how much gets exported to lower latitudes. The mean ice export can be calculated from the basin-averaged heat loss to the atmosphere. The heat exchange with the ocean, which would melt the ice, is only about 2 W m−2 on average and so is neglected here. The ice export F is related to the heat loss through the heat of fusion, which can be expressed as the area integral of the heat loss (7) as
The ice export varies strongly with model parameters, ranging from a low of about 1011 m3 month−1 to a high of about 4 × 1011 m3 month−1 (Fig. 11f). The export increases for increasing wind speed, increasing thermal conductivity, and increasing basin size. Both the model and theory show little sensitivity to ice strength. The ice export is well predicted by the theory at high wind speeds even though the ice volume and maximum thickness are not well predicted. This is because the heat loss, and thus ice formation, is low in regions of thick ice so the western basin, where the error occurs, is not very important for the ice export. These values for the mean export are within the range of estimates for the observed ice export through Fram Strait, which is approximately 2 × 1011 m3 month−1 (Haine et al. 2015). Clearly the model is very idealized and not intended to reproduce the observations in any detail, but it is reassuring that the model and theory are operating in realistic parameter space.
4. Discussion and summary
A theory has been developed to understand the large-scale mean ice thickness and Transpolar Drift in the Arctic Ocean. The solutions satisfy the leading-order momentum and thermodynamic balances and make use of a cavitating fluid ice rheology where the internal ice stress is proportional to the ice thickness. The leading-order terms in the momentum and continuity equations are used to identify two distinct regimes: thin ice in the eastern (upwind) Arctic and thick ice in the western (downwind) Arctic. The drift in the thin ice regime is determined by a balance between wind stress and ocean drag, while the thickness is then determined by the heat loss to the atmosphere, which is inversely proportional to the ice thickness. By contrast, in the thick ice regime the ice thickness is determined by a balance between wind stress and the ice thickness gradient while the ice drift is controlled by the heat loss to the atmosphere. This results in thin, fast moving ice in the eastern Arctic and thick, slow moving ice in the western Arctic. The thicker ice near Greenland is not due to stronger cooling but instead due to the wind. The export of ice through Fram Strait is a result of a three-way balance between the internal stress pushing out of the basin and opposing forces from the Coriolis and ocean drag terms. The basic predictions from the theory for the ice thickness, heat loss to the atmosphere, ice volume, and ice export compare well with an idealized coupled ocean/ice numerical model for the Arctic Ocean. The theory also qualitatively reproduces the observed ice thickness distribution and Transpolar Drift, at least away from the Beaufort Gyre where the theory does not apply. It is expected that the theory will fail close to Fram Strait, where the ocean velocity and lateral shear in the ice velocity, both neglected in the theory, are large. The Beaufort Gyre in the western Arctic is also not considered in the theory and, due to its cyclic nature, must be governed by different dynamics.
In addition to providing insight into the balances that determine the mean ice thickness and Transpolar Drift, the theory can be used to predict how ice in the Arctic might change with a changing climate. In general, it is expected that the surface temperature of the ice, which is a forcing parameter in the theory, will increase in a warming climate. It has also been predicted that the surface winds will increase with increasing atmospheric temperatures (Seo and Yang 2013; Mioduszewsku et al. 2018). The theory was applied to a square basin 1.8 × 106 m × 1.8 × 106 m with a uniform wind stress of strength −τ0 and surface ice temperature of Ts with P* = 2.75 × 104 N m−2, λ = 3 × 105 J kg−1, κ = 2.16 W m−1 °C−1, Cd = 0.001 m s−1, hbdy = 0. The critical ice thickness at the transition from the thin ice regime to the thick ice regime was chosen as hc = 1 m. Several quantities of climatic interest were calculated as a function of the wind strength and the surface temperature of the ice, since these are the parameters that are expected to change in a warming climate. The maximum ice thickness is shown in Fig. 12a. For weak winds, the maximum thickness is nearly independent of the surface ice temperature. This is because most of the basin is in the thick ice regime in which the ice thickness is controlled by the wind. However, for stronger winds or warmer surface ice temperature the maximum ice thickness begins to decrease. In this case the domain is dominated by the thin ice regime in which the ice thickness is determined by both the wind and the heat loss to the atmosphere. When the wind is strong and the atmosphere is warm, the ice is thin because the thin ice regime extends all the way to the western boundary. In these cases, neglect of a no-normal flow western boundary condition is probably important because the ice velocity is strongly downwind. The maximum ice thickness is found where the surface temperature is given by
A similar pattern is found for the volume of ice in the Arctic basin (Fig. 12b). At low winds the volume is controlled by wind but at higher winds the surface temperature becomes important. The maximum volume occurs at slightly weaker winds (or colder surface temperature) than the maximum thickness, but the overall behavior is similar. This may lead to a positive feedback where stronger surface winds result in thinner ice and, according to some recent studies (Seo and Yang 2013; Mioduszewsku et al. 2018), thinner ice can increase the surface winds due to enhanced heat loss through the ice.
Another quantity of interest for climate is the heat exchange with the atmosphere (Fig. 12c). Contrary to the ice thickness and volume, the heat exchange with the atmosphere increases monotonically with increasing wind stress and decreasing surface temperature. The heat loss depends in a very nonlinear way on the distribution of ice thickness, (7), and so it is not simply related to the mean ice thickness (or volume) in the Arctic basin. Based on the observed ice thickness distribution in Fig. 3, it appears that the real Arctic Ocean lies in the thin ice dominated regime and so, based on Fig. 12, increases in wind stress or surface temperature will result in decreased maximum ice thickness and ice volume. The heat exchange with the atmosphere can either increase, decrease, or remain constant depending on the relative changes in the wind stress and surface temperature. Recall that ice export is directly related to the total heat loss to the atmosphere, so while cases of strong wind may result in a decrease in ice volume, they also result in an increase in ice export.
The theory developed here is for the time-mean ice thickness and drift. Numerical calculations with a simple seasonal cycle indicate that the basic predictions from the theory apply to the mean of a seasonally forced model, although the steady theory slightly overpredicts the ice export. Seasonal forcing causes the upwind edge of the ice to oscillate between the boundary and some distance toward the central Arctic, which effectively changes the along-wind extent of the domain in the theory and may result in significant heat storage in the ocean. Variability on longer time scales is also of interest as large-scale wind patterns are known to vary on decadal time scales and it is expected that the Arctic will continue to warm.
MAS was supported by the National Science Foundation under Grant OPP-1822334. Comments and suggestions from Michael Steele, Gianluca Meneghello, and an anonymous reviewer helped to clarify the work.