The simulated climate of the Last Glacial Maximum and insights into the global marine carbon cycle

. The ocean’s ability to store large quantities of carbon, combined with the millennial longevity over which this reservoir is overturned, has implicated the ocean as a key driver of glacial-interglacial climates. However, the combination of processes that cause an accumulation of carbon within the ocean during glacial periods is still under debate. Here we present simulations of the Last Glacial Maximum (LGM) using the CSIRO Mk3L-COAL Earth System Model to test the contribution of physical and biogeochemical processes to ocean carbon storage. For the LGM simulation, we ﬁnd a signiﬁcant global 5 cooling of the surface ocean (3.2 ◦ C) and the expansion of both minimum and maximum sea ice cover broadly consistent with proxy reconstructions. The glacial ocean stores an additional 267 Pg C in the deep ocean relative to the Pre-Industrial (PI) simulation due to stronger Antarctic Bottom Water formation. However, 889 Pg C is lost from the upper ocean via equilibration with a lower atmospheric CO 2 concentration and a global decrease in export production, causing a net loss of carbon relative to the PI ocean. The LGM deep ocean also experiences an oxygenation ( > 100 mmol O 2 m − 3 ) and deepening of the calcite 10 saturation horizon (exceeds the ocean bottom) at odds with proxy reconstructions. With modiﬁcations to key biogeochemical processes, which include an increased export of organic matter due to a simulated release from iron limitation, a deepening of remineralisation and decreased inorganic carbon export driven by cooler temperatures, we ﬁnd that the carbon content of the glacial ocean can be sufﬁciently increased (317 Pg C) to explain the reduction in atmospheric and terrestrial carbon at the LGM (194 ± 2 and 330 ± 400 Pg C, respectively). Assuming an LGM-PI difference of 95 ppm p CO 2 , we ﬁnd that 55 ppm can be 15 attributed to the biological pump, 28 ppm to circulation changes, and the remaining 12 ppm to solubility. The biogeochemical modiﬁcations also improve model-proxy agreement in export production, carbonate chemistry and dissolved oxygen ﬁelds. Thus, we ﬁnd strong evidence that variations in the oceanic biological pump exert a primary control on the climate. ◦ × with 21 vertical levels. using both the full climate system model and the stand-alone ocean model. Two fully coupled model experiments were undertaken to simulate the Pre-Industrial (Cpl-PI) and Last Glacial Maximum 5 (Cpl-LGM) climates. The Cpl-PI climate was obtained by forcing the model with an atmospheric CO 2 concentration of 280 ppm and by prescribing 1950 CE values for the orbital parameters. This experiment was integrated for a total of 10,000 years (Phipps et al., 2013). The Cpl-LGM simulation followed the protocol developed by Phase III of the Palaeoclimate Modelling Intercomparison Project (PMIP3), with the exception that no changes were made to terrestrial topography, oceanic bathymetry or the positions of the coastlines. The closure of important oceanic connections due to sea level loss, such as the Bering Strait, 10 was not considered. The atmospheric CO 2 equivalent concentration was set to 167 ppm, providing a radiative forcing equivalent to the speciﬁed reductions in the atmospheric concentrations of CO 2 , CH 4 and N 2 O from 280 ppm/760 ppb/270 ppb for pre-industrial simulations to 185 ppm/350 ppb/200 ppb for LGM simulations. The orbital parameters were set to values for 21,000 years BP. The Cpl-LGM simulation was initialised from the state of the Cpl-PI simulation at the end of model year 100. The model was then integrated for a total of 3,900 model years, until it had reached quasi-equilibrium. Over this integration 15 the ocean experienced an increase in salinity by 0.5 psu due to increased evaporation, which reﬂected the coupling between a cooler, drier atmosphere and the ocean.


