Identifying the Mechanisms of DO-scale Oscillations in a Gcm: a Salt Oscillator Triggered by the Laurentide Ice Sheet

The driver mechanisms of Dansgaard-Oeschger (DO) events remain uncertain, in part because many climate models do not show similar oscillatory behaviour. Here we present results from glacial simulations of the HadCM3B coupled atmosphere-ocean-vegetation model that show stochastic, quasi-periodical DO-scale variability. This variability is driven by variations in the strength of the Atlantic Meridional Overturning Circulation in response to North Atlantic salinity uctuations. The mechanism represents a salt oscillator driven by the salinity gradient between the subtropical gyre and the Northern North Atlantic. Utilising a full set of model salinity diagnostics, we identify a complex ocean-atmosphere-sea-ice feedback mechanism that maintains this oscillator, driven by the interplay between surface freshwater uxes (tropical P-E balance and sea-ice), advection, and convection. The key trigger is the extent of the Laurentide ice sheet, which alters atmospheric and ocean circulation patterns, highlighting the sensitivity of the climate system to land-ice extent. This, in addition to the background climate state, pushes the climate beyond a tipping point and into an oscillatory mode on a timescale comparable to the DO events.


Introduction
The mechanisms that drive the Dansgaard-Oeschger (DO) events remain uncertain. The prevailing theory is that they are a response to changes in the strength of the Atlantic Meridional Overturning Circulation (AMOC; Broecker et al., 1985), which itself is strongly linked to salinity uctuations in the North Atlantic background conditions (see Li and Born, 2019). A key feature of the glacial climate was the presence of large continental ice sheets, enhanced wind strength and a more intense, less wobbly and shifted jet stream, which may have permitted the climate to support such oscillations and abrupt events (Armstrong et al., 2019a;Cook and Held, 1988; Li and Born, 2019; Merz et al., 2015;Ullman et al., 2014). Indeed, Wunsch (2006) concluded that DO events were reliant on the presence of continental ice sheets and their in uence on wind patterns. Similarly, Zhang et al. (2014) showed that only minor changes in the height of Laurentide ice sheet triggered DO events in the COSMOS model. Ice sheet height shifted wind patterns, which altered the atmosphere-ocean-ice system and impacted sea-ice/freshwater export from the Arctic and advection of salt. This increased AMOC strength and initiated an oscillation.
Problematically, understanding the drivers of DO events is hindered because many complex climate models do not simulate similar oscillations without additional forcing from freshwater hosing. This issue may be due to models lacking the correct sensitivity to boundary conditions rather than the necessary physics (Li and Born, 2019;Valdes, 2011). Saying that, several unforced modelling studies have demonstrated similar variability, however these were run with pre-industrial (PI) boundary conditions and events are typically short-lived and weak (Drijfhout et (Brown and Galbraith, 2016). They used a similar con guration of PI ice sheets and a glacial CO 2 concentration (180ppm), and the advection of salt between the STG and SPG was also shown to play an important role.
To our knowledge, the only fully glacial simulation that exhibits persistent and systematic variability has been performed using CESM1 (Peltier and Vettoretti, 2014;Vettoretti and Peltier, 2015;Vettoretti and Peltier, 2016;Vettoretti and Peltier, 2018). Here DO-type oscillations were linked to the 'salt oscillator' hypothesis, with transitions between a stadial and interstadial driven by instability in the water column within the Irminger Sea. Here, we add an additional model to CESM1 with results from the HadCM3B model , which demonstrates DO-like unforced millennial-scale oscillations in fully glacial simulations. We give a detailed explanation of the coupled ocean-atmosphere-sea-ice feedbacks that drive this oscillation utilising a full set of salinity diagnostics. The aim of this study is to 1) identify the feedbacks that maintain the oscillation and initiate the stadial/interstadial phases, and 2) identify what pushes the glacial climate system into an oscillatory state.

