Open Research Online A multi-variable box model approach to the soft tissue carbon pump

. The canonical question of which physical, chemical or biological mechanisms were responsible for oceanic uptake of atmospheric CO 2 during the last glacial is yet unan-swered. Insight from paleo-proxies has led to a multitude of hypotheses but none so far have been convincingly supported in three dimensional numerical modelling experiments. The processes that inﬂuence the CO 2 uptake and export production are inter-related and too complex to solve conceptually while complex numerical models are time consuming and expensive to run which severely limits the combinations of mechanisms that can be explored. Instead, an intermediate inverse box model approach of the soft tissue pump is used here in which the whole parameter space is explored. The glacial circulation and biological production states are derived from these using proxies of glacial export production and the need to draw down CO 2 into the ocean. We ﬁnd that circulation patterns which explain glacial observations include reduced Antarctic Bottom Water formation and high latitude upwelling and mixing of deep water and to a lesser extent reduced equatorial upwelling. The proposed mechanism of CO 2 uptake by an increase of eddies in the Southern Ocean, leading to a reduced residual circulation, is not supported. Regarding biological mechanisms, an increase in the nutrient in either the equatorial regions or the northern polar can atmospheric CO and satisfy of glacial Consistent with previous studies, CO is drawn more easily through increased productivity in the Antarctic region than the sub-Antarctic, but that observations of lower export production there. The glacial states are more sensitive to changes in the circulation and less sensitive to changes in nutrient utilization rates than the interglacial states.


Introduction
During the last 800 000 years atmospheric pCO 2 has varied in concert with Antarctic surface air temperature (EPICA Community Members, 2004;Lüthi et al., 2008;Petit et al., 1999;Siegenthaler et al., 2005). Determining the mechanisms behind the correlation remains a tantalizing question in paleoclimatology today. Kohfeld and Ridgwel (2009) provide an excellent review of the processes that are known or proposed to play an important role in lowering glacial atmospheric pCO 2 . They conclude from previous work that lower sea surface temperatures can explain approximately 26 ppmv (21-30 ppmv) of the interglacial to glacial pCO 2 decrease and that this decrease was countered by a reduction in the terrestrial biosphere and reduced ice volume which resulted in an increase of about 22 ppmv (12-36 ppmv) and 13 ppmv (11-17 ppmv) respectively. Other lesser understood mechanisms include changes in the ocean circulation, aeolian Fe fertilization, bacterial metabolic rate change, Si fertilization, and coral reef growth.
The observed pCO 2 -temperature correlation in the glacial record is much stronger in the Antarctic (AA) than in the Northern hemisphere Greenland records. This is partly because of the more complex and variable atmospheric circulation patterns that prevail over Greenland. However, the very good correspondence of temperature and pCO 2 in Antarctica suggest that the Southern Ocean (SO) played an important role in the glacial carbon cycle. Oceanic uptake of atmospheric CO 2 through air sea gas exchange is facilitated by the soft tissue pump, solubility pump and carbonate pump (Toggweiler et al., 2003a,b). Most of the lesser understood hypotheses of glacial oceanic carbon uptake relates foremost to the efficiency of the soft tissue biological pump where the efficiency is defined as the percentage of all nutrients exported from the surface layer that was exported as organic matter and is therefore affected by both the biology and the circulation of the ocean.
Published by Copernicus Publications on behalf of the European Geosciences Union. The zonally integrated transport in the SO can be divided into a deep poleward cell of which the downwelling branch represents Antarctic Bottom Water (AABW) formation and an upper, more equatorward cell referred to as the SO residual circulation (Karsten and Marshall, 2002;Marshall and Radko, 2003). Increased net uptake of atmospheric CO 2 (through reduced outgassing of CO 2 ) has been proposed to be a result of a weakening of the poleward cell though enhanced surface stratification or reduced winds (Sigman and Boyle, 2001;Sigman et al., 2004;Toggweiler et al., 2006;Francois et al., 1997) or a weakening of the residual circulation through reduced buoyancy or winds (Watson and Garabato, 2006;Fischer et al., 2010;Parekh et al., 2006). A reduction in either of these cells would imply less deep upwelling during the glacial for which evidence remains ambiguous (Anderson et al., 2009;De Pol-Holz et al., 2010). As noted, given a steady state circulation, the atmospheric CO 2 can also be reduced through increased biological production in the SO (Sarmiento and Toggweiler, 1984). The main hypothesis for lowering glacial pCO 2 as a consequence of increased SO biological export is that there was an enhanced supply of the limiting micro nutrient iron (Martin, 1990;Joos et al., 1991;Watson et al., 2000;Brovkin et al., 2007).
The role of the equatorial upwelling regions in the carbon cycle is still under debate. Matsumoto and Sarmiento (2002) suggested that, during the last glacial, smaller Si:N ratios in Antarctic plankton led to more Si leakage to equatorial regions where higher Si enhanced diatom production (at the expense of cocolithophores) weakened the carbonate pump. Sediment core data does not reveal the high opal fluxes needed to support the Si leakage hypothesis (Bradtmiller et al., 2006) but Pichevin et al. (2009) recently argued that the lower equatorial opal burial rate is due to a decrease in the silicon to carbon uptake ratio found in iron rich environments. The hypothesis is supported by clear evidence for an enhanced glacial iron supply to the Eastern Equatorial Pacific (McGee et al., 2007;Winckler et al., 2008). It should be noted that more organic matter export from the Equatorial Pacific does not necessarily imply a stronger soft tissue pump. Enhanced equatorial nutrient utilization could also imply a geographical shift in global export production (i.e., from the sub-tropics to the tropics) rather than an overall increase.
The traditional process through which solutions to the glacial carbon problem are sought can, conceptually at least, be seen as comprising of four steps. The first step is to measure or assimilate existing knowledge of the glacial state (i.e., proxies of the atmospheric and oceanic biochemistry and circulation). The second step is to form a hypothesis of a specific mechanism that may have led to a glacial atmospheric drawdown of CO 2 . As a next step, the hypothesis is often tested in a simple box model in which the general circulation is set as close as possible to observations and the sensitivity of the atmospheric pCO 2 to the specific mechanism is derived. If the result is promising, the final step is to verify the hypothesis in a more complex ocean or coupled model. Ideally not only the CO 2 uptake, but also the circulation and biology of the glacial simulation should be consistent with the observations in step one. All four of these steps are needed and each step has obvious advantages and drawbacks. The main problem with the box model approach is that it requires the simplification of the circulation and biology in ways that are rather ad hoc and the implications of these simplifications for the conclusions are not clear. The physical and biological processes in general circulation models are closer to that of the real world. However, they are expensive to run and, perhaps more importantly, it is often impossible to distinguish in a complex model between various mechanisms for CO 2 uptake. For example, increasing the SO winds stress has an effect on the AABW formation, sea ice, the residual circulation and the North Atlantic deep water formation and it is therefore not clear how the resulting atmospheric CO 2 concentration has come about.
We propose here a new type of box model that eliminates some of the common problems with box models associated with hidden sensitivities to input parameters. It is inexpensive to run so that all combinations of circulation and biological CO 2 uptake mechanisms can be explored in conjunction with each other or in isolation. It is akin to an inverse model in that we deduce the circulation and biological parameters for the interglacial and glacial from observations rather than assume any values a priori. The primary observations used to determine the interglacial states are observed phosphate distributions. For the glacial we isolate states that draw down the most CO 2 and whose export production in relation to the interglacial is consistent with the observational synthesis of Kohfeld et al. (2005). The bulk of the evidence suggests decreased export in the Antarctic region of the SO, and increased export in the Sub-Antarctic and the Equatorial upwelling regions. Because all our input parameters assume the whole range of realistic values, there are no hidden sensitivities other than those built into the geometry of the model. We explore two other model setups to ensure that our conclusions are not highly dependent on the geometry.
The purpose of this study is to identify which types of glacial circulation and productivity regimes are consistent with the observations and to relate these to current hypotheses for glacial oceanic uptake of atmospheric CO 2 . A detailed description of the model is provided in the next section as well as the link between preformed nutrient concentration and atmospheric pCO 2 . In Sect. 3 we present the most likely solutions of the glacial and interglacial, and explore the sensitivity of the pCO 2 and export production to the circulation and biology. The implication of the results for glacial carbon cycle hypotheses are discussed in Sect. 4 and overall conclusions are drawn in Sect. 5.

