Nonlinear Processes in Geophysics Albedo parametrization and reversibility of sea ice decay

The Arctic's sea ice cover has been receding rapidly in recent years, and global climate models typically predict a further decline over the next century. It is an open question whether a possible loss of Arctic sea ice is reversible. We study the stability of Arctic model sea ice in a conceptual, two-dimensional energy-based regular net- work model of the ice-ocean layer that considers ARM's longwave radiative budget data and SHEBA albedo mea- surements. Seasonal ice cover, perennial ice and perennial open water are asymptotic states accessible by the model. We show that the shape of albedo parameterization near the melting temperature differentiates between reversible contin- uous sea ice decrease under atmospheric forcing and hystere- sis behavior. Fixed points induced solely by the surface en- ergy budget are essential for understanding the interaction of surface energy with the radiative forcing and the under- lying body of ice/water, particularly close to a bifurcation point. Future studies will explore ice edge stability and re- versibility in this lattice model, generalized to a latitudinal transect with spatiotemporal lateral atmospheric heat trans- fer and high spatial resolution.


Introduction
The future state of the Arctic's sea ice cover is subject to a multi-faceted discussion in the literature.Results of general circulation models (GCMs) are scrutinized regarding the ice extent under various forcing scenarios (Fourth Assessment Report of the Intergovernmental Panel on Climate Change; IPCC, 2007).Single column models (ODE, by e.g., Maykut and Untersteiner, 1971;Thorndike, 1992;Eisenman and Wettlaufer, 2009;Notz, 2009;Müller-Stoffels and Wackerbauer, 2011a) and energy balance models (EBM, by e.g., Sellers, 1969;North, 1984) are used to discuss the stability of Arctic sea ice in terms of few interacting climate components.While GCMs contain the most detail of physical features found in the climate system and describe large-scale interactions between oceans, landmasses and atmosphere reasonably well, their complexity and highdimensional parameter space require hierarchical modeling in order to understand complex climate processes and feedbacks (Bony et al., 2006).
Arctic sea ice has a relatively high albedo in the visible shortwave spectrum, while open water has a low albedo.When sea ice melts, the surface turns from a very good reflector of sunlight to a very good absorber; more open ocean absorbs more solar energy to melt more ice.This ice albedo feedback is one of several known major feedbacks of similar magnitude in the Arctic (Bony et al., 2006) and its understanding is crucial to the understanding of the stability of sea ice.
A variety of albedo parameterizations are used in the literature, ranging from steplike to more gradual albedo changes around the melting point.Albedo parameterizations in GCMs can depend on snow depth, ice thickness, and the frequency band of shortwave radiation, and typically exhibit a sharp step in some neighborhood of the phase transition (e.g., Bony et al., 2006;Eisenman et al., 2007;DeWeaver et al., 2008;Shine and Henderson-Sellers, 1985;Winton, 2008;Liu et al., 2007).Small tuning of the ice albedo value severely changes the outcome of GCM's ice extent predictions (Eisenman et al., 2007).The IPCC (2007) reports large intra-and inter-model differences in ice extent and ice volume for the coming century, but these GCMs do not only differ by the albedo parameterization.In EBM models the form of the albedo parametrization influences the existence of a small ice cap instability (North, 1984).In ODE models the loss of Arctic sea ice is found to be reversible or nonreversible (Thorndike, 1992;Eisenman and Wettlaufer, 2009;Notz, 2009;Müller-Stoffels and Wackerbauer, 2011a).
Published by Copernicus Publications on behalf of the European Geosciences Union & the American Geophysical Union.
ODE-type ice albedo feedback models differ in forcing components and albedo parameterization, and are realized in energy (Müller-Stoffels and Wackerbauer, 2011a), enthalpy (Eisenman and Wettlaufer, 2009), or temperature space (Thorndike, 1992;Notz, 2009).Parameterizations for longwave radiative balances are similar among models insofar as radiative surface heat loss increases with temperature; they are assumed to be seasonally dependent (e.g., Thorndike, 1992;Eisenman and Wettlaufer, 2009) or only dependent on temperature (e.g., Notz, 2009;Müller-Stoffels and Wackerbauer, 2011a).Shortwave forcing was modeled as constant (Notz, 2009), binary on-off (Thorndike, 1992), monthly (Eisenman and Wettlaufer, 2009) and hourly (Müller-Stoffels and Wackerbauer, 2011a) surface heat influx.These differences in shortwave forcing did not result in qualitatively different model dynamics and hysteresis behavior for the model used in this paper.Thus, we hypothesize that the albedo parametrization is key to the model's outcome regarding stability and reversibility of Arctic sea ice, since albedo parameterizations differ across models, including smooth (Müller-Stoffels and Wackerbauer, 2011a) and sharp (Thorndike, 1992) step-like albedo behavior around the melting point and wide linear albedo parameterizations (Notz, 2009).In Sect. 2 we briefly describe the regular network model for the ice albedo feedback (Müller-Stoffels and Wackerbauer, 2011a), adjusted for ARM's longwave radiative budget data (Barrow, AK) and compared with albedo data from the SHEBA experiment.In Sect. 3 we discuss the relevance of fixed points induced by the annual average surface energy budget and the instantaneous surface energy budget for identifying bistability in parameter space as well as for understanding the system dynamics close to bifurcation points.The effect of albedo parameterization on the robustness of hysteresis is explored in Sect. 4.

Model and albedo parameterization
In the following we briefly describe the regular network model by Müller-Stoffels and Wackerbauer (2011a) for an Arctic sea ice-ocean layer under atmospheric and oceanic forcing.Lattice (regular network) models are commonly used to describe nonlinear spatiotemporal dynamics (Cross and Hohenberg, 1993;Stahlke and Wackerbauer, 2009).We review several albedo parameterizations from the literature that are used in conceptual models.The choice of parameterizations for longwave-radiative flux density and for albedo is motivated from observational data.