The HadCM3B Model
The Hadley Centre Coupled Model 3 Bristol (HadCM3B) is a coupled climate model and version of the more commonly known HadCM3, which was originally developed at the Hadley Centre/UK Meteorological O ce and has been adapted at the University of Bristol. The speci c version of the model we use is HadCM3B-M2.1aD which is outlined in detail in Valdes et al. (2017) and is only brie y discussed here. HadCM3B been shown to produce an accurate representation of climate and remains competitive with other more modern models used in CMIP5.
The atmosphere model (Pope et al., 2000) has a resolution of 3.75° × 2.75° and 19 vertical levels with a 30-minute timestep. It utilises a number of boundary conditions, including the land-sea mask, orography, and a number of soil and vegetation parameters. The ocean model (Gordon et al., 2000) has a resolution of 1.25° × 1.25°, 20 vertical levels and a 1-hour timestep. It utilises a rigid lid approach, so there is no variation in ocean volume. HadCM3B incorporates the MOSES 2 land surface scheme (Cox et al., 1999) which simulates physiological processes and uxes of energy and water. It also incorporates the fractional coverage of nine different surface types, simulated by the dynamic global vegetation model (DGVM) TRIFFID.
The model does not include an interactive land-ice model or carbon and methane cycle so these boundary conditions have been imposed (see below). The river runoff is instantaneously transferred to out ow grid points following a prescribed runoff map. Sea ice is simulated using a zero layer model (Semtner, 1976) calculated on top of the ocean with movement controlled by upper-ocean currents. Ice is formed predominantly in leads and by snowfall, with continuous removal from the base through the year and from the surface during summer melts. Sea-ice salinity is assumed to be constant with the ux into the ocean depending on ice formation or melt.

Experimental Set-up
We performed a set of glacial simulations with 30kyr BP (before present) boundary conditions (Table 1), which have previously been shown to exhibit millennial-scale variability in Armstrong et al. (2019b). The boundary conditions include prescribed orbital parameters which are very well constrained and taken from Berger et al. (1998). The concentration of atmospheric CO 2 is from the Vostok Ice core (Loulergue et al., 2008;Petit et al., 1999), whilst N 2 O and CH 4 concentrations are from the EPICA Dome C ice core (Spahni et al., 2005).
HadCM3B does not have an ice-sheet model, so the elevation and extent of the Greenland, Laurentide, Fennoscandian and Antarctic ice sheets have been prescribed. For this, we utilise the ICE-5G model (Peltier, 2004) which simulates the evolution of ice thickness, extent, and isostatic rebound. This is used in HadCM3B to calculate orographic height, bathymetry, ice fraction, and the land-sea mask. However, ICE-5G only simulates ice up to the last glacial maximum (LGM; 21kyr BP), so we calculate the pre-LGM ice sheet from the equivalent ice volume/sea level during the deglaciation. For instance, the sea level depression at 30 kyr BP is compared to the sea level during the deglaciation (since the LGM), where they are the same the ice extent is inferred to be the same as that simulated by Ice-5G but for 30kyr BP. This is the same approach used in Armstrong et al. (2019b) and Davies-Barnard et al. (2017). There are limitations as ice sheets show different structures during decay and growth phases, but it provides a better approximation of the ice area than other trialled methodologies, such as in Singarayer and Valdes (2010).
The boundary conditions have been incorporated into a 6000-year base simulation, labelled 30kyrBASE. A second 4000-year simulation, labelled 30kyrST, was carried out which outputs a full set of salinity diagnostics. This reveals the relative importance of the salinity tendencies in driving overall salinity change, which is crucial in understanding AMOC strength (see Section 4). The key tendencies that are relevant to this study are; advection, diffusion, Gent and McWilliams isopycnal diffusion (Gent and Mcwilliams, 1990; which is a parameterisation of isopycnal eddies), convection, surface uxes which includes the P-E balance and river runoff, sea-ice processes which is associated with brine rejection from sea-ice formation and meltwater, and mixed layer physics. The tendencies are plotted as a rate of change of salinity (PSU/yr), and when combined, are equivalent to the change in overall salinity between each year. These are discussed in detail in Section 4. Both of these base simulations have been initiated from a spin up experiment that has been run for more than 6000 years and therefore very long term trends are small.
This 30kyr BP set-up was originally performed as part of an ensemble study of simulations covering 0 to 60kyr BP outlined in Armstrong et al. (2019b). The 30kyr BP experiment from this study exhibits millennial variability, however a 28kyr BP experiment does not show comparable variability despite similar boundary conditions (Table 1). To investigate what pushes the climate system into an oscillatory state at 30kyr BP, we performed a range of additional 30kyr BP experiments but with perturbed boundary conditions, altered to that used for the 28kyr experiment (Table 2, Figure 1). These are 3000-year simulations with 28kyr; orbital parameters (30kyrORB), greenhouse gases (30kyrGHG) and ocean bathymetry (30kyrBATHY; Figure 1a).
To investigate the role of ice-sheet con guration, we performed a further ve simulations with 28kyr BP boundary conditions but with varying 30kyr BP ice sheet con guration (regional height and extent) in order to test if this initiates an oscillation. The difference in fractional extent and height of the ice sheets between the 28kyr and 30kyr BP experiments (Figure 1b,c) re ects the ice-sheet methodology that we use. We ran ve 3000-year simulations with 28kyr BP boundary conditions with 1) a 30kyr Fennoscandian ice sheet (28kyrFEN), 2) Fennoscandian and east Laurentide ice sheets (28kyrELAU), 3) Fennoscandian, east and west Laurentide ice sheets (28kyrEWLAU), 4) Fennoscandian and the entire Laurentide ice sheets (28kyrLAU), and nally 5) 50% of all the ice sheets (28kyr50). The speci c regions are outlined in Figure  1c. For all of the perturbed boundary simulations, the analysis excludes the initial 1000 years of the simulation to permit a model spin-up.

