# Assessing the sensitivity of the North Atlantic Ocean circulation to freshwater perturbation in various glacial climate states

Cédric van Meerbeeck, Didier M. Roche, Hans Renssen

► To cite this version:

Cédric van Meerbeeck, Didier M. Roche, Hans Renssen. Assessing the sensitivity of the North Atlantic Ocean circulation to freshwater perturbation in various glacial climate states. Climate Dynamics, 2011, 37 (9-10), pp.1909-1927. 10.1007/s00382-011-1043-z. hal-03203455

HAL Id: hal-03203455

https://hal.science/hal-03203455

Submitted on 21 Apr 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Assessing the sensitivity of the North Atlantic Ocean circulation to freshwater perturbation in various glacial climate states

Cédric J. Van Meerbeeck · Didier M. Roche ·
Hans Renssen

Received: 26 March 2010 / Accepted: 2 March 2011 / Published online: 11 March 2011
© The Author(s) 2011. This article is published with open access at Springerlink.com

**Abstract** A striking characteristic of glacial climate in the North Atlantic region is the recurrence of abrupt shifts between cold stadials and mild interstadials. These shifts have been associated with abrupt changes in Atlantic Meridional Overturning Circulation (AMOC) mode, possibly in response to glacial meltwater perturbations. However, it is poorly understood why they were more clearly expressed during Marine Isotope Stage 3 (MIS3, ~60–27 ka BP) than during Termination 1 (T1, ~18–10 ka BP) and especially around the Last Glacial Maximum (LGM, ~23–19 ka BP). One clue may reside in varying climate forcings, making MIS3 and T1 generally milder than LGM. To investigate this idea, we evaluate in a climate model how ice sheet size, atmospheric greenhouse gas concentration and orbital insolation changes between 56 ka BP (=56k), 21k and 12.5k affect the glacial AMOC response to additional freshwater forcing. We have performed three ensemble simulations with the earth system model LOVECLIM using those forcings. We find that the AMOC mode in the mild glacial climate type (56k and 12.5k), with deep convection in the Labrador Sea and the Nordic Seas, is more sensitive to a constant 0.15 Sv freshwater forcing than in the cold type (21k), with deep convection mainly south of Greenland and Iceland. The

initial AMOC weakening in response to freshwater forcing is larger in the mild type due to an early shutdown of Labrador Sea deep convection, which is completely absent in the 21k simulation. This causes a larger fraction of the freshwater anomaly to remain at surface in the mild type compared to the cold type. After 200 years, a weak AMOC is established in both climate types, as further freshening is compensated by an anomalous salt advection from the (sub-)tropical North Atlantic. However, the slightly fresher sea surface in the mild type facilitates further weakening of the AMOC, which occurs when a surface buoyancy threshold (–0.6 kg m⁻³ surface density anomaly to the 56k reference state) is stochastically crossed in the Nordic Seas. While described details are model-specific, our results imply that a more northern location of deep convection sites during milder glacial times may have amplified frequency and amplitude of abrupt climate shifts.

**Keywords** Atlantic meridional overturning circulation · Marine isotope stage 3 · Last glacial maximum · North Atlantic Ocean · Climate modelling

## 1 Introduction

Climate in Greenland witnessed 25 abrupt, millennial-scale shifts between cold stadials and mild interstadials during the Last Glacial period (~110–10 ka BP; Johnsen et al. 1992; Dansgaard et al. 1993; NorthGRIP Members 2004), of which those during the interval 70–10 ka BP are shown in Fig. 1. In parallel, abrupt shifts in ocean circulation and climate—that appear to be related to the Greenland climate shifts—affected the North Atlantic region (e.g. Bond et al. 1993; van Kreveld et al. 2000; Voelker 2002; Sánchez Goñi et al. 2002; Wohlfarth et al. 2008). However, the

<footnote>C. J. Van Meerbeeck (✉) · D. M. Roche · H. Renssen
Section Climate Change and Landscape Dynamics, Department of Earth Sciences, VU University Amsterdam, De Boelelaan 1085, 1081 HV Amsterdam, The Netherlands
e-mail: cedric.van.meerbeeck@falw.vu.nl</footnote>

<footnote>D. M. Roche
Laboratoire des Sciences du Climat et de l'Environnement (LSCE/IPSL), Laboratoire CEA/INSU-CNRS/UVSQ, C.E. de Saclay, Orme des Merisiers Bat. 701, 91190 Gif sur Yvette Cedex, France</footnote>
**Fig. 1** The NorthGRIP ice core record of δ¹⁸O for the last 70,000 years as a proxy for central-Greenland temperature (adapted from NorthGRIP Members 2004). Approximate time spans of the Last Glacial Termination (LGT), Last Glacial Maximum (LGM) and Marine Isotope Stage 3 (MIS 3) are shaded grey

<figure><img src="image_2.png" /></figure>

cause and climatic expression of these shifts are at present not well understood (Clement and Peterson 2008). Furthermore, it is not well established why the observed abrupt shifts were more clearly expressed during Marine Isotope Stage 3 (MIS 3, ~60–27 ka BP; e.g. Huber et al. 2006) and glacial Termination 1 (T1, ~18–10 ka BP; e.g. Bard et al. 2000; McManus et al. 2004; Steffensen et al. 2008), than during the Last Glacial Maximum (LGM, ~23–19 ka BP as defined by the EPILOG group in Mix 2003).

Broecker et al. (1985) suggested that fluctuations in the strength of the Atlantic Thermohaline Circulation (THC) caused the abrupt climate shifts, with cold stadials resulting from a weak or collapsed THC and mild interstadials from a subsequent resumption of the THC. Broecker et al. (1985) followed the model of THC bi-stability of Stommel (1961) to further develop the theory of Rooth (1982) who suggested that freshwater perturbation of the high-latitude North Atlantic during T1 may have caused the Younger Dryas stadial (YD, ~12.7–10.5 ka BP).

Numerous palaeoclimate modelling studies support THC-driven millennial-scale abrupt glacial climate change through shifts in the Atlantic Meridional Overturning Circulation (AMOC) between different modes and strengths (e.g. Sakai and Peltier 1997; Ganopolski and Rahmstorf 2001; Schmittner and Clement 2002; Knutti et al. 2004; Flückiger et al. 2006, Kageyama et al. 2010). At times of vigorous overturning, the shallow component of the Atlantic Meridional Overturning Circulation (AMOC) is characterised by a net northward water flow, advecting heat and salt northwards. At northern high latitudes, deep waters are formed by means of deep convection. Subsequently, these deep waters flow back south as part of the THC.

Deep convection occurs mainly in two forms, which at present-day result in North Atlantic Deep Water (NADW) formation in the North Atlantic Ocean, the Nordic Seas and the Labrador Sea (Dickson and Brown 1994). Firstly, open ocean convection occurs when relatively warm and saline Atlantic surface waters encounter much colder and fresher Arctic surface waters. At a certain point, the Atlantic water masses have lost enough heat to the atmosphere to obtain

equal density to that of Arctic waters. At that point, these masses mix to create a denser solution than their end components. Consequently, deep mixing occurs in the water column. Secondly, deep convection can also take place under sea-ice cover as a result of sea-ice formation, where brine rejection during the freezing increases density of the underlying surface waters up to the point that deep convection takes place (Kuhlbrodt et al. 2007).

With water column temperatures mostly near freezing in high northern latitudes at glacial times (Adkins et al. 2002; Paul and Schäfer-Neth 2003; Meland et al. 2005; MARGO Project Members 2009), surface water density could not rise much by heat loss to the atmosphere and sink. Therefore, deep convection in the North Atlantic Ocean, the Labrador Sea and the Nordic Seas in a glacial climate may have been particularly sensitive to differences in sea surface salinity (SSS) (Stommel 1961; Rahmstorf 2002). In a relatively simple coupled climate model (CLIMBER2), the glacial AMOC was found to be sensitive to small SSS changes, leading to a shallower AMOC and a southward shift in convection in an LGM state compared to modern climate (Ganopolski et al. 1998; Ganopolski and Rahmstorf 2001; Rahmstorf 2002). In the same model, Ganopolski and Rahmstorf (2002) found that a small periodical freshwater forcing and an additional freshwater flux in form of white noise were sufficient to stochastically shift the AMOC mode between a weak and a strong glacial mode. Similarly, stochastic AMOC mode shifts have also been found by amongst others Schulz et al. (2007) and Jongma et al. (2007) in Holocene simulations with a more comprehensive coupled model (LOVECLIM). Since the ocean component of the latter model is fully three-dimensional, Schulz et al. (2007) and Jongma et al. (2007) could distinguish that these mode shifts are characterised either by presence or absence of Labrador Sea Water (LSW) formation. With respect to simulated LGM climate, Roche et al. (2007) found no significant change in AMOC strength between LGM and modern climate in LOVECLIM. Nevertheless, clear differences in the location of deep convection occurred in that LGM state, with no LSW
formation—consistent with palaeodata (Hillaire-Marcel et al. 2001)—and only intermittent Nordic Seas Deep Water formation, but an enhanced convection cell south of Greenland and Iceland. Van Meerbeeck et al. (2009) showed that the AMOC strength in early MIS 3 simulations with LOVECLIM did not vary much from LGM strength either. However, deep convection was mainly shifted to the Nordic Seas and the Labrador Sea in the MIS 3 simulations.

While it has been shown for the Holocene that slight
freshwater forcing changes may lead to AMOC mode shifts
(Puinel-Cottet et al. 2004; Renssen et al. 2005; Jongma
et al. 2007), it is at present unclear how glacial background
climate controls the NADW formation. More specifically,
we do not know the impact of climate forcing factors on the
location and intensity of NADW production in glacial
times. To improve our understanding of this issue, we
therefore simulate the effect of changing glacial climate
forcings on AMOC weakening in response to freshwater
forcing in LOVECLIM.

