Journal topic
Clim. Past, 14, 887–900, 2018
https://doi.org/10.5194/cp-14-887-2018
Clim. Past, 14, 887–900, 2018
https://doi.org/10.5194/cp-14-887-2018

Research article 22 Jun 2018

Research article | 22 Jun 2018

# Arctic warming induced by the Laurentide Ice Sheet topography

Arctic warming induced by the Laurentide Ice Sheet topography
Johan Liakka1 and Marcus Lofverstrom2 Johan Liakka and Marcus Lofverstrom
• 1Nansen Environmental and Remote Sensing Center, Bjerknes Centre for Climate Research, Thormøhlensgate 47, Bergen 5006, Norway
• 2National Center for Atmospheric Research, 3090 Center Green Dr., 80301, Boulder, Colorado, USA

Correspondence: Johan Liakka (johan.liakka@nersc.no)

Abstract

It is well known that ice sheet–climate feedbacks are essential for realistically simulating the spatiotemporal evolution of continental ice sheets over glacial–interglacial cycles. However, many of these feedbacks are dependent on the ice sheet thickness, which is poorly constrained by proxy data records. For example, height estimates of the Laurentide Ice Sheet (LIS) topography at the Last Glacial Maximum (LGM;  21 000 years ago) vary by more than 1 km among different ice sheet reconstructions. In order to better constrain the LIS elevation it is therefore important to understand how the mean climate is influenced by elevation discrepancies of this magnitude. Here we use an atmospheric circulation model coupled to a slab-ocean model to analyze the LGM surface temperature response to a broad range of LIS elevations (from 0 to over 4 km). We find that raising the LIS topography induces a widespread surface warming in the Arctic region, amounting to approximately 1.5 C per km of elevation increase, or about 6.5 C for the highest LIS. The warming is attributed to an increased poleward energy flux by atmospheric stationary waves, amplified by surface albedo and water vapor feedbacks, which account for about two-thirds of the total temperature response. These results suggest a strong feedback between continental-scale ice sheets and the Arctic temperatures that may help constrain LIS elevation estimates for the LGM and explain differences in ice distribution between the LGM and earlier glacial periods.

1 Introduction

The Last Glacial Maximum (LGM),  21 000 years before present (21 kyr BP), was the apex of the last glacial period when the global ice volume reached its maximum value . In comparison with earlier glacial cycles, the LGM climate conditions are relatively well documented by proxy data records. In addition to the well-established records of Earth's orbital configuration , atmospheric greenhouse gas concentrations , and sea surface temperatures , the LGM is the only glacial period for which the margins of the North American and Eurasian ice sheets – the Laurentide Ice Sheet (LIS) and the Fennoscandian Ice Sheet (FIS), respectively – can be reliably reconstructed from geological and geomorphological observations. These data records reveal that the LIS was by far the larger of the two, covering most of the North American continent poleward of  40 N, while the FIS was comparatively small and mostly confined to northern Europe . However, while the horizontal margins of the LGM ice sheets have been established, their thickness and vertical extent remain uncertain. This uncertainty is perhaps best reflected in elevation estimates of the LIS, which vary by more than 1 km among contemporary reconstructions (e.g., Peltier2004; Kleman et al.2013; Abe-Ouchi et al.2015). As a result, the Paleoclimate Modeling Intercomparison Project (PMIP4; Kageyama et al.2017) now encourages sensitivity experiments with three distinct ice sheet reconstructions that differ in height by several hundred meters .

Modeling studies have revealed that the continental ice sheets – in particular the LIS due to its larger size and location in the westerly mean flow – had a substantial impact on the atmospheric and oceanic circulation during the last glacial cycle. For example, it has been shown that the LIS topography may have altered both the strength and orientation of the midlatitude Atlantic jet stream and the associated storm tracks and precipitation patterns . Model experiments have also revealed that a higher LIS elevation helps strengthen the Atlantic Meridional Overturning Circulation (AMOC) and wind-driven North Atlantic gyre circulation, which typically results in a poleward shift of the sea ice edge and increased temperatures in the subpolar North Atlantic . Moreover, several studies have shown that the LIS topography can significantly alter the atmospheric stationary wave field (i.e., the zonally asymmetric component of the time-mean atmospheric circulation; Cook and Held1988; Kageyama and Valdes2000; Roe and Lindzen2001; Langen and Vinther2008; Colleoni et al.2016b; Liakka et al.2016; Löfverström et al.2014, 2016). The stationary waves can in turn influence the local temperature and precipitation anomalies that are important for the surface mass balance .

Modeling studies have also shown that the LIS topography can influence the meridional (Equator-to-pole) temperature profile, especially in the northern high latitudes . For example, found that high-end estimates of the LGM LIS elevation increase the mean Arctic surface temperature by several degrees Celsius compared to lower reconstructions. While the authors argued that this response may be attributed to reduced snow cover in Siberia as a result of changes in the atmospheric stationary wave field, they did not explore this narrative further. This result is however noteworthy, as it shows that the elevation of a midlatitude topographic barrier can substantially influence the (zonal) mean climate in high latitudes. This further illuminates a gap in our current understanding of atmosphere–topography interactions and potential feedbacks between continental-scale ice sheets and the temperature in glacial climates.

Understanding the origin, physics/dynamics, and implications of such feedbacks is paramount as they can potentially help constrain the range of possible LIS elevations at the LGM, and also explain differences in the ice sheet extent between the LGM and earlier glacial periods. For example, the penultimate glacial maximum (PGM;  140 kyr BP) had a somewhat reversed ice-volume distribution compared to the LGM, with a larger ice sheet in Eurasia and a comparatively smaller ice sheet in North America . There is also evidence of extensive Arctic ice shelves from the PGM; similar evidence has not been found from the LGM when the LIS was larger .

Here we explore feedbacks between the LIS elevation and the LGM surface temperature in a comprehensive atmospheric general circulation model (AGCM) coupled to a slab-ocean model. In agreement with , we find that raising the LIS topography (from 0 to over 4 km) yields a widespread surface warming in the Arctic, culminating at about 6.5 C for the highest LIS. A thorough analysis of the Arctic energy budget reveals that the LIS-induced surface warming is primarily explained by an increased meridional energy-flux convergence from atmospheric stationary waves, amplified by positive feedbacks from the surface albedo and atmospheric water vapor.

The paper is organized as follows. Section 2 describes the AGCM and the experimental design and Sect. 3 presents the results, which are further discussed in Sect. 4.

2 Model and experiments

We use the National Center for Atmospheric Research (NCAR) Community Atmospheric Model version 3 (CAM3; Collins et al.2006a) with T85 spectral resolution ( 1.4 horizontal resolution) and 26 hybrid pressure–sigma levels in the vertical. Land processes are treated by the Community Land Model 3 (CLM3; Oleson et al.2004). The planetary boundary conditions are prescribed as typical LGM conditions: the orbital parameters are set to appropriate values for 21 kyr BP , and the concentrations of CO2, CH4, and N2O are prescribed as 185 ppmv (parts per million by volume), 350 ppbv (parts per billion by volume), and 200 ppbv, respectively ; CFCs (chlorofluorocarbons) are set to zero. Aerosols, vegetation, and non-glaciated land areas attain their preindustrial (PI) configuration. We use a standard PI simulation as reference climate, which has been evaluated against observations in and , for example – the model captures the amplitude and spatial variations of the observed climate conditions (e.g., surface temperature, precipitation, and geopotential height anomalies) to a high degree.

Figure 1Annual-mean surface temperature (C) (colored shading) and eddy geopotential height at 500 hPa (gray contours) in (a) the PI and (b) the LIStopo0 simulations. Panels (c–g) show the influence of the LIS elevation on the surface temperature with respect to LIStopo0 in the (c) LIStopo0.25, (d) LIStopo0.5, (e) LIStopo0.75, (f) LIStopo1, and (g) LIStopo1.25 simulations, respectively. The contour interval of the eddy geopotential height is 30 m (zero contour is omitted) and negative values are dotted. The black contours outline the LIS and FIS extents in the LGM simulations.

