Model intercomparison studies of coupled carbon-climate simulations have the potential to improve our understanding of the processes explaining the pCO2 drawdown at the Last Glacial Maximum (LGM) and to identify related model biases. Models participating in the Paleoclimate Modelling Intercomparison Project (PMIP) now frequently include the carbon cycle. The ongoing PMIP-carbon project provides the first opportunity to conduct multimodel comparisons of simulated carbon content for the LGM time window. However, such a study remains challenging due to differing implementation of ocean boundary conditions (e.g. bathymetry and coastlines reflecting the low sea level) and to various associated adjustments of biogeochemical variables (i.e. alkalinity, nutrients, dissolved inorganic carbon). After assessing the ocean volume of PMIP models at the pre-industrial and LGM, we investigate the impact of these modelling choices on the simulated carbon at the global scale, using both PMIP-carbon model outputs and sensitivity tests with the iLOVECLIM model. We show that the carbon distribution in reservoirs is significantly affected by the choice of ocean boundary conditions in iLOVECLIM. In particular, our simulations demonstrate a ~250 GtC effect of an alkalinity adjustment on carbon sequestration in the ocean. Finally, we observe that PMIP-carbon models with a freely evolving CO2 and no additional glacial mechanisms do not simulate the pCO2 drawdown at the LGM (with concentrations as high as 313, 331 and 315 ppm), especially if they use a low ocean volume. Our findings suggest that great care should be taken on accounting for large bathymetry changes in models including the carbon cycle.

Our limited understanding of millennial-scale variability in the context of the last glacial period can be explained by the lack of a reliable modelling framework to study abrupt climate changes under realistic glacial backgrounds. In this article, we describe a new set of long-run Last Glacial Maximum experiments where such climate shifts were triggered by different snapshots of ice-sheet meltwater derived from the early stages of the last deglaciation. Depending on the location and the magnitude of the forcing, we observe three distinct dynamical regimes and highlight a subtle window of opportunity where the climate can sustain oscillations between cold and warm modes. We identify the European-Arctic and Nordic Seas regions as being most sensitive to meltwater discharge in the context of switching to a cold mode, compared to freshwater fluxes from the Laurentide ice sheets. These cold climates follow a consistent pattern in temperature, sea ice and convection, and are largely independent from freshwater release as a result of effective AMOC collapse. Warm modes, on the other hand, show more complexity in their response to the regional pattern of the meltwater input, and within them, we observe significant differences linked to the reorganisation of deep water formation sites and the subpolar gyre. Broadly, the main characteristics of the oscillations, obtained under full-glacial conditions with realistically low meltwater discharge, are comparable to δ18O records of the last glacial period, although our simplified experiment design prevents detailed conclusions from being drawn on whether these represent actual Dansgaard-Oeschger events.

North Pacific atmospheric and oceanic circulations are key missing pieces in our understanding of the reorganisation of the global climate system since the Last Glacial Maximum (LGM). Here, using a basin-wide compilation of planktic foraminiferal δ18O, we show that the North Pacific subpolar gyre extended ~3 degrees further south during the LGM, consistent with sea surface temperature and productivity proxy data. Analysis of an ensemble of climate models indicates that the expansion of the subpolar gyre was associated with a substantial gyre strengthening. These gyre circulation changes were driven by a southward shift in the mid-latitude westerlies and increased wind-stress from the polar easterlies. Using single-forcing model runs, we show these atmospheric circulation changes are a non-linear response to the combined topographic and albedo effects of the Laurentide Ice Sheet. Our reconstruction suggests the gyre boundary (and thus westerly winds) began to migrate northward at ~17-16 ka, during Heinrich Stadial 1.

The maximum extent of the last North American ice sheet is well constrained empirically, but it has proven to be challenging to simulate with coupled Climate-Ice Sheet models. Coupled Climate-Ice Sheet models are often too computationally expensive to sufficiently explore uncertainty in input parameters, and it is unlikely values calibrated to reproduce modern ice sheets will reproduce the known extent of the ice at the Last Glacial Maximum. To address this, we run a series of ensembles with a coupled Climate-Ice Sheet model (FAMOUS-ice), simulating the final stages of growth of the last North American Ice Sheets’ maximum extent. Using this large ensemble approach, we explore the influence of uncertain ice sheet, albedo, atmospheric, and oceanic parameters on the ice sheet extent. We find that albedo parameters determine the majority of uncertainty when simulating the Last Glacial Maximum North American Ice Sheets. Importantly, different albedo parameters are needed to produce a good match to the Last Glacial Maximum North American Ice Sheets than have previously been used to model the contemporary Greenland Ice Sheet, due to differences in cloud cover over ablation zones. Thus calibrating coupled climate-ice sheet models solely for present day strongly biases simulations of past and future climates different from today.