CO 2 drawdown due to particle ballasting and iron addition by glacial aeolian dust : an estimate based on the ocean carbon cycle model MPIOM / HAMOCC version 1 . 6 . 2 p 3

Despite intense efforts, the mechanisms that drive glacial–interglacial changes in atmospheric pCO2 are not fully understood. Here, we aim at quantifying the potential contribution of aeolian dust deposition changes to the atmospheric pCO2 drawdown during the Last Glacial Maximum (LGM). To this end, we use the ocean circulation and carbon cycle model MPIOM/HAMOCC, including a new parameterisation of particle ballasting that accounts for the acceleration of sinking organic soft tissue in the ocean by higher density biogenic calcite and opal particles, as well as mineral dust. Sensitivity experiments 5 with reconstructed LGM dust deposition rates indicate that the acceleration of detritus by mineral dust likely played a small role for atmospheric pCO2 variations during glacial–interglacial cycles – on the order of 5 ppmv, compared to the reconstructed ∼80 ppmv-rise of atmospheric pCO2 during the last deglaciation. The additional effect of the LGM dust deposition, namely the enhanced fertilisation by the iron that is associated with the glacial dust, played a more important role – leading to a pCO2drawdown by more than 8 ppmv in our LGM sensitivity experiments despite an underestimated iron-limitation in the used 10 model setup.


Introduction
According to ice core data (Lüthi et al., 2008;Köhler et al., 2017), the atmospheric CO 2 concentration during the Last Glacial Maximum (LGM, ∼21 thousand years ago) was about 80 ppmv lower than during the early Holocene.Model simulations show that the reduced atmospheric CO 2 concentration caused a global cooling, and, combined with orbital effects on Earth's climate, allowed the Laurentide, Cordilleran, British and Scandinavian ice sheets to build up (Abe-Ouchi et al., 2007;Heinemann et al., 2014;Ganopolski and Brovkin, 2017).It is still a matter of active discussion, however, which of the many potential climate and carbon cycle changes in response to the orbital forcing caused how much of the atmospheric pCO 2 drawdown (e.g., Brovkin et al., 2012;Chikamoto et al., 2012;Ganopolski and Brovkin, 2017).
Here we introduce a new parameterisation of the effect of higher density particles such as calcite, opal and dust on the sinking speed of detritus into the ocean circulation and carbon cycle model MPIOM/HAMOCC.Based on this new model setup, we estimate to what extent enhanced aeolian dust deposition contributed to the reconstructed atmospheric pCO 2 drawdown during the last glacial period due to 1) its effect on particle ballasting, and 2) the enhanced iron fertilization associated with the dust.Observationally constrained simulations of the global dust cycle have suggested that dust production and deposition rates were at least twice as large during the LGM compared to the Holocene and present (e.g., Werner et al., 2002;Mahowald et al., 2006;Albani et al., 2016).The increased aeolian dust emissions were due to both, enhanced desert and enhanced glacigenic dust production.The increased desert dust production was caused mainly by a reduction in terrestrial vegetation (leading to larger desert dust source areas) due to inverse CO 2 -fertilisation, temperature, precipitation and cloudiness changes (Mahowald et al., 2006) or by stronger tropospheric winds (Werner et al., 2002).Glacigenic dust, which is produced by the grinding of ice sheets and glaciers over the bedrock and sediments, and the subsequent transport of the produced fine sediment by meltwater beyond the ice margin (e.g., Bullard, 2013), contributed to the enhanced dust deposition, and globally, the magnitude of this glacigenic dust contribution was comparable to that of desert dust (Mahowald et al., 2006).
It is well known that enhanced glacigenic and desert dust deposition over the oceans can strengthen the biological carbon pump due to fertilisation with the bio-available iron contained in the dust (Martin, 1990;Kumar et al., 1995;Sigman and Boyle, 2000;Werner et al., 2002).This iron hypothesis is supported for example by mesoscale iron enrichment experiments which demonstrated that primary production in one third of the world ocean is limited by iron-availability (Boyd et al., 2007).
Various modelling studies testing the iron hypothesis have concluded that a range of 15-40 ppmv of the reconstructed 80 ppmv atmospheric pCO 2 difference can be explained by iron fertilisation due to dust deposition changes (Watson et al., 2000;Bopp et al., 2003;Brovkin et al., 2007;Lambert et al., 2015;Ganopolski and Brovkin, 2017).
Enhanced dust deposition can also have an effect on the sinking speed and remineralisation length scale of particulate organic carbon (POC or detritus).Detritus sinks because of its higher density compared to seawater (excess density).As aggregates including POC form, ballast minerals, i.e., silicate and calcite biominerals, as well as lithogenic minerals such as aeolian dust, which are more dense than organic tissue and seawater, provide a major weight fraction of sinking particles in the deep ocean (Honjo et al., 1982).Hence, ballast minerals contribute largely to the excess density of aggregates including detritus in the deep ocean.This led to the hypothesis that ballast minerals determine accelerated deep-water POC fluxes (ballasting hypothesis; e.g., Armstrong et al., 2002).Several deep-moored sediment trap observations of particle flux rates supported this ballasting hypothesis (see Boyd and Trull, 2007, for a review).While, e.g., Klaas and Archer (2002) and François et al. (2002) confirmed the association of calcite minerals with efficient transfer of POC to the deep ocean but suggested that the role of lithogenic minerals is negligible (see also Berelson, 2001), a more recent study indicated that neglecting lithogenic material in carbon cycle models would lead to an underestimation of the POC flux to the seafloor by 16-51 % globally (Dunne et al., 2007).Noticeably high sinking speeds off Northwest Africa are also in line with the hypothesis that the deposition of Saharan dust accelerates the sinking of detritus (Fischer and Karakaş, 2009).
Still, to the best of our knowledge, the hypothesis that particle ballasting changes due to variations of fluxes of lithogenic material into the ocean contributed to the reconstructed atmospheric pCO 2 changes during glacial-interglacial cycles (Ittekkot, 1993) has not yet been tested using ocean-biogeochemical models.This manuscript is organised as follows.In Section 2, we describe our configuration of the ocean circulation and carbon cycle model MPIOM/HAMOCC (Section 2.1), introduce a simple particle ballasting scheme that includes the effect of lithogenic minerals from aeolian sources (Section 2.2), and describe the spinup procedure of the model (Section 2.3).In Section 3, we de-2 Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-137Manuscript under review for journal Geosci.Model Dev. Discussion started: 16 August 2018 c Author(s) 2018.CC BY 4.0 License.scribe the necessary adjustment of HAMOCC parameters to the modified sinking speeds (Section 3.1) and compare the resulting reference simulation with the new particle ballasting to a reference simulation with the standard prescribed particle sinking speeds, and to observations (Section 3.2).In Section 4, to estimate the effect of dust deposition changes on glacial-interglacial atmospheric pCO 2 variations, sensitivity experiments with LGM instead of modern dust deposition rates are presented.We conclude with a summary and discussion of our results (Section 5).