A further subject of debate has been the timing of abrupt
glacial climate shifts in Greenland. In particular, stadial-
interstadial transitions exhibit an apparent periodicity of
1470 years—and multiples thereof—in the GISP2 ice core
record (Grootes and Stuiver 1997; Alley et al. 2001).
Several arguments have been proposed in favour of and
against such periodicity (e.g. Wunsch 2000; Alley et al.
2001; Schulz 2002; Ditlevsen et al. 2007). However, it has
been suggested and shown in coupled atmosphere-ocean
models that abrupt millennial-scale AMOC shifts may be
triggered by stochastic processes in a glacial climate
(Weaver et al. 1993). In this view, the AMOC is a highly
non-linear system which can be shifted from one mode to
another by crossing a buoyancy threshold—i.e. freshwater
perturbation—(Kuhlbrodt et al. 2007). By running transient
ensemble simulations of glacial AMOC strength shifts, we
aim to give a first estimation of the role of internal vari-
ability of the coupled atmosphere-ocean system in crossing
a buoyancy threshold beyond which the AMOC shuts
down.

Finally, combining the findings on AMOC mode shifts
under varying climate forcings and on buoyancy threshold
crossings, we search for an explanation of the higher fre-
quency and more pronounced expression of abrupt glacial
climate shifts during MIS 3 and T1 than during the LGM.

**2 Methods**

## 2.1 Model

We use the three-dimensional coupled earth system model
of intermediate complexity LOVECLIM 1.0 (Driesschaert

et al. 2007). In this study, we only use three coupled
components, namely the atmospheric ECBilt and oceanic
CLIO components, and the vegetation module VECODE.

ECBilt is a quasi-geostrophic, T21 spectral model, with three vertical levels (Opsteegh et al. 1998). Due to its spectral resolution corresponding to a horizontal ~5.6° × ~5.6° resolution, its surface topography is simplified. Its parameterisation scheme allows for fast computing and includes a linear longwave radiation scheme. ECBilt contains a full hydrological cycle, including a simple bucket model for soil moisture over continents, and computes synoptic variability associated with weather patterns. Precipitation falls in the form of snow when surface air temperatures fall below 0°C. CLIO is a primitive-equation three-dimensional, free-surface ocean general circulation model coupled to a thermodynamical and dynamical sea-ice model (Goosse and Fichefet 1999). CLIO has a realistic bathymetry, a 3° × 3° horizontal resolution and 20 vertical levels. The free-surface of the ocean allows introduced freshwater fluxes to change sea level (Tartinville et al. 2001). In order to bring precipitation amounts in ECBilt-CLIO closer to observations, a negative precipitation-flux correction is applied over the Atlantic and Arctic Oceans to correct for excess precipitation. This flux is reintroduced in the North Pacific. The climate sensitivity of LOVECLIM 1.0 to a doubling of atmospheric CO₂ concentration is 1.8°C, associated with a global radiative forcing of 3.8 Wm⁻² (Driesschaert 2005). The dynamic terrestrial vegetation model VECODE computes the surface fraction of each land grid cell covered by herbaceous plants, trees and desert fractions (Brovkin et al. 1997). The potential vegetation is calculated from annual mean surface air temperatures and annual precipitation sums from ECBilt. VECODE then feeds back to ECBilt through the surface albedo (Driesschaert 2005).

LOVECLIM produces a generally realistic modern cli-
mate and ocean circulation (Driesschaert 2005) and an
LGM climate generally consistent with data (Roche et al.
2007). Aside from a slightly stronger and shallower AMOC
than in modern climate, the LGM ocean circulation is
characterised by too warm sea surface temperatures (SST)
in the eastern North Atlantic at LGM compared to modern
climate (Roche et al. 2007; MARGO Project Members
2009). Although LOVECLIM is best suited for studies
focusing mainly on extra-tropical climate dynamics
(because vertical convection is parameterised in ECBilt),
the model simulates strength changes and latitudinal shifts
of the tropical rain belt due to AMOC shifts broadly con-
sistently with general circulation models (e.g. Zhang and
Delworth 2005). ECBilt can thus explicitly simulate
atmospheric teleconnections involved in transferring cli-
matic signals across different parts of the globe, though
details may not be well-resolved. It should be noted that the
freshwater correction flux in the model implies that the zero point of the freshwater budget of the ocean is not absolute. Furthermore, this zero point is also model dependent, resulting in different amounts of freshwater perturbation needed to shut down the AMOC (Rahmstorf et al. 2005). Nevertheless, the quasi-equilibrium AMOC response to freshwater forcing in terms of hysteresis shape in LOVECLIM is similar to other models, with or without flux corrections (Rahmstorf et al. 2005). Therefore, the AMOC response to freshwater perturbation appears at least qualitatively model-independent, which allows us to assess the sensitivity in different climate states.

## 2.2 Experimental design

In order to obtain our reference climate states for early-MIS 3 (~56 ka BP) and LGM (~21 ka BP), we spun up the model with constant climate forcings for 7000 and 5000 years to near-equilibrium. The YD state (12.5k-REF) was obtained by running the model forward in time from the LGM state to 12.5 ka BP with (1) transient atmospheric GHG concentrations, (2) orbital forcing and (3) the topography and albedo changes due to ice sheet forcing is interpolated with 50-year time steps from the ICE-5G VM2 reconstruction of Peltier (2004). The forcings are shown on Fig. 2. GHG concentrations were taken from Flückiger et al. (1999), Indermühle et al. (2000), Monnin et al. (2001) and Flückiger et al. (2004); orbital parameters from Berger (1978) and Berger and Loutre (1992); the LGM land-sea mask and ice sheets from Roche et al. (2007); and finally the MIS 3 ice sheet topography from Van Meerbeeck et al. (2009). The three reference states are called 56k-REF, 21k-REF and 12.5k-REF. Climatologies of the reference states presented in this paper are based on the last 100 years of the experiments.

The three transient simulation ensembles—consisting of five members each—were designed to represent rapid glacial warming transitions from stadials to interstadials. We started all 56k, 21k and 12.5k ensembles from their respective reference state, but with a slightly modified initial atmospheric state for each ensemble member. We then ran the model for 500 years with all forcings kept constant, the only forcing change from the reference states being additional freshwater forcing, which is not compensated elsewhere and thus changes sea level. After this, 1,000 years were further computed for all simulations, but now with time-varying atmospheric GHG concentrations and orbital forcing. The ice sheet topography was kept constant throughout the transient experiments. All forcing schemes (except the fixed LGM land-sea mask) are shown in Fig. 2.

Positive freshwater fluxes were equally distributed over the sea surface throughout the first 900 years, with 0.07 Sv

(1 Sv = 1 Sverdrup = $10^6 \text{ m}^3 \text{ s}^{-1}$) added to an area of the mid-latitude North Atlantic Ocean and 0.08 Sv to the Fennoscandian Ice Shelf (for area definitions, see Roche et al. 2010). This was done to slow down the AMOC to a ‘weak’ mode. From year 901 to year 1300, the respective freshwater fluxes were 0.21 Sv (North Atlantic) and 0.04 Sv (Fennoscandian Ice Shelf) and finally from year 1301 onwards –0.1 Sv in each area. The relatively large positive freshwater forcings values were imposed to reach an AMOC ‘off’ mode in all simulations. The simulated North Atlantic Ocean circulation, obtained by applying such forcing, mimics that of stadials perturbed by Heinrich events (Van Meerbeeck 2010). The subsequent negative freshwater flux would ensure a recovery of the AMOC to a strong ‘on’ mode with stronger overturning than in the reference state. As such, the strong AMOC ‘on’ mode should result in a North Atlantic Ocean circulation resembling that of interstadials (Van Meerbeeck 2010). Note that the global mean ocean salinity in the 56k simulations is ~1.5 psu (psu = practical salinity unit) lower than in the 21k and 12.5k simulations. However, sensitivity tests we performed, show that a global mean ocean salinity increase or decrease of 1 psu or 2 psu does not significantly influence the mean state of the AMOC. This is because salinity differences remain uniform across depth, latitude and longitude. Therefore, to consistently discuss salinity, we will use anomalies compared to reference states.

In conclusion, the above detailed experimental setup is necessary to carefully estimate location and intensity of NADW formation in different glacial background climate states. Such climate simulations may appear better fit for model-data comparisons, which has been undertaken with a 56k simulation with nearly identical forcings elsewhere (Van Meerbeeck 2010). However, presence or absence of NADW formation in the Labrador Sea in LOVECLIM greatly depends on the physical state of the glacial North Atlantic Ocean (Roche et al. 2007; Van Meerbeeck et al. 2009). Therefore, subtle differences in forcings will deeply affect the process and thereby the AMOC. We thus argue that a thorough analysis of processes at play during AMOC weakening and shutdown in glacial climates requires a detailed set of forcings.

# 3 Results and discussion

## 3.1 North Atlantic surface ocean conditions and overturning circulation in three glacial reference climates

At present-day, February is the coldest month at sea surface throughout most of the North Atlantic Ocean, the Labrador
<figure><img src="image_3.png" /><figcaption>**Fig. 2** Forcing schemes for the three ensemble simulations. The initial values of the forcings are the same as those of the reference climates each simulation is started from. **Left panel:** Forcing scenarios for **a** freshwater forcing ($Sv = 10^6 m^3 s^{-1}$); atmospheric **b** $CH_4$, $C_2O_4$ and $CO_2$ concentrations; and **c** 60°N February and August insolation anomalies to present-day ($W m^{-2}$). **Right panel:** Ice sheets additional topography compared to present-day (km) for **d** MIS 3, **e** LGM and **f** LGT</figcaption></figure>

Sea and the Nordic Seas. Due to the maximum winter cooling, water surface densities reach their annual maximum in this month, thus destabilising the water column and allowing for deep convection. We therefore focus on sea surface conditions in February.