The Oscillation In Hadcm3b
The timeseries for a range of atmospheric and oceanic variables (Figure 2), including the AMOC index (mean AMOC strength between 40° and 50°N at 800m, which re ects peak AMOC ow in HadCM3B), for the 30kyrBASE simulation demonstrates an oscillation with a dominant timescale of approximately 1500 years. The oscillation frequency is compatible with the DO events in the ice core record, which have a characteristic periodicity of between 1000 and 5000 years (Dansgaard et al., 1993), and similarly shows longer stadials than interstadials. The DO events are commonly linked to oscillating AMOC strength (Broecker et al., 1985). The average AMOC strength (at 30°N) at 30kyr is 13.12 ± 2.47 Sv, which compares  Klockmann (2020) showed the opposite to be the case in the MPI-ESM. The barotropic streamfunction averaged over the SPG gyre (not shown) indicates that the SPG strength is positively correlated to AMOC strength in HadCM3B. This overall increase in gyre strength during interstadials re ects the overall increase in transport in the NA. This is discussed further in Section 4.
The simulated oscillation does not exhibit a clear sawtooth pattern that is shown for some of the DO events in the ice core record, characterised by rapid warming followed by slower cooling (Seager and Battisti, 2007). Instead it shows a more symmetric rapid warming/rapid cooling shape. This lack of abrupt change may be due to inadequate representation of physical processes, such as the response time of sea ice and/or sensitivity to meltwater. Saying this, some of the shorter DO events in observations are more symmetric (Buizert and Schmittner, 2015) which are more consistent with the results presented here.