Methods
To quantify the effect of aeolian dust variations on atmospheric pCO 2 , we use the Max-Planck-Institute Ocean Model MPIOM and the embedded Hamburg Ocean Carbon Cycle model HAMOCC (Section 2.1).To account for particle ballasting effects, the default prescribed Martin-type sinking for detritus and the constant sinking speeds for opal, dust and calcite are replaced by a simple prognostic particle ballasting parameterisation, which assumes that detritus and the heavier particles form aggregates that sink at speeds computed from the density differences between these imaginary aggregates and the surrounding seawater (Section 2.2).Using the new particle ballasting parameterisation, we then perform sensitivity tests to quantify the effect of dust deposition changes on the atmosphere-ocean CO 2 fluxes, separating the effects of 1) particle ballasting by the dust, and 2) fertilisation by the bio-available iron associated with the dust.The presented sensitivity tests are equilibrium simulations with prescribed modern or LGM dust deposition fields, atmospheric climatological forcing representing modern conditions, and pre-industrial greenhouse gas concentrations (GHGs).The atmosphere-ocean CO 2 flux anomalies in the simulations with LGM dust are used to estimate the effect of dust changes on atmospheric pCO 2 variations during glacial cycles.We do not perform transient simulations of glacial cycles yet, because 1) dust deposition fields are currently only available for snapshots during the last deglaciation, and 2) the computational cost and required wall-time to transiently simulate a glacial cycle or inception is currently too large.