Model description
Our standard model has seven boxes of which six are in the upper 500 m (Fig. 1a). As is typical in box models (Watson and Garabato, 2006;Johnson et al., 2007), we assume a closed system that represents a zonally integrated ocean.
Combining the Atlantic and Indo-Pacific into one basin is not ideal but the assumptions we shall make about latitudinal distributions of nutrients and export and where water upwells and sinks, hold in both basins. The SO is divided into an Antarctic box and a sub-Antarctic box and the division is roughly at the polar front, assumed to be at 55 • S. The low latitudes extend from 40 • S to 40 • N and are intersected by a narrow box, from 10 • S to 10 • N, which explicitly represents the high nutrient equatorial upwelling region. A northern box starts at 40 • N and extends to the north pole. The volumes of the surface boxes are the same as in the real world between these specified zonal boundaries and for the deep box between 500 m and the bottom. We shall use the term SO to describe both the Antarctic and the Sub-Antarctic boxes. In order to relate the model to the real world, we refer to the sinking in the Antarctic box as AABW.
The fluxes between the boxes are fully described by five independent model parameters: (1) Antarctic deep water formation, T A , (2) northward surface Ekman transport between the two SO surface boxes, T T , (3) a southward transport representing the contribution of SO mesoscale eddies, T ME , (4) an equatorial upwelling flux, T EU , and (5) a northern upwelling flux, T NU . These variables are indicated in Fig. 1a and the range of values they take is listed in Table 1. The rest of the transport terms are quantities derived from continuity. Note that the southward eddy flux is never more than the northward Ekman transport so that the residual flow into the equatorial box is always positive. Upwelling in the SO occurs in the Antarctic box (south of the polar front) and is equal to the net flow out of the box from Antarctic sinking and the residual circulation. In the northern box the upwelling is balanced by local downwelling so that it can be seen as a mixing term that brings up high nutrient water to the surface similar to winter mixing in the real ocean. The upwelling in the equatorial box represents both wind driven upwelling and that due to turbulent mixing. There is no simple way to represent the complex dynamics of the equatorial region -we have chosen to return the upwelled water to the deep ocean in the northern box but in reality some downward mixing of the surface properties occurs locally and some of the upwelled water is subducted in the low latitudes. To ensure that our results are not highly dependent on equatorial upwelling or the northern sinking that goes along with it, the study is repeated in a model that does not have an equatorial box (Fig. 1b). In this case the equatorial box is treated the same way as the low latitude boxes. In another test of the effect of the geometry the upper layer depth was reduced from 500 m to 300 m.
The soft tissue biological pump in our model is simplified to a one-nutrient formulation, P , denoting phosphate (PO 4 ). The nutrients are advected by the circulation in the usual way. Export production from a surface box is linearly proportional to the nutrient concentration in the box. The proportionality constant, α, differs from box to box and encapsulates factors which limit the uptake of nutrients such as light or iron. It is referred to in the rest of the text as the nutrient utilization rate and is the fraction of the available nutrients that is exported from the surface as organic matter (rather than as preformed nutrients) each second. This constant is allowed to vary over a wide range of values (Table 1). The export production from the box is then simply, αP , multiplied by the volume of the box. The overall nutrient content is set to be the same as that derived from the World Ocean Atlas of 2005 (Garcia et al., 2006). The above model is mathematically formulated in the following seven equations, where T is a transport variable and V the box volume. The subscripts A, S, LS, E, LN, N, and D denote the Antarctic box, Sub-Antarctic box, Low latitude Southern box, Equatorial box, Low Latitude Northern box, Northern box, and Deep box respectively. Upwelling fluxes into the boxes have a U at the end, i.e., T AU is upwelling into the Antarctic box. The SO Ekman transport and mesoscale eddy transport are denoted by T T and T ME respectively and P SUM is the total amount of PO 4 molecules in the ocean. These equations are solved for 10 7 random combinations of the five circulation and five nutrient utilization input parameters (see Table 1 for range). Note that the six surface boxes give only five nutrient utilization parameters because we assume that the conditions for biological production in the northern and southern low latitude regions are similar and therefore set α to be the same in these boxes (i.e., α LS = α LN ).