The Oscillator Mechanism
In order to understand the mechanism that drives the oscillation we need to determine what governs the strength of the AMOC. Therefore we need to understand the movement of salt and how this varies through time. The 30kyrST ( Table 2) simulation includes a full set of salinity tendency diagnostics (PSU/yr); positive/negative values indicate where salt is added/removed via a process. Figure 4 shows the composite spatial patterns of the tendencies during AMax (interstadials), AMin (stadials) and the anomaly. To highlight key sites of spatial variability, Figure 5 shows the rst EOF and corresponding PC1 timeseries for each tendency (additional EOFs do not further clarify the mechanism and are not shown).
Within the NA, advection is the dominant tendency that transports salt northwards along the NAC from the saline tropics to the subpolar regions, driven by the salinity gradient. The total (vertical and horizontal) advective tendency is negative in the tropics as salt is removed via advection, and broadly positive at higher latitudes, particularly the eastern SPG, re ecting the direction of the NAC (Figure 4a-b).
The key area of variability is the eastern edge of the NAC and the subpolar gyre, with interstadials showing an increase in horizontal advection in the Nordic Seas and a decrease in the SPG (Figure 5a-b).
The diffusive tendencies include horizontal and vertical diffusion and the Gent McWilliams scheme (Gent and Mcwilliams, 1990), which parameterises isopycnal eddy mixing (Figure 4d-e, Figure 5c-d). These broadly counter advection as they remove the salinity gradient developed via advective processes. As such, variability in diffusion largely opposes advection. The combined advection and diffusion tendencies (Figure 4f) highlight the changing role of the SPG and GIN Seas. The composite anomaly indicates that during interstadials these tendencies are more positive in the GIN Seas and negative across the SPG, particularly the Eastern SPG. This negative anomaly indicates contraction of the SPG and a broader in ow of salty water along its eastern ank, with increased transport of salt into the GIN Seas.
Convection (Figure 4c, Figure 5e) is a comparatively small tendency and is primarily present in surface (0-300m) waters, yet it is crucial in forming deep water and driving overturning. During stadials, convection is limited to the Iceland Basin. During interstadials, convection increases across the Norwegian Sea, Iceland-Scotland ridge, the Irminger Sea and across the central SPG. This indicates a widespread increase in deep-water formation.
Surface tendencies are comprised of sea-ice processes and freshwater uxes (Figure 4g-h, Figure 5f-g) associated with the P-E balance. Note that the model does not have prognostic ice sheet calving. The tropics are more saline primarily due to the negative precipitation-evaporation (P-E) balance which adds salt; this is key in maintaining the STG/SPG salinity gradient. In the tropics, variability in the P-E balance re ects the changing position of the ITCZ. Interstadials show a negative anomaly due to increased freshening from a northward shifted ITCZ and vice versa. This tendency is largely positive at higher latitudes during interstadials due to increased evaporation associated with warmer temperatures.
The impact of sea ice on salinity is complex and its role differs depending on location. In the Arctic and eastern GIN Seas, sea-ice formation removes freshwater and so increases salinity in these regions. However, in areas of the GIN Seas and SPG sea ice freshens during both stadials and interstadials (Figure 4h). This is due to the seasonal nature of sea-ice; with build up during the winter from freezing of sea water and accumulation from precipitation, which melts in the summer ushing out both frozen sea water and accumulated freshwater. The net impact is a freshening effect. The role of sea ice diminishes across the SPG during interstadials due to retreat of the sea-ice edge. However, in the Nordic and Irminger seas, the tendency becomes more negative due to a decrease in sea ice cover and consequent release of freshwater. The changing regional roles of the tendencies are key in driving the oscillation. In order to identify what triggers the transitions between stadials and interstadials, Figure 6 shows the timeseries for regional tendencies and different climate variables. Only those tendencies important to each region are shown.