Configuration of MPIOM/HAMOCC
Our results are based on the ocean-only model configuration of the Max-Planck-Institute Earth System Model (MPI-ESM; Giorgetta et al., 2013) version 1.2.00p4, which contains MPIOM version 1.6.2p3(corresponding to revision 3997).MPI-ESM is a comprehensive Earth system model that consists of the atmosphere general circulation model ECHAM (Stevens et al., 2013), the land surface scheme JSBACH (Raddatz et al., 2007), and the ocean circulation model MPIOM (Jungclaus et al., 2013), which contains the ocean biogeochemistry model HAMOCC (Ilyina et al., 2013;Paulsen et al., 2017).MPI-ESM has been used successfully for millennial-scale atmosphere-ocean-carbon cycle simulations, including a transient simulation of the last millennium with fully interactive atmospheric CO 2 (Jungclaus et al., 2010).Since our aim is to isolate the effect of particle ballasting changes in the ocean during glacial-interglacial cycles on the atmosphere-ocean CO 2 fluxes, rather than realistically simulating past climate evolution, we use the ocean-only configuration of MPI-ESM, namely MPIOM/HAMOCC in stand-alone mode (Marsland et al., 2003;Maier-Reimer et al., 2005;Paulsen et al., 2017).However, we do consider variable atmospheric pCO 2 with a simple atmospheric CO 2 box model based on atmosphere-ocean CO 2 flux anomalies (Section 4.2).This way, atmosphere and land surface feedbacks to the pCO 2 variations -that would alter the ocean circulation and complicate the analysis of dust-related CO 2 flux anomalies -are avoided, and computational time is saved.
Another measure to save computational time is the use of the pre-defined, coarse, bi-polar ocean model grid setup GR30L40.
The grid poles are located over Greenland and Antarctica; the curvilinear grid has a grid-point spacing of approximately 3 • , and 40 levels with thicknesses ranging from 10-12 m in the upper ocean to 600 m in the deepest ocean.
The ocean model MPIOM is forced by a repeating 1-year-long climatology of daily mean surface boundary conditions, including 2 m air temperature, zonal and meridional windstress, windspeed, shortwave radiation, total cloud cover, precipitation and river runoff, representing modern conditions based on the ERA40 and ERA15 reanalyses (Röske, 2005(Röske, , 2006) ) as in Paulsen et al. (2017).
To avoid salinity trends due to imbalances of precipitation, evaporation and river runoff, the surface salinity is restored with a relaxation coefficient of 3.3•10 −7 s −1 to the annual salinity field of the Polar science center Hydrographic Climatology (PHC).
Note that biogeochemical changes do not affect the ocean circulation in the default (and our) configuration of MPIOM/HAMOCC.That means, for example, any effect of phytoplankton concentration changes on light absorption (Wetzel et al., 2006) is neglected in MPIOM.Differences in the ocean circulation between our simulations with the standard MPIOM/HAMOCC setup and the setup that includes particle ballasting are only due to the fact that, unfortunately, different numbers of CPUs were used for the computations (Fig. 1a: black lines cover green and brown lines, but not the dark grey line of the setup without ballasting).For binary-identical results, the parameter 'nprocx' must not be changed.Otherwise, the spatial and temporal variations of the ocean circulation would be identical in the spinup with the standard model and the simulations with ballasting.

Particle sinking and ballasting in HAMOCC
In the standard version of HAMOCC, particle ballasting is not accounted for.Calcite and opal shells sink at a prescribed, constant speed of 30 m day −1 ; dust sinks at a speed of only 0.05 m day −1 , resulting from the assumption of spherical dust particles with a diameter of 1 µm following Stokes' law (Stokes' drag is balanced by gravitational force minus buoyancy force); the sinking speed of detritus w det is constant in the euphotic zone, and increases with depth z below the euphotic zone (z > 100 m) following Martin et al. (1987): where w eu = 3.5 m day −1 is the sinking speed in the euphotic zone, the parameter b is set to 2, and the remineralisation rate for detritus λ det is set to 0.0025 day −1 , resulting in sinking speeds up to 80 m day −1 (Fig. 2).
To account for particle ballasting, we add a simple parameterisation of particle sinking speeds to HAMOCC.We assume that, below the euphotic zone, all particulate matter forms aggregates and sinks together at a common speed w ballast .As suggested by Gehlen et al. (2006), the common speed is proportional to the difference between the mean density of all particles ρ particles w 0 is the proportionality factor, and can be interpreted as the sinking speed of the aggregate in the absense of particle ballast.
The density difference ρ particles − ρ sw is also called excess density.The mean density of all particles ρ particles is given by the mass of all particles per unit volume of seawater, divided by the volume of all particles per unit volume of seawater.In contrast to Gehlen et al. (2006), here, dust also contributes to the ballasting of detritus: where c dust is the mass concentration of dust, ρ dust is the density of dust (2.50 g cm −3 ), Ψ b are the molar concentrations of the biogenic particle types b (namely: detritus, calcite, and opal), m b are their molar masses (32.7 g (mol C) −1 for detritus, 100.0 g (mol C) −1 for calcite, and 72.8 g (mol Si) −1 for opal) and ρ b are their densities (1.06 g cm −3 for detritus, 2.71 g cm −3 for calcite, and 2.10 g cm −3 for opal).
Monthly mean dust deposition rates are prescribed based on numerical reconstructions for modern and LGM boundary conditions, including glacigenic dust (Fig. 5a-b; Mahowald et al., 2006).The dust deposition fields not only affect the sinking speed in our ballasting parameterisation, but also the availability of iron.It is assumed that 3.5 % (weight) of the dust in the water column is iron, of which 1 % is soluble and bio-available.