The relation of preformed nutrients to atmospheric CO 2
During glacial periods CO 2 is drawn into the ocean through an increase in the surface to deep gradient of dissolved inorganic carbon (DIC). This gradient relies on the soft tissue biological pump in which carbon is being sequestrated into organic matter in the surface and later remineralized in the deep ocean. Nutrients that are advected to the deep ocean as a result of the physical circulation of the ocean are not available for biological uptake and are referred to as preformed nutrients. If the proportion of preformed nutrients in the deep ocean is large, the biological pump is inefficient, the vertical gradient of DIC weak, and the atmospheric CO 2 high. The strong positive relation between preformed nutrients and atmospheric CO 2 has been derived theoretically (Ito and Follows, 2005;Goodwin et al., 2008) and confirmed in ocean general circulation models (Ito and Follows, 2005;Marinov et al., 2006). In our model we describe the biology in terms of the macro nutrient PO 4 and then we use this theory (Ito and Follows, 2005) to relate the preformed nutrient concentration to atmospheric CO 2 . For completeness and to alert the reader of the assumptions made, we briefly summarise the theory here. For a more comprehensive derivation please refer to the original reference. Dissolved inorganic carbon in the ocean can be subdivided into preformed and regenerated carbon. The preformed concentration can in turn be written as the sum of a saturated component C sat and a preindustrial air-sea disequilibrium component C and the regenerated carbon concentration can be split into a soft tissue, C org , and calcite component, C calcite . Assuming a constant reservoir of carbon in the ocean and atmosphere, a small change in the pCO 2 of the atmosphere can be written in terms of the changes in the oceanic carbon components, i.e., Mδ pCO atm 2 + V δC sat + δ C + δC org + δC calcite = 0 (1) where M is the total moles of gas in the atmosphere and V is the volume of the ocean. For simplicity it is assumed that changes in the soft tissue pump are independent of changes in the carbonate pump or saturation so that δ C = δC calcite = 0). Using the Buffer factor, B (Bolin and Erikson, 1959), one can approximate variations in δ C sat by Eliminating C sat from Eqs. (1) and (2) gives an expression for the sensitivity of pCO atm 2 to C org δpCO atm We now define the efficiency of the biological pump, where P reg is the mean regenerated phosphate concentration and P 0 is the global mean phosphate concentration. For a 100% efficient pump P * will be 1. A perturbation in the organic carbon concentration can be written in terms of an equivalent change in the biological pump using a fixed Redfield ratio, δC org = R C:P δP reg = R C:P P 0 δP * Combining Eqs. (3), (4) and (6), a linear relation between pCO atm 2 and the efficiency of the biological pump as defined by Eq. (5) is found, In our model, the preformed nutrient concentration is defined as the fraction of nutrients that reach the deep ocean through the circulation multiplied by the deep ocean nutrient concentration. The biological pump efficiency is the fraction of nutrients that is exported to the deep ocean through biological productivity and ranges from near 0 to 100% for our 10 7 solutions (Fig. 2). The preformed nutrient concentration ranges from 0 -2.5 µmol Kg −1 and can be related by Eq. (7) to a pCO 2 of approximately 300. Given that the theory that relates pCO atm 2 to a change in preformed nutrients concentration is not designed for large perturbations the results that pertain to large glacial-interglacial changes should be interpreted with caution.  Ito and Follows (2005). A zero pCO 2 is assumed at the average preformed nutrient concentration of our 300 best interglacial solutions (see Sect. 3.2 for an explanation of how these are derived). The pCO 2 and preformed nutrients are also linearly related to the efficiency of the biological pump (right axis).

The interglacial states
It is common in carbon cycle box models to set the circulation as close as possible to observations or model output or to tune it to give realistic solutions (Keir, 1988;Toggweiler, 1999;Watson and Garabato, 2006). While these approaches have merit, it is not clear to which extent the conclusions of the studies depend on the chosen circulation states. In this study the whole parameter range is explored to ensure that our conclusions are not dependant on a chosen parameter. Initially 10 7 random solutions are calculated to cover all parameter combinations within the limitations of our model setup. From these we select the modern states according to the circulation and nutrient and organic matter export distributions. In the first step, solutions are required to have higher export in the Sub-Antarctic than in the Antarctic box, and higher export in the equatorial box, northern and sub-Antarctic boxes than in the low latitude boxes. The circulation of the modern state is required to have at least 20 Sv (1 Sv = 10 6 m 3 s −1 ) of northward Ekman transport and AABW formation between 5 and 25 Sv (Table 2). Of the solutions that fulfil these criteria we consider the 300 solutions that best fit the modern phosphate distribution as derived from the World Ocean Atlas 2005 (Garcia et al., 2006). The maximum deviation from PO 4 observations for any box of any of the 300 solutions is 0.3 µmol kg −1 ). Figure 3 (top right) shows that these constraints can be satisfied using any among a wide range of input parameter   (bottom) parameters of the interglacial states (right) and the glacial states (left) are shown here against the preformed nutrient concentration in the deep ocean. The top axis relates the difference in preformed nutrient between the glacial and interglacial states to changes in atmospheric pCO 2 (as discussed in Sect. 2.2). Best fit lines are drawn through the solutions for visual aid. They do not imply an independent linear correlation between the preformed nutrient concentrations and the input parameters.
sets. In this sense we effectively allow for model structural error using the plausibility approach of Holden et al. (2009) or Edwards et al. (2010). The Ekman transport and AABW cover the whole range of permitted values in the 300 solutions. To satisfy the observations of high surface nutrients, the northern and equatorial upwelling fluxes are towards the higher range of their permitted values. The eddy mass flux is somewhat on the low side so that the average residual circulation of the interglacial parameter sets is about 16 Sv, with a minimum of 4 Sv. The interglacial nutrient utilization rates have more confined patterns than the circulation parameters. Low utilization rates in the AA box, the low latitudes, and the equatorial box are necessary to simulate modern phosphate distributions (Fig. 3, bottom right). The wide range of parameters that fulfil the basic interglacial criteria suggest that the traditional approach of picking one best guess state on which to do a sensitivity analysis may be unwise. The conclusions drawn in this study are instead based on all these states.
The purpose of isolating the best modern states is not to determine the real modern circulation or nutrient utilization. Rather, the modern states are obtained to (1) confirm that the model can reproduce a realistic interglacial state, (2) indicate the sensitivity of pCO 2 and nutrient concentrations to the input parameters, given a modern state, and (3) provide reference states with which to compare the glacial states. The Table 2. Constraints used to derive the 300 interglacial and glacial. E Z denotes the export flux from box Z and is equal to P Z· α z· V Z . The operator () i indicates an average of the 300 interglacial states.