The three reference states represent quasi-equilibrium climates with earth surface conditions ranging from cold to mild glacial (Fig. 3). In response to larger and taller Northern Hemisphere ice sheets, lowest GHG concentrations and lowest summer insolation, 21k-REF shows the coldest February sea surface temperatures (SST) over the North Atlantic (Fig. 3c). This is consistent with data that indicate overall very cold glacial conditions during the LGM chronozone (see e.g. Kucera et al. 2005; MARGO Project Members 2009). Due to smaller ice sheets,

enhanced summer insolation and slightly higher atmospheric GHG concentrations, February SSTs in 56k-REF are 1–6°C higher between ~50°N and ~70°N (Fig. 3a), while even smaller ice sheets and higher $CO_2$ concentrations make 12.5k-REF the mildest, being still 0.5–2°C warmer than 56k-REF regionally (Fig. 3e). As a consequence of these SST differences, the winter sea-ice edge has the most northern position in 12.5k-REF (Fig. 4e), while the sea-ice cover is largest in 21k-REF, reaching 5–20° of latitude further south compared to the other two reference states (Fig. 4c).

Based on differences in simulated SSTs as well as the geographic distribution of deep convection sites and sea-ice edge, the three reference climate states can be classified in 2 glacial climate types. The first type consists of only the
<caption>**Fig. 3** Average February North Atlantic Ocean surface conditions in the reference states (**left panels**) and three 56k ensemble members in years 688–712 (**right panels**).</caption>

a February Sea Surface Temperatures (SSTs) in 56k-REF; c and e SST anomalies of 21k-REF, respectively 12.5k-REF to 56k-REF. During years 688–712, b 56k-stb is characterised by a weak Atlantic Meridional Overturning Circulation (AMOC), d 56k-shd by a shut down AMOC and f 56k-red by a weak AMOC diverging to shutdown. SST is shown in color (°C) in all panels, while the Sea Surface Salinity (SSS) anomalies to 56k-REF in the right panels are shown as black contours (psu)

<figure><img src="image_4.png" /></figure>

21k-REF, namely the cold glacial climate type (Fig. 3c). In this type, deep convection mostly takes place in the open ocean south of Greenland and Iceland (Fig. 4c). In addition, deep convection episodically occurs in proximity of the mean sea-ice edge in the Norwegian Sea or due to sea-ice formation between Greenland and Spitsbergen. Sea-ice covers the entire Labrador Sea in February and even most of it in August. Further east, the sea-ice edge progressively lies more northward, going from ~55°N south of Greenland, to Iceland’s south coast and finally to between 65 and 70°N in the Nordic Seas.

In the mild glacial type—consisting of 56k-REF and 12.5k-REF, surface waters are warmer than in the cold type (Fig. 3a, e). Open waters are thus found further north than in the cold type (Fig. 4a, e). Consequently, at high latitudes the surface waters can gain density by heat loss to the atmosphere, increasing open ocean convection in the Nordic Seas, and allowing it in the Labrador Sea.

Concomitantly, deep convection south of Greenland and Iceland is much reduced.

Although the shift in convection sites between the cold and the mild glacial climate type does not result in remarkable differences in AMOC strength—with an average overturning of around 15 Sv near 60°N in all three reference states (Fig. 5a, c, e), neither in northward oceanic heat flux—~0.33 10¹⁵W in the Atlantic at 30°S in all three—three notable differences in AMOC mode may be distinguished. First, the main overturning cell—which extends from the southern end of the Atlantic Ocean at 33°S to around 60°N—is deeper and overall slightly more vigorous southwards of 50°N in 21k-REF (Fig. 5c) than in 56k-REF and 12.5k-REF (Fig. 5a, e). This overturning cell transports surface and shallow water masses northwards roughly between 0 m and 1,000 m depth, and intermediate and deep water masses southwards mostly between 2,000 m and 4,000 m depth.
**Fig. 4** February southern sea-ice edge (15% ice concentration, *dark grey* contours) and deep convection sites (*black* contours) in the North Atlantic. **a, c, e** the reference states; **b, d, f** three 56k ensemble members in years 688–712. Unshaded grid cells represent land surface only in the model, whereas cells shaded grey contain at least a fraction of sea surface

<figure><img src="image_5.png" /></figure>

This cell exports NADW out of the Atlantic Ocean at a rate of 17.3 Sv versus 16.4 Sv in 56k-REF and 15.9 Sv in 12.5k-REF at 30°S. (NADW export is a common measure of the AMOC strength.) Second, because deep convection is absent in the Labrador Sea and reduced in the Nordic Seas in 21k-REF, the secondary overturning cell between 65°N and 80°N is weaker, producing 1.9 Sv of Nordic Seas Deep-Water (NSDW), compared to 2.4 Sv in 56k-REF and 2.8 Sv 12.5k-REF. Third, the counter-clockwise deepest overturning cell centred at ~4,500 m water depth is somewhat weaker in 21k-REF compared to 56k-REF, and further reduced compared to 12.5k-REF. As a result, the latter cell imports Antarctic Bottom Water (AABW) northwards into the Atlantic at a rate of 1.4 Sv in 21k-REF, compared to 3.7 Sv in 56k-REF and even 4.8 Sv in 12.5k-REF. Transporting more AABW into the Atlantic Ocean in 12.5k, the deepest overturning cell reaches further north and occupies an up to 500 m larger fraction of the water column in 12.5k-REF than 56k-REF (Fig. 5).

## 3.2 Dynamical AMOC response to freshwater forcing in three ensemble simulations

The temporal evolution of North Atlantic Deep Water export out of the Atlantic (at 30°S) in the ensembles of transient glacial simulations shows strong similarities between the three climate states (Fig. 6). Firstly, within less than 200 years from the start of each simulation, NADW export is approximately halved to ~7–10 Sv compared to the respective reference states (i.e. at time = 0). Secondly, from year 200 to 900, NADW export reduction is much more gradual, with only ~1.5 Sv, ~2 Sv and ~3 Sv reduction in the 21k, 12.5k and 56k ensembles, respectively. Thirdly, between years 1050 and 1350 all simulations have a shutdown AMOC, with only ~3 Sv NADW export remaining. Lastly, after year 1350, NADW export resumes in all simulations with a maximum increase rate of ~8 Sv/100 years around year 1500.

However, the cold and a mild glacial climate types differ in terms of AMOC response to freshwater forcing. It is
<figure><img src="image_6.png" /><figcaption>**Fig. 5** The AMOC stream function in the reference states (**left panels**) and of three 56k ensemble members in years 688–712 (**right panels**). x is latitude, y water depth. Contour alues denote the vertically integrated water flow (where positive means a clockwise motion)</figcaption></figure>

clear from Fig. 6a that after 100 years, NADW export is reduced by ~7 Sv on average in the 56k and 12.5k ensembles, compared to only ~4.5 Sv in the 21k ensemble. Also, the reduction of NADW export after 900 years is smallest for 21k and largest for 56k. Finally, despite the differences in average NADW export at year 900, a final reduction to 2.5 Sv is reached within ~150 years in the three ensembles. The 0.25 Sv freshwater forcing after year 900 is apparently sufficient to shut down the glacial AMOC in our model.

By looking at NADW export changes in the individual ensemble members (Fig. 6b, c, d), we can categorise the behaviour of the overturning strength between years 200 and 900 into two types, namely A and B (shown in Fig. 6b–d). In type A, the weak AMOC is characterised by a 7–10 Sv NADW export persisting until year 900—i.e. until stronger freshwater forcing is applied; in type B, the

weak AMOC veers to shutdown before year 900. Of type A, 4 are found in the 21k ensemble (Fig. 6c), 2 in the 12.5k ensemble (Fig. 6d) and 1 in the 56k ensemble (Fig. 6b). Type B can be further subdivided in those members with a full shutdown AMOC before year 900 and those without. A relatively early full shutdown (before year 900) is found in 3 members of the 56k ensemble, 2 members of the 12.5k ensemble, but in none of the 21k members. In conclusion, with mostly a type A AMOC mode in response to the 0.15 Sv freshwater forcing, the 21k ensemble or cold glacial type appears to be less sensitive than the mild glacial type.

Since Labrador Sea deep convection shuts down after 9–20 years in the 56k and 38–45 years in the 12.5k members, whereas LSW formation is absent in the 21k-REF state, a 2.5 Sv larger NADW export reduction in the first 100 years is noted in the mild type. Earlier shutdown
**Fig. 6** AMOC strength changes in terms of North Atlantic Deep Water (NADW) export out of the Atlantic Ocean at 30°S (in Sv) in response to the freshwater forcing (scenario depicted in Fig. 2a) for the three glacial ensemble simulations.

a Ensemble means; b 56k ensemble (the AMOC mode in the three highlighted members is examined in detail for years 688–712); c 21k ensemble; d 12.5k ensemble. A and B refer to weak AMOC type A (remaining at 7–10 Sv until year 900), respectively type B (diverging from 7 to 10 Sv towards shutdown before year 900)

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

in the Labrador Sea in 56k than 12.5k is partly explained by SSTs closer to the freezing point (Fig. 3a, c), implying that surface water density cannot increase much by radiative cooling. Therefore, advective buoyancy—by which we mean the freshwater forcing—appears to have greater effect in 56k than in 12.5k.

More importantly, the first site displaying a convection shutdown in the mild glacial climate type is the Labrador Sea mainly because (1) SSTs ~6°C closer to the freezing point than south of Iceland (see Fig. 3e), implying that surface water density cannot increase much by radiative cooling; and (2) by the smaller influx of highly saline Atlantic water masses from the sub-tropics than at the other convection sites. This is because the latter sites lie in the pathway of the warm and saline Gulf Stream - North Atlantic Current system.

