# Simulated millennial-scale climate variability driven by a convection-advection oscillator

Yvan M Romé, Ruza F Ivanovic, Lauren J Gregoire, Didier Swingedouw, Sam Sherriff-Tadano, Reyk Börner

► To cite this version:

Yvan M Romé, Ruza F Ivanovic, Lauren J Gregoire, Didier Swingedouw, Sam Sherriff-Tadano, et al. Simulated millennial-scale climate variability driven by a convection-advection oscillator. Climate Dynamics, 2025, 63 (3), pp.150. 10.1007/s00382-025-07630-x. hal-05372694

HAL Id: hal-05372694

https://hal.science/hal-05372694v1

Submitted on 19 Nov 2025

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
# Simulated millennial-scale climate variability driven by a convection–advection oscillator

Yvan M. Romé<fnref n="1" /> Ⓒ · Ruza F. Ivanovic<fnref n="1" /> · Lauren J. Gregoire<fnref n="1" /> · Didier Swingedouw<fnref n="2" /> · Sam Sherriff-Tadano<fnref n="3" /> · Reyk Börner<fnref n="4" />

Received: 10 September 2024 / Accepted: 7 February 2025 / Published online: 7 March 2025
© The Author(s) 2025

## Abstract

The last glacial period, between around 115 and 12 thousand years before present, exhibited strong millennial-scale climate variability. This includes abrupt transitions between cold and warm climates, known as Dansgaard–Oeschger (D–O) cycles. D–O cycles have been linked to switches in dynamical regimes of the Atlantic Overturning Meridional Circulation (AMOC), but the exact mechanisms behind abrupt climate changes and AMOC regime shifts remain poorly understood. This paper introduces the convection–advection oscillator mechanism to explain the millennial-scale oscillations observed in a set of HadCM3 general circulation model simulations forced with snapshots of deglacial meltwater history. The oscillator can be separated into two components acting on different time scales. The fast convection component responds to changes in vertical stratification in the North Atlantic by activating or deactivating deep water formation sites. The slow advection component regulates the accumulation and depletion of salinity in the North Atlantic. This oscillator mechanism is triggered under specific background conditions and freshwater release patterns. The freshwater perturbation causes an instability that triggers a global salt reorganisation, modifying the North Atlantic stratification. For a given forcing pattern, the system oscillates if the salt transport can lead to an alternating reactivation and deactivation of the AMOC. Otherwise, the climate settles in a warm or cold steady state. This mechanism expands existing theories of millennial-scale variability and provides a general framework for understanding abrupt climate change in general circulation models.

**Keywords** Millennial-scale variability · AMOC · General circulation model · Salt oscillator · North Atlantic stratification · Convection–advection oscillator

## 1 Introduction

### 1.1 Millennial-scale variability during the last glacial period

The Greenland ice core records of the last glacial period, stretching between about 115 and about 12 ka BP (thousand years before present), are dominated by millennial-scale variability (Wolff et al. 2010). Dansgaard–Oeschger (D–O) cycles (Dansgaard et al. 1993; Bond et al. 1993), recurring transitions between cold stadial and warm interstadial North Atlantic climates (Thomas et al. 2009; Lohmann and Ditlevsen 2019), are the most prominent example of such variability. During these transitions, Greenland underwent temperature changes of up to 16 °C (Huber et al. 2006; Kindler et al. 2014; Buizert et al. 2014) within a few decades (Rasmussen et al. 2014). North Atlantic abrupt climate changes have been widely associated with changes in

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

Yvan M. Romé
y.m.rome@leeds.ac.uk

Ruza F. Ivanovic
r.ivanovic@leeds.ac.uk

Lauren J. Gregoire
l.j.gregoire@leeds.ac.uk

<footnote n="1">¹ School of Earth and Environment, University of Leeds,
Woodhouse Lane, Leeds LS2 9JT, UK</footnote>

<footnote n="2">² CNRS, Bordeaux INP, EPOC, UMR 5805, University of Bordeaux, Pessac 33600, France</footnote>

<footnote n="3">³ Faculty of Science, University of the Ryukyus, 1 Senbaru,
Nishihara, Nakagami District, Okinawa 903-0129, Japan</footnote>

<footnote n="4">⁴ Department of Mathematics and Statistics, University of Reading, Whiteknights House, Reading RG6 6UR, UK</footnote>
the Atlantic Meridional Overturning Circulation (AMOC, Rahmstorf 2002; Clark et al. 2002). This is primarily because the AMOC is believed to have shifted between different modes during the last glacial period (Broecker et al. 1985; Gottschalk et al. 2015; Henry et al. 2016; Lynch-Stieglitz 2017; Weijer et al. 2019), disrupting the flux of heat and salt into the North Atlantic. Three dynamical regimes of the AMOC have been hypothesised by Rahmstorf (2002): a *warm* mode where a vigorous AMOC reaches up to 70°N, a *cold* mode where a shallower AMOC does not extend past the Greenland-Iceland-Scotland (GIS) ridge, and an *off* mode where the AMOC is collapsed (Böhm et al. 2015). Simulations of the climate of the next centuries (Weijer et al. 2020; Romanou et al. 2023) and observed fingerprints of potential decline in AMOC strength (Thornalley et al. 2018; Caesar et al. 2018; Boers 2021; Michel et al. 2022; Ditlevsen and Ditlevsen 2023) have raised concerns that the AMOC might undergo an abrupt weakening or collapse again with potentially dramatic consequences on the North Atlantic climate (van Westen et al. 2024).

However, the climate projections are highly model-dependent (Weijer et al. 2020; Bellomo et al. 2021) and direct observations are still too short to draw any statistically significant conclusion on potential shifts in the AMOC regime (Lobelle et al. 2020; Jackson et al. 2022). Reducing the uncertainty in future AMOC stability requires a deeper understanding of the dynamic processes driving the overturning circulation, which can be examined through climate modelling of the last glacial period. In particular, the Last Glacial Maximum (LGM)—about 21 ka BP—and marine isotope stage 3—about 60–30 ka BP—have proven insightful for studying AMOC stability (Weijer et al. 2019). Freshwater hosing experiments have historically been an effective way to trigger abrupt climate changes in most climate models (Manabe and Stouffer 1995; Kageyama et al. 2010; Hawkins et al. 2011; Kageyama et al. 2013; van Westen and Dijkstra 2023). When it comes to simulating spontaneous oscillations, however, General Circulation Models (GCMs) were for a long time regarded as too biased towards producing more stable climates (Valdes 2011; Mecking et al. 2017; Liu et al. 2017). This bias was addressed by making use of modern supercomputers to explore a comprehensive range of model inputs (including uncertain model parameters and climate forcings) and running models for multiple millennia. As a result, more GCMs have been found to exhibit millennial-scale variability in the past decade (e.g., Peltier and Vettoretti 2014; Brown and Galbraith 2016; Klockmann et al. 2018; Zhang et al. 2021; Armstrong et al. 2022; Kuniyoshi et al. 2022; Malmierca-Vallet et al. 2023).

Barker and Knorr (2021) coined the term ‘window of opportunity’ to define a model parameter regime where millennial-scale variability can occur. Not all models share the same windows of opportunity, i.e. the location and extent

of the window of opportunity in the parameter space varies across models. Most models agree, however, on the most relevant parameters that determine its definition; namely the atmospheric greenhouse gases, orbital parameters, ice sheet topography and freshwater release. Atmospheric CO₂ concentrations (Brown and Galbraith 2016; Zhang et al. 2017; Klockmann et al. 2018; Vettoretti et al. 2022; Malmierca-Vallet et al. 2024) and insolation (Zhang et al. 2021; Kuniyoshi et al. 2022) are essential to the creation of conditions facilitating abrupt climate change through their influence on the Atlantic energy balance, atmospheric and oceanic circulations and stratification. The ice sheet configuration is equally important, as the extensive ice sheets covering North America and northern Europe during the last glacial maximum (Hughes et al. 2016; Batchelor et al. 2019) act to prevent the occurrence of AMOC mode shifts (Zhang et al. 2014; Brown and Galbraith 2016; Klockmann et al. 2018). Glacial ice sheets tend to drive a less variable, more southern and potentially faster eddy-driven jet (Merz et al. 2015; Li and Born 2019). This results in stronger surface winds and wind stress curl that intensify the salt and heat transport to the deep water formation sites (Czaja 2009; Montoya et al. 2011; Sherriff-Tadano et al. 2018, 2021). The topography and albedo of the ice sheets also modify the inter-hemispheric heat fluxes (Roberts and Valdes 2017). In addition, glacial ice sheets can produce larger freshwater fluxes than the pre-industrial ice sheers, as observed in sea level records of the last glacial period (Fairbanks 1989; Lambeck et al. 2014). This is another critical input because freshwater can enhance ocean stratification by decreasing the surface layer density, reducing deep convection and thereby weakening the AMOC. Freshwater forcing is particularly efficient at slowing down the AMOC when targeted at regions of deep water formation, such as within the Nordic Seas (Smith and Gregory 2009; Roche et al. 2010), but the exact rate and location of freshwater release during the last glacial period discharge events remains challenging to constrain (Birchfield et al. 1994; Bethke et al. 2012).

## 1.2 Mechanisms of millennial-scale variability

With his pioneering box model, Stommel (1961) demonstrated that heat and salt fluxes between two interconnected basins could lead to two distinct stable regimes of ocean circulation. This hypothesis was influential in the establishment of theories for AMOC regime shifts in the context of millennial-scale variability. As an alternative to multistability, Broecker et al. (1990) introduced the salt oscillator as a mechanism for millennial-scale variability. In this theory, when the AMOC is weak, salt is accumulated in the Atlantic by evaporation and increases the surface density, leading to the reactivation of North Atlantic vertical convection and the overturning circulation. When the AMOC is strong,

<figure><img src="image_3.png" /></figure>
the overturning circulation exports salt out of the Atlantic basin and the surface density decreases, leading to the deactivation of North Atlantic deep water formation sites and a weakening of the overturning circulation. The validity of this mechanism depends on the sign of the net salt export out of the Atlantic Ocean (Huisman et al. 2010; Drijfhout et al. 2011), and an initial perturbation is required to start the oscillator. Millennial-scale variability can also arise from local instabilities. Spontaneous oscillations in a simple mixed-layer model were demonstrated by Welander (1982), and Cessi (1996) argued that a freshwater input could trigger this oscillator. Winton and Sarachik (1993) also explored the role of vertical instability of the water column to introduce the concept of deep-decoupling oscillations. Deep-decoupling oscillations rely on the different response times between the surface and the deep North Atlantic waters, and the building and collapse of the halocline in the North Atlantic deactivates or reactivates deep convection. Colin de Verdière (2007) proposed a model where the salt-oscillator theory described by Broecker et al. (1990) and the vertical mixed-layer oscillations of Welander (1982) work in concert to produce millennial-scale oscillations. Finally, stochastic resonance was suggested as another possible mechanism to trigger abrupt regime shifts (Alley et al. 2001; Ganopolski and Rahmstorf 2002; Cimatoribus et al. 2013). This theory relies on stochastic noise to amplify internal ocean variability and kick-start oscillations (Schulz et al. 2002; Timmermann et al. 2003).