Spinup of the standard model setup
Before starting our experiments with particle ballasting, we spin-up the standard version of MPIOM/HAMOCC without particle ballasting for over 10,000 years (Fig. 1).During the spin-up phase, the ocean model, which is initialised from three-dimensional PHC temperature and salinity fields, approaches an equilibrium state (thanks to the salinity-restoring at the surface).The simulated ocean circulation is characterised by deepwater formation in the Atlantic and Southern Ocean, and no deepwater formation in the Pacific; the Atlantic meridional overturning at 26 • N amounts to about 17 Sv (Sverdrup; 10 6 m 3 s −1 ; Fig. 1a).
The ocean carbon cycle is not so easily tamed, since it remains an open system in our model (no restoring).To reach a stable equilibrium of the carbon cycle in the water column and the sediment, the net fluxes of organic carbon, calcite and opal from the ocean column to the sediment -which are computed in HAMOCC -need to be balanced by the fluxes into the ocean of organic matter, calcite and silicate from weathering on land -which are prescribed as namelist parameters (called deltaorg, deltacal, and deltasil).The weathering fluxes are adjusted several times during the spinup phase, aiming at 1) balancing the sediment burial fluxes as diagnosed from the trends in the sediment inventory over several hundred years prior to each weathering flux adjustment, 2) a stable calcite export production at 90 m depth of about 1 Gt C yr −1 , and 3) a rain ratio (ratio of calcite export to detritus export at 90 m depth) of about 10 % (Fig. 1f).An equilibrium of the sediment, however, is not reached in our simulations.Note that the selection of the weathering fluxes itself affects the simulated biogeochemistry.For example, if the silicate weathering flux is increased, more opal can be produced, and calcite production is suppressed (see Eqs. 12-14 of Ilyina et al., 2013).do not intend to interpret the biogeochemistry of this run.However, the endpoint of this older simulation is used as a starting point for our modern control simulation without particle ballasting, using the default-configuration of MPIOM/HAMOCC in MPI-ESM 1.2.00p4, including a new parameterisation of diazotrophs (Paulsen et al., 2017) and Martin-type sinking (darker grey).In model year 6500, the ballasting parameterisation is switched on, and, to compensate for the reduced sinking speed of opal in the euphotic zone, the opal dissolution rate λ opal (namelist parameter dremopal) is reduced from 0.01 day −1 by factors of 2, 4, 6 and 8.

Results for modern dust deposition
To simulate a reasonably realistic carbon cycle with atmosphere-ocean CO 2 fluxes comparable to observations, and to create a mean state with particle ballasting that is very close to the mean state of the standard setup, it is important to adjust HAMOCC to the substantially modified particle sinking speeds that result from the particle ballasting parameterisation.In this section, we argue that this adjustment can be achieved by a reduction of the opal dissolution rate (Section 3.1) and then compare

Adjusting HAMOCC to the new particle sinking
The simulations with particle ballasting are initialized from year 6500 of the reference simulation without particle ballasting (Fig. 1b).Using a no-ballast reference speed w 0 (Eq.2) of 0.5 m day −1 leads to global mean sinking speeds w ballast of the virtual aggregates between about 5 m day −1 at the surface and about 120 m day −1 below 5 km depth (Fig. 2a).A larger w 0 , such as 3.5 m day −1 as used by Gehlen et al. (2006), would lead to much higher global mean sinking speeds.The advantage of using the small w 0 is that it yields global mean aggregate sinking speeds that are comparable to the depth-dependent, Martintype sinking speed of detritus that is prescribed in the standard control simulation (Fig. 2a), which allows us to avoid a re-tuning of the remineralisation rate of detritus.
The dissolution rate of opal, however, needs to be re-tuned, since the sinking speed of opal is reduced from previously 30 m day −1 to now only ∼5 m day −1 in the euphotic zone (Fig. 2a).Consequently, more opal is dissolved within the euphotic 10 zone, leading to higher silicate concentrations, a higher opal production rate (Fig. 1e), and a very low rain ratio (Fig. 1f) due to a completely suppressed calcite production.Without adjusting the river input (weathering) of calcite to the reduced calcite precipitation and burial, surface alkalinity and DIC rise (Fig. 1b-c).For a reduction of the opal dissolution rate by a factor of 6, the trends in alkalinity and DIC become small, and comparable to (or, actually smaller than) the trend in the control simulation without ballasting (Fig. 1b-c), and, even without re-adjusting the weathering fluxes, a stable carbon cycle simulation with particle ballasting is achieved that is similar to the control simulation without ballasting.