Then, after deep convection collapses in the Labrador Sea, the weaker overturning in 56k and 12.5k compared to 21k results in less drawdown of the advective buoyancy from the surface into deeper ocean layers. This slower freshwater removal at surface constitutes a positive convective buoyancy feedback to the advective weakening of the AMOC (see Fig. 9). The combined effect prevents deep convection in the Labrador Sea and the Nordic Seas from re-invigorating while the AMOC is weak.

Surface freshening due to the positive advective and convective buoyancy feedbacks in the Nordic Seas is shown in Fig. 8a, c, e. This figure shows that SSS gradually decreases in three 56k members until year ~200 when the weak AMOC is established. It also shows that surface buoyancy gain (density loss) is mostly due to SSS reduction, with cooling only having a weak, opposite effect. For the Labrador Sea convection site, a similar process is noted, though obviously at a faster pace there than in the Nordic Seas, with a 1.5 kg m⁻³ surface density (ρsurface) decrease compared to only ~0.4 kg m⁻³ in the Nordic Seas and North Atlantic sites where freshwater is actually applied. The simulated fast convection shutdown in the Labrador—constituting a positive buoyancy gain—thus is the main cause of a faster and larger AMOC weakening by year 200 in the mild type.

In the mild glacial type, the combined convective and advective increase in surface buoyancy accelerates weakening of deep convection in the Nordic Seas (Fig. 9). Here, the maximum overturning rate decreases to 1 Sv within 150 years from the start of two 56k members and within 550 and 710 years in two other members; only one member maintains a 2 Sv Nordic Seas overturning rate after 900 years. Similarly, in the 12.5k ensemble, 1 Sv Nordic Seas overturning is reached within 250 and 555 years in
<figure><img src="image_11.png" /></figure>

**Fig. 7** Upper ocean salinity changes in 56k ensemble members between the AMOC strong and weak mode. **a** February sea surface salinity in the western sector of the subtropical North Atlantic (10–30°N 50–80°W) in the 56k-stb (*dark grey*), 56k-red (*black*) and 56k-shd (*light grey*) ensemble members (*thick lines*, scale on left hand side). Sea surface salinities are expressed in psu. Shown here is the time evolution of SSS for years 1 through 1,100 of the respective simulations. The western sector of the subtropical North Atlantic was chosen to represent the surface waters transported to the North Atlantic and Nordic Seas convection sites by the wind driven sub-tropical and sub-polar gyres. Further shown in (a) is the time evolution of meridional salinity transport in the top 100 m of the

ocean through the mid-latitudes of the eastern North Atlantic in ensemble member 56k-stb. Here, meridional salinity transport (*thin line*, scale on right hand side) is represented by the magnitude of the meridional velocity (v) times the practical salinity (S) at 50°N and 20°W. The unit of this measure of transport is m s⁻¹ psu. Positive values of v signify northward water flow. **b, c, d** indicate the meridional salinity transport (colour scale) in the upper ocean through the North Atlantic in the strong (b) and weak (c) AMOC modes as well as the anomaly of the weak minus strong mode (d). In addition, the three figures show the vector of the total (meridional plus zonal (u)) water velocity in the upper 100 m to indicate strength (in m s⁻¹) and direction of flow

the two members with a marginally stable weak AMOC, but remains near 2 Sv in the other three members. The weakening of deep convection in the Nordic Seas, if occurring prior to year 900 in the simulation, contributes to the positive convective buoyancy feedback (Fig. 9).

By contrast, in the cold type, a slower rise in surface water buoyancy leads to a delayed advective AMOC shutdown (Fig. 6c). Two factors contribute to the slower buoyancy build-up (shown on the right panel of Fig. 9). Firstly, the contribution of the convective buoyancy feedback is smaller than in the mild glacial type. The positive feedback only comes from advective weakening of deep convection south of Greenland near the sea-ice edge (see Fig. 4b). Secondly, an anomalous positive northward SSS

advection along the Gulf Stream—North Atlantic Current system counteracts the advective buoyancy gain for as long as a weak AMOC is maintained. This explains the persistence of the weak AMOC mode (i.e. type A) in spite of continued freshwater forcing.

The above anomalous SSS advection is the effect of a strong coupling between AMOC and atmosphere. It originates in increased (sub-)tropical SSS compared to the reference state, because of a strong reduction in precipitation minus evaporation. This effect is illustrated on Fig. 7a, depicting the SSS evolution in three 56k ensemble members in the western sector of subtropical North Atlantic. As can be seen, SSS increases by about 0.6 psu until year 200. Once the weak AMOC is established, internal variability
governs the further SSS evolution until year 900 when Heinrich event-like fresh water forcing is applied. During this interval with weak AMOC, SSS increases in all three ensemble members by up to 0.2 psu in about 150 years just before the AMOC veers to shutdown (Fig. 7a). In the case of experiment 56k-shd this SSS increase occurs between year 200 and 500, in 56k-red between year 600 and 750 and in 56k-stb between year 900 and 1000. This suggests further drying (lower precipitation minus evaporation) over this subtropical Alantic Ocean sector. After the initial increase, SSS gradually decreases because during AMOC shutdown, absence of deep convection allows for more freshwater accumulation at the surface, which is recirculated to the sub-tropics within the sub-tropical gyre. As the SSS evolution of 56k-stb suggests (Fig. 7a), SSS in this region stays at 0.6 psu above the 56k-ref level when the AMOC remains in its weak mode. In 56k-stb this is disturbed by the start of Heinrich forcing at year 900.

It is this positive SSS anomaly that, through transport to the convection sites, prevents further freshening (buoyancy gain) in the North. The actual northward salinity transport in the upper ocean is shown in Fig. 7. The bottom time series of meridional water velocity (v) times salinity (S) in the top 100 m in Fig. 7a suggests that, as soon as deep convection in the Labrador Sea is shut down, salinity transport to the two remaining convection sites is enhanced and remains higher for as long as the AMOC remains in its weak mode. Fig. 7d depicts the change in v * S from the strong (Fig. 7b) to weak (Fig. 7c) AMOC modes in the 56k-stb simulation. Although the annual average absolute salt transport at the surface of the northeast Atlantic is southward in both the strong and weak AMOC modes (Fig. 7b, c), it is less intense in the weak mode, implying an anomalous northward salt transport component compared to the strong AMOC mode, directed towards the convection sites (Fig. 7d).

The discussed anomalous northward salt transport is consistent with the local surface conditions over the Nordic Seas convection site, shown in Fig. 8a, c, e. Especially in the 56k-red and 56k-stb ensemble members, the SSS remains nearly constant between year 200 and the start of the AMOC (indicated by the crossing of the grey bars), despite the continuous 0.08 Sv freshwater forcing added in this period to the Nordic Seas.

Progressive drying over the (sub-)tropical North Atlantic as the AMOC weakens is caused by a southward shift in tropical rain band. These results are consistent with data (e.g. Peterson and Haug 2006; Schmidt et al. 2006) and modelling studies, which associate a southward shift of the Intertropical Convergence Zone to an AMOC weakening or shutdown (e.g. Zhang and Delworth 2005; Menviel et al. 2008; Kageyama et al. 2009; Khodri et al. 2009). Ultimately, it is the effect of lower atmospheric temperatures in

the North Atlantic region as a weak AMOC transports about half the heat northwards (Menviel et al. 2008).

Finally, while the positive convective buoyancy feedback is active for both the mild and the cold glacial type, a slightly larger negative SSS anomaly in the North Atlantic region has built up before this time in the mild type. Consequently, the weak AMOC exports ~1 Sv less NADW in the mild type than in the cold type (Fig. 6b–d).

### 3.3 Density threshold and AMOC shutdown

Aside from maintaining a stable weak AMOC, the negative advective buoyancy feedback further manifests itself in both the cold and the mild glacial type by maintaining deep convection in the North Atlantic Ocean south of Iceland for as long as convection in the Nordic Seas is still active (Fig. 9). While the AMOC weakens to shutdown, the convective layer depth diminishes nearly synchronously at the latter convection sites, the decrease south of Iceland lagging ~50 years behind the decrease in the Nordic Seas (see Fig. 8b, d, f for the convective layer depth in the Nordic Seas). However, the different timings of the onset of further AMOC weakening to shutdown in individual members of each ensemble is not explained yet. Therefore, we analyse three ensemble members of 56k in more detail (see Fig. 6b). Each member was identically forced and only differs in initial atmospheric circulation. Despite the climatologically insignificant difference, we note a sudden divergence of 56k-shd from 56k-stb and 56k-red in terms of NADW export around year 275. At this time, the former diverges from a type B weak AMOC to a shutdown. Similarly, the type B weak AMOC in 56k-red diverges from 56k-stb—characterised by a type A weak AMOC—around year 600.

To explain the sudden divergence between ensemble members, the influence of internal variability on a non-linear system with a threshold needs to be invoked. When the February surface density threshold of −0.6 kg m⁻³ compared to the reference state is reached in the Nordic Seas, the weak AMOC veers to a shutdown. Below this value, surface water density can no longer surpass that of local deepest water masses. As a result, deep convection in the Nordic Seas weakens. At this point, the positive convective buoyancy feedback is enhanced and overwhelms the negative advective buoyancy feedback. Consequently, the AMOC will further weaken to shutdown. Since the weak AMOC is characterised by a Nordic Seas ρ surface anomaly of −0.4 kg m⁻³, this threshold can be crossed either by imposing a larger freshwater forcing—e.g. 0.25 Sv after year 900, or stochastically, in terms of interannual variability, by a relatively large negative SSS anomaly to overcome the remaining necessary −0.2 kg m⁻³ anomaly. Even though surface warming may
**Fig. 8** Relationship between February SST, SSS, surface water density ($\rho_{\text{surface}}$) and depth of convection in the Nordic Seas in the 56k-stb (top panels), 56k-red (middle panels) and 56k-shd (bottom panels) ensemble members. **a, c, e** show the time evolution of SST (dark grey), SSS (light grey) and $\rho_{\text{surface}}$ for years 1 through 1,100 of the respective simulations. **b, d, f** display the time evolution of the $\rho_{\text{surface}}$-anomaly to 56k-REF and of the convective layer depth. The horizontal and vertical grey bars cross where the $-0.6 \text{ kg m}^{-3} \rho_{\text{surface}}$ anomaly threshold is crossed beyond which the convective layer depth gradually decreases (and remains below 600 m at all times). The same threshold value is found for the Labrador Sea and North Atlantic convection sites