We couple the atmospheric model to a computationally efficient slab-ocean (mixed-layer) model in order to facilitate a high number of experiments. Although the ocean representation is motionless and therefore does not account for changes in circulation, it retains the thermodynamic feedback between the ocean and the atmosphere. Sea surface temperatures (SSTs) and sea ice are explicitly calculated from the energy balance in the ocean mixed layer , where the monthly oceanic heat flux convergence field (commonly referred to as the “q flux”) and annual mixed-layer depth are derived from a 50-year LGM time slice from the fully coupled TraCE-21ka (Transient Climate Evolution over the last 21 000 years; Liu et al.2009; He2011). The PI simulation uses the q flux and mixed-layer depth derived from an atmospheric simulation with PI boundary conditions and prescribed observed sea surface conditions (see also Löfverström et al.2014; Liakka et al.2016). All simulations are integrated for 60 years, of which the first 35 model years are regarded as spin-up, and the remaining 25 years are averaged to create the atmospheric climatological fields used in the analysis.

The ice sheets in North America and Eurasia are derived from the LGM reconstruction in , which is broadly similar to other contemporary reconstructions (see black dashed contours in Fig. 1 with  Peltier2004; Abe-Ouchi et al.2015); the maximum LIS elevation is approximately 3.3 km.

We conduct a total of six steady-state simulations with different heights of the LIS. In each simulation the LIS height (evaluated with respect to the PI topography) is multiplied by a uniform constant N, which takes on values between 0 and 1.25 in increments of 0.25. The LIS morphology and spatial extent therefore remains the same in all experiments, but the elevation is altered: N=0 represents present-day North American orography and LGM land albedo (glacial mask), N=1 is the “standard” LGM case with unscaled LIS topography (3.3 km maximum elevation; Kleman et al.2013), and N=1.25 has a maximum LIS elevation of  4.1 km, approximately similar to the ICE-5G reconstruction (Peltier2004). The LGM sensitivity experiments are referred to as LIStopoN, for which N is the topography scaling factor. As the objective of this study is to evaluate the importance of the LIS topography for the surface temperature, the height of the FIS remains constant in all experiments.

3 Results

## 3.1 Surface temperature and eddy geopotential height

Figure 1a, b show the annual-mean surface temperature from the PI and LIStopo0 simulations. The LIStopo0 surface temperature is substantially colder than PI, particularly in high latitudes. This cooling is explained by the high surface albedo from the ice sheets, a more extensive sea ice cover, lower concentrations of atmospheric greenhouse gases, differences in insolation, and ocean-circulation parameterization (q flux).

The annual-mean surface temperature response to the LIS topography is shown in Fig. 1c–g. This is evaluated as the difference in surface temperature with respect to the LIStopo0 simulation. The largest cooling is confined to the LIS area (as a result of the surface elevation change), and some smaller cold anomalies are found in the mid- and subpolar North Atlantic. The cooling east of the LIS topography is a typical downstream response to topographically induced stationary waves as cold air is advected from the ice sheet interior by the westerly mean flow (e.g., Roe and Lindzen2001; Liakka et al.2011). Elsewhere, raising the LIS topography results in warmer temperatures, particularly in non-glaciated high-latitude land areas and in the Arctic basin (Fig. 1c–g).

Figure 2The Arctic surface temperature anomaly (C) (area-weighted average between 70 and 90 N) with respect to LIStopo0 for the annual mean (ANN), boreal winter (DJF), and boreal summer (JJA) seasons.

The spatial variability in the topographically forced temperature response correlates reasonably well with changes in the stationary wave field (shown here by the 500 hPa eddy geopotential height; gray contours in Fig. 1c–g): the slight cooling in the subpolar North Atlantic is typically associated with a stationary trough, while the warming in Alaska and eastern Siberia coincides with a ridge (Fig. 1c–g).

In contrast to the eddy geopotential height field (which by definition cancels in the zonal mean), the zonal-mean temperature response to the LIS topography is significantly different from zero, particularly at high latitudes. The annual-mean surface temperature in the Arctic increases by approximately 1.5 C per kilometer of LIS elevation (Fig. 2). The high-latitude warming is present in all seasons (Fig. 2), but is strongest in boreal winter (December–January–February; DJF). The seasonal persistence of the Arctic warming suggests that it is driven by changes in both atmospheric dynamics and physics. While the effects of dynamics are typically more pronounced in winter, many radiative features (such as the surface albedo–temperature feedback) are more important in summer when the insolation is higher (note that insolation in boreal winter is negligible at these latitudes).

## 3.2 Role of atmospheric dynamics

Motivated by the apparent connection between the LIS height and Arctic temperatures, the following sections analyze how the atmospheric meridional heat flux changes with the height of the LIS topography.

Figure 3Vertically and zonally integrated (annual-mean) atmospheric northward flux (PW = 1015W) of (a) total energy (F), (b) latent energy (FL), and (c) dry-static energy (FDS) and dry-static energy contributions from the (d) mean circulation (FDSM), (e) stationary eddies (FDSS), and (f) transient eddies (FDST).

### 3.2.1 Basic theory

In steady state, changes in the atmospheric energy storage are zero for annual climatologies , implying that the atmospheric energy balance can be written as (e.g Serreze et al.2007)

$\begin{array}{}\text{(1)}& C=S-R,\end{array}$

where $C\equiv -\mathrm{\nabla }\cdot F$ is the convergence of vertically integrated (annual-mean) horizontal energy flux, R the net radiation at the top of the atmosphere (TOA), and S the energy balance at the surface (R and S are both positive downward). The surface energy balance is defined as the sum of the net radiation and turbulent fluxes. Similarly, the latent energy-flux convergence (${C}_{\mathrm{L}}\equiv -\mathrm{\nabla }\cdot {F}_{\mathrm{L}}$) is proportional to the difference between precipitation (P) and evaporation (E):

$\begin{array}{}\text{(2)}& {C}_{\mathrm{L}}={L}_{v}\left(P-E\right),\end{array}$

with the latent heat of evaporation ${L}_{v}=\mathrm{2.5}×{\mathrm{10}}^{\mathrm{6}}$ J kg−1. The dry-static energy-flux convergence (CDS) is defined as the residual of the total and latent energy fluxes1 (CCL), i.e.,

$\begin{array}{}\text{(3)}& {C}_{\mathrm{DS}}=S-R-{L}_{v}\left(P-E\right).\end{array}$

The implied (zonally and vertically integrated) northward energy flux at each latitude is obtained by integrating Eq. (1):

$\begin{array}{}\text{(4)}& F\left(\mathit{\varphi }\right)=-{a}^{\mathrm{2}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\varphi }}{\int }}C\left({\mathit{\varphi }}^{\prime },{\mathit{\lambda }}^{\prime }\right)\mathrm{cos}\left({\mathit{\varphi }}^{\prime }\right)\mathrm{d}{\mathit{\varphi }}^{\prime }\mathrm{d}{\mathit{\lambda }}^{\prime },\end{array}$

where a is Earth's radius, λ is the longitude, and ϕ is the latitude (both defined in radians). The equivalent northward fluxes of latent (FL) and dry-static (FDS) energy are obtained by substituting C with CL and CDS in Eq. (4).

Table 1The atmospheric flux convergence (W m−2) in the Arctic (area-weighted average > 70 N) separated into total energy (C), latent energy (CL), and dry-static energy (CDS). The latter is further decomposed into contributions from the mean circulation (CDSM), stationary eddies (CDSS), and transient eddies (CDST). δC shows the change in the total energy-flux convergence with respect to LIStopo0.

n/a: not applicable

The energy flux quantities can be further decomposed into the relative contributions from the zonal-mean circulation, stationary eddies, and transient eddies by using atmospheric state variables at model levels; see and Appendix A for details.

### 3.2.2 Meridional flux of atmospheric energy