Period
Essential constraint Final constraint Interglacial T T > 20 Sv The 300 states, fulfilling the 5 Sv < T A < 25 Sv essential constraints, whose Table 1. Glacial The 300 states, fulfilling the average efficiency of the biological pump of the 300 interglacial solutions is 35 ± 4% which corresponds well with the estimate of 36% from Ito and Follows (2005).

The glacial states
Observations of the glacial ocean circulation are ambiguous but it is clear that the deep oceans were filled with a southern source water to a more shallow depth than today, and pore fluid measurements suggest this southern source water was much saltier than that from the North Atlantic (Ninnemann and Charles, 2002;Adkins et al., 2002b). Proxies for export production indicate increased export in the Sub-Antarctic and Equatorial regions and reduced export in the Antarctic region (Kohfeld et al., 2005). It is uncertain by how much the export has changed in each region -to select our most likely glacial states, we limited the solutions to those that have at least a 20% increase in export in the Sub-Antarctic and Equatorial boxes and a 20% decrease in the Antarctic box ( Table 2). The resulting states are only weakly dependent on this number (i.e., the percentage rate). From these we select the 300 that have the lowest preformed nutrients, i.e., the states that are associated with the greatest draw down of CO 2 (Fig. 3, left). The average preformed nutrient concentration, biological pump efficiency, and pCO 2 change for these states are 1.6 µmol kg −1 , 84%, and −165 ppmv respectively. We choose the 300 states with the most severe drawdown of CO 2 as it illuminates clearly the parameters that are important for CO 2 uptake through the soft tissue pump, but our conclusions are qualitative and hold also for smaller changes in pCO 2 , say or the order of −100 ppmv (confirmed, but not shown).
The circulation parameter which undergoes the greatest change in the glacial is the Antarctic bottom water (13 Sv average decrease), followed closely by upwelling (i.e., mixing) in the Northern high latitude box (12 Sv average decrease). The northward Ekman transport is at least 13 Sv and Table 3. The mean and standard deviation of the change in pCO 2 due to a change in the input parameters for the 300 glacial and interglacial states, and for 300 randomly chosen states. For the circulation parameter the difference in pCO 2 due to a change in the circulation from 5 to 20 Sv is shown and for the nutrient utilization the difference in pCO 2 due to a change in α from 2 to 8e-10 s −1 .  (T T ) 5 ± 12 0 ± 5 −7 ± 16 Eddies (T ME ) −10 ± 8 11 ± 4 0 ± 12 Eq-upwell (T E ) 30 ± 8 16 ± 6 4 ± 20 North-upwell (T N ) 67 ± 7 25 ± 4 52 ± 18 Model Boxes Glacial Interglacial Random (ppm/6 × 10 10 s −1 ) (ppm/6 × 10 10 s −1 ) (ppm/6 × 10 10 s −1 ) on average 24 Sv. Interestingly, the southward eddy flux is never more than 20 Sv and always at least 9 Sv less than the Ekman transport (see Sect. 4.3 for discussion). Equatorial upwelling can take on any value in the glacial. The nutrient utilization rates in the low latitudes and northern boxes are not tightly constrained in the glacial states. However, the Antarctic nutrient utilization rate is always low and the Equatorial and Sub-Antarctic utilization rates are high. The implications of the results are discussed in Sect. 4. Note that the average change can be viewed as the most likely change and the range gives the possible solutions that are consistent with the observations. The glacial and interglacial states have so far been presented in terms of the input parameters and the preformed nutrients. The actual variables that are solved for are the phosphate concentrations. The average nutrient distribution for the 300 glacial states has a similar pattern to the interglacial although surface nutrients are reduced everywhere north of 40 • S (Fig. 4, top). The export production is higher by design in the Equatorial and Sub-Antarctic boxes and lower in the Antarctic box of our modelled glacial states (Fig. 4, bottom). In the low latitude boxes we also find increased export which is consistent with the observations but this may well be fortuitous because the observations are in upwelling regions which we do not resolve in the model. The glacial states have a lower export production in the Northern box which is again consistent with the few data points north of 50 • N in the Kohfeld et al. (2005) reconstruction. Fig. 4. The average surface nutrient concentrations (top) for the 300 interglacial states (solid red lines and diamonds) and the 300 glacial states (solid blue lines and stars). Also shown for reference is the observed modern nutrient distribution (dashed red lines and crosses). Export production (bottom) from each surface box for the interglacial and glacial states. The modern nutrient distribution were obtained from the World Ocean Atlas (Garcia et al., 2006).