<figure><img src="image_12.png" /></figure>

also enhance buoyancy, a sufficiently strong warming is not found in the weak AMOC mode. Therefore, and as diminished deep convection seems to be slightly delayed in the North Atlantic convection site, interannual SSS variability in the Nordic Seas controls the density (buoyancy) threshold crossing (see Fig. 8),

Interestingly, we find that the same density threshold anomaly to shift deep convection to shutdown in the Nordic Seas is also seen south of Iceland. Figs. 3b, d, f show that, as the AMOC further weakens to shutdown from its weak mode, the positive convective buoyancy feedback overwhelms the negative advective feedback. This is reflected in the larger SSS-anomaly build-up in 56k-shd—in which the AMOC has already reached shutdown (Figs. 5f, 6b), smaller in 56k-red—in which NADW export

decreased from ~7 Sv to ~5 Sv between years 600 and 700 (Fig. 6b), and still much smaller in 56k-stb.

The non-linearity of buoyancy changes can be inferred from differences in overturning stream function between 56k-stb, 56k-red and 56k-shd around year 700 (Fig. 5b, d, f). The most northerly overturning cell—maintained by convection in the Nordic Seas in the weak AMOC—is only active in 56k-stb (Fig. 5b). The AMOC in 56k-stb is further characterised by a ~50% weaker and ~15% shallower main overturning cell than in 56k-REF. As such, it typifies the weak AMOC seen in all members. In contrast, 56k-red—in which the AMOC was essentially the same as 56k-stb until year 600—is characterised by a collapsed northern overturning cell and a further 30% weaker, 50% shallower main overturning cell (Fig. 5d).
**Fig. 9** The dynamic AMOC response to +0.15 Sv freshwater forcing at the deep convection sites (circles) for the mild (left) and cold glacial climate types (right) in the reference states = year 0; at year 200 and at year 900. Circle thickness reflects overturning strength. We numbered processes chronologically. A positive, convective (grey boxes) and a negative, advective buoyancy feedback (white boxes) affect speed and amplitude of AMOC weakening. The shutdown of a convection site occurs gradually (Nordic Seas and North Atlantic) or abruptly (Labrador Sea) once the surface density ($\rho_{\text{surface}}$) falls below $-0.6 \text{ kg m}^{-3}$ compared to 56k-REF. This density threshold is crossed by a local SSS anomaly related to internal variability of the ocean–atmosphere coupled system

<figure><img src="image_13.png" /></figure>

What probability for the timing of surface density crossing leading to AMOC collapse can we expect in the cold and mild glacial climate types? We calculated this probability based on our ensemble sets performed with 21k and 56k forcings, respectively. The expected timing of density threshold crossing is earlier in the mild type than in the cold type. Simulated threshold crossings vary from year 281 up to year 938 in the 56k ensemble, so only once after larger, Heinrich event-like freshwater forcing was started at year 900 (see e.g. Fig. 8). In the 21k ensemble on the other hand, the threshold crossing occurred before year 900 (in year 805) in only one member. With (1) an interannual SSS variability of ~0.067 psu (=1$\sigma$ = 1 standard deviation) in

the weak AMOC mode in all 56k and 21k simulations and (2) SSS controlling surface density, expected Nordic Seas density threshold crossing timings may be estimated from the difference in Nordic Seas SSS between the weak AMOC mode and the SSS at the time of simulated density. The respective values are 33.7 psu and 33.5 psu in the 56k ensemble members 56k-stb, 56k-red and 56k-shd, as highlighted by the grey bar crossings on Fig. 8a, c, e. This difference of ~-0.2 psu for the 56k ensemble—compared to ~-0.3 psu for the 21k ensemble (not shown)—is smaller in 56k simulation than in 21k simulations due to a larger decrease in AMOC strength from strong to weak mode (see Figs. 6, 9). This SSS difference means the
density threshold would be crossed if Nordic Seas SSS reaches an negative anomaly of ~3σ in 56k simulations and >4σ in 21k simulations. Assuming a normal distribution of interannual SSS variability in the Nordic Seas, this would mean an expected threshold crossing at around year 500 in 56k simulations and after several thousands of years in 21k simulations, not inconsistent with the simulated timings. This earlier timing in the 56k ensemble explains the higher probability of AMOC weakening to shutdown before year 900 in the mild type.

## 3.4 Discussion

### 3.4.1 Implications

Previously, it was shown that the sensitivity of the AMOC to freshwater perturbation differs between present-day or pre-industrial climate states and LGM states in earth system models of various complexity (e.g. Ganopolski et al. 1998; Schmittner et al. 2002; Rahmstorf et al. 2005; Roche et al. 2007; Weber et al. 2007; Hu et al. 2008). In these models, the sensitivity is generally larger in LGM states, which is generally thought to be the result of much colder climate conditions and the closure of the Bering Strait (see Hu et al. 2008). In the LOVECLIM model, a larger sensitivity was previously found in an LGM state than in a warm interglacial climate state (Roche et al. 2007). In conjunction with these studies, our results imply that a range of sensitivities under glacial and interglacial climate states may exist. We show that, at least in our model, the AMOC could be the most sensitive to freshwater perturbation in a mild glacial climate type.

If the hypothesis of Broecker et al. (1985) on the cause of rapid shifts between stadials and interstadials during the Last Glacial holds true, and if our model simulations represent the physics of the observed climatic changes well, we find an exciting new hypothesis that could help explain the more frequent recurrence of these shifts during mild times, for instance MIS 3 and during Termination 1, than during cold ones, for instance around the LGM. As sketched in Fig. 9, in mild glacial climates, the northerly position of convection sites—the Labrador Sea and Nordic Seas in our model—make the AMOC more sensitive to freshwater fluxes. Conversely, when climate is colder, i.e. at LGM, the AMOC sensitivity to freshwater forcing is reduced as convection sites are shifted southward, i.e. to the North Atlantic Ocean in our model. Similar results, with a southward shift of convection in LGM versus MIS 3 glacial climates have been found in other models (e.g. Merkel et al. 2010), further supporting the hypothesis.

The AMOC resumption may be governed by a similar, but reversed process sequence than the one shown in Fig. 9, again involving advective and convective buoyancy

feedbacks (Renold et al. 2010). The implication is that the latitudinal position of the sea-ice edge and the intimately linked location of the convection sites may control the frequency and amplitude of abrupt shifts in glacial climates, being more frequent and larger in a mild than a cold glacial climate. The new hypothesis might explain the findings of McManus et al. (1999) that abrupt, millennial-scale climate change in the North Atlantic region was most clearly expressed at times of intermediate Northern Hemisphere ice sheets during the last 500 ka and (nearly) absent from glacial maxima and interglacials.

### 3.4.2 Remaining uncertainties

Several caveats may limit the validity of our hypothesis that a northward shift in deep convection sites in response to changing glacial climate forcings between LGM and MIS 3 or T1 may explain more frequent abrupt climate shifts between stadials and interstadials in the latter two. We summarise six important caveats now.

Firstly, whereas absence of deep convection in an ice-covered Labrador Sea in our LGM background climate is consistent with data (e.g. Paul and Schäfer-Neth 2003; de Vernal et al. 2005), deep convection and absence of sea-ice in our Younger Dryas reference climate is not (Sarnthein et al. 1994; Hillaire-Marcel and Bilodeau 2000). Furthermore, to our knowledge, no data on LSW formation during MIS 3 is found in literature, except for a time slice at ~30 ka BP when no Labrador Sea Water formation took place as suggested by δ¹³C profiles of the North Atlantic Ocean (Sarnthein et al. 1994). However, our YD state was not realistically forced, since we did not apply a baseline meltwater flux from the melting ice sheets during the deglaciation. We argue that doing so would effectively shutdown deep convection in the Labrador Sea in our model, as shown by Renssen et al. (2009) for the early-Holocene. Similarly, Cottet-Puinel et al. (2004) found that the absence of LSW until 7 ka BP in their model could only be obtained with a baseline meltwater flux from the Laurentide Ice Sheet, highlighting the strong local sensitivity of deep convection to freshwater forcing. Regarding MIS 3, we have no reason to believe that LSW formation was prevented throughout the entire interval, since no major deglaciation occurred then (Peltier 2004). That is unless an ice shelf permanently and extensively covered the Labrador Sea across stadials and interstadials. On the contrary, Hulbe et al. (2004) suggest that an ice-shelf collapse may have triggered part of the massive ice-berg release to the North Atlantic Ocean associated with Heinrich events. This implies that at several occasions during MIS 3—at least during stadials coinciding with a Heinrich event—an ice shelf may have built up and largely disappeared after the Heinrich event. As Heinrich events during
MIS 3 were followed by particularly long and warm interstadials (Bond et al. 1993), it is conceivable that LSW may have occurred at least in the early phase of those interstadials. LSW formation may, in fact, have slightly amplified the warming into those interstadials, as suggested by our model where a shutdown in LSW formation decreases annual mean surface air temperatures in Central Greenland by 0.6°C—following the warming signal over the Labrador Sea. In any case, there are indications that a southward shift in deep convection sites did take place between interstadials of MIS 3 or T1 and the LGM (Rasmussen et al. 1997). This suggests that on/off switches in Nordic Seas convection may also contribute to the large amplitude and, though in our model less so than Labrador Sea convection switches, abruptness of glacial climate shifts.