Figure 3a–c show the implied atmospheric northward energy fluxes (F, FL, and FDS) from our simulations. There is a slight increase in total energy flux (F) in the LGM simulations with respect to PI, and the Northern Hemisphere (NH) peak value is shifted slightly equatorward (Fig. 3a). The LGM change in energy flux is largely represented by a comparable shift in dry-static energy (FDS) (Fig. 3c), while the latent energy flux (FL) shows an overall reduction in the NH midlatitudes (Fig. 3b).

The increase in total energy transport at the LGM (with respect to PI) is in agreement with previous results from both coupled atmosphere–ocean models and from atmosphere models forced by prescribed sea surface conditions . Similarly, the increase (decrease) in dry-static (latent) energy flux is a typical LGM response in fully coupled models and is partially explained by a weaker hydrological cycle in colder climates .

To investigate the role of the LIS topography in more detail, we analyze how different atmospheric circulation regimes influence the meridional energy flux. Figure 3d–f show how the meridional flux of dry-static energy is influenced by the time- and zonal-mean atmospheric circulation (FDSM), stationary eddies (FDSS), and transient eddies (FDST) (see Appendix A for details).

The majority of the increase in low-latitude dry-static energy flux is attributed to changes in the mean circulation (Fig. 3d). These changes are however not directly attributed to the LIS topography, as all LGM simulations show a similar response. The LIS topography is found to be more important for the meridional flux of dry-static energy from stationary (FDSS) and transient eddies (FDST). In LIStopo0, FDSS is roughly similar to the PI, whereas the peak FDST is somewhat higher and shifted equatorward (Fig. 3e, f). Raising the LIS elevation yields a gradual increase (decrease) in stationary (transient) dry-static energy flux in the NH extratropics (Fig. 3e, f), in broad agreement with the fully coupled simulations in and . Here we demonstrate that these changes can be primarily attributed to the LIS topography, as all other boundary conditions remain unchanged in our LGM sensitivity simulations.

### 3.2.3 Energy-flux convergence in the Arctic

While the meridional energy flux (as calculated by Eq. 4) is valuable for identifying structural changes of the large-scale atmospheric circulation, it reveals limited information on how the mean climate responds. For that purpose, the energy-flux convergence (meridional derivative) is a more useful metric than the flux itself (see Eq. 1).

Table 1 shows the horizontal atmospheric energy-flux convergence in the Arctic polar cap (area-weighted average of all grid points poleward of 70 N). The energy-flux convergence in the PI simulation amounts to 102 W m−2, which is in close agreement with estimates from atmospheric reanalysis data (100 to 103 W m−2Serreze et al.2007). The total energy-flux convergence (C) in LIStopo0 is reduced by 6 W m−2 compared to the PI, which is primarily explained by a decrease in the latent energy-flux convergence (CL) (Table 1). Raising the LIS topography yields a gradual increase in the total energy-flux convergence (C) by an average rate of  1.5 W m−2 km−1, resulting in values similar to PI for the highest LIS. This increase stems from an enhanced contribution from the dry-static energy-flux convergence (CDS), in particular from stationary waves (CDSS; Fig. 3 and Table 1); the latent energy-flux convergence (CL) is approximately the same in all LGM simulations. For the highest LIS, the stationary wave contribution to the total Arctic energy-flux convergence even dominates over the contribution from transient eddies (CDST; Table 1). The reduction in transient eddy activity is characteristic for a reduced storminess at the LGM; see , , and for further discussions on this topic.

Figure 4Estimated annual SMB (m yr−1) in the Arctic region from the LGM simulations.

## 3.3 Other feedbacks

The results in Fig. 3 and Table 1 demonstrate that the LIS topography has a dominant influence on the LGM stationary wave field (in agreement with Cook and Held1988; Kageyama and Valdes2000; Löfverström et al.2014). Although the increased energy-flux convergence from the LIS-induced stationary waves is largely compensated for by a comparable reduction from transient eddies, the net effect is a positive contribution (warming) to the Arctic energy balance (δC in Table 1). All else being equal, we can crudely estimate the influence of the LIS-induced atmospheric circulation on the Arctic temperature by assuming a typical value of the surface temperature (Planck) feedback parameter (λT=3.2 W m−2 K−1Flato et al.2013). For the highest LIS, the contribution from the atmospheric circulation is roughly $\mathit{\delta }C/{\mathit{\lambda }}_{\mathrm{T}}\approx \mathrm{2}$C, which is significantly lower than the  6.5 C seen in Fig. 2. This difference suggests that other (positive) feedbacks are important for the LIS-induced warming as well.

Table 2Changes in TOA net radiation (W m−2) with respect to LIStopo0 (δR), and estimated contributions to δR from changes in surface albedo (δRα), water vapor and lapse rate (δRwv+lr), cloudiness (δRcld), and surface temperature (Planck feedback; δRT).

As seen in Eq. (1), changes in the total energy-flux convergence (C) reflect an imbalance between the energy fluxes at the surface (S) and at the TOA (R). For a climate in balance, S is close to zero over land and represented by the horizontal heat flux divergence in the ocean (i.e., $S=-{C}_{\mathrm{ocean}}$Trenberth et al.2001; Serreze et al.2007). As we use a slab-ocean model, the representation of the ocean heat flux is identical in our simulations, implying that the surface energy balance (S) over the ocean remains unchanged when changing the LIS elevation. Any LIS-induced change in the atmospheric energy-flux convergence should therefore be compensated for by an equivalent change in the TOA net radiation balance:

$\begin{array}{}\text{(5)}& \mathit{\delta }C+\mathit{\delta }R=\mathrm{0}.\end{array}$

The influence of the LIS elevation on the Arctic TOA net radiation is shown in the left-most column of Table 2. As suggested in Eq. (5), these values have a magnitude similar to the atmospheric energy-flux convergence (δC in Table 1). The relatively small differences stem from the fact that the (annual-mean) Arctic surface energy balance is not completely identical in all simulations, but varies around 9.1 ± 0.2 W m−2 as a result of slow processes in the land model. These small inconsistencies are however not important for the interpretation of our results and conclusions.

To evaluate the contributions from individual feedbacks to the LIS-induced Arctic warming, the TOA net radiation change (δR) is separated into radiative contributions from changes in surface albedo (δRα), water vapor and lapse rate (δRwv+lr), total cloudiness (δRcld), and surface temperature (i.e., the Planck feedback: δRT). The estimated strengths of these feedbacks are obtained by separately calculating the contributions to the shortwave (SW) and longwave (LW) parts of the radiative spectrum. The SW decomposition is carried out using the approximative partial radiative perturbation (APRP) method, which uses a simplified model to separate changes in net incoming SW radiation into contributions from surface albedo, clouds, non-cloud (clear-sky) SW, and a residual term ; the LW decomposition is provided in Appendix B. While the surface albedo (δRα) and Planck (δRT) feedbacks, respectively, influence the SW and LW radiation separately, the other terms (δRwv+lr and δRcld) are assumed to contribute to both (see Appendix B). Note that the residual term from the APRP method is omitted here as it is not relevant for the discussion. This term is also comparably small so that the total TOA radiation change is approximately equal to the sum of all individual feedbacks, i.e.,

$\begin{array}{}\text{(6)}& \mathit{\delta }C+\mathit{\delta }{R}_{\mathit{\alpha }}+\mathit{\delta }{R}_{\mathrm{wv}+\mathrm{lr}}+\mathit{\delta }{R}_{\mathrm{cld}}+\mathit{\delta }{R}_{\mathrm{T}}\approx \mathrm{0}.\end{array}$