Model
A thin vertical ice-ocean layer is modeled by a twodimensional regular network of diffusively interacting iceocean cells (Müller-Stoffels and Wackerbauer, 2011a).A conceptual picture of the model domain is given in Fig. 1.
Fig. 1.The thin vertical ice-ocean layer is modeled by a regular network of diffusively interacting ice-ocean cells.The state of each cell is given by its energy content E i,j , which also defines the phase of the cell (ice, water).Each cell can exchange energy with its nearest neighbors via state-dependent nonlinear diffusion.Atmospheric heat fluxes are entering at the upper surface (top row), and oceanic heat fluxes are entering at the bottom row of the ice-ocean layer.The inset shows the steplike phase transition between ice (φ = 0) and water (φ = 1) and the gradual release of latent heat L as the energy of a cell increases (shaded area).
The state of each cell (i,j ) is given by its energy content E i,j , which also governs the liquid-solid phase transition.The model equation for the total energy (sensible and latent heat) E n+1 i,j of cell (i,j ) at discrete time n + 1 is The phase of each cell (i,j ) determines its physical properties such as surface albedo or thermal diffusivity.The phase transition is modeled by a local energy-dependent phase factor (inset Fig. 1), with melting energy E m and volumetric latent heat L. E m = C I h 3 T 0 represents the energy content of a cell of ice at temperature T 0 = 273K with C I the volumetric heat capacity of ice.The sharpness of the phase transition is determined by c.
The phase factor φ controls the phase dependent parameter X via a switch function according to where the subscripts I and W refer to ice and water respectively.The parameter X represents the surface albedo α, the thermal conductivity k, the volumetric heat capacity C, or the release of latent heat L.
The thermal diffusivity D i,j = k i,j /C i,j for a single cell follows Eq. ( 3).We use the molecular thermal diffusivity of ice or water as given in Table 1.Thermal exchange between two adjacent cells is calculated from the respective individual thermal diffusivities D i 1 ,j and D i 2 ,j as for two vertically adjacent cells (i 1 ,j ) and (i 2 ,j ).In this study we assume spatial homogeneity in the horizontal direction, which leads to D i,[j 1 ,j 2 ] = D i,j with j = j 1 = j 2 .
The absorbed solar shortwave radiation R SW depends on the surface albedo α and the actual downwelling shortwave radiation R ↓ SW (e.g., Campbell and Norman, 1998;Müller-Stoffels and Wackerbauer, 2011a), S p0 is the shortwave flux density received at the top of the atmosphere on a surface perpendicular to the incoming flux, and depends on the zenith angle.The factor in square brackets describes the transmission losses in the atmosphere for direct (first term) and diffuse (second term) radiation, with τ m the optical thickness.The absorbed shortwave radiation is both a function of time and local surface properties.
Longwave radiative forcing, R LW : one contribution to R LW is due to the longwave radiative budget between the atmosphere and the surface of the ice-ocean layer according to Stefan-Boltzmann's equation, ( a − s )σ T 4 , with Stefan-Boltzmann constant σ and temperature T in Kelvin.The emissivity s for the upwelling thermal radiation is rather constant, with s = 0.96 being a good approximation for sea ice's and ocean's surface emissivity (Campbell and Norman, 1998).The atmosphere's emissivity a depends on various parameters such as temperature, water vapor pressure, and gas concentration (e.g., Satterlund, 1979;Prata, 1996;Crawford and Duchon, 1999).Not all parameterizations for a given in the literature remain physical (i.e., a < 1) across reasonable model conditions.After exploring several parameterizations (Satterlund, 1979;Prata, 1996;Crawford and Duchon, 1999), we chose the parameterization by Satterlund (1979)  ) .
The vapor pressure e 0 is assumed to be 80 % of the saturation vapor pressure, since an average relative humidity of 80 % in the Arctic is comparable to annually averaged meteorological data for Barrow, AK (ARM, 2010).The saturation vapor pressure is calculated from the temperature (in Kelvin) using the World Meteorological Organization's recommended parameterization by Goff (1957).With these assumptions the exponent in Eq. ( 6) can be well approximated as −e T /2016 0 ≈ −7.7 × 10 −5 T 2 + 0.028T − 3.1.The other contribution to the longwave radiative forcing R LW is due to the thermal contact of the Arctic with lower latitudes via ocean and atmospheric heat transfer.This plays an important role, especially during times when temperature gradients between lower latitudes and the Arctic are large.This can cause strong local temperature inversions, especially in the winter, which complicates calculation of a radiation balance via just one (surface) temperature.To account for the actual transfer of sensible heat northward we assume a linear temperature relationship and fit the longwave radiation budget R LW to longwave radiation data provided by ARM for Barrow, AK (ARM, 2010), following and compares it to daily averages for the 14 yr record of longwave radiation data for Barrow, AK under all sky and clear sky conditions.For comparison we also show the longwave radiation budgets used by Thorndike (1992), Eisenman andWettlaufer (2009), andNotz (2009) in their respective conceptual models.Thorndike (1992) and Eisenman and Wettlaufer (2009) use an irradiance model for the atmosphere to derive the longwave budget.All parameterizations assume that the LW budget is a monotonously decreasing function of temperature.This is not a trivial assumption and does not simply emerge in Stefan-Boltzmann budgets ( a − s )σ T 4 from the parameterizations of a given in the literature (e.g., Satterlund, 1979;Prata, 1996;Crawford and Duchon, 1999), although it is consistent with the data record for Barrow, AK.Both Thorndike (1992) and Eisenman and Wettlaufer (2009) include a direct time-dependence into their longwave budget parameterization, based on observations.Interestingly, our parametrization (Eq.7) and Eisenman and Wettlaufer's August and March parameterization intersect at reasonable seasonal temperatures.
Our system exhibits slightly exaggerated high and low temperature values because several significant dampening factors are not included, as e.g., a connection to a stable temperature reservoir (deep ocean), evaporation effects (cooling due to absorption of latent heat of evaporation), and enhanced heat transfer in the ocean due to eddy thermal diffusivities.