Sensitivity of preformed nutrients to circulation and nutrient utilization
In the previous section we determined which solutions are consistent with the proxies of glacial export production and also give significantly decreased preformed nutrient concentrations in the deep ocean. This approach identifies what the possible combinations of circulation and nutrient utilization www.clim-past.net/6/827/2010/ Clim. Past, 6, 827-841, 2010 The response of atmospheric pCO 2 (solid black), Antarctic export production (from box A, solid dark blue), Sub-Antarctic export production (from box S, dashed magneta), and Equatorial export production (from box E, dashed light blue) to the five circulation input parameters, T A , T T , T ME , T EU , and T NU . The export fluxes (left axis) and pCO 2 (right axis) are expressed as the difference from that of the best fit interglacial change.
rates are that we need to look for in mechanisms to explain the glacial carbon cycle. It does not illustrate which of these parameters are responsible for the CO 2 drawdown and which are only necessary to explain the proxy data of export production. For instance, it is likely that the reduced nutrient utilization rate in the glacial in the Antarctic box (Fig. 3) is necessary to reduce export there (as forced by the observational constraint) and not to reduce preformed nutrients.
Here we look at the sensitivity of the pCO 2 and the export production to the input parameters of the best fit interglacial state (Figs. 5 and 6). The export flux and pCO 2 are expressed as the difference from the export and pCO 2 of the interglacial state. The sensitivities to input parameters are dependent on the chosen reference state so we have determined the sensitivity of pCO 2 to input parameters for all 300 interglacial states, for the 300 glacial states, and for 300 randomly chosen states from the complete set of solutions (Table 3). In general the glacial states are more sensitive to changes in the circulation parameters (especially those relating to surface-to-deep fluxes) and less sensitive to changes in nutrient utilization rates than the interglacial states. This is likely because the circulation is weaker in the glacial than the interglacial states so that a change in one circulation parameter has a bigger overall effect on the ventilation rate of the deep ocean. Similarly, the nutrient utilization rates are very high in the glacial states so that a change in nutrient utilization in one box has a smaller effect on preformed nutrients than the equivalent change in the interglacial states where average utilization rates are much lower. From Sect. 3.2 it is clear that the AABW formation, the northern upwelling, and the equatorial upwelling are all less in the glacial states than in the interglacial states. Figure 5 shows that pCO 2 is highly sensitive to the former two but indicates no great sensitivity of pCO 2 to equatorial upwelling. However, if one considers all of the 300 interglacial and glacial states, the picture that emerges is different. The mean change in pCO 2 due to a 15 Sv change in the equatorial upwelling is 16 ppmv in the interglacial and 30 ppmv in the glacial, almost half that of the AABW and northern upwelling ( Table 3). The example indicates how misleading it can be to consider the importance of a carbon uptake mechanism only in a best guess fixed state (and this may very well be true for a control state in a numerical model too). The pCO 2 has a weak negative relation to the eddy flux and is largely independent of the Ekman transport (Fig. 5, Table 3). The atmospheric pCO 2 is reduced when uptake of nutrients is more efficient in any of the surface boxes in our model, but the sensitivity is highest north of 40 • S (Fig. 6, Table 3). The most severe glacial states (i.e., those glacial states with lowest pCO 2 ) have high nutrient utilization rates in all boxes except the Antarctic box (Fig. 3). Here an increased utilization rate implies a higher export production (Fig. 6, top left) and that is not consistent with the observations. The surface transport is always from the south to the north so that biological productivity in the northern box does not influence the nutrient distributions or export in the other surface boxes. While it is clear that increased nutrient utilization in the low latitudes and northern boxes helps to draw down CO 2 (Table 3), the utilization parameters for these regions are not tightly grouped in the glacial states (Fig. 3). The lack of observations of export production here means that not enough information is available to indicate whether an increase in the nutrient utilization rate occurred.

The representation of the equatorial region
We have chosen to describe the equatorial region separately in our control model because (a) this region behaves differently from the sub-tropics, (b) there are a high number of export production observations here, and (c) various hypotheses for glacial CO 2 uptake revolve around this region.
Including an equatorial upwelling flux is necessary to produce the high nutrient concentrations that are observed here. We have assumed, perhaps more for simplicity that any other reason, that the upwelled water all sinks to the deep ocean in the northern box. In the reality the dynamics in the equatorial region are complex and some of this upwelled water sinks out of the surface ocean in the sub-tropical gyres. Upwelling in the Equatorial regions is also sometimes modelled as a mixing (two-way) flux but in reality is it unlikely that there is a downwards transport of mass. The location of the ultimate sinking of the upwelling water does not affect the nutrients in the equatorial region, but it does affect the dynamics in the boxes to its north. To ensure that the conclusions drawn from our results are not dependent on sensitive equatorial dynamics, we have repeated the analysis in a model in which the equatorial box has been removed (Fig. 1b). The low latitude boxes are combined into one box. The concentration to which the low latitude phosphate is restored is the observed average phosphate concentration between 40 • S and 40 • N, excluding the equatorial regions between 10 • S and 10 • N. Naturally the glacial constraint of higher equatorial export production is also not applied.
The circulation and biological parameters that produce the best glacial and interglacial states (Fig. 7) are remarkably similar to those of the standard model that includes an equatorial box (Fig. 3). The only clearly noticeable difference  is that the interglacial solutions are more tightly grouped and that the nutrient utilization rates in the sub-Antarctic and northern boxes are lower when no equatorial box exists. These are the boxes where the highest biological export production occurs. To compensate for the lack of upwelling of deep nutrients to the surface (without an equatorial box), the biological export from these two boxes are reduced.

Discussion
The main input parameters that are directly responsible for reduced atmospheric pCO 2 in the glacial states of the model are AABW formation, northern upwelling, and the nutrient utilization rates in the equatorial and low latitude boxes (Figs. 3, 5 and 6, Table 3). Perhaps surprisingly, the SO eddy fluxes and SO nutrient utilization rates are of secondary or no importance. We discuss the implication of each of our main glacial input parameter changes as well as those that are not as important as expected.