Introduction
The late Pleistocene is characterised by a sawtooth-like cycling between cool glacial and warm interglacial states (Emiliani, 1966;Shackleton, 1967). Global temperatures and atmospheric CO 2 are strongly correlated across these climate cycles with approximately 80-100 ppm of change corresponding to global-mean temperature variations of 3-4 • C (Grootes and Stuiver, 1997;Jouzel et al., 1987;Parrenin et al., 2013;Petit et al., 1999). This correlation was extended to causation by evidence that 5 increases in atmospheric CO 2 preceded decreases in global ice mass during glacial terminations (Sowers et al., 1991;Broecker and Henderson, 1998). The large potential of the ocean to store carbon, coupled with evidence that the terrestrial reservoir was diminished under glacial conditions (Shackleton, 1977), has implicated the ocean as the prime candidate for driving glacialinterglacial changes in atmospheric CO 2 (Broecker, 1982;Skinner et al., 2015;Wilson et al., 2015). However, identifying the combination of mechanisms that drove a flux of carbon into the ocean during glacial periods remains a fundamental and largely 10 unresolved problem.
If we first consider only physical changes, a net influx of CO 2 caused by cooling is a feature of the glacial ocean. However, the increase in solubility attributed to cooling is partially offset by increased salinity due to a lower sea level (Siddall et al., 2003), so that the contribution of solubility changes is estimated at ∼13 ppm of the total 80-100 ppm CO 2 drawdown 15 (Sigman and Boyle, 2000;Kohfeld and Ridgwell, 2009). Therefore, other physical changes associated with a glacial climate, notably changes to the large-scale circulation and sea ice fields (Stephens and Keeling, 2000), may make a considerable contribution to carbon sequestration. Proxy evidence (Curry and Oppo, 2005;Duplessy et al., 1988;McManus et al., 2004;Oliver et al., 2010;Skinner et al., 2010) and model experiments (Brovkin et al., 2007;Hain et al., 2010;Menviel et al., 2012;Watson and Naveira Garabato, 2006) of the glacial climate have shown that a greater proportion of the deep ocean was dominated 20 by southern source waters (see Adkins, 2013, for a review). The existence of this glacial-type circulation has been connected to an expanded sea ice field (Ferrari et al., 2014). An expansion of southern source waters throughout the deep ocean and an expanded sea ice field are now considered to be primary candidates for carbon sequestration during glacial climates (Adkins, 2013). 25 However, the most promising explanations of the glacial decline in atmospheric CO 2 involve changes to ocean productivity in concert with reorganisations of the global overturning circulation (Hain et al., 2010;Gottschalk et al., 2016). An increased glacial productivity, first proposed by Broecker (1982) and explored by Archer et al. (2000) by increasing the global nutrient inventory, is now an established feature of the sub-Antarctic zone due to enhanced aeolian deposition of iron (Martinez-Garcia et al., 2014). The Southern Ocean exerts a strong control on atmospheric CO 2 through its direct connection deep waters (Marinov et al., 2006), and enhanced productivity in this zone is thus a prime candidate for explaining a fraction of the glacialinterglacial CO 2 difference.
In numerous other regions, however, productivity appears to have been reduced during glacial climates. The affected regions include waters south of the Antarctic Polar Front (Francois et al., 1997;Jaccard et al., 2013), the North Pacific (Crusius et al., 2004;Jaccard et al., 2005;Kohfeld and Chase, 2011;Ortiz et al., 2004), tropical Indian Ocean (Singh et al., 2011) and the Equatorial Pacific (Costa et al., 2016;Herguera, 2000;Loubere et al., 2007). A weaker export production in these regions would have offset a strengthened biological pump in the sub-Antarctic, thereby weakening the ability of the ocean to store carbon during glacial conditions. Whether the strengthening of the biological pump in the glacial sub-Antarctic was sufficient 5 to increase the carbon content of the ocean despite losses in productivity in other regions therefore requires further testing.
This has led some authors to look for alternative biogeochemical mechanisms. A notable example is the application of temperature-dependent remineralisation to global fluxes of organics into the interior ocean. The positive relationship between microbial metabolism and temperature has been known for some time (Eppley, 1972), but it was only recently that this rela- 10 tionship was applied to a glacial setting. By prescribing a global cooling of 5 • C, Matsumoto (2007) reduced atmospheric CO 2 by ∼35 ppm. Further research has shown that even a slight deepening in the remineralisation profile can cause large changes in oceanic carbon content by reducing surface ocean inorganic carbon concentrations, which in turn strengthens the air-sea influx of carbon (Menviel et al., 2012). 15 Another biogeochemical mechanism that is proposed to increase oceanic carbon storage is an altered calcium carbonate to organic carbon (CaCO 3 :C org ) export production ratio (Sigman et al., 1998;Archer et al., 2000). A global decrease in the CaCO 3 :C org ratio would enhance carbon storage by reducing the pCO 2 of surface waters (see Sigman and Boyle, 2000, for a review). A decrease in CaCO 3 production could be caused by cooling (Stoll et al., 2002) and/or an increase in silicate delivery to the lower latitude oceans that simulated organic carbon production (Matsumoto et al., 2002). 20 Numerous physical and biogeochemical changes have been associated with a glacial ocean and all have been identified in some respect as important drivers of glacial-interglacial climate cycles. Now, recent insights into the distributions of dissolved oxygen  and carbonate species (Yu et al., 2014) at the Last Glacial Maximum (LGM; ∼21,000 years BCE) provide new constraints to identify which combination of physical and biogeochemical changes could have realisti- 25 cally sequestered carbon within the ocean at this time. Here, we use an Earth System Model with attached biogeochemistry, CSIRO Mk3L-COAL, to test current theories against these new insights. Using our simulated LGM ocean state, we quantify the contribution of physical and biogeochemical changes to the estimated increased of 520 ± 400 Pg of carbon within the oceanic reservoir at the LGM (Ciais et al., 2011), and demonstrate the importance of marine biogeochemistry to global climate. 30