An effective framework for comparing complex climate model simulations with theoretical models is to conceptualise millennial-scale variability mechanisms as a combination of slow and fast components (Roberts and Saha 2017; Vettoretti et al. 2022). Slow components modify the atmospheric and oceanic conditions in the North Atlantic, while fast components consist of a combination of thresholds and positive feedbacks that can rapidly alter deep water formation activity. Although the scope and definition of the salt oscillator has varied since its introduction (Broecker et al. 1990), it remains a central concept to describe the slow components (e.g. Vettoretti and Peltier 2018; Klockmann et al. 2020; Armstrong et al. 2022). Vettoretti and Peltier (2018) showed that southwards transport of Arctic sea ice was responsible for the freshening of the North Atlantic prior to a cooling event, and the northwards gyre transportation of salt into the North Atlantic was responsible for the salinification of the North Atlantic prior to a warming event. In Klockmann et al. (2020)'s simulations, phases of accumulation and depletion of the subtropical Atlantic salinity affect the density gradient across the subpolar gyre. Similar to Broecker et al. (1990), Armstrong et al. (2022)'s analysis points towards changes in Atlantic precipitation and evaporation patterns due to shifts of the ITCZ. Additionally, they discuss the influence

of sea ice transport on the meridional salinity gradients in the North Atlantic. In Kuniyoshi et al. (2022)'s study, sea ice transport is also the main driver of the mechanisms, although thermal changes of the North Atlantic water column are more impactful than changes in surface salinity gradient. Thermal control of millennial-scale oscillations was equally considered in previous works (e.g., Oka et al. 2012; Dokken et al. 2013; Brown and Galbraith 2016). Additionally, sources of salt and heat outside of the North Atlantic, such as the Southern Ocean (Knorr and Lohmann 2003; Banderas et al. 2015; Buizert and Schmittner 2015; Thompson et al. 2019; Oka et al. 2021) or the Indian Ocean (Nuber et al. 2023), were also identified as potential drivers of slow changes of the North Atlantic water masses. This is the case, for instance, in the mechanism by Vettoretti et al. (2022) where the density of the AABW modifies the Atlantic vertical stratification. Finally, Weijer and Dijkstra (2003) argues that a salinity anomaly that propagates along the different ocean basins could set the pacing of the oscillations.

The fast dynamics involve positive feedback that can accelerate the reactivation and deactivation of the deep water formation sites when one or several North Atlantic climate thresholds. The salt-advection feedback (Drijfhout et al. 2015; Weijer et al. 2019) is a central concept here: a stronger AMOC advects more salt to the deep water formation sites, enhancing deep convection which further amplifies the AMOC. Conversely, a decrease in AMOC strength further weakens the AMOC by advecting less salt northwards, reducing surface density at the deep convection sites. This is at stake, for instance, in Kuniyoshi et al. (2022) and Vettoretti and Peltier (2018)'s mechanism, where the positive feedback accelerates the initial changes in North Atlantic deep water formation when the vertical stratification in the region passes a certain threshold. In other instances, the feedback can amplify the impact of exceptional events, such as the opening of polynyas to release abruptly the heat accumulated under the surface of the ice (a phenomenon referred to as subsurface warming) in Vettoretti and Peltier (2016). Both the threshold and exceptional events approaches were considered in Klockmann et al. (2020), building on the theory introduced by Li and Born (2019), to describe a positive feedback between the winds and the subpolar gyre leading to shifts in the North Atlantic convection. This feedback can be activated either by a change in the cross-gyre density gradient or by stochastic wind forcing. In Armstrong et al. (2022), the meridional density gradient in the North Atlantic acts as the threshold, and both the wind-driven and the salt advection feedbacks accelerate the changes in North Atlantic convection. Kuniyoshi et al. (2022) describes a sea-ice positive feedback that cools down the North Atlantic during an interstadial to stadial transition. Finally, external forcing such as the freshwater input used in Brown and Galbraith
(2016) and the change of vertical diffusivity in Peltier and Vettoretti (2014) can act as the initial trigger of the fast components.

## 1.3 Objectives of the paper

The Hadley Centre Climate Model 3, also known as HadCM3 (Valdes et al. 2017), is one of the models that has managed to simulate abrupt climate changes under realistic glacial conditions. Signs of AMOC transitions were previously reported in response to freshwater forcing (Matero et al. 2017; Ivanovic et al. 2018), but only recently did two independent studies manage to observe oscillations of the North Atlantic climate: Romé et al. (2022) triggered oscillations using deglacial meltwater patterns in a last glacial maximum background state and Armstrong et al. (2022) described oscillations that occur for a specific snapshot of marine isotope stage 3 conditions (30 ka BP boundary conditions). In this work, we aim to analyse the processes at stake in the LGM simulations presented in Romé et al. (2022) and explore how they connect to other mechanisms derived from HadCM3 and other model simulations. We do not try to directly identify the mechanisms behind D–O cycles, as the experimental design prevents a direct comparison. After presenting the experimental design of the simulations (Sect. 2) and summarising the results and terminology introduced by Romé et al. (2022) (Sect. 3), we separately examine the changes in stratification at the North Atlantic convection sites (Sect. 4.1) as well as the global salt exchanges (Sect. 4.2). Section 5 brings these insights together to introduce the convection–advection oscillator mechanism. Section 6 analyses how our mechanism applies to non-oscillating simulations and could help characterise AMOC stability in general circulation models. Finally, Sect. 7 explores the limitations of this theory and how it connects to previous mechanisms for millennial-scale variability.

# 2 Methods

## 2.1 The HadCM3 general circulation model

All the simulations in this study used the BRIDGE (Bristol Research Initiative for the Dynamic Global Environment group) version of the HadCM3 (HadCM3B) atmosphere-ocean general circulation model developed by Valdes et al. (2017). The atmosphere model, described by Pope et al. (2000), consists of a 2.5° x 3.75° grid on 19 elevation levels. The ocean model, described by Gordon et al. (2000), is a 1.25° x 1.25° grid on 20 depth layers. In order to verify the equation of state (Bryan and Cox 1972; Fofonoff and Millard 1983; Fofonoff 1985), the ocean model uses a rigid surface lid (Gordon et al. 2000), meaning that the ocean volume

stays constant throughout the simulations. In addition to the ocean and atmosphere components, the model includes the MOSES 2.1 land model (Cox et al. 1999) and the TRIFFID dynamic vegetation model (Cox 2001). The land-sea mask is matched to the chosen palaeo ice sheet reconstruction. The Gibraltar Strait is closed during the LGM—as it is for the present, due to the spatial resolution of the model—and Mediterranean-Atlantic salt and heat exchanges are parametrised by a diffusive pipe (Ivanovic et al. 2013, 2014). HadCM3B was optimised for multi-millennial palaeoclimate runs (Valdes et al. 2017) and has been used previously for glacial simulations with an imposed meltwater forcing (e.g., Ivanovic et al. 2018; Romé et al. 2022).

## 2.2 Global mean salinity conservation

HadCM3 uses a rigid-lid ocean model to compute water fluxes. The rigid-lid model implies that the volume of the ocean remains constant throughout the simulations, and water fluxes, whether they are freshwater discharge, precipitation or evaporation, are implemented in the model as a change of the surface layer salinity rather than a physical input of water. This also means that the land-sea mask is not updated in response to the freshwater fluxes in such experiments.

This method does not guarantee the closing of the hydrological budgets causing shifts in the global mean salinity. Global mean salinity drifts of the order of 0.25 g.kg⁻¹.kyear⁻¹ can spontaneously occur in HadCM3 even without freshwater input (Dentith et al. 2019). In other cases, global salinity increase can also be observed if the model configuration is biased towards higher evaporation than precipitation, leading to salt accumulation in inland seas and the tropical oceans. On top of these background shifts in global salinity, an input of meltwater of around 0.1 Sv, the order of magnitude of the discharge used in this study, can drive a significant decrease in salinity comparable to that observed during the last deglaciation over 10.000 years of simulations.

The aim of this set of simulations is not to simulate a deglaciation and the continued freshening of the global ocean during the long integration should be avoided. To tackle this issue, we use a global salinity conservation algorithm—this corresponds to the VFLUX method described by Dentith et al. (2019). The VFLUX algorithm sets a global mean salinity target (shortened salinity target in the following) for the entire run. At the end of each time step, the model corrects the salinity uniformly throughout each depth so that the corrected global mean salinity matches the value of the salinity target. The global mean salinity target is calculated from the 1950 CE (Common Era) global mean salinity, adjusted to account for the differences in the ice sheets volume between the present and a given point in

<figure><img src="image_4.png" /></figure>
time. Note that the change in salinity target only covers ice volume changes and does not take into account additional freshwater input. The new salinity target is implemented at the start of each simulation, and a transient regime where the model adjusts to the new conditions can be observed in the initial millennium (Romé et al. 2022).

## 2.3 Experimental design

The simulations in this paper, presented in Table 1, are branched from the *CTRL_LGM* (or *control*) simulation at year 0. The *control* simulation is a last glacial maximum run following the 21 ka BP PMIP4 protocol for the boundary conditions (greenhouse gas concentrations, orbital parameters and solar constant, Kageyama et al. 2017). We used the GLAC-1D (Tarasov and Peltier 2002; Tarasov et al. 2012; Briggs et al. 2014; Ivanovic et al. 2016) ice sheet reconstruction and the corresponding land-sea mask, orography and bathymetry. *CTRL_LGM* extends a previous HadCM3 last glacial maximum simulation (Davies-Barnard et al. 2017), and starts after an initial 3,500 years spin-up with the updated boundary conditions. *CTRL_LGM*'s start year will be the year 0 in our time series.

The main focus of this paper will be the 20.7k simulation described in Romé et al. (2022). The 20.7k simulation differs from *CTRL_LGM* only by the addition of a **constant** meltwater discharge corresponding to the 20.7 ka BP slice of GLAC-1D meltwater history derived by the protocol introduced in Romé et al. (2022). For each time step of the GLAC-1D ice sheet reconstruction, freshwater discharge was computed from changes in ice sheet elevation at each grid cell and routed to an ocean grid cell using an offline drainage map of the last glacial maximum (Ivanovic et al. 2016). Surface freshwater saturation (i.e. a grid cell salinity of 0 g.kg⁻¹) can be reached in some grid cells when large freshwater fluxes are recorded. To avoid this situation, the meltwater is spread across a wider (proximal ocean) region

around the source grid box. Spreading regions are in the order of 10⁶km², about a hundred ocean grid cells, an area comparable to a typical Heinrich event (Roberts et al. 2014). From this meltwater forcing history, we derive the freshwater pattern corresponding to the time slice we are interested in, in this case 20.7 ka BP. This pattern will be called the meltwater snapshot and used as a constant forcing for the 10,000 integration years. Independently from the freshwater discharge, we change the global mean salinity target to match the meltwater snapshot, in this case 35.8255 g.kg⁻¹ for 20.7 ka BP.

To get a more precise picture of the ocean dynamics in our simulations, we repeated the 20.7k simulation for 4000 years, saving an extra set of salinity diagnostics, or salinity tendencies, in the model output. This simulation is called 20.7k_tdc. Armstrong et al. (2022) introduced the HadCM3 tendencies diagnosis and included a comprehensive description of the implementation. We summarise the essential information in Section S6.

We included two additional simulations from Romé et al. (2022) for comparison. The 21k simulation used the meltwater snapshot corresponding to 21 ka BP, a global salinity target of 35.8334 g.kg⁻¹ and was run for 4000 years. The 18.2k simulation used the meltwater snapshot corresponding to 18.2 ka BP, a global salinity target of 35.7348 g.kg⁻¹ and was run for 10,000 years. Finally, we investigated the impact of the global salinity target by running a simulation with the salinity correction algorithm turned off. This simulation, called 20.7k_no_gst, used the same meltwater snapshot as 20.7k and was run for 3000 years.

## 3 Anatomy of the simulations