Reduction in AABW and northern upwelling
Our model results suggest a glacial decrease in both AABW formation and mixing-driven upwelling of deep water in northern high latitudes. The importance of SO deep ventilation for atmospheric CO 2 has also been found by previous box model studies where the ventilation is usually viewed as a surface-to-deep mixing flux rather than a deep water formation transport (Köhler et al., 2005). Numerical model simulations of the ocean circulation during the Last Glacial Maximum give ambiguous results (Otto-Bliesner et al., 2007) but in general support the observation that the boundary between North Atlantic and Southern Source deep water was shallower than today. Observations of δ 13 C, anoxic conditions in the deep SO, and reduced δ 15 N ratios in the glacial SO suggest that the ventilation was weaker there and the surface ocean more stratified (Francois et al., 1997;Hodell et al., 2003;Toggweiler et al., 2006). Sigman et al. (2004) proposed that the increased glacial stratification in the SO and North Pacific is a result of a drop in the mean ocean temperature. At very low temperatures the seawater density becomes almost insensitive to temperature changes and is mostly affected by salinity. In cold climates, heat loss becomes unable to destabilize polar haloclines, leading to a reduction in deep water formation. De Boer et al. (2007) tested the effect of cold water induced stratification in an ocean general circulation model and found that it was especially efficient in reducing convection in strong halocline regions such as the North Pacific and the SO. Such an effect could thus have been responsible for reducing AABW and northern mixing-driven upwelling of deep water in the glacial ocean, as suggested by our model.
Another proposed mechanism for reducing AABW is reduced and northward shifted winds . Note that the winds would not affect preformed nutrients through a reduced northward Ekman transport and associated upwelling but rather through the AA stratification that results from less wind-enhanced winter mixing and a consequent reduction in AABW formation. However, numerical modelling study has so far found no or a very small atmospheric CO 2 reduction from reduced or equatorward shifted winds (Tschumi et al., 2008;Menviel et al., 2008). Moreover, Sime et al. (2010) made an extensive observations and numerical model comparison which suggested that the SO winds were neither stronger nor shifted northward during the last glacial.
An alternative explanation for reduced AABW formation during the last glacial that may have been active in conjunction with mean ocean temperature stratification is that destabilization of the water column was rare because of the extremely salty deep water (Adkins et al., 2002a). It is possible that the salty deep water formed through brine rejection in the Southern Ocean at a small rate. The high vertical stratification could have reduced the vertical mixing that might otherwise have caused the deep ocean water to decrease in density and rise to the surface (Watson and Garabato, 2006). The combination of very cold temperatures, salty deep water formed by brine rejection, and reduced deep mixing can explain a deep ocean filled with southern sourced water that nevertheless has very low rates of ventilation.