Ballasting effects and comparison to observations
For pre-industrial GHGs (pCO 2 =278 ppmv), modern ocean circulation, and modern dust deposition (Fig. 5a), the simulated atmosphere-ocean CO 2 fluxes in the standard setup and the setup with particle ballasting are hardly distinguishable (Fig. 3b-c).
There are, however, clear differences compared to the observed modern CO 2 fluxes (Fig. 3a).Most notably, the observed flux of CO 2 from the atmosphere into the ocean between 30 • S and 60 • S in the South Atlantic, Indian Ocean and West Pacific is underestimated, due to an overestimation of the ocean surface pCO 2 , as previously reported for HAMOCC within the coupled MPI-ESM (Ilyina et al., 2013).Also note that the CO 2 flux data (Takahashi et al., 2009) is based on modern observations, while in the model simulations, pre-industrial GHGs are prescribed.
Similarly, the surface dissolved phosphate concentrations are hardly distinguishable in the simulations with and without ballasting.Compared to observations, the simulated phosphate concentrations in the Southern Ocean, Indian Ocean, Equatorial Pacific, and in the North Pacific are underestimated (Fig. 3d-f), and these biases are again similar to previous HAMOCC studies (Ilyina et al., 2013).
The global primary production in the simulation with ballasting amounts to about 39 Gt C yr −1 compared to 44 Gt C yr −1 in the simulation without particle ballasting (Fig. 1d), not including the primary production by diazotrophs, which amounts to additional 5 Gt C yr −1 and 4 Gt C yr −1 for the simulations with ballasting and the standard setup, respectively.The global organic matter export is slightly reduced (7.0Gt C yr −1 compared to 7.5 Gt C yr −1 without ballasting), and the global calcite production and export are only moderately larger (1 Gt C yr −1 compared to 0.85 Gt C yr −1 without ballasting, not shown), leading to an elevated rain ratio of about 13 % compared to 11 % in the simulation without ballasting (Fig. 1f).
The simulated carbon export in 100 m depth (export production) in the standard model setup and the setup with ballasting are very similar to each other, with the area of high export production in the Equatorial Pacific extending further into the subtropics in the standard setup (Fig. 3h-i).Compared to estimates based on sediment trap data (Honjo et al., 2008) and satellite data (Henson et al., 2012), the simulated export production is underestimated in the North Atlantic, and overestimated in the Equatorial Pacific (Fig. 3g-i).
Only a fraction of this export production from the euphotic zone reaches the deep ocean.This fraction, also called transfer efficiency, is a measure for the efficiency of the biological pump.Compared to estimates based on sediment trap data (Honjo et al., 2008), the simulated transfer efficiency through the mesopelagic zone (here computed as the simulated POC flux at 960 m depth divided by the POC flux at 90 m depth) is too high for most of the sediment trap locations.
The mostly enhanced transfer efficiency in our simulation with ballasting compared to the simulation without ballasting is in line with the higher global mean sinking speed between 100 m and 1000 m compared to the applied Martin-type sinking speed profile (Fig. 2a).However, as far as we understand, even qualitatively, the global pattern of transfer efficiencies is uncertain: On the one hand, Henson et al. (2012) suggested based on satellite data that, in large areas of the tropics and subtropics except in the eastern Equatorial Pacific cold tongue, over 30 % of the export production in 100 m depth reach the deep ocean below In our simulation with particle ballasting, the Equatorial Pacific cold tongue is characterised by high transfer efficiencies, and the subtropics by lower efficiencies, which better fits to the results of Weber et al. (2016) than to those of Henson et al. (2012).
However, our results do not support the high transfer efficiencies at high latitudes suggested by Weber et al. (2016), despite the fact that opal and dust particles lead to an acceleration of detritus in the Arctic and Southern Ocean in our simulations with particle ballasting (Fig. 2d-e).The relatively high transfer efficiencies in the mid-latitudes in our simulation with ballasting are in line with high particle sinking speeds mostly due to ballasting by calcite particles (compare Fig. 2c to Fig. 3l).
Dust substantially affects sinking speeds and the transfer efficiencies in the Southern Ocean and South Atlantic off Patagonia, as well as in the Equatorial Atlantic due to Saharan dust deposition, in the North Atlantic, and in the Arctic Ocean.Note, however, that, in particular in the Arctic Ocean, where POC concentrations are low, even small quantities of dust can lead to high sinking speeds, since our simple particle ballasting parameterisation only takes into account the density of the dust (in the limit of zero POC), and not its diameter.This leads to high sinking speeds in the model, while in reality, when POC concentrations are so low that particle collisions and aggregate formation are unlikely, the small and mostly individual dust particles sink very slowly according to Stoke's law.While this illustrates the limitation of the simple particle ballasting parameterisation, the effect of the erroneously high sinking speeds in the Arctic on the carbon export is small: the integrated carbon export at 90 m depth north of 80 • N amounts to only 0.01 Gt C yr −1 , or just over 0.1 % of the global export of about 7 Gt C yr −1 .
In summary, although the simulated sinking speeds substantially differ between the simulations with and without particle ballasting, the biases of the setup with particle ballasting compared to observations remain similar to the biases of the standard setup.But the setup with particle ballasting now enables us to estimate the effect of glacial dust deposition on atmospheric pCO 2 .
4 Sensitivity to LGM dust deposition

Experimental setup
Starting from the simulation with particle ballasting described in the previous section (with the opal dissolution rate reduced by a factor of 6, modern dust deposition, modern ocean circulation and pre-industrial GHGs), 3 experiments are performed to estimate the sensitivity of the ocean carbon cycle to the reconstructed LGM dust deposition.In all 3 sensitivity experiments, the prescribed monthly-mean modern dust deposition fields are replaced instantaneously by monthly-mean LGM dust deposition fields (see Fig. 4 for timeseries of the sensitivity experiments, and Fig. 5a-c for maps of the annual mean dust deposition rates and the LGM-modern dust deposition difference).In the first experiment, the LGM dust deposition rates are only used for the computation of dust concentrations in the water column, isolating the ballasting effect of the LGM dust deposition (solid blue lines in Fig. 4).In the second experiment, the LGM dust deposition rates are only used for the computation of the bio-available iron concentrations, isolating the iron fertilisation effect associated with the dust (orange).In the third experiment, the LGM dust deposition is used for both, particle ballasting and iron fertilisation (pink).