Romé et al. (2022) introduced a method for clustering the meltwater snapshot simulations between different convection modes characterising the North Atlantic deep water

<table><caption>**Table 1** Experiments summary</caption><thead><tr><th>Simulation name</th><th>Meltwater snapshot<br>Meltwater flux (Sv)</th><th>Integration length<br>(years)</th><th>Salinity Target<br>(g.kg<sup>-1</sup>)</th><th>Tendencies<br>(included)</th></tr></thead><tbody><tr><td>CTRL_LGM</td><td>None</td><td>4000</td><td>35.8334</td><td>No</td></tr><tr><td>20.7k</td><td>20.7 ka BP (0.084)</td><td>10,000</td><td>35.8225</td><td>No</td></tr><tr><td>21k</td><td>21.0 ka BP (0.054)</td><td>4000</td><td>35.8334</td><td>No</td></tr><tr><td>18.2k</td><td>18.2 ka BP (0.109)</td><td>10,000</td><td>35.7348</td><td>No</td></tr><tr><td>20.7k_tdc</td><td>20.7 ka BP (0.084)</td><td>4000</td><td>35.8225</td><td>Yes</td></tr><tr><td>20.7k_no_gst</td><td>20.7 ka BP (0.084)</td><td>3000</td><td>None</td><td>No</td></tr></tbody></table>

All experiments were designed with LGM boundary conditions, using the LGM GLAC-1D ice sheet extent and associated geographies. The meltwater entries indicate the meltwater snapshot used (see Sect. 2.3) and the total meltwater flux in brackets. The meltwater flux only corresponds to the freshwater from the ice sheet melting and does not include the other freshwater fluxes (e.g., precipitation, sea ice, etc.). The Salinity Target (g.kg⁻¹) entries are the global salinity target used for the salinity correction described in Sect. 2.2. Tendencies indicates if the salinity tendencies diagnostic described in S6 were included
formation layouts. An oscillating simulation is character-
ised by the periodic alternation of ‘warm’ modes—vigorous
AMOC with relatively warm North Atlantic temperatures—
and ‘cold’ modes—weak AMOC with relatively cold North
Atlantic temperatures. The 20.7k simulation is an example
of an oscillating simulation where this sequence can be
observed in Fig. 1a. If instead of sustained quasi-oscillations
the AMOC settles in a warm or cold steady state, the simu-
lation is referred to as a *warm* or *cold* simulation, respec-
tively. The aim of Romé et al. (2022) was not to reproduce
D–O cycles and we will therefore avoid using the commonly
applied ‘stadial’/‘interstadial’ terminology. Nonetheless, it
can be helpful to think of warm and cold modes in the light
of interstadials and stadial states.

To characterise a simulation, we can look at its path in the
state space projected onto the Mixed Layer Depth (MLD)—
an indicator of the mixing and convection of the upper water
column—in the GIN Seas against the MLD in the Irminger
Sea (as plotted in Fig. 1b). A detail of the zones is given
in Figure S1b. In this state space projection (Fig. 1b), we
manually define five distinct phases of the oscillatory cycle.
The cold phase, where deep convection is weak in the North

Atlantic; the meridional *warm* phase (or simply meridional phase) with a main convection region in the GIN Seas and an intermediate convection spot in the Irminger Sea; the *zonal* *warm* phase (or simply *zonal* phase) dominated by convection in the Irminger Sea and low convection in the GIN Seas; the *warming* phase, a transition from the *cold* phase to the meridional phase and the *cooling* phase, a transition from the *zonal* phase to the *cold* phase. In this paper’s simulations, the transition from the meridional to the *zonal* phase is too fast to be considered a distinct phase. Each phase was assigned a colour, and coloured bars will be shown in time series plots of this article in the style of Fig. 1a.

The oscillating 20.7k simulation evolves along a well-
defined cycle with a periodicity of about 1540 years (see
Section S3). After an initial drop to a *cold* phase, it starts
by *warming* from the *cold* phase to the *meridional* phase,
switching to the *zonal* phase after a few hundred years before
*cooling* back to the *cold* phase. The surface temperature in
Greenland typically varies by 8°C over the course of an
oscillation (Romé et al. 2022). It takes a couple of cycles
for the pattern to become more regular as the climate system
has completely adjusted to the freshwater perturbation and

**Fig. 1** Anatomy of the oscillations. **a** Time-series of AMOC index (maximum overturning circulation at 26.5° N) in the 20.7k simulation. The background colour bars indicate the phases with *cold* in dark blue, *warming* in pink, the meridional in red, *zonal* in orange and *cooling* in light blue. Solid lines indicate the 30-year running mean and transparent lines the annual mean. **b** State space visualisation (i.e. scatter-plot) of the annual mean mixed layer depth in the Greenland Iceland Nordic (GIN) Seas against the mixed layer depth in the Irminger Sea (the zones are defined in Figure S1b). Each point corresponds to a year of the 20.7k simulation time series after being processed by a low-pass filter with a cut-off frequency of 2 10⁻³ year⁻¹—the algorithm was described in Section S3. Shaded boxes define the phases shown in panel **a**. In the small globes, the composite mean winter (DJF) mixed layer depth is plotted for the *cold* and the two *warm* phases

<figure><img src="image_5.png" /></figure>
the change of salinity target, and we will therefore refer to
the first ~4000 years as the *transient* regime.

The AMOC time series and state space plots of the 18.2k, 21k, 20.7k_no_gst and 20.7k_tdc simulations are shown in Figure S2. Across all simulations, the modes are associated with the same convective dynamics as the phases of the oscillators, although they span slightly different areas on the two-dimensional MLD state space projection due to the different freshwater release patterns. Note that in the general case, we call a *warming* mode a transition between a *cold* mode to any *warm* mode and a *cooling* mode from any *warm* mode to a *cold* mode. The *CTRL_LGM* simulation is a *warm* simulation at a steady state in a *zonal* mode. All meltwater release simulations respond to the initial melt-water release with a multi-centennial drop of the AMOC. The 21k simulation then recovers into a *zonal* steady state and thus is classified into the *warm* simulations, similarly to the 20.7k_no_gst simulation that follows one oscillatory cycle before stabilising in a *zonal* steady state. On the other hand, the 18.2k simulation mostly remains in a *cold* steady state after the initial perturbation, albeit interrupted by short-lived excitations into the *meridional* mode. Because

of the irregular and transient nature of these excitations, this
simulation is categorised as a *cold* simulation. Finally, even
though the *oscillating* 20.7k and 20.7k_tdc simulations are
effectively the same experiment, the noise introduced when
the model was restarted to initialise the salinity tendencies
diagnostics (a process called ‘reconfiguration’) results in
slight differences in the periodicity and shape of the oscil-
lations (Figure S2a, e, f, j). This is due to the chaotic nature
of the climate system, and we consider the two simulations
dynamically equivalent.

4 Results

**4.1 Stratification changes at the North Atlantic deep water formation sites**

The strength of the AMOC depends on the intensity of
the convection at the North Atlantic deep water forma-
tion sites, which in turn depends on the vertical stratifica-
tion of the water column. The stratification can be inferred
from the density profile in Fig. 2p-r, plotted here for the

<figure><img src="image_6.png" /><figcaption>**Fig. 2** North Atlantic dynamics at deep water formation sites. Plotted for the 20.7k simulation in the GIN Seas (left), the Irminger Sea (middle) and the Iceland Basin (right). The geographical zones are defined in Figure S1, and the background shading corresponds to the phases defined in Fig. 1. **a–c** Early spring (March) and early autumn (September) sea ice concentration. **d–f** Annual mean mixed layer depth. **g–i** Annual mean anti-clockwise barotropic stream function</figcaption></figure>

indicating gyre circulation strength. **j–l** Annual mean temperature in the upper, intermediate and deep depth slices. **m–o** Annual mean salinity in the upper, intermediate and deep depth slices. **p–r** Annual mean density in the upper, intermediate and deep depth slices. In panels **a–i**, solid lines indicate the 30-year running mean and transparent lines the annual mean
20.7k simulation between the *upper* waters (above 550 ms deep) and the *intermediate* waters (550 to 1800 meters deep), and the *deep* waters (below 1800 meters deep). The depth slices were chosen so that the *intermediate* waters always lie outside of the mixed layer boundaries while still feeling the influence of the upper waters' changes in temperature and salinity. In this section, we focused our analysis on the GIN, Irminger and Iceland basins (Figure S1b). We only included the last two cycles of the simulation to ensure minimum interference with the transient regime at the start of the simulations. A comprehensive mapping of the North Atlantic and Arctic ocean profiles is included in Figure S4. We consider the conclusions of the following analysis to be the same between the Arctic basin and the Nordic Seas, and the Labrador Sea and the Irminger Sea.

During the *cold* phase, **high stratification** prevents deep convection in the North Atlantic. The deep water formation (diagnosed through mixed layer depth, Fig. 2d–f) and the North Atlantic subpolar gyre strength (Fig. 2g, h) are at their weakest. Winter sea ice extends down to the Iceland Basin, and the Nordic seas are covered with ice all year long (Fig. 2a–c). The low density of the upper North Atlantic (Fig. 2p–r) is due to a deficit in salinity (Fig. 2m, n), as well as higher temperatures in the GIN and Irminger seas (Fig. 2j, k). This high stratification tends to decrease throughout the *cold* phase, and in particular in the GIN Seas through the combined effect of subsurface warming (warming under the mixed layer, most evident in intermediate waters temperatures) over the entire water column (Fig. 2j), a salinity increase of the upper waters and a salinity decrease of the intermediate and deep waters (Fig. 2m).

During the **warming** phase, the **stratification decreases** abruptly. The salinification of the upper GIN Seas gains momentum while the intermediate and deep waters keep freshening (Fig. 2m). As deep convection recovers (Fig. 2d) and sea ice retreats (Fig. 2a), the heat accumulated under the ice is evacuated to the atmosphere (Fig. 2j) further accelerating the stratification decrease (Fig. 2p). Similar effects are observed in the Irminger Sea (Fig. 2q) and in the Iceland Basin (Fig. 2r), although the changes in the density profile are less abrupt.

At the end of the warming phase, deep water formation resumes in the North Atlantic (Fig. 2d–f), marking the start of the **warm** phases characterised by **low stratification**. At the onset of the **meridional** phase, the subpolar gyre is reactivated (Fig. 2e) and sea ice cover reaches a minimum in the GIN Seas and the Iceland Basin (Fig. 2a, c). The density gradient reduces between the upper and deep waters (Fig. 2p–r), mostly driven by a drop in the upper waters salinity (Fig. 2m–o), salinity gradient, resulting in a minimum of stratification in the entire North Atlantic. In the GIN Seas, this minimum only lasts for a few hundred years before upper salinity starts to decrease again (Fig. 2m, p).

This decline eventually leads the deep convection to shut down (Fig. 2d), preceding a sharp recovery of the subpolar gyre strength and of the Irminger Sea deep water formation, marking the start of the **zonal** phase. During the **zonal** phase, the slow decrease of upper waters salinity is compensated for by a reduction of the temperatures (Fig. 2j–o), and the density profiles are mostly unaltered (Fig. 2q–s). Sea ice extent, however, is consistently increasing both in summer (Fig. 2a) and winter (Fig. 4.1b, c). Disruptions appear at the end of the **zonal** phase, when the salinity decrease in upper waters intensifies (Fig. 2o) and is no longer compensated for by subsurface heat evacuation (Fig. 2l).