Increase in low latitude and equatorial export production
As noted in our model description, box models are not ideal to simulate the complex dynamics of the equatorial regions. However, a consistent result on all our model setups is that the biological pump is sensitive to the nutrient utilization rate in our low latitude and equatorial boxes and it is therefore appropriate to put the result in context of previous work. One proposed explanation for glacial oceanic uptake of CO 2 is an increased dust supply to the Eastern Equatorial Pacific that invigorated the soft tissue pump (McGee et al., 2007;Winckler et al., 2008). This mechanism has been criticized due to observations of reduced opal accumulation rates in the region (Bradtmiller et al., 2006) but Pichevin et al. (2009 recently argued that the lower opal export can be explained by a decrease in the silicon to carbon uptake ratio in an iron rich environment. Our model results supports the concept of a global soft tissue pump being enhanced by a higher utilization of equatorial nutrients. However, it is not obvious that enhanced equatorial export would significantly decrease atmospheric CO 2 in the real ocean. The key to a global reduction in atmospheric CO 2 through the soft tissue pump lies in the concentration of preformed nutrients in the deep ocean. Uptake of nutrients at the equator will only reduce deep ocean preformed nutrients if those nutrients would otherwise escape to a deep water formation region and be advected to the deep ocean in the physical circulation. In the real ocean this is unlikely to happen as "escaped" nutrients have to pass through the low latitudes where the utilization is so high that surface nutrients are almost completely depleted. In our model the low latitude nutrients are not depleted because the nutrient concentration is averaged over a 500 m deep upper layer. Hence nutrients leaving the equatorial region in our model interglacial states can reach the northern region unscathed and sink through advection there. This is indeed what happens as we find the nutrient concentration in boxes north of the equator much reduced for the glacial states (Fig. 4). It is likely that the importance of the low latitude nutrient utilization in our glacial states is similarly exaggerated due to the fact that our surface ocean is 500 m deep. In the real ocean, the nutrients in the interglacial surface ocean are almost completely depleted. A large increase of the drawdown of low latitude nutrients, as suggested by our model for the glacial, is therefore unlikely. In spite of the fact that biological production is limited to the near surface of the ocean, say the top 100 m, we have chosen a deeper level of 500 m because the transport of Ekman and eddy flow and the resulting northward residual flow does not occur in the euphotic zone only. Indeed, it is not possible to drive a realistic northward transport through surface boxes of 100 m and maintain high gradients in surface nutrient concentrations.
To confirm that the high sensitivity of pCO 2 to equatorial and low latitude nutrient utilization rates in our model is exaggerated because of the assumption that the top 500 m is biologically active, we reproduced our study in a model with an upper layer of 300 m (Fig. 8). The meridional gradients in nutrient concentrations are very high at the ocean surface so that the solutions for the interglacial state in the shallow model are, as expected, not as close to the observations as in the standard model. The preformed nutrient concentrations are lower in the interglacial and it is more difficult for the ocean to take up CO 2 . As expected, the difference between the nutrient utilization in the glacial and interglacial states is much less pronounced than in the standard model (Fig. 8,  bottom).
We have argued that our model implied importance of the nutrient utilization rates in the Equatorial region is likely to be exaggerated because in the real ocean a decrease in nutrients that are taken up in biological production at the Equator would lead to an increase in the export production of the low latitude regions to which the nutrients escaped. However, the effect of equatorial production could still be non-negligible because, although nutrients in the surface low latitudes are almost completely depleted, nutrients can travel from the equator to the north in the sub-surface below the euphotic zone.

The Southern Ocean residual circulation
One of the leading hypotheses for glacial oceanic uptake of carbon is that of a reduced residual circulation in the Southern Ocean (Keeling and Visbeck, 2001;Watson and Garabato, 2006;Fischer et al., 2010). The residual circulation is the sum of the northward Ekman transport and the southward eddy fluxes and is controlled by the surface buoyancy forcing (Marshall and Radko, 2006). Reduced surface buoyancy input would lead to a reduced residual circulation and therefore reduced upwelling of nutrients and DIC to the surface. The residence time of the nutrients in the surface is also expected to increase, enhancing the chances that they will be consumed in production rather than sink as preformed nutrients. The theory is reasonable but the processes involved in linking the residual circulation with the carbon cycle are 838 A. M. de Boer et al.: A multi-variable box model approach to the soft tissue carbon pump complex and so far only one general circulation study that we know of has attempted to model it (Parekh et al., 2006). While a weaker upwelling of DIC rich water would limit the outgassing of CO 2 in the SO, it is possible that the decreased supply of nutrients to the surface will reduce the biological uptake of DIC by a similar amount so that the net effect on atmospheric CO 2 is zero. Another complication is the effect of the reduced nutrients on the productivity in the SO. If the water masses in the upper and lower branch of the overturning in the ACC are mostly separated and the upper branch is associated with the residual circulation, then one would expect that the reduced surface nutrients would affect the Sub-Antarctic region north of the polar front more than the Antarctic region to the south. The distinction is important because a reduction in nutrients to the south would reduce the deep ocean preformed nutrient transported by AABW while at the north it would only affect the export in the Sub-Antarctic region and probably not the CO 2 . The residence time of the nutrients in the SO may affect their chance of being biologically consumed rather than becoming preformed nutrients but again, one would expect this to be important in the Antarctic region where observations show no indication of increased export. To complicate things further, the strength of deep mixing in the Antarctic regions and AABW formation is not necessarily related to the residual circulation but would modulate its influence. Parekh et al. (2006) tested the affect of the residual circulation on atmospheric CO 2 in a general circulation model coupled to a carbon cycle model. They found a drop in atmospheric pCO 2 of about 3 ppmv per Sverdrup reduction in the residual circulation. However, they changed the residual circulation by adjusting the winds in the SO from 50% to 150% of modern winds. An increase in the SO winds would indeed increase the residual circulation if the buoyancy flux is not fixed, but it would probably also decrease the stratification further south and increase AABW formation (De Boer et al., 2008). Our results suggest that their reduction in atmospheric CO 2 for weaker winds was a consequence of reduced AABW formation in their model rather than a weaker residual circulation. Also, a reduction in the residual circulation due to weakened winds will have a different effect than the more likely case in which the winds are the same but the SO eddies increased. In the former case there will be little mixing of nutrients across the ACC while in the latter case the strong winds and eddies will result in more cross frontal mixing of nutrients.
In our model we explore the effect of the residual circulation beyond that of conceptual deductions and furthermore, unlike what is possible in GCMs, we isolate the various processes that play a role in the region. Thus, we change the residual circulation by changing the eddy fluxes or the Ekman transport without changing the AABW formation and we increase or decrease the nutrient utilization rate in both the Antarctic and Sub-Antarctic boxes in conjunction with the residual circulation changes and separately. What we find is that a strongly reduced residual circulation (strong southward eddy fluxes countering the Ekman flow) is not necessary to explain the glacial state (Fig. 3). In fact, the residual circulation is 18 ± 4 Sv in the glacial states as compared to 16 ± 6 Sv in the interglacial states. The range of the residual circulation is 9 to 29 Sv for the glacial and 4 to 29 Sv for the interglacial indicating that in none of the 300 glacial states is the residual circulation less than 9 Sv. The result is interesting because the pCO 2 shows a weak but clear negative correlation with the eddy flux (Table 3) so that one would expect higher eddy transport in the glacial. The curiosity is explained when one drops the requirement for higher export production in the sub-Antarctic box during the glacial. In this case the southward eddy flux still does not increase much during the glacial (remains 9 Sv on average), but the northward Ekman transport drops from 25 Sv to 9 Sv (and the soft tissue pump efficiency shoots to 95%). The observation of enhanced sub-AA production in the glacial suggests that there was a strong northward Ekman transport to supply nutrients for the export production.

The Southern Ocean nutrient utilization and the biochemical divide
It has been proposed that the polar front divides the Southern Ocean into a Sub-Antarctic region where increased biological productivity results in increased export but not oceanic uptake of CO 2 and an Antarctic region where it does take up CO 2 (Marinov et al., 2006). We also find in our model that increased nutrient utilization is slightly more effective at reducing atmospheric CO 2 in the Antarctic box than in the Sub-Antarctic box in our interglacial states (Fig. 6, Table 3). However, increased Antarctic nutrient uptake in our model leads to a large increase in export production there and a decrease in export from the Sub-Antarctic while the observations point clearly to the opposite effect. On the other hand, the observations are not violated by an increased export from the Sub-Antarctic region but the effect on our model is small. The average nutrient utilization rates in the Sub-AA box, α S , for the glacial and interglacial states are about 8 × 10 10 s −1 and 5 × 10 10 s −1 respectively. According to our sensitivity analysis for the interglacial states (which have a larger sensitivity to α S changes than the glacial states) a 3 × 10 10 s −1 increase in the nutrient utilization rate would lead to about an 8 ppmv decrease in atmospheric pCO 2 (Table 3). That is at the lower end of the estimate of the effect of Fe fertilization on atmospheric pCO 2 given in the review of Kohfeld and Ridgwell (2009). Note that this quantitative comparison is made for interest sake only and should be interpreted with extreme care because of the simplicity of the model. Whilst the model results imply that the effect of changes in export in the Sub-Antarctic region on atmospheric CO 2 is small it may have been responsible for some of the smaller variations in the glacial CO 2 record. In the glacial states the nutrient utilization is much decreased in both boxes and of a similar magnitude. Overall our model suggests that one must look for changes in the circulation rather than the biology to explain the strong link between the SO temperature and atmospheric CO 2 in the ice core record.

Conclusions
We have formulated a simple box model that represents the major circulation variables of the real ocean and explores the effect of nutrient utilization in each of the major biological production regimes. Due to the simplicity of the model we could explore the full parameters space which includes five circulation and five nutrient utilization rate parameters to produce 10 7 solutions. Of these solutions, we determined the subsets which best fit the current observations and those that fulfil the glacial observations. In addition, we explored the sensitivity of deep ocean preformed nutrients (and by deduction the atmospheric pCO 2 ) and export production to each input parameter separately. The strong link between temperature and CO 2 in the Vostok and Epica Dome C ice core records (Petit et al., 1999;Siegenthaler et al., 2005;Lüthi et al., 2008) suggest that the SO must be a player in the glacial carbon cycle. Our model suggests the SO-CO 2 link is through a change in AABW formation or SO deep mixing and not through changes in the residual circulation or nutrient utilization rates directly. Such a reduction in AABW can be brought on by increased surface stratification that inhibits convection and deep mixing or increased deep stratification that reduces vertical mixing and upwelling of AABW. The surface stratification would have been enhanced due to the near freezing water temperature at which the density is not very sensitive to temperature and thus haloclines are not easily overturned in winter. This type of stratification was probably also active in the North Pacific and to a lesser extent the North Atlantic and may have contributed to further uptake of CO 2 through reduced mixingdriven upwelling as suggested by our model. A decrease in SO winds would also lead to increased stratification although numerical modelling studies by Tschumi et al. (2008) and Menviel et al. (2008) have not found only a small or no atmospheric CO 2 decrease from reduced or equatorward-shifted SO winds.
A reduced residual circulation through increased southward eddy fluxes is one of the leading hypotheses for glacial CO 2 uptake. In our model the atmospheric CO 2 is only weakly sensitive to the SO residual circulation. All the glacial solutions have a residual circulation of at least 9 Sv and it can be as high as 30 Sv. Further analysis indicated that this is because of the observational constraint of higher export production in the sub-AA region. The implication is that the glacial sub-AA export production is not only controlled by the nutrient utilization rate, but is limited by the supply of macro nutrients. The reason why the residual circulation in itself affects the atmospheric CO 2 only weakly can be understood in terms of the productive and unproductive circulation conceptualization of Toggweiler et al. (2003bToggweiler et al. ( , 2006. The productive circulation represents the North Atlantic overturning cell and the unproductive circulation the AABW cell. Nutrients that upwell in the SO and are transported northward in the Atlantic or Pacific will be consumed before they reach the North Atlantic deep water formation region and thus all the associated DIC that was outgassed at their surfacing is subsequently taken up again. The circulation therefore has no net effect on atmospheric CO 2 . In the unproductive southern circulation, upwelled nutrients are not depleted by the time they reach the AABW formation region and this branch is thus less productive in terms of the biological pump. A reduction in the residual circulation will indeed reduce upwelling, but it will reduce only the productive circulation to which atmospheric CO 2 is less sensitive. It may affect the nutrient supply and utilization in the SO though. The biology of the SO does not emerge as a directly important player in the biological pump. The deep ocean preformed nutrient concentration is more sensitive to changes in the nutrient utilization north of 40 • S. We do find a weak biochemical divide (similar to Marinov et al., 2006) in the sense that nutrient uptake in the sub-AA box draws down a bit less CO 2 than in the AA box. It appears irrelevant though because an increase in the nutrient utilization rate in the AA box leads to higher export production there and that is counter to the observations. The strong link between atmospheric CO 2 and temperature in Antarctic ice core records has been put forward as indicating that the SO must have played an important role in the glacial carbon cycle. Our study suggests that this link is through changes in deep mixing and AABW formation in the Antarctic region and not through changes in biology of the SO or the residual circulation. Sub-Antarctic productivity enhanced by a supply of iron may have accounted for smaller variations in the glacial CO 2 record.
The equatorial region emerges from our model results as an important region for glacial CO 2 uptake. While there may be some real sensitivity here, one should keep in mind that the dynamics are not well resolved and the sensitivity of the equatorial and low latitude regions are exaggerated due to the fact that the biologically active upper layer in the model is 500m deep. The low latitude nutrients are therefore not depleted and nutrients can easily escape from the equatorial region to the northern region where it can sink as preformed nutrients. We confirmed that deep ocean preformed nutrients are less sensitive to increased production in the equatorial and low latitude regions in a shallower upper layer version of the model. The model is simple by construction and has by design a few limitations. Changes in preformed nutrients are directly related to changes in atmospheric pCO 2 , using the theoretical arguments of Ito and Follows (2005). The theory is in essence derived from small changes to the system and not full glacial-interglacial changes. Also, the theory (as most box and numerical models) assumes a fix C:P ratio. Decoupling of the C:P and C:N ratio in organic matter could lead to an enhanced drawdown of atmospheric CO 2 (Tagliabue et al., 2009). With regard to the geometry, a potentially important simplification in the models is the 500 m depth upper layer in which nutrient utilization occurs. In the real ocean biological production does not occur below about 100 m. Export from our upper layer is thus rather the amount of organic matter that escapes the top 500 m rather than export from euphotic zone. Another geometrical simplification is the combination of Indo-Pacific and Atlantic basins into one. In a future setup the explicit representation of these basins are planned as well as the division of the deep ocean in a mid-depth and bottom cell that are fed from the North Atlantic and SO respectively. These additions would aid with the interpretation of the results as it relates to the northern high latitude regions. In reality the northern North Pacific is richer in nutrients than the North Atlantic and water does not sink there to the same depths.
The simplicity of the model means that the conclusions drawn from it are qualitative in nature and, as with all models, should be interpreted within its limitations. However, unlike General Circulation Models and previous box models, our approach to cover the whole of parameter space means that there are no tuned parameters in the model. Rather, all input parameters are deduced from observations or are variables that are explored as part of the solution. Glacial pCO 2 hypothesis are tested in all the likely glacial and interglacial states rather than one best guess and the results indicate that the glacial ocean was probably more sensitive to circulation changes and less to nutrient utilization rates than the interglacial ocean. The use of this type of multistate statistical approach can add a valuable dimension to future box studies of the carbon cycle.