Atmospheric pCO 2 feedback
To quantify the roles of particle ballasting and iron fertilisation by dust deposition changes for atmospheric pCO 2 , the atmosphereocean CO 2 flux anomalies relative to the reference simulation with modern dust deposition (black lines in Fig. 4) are diagnosed at the end of each simulated year, and the prescribed atmospheric pCO 2 is updated accordingly.For an accumulated net atmosphere-ocean CO 2 flux anomaly of 2.1 Gt C into the ocean, the atmospheric CO 2 concentration is reduced by 1 ppmv.
Note that only the surface ocean biogeochemistry responds to this atmospheric pCO 2 feedback.There is no effect on the atmospheric radiative transfer or climate; the prescribed climatological atmospheric forcing is unchanged.
If the atmospheric CO 2 concentration was not adjusted to the flux anomalies, a for example positive anomalous ocean CO 2 uptake would not result in a reduced atmospheric pCO 2 , the atmosphere would behave like a CO 2 reservoir in our experiments, and the subsequent ocean CO 2 uptake would be overestimated.This is illustrated in a fourth experiment, in which LGM dust is prescribed for particle ballasting without applying the simple atmospheric pCO 2 box model described above (dashed blue line in Fig. 4a-b).
The long-term net ocean CO 2 uptake in our reference simulations (Fig. 4b) would lead to a negative CO 2 trend in our simple CO 2 box model, if absolute atmosphere-ocean CO 2 fluxes were used; we avoid this effect by computing atmospheric CO 2 from the flux anomalies relative to the reference simulation in each year.Moreover, this approach, compared to, e.g., the alternative of subtracting the mean flux imbalance of the reference run, avoids additional atmospheric pCO 2 variability due to, e.g., internal variability of the ocean circulation.

Role of LGM dust as ballast
The largest ballasting effect of the LGM dust in our sensitivity experiments occurs in the Southern Ocean off Patagonia.
Here, the higher dust deposition rates lead to about 5 to 10 m day −1 faster particle sinking in 100 m depth (Fig. 5e), which is equivalent to about a doubling of the sinking speed (compared to Fig. 2b).This faster sinking locally leads to a higher export production (Fig. 5f), to a reduced surface dissolved inorganic carbon concentration (Fig. 5i), and consequently to a positive downward net atmosphere-ocean CO 2 flux anomaly (Fig. 5d), which, accumulated over the first ∼3000 years, causes an atmospheric pCO 2 drawdown by about 4.5 ppmv (Fig. 4a).This anomalous ocean CO 2 uptake occurs despite a reduction of the globally integrated export production (Fig. 4c), as opposed to the locally enhanced export production off Patagonia (Fig. 5f).The smaller globally integrated export production is due to a reduced primary production over wide areas of the Atlantic, Pacific, and Indian Ocean, as well as parts of the Southern Ocean, just north of the area off Patagonia with enhanced organic matter export (Fig. 5g and f).This reduced primary production is consistent with a stronger nitrate depletion in surface waters (not shown) due to increased particle sinking speeds (Fig. 5e).However, at the same time, the production and export of calcite shells is reduced in the sensitivity simulation (Fig. 5h).Due to higher surface silicate concentrations (not shown), this calcite export reduction is disproportionately large (compared to the POC export reduction), which leads to a reduced rain LGM -modern dust deposition difference over the ocean (kg m -2 yr -1 ) LGM dust for ballasting LGM dust for iron LGM dust deposition (kg m -2 yr -1 ) modern dust deposition (kg m -2 yr - ratio (Fig. 4d).This weakening of the calcium carbonate counter pump overcompensates the effect of the global reduction of the organic matter export on the atmosphere-ocean CO 2 fluxes, leading to a net ocean CO 2 uptake in the first ∼3000 years of our sensitivity experiment.

Role of LGM dust as a source of iron
Iron fertilisation, due to the iron contained in the dust (Section 2.2), also leads to a net CO 2 flux anomaly into the ocean.That anomalous CO 2 flux is equivalent to an atmospheric pCO 2 drawdown by about 8 ppmv within 4500 years (Fig. 4a).
However, as evident from the lack of CO 2 flux anomalies in the Southern Ocean (Fig. 5j), this CO 2 drawdown in our sensitivity experiment is not due to the expected fertilisation of phytoplankton growth by enhanced dust deposition in the Southern Ocean.Hence, the mechanism at work in our iron sensitivity experiment differs from the frequently studied classical iron hypothesis (Martin, 1990;Kumar et al., 1995;Watson et al., 2000;Werner et al., 2002;Bopp et al., 2003;Brovkin et al., 2007;Ganopolski and Brovkin, 2017).The reason for this difference is, that, in contrast to observations (see Boyd et al., 2007, for a review), non-diazotroph phytoplankton growth in our reference simulations (for the standard setup and the setup with particle ballasting) is limited by nitrate everywhere, and nowhere by iron.Since the diazotroph growth rate, unlike the non-diazotroph phytoplankton growth rate in our model, is not limited by the least available nutrient (Equation 2 of Ilyina et al., 2013), but is proportional to the product of the limiting functions for temperature, iron, phosphate and light (Equations 3 and 8 of Paulsen et al., 2017), the addition of iron can always lead to enhanced diazotroph growth.In fact, we find enhanced diazotroph growth in particular in the tropical Pacific (Fig. 5k), which locally leads to a net positive downward atmosphereocean CO 2 flux anomaly (a reduced net CO 2 outgassing in the tropical Pacific; Fig. 5j).This CO 2 flux anomaly is caused by 1) an enhanced organic matter export due to the enhanced diazotroph production (Fig. 4c), and 2) the associated reduction of the rain ratio in the tropical Pacific leading to a weakening of the calcium carbonate counter pump (Figs.4d and 5l).