The radiative contribution from each feedback is shown in Table 2. It is evident that the surface albedo (δRα) and water vapor (δRwv+lr) are the most important feedbacks (contributing to increase the Arctic temperature) when raising the LIS elevation; cloud feedbacks (δRcld) are virtually negligible. For the highest LIS elevations, the combined influence of changes in surface albedo and water vapor yield a positive radiative contribution of about 14 W m−2, thus exceeding the contribution from the energy-flux convergence (δC in Table 1) by approximately a factor of 2. Of these two feedbacks, the water vapor feedback is overall about twice as large as the surface albedo feedback, mediated by a  30 % increase in the total precipitable water content between LIStopo0 and LIStopo1.25 (not shown). To compensate for the warming contributions from the atmospheric energy-flux convergence, and the water vapor and albedo feedbacks, there is an increased outgoing LW radiation from the surface as a result of higher temperatures (Planck feedback; δRT), amounting to about 21 W m−2 for the highest LIS reconstruction (Table 2).

4 Discussion and concluding remarks

Here we investigate how the LIS topography influences the Arctic surface temperature, using a comprehensive AGCM coupled to a slab-ocean model. Our results show that increasing the LIS elevation (from 0 to over 4 km), while keeping all other boundary conditions fixed at their LGM configuration, results in an annual-mean Arctic warming in excess of 6.5 C. This warming is primarily attributed to a net increase in the atmospheric energy-flux convergence in high latitudes, which is further reinforced by positive feedbacks from a reduced surface albedo and a higher atmospheric water vapor content.

The correlation between Arctic temperatures and the LIS elevation suggests that LGM LIS may have helped reduce the Equator-to-pole temperature gradient. This is also supported by annual-mean surface mass balance (SMB) estimates (Fig. 4), evaluated as the difference between accumulation (precipitation) and ablation using the positive degree-day (PDD) approach from the ice sheet model SICOPOLIS . Figure 4 shows that most (non-glaciated) Arctic land areas change from a positive to negative SMB when raising the LIS elevation, suggesting that the presence of the LIS topography may have helped keep Alaska and Siberia largely ice free at the LGM (in agreement with Roe and Lindzen2001; Liakka et al.2016; Löfverström and Liakka2016). Furthermore, areas with positive SMBs are found in Siberia in all simulations except LIStopo1.25, which suggests that the maximum LIS elevation at the LGM may have been higher than our default LIS reconstruction (LIStopo1; 3.3 km). It is important to stress that this result likely is model dependent (e.g., the global-mean surface temperature at the LGM is found to be between 3.1 and 5.8 C cooler than PI in the PMIP2 models; Braconnot et al.2007). In addition, the ablation calculations presented here should only be considered a crude first-order estimate, as the PDD model relies on the assumption that the annual melt potential is a function of the monthly-mean surface temperature. It is therefore important to assess the impact of LIS elevation on the SMB using more realistic ablation representations in other models (e.g., SMB parameterizations that are based on the surface energy balance) before we can use this information to constrain the range of possible LIS elevations. With the result presented here, however, we hope to encourage such experiments in the future.

The feedback between continental-scale ice sheets and meridional temperature distribution presented here may also provide a better understanding of glacial environments beyond the LGM. For example, showed evidence of a thick and partially grounded Arctic ice shelf during the PGM, when the LIS is believed to have been considerably smaller than its LGM size . Here we obtain positive SMB values across most of the Arctic coastal areas for the lowest LIS topographies (N≤0.5), while higher LIS elevations yield a negative SMB in these areas (Fig. 4). Hence, these results suggest that the smaller LIS size at the PGM may have been a contributing factor to the formation of an extensive Arctic ice shelf.

The main limitation of this study is that we use a slab-ocean model and thus neglect potential changes in the ocean circulation. However, these changes are not expected to influence the first-order conclusions from this study. There are two main reasons for this. First (i), the primary source of the Arctic warming found here, i.e., the increased energy flux by stationary waves, is also featured in many LGM experiments with fully coupled models (e.g., Li and Battisti2008; Murakami et al.2008). Second (ii), the direct impact of the ocean circulation on the Arctic temperature is only important in regions that were mostly free of sea ice at the LGM. Proxy data from the LGM show that essentially only the subpolar North Atlantic sector, a region strongly influenced by AMOC variability, was characterized by seasonally ice-free conditions poleward of 70 N . However, there is strong evidence from several fully coupled models that raising the LIS elevation yields a stronger AMOC , which is expected to amplify the Arctic warming signal.

Put in perspective, the LIS-induced Arctic energy-flux convergence found here ( 6.5 W m−2) even exceeds the radiative forcing from a doubling of atmospheric CO2 ( 4 W m−2; e.g., Hansen et al.1997), emphasizing the importance of the LIS topography for the LGM climate. A similar influence on the stationary waves has not been found for the FIS or the smaller (pre- and post-LGM) configurations of LIS . Hence, it is possible that the stationary-wave-induced energy flux, and thus also the associated temperature feedback, is only important when the continental ice sheets are sufficiently large to interact with the westerly mean flow (some evidence of this is shown in Löfverström et al.2014; Löfverström and Lora2017). To explore the limits of this feedback with respect to different ice sheet configurations and atmospheric mean states is beyond the scope of this study, but we hope that results from the PMIP4 experiments – in particular the sensitivity experiments with different ice sheet reconstructions – will help illuminate some of these issues.

Data availability
Data availability.

The model output files can be obtained from the first author (johan.liakka@nersc.no) upon request.

Appendix A: Contributions from the circulation to the meridional energy flux

Here we estimate the (zonally and vertically integrated) northward energy flux from the time- and zonal-mean circulation, as well as stationary and transient eddies. We refer to for a more comprehensive review of this topic. The total meridional energy flux (Eq. 4) can be expressed in atmospheric state variables as

$\begin{array}{}\text{(A1)}& F\left(\mathit{\varphi }\right)=\mathrm{2}\mathit{\pi }a\mathrm{cos}\left(\mathit{\varphi }\right)\underset{\mathrm{0}}{\overset{{p}_{\mathrm{s}}}{\int }}\left(\left[\stackrel{\mathrm{‾}}{v{h}_{\mathrm{DS}}}\right]+\left[\stackrel{\mathrm{‾}}{v{h}_{L}}\right]\right)\frac{\mathrm{d}p}{g},\end{array}$

where p is the pressure, ps the surface pressure, v the meridional wind, and g=9.8 m s−2 the gravitational acceleration. Here, ${h}_{\mathrm{DS}}\equiv {c}_{\mathrm{p}}T+gz$ represents the dry-static energy (per unit mass), defined as the sum of the internal and potential energy, where T is the temperature, cp=1004 J kg−1 K−1 the specific heat capacity, and z the geopotential height. The latent energy is given by hLLvq, where q is the specific humidity, and ${L}_{v}=\mathrm{2.5}×{\mathrm{10}}^{\mathrm{6}}$ J kg−1 is the latent heat of evaporation. Over bars denote time mean and square brackets zonal mean. The dry-static energy flux can be separated into the relative contributions from different components of the atmospheric circulation as

$\begin{array}{}\text{(A2)}& \left[\stackrel{\mathrm{‾}}{v{h}_{\mathrm{DS}}}\right]=\left[\stackrel{\mathrm{‾}}{v}\right]\left[\stackrel{\mathrm{‾}}{{h}_{\mathrm{DS}}}\right]+\left[\stackrel{\mathrm{‾}}{{v}^{\ast }}\stackrel{\mathrm{‾}}{{h}_{\mathrm{DS}}^{\ast }}\right]+\left[\stackrel{\mathrm{‾}}{{v}^{\prime }{h}_{\mathrm{DS}}^{\prime }}\right],\end{array}$

where the terms on the right-hand side represent the zonal-mean circulation, stationary eddies, and transient eddies, respectively; asterisks and primes represent deviations from the zonal- and time-mean states. The equivalent contributions from each circulation regime to the dry-static energy flux is given by