Stadial period
During stadials, temperatures are cold, sea ice extent is widespread, and the NAC and AMOC are weak.
Salinity and ocean temperature in the northern NA are at a minimum due to reduced transport of warm, saline waters from the tropics. Advection is at a minimum in the Nordic Seas but remains high in the Irminger Sea and across the SPG, indicating that advection is predominantly into this region during stadials as indicated by an enhanced Irminger current ( Figure 3). The transport of salt is too small to maintain substantial deep-water formation and the increased density gradient indicates a more strati ed water column and so convection is at a minimum.
In the Nordic Seas, colder conditions accelerate sea-ice formation, which adds salt and acts as a negative feedback, countering the decline in advection. In the SPG and Irminger Sea, the seasonal nature of sea ice means it freshens these regions. Across the south and central SPG, sea ice acts to freshen more during the stadial as sea ice is more widespread. However, in the EGC and Irminger Sea its role is less negative as colder conditions reduce seasonal sea ice uctuations.
The tropics exhibit an opposing salinity anomaly due to a decrease in advection out of the tropics and an increase in surface salinity ux. The latter re ects the southward shift of the ITCZ due to the increased meridional temperature gradient. This enhances the negative P-E balance (Figure 3d) and so the tendency is more positive. Due to reduced advection, salt is kept within the STG and advected into intermediate depths, where in contrast to surface waters, subtropical gyre circulation has strengthened (not shown). These factors cause an accumulation of salt in the tropics and an increase in the salinity gradient between the two regions during the stadial (Figure 7). The build-up and collapse of this salinity gradient is crucial in maintaining the oscillation. Note that across the Northern North Atlantic, salinity and AMOC strength increases steadily throughout the stadial (Figure 7). In the Nordic Seas, sea ice growth may play an important role in driving these increases.

Transition into an Interstadial
The transition into an interstadial is characterised by an increase in both AMOC strength and NA salinity/temperature. Advection decreases within the central SPG and Irminger Seas and increases into the Nordic Seas, indicating a shift in circulation. This shift may be a response to increasing temperatures, so water needs to move further north to lose enough heat to sink. Convection increases in all these regions however, despite reduced advection. In the Nordic Seas sea ice volume decreases to zero, this freshens the region, which acts as a negative feedback to increased advection.
In the Nordic Seas the increase in convection and temperature and decrease in sea ice cover and the density gradient, lead AMOC change by approximately 50 years (Figure 8) indicating this region is key in triggering the interstadial. We do not see a similar lead for any of the other regions. The driver of this may be a wind-driven positive feedback in response to a decrease in sea-ice cover and an increase in sea-level pressure (SLP) across the Nordic Sea ( Figure 6a). As regional salinity increases, circulation and temperatures increase, and sea ice cover decreases. Consequently, a negative SLP anomaly develops (Figure 6a), which occurs with the same lead-time as convection, temperature and sea ice changes ( Figure 8). This strengthens the Icelandic low and increases regional wind stress curl (Figure 6a). This wind-driven feedback may be the trigger that invigorates convection and circulation in the Nordic Seas, a similar mechanism to that proposed by Klockmann et al. (2020) in the MPI-ESM. This in turn accelerates the NAC, which ushes the anomalously saline waters from the STG into the SPG and Nordic Seas. This causes a rapid increase in salinity and active deep-water formation and the AMOC accelerates. Consequently, there is a rapid decline in the salinity gradient between these regions during the transition (Figure 7).
In contrast, temperatures decrease in the Greenland Sea and EGC during the onset of the interstadial. This is due to an increase in Arctic in ow along the EGC due to enhanced circulation throughout the GIN Seas, which drags cold water from the Arctic into both regions. In the Greenland Sea, this results in a decrease in salt added via advection but a more positive sea ice tendency, which contributes to the positive overall salinity change in the basin. This leads AMOC change by approximately 50 years and highlights the increase in circulation in the region prior to AMOC change.

Interstadial Period
The interstadial events are comparatively short-lived, characterised by elevated salinity, warm temperatures, a limited sea-ice extent, and a strong AMOC. Advection is high into the Nordic Seas but at a minimum into the Irminger and central SPG. Despite this, salinity and convection remain high across all regions indicating widespread deep-water formation. In the Irminger Sea, reduced advection is replaced by an increase in diffusive processes, responsible for the increase in regional salinity that maintains convection. In the central SPG however, salinity increases due to a decrease in the freshening role of surface and sea-ice processes due to the retreat of the sea-ice margin.
In the tropics, salinity is at a minimum due to a decrease in surface uxes and an increase in advection out of the basin. The former re ects a northward shift in the ITCZ, which increases precipitation and the P-E balance. As such, the salinity gradient between the Northern NA and STG decreases to a minimum (Figure 7).