Summary and discussion
We estimate the potential effect of aeolian dust deposition changes on glacial-interglacial atmospheric pCO 2 variations by computing the sensitivity of the atmosphere-ocean CO 2 fluxes to LGM dust deposition (Section 4).To this end, we use the ocean circulation and carbon cycle model MPIOM/HAMOCC combined with a new particle ballasting parameterisation (Section 2).
We find that the LGM dust deposition as estimated by Mahowald et al. (2006) leads to faster sinking and enhanced organic matter export off Patagonia, and to a reduction of the calcium carbonate counter pump globally, which results in an anomalous (glacial) ocean CO 2 uptake that is equivalent to an atmospheric pCO 2 drawdown by 4.5 ppmv within 3000 years (Section 4.3).
The iron input associated with the LGM dust deposition causes a fertilisation of diazotroph growth in the tropical Pacific, which leads to an additional atmospheric pCO 2 drawdown by more than 8 ppmv within the 4000-year-long sensitivity experiment (Section 4.4).The combination of these effects in response to the LGM dust deposition yields an atmospheric pCO 2 drawdown by about 12 ppmv (pink line in Fig. 4a).According to previous modelling studies, iron fertilisation in the Southern Ocean played a large role, leading to an atmospheric pCO 2 drawdown during the LGM by 15-40 ppmv (Watson et al., 2000;Bopp et al., 2003;Brovkin et al., 2007;Lambert et al., 2015;Ganopolski and Brovkin, 2017).However, in our simulations, due to the abundance of iron already for the applied modern dust deposition rates, non-diazotroph phytoplankton production is nowhere iron limited, and the additional iron input associated with the enhanced LGM dust deposition does not lead to a strengthening of the biological pump in the Southern Ocean.It has thus been assumed that the modern dust deposition from Mahowald et al. (2006) Mahowald et al. (2006) were replaced by an older dust deposition product (Mahowald et al., 2005) that is characterised by lower dust deposition rates over large areas of the ocean, including the South Atlantic, leading to iron limitation in the Southern Ocean in HAMOCC (Stemmler, personal communication).We expect that the addition of LGM dust anomalies using a version of MPIOM/HAMOCC that is newer than revision 4283 will lead to an additional atmospheric pCO 2 drawdown in response to the LGM dust/iron deposition in the Southern Ocean.Therefore, our results likely represent an estimate of the lower bound of the iron effect.
On the other hand, mesocosm dust deposition experiments in the Mediterranean Sea (Wagener et al., 2010), and their numerical simulation (Ye et al., 2011), as well as more recent simulations of the global iron cycle (Ye and Völker, 2017) indicated that dust can also be a sink for dissolved bio-available iron, because the dust particles provide surfaces for adsorption of dissolved bio-available iron, leading to adsorptive scavenging.By not taking adsorptive scavenging by dust particles into account, we potentially overestimate the effect of dust deposition changes on iron fertilisation.
Aggregate porosity and the degree of aggregate repackaging likely affects sinking speeds, too, as found for example in mesocosm experiments (Bach et al., 2016).In particular, opal ballast during diatom blooms did not lead to higher sinking speeds in the 25 m deep mesocosms, potentially due to the high porosities of the rather fresh aggregates.It is left for future studies to take into account the effect of porosity on the particle sinking speeds in HAMOCC.To avoid the computational expense of explicitly simulating aggregation and the porosities, the fraction of repackaged detritus (by zooplankton) to freshly produced detritus (e.g., by diatoms) could be computed and used to extend the parameterisation of particle sinking speeds.
Particle ballasting may play a larger role for glacial-interglacial pCO 2 variations than suggested by the dust-sensitivity alone.Variable ratios of opal/POC or PIC/POC can also affect particle sinking speeds via ballasting, and the ratios may have changed in response to, for example, variable temperature, variable input of weathering fluxes from land, sea-level or ocean circulation changes.Transient simulations with the coupled MPI-ESM are envisaged to quantify the role of these potential additional particle ballasting feedbacks for glacial-interglacial CO 2 variability.