$\begin{array}{}\text{(A3)}& & {F}_{\mathrm{DSM}}=\mathrm{2}\mathit{\pi }a\mathrm{cos}\left(\mathit{\varphi }\right)\underset{\mathrm{0}}{\overset{{p}_{\mathrm{s}}}{\int }}\left[\stackrel{\mathrm{‾}}{v}\right]\left[\stackrel{\mathrm{‾}}{{h}_{\mathrm{DS}}}\right]\frac{\mathrm{d}p}{g},\text{(A4)}& & {F}_{\mathrm{DSS}}=\mathrm{2}\mathit{\pi }a\mathrm{cos}\left(\mathit{\varphi }\right)\underset{\mathrm{0}}{\overset{{p}_{\mathrm{s}}}{\int }}\left[\stackrel{\mathrm{‾}}{{v}^{\ast }}\stackrel{\mathrm{‾}}{{h}_{\mathrm{DS}}^{\ast }}\right]\frac{\mathrm{d}p}{g},\text{(A5)}& & {F}_{\mathrm{DST}}={F}_{\mathrm{DS}}-{F}_{DSM}-{F}_{\mathrm{DSS}}.\end{array}$

Hence, the contribution from transient eddies is determined as the residual of the total, zonal mean, and stationary components of the circulation. The equivalent latent energy flux can be obtained from Eqs. (A2) to (A3) by substituting hDS with hL.

The corresponding horizontal energy-flux convergence is calculated by differentiating Eqs. (A3) to (A3) with respect to the latitude (ϕ), and dividing by −2πa2cos (ϕ), i.e.,

$\begin{array}{ll}& {C}_{\mathrm{DSM}}\equiv -\frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }{a}^{\mathrm{2}}\mathrm{cos}\left(\mathit{\varphi }\right)}\frac{\partial {F}_{\mathrm{DSM}}}{\partial \mathit{\varphi }}\\ \text{(A6)}& & \phantom{\rule{1em}{0ex}}=-{a}^{-\mathrm{1}}\underset{\mathrm{0}}{\overset{{p}_{\mathrm{s}}}{\int }}\frac{\partial }{\partial \mathit{\varphi }}\left[\stackrel{\mathrm{‾}}{v}\right]\left[\stackrel{\mathrm{‾}}{{h}_{\mathrm{DS}}}\right]\frac{\mathrm{d}p}{g},& {C}_{\mathrm{DSS}}\equiv -\frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }{a}^{\mathrm{2}}\mathrm{cos}\left(\mathit{\varphi }\right)}\frac{\partial {F}_{\mathrm{DSS}}}{\partial \mathit{\varphi }}\\ \text{(A7)}& & \phantom{\rule{1em}{0ex}}=-{a}^{-\mathrm{1}}\underset{\mathrm{0}}{\overset{{p}_{\mathrm{s}}}{\int }}\frac{\partial }{\partial \mathit{\varphi }}\left[\stackrel{\mathrm{‾}}{{v}^{\ast }}\stackrel{\mathrm{‾}}{{h}_{\mathrm{DS}}^{\ast }}\right]\frac{\mathrm{d}p}{g},\text{(A8)}& & {C}_{\mathrm{DST}}={C}_{\mathrm{DS}}-{C}_{\mathrm{DSM}}-{C}_{\mathrm{DST}},\end{array}$

where CDS is defined in Eq. (3). To calculate the integrals in Eqs. (A3) and (A3), the heat flux quantities are first interpolated from the model's 26 hybrid sigma–pressure levels to 20 equally thick pressure layers (layer midpoints range from 25 to 975 hPa), followed by a numerical integration from the top layer to the surface pressure.

Appendix B: Disentangling longwave feedbacks on the TOA net radiation balance

To estimate the TOA longwave (LW) contributions from changes in the surface temperature ($\mathit{\delta }{R}_{\mathrm{T}}^{\mathrm{LW}}$), water vapor and lapse rate ($\mathit{\delta }{R}_{\mathrm{wv}+\mathrm{lr}}^{\mathrm{LW}}$), and clouds ($\mathit{\delta }{R}_{\mathrm{cld}}^{\mathrm{LW}}$) between two simulations (subscripts 1 and 0), we use the following equations (positive flux downward):

$\begin{array}{}\text{(B1)}& & \mathit{\delta }{R}_{\mathrm{T}}^{\mathrm{LW}}={R}_{s,\mathrm{1}}^{\mathrm{LW}}-{R}_{s,\mathrm{0}}^{\mathrm{LW}},\text{(B2)}& & \mathit{\delta }{R}_{\mathrm{wv}+\mathrm{lr}}^{\mathrm{LW}}=\left({R}_{c,\mathrm{1}}^{\mathrm{LW}}-{R}_{s,\mathrm{1}}^{\mathrm{LW}}\right)-\left({R}_{c,\mathrm{0}}^{\mathrm{LW}}-{R}_{s,\mathrm{0}}^{\mathrm{LW}}\right),\text{(B3)}& & \mathit{\delta }{R}_{\mathrm{cld}}^{\mathrm{LW}}=\left({R}_{\mathrm{1}}^{\mathrm{LW}}-{R}_{c,\mathrm{1}}^{\mathrm{LW}}\right)-\left({R}_{\mathrm{0}}^{\mathrm{LW}}-{R}_{c,\mathrm{0}}^{\mathrm{LW}}\right).\end{array}$

The subscripts “s” and “c” represent surface and clear-sky fluxes, respectively. The clear-sky and total fluxes are taken directly from the model output, and the surface fluxes are computed from the surface temperature using Stefan Boltzmann's law for black body radiation (surface flux proportional to the fourth power of temperature).

An important property of Eqs. (B1) to (B1) is that the individual contributions add up to total LW change: $\mathit{\delta }{R}^{\mathrm{LW}}=\mathit{\delta }{R}_{\mathrm{T}}^{\mathrm{LW}}+\mathit{\delta }{R}_{\mathrm{wv}+\mathrm{lr}}^{\mathrm{LW}}+\mathit{\delta }{R}_{\mathrm{cld}}^{\mathrm{LW}}={R}_{\mathrm{1}}^{\mathrm{LW}}-{R}_{\mathrm{0}}^{\mathrm{LW}}$. Hence, an alternative way to interpret the individual contributions presented here is to first subtract the surface temperature and cloud contributions from the total LW change. The residual term then contains all other LW contributions, including aerosols, greenhouse gases, water vapor, and lapse rate. However, because the aerosol and greenhouse gas concentrations are identical in the LGM simulations, this term only reflects changes in water vapor and lapse rate.

Finally, the combined SW and LW contributions to the TOA net radiation changes are evaluated as

$\begin{array}{}\text{(B4)}& & \mathit{\delta }R=\mathit{\delta }{R}^{\mathrm{SW}}+\mathit{\delta }{R}^{\mathrm{LW}},\text{(B5)}& & \mathit{\delta }{R}_{\mathit{\alpha }}=\mathit{\delta }{R}_{\mathit{\alpha }}^{\mathrm{SW}},\text{(B6)}& & \mathit{\delta }{R}_{\mathrm{wv}+\mathrm{lr}}=\mathit{\delta }{R}_{\mathrm{clr}}^{\mathrm{SW}}+\mathit{\delta }{R}_{\mathrm{wv}+\mathrm{lr}}^{\mathrm{LW}},\text{(B7)}& & \mathit{\delta }{R}_{\mathrm{cld}}=\mathit{\delta }{R}_{\mathrm{cld}}^{\mathrm{SW}}+\mathit{\delta }{R}_{\mathrm{cld}}^{\mathrm{LW}},\text{(B8)}& & \mathit{\delta }{R}_{\mathrm{T}}=\mathit{\delta }{R}_{\mathrm{T}}^{\mathrm{LW}},\end{array}$

where δRSW is the total TOA SW change, and the quantities $\mathit{\delta }{R}_{\mathit{\alpha }}^{\mathrm{SW}}$, $\mathit{\delta }{R}_{\mathrm{clr}}^{\mathrm{SW}}$, and $\mathit{\delta }{R}_{\mathrm{cld}}^{\mathrm{SW}}$ represent the TOA SW contributions from the surface albedo, clear-sky atmosphere (mainly due to changes in the SW absorption by water vapor), and clouds derived from the APRP method .