This triggers the start of the **cooling** phase where **stratification increases** abruptly. The subpolar gyre shrinks (Fig. 2e), subsurface heat accumulation restarts (Fig. 2d–f) and the salinity decreases further in the upper waters (Fig. 2m–o). The intermediate waters density reduction does not match the pace of the upper waters (Fig. 2p–r), leading to a collapse of deep convection in the North Atlantic (Fig. 2d–f) and a return to a **cold** phase.

In this section, we showed that stratification changes are due to the combined effect of salinity fluxes and subsurface heat accumulation at the deep water formation sites. While the thermal effect can be linked to changes in sea ice cover, there is no obvious explanation for the vertical salinity changes through dynamical processes within the North Atlantic. A deeper analysis of the global salinity transfers is therefore required.

## 4.2 Global salt exchanges regulating North Atlantic salinity

Fluxes of salt in and out of the cluster formed by the North Atlantic and Arctic basins, referred to as the *North Atlantic waters* in the rest of the analysis, can disrupt the stratification of the deep water formation sites. Figure 3a shows that in the 20.7k simulation, these two regions experience mean salinity shifts larger than in any other basins (Figure S5, the zones used in this paragraph are defined in Figure S1a). Between the *cold* and *warm* phases, the mean salinity changes by around 0.2 g.kg⁻¹ in the North Atlantic and 0.35 g.kg⁻¹ in the Arctic, reaching up to 0.75 g.kg⁻¹ in the North Atlantic and 1.1 g.kg⁻¹ in the Arctic if we restrict the domain to the first 550 ms of the water column (Figure S5). Such variations of mean salinity only take small changes in the salinity budgets (Fig. 3b) because of the relatively small volume of these two basins. We observe a strong correlation between the AMOC time series and the salinity content in the North Atlantic waters in Fig. 3c. This implies that the AMOC is the main driver behind salinity changes in the domain. The North Atlantic waters' salinity reaches its minimum roughly at the same time as the AMOC index (~750 years after its maximum, Fig. 3c), but the North Atlantic
**Fig. 3** Inter-basin salinity oscillations. Time-series of the 20.7k simulation anomaly relative to the *CTRL_LGM* simulation for a mean salinity and b salt content. Background shading indicates the phases defined in Fig. 1. Cross-correlation between the normalised salinity budget anomalies and AMOC index (see Fig. 1a) in c the North Atlantic waters, d the Subtropical Atlantic waters and e the Deep waters domain. The lag phase in the cross-correlations is calculated with the 1540-year periodicity derived in Section S3. For readability purposes, only the signals from the upper and deep Arctic, North Atlantic, subtropical Atlantic, Southern Ocean and Pacific are plotted in the foreground in c–e, with the other cross-correlations visible in the background

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

waters' salinity variations are more muted than the AMOC signal during the cold mode.

The subtropical Atlantic mean salinity is ahead of the North Atlantic signal by ~500 years (Fig. 3a, d). The correlation and anti-correlation coefficients are around 0.8. and –0.8., respectively, due to less sharp phase transitions in the subtropical Atlantic salinity signal. Despite large changes in the salt budget (Fig. 3b), the mean salinity in the subtropical Atlantic varies by less than 0.1 g.kg⁻¹ throughout the simulation (Fig. 3a). The subtropical Atlantic signal is echoed in the deep North Atlantic and the upper South Atlantic, forming another domain that we call the *Subtropical waters*. The upper layer of the Pacific basin follows a similar trend but is delayed by ~200 years compared to the *Subtropical waters*. To conserve global salt content (see Sect. 2.2), the large changes of salt content in the *Subtropical waters* have to be compensated for in other basins.

We find changes proportional but of opposite sign to the salinity in the *Subtropical waters* in the Pacific and, to a lower extent, the Southern oceans (Fig. 3b). The salt content anomaly in these basins is dominated by the intermediate and deep water signal (below 550 ms, Figure S6). This domain, which we call the *Deep waters*, shows less spatial coherence as changes in the deep waters (below 1800 ms) experience a delay of ~250 years compared to changes in the intermediate waters (e.g. Pacific signal in Fig. 4.2e). Overall, the *deep waters* signal lags behind the AMOC signal by ~250 years, and its opposite signal leads the AMOC by ~500 years.

In summary, the global salt fluxes resemble a relaxation oscillator where salt is slowly transported from the deep Pacific reservoir to the upper waters of the Atlantic during the *cold* phases, and exported from the Atlantic to the Pacific during the *warm* phases (Fig. 3b). Salt anomalies in the North Atlantic are remarkably coupled to changes in the AMOC strength and follow the *Subtropical* and *Deep* waters salt exchanges with a 500-year delay. The global salt redistribution influences the North Atlantic stratification changes described in Sect. 4.1, and both components come together to explain the mechanisms for the AMOC oscillations that will be described in Sect. 5.

# 5 The convection–advection oscillator mechanism

In this section, we introduce the **convection–advection oscillator** mechanism to explain the millennial-scale variability simulated in 20.7k. This mechanism, summarised in Fig. 4, is based on the fast/slow components coupling framework introduced in (Roberts and Saha 2017). The North Atlantic deep water formation activation in response to stratification thresholds makes up the **fast convection** component, and the large scale reorganisation of global salt content changing buoyancy fluxes in the North Atlantic makes up the **slow advection** component. During *cold* phases, the North Atlantic stratification is high and the AMOC is weak (Fig. 2). Salt accumulation
<figure><img src="image_8.png" /><figcaption>**Fig. 4** The oscillator mechanism. Millennial-scale variability mechanism observed in 20.7k. The colours correspond to the phases defined in Fig. 1 to the regions defined in Figure S1</figcaption></figure>

in the subtropical Atlantic (Fig. 3b) leaks into the North Atlantic (Fig. 2p) while the accumulation of subsurface heat decreases the deep water formation sites' stratification (Fig. 2p). In about 500 years, the gradual drop of stratification passes a threshold where the fast convection component is reactivated. There is a positive feedback loop between the strengthening of North Atlantic deep convection and salt advection into the Upper North Atlantic that leads to the abrupt acceleration of the AMOC. After a transition to a *warm* phase, the North Atlantic stratification is low and the AMOC is strong (Fig. 2o, r). The overturning circulation flushes the salt excess to the deep oceans, and the Atlantic becomes progressively salt depleted (Fig. 3b), reducing its upper water density. This slow decrease lasts for 500–1000 years and the southwards transport of sea ice in the North Atlantic (Fig. 2a–c) progressively resumes the subsurface warming (Fig. 2k, l). A more abrupt drop is observed at the start of the *cooling* period when a new threshold is crossed and the fast convection component is reactivated. The positive feedback loop between the decrease of North Atlantic convection and the slowdown of salt import into the North Atlantic

leads to the deactivation of the AMOC and the system returns to a *cold* phase.

The **meltwater discharge** creates a fresh surface layer in the North Atlantic that reduces the upper salinity and can lead to deep water formation deactivation when targeted at deep convection sites [see Section 4 by Romé et al. (2022)]. The decrease of upper water salinity triggers an initial shift from the *control zonal* phase to a *cold* phase before the deep waters have time to adapt to the perturbation. This transition happens in the first few hundred years and kickstart the oscillator. It takes two to three cycles until the end of the transient regime where the influence of the freshwater perturbation is observed, for instance, in the deep density profiles in Figure S4). For this reason, we focus the following analysis on the last three cycles of the *oscillating* simulation. In the rest of the run, the constant meltwater input acts to decrease the upper North Atlantic density.

The rest of this analysis relies on the more detailed information gathered from the salinity tendency diagnostics presented in Section S6. In order to reactivate North Atlantic convection, the freshwater forcing needs to be compensated for by salinity fluxes into the North Atlantic basin: this is

<figure><img src="image_9.png" /></figure>
the role of the **slow advection** component. The subtropical Atlantic is a net evaporation basin and tends to accumulate salt in the absence of the AMOC advective flux (Fig. 3b, Figure S9). This excess of salt is eventually transferred to the upper North Atlantic by gyre circulation and diffusion (Figure S8), increasing the upper North Atlantic density (Fig. 2p–r, Figure S12). When the AMOC is reactivated, the salt excess is progressively evacuated from the upper North Atlantic to the deep ocean through convection (Fig. 3e, Figure S4), and then redistributed along the oceanic pathway. The global salinity transfer can be seen as a salinity anomaly slowly moving between the subtropical Atlantic and the Pacific (Weijer and Dijkstra 2003), and it takes more than a thousand years for the waters for the global salinity exchanges to reach a new equilibrium state (Fig. 3b).

The combination of salt fluxes entering the North Atlan-
tic, subsurface warming and sea ice transport modify strati-
fication at the deep water formation sites and trigger the
**fast convection component**. During the *cold* phases, the
insulating effect of sea ice described by Marcott et al. (2011)
led to a warming and density drop of the entire water column
(Figure S12). This effect is the most prominent at high lati-
tudes, including the GIN Seas, where the sea ice is perennial
(Fig. 2a, p). In the upper North Atlantic waters, the sea ice
effect is superseded by the salinity input from the subtropical
Atlantic (Figure S12a–c). In the GIN Seas, a vertical strati-
fication is eventually crossed leading to the reactivation of
deep convection in this region. The salt advection feedback
(Drijfhout et al. 2015; Weijer et al. 2019) accelerate this tran-
sition into a *warming* phase.

During the *warm* phases, the decreasing salt availabil-
ity in the subtropical Atlantic reduces upper North Atlan-
tic water density (Figure S8), and the slow evacuation of
subsurface heat accumulated in the Irminger and Iceland
basins during the cold phase increases deep waters density
(Figure S13). Despite the gradually expanding sea ice in the
North Atlantic, the absence of a clear salinity signal in the
region that would indicate the brine rejection (Figure S8)
indicate that the sea ice is more likely to be transported in
rather than formed in the regions. Instead, further reduc-
tion of upper water density (Fig. 2p–r) can be linked to an
increase of temperatures resulting from the insulating effect
of sea ice (Figure S13d–h). As the stratification increases,
it crosses a threshold where the positive feedback leads to a
sharp decline of the AMOC.

6 Sensitivity of the convection–advection
mechanism to climate conditions

To gain insight into the range of climate conditions under
which the convection–advection oscillator mechanism
(Sect. 5) is active, we extend our analysis to the additional

glacial meltwater experiments. A comparison of the dynam-
ics of the main simulations in this article is given in Fig. 5.

The magnitude of the freshwater perturbation is criti-
cal for sustaining the oscillations. In the 21k simulation, a
weaker freshwater forcing compared to the 20.7k simula-
tion does not significantly alter the global salt distribution
(Fig. 5i). The AMOC index only drops by 5 Sv (Fig. 5c)
and the salinity anomaly is dispersed before a stratification
threshold is crossed in the North Atlantic (Fig. 5l). The 21k
simulation eventually recovers to a zonal steady state, the
sea ice edge settles after another 500 years (Fig. 5f), and
the simulation shows no further signs of instability. The
warm steady state simulations, such as the 21k simulation,
behave like an over-damped oscillator: the applied forcing
transiently kicks the system out of the initial state, but the
system returns to a warm state after less than a thousand
years. This suggests that the zonal state retains its stability
and that the simulation remains in its basin of attraction fol-
lowing the onset of the forcing.