Transition into a Stadial
The transition into a stadial is triggered by the collapse of the salinity gradient between the STG and northern NA during the interstadial (Figure 7). This collapse is a response to increasing freshening from surface uxes and advection of the salinity anomaly out of the region. Figure 7 shows that the rapid decline in AMOC strength occurs when the salinity gradient falls from a high of between 1.25 and 1.30 PSU (gkg -1 ) during the stadial to between 0.48 and 0.55 PSU (gkg -1 ). Between these critical values, the salt required to drive convection and deep-water formation cannot be sustained and so the AMOC rapidly diminishes. This is comparable to the MPI-ESM (Klockmann et al., 2020), in which a gradient of approximately 0.5 gkg-1 was shown to trigger an interstadial.
In response, salt transported via advection shifts from the Nordic Seas to the Irminger Sea/central SPG, and convection decreases across all regions. In the Nordic Seas, sea ice formation accelerates as temperatures decrease, which shifts the Icelandic low southward causing an increase in regional SLP and decrease in wind stress curl (Figure 6a). Circulation in the GIN Seas weakens, which reduces the EGC and drawdown of water from the Arctic, so the Greenland Sea and EGC experience a comparative warming. In the tropics, advection out of the region decreases, and the ITCZ is once again shifted southward. This southward shift decreases the P-E balance and salinity begins to accumulate. The oscillation can then start again.

The Role Of Continental Ice Sheets
In order to determine what pushes the AMOC into an oscillating mode, we performed a range of experiments with perturbed glacial boundary conditions (see Methods; Table 2). Altered greenhouse gases, orbital parameters, land-sea mask and bathymetry had no impact on the oscillation. Next, we successively decreased the height and extent of the continental ice sheets within a 28kyr non-oscillating simulation to match the con guration of the 30kyr simulation (see Methods; Table 2). An oscillation in the AMOC is initiated following a decrease in the height and extent of the West and East Laurentide ice sheet in the 28kyrEWLAU and 28kyrLAU simulations (Figure 9c-d). This initiation indicates that AMOC oscillations are sensitive to the con guration of the Laurentide ice sheet in HadCM3B, the height and extent of which triggers the region to be able to oscillate.
The mechanism for this appears to be linked to a northward shift in the position of the jet stream following a decrease in the height and extent of the Laurentide ice sheet (Figure 10a). Due to the reduced topographic forcing, there is an increase in cold continental air drawn in from over the ice sheet, decreasing SATs over much of the SPG (Figure 10b). This enhances sea ice growth in the SPG, which has a freshening regional effect due to its seasonal nature as discussed in Section 4 and shown in Figures 4h   and 6f. Consequently, there is a decrease in salinity and density (Figure 10c) despite the colder temperatures. This spins down the SPG (Figure 10d), reducing deepwater formation and northward heat transport. This weakening is indicated by a decrease in average AMOC strength from 8.98 Sv for the nonoscillating simulations (28kyrFEN, 28kyrELAU and 28kyr50) to 8.16 Sv for the oscillating simulations (28kyrEWLAU and 28kyrLAU). There is also an increase in temperature across the STG/SPG boundary at approximately 40°N ( Figure 10b) and a northward expansion and slowdown of the STG (Figure 10d). This expansion may permit a greater salinity anomaly to develop in the STG. Czaja (2009) similarly showed that a more northerly jet expanded/contracted the subtropical and subpolar gyres.
The result of this is reduced advection out of the STG and so a salinity anomaly can build up in the region whilst salinity decreases in the SPG, as shown in the regional timeseries (Figure 9c-d). Furthermore, the increased meridional temperature gradient due to SPG cooling shifts the ITCZ to an average position further south in the oscillating experiments (not shown), which reduces the P-E balance and contributes to developing a more saline tropics. This build up of salinity initiates the set of feedbacks associated with the stadial identi ed in Section 4. Again, the salinity gradient between the two regions is shown to be key in determining the timing of the interstadial/stadial transition.