Competing interests
Competing interests.

The authors declare that they have no conflict of interest.

Acknowledgements
Acknowledgements.

We thank Johan Kleman for providing the LGM ice sheet reconstruction. The computational resources for the numerical simulations were provided by the Centre for Scientific Computing in Frankfurt, Germany. We thank the editor Uwe Mikolajewicz, David Ullman, and the one anonymous reviewer for providing insightful comments on the paper.

Edited by: Uwe Mikolajewicz
Reviewed by: David Ullman and one anonymous referee

References

Abe-Ouchi, A., Saito, F., Kageyama, M., Braconnot, P., Harrison, S. P., Lambeck, K., Otto-Bliesner, B. L., Peltier, W. R., Tarasov, L., Peterschmitt, J.-Y., and Takahashi, K.: Ice-sheet configuration in the CMIP5/PMIP3 Last Glacial Maximum experiments, Geosci. Model Dev., 8, 3621–3637, https://doi.org/10.5194/gmd-8-3621-2015, 2015. a, b

Alexeev, V., Langen, P. L., and Bates, J. R.: Polar amplification of surface warming on an aquaplanet in “ghost forcing” experiments without sea ice feedbacks, Clim. Dynam., 24, 655–666, 2005. a

Berger, A. and Loutre, M.: Astronomical theory of climate change, in: Journal de Physique IV (Proceedings), EDP Sciences, vol. 121, 1–35, 2004. a, b

Bitz, C. M., Shell, K., Gent, P., Bailey, D., Danabasoglu, G., Armour, K., Holland, M., and Kiehl, J.: Climate sensitivity of the community climate system model, version 4, J. Clim., 25, 3053–3070, 2012. a

Braconnot, P., Otto-Bliesner, B., Harrison, S., Joussaume, S., Peterchmitt, J.-Y., Abe-Ouchi, A., Crucifix, M., Driesschaert, E., Fichefet, Th., Hewitt, C. D., Kageyama, M., Kitoh, A., Laîné, A., Loutre, M.-F., Marti, O., Merkel, U., Ramstein, G., Valdes, P., Weber, S. L., Yu, Y., and Zhao, Y.: Results of PMIP2 coupled simulations of the Mid-Holocene and Last Glacial Maximum – Part 1: experiments and large-scale features, Clim. Past, 3, 261–277, https://doi.org/10.5194/cp-3-261-2007, 2007. a

Braithwaite, R. J. and Olesen, O. B.: Calculation of glacier ablation from air temperature, West Greenland, Glac. Quat. G., 6, 219–233, 1989. a

Calov, R. and Greve, R.: A semi-analytical solution for the positive degree-day model with stochastic temperature variations, J. Glaciol., 51, 173–175, 2005. a

Clark, P. U. and Mix, A. C.: Ice sheets and sea level of the Last Glacial Maximum, Quaternary Sci. Rev., 21, 1–7, 2002. a

Colleoni, F., Kirchner, N., Niessen, F., Quiquet, A., and Liakka, J.: An East Siberian ice shelf during the Late Pleistocene glaciations: Numerical reconstructions, Quaternary Sci. Rev., 147, 148–163, https://doi.org/10.1016/j.quascirev.2015.12.023, 2016a. a

Colleoni, F., Wekerle, C., Näslund, J.-O., Brandefelt, J., and Masina, S.: Constraint on the penultimate glacial maximum Northern Hemisphere ice topography (140 kyrs BP), Quaternary Sci. Rev., 137, 97–112, 2016b. a, b, c

Collins, W. D., Bitz, C. M., Blackmon, M. L., Bonan, G. B., Bretherton, C. S., Carton, J. A., Chang, P., Doney, S. C., Hack, J. J., Henderson, T. B., Kiehl, J. T., Large, W. G., McKenna, D. S., Santer, B. D., and Smith, R. D.: The community climate system model version 3 (CCSM3), J. Clim., 19, 2122–2143, 2006a. a

Collins, W. D., Rasch, P. J., Boville, B. A., Hack, J. J., McCaa, J. R., Williamson, D. L., Briegleb, B. P., Bitz, C. M., Lin, S.-J., and Zhang, M.: The formulation and atmospheric simulation of the Community Atmosphere Model: CAM3, J. Clim., 19, 2144–2161, https://doi.org/10.1175/jcli3760.1, 2006b. a

Cook, K. H. and Held, I. M.: Stationary waves of the ice age climate, J. Clim., 1, 807–819, 1988. a, b

Donohoe, A. and Battisti, D. S.: Causes of reduced North Atlantic storm activity in a CAM3 simulation of the Last Glacial Maximum, J. Clim., 22, 4793–4808, 2009. a

Dyke, A., Andrews, J., Clark, P., England, J., Miller, G., Shaw, J., and Veillette, J.: The Laurentide and Innuitian ice sheets during the Last Glacial Maximum, Quaternary Sci. Rev., 21, 9–31, https://doi.org/10.1016/s0277-3791(01)00095-6, 2002. a

Eisenman, I., Bitz, C. M., and Tziperman, E.: Rain driven by receding ice sheets as a cause of past climate change, Paleoceanogr. Paleoclimatol., 24, PA4209, https://doi.org/10.1029/2009PA001778, 2009. a, b, c

Flato, G., Marotzke, J., Abiodun, B., Braconnot, P., Chou, S. C., Collins, W. J., Cox, P., Driouech, F., Emori, S., Eyring, V., Forest, C., Gleckler, P., Guilyardi, É., Jakob, C., Kattsov, V., Reason, C., and Rummukainen, M.: Evaluation of climate models, in: climate change 2013: the physical science basis. Contribution of working group I to the fifth assessment report of the intergovernmental panel on climate change, Climate Change, 5, 741–866, 2013. a

Fyke, J., Sergienko, O., Löfverström, M., Price, S., and Lenaerts, J.: An overview of interactions and feedbacks between ice sheets and the Earth system, Rev. Geophys., 56, https://doi.org/10.1029/2018RG000600, 2018. a

Gong, X., Zhang, X., Lohmann, G., Wei, W., Zhang, X., and Pfeiffer, M.: Higher Laurentide and Greenland ice sheets strengthen the North Atlantic ocean circulation, Clim. Dynam., 45, 139–150, https://doi.org/10.1007/s00382-015-2502-8, 2015. a, b

Gregoire, L. J., Ivanovic, R. F., Maycock, A. C., Valdes, P. J., and Stevenson, S.: Holocene lowering of the Laurentide ice sheet affects North Atlantic gyre circulation and climate, Clim. Dynam., https://doi.org/10.1007/s00382-018-4111-9, 2018. a, b, c

Greve, R.: Application of a polythermal three-dimensional ice sheet model to the Greenland ice sheet: response to steady-state and transient climate scenarios, J. Clim., 10, 901–918, 1997. a

Hall, N. M. J., Dong, B., and Valdes, P. J.: Atmospheric equilibrium, instability and energy transport at the last glacial maximum, Clim. Dynam., 12, 497–511, https://doi.org/10.1007/s003820050123, 1996. a

Hansen, J., Sato, M., and Ruedy, R.: Radiative forcing and climate response, J. Geophys. Res.-Atmos., 102, 6831–6864, 1997. a

He, F.: Simulating transient climate evolution of the last deglaciation with CCSM3, Ph.D. thesis, University of Wisconsin-Madison, 2011. a

Held, I. M. and Soden, B. J.: Robust responses of the hydrological cycle to global warming, J. Clim., 19, 5686–5699, 2006. a

Herrington, A. R. and Poulsen, C. J.: Terminating the Last Interglacial: The role of ice sheet–climate feedbacks in a GCM asynchronously coupled to an ice sheet model, J. Clim., 25, 1871–1882, 2011. a

Hewitt, C., Stouffer, R., Broccoli, A., Mitchell, J., and Valdes, P. J.: The effect of ocean dynamics in a coupled GCM simulation of the Last Glacial Maximum, Clim. Dynam., 20, 203–218, 2003. a