In the 18.2k simulation, a strong meltwater input triggers a near-collapse of the AMOC and high North Atlantic stratification (Fig. 5k) prevents the recovery of deep water formation. It takes about 1000 years for the global salinity distribution to reach a new equilibrium (Fig. 5h), without any stratification threshold being crossed that could reinvigorate the AMOC. In the case of the 18.2k simulation, there are short-lived excursions from a stable cold state to a transient meridional state, reminiscent of excitable behaviour. The cold simulations are therefore defined as cold excitable steady state simulation with the cold mode acting as an attractor. As an alternative to multistable, oscillatory or resonant dynamics, Riechers et al. (2024) conducted a study where an excitable monostable model can capture the shape of D–O cycles. In our case, the 18.2k simulation does not resemble D–O variability because of short warm modes and the absence of slow cooling. The AMOC recoveries result from the warming of the North Atlantic water column visible in Fig. 5n, fuelled by sea ice accumulation over the deep water formation sites (Fig. 5e). The recoveries are strictly confined to the GIN Seas, and the stratification of the Irminger Sea remains too high to transition into a zonal steady state. In other words, the salinity import from the subtropical Atlantic is not strong enough to compensate for the freshwater forcing and, therefore, to trigger the positive feedback that would reactivate the AMOC. Although the excursions into a meridional mode occur at what seems to be random intervals, their frequency increases over time due to continuous subsurface heat accumulation, indicating that the asymptotic state may not have been reached.

In between the *warm steady state* and the *cold excit-
able steady state*, the *oscillating* simulations are obtained
in a small window of opportunity where the freshwater
forcing is strong enough to trigger a *cold mode*, but not
<figure><img src="image_10.png" /><figcaption>**Fig. 5** Comparison to the non-oscillating simulations. 30-year running means (solid lines) and annual mean (transparent lines) time series for **a–c** maximum AMOC strength, **d–f** winter sea ice concentration, **g–i** salinity budgets in global ocean basins (see Fig. 3a), **j–l** N2 stratification index (Li et al. 2020) in the deep water formation sites (low negative values indicate high stratification) and **m–o** ocean</figcaption></figure>

potential temperature profiles in the North Atlantic for three meltwater snapshot simulations (see Table 1). The abrupt increase of Pacific salinity in panels **j** and **k** correspond to a change of the basin volume introduced by the bathymetry smoothing algorithm in the Philippines described by Romé et al. (2022)

too strong to prevent a subsequent transition to one of the *warm* modes. In this window of opportunity, the solution is attracted to a chaotic limit cycle (closed curve in the state space) around which it evolves quasi-periodically. Possibly, the *warm* mode undergoes a supercritical Hopf bifurcation as the freshwater input is increased at the relevant locations, giving rise to stable oscillations (Boers et al. 2022). Finally, small changes in meltwater discharge can modify the shape and periodicity of the oscillations described by Romé et al. (2022), but do not change the underlying mechanism.

# 7 The scope of application of the convection–advection oscillator mechanism

The *convection–advection* oscillator mechanism presented here was identified in the 20.7k simulation introduced in Romé et al. (2022), and its relevance may depend on the simulation inputs, the model settings and the choice of the model. The comparison to Armstrong et al. (2022)'s simulations, using the same version of the model, indicates

<figure><img src="image_11.png" /></figure>
that several mechanisms for millennial-scale variability
can coexist within one model. Armstrong et al. (2022) use
the 30 ka BP background climate in which the Irminger
Sea and European coast convection sites are still active
during the cold modes, and they can trigger spontaneous
oscillations in the absence of additional freshwater. The
*convection–advection* oscillator should be seen as a new
occurrence of a *fast-slow* oscillator (see Figure S15) based
on the concept of a relaxation oscillator (Pikovsky et al.
2001; Crucifix 2012; Roberts and Saha 2017). Although
its specific interplay of fast and slow processes is new, the
building blocks of the mechanism have been identified in
former studies. For example, Klockmann et al. (2020) and
Vettoretti and Peltier (2018) highlighted the role of the
subtropical Atlantic salt accumulation and depletion. The
difference of time response between the upper and lower
waters in the North Atlantic as a driver of instability is
central to the deep-decoupling oscillations theory (Winton
and Sarachik 1993), and observed in Vettoretti and Peltier
(2018) and Kuniyoshi et al. (2022). The hypothesis that
salt oscillations can reach further than the Atlantic basin
was introduced in Weijer and Dijkstra (2003). The insu-
lating effect of sea ice is similar to the one described by
Kuniyoshi et al. (2022). Finally, the response of the North
Atlantic convection to changes in stratification fits well
with the description in Vettoretti and Peltier (2018).

On the other hand, our *convection–advection* oscillator
does not include previously proposed *slow* components such
as changes of Southern Ocean density (e.g., Oka et al. 2021;
Vettoretti et al. 2022). These changes can modify the Atlan-
tic stratification in models with a strong sensitivity of the
AABW to the Southern Ocean dynamics, but in the 20.7k
simulation the small density increase and decrease due to
the shrinking and expansion of the AABW during the cold
and warm phases, respectively, are too weak to significantly
impact the North Atlantic convection (Figure S14). Arctic
sea ice transport (e.g., Vettoretti and Peltier 2018) is at play
in our simulations, but sea ice formation and melting appear
to play a lesser role than inter-basin salt transport in regards
to the North Atlantic salinity changes, contrary to Vettoretti
and Peltier (2018) or Armstrong et al. (2022). Armstrong
et al. (2022) also identify large changes of precipitation and
evaporation in the Atlantic as a driver of the salt oscilla-
tor. This tendency is not dominant in our simulations (Fig.
S8). The *fast* components controlling the reactivation of
deep water formation sites can also be the result of differ-
ent physical processes in different experiments. Klockmann
et al. (2020) and Armstrong et al. (2022) discuss the role of
wind driven-feedback in the abrupt re(de)activation of deep
water formation. There is no clear evidence of a similar phe-
nomenon happening in our simulations. Additionally, both
Vettoretti and Peltier (2016) and Kuniyoshi et al. (2022),
associate the overshoot phase with polynyas opening in the

sea-ice-covered Irminger Sea, a phenomenon that is not
observed in any of our simulations. Finally, the regularity
of the oscillations indicates a weak role of stochastic forcing
in 20.7k (Figure S3), but may play a major role in 18.2k.

Another peculiarity of our oscillating simulations is
the existence of two warm states and the transition from a
*meridional* to a *zonal* phase. Cheng et al. (2011) simulated a
similar two-step warming phase, but in their case the abrupt
recovery of deep water formation was observed first in the
Labrador Sea and then in the GIN Seas. In this simulation,
the transition is smooth and the *meridional* phase is akin to
an **overshoot phase** of the *zonal* phase, a phenomenon regu-
larly simulated in GCM millennial-scale variability (Peltier
and Vettoretti 2014; Kuniyoshi et al. 2022; Vettoretti et al.
2022, e.g.,). It happens when deep water formation shifts
from the GIN Seas to the Irminger Sea due to local reor-
ganisation of convection and gyre circulation (Figure S8).
The *meridional* phase witnes a drop in the advective flux of
salt into the North Atlantic (Figure S8b, d) that leads to a
negative salinity import into the basin (Figure S8c), making
this phase unstable. While the transition from a *meridional*
to a *zonal* phase is not an essential part of the overarching
convection–advection mechanism, it may play a more central
role in other simulations using different experiment configu-
rations. It is possible that, for instance, longer *meridional*
phases can be sustained in simulations where the salinity
tendency inversion in the North Atlantic is not obtained,
such as in the 17.8k simulation in Romé et al. (2022).

The combination of fast and slow components involved in
millennial-scale variability implies a multiplicity of possible
mechanisms, probably leading to different shapes, duration
and regularity of the oscillations. In 20.7k, the pace of the
simulations is set by two limiting processes. During the *cold*
phase, it takes about 500 years for a sufficient amount of salt
to leak from the subtropical Atlantic into the North Atlantic
in order to cross a North Atlantic stratification threshold. In
the *warm* phases, the salt anomaly takes more than a thou-
sand years to travel across the entire track of the overturn-
ing circulation to reach the Subtropical Atlantic. These two
processes are dependent on model parametrisation and input
fields. Other factors impacting the emergence of sustained
oscillations could be associated with slow components such
as the rate of Arctic sea ice transport (Vettoretti and Peltier
2018) or the bipolar see-saw (Crowley 1992; Stocker 1998).