Discussion And Summary
This study presents results from a set of glacial simulations using the HadCM3B coupled climate model with 30kyr BP boundary conditions that demonstrates stochastic oscillations on a 1500-year timescale, comparable to the DO events observed in the ice core record (Dansgaard et al., 1993). To our knowledge, this is only the second GCM to demonstrate such variability utilising fully glacial boundary conditions, the rst being CESM1 (Peltier and Vettoretti, 2014; Vettoretti and Peltier, 2016). We utilise a set of salinity diagnostics to provide a comprehensive mechanism for the oscillation, including the triggers between the phases, and what conditions the North Atlantic to be able to oscillate. In summary: The oscillation represents a complex coupled ocean-atmosphere-sea-ice feedback system within the North Atlantic (NA).
This system is strongly linked to the build-up and collapse of the salinity gradient between the STG and northern NA. This determines the level of advection between the two regions, which governs AMOC strength and consequent climatic impacts.
Surface freshwater uxes and sea ice processes are crucial in developing this salinity gradient. In the tropics, the position of the ITCZ regulates the P-E balance, which acts as a positive feedback in increasing/decreasing regional salinity. Sea-ice feedbacks are location dependant. In the Nordic Seas, sea-ice acts as a negative feedback and may contribute in pushing the oscillation into an opposing phase.
The key site of advection and deep-water formation moves from the SPG during stadials to the Nordic Sea during interstadials.
The interstadial may be initiated by an atmospheric forcing; this is a wind-driven atmospheric feedback in the Nordic Seas in response to an increase in regional temperature, a decrease in sea ice cover and a consequent increase in SLP, which invigorates wind stress and subsequently convection.
A stadial is initiated by an ocean forcing; this involves the collapse of the salinity gradient between the STG and Northern NA. This collapse weakens advection into the Nordic Seas and diminishes deep-water formation.
This oscillation of the AMOC is sensitive to the con guration of the Laurentide ice-sheet. A reduction in ice-sheet height shifts the jet stream northward, causing regional cooling and an increase in seasonal sea ice concentration in the SPG. This freshens the region and reduces deep-water formation, which restricts and weakens the SPG. In contrast, the STG expands northward.
Weaker STG/SPG transport and a shift in the ITCZ permits a salinity anomaly to develop in the STG. This initiates the rst interstadial and associated feedbacks that maintain the oscillation.
The oscillation therefore re ects a salt oscillator mechanism in the NA, similar to that initially proposed by Broecker et al. (1990), and also identi ed in CESM1 (Peltier and Vettoretti, 2014;Vettoretti and Peltier, 2016;Vettoretti and Peltier, 2018). The salinity gradient between the STG and Northern North Atlantic . Sea ice evidently plays an important role in conditioning (i.e. freshening and weakening) the SPG and in triggering the interstadials. However the tropical P-E balance, which re ects the position of the ITCZ, is also of critical importance for maintaining the oscillation. Evidence for such ITCZ shifts on DO time scales are apparent in a wide range of observational studies using tropical speleothem and sediment records (Deplazes et al., 2013;Kanner et al., 2012;Peterson et al., 2000;Dupont et al., 2010;Wang et al., 2004), and a number of modelling studies (e.g. Brown and Galbraith, 2016;Chiang and Bitz, 2005;Zhang and Delworth, 2005). Indeed Vellinga and Wu (2004) and Armstrong et al. (2017) showed that the position of the ITCZ was important in driving centennial AMOC variability in HadCM3, these results also highlight its important role on millennial timescales.
The mechanism that initiates the oscillation differs to those previously proposed. Zhang et al. (2014) increased the height of the Laurentide ice sheet, which increased North Atlantic SATs, accelerating SPG circulation and instigating an oscillation. In contrast Madonna et al. (2017) showed that an equator-ward shift in the jet stream freshens and weakens the SPG via Greenland blocking, which reduces heat transport and initiates an oscillation. We show that a pole-ward shift of the jet weakens both Atlantic gyres and permits a salinity anomaly to develop to drive the interstadial. Therefore it is clear that the impact of the Laurentide ice sheet on the jet stream and the role this plays in driving AMOC oscillations remains unclear within models. The mechanism we identify is evidently linked to the seasonal and freshening role of sea ice in the SPG, which is likely to be a highly model-dependent characteristic of this region, in contrast to other models. It is worth noting that wind stress is directly applied to the surface current in HadCM3B through sea-ice cover. Therefore the increase in sea-ice cover does not decouple the surface currents from the wind forcing, such as in the real world . Therefore increased sea-ice cover does not contribute to weakening gyre strength, so density change is the leading factor.  Table   2). These simulations did not demonstrate an oscillation, which indicates the important role of the background climate state. This validates an idea discussed by Li and Born (2019) of a 'sweet spot' in model simulations, where the system is more prone to spontaneous oscillatory behaviour within a speci c range of conditions. It might be that the con guration of the ice sheet pushes the system across a transition point, but the background climate is important in permitting a salt oscillator mode to develop, at which point feedbacks within the system can maintain the oscillation.  The boundary conditions for the 28kyr and 30kyr simulations and anomaly between the two, a) bathymetry (m), b) ice sheet fractional extent, and c) orography which refers to ice-sheet height (m). The speci c regions used in the altered ice sheet simulations are identi ed in the boxes; Fennoscandian (green), East Laurentide (red) West Laurentide (blue) and entire Laurentide (purple).