Jakobsson, M., Nilsson, J., Anderson, L., Backman, J., Björk, G., Cronin, T. M., Kirchner, N., Koshurnikov, A., Mayer, L., Noormets, R., O'Regan, M., Stranne, C., Ananiev, R., Macho, N. B., Cherniykh, D., Coxall, H., Eriksson, B., Flodén, T., Gemery, L., Örjan Gustafsson, Jerram, K., Johansson, C., Khortov, A., Mohammad, R., and Semiletov, I.: Evidence for an ice shelf covering the central Arctic Ocean during the penultimate glaciation, Nat. Commun., 7, 10365, https://doi.org/10.1038/ncomms10365, 2016. a, b

Justino, F., Timmermann, A., Merkel, U., and Souza, E. P.: Synoptic Reorganization of Atmospheric Flow during the Last Glacial Maximum, J. Clim., 18, 2826–2846, https://doi.org/10.1175/jcli3403.1, 2005. a

Justino, F., Timmermann, A., Merkel, U., and Peltier, W. R.: An Initial Intercomparison of Atmospheric and Oceanic Climatology for the ICE-5G and ICE-4G Models of LGM Paleotopography, J. Clim., 19, 3–14, https://doi.org/10.1175/jcli3603.1, 2006. a, b

Kageyama, M. and Valdes, P. J.: Impact of the North American ice-sheet orography on the Last Glacial Maximum eddies and snowfall, Geophys. Res. Lett., 27, 1515–1518, https://doi.org/10.1029/1999gl011274, 2000. a, b, c

Kageyama, M., Albani, S., Braconnot, P., Harrison, S. P., Hopcroft, P. O., Ivanovic, R. F., Lambert, F., Marti, O., Peltier, W. R., Peterschmitt, J.-Y., Roche, D. M., Tarasov, L., Zhang, X., Brady, E. C., Haywood, A. M., LeGrande, A. N., Lunt, D. J., Mahowald, N. M., Mikolajewicz, U., Nisancioglu, K. H., Otto-Bliesner, B. L., Renssen, H., Tomas, R. A., Zhang, Q., Abe-Ouchi, A., Bartlein, P. J., Cao, J., Li, Q., Lohmann, G., Ohgaito, R., Shi, X., Volodin, E., Yoshida, K., Zhang, X., and Zheng, W.: The PMIP4 contribution to CMIP6 – Part 4: Scientific objectives and experimental design of the PMIP4-CMIP6 Last Glacial Maximum experiments and PMIP4 sensitivity experiments, Geosci. Model Dev., 10, 4035–4055, https://doi.org/10.5194/gmd-10-4035-2017, 2017. a, b, c

Kleman, J., Fastook, J., Ebert, K., Nilsson, J., and Caballero, R.: Pre-LGM Northern Hemisphere ice sheet topography, Clim. Past, 9, 2365–2378, https://doi.org/10.5194/cp-9-2365-2013, 2013. a, b, c, d, e

Klockmann, M., Mikolajewicz, U., and Marotzke, J.: The effect of greenhouse gas concentrations and ice sheets on the glacial AMOC in a coupled climate model, Clim. Past, 12, 1829–1846, https://doi.org/10.5194/cp-12-1829-2016, 2016. a, b

Lambeck, K., Rouby, H., Purcell, A., Sun, Y., and Sambridge, M.: Sea level and global ice volumes from the Last Glacial Maximum to the Holocene, P. Natl. Acad. Sci. USA, 111, 15296–15303, 2014. a

Langen, P. L. and Vinther, B. M.: Response in atmospheric circulation and sources of Greenland precipitation to glacial boundary conditions, Clim. Dynam., 32, 1035–1054, https://doi.org/10.1007/s00382-008-0438-y, 2008. a

Li, C. and Battisti, D. S.: Reduced Atlantic Storminess during Last Glacial Maximum: Evidence from a Coupled Climate Model, J. Clim., 21, 3561–3579, https://doi.org/10.1175/2007jcli2166.1, 2008. a, b, c, d, e, f

Liakka, J.: Interactions between topographically and thermally forced stationary waves: implications for ice-sheet evolution, Tellus A, 64, 11088, https://doi.org/10.3402/tellusa.v64i0.11088, 2012. a

Liakka, J. and Nilsson, J.: The impact of topographically forced stationary waves on local ice-sheet climate, J. Glaciol., 56, 534–544, https://doi.org/10.3189/002214310792447824, 2010. a

Liakka, J., Nilsson, J., and Löfverström, M.: Interactions between stationary waves and ice sheets: linear versus nonlinear atmospheric response, Clim. Dynam., 38, 1249–1262, https://doi.org/10.1007/s00382-011-1004-6, 2011. a, b

Liakka, J., Löfverström, M., and Colleoni, F.: The impact of the North American glacial topography on the evolution of the Eurasian ice sheet over the last glacial cycle, Clim. Past, 12, 1225–1241, https://doi.org/10.5194/cp-12-1225-2016, 2016. a, b, c, d, e, f

Lindeman, M. and Oerlemans, J.: Northern hemisphere ice sheets and planetary waves: A strong feedback mechanism, J. Climatol., 7, 109–117, https://doi.org/10.1002/joc.3370070202, 1987. a

Liu, Z., Otto-Bliesner, B. L., He, F., Brady, E. C., Tomas, R., Clark, P. U., Carlson, A. E., Lynch-Stieglitz, J., Curry, W., Brook, E., Erickson, D., Jacob, R., Kutzbach, J., and Cheng, J.: Transient Simulation of Last Deglaciation with a New Mechanism for Bolling-Allerod Warming, Science, 325, 310–314, https://doi.org/10.1126/science.1171041, 2009. a

Löfverström, M. and Liakka, J.: On the limited ice intrusion in Alaska at the LGM, Geophys. Res. Lett., 43, 11030–11038, https://doi.org/10.1002/2016gl071012, 2016. a, b

Lofverstrom, M. and Liakka, J.: The influence of atmospheric grid resolution in a climate model-forced ice sheet simulation, The Cryosphere, 12, 1499–1510, https://doi.org/10.5194/tc-12-1499-2018, 2018. a

Löfverström, M. and Lora, J. M.: Abrupt regime shifts in the North Atlantic atmospheric circulation over the last deglaciation, Geophys. Res. Lett., 44, 8047–8055, https://doi.org/10.1002/2017GL074274, 2017. a, b

Löfverström, M., Caballero, R., Nilsson, J., and Kleman, J.: Evolution of the large-scale atmospheric circulation in response to changing ice sheets over the last glacial cycle, Clim. Past, 10, 1453–1471, https://doi.org/10.5194/cp-10-1453-2014, 2014. a, b, c, d, e, f, g

Löfverström, M., Liakka, J., and Kleman, J.: The North American Cordillera – An impediment to growing the continent-wide Laurentide ice sheet, J. Clim., 28, 9433–9450, https://doi.org/10.1175/jcli-d-15-0044.1, 2015. a

Löfverström, M., Caballero, R., Nilsson, J., and Messori, G.: Stationary Wave Reflection as a Mechanism for Zonalizing the Atlantic Winter Jet at the LGM, J. Atmos. Sci., 73, 3329–3342, https://doi.org/10.1175/jas-d-15-0295.1, 2016. a, b, c

