Latitude of ice margin as a function of a non-dimensional total solar irradiance in the diffusive energy balance climate model described by North 1975, for different values of the non-dimensional diffusion . Stable states are indicated by a thicker line.
When we were first starting out as graduate students, Max Suarez and I became interested in ice age theories and found it very helpful as a starting point to think about energy balance models for the latitudinal structure of the surface temperature. At about the same time, Jerry North had simplified this kind of model to its bare essence: linear diffusion on the sphere with constant diffusivity, outgoing infrared flux that is a linear function of surface temperature, and absorbed solar flux equal to a specified function of latitude multiplied by a co-albedo that is itself a function of temperature to capture the different planetary albedos for ice-free and ice-covered areas. Playing with this kind of “toy” model is valuable pedagogically – I certainly learned a lot by building and elaborating this kind of model — and can even lead to some nuggets of insight about the climate system.
Using the same notation as in post #36,
.is the absorbed solar flux, where is the co-albedo. We set , where is the global mean incident flux = total solar irradiance/4 so that averages to unity over the sphere. The form , with the second Legendre polynomial is a pretty good fit to the annual mean incident solar flux as a function of latitude . The seasonal cycle is ignored here. is the outgoing infrared flux linearized about the reference temperature , is a heat capacity, and is a kinematic diffusivity. We can just set (and think of it as the freezing temperature when we set the albedo). We can define a non-dimensional temperature, , diffusivity , and mean incident solar flux or total solar irradiance, . Finally, we can use the simplest possible albedo formulation that provides some ice-albedo feedback: for and for . I use for the results described here.
Our equation for steady state solutions independent of longitude, writing out the divergence of the diffusive flux in spherical coordinates, is now
(Aug 12: Fixed a couple of typos in the last few days; hopefully OK now.) If both the diffusivity and the albedo are chosen to be spatially uniform one can solve this equation analytically for this specific choice of because is then a constant plus a term proportional to the second Legendre polynomial which is an eigenfunction of the Laplacian. Comparing the result with with that for non-zero , one finds that the presence of diffusion reduces the equator-to-pole temperature gradient by the factor . One needs this reduction to be about a factor of 2-3 to get a reasonable temperature gradient, which translates into a value of of 0.2-0.3.
We look for solutions that are symmetric about the equator and have temperatures below freezing poleward of a given latitude, the” icecap edge”, and above freezing equatorward of this latitude. The resulting ice edge as a function of for different values of is shown in the figure at the top of the post. Ice-free states are indicated by the horizontal line at and ice-covered “snowball Earth” states by the horizontal line at . Partially glaciated steady states also exist, some of which are unstable. The branches on which the ice edge moves equatorward with increasing are unstable, not surprisingly. There is a small ice cap instability, with ice caps smaller than this threshold receding unstably to the ice-free state. And there is a large ice cap instability beyond which the ice grows unstably due to a runaway albedo feedback, until one reaches the snowball state. This kind of model attracted considerable attention because of the rather small range of solar flux for which partially glaciated states exist and the proximity of these state to the large ice cap instability, for plausible values of the diffusivity. As a point of comparison, Voigt and Marotske 2010, using a modern comprehensive coupled atmosphere-ocean GCM, find that a reduction of 6-9% in the solar flux is sufficient to generate the large-icecap instability. This is obviously an interesting number.
The small ice cap instability in this simple model captures the basic idea that an icecap has to have a certain size to protect itself from diffusion of warm air from lower latitudes. The critical size increases with increasing diffusivity. But this small ice cap instability (unlike the large ice cap instability) turns out to be fragile to modifications to the model such as smoothing the temperature dependence of the albedo or adding a seasonal cycle or adding some noise. To decide if the Arctic ice possesses a small icecap instability requires much more realistic atmospheric and sea ice models. But this simple model does get you thinking about the importance of heat transport from lower latitudes for this issue.
In post #36 there are some plots indicating that the effective diffusivity for heat in the atmosphere should be thought of as having a maximum in midlatitudes. To mimic this schematically, I have set for and elsewhere. The solution is shown by the black line in the figure below. (The results for uniform values of wtih are copied over from the figure at the top for comparison.)
A new instability has been created, with no stable ice caps ending within the region in which the diffusivity has been given the larger value of 0.4. Comparison with the uniform case indicates that it is not simply the magnitude of the diffusivity that creates this instability but rather its horizontal structure. David Linder, Max, and I touched on this kind of behavior in an old paper Held et al 1981, Some related results are discussed in Rose and Marshall, 2009, coming from the direction of trying to include the effects of ocean heat transport in an energy balance model. This kind of stability diagram could generate some interesting hysteresis loops from a time-dependent parameter like the obliquity in Milankovitch ice age theories.
Irrespective of any imagined relevance for the ice age problem, I am interested in closing the gap between this kind of diffusive model and GCMs. In particular, it would be interesting to take an aqua planet atmospheric GCM over a slab ocean with some heat capacity but no horizontal oceanic heat fluxes, no seasonal cycle, no sea ice, but with a simple specified surface albedo as a function of temperature, and map out its behavior as a function of the incident solar flux (with and without cloud feedback). You could then try to mimic this idealized but still chaotic and turbulent GCM’s behavior with steady state energy balance models incorporating theories for the atmospheric heat flux. Points of interest include how an advancing ice edge affects the effective diffusivity and how best to represent cross-equatorial influences.
[I think it is an excellent project for beginning students to generate these solutions on their own. You can reintroduce the time-dependence and integrate forward in time, but this will only give you the stable states. For this simple case of a step function albedo, there is an easy way of getting all of the states shown: specify the ice edge and, therefore, the albedo, then solve the boundary value problem directly for the steady state by inverting the tridiagonal matrix that you get from simple finite-differencing. The solutions will not have the point coincide with the ice edge, so you then need to iterate and find the value of that gives you consistency between temperature and albedo.]
[The views expressed on this blog are in no sense official positions of the Geophysical Fluid Dynamics Laboratory, the National Oceanic and Atmospheric Administration, or the Department of Commerce.]