Secondly, our hypothesis stands only if the cold temperatures of MIS 3 stadials are due to perturbations of a generally milder background climate, more closely resembling interstadials, as appears to be the case in our model for MIS 3 (Van Meerbeeck et al. 2009). This is in agreement with climate reconstructions suggesting an overall mild MIS 3 climate in Europe (e.g. van der Hammen et al. 1967; van Huissteden and Vandenberghe 1988), although those mild palaeotemperature records are no proof of stability of mild conditons. Yet, we currently do not know whether stadials or interstadials or neither of these climate states were stable under MIS 3 climate forcings, with climate models offering different view points. For instance, Ganopolski and Rahmstorf (2001) find a stadial-like climate to be in equilibrium with glacial forcings. In their simulations, a shallower and slightly weaker AMOC is found in their glacial reference climate—generally consistent with data for LGM, along with a mid-latitudinal location of deep convection. Then, with a small negative freshwater forcing to the mid-latitudes of the North Atlantic Ocean, deep convection was shifted to Polar regions, thereby abruptly shifting climate to interstadial-like conditions (Ganopolski and Rahmstorf 2001). Along with the palaeoceanographic work of Sarnthein et al. (1994) and modelling work of e.g. Manabe and Stouffer (1988) and many others (e.g. Schmittner et al. 2002), this has led us to consider that that the glacial AMOC may be bi-stable between a stadial and an interstadial mode, with the possibility of a third, 'Heinrich' shut down mode during the longest stadials (Rahmstorf 2002). Furthermore, Rahmstorf et al. (2005) showed that the bi-stability of the AMOC under a certain freshwater forcing is model-dependent in earth system models of intermediate complexity—a hierarchy of models to which both Ganopolski and Rahmstorf (2001)'s and our model belong. An additional complication is that the glacial reference state of Ganopolski and Rahmstorf (2001)—and many other

modelling studies (e.g. Weber et al. 2007)—was setup with climate forcings similar to our LGM reference state, and might therefore not be representative of MIS 3 background climate—for which we find a stable and strong AMOC mode in our model. Finally, AMOC bi-stability seems to be strongly reduced or inexistent in coupled atmosphere-ocean general circulation models (Stouffer et al. 2006). This casts some further doubt on the possible stability of a weak AMOC mode during stadials and a strong AMOC mode during interstadials—which may both merely be transient features of glacial climate. Past studies have, therefore, thus far not been able to provide an unequivocal solution to this problem.

Thirdly, it is presently not known to what extent the suggested AMOC strength and mode shifts between stadials and interstadials were wind-driven or resulted from a perturbation of the thermohaline forced component of the AMOC (Wunsch 2006). Put into the context of this paper, even if AMOC shifts due to varying meltwater fluxes were facilitated in a mild compared to a cold glacial climate, this control may have been subordinate to changes in atmospheric circulation. This doubt is justified, since Kuhlbrodt et al. (2007) identified the main drivers of the AMOC—in terms of energy transfer—to be mid- and high-latitude winds in the Southern Hemisphere as well as vertical and horizontal diffusion in the water column.

Fourthly, even if thermohaline-controlled, the AMOC shifts associated with stadial-interstadial transitions may not have been triggered by variable glacial meltwater flow from Northern Hemisphere ice sheets (e.g. Knorr and Lohmann 2003). If a relatively strong meltwater flux were required to substantially reduce overturning strength during stadials—as seen in our and many other models (Rahmstorf et al. 2005; Stouffer et al. 2006), a sufficiently large meltwater source on the right time scale has yet to be found for MIS 3 (Clement and Peterson 2008). Indications are that stadial meltwater from the Eurasian Ice Sheets may have contributed to this flux (Rasmussen et al. 1997; Lekens et al. 2006), however no quantification has been undertaken. Also, whereas in models large freshwater forcing—said to represent meltwater release during the Younger Dryas cold event (~12,700–11,500 ka BP) or during Heinrich events—triggers an AMOC shutdown, a growing body of evidence suggest AMOC reduction or shutdown took place prior to Heinrich events (Hemming 2004). Moreover, the iceberg armadas associated with the massive ice rafted debris deposition during Heinrich events may have been triggered by an AMOC shutdown (e.g. Hulbe et al. 2004; Flückiger et al. 2006), with its meltwater release preventing an AMOC resumption. However, so far, error margins in chronologies have precluded unequivocal conclusions on this matter (Hemming 2004). In any case, with an AMOC shutdown achieved in several of our
simulations before the stronger freshwater forcing that represents the meltwater release by a Heinrich event, our results are at least partly consistent with Hulbe et al. (2004).

Fifthly, our understanding of the stochastic process leading to AMOC shutdown is presently poor (e.g. Ditlevsen et al. 2007). Three possibilities are under debate: (1) AMOC shifts were inherent features of the glacial ocean circulation through internal variability (Sakai and Peltier 1997; Ditlevsen et al. 2007) or (2) they were externally forced (Alley et al. 2001; Ganopolski and Rahmstorf 2001), or (3) both (Timmermann et al. 2003). In a simple coupled climate model, the latter authors found that an initial meltwater pulse such as a Heinrich event freshens the ocean could bring the AMOC to a bifurcation point (Hopf bifurcation), thus exciting the AMOC into a mode of Dansgaard-Oeschger type variability. Small noise arising from internal variability may then sustain this variability. While our results preclude neither explanation, they indicate that, at least in our model, the internal variability needed to cross a surface density threshold in the convection sites may involve the atmosphere, not just the ocean. Specifically, different initial atmospheric conditions between different ensemble members produce different weather patterns in an otherwise dynamically identical climate state. Further investigation is required to answer the questions whether and which weather patterns bring about SSS anomalies large enough to cross the density threshold. Although an in-depth analysis of the atmosphere-ocean coupling lies beyond the scope of this paper, we can provide a way to clarify this problem. Using LOVECLIM, the presented simulation could be repeated, only with a fixed annual cycle for the atmospheric conditions. If a weak AMOC is maintained in the ocean only experiments, a proof is given that simulated weather patterns are responsible for the timing of the threshold crossing in our model, and thus further AMOC weakening to shutdown.

Sixthly, although Fig. 6 gives clear indications of stronger AMOC sensitivity to a relatively large positive freshwater flux (0.15 Sv during years 1–900) in the mild glacial type than the cold type, five ensemble members for each ensemble are not sufficient to statistically underpin this finding. Therefore, with a large standard deviation of overturning strength response time for each ensemble, a thorough assessment of sensitivity differences would require several tens of members per climate type. Such an investigation remains to be done. Nevertheless, the process sequences and mechanisms discussed in this paper are consistent with our results and, to our knowledge, with data. We thus assume that 5 ensemble members are sufficient for a first assessment of AMOC sensitivity to additional freshwater fluxes in different glacial climate states.

## 4 Conclusions

By perturbing three different glacial reference climate states in the LOVECLIM model with one freshwater scenario, we find different Atlantic Meridional Overturning Circulation (AMOC) responses to this forcing.

We distinguish two glacial climate types based on sea surface conditions in the North Atlantic region and on AMOC mode in the reference states. (1) The cold glacial type consists of our Last Glacial Maximum (21 ka BP) simulations, with deep convection mainly found south of Greenland and Iceland. (2) The mild glacial climate type includes our early-Marine Isotope Stage 3 (56 ka BP) and our Late-Glacial (12.5 ka BP) simulations. A milder sea surface allows deep convection further north, namely in the Labrador Sea and the Nordic Seas.

Labrador Sea convection in the reference climates of the mild type makes the AMOC more sensitive to an imposed additional 0.15 Sv freshwater flux in the North Atlantic and Nordic Seas than the cold type. The mild type attains (1) a faster and larger reduction in AMOC strength on average; and (2) an abrupt shutdown of Labrador Sea convection within less than 50 years time from the start of the freshwater forcing.

After shutdown of Labrador Sea convection, a positive convective buoyancy feedback accelerates accumulation of the freshwater anomaly—thus lower density—at surface. The buoyancy build-up slows down deep convection at the other convection sites, leading to further AMOC weakening to about half its initial strength within 200 years. However, at this time a negative advective buoyancy feedback—northward transport of anomalously saline subtropical Atlantic surface waters along the AMOC—prevents further AMOC weakening in all glacial states.

Being close to a density threshold, the weak AMOC mode can be brought to full shutdown stochastically by a relatively strong SSS anomaly in the Nordic Seas. In the mild type, such a shutdown is often—but not always—achieved within 600 to 900 years. That is, before applying a larger freshwater forcing representing the meltwater released during Heinrich events. In contrast, slightly higher local SSS anomalies and thus densities (to the reference states) in the cold type reduces the probability of achieving a full shutdown, since the positive convective buoyancy feedback is weaker.

In conclusion, our results imply that MIS 3 (and possibly Late-Glacial) climate forcings may have made background climate more prone to abrupt climate shifts. Local-scale interannual variability of the North Atlantic Ocean surface circulation may then have controlled the recurrence of cold stadial intervals.
**Acknowledgments** This work is a contribution to the RESOLU-
TION-project (ESF EUROCORES on EuroCLIMATE program).
C.J.V.M. was sponsored by the Netherlands Organisation for Scien-
tific Research (NWO), under contract number 855.01.085. D.M.R.
was sponsored by the NWO under contract number 854.00.024 and
supported by INSU/CNRS.

**Open Access** This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.

References

Adkins J, McIntyre K, Shrag D (2002) The salinity, temperature, and δ¹⁸O of the Glacial Deep Ocean. Science 298:1769-1773

Alley RB, Anandakrishnan S, Jung P (2001) Stochastic resonance in the North Atlantic. Paleooceanography 16:190-198

Bard E, Rostek F, Turon J-L, Gendreau S (2000) Hydrological impact of Heinrich events in the subtropical Northeast Atlantic. Science 289:1321-1324

Berger AL (1978) Long-term variations of daily insolation and quaternary climatic changes. J Atmos Sci 35:2363-2367