Figure 1 .
Figure 1.Selection of the weathering rates during the spinup of MPIOM/HAMOCC, and adjustment of the opal remineralisation rate to modified sinking speeds.a, Strength of the Atlantic Meridional Overturning Circulation (AMOC) at 26 • N; b, surface alkalinity; c, surface dissolved inorganic carbon (DIC); d, global primary production; e, global opal production; and f, ratio of global calcite to organic carbon export at 90 m depth (rain ratio).The initial 5000-year long simulation was done with an older MPIOM/HAMOCC setup; we

Figure 2 .
Figure 2. Particle sinking speed w ballast for a, global mean profiles and b-e, at 150 m depth.The contributions of the different ballast types (calcite, opal, and dust) to the resulting, applied sinking speed (black) are estimated by re-computing the sinking speed assuming that only one type of ballast is present.The dashed grey line in panel a shows the Martin-type sinking speed applied to POC in every watercolumn in the simulation without ballasting.Note the logarithmic depth scale, the euphotic zone is indicated by grey shading.

Figure 3 .
Figure 3.Comparison of simulations with particle ballasting and without particle ballasting (standard setup with constant / Martin-type particle sinking) to observations.a-c, net atmosphere-to-ocean CO2 flux from observations (Takahashi et al., 2009) compared to our simulations; d-f, surface dissolved phosphate from the World Ocean Atlas (Garcia et al., 2013) compared to our MPIOM/HAMOCC simulations; g-i, carbon export at 100 m depth as derived from satellite data (Henson et al., 2012, filled circles) and from sediment trap data (Honjo et al., 2008, crosses) compared to the export at 90 m depth in our simulations; and j-l, transfer efficiency (carbon flux at 1000 m depth divided bythe carbon flux at 100 m depth) as derived from sediment trap data(Honjo et al., 2008) compared to the 90 m to 960 m transfer efficiency in our simulations.For all simulations in this figure, modern dust deposition rates were prescribed; shown are averages computed over the model years 12900 to 12999.

Figure 4 .
Figure 4. Effects of switching to LGM dust deposition for particle ballasting (blue), for the computation of iron fertilisation associated with the dust (orange), and for both (pink) compared to the control simulation with ballasting and iron fertilisation based on modern dust deposition (black).a, estimated atmospheric pCO2 change based on the annual mean atmosphere-ocean CO2 flux anomalies relative to the control run (black); b, global mean atmosphere-ocean CO2 fluxes; positive values indicate a net flux from the atmosphere into the ocean;solid lines are 50-year running means; the dots show annual means of the CO2 fluxes in the simulation with the standard setup without particle ballasting to illustrate the large interannual variability (e.g., due to ocean circulation variability) compared to the dust-induced anomalies; the blue dashed line is the CO2 flux in a simulation where the atmospheric pCO2 was not updated every year according to the atmosphere-ocean CO2 flux anomalies of the previous year; c, global organic matter export in 90 m depth; d, global ratio of particulate inorganic carbon (PIC, calcite) to particulate organic carbon (POC) within the particle export in 90 m depth.

Figure 5 .
Figure 5. Effects of switching from modern to LGM dust deposition.a, Modern dust deposition, b, LGM dust deposition, and c, LGMmodern dust deposition difference according toMahowald et al. (2006); d-i, LGM instead of modern dust deposition rates applied to particle ballasting leading to changes of d, the atmosphere-ocean CO2 flux (positive into the ocean), e, the common sinking speed w ballast in 100 m depth, f, the organic matter (OM) export in 90 m depth, g, the vertically integrated primary production, h, the calcite export in 90 m depth, and i, the surface dissolved inorganic carbon (DIC) concentration; and j-l, LGM instead of modern dust deposition rates applied to the computation of dust-related iron input leading to changes of j, the atmosphere-ocean CO2 flux, k, the concentration of diazotrophs, and l, the rain ratio (PIC/POC ratio of the export flux at 90 m depth).The changes are computed from the first 100 years of the sensitivity simulations relative to the control simulation with modern dust deposition (blue or orange minus black in Fig.4).
Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-137Manuscript under review for journal Geosci.Model Dev. Discussion started: 16 August 2018 c Author(s) 2018.CC BY 4.0 License.are provided as supplementary material.To build the model code used in this study, one would need to obtain MPI-ESM via the MPI-ESM-Forum and manually include the code modifications at the described locations.The data to reproduce the figures in this manuscript is available upon request.Stemmler and Katharina Six for sharing their HAMOCC expertise.The model simulations were performed at the German Climate Computing Centre (DKRZ).This work was supported by the German Federal Ministry of Education and Research (BMBF) as Research for Sustainable Development (FONA, www.fona.de)through the PalMod project (FKZ: 01LP1505D).16 Geosci.Model Dev.Discuss., https://doi.org/10.5194/gmd-2018-137Manuscript under review for journal Geosci.Model Dev. Discussion started: 16 August 2018 c Author(s) 2018.CC BY 4.0 License.
is an overestimate.Hence, since MPIOM revision 4283, which is included in the recent CMIP6-release of MPI-ESM (MPI-ESM version 1.2.01,MPIOM version 1.6.3,revision 4639), the dust deposition fields by