Albedo parameterization
The (surface) albedo describes the fraction of shortwave radiation reflected by the surface, and refers to a spectrally integrated bulk albedo herein.In the sea ice-ocean system the albedo can range from 0.9 for fresh snow to 0.07 for open ocean (Perovich et al., 2002).Values for ice lie within this interval, with young ice's albedo around 0.7, and melting ice and melt ponds down to 0.4 (Perovich et al., 2002).Albedo measurements at a fixed location do not support a clear temperature dependence, especially as snow or ice albedo display aging effects due to surface freeze-thaw cycles and snow metamorphosis.Daily average albedos along a SHEBA transect, however, do reveal a relatively sharp transition between ] Fig. 2. Longwave radiative energy budget R LW versus surface temperature T , with negative values signifying a heat loss of the surface to the atmosphere.The light grey dots represent daily averages for the 14 yr record of ARM's longwave radiative budget data for Barrow, AK (ARM, 2010) for all sky (cloudy and clear) conditions.The black dots mark the subset of data points for clear sky conditions, assuming that the sky is free of clouds if the skyward pointing IR thermometer reads values below 220 K.These scattered clear sky data are filtered with a 50 data point wide moving average (solid blue line).The longwave budget (solid red line) from Eq. ( 7) is fitted to these filtered clear sky data with a higher weight to the regions of higher data density, and yields a high values for ice (Perovich, 1998) and low values for open water and leads (Fig. 3).At temperatures below but close to the freezing point the fluctuations in the surface albedo grow as the average albedo decreases.The average albedo for a lead (Paulson, 1998) is much lower than that for ice at the same temperature (Andreas et al., 1998).
In conceptual models different approaches have been taken in parameterizing the albedo (Fig. 3).Notz (2009) follows Sellers (1969) suggestion of a monotonously decreasing linear function between −43 • C and 10 • C, which results in a very gradual transition of the surface albedo over a wide temperature range (ALW, Fig. 3).Eisenman and Wettlaufer (2009) parameterize the albedo as a smooth step-like function with a relatively steep transition around the freezing temperature (AST, Fig. 3).Thorndike (1992)  The albedo parameterization in Eq. (3) (Müller-Stoffels and Wackerbauer, 2011a) is routinely used in this paper; it follows a smooth step like energy dependence, which appears as a sharp transition in temperature space (ASE, Fig. 3).Two quadratic parameterizations are introduced as intermediate cases between the step-like and the linear transition.One (AQW) is for the wide temperature band (−43 • C to 10 • C) in which Sellers (1969) varies the albedo linearly, and the other (AQN) describes the same transition on a narrower (−18 • C to 0 • C) temperature band (Fig. 3).All these albedo parameterizations were adjusted to fit the SHEBA data for cold ice (α I = 0.85) and open water (α W = 0.07).