Figure 2
Timeseries of normalised ocean (top panel) and atmosphere (bottom panel) variables and the AMOC index (de ned as the mean AMOC strength between 40° and 50°N at 800m) for the 30kyrBASE simulation. A 100-year low pass Lanczos lter has been applied to remove high frequency variability.

Figure 3
Composite spatial climate variables for the 30kyrBASE simulation for the long-term mean (Ltm; left panels), and the anomalies relative to this baseline during periods of strong AMOC (AMax-Ltm; middle panels) and weak AMOC (AMin-Ltm; right panels). AMax (interstadials) and AMin (stadials) are de ned as when the AMOC index is one standard deviation above and below the mean. Ocean variables are averaged over 0-666m water depths.

Figure 4
Composite spatial patterns for the salinity tendencies (PSU/yr) averaged over 0-666m for the 30kyrST simulation during periods of AMax (interstadials; left panels), AMin (stadials; middle panels) and the difference between the phases (right panels). Convection is averaged over 0-300m to highlight changes more clearly. Positive/negative values indicate where salt is added/removed. AMax (interstadials) and AMin (stadials) are de ned as when the AMOC index is one standard deviation above and below the mean.

Figure 5
The rst EOF and the corresponding PC1 timeseries (PSU/yr) for the salinity tendencies averaged over 0-666m from the 30kyrST simulation. Convection is averaged over 0-300m to highlight variability more clearly. The amount of variance explained by EOF1 is stated. Note the change in scale.
are shown on the map. A 100-year low pass Lanczos lter has been applied to remove high frequency variability.  Lead lag analysis for the Nordic Seas of the salinity tendencies and a range of ocean and atmosphere variables with the AMOC index (de ned as the mean AMOC strength between 40° and 50°N at 800m) for the 30kyrST simulation. Tendencies and variables have been averaged over 0-666m water depth except for convection which is averaged over 0-300m. The data has been ltered with a low pass Lanczos lter prior to correlation analysis.

Figure 9
The same as Figure 7 but for the simulations with altered ice sheets (see Table 2).