Berger A, Loutre M (1992) Astronomical solutions for paleoclimate studies over the last 3 millions years. Earth Planet Sci Lett 111:369-382

Bond G, Broecker W, Johnsen S, McManus J, Labeyrie L, Jouzel J, Bonani G (1993) Correlations between climate records from North Atlantic sediments and Greenland ice. Nature 365:143-147

Broecker WS, Peteet DM, Rind D (1985) Does the ocean-atmosphere system have more than one stable mode of operation? Nature 315:21-26. doi:10.1038/315021a0

Brovkin V, Ganopolski A, Svirezhev Y (1997) A continuous climate-vegetation classification for use in climate-biosphere studies. Ecol Model 101:251-261

Clement AC, Peterson LC (2008) Mechanisms of abrupt climate change of the last glacial period. Rev Geophys 46:RG4002

Dansgaard W, Johnsen SJ, Clausen HB, Dahl-Jensen D, Gundestrup NS, Hammer CU, Hvidberg CS, Steffensen JP, Sveinbjörnsdottir AE, Jouzel J, Bond G (1993) Evidence for general instability of past climate from a 250-kyr ice-core record. Nature 364:218-220

de Vernal A, Eynaud F, Henry M, Hillaire-Marcel C, Londeix L, Mangin S, Matthiessen J, Marret F, Radi T, Rochon A, Solignac S, Turon J-L (2005) Reconstruction of sea-surface conditions at middle to high latitudes of the Northern Hemisphere during the last glacial maximum (LGM) based on dinoflagellate cyst assemblages. Quat Sci Rev. doi:10.1016/j.quascirev.2004.06.04

Dickson RR, Brown J (1994) The production of North Atlantic Deep Water: sources, rates, and pathways. J Geophys Res 99(C6):12319-12341

Ditlevsen PD, Andersen KK, Svensson A (2007) The DO-climate events are probably noise induced: statistical investigation of the claimed 1470 years cycle. Clim Past 3: 129-143. URL: www.clim-past.net/3/129/2007/

Driesschaert E (2005) Climate change over the next millennia using LOVECLIM, a new Earth system model including the polar ice sheets. Dissertation, Université Catholique de Louvain

Driesschaert E, Fichefet T, Goosse H, Huybrechts P, Janssens I, Mouchet A, Munhoven G, Brovkin V, Weber SL (2007) Modeling the influence of Greenland ice sheet melting on the Atlantic meridional overturning circulation during the next millennia. Geophys Res Lett 34:L10707. doi:10.1029/2007GL029516

Flückiger J, Dällenbach A, Blunier T, Stauffer B, Stocker T, Raynaud D, Barnola J-M (1999) Variations in atmospheric N₂O concentration during abrupt climatic changes. Science 285:227-230. doi:10.1126/science.285.5425.227

Flückiger J, Blunier T, Stauffer B, Chappellaz J, Spahni R, Kawamura K, Schwander J, Stocker TF, Dahl-Jensen D (2004) N₂O and CH₄ variations during the last glacial epoch: Insight into global processes. Global Biogeochem Cycl 18:GB1020. doi:10.1029/2003GB002122

Flückiger J, Knutti R, White WC (2006) Oceanic processes as potential trigger and amplifying mechanisms for Heinrich events. Paleooceanography 21:PA2014. doi:10.1029/2005PA001204

Ganopolski A, Rahmstorf S (2001) Rapid changes of glacial climate simulated in a coupled model. Nature 409:153-158

Ganopolski A, Rahmstorf S (2002) Abrupt glacial climate changes due to stochastic resonance. Phys Rev Lett. doi: 10.1103/Phys RevLett.88.038501

Ganopolski A, Rahmstorf S, Petoukhov V, Claussen M (1998) Simulation of modern and glacial climates with a coupled global model of intermediate complexity. Nature 391:351-356

Goosse H, Fichefet T (1999) Importance of ice-ocean interactions for the global ocean circulation: a model study. J Geophys Res 104(C10):23337-23355. doi:10.1029/1999JC900215

Grootes PM, Stuiver M (1997) Oxygen 18/16 variability in Greenland snow and ice with 10-3 to 105-year time resolution. J Geophys Res 102(C12):26455-26470

Hemming SR (2004) Heinrich events: massive late Pleistocene detritus layers of the North Atlantic and their global climate imprint. Rev Geophys 42:RG1005. doi:10.1029/2003RG000128

Hillaire-Marcel C, Bilodeau G (2000) Instabilities in the Labrador Sea water mass structure during the last climatic cycle. Can J Earth Sci 37:795-809

Hillaire-Marcel C, de Vernal A, Bilodeau G, Weaver AJ (2001) Absence of deep-water formation in the Labrador Sea during the last interglacial period. Nature 410:1073-1077

Hu A, Otto-Bliesner BL, Meehl GA, Han W, Morrill C, Brady EC, Briegleb B (2008) Response of thermohaline circulation to freshwater forcing under present-day and LGM conditions. J Clim 21:2239-2258. doi:10.1175/2007JCLI11985.1

Huber C, Leuenberger M, Spahni R, Flückiger J, Schwander J, Stocker TF, Johnsen SJ, Landais A, Jouzel J (2006) Isotope calibrated Greenland temperature record over Marine Isotope Stage 3 and its relation to CH4. Earth Planet Sci Lett 243:504-519

Hulbe CL, MacAyeal DR, Denton GH, Kleman J, Lowell TV (2004) Catastrophic ice shelf breakup as the source of Heinrich event icebergs. Paleooceanography 19:PA1004. doi:10.1029/2003PA000890

Indermühle A, Monnin E, Stauffer B, Stocker TF, Wahlen M (2000) Atmospheric CO₂ concentration from 60 to 20 kyr BP from the Taylor Dome ice core, Antarctica. Geophys Res Lett 27(5): 735-738

Johnsen SJ, Clausen HB, Dansgaard W, Fuhrer K, Gundestrup N, Hammer CU, Iversen P, Jouzel J, Stauffer B, Steffensen JP (1992) Irregular glacial interstadials recorded in a new Greenland ice core. Nature 359(6393):311-313. doi:10.1038/359311a0

Jongma JI, Prange M, Renssen H, Schulz M (2007) Amplification of Holocene multicentennial climate forcing by mode transitions in North Atlantic overturning circulation. Geophys Res Lett 34:L15706. doi:10.1029/2007GL030642

Kageyama M, Mignot J, Swingedouw D, Marzin C, Alkama R, Marti O (2009) Glacial climate sensitivity to different states of the Atlantic Meridional Overturning Circulation: results from the IPSL model. Clim Past 5:551-570

Kageyama M, Paul A, Roche DM, Van Meerbeeck CJ (2010) Modelling millennial-scale variability during glacial epochs: a
review. Quaternary Science Reviews 29(21–22):2931–2956. doi: 10.1016/j.quascirev.2010.05.029

Khedri M, Kageyama M, Roche DM (2009) Sensitivity of South American tropical climate to last glacial maximum boundary conditions: focus on teleconnections with tropics and extratrop- ics. In: Vimeux F, Sylvestre S, Khodri M (eds) Past climate variability in South America and surrounding regions—from the last glacial maximum to the holocene, vol 14. Springer, Berlin, pp 213–238

Knorr G, Lohmann G (2003) Southern Ocean origin for the resumption of Atlantic thermohaline circulation during deglaci- ation. Nature 424:532–536. doi:10.1038/nature01855

Knutti R, Flückiger J, Stocker TF, Timmermann A (2004) Strong hemispheric coupling of glacial climate through freshwater discharge and ocean circulation. Nature 430:851–856

Kucera M, Weinelt M, Kiefer T, Pflaumann U, Hayes A, Weinelt M, Chen M-T, Mix A, Barrows T, Cortijo E, Duprat J, Juggins S, Waelbroeck C (2005) Reconstruction of sea-surface tempera- tures from assemblages of planktonic foraminifera: multi- technique approach based on geographically constrained cali- bration data sets and its application to glacial Atlantic and Pacific Oceans. Quat Sci Rev 24:951–998. doi:10.1016/j.quascirev. 2004.07.014

Kuhlbrodt T, Griesel A, Montoya M, Levermann A, Hofmann M, Rahmstorf S (2007) On the driving processes of the Atlantic Meridional Overturning Circulation. Rev Geophys 45:RG2001. doi:10.1029/2004RG000166

Lekens WAH, Sejrup HP, Haflidason H, Knies J, Richter T (2006) Meltwater and ice rafting in the southern Norwegian Sea between 20 and 40 calendar kyr BP: implications for Fenno- scandian Heinrich events. Paleooceanography 21:PA3013. doi: 10.1029/2005PA001228

Manabe S, Stouffer RJ (1988) Two stable equilibria of a coupled ocean-atmosphere model. J Clim 1:841–866

MARGO Project Members (2009) Constraints on the magnitude and patterns of ocean cooling at the Last Glacial Maximum. Nat Geosci 2:127–132. doi:10.1038/NGEO411

McManus JF, Oppo DW, Cullen JL (1999) A 0.5-million-year record of millennial-scale climate variability in the North Atlantic. Science 283(5404):971. doi:10.1126/science.283.5404.971

McManus JF, Francois R, Gherardi J-M, Keigwin LD, Brown-Leger S (2004) Collapse and rapid resumption of Atlantic meridional circulation linked to deglacial climate changes. Nature 428:834–837. doi:10.1038/nature02494

Meland MY, Jansen E, Elderfield H (2005) Constraints on SST estimates for the northern North Atlantic/Nordic Seas during the LGM. Quat Sci Rev 24:835–852. doi:10.1016/j.quascirev.2004.05.011

Menviel L, Timmermann A, Mouchet A, Timm O (2008) Meridional reorganizations of marine and terrestrial productivity during Heinrich events. Paleooceanography 23:PA1203. doi:10.1029/ 2007PA001445