System stability, dynamics, and hysteresis
If the ice-ocean system would react to changes in surface forcing instantaneously, the system would precisely follow the instantaneous fixed points resulting from the surface drive.Technically, this would correspond to a (single row) ice ocean layer with atmospheric drive, where the instantaneous fixed points are determined from R LW + R SW = 0. Since the mass of ice and water can only take on or release a certain amount of heat per time step determined by their thermal diffusivities, the system lags behind these instantaneous fixed points.During Arctic night, i.e., the time in winter when the sun remains below the horizon, the surface energy budget takes its minimum value R LW (Fig. 4a), and only one fixed point exists for any b.In summer when shortwave radiation dominates the energy budget, one or three instantaneous fixed points exist depending on b.For a frozen cell the difference between minimum and maximum surface forcing is given by (1 , and for a liquid cell the corresponding difference is The ice albedo feedback is clearly visible as an increase in energy uptake of liquid cells, E − (E m + L/2) > 0.
The Arctic Ocean can be assumed to be in a dynamic equilibrium.Over the years annual forcing, temperature and ice extent cycles are more or less similar, if we disregard the recently observed trends and the paleoclimate record.Thus, information about system stability and system dynamics can be gleaned from annual averages.Figure 4a shows the annual average surface forcing for b = 0 as a function of surface cell energy.If we assume that the ice-ocean layer is driven by this annual average forcing, three fixed points exist: a stable fixed point below the freezing point, an unstable fixed point within the phase transition regime, and a stable fixed point above the freezing point.The bifurcation parameter b manipulates the existence and exact location of the fixed points.The height of the step between the average surface forcing on ice and on water (i.e., around E = E m + L/2) is directly related to the albedo difference between ice and water, and can be approximated as In the following we compare the (dynamic) equilibria observed in simulations of Eq. ( 1) with the fixed points resulting from the average surface forcing.By (dynamic) equilibrium we mean that the system's annual cycle reaches its asymptotic value within numerical accuracy.
From the cubic shape of the annual average surface forcing in Fig. 4a we expect hysteresis of the average cell energy E as a response to atmospheric forcing b.We started simulations in either the open water regime or the ice covered regime and let the system equilibrate under various atmospheric forcings.Figure 4b reveals a single stable fixed point corresponding to perennial ice cover for parameters below b = 28 W m −2 .In a narrow range of forcing between b = 28 and b = 30 W m −2 , bistability between perennial ice and seasonal ice cover exists.For seasonal ice cover all ice in the system has melted (E ij > E m + L/2 ∀ i,j ) for at least one day per annual cycle, and the surface row of cells is frozen (E 1j < E m +L/2 ∀ j ) for at least one day.In both graphs the energy scale is centered around the phase transition (Eq.2).Energies within −0.5L and 0.5L (grey area) correspond to the melting point temperature, i.e., they represent the absorption/release of latent heat of that cell.perennial open water and perennial ice exists, while only the perennial open water state is stable above that forcing regime.
Comparing the hysteresis curve with the fixed points from the annual average surface forcing (black line) in Fig. 4b reveals that the transition from bistability to a single stable perennial ice state at about b = 62 is well predicted by the annual average surface drive, while the width of the forcing interval for bistability is overestimated.
Figure 5 shows the existence and location of the instantaneous fixed points (black graph) for a typical atmospheric forcing ( b = 40 W m −2 ) in the parameter regime of bistability.The time of year during which three instantaneous fixed points exist is surprisingly small and confined to spring and fall.For most of the year only a single fixed point exists with a cell energy below the phase transition for the winter months and above the phase transition for the summer months.
The duration of the existence of three instantaneous fixed points is directly tied to the difference in ice albedo α I and ocean albedo α W .The larger the ice albedo, the longer into the summer the lower fixed point remains in existence (Fig. 4a).Or vice-versa, the larger the ocean albedo the later in the summer the higher fixed point comes to exist.The opposite is the case in the fall.In the extreme case of α I = α W the step in the surface forcing around E = E m +L/2 vanishes, and a single fixed point exists that varies with solar radiation.Variation of α I − α W = 0 gives insight into the dependence of bistability on instantaneous fixed points.For various ice albedos α I ∈ {0.6,0.7,0.85,0.95}and a fixed ocean albedo α W = 0.07 we numerically determine the bifurcation point b ↓ for which simulations starting from open water remain in the high energy state of seasonal ice.The relationship between b ↓ and α I is linear and follows the numerically determined equation b ↓ = 34.67αI − 1.267, which is consistent with our previous results in Müller-Stoffels and Wackerbauer (2011a).For the respective bifurcation points b ↓ we find that the duration of the existence of the upper stable fixed point is at least 150 days.While this value is likely to be strongly model dependent, the finding points towards a critical minimum existence time of a second stable fixed point for bistability to develop.
Figure 5 also shows how typical annual system trajectories are tracking the instantaneous fixed points within the parameter regime of bistability ( b = 40 W m −2 ).The surface cell energy for the equilibrated system, starting from open water initial conditions (red), decreases during winter but does not reach the stable fixed point for ice cover nor is it sufficiently cooled for the cell to freeze before the onset of Arctic day.During that first part of the year the surface cell is heated from below, indicated by a larger average cell energy than surface cell energy, but the surface energy budget is still sufficiently negative to cool the surface cell.Solar shortwave radiation eventually induces a time interval with three instantaneous fixed points; the surface cell trajectory moves away from the unstable steady state around E = E m + L/2 to approach the instantaneous stable fixed point in the open water energy regime.The surface energy increases rapidly due to the low albedo of open water, but does not truly follow the trajectory of the fixed point due to cooling from below the surface and due to the rapid change of the open water fixed point to higher and higher energy values.The surface cell energy reaches its maximum value slightly after the time of maximum solar radiation since the stable fixed point is still located at higher energies than the surface cell energy.Cooling from below together with a rapid decline of solar shortwave radiation reduces the surface cell energy to enter the time period of a stable single instantaneous fixed point, located in the frozen regime.The surface cell energy for the equilibrated system for an initially frozen ice ocean layer (blue) is very close to that of the fixed point in the freezing regime.During the part of summer when that fixed point no longer exists, the surface cell energy stays rather unchanged, since R SW is low due to the high ice albedo and since the open water fixed point is far from the melting point and only weakly attractive (low slope of R LW for high energies).
The surface energy dynamics near the onset and near the ending of the parameter regime of bistability is considered in more detail, focusing on the seasonal progression of the surface cell energy (Fig. 6a-c) and the interaction of surface cell energy with the surface energy budget R LW + R SW (Fig. 6df).Note that the surface energy budget differs from the rate of change of surface cell energy since energy exchange between surface layer and second layer is not included; Fig. 6d-f are therefore not phase space representations in the usual sense.We find that the transition from an open ocean to perennial ice is distinctly different than the transition from perennial ice to open ocean.
Figure 6a shows consecutive annual cycles of surface cell energy for a simulation with open water initial condition and atmospheric forcing parameter b = 27 W m −2 slightly before the onset of bistability (Fig. 4b).During the transition from initially open water to perennial ice the surface cell energy and especially its summer maximum diminishes with each new annual cycle.During the winter months when the surface cells enter the freezing regime, a latent heat signal becomes visible in November/December (day 300 to 365) lasting into the early days of the following year.This release of latent heat manifests itself as a local variation (squiggle) of the surface cell energy, where each consecutive squiggle within an annual cycle denotes a further row of cells freezing.The last annual cycle to leave the freezing regime during summer (seasonal ice as a transient state) exhibits fewer squiggles of latent heat release, since the absorbed solar shortwave radiation only resulted in surface melt.The following year even the surface remains frozen throughout the summer; the lack of a sharp increase in surface energy during summer is direct evidence that the albedo remained high.The annual cycle of instantaneous fixed points in Fig. 6a (black) consists of a stable fixed point in the frozen water energy regime that exists throughout the year.An additional stable fixed point in the open water energy regime together with an unstable fixed point within the phase transition regime exists during most part of the Arctic day.For this atmospheric forcing parameter the upper stable fixed point only occurs in conjunction with the lower stable fixed point.
Figure 6d shows the typical system trajectories from Fig. 6a in the space spanned by the surface cell energy and the surface energy budget.While the system experiences (transient) seasonal ice the trajectory spirals around the stable open water fixed point that results from the average annual surface drive.During the winter months when there is no shortwave radiation, the system trajectory follows the minimum surface budget (R LW ) curve along the horizontal straight line towards the lower stable fixed point.In the phase transition regime (shaded gray) the surface temperature is constant at the freezing point, and the system looses heat at a constant rate.With the onset of shortwave radiation in spring, the trajectory departs from the horizontal line upwards and to the right.The surface energy budget is still negative, but the heat loss is sufficiently reduced from the lower cells warming the surface from below.The trajectory then moves close to the unstable fixed point from the annual average surface drive.Now the direction of the trajectory depends on how much energy is entering the surface.If the incoming energy can melt the surface and thus clearly reduce the albedo, the absorbed shortwave radiation is increasing drastically to warm the surface cells; the trajectory moves quickly to the right.If the surface does not melt however, the albedo stays high and the absorption of solar shortwave radiation is small; the trajectory approaches the R LW curve at lower cell energies during fall, and moves towards the stable frozen water fixed point (resulting from the annual average surface drive) to eventually spiral around that state of perennial ice.A video illustrating the motion of the trajectory together with the instantaneous fixed points arising from the time varying surface energy budget is supplied at Müller-Stoffels and Wackerbauer (2011b).
Figure 6b shows consecutive annual cycles of surface cell energy for a simulation with open water initial condition and forcing parameter b = 28.5 W m −2 slightly above the onset of bistability (Fig. 4b).The system converges fast to a stable seasonal ice cover state denoted by the drop of the surface cell energy below 0L in the winter months and the sudden increase in energy uptake during surface melt at day 130.The corresponding system trajectory in Fig. 6e shows a stable closed path around the open water fixed point (from the annual average surface budget).During spring as the trajectory departs from the R LW curve, the system is heated from below before the surface melts and takes on large amounts of shortwave energy during the summer months, similar to Fig. 6d.During fall the trajectory reaches the R LW curve at high enough surface energies to stay on this side of the instantaneous unstable fixed point (Video, Müller-Stoffels and Wackerbauer, 2011b).The coalescence of the lower stable and the unstable instantaneous fixed point (Fig. 6b) is not an indicator for the onset of bistability as simulations with various ice albedo values show.The duration of the existence of the upper stable fixed point however is slightly longer than 150 days.
Figure 6c shows consecutive annual cycles of surface cell energy during a transition from perennial ice cover to open ocean for a forcing parameter b = 63 W m −2 slightly above the parameter regime of bistability.Because of the much larger forcing b, this transition occurs much more rapidly than the reverse transition in Fig. 6a, which is evident in the wider spacing of the annual cycles, especially during the winter months.Once the ice cover begins to exhibit surface melt it declines rapidly and does not refreeze; a latent heat signal is absent.The annual cycle of instantaneous fixed points still shows a time period when only the stable fixed point in the frozen water energy regime exists, but this period is much shorter than in the two cases described above.As soon as the surface cell energy reaches values above E m + L/2, it increases rapidly and moves fast towards the instantaneous fixed point in the open water regime.Figure 6f reveals that the fixed points resulting from the annual average forcing can describe the system dynamics.The low energy stable fixed point and the unstable fixed point are close together.Heating from below (when R LW + R SW < 0) followed by surface heating from above (when R LW + R SW > 0) moves the trajectory to higher energies past the unstable fixed point.The system heats rapidly due to surface melt and stabilizes on a trajectory outside the freezing regime to circle around the stable open water fixed point.
Figure 6d-f reveals that the fixed points resulting from the annual average forcing do approximate the trajectories in the space spanned by the surface cell energy and the surface energy budget.Close to a bifurcation point, however, the annual average information can fail to predict the transition correctly; the onset of bistability in Fig. is shifted whereas the ending of bistability is well predicted.In such parameter regimes the dynamics is more sensitive and instantaneous fixed points become more relevant to characterize the system dynamics, together with perturbations initiated by heating/cooling from below the surface.for the parameter b 27 W m −2 the average surface budget predicts bistability and consequently an open water asymptotic state, in contrast to the simulation results in Fig. 6a and d.Surface cells cool enough during winters such that the surface cell energy is below the instantaneous unstable fixed point when it comes back to exist in the spring.As the system looses more and more energy in subsequent winters the heating from below in the spring is no longer sufficient to increase the surface cell energy sufficiently far beyond the unstable fixed point to initiate surface melting.The high surface albedo together with the diminishing solar radiation in fall yield sufficient cooling of the surface cells to reduce the surface cell energy quickly to a value below that of the unstable instantaneous fixed point.At this point the system approaches the lower stable fixed point asymptotically.Empirical data support the asymmetry between the predictability of the onset and the ending of bistability (Fig. 4b) by the annual average surface budget.As reported earlier in this section, simulations starting from open water remain in a high energy state of seasonal ice as long as the open water instantaneous fixed point existed for at least 150 days within a year.This criterion is not fulfilled for the parameter range between the predicted onset by the annual average drive ( b ≈ −8 W m −2 ) and the actual onset ( b ↓ ≈ 18.5 W m −2 ).Simulations with initially frozen water remain a low energy state of perennial ice until the end of the predicted parameter regime for bistability ( b ≈ 62 W m −2 ).Throughout this parameter regime the instantaneous fixed point in the frozen water energy range existed for at least 150 days within a year (Fig. 6a-c).

Existence of hysteresis under albedo parameterization
Compared to data taken in a spatially confined area of the Arctic, step-like albedo parameterizations with temperature (or energy/enthalpy) are more suitable to describe albedo changes around the melting point (Fig. 3).This does not mean that the linear albedo parameterization on a wide temperature interval is unreasonable, especially since this parameterization was initially published for modeling climate data in large regions (latitudinal band model, Sellers, 1969).
In such a case it is conceivable that the average albedo across a latitudinal band changes with the average temperature of such a band in a wider temperature interval and with a smoother behavior than a step-like function.Satellite imagery for surface albedo and respective surface temperature would reveal more insight into the dependence of albedo on temperature, but is not part of this paper (Agarwal et al., 2011).We demonstrate in a simple thought experiment that the smoother approach of parameterizing the albedo is valid in ODE models of the Arctic and that the step-like (local) parameterization is linked to the smoother (nonlocal) parameterization.
We assume the Earth to be a perfect sphere with radius R, the Arctic Ocean to be the surface region between two given latitudes θ S and θ N , and the surface temperature T (θ,t) within this region to be a function of inclination angle θ and time t.For simplicity the time variation of surface temperature is independent on the location, and described as ODE models resolve time, and usually include averages of spatially varying physical quantities.A step-like nonlocal albedo parameterization could be obtained from calculating the average temperature of a latitudinal band and using this temperature to determine the corresponding albedo.Averaging the local albedos over the various latitudinal bands between θ S and θ N , however, yields a wider albedo parameterization, Spatially averaged surface albedo α(T ) [θ S ,θ N ] (Eq.9) versus spatially averaged surface temperature with α(T ) [θ S ,θ N ] a function of T 2 (t).For simplicity we consider a linear South to North surface temperature gradient between latitude θ S = 60 • N and θ N = 90 • N in Eq. ( 8), and an albedo α(T ) that can locally be described as a stepfunction with surface temperature T ( • C), with T 0 a scaling temperature used to control the sharpness of the step.The spatially averaged surface albedo is obtained from numerical integration of Eq. ( 9) under the assumptions of Eqs. ( 10) and ( 11).An analytic solution is accessible in this particular case, but not very illustrative.
Figure 7 shows the spatially averaged surface albedo as a function of spatially averaged surface temperature for two linear South to North surface temperature gradients.With increasing surface temperature gradient the spatially averaged surface albedo becomes wider and looses the step-like behavior, even though the local albedos obey a step-like temperature dependence (Eq.11).The actual shape of a spatially averaged albedo depends strongly on the spatial temperature distribution in the Arctic; the linear temperature gradient was assumed to illustrate scaling properties and possible difficulties when the extended system is interpreted from an ODE approach.
In the following we show that albedo parameterization has a significant effect on the properties of stable states and on the existence of bistability and hysteresis.The albedo parameterizations introduced in Fig. 3 cover the range from steplike to wide linear temperature dependence around the melting point with a fixed albedo for cold ice α I = 0.85 and open water α W = 0.07. Figure 8 shows the annual average surface energy balance, R LW +R SW = R LW +[1−α(T )] R ↓ SW , for all albedo parameterizations.The parameter regime of bistability in the annual average surface budget is largest for steplike albedo parameterizations (ASE, Eq. 3; AST).The corresponding cubic shape is progressively less pronounced for more gradual transitions (AQN, AQW, ALW) between ice albedo and open ocean albedo.For the linear albedo parameterization on a wide temperature range (ALW) the cubic feature in the surface energy budget is almost lost.For a slightly smaller difference between α I and α W this cubic feature would disappear, and the annual average surface budget would lack a parameter regime of bistability for that ALW parameterization.
Figure 9 shows the influence of albedo parameterization on hysteresis properties in the average cell energy E as a response to atmospheric forcing b.We find that hysteresis between perennial ice and seasonal ice (or open water) disappears for sufficiently wide/gradual albedo parameterization around the melting temperature.The step-like albedo in temperature space (AST, Fig. 3 parameterization in energy space (ASE, Fig. 4b), the onset of bistability is at a lower b for AST, and the parameter regime of bistability is shortened by about 40 %.Simulations that started with an open ocean initial condition (Fig. 9a, red circles) in the bistable region converge to seasonally ice covered solutions for a wider range of parameters than in the ASE case.The quadratic albedo parameterization on the narrow temperature interval (AQN) softens the step-like albedo change around the melting temperature; the albedo starts to deviate from α I = 0.85 around −15 • C (Fig. 3).This shifts the local minimum in the annual average surface energy budget from the neighborhood of the melting point (AST, ASE) to the neighborhood of −15 • C (AQN, Fig. 8) and reduces the height of the cubic part in the annual average surface energy budget to about two thirds of the height for the step-like functions.Bistability exists on a narrow b-interval of width 6 W m −2 , with the onset of bistability at even smaller atmospheric forcing, b = 6 W m −2 (Fig. 9b).Once the quadratic albedo parameterization is flattened (AQW, Fig. 3) by widening the temperature interval according to Sellers (1969), the local minimum in the annual average surface energy budget moves to lower temperatures, and the cubic structure is further weakened (Fig. 8).Hysteresis between perennial ice and seasonal ice (or open water) is lost.With increasing atmospheric forcing the system transitions rather steadily from perennial ice, to perennial ice with surface melt, to seasonal ice, and then to perennial open ocean (Fig. 9c).This steady transition is interrupted by two very narrow bistable regimes, in both of which the stable states are close together (insert, Fig. 9c).One hysteresis exists between perennial ice and perennial ice with surface melt, the other hysteresis exists between two stable states of perennial ice with surface melt but slightly different energies.For the linear albedo parameterization (ALW, Fig. 3) on the wide temperature interval by Sellers (1969) the cubic structure in the annual average surface energy budget is almost lost (Fig. 8).The transition from perennial ice to open water is very similar to the AQW case, with again no hysteresis between perennial ice and seasonal ice or open water (Fig. 9d).The insert reveals a very narrow region of bistability between two slightly different energies, i.e., perennial ice and perennial ice with surface melt.The transition between ice and ocean albedo can also be made more gradual by reducing the difference between α I and α W and considering a fixed temperature interval for the albedo transition around the melting point.Following Sellers (1969), we keep α I = 0.85 and replace α W with a "nonice" albedo value of 0.37, which is rather high as it incorporates both ocean and land.Under these conditions, the width of the bistability regime is reduced in comparison to the regimes in Fig. 9a and b, and the linear albedo parameterization does not yield any bistability.These results are consistent with the studies by Eisenman (2011) and Abbot et al. (2011) which indicate that bistability disappears or is weakened for smaller differences α I − α W .They are also consistent with Notz (2009), who uses these albedo values to define the slope in the albedo transition and to conclude that the ice pack is stable under perturbations.

Conclusions
The spatially extended conceptual dynamical systems model by Müller-Stoffels and Wackerbauer (2011a) addresses the stability of Arctic sea ice and feedback processes based on the thermodynamic interaction of nonlinearly coupled ice ocean cells with atmospheric and oceanic nonlinear spatiotemporal forcing.We adjust several model components to include observations from the Arctic region.The longwave radiative budget considers clear sky data provided by ARM, and variations in albedo parameterization were inspired by data provided through the SHEBA project and the work by Sellers (1969), Thorndike (1992), Eisenman and Wettlaufer (2009) and Notz (2009).These adjustments follow our core model development paradigm to base the change in physical variables and parameters on the change in system state and not on explicit time-dependences in physical quantities that are based on todays or past observation.This is particularly relevant when investigating system dynamics and feedback processes under a broader range of climatic conditions than accessible by reliable data.
The complexity of the dynamical system is evident in a seemingly minor change in the shape of a model component that significantly changes the outcome of certain model simulations.Adjusting the longwave radiative budget to observations provided by ARM from its site in Barrow, Alaska, results in seasonal ice being a stable state, where previously only perennial ice cover and perennially open ocean were stable solutions (Müller-Stoffels and Wackerbauer, 2011a).This stable annual cycle of seasonal ice is a special case of the higher energy stable state and becomes the open ocean stable state for larger atmospheric forcing.
The annual average surface forcing provides an approximate tool to locate bistability in parameter space.For bistability to exist in our model, it is necessary that the fixed points induced solely by the annual average surface forcing exhibit a parameter regime of bistability.Further analysis of the (instantaneous) fixed points induced solely by the instantaneous surface drive reveals that in much of the parameter regime of bistability only one of the two stable (instantaneous) fixed points exists for most of the annual cycle.The system can be considered to be of cyclically varying stability, possibly displaying time intervals with the existence of a single stable state in the frozen energy regime, a single stable state in the open water energy regime, or with coexistence of both of these stable states.We find that the instantaneous fixed points are essential for the understanding of the system dynamics, particularly in the neighborhood of a bifurcation point.They provide insight into the interaction of the surface energy with the radiative forcing and the underlying body of ice/water.Surface cell energies are tracking the instantaneous fixed points, but are rarely close due to the inertia of the underlying system, which releases/stores energy at a limited rate.
There is an ongoing discussion in the literature about the reversibility of sea ice loss and the existence of a tipping point (Walker, 2006;Lenton et al., 2008;Eisenman and Wettlaufer, 2009;Notz, 2009).The albedo parameterization (among others) is distinctly different between ODE models exhibiting hysteresis and tipping point behavior and models exhibiting a smooth continuous sea ice decrease under increasing atmospheric forcing (Sellers, 1969;Thorndike, 1992;Eisenman and Wettlaufer, 2009;Notz, 2009;Eisenman, 2011).We choose ice and ocean albedos that are commonly used in GCMs (Liu et al., 2007) and manipulate the shape of albedo parametrization around the melting point.We show that albedo parametrization is key to distinct transition behavior.The hysteresis between perennial ice and seasonal ice is less pronounced for albedo parameterizations with a more gradual transition between ice and ocean albedo around the melting temperature.For the linear albedo parameterization on the wide temperature interval by Sellers (1969) the transition from perennial ice to open water is smooth with increasing atmospheric forcing and stable under perturbations.
An interesting question is of course, which albedo parameterization to pick in a conceptual model.Local measurements (SHEBA) reveal that while the albedo starts to diminish slightly below the freezing point, step-like parameterizations capture the transition better than wide gradual transitions.Spatial averaging, however, reconciles the local steplike albedo behavior with the more gradual albedo transition, indicating that gradual parameterizations can provide a bulk albedo for the model Arctic as a whole.Consequently, the albedo parametrization is model dependent.If only local parameters enter into an ODE model then the outcome is local.For a regional ODE model the (spatially averaged) albedo parameterization depends on e.g., temperature gradients in the region, and the actual form of the local albedo parameterization.Further studies are necessary to better represent albedo parameterization in conceptual models and models per se, particularly, since the existence of hysteresis is crucially dependent on the albedo parameterization (this paper), and since small changes to the ice albedo severely change the outcome of GCM's ice extent predictions (Eisenman et al., 2007).
Many GCMs produce a seasonal ice state under increased greenhouse gas emission scenarios (IPCC, 2007).The initial conditions of these GCMs usually allow for some region of open water at lower latitudes.Latitudinal differences in shortwave fluxes generate an inherent spatial inhomogeneity that induce lateral fluxes.This is a significant difference to ODE models and cannot be overlooked when translating ODE model to GCM results.This paper explores a spatially confined model domain at latitude 80 • N. Insolation values for other latitudes would shift parameter ranges for bistability and seasonal ice.In ODE models the common definition of seasonal ice is that the whole model domain will be ice covered for some time every year.The definition of seasonal ice in a GCM is that there exists a region within the model domain that will be ice covered for some time every year.Inherent in the difference is the significant limitation of ODE models regarding spatial resolution.Our findings that (1) smoother albedo parametrizations have significant impact on the asymptotic behavior of ODE models, and that (2) spatial scaling gives a smoother albedo parametrization, are significant in understanding and interpreting ODE model results in a GCM context.
The lattice structure of our model enables spatial inhomogeneities in the ice-ocean layer to be implemented rather straightforward.We adjust the lattice to describe a latitudinal transect between 60 • and 90 • N (with one-tenth of a degree resolution) and include lateral atmospheric heat transfer via diffusion of surface temperatures and enhanced heat transfer in the ocean mixed layer via vertical and horizontal eddy thermal diffusivities.Preliminary results show that bistability in the annual average ice edge location is strongly influenced by the albedo parametrization, and that a locally step-like albedo parameterization results in a wider spatially averaged albedo along the transect.
with b = b 0 + b. a and b 0 are empirical constants, and b is a bifurcation parameter.Figure 2 shows the temperature dependence of this longwave radiation budget R LW for a = −1.05W m −2 K −1 , b 0 = 301 W m −2 , and b = 0 W m −2 , Fig. 2.Longwave radiative energy budget R LW versus surface temperature T , with negative values signifying a heat loss of the surface to the atmosphere.The light grey dots represent daily averages for the 14 yr record of ARM's longwave radiative budget data for Barrow, AK (ARM, 2010) for all sky (cloudy and clear) conditions.The black dots mark the subset of data points for clear sky conditions, assuming that the sky is free of clouds if the skyward pointing IR thermometer reads values below 220 K.These scattered clear sky data are filtered with a 50 data point wide moving average (solid blue line).The longwave budget (solid red line) from Eq. (7) is fitted to these filtered clear sky data with a higher weight to the regions of higher data density, and yields a = −1.05W m −2 K −1 and b = b 0 = 301 W m −2 .Longwave radiative energy budgets used in other conceptual models are plotted for comparison: parameterization by Notz (2009) (dashed grey line; increased by 100 W m −2 to fit the figure); Thorndike (1992) (dash-dotted lines) differentiates between a summer (red) and winter (blue) parameterization; and Eisenman and Wettlaufer (2009) (dashed lines) use monthly varying parameterizations that expand between a maximum (August, red) and minimum curve (March, blue).All parameterizations are monotonously decreasing with T , consistent with the clear sky 50 point running average.

Fig. 3 .
Fig. 3. Observed and parameterized temperature dependence in surface albedo: The data points correspond to daily average albedo along a transect in the SHEBA experiment versus daily average temperature (Perovich's ice SHEBA data (blue w/ error bars); Paulson's lead SHEBA data (red w/ error bars); and Andreas et al.'s temperature SHEBA data).The albedo parameterizations were adjusted to fit the data for cold ice (α I = 0.85) and open water (α W = 0.07).The parameterizations from conceptual models in the literature include a linear temperature dependence on a wide interval (ALW, grey line) suggested by Sellers (1969), a smooth step like temperature dependence (AST, purple line) suggested by Eisenman and Wettlaufer (2009), and a smooth step like energy dependence (Eq.3, ASE, orange line) suggested by Müller-Stoffels and Wackerbauer (2011a).For illustration of intermediate parameterizations a quadratic parameterization over Seller's wide temperature interval (−43 • C ≤ T ≤ 10 • C, AQW, blue dashed line), and a quadratic parameterization over a narrow temperature interval (−18 • C ≤ T ≤ 0 • C, AQN, grey dash-dotted line) are introduced.

Fig. 4 .
Fig. 4. (a) Surface energy budget R SW + R LW versus surface cell energy for b = 0, α I = 0.85, α W = 0.07, and the albedo parameterization in Eq. (3): annual maximum values (upper blue), annual minimum (lower blue), and annual averages (red).The zeros of the surface drives are marked by the horizontal line; the symbols indicate the fixed points resulting from the annual average surface drive.(b) Hysteresis of the average cell energy E as a response to varying atmospheric forcing b.Simulations with an initially ice covered state (blue dots) started with a uniformly frozen ice-ocean layer [E − (E m + L/2) = −0.61Lfor all cells], and simulations with an initially open water state (red circles) started from a uniform water layer [E − (E m + L/2) = 0.85L for all cells].The ice cover is seasonal within a narrow range of forcing b (shaded blue), if started with an open water initial condition.For comparison, the zeros of the annual average surface energy budget R SW + R LW are plotted as a function of b (black graph).In both graphs the energy scale is centered around the phase transition (Eq.2).Energies within −0.5L and 0.5L (grey area) correspond to the melting point temperature, i.e., they represent the absorption/release of latent heat of that cell.

Fig. 5 .
Fig. 5.Annual cycle of the instantaneous fixed points resulting from the surface energy budget (black curve) together with typical system trajectories for b = 40 W m −2 within the parameter regime of bistability (Fig.4b): Trajectory of average cell energy (dark red) and surface cell energy (light red) for the equilibrated system starting from open water initial conditions; trajectory of average cell energy (dark blue) and surface cell energy (light blue) for the equilibrated system for an initially frozen ice ocean layer.The vertical red line denotes the time of year when the system receives shortwave radiation (Arctic day).

Fig. 6 .
Fig. 6.Top row: annual cycle of the instantaneous fixed points (black curve) resulting from the surface energy budget together with typical system trajectories (blue) for various atmospheric forcing b near the onset and near the end of the parameter regime of bistability (Fig. 4b): (a) transition from initially open water to perennial ice ( b = 27 W m −2 ), (b) transition from initially open water to stable seasonal ice ( b = 28.5 W m −2 ), and c) transition from an initially frozen ice ocean layer to perennial open water ( b = 63 W m −2 ).The arrow marks the direction of consecutive annual cycles.The vertical red line denotes the time of year when the surface receives solar shortwave radiation (Arctic day).Bottom row: Corresponding panels that show the trajectories of the surface cells (red dots) in the surface cell energy vs surface energy budget space.Consecutive dots cover the same time interval.The black graphs in each panel show the annual maximum, annual average, and annual minimum surface energy budget (from top to bottom) as a function of surface cell energy.The stable (full circle) and unstable fixed points (open circle) resulting from the annual average drive are marked.A video illustrating the motion of the trajectory together with the instantaneous fixed points arising from the time varying surface energy budget is supplied at Müller-Stoffels and Wackerbauer (2011b).

Fig. 9 .
Fig. 9. Asymptotic states of the average cell energy E as a response to atmospheric forcing b for the albedo parameterizations in Fig. 3: (a) AST, (b) AQN, (c) AQW, and (d) ALW.Simulations with an initially ice covered state (blue dots) started with a uniformly frozen iceocean layer, and simulations with an initially open water state started from a uniform water layer as in Fig.4b.The ice cover is seasonal for the data points within the blue shaded area, and surface melt exists for the data points within the dark grey shaded area.Convergence is reached within numerical accuracy, i.e., doubling the simulation time does not change the values of E .
Müller-Stoffels and Wackerbauer (2011a) (2011a)for model details, and Table1for the physical and numerical parameters used in this study.
main is forced by absorbed solar shortwave radiation R SW , longwave radiative exchange R LW between atmosphere and ice/ocean surface, and an oceanic heat flux density F o that is assumed to be a constant in this study.The factor h 4 ensures correct dimension with h being the cell size.δ ij represents the Kronecker delta.

Table 1 .
Physical constants and system parameters used in all simulations unless stated otherwise.