Margo Project Members, Waelbroeck, C., Paul, A., Kucera, M., Rosell-Melé, A., Weinelt, M., Schneider, R., Mix, A. C., Abelmann, A., Armand, L., Bard, E., Barker, S., Barrows, T. T., Benway, H., Cacho, I., Chen, M.-T., Cortijo, E., Crosta, X., de Vernal, A., Dokken, T., Duprat, J., Elderfield, H., Eynaud, F., Gersonde, R., Hayes, A., Henry, M., Hillaire-Marcel, C., Huang, C.-C., Jansen, E., Juggins, S., Kallel, N., Kiefer, T., Kienast, M., Labeyrie, L., Leclaire, H., Londeix, L., Mangin, S., Matthiessen, J., Marret, F., Meland, M., Morey, A. E., Mulitza, S., Pflaumann, U., Pisias, N. G., Radi, T., Rochon, A., Rohling, E. J., Sbaffi, L., Schäfer-Neth, C., Solignac, S., Spero, H., Tachikawa, K., and Turon, J.-L.: Constraints on the magnitude and patterns of ocean cooling at the Last Glacial Maximum, Nat. Geosci., 2, 127–132, https://doi.org/10.1038/ngeo411, 2009. a, b

Murakami, S., Ohgaito, R., Abe-Ouchi, A., Crucifix, M., and Otto-Bliesner, B. L.: Global-Scale Energy and Freshwater Balance in Glacial Climate: A Comparison of Three PMIP2 LGM Simulations, J. Clim., 21, 5008–5033, https://doi.org/10.1175/2008jcli2104.1, 2008. a, b, c, d

Niessen, F., Hong, J. K., Hegewald, A., Matthiessen, J., Stein, R., Kim, H., Kim, S., Jensen, L., Jokat, W., Nam, S.-I., and Kang, S.-H.: Repeated Pleistocene glaciation of the East Siberian continental margin, Nat. Geosci., 6, 842–846, https://doi.org/10.1038/ngeo1904, 2013. a

Nilsson, J., Jakobsson, M., Borstad, C., Kirchner, N., Björk, G., Pierrehumbert, R. T., and Stranne, C.: Ice-shelf damming in the glacial Arctic Ocean: dynamical regimes of a basin-covering kilometre-thick ice shelf, The Cryosphere, 11, 1745–1765, https://doi.org/10.5194/tc-11-1745-2017, 2017. a

Oleson, K., Dai, Y., Bonan, G., Bosilovich, M., Dickinson, R., Dirmeyer, P., Hoffman, F., Houser, P., Levis, S., Niu, G.-Y., Thornton, P., Vertenstein, M., Yang, Z.-L., and Zeng, X.: Technical Description of the Community Land Model (CLM), NCAR Tech. Note NCAR/TN-461+STR, Natl. Cent. for Atmos. Res., Boulder, Colorado, 2004. a

Pausata, F. S. R., Li, C., Wettstein, J. J., Kageyama, M., and Nisancioglu, K. H.: The key role of topography in altering North Atlantic atmospheric circulation during the last glacial period, Clim. Past, 7, 1089–1101, https://doi.org/10.5194/cp-7-1089-2011, 2011. a, b

Peixoto, J. P. and Oort, A. H.: Physics of climate, New York, NY (United States), American Institute of Physics, 1992. a, b, c

Peltier, W.: Global glacial isostasy and the surface of the ice-age Earth: the ICE-5G (VM2) model and GRACE, Annu. Rev. Earth Pl. Sc., 32, 111–149, https://doi.org/10.1146/annurev.earth.32.082503.144359, 2004. a, b, c

Peltier, W. and Fairbanks, R. G.: Global glacial ice volume and Last Glacial Maximum duration from an extended Barbados sea level record, Quaternary Sci. Rev., 25, 3322–3337, 2006. a

Petit, J.-R., Jouzel, J., Raynaud, D., Barkov, N. I., Barnola, J.-M., Basile, I., Bender, M., Chappellaz, J., Davis, M., Delaygue, G., Delmotte, M., Kotlyakov, V. M., Legrand, M., Lipenkov, V. Y., Lorius, C., Pépin, L., Ritz, C., Saltzman, E., and Stievenard, M.: Climate and atmospheric history of the past 420,000 years from the Vostok ice core, Antarctica, Nature, 399, 429–436, 1999. a, b

Reeh, N.: Parameterization of melt rate and surface temperature on the Greenland ice sheet, Polarforschung, 59, 113–128, 1991. a

Rivière, G., Berthou, S., Lapeyre, G., and Kageyama, M.: On the reduced North Atlantic storminess during the last glacial period: the role of topography in shaping synoptic eddies, J. Clim., 31, 1637–1652, 2018.  a

Roe, G. H. and Lindzen, R. S.: The mutual interaction between continental-scale ice sheets and atmospheric stationary waves, J. Clim., 14, 1450–1465, https://doi.org/10.1175/1520-0442(2001)014<1450:TMIBCS>2.0.CO;2, 2001. a, b, c, d

Serreze, M. C., Barrett, A. P., Slater, A. G., Steele, M., Zhang, J., and Trenberth, K. E.: The large-scale energy budget of the Arctic, J. Geophys. Res.-Atmos., 112, D11122, https://doi.org/10.1029/2006JD008230, 2007. a, b, c, d

Shin, S.-I., Liu, Z., Otto-Bliesner, B., Brady, E., Kutzbach, J., and Harrison, S.: A simulation of the Last Glacial Maximum climate using the NCAR-CCSM, Clim. Dynam., 20, 127–151, 2003. a

Spahni, R., Chappellaz, J., Stocker, T. F., Loulergue, L., Hausammann, G., Kawamura, K., Flückiger, J., Schwander, J., Raynaud, D., Masson-Delmotte, V., and Jouzel, J.: Atmospheric Methane and Nitrous Oxide of the Late Pleistocene from Antarctic Ice Cores, Science, 310, 1317–1321, https://doi.org/10.1126/science.1120132, 2005. a, b

Svendsen, J. I., Alexanderson, H., Astakhov, V. I., Demidov, I., Dowdeswell, J. A., Funder, S., Gataullin, V., Henriksen, M., Hjort, C., Houmark-Nielsen, M., Hubberten, H. W., Ingólfsson, Ó., Jakobsson, M., Kjær, K. H., Larsen, E., Lokrantz, H., Lunkka, J. P., Lyså, A., Mangerud, J., Matiouchkov, A., Murray, A., Möller, P., Niessen, F., Nikolskaya, O., Polyak, L., Saarnisto, M., Siegert, C., Siegert, M. J., Spielhagen, R. F., and Stein, R.: Late Quaternary ice sheet history of northern Eurasia, Quaternary Sci. Rev., 23, 1229–1271, https://doi.org/10.1016/s0277-3791(03)00342-1, 2004. a, b

Taylor, K., Crucifix, M., Braconnot, P., Hewitt, C., Doutriaux, C., Broccoli, A., Mitchell, J., and Webb, M.: Estimating shortwave radiative forcing and response in climate models, J. Clim., 20, 2530–2543, 2007. a, b

Trenberth, K. E., Caron, J. M., and Stepaniak, D. P.: The atmospheric energy budget and implications for surface fluxes and ocean heat transports, Clim. Dynam., 17, 259–276, 2001. a, b

Ullman, D., LeGrande, A., Carlson, A., Anslow, F., and Licciardi, J.: Assessing the impact of Laurentide Ice Sheet topography on glacial climate, Clim. Past, 10, 487–507, https://doi.org/10.5194/cp-10-487-2014, 2014. a, b, c, d

Wekerle, C., Colleoni, F., Näslund, J.-O., Brandefelt, J., and Masina, S.: Numerical reconstructions of the penultimate glacial maximum Northern Hemisphere ice sheets: sensitivity to climate forcing and model parameters, J. Glaciol., 62, 607–622, https://doi.org/10.1017/jog.2016.45, 2016. a, b

Zhang, X., Lohmann, G., Knorr, G., and Purcell, C.: Abrupt glacial climate shifts controlled by ice sheet changes, Nature, 512, 290–294, 2014. a, b, c

Zhu, J., Liu, Z., Zhang, X., Eisenman, I., and Liu, W.: Linear weakening of the AMOC in response to receding glacial ice sheets in CCSM3, Geophys. Res. Lett., 41, 6252–6258, 2014. a, b

The kinetic energy is neglected as it is much smaller (typically 2 orders of magnitude) than the dry-static and latent energy contributions.