Merkel U, Prange M, Schulz M (2010) ENSO variability and teleconnections during glacial climates. Quat Sci Rev 29(1–2):86–100. doi:10.1016/j.quascirev.2009.11.006

Mix A (2003) Chilled out in the ice-age Atlantic. Nature 425:32–33

Monnin E, Indermuehle A, Daellenbach A, Flueckiger J, Stauffer B, Stocker T, Raynaud D, Barnola J-M (2001) Atmospheric CO₂ concentrations over the last glacial termination. Science 291:112–114

NorthGRIP-Members (2004) High-resolution record of Northern Hemisphere climate extending into the last interglacial period. Nature 431:147–151

Opsteegh J, Haarsma R, Selten F, Kattenberg A (1998) ECBILT: a dynamic alternative to mixed boundary conditions in ocean models. Tellus 50:348–367

Paul A, Schäfer-Neth C (2003) Modeling the water masses of the Atlantic Ocean at the Last Glacial Maximum. Paleooceanography 18:1058. doi:10.1029/2002PA000783

Peltier W (2004) Global glacial isostasy and the surface of the ice-age earth: the ICE-5G (VM2) Model and GRACE. Annu Rev Earth Planet Sci 32:111–149. doi:10.1146/annurev.earth.32.082503. 144359

Peterson LC, Haug GH (2006) Variability in the mean latitude of the Atlantic Intertropical Convergence Zone as recorded by riverine input of sediments to the Cariaco Basin (Venezuela). Palaeoge- ogr Palaeocl 234(1):97–113. doi:10.1016/j.palaeo.2005.10.021

Puinel-Cottet M, Weaver AJ, Hillaire-Marcel C, de Vernal A, Clark PU, Eby M (2004) Variation of Labrador Sea water formation over the Last Glacial cycle in a climate model of intermediate complexity. Quat Sci Rev 23:449–465. doi:10.1016/S0277-3791 (03)00123-9

Rahmstorf S (2002) Ocean circulation and climate during the past 120,000 years. Nature 419:207–214

Rahmstorf S, Crucifix M, Ganopolski A, Goosse H, Kamenkovich I, Knutti R, Lohmann G, Marsh R, Mysak LA, Wang Z, Weaver AJ (2005) Thermohaline circulation hysteresis: a model inter- comparison. Geophys Res Lett 32:L23605. doi:10.1029/ 2005GL023655

Rasmussen TL, van Weering TCE, Labeyrie L (1997) Climatic instability, ice sheets and ocean dynamics at high northern latitudes during the last glacial period (58–10 ka BP). Quat Sci Rev 16:71–80

Renold M, Raible CC, Yoshimori M, Stocker TF (2010) Simulated response of the North Atlantic meridional overturning circula- tion—slow basin-wide advection and abrupt local convection. Quat Sci Rev 29(1–2):101–112. doi:10.1016/j.quascirev. 2009.11.005

Renssen H, Goosse H, Fichefet T (2005) Contrasting trends in North Atlantic deep-water formation in the Labrador Sea and Nordic Seas during the Holocene. Geophys Res Lett 32:08711. doi: 10.1029/2005GL022462

Roche DM, Dokken TM, Goosse H, Renssen H, Weber SL (2007) Climate of the Last Glacial Maximum: sensitivity studies and model-data comparison with the LOVECLIM coupled model. Clim Past 3:205–224. http://www.clim-past.net/3/205/2007/

Roche DM, Wiersma AP, Renssen H (2010) A systematic study of the impact of freshwater pulses with respect to different geograph- ical locations. Clim Dyn 34:997–1013. doi:10.1007/s00382- 009-0578-8

Rooth CGH (1982) Hydrology and ocean circulation. Prog Oceanogr 11:131–149. doi:10.1016/0079-6611(82)90006-4

Sakai K, Peltier WR (1997) Dansgaard-Oeschger oscillations in a coupled atmosphere-ocean climate model. J Clim 10:949–970

Sánchez Goñi MF, Cacho I, Turon J-L, Guiot J, Sierro FJ, Peypouquet J-P, Grimalt JO, Shackleton NJ (2002) Synchronicity between marine and terrestrial responses to millennial scale climatic variability during the last glacial period in the Mediterranean region. Clim Dyn 19:95–105. doi:10.1007/s00382-001-0212-x

Sarnthein M, Winn K, Jung SJJ, Duplessy J-C, Labeyrie L, Erlenkeuser H, Ganssen G (1994) Changes in east Atlantic deepwater circulation over the last 30, 000 years: Eight time slice reconstructions. Paleooceanography 9(2):209–267. doi: 10.1029/93PA03301

Schmidt MW, Vautravers MJ, Spero HJ (2006) Rapid subtropical North Atlantic salinity oscillations across Dansgaard-Oeschger cycles. Nature 443:561–564. doi:10.1038/nature05121

Schmittner A, Clement AC (2002) Sensitivity of the thermohaline circulation to tropical and high latitude freshwater forcing during the last glacial-interglacial cycle. Paleooceanography 17(2):1017. doi:10.1029/2000PA000591
Schmittner A, Yoshimori M, Weaver AJ (2002) Instability of glacial climate in a model of the ocean-atmosphere-cryosphere system. Science 295:1489–1493

Schulz M (2002) On the 1470-year pacing of Dansgaard-Oeschger warm events. Paleooceanography 17(2):1014. doi:10.1029/2000 PA000571

Schulz M, Prange M, Klocke A (2007) Low-frequency oscillations of the Atlantic meridional overturning circulation in a coupled climate model. Clim Past 3:97–107. http://www.clim-past. net/3/97/2007/

Steffensen JP, Andersen KK, Bigler M, Clausen HB, Dahl-Jensen D, Fischer H, Goto-Azuma K, Hansson M, Johnsen SJ, Jouzel J, Masson-Delmotte V, Popp T, Rasmussen SO, Röthlisberger R, Ruth U, Stauffer B, Siggaard-Andersen M-L, Sveinbjörnsdóttir AE, Svensson A, White JWC (2008) High-resolution Greenland Ice Core Data show abrupt climate change happens in few years. Science 321:680–684. doi:10.1126/science.1157707

Stommel H (1961) Thermohaline convection with two stable regimes of flow. Tellus 13:224–230

Stouffer RJ, Yin J, Gregory JM, Dixon KW, Spelman MJ, Hurlin W, Weaver AJ, Eby M, Flato GM, Hasumi H, Hu A, Jungclaus JH, Kamenkovich IV, Levermann A, Montoya M, Murakami S, Nawrath S, Oka A, Peltier WR, Robitaille DY, Sokolov A, Vettoretti G, Weber SL (2006) Investigating the causes of the response of the thermohaline circulation to past and future climate changes. J Clim 19:1365–1387

Tartinville B, Campin J, Fichefet T, Goosse H (2001) Realistic representation of the surface freshwater flux in an ice-ocean general circulation model. Ocean Model 3:95–108

Timmermann A, Gildor H, Schulz M, Tziperman E (2003) Coherent resonant millennial-scale climate oscillations triggered by massive meltwater pulses. J Clim 16:2569–2585

van der Hammen T, Maarleveld GC, Vogel JC, Zagwijn W (1967) Stratigraphy, climatic succession and radiocarbon dating of the last glacial in the Netherlands. Geol Mijnb 46:79–95

van Huissteden J, Vandenberghe J (1988) Changing fluvial style of periglacial lowland rivers during the Weichselian Pleniglacial in the eastern Netherlands. Z Geomorph Suppl 71:131–146

van Kreveld S, Sarnthein M, Erlenkeuser H, Grootes P, Jung S, Nadeau MJ, Pflaumann U, Völker A (2000) Potential links between surging ice sheets, circulation changes, and the Dansgaard-Oeschger cycles in the Irminger Sea, 60–18 kyr. Paleooceanography 14(4):425–442

Van Meerbeeck CJ (2010) Modelling glacial abrupt climate change—on the expression of millennial-scale events over Europe and the North Atlantic region during Marine Isotope Stage 3. Dissertation, VU University Amsterdam. ISBN: 978-90-664-396-5

Van Meerbeeck CJ, Renssen H, Roche DM (2009) How did Marine Isotope Stage 3 and Last Glacial Maximum climates differ? Perspectives from equilibrium simulations. Clim Past 5:33–51. http://www.clim-past.net/5/33/2009/

Voelker AHL (2002) Global distribution of centennial-scale records for Marine Isotope Stage (MIS) 3: a database. Quat Sci Rev 21(10):1185–1212. doi:10.1016/S0277-3791(01)00139-1

Weaver AJ, Marotzke J, Cummins PF, Sarachik ES (1993) Stability and variability of the thermohaline circulation. J Phys Ocean 23:39–60

Weber S, Drijfhout S, Abe-Ouchi A, Crucifix M, Eby M, Ganopolski A, Murakami S, Otto-Bliesner B, Peltier W (2007) The modern and glacial overturning circulation in the Atlantic ocean in PMIP coupled model simulations. Clim Past 2:923–949. doi:1814- 9359/cpd/2006-2-923

Wohlfarth B, Veres D, Ampel L, Lacourse T, Blaauw M, Preusser F, Andrieu-Ponel V, Kéravis D, Lallier-Vergès E, Björck S, Davies SM, de Beaulieu J-L, Risberg J, Hormes A, Kasper HU, Possnert G, Reille M, Thouveny N, Zander A (2008) Rapid ecosystem response to abrupt climate changes during the last glacial period in western Europe, 40–16 ka. Geology 36(5):407–410. doi: 10.1130/G24600A.1

Wunsch C (2000) On sharp spectral lines in the climate record and the millennial peak. Paleooceanography 15:417–424. doi:10.1029/ 1999PA000468

Wunsch C (2006) Abrupt climate change. An alternative view. Quat Res 65(2):191–203

Zhang R, Delworth TL (2005) Simulated tropical response to a substantial weakening of the Atlantic thermohaline circulation. J Clim 18:1853–1860