Model and experiments
The model simulations were performed using the CSIRO Mk3L climate system model version 1.2 (Phipps et al., 20111.2 (Phipps et al., , 2012, which includes components that describe the atmosphere, land, sea ice and ocean. The horizontal resolution of the atmosphere, land and sea ice models are 5.6 • ×3.2 • in the longitudinal and latitudinal dimensions, respectively, with 18 vertical levels. The ocean model has a horizontal resolution of 2.8 • ×1.6 • with 21 vertical levels. For this study, we conduct simulations using both the full climate system model and the stand-alone ocean model. Two fully coupled model experiments were undertaken to simulate the Pre-Industrial (Cpl-PI) and Last Glacial Maximum 5 (Cpl-LGM) climates. The Cpl-PI climate was obtained by forcing the model with an atmospheric CO 2 concentration of 280 ppm and by prescribing 1950 CE values for the orbital parameters. This experiment was integrated for a total of 10,000 years (Phipps et al., 2013). The Cpl-LGM simulation followed the protocol developed by Phase III of the Palaeoclimate Modelling Intercomparison Project (PMIP3), with the exception that no changes were made to terrestrial topography, oceanic bathymetry or the positions of the coastlines. The closure of important oceanic connections due to sea level loss, such as the Bering Strait, 10 was not considered. The atmospheric CO 2 equivalent concentration was set to 167 ppm, providing a radiative forcing equivalent to the specified reductions in the atmospheric concentrations of CO 2 , CH 4 and N 2 O from 280 ppm/760 ppb/270 ppb for pre-industrial simulations to 185 ppm/350 ppb/200 ppb for LGM simulations. The orbital parameters were set to values for 21,000 years BP. The Cpl-LGM simulation was initialised from the state of the Cpl-PI simulation at the end of model year 100.
The model was then integrated for a total of 3,900 model years, until it had reached quasi-equilibrium. Over this integration 15 the ocean experienced an increase in salinity by 0.5 psu due to increased evaporation, which reflected the coupling between a cooler, drier atmosphere and the ocean.
With the climate state provided by the Cpl-LGM experiment, nine additional ocean biogeochemical simulations were made with different physical and biogeochemical conditions to explore the effect on the carbon cycle (Table 1). These experiments 20 utilised Mk3L-COAL (Carbon-Ocean-Atmosphere-Land), an enhanced version of the Mk3L climate system model which includes biogeochemical modules embedded within the ocean, atmosphere and terrestrial models. All experiments were forced by key boundary conditions (wind stresses, temperature, salinity, incident radiation, sea ice), which were obtained as monthly averages over the final 50 years of the fully coupled model experiments. The heat and freshwater fluxes into the ocean were determined by relaxing the SST and SSS towards the prescribed fields using a 20 day timescale. An additional 0.5 psu was 25 added to the salinity field of the LGM experiments to ensure that the ocean was 1 psu more saline than the PI. For a description of the ocean biogeochemistry the reader is directed towards Appendix A of Matear and Lenton (2014)  LGM experiments in which the equations controlling biogeochemical cycling were altered. It should be made clear that these experiments did not explicitly simulate the biogeochemical changes caused by an altered climate in any mechanistic sense. However, the prescription of the following changes allowed us to undertake a theoretical investigation into their capacity to sequester carbon at the LGM. The experiments were as follows: The scaling factor (S O npp ) was increased by a factor of 10 (Eq. (1)) to increase the export of Particulate Organic Carbon (POC) from the surface ocean, and therefore strengthen the biological carbon pump. Increasing POC export in the LGM ocean was motivated by an enhanced delivery of iron to the surface ocean via aeolian dust at the LGM (Martin, 1990;Martinez-Garcia et al., 2014).
where, V max = The temperature-depedent maximum growth of phytoplankton (Eppley, 1972) P k = The half-saturation constant for nutrient-limited growth, set to 0.1 mmol PO 4 m −3 F (I) = The productivity versus irradience equation for determining light-limited growth (Clementson et al., 1998) 20 O-LGM BGC rem . The POC remineralisation depth was increased by changing the power law exponent (b) in Equation (2) from -0.9 to -0.7, which replicated a bulk shift of POC from the upper to the deep ocean. The motivation for increasing the amount of POC that reaches deeper levels is the expectation that a cooler ocean would reduce the rate of bacterial remineralisation of POC in the upper ocean (Rivkin and Legendre, 2001;Matsumoto, 2007). This change increased the simulated POC that reaches the 1,000 m depth level from 12.5 % to 20 %.

Remin
O-LGM BGC pic . Export production of Particulate Inorganic Carbon (PIC) was turned off by setting the rain ratio (R P IC ) of PIC:POC to zero in Equation (3). The motivation for reducing PIC export in the glacial ocean is related to the positive linear relationship between calcification and temperature (Stoll et al., 2002), and an increased silicate supply to lower latitudes that potentially favoured non-calcifying producers (Matsumoto et al., 2002).
O-LGM BGC all . All three modifications to ocean biogeochemistry were employed. All ocean-only simulations were integrated for 10,000 years to ensure that the ocean carbon cycle reached a steady state.
To assess whether the behaviour of the biogeochemical tracers within the coupled model differed from those in the oceanonly model, we ran the coupled model with online ocean biogeochemistry for a further 1,000 years using the steady-state 5 biogeochemical fields from the ocean-only experiments. This assessment was made using both the PI and LGM climates. For key diagnostics, such as the meridional overturning circulation, ocean carbon content and global export production, the behaviour of the ocean-only simulation differed by less than 1 % from the coupled simulations. Given the computational speed of the ocean-only model, these experiments provide an ideal platform to test the sensitivity of the ocean biogeochemical fields to the parameterisations used in the biogeochemical model. 10

Results and discussion
In the following, we first discuss the simulated physical changes to the ocean observed between the Cpl-PI and Cpl-LGM simulations. Second, we discuss how the ocean biogeochemical fields differed between the O-PI and O-LGM simulations, which were forced with the output of the coupled simulations. Finally, we explore how modifying biogeochemical parameterisations 15 alters the biogeochemistry, including changes to carbon storage, export production, carbonate chemistry and dissolved oxygen. The simulated change in SST between the Cpl-PI and the Cpl-LGM simulations shows a similar magnitude and spatial struc-20 ture to proxy reconstructions and prior modelling studies, with greatest cooling in the equatorial oceans, high latitudes and eastern boundary currents, and the least cooling in the subtropics and western boundary current regions ( Fig. 1; Table 2). The global SST mean of the Cpl-LGM was 3.2 • C cooler than the Cpl-PI. This change falls within the range of estimates (∼2-4 • C) produced by other climate models (Alder and Hostetler, 2015;Annan and Hargreaves, 2013;Braconnot et al., 2007), but sits towards the cooler limits of previous multiproxy SST reconstructions that estimate a change of 2 ± 1.8 • C (Ballantyne et al., 25 2005;Waelbroeck et al., 2009). However, a recent reanalysis of the proxy data presented by Waelbroeck et al. (2009) showed past estimates may have underestimated cooling by as much as 50 % (Ho and Laepple, 2015). This finding reconciles some disagreement between climate models and palaeoproxies, and places our simulated cooling of 3.2 • C well within the bounds of uncertainty in reconstructions. 30 Regionally, the greatest cooling took place in the high latitudes and in the Equatorial Pacific, where temperatures were in excess of 4 • C cooler than the Cpl-PI climate. Meanwhile, the Western Pacific Warm Pool, subtropical gyres and western boundary currents cooled less (0.5-3.0 • C). Again, proxy (Waelbroeck et al., 2009) and climate modelling (Annan and Hargreaves, 2013) are consistent with both the magnitude and spatial pattern of cooling. A notable example of model-data agreement is in the Pacific sector of the Southern Ocean, where SSTs were up to and in excess of 4 • C cooler at the LGM (Benz et al., 2016).
Enhanced cooling in the high latitudes and in the eastern boundary currents generated strong zonal and meridional temperature 5 gradients relative to Cpl-PI SST. There is a consistent regional pattern to SST cooling in the LGM emerging from proxy and model simulations (Annan and Hargreaves, 2013;Braconnot et al., 2007) that is broadly consistent with our simulated cooling.
Where there is still large uncertainty in SST change at the LGM is in the tropical ocean (see Annan and Hargreaves, 2015, for a review). The Cpl-LGM cooling of 3.3 • C across the tropical ocean (15 • S -15 • N) is greater than other simulations (Annan 10 and Hargreaves, 2013;Braconnot et al., 2007), but falls well within the −5.1 to −2.17 • C estimated by Ho and Laepple (2015).
Regionally, climate models (Otto-Bliesner et al., 2009) and proxies (Waelbroeck et al., 2009) agree that cooling in the tropical Atlantic Ocean probably exceeded cooling in the tropical Pacific and Indian Oceans by roughly 1 • C. In contrast, the tropical Pacific Ocean cooled by 2 • C more than the Tropical Atlantic and Indian Oceans in the Cpl-LGM simulation. Although SSTs in the east equatorial Pacific have been reported as 1.5-3.0 • C cooler than the PI (Dubois et al., 2014;Kucera et al., 2005), 15 the simulated cooling over much of the tropical Pacific appears excessive compared to previous simulations (Braconnot et al., 2007).

Sea ice extent
The Cpl-PI sea ice extent is consistent with estimates made using satellite measurements during the 1979-1987 period (Glo-20 ersen et al., 1993, Table 2). These measurements represent the first global estimates of sea ice coverage, and although there is evidence that sea ice has declined by 20 % since the 1950's (Curran et al., 2003), the strong agreement between the Cpl-PI sea ice fields and the observations of Gloersen et al. (1993) provide a good benchmark for assessing LGM sea ice changes.
Associated with cooler SSTs, sea ice coverage (fractional sea ice area ≥ 15 %) was greatly expanded in the Cpl-LGM for 25 both hemispheres relative to the Cpl-PI (Fig. 2, Table 2). In the Southern Hemisphere, total sea ice coverage increased by ∼120 % and ∼225 % at its seasonal maximum and minimum, respectively, relative to the Cpl-PI. In the Northern Hemisphere, total sea ice coverage increased by ∼145 % and ∼105 % at its seasonal maximum and minimum, respectively, relative to the Cpl-PI.
These increases correspond to equatorward expansions of the sea ice field of between 5-10 • around the Southern Ocean, and in excess of 15 • in both the North Atlantic and Pacific Oceans. 30 The simulated expansion of sea ice around much of the Southern Ocean agrees well with proxy reconstructions. Maximum sea ice extent reached as far north as 47 • S in both the Atlantic and Indian sectors (Gersonde et al., 2005) and as far north as 55 • S in the Pacific sector of the Southern Ocean (Benz et al., 2016;Gersonde et al., 2005, Fig. 2). The magnitude of growth in the Atlantic and Indian sectors has been tested and largely supported by subsequent studies (Collins et al., 2012;Xiao et al., 2016), and is consistent with our Cpl-LGM sea ice field. In the Pacific sector, however, the simulated maximum sea ice edge extends well equatorward of the 55 • S boundary that has been defined by Benz et al. (2016) (Fig. 2). By comparing the coverage of sea ice in the Southern Hemisphere of the Cpl-LGM (∼46×10 6 km 2 ) with that estimated by Gersonde et al. (2005) km 2 ), we can attribute the simulated excess of sea ice in the glacial Southern Ocean to a possible overestimate in the Pacific 5 sector.
The Cpl-LGM sea ice was broadly consistent with reconstructions in the North Atlantic, with the exception that too much ice covered the Nordic Seas. The central and eastern parts of the subpolar North Atlantic, including the Nordic Seas, are thought to have been at least seasonally ice-free (Pflaumann et al., 2003;De Vernal et al., 2005). The Cpl-LGM sea ice field 10 showed strong, year-round cover in these regions. However, good model-proxy agreement was produced in other parts of the North Atlantic. Perennial sea ice cover was present in the Greenland Sea and Fram Strait during the LGM (Müller et al., 2009;Telesiński et al., 2014). There is also evidence that winter sea ice reached south of Iceland to fill much of the Labrador Sea (Pflaumann et al., 2003) and extended along the eastern Canadian margin (De Vernal et al., 2005). These features were produced in the Cpl-LGM simulation. 15 In the North Pacific, the Cpl-LGM sea ice field expanded across a large area. Proxy reconstructions indicate the presence of strong cover in the Okhotsk Sea (Sakamoto et al., 2005;Nürnberg and Tiedemann, 2004), the Japan Sea (Ikehara, 2003) and the western Bering Sea (Riethdorf et al., 2013) during the LGM, with seasonally ice-free conditions in the central west (Jaccard et al., 2005). Thus, like the North Atlantic, sea ice presence was likely more extensive in the western margins of the 20 North Pacific, and this pattern is replicated by other climate models (Fig. 2). In the Cpl-LGM simulation, more intense sea ice cover developed in the west, but year-round cover also developed over the central North Pacific and contrasts directly with the findings of Jaccard et al. (2005), who argued for ice-free conditions during the summer. Thus, our simulated sea ice extent in the North Pacific may be an exaggeration.

Meridional overturning circulation
The changes observed in the surface ocean within the Cpl-LGM climate were accompanied by changes in the global meridional overturning circulation ( Fig. 3; Table 2). The rate of Antarctic Bottom Water (AABW) formation in the Southern Ocean doubled between the Cpl-PI and Cpl-LGM experiments, increasing from from 10.3 to 20.2 Sv. An increase in surface density of 0.9 kg m −3 between 60 • S to 40 • S drove this intensification, and also strengthened transport by the Antarctic Circumpolar 30 Current (ACC) from 140 to 309 Sv. Meanwhile, the Atlantic Meridional Overturning Circulation (AMOC) weakened from 15.6 to 11.4 Sv. The weakened glacial AMOC was also associated with a shoaling of its lower boundary from approximately 3,000 to 1,500 metres. As a result, much of the Atlantic Ocean below 1,500 metres was dominated by AABW as part of the lower overturning cell.
These changes in the lower and upper overturning cells were conducive to the development of a global overturning circulation dominated by a denser AABW and a shallower AMOC. These results are supported by numerous palaeonutrient tracers showing an increased presence of southern source waters within the deep North Atlantic Ocean during glacial periods (Curry and Oppo, 2005;Duplessy et al., 1988;Keigwin, 2004;Marchitto and Broecker, 2006;Oliver et al., 2010;Skinner et al., 2010).

5
The prevailing interpretation of the palaeonutrient tracers is that the maximum depth of the AMOC was displaced above about 2,000 metres, and that this shoaling facilitated the development of a saltier, more stratified glacial deep ocean (Adkins, 2013). Ferrari et al. (2014) have linked these changed to the expansion of sea ice in the Southern Ocean, which caused a greater proportion of Circumpolar Deep Water to rise into a zone of negative buoyancy flux and thereby produce greater quantities of denser AABW. 10 However, contradictory changes in the glacial overturning circulation have been simulated in other climate system models. The rates of AABW formation tend to increase under LGM conditions , but responses of the AMOC among models are highly variable. In earlier experiments as part of the PMIP2 project, the AMOC response ranged between 40 % above and below the PI rate of overturning . Our weakened (∼35 %) and shallower (∼50 %) glacial 15 AMOC is therefore consistent with the lower bounds of these simulations, as are our rates of AABW formation (Table 2). More recent LGM simulations as part of the PMIP3 project, however, developed stronger and deeper glacial AMOCs (Muglia and Schmittner, 2015). Furthermore, a recent reconstruction of Southern Ocean circulation indicates that extreme intensifications of AABW formation and ACC transport at the LGM are unlikely (Lynch-Stieglitz et al., 2016). Consequently, these results challenge our simulated changes in meridional overturning, as well as the prevailing interpretation of palaeonutrient evidence. 20 Despite inconsistencies between climate model simulations, palaeonutrient reconstructions continue to support the existence of a shallower AMOC overlying southern source waters during glacial periods. Variations in the isotopic signature of Neodynium, for instance, indicate that AABW was more dominant in the deep ocean during the LGM, and that its mixing with a glacial form of NADW was more intense (Howe et al., 2016). These and other authors (Burckel et al., 2016) find further support for 25 the presence of a shallower AMOC above 2,500 metres. Moreover, simulated distributions of carbon isotopes across a range of idealised circulations have shown that a shallower AMOC is necessary to ensure good model-proxy agreement at the LGM (Menviel et al., 2016). Importantly, our Cpl-LGM simulation developed an increased presence of AABW throughout the global deep ocean and the development of a shallower AMOC.

LGM climate: biogeochemical fields
The physical changes in the ocean between the Cpl-PI and the Cpl-LGM, as described above, caused significant changes in ocean biogeochemistry within the ocean-only simulations (Table 3). To assist in the discussion of the large-scale biogeochemical changes we divide the upper and the deep ocean based on the 2,000 m depth. This approach also allows for more clearly distinguishing between changes to the global overturning circulation, air-sea exchange and biological processes on the biogeochemical fields.

Carbon
For experiment O-LGM, the carbon content of the ocean was 622 Pg C less than the O-PI experiment ( Fig. 4; Table 3). The net 5 loss of carbon reflected the combined effect of physical changes to the ocean, which include an increase in solubility due to cooling, an expanded sea ice field, an altered overturning circulation, and the tendency for outgassing caused by a lower pCO 2 .
The physical changes were sufficient in combination to increase carbon in the deep ocean by 267 Pg C. However, lowering the atmospheric pCO 2 to 185 ppm drove a large amount of carbon out of the ocean, causing a loss of 889 Pg C from waters in the upper 2,000 metres. The combined effect of our simulated LGM physical state could therefore not overcome the global 10 equilibration with a lower atmospheric pCO 2 concentration.
The increase in carbon content in the glacial deep ocean did suggest that despite the net loss caused by outgassing, the glacial ocean was indeed conducive to storing carbon. The loss of carbon from the ocean of experiment O-PI LGM CO2 demonstrated this ( Fig. 4; Table 3). The ocean carbon content of O-PI LGM CO2 evolved to be 1,290 Pg less than the O-PI experiment, which placed 15 the O-LGM carbon content as 668 Pg greater than O-PI LGM CO2 . This confirmed that the glacial ocean had a greater ability to store carbon than the PI ocean. Therefore, we primarily implicate circulation and secondarily solubility changes as the physical drivers of oceanic carbon 30 storage at the LGM, while presenting sea ice expansion as a mechanism for reducing ocean carbon storage. Regarding solubility, previous work has constrained the effect of glacial solubility on atmospheric carbon to approximately 15 ppm (Kohfeld and Ridgwell, 2009;Sigman and Boyle, 2000). By taking the 1,290 Pg C difference between the O-PI and O-PI LGM CO2 experiments and assuming an a priori addition of 520 Pg C in a glacial ocean (Ciais et al., 2011), we can attribute approximately 18 ppm to our simulated solubility changes (Table 4). Our results are therefore roughly consistent with other estimates. Clearly though, the tendency for carbon to be released to the atmosphere would have been reduced by storing more carbon in the deep ocean.
We therefore note that circulation changes and cooling would have had a complementary effect on carbon sequestration, and magnified their individual effects.

5
Regarding sea ice, a prevailing view is that sea ice expansion would enhance oceanic carbon storage by limiting air-sea gas exchange (Stephens and Keeling, 2000). This theory largely focuses on the restriction of outgassing from carbon-rich deep waters that upwell in the Southern Ocean. However, this neglects responses in the Northern Hemisphere and in the biological pump. In a seminal paper, Marinov et al. (2006) showed that export production and circulation in the Antarctic zone has a strong effect on atmospheric CO 2 . We found that light limitation of highly productive regions in both hemispheres caused by 10 increased sea ice cover led to a global loss of 160 Pg C, equivalent to 8 ppm pCO 2 . Similar responses have been simulated in other models that consider the impact of sea ice on biological production and do not consider temperature and salinity changes associated with sea ice growth (Kurahashi-Nakamura et al., 2007;Sun and Matsumoto, 2010).
Our results demonstrated that although the storage of carbon was enhanced in the glacial ocean (668 Pg C) due to physi-15 cal changes, namely due to the overturning circulation, they could not overcome the loss of carbon (1,290 Pg C) caused by equilibration with a lower atmospheric pCO 2 .

Nutrients and export production
Like carbon, phosphate (PO 4 ) concentrations in experiment O-LGM were redistributed from the upper to the deep ocean 20 (Fig. 5). This change was driven by a strengthened AABW formation and is consistent with proxy reconstructions. Cadmium and δ 13 C measurements from the Atlantic Ocean show increased nutrient concentrations in the deep ocean, but reduced levels above 2,000 m at the LGM (Boyle, 1992;Gebbie, 2014;Marchitto and Broecker, 2006;Tagliabue et al., 2009).
A direct consequence of the redistribution of PO 4 was the reduction in the production of particulate organic matter across 25 many regions of the O-LGM ocean (Fig. 5). With the exception of the South Pacific and isolated areas in the subtropics, export production in the O-LGM experiment decreased relative to the O-PI experiment, so that global export production was 56 % of O-PI. The global reduction was also illustrated by a decrease in regenerated carbon (C org ), which indicated that the biological carbon pump was weakened (Table 3). The reduction in export production and regenerated carbon for the O-LGM experiment is significant when compared with other studies that argue for a more efficient glacial biological pump than that of the Holocene 30 (Galbraith and Jaccard, 2015;Schmittner and Somes, 2016). The weakened efficiency of our simulated biological pump can be attributed, in part, to a large decrease in export production from the sub-Antarctic zone. This feature is in direct conflict with palaeoproductivity proxies in the Atlantic and Indian sectors of the sub-Antarctic Ocean (Anderson et al., 2002, 2014Chase et al., 2001;Jaccard et al., 2013;Nürnberg et al., 1997), and some parts of the Pacific sector (Bradtmiller et al., 2009;Lamy et al., 2014). Outside of the Southern Ocean, the reduction in export production in the O-LGM experiment is largely consistent with palaeoproductivity evidence (see Introduction).

Carbonate chemistry
The loss of phosphate from the upper ocean and its increase at depth was mirrored by changes in alkalinity and salinity, so that 5 a more alkaline and saline signature of AABW, relative to NADW, dominated the deep ocean in experiment O-LGM. Alkalinity and salinity decreased by 66 mmol Eq m −3 and 0.37 psu in the surface ocean and increased by 147 mmol Eq m −3 and 2.71 psu in the deep ocean (Fig. 6).
Because the majority of LGM-PI change in salinity occurred in the deep ocean, the changes in carbonate chemistry across 10 the surface ocean were small. Little change between experiments O-PI and O-LGM was found in the aragonite saturation state (Ω ar ), which is a unitless index indicating under-and super-saturation at values below and above one (Fig. 7). Surface  However, the magnitude of increase in alkalinity in the glacial deep ocean, which was not accompanied by a stoichiometrically matched increase in carbon, caused unrealistic increases in the calcite saturation horizon (Ω ca = 1). The average position of the calcite saturation horizon increased from 2666 metres in O-PI to 3208 metres in O-LGM. While this increase may seem modest, the increase exceeded the ocean floor over the majority of the ocean, so that seawater was completely saturated for calcite outside of the eastern tropical Pacific (Fig. 8). There is good evidence that the carbonate chemistry of the LGM ocean 25 was not appreciably different to the Late Holocene (Yu et al., 2014), and this information places our simulated changes in deep ocean Ω ca at the LGM as unrealistic.  (Garcia, 2005). The combination of cooler SSTs, an enhanced subduction of AABW and the reduction in export production in experiment O-LGM dramatically increased the oxygen content in both the upper and deep ocean by ∼80 and ∼120 mmol m −3 , respectively, which constituted a global increase of 55 % ( Fig. 9; Table 3).

Dissolved oxygen
The increase in dissolved oxygen in O-LGM was considerable, but agreed well with proxy reconstructions for the upper ocean. The oxygen-poor intermediate waters of the western North Pacific (Ishizaki et al., 2009;Shibahara et al., 2007), eastern North Pacific (Cannariato andKennett, 1999;Cartapanis et al., 2011;Chang et al., 2014;Dean, 2007;Nameroff et al., 2004;Pride et al., 1999;Ohkushi et al., 2013;van Geen et al., 2003), eastern South Pacific (Martinez et al., 2006;Muratli et al., 5 2010; Salvatteci et al., 2016), Equatorial Pacific (Leduc et al., 2010 and Indian Ocean (Reichart et al., 1998;Suthhof et al., 2001;van der Weijden et al., 2006) were better oxygenated at the LGM relative to the PI climate. An important consequence of oxygenating the upper ocean is a reduction in the strength of denitrification in the these regions. Sedimentary δ 15 N records suggest that global aggregate rates of water column denitrification rates over the past 200,000 years were lower during glacial periods and higher during interglacial periods (Galbraith et al., 2004), and this is consistent with the simulated oxygenation of

Importance of ocean biogeochemistry for climate
The carbon content, export production field, carbonate chemistry, and deep ocean oxygen content of experiment O-LGM are outstanding in their disagreement with proxy evidence. Notably, 622 Pg C was lost from experiment O-LGM relative to O-PI. 25 The standard LGM simulation was therefore unable to explain the glacial-interglacial drawdown of atmospheric CO 2 , despite the existence a physical ocean state within realistic bounds. If we are to reconcile the biogeochemistry of the glacial ocean with that inferred from proxy evidence, we must therefore consider altering ocean biogeochemistry. 30 Three plausible modifications to ocean biogeochemistry (see methods) were considered: (1) increased POC export production,

Reconciling the carbon budget
(2) increased depth of POC remineralisation, and (3) reduced PIC export. In the following we step through the changes to carbon content caused by each modification, and the reader is directed to Table 3 for reference.
(1) Experiment O-LGM BGC poc . Although the scaling factor controlling the export production of organic matter was increased 10-fold, the actual increase in POC export production averaged over the global ocean was more modest at roughly 30 %.
Because most of the ocean became phosphate limited as greater quantities of nutrients were redistributed into the deep ocean, the increase in export production in experiment O-LGM BGC poc was only felt in those regions where PO 4 was not limiting. The 5 sub-Antarctic zone of the Southern Ocean experienced the greatest increase in export production (∼250 %), followed by a few small regions along the Chilean margin and in the Northwest Pacific (Fig. 10). These regional responses caused the global net export production rate to increase from 4.48 to 5.96 Pg C yr −1 . Although this rate of POC export production was still lower than the O-PI experiment of 8.02 Pg C yr −1 , this increased carbon content by 179 Pg C.

10
(2) Experiment O-LGM BGC rem . The shift of organic matter to depth was associated with a global reduction in POC export production of ∼1.2 Pg C yr −1 as remineralisation released PO 4 and carbon further from the photic zone. Despite the reduction in the biological pump, the bulk transfer of POC to depth generated an increase in ocean carbon storage of 167 Pg C. This magnitude of increase places our glacial ocean state within the plausible bounds required to offset the loss of atmospheric and terrestrial carbon reported by Ciais et al. (2011) of ∼520 ± 400 Pg C at the LGM (Table 4). By assuming a glacial-interglacial difference in atmospheric CO 2 of 95 ppm and applying this to our changes in carbon content, we attribute 25 roughly 40 ppm to changes in ocean physics and 55 ppm to changes in the biological pump. Within the physical changes, 28 ppm is attributed to the reorganisation of the global overturning circulation, while 12 ppm can be attributed to changes in surface properties, including sea ice expansion, cooling and salinification.

30
Of the three biogeochemical modifications applied to the LGM ocean, only two had an effect on POC export, as the amount of PIC exported from the photic zone has no influence on the amount of POC export. Deepening the remineralisation of POC (O-LGM BGC rem ) shifted a greater fraction of regenerated PO 4 into the deep ocean, which resulted in a global reduction of export production. Increasing the scaling factor (O-LGM BGC poc ), however, caused an increase in global export production from 4.48 to 5.96 Pg C yr −1 . Most of this increase occurred in the Southern Ocean, particularly the sub-Antarctic zone, and in a few isolated pockets in the Northwest Pacific and North Atlantic where excess nutrients were available (Fig. 10).
The increase in the scaling factor dominated the change in export production produced when combining all three biogeochemical modifications (O-LGM BGC all ). The strong increase in export production observed in the sub-Antarctic was clearly 5 replicated within this experiment and reconciles our simulated export production field with current evidence of productivity at the LGM. In the Southern Ocean, the Atlantic and Indian sectors of the sub-Antarctic zone experienced a greater flux of organics (Anderson et al., 2002Chase et al., 2001;Jaccard et al., 2013;Nürnberg et al., 1997). Whether this was also the case for the Pacific sector remains under debate, with some evidence for increase (Bradtmiller et al., 2009;Lamy et al., 2014) and some for no change (Bostock et al., 2013;Chase et al., 2003). Meanwhile, it is widely accepted that waters south of

15
In experiment O-LGM BGC all , net export production remained less than the O-PI experiment by 3.16 Pg C yr −1 despite the application of biogeochemical modifications. The weakened carbon transfer to the interior ocean was also observed in the regenerated carbon content of the ocean (C org ), which was 646 Pg less than the O-PI experiment ( Table 3). The net decline in export production observed in this study was dominated by the decline in tropical and subtropical waters. Many palaeoproductivity studies located outside of the sub-Antarctic zone have found weakened productivity at the LGM (Chang et al., 2014, 20 2015; Costa et al., 2016;Crusius et al., 2004;Kohfeld et al., 2005;Kohfeld and Chase, 2011;Jaccard et al., 2005;McKay et al., 2015;Ortiz et al., 2004;Riethdorf et al., 2013;Salvatteci et al., 2016;Singh et al., 2011;Thomas et al., 1995). Additionally, an enhanced utilisation of available nutrients in the sub-Antarctic zone (Martinez-Garcia et al., 2014) would reduce the nutrient content of intermediate waters formed in the Southern Ocean and would thus reduce the delivery of nutrients to lower latitudes (Sarmiento et al., 2004). This mechanism coupled with cooler temperatures caused reductions in export production across 25 much of the mid and lower latitude oceans in experiment O-LGM BGC all , which maintains the qualitative agreement between the simulated export production field and proxy observations. However, the global weakening of the biological pump in our simulations is contrary to proxy and model-based evidence for a strengthened biological pump (> C org ) at the LGM (Galbraith and Jaccard, 2015;Schmittner and Somes, 2016). Hence, while we present good spatial agreement between O-LGM BGC all and palaeoproductivity proxies at the LGM, which was essential to increasing the carbon content of the ocean, we note that addi-30 tional increases in the export production field may be valid.

Reconciling carbonate chemistry
There is good evidence that the carbonate chemistry of the glacial deep ocean was not appreciably different to the Late Holocene (Yu et al., 2014). Because much of the ocean was saturated for calcite in the O-LGM experiment, additional processes were required to shoal the depth at which calcite becomes unsaturated (Ω ca = 1), and thereby reconcile proxy evidence.
One mechanism to reduce deep ocean Ω ca would be to reduce continental inputs of alkalinity at the LGM. However, the presence of glaciers, drier atmospheric conditions and the exposure of continental shelves due to lower sea level would have 5 increased the supply of carbonates to the ocean (Gibbs and Kump, 1994;Riebe et al., 2004), thereby increasing ocean alkalinity and further deepening the carbonate saturation horizon. This mechanism has been largely refuted as having a significant effect on the glacial-interglacial difference in the carbon budget (Brovkin et al., 2007;Foster and Vance, 2006;Jones et al., 2002), and can therefore be ignored. 10 The individual biogeochemical modifications were also insufficient to reduce Ω ca to be consistent with palaeo evidence.
However, combining all three modifications in experiment O-LGM BGC all reduced Ω ca significantly (Fig. 11), and produced a globally-averaged calcite saturation horizon at 2839 metres. Remarkably, our simulated PI to LGM changes in Ω ca are also consistent with palaeoproxy reconstructions in a regional sense. The calcite saturation horizon in the Pacific Ocean was deeper in experiment O-LGM BGC all relative to O-PI, but was shallower in the Atlantic Ocean and within the Atlantic and Indian sectors 15 of the sub-Antarctic zone. Similar changes are seen in the proxy record, with a deepening of less than 1,000 metres in the North Pacific and Southern Ocean at the LGM (Anderson et al., 2002;Catubig et al., 1998), and a shoaling in the Atlantic Ocean (Anderson et al., 2002).
An important caveat of this study, which cannot be ignored, is the exclusion of calcium carbonate (CaCO 3 ) burial within 20 ocean sediments. Because this process is not included in the model, it is highly likely that the deepening of the calcite saturation horizon that occurred in the standard LGM experiment (O-LGM) was too extreme. CaCO 3 burial lowers the alkalinity of the glacial ocean, and is therefore a negative feedback mechanism that modulates changes in carbonate chemistry (see Sigman et al., 2010, for a review). If calcite saturation increases through the water column, as found in O-LGM, the burial of CaCO 3 would increase and subsequently mitigate the rise in whole-ocean alkalinity, which would in turn the rise in calcite saturation 25 and the ability of the ocean to store carbon. By not taking this process into account, both the deepening of the calcite saturation horizon and the storage of carbon were overestimated in experiment O-LGM. However, the same reasoning can be applied to the experiments with biogeochemical modifications. If the burial of CaCO 3 was included in these experiments, the shoaling of the calcite saturation horizon would have been somewhat mitigated by 30 decreased CaCO 3 burial that increased ocean alkalinity (Archer and Maier-Reimer, 1994). Consequently, the shoaling that was observed in these experiments was likely exaggerated, just as the deepening observed in experiment O-LGM was exaggerated.
Again, this can be applied to changes in the carbon content of the ocean, as a shallower calcite saturation horizon would have increased whole-ocean alkalinity and thereby increased carbon storage. This effect would have been particularly important for experiment O-LGM BGC pic , where inorganic carbon export was eliminated. If whole-ocean alkalinity was able to respond  Fig. 12; Table 3), and therefore did not compromise the agreement between simulated and proxy oxygen reconstructions. 10 Modifying ocean biogeochemistry did, however, have a large effect on the oxygen concentrations of the deep ocean (Fig. 13).
Increasing export production (O-LGM BGC poc ) and deepening the remineralisation depth (O-LGM BGC rem ) both reduced oxygen concentrations by 28 and 13 mmol m −3 , respectively. The combination of these modifications (O-LGM BGC all ) amplified their individual effects, so that deep ocean oxygen was reduced by 63 mmol m −3 relative to O-LGM. The increased sensitivity of deep ocean oxygen to the combination of increased export production and a deeper remineralisation depth was also observed in 15 the increase in the quantity of regenerated nutrients (C org ) that resulted (Table 3). A greater proportion of regenerated nutrients relative to preformed nutrients at the LGM has been identified as a key driver of interior ocean deoxygenation (Jaccard and Galbraith, 2012;Sigman et al., 2010) LGM simulation (see sections 3.1.2 and 3.1.3), and we therefore suggest that a sluggish circulation is necessary for reducing deep ocean oxygen. Key biogeochemical mechanisms include a further increase in global export production, and/or a different spatial pattern of export production, and/or increasing the injection of organic matter to deep water via further lengthening the remineralisation profile. The weakened biological pump in our LGM simulations contrasts with other studies (Galbraith and combination of a more sluggish deep ocean circulation with an enhanced export of organics would significantly reduce the oxygen content of the deep ocean, while potentially further increasing carbon storage.
Thomas, E., Booth, L., Maslin, M., and Shackleton, N. J.: Northeastern Atlantic benthic foraminifera and implications of productivity during the last 40,000 years., Paleoceanography, 10, 545-562, 1995.   LGM experiment, and is broadly consistent with the results of other PMIP3 models. In panel (d), the coloured markers represent locations were winter sea ice was deemed to have been present (blue) and absent (red) at the LGM according to Gersonde et al. (2005). Table 2. Changes in sea surface temperature (SST), sea ice extent and large scale circulation between the LGM and PI simulations. SST changes are compared with proxy reconstructions of SST generated by Waelbroeck et al. (2009) and Ho and Laepple (2015) whom use different proxies for their reconstructions to produce the differences depicted here and discussed in the text. The values of AABW and NADW formation provided by the models from the Palaeoclimate Modelling Intercomparison Project Phase II (PMIP2) are presented by Otto-Bliesner et al. (2007). The transport rate of the Antarctic Circumpolar Current (ACC) for the PMIP2 models is presented by Lynch-Stieglitz et al. (2016). The estimates of NADW production provided by the PMIP3 models were taken from Muglia and Schmittner (2015).     LGM (solubility) -unmodified BGC −941 +21 LGM (sea ice) -unmodified BGC −1450 −9 LGM (circulation) -unmodified BGC −857 +28 LGM (all physical changes) -unmodified BGC −622 +40 LGM (biological pump) -modified BGC −351 +55 a Estimate of increase in ocean carbon content during the LGM made by Ciais et al. (2011), whereby atmospheric carbon was reduced by 194 ± 2 Pg C and terrestrial carbon was reduced by 330 ± 400 Pg C. b Assumes all three biological modifications that were postulated (see Table 1         It should be noted that the export production field of particulate organic carbon for experiment O-LGM BGC pic , whereby particulate inorganic carbon was set to zero, did not differ from unmodified experiment O-LGM and is therefore not shown. For this comparison, the reader is directed to Figure 5.  that the oxygen field for experiment O-LGM BGC pic , whereby particulate inorganic carbon was set to zero, did not differ from the unmodified glacial experiment O-LGM and can therefore be used here as a reference to that simulation. It should be noted that the oxygen field for experiment O-LGM BGC pic , whereby particulate inorganic carbon was set to zero, did not differ from the unmodified glacial experiment O-LGM and can therefore be used here as a reference to that simulation.