In Sect. 5, we showed that the convection–advection
mechanism requires the propagation of a salinity anomaly
between the ocean basins. One of the main parameters that
can impact global salt distribution is the choice of the salin-
ity conservation scheme. In our simulations, the ice sheet
extent does not vary with the meltwater discharge, and there-
fore requires the use of a global salinity conservation algo-
rithm to prevent salinity drifts (see Sect. 2.2 Dentith et al.
2019). In order to test the effect of the salinity conservation
scheme on the convection–advection mechanism, we ran
a simulation without any salinity conservation: the 20.7k_
no_gst simulation. This simulation follows the evolution of
the 20.7k simulation for the first 500 years (Figure S16a, b)
before the salt distribution becomes unstable at the end of
the cold phase: the salinity of all the oceans then increases
in every basin during the zonal phase (Figure S16c, d). In
this configuration, the tropical Atlantic salt content is in
excess, and the North Atlantic upper density remains high,
maintaining active deep-water formation (Figure S16f),
pushing the system away from the window of opportunity.
Salt accumulation may seem counter-intuitive in a fresh-
water hosing experiment, but the rigid lid parametrisation
prevents the model from fully resolving the water cycle
(Demory et al. 2014), which can lead to these unrealistic
drifts. In addition to its role in closing the water budget,
the applied salinity correction is a way to emulate the sta-
bilising power of dynamical ice sheets on global salinity.
This phenomenon could play an active role in the occur-
rence of millennial-scale variability: Wickert et al. (2023),
for instance, proposed a mechanism where the storage and
release of freshwater in the Laurentide ice sheet during the
stadials and interstadials, respectively, can provide a nega-
tive feedback that impedes the deglaciation of the ice sheet.
On the other hand, there are certain times during the last
glacial period when the meltwater released by the desta-
bilisation of the ice sheets, such as during Heinrich events
(Menviel et al. 2014, 2020) or saddle collapse (Gregoire
et al. 2012; Ivanovic et al. 2017), has acted as a driver for
abrupt environmental change. This makes the case for fur-
ther coupled-climate ice sheet experiments (e.g., Alvarez-
Solas et al. 2019; Schannwell et al. 2024; Mikolajewicz et al.
2024 in order to include prognostic feedbacks from the ice
sheets and better constrain the freshwater forcing during the
last glacial period and deglaciation.

**8 Conclusion**

We presented a mechanism explaining the millennial-scale
variability observed in a last glacial maximum HadCM3
general circulation model simulation, forced with a con-
stant pattern of deglacial meltwater. The simulation dis-
plays a periodic evolution along a limit cycle through a cold
phase—where the AMOC is weak and there is no deep water
formation in the North Atlantic, a meridional warm phase—
where the AMOC is strong and convection predominates in
the Nordic seas, and a zonal warm phase—where the AMOC
is strong and convection shifts primarily to the Irminger Sea.

The convection–advection oscillator is composed of a
slow advection component and a fast convection compo-
nent, with the AMOC coupling the two. The slow advec-
tion component transports a salt anomaly between the

subtropical Atlantic and mainly the deep Pacific. The *fast
convection* component responds to thresholds in North
Atlantic stratification by reactivating or deactivating deep
water formation sites, amplified by the positive salt-advec-
tion feedback. During *cold* phases, the salt accumulates in
the subtropical Atlantic and slowly leaks into the North
Atlantic to increase its upper-ocean density. In addition,
progressive subsurface warming decreases the deep den-
sity in sea ice-covered regions. This leads to a decrease
in the stratification of the Nordic Seas until reaching a
threshold where convection rapidly reactivates. The result-
ing *meridional* phase acts as an overshoot mode before
deep water formation shifts to a predominantly *zonal* con-
figuration. During the *warm* phases, the Atlantic becomes
depleted in salt because of the slow circulation of the
subtropical salt anomaly along the deep ocean pathway.
In the North Atlantic, the subsurface heat accumulated
during *cold* phases is progressively being evacuated, while
the upper waters remain warm due to southwards sea ice
transport into the North Atlantic. The salt and tempera-
ture changes lead to a slow increase in the stratification of
the Irminger Sea and Iceland Basin. The *cold* mode ends
when deep convection reaches a threshold of deactivation,
leading to a rapid weakening of the AMOC driven by the
salt-advection feedback.

Starting from the AMOC control state, the model is
forced with the respective deglacial freshwater pattern at
the beginning of each simulation. This perturbation has to
be strong enough to deactivate the North Atlantic deep water
formation sites, but not too strong to allow for the crossing of
a threshold from a *cold* mode back into a *warm* mode within
the span of the global salt redistribution. In the absence of
salinity correction, the global salinity anomaly between the
Pacific and Subtropical Atlantic is superseded by a global
salinity drift, leading to an abnormally high upper North
Atlantic density and moving the system out of the window
of opportunity.

The simulations analysed in the paper demonstrate the
concept of a window of opportunity, associated with melt-
water discharge. Depending on the forcing, we observe a
*warm steady state* simulation where the system sustains
North Atlantic deep water formation despite the freshwater
input, a *limit cycle* simulation where the system follows the
convection–advection oscillator mechanism, and an *excit-
able cold steady state* simulation, where a persistence in the
weak AMOC state is interrupted by short-lived recoveries of
North Atlantic convection. Further sensitivity experiments
involving a more comprehensive range of boundary condi-
tions (such as atmospheric greenhouse gas concentrations,
orbital parameters, ice sheet geometry, salinity conservation
scheme), climate forcing and physical parameter settings are
needed to categorically constrain the regions in parameter
space where these different behaviours occur.
While model and boundary condition dependencies
influence the results in any study of simulated abrupt cli-
mate change, we are at an exciting juncture where com-
mon patterns between the different mechanisms presented
across the literature are starting to crystallise. This study
evidences how concepts of dynamical systems theory can
guide the analysis of model simulations towards a uni-
fied understanding of millennial-scale climate variability.
Finally, in the absence of dynamical ice sheets, salin-
ity conservation is a necessary condition for obtaining
millennial-scale oscillations in our simulations, and the
introduction of climate-ice sheet feedback would repre-
sent a significant step towards simulations of more realistic
millennial-scale variability.

**Supplementary Information** The online version contains supplementary material available at https://doi.org/10.1007/s00382-025-07630-x.

**Acknowledgements** YMR was supported by the Leeds-York-Hull Natural Environment Research Council (NERC) Doctoral Training Partnership (DTP) Panorama under Grant NE/S007458/1. LJG was funded by UKRI Future Leaders Fellowship under Grant agreement no. MR/S016961/1. DS was supported by the TipESM project funded by the European Union's Horizon Europe Climate, Energy and Mobility research and innovation programme under Grant agreement no. 101137673. RB gratefully acknowledges funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie Grant agreement no. 956170 (CriticalEarth). Simulations were run on the high-performance facilities ARC3 and ARC4 at the University of Leeds. Thanks are due to Dr Edward Armstrong from the University of Helsinki (previously at the University of Bristol) and Prof. Paul Valdes from the University of Bristol for providing the salinity tendencies scripts and helping with their interpretations. Thanks to Josselin Romé for his help on the figure design.

**Author Contributions** The authors are listed as per their estimated
contribution to this article. YMR wrote the manuscript, ran the simu-
lations, conducted the analysis and produced the figures. RI helped run-
ning the simulations, supervised the analysis and proofread the article.
LJG supervised the analysis and proofread the article. DS provided
the oceanic expertise in the definition of the mechanism and proofread
the article. SST helped with the description of the mechanism and
proofread the manuscript. RB helped with the theoretical framing of
the analysis and proofread the manuscript.

**Funding** Funding was provided by Natural Environment Research Council (NE/S007458/1), UK Research and Innovation (MR/ S016961/1), HORIZON EUROPE Climate, Energy and Mobility (101137673), HORIZON EUROPE Marie Skłodowska-Curie Actions (956170).

**Data availability** The data used to produce the Figures of this article are available on the CEDA archive under Romé et al. (2024). Additional data can be requested to the corresponding author.

**Code availability** The code used to produce the Figures of this article is available on Zenodo under Romé (2024).

**Declarations**

**Conflict of interest** The authors declare no conflict of interest or con-
flict of interest.

**Prior** A prior version of this manuscript was published in YMR's PhD thesis and can be accessed at https://etheses.whiterose.ac.uk/34761/.

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

**References**

Alley RB, Anandakrishnan S, Jung P (2001) Stochastic resonance in the North Atlantic. Paleoceanography 16(2):190–198. https://doi.org/10.1029/2000PA000518

Alvarez-Solas J, Banderas R, Robinson A et al (2019) Ocean-driven millennial-scale variability of the Eurasian ice sheet during the last glacial period simulated with a hybrid ice-sheet-shelf model. Clim Past 15(3):957–979. https://doi.org/10.5194/cp-15-957-2019

Armstrong E, Izumi K, Valdes P (2022) Identifying the mechanisms of DO-scale oscillations in a GCM: a salt oscillator triggered by the Laurentide ice sheet. Clim Dyn. https://doi.org/10.1007/s00382-022-06564-y

Banderas R, Alvarez-Solas J, Robinson A et al (2015) An interhemispheric mechanism for glacial abrupt climate change. Clim Dyn 44(9):2897–2908. https://doi.org/10.1007/s00382-014-2211-8

Barker S, Knorr G (2021) Millennial scale feedbacks determine the shape and rapidity of glacial termination. Nat Commun 12(1):2273. https://doi.org/10.1038/s41467-021-22388-6

Batchelor CL, Margold M, Krapp M et al (2019) The configuration of Northern Hemisphere ice sheets through the Quaternary. Nat Commun 10(1):1–10. https://doi.org/10.1038/s41467-019-11601-2

Bellomo K, Angeloni M, Corti S et al (2021) Future climate change shaped by inter-model differences in Atlantic meridional overturning circulation response. Nat Commun 12(1):3659. https://doi.org/10.1038/s41467-021-24015-w

Bethke I, Li C, Nisancioglu KH (2012) Can we use ice sheet reconstructions to constrain meltwater for deglacial simulations? Pale-oceanography. https://doi.org/10.1029/2011PA002258

Birchfield EG, Wang H, Rich JJ (1994) Century/millennium internal climate oscillations in an ocean-atmosphere-continental ice sheet model. J Geophys Res Oceans 99(C6):12459–12470. https://doi.org/10.1029/94JC00523. (important-abstract)

Boers N (2021) Observation-based early-warning signals for a col-
lapse of the Atlantic Meridional Overturning Circulation.
Nat Clim Change 11(8):680–688. https://doi.org/10.1038/
s41558-021-01097-4

Boers N, Ghil M, Stocker TF (2022) Theoretical and paleoclimatic evidence for abrupt transitions in the Earth system. Environ Res Lett 17(9):093006. https://doi.org/10.1088/1748-9326/ac8944

Böhm E, Lippold J, Gutjahr M et al (2015) Strong and deep Atlantic meridional overturning circulation during the last glacial cycle. Nature 517(7532):73–76. https://doi.org/10.1038/nature14059
Bond G, Broecker W, Johnsen S et al (1993) Correlations between climate records from North Atlantic sediments and Greenland ice. Nature 365(6442):143-147. https://doi.org/10.1038/365143a0

Briggs RD, Pollard D, Tarasov L (2014) A data-constrained large ensemble analysis of Antarctic evolution since the Eemian. Quat Sci Rev 103:91-115. https://doi.org/10.1016/j.quascirev.2014.09.003

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

Broecker WS, Bond G, Klas M et al (1990) A salt oscillator in the glacial Atlantic? 1. The concept. Paleooceanography 5(4):469-477. https://doi.org/10.1029/PA005i004p00469

Brown N, Galbraith E (2016) Hosed vs. unhosed: interruptions of the Atlantic Meridional Overturning Circulation in a global coupled model, with and without freshwater forcing. Clim Past 12:1663-1679. https://doi.org/10.5194/cp-12-1663-2016

Bryan K, Cox MD (1972) An approximate equation of state for numerical models of ocean circulation. J Phys Oceanogr 2(4):510-514. https://doi.org/10.1175/1520-0485(1972)002<0510:AAEOSF>2.0.CO;2

Buizert C, Schmittner A (2015) Southern Ocean control of glacial AMOC stability and Dansgaard-Oeschger interstadial duration. Paleooceanography 30(12):1595-1612. https://doi.org/10.1002/2015PA002795

Buizert C, Gkinis V, Severinghaus JP et al (2014) Greenland temperature response to climate forcing during the last deglaciation. Science 345(6201):1177-1180. https://doi.org/10.1126/science.1254961

Caesar L, Rahmstorf S, Robinson A et al (2018) Observed fingerprint of a weakening Atlantic Ocean overturning circulation. Nature 556(7700):191-196. https://doi.org/10.1038/s41586-018-0006-5

Cessi P (1996) Convective adjustment and thermohaline excitability. J Phys Oceanogr 26(4):481-491. https://doi.org/10.1175/1520-0485(1996)026<0481:CAATE>2.0.CO;2

Cheng J, Liu Z, He F et al (2011) Simulated two-stage recovery of Atlantic meridional overturning circulation during the last deglaciation. Abrupt Clim Change Mech Patterns Impacts 193:75. https://doi.org/10.1029/2010gm001014

Cimatoribus AA, Drijfhout SS, Livina V et al (2013) Dansgaard-Oeschger events: bifurcation points in the climate system. Clim Past 9(1):323-333. https://doi.org/10.5194/cp-9-323-2013

Clark PU, Pisias NG, Stocker TF et al (2002) The role of the thermohaline circulation in abrupt climate change. Nature 415(6874):863-869. https://doi.org/10.1038/415863a

Colin de Verdière A (2007) A simple model of millennial oscillations of the thermohaline circulation. J Phys Oceanogr 37(5):1142-1155. https://doi.org/10.1175/JPO3056.1

Cox P (2001) Description of the TRIFFID dynamic global vegetation model. Hadley Centre technical note 24

Cox PM, Betts RA, Bunton CB et al (1999) The impact of new land surface physics on the GCM simulation of climate and climate sensitivity. Clim Dyn 15(3):183-203. https://doi.org/10.1007/s003820050276

Crowley TJ (1992) North Atlantic Deep Water cools the southern hemisphere. Paleooceanography 7(4):489-497. https://doi.org/10.1029/92PA01058

Crucifix M (2012) Oscillators and relaxation phenomena in Pleistocene climate theory. Philos Trans R Soc A Math Phys Eng Sci 370(1962):1140-1165. https://doi.org/10.1098/rsta.2011.0315

Czaja A (2009) Atmospheric Control on the Thermohaline Circulation. J Phys Oceanogr 39(1):234-247. https://doi.org/10.1175/2008JPO3897.1

Dansgaard W, Johnsen SJ, Clausen HB et al (1993) Evidence for general instability of past climate from a 250-kyr ice-core record. Nature 364(6434):218-220. https://doi.org/10.1038/364218a0

Davies-Barnard T, Ridgwell A, Singarayer J et al (2017) Quantifying the influence of the terrestrial biosphere on glacial-interglacial climate dynamics. Clim Past 13(10):1381-1401. https://doi.org/10.5194/cp-13-1381-2017

Demory ME, Vidale PL, Roberts MJ et al (2014) The role of horizontal resolution in simulating drivers of the global hydrological cycle. Clim Dyn 42(7):2201-2225. https://doi.org/10.1007/s00382-013-1924-4

Dentith JE, Ivanovic RF, Gregoire LJ et al (2019) Ocean circulation drifts in multi-millennial climate simulations: the role of salinity corrections and climate feedbacks. Clim Dyn 52(3):1761-1781. https://doi.org/10.1007/s00382-018-4243-y

Ditlevsen P, Ditlevsen S (2023) Warning of a forthcoming collapse of the Atlantic meridional overturning circulation. Nat Commun 14(1):4254. https://doi.org/10.1038/s41467-023-39810-w

Dokken TM, Nisancioglu KH, Li C et al (2013) Dansgaard-Oeschger cycles: interactions between ocean and sea ice intrinsic to the Nordic seas. Paleooceanography 28(3):491-502. https://doi.org/10.1002/palo.20042

Drijfhout SS, Weber SL, van der Waluw E (2011) The stability of the MOC as diagnosed from model projections for pre-industrial, present and future climates. Clim Dyn 37(7):1575-1586. https://doi.org/10.1007/s00382-010-0930-z

Drijfhout S, Bathiany S, Beaulieu C et al (2015) Catalogue of abrupt shifts in intergovernmental panel on climate change climate models. Proc Natl Acad Sci 112(43):E5777-E5786. https://doi.org/10.1073/pnas.1511451112

Fairbanks RG (1989) A 17,000-year glacio-eustatic sea level record: influence of glacial melting rates on the Younger Dryas event and deep-ocean circulation. Nature 342:637-642. https://doi.org/10.1038/342637a0

Fofonoff NP (1985) Physical properties of seawater: a new salinity scale and equation of state for seawater. J Geophys Res Oceans 90(C2):3332-3342. https://doi.org/10.1029/JC090iC02p03332

Fofonoff NP, Millard RC Jr (1983) Algorithms for the computation of fundamental properties of seawater. UNESCO Tech Pap Mar Sci. https://doi.org/10.25607/OBP-1450

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

Gordon C, Cooper C, Senior CA et al (2000) The simulation of SST, sea ice extents and ocean heat transports in a version of the Hadley Centre coupled model without flux adjustments. Clim Dyn 16(2):147-168. https://doi.org/10.1007/s003820050010

Gottschalk J, Skinner LC, Misra S et al (2015) Abrupt changes in the southern extent of North Atlantic Deep Water during Dansgaard-Oeschger events. Nat Geosci 8(12):950-954. https://doi.org/10.1038/ngeo2558

Gregoire LJ, Payne AJ, Valdes PJ (2012) Deglacial rapid sea level rises caused by ice-sheet saddle collapses. Nature 487(7406):219-222. https://doi.org/10.1038/nature11257

Hawkins E, Smith RS, Allison LC et al (2011) Bistability of the Atlantic overturning circulation in a global climate model and links to ocean freshwater transport. Geophys Res Lett. https://doi.org/10.1029/2011GL047208

Henry LG, McManus JF, Curry WB et al (2016) North Atlantic ocean circulation and abrupt climate change during the last deglaciation. Science 353(6298):470-474. https://doi.org/10.1126/science.aaf5529

Huber C, Leuenberger M, Spahni R et al (2006) Isotope calibrated Greenland temperature record over Marine Isotope Stage 3 and its relation to CH4. Earth Planet Sci Lett 243(3):504-519. https://doi.org/10.1016/j.epsl.2006.01.002

Hughes ALC, Gyllencreutz R, Lohne ØS et al (2016) The last Eurasian ice sheets—a chronological database and time-slice
reconstruction, DATED-1. Boreas 45(1):1–45. https://doi.org/
10.1111/bor.12142

Huisman SE, Toom Md, Dijkstra HA et al (2010) An indicator of the
multiple equilibria regime of the Atlantic meridional overturn-
ing circulation. J Phys Oceanogr 40(3):551–567. https://doi.
org/10.1175/2009JPO4215.1

Ivanovic RF, Valdes PJ, Flecker R et al (2013) The parameterisation of Mediterranean-Atlantic water exchange in the Hadley Centre model HadCM3, and its effect on modelled North Atlantic climate. Ocean Model 62:11–16. https://doi.org/10.1016/j.ocemod.2012.11.002

Ivanovic RF, Valdes PJ, Gregoire L et al (2014) Sensitivity of modern climate to the presence, strength and salinity of Mediterranean-Atlantic exchange in a global general circulation model. Clim Dyn 42(3):859–877. https://doi.org/10.1007/s00382-013-1680-5

Ivanovic RF, Gregoire LJ, Kageyama M et al (2016) Transient climate simulations of the deglaciation 21–9 thousand years before present (version 1) - PMIP4 Core experiment design and boundary conditions. Geosci Model Dev 9(7):2563–2587. https://doi.org/10.5194/gmd-9-2563-2016

Ivanovic RF, Gregoire LJ, Wickert AD et al (2017) Collapse of the North American ice saddle 14,500 years ago caused widespread cooling and reduced ocean overturning circulation. Geophys Res Lett 44(1):383–392. https://doi.org/10.1002/2016GL071849. (important - Read)

Ivanovic RF, Gregoire LJ, Burke A et al (2018) Acceleration of North-
ern Ice sheet melt induces AMOC Slowdown and northern cool-
ing in simulations of the early last deglaciation. Paleooceanogr
Paleoclimatol 33(7):807–824. https://doi.org/10.1029/2017P
A003308

Jackson LC, Biastoch A, Buckley MW et al (2022) The evolution of the North Atlantic Meridional Overturning Circulation since 1980. Nat Rev Earth Environ 3(4):241–254. https://doi.org/10.1038/s43017-022-00263-2

Kageyama M, Paul A, Roche DM et al (2010) Modelling glacial climatic millennial-scale variability related to changes in the Atlantic meridional overturning circulation: a review. Quat Sci Rev 29(21):2931–2956. https://doi.org/10.1016/j.quascirev.2010.05.029

Kageyama M, Merkel U, Otto-Bliesner B et al (2013) Climatic impacts of fresh water hosing under Last Glacial Maximum conditions: a multi-model study. Clim Past 9(2):935–953. https://doi.org/10.5194/cp-9-935-2013

Kageyama M, Albani S, Braconnot P et al (2017) The PMIP4 contribu-
tion to CMIP6 - Part 4: Scientific objectives and experimental
design of the PMIP4-CMIP6 Last Glacial Maximum experi-
ments and PMIP4 sensitivity experiments. Geosci Model Dev
10(11):4035–4055. https://doi.org/10.5194/gmd-10-4035-2017

Kindler P, Guillevic M, Baumgartner M et al (2014) Temperature reconstruction from 10 to 120 kyr b2k from the NGRIP ice core. Clim Past 10(2):887–902. https://doi.org/10.5194/cp-10-887-2014

Klockmann M, Mikolajewicz U, Marotzke J (2018) Two AMOC states in response to decreasing greenhouse gas concentrations in the coupled climate model MPI-ESM. J Clim 31(19):7969–7984. https://doi.org/10.1175/JCLI-D-17-0859.1

Klockmann M, Mikolajewicz U, Kleppin H et al (2020) Coupling of the subpolar gyre and the overturning circulation during abrupt glacial climate transitions. Geophys Res Lett 47(21):e2020GL090361. https://doi.org/10.1029/2020GL090361

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

Kuniyoshi Y, Abe-Ouchi A, Sherriff-Tadano S et al (2022) Effect of climatic precession on Dansgaard-Oeschger-like oscillations.

Geophys Res Lett 49(6):e2021GL095695. https://doi.org/10.
1029/2021GL095695

Lambeck K, Rouby H, Purcell A et al (2014) Sea level and global ice volumes from the last glacial maximum to the Holocene. Proc Natl Acad Sci 111(43):15296–15303. https://doi.org/10.1073/pnas.1411762111

Li C, Born A (2019) Coupled atmosphere-ice-ocean dynamics in Dansgaard-Oeschger events. Quat Sci Rev 203:1–20. https://doi.org/10.1016/j.quascirev.2018.10.031

Li G, Cheng L, Zhu J et al (2020) Increasing ocean stratification over the past half-century. Nat Clim Change 10(12):1116–1123. https://doi.org/10.1038/s41558-020-00918-2

Liu W, Xie SP, Liu Z et al (2017) Overlooked possibility of a col-
lapsed Atlantic meridional overturning circulation in warming
climate. Sci Adv 3(1):e160166. https://doi.org/10.1126/sciadv.
1601666. (high)

Lobelle D, Beaulieu C, Livina V et al (2020) Detectability of an AMOC decline in current and projected climate changes. Geophys Res Lett 47(20):e2020GL089974. https://doi.org/10.1029/2020GL089974

Lohmann J, Ditlevsen PD (2019) Objective extraction and analysis of statistical features of Dansgaard-Oeschger events. Clim Past 15(5):1771–1792. https://doi.org/10.5194/cp-15-1771-2019

Lynch-Stieglitz J (2017) The Atlantic meridional overturning circula-
tion and abrupt climate change. Annu Rev Mar Sci 9(1):83–104.
https://doi.org/10.1146/annurev-marine-010816-060415

Malmierca-Vallet I, Sime LC, the D-O Community Members (2023) Dansgaard-Oeschger events in climate models: review and baseline Marine Isotope Stage 3 (MIS3) protocol. Clim Past 19(5):915–942. https://doi.org/10.5194/cp-19-915-2023

Malmierca-Vallet I, Sime LC, Valdes PJ et al (2024) The impact of CO2 and climate state on whether Dansgaard-Oeschger type oscillations occur in climate models. Geophys Res Lett 51(13):e2024GL110068. https://doi.org/10.1029/2024GL110068

Manabe S, Stouffer RJ (1995) Simulation of abrupt climate change induced by freshwater input to the North Atlantic Ocean. Nature 378(6553):165–167. https://doi.org/10.1038/378165a0

Marcott SA, Clark PU, Padman L et al (2011) Ice-shelf collapse from subsurface warming as a trigger for Heinrich events. Proc Natl Acad Sci USA 108(33):13415. https://doi.org/10.1073/PNAS.1104772108

Matero ISO, Gregoire LJ, Ivanovic RF et al (2017) The 8.2 ka cooling event caused by Laurentide ice saddle collapse. Earth Planet Sci Lett 473:205–214. https://doi.org/10.1016/j.epsl.2017.06.011

Mecking JV, Drijfhout SS, Jackson LC et al (2017) The effect of model bias on Atlantic freshwater transport and implications for AMOC bi-stability. Tellus A Dyn Meteorol Oceanogr 69(1):1299910. https://doi.org/10.1080/16000870.2017.1299910

Menviel L, Timmermann A, Friedrich T et al (2014) Hindcasting the continuum of Dansgaard-Oeschger variability: mechanisms, patterns and timing. Clim Past 10(1):63–77. https://doi.org/10.5194/cp-10-63-2014

Menviel LC, Skinner LC, Tarasov L et al (2020) An ice-climate oscillatory framework for Dansgaard-Oeschger cycles. Nat Rev Earth Environ 1(12):677–693. https://doi.org/10.1038/s43017-020-00106-y

Merz N, Raible CC, Woollings T (2015) North Atlantic eddy-driven jet in interglacial and glacial winter climates. J Clim 28(10):3977–3997. https://doi.org/10.1175/JCLI-D-14-00525.1

Michel SLL, Swingedouw D, Ortega P et al (2022) Early warning signal for a tipping point suggested by a millennial Atlantic Multidecadal Variability reconstruction. Nat Commun 13(1):5176. https://doi.org/10.1038/s41467-022-32704-3

Mikolajewicz U, Kapsch ML, Schannwell C et al (2024) Deglaciation and abrupt events in a coupled comprehensive
atmosphere-ocean-ice sheet-solid earth model. Clim Past Discuss. https://doi.org/10.5194/cp-2024-55

Montoya M, Born A, Levermann A (2011) Reversed North Atlantic gyre dynamics in present and glacial climates. Clim Dyn 36(5):1107–1118. https://doi.org/10.1007/s00382-009-0729-y

Nuber S, Rae JWB, Zhang X et al (2023) Indian Ocean salinity build-up primes deglacial ocean circulation recovery. Nature 617(7960):306–311. https://doi.org/10.1038/s41586-023-05866-3

Oka A, Hasumi H, Abe-Ouchi A (2012) The thermal threshold of the Atlantic meridional overturning circulation and its control by wind stress forcing during glacial climate. Geophys Res Lett. https://doi.org/10.1029/2012GL051421

Oka A, Abe-Ouchi A, Sherriff-Tadano S et al (2021) Glacial mode shift of the Atlantic meridional overturning circulation by warming over the Southern Ocean. Commun Earth Environ 2(1):1–8. https://doi.org/10.1038/s43247-021-00226-3

Peltier WR, Vettoretti G (2014) Dansgaard-Oeschger oscillations predicted in a comprehensive model of glacial climate: a “kicked” salt oscillator in the Atlantic. Geophys Res Lett 41(20):7306–7313. https://doi.org/10.1002/2014GL061413

Pikovsky A, Rosenblum M, Kurths J (2001) Synchronization: a universal concept in nonlinear sciences. Cambridge Nonlinear Science Series, Cambridge University Press, Cambridge

Pope VD, Gallani ML, Rowntree PR et al (2000) The impact of new physical parametrizations in the Hadley Centre climate model: HadAM3. Clim Dyn 16(2):123–146. https://doi.org/10.1007/s003820050009

Rahmstorf S (2002) Ocean circulation and climate during the past 120,000 years. Nature 419(6903):207–214. https://doi.org/10.1038/nature01090

Rasmussen SO, Bigler M, Blockley SP et al (2014) A stratigraphic framework for abrupt climatic changes during the Last Glacial period based on three synchronized Greenland ice-core records: refining and extending the INTIMATE event stratigraphy. Quat Sci Rev 106:14–28. https://doi.org/10.1016/j.quascirev.2014.09.007

Riechers K, Gottwald G, Boers N (2024) Glacial abrupt climate change as a multiscale phenomenon resulting from monostable excitable dynamics. J Clim 37(8):2741–2763. https://doi.org/10.1175/JCLI-D-23-0308.1

Roberts A, Saha R (2017) Relaxation oscillations in an idealized ocean circulation model. Clim Dyn 48(7):2123–2134. https://doi.org/10.1007/s00382-016-3195-3

Roberts WHG, Valdes PJ (2017) Green mountains and white plains: the effect of northern hemisphere ice sheets on the global energy budget. J Clim 30(10):3887–3905. https://doi.org/10.1175/JCLI-D-15-0846.1

Roberts WHG, Valdes PJ, Payne AJ (2014) A new constraint on the size of Heinrich Events from an iceberg/sediment model. Earth Planet Sci Lett 386:1–9. https://doi.org/10.1016/j.epsl.2013.10.020

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

Romanou A, Rind D, Jonas J et al (2023) Stochastic bifurcation of the north Atlantic circulation under a midrange future climate scenario with the NASA-GISS ModelE. J Clim 36(18):6141–6161. https://doi.org/10.1175/JCLI-D-22-0536.1

Romé Y (2025) Olnavy/ROME2025 cd mechanism: v.1.1 - Release for publication in Climate Dynamics. https://doi.org/10.5281/zenodo.14925907. https://zenodo.org/records/14925907. Accessed 25 Feb 2025

Romé YM, Ivanovic RF, Gregoire LJ et al (2022) Millennial-scale climate oscillations triggered by deglacial meltwater discharge in last glacial maximum simulations. Paleooceanogr Paleoclimatol 37(10):e2022PA004451. https://doi.org/10.1029/2022PA004451

Romé Y, Ivanovic R, Gregoire L (2024) The convection–advection oscillator identified in HadCM3 LGM simulations. https:// catalogueda.ac.uk/uuid/1f4ebb2944ec43a39ce6c69a8f1942fb/. Accessed 25 Feb 2025

Schannwell C, Mikolajewicz U, Kapsch ML et al (2024) A mechanism for reconciling the synchronisation of Heinrich events and Dansgaard-Oeschger cycles. Nat Commun 15(1):2961. https://doi.org/10.1038/s41467-024-47141-7

Schulz M, Paul A, Timmermann A (2002) Relaxation oscillators in concert: a framework for climate change at millennial timescales during the late Pleistocene. Geophys Res Lett 29(24):46–46.4. https://doi.org/10.1029/2002GL016144

Sherriff-Tadano S, Abe-Ouchi A, Yoshimori M et al (2018) Influence of glacial ice sheets on the Atlantic meridional overturning circulation through surface wind change. Clim Dyn 50(7):2881–2903. https://doi.org/10.1007/s00382-017-3780-0

Sherriff-Tadano S, Abe-Ouchi A, Oka A (2021) Impact of mid-glacial ice sheets on deep ocean circulation and global climate. Clim Past 17(1):95–110. https://doi.org/10.5194/cp-17-95-2021

Smith RS, Gregory JM (2009) A study of the sensitivity of ocean overturning circulation and climate to freshwater input in different regions of the North Atlantic. Geophys Res Lett. https://doi.org/10.1029/2009GL038607

Stocker TF (1998) The seesaw effect. Science 282(5386):61–62. https://doi.org/10.1126/science.282.5386.61

Stommel H (1961) Thermohaline convection with two stable regimes of flow. Tellus 13(2):224–230. https://doi.org/10.1111/j.2153-3490.1961.tb00079.x

Tarasov L, Peltier RW (2002) Greenland glacial history and local geodynamic consequences. Geophys J Int 150(1):198–229. https://doi.org/10.1046/j.1365-246X.2002.01702.x

Tarasov L, Dyke AS, Neal RM et al (2012) A data-calibrated distribution of deglacial chronologies for the North American ice complex from glaciological modeling. Earth Planet Sci Lett 315–316:30–40. https://doi.org/10.1016/j.epsl.2011.09.010

Thomas ER, Wolff EW, Mulvaney R et al (2009) Anatomy of a Dansgaard-Oeschger warming transition: high-resolution analysis of the North Greenland Ice Core Project ice core. J Geophys Res Atmos. https://doi.org/10.1029/2008JD011215

Thompson AF, Hines SK, Adkins JF (2019) A Southern Ocean mechanism for the interhemispheric coupling and phasing of the bipolar seesaw. J Clim 32(14):4347–4365. https://doi.org/10.1175/JCLI-D-18-0621.1

Thornalley DJR, Oppo DW, Ortega P et al (2018) Anomalously weak Labrador Sea convection and Atlantic overturning during the past 150 years. Nature 556(7700):227–230. https://doi.org/10.1038/s41586-018-0007-4

Timmermann A, Gildor H, Schulz M et al (2003) Coherent resonant millennial-scale climate oscillations triggered by massive meltwater pulses. J Clim 16(15):2569–2585. https://doi.org/10.1175/1520-0442(2003)016<2569:CRMCOT>2.0.CO;2

Valdes P (2011) Built for stability. Nat Geosci 4(7):414–416. https://doi.org/10.1038/ngeo1200

Valdes PJ, Armstrong E, Badger MPS et al (2017) The BRIDGE HadCM3 family of climate models: HadCM3@Bristol v1.0. Geosci Model Dev Discuss. https://doi.org/10.5194/gmd-2017-16

van Westen RM, Dijkstra HA (2023) Asymmetry of AMOC hysteresis in a state-of-the-art global climate model. Geophys Res Lett 50(22):e2023GL106088. https://doi.org/10.1029/2023GL106088
van Westen RM, Kliphuis M, Dijkstra HA (2024) Physics-based early warning signal shows that AMOC is on tipping course. *Sci Adv* 10(6):eadk1189. https://doi.org/10.1126/sciadv.adk1189

Vettoretti G, Peltier WR (2016) Thermohaline instability and the formation of glacial North Atlantic super polynyas at the onset of Dansgaard–Oeschger warming events. *Geophys Res Lett* 43(10):5336–5344. https://doi.org/10.1002/2016GL068891

Vettoretti G, Peltier WR (2018) Fast physics and slow physics in the nonlinear Dansgaard–Oeschger relaxation oscillation. *J Clim* 31(9):3423–3449. https://doi.org/10.1175/JCLI-D-17-0559.1

Vettoretti G, Ditlevsen P, Jochum M et al (2022) Atmospheric CO2 control of spontaneous millennial-scale ice age climate oscillations. *Nat Geosci* 15(4):300–306. https://doi.org/10.1038/s41561-022-00920-7

Weijer W, Cheng W, Drijfhout SS et al (2019) Stability of the Atlantic meridional overturning circulation: a review and synthesis. *J Geophys Res Oceans* 124(8):5336–5375. https://doi.org/10.1029/2019JC015083

Weijer W, Cheng W, Garuba OA et al (2020) CMIP6 models predict significant 21st century decline of the Atlantic meridional overturning circulation. *Geophys Res Lett* 47(12):e2019GL086075. https://doi.org/10.1029/2019GL086075

Weijer W, Dijkstra HA (2003) Multiple oscillatory modes of the global ocean circulation. *J Phys Oceanogr* 33(11):2197–2213. https://doi.org/10.1175/1520-0485(2003)033<2197:MOMOTG>2.0.CO;2

Welander P (1982) A simple heat-salt oscillator. *Dyn Atmos Oceans* 6(4):233–242. https://doi.org/10.1016/0377-0265(82)90030-6

Wickert AD, Williams C, Gregoire LJ et al (2023) Marine-calibrated chronology of southern laurentide ice sheet advance and retreat:

~2,000-year cycles paced by meltwater-climate feedback. *Geo-phys Res Lett* 50(10):e2022GL100391. https://doi.org/10.1029/2022GL100391

Winton M, Sarachik ES (1993) Thermohaline Oscillations Induced by Strong Steady Salinity Forcing of Ocean General Circulation Models. *J Phys Oceanogr* 23(7):1389–1410. https://doi.org/10.1175/1520-0485(1993)023<1389:TOIBSS>2.0.CO;2

Wolff EW, Chappellaz J, Blunier T et al (2010) Millennial-scale variability during the last glacial: The ice core record. *Quat Sci Rev* 29(21):2828–2838. https://doi.org/10.1016/j.quascirev.2009.10.013

Zhang X, Lohmann G, Knorr G et al (2014) Abrupt glacial climate shifts controlled by ice sheet changes. *Nature* 512(7514):290–294. https://doi.org/10.1038/nature13592

Zhang X, Knorr G, Lohmann G et al (2017) Abrupt North Atlantic circulation changes in response to gradual CO₂ forcing in a glacial climate state. *Nat Geosci* 10(7):518–523. https://doi.org/10.1038/ngeo2974

Zhang X, Barker S, Knorr G et al (2021) Direct astronomical influence on abrupt climate variability. *Nat Geosci* 14(11):819–826. https://doi.org/10.1038/s41561-021-00846-6

**Publisher's Note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.